Markov Chain Monte Carlo（1）

最近因为种种原因不得已去钻研一下一些基础的机器学习的知识。当初统计和矩阵论的知识都忘得差不多了。好在我爱自然语言处理这里的文章写得都比较清楚，剩下不懂的自己写写算算和朋友讨论讨论，也有一些收获。大部分都是转载，附加了一些学习过程中自己的笔记（免得将来忘了又得从头思考）

随机模拟(或者统计模拟)方法有一个很酷的别名是蒙特卡罗方法(Monte Carlo Simulation)。统计模拟中有一个重要的问题就是给定一个概率分布p(x)，我们如何在计算机中生成它的样本。一般而言均匀分布Uniform(0,1)的样本是相对容易生成的。 通过线性同余发生器可以生成伪随机数，我们用确定性算法生成[0,1]之间的伪随机数序列后，这些序列的各种统计指标和均匀分布Uniform(0,1) 的理论计算结果非常接近。这样的伪随机序列就有比较好的统计性质，可以被当成真实的随机数使用。

而我们常见的概率分布，无论是连续的还是离散的分布，都可以基于Uniform(0,1)的样本生成。例如正态分布可以通过著名的 Box-Muller 变换得到。不过我们并不是总是这么幸运的，当p(x)的形式很复杂，或者 p(x)是个高维的分布的时候，样本的生成就可能很困难了。此时就需要使用一些更加复杂的随机模拟的方法来生成样本。而本节中将要重点介绍的 MCMC(Markov Chain Monte Carlo) 和 Gibbs Sampling算法就是最常用的一种，这两个方法在现代贝叶斯分析中被广泛使用。要了解这两个算法，我们首先要对马氏链的平稳分布的性质有基本的认识。

马氏链及其平稳分布

马氏链的数学定义很简单，也就是状态转移的概率只依赖于前一个状态。

P(Xt+1=x|Xt,Xt−1,⋯)=P(Xt+1=x|Xt)

假设有初始分布π0=(π1,π2,...πm)和状态转移矩阵P=(p11p21⋮pm1p12p22⋮pm2⋯⋯⋱⋯p1mp2m⋮pmm) ，那么有

马氏链定理： 如果一个非周期马氏链具有转移概率矩阵P,且它的任何两个状态是连通的，那么 limn→∞Pnij存在且与i无关（这里Pnij表示矩阵P的n次幂的i行j列的元素），则limn→∞Pnij=π(j), 并且我们有：

πP=π，其中π=limn→∞π0Pn且π与π0无关；记π=[π(0),π(1),⋯π(m)]，π称为马氏链的平稳分布。

π是方程πP=π的唯一非负解。

π(j)=∑mi=0π(i)Pij

limn→∞P=⎡⎣⎢⎢⎢⎢⎢π(0)π(0)⋮π(0)π(1)π(1)⋮π(1)⋯⋯⋱⋯π(m)π(m)⋮π(m)⎤⎦⎥⎥⎥⎥⎥

这个马氏链的收敛定理非常重要，所有的 MCMC(Markov Chain Monte Carlo) 方法都是以这个定理作为理论基础的。 定理的证明相对复杂，一般的随机过程课本中也不给证明，所以我们就不用纠结它的证明了，直接用这个定理的结论就好了。就定义3做一些说明。

由矩阵乘法，πn(j)=∑ni=0πn−1(i)pij，对等式两侧n取极限有：

limn→∞πn(j)=∑ni=0limn→∞πn−1(i)pij

由于limn→∞πn(j)=π(j)，故：π(j)=∑mi=0π(i)Pij