R语言实现混合正态分布EM最大期望估计法.docx

R 语言实现：混合正态分布 EM 最大期望估计法 因为近期在分析数据时用到了 EM 最大期望估计法这个算法，在参数估计中也用到的比较多。然而，发现国内在R 软件上实现高斯混合分布的 EM 的实例并不多，大多数是关于1 到 2 个高斯混合分布的实现，不易于推广，因此这里分享一下自己编写的k 个高斯混合分布的 EM 算法实现请大神们多多指教。并结合EMCluster 包对结果进行验算。 本文使用的密度函数为下面格式： 对应的函数原型为 em.norm(x,means,covariances,mix.prop) x 为原数据，means 为初始均值，covariances 为数据的协方差矩阵，mix.prop 为混合参数初始值。 使用的数据为 MASS 包里面的 synth.te 数据的前两列x - synth.te[,-3] 首先安装需要的包，并读取原数据。 install.packages(MASS) library(MASS) install.packages(EMCluster) library(EMCluster) install.packages(ggplot2) library(ggplot2) Y=synth.te[,c(1:2)] qplot(x=xs, y=ys, data=Y) 然后绘制相应的变量相关图： 从图上我们可以大概估计出初始的平均点为(-0.7,0.4) (-0.3,0.8)(0.5,0.6) 当然 为了试验的严谨性，我可以从两个初始均值点的情况开始估计首先输入初始参数： mustart = rbind(c(-0.5,0.3),c(0.4,0.5)) covstart = list(cov(Y), cov(Y)) probs = c(.01, .99) 然后编写 em.norm 函数，注意其中的 clusters 值需要根据不同的初始参数进行修改， em.norm = function(X,mustart,covstart,probs){ params = list(mu=mustart, var=covstart, probs=probs) clusters = 2 tol=.00001 maxits=100 showits=T