EM算法 - 1 - 介绍

时间：2016-05-12 11:28:38 阅读：273 评论：0 收藏：0 [点我收藏+]

标签：

声明：

1，本篇为个人对《2012.李航.统计学习方法.pdf》的学习总结，不得用作商用，欢迎转载，但请注明出处（即：本帖地址）。

2，由于本人在学习初始时有很多数学知识都已忘记，所以为了弄懂其中的内容查阅了很多资料，所以里面应该会有引用其他帖子的小部分内容，如果原作者看到可以私信我，我会将您的帖子的地址付到下面。

3，如果有内容错误或不准确欢迎大家指正。

4，如果能帮到你，那真是太好了。

在开始之前需要准备些预备知识，如果你已经全部掌握那就直接向下翻到EM算法这一章吧。

极大似然估计

由于EM算法和极大似然估计有斩不断的联系，所以我们先复习下这个。

多数情况下，我们是根据条件/概率分布来推算结果，而极大似然估计就是反过来，即：已经知道了结果，然后寻求使该结果出现的可能性最大的条件/概率分布。

举个例子。

如果其他条件都固定的话“抽烟者发生肺癌的概率是不抽烟者的5倍”，那么如果我已经知道有个人是肺癌，我想问这个人抽烟还是不抽烟时你怎么判断？我想你会说这个人抽烟，因为这个条件是其患肺癌的最大可能。

代入些数学用语描述上面的例子就是：

输入：

1，观测结果集合T = {此人患肺癌}

2，已知是否患肺癌满足这个分布：抽烟者发生肺癌的概率是不抽烟者的5倍（如：知道结果符合高斯分布f(x | u,σ²)，当然了，我们不知道该分部的参数，不然还求什么极大似然....）

求：

这个分布的参数θ（如：知道结果符合高斯分布，但却不知道高斯分布的参数u,σ²，目的是求参数）

我们把上面的内容再数学化些：

输入：

1，样本集合T={x1,x2, ..., xn}，且样本集合中的元素相互独立。

2，该样本集合服从某分布：p(x| θ)（假设符合高斯分布f(x | u,σ²)，那θ=( u,σ²)。

求解：

得到该样本集合T的极大似然估计θ（即：当θ是多少时，以θ为参数的分布最能拟合样本集合。假设样本集合符合高斯分布，那就是求u,σ²）。

到此，极大似然估计是怎么回事应该清楚了，下面看看如何求。

极大似然估计，极大似然估计，既然是求极大，那当然离不开求导/求偏导了，为此我们就构建个关于θ的函数吧。

首先，既然是求最可能得到该样本集合的θ，那我们得把得到该样本集合中某个元素xi的概率表示出来，记为：P(xi; θ)，其中i = 1, 2, ..., n，θ服从某分布(别问我服从什么分布，这里是为了讲解，应用时自己根据情况判断)。

既然样本集合中的元素相互独立，那样本集合的联合概率就是：

技术分享

但上面的式子写起来好麻烦啊，那么我们想办法用个方式表示上面的式子吧。嗯....用L(x1, x2, ..., xn; θ)怎么样？....每次都要写x1, x2, ..., xn好麻烦。这个只是求θ，那就用L(θ)表示吧，于是：

技术分享

上面的L(θ)就是似然函数。

不过上面的L(θ)是连乘的，有时为了便于分析，会对L(θ)取对数，将其变成连加的：

技术分享

PS：其实对L(θ)取对数将其变成连加的H(θ)还有一个原因：通常L(θ)中每个p(xi; θ)都很小，许多很小的数字相乘起来在计算机里很容易造成浮点数下溢，所以对其取对数将其变成连加的。

似然函数定义好了，那我们的最终目的也用纯数学语言再次表示下吧，即，我们需要求最大似人估计量：

θ*= argmaxL(θ) 或θ* = argmaxH(θ)

最后，总结下求最大似然估计值的一般步骤：

1，写出似然函数；

2，对似然函数取对数，并整理；

3，求导数/偏导数，令其为0，得到似然方程；

4，解似然方程，得到的参数即为所求。

边缘概率与联合概率的关系

技术分享

如上图所示：

X={x1=1,x2=2}，Y={y1=1, y2=2, y3=3}

在Y确定为1的情况下

P.₁= P{Y=1} = P{X=1,Y=1} + P{X=2, Y=1} = P{X=x1, Y=1} + P{X=x2, Y=1} = P{X=xi,Y=1}

其他同理。

条件概率公式

