EM（expectation maximization）算法的理解和证明_聚类,em算法,数理统计证明

本文记录个人对EM算法的理解

首先为什么要使用EM算法，他适用于哪些场景呢？ 这个算法是在数据已知的前提下讨论参数的合理性。
换一句话说， 他是可以用似然函数建模的。等下会记录如何用数学语言描述这个算法。但是首先先使用经典的投硬币的问题将这个抽象的算法形象化。

什么情况下我们要使用 EM算法呢
假设我们现在有2枚硬币，分别记为硬币A 和硬币B，如果我拿起硬币A 投了5次得到了以下的结果：A:[+,+,-,-,+] 。这里“+”表示硬币正面朝上，“-”表示反面朝上。之后我拿起硬币B我也投5次得到了以下的结果B：[-,-,+,+,+]。 那么现在想要求到A出现正面的概率和B出现正面的概率是很简单的。

上述公式中N(A)所表示的是A抛的总次数，N(A)=+ 意思是在这些次数中出现正面的次数。
同理B出现正面的概率也可以用同样的方式表示。可以记为：

上述两个式子相对比较好理解。但是现在请设想这样一个情况。
我从AB两枚硬币中随机取一枚来投掷，并记录下投掷结果。经过N次这样的投掷之后，我也可以得到一组记录着正反面的数据。但是如果只给你这样一组数据，告诉你这一组数据是随机投掷AB两枚硬币产生的，有没有可能求出硬币A出现正面的概率和硬币B出现正面的概率呢？？这个时候就可以考虑使用EM算法了！

EM 算法的一般思路

这个问题和最开始的那一个投硬币的问题不同在哪里？
其实我们仔细的看看，会发现这一个问题中，我们不知道每一次投掷硬币的结果到底是A产生的还是B产生的，要是知道了就和最开始的那一个投硬币问题一样了。但是问题是我们不知道，而且它潜移默化地影响着我们最后的计算结果。 这种我们无法从数据中直接观测到的，但是又影响着最后模型输出的变量，我们叫它 ‘隐变量 / 潜变量’（Latent variable）

那这个时候我们就应该想方设法的把这个隐变量给表示出来呀！
假设我们现在有一个由投掷AB两枚硬币产生的一组结果 T：[+,+,-,-+,-]. 我们可以用另一个概率分布Q 来表示每一个结果是由A或者B产生的概率，用 Z 来表示那个隐藏变量（那么我们先明确隐藏变量 Z 其实就只有两个值呀！除了A就是B，不可能有C的情况出现，对吧！）。那么现在，我们就可以假设在 T 这一组结果中由A产生的概率是 Q(Z=A) ,同理也可以假设由B产生的概率是 Q(Z=B).那么这个时候求得A和B出现正面的概率也就不是难事了呀~

假如，我初始化Q(Z=A)=0.4, Q(Z=B)=0.6，那么在T中出现的3个正面朝上的结果中，就有
3 x 0.4 = 1.2 个是由硬币A产生的，3 x 0.6 = 1.8个是硬币B产生的。那同样的，在T中所有反面结果中，由A这个硬币产生的概率就是 3 x 0.4 = 1.2, 由B产生的概率是 3 x 0.6 = 1.8
那么这样一来 这个问题是不是就变成了最开始的那个最简单的投硬币问题了呢？
所以A出现正面和B出现正面的概率就可以这样来计算了（我们在这里把这两个概率用θ来表示吧~）：


（这个概率看上去好像很符合我们得到的数据，但是这是因为为了方便计算我举了一个很简单的例子啦~。真实情况比这个要复杂哦！）
得到这个之后我们继续！现在这个问题已经和最开始的那个投硬币的问题一样的，知道数据，知道两枚硬币出现正反面的概率了！但是这个值不一定准确！！！ 因为所有的前提都是我们随机初始化的，现在我们要做的就是要来优化我们的θ，争取找到最好的那一个θ。可以想想爬山的例子来帮助理解，你想去山顶，将你随机置于山上一个地方，然后一点一点向上爬的感觉。具体怎么做其实很简单。
既然现在我们已经知道了AB出现正反面的概率了，那么我们可以用这个概率去更新我们的Q，然后再用新的Q又更新θ，一直执行这个过程直到收敛！这里画个图可能会更加的直观一点儿~。这个更新的过程就类似图中所示 （图有点佛系。。。。）
这个过程一直进行直到收敛。收敛的条件可以是我们已经找到了最好的那个θ了，自然可以退出循环了。或者是达到了最大的迭代次数，也可以作为收敛条件。

所以正如算法的名字所述，整个算法的步骤就只有两步Expectation(E-step) 和 maximization（M-step）。从上文的描述中，可以我们在更新Q的那个过程其实就是在做Expectation。更新θ的那一步就是maximization，因为我们想要最好的θ。

