PyTorch深度学习实战 | 高斯混合模型聚类原理分析

2023-04-29

函数组件 plot

01、问题描述为理解高斯混合模型解决聚类问题的原理，本实例采用三个一元高斯函数混合构成原始数据，再采用GMM来聚类。1)数据三个一元高斯组件函数可以采用均值和协方差表示如表1所示：▍表1三个一元高斯组件函数的均值和协方差每个高斯组件函数分配不同的权重，其中1号组件权重为30%,2号组件权重为50%，

01、问题描述

为理解高斯混合模型解决聚类问题的原理，本实例采用三个一元高斯函数混合构成原始数据，再采用GMM来聚类。

1) 数据

三个一元高斯组件函数可以采用均值和协方差表示如表1所示：

▍表1 三个一元高斯组件函数的均值和协方差

每个高斯组件函数分配不同的权重，其中1号组件权重为30%, 2号组件权重为50%，3号组件权重为20%，随机生成1000个样本数据。

2) 可视化

为了理解三个高斯组件函数是如何混合的，可以将三个一元高斯函数显示在二维坐标中，显示三个高斯组件函数的钟形图。然后，三个组件按照权重比率混合，显示三个组件函数混合后的图形。

3) 聚类

为了找到混合后的数据属于哪一个组件，可以采用聚类的方法来对数据分类。聚类后给每个数据分配1，2或者3其中的一个标签，回顾在混合三个高斯函数时的顺序，对于1000个样本数据，是否对应前300个属于1号组件，正确标签应该为1，中间500个属于2号组件，正确标签应该为2，最后200个属于3号组件，正确标签应该为3，查看聚类后得到分类标签的准确率。

02、实例分析参考解决方案

数据生成MATLAB/Octave参考代码：


mu1=[-1];
mu2=[0];
mu3=[3];
sigma1=[2.25];
sigma2=[1];
sigma3=[.25];

每个高斯组件函数分配不同的权重，其中1号组件权重为30%, 2号组件权重为50%，3号组件权重为20%，随机生成1000个样本数据，MATLAB代码如下所示：


weight1=[.3];
weight2=[.5];
weight3=[.2];
component_1=mvnrnd(mu1,sigma1,300);
component_2=mvnrnd(mu2,sigma2,500);
component_3=mvnrnd(mu3,sigma3,200);
X=[component_1;component_2;component_3];

三个一元高斯函数显示在二维坐标中，MATLAB代码如下所示：


gd1=exp(-0.5*((component_1-mu1)/sigma1).^2)/(sigma1*sqrt(2*pi));
gd2=exp(-0.5*((component_2-mu2)/sigma2).^2)/(sigma2*sqrt(2*pi));
gd3=exp(-0.5*((component_3-mu3)/sigma3).^2)/(sigma3*sqrt(2*pi));
figure;
plot(component_1,gd1,'.');hold on;
plot(component_2,gd2,'.');hold on;
plot(component_3,gd3,'.');
title('Bell cureves of three components');
xlabel('Randomly produced numbers');ylabel('Gauss distribution');

运行以上代码后，可看到三个组件函数的钟形图如图1所示。

▍图1 三个一元高斯函数的钟形图

三个组件按照权重比率混合，MATLAB代码如下所示：


gm1=gmdistribution.fit(X,3);
a=pdf(gm1,X);
figure;plot(X,a,'.');
title('Curve of Gaussian mixture distribution');
xlabel('Randomly produced numbers');
ylabel('Gauss distribution');

运行以上代码，获得三个组件混合后的图形如图2所示。

▍图2 三个一元高斯函数混合后的图形

为了找到混合后的数据属于哪一个组件，可以采用聚类的方法来对数据分类，MATLAB实现代码如下：

idx=cluster(gm1,X);

聚类后给每个数据分配1，2或者3其中的一个标签，回顾在混合三个高斯函数时的顺序，对于1000个样本数据，前300个属于1号组件，正确标签应该为1，中间500个属于2号组件，正确标签应该为2，最后200个属于3号组件，正确标签应该为3，聚类结果后得到分类标签的准确率可以采用如下代码来查看：


figure;
hold on;
for i=1:1000
ifidx(i)==1
plot(X(i),0,'r*');
elseifidx(i)==2
plot(X(i),0,'b+');
else
plot(X(i),0,'go');
    end
end
title('Plot illustrating the cluster assignment');
xlabel('Randomly produced numbers');
ylim([-0.1 0.1]);