从分布函数估计力矩

Question

我有一个非正态分布函数，我想计算它的矩（均值、方差、偏度和峰度）。我知道的包是 e1071 和 moments，它们计算离散值向量的力矩。是否有估计连续分布函数矩的包？例如假设我有这个分布：

tmp <- runif(100,-10,10)
PDFdata <- density(tmp , kernel = "gaussian")
PDF <- splinefun(PDFdata$x , PDFdata$y)

现在我要计算：

Answer 1

您可以使用函数 integrate 来计算力矩。
您只需要定义一个被积函数，我在下面的代码中将其称为 int。

首先，我将通过设置 RNG 种子使结果可重现。

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")

Answer 2

如果您根据数据估算密度，最好使用数据中的经验矩来估算分布矩。如果您只是将其用作函数示例，则可以使用 stats 包中的 integrate 函数。例如，

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

事实上，对于这个数据集，PDF 函数不是密度：

plot(PDF, from = -100, to = 100)

所以我们强制它在密度估计范围之外为零：

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

从分布函数估计力矩

Estimating moments from a distribution function

r

distribution

variance

skew

kurtosis