采样理论深度基础

1 引言

采样（Sampling）是现代机器学习和统计推断的基石技术。从贝叶斯后验推断到扩散模型生成，采样方法无处不在。本章系统介绍采样理论的深度基础，包括蒙特卡洛方法的理论基础、方差缩减技术、以及粒子方法的现代进展。

---

2 蒙特卡洛方法理论基础

2.1 基本原理

蒙特卡洛方法（Monte Carlo Method）的核心思想是通过随机采样来估计数学期望。¹

问题设定：计算期望 ${\mathbb{E}}_p[f(X)]=\int f(x)p(x)dx$

蒙特卡洛估计：

$${\hat{I}}_N=\frac{1}{N}\sum _{i=1}^{N}f(x^{(i)}),x^{(i)}\sim p(x)$$

2.2 收敛性分析

大数定律

设 $X_1,X_2,\dots$ 为独立同分布随机变量，$\mathbb{E}[X_i]=\mu$，则：

$$\frac{1}{N}\sum _{i=1}^N{X}_i\stackrel{a.s.}{\to }\mu \text{当 }N\to \infty$$

意义：随着样本数增加，蒙特卡洛估计几乎必然收敛到真实值。

中心极限定理

$$\sqrt{N}({\hat{I}}_N-I)\stackrel{d}{\to }\mathcal{N}(0,{\sigma }^2)$$

其中 ${\sigma }^2={\text{Var}}_p(f(X))$ 为方差。

置信区间：

$${\hat{I}}_N\pm z_{\alpha /2}\cdot \frac{\sigma }{\sqrt{N}}$$

2.3 误差分析

蒙特卡洛方法的标准误差为：

$$\text{SE}=\frac{\sigma }{\sqrt{N}}$$

关键洞察：

• 误差与 $N^{-1/2}$ 成正比
• 精度翻倍需要4倍样本量
• 与问题维度无关（这是蒙特卡洛相比数值积分的优势）

2.4 偏差-方差分解

蒙特卡洛估计的均方误差：

$$\text{MSE}({\hat{I}}_N)=\mathbb{E}[({\hat{I}}_N-I{)}^2]={\text{Bias}}^2+\text{Var}({\hat{I}}_N)$$

对于简单蒙特卡洛估计：

• 无偏性：$\text{Bias}=0$
• 方差：$\text{Var}({\hat{I}}_N)={\sigma }^2/N$

---

3 方差缩减技术

方差缩减是提高蒙特卡洛效率的核心技术。² 基本思想是构造控制变量，在保持无偏性的同时降低方差。

3.1 重要性采样（Importance Sampling）

基本原理

从提议分布 $q(x)$ 采样，重加权以修正偏差：

$${\mathbb{E}}_p[f(X)]=\int f(x)p(x)dx=\int f(x)\frac{p(x)}{q(x)}q(x)dx$$

重要性采样估计：

$${\hat{I}}_{IS}=\frac{1}{N}\sum _{i=1}^{N}f(x^{(i)})w(x^{(i)}),x^{(i)}\sim q,w(x)=\frac{p(x)}{q(x)}$$

权重归一化

实践中常使用归一化权重：

$$\tilde{w}(x^{(i)})=\frac{w(x^{(i)})}{{\sum }_{j=1}^{N}w(x^{(j)})}$$ $${\hat{I}}_{IS}=\sum _{i=1}^N\tilde{w}(x^{(i)})f(x^{(i)})$$

方差分析

重要性采样的方差：

$${\text{Var}}_q({\hat{I}}_{IS})={\mathbb{E}}_q\left[{\left(\frac{p(X)}{q(X)}f(X)-I\right)}^2\right]$$

最优提议分布（重要性采样方差最小）：

$$q^{\star }(x)\propto |f(x)|p(x)$$

这正是方差为零的理论最优（若 $f(x)\ge 0$）。

自适应重要性采样

import numpy as np
 
