[机器学习] 小傻学EM：嚼烂EM算法

EM 算法

简介

最大期望算法（Expectation-maximization algorithm，又译为期望最大化算法），是在概率模型中寻找参数最大似然估计或者最大后验估计的算法，其中概率模型依赖于无法观测的隐性变量。
计算步骤：
输入：观测变量Y，隐变量数据Z，联合分布$P(Y,Z|\theta )$，条件概率分布$P(Y|Z,\theta )$
输入：模型参数$\theta$

1. 选择模型参数的初始值${\theta }_0$
2. E步：记${\theta }_i$为第i次迭代后的模型参数值，在第i+1次迭代的E步，计算期望：
   $$\begin{gathered} Q(\theta ,{\theta }_i)=E_{zP(Z|Y,{\theta }_i)}[\log P(Y,Z|\theta )] \\ =\sum _{Z}P(Z|Y,{\theta }_i)\log P(Y,Z|\theta ) \end{gathered}$$
3. M步：求使得$Q(\theta ,{\theta }_i)$极大化的$\theta$，确定第i+1次迭代第参数的估计值${\theta }_{i+1}$:
   $${\theta }_{i+1}=arg\max \theta Q(\theta ,{\theta }_i)$$
4. 重复2～3，直至收敛

Jenson（琴生）不等式：

若f为凸函数，有：
$$f(t{x}_1+(1-t)x_2)\le tf(x_1)+(1-t)f(x_2)$$
其中$t\in [0,1]$,同理，若f为凹函数，只需将上式中的$\le 变为\ge 即可。$
将上式中的t推广到n同样成立，即：
$$f(t_1{x}_1+t_2{x}_2+...+t_n{x}_n)\le t_{1}f(x_1)+t_{2}f(x_2)+...+t_{n}f(x_n)$$
其中$t_1,t_2,...,t_n\in [0,1],t_1+t_2+...+t_n=1$。在概率论中常以以下形式出现：
$$\psi (E(x))\le E(\psi (x))$$
X为随机变量，$\psi$是凸函数，E[X]表示X的期望。

Jensen不等式

设$f$是定义域为$R$的函数

• 如果对于任意的$x，f(x)$的二次导数$f^{\prime \prime }(x)\ge 0$，那么$f$是凸函数。
• 当$X$是向量时，若其$Hessian$矩阵$H$是半正定的，那么$f$是凸函数。

如果$f$是凸函数，$X$是随机变量，那么：$E[f(X)]\ge f(E[X])$，通俗的说法是函数值的期望大于等于期望的函数值。

log是凹函数，所以E[log(x)]<=log(E[X])
$$\begin{gathered} x均匀分布g(x)=\frac{1}{b-a} \\ EX=\int xg(x)dx=\frac{1}{b-a}\int xdx=\frac{a+b}{2} \\ f(EX)=f(\frac{a+b}{2}) \\ E(f(x))=\int f(x)g(x)dx=\frac{1}{b-a}\int f(x)dx \end{gathered}$$

原期望不等式等价于：${\sum }_i^n{p}_{i}f(x_i)\ge f({\sum }_i^n{p}_i{x}_i)$

凸函数的性质：$p_{1}f(x_1)+p_{2}f(x_2)\ge f(p_1{x}_1+p_2{x}_2);p_1+p_2=1$

归纳法，假设n=k时成立，${\sum }_i^k{p}_i=1$
$$\begin{gathered} \sum _i^{k+1}{p}_i{x}_i=p_{k+1}f(x_{k+1})+\sum _i^k{p}_{i}f(x_i) \\ =p_{k+1}f(x_{k+1})+Z_k\sum _i^k\frac{p_i}{Z_k}f(x_i);Z_k=\sum _i^k{p}_i \\ \ge p_{k+1}f(x_{k+1})+Z_{k}f(\sum _i^k\frac{p_i}{Z_k}{x}_i) \\ \ge f(p_{k+1}{x}_{k+1}+Z_k\sum _i^k\frac{p_i}{Z_k}{x}_i) \\ =f(\sum _i^{k+1}{p}_i{x}_i) \end{gathered}$$

已知：（1）样本服从分布的模型， （2）观测到的样本 求解：模型的参数

