马尔可夫链蒙特卡洛（MCMC）：为复杂概率分布采样

有些概率分布复杂到您无法直接与之打交道。

在对真实世界数据建模时，数学常常在变得有用之前就已经崩溃。很多时候，积分在纸面上看起来可控，但一旦加入几个潜变量，就立刻变得难以处理。这在贝叶斯推断中尤为常见：后验分布将您的先验与观测数据结合——而结果往往无法用一个简单公式概括。

马尔可夫链蒙特卡洛（MCMC）的基本思想是：与其直接做数学推导，不如通过模拟来探索分布。它抽取能够反映分布形状的样本，而无需完整计算该分布。

本文将介绍 MCMC 背后的核心概念，讲解最常见的算法，并展示如何在 Python 中将其用于实际问题。

您是否需要回顾一下 Python 数学？请阅读我们的 揭秘深度学习中的数学概念 系列博文，看看如何在 Numpy 中应用数学。

什么是马尔可夫链蒙特卡洛？

马尔可夫链蒙特卡洛（MCMC）是一类算法，用于从概率分布中生成样本——即使这些分布复杂到无法直接处理。

这个名称可以拆成两部分。马尔可夫链控制算法如何在可能的状态间移动。每一步只依赖于当前所在位置，而不依赖到达这里的完整历史。蒙特卡洛部分表示您使用随机采样来估计感兴趣的量。

两者结合后，MCMC 会构建一条随机样本链，随着时间推移，反映目标分布的形状。它首先是一种抽样技术。您并未精确求解数学问题，而是通过模拟进行近似。

为何需要马尔可夫链蒙特卡洛

真实世界的数据分布远没有教科书中的那么整洁。

在贝叶斯推断中，您常常需要计算后验分布——即在看到数据后对模型参数更新的概率。公式在纸面上很简单：先验乘以似然，再除以边际似然。最后一项需要对所有可能的参数值进行积分。在高维情况下，这个积分几乎不可能算出来。

随着模型增大，问题只会更糟。增加更多参数或潜变量时，精确计算会走向死胡同。以下常见情境中您都会遇到这种情况：

• 贝叶斯分层模型： 参数彼此依赖，形成嵌套分布，难以整洁地积分
• 高维后验： 参数数量庞大，使得基于网格的近似在计算上不可行
• 非共轭先验： 先验与似然无法组合成可识别的分布

在这些场景下，MCMC 是很好的变通方法。它不去计算分布，而是从分布中抽取样本。这些样本已经包含您所需的一切，而无需真正解出积分。

MCMC 背后的两个思想

MCMC 结合了两个各自简单却结合后威力巨大的思想。让我来依次说明。

马尔可夫链

马尔可夫链是一连串状态，每一步只依赖于您当前所处的位置。

过去到过哪里并不重要。只有当前状态决定下一步去向。这种“无记忆”特性——正式称为马尔可夫性——使数学更易处理，也让算法切实可行。

链条一次一步地穿越状态，在适当设置下，最终会收敛到一个平稳分布——即处于任一状态的概率不再变化的稳定模式。MCMC 正是为处理这种稳定分布而设计。

蒙特卡洛方法

蒙特卡洛方法通过随机采样来估计难以直接计算的量。

思路是对某个分布抽取足够多的随机样本，然后仅根据样本来估计其均值、方差或其他属性。样本越多，估计就越接近真实值。

单独使用时，蒙特卡洛方法要求您能够直接从分布采样——这恰恰是我们要解决的问题。马尔可夫链就负责这一部分。

马尔可夫链蒙特卡洛如何工作

MCMC 是一个循环，每一步都有一个简单决策。

1. 从初始状态开始： 在参数空间中选择一个起点，即使只是一个粗略猜测
2. 提出新状态： 使用提议机制建议下一步移动到哪里
3. 决定是否接受： 将拟议状态与当前状态比较，并根据它们在目标分布下的相对概率来接受或拒绝
4. 重复： 将此循环运行成千上万次，并在此过程中收集样本
5. 使用这些样本： 根据收集的样本估计均值、可信区间或其他感兴趣的量

接受/拒绝这一步就是“魔法”所在。

