I have a non-normal distribution function which I want to calculate it's moments (mean, variance, skewness, and kurtosis).
The packages I know are e1071 and moments which calculate the moments for a vector of discrete values. Is there a package that estimates moments for a continuous distribution function?
As an example assume I have this distribution:
tmp <- runif(100,-10,10)
PDFdata <- density(tmp , kernel = "gaussian")
PDF <- splinefun(PDFdata$x , PDFdata$y)
Now I want to calculate:
You can use function integrate to compute the moments.
All you need is to define an integrand, which I call int in the code below.
First, I will make the results reproducible by setting the RNG seed.
set.seed(6126) # Make the results reproducible
tmp <- runif(100,-10,10)
PDFdata <- density(tmp , kernel = "gaussian")
PDF <- splinefun(PDFdata$x , PDFdata$y)
int <- function(x, c = 0, g = PDF, n = 1) (x - c)^n * g(x)
integrate(int, lower = -15, upper = 15, n = 1)
#-0.3971095 with absolute error < 7.8e-06
integrate(int, lower = -15, upper = 15, n = 2)
#35.76295 with absolute error < 0.0012
plot(PDFdata)
curve(PDF, -15, 15, add = TRUE, col = "blue", lty = "dashed")
If you are estimating your density from data, you're better off using empirical moments from the data to estimate the moments of the distribution. If you just used this as an example of a function, then you could use the integrate function from the stats package. For example,
set.seed(123)
tmp <- runif(100,-10,10)
PDFdata <- density(tmp , kernel = "gaussian")
PDF <- splinefun(PDFdata$x , PDFdata$y)
mean(tmp)
[1] -0.02882012
integrate(function(x) x*PDF(x), -Inf, Inf)
Error in integrate(function(x) x * PDF(x), -Inf, Inf) :
the integral is probably divergent
and indeed, for this dataset that PDF function is not a density:
plot(PDF, from = -100, to = 100)
So we force it to zero outside the range of the density estimate:
PDF2 <- function(x) ifelse(x < min(PDFdata$x) | x > max(PDFdata$x), 0,
PDF(x))
integrate(function(x) x*PDF2(x), -Inf, Inf)
-0.02900585 with absolute error < 8.2e-05