总的来说：极大似然估计就是用来估计模型参数的统计学方法
1. EM算法的导入，为什么要有EM算法？
概率模型有时既有观测变量，又有隐变量；若用概率模型的变量都是观察变量，
对于给定数据可使用极大似然估计或贝叶斯估计来计算模型参数。
但当含有隐变量时，就不能简单地使用这些方法，EM算法就是含有隐变量的概率模型参数的极大似然估计法。
《统计学习的方法》例9.1 （三银币模型）：
假设有三枚银币，分别计作A，B，C。这些硬币正面出现的概率分别是$\pi ,p和q$。进行如下掷硬币试验：先掷银币A，根据其结果选出银币B或者银币C，正面选银币B，反面选银币C；然后掷选出的硬币。掷银币的结果，出现正面计作1，出现反面计作0；独立重复n次试验（这里n=10），观测结果如下：
$$1,1,0,1,0,0,1,0,1,1$$
假设只能观测到掷硬币的结果，不能观测掷硬币的过程。问如何估计三银币正面出现的概率，即三硬币模型的参数。
解：
对于每一次试验可进行如下建模：
$$\begin{gathered} P(y|\theta )=\sum _{z}P(y,z|\theta )=\sum _{z}P(z|\theta )P(y|z,\theta ) \\ =P(z=1|\theta )P(y|z=1,\theta )+P(z=0|\theta )P(y|z=0,\theta ) \\ =\pi p^y(1-p{)}^{1-y}+(1-\pi )q^y(1-q{)}^{1-y} \end{gathered}$$
式中：
随机变量y为观测变量，表示每一次试验观测的结果是1或者0；
随机变量z是隐变量，表示未观测到的掷银币A的结果；$\theta =(\pi ,p,q)$是模型参数；
将观测数据表示为$Y=(Y_1,Y_2,...,Y_n{)}^T$；
未观测数据表示为$Z=(Z_1,Z_2,...,z_n{)}^T$；
则观测数据的似然函数为：
$$\begin{gathered} P(Y|\theta )=\sum _{Z}P(Z|\theta )P(Y|Z,\theta ) \\ =\prod _{j=1}^{n}P(y_i|\theta ) \\ =\prod _{j=1}^n[\pi p^{y_i}(1-p{)}^{1-y_i}+(1-\pi )q^{y_i}(1-q{)}^{1-y_i})]\text{\hspace{1em}(1)} \end{gathered}$$
考虑求模型参数$\theta =(\pi ,p,q)$的极大似然估计，使用对数似然函数来进行参数估计可得：
$$\begin{gathered} \hat{\theta }==\arg \underset{\theta }{\max }logP(Y|\theta ) \\ =argma{x}_{\theta }log\prod _{j=1}^n[\pi p^{y_i}(1-p{)}^{1-y_i}+(1-\pi )q^{y_i}(1-q{)}^{1-y_i})](这个地方把log提到后边，连乘变为连加) \\ =argma{x}_{\theta }\sum _{j=1}^{n}log[\pi p^{y_i}(1-p{)}^{1-y_i}+(1-\pi )q^{y_i}(1-q{)}^{1-y_i})] \end{gathered}$$
我们面对一个函数隐变量的概率模型，目标是极大化观测数据Y关于参数$\theta$的对数似然函数，即极大化式（1）的对数似然表示：
$$L(\theta )=logP(Y|\theta )=log\sum _{Z}P(Y,Z|\theta )=log\sum _{Z}P(Z|\theta )P(Y|Z,\theta )\text{\hspace{1em}(2)}$$
这一极大化的主要困难是上式中有未观测数据Z并包含和（Z为离散型时）或这积分（Z为连续型时）的对数。EM算法采用的是通过迭代逐步近似极大化$L(\theta )$:
假设这第$i次迭代后\theta 的估计值为{\theta }^{(i)},我们希望新的估计值\theta 能使L(\theta )增加，即L(\theta )>L({\theta }^{(i)})，并逐步达到极大值，为此，我们考虑两者的差：$
$$\begin{gathered} L(\theta )-L({\theta }^i)=\log \sum _{Z}P(X|Z,\theta )P(Z;\theta )-\log P(X;{\theta }^i) \\ =\log \sum _{Z}P(Z|X,{\theta }^i)\frac{P(X|Z;\theta )P(Z;\theta )}{P(Z|X,{\theta }^i)}-\log P(X;{\theta }^i) \\ \ge \sum _{Z}P(Z|X,{\theta }^i)\log \frac{P(X,Z;\theta )}{P(Z|X,{\theta }^i)}-\log P(X;{\theta }^i);(jensen) \\ 因为：\sum _{Z}P(Z|X,{\theta }^i)=1 \\ =\sum _{Z}P(Z|X,{\theta }^i)\log \frac{P(X,Z;\theta )}{P(Z|X,{\theta }^i)}-\sum _{Z}P(Z|X,{\theta }^i)\log P(X;{\theta }^i) \\ =\sum _{Z}P(Z|X,{\theta }^i)(\log \frac{P(X,Z;\theta )}{P(Z|X,{\theta }^i)}-\log P(X;{\theta }^i))(log相减，这个地方直接除过去就好，log函数性质) \\ =\sum _{Z}P(Z|X,{\theta }^i)\log \frac{P(X,Z;\theta )}{P(Z|X,{\theta }^i)P(X;{\theta }^i)}=M \end{gathered}$$
这里我们令: $B(\theta ,{\theta }^i)=L({\theta }^i)+M$, 则有
$$L(\theta )\ge B(\theta ,{\theta }^i)$$
即函数 $B(\theta ,{\theta }^i)是函数L(\theta )的一个下界，此时若设{\theta }^{i+1}达到最大，即$
$$B({\theta }^{i+1},{\theta }^i)\ge B(\theta ,{\theta }^i)$$
又 $L({\theta }^i)=B({\theta }^i,{\theta }^i)$，可进一步推得：
$$\begin{gathered} L({\theta }^{i+1})\ge B({\theta }^{i+1},{\theta }^i)\ge B({\theta }^i,{\theta }^i)=L({\theta }^i)，即： \\ L({\theta }^{i+1})\ge L({\theta }^i) \end{gathered}$$
因此，任何可以使$B(\theta ,{\theta }^i)增大的\theta 也可以使L(\theta )增大，于是问题转化为求能使B(\theta ,{\theta }^i)达到极大的{\theta }^{i+1}$，即：
$$\begin{gathered} {\theta }^{i+1}=\arg \underset{\theta }{\max }B(\theta ,{\theta }^i); \\ =\arg \underset{\theta }{\max }(L({\theta }^i)+M) \\ =\arg \underset{\theta }{\max }(\sum _{Z}p(Z|X,{\theta }^i)\log \frac{p(X,Z;\theta )}{p(Z|X,{\theta }^i)p(X;{\theta }^i)}) \\ 去掉对\theta 无影响的项 \\ =\sum _{z}p(Z|X,{\theta }^i)\log p(X,Z;\theta ) \\ =\arg \underset{\theta }{\max }Q(\theta ,{\theta }^i) \\ =E_{p(Z|X,{\theta }^i)}\log p(X,Z;\theta ) \end{gathered}$$
到此即完成连EM算法的一次迭代，求出的${\theta }^{i+1}作为下一次迭代的初始{\theta }^i。综上，可总结出EM算法的‘E步’和‘M步’分别为$：
E步：计算完全数据的对数似然函数$logP(Y,Z|\theta )关于在给定观测数据Y和当前参数\theta 下对未观测数据Z对条件概率分布P(Z|Y,{\theta }^i)的期望，也即Q函数$:
$Q(\theta ,{\theta }^i)=E_Z[logP(Y,Z|\theta )|Y,{\theta }^i]={\sum }_{Z}P(Z|Y,{\theta }^i)logP(Y,Z|\theta )$
M 步：求使得Q函数达到极大的${\theta }^{i+1}$