通过更频繁地接受更好的状态、而较少接受较差的状态，链条会自然聚向高概率区域——而无需完整计算整个分布。

早期样本受起点影响较大，因此会被丢弃。经过足够多次迭代，链会“忘记”起点，剩下的样本就能反映目标分布的真实形状。

目标分布与平稳分布

MCMC 的宗旨是从您无法直接采样的目标分布中生成样本。

目标分布是您想要了解的一切——在贝叶斯推断中通常是后验分布。您知道它的形状，差一个归一化常数，但这个常数无法直接算出。MCMC 并不需要它。

每种 MCMC 算法都被设计为其马尔可夫链的平稳分布与目标分布一致。平稳分布就是链在足够多步后所稳定到的分布。

让链持续运行，它就会开始产生与目标分布抽样无异的样本。积分由此被绕开。

常见的马尔可夫链蒙特卡洛算法

实践中常见几类 MCMC 算法。它们都遵循同一个核心循环，但在如何提出新状态以及如何利用目标分布信息方面有所不同。

Metropolis 算法

Metropolis 是最简单、也是开创先河的 MCMC 算法。

每一步，它通过在当前状态上添加随机噪声来提出新状态。如果拟议状态在目标分布下具有更高概率，则总是接受；若概率更低，则按两者概率之比成比例地接受，否则链保持不动。

这种接收/拒绝机制意味着链会在高概率区域停留更久，而无需计算完整分布。

Metropolis 使用对称的提议分布，即各个方向上的步伐被提出的可能性相同。随着模型的增大，它往往会失效。

Metropolis-Hastings 算法

Metropolis-Hastings（MH）将 Metropolis 推广到允许非对称的提议分布。

MH 会调整接受概率，以考虑某些提议比其他提议更可能。您可以将提议分布调到更贴合目标分布的形状，从而更好地探索并更快收敛。

大多数现代 MCMC 方法都是对 MH 的扩展或基于相同原则构建。因此，只要理解了 Metropolis-Hastings，您就掌握了这一领域的基础。

Gibbs 采样

Gibbs 采样一次只更新一个变量，而不是为所有参数同时提出一个新状态。

每一步，它从各变量的条件分布中采样——即在其他变量当前取值给定的条件下该变量的分布。循环遍历所有变量后，即完成一次完整迭代。

这完全避免了接受/拒绝步骤，因为每个条件抽样都会被接受。当整体联合分布难以采样、但条件分布是可处理的情况下（分层贝叶斯模型中常见），它就很方便。

哈密顿蒙特卡洛

哈密顿蒙特卡洛（HMC）是让现代贝叶斯推断在规模上变得可行的首个算法。

HMC 不再随机提出新状态，而是利用目标分布的梯度信息，提出距离当前位置较远但仍可能被接受的状态。它在参数空间中的移动远胜于随机游走方法。拒绝的提议更少，对高维分布的探索更好。

像 Metropolis 这样的随机游走方法在参数数量增加时难以扩展。HMC 在很大程度上没有这个问题。

HMC 是 Stan（最广泛使用的概率编程平台之一）的引擎。No-U-Turn 采样器（NUTS）是 PyMC 中使用的 HMC 自适应扩展，免除了手动调节步长与步数的需要。

MCMC 在贝叶斯统计中的应用

如果说 MCMC 在哪个领域影响最大，那就是贝叶斯推断。

贝叶斯统计以后验分布为中心，即在看到数据后模型参数的更新概率。计算它意味着先验乘以似然并归一化。归一化这一步需要的积分几乎从不易解。

MCMC 完全移除了这一步。您只需在任意一点评估未归一化的后验，其余交给链来完成。

来看一个简单例子：估计硬币的偏倚。您以“硬币大概率是公平的”为先验，随后观察一系列抛掷。对于简单硬币模型，后验有封闭形式；若加入分层结构，即同时估计上百枚硬币的偏倚，就变得无法计算。

使用 MCMC，运行链、收集来自后验的样本，再用这些样本计算所需量。

理解预热（Burn-In）、收敛与混合（Mixing）

这三个概念常让初学 MCMC 的数据科学家困惑。若您理解有误，可能拿到结果却不知道其不可靠的原因。

预热（Burn-in）

