K-Means术语
- 簇: 所有数据点点集合,簇中的对象是相似的。
- 质心: 簇中所有点的中心(计算所有点的均值而来).
- SSE: Sum of Sqared Error(平方误差和), SSE 值越小,表示越接近它们的质心. 由于对误差取了平方,因此更加注重那么远离中心的点.
K-均值聚类算法
类别未知,发现数据集的内在联系
工作流程
- 首先, 随机确定 K 个初始点作为质心(不是数据中的点).
- 然后将数据集中的每个点分配到一个簇中, 具体来讲, 就是为每个点找到距其最近的质心, 并将其分配该质心所对应的簇. 这一步完成之后, 每个簇的质心更新为该簇说有点的平均值.
上述过程的 伪代码 如下:
- 创建 k 个点作为起始质心(通常是随机选择)
- 当任意一个点的簇分配结果发生改变时
- 对数据集中的每个数据点
- 对每个质心
- 计算质心与数据点之间的距离
- 将数据点分配到距其最近的簇
- 对每个质心
- 对每一个簇, 计算簇中所有点的均值并将均值作为质心
算法函数
def loadDataSet(fileName):
dataMat = []
fr = open(fileName)
for line in fr.readlines():
curLine = line.strip().split('\t')
fltLine = list(map(float,curLine)) # python3 map 返回map对象,要先转换为list
dataMat.append(fltLine)
return dataMat
# 计算两个向量之间的欧式距离
def distEclud(vecA, vecB):
return sqrt(sum(power(vecA - vecB, 2)))
from numpy import *
def randCent(dataSet, k):
"""
为给定数据集构建一个包含 k 个随机质心的集合。随机质心必须要在整个数据集
的边界之内,这可以通过找到数据集每一维的最小和最大值来完成。
然后生成 0~1.0 之间的随机数并通过取值范围和最小值,
以便确保随机点在数据的边界之内。
"""
n = shape(dataSet)[1] # 列的数量
centroids = mat(zeros((k,n))) # 创建k个质心矩阵
for j in range(n): # 创建随机簇质心,并且在每一维的边界内
minJ = min(dataSet[:,j]) # 最小值
rangeJ = float(max(dataSet[:,j]) - minJ) # 范围 = 最大值 - 最小值
centroids[:,j] = mat(minJ + rangeJ * random.rand(k,1)) # 随机生成
return centroids
dataMat = mat(loadDataSet('./testSet.txt'))
randCent(dataMat, 2)
matrix([[ 1.14709253, -3.70837142],
[ 1.38629699, 2.68556404]])
def kMeans(dataSet, k, distMeas=distEclud, createCent=randCent):
"""
k-means 聚类算法
该算法会创建k个质心,然后将每个点分配到最近的质心,再重新计算质心。
这个过程重复数次,直到数据点的簇分配结果不再改变位置。
运行结果(多次运行结果可能会不一样,可以试试,原因为随机质心的影响,但总的结果是对的, 因为数据足够相似,也可能会陷入局部最小值)
"""
m = shape(dataSet)[0] # 行数
clusterAssment = mat(zeros((m, 2))) # 创建一个与 dataSet 行数一样,但是有两列的矩阵,用来保存簇分配结果
centroids = createCent(dataSet, k) # 创建质心,随机k个质心
clusterChanged = True
while clusterChanged:
clusterChanged = False
for i in range(m): # 循环每一个数据点并分配到最近的质心中去
minDist = inf; minIndex = -1
for j in range(k):
distJI = distMeas(centroids[j,:],dataSet[i,:]) # 计算数据点到质心的距离
if distJI < minDist: # 如果距离比 minDist(最小距离)还小,更新 minDist(最小距离)和最小质心的 index(索引)
minDist = distJI; minIndex = j
if clusterAssment[i, 0] != minIndex: # 簇分配结果改变
clusterChanged = True # 簇改变
clusterAssment[i, :] = minIndex,minDist**2 # 更新簇分配结果为最小质心的 index(索引),minDist(最小距离)的平方
print(centroids)
for cent in range(k): # 更新质心
ptsInClust = dataSet[nonzero(clusterAssment[:, 0].A==cent)[0]] # 获取该簇中的所有点
centroids[cent,:] = mean(ptsInClust, axis=0) # 将质心修改为簇中所有点的平均值,mean 就是求平均值的
return centroids, clusterAssment
myCentroids, clusteAssing = kMeans(dataMat, 4)
[[-3.65984336 -4.13676524]
[-4.36183554 -3.94553196]
[-3.08009682 2.86472852]
[-4.79304975 -0.87171595]]
[[ 0.75862959 -3.06925659]
[-4.31794614 -3.15463957]
[ 0.30731902 2.68529874]
[-3.92038367 -2.057115 ]]
[[ 2.65077367 -2.79019029]
[-3.7122346 -3.4052254 ]
[ 0.25122803 3.04222063]
[-3.30723845 -1.69218527]]
[[ 2.72102136 -2.61215086]
[-3.5980785 -3.32781167]
[ 0.23721978 3.13048439]
[-3.3514497 -1.0952283 ]]
[[ 2.8692781 -2.54779119]
[-3.458902 -3.12120435]
[ 0.38075386 3.12396831]
[-3.295564 0.22418486]]
[[ 2.8692781 -2.54779119]
[-3.38237045 -2.9473363 ]
[ 1.21336621 3.14825539]
[-3.17006745 2.60393509]]
[[ 2.80293085 -2.7315146 ]
[-3.38237045 -2.9473363 ]
[ 2.31553173 3.07737886]
[-2.64677572 2.78993217]]
[[ 2.80293085 -2.7315146 ]
[-3.38237045 -2.9473363 ]
[ 2.6265299 3.10868015]
[-2.46154315 2.78737555]]
myCentroids
matrix([[ 2.80293085, -2.7315146 ],
[-3.38237045, -2.9473363 ],
[ 2.6265299 , 3.10868015],
[-2.46154315, 2.78737555]])
import matplotlib.pyplot as plt
def showCluster(dataSet, k, centroids, clusterAssment):
numSamples, dim = dataSet.shape
if dim != 2:
print("Sorry! I can not draw because the dimension of your data is not 2!")
return 1
mark = ['or', 'ob', 'og', 'ok', '^r', '+r', 'sr', 'dr', '<r', 'pr']
if k > len(mark):
print("Sorry! Your k is too large! please contact Zouxy")
return 1
# draw all samples
for i in range(numSamples):
markIndex = int(clusterAssment[i, 0])
plt.plot(dataSet[i, 0], dataSet[i, 1], mark[markIndex])
mark = ['Dr', 'Db', 'Dg', 'Dk', '^b', '+b', 'sb', 'db', '<b', 'pb']
# draw the centroids
for i in range(k):
plt.plot(centroids[i, 0], centroids[i, 1], '+', markersize = 12)
plt.show()
showCluster(dataMat, 4, myCentroids, clusteAssing)
K-Means的缺陷
在 kMeans 的函数测试中,可能偶尔会陷入局部最小值(局部最优的结果,但不是全局最优的结果).
为了解决这个问题,我们可以对生成的簇进行后处理,一种方法是将具有最大SSE值的簇划分成两个簇。具体实现时可以将最大簇包含的点过滤出来并在这些点上运行K-均值算法,令k设为2。
为了保持簇总数不变,可以将某两个簇进行合并。从上图中很明显就可以看出,应该将上图下部两个出错的簇质心进行合并。那么问题来了,我们可以很容易对二维数据上的聚类进行可视化, 但是如果遇到40维的数据应该如何去做?
有两种可以量化的办法:合并最近的质心,或者合并两个使得SSE增幅最小的质心。 第一种思路通过计算所有质心之间的距离, 然后合并距离最近的两个点来实现。第二种方法需要合并两个簇然后计算总SSE值。必须在所有可能的两个簇上重复上述处理过程,直到找到合并最佳的两个簇为止。
因为上述后处理过程实在是有些繁琐,所以有更厉害的大佬提出了另一个称之为二分K-均值(bisecting K-Means)的算法.
二分 K-Means 聚类算法
该算法首先将所有点作为一个簇,然后将该簇一分为二。
之后选择其中一个簇继续进行划分,选择哪一个簇进行划分取决于对其划分时候可以最大程度降低 SSE(平方和误差)的值。
上述基于 SSE 的划分过程不断重复,直到得到用户指定的簇数目为止。
伪代码描述如下:
- 将所有点看成一个簇
- 当簇数目小于 k 时
- 对于每一个簇
- 计算总误差
- 在给定的簇上面进行 KMeans 聚类(k=2)
- 计算将该簇一分为二之后的总误差
- 选择使得误差最小的那个簇进行划分操作
def biKMeans(dataMat, k, distMeas=distEclud):
m = shape(dataMat)[0]
clusterAssment = mat(zeros((m, 2))) # 保存每个数据点的簇分配结果和平方误差
centroid0 = mean(dataMat, axis=0).tolist()[0] # 质心初始化为所有数据点的均值
centList = [centroid0] # 初始化只有 1 个质心的 list
for j in range(m): # 计算所有数据点到初始质心的距离平方误差
clusterAssment[j, 1] = distMeas(mat(centroid0), dataMat[j, :])**2
while (len(centList) < k): # 当质心数量小于 k 时
lowestSSE = inf
for i in range(len(centList)): # 对每一个质心
ptsInCurrCluster = dataMat[nonzero(
clusterAssment[:, 0].A == i)[0], :] # 获取当前簇 i 下的所有数据点
centroidMat, splitClustAss = kMeans(
ptsInCurrCluster, 2, distMeas) # 将当前簇 i 进行二分 kMeans 处理
sseSplit = sum(splitClustAss[:, 1]) # 将二分 kMeans 结果中的平方和的距离进行求和
sseNotSplit = sum(
clusterAssment[nonzero(clusterAssment[:, 0].A != i)[0],
1]) # 将未参与二分 kMeans 分配结果中的平方和的距离进行求和
print("sseSplit, and notSplit: ", sseSplit, sseNotSplit)
if (sseSplit + sseNotSplit) < lowestSSE:
bestCentToSplit = i
bestNewCents = centroidMat
bestClustAss = splitClustAss.copy()
lowestSSE = sseSplit + sseNotSplit
# 找出最好的簇分配结果
bestClustAss[nonzero(bestClustAss[:, 0].A == 1)[0], 0] = len(
centList) # 调用二分 kMeans 的结果,默认簇是 0,1. 当然也可以改成其它的数字
bestClustAss[nonzero(bestClustAss[:, 0].A == 0)[0],
0] = bestCentToSplit # 更新为最佳质心
print('the bestCentToSplit is: ', bestCentToSplit)
print('the len of bestClustAss is: ', len(bestClustAss))
# 更新质心列表
centList[bestCentToSplit] = bestNewCents[0, :].tolist()[
0] # 更新原质心 list 中的第 i 个质心为使用二分 kMeans 后 bestNewCents 的第一个质心
centList.append(
bestNewCents[1, :].tolist()[0]) # 添加 bestNewCents 的第二个质心
clusterAssment[nonzero(clusterAssment[:, 0].A == bestCentToSplit)[
0], :] = bestClustAss # 重新分配最好簇下的数据(质心)以及SSE
return mat(centList), clusterAssment
dataMat = mat(loadDataSet('./testSet2.txt'))
centList, myNewAssments = biKMeans(dataMat, 3)
[[-3.58409225 3.86943555]
[ 0.29051049 2.60774022]]
[[-3.14707283 3.35698672]
[ 1.12342729 0.311763 ]]
[[-2.94737575 3.3263781 ]
[ 1.23710375 0.17480612]]
sseSplit, and notSplit: 813.9199514542203 0.0
the bestCentToSplit is: 0
the len of bestClustAss is: 60
[[-3.12034252 1.54045889]
[-3.81393349 3.31713546]]
[[-2.3134116 2.311198 ]
[-3.15869713 3.66477147]]
[[-2.427775 2.50627357]
[-3.22716077 3.76797285]]
[[-2.56363745 2.70664609]
[-3.41638922 4.08382833]]
[[-2.58990523 2.82015346]
[-3.61124957 4.26650957]]
[[-2.56458833 2.9616746 ]
[-4.095738 4.4204886 ]]
sseSplit, and notSplit: 19.867822980685162 805.2680903943892
[[ 4.01124567 -2.89711003]
[ 2.51935237 3.88501664]]
[[-0.45965615 -2.7782156 ]
[ 2.93386365 3.12782785]]
sseSplit, and notSplit: 54.43238789729306 8.651861059831095
the bestCentToSplit is: 1
the len of bestClustAss is: 40
centList
matrix([[-2.94737575, 3.3263781 ],
[-0.45965615, -2.7782156 ],
[ 2.93386365, 3.12782785]])
showCluster(dataMat, 3, centList, myNewAssments)