EM算法

EM是由于样本中存在隐变量，样本分布的参数又未知，则需要用EM算法来进行计算。

三硬币问题

例：（三硬币模型）三枚硬币，分别记作A,B,C.这些硬币出现正面的概率分别为$\pi ,p,q$，进行如下试验：先抛A硬币，A为正面选B，A为反面选C；然后抛选出的B或C硬币，正面记作1，反面记作0；独立重复n次试验（n=10），观测结果为$1,1,0,1,0,0,1,0,1,1$,试求出$\pi ,p,q$.

符号定义

观测值$Y=(Y_1,Y_2,...,Y_n{)}^T=(1,1,0,1,0,0,1,0,1,1{)}^T$，
未观测值$Z=(Z_1,Z_2,...,Z_n{)}^T$表示硬币A的正反面，取值为0或1，Y和Z为随机变量。
$\theta =(\pi ,p,q)$为模型参数
则观测数据的似然函数为$L(\theta )$,是产生这组观测值的概率,$$L(\theta )=P(Y|\theta )=\sum _{Z}P(Z|\theta )P(Y|Z,\theta )={\Pi }_{j=1}^n[\pi p^{y_j}(1-p{)}^{1-y_j}+(1-\pi )q^{y_j}(1-q{)}^{1-y_j}]$$
其对数似然为$LL(\theta )=log{\sum }_{j=1}^n[\pi p^{y_j}(1-p{)}^{1-y_j}+(1-\pi )q^{y_j}(1-q{)}^{1-y_j}]$
需要求使$L(\theta )$最大的$\theta$,即$\hat{\theta }=arg\underset{\theta }{\max }logL(\theta )$

但是这式子无法求解，所以需要EM算法来迭代求解：

EM算法

E步：

计算完全数据的对数似然函数对于Z的期望，即在模型$\theta$的结果下，产生Y观测值，隐变量是Z的概率的均值？？？
$$\begin{array}{rl}Q(\theta ,{\theta }^{(i)})& =E_{Z|Y,{\theta }^{(i)}}[log(P(Y,Z|\theta ))]\\ & =\sum _{j=1}^n\sum _{z_j}[p(z_j|y_j,{\theta }^{(i)})log(p(y_j,z_j|\theta ))]\\ & =\sum _{j=1}^n[p(z_j=1|y_j,{\theta }^{(i)})log(p(y_j,z_j=1|\theta ))]+[p(z_j=0|y_j,{\theta }^{(i)})log(p(y_j,z_j=0|\theta ))]\end{array}$$

$\mathbf{l}\mathbf{o}\mathbf{g}(\mathbf{P}(\mathbf{Y},\mathbf{Z}|\mathbf{\theta }))={\sum }_{j=1}^{n}log(p(y_j,z_j|\theta ))$是一个联合概率，为完全数据的对数似然函数，表示在参数$\theta$的条件下，产生这组$Y,Z$的观测值的概率
$p(y_j,z_j=1|\theta )=p(y_j|z_j=1,\theta )p(z_j=1|\theta )=\pi p^{y_j}(1-p{)}^{1-y_j}$
$p(y_j,z_j=0|\theta )=p(y_j|z_j=0,\theta )p(z_j=0|\theta )=(1-\pi )q^{y_j}(1-q{)}^{1-y_j}$

