EM算法
引言
EM 算法是一种选代算法，1977 年 Dempster 等人总结提出，用于含有隐变量(hidden variable)的概率模型参数的极大似然估计，或极大后验概率估计EM算法的每次选代由两步组成:E步，求期望 (expectation);M步，求极大(maximization)所以这一算法称为期望极大算法(expectation maximizationalgorithm)，简称EM算法。
EM算法例子及解法
三硬币模型: 假设有3 硬币，分别作
A，B，CA，B，CA，B，C这些硬币正面出现的概率分别是π\piπ，ppp和qqq进行如下硬币试验:先硬币AAA，根据其结果选出硬币BBB或硬币CCC，正面选硬币BBB，反面选硬币CCC;然后掷选出的硬币，掷硬币的结果，出现正面记作 1，出现反面记作 0;独立地重复$n $次试验(这里，n=10)，观测结果如下:1,1,0,1,0,0,1,0,1,11,1,0,1,0,0,1,0,1,11,1,0,1,0,0,1,0,1,1
假设只能观测到掷硬币的结果，不能观测掷硬币的过程。问如何估计三硬币正面出现的概率，即三硬币模型的参数。
模型表达是为：
P(y∣θ)=∑zP(y,z∣θ)=∑zP(z∣θ)P(y∣z,θ)=πpy(1−p)y+(1−π)qy(1−q)1−y\begin{align} P(y|\theta) &= \sum_{z}P(y,z|\theta)=\sum_{z}P(z|\theta)P(y|z,\theta) \nonumber\\ &=\pi p^y (1-p)^y+(1-\pi)q^y(1-q)^{1-y} \nonumber \end{align}P(y∣θ)=z∑P(y,z∣θ)=z∑P(z∣θ)P(y∣z,θ)=πpy(1−p)y+(1−π)qy(1−q)1−y
这里，随机变量
yyy是观测变量，表示一次试验观测的结果是 1或0;随机变量是隐变量，表示未观测到的掷硬币AAA 的结果;θ=(π,p,q)\theta=(\pi,p,q)θ=(π,p,q)是模型参这一模型是以上数据的生成模型，注意，随机变量yyy的数据可以观测，随机变量zzz的数据不可观测。
将观察数据表示为
Y=(Y1,Y1,...,Y1)TY=(Y_1,Y_1,...,Y_1)^TY=(Y1,Y1,...,Y1)T,未观察数据表示为Z=(Z1,Z2,...,Zn)TZ=(Z_1,Z_2,...,Z_n)^TZ=(Z1,Z2,...,Zn)T,则观察数据的似然函数为P(Y∣θ)=∑zP(Z∣θ)P(Y∣Z,θ)P(Y|\theta)=\sum_zP(Z|\theta)P(Y|Z,\theta)P(Y∣θ)=z∑P(Z∣θ)P(Y∣Z,θ)
即
P(Y∣θ)=∏j=1n[πpyj(1−p)yj+(1−π)qyj(1−q)1−yj]P(Y|\theta)=\prod_{j=1}^{n}[\pi p^{y_j} (1-p)^{y_j}+(1-\pi)q^{y_j}(1-q)^{1-{y_j}}]P(Y∣θ)=j=1∏n[πpyj(1−p)yj+(1−π)qyj(1−q)1−yj]
考虑求模型参数
θ=(π,p,q)\theta=(\pi,p,q)θ=(π,p,q)的极大似然估计，即θ^=argmaxθlogP(Y∣θ)\hat{\theta}=arg \max\limits_{\theta}logP(Y|\theta)θ^=argθmaxlogP(Y∣θ)
这个问题没有解析解，只有通过迭代的方法求解。EM算法就是可以用于求解这个问题的一种迭代算法下面给出针对以上问题的EM算法其推导过程省略。
EM算法首先选取参数的初值，记作
θ(0)=(π(0),p(0),q(0))\theta^{(0)}=(\pi^{(0)},p^{(0)},q^{(0)})θ(0)=(π(0),p(0),q(0)),然后通过下面的迭代参数的估计值，直至收敛为止。第iii次迭代参数的估计值为θ(0)=(π(i),p(i),q(i))\theta^{(0)}=(\pi^{(i)},p^{(i)},q^{(i)})θ(0)=(π(i),p(i),q(i))。EM 算法的第i+1i+1i+1次迭代如下。
E步：计算在模型参数
π(i),p(i),q(i)\pi^{(i)},p^{(i)},q^{(i)}π(i),p(i),q(i)下观测数据yiy_iyi来自掷硬币B的概率μi+1=π(p(i))yj(1−(p(i)))yjπ(p(i))yj(1−(p(i)))yj+(1−π)(q(i))yj(1−(q(i)))1−yj\mu^{i+1}=\frac{\pi (p^{(i)})^{y_j} (1- (p^{(i)}))^{y_j}}{\pi (p^{(i)})^{y_j} (1- (p^{(i)}))^{y_j}+(1-\pi)(q^{(i)})^{y_j}(1-(q^{(i)}))^{1-{y_j}}}μi+1=π(p(i))yj(1−(p(i)))yj+(1−π)(q(i))yj(1−(q(i)))1−yjπ(p(i))yj(1−(p(i)))yj
M步：计算模型参数的新估计值
π(i+1)=1n∑j=1nμj(i+1)\pi^{(i+1)}=\frac{1}{n} \sum_{j=1}^{n} \mu_j^{(i+1)}π(i+1)=n1j=1∑nμj(i+1)p(i+1)=∑j=1nμj(i+1)yj∑j=1nμj(i+1)p^{(i+1)}=\frac{ \sum_{j=1}^{n} \mu_j^{(i+1)}y_j}{ \sum_{j=1}^{n} \mu_j^{(i+1)}}p(i+1)=∑j=1nμj(i+1)∑j=1nμj(i+1)yjq(i+1)=∑j=1n(1−μj(i+1))yj∑j=1n(1−μj(i+1))q^{(i+1)}=\frac{ \sum_{j=1}^{n} (1-\mu_j^{(i+1)})y_j}{ \sum_{j=1}^{n} (1-\mu_j^{(i+1)})}q(i+1)=∑j=1n(1−μj(i+1))∑j=1n(1−μj(i+1))yj
进行数字计算，假设模型参数的初值为
π(0)=0.5,p(0)=0.5,q(0)=0.5\pi^{(0)}=0.5,p^{(0)}=0.5,q^{(0)}=0.5π(0)=0.5,p(0)=0.5,q(0)=0.5
对
yj=1y_j=1yj=1与yj=0y_j=0yj=0均有μj(1)=0.5\mu_j^{(1)}=0.5μj(1)=0.5。
根据M步计算得到
π(0)=0.5,p(0)=0.6,q(0)=0.6\pi^{(0)}=0.5,p^{(0)}=0.6,q^{(0)}=0.6π(0)=0.5,p(0)=0.6,q(0)=0.6
根据E步，可以得到
μj(2)=0.5,j=1,2,...,10\mu_j^{(2)}=0.5,j=1,2,...,10μj(2)=0.5,j=1,2,...,10
继续迭代，得
π(0)=0.5,p(0)=0.6,q(0)=0.6\pi^{(0)}=0.5,p^{(0)}=0.6,q^{(0)}=0.6π(0)=0.5,p(0)=0.6,q(0)=0.6
于是得到模型参数
θ\thetaθ的极大似然估计：π^=0.5,p^=0.6,q^=0.6\hat{\pi}=0.5,\hat p=0.6,\hat q=0.6π^=0.5,p^=0.6,q^=0.6π=0.5\pi=0.5π=0.5表示硬币A是匀称的，这一结果容易理解。
如果取初始值
π(0)=0.4,p(0)=0.6,q(0)=0.7\pi^{(0)}=0.4,p^{(0)}=0.6,q^{(0)}=0.7π(0)=0.4,p(0)=0.6,q(0)=0.7,那么得到的模型参数的极大似然估计π^=0.4064,p^=0.5368,q^=0.6432\hat{\pi}=0.4064,\hat p=0.5368,\hat q=0.6432π^=0.4064,p^=0.5368,q^=0.6432。这就是说，EM算法与初始值的选择有关，选择不同的初值可能得到不同的参数估计值。
EM算法步骤和说明
一般地，用Y表示观测随机变量的数据，Z表示隐随机变量的数据Y和Z连在一起称为完全数据 (complete-data)，观测数据Y又称为不完全数据(incomplete-data)。假设给定观测数据Y，其概率分布是
P(Y∣θ)P(Y|\theta)P(Y∣θ)，其中是需要估计的模型参数，那么不完全数据Y的似然函数是P(Y∣θ)P(Y|\theta)P(Y∣θ)，对数似然函数L(0)=logP(Y∣θ)L(0)=logP(Y|\theta)L(0)=logP(Y∣θ);假设Y和Z的联合概率分布是P(YZ∣θ)P(YZ|\theta)P(YZ∣θ)，那完全的对数似然函数是logP(Y,Z∣θ)logP(Y,Z|\theta)logP(Y,Z∣θ)。
EM算法通过选代求
L(θ)=logP(Y∣θ)L(\theta)=logP(Y|\theta)L(θ)=logP(Y∣θ)的极大似然估计每次代包含两步:E 步，求期望:M步，求极大化，下面来介绍 EM 算法。
输入:观测变量数据
YYY，隐变量数据ZZZ，联合分布P(Y,Z∣θ)P(Y,Z|\theta)P(Y,Z∣θ)，条件分布P(Z∣Y,θ)P(Z|Y,\theta)P(Z∣Y,θ);
输出:模型参数
θ\thetaθ。
(1)选择参数的初值
θ(0)\theta^{(0)}θ(0)，开始迭代;
(2)E步:记
θ(i)\theta^{(i)}θ(i)为第iii次代参数的估计值，在第i+1i+1i+1次选代的E步，计算Q(θ,θ(i))=EZ[logP(Y,Z∣θ)∣Y,θ(i)]=∑ZlogP(Y,Z∣θ)P(Z∣Y,θ(i))\begin{align} Q(\theta,\theta^{(i)}) &=E_Z[logP(Y,Z|\theta)|Y,\theta^{(i)}] \nonumber\\ &=\sum_{Z}logP(Y,Z|\theta)P(Z|Y,\theta^{(i)}) \nonumber\\ \end{align}Q(θ,θ(i))=EZ[logP(Y,Z∣θ)∣Y,θ(i)]=Z∑logP(Y,Z∣θ)P(Z∣Y,θ(i))
这里，
P(Z∣Y,θ)P(Z|Y,\theta)P(Z∣Y,θ)在给定观测数据YYY和当前的参数估计θ(i)\theta^{(i)}θ(i)下隐变量数据ZZZ的条件概率分布；
(3)M步：求使
Q(θ,θ(i))Q(\theta,\theta^{(i)})Q(θ,θ(i))极大化的θ\thetaθ,确定第i+1i+1i+1次迭代的参数的估计值θ(i+1)\theta^{(i+1)}θ(i+1)θ(i+1)=argmaxθQ(θ,θ(i))\theta^{(i+1)}=\arg \max \limits_{\theta}Q(\theta,\theta^{(i)})θ(i+1)=argθmaxQ(θ,θ(i))
(4)重复第（2）步和第（3）步，直到收敛。
函数
Q(θ,θ(i))Q(\theta,\theta^{(i)})Q(θ,θ(i))是EM算法的核心，称为Q函数（Q function）。
下面关于 EM 算法作几点说明:
步骤(1)参数的初值可以任意选择，但需注意 EM算法对初值是敏感的步骤(2)E步求
Q(θ,θ(i))Q(\theta,\theta^{(i)})Q(θ,θ(i))。Q函数式中Z是未观察数据，Y是观测数据注意，Q(θ,θ(i))Q(\theta,\theta^{(i)})Q(θ,θ(i))的第1个变元表示要极大化的参数，第2个变元表示参数的当前估计值。每次迭代实际在求QQQ函数及其极大。
步骤(3)M步求
Q(θ,θ(i))Q(\theta,\theta^{(i)})Q(θ,θ(i))的极大化，得到θ(i+1)\theta^{(i+1)}θ(i+1)，完成一次代θ(i)−>θ(i+1)\theta^{(i)}->\theta^{(i+1)}θ(i)−>θ(i+1)。后面将证明每次迭代使似然函数增大或达到局部极值。
步骤(4)给出停止迭代的条件，一般是对较小的正数
ϵ1,ϵ2\epsilon_1,\epsilon_2ϵ1,ϵ2，若满足∣∣θ(i+1)−θ(i)∣∣<ϵ1或∣∣Q(θ(i+1),θ(i))−Q(θ(i),θ(i))∣∣<ϵ2|| \theta^{(i+1)}- \theta^{(i)}||<\epsilon_1 或||Q(\theta^{(i+1)},\theta^{(i)})- Q(\theta^{(i)},\theta^{(i)})||<\epsilon_2∣∣θ(i+1)−θ(i)∣∣<ϵ1或∣∣Q(θ(i+1),θ(i))−Q(θ(i),θ(i))∣∣<ϵ2
则停止迭代
本文转载自: https://blog.csdn.net/weixin_42491648/article/details/132031864
版权归原作者 菜菜的小粉猪 所有， 如有侵权，请联系我们删除。
版权归原作者 菜菜的小粉猪 所有， 如有侵权，请联系我们删除。