用稍微数学一点的语言描述一下这个过程

在前文的时候我提到，类似这种知道数据讨论参数的情况，我们可以用最大似然函数对其进行建模，假设是针对参数θ的最大似然函数图像：

因为我们是在讨论参数θ的合理性，想要找到一个最佳的θ，所以这个最大似然函数的一定是关于θ的。这样的函数存在一个问题，他可能不是凸函数，在优化问题上，我们知道凸函数它一定是有一个全局最优解的，但是如果是非凸，那么可能存在多个局部最优解，不是很好进行优化。对于图中的这个函数如何来找到最佳的θ呢? 还是以爬山的那个场景来类比。我既然不知道哪个θ是最好的，我可以随机的初始化一个θ呀，就像随机把人放在山的某个位置一样。就像下图所示，我随机的初始化了一个θ1：
现在的问题是怎么让这个θ向上爬呢？既然我们知道凸函数一定有一个全局最优解，那为什么不利用起来呢？现在我创建一个简单的凸函数g（图中绿色表示的曲线）, 使这个函数g的交于θ1，并且保证g(θ1)的极大值要小于L(θ1)的极小值，其实就是要保证创建的函数g 要在L 下方。上述过程用图像画出来差不多，如下图所示：

那么凸函数g的最大值是很容易求的。这个时候我们可以把这个凸函数g的全局最优解所对应的那个θ值作为新的θ的值，我将其记为θ2。这个时候可以发现θ已经向上爬了。我们可以对上述步骤进行迭代，在针对θ2创建一个凸函数g2并得到θ3的值，我们可以得到如下图像：

上述步骤一直迭代进行，直到找到一个θ值能最大化这个最大似然函数。也就是找到‘山顶’为止。

数学推导EM算法
在理解EM的数学推导过程之前，首先要理解一个很重要的概念那就是Jensen’s inequality（琴生不等式）。
他可以被定义成：
这里的 f 是一个凸函数， $\lambda$ 是指的一种概率分布。那应该如何来理解这一个不等式呢？我们可以先假设目前就只有两个$\lambda$ 和两个$x_j$的情况就比较好理解了。
我们设${\lambda }_1=t$，那么就应该是${\lambda }_2=1-t$，因为一个概率分布中所以的概率相加应该是1嘛。 同时也假设我们有两个变量$x_1$和$x_2$然后将我们所假设好的变量带回到上述的不等式中我们可以得到以下式子：

$f((1-t)x_1+t{x}_2)⩾(1-t)f(x_1)+tf(x_2)$

这样乍一看还是觉得没有感觉。但是至少要意识到，不等式的左边是一条曲线，而右边是一条直线。那么我们要是把图画出来就会清晰很多了。假设这里的这个凸函数是长图里这个样子的：
这样子就很容易理解这个大于等于号是怎么来的了，如果在$x_1$和$x_2$之间随意取一点$x_T$那么这个点所对应的值永远小于那个凸函数,取到段直线两个端点的时候等号成立，正如下图所示：

那么上面这一切和EM算法有什么关系吗？
假设我们有一组独立的数据集{x1,x2,…},我们希望使用一个模型来拟合它，我们可以将这个模型写成下面这个似然函数，这里的 X 就像是我们在第二个投硬币例子里的投掷结果，θ就像是硬币A和B 出现正面的概率。

为什么要取log呢？我的理解是为了方便数学运算，而且取log是不会影响寻找最大值的。
但是上述公式没有把隐变量对结果的影响表达出来，所以我们想办法把隐变量塞进我们的模型中。为了方便理解我将公式拆开来记录。正如前文提到，对于一个投掷结果，它可能是硬币A产生的，也有可能是B产生的，所以我们应该把这个两种情况的概率加在一起:

这里的大Z指的就是隐变量有多少种情况。对于投硬币那个情况而言，就是大Z等于2。所以新的模型可以表达成：

做到这里再来看看当前这个模型离Jensen‘s inequality 还差点啥？？ 是不是只差一个概率分布了？？那在我们已知的条件中还有什么是概率分布呢？？结合抛硬币的例子来看，只有隐变量z的分布Q了呀！！那为了保持等式的值不变，我在等式的右边同时除以和乘以一个$Q(z_j)$得到以下这个等式：

现在上述这个等式就完全和Jensen’s inequality 的形式一致了 (log函数本身就是一个凸函数哦！) 。所以就可以表示成Jensen‘s inequality的形式：

