EM算法

---

简介

期望最大化(Expectation-Maximization)算法最初是由Ceppellini等人1950年提出的，目前引用的较多的是1977年Dempster等人的工作。其主要用于从不完整的数据中计算最大似然估计，是一个在含有隐变量的模型中常用的算法，比如最大似然估计（MLE）和最大后验估计（MAP）。在GMM、HMM、PCFG、机器翻译对齐模型以及K-Means聚类方法中均有它的影子。

---

预备知识

为了理解便利，复习以下最大似然，Jensen不等式相关概念。

最大似然

假设通过抽样调查学校男生的身高分布，我们随机挑选100个男生统计他们的身高。假设身高服从高斯分布$p(x|\mu,\sigma)$，但是这个分布的参数（均值$\mu$和方差$\sigma^{2}$）事前是不知道的。

那么这两个参数该如何确定呢？这个问题先放一放。我们先看另一个问题，取到这100个男生身高的概率是多少？

首先，假设每个抽取到的男生身高是独立地从高斯分布$p(x|\mu,\sigma)$中抽取出来的（记$\theta = [\mu, \sigma]^T$）,则高斯分布为$p(x|\theta)$），那么抽到这100个男生的概率则为：

$L(\theta)=L(x_1, \cdots,x_n;\theta)=\prod\limits_{i=1}^{n}{p(x_i;\theta)},\theta \in \Theta$

其中$x_i$表示第$i$个男生的身高，$\Theta$为模型的参数空间。可以发现$\{x_1,\cdots,x_{100}\}$是已知的，而$\theta$是未知的，所以$L(\theta)$是$\theta$的函数。

这个函数反映的是在不同的参数$\theta$取值下，取得当前这个样本集的可能性。我们称这个函数为参数$\theta$相对于样本集X的似然函数(likehood function)。

我们让抽到这100个男生身高出现的概率最大，即需要找到一个参数$\theta$，使似然函数$L(\theta)$最大（符号描述为$\hat{\theta}=\arg\max l(\theta)$）。为了便于求解对似然函数取对数：

$H(\theta)=\ln L(\theta)=\ln \prod\limits_{i=1}^{n}{p(x_i;\theta)}=\sum\limits_{i=1}^{n}{\ln p(x_i; \theta)}$

要使$\theta$的似然函数$L(\theta)$极大化，也就是通过函数求导令导数为0，构建关于$\theta$的方程组，最后对方程组进行求解（前提是$L(\theta)$连续可微）。当$\theta$是包含多个参数的向量时，求$L(\theta)的梯度$。

求最大似然函数估计值的一般步骤：
1. 定义分布（分布的函数形式），写出似然函数；
2. 对似然函数取对数，并整理；
3. 求导数，令导数为0，得到似然方程；
4. 解似然方程，得到的参数即为所求。

Jensen不等式

设$f$是定义域为实数的函数，如果对于所有的实数$x$，$f^{''}(x) \ge 0$，那么$f$是凸函数。当$x$是向量时，如果其Hessian矩阵$H$是半正定的（$H \geq 0$），那么$f$是凸函数。如果$f^{''}(x)>0$或者$H>0$，那么称$f$是严格凸函数。如果$f$是凸函数，$X$是随机变量，那么$E[f(X)] \geq f(EX)$。特别地，如果$f$是严格凸函数，那么$E[f(X)]=f(EX)$，当且仅当$p(x=E[X])=1$，也就是说X是常量。这里我们将$f(E[X])$简写为$f(EX)$。如果用图表示会很清晰：

图中，实线$f$是凸函数，$X$是随机变量，有$0.5$的概率是$a$，有$0.5$的概率是$b$。（就像掷硬币一样）。$X$的期望值就是$a$和$b$的中值了，图中可以看到$E[f(X)] \geq f(EX)$成立。当$f$是（严格）凹函数当且仅当$-f$是（严格）凸函数。Jensen不等式应用于凹函数时，不等号方向反向，也就是$E[f(X)] \leq f(EX)$。

---

EM算法的实例描述

接着上面学生身高的例子，我们将问题再复杂一点。现在我们随机挑选100个男生100个女生，我们将他们混在一起，不知道每个样例是男生还是女生，这是我们需要估计：1. 男生身高正态分布的参数; 2. 女生身高正态分布的参数; 3. 样例是男是女的概率分布。

当估计两个正态分布的参数时，需要将样例划分为男生身高集合和女生身高集合，然后在男生集合和女生集合分别利用最大似然估计来估计参数。

当在划分样例到男生集合和女生集合时，我们需要将样例身高分别代入到分布中，比较在两个分布中该样例出现的概率大小以确定样例是男是女。

