读书笔记 - 机器学习实战 - 10 k-均值聚类

10 $k$-均值聚类（Grouping unlabeled items using k-means clustering）

聚类（clustering）：一种无监督学习（unsupervised learning）算法，自动生成相似样本簇（cluster）。

$k$-均值（$k$-means）：生成$k$个簇，各簇中心为簇内样本均值。

聚类也称非监督分类（unsupervised classification），其输出与分类（classification）相同，但没有预定义类别。

聚类分析将相似（similar）样本归类到同一簇中，非相似（dissimilar）样本归类到不同簇中。相似性（similarity）取决于相似性测度（similarity measurement）。

import matplotlib.pyplot as plt
import numpy as np

random_seed = 13
np.random.seed(seed=random_seed)

10.1 $k$-均值聚类算法

$k$-均值聚类（$k$-means clustering）：

优点：易于使用

缺点：收敛于局部极小值（converge at local minima）；大规模数据集上运行很慢

适用范围：数值

给定数据集，$k$-均值从中寻找$k$个簇，$k$为用户定义的超参数。各簇用簇中心（centroid）表示。

$k$-均值伪代码：

（随机）创建k个初始中心
当任意样本点的簇类别改变时：
    遍历数据集：
        遍历簇中心：
            计算簇中心与样本点距离
        将样本点分配距离最近的簇类别
    遍历所有簇：
        计算各簇均值，将簇均值分配给簇中心

$k$-均值聚类步骤：

1. 收集数据
   
2. 准备：数值型数据，标称值需映射为二进制数值
   
3. 分析：
   
4. 训练：
   
5. 测试：应用聚类算法、检查结果，测量定量误差
   
6. 使用：
   

# Listing 10.1 k-means support functions

def loadDataSet(fileName):
    
    dataMat = []
    with open(fileName, "r") as fr:
        for line in fr.readlines():
            curLine = line.strip().split("\t")
            fltLine = list(map(float, curLine))
            dataMat.append(fltLine)
            
    return dataMat


def distEclud(vecA, vecB):
    
    return np.sqrt(np.sum(np.power(vecA - vecB, 2)))


def randCent(dataSet, k):
    
    n = np.shape(dataSet)[1]
    centroids = np.matrix(np.zeros((k, n)))
    
    for j in range(n):
        # create cluster centroids
        minJ = np.min(dataSet[:, j])
        rangeJ = float(np.max(dataSet[:, j]) - minJ)
        centroids[:, j] = minJ + rangeJ * np.random.rand(k, 1)
        
    return centroids

datMat = np.matrix(loadDataSet("./testSet.txt"))
# 最大、最小值
print(np.min(datMat[:, 0]))
print(np.min(datMat[:, 1]))
print(np.max(datMat[:, 0]))
print(np.max(datMat[:, 1]))

# 随机中心
print(randCent(datMat, 2))

# 距离测度
print(distEclud(datMat[0], datMat[1]))

fig = plt.figure(figsize=(8, 6))
ax = fig.add_subplot(111)
ax.scatter(datMat[:, 0].A.ravel(),
           datMat[:, 1].A.ravel(),
           label="samples")
plt.legend(loc="best")
plt.show()

-5.379713
-4.232586
4.838138
5.1904
[[ 2.56673435  3.53457907]
 [-2.95255221  4.86765517]]
5.184632816681332

# Listing 10.2 the k-means clustering algorithm
def kMeans(dataSet, k, distMeas=distEclud, createCent=randCent):
    
    m = np.shape(dataSet)[0]
    # cluster assignment, col_0: cluster index, col_1: dist^2 
    clusterAssment = np.matrix(np.zeros((m, 2)))
    centroids = randCent(dataSet, k)
    clusterChanged = True
    
    while clusterChanged:
        clusterChanged = False
        
        for i in range(m):
            minDist = np.inf
            minIndex = -1
            
            # find the closest centroid
            for j in range(k):
                distJI = distMeas(centroids[j, :], dataSet[i, :])
                if distJI < minDist:
                    minDist = distJI
                    minIndex = j
                    
            if clusterAssment[i, 0] != minIndex:
                clusterChanged = True
                clusterAssment[i, :] = minIndex, minDist ** 2
                
        print(centroids)
        
        # update centroid location
        for cent in range(k):
            ptsInClust = dataSet[np.nonzero(clusterAssment[:, 0].A == cent)[0]]
            centroids[cent, :] = np.mean(ptsInClust, axis=0)
            
    return centroids, clusterAssment

myCentroids, clustAssing =  kMeans(datMat, 4)

fig = plt.figure(figsize=(8, 6))
ax = fig.add_subplot(111)
ax.scatter(datMat[:, 0].A.ravel(),
           datMat[:, 1].A.ravel(),
           c=clustAssing[:, 0].A.ravel(),
           label="samples")
ax.scatter(myCentroids[:, 0].A.ravel(),
           myCentroids[:, 1].A.ravel(),
           s=80,
           marker="x",
           c="red",
           label="centroids")
plt.legend(loc="best")
plt.show()

[[ 4.55818026  1.81332757]
 [-0.74643615  2.57098167]
 [ 0.84339214 -3.90243732]
 [ 2.54450137 -1.42030081]]
