机器学习实战（十）——利用K-均值聚类算法对未标注数据分组

机器学习实战（十）——利用K-均值聚类算法对未标注数据分组

聚类是一种无监督的学习，即对于训练数据来说并没有已知标签，我们是根据样本之间的相似程度将其划分为多个类。

一、K-均值聚类算法

K-均值算法可以通过计算样本之间的相似程度，将训练集划分为K个类别。这里类别数K是用户指定的。每一个类别通过其质心（类别中所有点的中心）来描述。其特点是思路和算法都很简单，但是比较容易收敛到局部最小值，并在在大规模的数据集上收敛速度较慢

用伪代码来描述可以更好地理解该算法流程：

创建k个点作为起始质心（通常是随机选择）
当任意一个点的类别分配结果发生改变时：
	遍历数据集中的每一个数据点：
		遍历每一个质心：
			计算当前数据点和这个质心的距离
		将当前数据点分配到最近的质心，即最近的类别之中
	遍历完每一个数据点分类后重新计算每一个类别的质心（取类别中所有点的均值）

K-均值算法所需要的支持函数为：

from numpy import *

def loadDataSet(fileName):
    dataMat = []
    fr = open(fileName)
    for line in fr.readlines():
        curLine = line.strip().split('\t')
        fltLine = list(map(float,curLine))  # 这里原文错了
        dataMat.append(fltLine)
    return dataMat

def distEclud(vecA,vecB):
    disTance = sqrt(sum(power(vecA-vecB,2)))
    return disTance

# 用于初始化聚类中心的函数
def randCent(dataSet,k):
    n = shape(dataSet)[1]
    centroids = mat(zeros((k,n)))  # 每一行都是一个质心的坐标，有n维度
    for j in range(n):
        minJ = min(dataSet[:,j])
        rangeJ = float(max(dataSet[:,j]) - minJ)
        centroids[:,j] = minJ + rangeJ * random.rand(k,1)
        # 随机构建k个聚类中心，每一个中心的每一个坐标都不会超过已有样本的范围
    return centroids

注意：如果直接按照原文的代码，会报错，错误为：

TypeError: unsupported operand type(s) for -: 'map' and 'map'

因为在python3中，map函数返回的类型是该对象的内存地址，不是我们想要的列表类型，因此需要加上list，即：

fltLine = list(map(float,curLine))

接下来为具体的K-均值算法函数：

def kMeans(dataSet, k, distMeas = distEclud, createCent = randCent):
    # 与之前的类似，后两个参数为函数的引用，当想用另外的方式计算距离或者初始化
    # 就改这两个参数即可
    m = shape(dataSet)[0]
    clusterAssment = mat(zeros((m,2)))  # 用来存储样本点所属的类和距离
    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):
                distIJ = distMeas(centroids[j,:],dataSet[i,:])  # 计算距离
                if distIJ < minDist:  # 更新距离
                    minDist = distIJ
                    minIndex = j
            if clusterAssment[i,0] != minIndex:  # 如果不相等则说明发生了更新
                clusterChanged = True
                clusterAssment[i,:] = minIndex,minDist**2
        print(centroids)
        for cent in range(k):  # 更新质心的位置
            # 先找出属于当前类的样本点
            ptsInClust = dataSet[nonzero(clusterAssment[:,0].A == cent)[0]]
            centroids[cent,:] = mean(ptsInClust,axis = 0)  # 按照平均值来更新
    return centroids,clusterAssment

二、使用后处理来提高聚类性能

由于上述K-均值的过程中聚类数目k是需要用户指定的，那么问题就是什么样的k才是最好的，以及如何评价。

思路就在于之前分类的时候存储了每一个样本点到其对应质心的距离。一种用于度量聚类效果的指标为SSE（误差平方和），就是类别种各个样本点距离类别质心的距离之和。

降低SSE的方法最简单的当然是增加划分的类别数，例如可以将SSE最大的类划分为多个类。但为了在保持类别数不变的情况下降低SSE，可以划分为较多个类后再对类进行合并直到满足数目要求，例如将最近的质心合并，或者合并两个使得SSE增幅最小的质心。

三、二分 K-均值算法

该算法提出是为了克服K-均值收敛于局部最小值的问题。具体的思路为：

首先将所有点作为一个类别，然后将该类别按照K-均值的分类方法分成两个类

之后再选在其中一个类别继续划分，直到满足指定的类别数目为止

而如何选择划分的类别，取决于对哪个类别划分能够最大程度的降低整体SSE的值

用伪代码的形式加以理解：

将所有点看成一个类别
当类别数目小于k时：
	对于每一个类别：
		计算总误差
		对当前类进行K-均值聚类（k=2)
		计算划分后的总误差
	选择使得误差减少最多的那个类进行划分

具体代码为：

