基于R的密度估计：如何生成连续线条的密度曲线？

你的代码用的是均匀核（矩形窗）做密度估计，这种核的输出是分段常数的，所以画出来的线条是阶梯状、不连续的。要得到连续线条，你需要换成连续核函数，比如高斯核、Epanechnikov核这类。

下面是用高斯核修改后的代码，能生成连续平滑的密度曲线：

J = 512
Y = faithful$waiting
y = seq(min(faithful$waiting), max(faithful$waiting), length=J)

h = 5  # 带宽，可根据需求调整
n = length(Y)
d = rep(NA, J)

for(j in 1:J) {
  # 用高斯核计算每个y[j]处的密度估计
  kernel_vals = dnorm(y[j], mean = Y, sd = h)
  d[j] = sum(kernel_vals) / n
}

plot(y, d, type = "l", lwd = 2, xlab = "Waiting time", ylab = "Density")

关键修改说明： #

• 把原来的均匀核区间判断替换成了高斯核函数dnorm()，它会为每个样本点计算在y[j]处的连续核贡献值，最终得到的密度估计是连续平滑的。
• 高斯核的密度估计公式是：$\hat{f}(y) = \frac{1}{n} \sum_{i=1}^n \frac{1}{h\sqrt{2\pi}} \exp\left(-\frac{(y - Y_i)^2}{2h2}\right)$，代码里dnorm(y[j], mean=Y, sd=h)已经包含了$\frac{1}{h\sqrt{2\pi}}$的部分，所以直接求和除以样本量n即可。

如果你想尝试其他连续核，比如Epanechnikov核，也可以替换循环内的计算：

# Epanechnikov核的实现
for(j in 1:J) {
  u = (y[j] - Y) / h
  kernel_vals = ifelse(abs(u) <= 1, 0.75*(1 - u^2)/h, 0)
  d[j] = sum(kernel_vals) / n
}

另外，R里也有现成的密度估计函数density()，可以直接用它生成连续曲线，用来验证自己实现的结果：

plot(density(faithful$waiting, bw = 5), lwd = 2, col = "red")
lines(y, d, lwd = 2, col = "blue")  # 和自己实现的高斯核结果对比

内容的提问来源于stack exchange，提问作者Emiliano Rinato