EM Algorithm

EM 算法是极大似然法的推广，用于解决存在隐变量（hidden variables / latent factors）的参数估计问题。

1 EM 算法

1.1 问题描述

极大似然法是最常用的参数估计方法之一。设观测变量为 $x$，模型参数为 $\theta$，则极大似然法通过最大化对数似然 $\log p(x;\theta )$ 来求解最优的参数 $\theta$. 然而在一些问题中，观测变量 $x$ 依赖于隐变量 $z$，此时根据全概率公式有： $$p(x;\theta )=\int p(x,z;\theta )\mathrm{d}z$$ 于是对数似然为： $$L(\theta )=\log p(x;\theta )=\log (\int p(x,z;\theta )\mathrm{d}z)$$ 如果仍然使用极大似然法，我们会发现 $L(\theta )$ 的导数将变得非常复杂，要优化的参数之间无法分离，导致无法写出封闭形式的解。这时就需要用到 EM 算法了。

1.2 理论推导

为了解决上述问题，我们为隐变量 $z$ 引入一个概率分布 $q(z)$，则有： $$L(\theta )=\log (\int p(x,z;\theta )\mathrm{d}z)=\log (\int q(z)\frac{p(x,z;\theta )}{q(z)}\mathrm{d}z)\ge \int q(z)\log \frac{p(x,z;\theta )}{q(z)}\mathrm{d}z=J(\theta ,q)$$ 其中应用了 Jensen 不等式。于是我们找到了“一族” $L(\theta )$ 的下界 $J(\theta ,q)$，称为证据下界 (Evidence Lower BOund, ELBO)。“一族”的意思是通过改变 $q$，我们能够得到不同的下界函数。特别地，根据 Jensen 不等式的取等条件可知，当 $q(z)=p(z|x;\theta )$ 时等号成立，我们可以代入 $q^{\ast }(z)=p(z|x;\theta )$ 来验证这一点： $$J(\theta ,q^{\ast })=\int p(z|x;\theta )\log \frac{p(x,z;\theta )}{p(z|x;\theta )}\mathrm{d}z=\int p(z|x;\theta )\log p(x;\theta )\mathrm{d}z=\log p(x;\theta )=L(\theta )$$

不等式能取到等号说明 $J(\theta ,q)$ 不仅是 $L(\theta )$ 的下界，还是紧的下界，所以我们可以通过优化 $J(\theta ,q)$ 达到优化 $L(\theta )$ 的目的，即： $$\underset{\theta }{max}L(\theta )⟺\underset{\theta ,q}{max}J(\theta ,q)$$ 代价是原本只需要优化一个变量 $\theta$，现在要优化两个变量 $\theta$ 和 $q$. 同时优化二者比较困难，因此 EM 算法采用交替迭代的方式优化。具体而言，首先初始化 ${\theta }^{(0)}$，然后交替执行 E-step 和 M-step. 在第 $t$ 轮中：

• E-step. 固定 $\theta$ 优化 $q$： $$q^{(t+1)}=\arg \underset{q}{max}J({\theta }^{(t)},q)$$ 根据前文的讨论，这一步的最优解就是取 $q^{(t+1)}(z)=p(z|x;{\theta }^{(t)})$，以使得 $J({\theta }^{(t)},q^{(t+1)})=L({\theta }^{(t)})$.

• M-step. 固定 $q$ 优化 $\theta$： $${\theta }^{(t+1)}=\arg \underset{\theta }{max}J(\theta ,q^{(t+1)})$$ 注意到该优化目标可以拆出一个与 $\theta$ 无关的常数项： $$\begin{array}{rl}J(\theta ,q^{(t+1)})& =\int q^{(t+1)}(z)\log p(x,z;\theta )\mathrm{d}z-\int q^{(t+1)}(z)\log q^{(t+1)}(z)\mathrm{d}z\\ & =\underset{Q(\theta ,{\theta }^{(t)})}{\underbrace{\int p(z|x;{\theta }^{(t)})\log p(x,z;\theta )\mathrm{d}z}}+\underset{\text{constant}}{\underbrace{H(q^{(t+1)}(z))}}\end{array}$$ 因此这一步等价于优化 $Q(\theta ,{\theta }^{(t)})$： $${\theta }^{(t+1)}=\arg \underset{\theta }{max}Q(\theta ,{\theta }^{(t)})$$