E 步，导出Q函数

Q ( θ , θ i ) = ∑ Z P ( Z ∣ Y , θ i ) l o g P ( Y , Z ∣ θ ) = ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) l o g [ ∏ j = 1 N P ( Z j ∣ y j , θ ] } = ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) [ ∑ j = 1 N l o g P ( Z j ∣ y j , θ ] } = ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) [ l o g P ( y 1 , z 1 ∣ θ ) + ∑ j = 2 N l o g P ( Z j ∣ y j , θ ] } = ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) l o g P ( y 1 , z 1 ∣ θ ) + ∏ j = 1 N P ( Z j ∣ y j , θ i ) ∑ j = 2 N l o g P ( Z j ∣ y j , θ ] } = ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) l o g P ( y 1 , z 1 ∣ θ ) } + ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) [ ∑ j = 2 N l o g P ( Z j ∣ y j , θ ] } = ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) l o g P ( y 1 , z 1 ∣ θ ) } + ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) [ ∑ j = 2 N l o g P ( Z 2 ∣ y 2 , θ ] } + . . . + ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) [ ∑ j = 2 N l o g P ( Z N ∣ y N , θ ] } 单 独 考 察 ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) l o g P ( y 1 , z 1 ∣ θ ) } = ∑ z 1 , z 2 , . . . , z N { ∏ j = 1 N P ( Z j ∣ y j , θ i ) P ( z 1 ∣ y 1 , θ i ) l o g P ( y 1 , z 1 ∣ θ ) } Q(\theta,\theta^i) = \sum _Z P(Z|Y, \theta^i)logP(Y, Z|\theta)\\ =\sum_{z_1, z_2, ... , z_N}\{\prod_{j=1}^NP(Z_j|y_j, \theta^i)\color{red}log[\prod_{j=1}^N\color{l}P(Z_j|y_j, \theta]\}\\ =\sum_{z_1, z_2, ... , z_N}\{\prod_{j=1}^NP(Z_j|y_j, \theta^i)\color{red}[\sum_{j=1}^Nlog\color{rd}P(Z_j|y_j, \theta]\}\\ =\sum_{z_1, z_2, ... , z_N}\{\prod_{j=1}^NP(Z_j|y_j, \theta^i)\color{red}[logP(y_1,z_1|\theta) + \sum_{j=2}^NlogP(Z_j|y_j, \theta]\}\\ \color{black}=\sum_{z_1, z_2, ... , z_N}\{\prod_{j=1}^NP(Z_j|y_j, \theta^i)\color{red}logP(y_1,z_1|\theta) +\color{l} \prod_{j=1}^NP(Z_j|y_j, \theta^i)\color{red}\sum_{j=2}^NlogP(Z_j|y_j, \theta]\}\\ \color{black}=\sum_{z_1, z_2, ... , z_N}\{\prod_{j=1}^NP(Z_j|y_j, \theta^i)\color{rd}logP(y_1,z_1|\theta) \}+\sum_{z_1, z_2, ... , z_N}\{\color{l} \prod_{j=1}^NP(Z_j|y_j, \theta^i)[\sum_{j=2}^NlogP(Z_j|y_j, \theta]\}\\ =\sum_{z_1, z_2, ... , z_N}\{\prod_{j=1}^NP(Z_j|y_j, \theta^i)\color{rd}logP(y_1,z_1|\theta) \}+\color{red}\sum_{z_1, z_2, ... , z_N}\{\prod_{j=1}^NP(Z_j|y_j, \theta^i)[\sum_{j=2}^NlogP(Z_2|y_2, \theta]\} + ... + \sum_{z_1, z_2, ... , z_N}\{\prod_{j=1}^NP(Z_j|y_j, \theta^i)[\sum_{j=2}^NlogP(Z_N|y_N, \theta]\}\\ \color{black}单独考察\sum_{z_1, z_2, ... , z_N}\{\prod_{j=1}^NP(Z_j|y_j, \theta^i)\color{rd}logP(y_1,z_1|\theta) \}\\ =\sum_{z_1, z_2, ... , z_N}\{\prod_{j=1}^NP(Z_j|y_j, \theta^i)\color{red}P(z_1|y_1, \theta^i)\color{rd}logP(y_1,z_1|\theta) \}\\ Q(θ,θi)=Z∑​P(Z∣Y,θi)logP(Y,Z∣θ)=z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)log[j=1∏N​P(Zj​∣yj​,θ]}=z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)[j=1∑N​logP(Zj​∣yj​,θ]}=z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)[logP(y1​,z1​∣θ)+j=2∑N​logP(Zj​∣yj​,θ]}=z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)logP(y1​,z1​∣θ)+j=1∏N​P(Zj​∣yj​,θi)j=2∑N​logP(Zj​∣yj​,θ]}=z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)logP(y1​,z1​∣θ)}+z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)[j=2∑N​logP(Zj​∣yj​,θ]}=z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)logP(y1​,z1​∣θ)}+z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)[j=2∑N​logP(Z2​∣y2​,θ]}+...+z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)[j=2∑N​logP(ZN​∣yN​,θ]}单独考察z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)logP(y1​,z1​∣θ)}=z1​,z2​,...,zN​∑​{j=1∏N​P(Zj​∣yj​,θi)P(z1​∣y1​,θi)logP(y1​,z1​∣θ)}

推不动了，什么时候有空在推吧～