EM算法训练GMM的Matlab实现过程(总结)
2013-11-05 14:06
471 查看
最近看到论文中很多地方提到EM算法,之前对EM算法只是大概知道是一个参数优化算法,而不知道具体的过程,通过阅读相关的资料,大概了解了其推导过程以及实现过程。
GMM模型就是由若干个高斯分量相互组成的,通过混合的高斯模型来逼近样本的真实分布。
GMM模型估计包括三个参数:混合权重,每个高斯函数的均值以及方差,他们的递推公式如下:
权重的递推公式如下:
![](https://oscdn.geek-share.com/Uploads/Images/Content/202011/10/39af88f338edf24d54883c1cffb9db73.gif)
均值和方差的递推公式如下:
![](https://oscdn.geek-share.com/Uploads/Images/Content/202011/10/1e86dd23db0bf8da9a7107d1957bb746.gif)
![](https://oscdn.geek-share.com/Uploads/Images/Content/202011/10/999594596bafe4a4a2ea99401fa2b850.gif)
![](https://oscdn.geek-share.com/Uploads/Images/Content/202011/10/92bcc44c567294037c03dec973870589.gif)
其中M为混合高斯数,n为训练的样本数
假设现在有训练样本data集合,每一列为一个样本,行数代表样本的特征维数,采用Matlab实现EM算法的训练过程如下:
[plain]
view plaincopyprint?
%演示EM训练算法的实现过程
clc;
clear all;
load data;
[dim,Num]=size(data);
max_iter=10;%最大迭代次数
min_improve=1e-4;% 提升的精度
Ngauss=3;%混合高斯函数个数
Pw=zeros(1,Ngauss);%保存权重
mu= zeros(dim,Ngauss);%保存每个高斯分类的均值,每一列为一个高斯分量
sigma= zeros(dim,dim,Ngauss);%保存高斯分类的协方差矩阵
fprintf('采用K均值算法对各个高斯分量进行初始化\n');
[cost,cm,cv,cc,cs,map] = vq_flat(data, Ngauss);%聚类过程 map:样本所对应的聚类中心
mu=cm;%均值初始化
for j=1:Ngauss
gauss_labels=find(map==j);%找出每个类对应的标签
Pw(j)= length(gauss_labels)/length(map);%类别为1的样本个数占总样本的个数
sigma(:,:,j) = diag(std(data(:,gauss_labels),0,2)); %求行向量的方差,只取对角线,其他特征独立,并将其赋值给对角线
end
last_loglik = -Inf;%上次的概率
% 采用EM算法估计GMM的各个参数
if Ngauss==1,%一个高斯函数不需要用EM进行估计
sigma(:,:,1) = sqrtm(cov(data',1));
mu(:,1) = mean(data,2);
else
sigma_i = squeeze(sigma(:,:,:));
iter= 0;
for iter = 1:max_iter
%E 步骤
%求每一样样本对应于GMM函数的输出以及每个高斯分量的输出,
sigma_old=sigma_i;
%E步骤。。。。。
for i=1:Ngauss
P(:,i)= Pw(i) * p_single(data, squeeze(mu(:,i)), squeeze(sigma_i(:,:,i)));%每一个样本对应每一个高斯分量的输出
end
s=sum(P,2);%
for j=1:Num
P(j,:)=P(j,:)/s(j);
end
%%%Max步骤
Pw(1:Ngauss) = 1/Num*sum(P);%权重的估计
%均值的估计
for i=1:Ngauss
sum1=0;
for j=1:Num
sum1=sum1+P(j,i).*data(:,j);
end
mu(:,i)=sum1./sum(P(:,i));
end
%方差估计按照公式类似
%sigma_i
if((sum(sum(sum(abs(sigma_i- sigma_old))))<min_improve))
break;
end
end
end
子函数:
[plain]
view plaincopyprint?
function p = p_single(x, mu, sigma)
%返回高斯函数的值
[dim,N]=size(x);
p=zeros(1,N);
for i=1:N
p(i)= 1/(2*pi*abs(det(sigma)))^(length(mu)/2)*exp(-0.5*(x(:,i)-mu)'*inv(sigma)*(x(:,i)-mu));
end
注明:鉴于大家都要求vq_flat代码,这里就不一一发送到邮箱了,提供下载地http://download.csdn.net/detail/xiaoding133/5501211
GMM模型就是由若干个高斯分量相互组成的,通过混合的高斯模型来逼近样本的真实分布。
GMM模型估计包括三个参数:混合权重,每个高斯函数的均值以及方差,他们的递推公式如下:
权重的递推公式如下:
![](https://oscdn.geek-share.com/Uploads/Images/Content/202011/10/39af88f338edf24d54883c1cffb9db73.gif)
均值和方差的递推公式如下:
![](https://oscdn.geek-share.com/Uploads/Images/Content/202011/10/1e86dd23db0bf8da9a7107d1957bb746.gif)
![](https://oscdn.geek-share.com/Uploads/Images/Content/202011/10/999594596bafe4a4a2ea99401fa2b850.gif)
![](https://oscdn.geek-share.com/Uploads/Images/Content/202011/10/92bcc44c567294037c03dec973870589.gif)
其中M为混合高斯数,n为训练的样本数
假设现在有训练样本data集合,每一列为一个样本,行数代表样本的特征维数,采用Matlab实现EM算法的训练过程如下:
[plain]
view plaincopyprint?
%演示EM训练算法的实现过程
clc;
clear all;
load data;
[dim,Num]=size(data);
max_iter=10;%最大迭代次数
min_improve=1e-4;% 提升的精度
Ngauss=3;%混合高斯函数个数
Pw=zeros(1,Ngauss);%保存权重
mu= zeros(dim,Ngauss);%保存每个高斯分类的均值,每一列为一个高斯分量
sigma= zeros(dim,dim,Ngauss);%保存高斯分类的协方差矩阵
fprintf('采用K均值算法对各个高斯分量进行初始化\n');
[cost,cm,cv,cc,cs,map] = vq_flat(data, Ngauss);%聚类过程 map:样本所对应的聚类中心
mu=cm;%均值初始化
for j=1:Ngauss
gauss_labels=find(map==j);%找出每个类对应的标签
Pw(j)= length(gauss_labels)/length(map);%类别为1的样本个数占总样本的个数
sigma(:,:,j) = diag(std(data(:,gauss_labels),0,2)); %求行向量的方差,只取对角线,其他特征独立,并将其赋值给对角线
end
last_loglik = -Inf;%上次的概率
% 采用EM算法估计GMM的各个参数
if Ngauss==1,%一个高斯函数不需要用EM进行估计
sigma(:,:,1) = sqrtm(cov(data',1));
mu(:,1) = mean(data,2);
else
sigma_i = squeeze(sigma(:,:,:));
iter= 0;
for iter = 1:max_iter
%E 步骤
%求每一样样本对应于GMM函数的输出以及每个高斯分量的输出,
sigma_old=sigma_i;
%E步骤。。。。。
for i=1:Ngauss
P(:,i)= Pw(i) * p_single(data, squeeze(mu(:,i)), squeeze(sigma_i(:,:,i)));%每一个样本对应每一个高斯分量的输出
end
s=sum(P,2);%
for j=1:Num
P(j,:)=P(j,:)/s(j);
end
%%%Max步骤
Pw(1:Ngauss) = 1/Num*sum(P);%权重的估计
%均值的估计
for i=1:Ngauss
sum1=0;
for j=1:Num
sum1=sum1+P(j,i).*data(:,j);
end
mu(:,i)=sum1./sum(P(:,i));
end
%方差估计按照公式类似
%sigma_i
if((sum(sum(sum(abs(sigma_i- sigma_old))))<min_improve))
break;
end
end
end
%演示EM训练算法的实现过程 clc; clear all; load data; [dim,Num]=size(data); max_iter=10;%最大迭代次数 min_improve=1e-4;% 提升的精度 Ngauss=3;%混合高斯函数个数 Pw=zeros(1,Ngauss);%保存权重 mu= zeros(dim,Ngauss);%保存每个高斯分类的均值,每一列为一个高斯分量 sigma= zeros(dim,dim,Ngauss);%保存高斯分类的协方差矩阵 fprintf('采用K均值算法对各个高斯分量进行初始化\n'); [cost,cm,cv,cc,cs,map] = vq_flat(data, Ngauss);%聚类过程 map:样本所对应的聚类中心 mu=cm;%均值初始化 for j=1:Ngauss gauss_labels=find(map==j);%找出每个类对应的标签 Pw(j)= length(gauss_labels)/length(map);%类别为1的样本个数占总样本的个数 sigma(:,:,j) = diag(std(data(:,gauss_labels),0,2)); %求行向量的方差,只取对角线,其他特征独立,并将其赋值给对角线 end last_loglik = -Inf;%上次的概率 % 采用EM算法估计GMM的各个参数 if Ngauss==1,%一个高斯函数不需要用EM进行估计 sigma(:,:,1) = sqrtm(cov(data',1)); mu(:,1) = mean(data,2); else sigma_i = squeeze(sigma(:,:,:)); iter= 0; for iter = 1:max_iter %E 步骤 %求每一样样本对应于GMM函数的输出以及每个高斯分量的输出, sigma_old=sigma_i; %E步骤。。。。。 for i=1:Ngauss P(:,i)= Pw(i) * p_single(data, squeeze(mu(:,i)), squeeze(sigma_i(:,:,i)));%每一个样本对应每一个高斯分量的输出 end s=sum(P,2);% for j=1:Num P(j,:)=P(j,:)/s(j); end %%%Max步骤 Pw(1:Ngauss) = 1/Num*sum(P);%权重的估计 %均值的估计 for i=1:Ngauss sum1=0; for j=1:Num sum1=sum1+P(j,i).*data(:,j); end mu(:,i)=sum1./sum(P(:,i)); end %方差估计按照公式类似 %sigma_i if((sum(sum(sum(abs(sigma_i- sigma_old))))<min_improve)) break; end end end
子函数:
[plain]
view plaincopyprint?
function p = p_single(x, mu, sigma)
%返回高斯函数的值
[dim,N]=size(x);
p=zeros(1,N);
for i=1:N
p(i)= 1/(2*pi*abs(det(sigma)))^(length(mu)/2)*exp(-0.5*(x(:,i)-mu)'*inv(sigma)*(x(:,i)-mu));
end
function p = p_single(x, mu, sigma) %返回高斯函数的值 [dim,N]=size(x); p=zeros(1,N); for i=1:N p(i)= 1/(2*pi*abs(det(sigma)))^(length(mu)/2)*exp(-0.5*(x(:,i)-mu)'*inv(sigma)*(x(:,i)-mu)); end
注明:鉴于大家都要求vq_flat代码,这里就不一一发送到邮箱了,提供下载地http://download.csdn.net/detail/xiaoding133/5501211
相关文章推荐
- EM算法训练GMM的Matlab实现过程(总结)
- 使用EM算法估计GMM参数的原理及matlab实现
- Kmeans和GMM参数学习的EM算法原理和Matlab实现
- 如何实现在Oracle中应用存储过程调用MatLab函数(3)
- GMM的EM算法实现
- PCA降维算法总结以及matlab实现PCA(个人的一点理解)
- 使用ganglia 实现监控 hadoop 和 hbase(详细过程总结)
- 实现google suggest技术过程总结
- GMM的EM算法实现
- EM方法解高斯混合模型(GMM)Matlab实现
- 用户增长量统计项目实现过程踩坑总结 推荐
- EM算法与GMM的训练应用
- PCA降维算法总结以及matlab实现PCA(个人的一点理解)
- GMM混合高斯模型的EM算法及Python实现
- Matlab使用过程中内存不足问题的总结
- GMM的EM算法实现
- 浅谈Matlab中imread函数读取图像的实现过程
- 吐血总结:HOG特征+SVM训练过程
- 多级树集合分裂(SPIHT)算法的过程详解与Matlab实现(6)解码过程——主程序
- 百度地图实现定位的过程总结