从两个例子理解EM算法

本文是作者对EM算法学习的笔记，从EM算法出发介绍EM算法，为了更好理解，用两个应用EM算法求解的例子进一步解释EM的应用。

---

EM算法

EM算法引入

EM算法，指的是最大期望算法（Expectation Maximization Algorithm，期望最大化算法），是一种迭代算法，在统计学中被用于寻找，依赖于不可观察的隐性变量的概率模型中，参数的最大似然估计。基本思想是首先随机取一个值去初始化待估计的参数值，然后不断迭代寻找更优的参数使得其似然函数比原来的似然函数大。

• EM算法当做最大似然估计的拓展，解决难以给出解析解（模型中存在隐变量）的最大似然估计（MLE）问题
• 在算法中加入隐变量的思想可以类比为几何题中加入一条辅助线的做法。

假定有训练集{x(1),x(2),...x(m)$x^{(1)},x^{(2)},...x^{(m)}$}，包含m$m$个独立样本，希望从中找到该组数据的模型p(x,z)$p(x,z)$的参数。
对数似然函数表达如下：

在表达式中因为存在隐变量，直接找到参数估计比较困难，所以我们通过EM算法迭代求解下界的最大值，直到收敛。

我们通过以下的图片来解释这一过程：

图片上的紫色部分是我们的目标模型p(x|θ)$p(x|\theta )$曲线，该模型比较复杂，难以直接求解其解析解，为了消除隐变量z$z$带来的影响，我们可以得到一个不包含的z$z$的模型r(x|θ)$r(x|\theta )$(该函数是我们自己选定的，因此最大值可求解)， 同时满足条件r(x|θ)≤p(x|θ)$r(x|\theta )\le p(x|\theta )$。

• 我们先取一个θ1${\theta }_1$，使得r(x|θ1)=p(x|θ1)$r(x|{\theta }_1)=p(x|{\theta }_1)$(如绿线所示)，然后再对此时的r$r$求其最大值，得到极值点θ2${\theta }_2$,实现参数的更新。
• 不断重复以上过程，在更新过程中始终满足r≤p$r\le p$直到收敛。

从以上过程来看，EM算法的核心就是如何找到这个r$r$，即p$p$的下界函数。

这个下界函数有多种方法理解，我们从Jensen不等式的角度来理解。

上述等号成立的条件是
p(x(i),z(i);θ)Qi(z(i))=c$\frac{p(x^{(i)},z^{(i)};\theta )}{Q_i(z^{(i)})}=c$, ∑zQi(z(i))=1$\sum _z{Q}_i(z^{(i)})=1$,所以：

最终框架如下：

EM推导高斯混合模型

高斯混合模型GMM

设有随机变量X$X$, 则高斯混合模型可以用p(x)=∑Kπk(x|μk,Σk)$p(x)=\sum _K{\pi }_k\mathcal{N}(x|{\mu }_k,{\mathrm{\Sigma }}_k)$，其中(x|μk,Σk)$\mathcal{N}(x|{\mu }_k,{\mathrm{\Sigma }}_k)$表示混合模型中的第k$k$个分量πk${\pi }_k$表示混合系数，满足
∑kπk=1，0≤πk≤1$\sum _k{\pi }_k=1，0\le {\pi }_k\le 1$。

我们知道高斯函数的概率分布为f(x)=1(√2π)σexp(−(x−μ)22σ2)$f(x)=\frac{1}{\sqrt{(}2\pi )\sigma }exp(-\frac{(x-\mu {)}^2}{2{\sigma }^2})$, 在混合高斯分布中待估计变量就包括了μ,σ,π$\mu ,\sigma ,\pi$。

对数似然函数为lμ,Σ,π=∑Ni=1log(∑Kk=1)πk(xi|μk,Σk))$l_{\mu ,\mathrm{\Sigma },\pi }=\sum _{i=1}^{N}log(\sum _{k=1}^K){\pi }_k\mathcal{N}(x_i|{\mu }_k,{\mathrm{\Sigma }}_k))$