$\mathbf{P}({\mathbf{z}}_{\mathbf{j}}|{\mathbf{y}}_{\mathbf{j}},{\mathbf{\theta }}^{(\mathbf{i})})$是一个后验概率，表示在观测数据$y_j$和当前参数${\theta }^{(i)}$下隐变量$z_j$的条件概率分布，${\mu }^{(i+1)}$是一个定值，相当于常数
如：计算在模型参数${\pi }^{(i)},p^{(i)},q^{(i)}$下观测数据$y_j$是来自B硬币的条件概率，然后进行加权求和就为选中硬币B的期望
$$\begin{array}{rl}P(z_j=1|y_j,{\theta }^{(i)})& ={\mu }^{(i+1)}\\ & =\frac{p(y_j,z_j=1)}{p(y_j)}\\ & =\frac{p(y_j|z_j=1)p(z_j=1)}{{\sum }_{z_j}p(y_j,z_j)}\\ & =\frac{p(y_j|z_j=1)p(z_j=1)}{p(y_j|z_j=1)+p(y_j|z_j=0)}\\ & =\frac{{\pi }^{(i)}(p^{(i)}{)}^{y_j}(1-p^{(i)}{)}^{1-y_j}}{{\pi }^{(i)}(p^{(i)}{)}^{y_j}(1-p^{(i)}{)}^{1-y_j}+(1-{\pi }^{(i)})(q^{(i)}{)}^{y_j}(1-q^{(i)}{)}^{1-y_j}}\\ & \\ P(z_j=0|y_j,{\theta }^{(i)})& =1-{\mu }^{(i+1)}\end{array}$$
代入$Q$有：

$$\begin{array}{rl}Q(\theta ,{\theta }^{(i)})& =E_{Z|Y,{\theta }^{(i)}}[log(P(Y,Z|\theta ))]\\ & =\sum _{j=1}^n\sum _{z_j}[p(z_j|y_j,{\theta }^{(i)})log(p(y_j,z_j|\theta ))]\\ & =\sum _{j=1}^n[p(z_j=1|y_j,{\theta }^{(i)})log(p(y_j,z_j=1|\theta ))]+[p(z_j=0|y_j,{\theta }^{(i)})log(p(y_j,z_j=0|\theta ))]\\ & =\sum _{j=1}^n{\mu }_j^{(i+1)}log(\pi p^{y_j}(1-p{)}^{1-y_j})+(1-{\mu }_j^{(i+1)})log((1-\pi )q^{y_j}(1-q{)}^{1-y_j})\end{array}$$

M步：

求使$Q(\theta )$极大化的$\theta$
用$Q$函数对参数求偏导，$\frac{\partial Q}{\partial \pi },\frac{\partial Q}{\partial p},\frac{\partial Q}{\partial q}$,再令其为0即可求得参数（${\mu }^{(i+1)}$是一个定值，相当于常数）
(1)对$\pi$求偏导和估计
∂ Q ∂ π = ∑ j = 1 n { μ j ( i + 1 ) ∗ 1 p y j ( 1 − p ) 1 − y j π ∗ p y j ( 1 − p ) 1 − y j − ( 1 − μ j ( i + 1 ) ) ∗ 1 q y j ( 1 − q ) 1 − y j ( 1 − π ) ∗ q y j ( 1 − q ) 1 − y j } = ∑ j = 1 n { μ j ( i + 1 ) ∗ 1 π − ( 1 − μ j ( i + 1 ) ) ∗ 1 1 − π } = ∑ j = 1 n μ j ( j + 1 ) − π π ( 1 − π ) = 0 则 π ( i + 1 ) = 1 n ∑ j = 1 n μ j ( i + 1 ) \begin{aligned} {{\partial Q}\over{\partial \pi}}&=\sum^n_{j=1} \{\mu ^{(i+1)}_j*{1\over {p^{y_j}(1-p)^{1-y_j}\pi}}*p^{y_j}(1-p)^{1-y_j}-(1-\mu^{(i+1)}_j)*{1\over{q^{y_j}(1-q)^{1-y_j}(1-\pi)}}*q^{y_j}(1-q)^{1-y_j}\}\\ &=\sum^n_{j=1}\{\mu ^{(i+1)}_j*{1\over \pi}-(1-\mu^{(i+1)}_j)*{1\over {1-\pi}}\}\\ &=\sum^n_{j=1}{{\mu^{(j+1)}_j}-\pi\over{\pi(1-\pi)}}\\ &=0 \\ 则\pi^{(i+1)}&={1\over n}\sum^n_{j=1}\mu^{(i+1)}_j \end{aligned} ∂π∂Q​则π(i+1)​=j=1∑n​{μj(i+1)​∗pyj​(1−p)1−yj​π1​∗pyj​(1−p)1−yj​−(1−μj(i+1)​)∗qyj​(1−q)1−yj​(1−π)1​∗qyj​(1−q)1−yj​}=j=1∑n​{μj(i+1)​∗π1​−(1−μj(i+1)​)∗1−π1​}=j=1∑n​π(1−π)μj(j+1)​−π​=0=n1​j=1∑n​μj(i+1)​​
(2)对$p$求偏导和估计

