本文目的
最近这几天一直在研究如何评估Kmeans聚类算法中的最优K值。主要理论依据是《数据挖掘导论》8.5.5节中介绍的SSE和Silhouette Coefficient系数的方法评估最优K。现在记录整个实验过程，作为备忘。不过，体验过程中，由于R软件使用的还不太熟练，实现过程中有些地方可能不准确，还请大牛指点。
实验步骤概述
- 下载实验数据，点击这里。
- 取k值范围，计算出SSE，并绘制出曲线图，观察规律
- 取步骤2同样的范围，计算出Silhouette Coefficient，并绘制出曲线图，观察规律
- 根据步骤2，3观察的规律，评估出最优K值
- 验证最优聚类
步骤1：下载实验数据
按照上面给出的链接，下载解压后的数据如下所示：
数据代表的是一些文本的特征矩阵，第一行是代表的词语，出于保密原因，词语被转义了。每一行代表一个文本。第二行开始的第一列是文本的ID。所以，上图中红框中的部分才是文本的特征矩阵。矩阵中每个元素的值代表所在文本当前词语的词频。
步骤2：计算SSE
主要采用R软件计算不同K值的SSE，R脚本如下：
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
|
# 开始与结果边界
begin = 1;
length = 15;
count = 50;
end = begin + length - 1;
# 结果容器
result = c();
result[begin:end] = 0;
# 遍历计算kmeans的SSE
qc =
read
.table(
"d:/question_cluster.txt"
, header=T);
for
(i
in
begin:end) {
# 计算SSE
tmp = c();
tmp[1:count] = 0;
for
(j
in
1:count) {
kcluster = kmeans(qc, i);
tmp[j] = kcluster$tot.withinss;
}
result[i] = mean(tmp);
}
# 绘制结果
plot(result,
type
=
"o"
, xlab=
"Number of Cluster"
, ylab=
"Sum of Squer Error"
);
|
此脚本按照K从1到15，计算不同的聚类的SSE，由于kmeans算法中的随机因数，每次结果都不一样，为了减少时间结果的偶然性，对于每个k值，都重复运行50次，求出平均的SSE，最后绘制出SSE曲线，如下所示：
步骤3：计算Silhouette Coefficient
仍然采用R脚本计算，脚本如下：
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
|
# 开始与结果边界
begin = 2;
length = 15;
count = 50;
end = begin + length - 1;
# 结果容器
result = c();
result[begin:end] = -1;
# 遍历计算kmeans的SSE
library(cluster);
qc =
read
.table(
"d:/question_cluster.txt"
, header=T);
for
(i
in
begin:end) {
# Silhouette coefficient
tmp = c();
tmp[1:count] = 0;
for
(j
in
1:count) {
kcluster = pam(qc, i);
tmp[j] = kcluster$silinfo$avg.width;
}
result[i] = mean(tmp);
}
# 绘制结果
plot(result,
type
=
"o"
, xlab=
"Number of Cluster"
, ylab=
"Silhouette Coefficient"
);
|
K从2到15（k=1时无法计算），重复执行50次，得到结果如下：
步骤4：估算K值
SSE图中没有什么特点，但是充Silhouette Coefficient图中可以明显看到K=8与K=9之间有一个巨大的深沟，根据Silhouette Coefficient的定义，值较大时的K较优，所以估算最优K=8为最优点。但是，书上给出的示意图是可以看到SSE中对应的地方应该改也有一个比较大的特征点，可能是kmeans的实现算法和pam具有一定的差异造成的。
步骤5：验证最优K
采用如下代码计算8个聚类的kmeans
|
1
2
|
qc =
read
.table(
"d:/question_cluster.txt"
, header=T);
kc = kmeans(qc, centers=8)
|
使用下面的代码将实验数据从N维空间降低到2维平面，并绘制图像
|
1
2
3
|
eucQcd = dist(qc,method=
"euclidean"
)
mds = cmdscale(eucQcd)
plot(mds, col=kc$cluster, main=”8 clusters”)
|
输出图像如下所示
实验感受
上述面描述的实验过程，大致与《数据挖掘导论》8.5.5节中的一致。通过实验，发现的确可以通过Silhouette Coefficient评估最优K值，但是如何将此评估过程自动化，这还有待进一步研究。最后，也希望高人可以指出kmeans的最佳实践。
参考资料
《数据挖掘导论》（by Pang-Ning Tan） 8.5.5节