这就成了一个先有鸡还是先有蛋的问题了，即要估计两个分布参数，需要先判定所有样例是男是女。要知道是男是女，先要知道分布的参数。EM算法针对这样的问题给出了解决方案。

对于上面的例子，首先随便猜两个正态分布的参数，每个样例就可以进行男女判定，这是Expectation步骤。每个样例有了性别归属后，根据之前说的极大似然，分别用男生身高样例估计男生身高正态分布参数，女生身高样例估计女生身高正态分布参数，这个是Maximization步骤。然后，当我们更新了这两个分布的时候，身高样例是男是女的概率又变了，那么我们就再需要调整E步……如此往复，直到参数基本不再发生变化为止。下面给出EM算法的推导过程。

---

EM算法推导

假设：

1. $\{x^{(1)}, x^{(2)}, \cdots, x^{(m)}\}$为样例集合(上标表示第$i$个样例)，且样例间相互独立, 每个$x^{(i)}为一个随机变量$。

2. $\{z^{(1)}, z^{(2)}, \cdots, z^{(m)}\}$,其中$z^{(i)}$表示第$i$个样例的隐含类别，也是一个随机变量。

3. $Q$为样例的隐含类别分布，$Q$满足的条件是$\sum\limits_{z^{(i)}}{Q{(z^{(i)})}}=1, Q{(z^{(i)})}\geq 0$（如果$z^{(i)}$是连续性的，那么$Q$是概率密度函数，需要将求和符号换做积分符号。）

由上述描述可得：

$$\begin{equation} \begin{array}{rcl} \sum\limits_{i}{\log{p(x^{(i)};\theta)}} & = & \sum\limits_{i}{\log\sum\limits_{z^{(i)}}{p(x^{(i)},z^{(i)};\theta)}} \\ & = & \sum\limits_{i}{\log\sum\limits_{z^{(i)}}{Q(z^{(i)})\frac{p(x^{(i)},z^{(i)};\theta)}{Q(z^{(i)})}}}\\ & \geq & \sum\limits_{i}{\sum\limits_{z^{(i)}}{Q(z^{(i)})\log{\frac{p(x^{(i)},z^{(i)};\theta)}{Q(z^{(i)})}}}} \end{array} \end{equation}$$
从等式到不等式利用了jensen不等式:

a. $log(x)$是凹函数（二阶导$-\frac{1}{x^2}<0$），则有$E[\log(X)] \leq \log(EX)$

b. 令$p(k) = Q(z^{(i)})$，$g(z^{(i)})=\frac{p(x^{(i)},z^{(i)};\theta)}{Q(z^{(i)})}$
c. 则

$$\begin{equation}\begin{array}{rcl}\log\sum\limits_{z^{(i)}}{Q(z^{(i)})\frac{p(x^{(i)},z^{(i)};\theta)}{Q(z^{(i)})}}&=&\log\sum\limits_{z^{(i)}}{g(z^{(i)})p(k)}\\&=&log(E(g(z^{(i)})))\end{array}\end{equation}$$
d. 由 $$\begin{equation}\begin{array}{rcl}log(E(g(z^{(i)}))) & \geq & E[\log(g(z^{(i)}))]\\&=&\sum\limits_{z^{(i)}}{Q(z^{(i)})\log{\frac{p(x^{(i)},z^{(i)};\theta)}{Q(z^{(i)})}}}\end{array}\end{equation}$$

e. 故得$\log\sum\limits_{z^{(i)}}{Q(z^{(i)})\frac{p(x^{(i)},z^{(i)};\theta)}{Q(z^{(i)})}} \geq \sum\limits_{z^{(i)}}{Q(z^{(i)})\log{\frac{p(x^{(i)},z^{(i)};\theta)}{Q(z^{(i)})}}}$

现将不等式改写为$L(\theta)\geq J(z,Q)$，那么可以通过不断的最大化下界$J(z,Q)$，来使得$L(\theta)$不断提高，最终达到它的最大值。

首先固定$\theta$（对应于上例就是确定两个分布的参数），调整$Q(z)$使下界$J(z,Q)$上升至与$L(\theta)$在此点$\theta$处相等（绿色曲线到蓝色曲线），这个过程对应Expectation步骤。

然后固定$Q(z)$，调整θ使下界$J(z,Q)$达到最大值（$\theta_{t}$到$\theta_{t+1}$），这个过程对应于Maximization步骤。

然后再固定$\theta$，调整$Q(z)\cdots$直到收敛到似然函数$L(\theta)$的最大值处的$\theta^{*}$。

