EM算法的改进.docx
《EM算法的改进.docx》由会员分享，可在线阅读，更多相关《EM算法的改进.docx（13页珍藏版）》请在冰点文库上搜索。
EM算法的改进
EM算法的改进
3.1算法介绍
3.1.1算法起步
EM方法是在十九世纪七十年代被提出的求因子极大似然估算的一种方法，这个方法适应于大多数情况，主要是从不完整数据集中对参数进行一种估算。
这种方法可以很好地实施到解决存有缺陷的数据中。
我们有一种很形象的描述来讲解这个算法。
比如两个人盛一锅不多的米饭，他们都希望两个人盛的一样多，因为他们谁都不想吃亏，这样显然也没必要找一个称来称精确的两份，可行的方法是先给每个人碗里随意盛一些，然后再观察哪
个人碗里的米饭多盛出来一些给另一个人，然后再观察，不断的将米饭多的那个
人碗里的米饭盛出来给少的那个人，这个过程不断地持续下去，用计算机语言描述就是迭代过程，直到大家看不出两个人碗里的米饭在分量上有什么差别为止。
EM算法的过程如下：
如果我们猜想得到A和B两个数据因子，这个时候我们并不知道他们的任何信息，但是他们又是相互关联的，知道A就可以推出B,相反也是一样的。
这个时候我们给定一个数字为A的值，这样我们就可以计算出来B的值，就这样不停地循环计算一直到A和B的值相差一个给定数为止。
3.1.2算法介绍
EM算法在数据挖掘中称得上是一种著名的算法，主要有以下两个过程：
E步骤：
eatimatetheexpectedvalues
M步骤：
re-estimateparameters
这个方法的一个重要作用是在于对数据因子的估算上。
可是因为算法的运算素的
相当慢，也就没有被很多人运用。
EMW法和K-means方法相结合，主要是希望找到一个极大似然估计，通过当前的设想估算来找到未知参数的期望。
另外，可以通过未知参数的期望值反过来得到极大似然估计值。
下面的图很好的介绍了最大似然算法的原理：
图3.1正态分布生成数据
这里的计算步骤如下：
Step1:
计算隐藏变量Zj的期望值E（Zj）,并且我们假定这里E（Zj）代表了参数Xi是通过第j个分布生成的数据：
E（Zj）
P（xXI
2
p（xXiI
n1
j）
n）
（6）
2—2"（Xi
（Xi
Step2：
我们可以计算其余的一个极大似然估计’（1,2）,如果能够假定由隐
藏变量Zj取得的值是第一部计算得到的E（Zj），那么我们用新的假设
’（1,2）来代替以前的假设（1,2），重复的进行迭代。
那么我们得到极
大似然估计是：
m
E（Zj）Xi
i1
jm
E（Zj）⑺
i1
这个表达式和（6）的样本均值有些类似，它是针对于单个正态分布的。
现有的式子计算了加权之后的因子数据的均值，权重是根据前面的计算。
这里我写出了主程序算法代码：
%演示EM训练算法的实现过程
[dim,Num]=size（data）;
max_iter=10;%最大迭代次数
min_improve=1e-4;
Ngauss=3;%混合高斯函数个数
Pw=zeros（1,Ngauss）;啾存权重
mu=zeros（dim,Ngauss）;%S：
存每个高斯分类的均值
sigma=zeros（dim,dim,Ngauss）;%®存协方差矩阵
fprintf（'对各个高斯分量进行初始化\n'）;
[ctcm,cv,cc,cs,mp]=qv_ft（data,Ngs）;
mu=cm;%均值初始化
forj=1:
Ngauss
gauss_labels=find（map==j）;%tfe出每个类对应的标签
Pw（j）=length（gsls）/length（mp）;%分类后是1的样本数的个数
sigma（:
:
j）=diag（std（data（:
gslbs）,0,2））;%^行向量白勺方差end
last_loglik=-Inf;%上次的概率
ifNgauss==1,%一个高斯函数不需要用EM进行估计
sa（:
:
1）=sqrtm（cov（data',1））;mu（:
1）=mean（data,2）;
else
Sa_j=se（sa（:
:
:
））;ir=0;
foriter=1:
max_iter%求样本对应函数的输出及每个高斯分量的输出
sigma_old=sigma_i;
fori=1:
Ngauss
P（:
i）=Pw（i）*p_se（data,se（mu（:
i））,squeeze（sigma_i（:
:
i）））;%每一个样本对应每一个高斯分量的输出ends=sum（P,2）;%
forj=1:
Num
P（j,：
）=P（j,：
）/s（j）；end
Pw（1:
Ngauss）=1/Num*sum（P）;%M重的估计
%均值的估计
fori=1:
Ngausssum1=0;
forj=1:
Num
sum1=sum1+P（j,i）.*data（:
j）;
end
mu（:
i）=sum1./sum（P（:
i））;
end%方差估计按照公式类似
if（（sum（sum（sum（abs（s_j-sa_d））））endendend
3.2算法公式
3.2.1理论公式
我们给定训练样本是相互独立的,
是x
（1）,x
（2）,,x（m）,这样我们推出了p（x,z）
的模拟函数为:
m
（）logp（x;）
i1
logp（x,z;）
z
（8）
在这里我们把每一类别的每种可能情况
z解出复合的分布概率之和。
我们可
以转换成先求解隐臧变量
z的概率，通过它求解所需似然估计就容易多了。
logp（x（i）;
i
log
p（x（i）;z⑴；）z⑴
（9）
log
z（i）
Q^ogW^
（10）
Qi（z⑴）log
z（i）
p（x⑴,z（i）;）
Qi（z⑴）
（11）
我们也可从中得到
Ep（x（i）,z（i）;
）/Qi（z（i））
c（z（i））畸w
（12）
此时，我们可以用一个图生动表现出算法抬高下界的过程:
Jensen不等式：
接着我们结合针对凹函数的
（i）
/（i）（i）
p（x,z;
Q（z⑴）
）EZ⑴~Qi
p（x（i）,z（i）;
Qi（z⑴）
（13）
在这里我们定义等式为一个常数，则有
/（i）（i）
（14）
p（x_,z_；_c
（i）c
Qi（z（i））
接着有
（i）⑴.
（15）
⑴）p（x,z；）（z（i）|x（i）.）Qi（z），⑴、p（z1x；）
p（x,z;）
z
这样我们能够写出EM算法的主要过程如下:
这是EM算法似然函数图像:
图3.3似然估计函数
循环重复下面程序知道它到达某个给定值：
（1）E步：
对于每个i,计算
Qi（z⑴）：
p（z（i）|x（i）;）（16）
（2）M步：
计算
（i）p（x⑴,z⑴;）
:
argmaxQi（z'，）log——（17）
iz（i）Qi（z（i））
这就是一般的EM的步骤。
3.2.2部分c语言算法
int_nVec;//待分类的向量的个数
double**_pplfVector;//待分类的向量
double**_pplfZ;//第i个向量属于第j类的概率
double**_pplfU;//对于每一个分类的均值向量
double*_plfPi;//对于每一个分类的先验概率
double***_ppplfDelta;//对于每一个分类的协方差矩阵
double_lfL;//保存他们的函数值，他们的值通过E-Step来求解
voidCluster（double**pplfVector）;〃求初始聚类中心
voidInitCluster（）;//E-Step
voidExpectation。
；//M-Step
voidMaximization。
；//矩阵求行列式//参数pplfMat为待求矩阵
//返回值为行列式的值
doubleDeterminant（double**pplfMat）;〃矩阵求逆
//返回值double型指针的指针，指向新开辟空间，其保存有原矩阵的逆
float**Ie（float**ppt）;//多维高斯分布密度函数
//paramplfVec多维空间中的点坐标，即要求该点上的概率密度
//高斯分布的均值向量//parampplfDelta高斯分布的协方差矩阵
//返回值多维空间中对应点的概率密度
floatGas（float*pc,float*pu,float**ppa）;
voidCluster（double**pplfVector）{
_pplfVector=pplfVector;//随机产生聚类中心
InitCluster（）;//下面开始E-M过程
Expectation（）;Maximization（）;
intnCount=0;//记录收敛的次数，初始化为0
doublelfL_bf;//前一次E-Step计算出的对数似然函数值
do
{lfL_bf=_lfL；
nCount++;//每循环一次，计数器加1
Etn（）;//E-M过程
Maxon（）;}
while（fabs（_lfL-lfL_bf）>ce）;//如果函数未达到要求，那么接着这个步骤；
//若收敛，则退出循环
}
用以上程序可进行自动划分分割区域，以下是通过EM算法实现图像分割的
实例：
图3.4自动划分图像图像区域（转）
其中左边的是原图，右边的是当M=2时的分割图像。
当然对于不同的M值会产生不同的划分。
这些都是EM算法在很多领域中的应用，因此说EM算法具有很好的操作性，进一步研究EM算法并提高它的算法效率是必要的。
3.3改进算法
3.3.1改进动机
传统的EM算法是基于最大似然函数的迭代过程，它的优点就是简单稳定，但是简单相对的就是算法迭代的次数多，速度慢，这就造成了时间复杂度的增大，对于很多问题来说都不具有吸引力，其他的就是方法会得到部分最好解，这对精度要求稍高的问题来说也不太现实。
因此很多人选择使用K-均值算法。
本文希
望引入cs算法解决这些问题，因为Levy飞行具有跳跃性强，不会陷入局部最优的特点，cs搜索模型速度快，可以克服EM模型的几个缺点。
因此有很大的可行性。
3.3.2算法改进
我们将cs搜索算法应用到EM中，以下是改进算法的步骤：
Step1:
最大似然的参数有多种取值，我们给定一个很大范围，参数在此范围内取得，将这一取值范围看成是cs算法中的不同的区域坐标，因此引入cs算法；
Step2:
从m个范围值中选择参数，选择参数的过程依据cs算法的过程：
第t次迭
代的第i个参数，通过Levy飞行，产生下一次迭代的初始参数s（t1）,有
s（t1）s（t）Levy（）（18）
其中是点对点乘法，Levy（）就是一个步长大小服从Levy分布的随机游走，可以表示为
Levy~ut,（13）（19）
此处，具体的步长大小是利用Mantegna法则来实现的。
Step3:
设置一个限度迭代次数，用每次选定的参数进行重复期望的计算，到达指定迭代次数算法如果还未收敛，利用cs搜索算法进行下一次参数搜索；
Step4:
循环上面过程，利用Step1的过程继续产生不同的参数，接着进行Step2,不断的迭代计算参数，当我们最后计算的似然函数参数差距小于我们给定的值，停止迭代，进入Step5;
Step5:
在计算完成后，得到了我们需要的解。
3.3.3改进算法程序
importJAVA.clusterers.EM;
importJAVA.core;
importJAVA.core.converters.*;
importJAVA.clusterers.*;
importjava.io.*;
publicclassEM{
publiclEM（）{}
publicvoidMain（）
throwsException{
EMm_cluster=newEM（）;
FileinputFile=newFile（"MyJAVA"）;
ArffLoaderarf=newArffLoader（）;
arf.setFile（inputFile）;
InstancesinstancesTrain=arf.getDataSet（）;
inle=newFile（"JAVA"）;
arf.setFile（inputFile）;
Inesinst=arf.get（）;
m_cluster.ber（insTn）;
System.out.println（"Thenumberofcluster:
"+mcl.nOCl（））;
intnum=m_cluster.numberOfClusters（）;
System.out.println（""）;
System.out.println（m_cluster.toString（））;
System.out.println（""）;
double口predict=m_cluster.clusterPriors（）;
for（inti=0;i"+pt[i]）;}}
publicstaticvoidmain（String口args）
throwsException
{lEMa=newEM（）;a.Main（）;}}
3.3.4改进算法的优点
（1）由于cs算法有不会陷入局部最优的特点，因此算法参数可以不断地跳出局部最优限制，EMT法这样就不易得到部分最好的解；
（2）改进后的算法提高了参数计算的准确性，改善了算法效果；
（3）在一定程度上降低了算法的迭代次数；
（4）在运算速度和算法适应性方面都有了很大的改进。
3.4本章小结
本章先介绍了一种著名的中的算法-EM算法，对算法的程序也进行了部分介绍，并且分析了此算法的优缺点，针对它的缺点我们加入了cs算法的Levy行走模型，解决了初始方法计算次数太多，容易得到部分最好解的缺点。
通过改进后的算法和原始算法进行比较发现改进算法的优势。
升级会员 