$$\begin{array}{rl}\frac{\partial Q}{\partial p}& =\sum _{j=1}^n{\mu }_j^{(i+1)}\ast \frac{1}{\pi p^{y_j}(1-p{)}^{1-y_j}}\ast [\pi y_j\frac{{(1-p)}^{1-y_j}}{p^{1-y_j}}-\pi (1-y_j)\frac{p^{y_j}}{(1-p{)}^{y_j}}]\\ & =\sum _{j=1}^n{\mu }_j^{(i+1)}\ast \frac{1}{\pi p^{y_j}(1-p{)}^{1-y_j}}\ast \pi \ast \frac{y_j\ast (1-p{)}^{1-y_j}\ast (1-p{)}^{y_j}-(1-y_j)\ast p^{y_j}\ast p^{1-y_j}}{p^{1-y_j}\ast (1-p{)}^{y_j}}\\ & =\sum _{j=1}^n{\mu }_j^{(i+1)}\ast \frac{y_j(1-p)-(1-y_j)p}{p\ast (1-p)}\\ & =\sum _{j=1}^n{\mu }_j^{(i+1)}\ast \frac{y_j-p}{p\ast (1-p)}\\ & =0\\ 则p^{(i+1)}& =\frac{{\sum }_{j=1}^n{\mu }_j^{(i+1)}{y}_j}{{\sum }_{j=1}^n{\mu }_j^{(i+1)}}\end{array}$$

(3)对$q$求偏导和估计
$$\begin{array}{rl}\frac{\partial Q}{\partial q}& =\sum _{j=1}^n(1-{\mu }_j^{(i+1)})\frac{1}{(1-\pi )q^{y_j}(1-q{)}^{1-y_j}}\ast [(1-\pi )y_j\frac{(1-q{)}^{1-y_j}}{q^{1-y_j}}-(1-\pi )(1-y_j)\frac{q^{y_j}}{(1-q{)}^{y_j}}]\\ & =\sum _{j=1}^n(1-{\mu }_j^{(i+1)})\frac{1}{(1-\pi )q^{y_j}(1-q{)}^{1-y_j}}\ast (1-\pi )\ast \frac{y_j(1-q{)}^{1-y_j}\ast (1-q{)}^{y_j}-(1-y_j)q^{1-y_j}\ast q^{y_j}}{(1-q{)}^{y_j}\ast q^{1-y_j}}\\ & =\sum _{j=1}^n(1-{\mu }_j^{(i+1)})\frac{y_j(1-q)-(1-y_j)q}{(1-q)\ast q}\\ & =\sum _{j=1}^n(1-{\mu }^{(i+1)})\frac{y_j-q}{q(1-q)}\\ & =0\\ 则q^{(i+1)}& =\frac{{\sum }_{j=1}^n(1-{\mu }_j^{(i+1)})y_j}{{\sum }_{j=1}^n(1-{\mu }_j^{(i+1)})}\end{array}$$

进行数字计算，假定模型参数${\pi }^{(0)}=0.5,p^{(0)}=0.5,q^{(0)}=0.5$

则对于$y_j=1$和$y_j=0$，$E$步中都有${\mu }_j^{(1)}=\frac{{\pi }^{(i)}(p^{(i)}{)}^{y_j}(1-p^{(i)}{)}^{1-y_j}}{{\pi }^{(i)}(p^{(i)}{)}^{y_j}(1-p^{(i)}{)}^{1-y_j}+(1-{\pi }^{(i)})(q^{(i)}{)}^{y_j}(1-q^{(i)}{)}^{1-y_j}}=0.5$

则根据$M$步可以得到${\pi }^{(1)}=0.5,p^{(1)}=0.6,q^{(1)}=0.6$已经收敛

不懂就问

？为何会收敛
？期望为什么是条件概率与联合概率的求乘积
？$\mu$为何是定值