随机模拟:马尔科夫链蒙特卡洛采样MCMC与EM算法「2.3」

最近学习了机器学习中的马尔科夫链蒙特卡洛(Markov Chain Monte Carlo, 简称MCMC) 相关的知识。

主要内容包括：

【1】蒙特卡洛原则，及其应用于采样的必要性（已经发布在头条）

【2】用于求解最大似然、近似推断、期望问题的经典采样算法：Metropolis-Hastings，Rejection，Importan，Metropolis和Gibbs算法。（本文属于此部分）

【3】马尔可夫链各个性质在蒙特卡洛采样问题中的应用，包括同质性，平移不变性

—————【2】—————

上一篇【2.2】中详细讨论了EM优化算法的推导和性质，EM算法通过不断提高下界来逼近最大似然，其中E步求Q(z)=p(z|x,θ)，此时的θ是上一个M步已经固定的，求得此θ对应的最大下界的Q(z)，更新Q(z)。M步固定Q(z)，求得此时使得下界最大的θ，更新θ，如此迭代，直到收敛到局部最优，得到最终估计值。

M步中对L最大化求新θ，等价于求使得Q(θ，θold)最大的θ值。因此，如果此后验概率p(z|x,θ)不能用分析式直接求得，可以用采样的方法近似之。Q(θ，θold)写成积分形式为 integral(p(z|x,θold)ln(p(z,x|θ)))dz。

这里，p(z|x,θ)成为我们采样的目标函数，分析式未知。

【蒙特卡洛期望最大化算法Monte Carlo EM algorithm】

MCEM算法可以解决此问题，该算法从当前p(z|x,θold)中取L个样本，根据积分的定义，Q(θ，θold)可由下式的样本逼近：

使用此算法计算Q，完成EM算法的运算。

以上就是EM算法和蒙特卡洛采样的关系。在上上篇文章讨论重要性采样算法时提到，其局限性在于假设分布和目标分布的差异不能太大。而MCMC采样可以克服此问题。

【马尔科夫链蒙特卡洛采样算法Markov Chain Monte Carlo】

MCMC是使用马尔科夫链生成样本的方法，这个链设计为在最重要的区域花费最多的时间，样本x的生成模拟从目标分布p(x)中产生的过程。

对比MCMC和重要性采样（见上上篇文章【2采样算法】），MCMC采用可变的假设分布，即使用q(x'|x)而不是q(x)，x'是新的要采样的状态，x是前一个样本。每产生一个样本，原假设分布都会更新为一个新的假设分布。

MCMC算法不是一个特定的算法，而是使用上述思想的一类算法。

从最简单的一个成员开始：

【Metropolis-Hasting算法】

1、初始化x0

2、循环N-1次：对i=0，1，一直到N-1

从【0，1】均匀分布中随机采样小数u；

从假设分布q(x'|xi)中随机采样x'；

判断：若u< min {1, p_(x')*q(xi|x')/(p_(xi)*q(x'|xi))},x(i+1)=x'接受x'为新样本，

否则x(i+1)=x(i)

例如，令每一个假设分布为均值等于x的高斯分布，【1蒙特卡洛原则】中已经讨论过p_(x)也易得，因此可以很方便对p(x)采样。

如图所示，每次取一个样本修正一次假设分布，避免了重要性采样算法的误差。马尔科夫链性质体现在假设分布的设置。

【废弃样本Burn-In period】

在MCMC算法中，初次取样本时，可能取在上图中的长尾之中，导致经过很多次迭代才把假设分布逐渐逼近到目标分布去，因此可能需要在采样完成后将前面这些废弃样本删除，例如前1000个样本可能都是无效的。

本文讨论了EM算法中应用蒙特卡洛采样，以及一种马尔可夫链蒙特卡洛采样算法-Metropolis hasting算法，之后将继续讨论其他MCMC算法。