127 lines
7.2 KiB
Markdown
127 lines
7.2 KiB
Markdown
|
# bootstrap_ci
|
|||
|
|
|
|||
|
|
| 项 | 内容 |
|
|||
|
|
|---|---|
|
|||
|
|
| 模块 | `mcm_toolkit.stats.bootstrap` |
|
|||
|
|
| 档位 | S |
|
|||
|
|
| 参考来源 | Efron (1979)；Efron & Tibshirani (1993) *An Introduction to the Bootstrap*；
|
|||
|
|
| | `scipy.stats.bootstrap`（只做对照，没有调用） |
|
|||
|
|
| 状态 | 实现已落地（2026-10-03），§10 有 4 个口径待拍板 |
|
|||
|
|
|
|||
|
|
## 1. 数学定义
|
|||
|
|
|
|||
|
|
给定样本 `x_1..x_n` 与统计量 `T(·)`：从 `x` 中**有放回**抽 `n` 个得到 `x*`，重复 `B` 次，
|
|||
|
|
得 `T*_1..T*_B`。百分位法区间
|
|||
|
|
|
|||
|
|
```text
|
|||
|
|
estimate = T(x) # 用**原样本**，不是重采样分布的均值
|
|||
|
|
alpha = 1 - confidence
|
|||
|
|
low = T*_(alpha/2) # 经验分布的 alpha/2 分位
|
|||
|
|
high = T*_(1 - alpha/2)
|
|||
|
|
```
|
|||
|
|
|
|||
|
|
输入是一个一维样本，输出是「点估计 + 区间」两个标量 —— 这是**非参数**区间，
|
|||
|
|
不假设正态、也不需要 `T` 的标准误解析式。
|
|||
|
|
|
|||
|
|
## 2. 适用场景
|
|||
|
|
|
|||
|
|
**用**：`T` 没有解析分布（中位数、极差、自定义指标）；样本有偏或分布未知；
|
|||
|
|
想给区间但不想假设正态（`mean ± z·s/√n` 的替代品）。
|
|||
|
|
|
|||
|
|
**不用**：`n < 10` —— 经验分布太粗，`x` 能重采样出的不同样本极少；
|
|||
|
|
数据有**强时间相关性**（股票、气温），此时必须用 block bootstrap；
|
|||
|
|
要估的是**极值**（最大值本身就是上界，见 §5）。
|
|||
|
|
|
|||
|
|
## 3. 参数与返回值
|
|||
|
|
|
|||
|
|
| 参数 | 类型 | 默认 | 含义 | 取值影响 |
|
|||
|
|
|---|---|---|---|---|
|
|||
|
|
| `x` | Any | — | 样本；走 `core.as_array(x, numeric_only=True, ndim=1)` | DataFrame 只取数值列，再摊平 |
|
|||
|
|
| `statistic` | Callable | `np.mean` | 任意 `f(sample) -> float` | 换 `np.median` 就得到中位数的区间 |
|
|||
|
|
| `confidence` | float | 0.95 | 置信水平 | 越大区间越宽 |
|
|||
|
|
| `n_resamples` | int | 10000 | 重采样次数 `B` | 越大端点越稳；10000 通常够 |
|
|||
|
|
| `method` | str | `"percentile"` | 区间算法 | `"bca"` 尚未实现，传了抛 `NotImplementedError` |
|
|||
|
|
| `random_state` | int \| None | `None` | 随机种子 | `None` 时复用 `core.get_rng` 的模块级生成器且会推进 |
|
|||
|
|
|
|||
|
|
返回 `Result`：
|
|||
|
|
|
|||
|
|
| 字段 | 内容 |
|
|||
|
|
|---|---|
|
|||
|
|
| `params` | `estimate` / `low` / `high` / `confidence` / `n_resamples` |
|
|||
|
|
| `metrics` | `statistic_name`（`np.mean` → `"mean"`；lambda → `"custom"`）|
|
|||
|
|
| `data` | 不使用 |
|
|||
|
|
|
|||
|
|
## 4. 复杂度
|
|||
|
|
|
|||
|
|
- 时间：`O(B · n)`（`B` 次重采样，每次都要在 `n` 个值上算一次 `T`）
|
|||
|
|
- 空间：`O(B)` —— **按块处理**，每批 `batch_size = 1000` 行，绝不落地 `B × n` 的整矩阵
|
|||
|
|
- **什么时候会跑不完**：
|
|||
|
|
- `n = 1e4, B = 1e4` → `1e8` 次操作，秒级，正常；
|
|||
|
|
- `n = 1e6, B = 1e4` → `1e10` 次操作，**跑不完**；而且 1000×1e6 的单批索引矩阵就要 8 GB，
|
|||
|
|
此时应改用解析近似，或把 `B` 降到 1e3 以下、换更便宜的 `T`。
|
|||
|
|
|
|||
|
|
## 5. 数值稳定性
|
|||
|
|
|
|||
|
|
- **统计量有上界时，百分位区间会「塌」到上界。** 极差 `max - min` 的重采样值**永远不可能超过**
|
|||
|
|
原样本极差，`B` 个值里有相当比例正好等于它（抽到原样本极值点的那些），
|
|||
|
|
于是 `high == estimate` 精确成立。这不是 bug，是 bootstrap 的固有性质 ——
|
|||
|
|
区间端点是**样本**的经验分位数，不能给出样本之外的结论。要估极值请用极值理论。
|
|||
|
|
- 分块只是省内存，不改结果：每批都从同一个 `rng` 按顺序抽索引，等价于一次抽 `B × n`。
|
|||
|
|
- 随机数**只来自 `core.get_rng`**（决策 24）：给了种子可复现，不给则跟随 `set_seed`。
|
|||
|
|
- **NaN 目前不处理**：`x` 里的 NaN 会直接进 `T`，`np.mean` 于是给出 `nan` 区间。见 §10。
|
|||
|
|
|
|||
|
|
## 6. 常见错误
|
|||
|
|
|
|||
|
|
1. 不写 seed → 每次跑出来的区间都不一样，论文里的数字对不上 → 传 `random_state` 或先 `set_seed`。
|
|||
|
|
2. 以为 bootstrap 能得出**比样本更极端**的结论 → 区间端点就是样本的经验分位数，
|
|||
|
|
"95% 置信上限 173.2cm" 不可能超过样本最大值。
|
|||
|
|
3. 拿它估最大值 / 最小值 → 区间会塌成一个点（原因同 §5），该换极值理论。
|
|||
|
|
4. `n < 10` 硬跑 → 能抽出的不同重采样样本太少，区间没有统计意义。
|
|||
|
|
5. 把 `n_resamples` 当成"样本量"写进论文 → 论文里要同时写清 `n`、`B`、seed、method。
|
|||
|
|
|
|||
|
|
## 7. 最小示例
|
|||
|
|
|
|||
|
|
```python
|
|||
|
|
import numpy as np
|
|||
|
|
|
|||
|
|
from mcm_toolkit.stats import bootstrap_ci
|
|||
|
|
|
|||
|
|
r = bootstrap_ci(heights, statistic=np.median, n_resamples=10000, random_state=42)
|
|||
|
|
print(r.params["low"], r.params["estimate"], r.params["high"])
|
|||
|
|
print(r.summary()) # 可直接贴进论文的 Markdown 表
|
|||
|
|
```
|
|||
|
|
|
|||
|
|
## 8. 单元测试要点
|
|||
|
|
|
|||
|
|
| 测什么 | 期望 | 对应测试 |
|
|||
|
|
|---|---|---|
|
|||
|
|
| 大样本正态，`statistic=np.mean` | 区间半宽 ≈ `2 × 1.96 × se`，容差 10% | `test_bootstrap_ci_matches_normal_theory_on_large_n` |
|
|||
|
|
| 任意统计量（中位数，对数正态样本） | 区间包含样本中位数 | `test_bootstrap_ci_works_for_arbitrary_statistic` |
|
|||
|
|
| 同种子两次 | `low` 完全相等 | `test_bootstrap_ci_is_reproducible` |
|
|||
|
|
| 不给种子 | 跟随 `set_seed`（决策 24） | `test_bootstrap_ci_follows_get_rng_when_no_seed_is_given` |
|
|||
|
|
| `estimate` | 等于原样本上的 `T(x)` | `test_bootstrap_ci_reports_the_original_estimate` |
|
|||
|
|
| `confidence` 0.80 → 0.95 | 区间单调变宽 | `test_bootstrap_ci_interval_widens_with_confidence` |
|
|||
|
|
| 极差作为统计量 | `high <= estimate`（原子点塌陷） | `test_bootstrap_ci_works_for_a_custom_lambda` |
|
|||
|
|
| 五种容器 | DataFrame/Series/ndarray/list/tuple 结果一致 | `test_bootstrap_ci_accepts_the_documented_containers` |
|
|||
|
|
| `method` 非法 / `"bca"` | `ValueError` / `NotImplementedError` | `test_bootstrap_ci_rejects_unknown_method` |
|
|||
|
|
|
|||
|
|
## 9. MCM/ICM 用法
|
|||
|
|
|
|||
|
|
- **哪类题**：任何要"给结论 + 给不确定性"的题 —— 评价类（排名稳定性）、预测类（区间预测）。
|
|||
|
|
- **论文里怎么表述**："采用非参数 bootstrap（`B = 10000`，`seed = 42`）估计中位数的 95% 置信区间。"
|
|||
|
|
- **图**：预测曲线 + 阴影置信带（`viz.plot` / `viz.save_figure`）。
|
|||
|
|
- **踩坑提醒**：必须报 `B` 和 seed；说清用的是 percentile 还是 BCa（目前只有 percentile）；
|
|||
|
|
被问到"为什么不用 ±1.96σ"时就答"统计量是中位数，没有解析标准误"。
|
|||
|
|
|
|||
|
|
## 10. 待拍板（未定稿，Codex 不代选）
|
|||
|
|
|
|||
|
|
| # | 口径 | 现状 | 可选 |
|
|||
|
|
|---|---|---|---|
|
|||
|
|
| 1 | `x` 里的 NaN | 静默传进 `T` → `nan` 区间 | 像 `confidence_interval` 一样丢掉 / 直接报错 / 维持现状 |
|
|||
|
|
| 2 | `n_resamples <= 0` | 走到 `np.percentile` 抛 `IndexError: index -1 is out of bounds` | 入口 `ValueError（必须 ≥ 1）` |
|
|||
|
|
| 3 | `confidence` 越界 | 由 numpy 抛 `ValueError: Percentiles must be in the range [0, 100]` | 入口 `ValueError + 可读文案` |
|
|||
|
|
| 4 | `method="bca"` | docstring 承诺过，实测 `NotImplementedError` | 补 BCa 实现（S 档，约 25 行）/ 从 docstring 删掉 |
|
|||
|
|
| 5 | `bootstrap_statistic` | 仍是 `NotImplementedError`，而 `bootstrap_ci` 把重采样循环内联了 | 抽出去让 `bootstrap_ci` 调它（易错点 16 的变体） |
|
|||
|
|
|
|||
|
|
> 口径 5 尤其值得做：`bootstrap_statistic` 的 docstring 与 `bootstrap_ci` 里那段循环说的是同一件事，
|
|||
|
|
> 现在等于把同一段算法写了两遍，只实现了一遍。
|