关于EM算法的理解
EM算法(Expectation-Maximization algorithm)是一种迭代算法，用于含有隐变量的概率模型参数的最大似然估计，或最大后验概率估计。与其他算法相比最大的区别在于，EM算法是个一般的方法，没有具体的模型。
最开始研读EM算法的时候，总是因为混合高斯模型的原因，没有理解EM算法的“一般”之处体现在哪里。后来才明白，混合高斯模型本身就是EM算法的一个重要应用。
我用一个例子来解释EM算法的一般之处。已知学校男女生身高各自独立的服从一个高斯分布，随机抽查了学校200个学生的身高，但是没有记录学生的性别。现在的问题是，如何通过这200个数据样本来估计男女生身高所对应的高斯分布的参数。需要解释的，每个样本记录是来自男生还是女生就是EM算法当中所谓的隐变量，对于两个高斯分布参数的估计，就是EM算法所求的结果。EM算法的一般之处体现在，男女生各自所对应的概率分布不一定非要是高斯分布，可以是其他任何可能的带有参数的分布，而EM算法的目的就是对这两个带有参数的分布进行参数估计。
最大似然估计
EM算法中，使用的参数估计的方法是最大似然估计。这部分内容主要参照《概率论与数理统计教程》。
定义
设总体的概率函数$p(x;\theta),\theta \in \Theta$，其中$\theta$是一个位置参数或几个位置参数组成的参数向量，$\Theta$是参数$\theta$可能的取值空间，$x_1,…,x_n$是来自该总体的样本，将样本的联合概率函数看成$\theta$的函数，用$L(\theta;x_1,…,x_n)$表示，简记为$L(\theta)$，
$$L(\theta)=L(\theta;x_1,…,x_n)=p(x_1;\theta)…p(x_n;\theta)=\prod_{i=1}^n p(x_i;\theta)$$
$L(\theta)$称为样本的似然函数。如果某估计量$\hat{\theta}=\hat{\theta}(x_1,…,x_n)$满足
$$L(\hat{\theta})=\max_{\theta \in \Theta}L(\theta)$$
则称$\hat{\theta}$是$\theta$的最大似然估计(Maximum Likelihook Estimate)。
对数似然函数
事实上，由于$\ln x$是$x$的单调函数，因此，在对数似然函数$\ln L(\theta)$上达到最大与使似然函数$L(\theta)$达到最大是等价的。在具体的计算中，人们通常更习惯于使用对数似然函数$\ln L(\theta)$除法来寻找$\theta$的最大似然估计。
Jensen不等式
Jensen不等式用于EM算法的推导过程，以下只给出其表述，不涉及证明。
如果$f$是凸函数，根据Jensen不等式得
$$f(\sum_{i=1}^n \lambda_ix_i)\leq \sum_{i=1}^n \lambda_if(x_i)\quad \lambda_i\geq 0\quad\sum_{i=1}^n\lambda_i=1$$
如果将$x_i$看做是随机变量，将$\lambda_i$看成是随机变量$x_i$出现的概率，那么$E[f(x)] \geq f(E(x))$。
特别的，如果$f$是严格凸函数，那么当且仅当$x$为常量时，等号成立。
对应的，如果$f$是凹函数，改变对应不等号的方向，也就是$E[f(x)] \leq f(E(x))$。
另外，关于Jensen不等式的图形化表述可以参照：EM-Algorithm。
EM算法
为了更好的表述，现在重新回到前面估计学生身高分布参数的例子中。但是与之前稍许有些不同的是，为了保证EM算法的一般性，我并不准备假设男女学生的身高服从高斯分布，而是假设他们服从某个含有参数$\theta$的任意分布。
一个比较完整的叙述请参考高斯混合模型及EM。
问题的描述
设采集学校男女生的身高样本$X={x_1,x_2,\cdots,x_N}$（但是并没有记录样本来自男生还是女生），用二元组$(x_i,z_i)$来表示每一项数据：$x_i$表示第$i$个学生的身高，$z_i$表示第$i$个学生究竟是男生（用1表示）还是女生（用0表示）。
同时假设，男女生的身高服从某个含有参数$\theta$的任意分布。目标是，通过样本数据求出参数$\theta$的最大似然估计。
过程的推导
根据联合分布的概率和，得到：$p(x_i;\theta) = \sum\limits_{z_i}p(x_i,z_i;\theta)$
利用对数极大似然估计函数，得到
\begin{align}
H(\theta)&= \sum_{i=1}^N\ln p(x_i;\theta) \\
&= \sum_{i=1}^N \ln \sum_{z_i} p(x_i,z_i;\theta)
\end{align}
是和极大似然估计中一样，我们需要使得$H(\theta)$达到最大值，最常规的做法就是对$H(\theta)$进行求导，令其导数为0。但是这和之前不含隐变量的函数不同，这个函数对$\theta$比较复杂。
所以，为了能解决这个为题，可以尝试做一些改变：令$Q_i$表示变量$z_i$的某种分布，且$Q_i$满足$\sum\limits_{z_i}Q_i(z_i)=1$，$Q_i(z_i) \ge 0$,利用这个$Q_i$，我们可以将关于$H(\theta)$等式改变为：
$$H(\theta) = \sum_{i=1}^N\ln \sum_{z_i}Q_i(z_i) \cdot \frac{p(x_i,z_i;\theta)}{Q_i(z_i)}$$
由于对数函数$\ln x$是个凹函数，故此处可以直接应用Jensen不等式：
\begin{align}
H(\theta) &= \sum_{i=1}^N\ln \sum_{z_i}Q_i(z_i) \cdot \frac{p(x_i,z_i;\theta)}{Q_i(z_i)} \\
&\ge \sum_{i=1}^N \sum_{z_i} Q_i(z_i) \ln \frac{p(x_i,z_i;\theta)}{Q_i(z_i)}
\end{align}
但是，到了这里你也许会问，我们明明要求的是$H(\theta)$的最大值，可是这里怎么会是一个”大于等于”呢？的确，目前我们所得到的确实是一个下界，而且这个下界和$Q_i$以及$p(x_i,z_i;\theta)$有关。但是在EM算法中，我们能够通过迭代，不断提高这个下界，来逐渐逼近$H(\theta)$。当不等式变为等式时，迭代结束。
关于逼近过程的图像化描述可以参照资料EM算法。
根据Jensen不等式等号成立的条件可得：
$$\frac{p(x_i,z_i;\theta)}{Q_i(z_i)}=c$$
又因为$\sum\limits_{z_i}Q_i(z_i)=1$，因此可以得$p(x_i;\theta)=c$。
重新将其带回，可以得到：
\begin{align}
Q_i(z_i) &= \frac{p(x_i,z_i;\theta)}{c}\\
&= \frac{p(x_i,z_i;\theta)}{p(x_i,;\theta)}\\
&= p(z_i\vert x_i;\theta)
\end{align}
这样一来，如果给定$x_i,\theta$，计算出$Q_i(z_i)$；知道$Q_i(z_i)$，又可以重新估计参数$\theta$。因此，如果我们开始时，赋给$\theta$一个初始值，通过这个过程就可以反复迭代直至收敛，这也就是EM算法的思想了。
实际上，总结EM算法的过程就是两步，E&M。
E-Step
对于每个$i$，估计隐变量$z_i$的分布：
$$Q_i(z_i)=p(z_i\vert x_i;\theta_t)$$
M-Step
重新估计参数$\theta_{t+1}$：
$$\theta_{t+1}=\arg\max_\theta \sum_{i=1}^N \sum_{z_i} Q_i(z_i) \ln \frac{p(x_i,z_i;\theta_t)}{Q_i(z_i)}$$
对数似然函数
当对数似然函数值或者参数收敛后，迭代停止
$$H(\theta)=\sum_{i=1}^N \sum_{z_i} Q_i(z_i) \ln \frac{p(x_i,z_i;\theta)}{Q_i(z_i)}$$
另外，关于EM算法收敛的证明参照The EM Algorithm。
高斯混合模型与EM算法
EM算法的一个重要应用就是对高斯混合模型进行参数估计。
高斯混合模型的定义
该部分的内容主要参考《统计学习方法》。
高斯混合模型是指具有如下形式的概率分布模型：
$$p(x|\theta)=\sum_{k=1}^K \alpha_k \phi(x|\theta_k)$$
其中，$\alpha_k$是系数，并且满足
$$\sum_{k=1}^K \alpha_k=1,\alpha_k \geq 0$$
$\phi(x|\theta_k)$是高斯分布函数，$\theta_k=(\mu_k,\sigma_k^2)$，并且有
$$\phi(x|\theta_k)=\frac{1}{\sqrt{2\pi\sigma_k^2}}\exp\left(-\frac{(x-\mu_k)^2}{2\sigma_k^2}\right)$$
称为第$k$个分模型。
形式化表述
这里直接给出混合高斯模型对应的E-Step和M-Step，此部分内容参考Mixtures of Gaussians and the EM algorithm
。
E-Step
对于每个$i$，估计隐变量$w_i^{(j)}$的分布：
\begin{align}
w_j^{(i)}&=p(z^{(i)}=j|x^{(i)};\alpha,\mu,\Sigma)\\
&=\frac{\alpha_j p(z^{(i)}=j|x^{(i)};\mu,\Sigma)}{\sum_{t=1}^K \alpha_t p(z^{(i)}=j|x^{(i)};\mu,\Sigma)}
\end{align}
M-Step
重新估计参数$\alpha_j,\mu_j,\Sigma_j$：
\begin{align}
\alpha_j&=\frac{1}{N} \sum_{i=1}^N w_j^{(i)}\\
\mu_j&=\frac{\sum_{i=1}^N w_j^{(i)}x^{(i)}}{\sum_{i=1}^N w_j^{(i)}}\\
\Sigma_j&=\frac{\sum_{i=1}^N w_j^{(i)}(x^{(i)}-\mu_j)(x^{(i)}-\mu_j)^T}{\sum_{i=1}^N w_j^{(i)}}\\
\end{align}
对数似然函数
\begin{align}
\ln l(\alpha,\mu,\Sigma)&=\sum_{i=1}^N \ln p(x^{(i)};\alpha,\mu,\Sigma)\\
&=\sum_{i=1}^N \ln \left( \sum_{j=1}^K \alpha_j p(x^{(i)}|z^{(i)};\mu,\Sigma) \right)
\end{align}
参考资料
[1] EM-Algorithm
[2] Mixtures of Gaussians and the EM algorithm