【课堂笔记】EM算法

背景

  EM算法（期望最大化算法，Expectation-Maximization Algorithm）是一种迭代优化算法，用于在含有隐变量的概率模型中估计最大似然参数。
  这是概括性的定义，下面我会解释其中的名词并用具体例子来引入EM算法。

极大似然估计

  先复习一下极大似然函数估计，我们假设数据满足某个分布（例如正态分布$N(\mu ,{\sigma }^2)$），但我们不知道其中的参数（$\mu ,\sigma$），于是我们需要从已知的数据中去拟合或估计出这些参数。
  进行极大似然估计的一般过程为（以正态分布为例）：
（1）确定概率模型
$$p(x_i;\mu ,{\sigma }^2)=\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp \left(-\frac{(x_i-\mu {)}^2}{2{\sigma }^2}\right)$$
（2）确定似然函数，并取负对数得到负对数似然
$$\begin{gathered} L(\mu ,{\sigma }^2)=\prod _{i=1}^{n}p(x_i;\mu ,{\sigma }^2)=\prod _{i=1}^n\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp \left(-\frac{(x_i-\mu {)}^2}{2{\sigma }^2}\right) \\ l(\mu ,{\sigma }^2)=-\ln L(\mu ,{\sigma }^2)=\frac{n}{2}\ln (2\pi {\sigma }^2)+\frac{1}{2{\sigma }^2}\sum _{i=1}^n(x_i-\mu {)}^2 \end{gathered}$$
（3）计算极值点，$\hat{\theta }=(\mu ,{\sigma }^2)$是待估计的参数
$$\hat{\theta }=\arg \underset{\mu ,{\sigma }^2}{\min }l(\mu ,{\sigma }^2)$$
  如果方程简单，可以解析求解，即导数为0，得到似然方程。
  如果方程复杂，可能需要数值方法（如梯度下降、牛顿法）。

隐变量

  什么是隐变量？考虑以下情景：
  现在有一个蛋糕集，里面有巧克力蛋糕和草莓蛋糕，用$k_i=0,1$来表示。假设它们分别满足分布$N({\mu }_1,{\sigma }_1^2)$和$N({\mu }_2,{\sigma }_2^2)$，然而，我们并不知道某个蛋糕具体是巧克力做的还是草莓做的。换句话说，我们需要估计出所有的$k_i$以及参数${\mu }_1,{\sigma }_1,{\mu }_2,{\sigma }_2$。这里$k_i$就是隐变量，即无法直接观测的变量。
  这里就产生了一个“鸡生蛋，蛋生鸡”问题：
（1）要想估计模型参数${\mu }_1,{\sigma }_1,{\mu }_2,{\sigma }_2$，需要知道每个样本的类别$k_i$
（2）要想确认样本的类别$k_i$需要事先知道模型的参数。

高斯混合模型

  事实上，我们可以用高斯混合模型来同时表示两个正态分布模型：
$$p(x_i|{\mu }_1,{\sigma }_1^2,{\mu }_2,{\sigma }_2^2,{\pi }_0,{\pi }_1)={\pi }_0\cdot \frac{1}{\sqrt{2\pi {\sigma }_1^2}}\exp \left(-\frac{(x_i-{\mu }_1{)}^2}{2{\sigma }_1^2}\right)+{\pi }_1\cdot \frac{1}{\sqrt{2\pi {\sigma }_2^2}}\exp \left(-\frac{(x_i-{\mu }_2{)}^2}{2{\sigma }_2^2}\right)$$
  这里引入了两个新的参数${\pi }_0,{\pi }_1$满足${\pi }_0+{\pi }_1=1$，分别表示巧克力蛋糕和草莓蛋糕的占比。注意由于我们不知道样本的类别，这两个参数也是需要估计的。
  看起来这样干很美好，只要把刚才最大似然估计的步骤来一遍就好了。然而在计算负对数似然的时候你就会卡住了。先前取对数的时候可以很容易把$\log \prod p(x_i)$ 变成$\sum \log p(x_i)$，然后$\log$和$p(x_i)$中的$\exp$消掉了。但这里不行。然后后面不管是求导还是梯度下降都完蛋了。

EM算法

  EM算法通过另一个思路解决了这个问题，那就是互相猜！
  这里先讲一下要算点啥，后面再讲来路。
