R语言入门--第十四节（聚类分析）
聚类分析（Cluster analysis）是一种数据归约技术，旨在揭露一个数据集中观测值的子集；关于“类”的概念是指由若干个观测值组成的群组。所谓聚类分析就是把大量的观测值归约为若干个类，群组内的观测值间的相似度比群间的观测相似度高。常用的聚类方法有层次聚类、划分聚类，下面一起来学习吧~
层次聚类 hierarchical agglomerative clustering
1、基本概念
- 特征： 起初定义每一个观测值属于一类；聚类就是每一次把两类聚成新的一类，直到所有的类聚成单个类为止。
- 适用范围：对于小样本来说很实用（150个观测以下）
- 聚类算法：常见的算法有 单联动(single)、全联动(complete)、平均联动(average)、质心(centroid)，Ward法。常用的是平均联动
2、算法步骤
（1）定义每一个观测值为一类；
（2）计算每一类和其它各类的距离；
- 不同算法都需要首先计算被聚类的实体之间的距离；
- 每个观测值之间常用的距离量度是欧几里得距离（针对连续型变量），例见p345；
- 由于观测的变量性质不同，对距离的贡献可能差异很大，因此，需要对各变量做归一化处理（均值为0，标准差为1）
（3）把“距离”最短的两类合并成一类，这样类的个数就减少一个；
- 这里的“距离”应该不是指第二步直接的欧几里得距离，而是涉及的算法距离。
- 不同算法主要区别就是两类之间的距离定义不同。（详见p346）
（4）重复步骤1和步骤2，直到包含所有观测值的类合并成单个类为止。
3、实验示例
基于5种营养标准含量（变量）的27种食物（观测）进行层次聚类分析，探索不同食物的相同点与不同点，并分成有意义的类。此处层次聚类算法以平均联动(average)为例。
注意！在选择变量时，一定要慎重。要选择可能对识别和理解数据中不同观测值分组有重要影响的变量。~以免“烂泥扶不上墙”~
（1）数据预处理--归一化
par(ask=TRUE)
opar <- par(no.readonly=FALSE)
data(nutrient, package="flexclust") #加载数据包
head(nutrient, 10)
#可以看出energy这一列的数值差异很大，十分影响后几个变量的贡献，因此需要归一化。
nutrient.scaled <- scale(nutrient)
head(nutrient.scaled, 10)
#这样就舒服很多了
原始nutrient部分数据截图
（2）计算欧几里得距离
d <- dist(nutrient.scaled)
#好像dist还是专门的一种储存格式
as.matrix(d)[1:4,1:4]
#局部查看一下
as.matrix(d)[1:4,1:4]
（3）平均联动层次聚类分析
fit.average <- hclust(d, method="average")
plot(fit.average, hang=-1, cex=.8, main="Average Linkage Clustering")
#绘树状图图，聚类可视化。
聚类树状图
-
hang=-1设定展示观测值标签的位置； - 树状图高度刻度代表了该高度（距离）类之间合并判定值--
如图注释：（1）首先是beef braised与smoked ham的最先合并（它们最矮--距离最近）；
（2）然后是pork roast 和 pork simmered 的合并（重新计算距离后的最相近的点，下同）；
（3）其次为chicken canned 和 tuna canned的合并；
（4）再然后是beef braised/smoked ham的类与pork roast/pork simmered的类的两小类的合并；
.............
最后是sardines canned 与一个大类的合并。
上图已经足够我们理解27种食物基于营养成分的相似性与差异性，也可进行深入的分析：将27种食物分为较少的、有意义的类。基于此，做如下分析步骤：
（1）确定聚类个数
NbClust包提供了众多的指数来确定在一个聚类分析里类的最佳数目。
library(NbClust)
nc <- NbClust(nutrient.scaled, distance="euclidean",
min.nc=2, max.nc=15, method="average")
#直接返回的结论推荐为2
- 提供的数据仍为归一化后的数据；
- distance="euclidean"交代距离量度为欧几里得方法；
- min.nc与max.nc设置可考虑的最小和最大聚类数;
- method交代层次聚类算法。
NbClust结论
par(opar)
table(nc$Best.n[1,])
#分别各有四个评判准则赞同聚类个数为2、3、5、15
barplot(table(nc$Best.n[1,]),
xlab="Numer of Clusters", ylab="Number of Criteria",
main="Number of Clusters Chosen by 26 Criteria")
#将上述结果绘制成柱状图。
26个评判准则的推荐聚类个数
最终的决策聚类方案可以均试一试，选择一个最具有解释意义的聚类方案。下面以聚类数为5，演示一下。（后来自己做了一下聚类数为2的结果，可能理解了教材为什么不用2的原因了；就是树状图最左边的那个食物单分成了一类...）
（2）获取最终的聚类方案
- 利用cutree()函数把树状图分成5类（砍树.......）
clusters <- cutree(fit.average, k=5)
- 查看5种小类的观测数目
table(clusters)
每类观测数
- 观察5种类的变量平均水平，分别为原水平、归一化后的水平
aggregate(nutrient, by=list(cluster=clusters), median)
aggregate(as.data.frame(nutrient.scaled), by=list(cluster=clusters),
median)
每类对应变量均值（原数据）
- 可视化聚类方案
plot(fit.average, hang=-1, cex=.8,
main="Average Linkage Clustering\n5 Cluster Solution")
rect.hclust(fit.average, k=5)
5类聚类树状图
由上图，尝试解释每类变量的含义：
- sardines canned 单独生成一类，因为钙的含量多；
- beef heart 单独生成一类，因为富含蛋白质和铁；
- clams类是低蛋白和高铁的；
- 另一较小类都是高能量和高脂肪的食物；
- 最大的一类均含有相对较低的铁。
划分聚类分析 partitioning clustering
1、基本概念
- 特征：首先指定类的个数K，然后观测值被随机分成K类，再重新形成聚合的类。
- 适用范围：大样本观测（成百上千）
- 算法：K均值、基于中心点的划分（PAM）
实验示例为178种意大利葡萄酒（观测）的13种化学成分含量数据（变量）。其实是3类葡萄酒，但是我们假装不知道，看能否分析出来，并且比较两种算法（K均值与PAM）的结果差异。
2、K均值聚类
K均值聚类为最常见的划分方法。
2.1、算法步骤
（1）选择K个中心点（随机选择K个观测），K数值就是我们预期的聚类数。
（2）把每个数据点分配给离它最近的中心点；第一次中心点是随机选择的，但也可以设置参数，选择最优的初始值。
（3）重新计算每类中的点到该类中心点距离的平均值；此时的中心点应该为每一类的均值中心点，对异常值敏感（之后都是如此）
（4）分配每个数据到它最近的中心点；
（5）重复步骤3、4，直到所有的观测值不在被分配或是达到最大的迭代次数（默认10次）
2.2、实验演算
（1）数据预处理：去除第一列干扰数据，并归一化数据。
data(wine, package="rattle")
head(wine)
df <- scale(wine[-1])
（2）确定待提取的聚类个数，同样可用NbClust包判断（顺序与层次聚类分析不同，如前所述，层次聚类分析在最后才确定聚类个数）
library(NbClust)
set.seed(1234)
nc <- NbClust(df, min.nc=2, max.nc=15, method="kmeans")
par(opar)
table(nc$Best.n[1,])
barplot(table(nc$Best.n[1,]),
xlab="Numer of Clusters", ylab="Number of Criteria",
main="Number of Clusters Chosen by 26 Criteria")
NbClust()直接结论推荐3个
26个评判准则的推荐聚类个数
因为K值分析在开始时，要随机选择K个中心点，每次可能获得稍有差异的方案。使用
set.seed()可以保证结果可复制。
此外，在K均值聚类中，类中总的平方值对聚类数量的曲线可能对确定聚类数量是有帮助的，详见p351
（3）K均值聚类分析
set.seed(1234)
fit.km <- kmeans(df, 3, nstart=25)
由于聚类方法对初始中心值的选择很敏感，
nstart=参数允许尝试多种初始配置并输出最好的一个。比如上例，会生成25个初始配置。
fit.km$size
#三类所包含的样本数
fit.km$centers
#各类中心值（归一化后的数据）
aggregate(wine[-1], by=list(cluster=fit.km$cluster), mean)
#还原数据
各类对应变量中心值
（4）最后将聚类结果与原始数据标准结果（第一列数据）进行比对，看看分析质量如何。
- 通过flexclust包的兰德指数评价，变化范围从-1（完全不同意）到1（完全同意）
ct.km <- table(wine$Type, fit.km$cluster)
ct.km
library(flexclust)
randIndex(ct.km)
兰德指数接近0.9，看来K均值聚类算法还不错~
3、基于中心点的划分（PAM）
K均值法对均值异常敏感，相比来说，PAM为更稳健的方法。
3.1、算法步骤
（1）随机选择K个观测（每个都称为中心点）；
（2）计算观测值到各个中心的距离；
（3）把每个观测值分配到最近的中心点；
（4）计算每个中心点到每个观测值的距离的总和（总成本）；
（5）选择一个该类中不是中心的点，并和中心点互换；
（6）重新把每个点分配到距它最近的中心点；
（7）再次计算总成本；
（8）若新的总成本比步骤4计算的总成本少，就把新的点作为中心点；
（9）重复步骤5-8，直到中心点不变。
3.2、实验演算
- 基于上述2.2，判断聚类数为3后的基础上，进行分析
library(cluster)
set.seed(1234)
fit.pam <- pam(wine[-1], k=3, stand=TRUE)
# 进行PAM分析
fit.pam$medoids
# 查看各类中心值
fit.pam$id.med
# 查看各类观测数
clusplot(fit.pam, main="Bivariate Cluster Plot")
# 聚类可视化
PAM算法的三组聚类图
- 最后计算下兰德指数，0.7左右，看来不比K均值法~
ct.pam <- table(wine$Type, fit.pam$clustering)
ct.pam
randIndex(ct.pam)
由于聚类分析有时也能分析出原本肯定不存在的类，因此需要验证（立方聚类规则，见p356）。此外，要想对自己的聚类结果更自信，可以尝试不同的聚类方法，若同一类持续复原，那么就增加了把握~~~
以上就是聚类分析的R语言学习，参考教材《R语言实战（第2版）》
