Bisecting K-Means
什么是二分K-Means
二分K-Means其实就是基于K-Means改进的算法,他的主要核心还是在于K-Means算法中,只不过它的算法思想是先从一个总簇,不断通过二分裂,直到分裂成k个簇则停止。在K-Means博文当中,我们知道经过算法后,返回了2个参数:
- centroids:
返回的是k质心的坐标矩阵 - clusterAssment
返回的是m*2的簇矩阵,其中m是dataSet里数据点的总个数。第0列储存的是数据点被归类到哪个质心簇的索引(对应centroids矩阵), 第1列储存的是与对应质心的误差平方。
而正因为有了与质心的误差平方,通过求和算出误差平方和(SSE),从何判断聚类的整体效果。SSE可以理解为物理热力学上的熵,熵的物理意义是体系混乱程度的度量。
二分K-Means原理
这里的二分裂并不是1 -> 2 -> 4 -> 8,而是1 -> 2 -> 1+2 -> 2+2 > 3+2 -> …
二分算法的核心:
- 初始化整个数据集为1个大簇
- 通过尝试分裂原有的质心点,找出使得SSE(熵)变小的质心点x
- 采取质心点 x 分裂,其他质心点不变的方案
- 将 x 点分裂后2个簇的数据, 更新至将原本 x 点的簇数据
- 更新质心点和总簇
- 直到分裂出 k 个质心
算法优缺点
优点:相对K-Means算法,不会掉落局部极值点,使得更精准了。
缺点:速度更慢了一点。
代码实现
实现方式写了两种
- 代码一
每次计算质心矩阵里的所有质心分裂后的总SSE,选取总SSE最小的方案进行分裂。
优点:可以保证总的SSE是最小的。缺点:因为每次都要对质心进行K-Means,可能速度上会慢一点 - 代码二
每次计算质心矩阵里的所有质心的SSE,选取最大SSE的那簇进行分裂。
优点:相比与代码一的方案,速度可能会快一些,因为并不需要先K-Means再计算SSE,而是直接计算summation值,省去了K-Means的那段时间消耗。
缺点:不能保证总的SSE是最小的
代码(一):找到最小的SSE
from numpy import mat, shape, zeros, mean, sqrt, power, sum, inf, nonzero
from numpy.random import rand
from numpy import array
def distanceEuclid(vecA, vecB):
# 用欧几里得坐标算距离
return sqrt(sum(power(vecA - vecB, 2)))
def randCenter(dataSet, k):
# 生成k个随机质心
# 获取数据的总列数
n = shape(dataSet)[1]
# 初始化一个质心矩阵 k * n k行n列的0矩阵
centroids = mat(zeros((k, n)))
for j in range(n):
minJ = min(dataSet[:, j]) # 计算数据集第j列的最小值
maxJ = max(dataSet[:, j]) # 计算数据集第j列的最大值
rangeJ = float(maxJ - minJ) # 最大值-最小值 计算最大差值
# 生成质心矩阵第j列的值 通过均匀分布的随机数
centroids[:, j] = minJ + rangeJ * rand(k, 1) # 值的范围∈[minJ, maxJ]
return centroids
def k_means(dataSet, k, distMeas=distanceEuclid, createCent=randCenter):
# 将数据集初始化为mat矩阵对象
dataSet = mat(dataSet)
# 获取数据集的总点数
m = shape(dataSet)[0]
# 生成数据集点数的簇矩阵并初始化为0 第一列保存该点被归类的质心索引,第二列保存该点与质心的距离
clusterAssment = mat(zeros((m, 2)))
# 初始化k个质心矩阵
centroids = createCent(dataSet, k)
clusterChanged = True
while clusterChanged:
clusterChanged = False
for i in range(m): # 分别取出全部点 并逐一与所有质心进行距离计算 把距离第j个质心最近时,将该点归类为j
minDist = inf # 先初始化最小值为无限大
minIndex = -1
for j in range(k): # 第i点逐一与所有质心进行距离计算
# 算出第i点与第j个质心的距离
distJI = distMeas(centroids[j, :], dataSet[i, :])
# 与上次的minDist判断 如果更小则更新质心
if distJI < minDist:
minDist = distJI # 质心距离
minIndex = j # 质心索引
# 如果其中第i个点被归类的质心与上次的不一样 则继续下次的计算,否则如果全部点的质心都没有发生变化则认为收敛
if clusterAssment[i, 0] != minIndex:
clusterChanged = True
# 并更新簇矩阵第i个点的记录[质心索引, 与该质心距离的平方]
clusterAssment[i, :] = minIndex, minDist ** 2
# 重新初始化质心
for cent in range(k): # 遍历k个质心
ptsInClust = dataSet[nonzero(clusterAssment[:, 0].A == cent)[0]]
# 被分类到当前第cent个点的所有点的各维度按列求平均值 生成1 * n列的矩阵并更新到质心矩阵中
centroids[cent, :] = mean(ptsInClust, axis=0)
# 返回 质心点 簇矩阵[质心索引, 与该质心距离的平方]
return centroids, clusterAssment
def bisecting_k_means(dataSet, k, distMeas=distanceEuclid):
dataSet = mat(dataSet)
m = shape(dataSet)[0]
# 初始化一个簇
clusterAssment = mat(zeros((m, 2)))
# 计算数据集每列的平均值, 生成初始质心
centroid0 = list(mean(dataSet, axis=0).tolist())[0]
# 质心集合
centList = [centroid0]
# 将数据集所有的点与初始质心进行距离计算
for j in range(m):
clusterAssment[j, 1] = distMeas(mat(centroid0), dataSet[j, :]) ** 2
while len(centList) < k:
bestClustAss, bestNewCents = zeros(1), zeros(1)
bestCentToSplit = -1
lowestSSE = inf
# 遍历质心集合里的全部质心, 找到第i个质心分裂后的SSE与其他质心SSE之和最小的 i点
for i in range(len(centList)):
# 取出第i个质心的数据簇
ptsIncurrCluster = dataSet[nonzero(clusterAssment[:, 0].A == i)[0], :]
# 将第i个质心的数据簇再通过k-means一分为二
centroidMat, splitClustAss = k_means(ptsIncurrCluster, 2, distMeas)
# 计算一分为二的两个子簇后 总的误差平方和SSE
sseSplit = sum(splitClustAss[:, 1])
# 计算未一分为二之前,除了第i个质心以外的数据簇误差平方和
sseNotSplit = sum(clusterAssment[nonzero(clusterAssment[:, 0].A != i)[0], 1])
# 如果第i个点分裂后SSE + 除了i点的SSE和是最小的 则更新最好的
print(f'center {i} sseSplit, notSplit and summation: ', sseSplit, sseNotSplit, sseSplit+sseNotSplit)
if (sseSplit + sseNotSplit) < lowestSSE:
# 最好的分裂点为i
bestCentToSplit = i
# 最好的分裂后的2个质心
bestNewCents = centroidMat
# 分裂后的簇
bestClustAss = splitClustAss.copy()
# 最小的SSE
lowestSSE = sseSplit + sseNotSplit
# 更新总的簇矩阵 即1分2, 2中选出最小SSE分=1+2, 3中选出最小SSE=2+2, ...
# 将分裂后的分裂点1 更新为 质心列表长度, 因为簇矩阵第一列是质心索引由0开始递增,所以直接赋值长度正好能+1
# 所以将分裂后的第二个质心定义为新增的质心
bestClustAss[nonzero(bestClustAss[:, 0].A == 1)[0], 0] = len(centList)
# 将分裂后的第一个质心定义为更新为最大SSE值质心X
bestClustAss[nonzero(bestClustAss[:, 0].A == 0)[0], 0] = bestCentToSplit
print('the bestCentToSplit is :', bestCentToSplit)
print('the len of bestClustAss is :', len(bestClustAss))
print('=========')
print('the SSE summation of clusterAssment is :', sum(clusterAssment[:, 1]))
print('=========')
# 因为分裂后的第一个质心定义为更新为最大SSE值质心X, 所以放回质心X的位置
centList[bestCentToSplit] = bestNewCents[0, :]
# 因为将分裂后的第二个质心定义为新增的质心, 所以放到质心矩阵的最后
centList.append(bestNewCents[1, :])
# 更新变化后簇矩阵 属于第i个质心的数据点
clusterAssment[nonzero(clusterAssment[:, 0].A == bestCentToSplit)[0], :] = bestClustAss
return mat(array(centList)), clusterAssment
测试算法结果(测试例子中k设置为4)
center 0 sseSplit, notSplit and summation: 792.9168565373268 0.0 792.9168565373268
the bestCentToSplit is : 0
the len of bestClustAss is : 80
=========
the SSE summation of clusterAssment is : 1465.5800234838164
=========
center 0 sseSplit, notSplit and summation: 83.5874695564185 326.2840752011824 409.87154475760093
center 1 sseSplit, notSplit and summation: 66.36683512000786 466.63278133614426 532.9996164561521
the bestCentToSplit is : 0
the len of bestClustAss is : 40
=========
the SSE summation of clusterAssment is : 792.9168565373268
=========
center 0 sseSplit, notSplit and summation: 31.62960698022571 358.8853180661336 390.51492504635934
center 1 sseSplit, notSplit and summation: 66.36683512000786 83.5874695564185 149.95430467642637
center 2 sseSplit, notSplit and summation: 18.398749985455712 377.2703018926498 395.6690518781055
the bestCentToSplit is : 1
the len of bestClustAss is : 40
=========
the SSE summation of clusterAssment is : 409.87154475760093
=========
代码(二):每次分裂SSE最大的质心
from numpy import mat, shape, zeros, mean, sqrt, power, sum, inf, nonzero
from numpy.random import rand
from numpy import array
def distanceEuclid(vecA, vecB):
# 用欧几里得坐标算距离
return sqrt(sum(power(vecA - vecB, 2)))
def randCenter(dataSet, k):
# 生成k个随机质心
# 获取数据的总列数
n = shape(dataSet)[1]
# 初始化一个质心矩阵 k * n k行n列的0矩阵
centroids = mat(zeros((k, n)))
for j in range(n):
minJ = min(dataSet[:, j]) # 计算数据集第j列的最小值
maxJ = max(dataSet[:, j]) # 计算数据集第j列的最大值
rangeJ = float(maxJ - minJ) # 最大值-最小值 计算最大差值
# 生成质心矩阵第j列的值 通过均匀分布的随机数
centroids[:, j] = minJ + rangeJ * rand(k, 1) # 值的范围∈[minJ, maxJ]
return centroids
def k_means(dataSet, k, distMeas=distanceEuclid, createCent=randCenter):
# 将数据集初始化为mat矩阵对象
dataSet = mat(dataSet)
# 获取数据集的总点数
m = shape(dataSet)[0]
# 生成数据集点数的簇矩阵并初始化为0 第一列保存该点被归类的质心索引,第二列保存该点与质心的距离
clusterAssment = mat(zeros((m, 2)))
# 初始化k个质心矩阵
centroids = createCent(dataSet, k)
clusterChanged = True
while clusterChanged:
clusterChanged = False
for i in range(m): # 分别取出全部点 并逐一与所有质心进行距离计算 把距离第j个质心最近时,将该点归类为j
minDist = inf # 先初始化最小值为无限大
minIndex = -1
for j in range(k): # 第i点逐一与所有质心进行距离计算
# 算出第i点与第j个质心的距离
distJI = distMeas(centroids[j, :], dataSet[i, :])
# 与上次的minDist判断 如果更小则更新质心
if distJI < minDist:
minDist = distJI # 质心距离
minIndex = j # 质心索引
# 如果其中第i个点被归类的质心与上次的不一样 则继续下次的计算,否则如果全部点的质心都没有发生变化则认为收敛
if clusterAssment[i, 0] != minIndex:
clusterChanged = True
# 并更新簇矩阵第i个点的记录[质心索引, 与该质心距离的平方]
clusterAssment[i, :] = minIndex, minDist ** 2
# 重新初始化质心
for cent in range(k): # 遍历k个质心
ptsInClust = dataSet[nonzero(clusterAssment[:, 0].A == cent)[0]]
# 被分类到当前第cent个点的所有点的各维度按列求平均值 生成1 * n列的矩阵并更新到质心矩阵中
centroids[cent, :] = mean(ptsInClust, axis=0)
# 返回 质心点 簇矩阵[质心索引, 与该质心距离的平方]
return centroids, clusterAssment
def bisecting_k_means(dataSet, k, distMeas=distanceEuclid):
dataSet = mat(dataSet)
m = shape(dataSet)[0]
# 初始化一个簇
clusterAssment = mat(zeros((m, 2)))
# 计算数据集每列的平均值, 生成初始质心
centroid0 = list(mean(dataSet, axis=0).tolist())[0]
# 质心集合
centList = [centroid0]
# 将数据集所有的点与初始质心进行距离计算
for j in range(m):
clusterAssment[j, 1] = distMeas(mat(centroid0), dataSet[j, :]) ** 2
while len(centList) < k:
biggestToSplit = 0
biggestDist = 0
for i in range(len(centList)):
# 计算每个簇的总SSE值
clusterDist = sum(clusterAssment[nonzero(clusterAssment[:, 0] == i)[0], 1])
print(f'center {i} SSE: {clusterDist}')
# 得到最大SSE值得簇得索引
if clusterDist > biggestDist:
biggestDist = clusterDist
biggestToSplit = i
else:
print('============')
print(f'center {biggestToSplit} is biggest:', biggestDist)
print('the SSE summation of clusterAssment is :', sum(clusterAssment[:, 1]))
print('============')
# 取出最大SSE值质心X的数据簇
ptsIncurrCluster = dataSet[nonzero(clusterAssment[:, 0].A == biggestToSplit)[0], :]
# 将质心X的数据簇再通过k-means一分为二
centroidMat, splitClustAss = k_means(ptsIncurrCluster, 2, distMeas)
# 将分裂后的分裂点1 更新为 质心列表长度, 因为簇矩阵第一列是质心索引由0开始递增,所以直接赋值长度正好能+1
# 所以将分裂后的第二个质心定义为新增的质心
splitClustAss[nonzero(splitClustAss[:, 0].A == 1)[0], 0] = len(centList)
# 将分裂后的第一个质心定义为更新为最大SSE值质心X
splitClustAss[nonzero(splitClustAss[:, 0].A == 0)[0], 0] = biggestToSplit
# 因为分裂后的第一个质心定义为更新为最大SSE值质心X, 所以放回质心X的位置
centList[biggestToSplit] = centroidMat[0, :]
# 因为将分裂后的第二个质心定义为新增的质心, 所以放到质心矩阵的最后
centList.append(centroidMat[1, :])
# 更新变化后簇矩阵 属于第i个质心的数据点
clusterAssment[nonzero(clusterAssment[:, 0].A == biggestToSplit)[0], :] = splitClustAss
return mat(array(centList)), clusterAssment
测试算法结果:(测试例子中k设置为7)
center 0 SSE: 1465.5800234838164
============
center 0 is biggest: 1465.5800234838164
the SSE summation of clusterAssment is : 1465.5800234838164
============
center 0 SSE: 326.2840752011824
center 1 SSE: 466.63278133614426
============
center 1 is biggest: 466.63278133614426
the SSE summation of clusterAssment is : 792.9168565373268
============
center 0 SSE: 326.2840752011824
center 1 SSE: 50.986226691467344
center 2 SSE: 32.60124286495115
============
center 0 is biggest: 326.2840752011824
the SSE summation of clusterAssment is : 409.87154475760093
============
center 0 SSE: 30.0483491615835
center 1 SSE: 50.986226691467344
center 2 SSE: 32.60124286495115
center 3 SSE: 36.318485958424354
============
center 1 is biggest: 50.986226691467344
the SSE summation of clusterAssment is : 149.95430467642635
============
center 0 SSE: 30.0483491615835
center 1 SSE: 1.8848153078419991
center 2 SSE: 32.60124286495115
center 3 SSE: 36.318485958424354
center 4 SSE: 26.45378625997894
============
center 3 is biggest: 36.318485958424354
the SSE summation of clusterAssment is : 127.30667955277994
============
center 0 SSE: 30.0483491615835
center 1 SSE: 1.8848153078419991
center 2 SSE: 32.60124286495115
center 3 SSE: 15.838658007214358
center 4 SSE: 26.45378625997894
center 5 SSE: 5.980214281567667
============
center 2 is biggest: 32.60124286495115
the SSE summation of clusterAssment is : 112.8070658831376
============