P(A|B)= P(AB)/P(B)

Jensen不等式

我们常见到的Jensen不等式的版本是：

当φ是凸函数时(φ’’>0)：φ( E(x) )<= E(φ(x) )

不过这里用到的是下面这个版本：

当φ是凸函数时

技术分享

EM算法

EM算法与“鸡生蛋蛋生鸡”问题

EM算法是为了解决“最大似然估计”中更复杂的情形而存在的。

这里“极大似然估计中更复杂的情形”是什么情形呢？

我们知道极大似然估计是求解实现结果的最佳参数θ，但极大似然估计需要面临的概率分布只有一个或者知道结果是通过哪个概率分布实现的，只不过你不知道这个概率分布的参数。而如果概率分布有多个呢或者你不知道结果是通过哪个概率分布实现的？于是别说去确定“这些概率分布”的最佳参数了，我们连最终结果是根据哪个概率分布得出来的都不知道，这就是EM算法要面临的情况了。

用现实生活中的例子来说明的话就是：

极大似然估计面临的情况：

一个快递员给你送货。若他到你家只有一条路（结果的实现依赖一个概率分布），但却不知道这条路今天修不修路（不知道该概率分布的参数），修路的话今天快递员就没法送货，若结果是快递员送到货了，那这条路修路了没？答案很明显：没修路。

EM算法面临的情况：

还是上面的快递员送货，不过我把条件改成：快递员到你家的路有N条（结果的实现依赖多个概率分布），但快递员只会选择一条路，即，今天他不会选择第二条路，若他选择的路修路，那他就不给你送货了，即使这会而让你暴跳如雷。问：如果今天快递员送到货了，则他选择的哪条路？那条路修路了吗？对于这个，因为你不知道他选择的哪条路（他把货送到就走了，根本不给你问他话的时间），所以你唯一能做的就是估计出这N条路被他选择的N个概率（即：每个概率分布的权值），然后在根据极大似然估计来得出：这条路没修路（求出每个概率分布的参数）。

一句话总结区别就是：最大似然估计和EM算法都是根据实现结果求解概率分布的最佳参数θ，但最大似然估计中知道每个结果对应哪个概率分布（我知道哪个概率分布实现了这个结果），而EM算法面临的问题是：我不知道哪个概率分布实现了该结果。

于是我们面临了这样一个困境：

一个样本集合的元素是从多个概率分布中得出的，其中，某个元素只可能是一个概率分布的结果，但我们不知道这个元素从哪个概率分布中得出。

所以，我们不仅要求出多个概率分布的最佳参数，而且要求出多个概率分布在模型中的比率。

而如果我们知道样本集合中的一个个元素是从哪个概率分布得出的，那直接运用极大似然估计就可求出多个概率分布的最佳参数，进而求出概率分布的比率；

反过来，如果我们知道这些概率分布的最佳参数，那我们也可以知道样本集合中的一个个元素分别对应哪个概率分布，也能进而求出概率分布的比率；

于是，就陷入“鸡生蛋蛋生鸡”的死循环了。

而EM算法就是解决这样的问题。

EM算法的求解思想

在说明EM算法的求解思想前，我们先总结下上面的内容。

一般的用Y表示观测到的随机变量的数据，Z表示隐随机变量的数据(因为我们观测不到结果是从哪个概率分布中得出的，所以将这个叫做隐变量)。于是Y和Z连在一起被称为完全数据，仅Y一个被称为不完全数据。

这时有没有发现EM算法面临的问题主要就是：有个隐变量数据Z。而如果Z已知的话，那问题就可用极大似然估计求解了。

于是乎，怎么把Z变成已知的？

我先举个日常生活的例子。

结果：大厨把锅里的菜平均分配到两个碟子里

难题：如果只有一个碟子乘菜那就什么都不用说了，但问题是有2个碟子，而因为根本无法估计一个碟子里应该乘多少菜，所以无法一次性把菜完全平均分配。

解法：大厨先把锅里的菜一股脑倒进两个碟子里，然后看看哪个碟子里的菜多，就把这个碟子中的菜往另一个碟子中匀匀，之后重复多次匀匀的过程，直到两个碟子中菜的量大致一样。

上面的例子中，平均分配这个结果是“观测数据”，为实现平均分配而给每个盘子分配多少菜是“待求参数θ”，分配菜的手感就是“概率分布”。于是若只有一个盘子，那概率分布就确定了（“把锅里的菜全部倒到一个盘子”这样的手感是个人都有吧），而因为有两个盘子，所以“给一个盘子到多少菜才好”的手感就有些模糊不定，不过我们可以采用上面的解法来实现最终目标。

