28 发生率的置信区间
log-Wald 与精确 Poisson：推导、比较与交互计算
log-Wald 和精确 Poisson 都是常用的单个发生率区间估计方法。 前者使用正态近似，后者直接利用 Poisson 分布。它们有共同的模型前提，不能无条件用于所有称为“率”的指标。
本文以事件数 \(a=9\)、总人时 \(T=4028.16\) 人年为初始例子。原计算器用 \(N\) 表示人时；本文改用 \(T\)，避免与人数混淆。
在终端运行 quarto render incidence-rate-ci.qmd，或在 RStudio / Positron 中点击 Render。交互计算器在生成的 HTML 中运行；不需要 R、Shiny 服务器或 Crosstalk。公式由 Quarto 排版，滑块与计算逻辑直接写在本文中。本文设置不执行代码，后面的 R 代码用于自行复核。
28.1 先动手：两种方法算出的区间有何不同？
拖动事件数和置信水平滑块，或直接修改输入。区间线段共用同一横轴；点表示相同的发生率估计 \(a/T\)。图会自动调整横轴范围，跨不同输入比较时请同时看数值。
| 方法 | 下限 | 上限 | 区间宽度 |
|---|
查看当前输入对应的中间计算值
28.2 模型前提与所需数据
令随机事件数为 \(A\)，实际观察到的值为 \(a\)。假定
\[ A\sim\operatorname{Poisson}(rT),\qquad \mu=rT,\qquad \widehat r=\frac{a}{T}. \]
其中 \(r\) 是目标发生率，\(T\) 是已知总人时。需要的数据只有观察事件数、总人时和置信水平；在此模型下不需要个体原始数据。
Poisson 模型要求事件计数的分布可以合理地用 Poisson 描述，满足 \(\operatorname{Var}(A)=E(A)=rT\)。独立事件及给定风险率下的计数过程是典型情形。风险率可以随时间变化，只要总计数仍符合 Poisson 分布，此处的 \(r\) 可以解释为观察期平均发生率；不能据此推断具体的年龄或季节效应。
下列情况需另外处理：
- 累计发病风险 / 比例：每人固定随访期内是否至少发生一次事件，常用二项分布区间；即使随访时长相同，也不能仅把人数当人时就认为精确 Poisson 适用。
- 聚集或相关事件：如同一儿童反复感染、家庭聚集、研究间异质性，可能需要稳健方差、负二项模型、混合模型或适合设计的 bootstrap。
- 加权、校正或模型预测的病例数：不能直接当作原始 Poisson 计数；需传播模型、抽样和其他输入的不确定性。
- 人时本身不确定：本文将分母视为已知，未估计分母误差。
对你的 RSV 工作，直接观察的病例数和队列总人时可按上述假设评估；国家负担模型的预测率及其不确定性区间应从完整模型的抽样结果计算。
28.3 log-Wald：delta 法为什么得到 \(1/\sqrt a\)？
28.3.1 从事件数的方差到发生率的方差
由于 \(T\) 是已知常数，
\[ \operatorname{Var}(\widehat r) =\operatorname{Var}\left(\frac{A}{T}\right) =\frac{\operatorname{Var}(A)}{T^2} =\frac{rT}{T^2} =\frac rT. \]
真实的 \(r\) 未知，用 \(\widehat r=a/T\) 代入，得到估计的方差和标准误：
\[ \widehat{\operatorname{Var}}(\widehat r)=\frac{a}{T^2}, \qquad \widehat{SE}(\widehat r)=\frac{\sqrt a}{T}. \]
这里 \(\operatorname{Var}(A)=rT\) 是模型等式；用观察值 \(a\) 替代未知均值 \(rT\) 才是代入估计。
28.3.2 一阶 Taylor 展开就是 delta 法的出发点
对可微函数 \(g\)，在真实参数 \(r\) 附近展开：
\[ g(\widehat r) =g(r)+g'(r)(\widehat r-r) +\frac12g''(r)(\widehat r-r)^2+\cdots. \]
忽略二阶及更高阶项：
\[ g(\widehat r)\approx g(r)+g'(r)(\widehat r-r). \]
右侧的 \(g(r)\) 和 \(g'(r)\) 是常数。利用 \(\operatorname{Var}(bX+c)=b^2\operatorname{Var}(X)\)，
\[ \boxed{\operatorname{Var}[g(\widehat r)] \approx [g'(r)]^2\operatorname{Var}(\widehat r).} \]
因此，delta 法把非线性变换的局部波动近似为其切线的波动，导数描述变换对误差的放大或缩小。它是近似方差的方法，不会单独保证变换后的分布服从正态。牛津大学的统计理论讲义给出了这一推导。
28.3.3 令 \(g(r)=\log r\)
当 \(r>0\) 时，\(g'(r)=1/r\)，所以
\[ \operatorname{Var}[\log(\widehat r)] \approx\left(\frac1r\right)^2\frac rT =\frac1{rT}. \]
再用 \(\widehat rT=a\) 代入：
\[ \boxed{\widehat{SE}[\log(\widehat r)]\approx\frac1{\sqrt a},\qquad a>0.} \]
为什么没有 \(T\)？在 Poisson 假设下，log 尺度的标准误对应相对不确定性，由事件数控制。在事件数不变时改变 \(T\)，会改变率及其绝对区间宽度，但不会改变相对区间宽度。
28.3.4 在 log 尺度构造区间，再指数变换
令 \(z=z_{1-\alpha/2}\)。使用 log 尺度的渐近正态近似：
\[ \log(\widehat r)\pm z\frac1{\sqrt a}. \]
因为指数函数单调递增，可对两个端点分别取指数：
\[ \boxed{ L_{\log W}=\frac aT\exp\left(-\frac z{\sqrt a}\right),\qquad U_{\log W}=\frac aT\exp\left(\frac z{\sqrt a}\right). } \]
区间在 log 尺度对称，在原始尺度呈乘法对称：\(L_{\log W}U_{\log W}=\widehat r^2\)。因此上端距离点估计更远，且下限始终为正。正的下限并不等于准确的覆盖率；事件很少时正态近似仍可能不可靠。
当 \(a=0\) 时，\(\log(\widehat r)\) 与 \(1/\sqrt a\) 均无定义，不能使用这条公式。随意加 0.5 会改变方法，并不是原始 log-Wald。
28.4 精确 Poisson：不依赖 log 正态近似
28.4.1 先求事件均值 \(\mu=rT\) 的区间
中央等尾的 Garwood 精确区间通过反演 Poisson 尾概率得到。对 \(a>0\)，下端点满足
\[ P_{\mu_L}(A\ge a)=\frac\alpha2; \]
上端点满足
\[ P_{\mu_U}(A\le a)=\frac\alpha2. \]
利用 Poisson 尾概率与 Gamma / 卡方分布的关系，可写为
\[ \mu_L=\frac12\chi^2_{2a,\,\alpha/2},\qquad \mu_U=\frac12\chi^2_{2(a+1),\,1-\alpha/2}. \]
\(\chi^2_{\nu,p}\) 指自由度为 \(\nu\) 的卡方分布的左尾 \(p\) 分位数。等价地，\(\mu_L\) 和 \(\mu_U\) 分别是 \(\operatorname{Gamma}(a,\mathrm{scale}=1)\) 与 \(\operatorname{Gamma}(a+1,\mathrm{scale}=1)\) 的相应分位数。交互计算器使用这一等价表达式计算。
28.4.2 两端除以总人时
\[ \boxed{ L_{\mathrm{exact}}=\frac{\chi^2_{2a,\,\alpha/2}}{2T},\qquad U_{\mathrm{exact}}=\frac{\chi^2_{2(a+1),\,1-\alpha/2}}{2T}. } \]
这里的“精确”是指：在 Poisson 模型成立时，重复抽样覆盖率不低于名义置信水平，而不是每个参数值处都恰好为 95%。由于计数离散，Garwood 区间往往保守。关于离散分布精确区间的研究讨论了 Garwood 区间及其性质。
28.4.3 零事件怎么处理？
当 \(a=0\) 时，下限单独定义为 0；不能将自由度为 0 的卡方分位数直接代入一般公式。因为 \(P_\mu(A=0)=e^{-\mu}\)，双侧区间的上限为
\[ \boxed{L=0,\qquad U=\frac{-\log(\alpha/2)}T.} \]
双侧 95% 上限为 \(3.688879/T\)。如果求的是单侧 95% 上限，则为 \(-\log(0.05)/T=2.995732/T\)，约为常见的“rule of three”。这两个数对应不同的问题，不能混用。
28.5 怎么选？还有哪些方法？
| 方法 | 核心思路 | 主要特点 |
|---|---|---|
| log-Wald | log 尺度上做 Wald 区间 | 计算简单、下限为正；小事件数下依赖近似；零事件未定义 |
| 精确 Poisson / Garwood | 反演 Poisson 尾概率 | 支持零事件；Poisson 模型下保证至少名义覆盖率，常偏保守 |
| 普通 Wald | \(a/T\pm z\sqrt a/T\) | 原尺度对称；小事件数可能有负下限，零事件退化为 \([0,0]\) |
| Poisson Score | 反演采用假设均值方差的 score 检验 | 区间非负；仍是近似法，不保证每个参数值处达到名义覆盖率 |
| 似然比区间 | 反演似然比检验 | 可体现似然不对称；常用卡方近似，零事件等边界情形需谨慎 |
| mid-P | 精确尾概率扣除一半观察值概率 | 通常减轻保守性；可能低于名义覆盖率，不能当作同样的严格精确区间 |
| 平方根 / Anscombe | 方差稳定变换后构造区间 | 另一类近似；不同修正版公式和边界处理有所差异 |
| 贝叶斯 Gamma–Poisson | 由先验与似然得到后验分位数 | 得到的是可信区间 CrI；解释与频率学 CI 不同，需说明先验 |
statsmodels 官方文档提供精确、Score、Wald、mid-P、平方根及 Jeffreys 等实现和方法比较参考文献。
实用选择： 若是单个原始事件率，尤其事件数较少或为零，Garwood 精确区间是容易说明和复核的选择。log-Wald 可以用于事件数足够大、近似合理的情形。没有一个通用的事件数阈值能保证近似良好；方法应预先选定并说明，不能按哪个区间更窄挑选。
95% CI 的频率学解释是：如果重复同样的抽样与计算过程，所构造的区间中约 95% 会覆盖固定的真实参数；不能直接解释成“给定这次数据，真实参数有 95% 概率在这个区间内”。
28.6 用 R 复核
下面仅使用 R 的基础 stats 函数；代码展示，不会在渲染本文时执行。
a <- 9L
T <- 4028.16 # person-years
conf_level <- 0.95
multiplier <- 1000 # report per 1,000 person-years
alpha <- 1 - conf_level
z <- qnorm(1 - alpha / 2)
rate <- a / T
log_wald <- if (a > 0) {
rate * exp(c(-1, 1) * z / sqrt(a))
} else {
c(NA_real_, NA_real_)
}
# Garwood exact Poisson CI; also valid for a = 0
exact <- poisson.test(a, T = T, conf.level = conf_level)$conf.int
# Equivalent chi-square formula, with a separate zero-count lower limit
exact_chisq <- c(
if (a == 0) 0 else qchisq(alpha / 2, df = 2 * a) / (2 * T),
qchisq(1 - alpha / 2, df = 2 * (a + 1)) / (2 * T)
)
stopifnot(isTRUE(all.equal(as.numeric(exact), exact_chisq)))
data.frame(
method = c("log-Wald", "Exact Poisson (Garwood)"),
rate = rep(rate * multiplier, 2),
lower = c(log_wald[1], exact[1]) * multiplier,
upper = c(log_wald[2], exact[2]) * multiplier
)poisson.test() 的事件数、人时和置信水平参数见 R 官方文档。
28.7 参考来源与计算说明
- CUHK：C.I. Calculator—Single Incidence Rate。其 \(a=9,T=4028.16\) 示例结果与 log-Wald 公式近似吻合；这是根据示例的判断。本文直接实现所列公式，不保证逐位复现旧网页的舍入或近似算法。
- Oxford：Bias and variance approximations—the delta method。一阶 Taylor 方差近似。
- R：
poisson.test官方文档。精确 Poisson 检验与区间接口。 - statsmodels：
confint_poisson官方文档。常用 Poisson 区间及比较文献。 - Thulin & Zwanzig：Exact confidence intervals and hypothesis tests for parameters of discrete distributions。精确区间的覆盖率与 Garwood 区间性质。
交互计算器用正态分位数近似以及 Gamma CDF 的数值反演求端点，显示至最多 6 位小数；统计意义上的“精确 Poisson”不等于浮点运算没有误差。人时单位固定为人年；若输入人月，可将总人月除以 12 后输入。仅改变每 1 / 1,000 / 100,000 人年的显示单位，不会改变区间的相对宽度。