Jupyter AI

21 Metropolis-Hastings算法

📅 发表日期: 2024年8月15日

分类: 📊贝叶斯学习入门

👁️阅读: --

在马尔可夫链蒙特卡洛(MCMC)方法中,Metropolis-Hastings算法是一个重要的产生采样方法,广泛应用于贝叶斯学习和统计推断的领域。它是一种构建马尔可夫链的方法,通过从一个指定的概率分布中生成样本,来近似该分布。接下来,我们将详细探讨该算法的原理、步骤和应用案例。

算法原理

Metropolis-Hastings算法Gibbs采样的推广,适用于更为复杂的联合分布。该算法的主要思想是通过在目标分布和建议分布之间建立一种接受-拒绝机制,来有效生成所需样本。

目标分布

假设我们要从一个复杂的概率分布 p(x)p(x) 中采样。通常,目标分布为后验分布,通常是难以直接采样的。

提议分布

为了生成样本,我们引入一个可行的提议分布 q(xx)q(x'|x),其用于从当前位置 xx 中生成一个新样本 xx'。典型的提议分布可以是对称的,如正态分布。

接受率

在生成新样本后,我们计算接受率 αα,决定是否接受新样本。接受率的计算公式为:

α=min(1,p(x)q(xx)p(x)q(xx))\alpha = \min\left(1, \frac{p(x')q(x|x')}{p(x)q(x'|x)}\right)
  • α1α \geq 1,则无条件接受新样本。
  • 否则,以概率 αα 接受新样本。

算法步骤

  1. 选择初始点 x0x_0
  2. 迭代 nn 次:
    • 从提议分布 q(xx)q(x'|x) 中生成新样本 xx'
    • 计算接受率 αα
    • 生成一个均匀随机数 uUniform(0,1)u \sim \text{Uniform}(0, 1)
    • 如果 u<αu < α,则接受 xx', 令 xt+1=xx_{t+1} = x';否则,保留旧样本,令 xt+1=xx_{t+1} = x

算法示例

让我们通过一个具体的案例,使用Metropolis-Hastings算法从一个目标分布中采样。

目标分布示例

假设我们的目标是从一个一维的正态分布中采样,设其均值为 μ=0\mu=0,标准差为 σ=1\sigma=1。我们设定提议分布为当前点附近的正态分布,例如 q(xx)N(x,1)q(x'|x) \sim \mathcal{N}(x, 1)

Python实现

下面是使用Python实现Metropolis-Hastings算法的一个示例代码:

import numpy as np
import matplotlib.pyplot as plt

# 目标分布: 标准正态分布
def target_distribution(x):
    return np.exp(-0.5 * x**2) / np.sqrt(2 * np.pi)

# 提议分布: 正态
def proposal_distribution(x):
    return np.random.normal(x, 1)

# Metropolis-Hastings算法
def metropolis_hastings(iterations):
    samples = []
    x = 0  # 初始化
    for _ in range(iterations):
        x_new = proposal_distribution(x)
        acceptance_ratio = target_distribution(x_new) / target_distribution(x)
        
        if np.random.rand() < acceptance_ratio:
            x = x_new
        
        samples.append(x)
    return np.array(samples)

# 生成样本
n_samples = 10000
samples = metropolis_hastings(n_samples)

# 绘图
plt.figure(figsize=(12, 6))
plt.hist(samples, bins=30, density=True, alpha=0.6, color='g')

# 绘制目标分布
x = np.linspace(-4, 4, 1000)
plt.plot(x, target_distribution(x), 'r', lw=2)
plt.title('Metropolis-Hastings Sampling')
plt.xlabel('x')
plt.ylabel('Density')
plt.legend(['Target Distribution', 'Samples'])
plt.show()

这段代码首先定义了目标分布和提议分布。然后实施了Metropolis-Hastings算法,并对生成的样本进行了可视化。你可以看到,样本的分布逐渐趋近于目标分布。

结论

Metropolis-Hastings算法是MCMC的一个强大工具,能够在无法直接采样的情况下,从复杂的后验分布中生成样本。本节中,我们不仅讨论了算法的原理及其步骤,还提供了实际代码的示例,帮助你在实际应用中理解和实现这一算法。

在下一篇中,我们将探讨贝叶斯学习在实际中的应用案例,深入理解这一理论如何通过实际问题得以体现。