下面，我们把上面的做法对应到EM算法。

EM算法的思想是：

1，给θ自主规定个初值（既然我不知道想实现“两个碟子平均分配锅里的菜”的话每个碟子需要有多少菜，那我就先估计个值）；

2，根据给定观测数据和当前的参数θ，求未观测数据z的条件概率分布的期望（在上一步中，已经根据手感将菜倒进了两个碟子，然后这一步根据“两个碟子里都有菜”和“当前两个碟子都有多少菜”来判断自己倒菜的手感）；

3，上一步中z已经求出来了，于是根据极大似然估计求最优的θ’（手感已经有了，那就根据手感判断下盘子里应该有多少菜，然后把菜匀匀）；

4，因为第二步和第三步的结果可能不是最优的，所以重复第二步和第三步，直到收敛（重复多次匀匀的过程，直到两个碟子中菜的量大致一样）。

而上面的第二步被称作E步（求期望），第三步被称作M步（求极大化），于是EM算法就在不停的EM、EM、EM....，所以被叫做EM算法，你看，多形象（摊手）。

EM算法推导

该来的总会来的，最终我们还是要用数学语言来描述EM算法，那么下面我们看看EM算法的推导过程。

在此先将问题抽象：

已知模型为p(Y|θ)，Y=(y1, y2,..., yn)，求θ。

为了满足EM算法的情形，我们引入隐含变量Z=(z1, z2, ..., zn)。

解：

1，虽然我们面对的概率模型里含有隐变量，但我们的目标是不会变的，即，极大化观测数据（不完全数据）Y关于参数θ的对数似然函数，即，极大化：

技术分享

PS：上面式子的第二个等号利用了边缘概率与联合概率方面的知识；第三个等号利用了条件概率公式，于是P(y,z) = P(y|z) P(z)，最后把参数θ加上去，即，所有式子都添加θ，就成了9.12式。

2，我们开始对L(θ)进行迭代，使得新的估计值θ能使L(θ)增加，假设现在已经迭代了i次，于是就是使当前的L(θ)大于L(θ⁽ⁱ⁾)，即L(θ) >L(θ⁽ⁱ⁾)，就这样使L(θ)逐步达到最大值。

PS：这就是EM算法中的M步。

为此，我们考虑两者的差：

技术分享

利用Jensen不等式得到其下界：

技术分享

为了便于之后的书写，令

技术分享

则

L(θ)>= B(θ, θ⁽ⁱ⁾)

即函数B(θ, θ⁽ⁱ⁾)是L(θ)的一个下界。PS：L(θ⁽ⁱ⁾)= B(θ⁽ⁱ⁾, θ⁽ⁱ⁾)

因此，任何可以使B(θ, θ⁽ⁱ⁾)增大的θ也可以使L(θ)增大。

于是，为了使L(θ)有尽可能大的增长，选择θ⁽ⁱ⁺¹⁾使B(θ, θ⁽ⁱ⁾)达到极大，即：

技术分享

现在求θ⁽ⁱ⁺¹⁾的表达式

技术分享

3，为了计算上面的式子，我们就要计算（所以在算法中，这一步在上一步之前）

技术分享

而这个就是完全观测数据的对数似然函数logP(Y,Z | θ)关于在给定观测数据Y和当前参数θ⁽ⁱ⁾下对未观测数据Z的条件概率分布P(Z | Y, θ⁽ⁱ⁾)的期望，这个被称为Q函数，即：

技术分享

也就是说，要计算期望。

而这，就是EM算法中的E步。

4，不断迭代上面两步，直到收敛。

用图形解释EM算法的推导过程的话就是下面这样：

技术分享

如上图所示，图中上方曲线为L(θ)，下方曲线为B(θ, θ⁽ⁱ⁾)。在推导过程中已经知道B(θ,θ⁽ⁱ⁾)是对数似然函数L(θ)的下界，且两个函数在点θ = θ⁽ⁱ⁾处相等。之后EM算法找到下一个点θ⁽ⁱ⁺¹⁾使函数B(θ, θ⁽ⁱ⁾)极大化，也使函数Q(θ, θ⁽ⁱ⁾)极大化。这时L(θ) >= B(θ, θ⁽ⁱ⁾)，函数B(θ, θ⁽ⁱ⁾)的增加保证对数似然函数L(θ)在每次迭代中也是增加的。EM算法再点θ⁽ⁱ⁺¹⁾重新计算Q函数值，进行下次迭代。在这个过程中，对数似然函数L(θ)不断增大。