上述不等式在什么情况况下取等号呢？在$log\frac{x,z_j;\theta }{Q(z_j}$等于一个常数的时候等号成立。

高斯混合模型
接下来会有一些"比较恐怖"数学公式，但是只要细心去理解每一个参数表达的含义，其实是容易理解的。

我们同样设想这样一个场景，我们有一组数据，知道这一组数据中的每一个元素都有可能属于多种高斯分布的某一种，然后需要找出每一个元素属于某一种概率分布的概率。结合前面的抛硬币问题，其实就是把硬币A和硬币B换成了多个高斯分布了。
要是我们把图画出来它长得就像下图：

这里可以看到灰色圈圈表示的是数据的每一个元素，都有可能属于三种高斯分布的一种，比如用红色标注出来的那一个点，他就属于三种分布。所以这样的聚类（分类）方式被称之为 Soft Decision (软分类)。因为它并不是像K-means算法那样，样本不是被分为A就是被分为B，soft decision 阐述的是样本属于某一种分布的可能性，如果某一个样本属于某一种的分布的概率为98% 属于其于两种分布的概率之和才2%，那不用多说，该样本就是属于第一种分布。

那么带着这样的思想，应该如何将这样一个复杂的情况进行建模呢？？
首先和抛硬币的问题一样，这个问题也存在着一个隐变量：就是我不知道到当前的这个样本是属于哪一种分布的，我们可以用 $z_j$来表示这个隐变量，它用来表示每一种分布。
这样一来的话，我要想知道样本$x_i$和$z_j$同时出现的概率，就可以用一下联合概率分布来表示了：

这个时候来观察一下$z_j$，它是不是一个变量，它的值是不确定的，有可能是三种分布中的任意一种，那么换一个角度来思考，这个$z_j$本身就是一个随机变量啊！那么随机变量一定是符合某一种概率分布的呀！我们在这里用${\varphi }_j$来表示$z_j$是某一种分布的概率，那么就一定有${\sum }_{i=1}^k{\varphi }_j=1$。

看到这里我们停下来思考一下，我们需要优化和更新的参数到底是什么呢？我们想让三种高斯分布尽可能正确的来拟合我们现在的数据，说明我们需要对高斯分布进行更新。那什么高斯分布的关键参数是什么？或者什么参数可以决定高斯分布的形状呢？？那一定是期望和方差(标准差)了！！如果说我们假设我们数据集中某一个灰色的圈圈是属于某一种高斯分布的话，那么这个点和这个高斯分布的联合概率分布就应该是服从高斯分布的哎！！
可以表达成这样:

然后如果我们的高斯分布更新了，那势必在没有收敛的前提之下${\mu }_j$和${\sigma }_j$肯定发生了变化，因为分布形状变了嘛。高斯分布更新了那是不是他被选取到的可能性也会发生变化，也就意味着${\varphi }_j$发生了变化。

所以我们可以照着上面的思路整理出E-step 和 M-step了，就和投硬币那个问题的思路很像。
E-step： 我们的分布不变也就是${\mu }_j$和${\sigma }_j$不发生改变，${\varphi }_j$更新
M-step: ${\varphi }_j$不发生改变，更新${\mu }_j$和${\sigma }_j$
类比投硬币的那个问题就很容易理解。

接下来用数学推导：
再次声明一下，我们要更新的参数是 ${\varphi }_j$，${\mu }_j$和 ${\sigma }_j$
首先整个问题也可以表示成一个似然函数：

上述式子表达的是在给定的${\varphi }_j$，${\mu }_j$和 ${\sigma }_j$的情况下，产生$x_i$的概率是多少。
同样的这个公式没有阐述隐变量对模型的影响，前面已经定义了$z_j$来表示，
将隐变量加入到模型中我们可以得到：

仔细看看上面这个式子，它是不是很像琴生不等式的左边呀！
如果我想最大化这个似然函数，那我可以试着写出琴生不等式的右边并让左右两侧取等，
那么回想我们投硬币的例子，我们要构建琴生不等式的右边是不是又差了一个概率分布呀！
那么我们可不可以把每一个点属于三种高斯分布的概率表示出来呢？
用$w_j$来表示概率值，Q来表示分布，那么我么可以得到：

上述式子表达的就是在知道样本 $\varphi$,$\sigma$,$\mu$的前提下出现第j个高斯分布的概率。
接下来我们就来构建这个琴生不等式，并将它取等号！
那么我们又可以得到:
那这个时候你就会发现了等式右侧的$P(x_i|z_j;\mu ,\sigma )$实际上是一个高斯分布！
它符合高斯分布的概率密度公式！若将其展开我们可以的到以下这个看着十分复杂的式子：

这就是我们的最终的最大似然函数了！
现在要来更新${\varphi }_j$,${\mu }_j$,${\sigma }_j$
那就需要对这个三个参数求其偏导数，偏导反应出来的是这些参数的变化对整个函数的影响。
我们只需要关注包含这些包含参数的部分，其余不参与的部分不用考虑。
那么我们可以得到：

要注意的是在用拉格朗日找$\sigma$的时候不要忘记他的约束条件是所有$\varphi$之和为1.

求到上述三个重要参数的表达式之后，就可以按照E, M两个步骤进行更新了！