聚类练习：对地理数据应用二分k-均值算法聚类

编程作业：聚类

实验原理：

参考：

机器学习实战 (豆瓣) (douban.com)

(13条消息) 机器学习实战 第十章 利用K-均值聚类算法对未标注数据分组_无名的博客-CSDN博客

k-means聚类：

k-means聚类将相似的对象归到同一个簇中，每个簇的中心采用簇中所含值的均值计算而成。

优点：容易实现

缺点：可能收敛到局部最小值，在大规模数据上收敛较慢

适用数据类型：数值型数据

算法流程：

    创建k个点作为起始质心（随机选择）
    当任意一个点的簇分配结果发生改变时：
        对数据集中的每个数据点：
         对每个质心：
           计算数据点到质心的距离
       将数据点分配到距其最近的簇
        对每一个簇重新计算其均值作为质心

二分k-means聚类：

k-means的缺点是可能收敛于局部最优解，使用二分k-means聚类算法来解决这个问题。

可以使用SSE(sum of squared error误差平方和)来度量聚类的效果，SSE值越小表示数据点接近于它的质心，聚类效果越好。因为对误差取了平方，因此更加重视那些远离质心的点。

二分k-means聚类算法的思想是先将所有数据点当作一个簇，然后将该簇一分为二。之后在现有的簇中选择一个簇进行划分，选择哪个簇取决于划分哪个簇后能使SSE值最小。不断重复上述过程，直到达到用户要求的簇的个数。

伪代码：

        将所有数据点都看成一个簇
        当簇的数目小于k时：
            初始化lowestSSE = inf
            对于每一个簇：
                对该簇进行k-means聚类(k=2)
                计算聚类后的总误差
                如果小于lowestSSE,则保存聚类后的参数，更新lowestSEE
            选择划分后使得误差值最小的那个簇进行划分

任务：对地理数据应用二分k-均值算法聚类

（1）问题描述：

你的朋友Drew希望你带他去城里庆祝他的生日。由于其他一些朋友也会过来，所以需要你提供一个大家都可行的计划。Drew给了你希望去的69个地址和相应的经纬度。你要决定将这些地方进行聚类的最佳策略，这样可以安排交通工具抵达这些簇的质心，然后步行到每个簇内地址。

（2）具体实现：

1. 读入数据：

   datList = []
   #导入数据
   for line in open('places.txt').readlines():
       lineArr = line.split('\t')
       datList.append([float(lineArr[4]), float(lineArr[3])])
   datMat = mat(datList)
   

2. 计算两个向量之间的距离：

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

3. 随机初始化聚簇中心：

   #随机生成簇中心函数
   def randCent(dataSet, k):
       n = shape(dataSet)[1]
       centroids = mat(zeros((k,n)))
       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
   

4. 球面距离计算：

   #根据经纬度计算球面距离，vecA[0,：]表示A点经纬度
   def distSLC(vecA, vecB):
       a = sin(vecA[0,1]*pi/180) * sin(vecB[0,1]*pi/180)
       b = cos(vecA[0,1]*pi/180) * cos(vecB[0,1]*pi/180) * \
                         cos(pi * (vecB[0,0]-vecA[0,0]) /180)
       return arccos(a + b)*6371.0 
   

5. k-means聚类:

   K-均值是发现给定数据集的k个簇的算法。簇个数k时用户给定的，每个簇通过其质心（centroid），即簇中所有点的中心来描述。

   #Kmeans
   def kMeans(dataSet, k, distMeas=distEclud, createCent=randCent):
       m = shape(dataSet)[0]  #样本数
       clusterAssment = mat(zeros((m,2))) #获得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):
                   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
               #每轮结束后调整簇心
               for cent in range(k):
                   if clusterAssment[i, 0] != minIndex: clusterChanged = True
                   clusterAssment[i, :] = minIndex, minDist ** 2 #计算每列平均值
           return centroids, clusterAssment
   

