第八课.EM算法的合理性与算法推导

EM算法背景

在硬币投掷和班级身高问题中，引入了隐变量问题，模型变成了混合模型，我们不能直接利用极大似然估计求出模型参数，通常需要用迭代法逐步尝试，这个迭代法正是EM算法的基础。

在极大似然估计中，我们采用计算最值的方法，将使得$p(x|\theta )$取得最大值的参数$\theta$作为估计值，比如一类简单的概率模型：高斯分布，他只有观测变量$x$，我们可以直接利用观测变量，基于极大似然，估计模型的参数。

对于一些模型，还有隐变量$z$，这一类变量不能观测到，这使得在模型中$p(x|\theta )$无法仅仅利用观测变量$x$直接通过求导等解析方法求出估计值$\theta$，只能通过迭代，不断逼近真实解。这个迭代的算法被总结为EM算法。

EM算法迭代的合理性

首先，给出EM算法的参数迭代公式：
$${\theta }^{t+1}=argma{x}_{\theta }{\int }_{z}log[p(x,z|\theta )]p(z|x,{\theta }^t)dz$$
在第$t$轮迭代中，可以利用参数${\theta }^t$，估计出${\theta }^{t+1}$；

但要思考，为什么每一轮迭代都能使似然函数$p(x|\theta )$的值不断增大，即问题描述为，对于任意轮数$t$，通过：
$${\theta }^{t+1}=argma{x}_{\theta }{\int }_{z}log[p(x,z|\theta )]p(z|x,{\theta }^t)dz$$
实现${\theta }^t$到${\theta }^{t+1}$后，使得$log[p(x|{\theta }^t)]\le log[p(x|{\theta }^{t+1})]$；

证明如下，首先利用贝叶斯公式得到观测变量$x$和隐变量$z$之间的概率关系：
$$p(x|\theta )p(z|x,\theta )=p(x,z|\theta )\Rightarrow log(p(x|\theta ))+log(p(z|x,\theta ))=log(p(x,z|\theta ))$$
因此，将隐变量引入到对数似然函数中有：
$$log(p(x|\theta ))=log(p(x,z|\theta ))-log(p(z|x,\theta ))$$
对等式两边同时求关于$p(z|x,{\theta }^t)$的期望，即积分：
$${\int }_{z}p(z|x,{\theta }^t)log(p(x|\theta ))dz={\int }_{z}p(z|x,{\theta }^t)log(p(x,z|\theta ))dz-{\int }_{z}p(z|x,{\theta }^t)log(p(z|x,\theta ))dz$$
对左边进行化简：
$${\int }_{z}p(z|x,{\theta }^t)log(p(x|\theta ))dz=log(p(x|\theta )){\int }_{z}p(z|x,{\theta }^t)dz=log(p(x|\theta ))\cdot 1=log(p(x|\theta ))$$
因此有：
$$log(p(x|\theta ))={\int }_{z}p(z|x,{\theta }^t)log(p(x,z|\theta ))dz-{\int }_{z}p(z|x,{\theta }^t)log(p(z|x,\theta ))dz$$
对于$p(x|{\theta }^{t+1})\ge p(x|{\theta }^t)$转化为验证：
$${\int }_{z}p(z|x,{\theta }^t)log(p(x,z|{\theta }^{t+1}))dz-{\int }_{z}p(z|x,{\theta }^t)log(p(z|x,{\theta }^{t+1}))dz$$$$\ge {\int }_{z}p(z|x,{\theta }^t)log(p(x,z|{\theta }^t))dz-{\int }_{z}p(z|x,{\theta }^t)log(p(z|x,{\theta }^t))dz$$
拆解为两个部分：

• 不等式1：
  $${\int }_{z}p(z|x,{\theta }^t)log(p(x,z|{\theta }^{t+1}))dz\ge {\int }_{z}p(z|x,{\theta }^t)log(p(x,z|{\theta }^t))dz$$
• 不等式2：
  $${\int }_{z}p(z|x,{\theta }^t)log(p(z|x,{\theta }^{t+1}))dz\le {\int }_{z}p(z|x,{\theta }^t)log(p(z|x,{\theta }^t))dz$$

