EM-最大期望算法

Posted on   |   In   |     |   views

对于EM算法，一直都是云里雾里。所以，今天索性就下个决定，不搞懂它，决不罢休。通过今天的学习，加上之前的基础，EM算法终于算是理清思绪了。回头想想，果真是如何做事不下定决心，真的很难有结果。下面，打算将EM算法的整个推导过程总结一遍，达到理解并掌握的目的。

首先来看一张EM算法的聚类图，来自wikipedia，效果比较直观。

期望最大算法是一种从不完全数据或有数据丢失的数据集（存在隐含变量）中求解概率模型参数的最大似然估计方法。EM算法是机器学习十大算法之一，或许确实是因它在实际中的效果很好吧。下面先来说说它的定义。

一、定义

EM算法，全称Expectation Maximization Algorithm，译作最大期望化算法或期望最大算法，它是一种迭代算法，用于含有隐变量（hidden variable）的概率参数模型的最大似然估计或极大后验概率估计。

二、Jensen不等式

在完善EM算法之前，首先来了解下Jensen不等式，因为在EM算法的推导过程中会用到。

Jensen不等式在优化理论中大量用到，首先来回顾下凸函数和凹函数的定义。假设f是定义域为实数的函数，如果对于所有的x，f(x)的二阶导数大于等于0，那么f是凸函数。当x是向量时，如果hessian矩阵H是半正定（即H>=0），那么f是凸函数。如果，f(x)的二阶导数小于0或者H>0，那么f就是凹函数。

Jensen不等式描述如下：

  • 如果$f$是凸函数，$X$是随机变量，则$E[f(X)]>=f(E[X])$，当$f$是严格凸函数时，则$E[f(X)]>f(E[X])$；
  • 如果$f$是凹函数，$X$是随机变量，则$f(E[X])<=E[f(X)]$，当$f$是（严格）凹函数当且仅当$-f$是（严格）凸函数。

通过下面这张图，可以加深印象：

上图中，函数$f$是凸函数，$X$是随机变量，有0.5的概率是$a$，有0.5的概率是$b$。X的期望值$E(x)$就是a和b的中值$\frac{a+b}{2}$了，图中可以看到E[f(X)]>=f(E[X])成立。

三、EM思想

EM算法推导过程中，会使用到极大似然估计法估计参数，所以，首先给出一个求最大似然函数估计值的一般步骤：

  • （1）写出似然函数；
  • （2）对似然函数取对数，并整理；
  • （3）求导数，令导数为0，得到似然方程；
  • （4）解似然方程，得到的参数即为所求；

关于极大似然估计的实例，这里就不再提及了，下面介绍EM算法。

给定m个训练样本{$x^{(1)},…,x^{(m)}$},假设样本间相互独立，我们想要拟合模型$p(x ,z)$到数据的参数。根据分布，我们可以得到如下这个似然函数：

第一步是对极大似然函数取对数，第二步是对每个样本实例的每个可能的类别$z$求联合分布概率之和。然而，直接求这个参数$\theta$会比较困难，因为上式存在一个隐含随机变量（latent random variable）$z$。如果$z$是个已知的数，那么使用极大似然估计来估算会很容易。在这种$z$不确定的情形下，EM算法就派上用场了。

EM算法是常用的估计参数隐变量的利器。对于上述情况，由于存在隐含变量，不能直接最大化$l(\theta)$，所以只能不断地建立$l$的下界（E-step），再优化下界（M-step），依次迭代，直至算法收敛到局部最优解。这就是EM算法的核心思想，简单的归纳一下：

EM算法通过引入隐含变量,使用MLE（极大似然估计）进行迭代求解参数。通常引入隐含变量后会有两个参数，EM算法首先会固定其中的第一个参数，然后使用MLE计算第二个变量值；接着通过固定第二个变量，再使用MLE估测第一个变量值，依次迭代，直至收敛到局部最优解。

E-Step和M-Step。

  • E-Step：通过observed data和现有模型估计参数估计值 missing data；
  • M-Step：假设missing data已知的情况下，最大化似然函数。

是否会收敛：由于算法保证了每次迭代之后，似然函数都会增加，所以函数最终会收敛(最后有推到)。

四、EM推导

下面来推导EM算法：

对于每个实例$i$,用$Q_{i}$表示样本实例隐含变量$z$的某种分布，且$Q_i$满足条件（$\sum_zQ_i(z)=1,Q_i(z)>=0$） ,如果$Q_i$是连续性的，则$Q_i$表示概率密度函数，需要将求和符号换成积分符号。

对于上面的式子，做如下变换：

上面三个式子中，式（1）是根据联合概率密度下某个变量的边缘密度函数求解的（这里把z当作是随机变量）。对每一个样本i的所有可能类别$z$求等式右边的联合概率密度函数和，也就得到等式左边为随机变量x的边缘概率密度。由于对式（1）直接求导非常困难，所以将其分子分母都乘以一个相等的函数$Q_{i}(z^{(i)})$，转换为式（2）。而在式（2）变为式（3）的过程，采用的是上面提到的Jensen不等式。分析过程如下：

