马尔可夫链蒙特卡洛（MCMC）中的Metropolis-Hastings算法原理与采样过程

题目描述
Metropolis-Hastings（MH）算法是马尔可夫链蒙特卡洛（MCMC）方法的核心算法之一，用于从复杂概率分布中抽取随机样本。假设我们需要从一个目标分布 \(\pi(x)\) 中采样，但其形式复杂（如未归一化或高维），直接采样困难。MH算法通过构建一个马尔可夫链，使其平稳分布收敛到目标分布 \(\pi(x)\)，从而通过迭代生成样本。本题目将详细讲解MH算法的原理、接受概率计算及采样步骤。

解题过程

1. 问题定义与核心思想
   • 目标：从目标分布 \(\pi(x)\) 中采样，已知 \(\pi(x)\) 可能未归一化（即仅知 \(\tilde{\pi}(x) \propto \pi(x)\)）。
   • 核心思想：构造一条马尔可夫链，其状态转移概率满足细致平衡条件（Detailed Balance）：

\[ \pi(x) T(x' \mid x) = \pi(x') T(x \mid x') \]

 其中 $ T(x' \mid x) $ 为转移概率。MH算法通过设计一个提议分布（Proposal Distribution） $ q(x' \mid x) $ 和接受概率 $ A(x' \mid x) $ 来修正转移，使得修正后的转移核 $ T(x' \mid x) = q(x' \mid x) A(x' \mid x) $ 满足细致平衡条件。

2. 接受概率的推导
   • 细致平衡条件要求：

\[ \pi(x) q(x' \mid x) A(x' \mid x) = \pi(x') q(x \mid x') A(x \mid x') \]

• 为满足该条件，定义接受概率为：

\[ A(x' \mid x) = \min \left(1, \frac{\pi(x') q(x \mid x')}{\pi(x) q(x' \mid x)} \right) \]

 此设计确保：  
 - 若 $ \frac{\pi(x') q(x \mid x')}{\pi(x) q(x' \mid x)} \geq 1 $，则接受新状态 $ x' $（即 $ A(x' \mid x) = 1 $）；  
 - 否则以概率 $ \frac{\pi(x') q(x \mid x')}{\pi(x) q(x' \mid x)} $ 接受 $ x' $。

3. 算法步骤详解
   • 输入：目标分布 \(\pi(x)\)（可未归一化），提议分布 \(q(x' \mid x)\)，初始状态 \(x_0\)，采样次数 \(N\)。
   • 过程：
     1. 初始化 \(t = 0\)，当前状态 \(x^{(0)} = x_0\)。
     2. 对于 \(t = 1\) 到 \(N\)：
        a. 从提议分布 \(q(x' \mid x^{(t-1)})\) 生成候选状态 \(x'\)。
        b. 计算接受概率：

\[ \alpha = \min \left(1, \frac{\tilde{\pi}(x') q(x^{(t-1)} \mid x')}{\tilde{\pi}(x^{(t-1)}) q(x' \mid x^{(t-1)})} \right) \]

       其中 $ \tilde{\pi}(x) $ 为未归一化的目标分布（因归一化因子在比值中抵消）。  
    c. 从均匀分布 $ U(0,1) $ 生成随机数 $ u $。  
    d. 若 $ u \leq \alpha $，接受候选状态：$ x^{(t)} = x' $；否则拒绝：$ x^{(t)} = x^{(t-1)} $。  
 3. 输出样本序列 $ \{x^{(1)}, x^{(2)}, \dots, x^{(N)}\} $。  

4. 关键点说明

   • 提议分布的选择：常用对称分布（如高斯分布），此时 \(q(x' \mid x) = q(x \mid x')\)，接受概率简化为 \(\alpha = \min \left(1, \frac{\tilde{\pi}(x')}{\tilde{\pi}(x)} \right)\)（即Metropolis算法）。
   • 收敛性：生成的马尔可夫链需不可约、非周期且正常返，才能保证平稳分布为 \(\pi(x)\)。
   • 样本相关性：连续样本间可能存在自相关性，实际使用时需丢弃前期样本（Burn-in）或稀释采样（Thinning）。

5. 示例说明
   假设目标分布为 \(\pi(x) \propto e^{-x^2/2}\)（标准正态分布），提议分布为 \(q(x' \mid x) = \mathcal{N}(x' \mid x, 1)\)。

   • 从 \(x^{(0)} = 0\) 开始，生成候选 \(x' \sim \mathcal{N}(x^{(0)}, 1)\)。
   • 计算 \(\alpha = \min \left(1, \frac{e^{-(x')^2/2}}{e^{-(x^{(0)})^2/2}} \right)\)。
   • 以概率 \(\alpha\) 接受 \(x'\)，否则保留原状态。重复过程得到近似服从 \(\mathcal{N}(0,1)\) 的样本。

通过以上步骤，MH算法实现了对复杂分布的高效采样，仅需计算目标分布的比值，避免了归一化常数的求解。