如此经过一轮 E-step 和 M-step，我们有： $$L({\theta }^{(t+1)})\ge J({\theta }^{(t+1)},q^{(t+1)})\ge J({\theta }^{(t)},q^{(t+1)})=L({\theta }^{(t)})$$ 因此 $L(\theta )$ 的确得到了优化。下图可视化了这样一轮迭代过程：

摘取自PRML，符号略有不同

可以看见，红色线条是要优化的目标，蓝色线条是在当前 $$ 下通过 E-step 找到的下界，M-step 找到了蓝色线条的最高点 $$，然后下一轮的 E-step 将改变 $$ 使得下界在 $$ 处与红色线条相等，即得到绿色线条……如此迭代即可。

1.3 进一步分析

在上一节中，我们推出了 $$，那他俩之间究竟差了个什么呢？ $$ 啊！原来是引入的分布 $$ 与后验分布 $$ 之间的 KL 散度。之前推导过程中使用的 Jensen 不等式，其实对应着这里的 KL 散度非负。事实上，我们能避开 Jensen 不等式直接推导出 $$，过程如下所示：

$$ 观察上式，可以看见 $$ 被划分成了两部分——ELBO 和 KL 散度。基于这个角度，我们重新审视 EM 算法——在 E-step 中，我们固定 $$ 并取 $$，使得 KL 项为零，也就是让 $$；在 M-step 中，我们固定 $$ 优化 $$，进而优化了 $$；这时由于 $$ 变了，KL 散度不再是零，我们便继续下一轮迭代，如下图所示：

摘自PRML，符号略有不同

1.4 算法步骤

综上所述，EM 算法的步骤如下：

1. 随机初始化 $$；

2. E-step. 给定 $$，求隐变量的后验分布 $$ 并计算优化目标： $$

3. M-step. 优化 $$，得到 $$： $$

4. 迭代执行 2、3 步直至收敛。

1.5 收敛性

EM 算法一定收敛吗？这里的收敛其实包含两层意思：$$ 收敛和 $$ 收敛，前者并不能蕴含后者。

1. $$ 收敛：上文已经说明了 $$ 在优化过程中单调不减，又它显然有上界，根据单调有界定理，$$ 必然收敛。
2. $$ 收敛：见参考文献^[5]。

需要注意的是，虽然 EM 算法一定收敛，但是可能收敛到局部最优解，收敛结果取决于初始化的情况。

1.6 关于独立重复试验的注解

在上面的推导中，我们把所有的观测结果记作 $$，这使得结论具有普适性（或者，你也可以认为我们只对单个观测结果做了推导）。实践中我们往往做的是独立重复试验，有 $$ 个独立同分布的样本 $$，各样本对应各自的隐变量 $$. 基于独立性，有：

$$ 如果把它们直接代入 $$ 之中，会得到非常丑陋的结果： $$ 这是因为我们代入得太晚——在 Jensen 不等式把对数和求和的顺序交换后，再代入就很难化简式子了（当然如果你头够铁，也不是不能化简，见参考资料^[9][10]）。但如果我们在最开始就代入： $$ 继续推导，会发现只需要在相应地方加上求和 $$，并给 $$ 和 $$ 加上下标 $$ 即可，整个推导过程没有什么大的变化。最终，优化目标写作： $$

2 例子：双硬币问题

假设有两枚硬币 A,B，它们正面朝上的概率未知，分别记作 $$. 为了估计 $$​，考虑做如下试验：首先等概率地随机选取一枚硬币，然后抛掷这枚硬币 10 次。假设我们一共做了 5 次试验，得到的结果如下：

试验 id	选取的硬币	试验结果
1	B	5 H, 5 T
2	A	9 H, 1 T
3	A	8 H, 2 T
4	B	4 H, 6 T
5	A	7 H, 3 T

可以看见，综合 5 次试验，硬币 A 一共 24 次正面朝上、6 次反面朝上，硬币 B 一共 9 次正面朝上、11 次反面朝上，所以 $$ 的估计值为： $$ 值得注意的是，虽然这个结果直观上很容易理解，但本质上是极大似然估计。其对数似然为： $$ 注意 $$ 和 $$ 可以分开考虑，因此分别求导并令为零，有： $$ 现在考虑一个更有挑战性的问题——我们并不知道每一次试验选的是哪枚硬币，只知道投掷结果，即：

试验 id	选取的硬币	试验结果
1	$$	5 H, 5 T $$
2	$$	9 H, 1 T $$
3	$$	8 H, 2 T $$
4	$$	4 H, 6 T $$
5	$$	7 H, 3 T $$

