6.3 无偏性与偏差
6.2 留下两个数字：maxDepth = 1 时图像比参考解暗 22.28%，打开俄罗斯轮盘赌后亮度又偏了 0.40%。两个都是「不准」，但只有一个是错的。实时混合光追那篇用屏幕空间近似加时间累积把噪声抹平（见 实时光线追踪的混合方案
），代价是它没法告诉你剩下的误差里哪一部分再多跑几帧也消不掉；本节要把这两种误差分开：偏差是系统性的，加样本不会消；方差是随机的，加样本就会降。判据不能只看「和参考解差多少」，必须看「差值随样本数怎么变」。
两种误差的定义
设真值 $I$，估计量 $\langle I \rangle$（随机变量，取决于采样的随机数），$N$ 是样本数。
- 偏差 $b = E[\langle I \rangle] - I$。与 $N$ 无关。$b \ne 0$ 就是有偏，$N \to \infty$ 时估计量收敛到 $I + b$——一个固定的错值。
- 方差 $\mathrm{Var}[\langle I \rangle]$。无偏估计量的误差只来自这里，而它按 $1/N$ 下降，所以标准差按 $1/\sqrt{N}$ 下降。
「有偏」和「方差大」是两个正交的维度。一个估计量可以既无偏又方差大（本节的轮盘赌），也可以方差很小却有偏（本节的硬截断）。图像上看不出区别——两者都是「和参考解不一样」。
双板腔：一个能把两条曲线解析化的构型
6.1 的双板腔（两块无限大平行板，反照率 $a$、$b$，板间均匀发射 $E$）在这里派上用场：它的真值 $C = E/(1-ab)$ 可以算出来，而且两种误差都能写成闭式。
硬截断只统计前 $K$ 次到达，估计量是
$$ \langle C \rangle_K = \sum_{k=0}^{K-1}(ab)^k E, \qquad b_K = \frac{C - \langle C \rangle_K}{C} = (ab)^K $$
偏差与样本数完全无关，只由反照率乘积和截断阶数决定。
俄罗斯轮盘赌在每次到达之后（本节的 rrStart = 2）以概率 $p$ 存活；存活时把权重除以 $p$。单次到达的期望贡献因此不变：
$$ E[\text{到达}] = p\cdot\frac{w}{p} + (1-p)\cdot 0 = w $$
但权重的平方期望涨了：$E[\text{到达}^2] = p\cdot(w/p)^2 = w^2/p$。多弹一次就多乘一个 $1/p$ 因子，于是方差随 $p$ 减小而膨胀。当 $p < (ab)^2$ 时，权重平方的期望发散——期望仍然存在（无偏），方差却不存在。这个阈值是本节最可断言的一条：
$$ p_{crit} = (ab)^2 = 0.444427^2 = 0.197516 $$
#!/usr/bin/env python3
# file: lantern/sim/cg_trace_06_3_rr_walk.py
# tier: S
# verify: /opt/homebrew/bin/python3 sim/cg_trace_06_3_rr_walk.py
"""6.3 最小复现：把「截断」与「轮盘赌」两件事分开量。
双板腔 a = b = 170/255、E = 1，真值 C = E/(1 − ab) = 1.799942。
· 截断到 K 项：估计量少了尾部 Σ_{k≥K}(ab)^k·E，相对偏差 = (ab)^K，**不随样本数下降**。
· 轮盘赌存活 p：每次到达后以概率 p 存活并把权重除以 p。期望不变（无偏），
但方差随 p 减小而膨胀；p < (ab)² = 0.197516 时期望仍存在而方差**发散**。
两组数字都从同一条 SplitMix64 整数流里跑出来。
"""
import os
import sys
sys.path.insert(0, os.path.dirname(os.path.abspath(__file__)))
import cglib # noqa: E402
import cg_trace_06_3_unbiasedness_bias as ub # noqa: E402
ONE = cglib.ONE
WALKS = 200000
PS = [32768, 16384, 8192, 4096]
def truncation(ab, e, exact, kmax=16):
"""截断到 K 项（k = 0..K−1）的估计值与相对偏差。"""
rows = []
acc = 0
thr = ONE
for k in range(kmax):
acc += cglib.fmul(thr, e)
thr = cglib.fmul(thr, ab)
rows.append((k + 1, acc, (exact - acc) / exact))
return rows
def main():
a = b = ub.C8(170)
e = ONE
exact, ab, rr_rows = ub.plate_experiment(a, b, e, PS, WALKS, 20261016)
abf = ab / ONE
print("== 6.3 无偏性与偏差 ==")
print(f"双板腔 a = b = 170/255；ab = {abf:.6f}；真值 C = E/(1−ab) = {exact / ONE:.6f}")
print(f"{'截断 K 项':>9} {'估计值':>12} {'相对偏差':>14} {'(ab)^K':>14}")
for k, est, bias in truncation(ab, e, exact, 16):
if k in (1, 2, 4, 8, 16):
print(f"{k:>9} {est / ONE:>12.8f} {bias:>14.8f} {abf ** k:>14.8f}")
print(f"轮盘赌（每个 p 跑 {WALKS} 次独立游走，同一批种子）：")
print(f"{'存活概率 p':>10} {'经验均值':>12} {'经验标准差':>14} {'最大值':>14}")
for p, mean, std, mx in rr_rows:
print(f"{p / ONE:>10.4f} {mean:>12.6f} {std:>14.6f} {mx:>14.4f}")
print(f"方差发散阈值 (ab)² = {abf * abf:.6f}"
f"（p = 0.125 / 0.0625 已在阈值之下，方差不存在）")
if __name__ == "__main__":
main()
实跑输出：
== 6.3 无偏性与偏差 ==
双板腔 a = b = 170/255；ab = 0.444427；真值 C = E/(1−ab) = 1.799942
截断 K 项 估计值 相对偏差 (ab)^K
1 1.00000000 0.44442655 0.44442749
2 1.44442749 0.19751443 0.19751579
4 1.72970581 0.03902137 0.03901249
8 1.79711914 0.00156831 0.00152197
16 1.79975891 0.00010173 0.00000232
轮盘赌（每个 p 跑 200000 次独立游走，同一批种子）：
存活概率 p 经验均值 经验标准差 最大值
0.5000 1.800495 0.459160 4.6189
0.2500 1.798277 1.166232 180.5572
0.1250 1.791037 2.636715 352.1297
0.0625 1.774285 5.258522 1323.1475
方差发散阈值 (ab)² = 0.197516（p = 0.125 / 0.0625 已在阈值之下，方差不存在）
曲线一：截断偏差有地板
| 截断 $K$ 项 | 估计值 | 相对偏差 | 理论 $(ab)^K$ |
|---|---|---|---|
| 1 | 1.00000000 | 0.44442655 | 0.44442749 |
| 2 | 1.44442749 | 0.19751443 | 0.19751579 |
| 4 | 1.72970581 | 0.03902137 | 0.03901249 |
| 8 | 1.79711914 | 0.00156831 | 0.00152197 |
| 16 | 1.79975891 | 0.00010173 | 0.00000232 |
$K = 1, 2, 4$ 三行的实测与理论差在 0.1% 以内，$(ab)^K$ 这条等式成立。$K = 8$ 分开（实测 0.00156831 对理论 0.00152197），$K = 16$ 分得更开（0.00010173 对 0.00000232，差 44 倍）。
分开的原因和 6.1 里是同一个：Q16.16 的 $1 - ab$ 偏大，真值 117961/65536 本身就带 3.2e-5 的相对误差，而 $(ab)^{16}$ 只有 2.3e-6——理论曲线已经钻到定点算术的分辨率以下了。这条偏差曲线因此不是「无限下降」，而是降到一个地板就停住。
真正要记住的是它不随样本数变化：把每个 $p$ 的游走次数从 20 万加到 200 万，这些数字一位都不会动。这就是「有偏」的含义。
曲线二：轮盘赌标准差无上限
| 存活概率 $p$ | 经验均值 | 经验标准差 | 最大值 |
|---|---|---|---|
| 0.5000 | 1.800495 | 0.459160 | 4.6189 |
| 0.2500 | 1.798277 | 1.166232 | 180.5572 |
| 0.1250 | 1.791037 | 2.636715 | 352.1297 |
| 0.0625 | 1.774285 | 5.258522 | 1323.1475 |
四个 $p$ 的经验均值都在真值 1.799942 附近（最大偏离 1.43%），无偏性成立。但标准差从 0.459160 涨到 5.258522（11.5 倍），最大值从 4.6189 涨到 1323.1475（286 倍）。
最大值这一列的膨胀比标准差更值得看：$p = 0.0625$ 时，200000 次游走里出现了一条权重累积到 1323 的路径。它的贡献是「一次幸运的连续存活」，权重被 $1/0.0625 = 16$ 连乘放大了。这条路径在物理上完全合法——它就是「光子连续弹射了七八次才被吸收」——但它的权重远大于平均值，于是方差被少数几条这样的路径主导。
$p = 0.125$ 与 $p = 0.0625$ 都在阈值 0.197516 之下，它们的理论方差是无穷大。经验标准差仍是有限数，是因为 20 万次游走采不到那条无穷长的尾部；换一批种子，这两个数字会明显跳动，而 $p = 0.5$ 那一行几乎不动。
灯箱：两条误差曲线的形状
双板腔是解析的，但它是确定性的——只有到达次数在变，没有逐像素的噪声。把同一组对照搬到 6.2 的灯箱上（32×24，参考解取 maxDepth=64、RR p=0.5、768 spp，平均辐射亮度 33970.39、每条路径 2.7689 次求交）：
| 配置 | 每像素样本数 | 平均亮度 | 相对参考解 | 求交/路径 |
|---|---|---|---|---|
| 截断 maxDepth=4（RR 关） | 32 | 32554.27 | −4.168713% | 4.4482 |
| 截断 maxDepth=4（RR 关） | 128 | 32606.97 | −4.013558% | 4.4489 |
| 截断 maxDepth=4（RR 关） | 512 | 32527.82 | −4.246552% | 4.4503 |
| 轮盘赌 maxDepth=64（RR p=0.5） | 32 | 34297.80 | +0.963795% | 2.7784 |
| 轮盘赌 maxDepth=64（RR p=0.5） | 128 | 34049.96 | +0.234222% | 2.7703 |
| 轮盘赌 maxDepth=64（RR p=0.5） | 512 | 33927.56 | −0.126100% | 2.7695 |
两条曲线的形状完全不同：
- 截断组是一条水平线。 样本数从 32 加到 512（16 倍），误差在 −4.17%、−4.01%、−4.25% 之间来回，没有下降趋势。这 4% 就是 maxDepth=4 丢掉的那部分间接光。加样本只会让这条水平线更精确地停在 −4% 上。
- 轮盘赌组按 $1/\sqrt{N}$ 下降。 +0.963795% → +0.234222% → −0.126100%：样本数翻两番，误差大致减半（0.96 → 0.23 是 4 倍，0.23 → 0.13 是 1.9 倍），且符号会变——这正是随机波动的特征。样本数够多时它会继续往 0 收。
两者还顺带给出了代价对比：截断组每条路径 4.45 次求交，轮盘赌组 2.77 次。轮盘赌组不仅更准，还更快——因为 maxDepth=64 的硬截断在封闭盒子里要追完所有长路径，而轮盘赌用概率提前砍掉了它们，且不引入偏差。
判据：怎么区分「有偏」和「方差大」
从这两组曲线可以提炼出可操作的判据：
| 观察 | 结论 |
|---|---|
| 误差随样本数下降且符号变化 | 无偏，纯方差，加样本即可 |
| 误差是水平线、符号不变 | 有偏差，加样本无用，必须改算法 |
| 误差随参数（如 maxDepth）单调收缩到 0 | 参数截断型偏差，放宽参数可修 |
| 标准差随参数（如 $p$）爆炸、最大值远超均值 | 方差发散前兆，参数已越过阈值 |
6.2 里那个 −0.397915% 属于第一行：它随样本数下降、符号可正可负，是方差不是偏差。而 maxDepth=4 的 −4% 属于第二行：它不会因为多跑几小时而变小。
P 层：游走序列逐位复现
TypeScript 侧复用同一套 fmul / fdiv / rngQ16，把前 16 次游走、截断表与轮盘赌统计量全部重算一遍：
// file: lantern/core/test/cg_trace_06_3_rr_walk.test.ts
// tier: P
// verify: node --test 'core/test/cg_trace_06_3_rr_walk.test.ts'
// 6.3 的最小 P 复现：同一批整数游走序列，TypeScript 侧必须给出与 Python 侧相同的
// 截断表与轮盘赌统计量。「截断有偏」与「轮盘赌无偏但方差大」是两条独立曲线。
import assert from "node:assert/strict";
import { readFileSync } from "node:fs";
import { dirname, join } from "node:path";
import { test } from "node:test";
import { fileURLToPath } from "node:url";
import { ONE } from "../src/lantern/fixed.ts";
import {
c8, plateExact, plateExperiment, plateWalkSeries, truncatePlate, varianceThreshold,
} from "../src/lantern/cg_trace_06_3_unbiasedness_bias.ts";
const root = join(dirname(fileURLToPath(import.meta.url)), "..", "..");
const doc = JSON.parse(
readFileSync(join(root, "out", "cg_trace_06_3_walks.json"), "utf8"),
) as {
aQ16: number; bQ16: number; eQ16: number; rrStart: number; maxArrivals: number;
seed: number; count: number;
walks: Record<string, { values: number[]; arrivals: number[] }>;
};
const py = JSON.parse(
readFileSync(join(root, "out", "cg_trace_06_3_bias.json"), "utf8"),
) as {
twoPlate: {
varianceThresholdP: number; walksPerP: number;
truncation: { K: number; relBias: number }[];
roulette: { surviveQ16: number; meanQ16: number; stdQ16: number; maxQ16: number }[];
};
};
const a = c8(170);
const b = c8(170);
const e = ONE;
const { ab, exact } = plateExact(a, b, e);
const PS = [32768n, 16384n, 8192n, 4096n];
test("【核心】四个 p 下的前 16 次游走值与到达次数逐字节相同", () => {
assert.equal(Number(a), doc.aQ16);
assert.equal(Number(exact), 117961);
assert.equal(exact, 117961n);
for (const p of PS) {
const got = plateWalkSeries(a, b, e, p, doc.rrStart, doc.maxArrivals, doc.seed, doc.count);
assert.deepEqual(got.values.map(Number), doc.walks[String(p)].values);
assert.deepEqual(got.arrivals, doc.walks[String(p)].arrivals);
}
});
test("【核心】截断偏差逐位一致：K = 1/2/4/8/16 与 (ab)^K 分离", () => {
for (const row of py.twoPlate.truncation) {
const got = truncatePlate(ab, e, exact, row.K);
assert.equal(got.relBias.toFixed(8), row.relBias.toFixed(8));
}
const k1 = truncatePlate(ab, e, exact, 1).relBias;
const k16 = truncatePlate(ab, e, exact, 16).relBias;
assert.equal(k1.toFixed(8), "0.44442655");
// K = 16 时理论 (ab)^16 = 2.32e-6，实测停在 1.02e-4：定点噪声设了地板
assert.equal(k16.toFixed(8), "0.00010173");
assert.ok(k16 > 40 * Math.pow(Number(ab) / 65536, 16));
});
test("【核心】轮盘赌无偏但方差膨胀：均值贴 1.8，标准差涨 11.5 倍", () => {
const got = plateExperiment(a, b, e, PS, py.twoPlate.walksPerP, doc.seed);
for (let i = 0; i < PS.length; i += 1) {
const want = py.twoPlate.roulette[i];
assert.ok(Math.abs(got[i].mean - want.meanQ16) < 1e-3);
assert.ok(Math.abs(got[i].std - want.stdQ16) < 1e-3);
assert.ok(Math.abs(got[i].max - want.maxQ16) < 1e-4);
assert.ok(Math.abs(got[i].mean - 1.799942) < 0.036, "均值始终无偏");
}
assert.ok(got[3].std / got[0].std > 11);
assert.ok(got[3].max / got[0].max > 280);
});
test("方差发散阈值 (ab)² = 0.197516 落在 p = 0.25 与 p = 0.125 之间", () => {
const th = varianceThreshold(ab);
assert.equal(th.toFixed(6), py.twoPlate.varianceThresholdP.toFixed(6));
assert.ok(16384 / 65536 > th);
assert.ok(8192 / 65536 < th);
});
$ node --test 'core/test/cg_trace_06_3_rr_walk.test.ts'
✔ 【核心】四个 p 下的前 16 次游走值与到达次数逐字节相同
✔ 【核心】截断偏差逐位一致：K = 1/2/4/8/16 与 (ab)^K 分离
✔ 【核心】轮盘赌无偏但方差膨胀：均值贴 1.8，标准差涨 11.5 倍
✔ 方差发散阈值 (ab)² = 0.197516 落在 p = 0.25 与 p = 0.125 之间
最硬的一条是游走序列：四个 $p$ 各 16 次游走的值与到达次数都逐项相同（assert.deepEqual 比整数数组）。到达次数也要对，是因为它记录了拒绝采样与轮盘赌各消耗了几次 RNG——序列对而次数不对，说明某处的 rngQ16 调用次序被重排过，后续会全面错位。
本节的位置
到这里，6.2 的两个「不准」被彻底分开了：硬截断有偏、加样本不消；轮盘赌无偏、方差会涨且当 $p < (ab)^2$ 时发散。剩下的问题是怎么把方差压下来——同样的无偏估计量，换一个采样分布就能让方差差出几个数量级。下一章从 pdf 与权重的关系讲起，再给出余弦加权、光源采样与多重重要性采样三种降方差手段。
阅读导航：上一节：6.2 路径采样 · 下一节：7.1 概率密度与权重 。
继续阅读
探索更多技术文章
浏览归档，发现更多关于系统设计、工具链和工程实践的内容。