EM 推导过程

第一步：估算数据来自于哪个组分，即估计每一个组分生成的概率，对每个样本xi$x_i$，它由第k$k$个组份生成的概率可以记作：γ(i,k)=πk(xi|μk,Σk)∑jπj(xi|μj,Σj)$\gamma (i,k)=\frac{{\pi }_k\mathcal{N}(x_i|{\mu }_k,{\mathrm{\Sigma }}_k)}{\sum _j{\pi }_j\mathcal{N}(x_i|{\mu }_j,{\mathrm{\Sigma }}_j)}$

第二步：估计每个组份的参数

E-step: 在给定了样本和每个高斯分布的参数以及组份的分布函数的情况下

w(i)j=Qi(z(i)=j)=p(z(i))=j|x(i);ϕ,μ,Σ)$w_j^{(i)}=Q_i(z^{(i)}=j)=p(z^{(i)})=j|x^{(i)};\varphi ,\mu ,\mathrm{\Sigma })$

M-step：将多项式分布和高斯分布的参数带入：
∑mi=1∑z(i)Qi(z(i))logp(x(i),z(i);ϕ,μ,Σ)Qi(z(i))$\sum _{i=1}^m\sum _{z^{(i)}}{Q}_i(z^{(i)})log\frac{p(x^{(i)},z^{(i)};\varphi ,\mu ,\mathrm{\Sigma })}{Q_i(z^{(i)})}$
=∑mi=1∑kj=1Qi(z(i)=j)logp(x(i)|z(i)=j;ϕ,μ,Σ)p(z(i)=j;ϕ)Qi(z(i))$=\sum _{i=1}^m\sum _{j=1}^k{Q}_i(z^{(i)}=j)log\frac{p(x^{(i)}|z^{(i)}=j;\varphi ,\mu ,\mathrm{\Sigma })p(z^{(i)}=j;\varphi )}{Q_i(z^{(i)})}$
∑mi=1∑kj=1w(i)jlog1(2π)n2|Σj|(12)exp(−12(x(i)−μj)TΣ−1j(x(i)−μj))ϕjw(i)j$\sum _{i=1}^m\sum _{j=1}^k{w}_j^{(i)}log\frac{\frac{1}{(2\pi {)}^{\frac{n}{2}}|{\mathrm{\Sigma }}_j{|}^{(\frac{1}{2})}}exp(-\frac{1}{2}(x^{(i)}-{\mu }_j{)}^T{\mathrm{\Sigma }}_j^{-1}(x^{(i)}-{\mu }_j)){\varphi }_j}{w_j^{(i)}}$

分别对其中的未知参数求偏导数：

• 对均值求偏导
  ∇uj∑mi=1∑kj=1w(i)jlog1(2π)n2|Σj|(12)exp(−12(x(i)−μj)TΣ−1j(x(i)−μj))ϕjw(i)j${\mathrm{\nabla }}_{u_j}\sum _{i=1}^m\sum _{j=1}^k{w}_j^{(i)}log\frac{\frac{1}{(2\pi {)}^{\frac{n}{2}}|{\mathrm{\Sigma }}_j{|}^{(\frac{1}{2})}}exp(-\frac{1}{2}(x^{(i)}-{\mu }_j{)}^T{\mathrm{\Sigma }}_j^{-1}(x^{(i)}-{\mu }_j)){\varphi }_j}{w_j^{(i)}}$

=−∇uj∑mi=1∑kj=1w(i)j12(x(i)−μj)TΣ−1j(x(i)−μj)$=-{\mathrm{\nabla }}_{u_j}\sum _{i=1}^m\sum _{j=1}^k{w}_j^{(i)}\frac{1}{2}(x^{(i)}-{\mu }_j{)}^T{\mathrm{\Sigma }}_j^{-1}(x^{(i)}-{\mu }_j)$
=∑mi=1w(i)jΣ−1j(x(i)−μj)=0$=\sum _{i=1}^m{w}_j^{(i)}{\mathrm{\Sigma }}_j^{-1}(x^{(i)}-{\mu }_j)=0$