马尔可夫链一开始时，并不知道目标分布的高概率区域在哪里。

这些早期样本受起点影响，而非目标分布。预热就是丢弃它们的做法。您先让链运行若干步，丢掉这些样本，只保留链有足够时间找到良好起点之后的样本。

预热时长没有放之四海而皆准的规则。它取决于您的模型、起点，以及链的混合程度。实践中，通常通过迹线图可视化诊断，而不是事先固定一个数。

收敛

收敛意味着链已不再受起点影响，开始抽取能反映目标分布的样本。

未收敛的链会产生有偏样本。您据此计算的均值不会等于真实后验均值，而是反映链当时卡住的区域。

收敛需要事后通过诊断评估。用不同起点运行多条链并检查它们是否一致，是发现收敛失败最可靠的方法之一。

混合（Mixing）

链即便收敛，但混合不佳仍然成问题。

混合描述链对目标分布的探索程度。良好混合的链能自由移动，既访问高概率区域也访问低概率区域，且产生彼此近似独立的样本。混合不佳的链会在某一区域停留许多步才移动，产生高度相关的样本，无法代表整个分布。

混合不佳常在迹线图上表现为缓慢曲折的“河流”，而不是嘈杂的水平带。看到这种情况，说明您的采样器需要调参——要么改进提议分布，要么直接更换算法。

混合对比图

如何评估 MCMC 结果

接下来我将介绍四种评估 MCMC 的方法，并说明各自适用场景。

迹线图（Trace plots）

迹线图展示每次迭代时参数的采样值。它是运行 MCMC 后您首先应该查看的内容。

健康的迹线图看起来像围绕稳定均值的白噪声。您不应看到趋势、长时间的平直段或缓慢漂移。如果看到链在某一区域徘徊或长时间卡住，那是混合问题，样本不可靠。

迹线图可视化

自相关

MCMC 样本从不完全独立。每个样本都受前一个样本影响。自相关衡量样本跨迭代的相关强度。

高自相关意味着样本所含信息少于其数量所暗示的。两千个相关样本的信息量可能相当于两百个独立样本。大多数 MCMC 库都包含自相关图，方便您查看随样本间距增大，相关性下降的速度。

自相关图可视化

有效样本量

有效样本量（ESS）将自相关转化为一个实用指标：您的链等效于多少个独立样本。

如果抽取了 5,000 个样本，但 ESS 为 200，那么您实际拥有的统计效力等同于 200 次独立抽样。ESS 过低意味着需要更长的链、更好地调参，或两者兼有。多数实践者会在信任估计之前，力求每个参数的 ESS 至少达到数百。

收敛诊断

当您运行多条链时，可以正式检验它们是否收敛到相同分布。Gelman-Rubin 诊断，以R-hat 报告，比较链内方差与链间方差。

R-hat 接近 1.0 表示多条链一致，这是好迹象。超过 1.01 或 1.05（取决于库采用的阈值）表明链尚未收敛，需要更多迭代。PyMC 等现代库会自动计算 R-hat，并在过高时给出警告。

Python 中的马尔可夫链蒙特卡洛

Python 有几种用于 MCMC 的库，各自理念不同。

• PyMC： 对大多数数据科学家来说最易上手。您用 Python 语法定义模型，PyMC 在底层用 NUTS（HMC 的自适应版本）进行采样。它是 Python 中进行贝叶斯建模的首选。
• CmdStanPy / PyStan： 指向 Stan 的 Python 接口，Stan 是专用的概率编程语言。Stan 会在运行前将模型编译为 C++，因此速度快，适合复杂分层模型。代价是学习曲线更陡峭。
• NumPy： 没有专门的 MCMC 库，但您可以从零实现诸如 Metropolis-Hastings 的算法。先理解底层机制，再转向更高层工具，这样做很有帮助。

在大多数实际工作中，PyMC 往往是起点。本文也将使用它，因此如果您跟着操作，请先安装该库：

pip install pymc

为保持简单，我们用一个容易的例子：根据一系列抛掷结果来估计硬币的偏倚。

步骤 1：定义模型

import pymc as pm
import numpy as np

