模糊c均值聚类（FCM）及其在R中实现

先前所简介的层次聚类、划分聚类等，所产生的聚类簇都是非重叠的实体，即对于每个对象而言只有一个确定的分类，这种聚类方法也成硬聚类（hard-clustering）。

相比之下，模糊聚类（fuzzy clustering，或称软聚类，soft-clustering）是涉及对象之间的模糊界限时按一定要求对对象分类的数学方法，用于定量地确定对象间的亲疏关系，达到客观准确的分类目的（Bezdek, 1981）。聚类根据对象之间的相似度将对象分成组，使组内对象尽可能相似（或距离短），组间尽可能不同（距离大）。模糊聚类产生的聚类簇的界限不明显，一个对象具有属于多个类的可能性，在描述对象分类时更加谨慎，更贴合自然界的实际情况。在计算上，创建模糊边界也比在一个类中为一个对象确定边界容易得多。

模糊聚类类似于原子轨道和电子行为：电子不处于固定的位置，而是有可能出现在特定轨道壳中，如果将轨道壳视为“聚类簇”，将电子视为“对象”，那么这种可能性就是将对象分配到特定簇的概率。

模糊c均值聚类（FCM）方法概述

模糊c均值聚类（Fuzzy C-Means Clustering，FCM）是最广泛使用的模糊聚类算法之一，与k均值划分聚类较为相似，通过不断迭代使组内平方和最小化。二者最主要的区别在于，k均值划分聚类中每个对象只由一个聚类中心约束，模糊c均值聚类中对象与所有聚类中心都有关。

组内平方和计算如下：

其中u_ij^m通过下式计算：

u_ij是对象x_i属于簇c_j的程度，与对象x_i离c_j质心（聚类中心）的距离成反比；u_j是簇j的质心。模糊c均值聚类中，一个对象理论上可以被分配为所有组，并通过给定对象在各组中的成员值（membership value）描述对象属于某一类的程度（可以理解为概率）。成员值介于0和1之间，接近0代表了对象远离该聚类中心，接近1代表对象靠近该聚类中心，每个对象在各类中成员值的总和为1。

聚类中心计算为组内所有对象的质心，并根据对象的成员度进行加权。

m是模糊参数，理论上可以取大于或等于1的任何值，当接近于1时聚类结果将与k均值划分等硬聚类相似，值接近无限将导致完全模糊，通常使用m=2。

模糊c均值聚类过程可以总结如下：

（1）指定聚类数k；

（2）为每个对象任意分配一组初始成员值；

（3）使用上述公式计算每个聚类簇的质心；

（4）使用上述公式重新计算每个对象在各聚类簇中的成员值；

（5）重复（3）、（4）过程，直到成员值收敛（不再变动）或达到最大迭代次数为止。

由于在计算方式上与k均值划分聚类具有类似点，因此也可能存在相似的问题：聚类结果受初始k的数量影响显著；不同的初始化配置可能导致不同的结果；如果到达最大迭代次数时仍未收敛，则可能会陷入局部最小值。

R包cluster的模糊c均值聚类

接下来展示R语言执行模糊c均值聚类的方法。

示例数据和R代码的百度盘链接：

https://pan.baidu.com/s/1KjnCaP4Gfdys9lSrw5mUUA

cluster包提供了模糊c均值聚类方法，通过函数fanny()实现。

library(cluster)
 
#示例数据包含 15 个样本（对象），20 个变量
dat <- read.delim('data.txt', sep = '\t', row.names = 1, stringsAsFactors = FALSE, check.names = FALSE)
#推荐变量标准化，主要用于消除量纲差异，降低极端值比重
dat <- scale(dat)
 
#已知这些样本来自于 3 组数据，但不明确具体哪些样本是一组的，用来演示模糊 c 均值聚类过程
#模糊 c 均值聚类，详情 ?fanny
set.seed(123)
cm.fanny <- fanny(dat, k = 3, metric = 'euclidean', maxit = 100)
cm.fanny
 
#或者也可计算对象间距离，将距离作为模糊 c 均值聚类的输入，不使用原始数据
dis_euc <- vegan::vegdist(dat, method = 'euclidean')
set.seed(123)
cm.fanny <- fanny(dis_euc, k = 3, diss = TRUE, maxit = 100)
cm.fanny
 
#查看主要结果
names(cm.fanny)
 
head(cm.fanny$membership)    #对象属于某一类的程度，即成员值，行和为 1
cm.fanny$clustering    #最佳分类，即各对象最接近的聚类簇
#轮廓图
plot(silhouette(cm.fanny))

