原文：https://mp.weixin.qq.com/s/1aBlwX11cBxw0sxONaHJRQ
笔记：修改后代码
①
a<-scan()
7.90 39.77 8.49 12.94 19.27 11.05 2.04 13.29
7.68 50.37 11.35 13.30 19.25 14.59 2.75 14.87
9.42 27.93 8.20 8.14 16.17 9.42 1.55 9.76
9.16 27.98 9.01 9.32 15.99 9.10 1.82 11.35
10.06 28.64 10.52 10.05 16.18 8.39 1.96 10.81
provinceData<-matrix(a,5)
rownames(provinceData)<-c("辽宁","浙江",'河南','甘肃','青海')
d<-dist(provinceData,method="euclidean",diag=T,upper=F,p=2) ;d
HC<-hclust(d,method="complete")
library(showtext)
plot(HC,hang=-1) #hang是表明谱系图中各类所在的位置
showtext_begin()
re4<-rect.hclust(HC,k=2,border="green");
showtext_end()
②
a<-scan()
99.00 98.00 78.00 80.00
88.00 89.00 89.00 90.00
79.00 80.00 95.00 97.00
89.00 78.00 81.00 82.00
75.00 78.00 95.00 96.00
60.00 65.00 85.00 88.00
79.00 87.00 50.00 51.00
75.00 76.00 88.00 89.00
60.00 56.00 89.00 90.00
100.00 100.00 85.00 84.00
score.data<-matrix(a,10)
colnames(score.data)<-c('数学','物理','语文','政治')
corMat<-cor( score.data[,c(1,2,3,4)] )
d<- as.dist( 1-corMat^2 )
hc<-hclust( d,method="average" )
③
companyData<-scan()
80.00 85.00 75.00 90.00
85.00 85.00 90.00 90.00
85.00 85.00 86.00 60.00
90.00 90.00 75.00 90.00
99.00 98.00 78.00 80.00
88.00 89.00 89.00 90.00
79.00 80.00 95.00 97.00
89.00 78.00 81.00 82.00
75.00 78.00 95.00 96.00
60.00 65.00 85.00 88.00
79.00 87.00 50.00 51.00
75.00 76.00 88.00 89.00
60.00 56.00 89.00 90.00
100.00 100.00 85.00 84.00
61.00 64.00 89.00 60.00
companyData<-matrix(companyData,,4)
colnames(companyData)<-c('组织文化' ,'组织氛围' ,'领导角色','员工发展')
KM<-kmeans(companyData[,c(1,2,3,4)],3,nstart=3,algorithm="Hartigan-Wong") ;KM
1、引言
聚类分析是一种经典的多元数据分析方法，在各领域有着广泛的应用，其统计思想及数学过程在各种多元统计分析教材中都有详细介绍。该方法在软件中的实现也散布于各类教材及网上的博文中，笔者在本文中总结了聚类分析方法在三种常用软件：SPSS、Matlab与R中的实现，通过几个简单的数据实例介绍各种软件的具体分析过程，供读者参考，三种软件各有优劣，读者可根据自身情况加以选择使用。说明：下文所用数据在许多教材和网上博文中都可见到，原始出处已无从追溯，故不再说明来源。
2、系统聚类法
2.1 样本聚类
例1：为了研究辽宁等5省区某年城镇居民生活消费的分布规律，根据调查资料做类型划分。
2.1.1 在SPSS中的实现
在SPSS中通过点击菜单选项的方式分以下几步实现聚类过程：
-
打开数据文件，点击菜单Analyze —> Classify —> Hierarchical cluster。
-
将变量province输入Label cases by对应的框中，变量X1-X8输入variables框中，cluster选项框默认选择Cases，即对样品聚类。
-
单击Statistics，出现对话框，选择Agglomeration schedule输出聚类过程表和Proximity matrix关系矩阵。不事先确定分类数的话在Cluster membership下默认第一项，事先确定类数的话可选择第二项并输入相应类数，然后单击Continue。
-
单击Plots，出现对话框，选择Dendrogram系统树图。Icicle（冰柱图）选择All clusters，Orientation选择Vertical，单击Continue。
-
单击Method，出现对话框，此对话框设置样本点间距离和类间距离计算方法，以及是否对原始数据进行变换等。此处我们对样本点间距离使用平方欧氏距离，类间距离用Furthest neighbor计算，不做其他数据变换，单击Continue。
-
上一步结束后返回对话框，单击OK，此时出现结果输出页面。下面依次是主要输出结果，下图是距离矩阵，按前述设置计算的样本点间距离。
-
下面三个图表分别是输出的聚类过程表，冰柱图和谱系图。
结论：从冰柱图、谱系图看，本问题分成两类是比较合适的，即河南、青海、甘肃为一类；浙江、辽宁为一类。
2.1.2 在Matlab中的实现
Matlab中实现系统聚类主要有三个命令，分别是pdist(),linkage()和dendrogram()。
其中pdist()函数计算样品间距离，得到的结果为一距离向量，分别是(2,1),(3,1),…,(n,1),(3,2),…,(n,2),…,(n,n-1)对样品间的距离；
linkage()函数利用指定的方法，如最短距离，最长距离，平均距离等建立系统聚类树；
最后由den drogram()函数画出聚类树形图。代码如下：
%%系统聚类法
provinceData=[7.90 39.77 8.49 12.94 19.27 11.05 2.04 13.29;
7.68 50.37 11.35 13.30 19.25 14.59 2.75 14.87;
9.42 27.93 8.20 8.14 16.17 9.42 1.55 9.76;
9.16 27.98 9.01 9.32 15.99 9.10 1.82 11.35;
10.06 28.64 10.52 10.05 16.18 8.39 1.96 10.81];
%pdist()函数计算样品间距离，得到的结果为一距离向量
d1=pdist(provinceData);
%计算距离方法默认为欧氏距离，即`euclidean`,
%其他可选项有'cityblock'，'minkowski'等
D1=squareform(d1);%距离向量转换为矩阵形式
%linkage()函数利用指定的方法，如最短距离，最长距离，平均距离等建立系统聚类树
z1=linkage(d1,'complete');%疑问：这里是d1还是转化为矩阵的D1？
%根据指定类间距离计算方法计算类间距离，
%可选项有'singe','centroid','ward','average'等
varlabel={'辽宁','浙江','河南','甘肃','青海'};%在聚类图上给相应样品命名
%den drogram()函数画出聚类树形图
H1=dendrogram(z1,'labels',varlabel)Matlab还可输出不同分类数下样品所在类的结果，即
T1=cluster(z1,2) %其中2是分类数,输出每个样品的分类结果
find(T1==1) %找出属于某个类的样品编号，此处为输出第一类样品序号关于聚类效果的优劣，Matlab还提供了一个cophenet()命令，用以计算系统聚类树的Cophenetic相关系数，该系数可用来对比各种不同的距离计算和不同的类间距离定义下的聚类效果，该值越接近于1，表明聚类效果越好。很多情形下，不同距离的计算方法对该系数值影响较小。
2.1.3 在R中的实现
R中实现系统聚类主要命令，分别是dist()、hclust()，和Matlab中的pdist()、linkage()功能基本一致，最后由plot()函数画出聚类树形图。
provinceData<-read.table("E:/Datasource4R/provinceData.txt",header=T)
provinceData<-as.matrix( provinceData [ , c(-1) ] )
rownames (provinceData)<-c("辽宁","浙江",'河南','甘肃','青海')
d<- dist( provinceData,method="euclidean",diag=T,upper=F,p=2 )
HC<-hclust( d,method="complete" )
showtext_begin()
plot(HC,hang=-1) #hang是表明谱系图中各类所在的位置
re4<-rect.hclust( HC,k=2,border="green" );
showtext_end()
R的作图功能很强大，可以利用rect.hclust()函数在聚类图上画出事先指定的类数的分类情况，并用矩形框框起来，非常直观。另外，利用R的ape包还可以画其他类型的各种有意思的谱系图，见下图。横向谱系图、进化分枝图、扇形谱系图、unrooted谱系图
2.2 变量聚类
例2：为了研究学科之间的关系，根据学生成绩对科目作聚类分析。
2.2.1 在SPSS中的实现
变量系统聚类与样品系统聚类在SPSS中的实现步骤基本一致，除了在下述第二步中的选项由Cases变为Variables之外。在SPSS中通过点击菜单选项的方式分以下几步实现变量系统聚类过程：
1.打开数据文件，点击菜单Analyze —> Classify —> Hierarchical cluster。
2.将4个变量输入variables框中，在Cluster选项框中选择Variables，表明要做变量聚类。
3.单击Statistics，出现对话框，选择Agglomeration schedule和Proximity matrix。不事先确定分类数的话在Cluster membership下默认第一项，事先确定类数的话可选择第二项并输入相应类数，然后单击Continue。
4.单击Plots，出现对话框，选择Dendrogram。Icicle（冰柱图）选择All clusters，Orientation选择Horizontal（Vertical亦可），单击Continue。
5.单击Method，出现对话框，此对话框设置样本点间距离和类间距离计算方法，以及是否对原始数据进行变换等。此处我们对样本点间距离使用Pearson系数，类间距离用Between-groups linkage方法，不做其他数据变换，单击Continue。
6.上一步结束后返回对话框，单击OK，此时出现结果输出页面。下面依次介绍一下主要输出结果，图21是距离矩阵和聚类过程表，按前述设置计算的样本点间距离及并类过程。
7.下两图分别是输出的冰柱图和谱系图。
2.2.2 在Matlab中的实现
Matlab中实现变量系统聚类与样品系统聚类过程一致，下面直接给出代码及结果。代码如下：
%%变量系统聚类
scoreData=[99.00 98.00 78.00 80.00;
88.00 89.00 89.00 90.00;
79.00 80.00 95.00 97.00;
89.00 78.00 81.00 82.00;
75.00 78.00 95.00 96.00;
60.00 65.00 85.00 88.00;
79.00 87.00 50.00 51.00;
75.00 76.00 88.00 89.00;
60.00 56.00 89.00 90.00;
100.00 100.00 85.00 84.00]; %%1,2,3,4列分别对应数学、物理、语文、政治
d1=pdist( scoreData','correlation' )
varlabel ={'数学','物理','语文','政治'};
z1=linkage( d1,'average' )
H2 = dendrogram( z1,0,'orientation','right','labels',varlabel )
T1=cluster(z1,2) % 2是分类数 输出每个变量的分类结果
find(T1==2) %% 找出属于某个类的变量编号
2.2.3 在R中的实现
R中实现变量系统聚类也与样品系统聚类过程一致，不再详述。
score.data<-read.table("E:/Datasource4R/scoreData.txt",header=T)
corMat<-cor( score.data[,c(-1)] )
d<- as.dist( 1-corMat^2 )
hc<-hclust( d,method="average" )另外值得一提的是，利用R的corrplot包可以作出非常直观的相关系数矩阵图，作为分类的参考也很直观。有Matlab爱好者也把这种图的生成编成Matlab函数，并将代码放在网上供大家下载，感兴趣的读者可以自行搜索。
3、K均值聚类法
例3：为研究不同公司的运营特点，调查了15个公司的组织文化、组织氛围、领导角色和员工发展4方面的内容。现要将这15个公司按照其各自的特点进行分成4类，数据如下所示。
3.1 在SPSS中的实现
在SPSS中通过点击菜单选项的方式分以下几步实现k均值聚类过程：
-
打开数据文件，点击菜单Analyze —> Classify —> K-Means cluster。
-
将变量“公司”输入Label cases by对应的框中，其余4个变量输入variables框中，Number of cluster选项框填4（即事先确定好的类数）。
-
单击Iterate，出现对话框，填最大迭代次数，默认即可，然后单击ontinue。
-
单击Save，出现对话框，两个选项都选上，单击Continue。
-
单击Options，出现对话框，如下图选择。
-
上一步结束后返回对话框，单击OK，此时出现结果输出页面。下面依次是主要输出结果，下图是初始聚类中心和迭代历史。
-
下两图分别是输出的各个类的成员信息以及最终聚类中心等其他信息。
3.2 在Matlab中的实现
Matlab实现K均值聚类的命令如下：
companyData=[80.00 85.00 75.00 90.00; 85.00 85.00 90.00 90.00; 85.00 85.00 86.00 60.00;90.00 90.00 75.00 90.00; 99.00 98.00 78.00 80.00; 88.00 89.00 89.00 90.00; 79.00 80.00 95.00 97.00; 89.00 78.00 81.00 82.00; 75.00 78.00 95.00 96.00; 60.00 65.00 85.00 88.00; 79.00 87.00 50.00 51.00; 75.00 76.00 88.00 89.00; 60.00 56.00 89.00 90.00; 100.00 100.00 85.00 84.00; 61.00 64.00 89.00 60.00]; opts = statset('Display','final'); [Idx,Ctrs,SumD,D] = kmeans(companyData,3,'Replicates',6,'Options',opts); %调用Kmeans函数 %Idx N*1的向量,存储的是每个点的聚类标号 %Ctrs K*P的矩阵,存储的是K个聚类质心位置 %SumD 1*K的和向量,存储的是类间所有点与该类质心点距离之和 %D N*K的矩阵，存储的是每个点与所有质心的距离;
Matlab还提供了一个silhouette()函数用以绘制聚类轮廓图，可根据图判断样品的聚类合理性。
[S, H] = silhouette(companyData,Idx)
3.3 在R中的实现
R中实现K均值聚类的代码及输出结果如下：
companyData<- read.table("E:/Datasource4R/companyData.txt",header=T)
KM<-kmeans(companyData[,c(-1)], 3,nstart=3, algorithm="Hartigan-Wong")
KM
## K-means clustering with 3 clusters of sizes 10, 3, 2
##
## Cluster means:
## 组织文化 组织氛围 领导角色 员工发展
## 1 86.00000 85.90000 85.10000 88.80000
## 2 60.33333 61.66667 87.66667 79.33333
## 3 82.00000 86.00000 68.00000 55.50000
##
## Clustering vector:
## [1] 1 1 3 1 1 1 1 1 1 2 3 1 2 1 2
##
## Within cluster sum of squares by cluster:
## [1] 2139.4000 622.6667 708.5000
## (between_SS / total_SS = 60.8 %)
致谢原文：https://mp.weixin.qq.com/s/1aBlwX11cBxw0sxONaHJRQ
感谢上帝的带领，祈求上帝使我长进，在基督里得着怜悯和恩惠。
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
