R:二元(或多变量)核密度的概率/数值积分

Question

我正在使用包 ks 进行核密度估计。这是一个简单的例子:

n <- 70
x <- rnorm(n)

library(ks)
f_kde <- kde(x) 

我实际上对输入数据各自的超出概率感兴趣，这可以通过具有 f_kde 的 ks 轻松返回:

p_kde <- pkde(x, f_kde)

这是使用辛普森规则在 ks 中通过数值积分完成的。不幸的是，他们只针对 1d 案例实现了这一点。在双变量情况下，ks 中没有任何返回概率的方法的实现:

y <- rnorm(n)
f_kde <- kde(data.frame(x,y))
# does not work, but it's what I am looking for:
p_kde <- pkde(data.frane(x,y), f_kde) 

我找不到任何包或帮助在 stackoverflow 中搜索来解决 R 中的这个问题(存在一些针对 Python 的建议，但我想将其保留在 R 中)。任何一行代码或包推荐都会受到赞赏。尽管我对双变量情况最感兴趣，但任何有关多变量情况的想法也值得赞赏。

Answer 1

kde 允许 multidimensional kernel estimate ，所以我们可以使用kde来计算pkde。
为此，我们使用 eval.points 参数在足够小的 dx 和 dy 步长上计算 kde:这给了我们dx*dy 上的局部密度估计 正方形。
我们验证了估计值乘以正方形面积的总和几乎等于 1:

library(ks)
set.seed(1)
n <- 10000
x <- rnorm(n)
y <- rnorm(n)
xy <- cbind(x,y)

xmin <- -10
xmax <- 10
dx <- .1

ymin <- -10
ymax <- 10
dy <- .1

pts.x <- seq(xmin, xmax, dx)
pts.y <- seq(ymin, ymax, dy)
pts <- as.data.frame(expand.grid(x = pts.x, y = pts.y))
f_kde <- kde(xy,eval.points=pts)

pts$est <- f_kde$estimate

sum(pts$est)*dx*dy
[1] 0.9998778

您现在可以查询 pts 数据框以获取您选择的区域的累积概率:

library(data.table)
setDT(pts)
# cumulative density
pts[x < 1 & y < 2 , .(pkde=sum(est)*dx*dy)]
        pkde
1: 0.7951228

# average density around a point
tolerance <-.1
pts[pmin(abs(x-1))<tolerance & pmin(abs(y-2))<tolerance, .(kde = mean(est))]
          kde
1: 0.01465478