# 1 = heads, 0 = tails
observed_flips = np.array([1, 0, 1, 1, 0, 1, 1, 1, 0, 1])

with pm.Model() as coin_model:
    # Prior: we believe the coin is probably fair
    bias = pm.Beta("bias", alpha=2, beta=2)

    # Likelihood: observed flips given the bias
    flips = pm.Bernoulli("flips", p=bias, observed=observed_flips)

pm.Beta 先验表达了我们对硬币公平的弱信念。pm.Bernoulli 似然将模型与观测数据连接起来。

步骤 2：运行采样器

with coin_model:
    trace = pm.sample(2000, tune=1000, return_inferencedata=True)

运行采样器的输出

tune 设定预热步数——这些样本会被丢弃。sample 在调参后，每条链抽取 2000 个后验样本。

步骤 3：检查后验

import arviz as az

az.plot_trace(trace, var_names=["bias"])
az.summary(trace, var_names=["bias"])

模型迹线图与汇总结果

az.summary() 会给出每个参数的后验均值、标准差和 R-hat。若 R-hat 接近 1.0，说明链已收敛。az.plot_trace() 为每个参数并排绘制迹线与后验分布。

对于本数据集——10 次抛掷中 7 次正面——后验均值为 0.642，标准差为 0.124。这既反映了数据中的证据，也与公平硬币的先验保持接近。R-hat 为 1.00，ESS 远超 2000，说明链已收敛且样本可靠。

MCMC 常见错误

MCMC 很容易上手，但也很容易用错。以下是最常见的错误。

• 过早认为已收敛：采样器没有报错并不代表已收敛。在信任结果前务必检查 R-hat 和迹线图。链看起来稳定也可能只是卡在参数空间的某个区域。
• 忽视预热：早期样本偏向起点而非目标分布，需要丢弃。多数库通过 tune 参数自动处理，但要确认没有把这些样本误算进分析。
• 混合不佳： 移动过慢的链会产生高度相关的样本，其信息量低于数量所暗示。如果迹线图像缓慢漂移而非嘈杂水平带，需要调优提议分布或更换算法。从 Metropolis 切换到 NUTS 往往能立刻解决。
• 样本过少： 有效样本量（ESS）过低意味着估计不可靠，即便链“技术上”已收敛。若 ESS 仅为总抽样数的一小部分，请延长链或先改善混合再下结论。

MCMC 与其他采样方法对比

MCMC 不是近似分布的唯一方法。下面是它与替代方案的对比。

• 拒绝采样 从更简单的分布提出样本，再根据与目标匹配程度来接受。可行时它是精确的，但在高维中接受率会崩塌。
• 重要性采样 对来自提议分布的样本加权，以近似目标分布下的期望。在低维中速度快且效果好，但权重方差会随维度上升而爆炸。实践中，少数样本往往承担几乎全部权重，导致估计不可靠。
• 变分推断（VI） 拟合更简单的分布以近似目标分布，通过最小化某种散度来实现。VI 远比 MCMC 快，能扩展到大数据集，但会引入 MCMC 所避免的近似误差。拟合的分布可能错过 MCMC 能捕捉到的后验结构。

简而言之：

MCMC 与替代方法对比

当准确性比速度更重要时，MCMC 是正确的选择。若您需要扩展到大型数据集或进行实时推断，变分推断或许值得用准确性换取速度。

结语

MCMC 这种工具外表看上去吓人，但一旦理解其真实在做什么——构建一条逐步反映某个无法直接计算的分布形状的样本链——就很容易理解。

将其拆分为两部分来理解也更容易：马尔可夫链与蒙特卡洛方法。

它在贝叶斯统计中的作用不言而喻。否则无法触及的后验分布，一旦有了可靠的采样器便能迎刃而解。这就是为何 MCMC 位于 PyMC 和 Stan 等概率编程库的核心。

但在动手实现之前，您应先把直觉理顺。理解为何链需要预热、混合到底意味着什么，以及如何阅读迹线图。代码本身是最容易的部分，因为 Python 库将所有抽象都封装在简单的函数调用之后。

如果您想熟练掌握机器学习，请报名我们的 Machine Learning Scientist in Python 学习路径。85 小时的材料将助您在 2026 年为求职做好准备。