EM算法通用形式（ESL 8.5.2）

8.5.2的EM算法的描述比较晦涩，这里总结一下EM算法的general形式。EM算法主要用于最大化似然函数。某些似然函数要最大化可能比较困难，但如果可以引入某些隐变量（latent data），那么这个最大化问题可能会变得简单一些。在The Element of Statistical Learning这本书里面这叫做data augmentation。需要注意的是，一般机器学习里面（尤其深度学习）说的data augmentation指的是对原来的数据进行一定的变换（比如计算机视觉里面把图片旋转拉伸放缩）再加入到原始数据里面，这样数据量大了，一定程度上面可以减少泛化误差，或者增加训练样本的量。可以参见巨著Deep Learning里面7.4节。还要注意，EM算法本质上是一个优化算法。从机器学习算法的三大构成（模型，目标函数，优化方法）来看，假定的数据分布是模型，最大似然是目标函数，而EM算法是用以优化出结果的方法。用形式化的方法来说，如果假定的数据分布的具有参数θ，拟合的数据是xi，那么似然函数则为

∑ilog p(xi;θ).

很多人可能对于那个分号有点懵逼，其实就是表示给定参数θ下xi的概率。

8.5.2对EM算法通用形式的描述有点不太好懂，先从吴恩达的比较好懂的版本开始说起。理解EM算法需要一个数据工具，Jensen不等式。

Theorem 如果f是一个凸函数，X是一个随机变量，那么

E[f(X)]≥f(EX)

如果f是严格凸的，那么仅当X=EX的时候上式取等号。

证明这个不等式也不难，需要用到凸函数的性质。如果x,c在f的定义域上，那么

f(x)−f(c)≥f′(c)(x−c)。

即f(x)−f(c)−f′(c)(x−c)≥0。假设c=EX，那么等式两边对X的密度函数g(x)积分（如果是离散变量就是乘以对应的概率函数求和了，或者说对概率测度求积分），可以得到

∫f(x)g(x)dx−f(c)−f′(c)(∫xg(x)dx−c)≥0E[f(x)]−f(c)≥0

因为c=EX，第三项消掉。继续代入，就可以证明出Jensen不等式了。因为凸函数和凹函数是对称的，因此如果f是一个凹函数，那么E[f(x)]≤f(EX)

回到最大似然的问题，假设我们额外增加的latant variable隐变量是z，对于混合高斯模型来说，这个隐变量就是指定某个xi属于哪一个高斯分布。对于不同的xi，z会有不同的分布。我们记对于xi来说z的分布是Qi（\sum Q_i(z) = 1, Q_i(z)\geq 0，z是连续变量时要积分）。另外，对于某个数据点xi，假设z的取值是zi。

l(θ)=∑ilog p(xi;θ)=∑ilog∑zip(xi,zi;θ)=∑ilog∑ziQi(zi)p(xi,zi;θ)Qi(zi)≥∑i∑ziQi(zi)logp(xi,zi;θ)Qi(zi)

最后一行用了Jensen不等式，logEzi∼Qip(xi,zi;θ)Qi(zi)≥Ezi∼Qilogp(xi,zi;θ)Qi(zi)。我们把最后一行的公式写成J(Q,θ)，那么这就是原来似然函数的一个下界lower bound。如果给定θ0，存在Q使得J(Q,θ0)=l(θ0)，那么最大化J(Q,θ)相当于逐渐地增加l(θ)。假设θ′=argmaxJ(Q,θ)，那么l(θ0)=J(Q,θ0)≤J(Q,θ′)≤l(θ′)。

根据Jensen不等式，如果J(Q,θ)=l(θ)成立，那么p(xi,zi;θ)Qi(zi)=c是一个常数，但是Qi(zi)是一个概率分布，因此

∑zQi(z)=1c∑zp(xi,z;θ)=1

因此

Qi(zi)=p(xi,zi;θ)∑zp(xi,z;θ)=p(xi,zi;θ)p(xi;θ)=p(zi∣xi;θ)