def biKmeans(dataSet, k, distMeas = distEclud):
    m = shape(dataSet)[0]
    clusterAssment = mat(zeros((m,2)))
    centroid0 = 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):  # 当前类别数目未达到要求k就循环
        print("the len of centList is :",len(centList))
        lowestSSE = inf
        for i in range(len(centList)):  # 对当前已有的类
            # 找到属于当前类的样本
            ptsInCurrCluster = dataSet[nonzero(clusterAssment[:,0].A == i)[0],:]
            # 对当前类进行二分
            centroidMat,splitClustAss = kMeans(ptsInCurrCluster,2,distMeas)
            seeSplit = sum(splitClustAss[:,1])  # 对划分之后的当前的误差求和
            # 计算其他类别的总误差，可以和这个划分的求和来计算样本的总误差
            seeNoSplit = sum(clusterAssment[nonzero(clusterAssment[:,0].A != i)[0],1])
            print("seeSplit and seeNoSplit:",seeSplit,seeNoSplit)
            if (seeSplit + seeNoSplit) < lowestSSE:
                bestCentToSplit = i  # 记录最好的那个划分类别
                bestNewCents = centroidMat  # 记录划分后的两个新类的质心
                bestClustAss = splitClustAss.copy()  # 记录类内样本的划分情况
                lowestSSE = seeSplit + seeNoSplit
        # 找到了那个最适合的类来划分了
        # 以下两个只是对新划分出来的类更新数字，因此划分为01，但现在已有了其他类
        bestClustAss[nonzero(bestClustAss[:,0].A == 1)[0],0] = len(centList)  # 作为新的类
        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,:]  # 更新原来的类的质心
        centList.append(bestNewCents[1,:])  # 增加新类的质心
        # 原来属于划分类的那些拿出来，令它们等于新的划分结果，相当于对他们更新
        clusterAssment[nonzero(clusterAssment[:,0].A == bestCentToSplit)[0],:] = bestClustAss
    centList = list(map(lambda x: [x[0], x[1]], [matrix.tolist(i)[0] for i in centList]))
    return mat(centList),clusterAssment

注意：如果按照原文代码，将会有以下错误：

发生在：

return mat(centList),clusterAssment

错误为：

ValueError: matrix must be 2-dimensional

将centList打印出来发现：

[matrix([[2.93386365, 3.12782785]]), matrix([[-0.45965615, -2.7782156 ]]), matrix([[-2.94737575,  3.3263781 ]])]

就是说这个只有一维，但是矩阵要求是二维的。那么可以进行如下处理：

centList = list(map(lambda x: [x[0], x[1]], [matrix.tolist(i)[0] for i in centList]))

分析：

>>> A = [np.matrix([[-1.70351595,  0.27408125]]), np.matrix([[2.836743 , 3.6383055]]), np.matrix([[3.1604785 , 1.93671333]])]
>>> [np.matrix.tolist(i)[0] for i in A]
Out[10]: [[-1.70351595, 0.27408125], [2.836743, 3.6383055], [3.1604785, 1.93671333]]
>>> map(lambda x: [(x[0]), x[1]], [np.matrix.tolist(i)[0] for i in A])
Out[11]: <map at 0x22c6ba6e5c0>
>>> list(map(lambda x: [(x[0]), x[1]], [np.matrix.tolist(i)[0] for i in A]))
Out[12]: [[-1, 0.27408125], [2, 3.6383055], [3, 1.93671333]]

$labmda$是匿名函数，就是创建一个匿名函数，其输入参数为x，返回值为$[x[0],x[1]]$，然后map即是将这个匿名函数应用于第二个参数的每一个元素中。

四、示例：对地图上的点进行聚类

从网站获取经纬度的代码就跳过了，因为文件places.txt中的最后两列就是经纬度信息，直接利用这些进行聚类

都是画图的代码，就不做解释了

def distSLC(vecA, vecB):  # 根据经纬度算地球表面距离
    a = sin(vecA[0, 1] * math.pi / 180) * sin(vecB[0, 1] * math.pi / 180)
    b = cos(vecA[0, 1] * math.pi / 180) * cos(vecB[0, 1] * math.pi / 180) * cos(
        math.pi * (vecB[0, 0] - vecA[0, 0]) / 180)
    return arccos(a + b) * 6371.0


def clusterClubs(numClust=5):
    datList = []
    for line in open('D:\学习\大四上学习\Python机器学习代码实现\MLiA_SourceCode\machinelearninginaction\Ch10/places.txt').readlines():
        lineArr = line.split('\t')
        datList.append([float(lineArr[4]), float(lineArr[3])])
    datMat = mat(datList)
    myCentroids, clustAssing = biKmeans(datMat, numClust, distMeas=distSLC)
    fig = plt.figure()
    rect = [0.1, 0.1, 0.8, 0.8]
    scatterMarkers = ['s', 'o', '^', '8', 'p',
                      'd', 'v', 'h', '>', '<']
    axprops = dict(xticks=[], yticks=[])
    ax0 = fig.add_axes(rect, label='ax0', **axprops)
    imgP = plt.imread('D:\学习\大四上学习\Python机器学习代码实现\MLiA_SourceCode\machinelearninginaction\Ch10/Portland.png')
    ax0.imshow(imgP)
    ax1 = fig.add_axes(rect, label='ax1', frameon=False)
    for i in range(numClust):
        ptsInCurrCluster = datMat[nonzero(clustAssing[:, 0].A == i)[0], :]
        markerStyle = scatterMarkers[i % len(scatterMarkers)]
        ax1.scatter(ptsInCurrCluster[:, 0].flatten().A[0], ptsInCurrCluster[:, 1].flatten().A[0], marker=markerStyle,
                    s=90)
    print(myCentroids)
    ax1.scatter(myCentroids[:, 0].flatten().A[0], myCentroids[:, 1].flatten().A[0], marker='+', s=300)
    plt.show()