获取ECDF的导数

Question

获取ECDF的导数

3

是否可能区分ECDF？以以下示例中获得的为例。

set.seed(1)

a <- sort(rnorm(100))
b <- ecdf(a)

plot(b)

我希望对b进行求导，以获得其概率密度函数（PDF）。

- MaxPlank

1

由于ECDF是一个阶梯函数，它在观察值处不可导，并且在其他地方的导数为0。 - Mikko Marttila

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Zheyuan Li · Accepted Answer

n <- length(a)  ## `a` must be sorted in non-decreasing order already
plot(a, 1:n / n, type = "s")  ## "staircase" plot; not "line" plot

但我正在寻找b的导数

在基于样本的统计学中，对于连续随机变量的估计密度并不是通过对ECDF进行微分得到的，因为样本大小是有限的而且ECDF不可微。相反，我们直接估计密度。我猜您真正需要的是plot(density(a))。

几天后..

警告：以下仅为数值解决方案，没有统计依据！

我将其视为练习学习R软件包scam用于形状约束加性模型，这是Wood教授早期博士生Pya博士的子软件包。

逻辑如下：

使用scam::scam，对ECDF拟合单调递增的P样条（您必须指定想要的结节点数）; [请注意，单调性不是唯一的理论限制。需要对平滑的ECDF进行"剪裁"，即在两个边缘上剪裁：在0处的左边缘和在1处的右边缘。我目前使用weights来施加这种约束，通过在两个边缘处给予非常大的权重]
使用stats::splinefun，通过节点和结节点预测值拟合单调插值样条;
返回插值样条函数，该函数还可以计算1、2和3阶导数。

为什么我希望这能够工作：

随着样本大小增长，

ECDF收敛于CDF；
P-spline是一致的，因此平滑的ECDF将越来越无偏差；
平滑的ECDF的一阶导数将越来越无偏估计PDF。

请小心使用：

您必须自己选择结节点数；
导数未经过规范化处理，使曲线下面积为1；
结果可能相当不稳定，并且仅适用于大样本量。

函数参数：

x：样本向量；
n.knots：结节点数量；
n.cells：绘制导数函数时网格点的数量

您需要从CRAN安装scam软件包。

library(scam)

test <- function (x, n.knots, n.cells) {

  ## get ECDF
  n <- length(x)
  x <- sort(x)
  y <- 1:n / n
  dat <- data.frame(x = x, y = y)  ## make sure `scam` can find `x` and `y`

  ## fit a monotonically increasing P-spline for ECDF
  fit <- scam::scam(y ~ s(x, bs = "mpi", k = n.knots), data = dat,
                    weights = c(n, rep(1, n - 2), 10 * n))
  ## interior knots
  xk <- with(fit$smooth[[1]], knots[4:(length(knots) - 3)])
  ## spline values at interior knots
  yk <- predict(fit, newdata = data.frame(x = xk))
  ## reparametrization into a monotone interpolation spline
  f <- stats::splinefun(xk, yk, "hyman")

  par(mfrow = c(1, 2))

  plot(x, y, pch = 19, col = "gray")  ## ECDF
  lines(x, f(x), type = "l")          ## smoothed ECDF
  title(paste0("number of knots: ", n.knots,
               "\neffective degree of freedom: ", round(sum(fit$edf), 2)),
        cex.main = 0.8)

  xg <- seq(min(x), max(x), length = n.cells)
  plot(xg, f(xg, 1), type = "l")     ## density estimated by scam
  lines(stats::density(x), col = 2)  ## a proper density estimate by density

  ## return smooth ECDF function
  f
  }

## try large sample size
set.seed(1)
x <- rnorm(1000)
f <- test(x, n.knots = 20, n.cells = 100)

f 是由stats::splinefun返回的函数（请阅读 ?splinefun）。

一种简单的类似解决方案是在ECDF上进行插值样条，而不进行平滑处理。但这是一个非常糟糕的想法，因为我们没有一致性。

g <- splinefun(sort(x), 1:length(x) / length(x), method = "hyman")
curve(g(x, deriv = 1), from = -3, to = 3)

提醒：强烈建议使用stats::density进行直接密度估计。