那么每次迭代什么时候算是调整好了呢？当不等式变成等式时，说明我们调整后的概率能够等价$L(\theta)$，按照这个思路，我们要找到等式成立的条件。在Jensen不等式中说到，当自变量X是常数的时候，$E[\log(X)] = \log(EX)$等式成立。而在这里，$X=g(z^{(i)})=\frac{p(x^{(i)},z^{(i)};\theta)}{Q(z^{(i)})}=c$。

再推导下，由于$\sum\limits_{z}{Q(z^{(i)})}=1$，则得到分子的和等于c。

$p(x_{(i)},z^{(i)};\theta) = c\cdot Q({z^{(i)}})$

$\sum\limits_{z}{p(x_{(i)},z^{(i)};\theta)}=\sum\limits_{z}{c\cdot Q(z^{(i)})}=c\sum\limits_{z}{Q(z^{(i)})}=c$

则：

$$\begin{equation}\begin{array}{rcl}Q(z^{(i)})&=&\frac{p(x^{(i)}, z^{(i)};\theta)}{c}\\&=&\frac{p(x^{(i)}, z^{(i)};\theta)}{\sum\limits_{z}{p(x^{(i)},z^{(i)};\theta)}}\\&=&\frac{p(x^{(i)}, z^{(i)};\theta)}{p(x^{(i)};\theta)}\\&=& p(z^{(i)}|x^{(i)};\theta)\end{array}\end{equation}$$

至此，我们推出了在固定参数$\theta$后，使下界拉升的$Q(z)$的计算公式就是后验概率，解决了$Q(z)$如何选择的问题，这就是Expectation步骤，建立了$L(\theta)$的下界。

接下来的Maximization步骤，就是在给定$Q(z)$后，调整$\theta$，去极大化$L(θ)$的下界$J$（在固定$Q(z)$后，下界还可以调整的更大）。那么一般的EM算法的步骤如下：

---

EM算法的流程

1、初始化分布的参数$\theta$
2、E步骤：根据参数初始值或上一次迭代的模型参数来计算出隐性变量的后验概率，其实就是隐性变量的期望。作为隐藏变量的现估计值：$Q(z^{(i)})=p(z^{(i)}|x^{(i)};\theta)$
3、M步骤：将似然函数最大化以获得新的参数值：$\theta=\arg\max\limits_{\theta}{\sum\limits_{i}{\sum\limits_{z^{(i)}}{Q(z^{(i)})\log{\frac{p(x^{(i)},z^{(i)};\theta)}{Q(z^{(i)})}}}}}$
4、这个不断的迭代，就可以得到使似然函数L(θ)最大化的参数θ了。

---

EM的收敛性说明

感性的说，因为下界不断提高，所以极大似然估计单调增加，那么最终我们会到达最大似然估计的最大值。理性分析的话，就会得到下面的东西：

$$\begin{equation}\begin{array}{rcl}l(\theta^{(t+1)})&\geq &\sum\limits_i{\sum\limits_{z^{(i)}}{Q^{(t)}(z^{(i)})\log{\frac{p(x^{(i)},z^{(i)};\theta^{(t+1)})}{Q^{(t)}(z^{(i)})}}}}\\&\geq & \sum\limits_i{\sum\limits_{z^{(i)}}{Q^{(t)}(z^{(i)})\log{\frac{p(x^{(i)},z^{(i)};\theta^{(t)})}{Q^{(t)}(z^{(i)})}}}}\\&=&l(\theta^{(t)})\end{array}\end{equation}$$

参考资料：

http://blog.csdn.net/zouxy09/article/details/8537620
http://www.cnblogs.com/jerrylead/archive/2011/04/06/2006936.html
http://www.cnblogs.com/jerrylead/archive/2011/04/06/2006924.html
http://www.cnblogs.com/jerrylead/archive/2011/04/06/2006910.html
http://www.tuicool.com/articles/Av6NVzy
http://wenku.baidu.com/link?url=sp9YuoMsTErNUW6lHmnt2P5frnrMYbUQiLvuGo9R4qKDMgiYqUTPTOJIBNCjb5emAJc7WNhbyLFemCXX2F4lc2y90QjUVdtGcs7VexHbEei
http://luowei828.blog.163.com/blog/static/3103120420120142193960/
http://www.cnblogs.com/zhengyuhong/articles/3461639.html
http://www.52nlp.cn/tag/em%E7%AE%97%E6%B3%95
http://www.cnblogs.com/tornadomeet/archive/2012/07/14/2591685.html
http://blog.csdn.net/abcjennifer/article/details/8198352
http://blog.163.com/huai_jing@126/blog/static/17186198320119231094873/