[科普中国]-马尔科夫链蒙特卡洛方法

从理论上说，贝叶斯推断和分析是容易实施的，即对于任何先验分布，只需要计算所需后验分布的性质，如后验分布的矩(如后验均值、后验方差)、后验概率密度函数等，而这些计算木质上就是计算后验分布某一函数的高维积分。但在实践中，鉴于未知参数的后验分布多为高维、复杂的非常见分布，对这些高维积分进行计算十分困难，这一困难使得贝叶斯推断方法在实践中的应用受到很大的限制，在很长一段时间，贝叶斯推断主要用于处理简单低维的问题，以避免计算上的困难。MCMC方法突破了这一原木极为困难的计算问题，它通过模拟的方式对高维积分进行计算，进而使原本异常复杂的高维积分计算问题迎刃而解，使贝叶斯方法仅适用于解决简单低维问题的状况大有改观为贝叶斯方法的应用开辟了新的道路。1

MCMC方法的结构如图1所示。2

MCMC方法是使用马尔科夫链的蒙特卡洛积分，其基木思想是：构造一条Markov链，使其平稳分布为待估参数的后验分布，通过这条马尔科夫链产生后验分布的样木，并基于马尔科夫链达到平稳分布时的样木(有效样本)进行蒙特卡洛积分。

设 为某一空间，n为产生的总样木数，m为链条达到平稳时的样本数，则MCMC方法的基木思路可概括为：

(1)构造Markov链。构造一条Markov链，使其收敛到平稳分布 ；

(2)产生样本：由 中的某一点 出发，用(1)中的Markov链进行抽样模拟，产生点序列： ；

(3)蒙特卡洛积分：任一函数f(x)的期望估计为： 。1

MCMC中采样算法如图2所示。

在采用MCMC方法时，马尔科夫链转移核的构造至关重要，不同的转移核构造方法，将产生不同的MCMC方法，目前常用的MCMC方法主要有两种：Gibbs抽样和Metropolis-Hastings算法。1

Metropolis 算法是马尔科夫链蒙特卡罗的基石。它是由 Metropolos 等人在 1953年的一篇仅 4 页的文章中提出。Metropolis 算法用一个对称的建议分布 T(x,y)来产生一个潜在的转移点，然后根据特定的接受拒绝方法来决定是否转移到该潜在点。最初的 Metropolis 算法将 T 取为对称的函数，而 Metropolis-Hasting 方法将之推广到非对称的 T，其中 T 满足 T(x,y)>0 的充要条件是 T(y,x)>0。已知当前状态 ，该算法的迭代步骤如下：

第 1 步：从建议分布 抽取潜在转移点 y。

第 2 步：在[0,1]均匀分布中抽取一个随机数 ，并且更新 ，若满足 ，否则令 。对于 r(x,y)，Metropolis 和 Hastings 建议采用如下形式的：

Gibbs 算法是一种特殊的 Metropolis 算法。对于一个 d 维随机变量 ，Gibbs 算法在生成 的第 i 个分量时选择建议分布 T 为第 i 个分量基于其他所有分量的条件分布。即已知当前状态 ， Gibbs 算法通过如下的方法更新 d 个分量来更新 ：

第 1 步 ： 对 于 固 定 的 i（初值为1），从条件分件=(|)中抽样出。

第 2 步：令 i=i+1，重复第 1 步。3

MCMC方法依赖于模拟的收敛性，下面介绍三种常用的收敛性诊断方法。

这种方法的思路是选取多个不同的初值，同时产生多条马尔科夫链，经过一段时间后，若这几条链稳定下来，则说明算法收敛了。在实际操作中，可在同一个二维图中画出这些不同的马尔科夫链产生的后验样本值对迭代次数的散点图，如果经过若干次迭代后，这些散点图基木稳定，重合在一起，则可判断其算法收敛。

这种方法的思路是选取多个不同的初值，同时产生T条马尔科夫链，记第j条链的方差估计为，链内方差的均值为W，链间方差为B，则：。R的值趋近1，则表明MCMC模拟收敛，且比较稳定，通常R