如何估计 R 中零膨胀参数的密度?

机器算法验证 r 可能性 克德
2022-03-10 14:54:32

我有一个包含很多零的数据集,如下所示:

set.seed(1)
x <- c(rlnorm(100),rep(0,50))
hist(x,probability=TRUE,breaks = 25)

我想为它的密度画一条线,但该density()函数使用一个移动窗口来计算 x 的负值。

lines(density(x), col = 'grey')

有一个density(... from, to)论点,但这些似乎只是截断计算,而不是改变窗口,以便 0 处的密度与数据一致,如下图所示:

lines(density(x, from = 0), col = 'black')

(如果插值发生了变化,我希望黑线在 0 处的密度高于灰线)

这个函数是否有替代方法可以更好地计算零时的密度?

在此处输入图像描述

4个回答

密度在零处是无限的,因为它包含一个离散的尖峰。您需要使用零的比例来估计尖峰,然后假设它是平滑的,估计密度的正部分。KDE 会在左端引起问题,因为它会给负值施加一些权重。一种有用的方法是转换为对数,使用 KDE 估计密度,然后再转换回来。参见Wand, Marron & Ruppert (JASA 1991)作为参考。

以下 R 函数将执行转换后的密度:

logdensity <- function (x, bw = "SJ") 
{
    y <- log(x)
    g <- density(y, bw = bw, n = 1001)
    xgrid <- exp(g$x)
    g$y <- c(0, g$y/xgrid)
    g$x <- c(0, xgrid)
    return(g)
}

然后下面将给出你想要的情节:

set.seed(1)
x <- c(rlnorm(100),rep(0,50))
hist(x,probability=TRUE,breaks = 25)
fit <- logdensity(x[x>0]) # Only take density of positive part
lines(fit$x,fit$y*mean(x>0),col="red") # Scale density by proportion positive
abline(v=0,col="blue") # Add spike at zero.

在此处输入图像描述

我同意 Rob Hyndman 的观点,即您需要分别处理零。有几种方法可以处理有界支持的变量的核密度估计,包括“反射”、“重新归一化”和“线性组合”。这些似乎没有在 R 的density函数中实现,但在Benn Jann 的kdensStata 包中可用。

当您的数据具有逻辑下限(例如 0,但可能是其他值)并且您知道数据不会低于并且常规内核密度估计将值置于该界限以下(或者如果您有上限)时,另一种选择,或两者)是使用对数样条估计。R 的 logspline 包实现了这些,并且函数具有用于指定边界的参数,因此估计值将达到边界,但不会超出边界,并且仍然缩放为 1。

还有一些方法(oldlogspline函数)会考虑间隔审查,所以如果那些 0 不是精确的 0,而是四舍五入,以便您知道它们代表 0 和其他数字(例如检测限制)之间的值,那么您可以将该信息提供给拟合函数。

如果额外的 0 是真 0(未四舍五入),则估计尖峰或点质量是更好的方法,但也可以与对数样条估计相结合。

您可以尝试降低带宽(蓝线代表adjust=0.5), 在此处输入图像描述

但可能 KDE 并不是处理此类数据的最佳方法。