EM算法matlab和Java实现

EM算法具体过程看前一篇博客

一、matlab实现

1.matlab代码

close all;
clear;
clc;%%
M=3;          % 高斯数量
N=600;        % 数据样本总数
th=0.000001;  % 聚合阀值
K=2;          % 输出信号保留% 待生成数据的参数
a_real =[2/3;1/6;1/6];
mu_real=[3 4 6;5 3 7];
cov_real(:,:,1)=[5 0;0 0.2];
cov_real(:,:,2)=[0.1 0;0 0.1];
cov_real(:,:,3)=[0.1 0;0 0.1];
% 这里生成的数据全部符合标准
x=[ mvnrnd( mu_real(:,1) , cov_real(:,:,1) , round(N*a_real(1)) )' ,...mvnrnd( mu_real(:,2) , cov_real(:,:,2) , round(N*a_real(2)) )' ,...mvnrnd( mu_real(:,3) , cov_real(:,:,3) , round(N*a_real(3)) )' ];figure(1),plot(x(1,:),x(2,:),'.')%% EM
% 参数初始化
a=[1/3,1/3,1/3]; %各类的比例（权重）
mu=[1 2 3;       %均值初始化2 1 4];
cov(:,:,1)=[1 0; %协方差初始化0 1];
cov(:,:,2)=[1 0;0 1];
cov(:,:,3)=[1 0;0 1];t=inf;
count=0;
figure(2),hold on
while t>=tha_old  = a;mu_old = mu;cov_old= cov;      rznk_p=zeros(M,N);%生成M行N列零矩阵for cm=1:Mmu_cm=mu(:,cm);cov_cm=cov(:,:,cm);for cn=1:N  %计算Pi(x)p_cm=exp(-0.5*(x(:,cn)-mu_cm)'/cov_cm*(x(:,cn)-mu_cm));rznk_p(cm,cn)=p_cm;endrznk_p(cm,:)=rznk_p(cm,:)/sqrt(det(cov_cm));endrznk_p=rznk_p*(2*pi)^(-K/2);
%E step%开始求rznk 相当于Pr(i|Xt)rznk=zeros(M,N);%r(Zpikn=zeros(1,M);%r(Zpikn_sum=0;for cn=1:Nfor cm=1:M%计算p(x|*)概率分布pikn(1,cm)=a(cm)*rznk_p(cm,cn);
%           pikn_sum=pikn_sum+pikn(1,cm);endfor cm=1:Mrznk(cm,cn)=pikn(1,cm)/sum(pikn);endend%求rank结束
% M stepnk=zeros(1,M);for cm=1:Mfor cn=1:Nnk(1,cm)=nk(1,cm)+rznk(cm,cn);endenda=nk/N;rznk_sum_mu=zeros(M,1);% 求均值MU%nk(cm)  相当于ni%rznk_sum_mu  ni*Xtfor cm=1:Mrznk_sum_mu=0;for cn=1:Nrznk_sum_mu=rznk_sum_mu+rznk(cm,cn)*x(:,cn);endmu(:,cm)=rznk_sum_mu/nk(cm);end% 求协方差COV   for cm=1:Mrznk_sum_cov=zeros(K,K);for cn=1:N%求协方差（Hi-u)^2rznk_sum_cov=rznk_sum_cov+rznk(cm,cn)*(x(:,cn)-mu(:,cm))*(x(:,cn)-mu(:,cm))';endcov(:,:,cm)=rznk_sum_cov/nk(cm);endt=max([norm(a_old(:)-a(:))/norm(a_old(:));norm(mu_old(:)-mu(:))/norm(mu_old(:));norm(cov_old(:)-cov(:))/norm(cov_old(:))]);temp_f=sum(log(sum(pikn)));plot(count,temp_f,'r+')count=count+1;
end  %while hold off
f=sum(log(sum(pikn)));% 输出结果
a
mu
covfigure(3),
hold on
plot(x(1,:),x(2,:),'k.');
plot(mu_real(1,:),mu_real(2,:),'*c');
plot(mu(1,:),mu(2,:),'+r');
hold offfigure(4),
hold on
for i=1:N[max_temp,ind_temp]=max(rznk(:,i));if ind_temp==1plot(x(1,i),x(2,i),'k.');endif ind_temp==2plot(x(1,i),x(2,i),'b.');endif ind_temp==3plot(x(1,i),x(2,i),'r.');end
end%fcm聚类
[center, U, OBJ_FCN]=fcm(x',3);
figure(5),
hold on
for i=1:N[max_temp,ind_temp]=max(U(:,i));if ind_temp==1plot(x(1,i),x(2,i),'k.');endif ind_temp==2plot(x(1,i),x(2,i),'b.');endif ind_temp==3plot(x(1,i),x(2,i),'r.');end
endplot(center(:,1),center(:,2),'c*')hold off

2.结果显示

二、Java实现

编译工具为：eclipse
1.Java代码

package cn.sxt.oo2;/*** 一維情況下的EM算法實現*1、求期望（e-step)*2、期望最大化（估值）（M-step)*3、循環以上兩部直到收斂*/
public class MyEM {private static final double[] points={1.0,1.3,2.2,2.6,2.8,5.0,7.3,7.4,7.5,7.7,7.9};private static double[][] w;//权值private static double[] means = {7.7,2.3};//均值private static double[] variances= {1,1};//方差private static double[] probs = {0.5,0.5};//每个类的概率；这里默认选择k=2了；/*** 高斯分布计算公式，也就是先验概率*///p(x|c)的计算private static double gaussianPro(double point,double mean,double variance){double prob = 0.0;prob = (1/(Math.sqrt(2*Math.PI)*Math.sqrt(variance)))*Math.exp(-(point-mean)*(point-mean)/(2*variance));return prob;}/*** E-step的主要逻辑*/private static double[][] countPostprob(double[] means,double[] variances,double[] points,double[] probs){int clusterNum = means.length;int pointNum = points.length;double[][] postProbs = new double[clusterNum][pointNum];  double[] denominator = new double[pointNum];for(int m = 0;m <pointNum;m++){denominator[m] = 0.0;for(int n = 0;n<clusterNum;n++){denominator[m]+=(gaussianPro(points[m], means[n], variances[n])*probs[n]);}}for(int i = 0;i<clusterNum;i++){for(int j = 0;j<pointNum;j++){postProbs[i][j]=(gaussianPro(points[j], means[i], variances[i])*probs[i])/(denominator[j]);}}return postProbs;}private static void  eStep(){w = countPostprob(means, variances, points, probs);}/*** M-step的主要逻辑之一：由E-step得到的期望，重新计算均值*/private static double[] guessMean(double[][] w,double[] points){int wLength = w.length;double[] means = new double[w.length];double[] wi = new double[wLength];for (int m = 0; m < wLength; m++) {wi[m] = 0.0;for(int n = 0; n<points.length;n++){wi[m] += w[m][n];}}for(int i = 0;i<w.length;i++){means[i] = 0.0;for(int j = 0;j<points.length;j++){means[i]+=(w[i][j]*points[j]);}means[i] /= wi[i];}return means;}/*** M-step的主要逻辑之一：由E-step得到的期望，重新计算方差*/private static double[] guessVariance(double[][] w,double[] points){int wLength = w.length;double[] means = new double[w.length];double[] variances = new double[wLength];double[] wi = new double[wLength];for (int m = 0; m < wLength; m++) {wi[m] = 0.0;for(int n = 0; n<points.length;n++){wi[m] += w[m][n];}}means = guessMean(w, points);for(int i = 0;i<wLength;i++){variances[i] = 0.0;for(int j = 0;j<points.length;j++){variances[i] +=(w[i][j]*(points[j]-means[i])*(points[j]-means[i])); }variances[i] /= wi[i];}return variances;}/*** M-step的主要逻辑之一：由E-step得到的期望，重新计算概率**/private static double[] guessProb(double[][] w){int wLength = w.length;double[] probs = new double[wLength];for(int i = 0;i<wLength;i++){probs[i] = 0.0;for(int j = 0;j<w[i].length;j++){probs[i]+=w[i][j];}probs[i] /=w[i].length;}return probs;}private static void mStep(){means = guessMean(w, points);variances = guessVariance(w, points);probs = guessProb(w);}/*** 计算前后两次迭代的参数的差值* */private static double threshold(double[] bef_values,double[] values){double diff = 0.0;for(int i = 0 ; i < values.length;i++){diff += (values[i]-bef_values[i]);}return Math.abs(diff);}public static void main(String[] args)throws Exception{int k = 2;w = new double[k][points.length];double[] bef_means;double[] bef_var;do{bef_means = means;bef_var = variances;eStep();mStep();}while(threshold(bef_means, means)<0.01&&threshold(bef_var, variances)<0.01);for(double prob:probs)System.out.println(prob);   }
}

2.结果显示

EM算法matlab和Java实现相关推荐

em算法matlab图像应用,em算法matlab程序
EM 算法作业 EM 算法简单介绍及应用 EM 算法是当存在数据缺失问题时,极... Matlab 实现根据以上推导,可以很容易实现 EM 算法估计 GMM 参数.现... 题目:matlab 实现 ...
GMM高斯混合模型及EM算法(matlab实现)
单元 %绘制男女生身高的GMM clc clear all %男女生共取2000人,女生平均身高163,男声平均身高180 male=180+sqrt(10)*randn(1,1000); %产生均值 ...
朴素贝叶斯算法matlab实现以及EM算法
这周,继续学习了朴素贝叶斯算法的一部分知识,看了matlab的贝叶斯分类算法.采用草地潮湿原因模型的一个例子来求证贝叶斯概率以及条件概率.联合概率的分析,详见日志http://blog.sina.co ...
视觉机器学习20讲-MATLAB源码示例（7）-EM算法
视觉机器学习20讲-MATLAB源码示例(7)-EM算法 1. EM算法 2. Matlab仿真 3. 仿真结果 4. 小结 1. EM算法最大期望算法(Expectation-Maximizati ...
使用EM算法估计GMM参数的原理及matlab实现
相关数学概念协方差矩阵多维高斯分布其中k=n,即x的维度. GMM的原理 GMM,高斯混合模型,是一种聚类算法. 1.GMM概念: -将k个高斯模型混合在一起,每个点出现的概率是几个高斯混合的结 ...
em算法的java实现_EM算法 - Java教程 - 找一找教程网
1.背景 2.理论 2.1.Jensen不等式优化理论中,假设 \(f\) 是定义域为实数的函数,如果对于所有的实数 \(x\) ,且二阶导数\(f''(x)\geq 0\) ,则 \(f\) 是凸 ...
java em算法_python em算法的实现
''' 数据集:伪造数据集(两个高斯分布混合) 数据集长度:1000 ------------------------------ 运行结果: ---------------------------- ...
gmm的java实现_4. EM算法-高斯混合模型GMM详细代码实现
1. 前言 EM的前3篇博文分别从数学基础.EM通用算法原理.EM的高斯混合模型的角度介绍了EM算法.按照惯例,本文要对EM算法进行更进一步的探究.就是动手去实践她. 2. GMM实现我的实现逻辑基 ...
em算法的java实现_机器学习——python模拟EM算法
<统计学习方法> 李航著第九章 EM算法我是小白一个:本文代码转载地址文末有注释:除代码和部分注释外大部分自己书写,有问题请多指教模拟课本第一个例子,即用EM算法估计三个硬币模型的参 ...

EM算法matlab和Java实现

一、matlab实现

二、Java实现

EM算法matlab和Java实现相关推荐

最新文章

热门文章