6. 二分k-means聚类:

   二分K-均值算法（bisecting K-means），它为了克服k-均值算法收敛于局部最小值的问题。该算法首先将所有点作为一个簇，然后将该簇一分为二。之后选择一个簇继续划分，选择哪个簇进行划分取决于对其划分是否可以最大程度降低SSE的值。上述基于SSE的划分过程不断重复，直到得到指定的簇数目为止。

   代码如下：

   # k:簇个数   distMeas：距离生成器
   def biKmeans(dataSet, k, distMeas=distEclud):
       m = shape(dataSet)[0]  #数据集矩阵的行数
       clusterAssment = mat(zeros((m,2)))
       centroid0 = mean(dataSet, axis=0).tolist()[0]   # 创建一个初始簇
       centList =[centroid0] #create a list with one centroid
       for j in range(m):  # 计算每个样本点到初始簇的距离
           clusterAssment[j,1] = distMeas(mat(centroid0), dataSet[j,:])**2
       while (len(centList) < k):  # 迭代直到簇的数目等于k
           lowestSSE = inf
           for i in range(len(centList)):
               # 尝试划分每一簇
               ptsInCurrCluster = dataSet[nonzero(clusterAssment[:,0].A==i)[0],:]
               centroidMat, splitClustAss = kMeans(ptsInCurrCluster, 2, distMeas)
               sseSplit = sum(splitClustAss[:,1])
               # 剩余数据集的误差
               sseNotSplit = sum(clusterAssment[nonzero(clusterAssment[:,0].A!=i)[0],1])
               # 记录总误差最小的那个簇
               if (sseSplit + sseNotSplit) < lowestSSE:
                   bestCentToSplit = i
                   bestNewCents = centroidMat
                   bestClustAss = splitClustAss.copy()
                   lowestSSE = sseSplit + sseNotSplit
           # 因为二分均值，所以结果簇编号只能是0或1，修改结果簇编号为：原簇号和新增簇号
           bestClustAss[nonzero(bestClustAss[:,0].A == 1)[0],0] = len(centList)
           bestClustAss[nonzero(bestClustAss[:,0].A == 0)[0],0] = bestCentToSplit
           # 更新簇列表centlist和样本点分配簇结果矩阵clusterAssment
           centList[bestCentToSplit] = bestNewCents[0,:].tolist()[0]
           centList.append(bestNewCents[1,:].tolist()[0])
           clusterAssment[nonzero(clusterAssment[:,0].A == bestCentToSplit)[0],:]= bestClustAss
       return mat(centList), clusterAssment;
   
   

（3）实验结果：

1. k-means聚类：

   from matplotlib import pyplot as plt
   import Kmeansch10 as km
   
   data_X = km.loadDataSet("places.txt")
   centroids, clusterAssment = km.kMeans(data_X, 4)
   
   plt.figure()
   plt.scatter(data_X[:, 0].flatten().tolist(), data_X[:, 1].flatten().tolist(), c="b", marker="o")
   plt.scatter(centroids[:, 0].flatten().tolist(), centroids[:, 1].flatten().tolist(), c='r', marker="+")
   plt.show()
   

   发现，在经过多次测试后，会出现聚簇收敛到局部最小值，如下所示：

   

2. 二分k-means聚类：

   k-means的缺点是可能收敛于局部最优解，使用二分k-means聚类算法来解决这个问题。

   def clusterClubs(numClust=5):
       datList = []
       #导入数据
       for line in open('places.txt').readlines():
           lineArr = line.split('\t')
           datList.append([float(lineArr[4]), float(lineArr[3])])
       datMat = mat(datList)
       #采用二分k-均值算法进行聚类
       myCentroids, clustAssing = biKmeans(datMat, numClust, distMeas=distSLC)
       #定义画布，背景
       fig = plt.figure()
       rect=[0.0,0.0,1.0,1.0]
       #不同图形标识
       scatterMarkers=['s', 'o', '^', '8', 'p', \
                       'd', 'v', 'h', '>', '<']
       axprops = dict(xticks=[], yticks=[])
       ax0=fig.add_axes(rect, label='ax0', **axprops)
       #导入地图
       imgP = plt.imread('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)
       #采用‘+’表示簇中心
       ax1.scatter(myCentroids[:,0].flatten().A[0], myCentroids[:,1].flatten().A[0], marker='+', s=300)
       plt.show()
   

   结果如下：

   • K=2时：

     

     

   • K=3时：

     

     

   • K=4时：

     

     

   • K=5时：

     

     

实验小结：

聚类算法时一种无监督的学习方法。所谓无监督学习是指实现不知道要寻找的内容，即没有目标变量。聚类将数据点归到多个簇中，其中相似的数据点处于同一簇，而不相似的数据点处于不同簇中。

同时，学习到关于聚类中K选取方法的经验：

• 分析数据，人工判断；
• 基于变化的算法：即定义一个函数，随着K的改变，认为在正确的K时会产生极值；
• 基于变化的算法：即定义一个函数，随着K的改变，认为在正确的K时会产生极值；
• 基于结构的算法：即比较类内距离、类间距离以确定K；
• 基于一致性矩阵的算法：即认为在正确的K时，不同次聚类的结果会更加相似；
• 基于层次聚类：即基于合并或分裂的思想，在一定情况下停止从而获得K；
• 基于采样的算法：即对样本采样，分别做聚类；根据这些结果的相似性确定K；
• 使用Canopy Method算法进行初始划分；
• 使用BIC算法进行初始划分。