可得
μj=∑mi=1w(i)jx(i)∑mi=1w(i)j${\mu }_j=\frac{\sum _{i=1}^m{w}_j^{(i)}{x}^{(i)}}{\sum _{i=1}^m{w}_j^{(i)}}$

• 对方差求导：
  ∇Σj∑mi=1∑kj=1w(i)jlog1(2π)n2|Σj|(12)exp(−12(x(i)−μj)TΣ−1j(x(i)−μj))ϕjw(i)j${\mathrm{\nabla }}_{{\mathrm{\Sigma }}_j}\sum _{i=1}^m\sum _{j=1}^k{w}_j^{(i)}log\frac{\frac{1}{(2\pi {)}^{\frac{n}{2}}|{\mathrm{\Sigma }}_j{|}^{(\frac{1}{2})}}exp(-\frac{1}{2}(x^{(i)}-{\mu }_j{)}^T{\mathrm{\Sigma }}_j^{-1}(x^{(i)}-{\mu }_j)){\varphi }_j}{w_j^{(i)}}$

=∇Σj∑mi=1∑kj=1w(i)j(logΣ−12−12(x(i)−μj)TΣ−1j(x(i)−μj))$={\mathrm{\nabla }}_{{\mathrm{\Sigma }}_j}\sum _{i=1}^m\sum _{j=1}^k{w}_j^{(i)}(log{\mathrm{\Sigma }}^{-\frac{1}{2}}-\frac{1}{2}(x^{(i)}-{\mu }_j{)}^T{\mathrm{\Sigma }}_j^{-1}(x^{(i)}-{\mu }_j))$
=∑mi=1w(i)jΣ(−1)−∑mi=1w(i)j(x(i)−μj)(x(i)−μj)TΣ(−2)=0$=\sum _{i=1}^m{w}_j^{(i)}{\mathrm{\Sigma }}^{(-1)}-\sum _{i=1}^m{w}_j^{(i)}(x^{(i)}-{\mu }_j)(x^{(i)}-{\mu }_j{)}^T{\mathrm{\Sigma }}^{(-2)}=0$
可得
Σj=∑mi=1w(i)j(x(i)−μj)(x(i)−μj)T∑mi=1w(i)j${\mathrm{\Sigma }}_j=\frac{\sum _{i=1}^m{w}_j^{(i)}(x^{(i)}-{\mu }_j)(x^{(i)}-{\mu }_j{)}^T}{\sum _{i=1}^m{w}_j^{(i)}}$

• 对ϕ$\varphi$求偏导， 等式约束，用到拉格朗日乘子法， 删除常数项目得到：
  ∇ϕj∑mi=1∑kj=1w(i)jlog(ϕj)+β(∑kj=1ϕj−1)${\mathrm{\nabla }}_{{\varphi }_j}\sum _{i=1}^m\sum _{j=1}^k{w}_j^{(i)}log({\varphi }_j)+\beta (\sum _{j=1}^k{\varphi }_j-1)$
  =∑mi=1w(i)jϕj+β=0$=\sum _{i=1}^m\frac{w_j^{(i)}}{{\varphi }_j}+\beta =0$
  −β=∑mi=1∑kj=1w(i)j=m$-\beta =\sum _{i=1}^m\sum _{j=1}^k{w}_j^{(i)}=m$
  可得
  ϕj=1m∑mi=1w(i)j${\varphi }_j=\frac{1}{m}\sum _{i=1}^m{w}_j^{(i)}$

---

EM推导PLSA模型

详细过程可参考作者的另一篇博客plsaEM的详细推导