Search code examples
rinversecdfkolmogorov-smirnov

R - Inverse cumulative distribution method with given function


I have a given function (let's call it f(x)) and I used the Monte Carlo method to normalized it. I calculated the probability density function and I got the cumulative distribution function from integrating that.

f = function(x) ...

plot(f,xlim = c(0, 5), ylim = c(0, 1),main="f(x)")

mc.integral = function(f, n.iter = 1000, interval){
  x = runif(n.iter, interval[1], interval[2])
  y = f(x)
  mean(y)*(interval[2] - interval[1])
}

MC = mc.integral(f, interval = c(0, 8))
print(MC)

densityFunction <- function(x){
  return ((f(x)/MC)) 
}

distributionFunction <- function(x){
  return  (integrate(densityFunction,0,x)$value)
}

vd <- Vectorize(distributionFunction)
plot(vd,xlim = c(0, 8), ylim = c(0, 1),ylab = "y",main="E(f(x))")

Now my next task is to use the inverse transform method / inverse cumulative distribution method to generate samples and test it with the Kolmogorov-Smirnov Test, but I don't know how should I do in R.

Can you please give me some help?


Solution

  • Well, this thread shows us how to generate a sample using the inverse transform method:

    sample <- vd(runif(1000))
    
    > head(sample)
    [1] 0.28737403 0.59295499 0.30814305 0.27998306 0.07601228 0.52753327
    

    Therefore, generating 10 different random samples could be done with:

    sample <- list()
    for(i in 1:10){
      set.seed(i)
      sample[[i]] <- vd(runif(1000))
    }
    

    Afterwards, loop ks.test over the list:

    lapply(sample, function(x) ks.test(x, pnorm))
    

    will give you the output of a test vs. normality for each sample. Choose the size of your samples wisely, as most tests for normality are prone to be significant for large samples even with small differences (reference here).