不过，从图中可以推断出EM算法可能因为陷入局部最优值而找不到全局最优解。这点需要注意。

EM算法描述：

下面就把上面的推导汇总成EM算法吧。

输入：

观测变量数据Y，隐变量数据Z，联合分布P(Y,Z | θ)，条件分布P(Z | Y,θ)；

输出：

模型参数θ；

解：

1，选择参数的初值θ⁽⁰⁾；

PS：这个初值不管是经验也好猜也好，反正你给它一个初始值就行。当然了，在实际使用中这个初始值往往由其他算法的结果给出，当然随机给他分配一个符合定义域的值也可以。

2，E步：记θ⁽ⁱ⁾为第i次迭代参数θ的估计值，在第i+1次迭代的E步，计算

技术分享

这里，P(Z |Y, θ⁽ⁱ⁾)是在给定观测数据Y和当前的参数估计θ⁽ⁱ⁾下隐变量数据Z的条件概率分布；

3，M步：求使Q(θ, θ⁽ⁱ⁾)极大化的θ，确定第i+1次迭代的参数的估计值θ⁽ⁱ⁺¹⁾

技术分享

4，重复上面两步，直到收敛。

EM算法收敛性的证明这里跳过，感兴趣的请自行查阅资料。

Q函数：

Q函数就是EM算法第二步的Q(θ, θ⁽ⁱ⁾)，我将Q函数单独放在一章是因为Q函数是EM算法的核心。

如果你已经深刻理解了上面的内容，那么会明白，整个EM算法最难的部分就是构建Q函数。为什么？因为EM算法的两步E步和M步在实际应用中就是“构建Q函数”和“通过偏导求极大值”，而后者我想大家都会，于是如何构造前者就是我们需要掌握的技能了。

PS：当然，对于最常用的“高斯混合模型”已经有现成的公式可以套用了。

首先先看下Q函数的定义：

完全观测数据的对数似然函数logP(Y,Z| θ)关于在给定观测数据Y和当前参数θ⁽ⁱ⁾下对未观测数据Z的条件概率分布P(Z | Y, θ⁽ⁱ⁾)的期望称为Q函数，即：

技术分享

于是你看懂了吗？

反正我一开始没看懂，虽然根据定义知道是求期望，但为了弄懂第二个等号真的用了我一番功夫。

下面我解释下。

首先Q函数是求期望，这个不用多说。

然后，是求什么的期望？为了方面说明，我们把上面的式子的条件给去掉，即Q(θ, θ⁽ⁱ⁾) = E[logP(Y, Z)]。这下一目了然了，是求函数logP(Y,Z)的期望，而既然是求函数的期望，那就需要用到下面的知识了：

技术分享

这里使用第一个。

首先对于logP(Y,Z)，因为Y是观测变量数据，Z是隐变量数据，因此Z是未知量。所以，对应到上面第一个求法就是：x(变量)=Z，g(x) = g(z) = logP(Y, Z)，又因为Q函数的定义中已经告知P(Z| Y, θ⁽ⁱ⁾)是条件概率分布，所以P(Z | Y, θ⁽ⁱ⁾)就对应上面的分布律，不过为了方面说明，这里还是先将条件给去掉，即分布律为P(Z)。

这样一来，套用上面的公式，就有了：

Ez[logP(Y,Z)] = logP(Y,Z)P(Z)

上面为了方面理解而把表示条件的参数都给去掉了，下面按照Q函数的定义把条件给还原回来，于是就有了Q函数定义中的式子。

PS：如果上面还原参数的这步你还有点晕，那就这样想：

联合概率P(Y, Z)是在参数为θ的某分布中，所以为了用公式表述完整，我们将其写成P(Y,Z | θ)，

分布律P(Z)是在观测变量Y和在EM算法中一步步迭代时当前步骤中已知的“以θ⁽ⁱ⁾为参数的某分布”下的分布律，所以为了用公式表述完整，我们将其写成P(Z | Y, θ⁽ⁱ⁾)

至于Q函数为什么写成Q(θ, θ⁽ⁱ⁾)，这是因为它想表达：在当前迭代中我要找出一个新θ，使得“以新θ为参数的某分布”优于“用上一次迭代中已找出的θ⁽ⁱ⁾为参数的某分布”。