白话机器学习算法（十一） EM（附三个硬币模型代码）

http://blog.csdn.net/wangxin110000/article/details/22404983

EM有很多公式推导，这里还是先说下EM思想，前面说过EM是最大似然估计的思路，但是他的方法是用 参数与数据  去估计 参数，直到估计出的参数基本正确；

很多初学的对EM感到困惑，常常被一堆公式所迷惑，这样不利于深入理解EM，在此我也不想用一堆公式去推倒这个思想的数学证明；我只想给初学EM的读者一个感性的认识以及一些学习建议：

拿K-means来做个例子，我们要求的是聚类中心，我们是怎么做的呢？我们假设若一个聚类中心，然后通过这些聚类中心去得到某些东西，然后我们通过得到这些东西再去更新聚类中心；放到GMM模型中，就是我们假设某个参数，然后通过这个参数得到某个中间变量，再用这个中间变量去更新参数！

如果你还不能理解？ 

在概率中我们一个数据的期望是不是可以用成我们观测数据的均值去估计，如果两者一致，就叫无偏估计，好了，这是一个统计与估计相互验证的问题，GMM中就是这样一个思路：数据的统计与估计相互验证，我猜个参数，用参数去进行统计，再用统计值对参数进行估计！（再次强调我这里说的不是EM的证明，EM的证明指的的为什么我这么做最后能趋向于正确！是一个收敛性证明！）

 

如果你想从感性到理性逐渐理解：

建议你首先看看三个硬币模型，再去看下GMM的实际执行过程，再去看下EM的证明，结合起来，相信很快就能理解！

下面是EM算法在三个硬币模型中的代码

[cpp]  view plain copy

1. #include <iostream>  
2. #include <math.h>  
3. using namespace std;  
4.   
5. typedef struct theta  
6. {  
7.     double pi;  
8.     double p;  
9.     double q;  
10.     theta(double piInput,double pInput,double qInput):pi(piInput),p(pInput),q(qInput)  
11.     {  
12.     }  
13. }Theta;  
14. void threeCoins(int *pData,int nSize,int nNum,Theta &thetaStart);  
15.   
16. int main()  
17. {  
18.   
19.     int pData[10]={1,1,0,1,0,0,1,0,1,1};  
20.     Theta thetaStart(0.5,0.5,0.5);  
21.     threeCoins(pData,10,10,thetaStart);  
22.     getchar();  
23.     return 0;  
24.   
25.   
26.   
27.   
28. }  
29. //nNum是迭代次数  
30. void threeCoins(int *pData,int nSize,int nNum,Theta &thetaStart)  
31. {     
32.     double *u=new double[nSize];  
33.     memset(u,0,sizeof(double)*nSize);  
34.     Theta thetaTemp=thetaStart;  
35.     int iNum=0;  
36.     while(iNum<nNum)  
37.     {  
38.         for (int i=0;i<nSize;i++)  
39.         {  
40.             u[i]=thetaTemp.pi*pow(thetaTemp.p,pData[i])*pow(1-thetaTemp.p,pData[i])  
41.                 /(thetaTemp.pi*pow(thetaTemp.p,pData[i])*pow(1-thetaTemp.p,pData[i])+(1-thetaTemp.pi)*pow(thetaTemp.q,pData[i])*pow(1-thetaTemp.q,1-pData[i]));  
42.         }  
43.         double uSum=0;  
44.         for (int i=0;i<nSize;i++)  
45.         {  
46.             uSum+=u[i];  
47.         }  
48.         thetaTemp.pi=uSum/nSize;  
49.         double uSum2=0;  
50.         for (int i=0;i<nSize;i++)  
51.         {  
52.             uSum2+=u[i]*pData[i];  
53.         }  
54.         thetaTemp.p=uSum2/uSum;  
55.           
56.         double uSum3=0;  
57.         for (int i=0;i<nSize;i++)  
58.         {  
59.             uSum3+=pData[i];  
60.         }  
61.         thetaTemp.q=(uSum3-uSum2)/(nSize-uSum);  
62.         cout<<thetaTemp.pi<<'\t'<<thetaTemp.p<<'\t'<<thetaTemp.q<<endl;  
63.   
64.         iNum++;  
65.           
66.     }  
67.     delete []u;  
68. }