对于不等式1，先回顾${\theta }^{t+1}$的计算：
$${\theta }^{t+1}=argma{x}_{\theta }{\int }_{z}log[p(x,z|\theta )]p(z|x,{\theta }^t)dz$$
即，令$\theta ={\theta }^{t+1}$，可以使得${\int }_{z}log[p(x,z|\theta )]p(z|x,{\theta }^t)dz$取到最大值，所以必然有：
$${\int }_{z}log[p(x,z|{\theta }^{t+1})]p(z|x,{\theta }^t)dz\ge {\int }_{z}log[p(x,z|{\theta }^t)]p(z|x,{\theta }^t)dz$$
对于不等式2，需要变形一下：
$${\int }_{z}p(z|x,{\theta }^t)log(p(z|x,{\theta }^t))dz-{\int }_{z}p(z|x,{\theta }^t)log(p(z|x,{\theta }^{t+1}))dz={\int }_{z}p(z|x,{\theta }^t)log(\frac{p(z|x,{\theta }^t)}{p(z|x,{\theta }^{t+1})})dz$$

---

现引入KL散度，即相对熵。设$P(x)$和$Q(x)$是随机变量$x$上的两个概率分布，在离散和连续随机变量的情形下，相对熵定义为：
$$KL(P||Q)=\sum P(x)log\frac{P(x)}{Q(x)}$$
$$KL(P||Q)=\int P(x)log\frac{P(x)}{Q(x)}dx$$
KL散度用于衡量两个分布之间的距离，因此，KL散度具有非负性：$KL(P||Q)\ge 0$，在两个分布完全相同时，等号成立。

---

根据KL散度的定义，容易发现：
$${\int }_{z}p(z|x,{\theta }^t)log(\frac{p(z|x,{\theta }^t)}{p(z|x,{\theta }^{t+1})})dz=KL(p(z|x,{\theta }^t)||p(z|x,{\theta }^{t+1}))\ge 0$$
此时，等式1和等式2均成立，因此证实，每轮迭代${\theta }^t\to {\theta }^{t+1}$确实可以让似然函数的值有所上升，当似然函数值收敛时，待估计的参数也就趋近于真实值（不保证全局最优）。至此，我们可以安心使用EM算法，下面将解释EM算法参数迭代公式的推导。

EM算法推导

EM算法的E步和M步

在前面的内容中，直接给出了EM算法的迭代公式，并证明了其迭代的合理性，现在需要了解其公式的由来，首先定义一些符号：

• $x$：观测数据，observed data
• $z$：非观测数据，隐变量，latent variable
• $(x,z)$：完整数据，complete data
• $\theta$：待估计参数

对于参数迭代公式：
$${\theta }^{t+1}=argma{x}_{\theta }{\int }_{z}log[p(x,z|\theta )]p(z|x,{\theta }^t)dz$$
其实包含了E步和M步：

• E步为积分的计算，即计算$log[p(x,z|\theta )]$关于条件概率$p(z|x,{\theta }^t)$的期望：
  $$E_{z|x,{\theta }^t}[log(p(x,z|\theta ))]={\int }_{z}log[p(x,z|\theta )]p(z|x,{\theta }^t)dz$$
• M步则是获取令期望取得最大值的$\theta$，并将其作为下一轮迭代的已知量：
  $${\theta }^{t+1}=argma{x}_{\theta }{E}_{z|x,{\theta }^t}[log(p(x,z|\theta ))]$$

分析EM算法的由来

首先，最初的极大似然估计目标为：
$${\theta }_{mle}=argma{x}_{\theta }log(p(x|\theta ))$$
根据贝叶斯公式有：
$$p(x|\theta )p(z|x,\theta )=p(x,z|\theta )\Rightarrow p(x|\theta )=\frac{p(x,z|\theta )}{p(z|x,\theta )}$$
则有：
$$log(p(x|\theta ))=log(p(x,z|\theta ))-log(p(z|x,\theta ))$$
此处引入一个关于$z$的分布$q(z)$，我们暂时不用关系$q(z)$具体是什么物理意义，代入式子：
$$log(p(x|\theta ))=log(\frac{p(x,z|\theta )}{q(z)})-log(\frac{p(z|x,\theta )}{q(z)})$$
此时，对等式左右两边同时求关于分布$q(z)$的期望，即计算积分：