[[ 2.942346    2.80047613]
 [-1.6334182   3.03655888]
 [-1.76738661 -3.13626093]
 [ 3.15816942 -2.14680933]]
[[ 2.6265299   3.10868015]
 [-2.46154315  2.78737555]
 [-3.01169468 -3.01238673]
 [ 3.03713839 -2.62802833]]
[[ 2.6265299   3.10868015]
 [-2.46154315  2.78737555]
 [-3.38237045 -2.9473363 ]
 [ 2.80293085 -2.7315146 ]]

10.2 后处理（Improving cluster performance with postprocessing）

$k$-均值聚类收敛到局部极小值点，导致聚类结果不正确：

平方误差和（sum of squared error，SSE）

簇分裂（cluster-splitting）降低SSE：选取SSE最大的簇进行分裂（对该簇内的点进行2-均值聚类），然后再选取2个簇合并。

datMat3 = np.matrix(loadDataSet("testSet2.txt"))

fig = plt.figure(figsize=(8, 6))
ax = fig.add_subplot(111)
ax.scatter(datMat3[:, 0].A.ravel(),
           datMat3[:, 1].A.ravel(),
           label="samples")
plt.legend(loc="best")
plt.show()

10.3 二分$k$-均值（Bisecting $k$-means）

二分$k$-均值（bisecting $k$-means）伪代码：

初始化，将所有样本点聚类到同一簇
当簇总数小于k时：
    遍历簇：
        计算总误差
        对当前簇进行2-均值聚类，计算总误差
    选择使SSE最小的簇分裂

# Listing 10.3 the bisecting k-means clustering algorithm

def biKmeans(dataSet, k, distMeas=distEclud):
    
    m = np.shape(dataSet)[0]
    clusterAssment = np.matrix(np.zeros((m, 2)))
    
    # initially create one cluster
    centroid0 = np.mean(dataSet, axis=0).tolist()[0]
    centList = [centroid0]
    
    for j in range(m):
        clusterAssment[j, 1] = distMeas(np.matrix(centroid0), dataSet[j, :]) ** 2
        
    while (len(centList) < k):
        lowestSSE = np.inf
        for i in range(len(centList)):
            # try splitting every cluster
            ptsInCurrCluster = dataSet[np.nonzero(clusterAssment[:, 0].A == i)[0], :]
            centroidMat, splitClustAss = kMeans(ptsInCurrCluster, 2, distMeas)
            
            sseSplit = np.sum(splitClustAss[:, 1])
            sseNotSplit = np.sum(clusterAssment[np.nonzero(clusterAssment[:, 0].A != i)[0], 1])
            
            print("sseSplit: {}, and notSplit: {}".format(sseSplit, sseNotSplit))
            
            if (sseSplit + sseNotSplit) < lowestSSE:
                bestCentToSplit = i
                bestNewCents = centroidMat
                bestClustAss = splitClustAss.copy()
                lowestSSE = sseSplit + sseNotSplit
                
        # update the cluster assignments
        bestClustAss[np.nonzero(bestClustAss[:, 0].A == 1)[0], 0] = len(centList)
        bestClustAss[np.nonzero(bestClustAss[:, 0].A == 0)[0], 0] = bestCentToSplit
        print("the bestCentToSplit is: {}".format(bestCentToSplit))
        print("the len of bestClustAss is: {}".format(len(bestClustAss)))
        
        print(centList)
        centList[bestCentToSplit] = bestNewCents[0, :].A.ravel()
        centList.append(bestNewCents[1, :].A.ravel())
        clusterAssment[np.nonzero(clusterAssment[:, 0].A == bestCentToSplit)[0], :] = bestClustAss
        
    return np.matrix(centList), clusterAssment

centList, myNewAssments = biKmeans(datMat3, 3)

fig = plt.figure(figsize=(8, 6))
ax = fig.add_subplot(111)
ax.scatter(datMat3[:, 0].A.ravel(),
           datMat3[:, 1].A.ravel(),
           c=myNewAssments[:, 0].A.ravel(),
           label="samples")
ax.scatter(centList[:, 0].A.ravel(),
           centList[:, 1].A.ravel(),
           s=80,
           marker="x",
           c="red",
           label="centroids")
plt.legend(loc="best")
plt.show()

[[ 2.46075605  2.05403444]
 [-4.18161624  0.80409831]]
[[ 2.29222778  1.72074989]
 [-2.16222773  0.81998667]]
[[ 2.93386365  3.12782785]
 [-1.70351595  0.27408125]]
sseSplit: 656.6855191745456, and notSplit: 0.0
the bestCentToSplit is: 0
the len of bestClustAss is: 60
[[-0.15772275000000002, 1.2253301166666664]]
[[1.62610908 1.52989103]
 [1.49043987 1.11553377]]
[[2.95977168 3.26903847]
 [2.441611   0.444826  ]]
sseSplit: 1.3545754399028473, and notSplit: 656.6855191745456
[[ 0.42631551 -2.01274686]
 [-0.96733016  0.99114804]]
[[-0.45965615 -2.7782156 ]
 [-2.94737575  3.3263781 ]]
sseSplit: 225.54055003224968, and notSplit: 0.0
the bestCentToSplit is: 1
the len of bestClustAss is: 40
[array([2.93386365, 3.12782785]), array([-1.70351595,  0.27408125])]