因此，当Qi(zi)是给定xi的后验概率，J(Q,θ)=l(θ)。接下来只要最大化J(Q,θ)，那么就可以迭代新的θ′使得l(θ)≤l(θ′)。计算Qi(zi)=p(zi∣xi;θ)就是EM的E-step，之后最大化J(Q,θ)就是EM的M-step。所以本来只是最大化l(θ)的问题变成了J(Q,θ)的坐标上升方法。可能有人会发现Qi(zi)=p(zi∣xi;θ)，此处的θ是前一次迭代得到的θ，而后面优化J(Q,θ)得到的θ是本次迭代得到的θ。这个在下面的版本中会体现得更好。

接下来是ESL上面8.5.2的版本，稍微不一样的使得问题复杂的就是一开始上述的Qi就固定了。为了化简符号，我们以X代表要拟合的数据，Z是隐变量。因此要拟合的最大似然就是l(θ)=P(X;θ)。然而

P(Z∣X;θ)=P(Z,X;θ)P(X;θ)

因此

P(X;θ)=P(Z,X;θ)P(Z∣X;θ)

取对数，那么就是

logP(X;θ)=logP(Z,X;θ)P(Z∣X;θ)

如果对于上式取对P(Z∣X;θ)的期望，因为左边不涉及Z，那么可以得到

logP(X;θ)=∑ZP(Z∣X;θ)logP(Z,X;θ)P(Z∣X;θ)＝EZlogP(Z,X;θ)P(Z∣X;θ)

这个已经很像是上面的J(Q,θ)了。但是这边一开始就用P(X;θ)=P(Z,X;θ)P(Z∣X;θ)构造了l(θ)=J(Q,θ)的情形。

然而只有两式相等的情形是不够的，需要知道如何迭代使得l(θ)不断上升。 接下来，假设我们已经知道了在第t次迭代θt的值，对于上面取期望的操作，我们取的是P(Z∣X;θt)的期望，因此

logP(X;θ)=∑ZP(Z∣X;θt)logP(Z,X;θ)P(Z∣X;θ)=∑ZP(Z∣X;θt)logP(Z,X;θ)−∑ZP(Z∣X;θt)logP(Z∣X;θ)=J1(θt,θ)−J2(θt,θ)

接下来是证明只需要优化J1(θt,θ)当中的θ就相当于优化了l(θ)。原因是J2(θt,θ)≥J2(θt,θt)。证明这个也是用Jensen不等式。

J2(θt,θ)−J2(θt,θt)=∑ZP(Z∣X;θt)logP(Z∣X;θ)P(Z∣X;θt)≥log∑ZP(Z∣X;θt)P(Z∣X;θ)P(Z∣X;θt)＝log∑ZP(Z∣X;θ)=0

因此对于任意θ必有J2(θt,θ)≥J2(θt,θt)。因此如果θ′最大化了J1(θt,θ)，那么可以令θt+1=θ′使得l(θt)≤l(θt+1)。因此在ESL的书上算法8.2，E-step是计算P(Z∣X;θt)形成J1(θt,θ)的确切形式，M-step是优化J1(θt,θ)。而且这里可以看到最大化J1(θt,θ)是不必要的，只需要升高一点就可以了。

以上两个版本是对EM的不同解读。共同之处都是需要额外的隐变量使得问题简化。第一个版本是用Jensen不等式求出一个过当前解的下界函数，通过优化这个相对简单的函数来迭代解。而第二个版本则是直接把似然函数变成两个部分来看，第二部分必然减少但是因为负号而增加，因而只需要优化第一部分。相对来说第二版本，即ESL书上的版本，相对比较晦涩。

可能有人会发现在这个版本里面，在对数log里面分母的一项是θ，是参与M-step优化的，而第一个版本里面的是θt，是不参与M-step优化的。然而因为第二个版本只需要优化J1(θt,θ)，即log里面分子的一项，所以和第一个版本殊途同归。

ESL上面8.5.3节的本质和第一个版本是一回事。这里也间接解决了课本习题8.1和8.2。吴恩达版本也是习题8.7的一个说明。

熟悉非负矩阵分解算法的人可能一下子就觉得很熟悉，其实NMF的paper里面也说到用了EM-like的方法。而且NMF里面也是通过不断放缩找到过当前解的一个上界函数来优化以迭代，思路本质上是一样的。