def adaptive_importance_sampling(f, p, q_init, n_iter=10, n_samples=1000):
    """
    自适应重要性采样
    f: 被积函数
    p: 目标分布（未归一化）
    q_init: 初始提议分布
    """
    q = q_init
    
    for iteration in range(n_iter):
        # 从当前提议分布采样
        x = q.sample(n_samples)
        
        # 计算权重
        w = p(x) / q.pdf(x)
        
        # 估计积分
        I_est = np.mean(w * f(x))
        
        # 自适应更新提议分布（基于样本加权）
        # 这里使用简单的均值移动
        weighted_mean = np.sum(w * x) / np.sum(w)
        weighted_std = np.sqrt(np.sum(w * (x - weighted_mean)**2) / np.sum(w))
        
        q = GaussianDistribution(mean=weighted_mean, std=weighted_std)
        
    return I_est

存在的问题

1. 权重退化：当提议分布与目标分布差异大时
2. 方差爆炸：当 $p/q$ 的尾部差异大时

解决方案：多重重要性采样、减方差采样

3.2 拒绝采样（Rejection Sampling）

基本原理

通过拒绝不合符目标的样本，从复杂分布采样。

算法：

1. 从提议分布 $q(x)$ 采样 $x$
2. 计算接受率 $A=\frac{\tilde{p}(x)}{Mq(x)}$
3. 以概率 $A$ 接受 $x$，否则拒绝
4. 重复直到接受

其中 $M\ge {\sup }_x\frac{\tilde{p}(x)}{q(x)}$ 是接受常数。

收敛性

接受率 $\alpha =\frac{1}{M}$，平均需要 $M$ 次提议才能接受一次。

最优策略：选择 $M$ 接近 $\sup \frac{\tilde{p}}{q}$ 以最小化拒绝率。

自适应拒绝采样（Adaptive Rejection Sampling）

当 $\tilde{p}$ 对数凹时，可以使用自适应拒绝采样，自动调整包络函数。

def rejection_sampling(p_tilde, q, M, n_samples):
    """
    拒绝采样
    p_tilde: 未归一化的目标分布
    q: 提议分布
    M: 接受上界
    """
    samples = []
    
    while len(samples) < n_samples:
        x = q.sample()
        u = np.random.uniform(0, 1)
        
        if u < p_tilde(x) / (M * q.pdf(x)):
            samples.append(x)
    
    return np.array(samples)

3.3 分层采样（Stratified Sampling）

基本原理

将积分区域划分为互不相交的层，分别在每层采样：

$$I=\int f(x)p(x)dx=\sum _{k=1}^K{\int }_{{\mathcal{X}}_k}f(x)p(x)dx$$

方差缩减

设每层采样 $n_k$ 个样本，则：

$$\text{Var}\left({\hat{I}}_{SS}\right)=\sum _{k=1}^K\frac{{\sigma }_k^2}{n_k}$$

其中 ${\sigma }_k^2$ 是第 $k$ 层内的方差。

最优分配（Neyman分配）：

$$n_k\propto {\sigma }_k$$

3.4 控制变量法（Control Variates）

基本原理

利用已知期望的控制变量 $g(X)$ 来修正估计：

$${\hat{I}}_{CV}=\frac{1}{N}\sum _{i=1}^N\left(f(X^{(i)})-\lambda (g(X^{(i)})-\mathbb{E}[g(X)])\right)$$

最优系数

最优系数 ${\lambda }^{\star }$ 使方差最小：

$${\lambda }^{\star }=\frac{\text{Cov}(f(X),g(X))}{\text{Var}(g(X))}$$

方差缩减量

$$\text{Var}({\hat{I}}_{CV})=(1-{\rho }_{fg}^2)\cdot \text{Var}(\hat{I})$$

其中 ${\rho }_{fg}$ 是 $f$ 和 $g$ 的相关系数。

3.5 对偶变量法（Antithetic Variates）

基本原理

利用负相关的样本对：