其中 $$ 表示第 $$ 次试验中正面朝上的次数，$$ 表示选取的硬币，是问题中的隐变量。应用 EM 算法：

1. 随机初始化 $$；

2. E-step. 求隐变量的后验分布 $$，以第一次试验为例： $$ 为了书写方便，以下将 $$ 和 $$ 分别简记为 $$ 和 $$. 计算优化目标 $$： $$

3. M-step. 优化上述目标，求偏导并令为零： $$ 解得： $$ $$ 同理求解。

4. 迭代执行 2、3 步直至收敛。

3 应用：高斯混合模型 (GMM)

设有 $$ 个 $$ 维高斯分布： $$ 假设样本 $$ 是由如下过程产生的：首先随机选取一个高斯分布（第 $$ 个高斯分布被选中的概率是 $$，其中 $$），再从被选中的高斯分布中采样。试估计这 $$ 个高斯分布的参数 $$.

与双硬币问题非常相似，这个问题中我们不知道每次采样究竟选取的是哪个高斯分布，所以可以定义隐变量 $$ 表示第 $$ 个样本来自哪个高斯分布。使用 EM 算法求解：

E-step. 给定 $$，隐变量的后验分布为： $$

为书写方便，以下将其简记为 $$，注意这是一个常数。计算优化目标 $$： $$ M-step（大计算量警告）. 首先解 $$：结合条件 $$，去掉 $$ 中与 $$ 无关的部分，优化问题化作： $$ 引入 Lagrange 乘子： $$ 求偏导并令为零： $$ 解得： $$ 下面解 $$：去掉 $$ 中与 $$​ 无关的部分，优化目标化作： $$ 求导令为零： $$

解得： $$ 最后解 $$：去掉 $$ 中与 $$ 无关的部分，优化目标化作： $$ 求导令为零： $$ 解得： $$ 总结一下，用 EM 算法求解 GMM 模型的步骤如下：

1. 随机初始化模型参数 $$.

2. E-step. 计算隐变量的概率分布： $$

3. M-step. 更新参数： $$

4. 迭代执行第 2、3 步直至收敛。

4 一点注记

EM 算法在 E-step 要求我们能表达出隐变量的后验概率 $$，从上面两个例子可以看出，这个后验概率可以通过贝叶斯公式计算： $$ 所以我们并没有避开计算似然 $$. 因此要明晰的一点是，阻碍我们直接使用极大似然法的不是计算似然，而是优化似然，即解令对数似然偏导为零后得到的那个方程组。

这自然带来了一个问题，如果似然本身就是不可解 (intractable) 的，那么 EM 算法就做不下去了。例如，在 GMM 中，隐变量取值是离散的，表示样本来自于哪个高斯分布。倘若我们把高斯分布的数量不断变大乃至无穷，那么隐变量取值就变成连续的，后验概率的分母（似然）变成了一个积分，无法通过枚举计算。为了解决这个问题，人们提出了变分推断的概念，我们在后面的文章中学习。

参考资料

1. Do, Chuong B., and Serafim Batzoglou. What is the expectation maximization algorithm?. Nature biotechnology 26, no. 8 (2008): 897-899. ↩︎
2. Borman, Sean. The expectation maximization algorithm-a short tutorial. https://www.lri.fr/~sebag/COURS/EM_algorithm.pdf ↩︎
3. Tengyu Ma and Andrew Ng. CS229 Lecture notes: The EM algorithm. https://cs229.stanford.edu/notes2021fall/cs229-notes8.pdf ↩︎
4. Bishop, Christopher. Pattern recognition and machine learning. Springer google schola 2 (2006): 531-537. ↩︎
5. Wu, CF Jeff. On the convergence properties of the EM algorithm. The Annals of statistics (1983): 95-103. ↩︎
6. 李航.统计学习方法 ↩︎
7. 人人都懂EM算法 - August的文章 - 知乎 https://zhuanlan.zhihu.com/p/36331115 ↩︎
8. 【机器学习】EM——期望最大（非常详细） - 阿泽的文章 - 知乎 https://zhuanlan.zhihu.com/p/78311644 ↩︎
9. 高斯混合模型（GMM）推导及实现 - 永远在你身后的文章 - 知乎 https://zhuanlan.zhihu.com/p/85338773 ↩︎
10. 机器学习-白板推导系列(十一)-高斯混合模型GMM（Gaussian Mixture Model）https://www.bilibili.com/video/BV13b411w7Xj?p=3&share_source=copy_web&vd_source=a43b4442e295a96065c7ae919b4866d3 ↩︎