（1）设置初始参数$\theta =({\mu }_1,{\sigma }_1^2,{\mu }_2,{\sigma }_2^2,{\pi }_0,{\pi }_1)$
（2）E-Step：根据目前的参数${\mu }_1,{\sigma }_1^2,{\mu }_2,{\sigma }_2^2$，分别计算每个样本属于草莓蛋糕的概率和属于巧克力蛋糕的概率。
似然：
$$\begin{gathered} p(x_i|k_i=0)=\frac{1}{\sqrt{2\pi {\sigma }_1^2}}\exp \left(-\frac{(x_i-{\mu }_1{)}^2}{2{\sigma }_1^2}\right) \\ p(x_i|k_i=1)=\frac{1}{\sqrt{2\pi {\sigma }_2^2}}\exp \left(-\frac{(x_i-{\mu }_2{)}^2}{2{\sigma }_2^2}\right) \end{gathered}$$
先验概率：
$$P(k_i=0)={\pi }_0,P(k_i=1)={\pi }_1$$
边缘概率：
$$\begin{array}{rl}p(x_i)& =p(x_i|k_i=0)P(k_i=0)+p(x_i|k_i=1)P(k_i=1)\\ & ={\pi }_0\cdot \frac{1}{\sqrt{2\pi {\sigma }_1^2}}\exp \left(-\frac{(x_i-{\mu }_1{)}^2}{2{\sigma }_1^2}\right)+{\pi }_1\cdot \frac{1}{\sqrt{2\pi {\sigma }_2^2}}\exp \left(-\frac{(x_i-{\mu }_2{)}^2}{2{\sigma }_2^2}\right)\end{array}$$
然后通过贝叶斯定理计算后验概率：
$$P(k_i=k|x_i)=\frac{p(x_i|k_i=k)P(k_i=k)}{p(x_i)}$$
（3）M-Step：这一步最大化期望对数似然，并更新参数。E 步计算隐变量期望后得到的优化目标为：
$$Q(\theta ,{\theta }^{(t)})=\sum _{i=1}^n\sum _{k=0}^{1}P(k_i=k|x_i,{\theta }^{(t)})\ln [p(x_i,k_i=k|\theta )]$$
  其中$\theta$是待更新参数，${\theta }^{(t)}$是当前迭代（第$t$轮）的参数，来自上一次估计。$P(k_i=k|x_i,{\theta }^{(t)})$是E-Step得到的后验概率。$p(x_i,k_i=k|\theta )$是联合概率。
  这个式子的来源和含义一会儿再讲。
  然后通过求偏导的方式更新参数：
$$\begin{gathered} {\pi }_k=\frac{1}{n}\sum _{i=1}^{n}P(k_i=k|x_i,{\theta }^{(t)}) \\ {\mu }_k=\frac{{\sum }_{i=1}^{n}P(k_i=k|x_i,{\theta }^{(t)})x_i}{{\sum }_{i=1}^{n}P(k_i=k|x_i,{\theta }^{(t)})} \\ {\sigma }_k^2=\frac{{\sum }_{i=1}^{n}P(k_i=k|x_i,{\theta }^{(t)})(x_i-{\mu }_k{)}^2}{{\sum }_{i=1}^{n}P(k_i=k|x_i,{\theta }^{(t)})} \end{gathered}$$
（4）返回步骤（2），重复E-Step，M-Step循环直到稳定。

合理性分析

  问题来了，为什么这么做可以拟合出最优状态，并且能保证算法的收敛性呢？
  回忆一开始最大似然估计解决不了的式子：
$$\begin{gathered} p(x_i|\theta )={\pi }_0\cdot \frac{1}{\sqrt{2\pi {\sigma }_1^2}}\exp \left(-\frac{(x_i-{\mu }_1{)}^2}{2{\sigma }_1^2}\right)+{\pi }_1\cdot \frac{1}{\sqrt{2\pi {\sigma }_2^2}}\exp \left(-\frac{(x_i-{\mu }_2{)}^2}{2{\sigma }_2^2}\right) \\ l(\theta )=-\ln \prod _{i=1}^{n}p(x_i|\theta ) \end{gathered}$$
  这个似然的问题在于无法直接优化，是因为求导时涉及 log(sum) 的形式，梯度不好算。于是我们引入一个任意$\mathbf{k}$的概率分布，并使用琴生不等式变化一下：
$$\ln p(\mathbf{x}\mid \theta )=\ln \sum _{\mathbf{k}}q(\mathbf{k})\cdot \frac{p(\mathbf{x},\mathbf{k}\mid \theta )}{q(\mathbf{k})}\ge \sum _{\mathbf{k}}q(\mathbf{k})\ln \frac{p(\mathbf{x},\mathbf{k}\mid \theta )}{q(\mathbf{k})}:=\mathcal{L}(q,\theta )$$
  于是我们得到了一个下界$\mathcal{L}$，这个下界被称为 Evidence Lower Bound（ELBO）。
  现在我们来研究一下这个下界，将它分解一下：