$${\hat{I}}_{AV}=\frac{1}{2N}\sum _{i=1}^N\left(f(X^{(i)})+f(-X^{(i)})\right)$$

当 $f$ 是单调函数时，$X$ 和 $-X$ 负相关，从而降低方差。

3.6 方差缩减技术对比

方法	原理	适用场景	实现难度
重要性采样	改变采样分布	难以直接采样	中等
拒绝采样	接受-拒绝机制	有上界可计算	简单
分层采样	分区域控制	可分层的问题	简单
控制变量	相关性利用	有已知期望的关联变量	中等
对偶变量	负相关对	单调函数	简单

---

4 序贯重要性采样（Sequential Importance Sampling）

4.1 基本框架

序贯重要性采样（SIS）是重要性采样在序列问题上的扩展。³

问题设定：估计序列分布 $p(x_{1:t}\mid y_{1:t})$

递归形式：

$$p(x_{1:t}\mid y_{1:t})\propto p(y_t\mid x_t)\cdot p(x_t\mid x_{t-1})\cdot p(x_{1:t-1}\mid y_{1:t-1})$$

4.2 重要性权重递归

提议分布分解：

$$q(x_{1:t}\mid y_{1:t})=q(x_t\mid x_{t-1},y_t)\cdot q(x_{1:t-1}\mid y_{1:t-1})$$

权重递归：

$$w_t^{(i)}\propto w_{t-1}^{(i)}\cdot \frac{p(y_t\mid x_t^{(i)})\cdot p(x_t^{(i)}\mid x_{t-1}^{(i)})}{q(x_t^{(i)}\mid x_{t-1}^{(i)},y_t)}$$

4.3 粒子滤波器（SIR）

SIS的主要问题是权重退化。粒子滤波器通过重采样解决。

def particle_filter(y_obs, n_particles, p_transition, p_emission, q_proposal):
    """
    粒子滤波器（序贯重要重采样）
    y_obs: 观测序列
    n_particles: 粒子数
    """
    T = len(y_obs)
    particles = []  # 存储所有时刻的粒子
    weights = []
    
    # 初始化
    x = q_proposal.sample(n_particles)
    w = np.ones(n_particles) / n_particles
    particles.append(x)
    weights.append(w)
    
    for t in range(1, T):
        # 重要性采样
        x_new = q_proposal.sample(particles[-1], y_obs[t])
        
        # 计算权重
        w_new = weights[-1] * p_emission(y_obs[t], x_new) * p_transition(x_new, particles[-1])
        w_new /= np.sum(w_new)  # 归一化
        
        # 重采样
        indices = np.random.choice(n_particles, size=n_particles, p=w_new)
        x_resampled = x_new[indices]
        w_resampled = np.ones(n_particles) / n_particles
        
        particles.append(x_resampled)
        weights.append(w_resampled)
    
    return particles, weights

4.4 重采样技术

多重重采样方法：

1. 多项式采样：简单但方差大
2. 分层采样：方差更小
3. 系统采样：推荐使用
4. 残差采样：计算高效

def systematic_resample(weights, n_samples):
    """系统重采样"""
    # 分层
    u0 = np.random.uniform(0, 1/n_samples)
    u = (u0 + np.arange(n_samples)) / n_samples
    
    # 累积权重
    cumsum = np.cumsum(weights)
    
    # 选择
    indices = np.searchsorted(cumsum, u)
    return indices

---

5 蒙特卡洛方法的现代进展

5.1 拉丁超立方体采样

拉丁超立方体采样（Latin Hypercube Sampling, LHS）在高维空间中更均匀地覆盖采样区域。

性质：

• 边缘分布在每维上均匀
• 相关性可控制
• 方差通常比简单随机采样小