对于式子左边有：
$${\int }_{z}q(z)log(p(x|\theta ))dz=log(p(x|\theta )){\int }_{z}q(z)dz=log(p(x|\theta ))$$
对于式子的右边有：
$${\int }_{z}q(z)log(\frac{p(x,z|\theta )}{q(z)})dz-{\int }_{z}q(z)log(\frac{p(z|x,\theta )}{q(z)})dz={\int }_{z}q(z)log(\frac{p(x,z|\theta )}{q(z)})dz+{\int }_{z}q(z)log(\frac{q(z)}{p(z|x,\theta )})dz$$
第一项被称为"证据下界"（Evidence lower bound，ELBO）：
$$ELBO={\int }_{z}q(z)log(\frac{p(x,z|\theta )}{q(z)})dz$$
第二项是KL散度：
$$KL(q(z)||p(z|x,\theta ))={\int }_{z}q(z)log(\frac{q(z)}{p(z|x,\theta )})dz\ge 0$$
当且仅当$q(z)=p(z|x,\theta )$时，$KL(q(z)||p(z|x,\theta ))=0$；

至此，式子可写作：
$$log(p(x|\theta ))=ELBO+KL(q(z)||p(z|x,\theta ))\ge ELBO$$
从以上形式可以看出，为什么ELBO被称为证据下界。ELBO成为了$log(p(x|\theta ))$取值的下界，所以可以得到一种思路：我们应当使得ELBO达到最大，从而间接地让$log(p(x|\theta ))$达到最大，即使用ELBO等效地控制$log(p(x|\theta ))$。但目前存在的问题是KL散度会破坏上面单纯的思路，因为$log(p(x|\theta ))$最终是ELBO与KL散度的总和。

这时，需要利用KL散度的性质，在第$t$轮迭代时，已知${\theta }^t$，要估计${\theta }^{t+1}$，我们把KL散度中概率$p(z|x,\theta )$的变量$\theta$固定，令$\theta ={\theta }^t$，同时，令$q(z)=p(z|x,{\theta }^t)$，从而使每一轮迭代时，KL散度刚好为0，此时有$log(p(x|\theta ))=ELBO$，取ELBO的极大值，就等效于取$log(p(x|\theta ))$的极大值：
$${\theta }^{t+1}=argma{x}_{\theta }log(p(x|\theta ))=argma{x}_{\theta }[{\int }_{z}p(z|x,{\theta }^t)log(\frac{p(x,z|\theta )}{p(z|x,{\theta }^t)})dz]$$$$=argma{x}_{\theta }{\int }_z[p(z|x,{\theta }^t)log(p(x,z|\theta ))-p(z|x,{\theta }^t)p(z|x,{\theta }^t)]dz$$
注意到，积分内的第二项$p(z|x,{\theta }^t)p(z|x,{\theta }^t)$已经不包含变量$\theta$，因此对$\theta$求极大值时可以舍去：
$${\theta }^{t+1}=argma{x}_{\theta }{\int }_z[p(z|x,{\theta }^t)log(p(x,z|\theta ))-p(z|x,{\theta }^t)p(z|x,{\theta }^t)]dz$$$$=argma{x}_{\theta }{\int }_z[p(z|x,{\theta }^t)log(p(x,z|\theta ))]dz$$
可见，这就得到了EM算法的迭代公式。

---

额外说明：

关于 $log(p(x|\theta ))=ELBO+KL(q(z)||p(z|x,\theta ))\ge ELBO$，有的文档采用jensen不等式来表达，并将$p(z|x,{\theta }^t)$，即$q(z)$称为$Q$函数（质量函数），这些其实都是EM算法推导的另一种解读

---