原文链接：http://tecdat.cn/?p=4815

原文出处：拓端数据部落公众号

因为近期在分析数据时用到了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=.00001maxits=100showits=Trequire(mvtnorm)N = nrow(X)mu = params$muvar = params$varprobs = params$probsri = matrix(0, ncol=clusters, nrow=N)         ll = 0                                        it = 0                                         converged = FALSE                            if (showits)                                 cat(paste("Iterations of EM:", "\n"))while (!converged & it < maxits) { probsOld = probsllOld = llriOld = ri# Compute responsibilitiesfor (k in 1:clusters){ri[,k] = probs[k] * dmvnorm(X, mu[k,], sigma = var[[k]], log=F)}ri = ri/rowSums(ri)rk = colSums(ri)                             probs = rk/Nfor (k in 1:clusters){varmat = matrix(0, ncol=ncol(X), nrow=ncol(X))         for (i in 1:N){varmat = varmat + ri[i,k] * X[i,]%*%t(X[i,])}mu[k,] = (t(X) %*% ri[,k]) / rk[k]var[[k]] =  varmat/rk[k] - mu[k,]%*%t(mu[k,])ll[k] = -.5 * sum( ri[,k] * dmvnorm(X, mu[k,], sigma = var[[k]], log=T) )}ll = sum(ll)parmlistold =  c(llOld, probsOld)            parmlistcurrent = c(ll, probs)             it = it + 1if (showits & it == 1 | it%%5 == 0)         cat(paste(format(it), "...", "\n", sep = ""))converged = min(abs(parmlistold - parmlistcurrent)) <= tol}clust = which(round(ri)==1, arr.ind=T)       clust = clust[order(clust[,1]), 2]           out = list(probs=probs, mu=mu, var=var, resp=ri, cluster=clust, ll=ll)}

结果，可以用图像化来表示：

qplot(x=xs, y=ys, data=Y) ggplot(aes(x=xs, y=ys), data=Y) +geom_point(aes(color=factor(test$cluster)))

类似的其他情况这里不呈现了，另外r语言提供了EMCluster包可以比较方便的实现EM进行参数估计和结果的误差分析。

ret <- init.EM(Y, nclass = 2)em.aic(x=Y,emobj=list(pi = ret$pi, Mu = ret$Mu, LTSigma = ret$LTSigma))#计算结果的AIC

通过比较不同情况的AIC，我们可以筛选出适合的聚类数参数值。（欢迎转载，请注明出处。）

【大数据部落】R语言实现：混合正态分布EM最大期望估计法相关推荐

大数据之R语言速成与实战
什么是R语言? R语言由新西兰奥克兰大学的Ross Ihaka和Robert Gentleman两人共同发明.其词法和语法分别源自Scheme和S语言. R定义:一个能够自有有效的用于统计计算和绘图的 ...
经典书单、站点 —— 大数据/数据分析/R语言
1. 科普.入门 <大数据智能>,刘知远.崔安顺等著: 特色:系统,宏观和全面: 2. R 语言站点 http://langdawei.com/:R 语言数据采集与可视化:
数据分享|R语言逻辑回归、线性判别分析LDA、GAM、MARS、KNN、QDA、决策树、随机森林、SVM分类葡萄酒交叉验证ROC...
全文链接:http://tecdat.cn/?p=27384 在本文中,数据包含有关葡萄牙"Vinho Verde"葡萄酒的信息(点击文末"阅读原文"获取完整代 ...
数据分享|R语言主成分PCA、因子分析、聚类对地区经济研究分析重庆市经济指标...
原文链接:http://tecdat.cn/?p=27515 建立重庆市经济指标发展体系,以重庆市一小时经济圈作为样本,运用因子分析方法进行实证分析,在借鉴了相关评价理论和评价方法的基础上,本文提取 ...
R语言对数线性模型loglm函数_使用R语言进行混合线性模型(mixed linear model) 分析代码及详解...
1.混合线性模型简介混合线性模型,又名多层线性模型(Hierarchical linear model).它比较适合处理嵌套设计(nested)的实验和调查研究数据.此外,它还特别适合处理带有被试内 ...
大数据场景中语言虚拟机的应用和挑战
点击上方蓝字关注我们大数据场景中语言虚拟机的应用和挑战吴明瑜1,2, 陈海波1,2, 臧斌宇1,2 1 领域操作系统教育部工程研究中心,上海 200240 2 上海交通大学软件学院并行与分布式系统 ...
数据分享|R语言关联规则挖掘apriori算法挖掘评估汽车性能数据
全文链接:http://tecdat.cn/?p=32092 我们一般把一件事情发生,对另一件事情也会产生影响的关系叫做关联.而关联分析就是在大量数据中发现项集之间有趣的关联和相关联系(形如" ...
数据分享|R语言因子分析、相关性分析大学生兼职现状调查问卷数据可视化报告...
全文链接:http://tecdat.cn/?p=31765 随着大学的普及教育,大学生就业形势变得更加困难,很多学生都意识到这个问题(点击文末"阅读原文"获取完整代码数据). 相 ...
《数据科学R语言实践：面向计算推理与问题求解的案例研究法》一一2.4　探索所有男选手的跑步时间...
本节书摘来自华章计算机<数据科学R语言实践:面向计算推理与问题求解的案例研究法>一书中的第2章,第2.4节,作者:[美] 德博拉·诺兰(Deborah Nolan) 邓肯·坦普·朗(Dun ...
《数据科学R语言实践：面向计算推理与问题求解的案例研究法》一一2.5　为跨年度的个人参赛选手构造记录...
本节书摘来自华章计算机<数据科学R语言实践:面向计算推理与问题求解的案例研究法>一书中的第2章,第2.5节,作者:[美] 德博拉·诺兰(Deborah Nolan) 邓肯·坦普·朗(Dun ...

【大数据部落】R语言实现：混合正态分布EM最大期望估计法

原文链接：http://tecdat.cn/?p=4815

原文出处：拓端数据部落公众号

【大数据部落】R语言实现：混合正态分布EM最大期望估计法相关推荐

最新文章

热门文章