轮廓图显示对象分类为3组的情况良好。

将模糊c均值聚类和排序图（如PCoA）结合，观测对象之间的整体相似性或差异，以及识别的分类特征。

#以 PCoA 为例降维，并将聚类结果标注在排序图中
#以上述所得对象间欧几里得距离为例，计算 PCoA
pcoa <- cmdscale(dis_euc, k = (nrow(dat) - 1), eig = TRUE)
eig_prop <- 100*pcoa$eig/sum(pcoa$eig)
pcoa_site <- pcoa$point[ ,1:2]
plot(pcoa_site,
xlab = paste('PCoA axis1:', round(eig_prop[1], 2), '%'),
ylab = paste('PCoA axis2:', round(eig_prop[2], 2), '%'))

#标注各对象最接近的聚类簇
for (i in 1:3) {
gg <- pcoa_site[cm.fanny$clustering == i, ]
hpts <- chull(gg)
hpts <- c(hpts, hpts[1])
lines(gg[hpts, ], col = 1+1)
}

#星图展示了各对象的成员值
stars(cm.fanny$membership, location = pcoa_site, draw.segments = TRUE,
add = TRUE, scale = FALSE, col.segments = 2:4, len = 0.5)

星图展示了对象属于各聚类簇的成员值，面积越大代表成员值越高，并将各对象所属的最佳聚类簇用红线连接。数据集整体特征一目了然。

相比上述plot()作图，factoextra包提供了更便捷的方法，以下是一些参考。

#factoextra 包的可视化方案
library(factoextra)

#聚类和 PCA 结合，点的颜色和形状代表分组，分组椭圆展示为 95% 置信区间
fviz_cluster(cm.fanny, ellipse.type = 'norm', repel = TRUE, palette = 'jco',
ellipse.level = 0.95, ggtheme = theme_minimal(), legend = 'right')

#轮廓图
fviz_silhouette(cm.fanny, palette = 'jco', ggtheme = theme_minimal())

R包e1071的模糊c均值聚类

e1071包也是可用于模糊c均值聚类的常用R包，继续使用上文的示例展示e1071包中的方法。

做过表达谱时间趋势分析的同学们可能用过一个R包，Mfuzz，它就是以e1071包中的模糊c均值方法为基础进行时间表达模式聚类的。

library(e1071)

#读取数据
dat <- read.delim('data.txt', sep = '\t', row.names = 1, stringsAsFactors = FALSE, check.names = FALSE)
#推荐标准化，主要用于消除量纲差异，降低极端值比重
dat <- scale(dat)

#模糊 c 均值聚类，详情 ?cmeans
#期望分为 3 类，最大迭代 100 次
set.seed(123)
cm.cmeans <- cmeans(dat, centers = 3, iter.max = 100, dist = 'euclidean', method = 'cmeans')
cm.cmeans

#查看主要结果
summary(cm.cmeans)

head(cm.cmeans$membership) #对象属于某一类的程度，即成员值，行和为 1
cm.cmeans$cluster #最佳分类，即各对象最接近的聚类簇

可以结合一些其它包的可视化方法，观测对象分类详情。

#使用 corrplot 包描述对象归属的成员值
library(corrplot)

corrplot(cm.cmeans$membership, is.corr = FALSE)

该图代表了模糊c均值聚类计算的各对象在各类别中的归属程度，明确属于某一聚类簇的对象在该簇中具有较高的成员值。

同样地，再结合排序图观测各类别的边界，将聚类结果标注在图中。

#factoextra 包的可视化方案，和 PCA 结合
#点的颜色和形状代表分组，分组椭圆展示为 95% 置信区间
library(factoextra)
fviz_cluster(list(data = dat, cluster = cm.cmeans$cluster), palette = 'jco',
ellipse.type = 'norm', ellipse.level = 0.95, ggtheme = theme_minimal())

不同颜色和形状的点代表对象属于不同的类别，各组中心最大的点代表了聚类中心，显示该示例数据的样本分类还是很好的。

参考资料

https://www.datanovia.com/en/lessons/fuzzy-clustering-essentials/fuzzy-c-means-clustering-algorithm/

https://www.datanovia.com/en/lessons/fuzzy-clustering-essentials/cmeans-r-function-compute-fuzzy-clustering/

Bezdek J C. 1981. Pattern Recognition with Fuzzy Objective Function Algorithms. New York: Plenum Press.

友情链接