def latin_hypercube(n_samples, dim):
    """拉丁超立方体采样"""
    # 每维划分为n_samples个区间
    intervals = np.linspace(0, 1, n_samples + 1)
    
    # 每维随机选择区间
    samples = np.zeros((n_samples, dim))
    for d in range(dim):
        # 打乱区间顺序
        perm = np.random.permutation(n_samples)
        # 在每个区间内均匀采样
        lower = intervals[perm]
        upper = intervals[perm + 1]
        samples[:, d] = np.random.uniform(lower, upper)
    
    return samples

5.2 拟蒙特卡洛方法

拟蒙特卡洛（Quasi-Monte Carlo, QMC）使用确定性序列，具有更均匀的覆盖。

常用序列：

• Halton序列
• Sobol序列
• Niederreiter序列

低差异序列：这些序列的点分布比随机采样更均匀。

5.3 桥式采样（Bridge Sampling）

用于计算两个分布之间的归一化常数比：

$$B=\frac{\int f(x)p(x)dx}{\int f(x)q(x)dx}$$

估计量：

$$\hat{B}=\frac{\frac{1}{N_1}{\sum }_{i=1}^{N_1}{w}_p(x^{(i)})}{\frac{1}{N_2}{\sum }_{j=1}^{N_2}{w}_q(y^{(j)})}$$

其中 $w_p(x)=\frac{f(x)p(x)}{p_0(x)}$, $w_q(y)=\frac{f(x)q(x)}{q_0(x)}$，$p_0,q_0$ 是辅助分布。

---

6 与相关内容的联系

6.1 与变分推断的关系

采样方法和变分推断是近似后验推断的两大范式：

方面	采样方法	变分推断
精度	渐近无偏	有偏但方差小
计算	通常昂贵	可微分、可扩展
表达能力	灵活	受限于近似族

详见 variational-inference-deep-dive。

6.2 与MCMC的联系

SIS/粒子滤波是序贯MCMC的基础。完整的MCMC理论见：

• mcmc-methods — 基础MCMC方法
• mcmc-advanced — NUTS、切片采样

6.3 在贝叶斯深度学习中的应用

采样方法是贝叶斯神经网络的理论基础：

• bayesian-neural-networks — 贝叶斯神经网络基础
• bayesian-last-layer-deep-learning — BNN的最后一层贝叶斯化

---

7 实践指南

7.1 方法选择流程

开始
  │
  ├─ 能否直接采样？ ──是──> 简单蒙特卡洛
  │
  └─ 否
       │
       ├─ 有建议分布 ──是──> 重要性采样
       │                    (检查方差)
       │
       ├─ 可计算上界 ──是──> 拒绝采样
       │
       └─ 其他情况
             │
             ├─ 序列问题 ──> 粒子滤波器
             │
             └─ 静态问题 ──> MCMC方法

7.2 常见问题与解决方案

问题	原因	解决方案
权重退化	提议分布偏离目标	重采样、自适应提议
高拒绝率	上界M太大	改进提议分布
高方差	采样不均匀	方差缩减技术
维度灾难	空间稀疏性	降维、MCMC

7.3 软件实现推荐

• NumPyro：基于JAX的灵活概率编程
• PyMC：成熟的贝叶斯建模工具
• TensorFlow Probability：深度学习集成
• Stan：高效的MCMC实现

---

8 总结

本章系统介绍了采样理论的深度基础：

1. 蒙特卡洛方法的理论基础：大数定律、中心极限定理、收敛性分析
2. 方差缩减技术：重要性采样、拒绝采样、分层采样、控制变量法
3. 序贯重要性采样：粒子滤波器的理论基础和实现
4. 现代进展：QMC、拉丁超立方体、桥式采样

这些技术是现代机器学习中近似推断和生成建模的核心工具。

---

参考文献

Footnotes

1. Robert, C.P. & Casella, G. (2004). Monte Carlo Statistical Methods. Springer. ↩

2. Owen, A.B. (2013). Monte Carlo Theory, Methods and Examples. ↩

3. Doucet, A., Godsill, S., & Andrieu, C. (2000). On sequential Monte Carlo sampling methods for Bayesian filtering. Statistics and Computing, 10(3), 197-208. ↩