$$\begin{array}{rl}\mathcal{L}(q,\theta )& =\sum _{\mathbf{k}}q(\mathbf{k})\ln \frac{p(\mathbf{x},\mathbf{k}\mid \theta )}{q(\mathbf{k})}=\sum _{\mathbf{k}}q(\mathbf{k})\ln \frac{p(\mathbf{k}\mid \mathbf{x},\theta )p(\mathbf{x}\mid \theta )}{q(\mathbf{k})}\\ & =\sum _{k}q(\mathbf{k})\ln \frac{p(\mathbf{k}\mid \mathbf{x},\theta )}{q(\mathbf{k})}+\ln p(\mathbf{x}\mid \theta )\cdot \sum _{k}q(\mathbf{k})\end{array}$$
  注意到最后一项求和是$1$（概率分布的性质），于是：
$$\begin{array}{rl}\mathcal{L}(q,\theta )& =\sum _{k}q(\mathbf{k})\ln \frac{p(\mathbf{k}\mid \mathbf{x},\theta )}{q(\mathbf{k})}+\ln p(\mathbf{x}\mid \theta )\\ & =-\text{KL}(q(\mathbf{k})\Vert p(\mathbf{k}\mid \mathbf{x},\theta ))+\ln p(\mathbf{x}\mid \theta )\end{array}$$
  这里$\text{KL}$被称为KL散度，衡量你选的隐变量分布 $q(\mathbf{k})$离真实后验$p(\mathbf{k}\mid \mathbf{x},\theta )$ 有多远。没学过也不要紧，看上面的式子就可以了。
  在EM算法中，实际上干了以下两件事：
（1）设当前参数为${\theta }^{(t)}$，E步取$q(\mathbf{k})=p(\mathbf{k}|\mathbf{x},{\theta }^{(t)})$，此时$\mathcal{L}(q,\theta )$对$q$来说取到最大值，即KL散度为0。
（2）M步最大化$\theta$上$\mathcal{L}(q,\theta )$的值
$${\theta }^{(t+1)}=\arg \underset{\theta }{\max }\mathcal{L}(q,\theta )$$
  再看一下此时的$\mathcal{L}(q,\theta )$
$$\begin{array}{rl}\mathcal{L}(q,\theta )& =\sum _{k}p(k|\mathbf{x},{\theta }^{(t)})\ln \frac{p(\mathbf{x},k|\theta )}{p(k|\mathbf{x},{\theta }^{(t)})}\\ & =\sum _{k}p(k|\mathbf{x},{\theta }^{(t)})\ln (p(\mathbf{x},k|\theta ))-\sum _{k}p(k|\mathbf{x},{\theta }^{(t)})\ln (p(k|\mathbf{x},{\theta }^{(t)}))\\ & ={\mathbb{E}}_{k\sim p(k|\mathbf{x},{\theta }^{(t)})}[\ln p(\mathbf{x},k|\theta )]-{\mathbb{E}}_{k\sim p(k|\mathbf{x},{\theta }^{(t)})}[\ln p(k|\mathbf{x},{\theta }^{(t)})]\end{array}$$
  我们现在要优化$\mathcal{L}(q,\theta )$在$\theta$上的值，而后面那一项与$\theta$无关，于是我们只需要优化第一项，也就是前文直接给出的优化目标$Q(\theta ,{\theta }^{(t)})$
$$\begin{array}{rl}Q(\theta ,{\theta }^{(t)})& ={\mathbb{E}}_{k\sim p(k|\mathbf{x},{\theta }^{(t)})}[\ln p(\mathbf{x},k|\theta )]\\ & =\sum _{k}p(k|\mathbf{x},{\theta }^{(t)})\ln (p(\mathbf{x},k|\theta ))\end{array}$$
  再使用样本$x_i$和隐变量$k_i$的独立性，有：
$$\begin{gathered} p(\mathbf{x},\mathbf{k}\mid \theta )=\underset{i=1}{\prod ^n}p(x_i,k_i\mid \theta ) \\ p(\mathbf{k}\mid \mathbf{x},{\theta }^{(t)})=\underset{i=1}{\prod ^n}p(k_i\mid x_i,{\theta }^{(t)}) \end{gathered}$$
  于是：
$$Q(\theta ,{\theta }^{(t)})=\sum _{i=1}^n\sum _{k=0}^{1}p(k_i=k|x_i,{\theta }^{(t)})\ln [p(x_i,k_i=k|\theta )]$$
  这就是最终M步需要优化的目标。

相关好文章

【机器学习】EM算法详细推导和讲解
如何通俗理解EM算法