1、EM 算法及其在高斯混合模型中的应用 EM 算法在高斯混合模型中的应用1.定义对于一个随机信号生成器，假设他的模型参数为 ，我们能观测到的数据输出为 X，不能观测到的数据输出为 Y，且随机系统模型结构的概率密度函数为（1）(,|)pxy能够观测到的一部分数据输出数据 ， 模型的另一部分输出数据 12,.Nx未知,模型的参数 也未知。 EM 算法就是要求我们从观测数据 中估12,.Nx计出参数 。2.EM 算法的描述假设每一对随机系统的输出样本 对于不同的 n相互独立，这样当(,)nxy，x 和y都已知的情况下，概率 也已知。未观测的输出y 的概(,)p p率分布也属于待求参数 。根据独立性假设
2、有：（2）1(,|)(,|)Nnpxyx3.EM算法的基本思路基本问题是求解下面的方程的解：（3）argm(,|)pxy由于X是确定量，Y是未知的，因此即使给定了 ，也无法求得 的(,|)pxy值，因此我们只能退一步求：（4）r(|)其中 （5）(|),|(|),(|,)yYyYpxpxpxy表示考虑了未知数据y的所有可能的取值Y后对 求平均值。最后根据log函数的单调性得到（4）的等效形式：（6）argmlo(|)x对于(6)给出的最优化问题，考虑用下面的递推算法解决，即：先给定一个估值并计算 ，然后更新 得到 并且有k(|)kpxk1k（7）1lo(|)l(|)kp（8）g|og|(|,)
3、l(|,)|,()|,log|,(,)yYkkyYkkyk pxypxyxpB 其中，等号在 时成立，即：k（9）,log(|)kkxEM 算法及其在高斯混合模型中的应用 于是对 的递推算法（7）可通过 进行，步骤为：log(|)px(,)kB1） 令k=0，先给出估值 k2） 然后找出 满足 （10）1k1(,),kB3） k更新为k+1并返回步骤2)直到收敛 令 (11)argmx,kk处理后 (12)1argx(,)(|)|,)m|logar(|,)|,(|,)log(|,)gxl()r(,)kk kyYk kkyYkBpyxpPxpyxpyxyxC其中 (13)(,)(|,)log(|
4、)(|,)kkyYpxpyx4.EM算法与高斯混合模型在随机系统模型中，假设 是通道 的随机信号生成器的概率密度函数的m参数， 是选中通道 的概率。记为 。()pymma假设 个随机信号生成器和通道选择随机生成器是相互独立的，从通道M输出的数据 的概率是：x(14)(|)mpx不考虑通信信息，输出 的概率为：x(15)1(|)|Mma其中：是第 个通道随机信号生成器的参数。m：参数集合 。1,2.,mMa观测数据为一批随机产生的输出信号，并且每个输出都是相互独立的，而每个输出来自哪个通道不可测。于是系统模型参数估计问题就变为通过有限的输出样本 估计 个通道参数 .12,.Nx,(1,2)maM
5、应用（12）求解，其中 可以简化为：(,)kC（161 112(,)log)|log(|)(|,)(,MMNk kmnmnmn mnkCapxpxx ）EM 算法及其在高斯混合模型中的应用 其中：11(,)log()|,)MNk kmnmnCapx2, |(|,k kn这样我们把 和 分别放在两项里面，他们不相关，可以独立考虑。map在 中应用约束条件：(,)kC1Mma用拉格朗日乘子优化 得到：ma11(|,)Nk knnpx上式的含义是，选中 号通道的概率估计 是每个观测数据 来自于1manx通道的条件概率（根据上一次估值 估算）的平均。其中的 通过m (|,)kp下式得出。 1(|)(|
6、,)kkmnnMmpxpxa中的 的优化取决于分布函数的类型，对于 为高斯分2(,)kCm (|)mpx布时， ,m其中 是分布的均值， 是方差。再经过推导，有：m11(|,)Nk kmnnapx， 1|,(|,)kkNnnpx1/211|,|=(|)kknnmkmNnpx通道参数 得更新可以看作是对 的加权，加权系数m,n可以看成是根据上一次的参数估计 算出来得 率属于 通道的(|,)kpxknxm概率。最后，上面的EM算法可能收敛到局部极大点，因此需要选择多个参数 的初始值进行迭代计算，并选择使得 最大的解， 最大的解可由下(|)px(|)p式算出：EM 算法及其在高斯混合模型中的应用 1
7、21(|)(|)|,)|(|)nmNyYMmmnpxpxyapx5.EM算法在matlab中的实现 利用上面推导出的公式，我们以二个一维的高斯分布（ , ）来1N2验证EM算法的性能，首先用二个一维的高斯分布来建立一个高斯混合模型 。H假设： ,211()N:2(,)N:H其中 与 为混合系数，且有 ，我们要用EM算法估计混合系数1212和各一维高斯分布的均值和方差 。并将利用EM算法估计出21(,)的值与真实值做比较，就可以得到该算法的性能。首先我们取 的真实值为（0.4,0.6,1,2,0.25,0.36）2121(,)这样我们得到一个混合高斯分布，他的密度函数为 ,然后产生1000个服H
8、N从 的分布的观测样本点。接下来要做的就是对这1000个样本点用EM算法进行HN处理，来估计出一组 的值。2121(,)在使用EM算法时，要首先给 设定一组初值21(,)这里假设初值为=0.3， =0.7， 0.8， 1.8， 0.2， 0.2512122Matlab具体程序如下：Y=zeros(1,10000);for i=1:10000if rand(1)0.3Y(i)=normrnd(2,sqrt(0.36),1,1);elseY(i)=normrnd(1,sqrt(0.25),1,1);endend %高斯混合模型A=0.3 0.7; %设置参数 的初值M=0.8 1.8; %设置均值
9、 的初值S=0.2 0.25; %设置方差 的初值2for n=1:1000for j=1:2a3=0;a4=0;a5=0;for k=1:10000a1=0;EM 算法及其在高斯混合模型中的应用 for t=1:2a1=A(t)*1/sqrt(2*pi*S(t)*exp(-(Y(k)-M(t)2/(2*S(t)+a1; endf= A(j) * 1/sqrt(2*pi*S(j)*exp(-(Y(k)-M(j)2/(2*S(j); a2=f/a1; a3=a2+a3; %a3 对应公式 1(|,)Nknnpmxa4=a2*Y(k)+a4; %a4 对应公式 1(|,)knna5=a2*(Y(k)-M(j)2+a5; %a5 对应公式 121(|,)()Nkknmnpxend A(j)=a3/10000; %循环更新系数值M(j)=a4/a3; %循环更新均值值S(j)=a5/a3; %循环更新方差值endend运行程序，查看变量 A,M,S 的值，与真实值比较一下，就可以得到用 EM 算法估计的高斯混合模型的性能了。得到参数为 A= 0.3063 0.6937，M= 1.0093 2.0013，S= 0.2558 0.3632，而真实值为 A=0.4 0.6，M=1 2，S=0.25 0.36 ，我们可以从中看出用 EM 算法估计的高斯混合模型的相关参数与真实值是比较接近的。
