计算机图形学 20:随机样本怎样估计积分
第 19 篇把半球划成规则小格,累加入射光。若场景里多一次反射,积分又增加一组方向;继续铺规则网格,组合数量会迅速增长。蒙特卡洛方法改用随机样本估计积分,但随机位置不是任意挑选的:样本来自哪里、以什么概率出现,决定了每个样本应有多大权重。
先用一个答案已知的问题检查这个关系:在 [0,1] 上积分 x²,答案为 1/3。同一份程序比较均匀采样、重要性采样,以及故意遗漏半个区间的错误策略。第三种结果波动很小,却始终离正确答案很远。
概率密度不是一个点的概率
连续随机变量 X 的概率密度为 p(x) 时,落进区间 [a,b] 的概率是该区间内 p 的积分。单个点通常具有零概率;不能把 p(x)=2 解释为某点出现的概率是 200%。例如 [0,1/2] 上恒定密度 2,在其他位置为零,总概率仍为 1。
这也与第 03 篇的规则采样不同。中心点序列由网格确定,没有为每个位置分配抽样概率;把确定性求积中的步长随意改名为 PDF,并不能直接得到同一组统计结论。
设目标积分为 I=∫f(x)dx,随机样本按 p 产生。对于 p 非零的位置,把样本贡献写成 Y=f(X)/p(X),就有:
关键条件是:f 非零的区域必须被 p 覆盖,允许忽略零测集;还需要 f 可积。如果抽样概率为零的整个区间仍有贡献,约分后的积分范围就已经变了。除以 PDF 只能校正出现频率,无法补回永远不会出现的区域。
本篇 x 是无量纲变量,因此密度单位没有额外长度因素。到方向积分时,PDF 相对于立体角 dω 定义,单位为 sr⁻¹。把相对于面积的 PDF 直接塞进方向积分,会漏掉变量变换的几何因子。
多次采样为什么能降低波动
用 N 个独立同分布样本的贡献求平均:
每项期望为 I,所以平均的期望仍为 I。这叫无偏,不是说每次计算都恰好等于 I。一次结果可以在答案之上,也可以在答案之下。
若单样本贡献方差有限,独立性使协方差项消失:
标准差因而按 1/√N 下降。将典型误差缩小一半,通常需要四倍样本,而不是两倍。这里讨论的是概率模型下的波动尺度,不保证某次运行每增加一个样本就更准确;高维也不自动改变这个指数,但会影响单样本方差与计算成本。
相关样本需要保留协方差分析。复制同一个随机样本 N 次,平均值仍是原来的一个样本,不会因此获得 1/N 的方差。伪随机生成器是有限状态的确定程序,使用不同种子是复现实验的安排,不能当作数学上的独立性证明。
怎样选择更合适的抽样分布
对 f(x)=x²,均匀密度 p=1 产生贡献 Y=X²。E[Y²]=∫x⁴dx=1/5,所以单样本方差为 1/5−1/9=4/45。
较大的 x 对积分贡献更大,可以提高这些位置的采样密度。选择 p(x)=2x,其累积分布 F(x)=x²。若 U 在 [0,1] 均匀分布,令 X=√U,就得到这一密度。贡献变为 X²/(2X)=X/2;零点按连续极限处理。
新贡献的二阶矩为 ∫(x/2)²·2x dx=1/8,单样本方差为 1/8−1/9=1/72。相同 N 下的理论方差比为 (4/45)/(1/72)=6.4。比较对象是同样数量的样本,不是相同执行时间;计算平方根也有成本,本篇没有计时数据。
如果忘记除以 2X,只平均按新分布生成的 X²,其期望会变成 ∫x²·2x dx=1/2。更密集地采到大值,同时还给每个样本原来的权重,就改变了要估计的量。
重要性采样并不保证任何新 PDF 都降低方差。在贡献较大的区域把密度设得极小,会产生很大的 f/p 权重。理论上对非负 f 选择 p=f/I 可使每次贡献等于 I,方差为零;但构造归一化分布通常需要先知道 I,不能把这个形式直接当作通用求解器。
方差很小,为什么仍然算错
故意只在 [0,1/2] 上采样,令 X=U/2,该区间内 p=2。即使老老实实计算 Y=X²/2,期望也只有:
这套采样漏掉了 (1/2,1]。它对半区间积分是正确估计,对原来的完整积分则有 −7/24 的偏差。增加样本不会改变支持范围,只会让结果更集中地靠近 1/24。
单样本方差为 1/720,比前两种都小。若只挑一张看起来平滑的结果图,错误策略反而可能显得更可靠。判断估计质量要同时考虑偏差与方差:MSE=方差+偏差²。本例针对完整积分的均方误差最终趋向 49/576≈0.08507,无法趋零。
不要把结果乘 8 当作通用修复。这个倍率来自 x² 的特殊积分关系,换成任意函数后,未采区域仍然没有证据。正确修复是恢复所需的采样支持,或用另一个覆盖缺失区域的合法估计器组合起来。
256 个种子实际得到了什么
C++17 程序使用 std::mt19937,种子为 1 到 256。每个种子分别运行 N=64、256、1024、4096,固定样本数结束,不按当前误差提前停止。同一个种子下三种策略共享 U,便于配对比较;每个 N 重置生成器,因此较短运行是较长运行的前缀,不能把这四档当作互相独立的试验。
U 按 (double(rng())+0.5)/4294967296.0 产生,落在开区间 (0,1),避免端点形成 0/0。它只有有限个离散取值,是连续均匀变量的数值近似;前面的连续无偏推导并不意味着有限精度程序具有完全相同的实数分布。
对每个 N,将 256 个完整运行的估计值再统计均值与样本方差,后者分母为 255。这里的方差描述一次 N 样本估计在重复运行间的变化,不是把一次运行内部 N 个贡献的方差直接贴到表中。
| N | 均匀采样均值 | 均匀方差 | 重要性均值 | 重要性方差 | 漏域均值 |
|---|---|---|---|---|---|
| 64 | 0.33291608 | 0.00145868 | 0.33292752 | 0.000202289 | 0.04161451 |
| 256 | 0.33406302 | 0.000369125 | 0.33331751 | 0.0000579227 | 0.04175788 |
| 1024 | 0.33407419 | 0.0000962605 | 0.33340335 | 0.0000152789 | 0.04175927 |
| 4096 | 0.33371895 | 0.0000230530 | 0.33343119 | 0.00000364153 | 0.04171487 |
4096 样本下,观测方差比约为 6.33,接近理论的 6.4。有限次重复的方差本身也会波动;四档均匀采样的实测/理论方差比为 1.050、1.063、1.109、1.062,重要性采样为 0.932、1.068、1.126、1.074,并非精确相等。
漏域策略在 N=4096 时的方差只有约 3.60×10⁻⁷,但针对完整积分的 MSE 约为 0.0850417。重要性策略的 MSE 约为 3.64×10⁻⁶。低波动与低误差在这个反例中清楚分开。
程序保留全部 3072 次估计和12 行汇总。没有选取最接近解析值的种子,也没有删掉表现较差的运行。Moments 用在线更新计算均值与平方偏差和,另用序列 1、2、3 检查均值 2、样本方差 1。
回归检查要求均值距离各自解析期望不超过五个理论标准误,方差比在预设的 0.65 到 1.35 内。这是当前固定种子组的宽容差诊断,不是用一个通过标记证明所有随机运行都会成功。尤其漏域策略相对于 1/24 的检查通过,绝不能写成它通过了完整积分的正确性验收。
怎样接回光照积分
第 19 篇的反射积分中,把 fᵣLᵢcosθ 看作被积函数,方向按 p(ω) 抽样,每次贡献就是 fᵣLᵢcosθ/p(ω)。Lambert 的常量 BRDF 不意味着整项恒定,因为光照和可见性可能随方向变化。
如果正面半球光照恒定,按 cosθ/π 采样能够抵消余弦因子,得到常量贡献。这与本篇匹配函数形状的想法相同。不过复杂场景中的遮挡、光源和多次散射会留下波动,不能据这个简化例子许诺路径追踪没有噪声。
本篇只验证一维积分,没有生成间接光场景。下一篇将把抽样方向接到最近交点查询,并单独讨论路径终止;求交正确、抽样正确与终止无偏,需要各自的证据。
复跑与练习
1 | |
输出目录需已存在。monte_carlo.hpp 提供开区间伪随机值与在线矩统计,monte_carlo_check.cpp 实现三个积分策略、原始输出与数值检查。完整日志在 writing-plans/computer-graphics/evidence/20-cpu.txt。
练习一:若保持分布不变,把 N 从 256 增到 4096,理论方差与标准差各缩小多少?是否保证某个种子的绝对误差也按这个比例减少?
答案:样本数扩大 16 倍,方差缩小 16 倍,标准差缩小 4 倍。不保证单次绝对误差的比例,甚至不保证每次单调下降。
练习二:按 p=2x 生成 X,却用 X² 直接平均,极限是多少?缺了哪一步?
答案:期望为 ∫₀¹2x³dx=1/2。缺少除以采样 PDF 的权重。应平均 X²/(2X),而不是仅改变样本位置。
练习三:漏域策略的重复结果越来越集中,能否据此声称“已收敛到完整积分”?
答案:不能。它的期望为 1/24,目标为 1/3,偏差不随 N 消失。应先核查支持范围,再用完整目标衡量误差;只看跨次方差会遗漏系统错误。
系列导航与资料
前篇:19:表面颜色从哪里来。系列入口。下一篇:21:光怎样经过多次反射。
PBRT 第四版 §2.1.3–2.1.4 的式 2.7、2.10、2.11分别用于核对估计器、独立变量方差与样本方差;§2.2.2核对重要性采样的收益和限制。实际核验日期为 2026-09-20,研究记录见 research-20.md;三种 x² 策略的解析推导与数据为本系列有界实验。






