Eri*_*ric 7 r probability kernel-density density-plot probability-density
我density对我的数据有一个密度估计(使用函数)learningTime(见下图),我需要找到概率Pr(learningTime > c),即从给定数字c(红色垂直线)到曲线末端的密度曲线下面积.任何的想法?
密度估计曲线下的计算区域并不困难.这是一个可重复的例子.
假设我们有一些观察到的数据x,为简单起见,正常分布:
set.seed(0)
x <- rnorm(1000)
Run Code Online (Sandbox Code Playgroud)
我们进行密度估算(有一些定制,见?density):
d <- density.default(x, n = 512, cut = 3)
str(d)
# List of 7
# $ x : num [1:512] -3.91 -3.9 -3.88 -3.87 -3.85 ...
# $ y : num [1:512] 2.23e-05 2.74e-05 3.35e-05 4.07e-05 4.93e-05 ...
# ... truncated ...
Run Code Online (Sandbox Code Playgroud)
我们想要计算右边曲线下面积x = 1:
plot(d); abline(v = 1, col = 2)
Run Code Online (Sandbox Code Playgroud)
在数学上,这是估计密度曲线的数值积分[1, Inf].
估计的密度曲线以离散格式存储在d$x和中d$y:
xx <- d$x ## 512 evenly spaced points on [min(x) - 3 * d$bw, max(x) + 3 * d$bw]
dx <- xx[2L] - xx[1L] ## spacing / bin size
yy <- d$y ## 512 density values for `xx`
Run Code Online (Sandbox Code Playgroud)
数值积分有两种方法.
方法1:黎曼和
估计密度曲线下面积为:
C <- sum(yy) * dx ## sum(yy * dx)
# [1] 1.000976
Run Code Online (Sandbox Code Playgroud)
由于黎曼和只是一个近似值,因此它略微偏离1(总概率).我们称这个C值为"归一化常数".
数值积分[1, Inf]可近似为
p.unscaled <- sum(yy[xx >= 1]) * dx
# [1] 0.1691366
Run Code Online (Sandbox Code Playgroud)
应进一步按比例缩放以C进行适当的概率估计:
p.scaled <- p.unscaled / C
# [1] 0.1689718
Run Code Online (Sandbox Code Playgroud)
由于我们模拟的真密度x已知,我们可以将此估计值与真实值进行比较:
pnorm(x0, lower.tail = FALSE)
# [1] 0.1586553
Run Code Online (Sandbox Code Playgroud)
这是相当接近的.
方法2:梯形规则
我们对这种线性插值进行线性插值(xx, yy)并应用数值积分.
f <- approxfun(xx, yy)
C <- integrate(f, min(xx), max(xx))$value
p.unscaled <- integrate(f, 1, max(xx))$value
p.scaled <- p.unscaled / C
#[1] 0.1687369
Run Code Online (Sandbox Code Playgroud)
关于罗宾的回答
答案是合法的,但可能是作弊.OP的问题始于密度估计,但答案完全绕过了它.如果允许,为什么不简单地执行以下操作?
set.seed(0)
x <- rnorm(1000)
mean(x > 1)
#[1] 0.163
Run Code Online (Sandbox Code Playgroud)
ecdf()以 R 为基础的经验累积分布函数使它变得非常容易。使用???的例子...
#Reproducible sample data
set.seed(0)
x <- rnorm(1000)
#Create empirical cumulative distribution function from sample data
d_fun <- ecdf (x)
#Assume a value for the "red vertical line"
x0 <- 1
#Area under curve less than, equal to x0
d_fun(x0)
# [1] 0.837
#Area under curve greater than x0
1 - d_fun(x0)
# [1] 0.163
Run Code Online (Sandbox Code Playgroud)
关于 ??? 对我的回答的回应。他们的回答假设您只有密度估计曲线。我的回答假设您拥有原始数据,这适用于 OP 的问题,因为它们曾经density()用于获取密度估计曲线。