首先，把（1）式中的log函数体看成是一个整体，由于log(x)的二阶导数为$-\frac{1}{x^2}$,小于0，为凹函数。所以使用Jensen不等式时，应用第二条准则：f(E[X])>=E[f(x)]。

到这里，问题简化为如何求解随机变量的期望。还记得当年读大学的时候，概率论中的随机变量的期望计算方法么，不记得也没关系，下面这张图比较详细：

因此，结合上面的知识点，我们可以把(2)式当中的$Q_i(z^{(i)})$看成相应的概率$p_i$，把$\frac{p(x^{i},z^{(i)};\theta)}{Q_i(z^{(i)})}$看作是$z^{(i)}$的函数$g(z)$，类似地，根据期望公式$E(x)=\sum x*p(x)$可以得到：

其实就是$\frac{p(x^{i},z^{(i)};\theta)}{Q_i(z^{(i)})}$的期望，再根据凹函数对应的Jensen不等式性质：

因此便得到了公式（3）。OK，现在我们知道上面的式（2）和式（3）两个不等式可以写成：似然函数$L(θ)>=J(z,Q)$的形式（$z$为隐含变量），那么我们可以通过不断的最大化$J$的下界，来使得$L(θ)$不断提高，最终达到它的最大值。使用下图会比较形象：

这里来说说上图的内在含义：首先我们固定θ，调整Q(z)使下界J(z,Q)上升至与L(θ)在此点θ处相等（绿色曲线到蓝色曲线），然后固定Q(z)，调整θ使下界J(z,Q)达到最大值（θt到θt+1），然后再固定θ，调整Q(z)……直到收敛到似然函数L(θ)的最大值处的θ。

这里有两个问题：

  1. 什么时候下界J(z,Q)与L(θ)在此点θ处相等？
  2. 为什么一定会收敛？

首先来解释下第一个问题（1.什么时候下界J(z,Q)与L(θ)在此点θ处相等？）。在Jensen不等式中说到，当自变量X=E(X)时，即为常数的时候，等式成立。而在这里，为：

对该式做个变换，并对所有的$z$求和，得到

$$\sum_z{p(x^{i},z^{(i)};\theta)}=\sum_z{Q_i(z^{(i)})}c$$

因为前面提到$\sum_zQ_i(z)=1$（概率之和为1），所以可以推导出：

$$\sum_z{p(x^{i},z^{(i)};\theta)}=c$$

因此也就可以得到下面的式子：

至此，我们推出了在固定参数$\theta$后，使下界拉升的$Q(z)$的计算公式就是后验概率（条件概率），一并解决了$Q(z)$如何选择的问题。此步就是EM算法的E-step，目的是建立$L(\theta)$的下界。接下来的M-step，目的是在给定$Q(z)$后，调整$\theta$，从而极大化$L(\theta)$的下界$J$（在固定$Q(z)$后，下界还可以调整的更大）。到此，可以说是完美的展现了EM算法的E-step & M-step，完整的流程如下：

第一步，初始化分布参数$\theta$；
第二步，重复E-step 和 M-step直到收敛：

  • E步骤：根据参数的初始值或上一次迭代的模型参数来计算出的隐性变量的后验概率（条件概率），其实就是隐性变量的期望值。作为隐藏变量的现有估计值：

  • M步骤：最大化似然函数从而获得新的参数值：

通过不断的迭代，然后就可以得到使似然函数L(θ)最大化的参数θ了。

上面多次说到直至函数收敛，那么该怎么确保EM收敛呢？(②为什么一定会收敛？)，下面进入证明阶段.

假定$\theta^{(t)}$和$\theta^{(t+1)}$是EM第t次和t+1次迭代后的结果。如果我们证明了$l(\theta^{(t)})<=l(\theta^{(t+1)})$，也就是说极大似然估计单调增加，那么最终我们就会得到极大似然估计的最大值。下面来证明，选定$\theta^{(t)}$之后，我们得到E步：

这一步保证了在给定$\theta^{(t)}$时，Jensen不等式中的等式成立，也就是

然后进行M步，固定$Q_i^{(t)}(z^{(i)})$,并将$\theta^{(t)}$试作变量，对上面的式子求导，得到$\theta^{(t+1)}$,这样经过一些推导会有以下式子成立：

在公式（4）中，得到$\theta^{(t+1)}$,只是最大化$l(\theta^{(t)})$，也就是$l(\theta^{(t+1)})$的下界，并没有使等式成立，等式成立只有在固定$\theta$，并按E步得到$Q_i$时才能成立。

这样就证明了$l(\theta)$会单调增加。如果要判断收敛情况，可以这样来做：一种收敛方法是$l(\theta)$不再变化，还有一种就是变化幅度很小,即根据$l(\theta)^{(t+1)}-l(\theta)^{(t)}$的值来决定。

EM算法类似于坐标上生法（coordinate ascent）：E步：固定θ，优化Q；M步：固定Q，优化θ；交替将极值推向最大。

五、应用

六、References