1. 程式人生 > >深入理解K-Means聚類演算法

深入理解K-Means聚類演算法

概述

什麼是聚類分析

聚類分析是在資料中發現數據物件之間的關係,將資料進行分組,組內的相似性越大,組間的差別越大,則聚類效果越好。

不同的簇型別

聚類旨在發現有用的物件簇,在現實中我們用到很多的簇的型別,使用不同的簇型別劃分資料的結果是不同的,如下的幾種簇型別。

這裡寫圖片描述

明顯分離的

可以看到(a)中不同組中任意兩點之間的距離都大於組內任意兩點之間的距離,明顯分離的簇不一定是球形的,可以具有任意的形狀。

基於原型的

簇是物件的集合,其中每個物件到定義該簇的原型的距離比其他簇的原型距離更近,如(b)所示的原型即為中心點,在一個簇中的資料到其中心點比到另一個簇的中心點更近。這是一種常見的基於中心的簇

,最常用的K-Means就是這樣的一種簇型別。
這樣的簇趨向於球形。

基於密度的

簇是物件的密度區域,(d)所示的是基於密度的簇,當簇不規則或相互盤繞,並且有早上和離群點事,常常使用基於密度的簇定義。

關於更多的簇介紹參考《資料探勘導論》。

基本的聚類分析演算法

1. K均值:
基於原型的、劃分的距離技術,它試圖發現使用者指定個數(K)的簇。

2. 凝聚的層次距離:
思想是開始時,每個點都作為一個單點簇,然後,重複的合併兩個最靠近的簇,直到嘗試單個、包含所有點的簇。

3. DBSCAN:
一種基於密度的劃分距離的演算法,簇的個數有演算法自動的確定,低密度中的點被視為噪聲而忽略,因此其不產生完全聚類。

距離量度

不同的距離量度會對距離的結果產生影響,常見的距離量度如下所示:

這裡寫圖片描述

這裡寫圖片描述

K-Means演算法

下面介紹K均值演算法:

優點:易於實現
缺點:可能收斂於區域性最小值,在大規模資料收斂慢

演算法思想較為簡單如下所示:

選擇K個點作為初始質心  
repeat  
    將每個點指派到最近的質心,形成K個簇  
    重新計算每個簇的質心  
until 簇不發生變化或達到最大迭代次數  

這裡的重新計算每個簇的質心,如何計算的是根據目標函式得來的,因此在開始時我們要考慮距離度量和目標函式。

考慮歐幾里得距離的資料,使用誤差平方和(Sum of the Squared Error,SSE)

作為聚類的目標函式,兩次執行K均值產生的兩個不同的簇集,我們更喜歡SSE最小的那個。

這裡寫圖片描述

k表示k個聚類中心,ci表示第幾個中心,dist表示的是歐幾里得距離。
這裡有一個問題就是為什麼,我們更新質心是讓所有的點的平均值,這裡就是SSE所決定的。

這裡寫圖片描述

下面用Python進行實現

# dataSet樣本點,k 簇的個數
# disMeas距離量度,預設為歐幾里得距離
# createCent,初始點的選取
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
            # 第1列為所屬質心,第2列為距離
            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

重點理解一下:

  for cent in range(k):
      ptsInClust = dataSet[nonzero(clusterAssment[:,0].A==cent)[0]]
      centroids[cent,:] = mean(ptsInClust, axis=0) 

迴圈每一個質心,找到屬於當前質心的所有點,然後根據這些點去更新當前的質心。
nonzero()返回的是一個二維的陣列,其表示非0的元素位置。

>>> from numpy import *
>>> a=array([[1,0,0],[0,1,2],[2,0,0]])
>>> a
array([[1, 0, 0],
       [0, 1, 2],
       [2, 0, 0]])
>>> nonzero(a)
(array([0, 1, 1, 2]), array([0, 1, 2, 0]))

表示第[0,0],[1,1] … 位非零元素。第一個陣列為行,第二個陣列為列,兩者進行組合得到的。

ptsInClust = dataSet[nonzero(clusterAssment[:,0].A==cent)[0]]
因此首先先比較clusterAssment[:,0].A==cent的真假,如果為真則記錄了他所在的行,因此在用切片進行取值。

一些輔助的函式:

def loadDataSet(fileName):      #general function to parse tab -delimited floats
    dataMat = []                #assume last column is target value
    fr = open(fileName)
    for line in fr.readlines():
        curLine = line.strip().split('\t')
        fltLine = map(float,curLine) #map all elements to float()
        dataMat.append(fltLine)
    return dataMat

def distEclud(vecA, vecB):
    return sqrt(sum(power(vecA - vecB, 2))) #la.norm(vecA-vecB)

def randCent(dataSet, k):
    n = shape(dataSet)[1]
    centroids = mat(zeros((k,n)))#create centroid mat
    for j in range(n):#create random cluster centers, within bounds of each dimension
        minJ = min(dataSet[:,j]) 
        rangeJ = float(max(dataSet[:,j]) - minJ)
        centroids[:,j] = mat(minJ + rangeJ * random.rand(k,1))
    return centroids

執行和結果

將上述程式碼寫到kMeans.py中,然後開啟python互動端。

>>> from numpy import *
>>> import kMeans
>>> dat=mat(kMeans.loadDataSet('testSet.txt')) #讀入資料
>>> center,clust=kMeans.kMeans(dat,4)
[[ 0.90796996  5.05836784]
 [-2.88425582  0.01687006]
 [-3.3447423  -1.01730512]
 [-0.32810867  0.48063528]]
[[ 1.90508653  3.530091  ]
 [-3.00984169  2.66771831]
 [-3.38237045 -2.9473363 ]
 [ 2.22463036 -1.37361589]]
[[ 2.54391447  3.21299611]
 [-2.46154315  2.78737555]
 [-3.38237045 -2.9473363 ]
 [ 2.8692781  -2.54779119]]
[[ 2.6265299   3.10868015]
 [-2.46154315  2.78737555]
 [-3.38237045 -2.9473363 ]
 [ 2.80293085 -2.7315146 ]]
# 作圖
>>>kMeans(dat,center)

這裡寫圖片描述

繪圖的程式如下:

def draw(data,center):
    length=len(center)
    fig=plt.figure
    # 繪製原始資料的散點圖
    plt.scatter(data[:,0],data[:,1],s=25,alpha=0.4)
    # 繪製簇的質心點
    for i in range(length):
        plt.annotate('center',xy=(center[i,0],center[i,1]),xytext=\
        (center[i,0]+1,center[i,1]+1),arrowprops=dict(facecolor='red'))
    plt.show()

K-Means演算法的缺陷

k均值演算法非常簡單且使用廣泛,但是其有主要的兩個缺陷:
1. K值需要預先給定,屬於預先知識,很多情況下K值的估計是非常困難的,對於像計算全部微信使用者的交往圈這樣的場景就完全的沒辦法用K-Means進行。對於可以確定K值不會太大但不明確精確的K值的場景,可以進行迭代運算,然後找出Cost Function最小時所對應的K值,這個值往往能較好的描述有多少個簇類。
2. K-Means演算法對初始選取的聚類中心點是敏感的,不同的隨機種子點得到的聚類結果完全不同
3. K均值演算法並不是很所有的資料型別。它不能處理非球形簇、不同尺寸和不同密度的簇,銀冠指定足夠大的簇的個數是他通常可以發現純子簇。
4. 對離群點的資料進行聚類時,K均值也有問題,這種情況下,離群點檢測和刪除有很大的幫助。

下面對初始質心的選擇進行討論:

拙劣的初始質心

當初始質心是隨機的進行初始化的時候,K均值的每次執行將會產生不同的SSE,而且隨機的選擇初始質心結果可能很糟糕,可能只能得到區域性的最優解,而無法得到全域性的最優解。如下圖所示:

這裡寫圖片描述
可以看到程式迭代了4次終止,其得到了區域性的最優解,顯然我們可以看到其不是全域性最優的,我們仍然可以找到一個更小的SSE的聚類。

隨機初始化的侷限

你可能會想到:多次執行,每次使用一組不同的隨機初始質心,然後選擇一個具有最小的SSE的簇集。該策略非常的簡單,但是效果可能不是很好,這取決於資料集合尋找的簇的個數。

關於更多,參考《資料探勘導論》

K-Means優化演算法

為了克服K-Means演算法收斂於區域性最小值的問題,提出了一種二分K-均值(bisecting K-means)

bisecting K-means

演算法的虛擬碼如下:

將所有的點看成是一個簇
當簇小於數目k時
    對於每一個簇
        計算總誤差
        在給定的簇上進行K-均值聚類,k值為2
        計算將該簇劃分成兩個簇後總誤差
    選擇是的誤差最小的那個簇進行劃分

完整的Python程式碼如下:

def biKmeans(dataSet, k, distMeas=distEclud):
    m = shape(dataSet)[0]
    # 這裡第一列為類別,第二列為SSE
    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):
        lowestSSE = inf
        # 對於每一個質心,嘗試的進行劃分
        for i in range(len(centList)):
            # 得到屬於該質心的資料
            ptsInCurrCluster =\ dataSet[nonzero(clusterAssment[:,0].A==i)[0],:]
            # 對該質心劃分成兩類
            centroidMat, splitClustAss = kMeans(ptsInCurrCluster, 2, distMeas)
            # 計算該簇劃分後的SSE
            sseSplit = sum(splitClustAss[:,1])
            # 沒有參與劃分的簇的SSE
            sseNotSplit = sum(clusterAssment[nonzero(clusterAssment[:,0].A!=i)[0],1])
            print "sseSplit, and notSplit: ",sseSplit,sseNotSplit
            # 尋找最小的SSE進行劃分
            # 即對哪一個簇進行劃分後SSE最小
            if (sseSplit + sseNotSplit) < lowestSSE:
                bestCentToSplit = i
                bestNewCents = centroidMat
                bestClustAss = splitClustAss.copy()
                lowestSSE = sseSplit + sseNotSplit

        # 較難理解的部分
        bestClustAss[nonzero(bestClustAss[:,0].A == 1)[0],0] = len(centList) #change 1 to 3,4, or whatever
        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]#replace a centroid with two best centroids 
        centList.append(bestNewCents[1,:].tolist()[0])
        clusterAssment[nonzero(clusterAssment[:,0].A == bestCentToSplit)[0],:]= bestClustAss#reassign new clusters, and SSE
    return mat(centList), clusterAssment

下面對最後的程式碼進行解析:

      bestClustAss[nonzero(bestClustAss[:,0].A == 1)[0],0] = len(centList) #change 1 to 3,4, or whatever
      bestClustAss[nonzero(bestClustAss[:,0].A == 0)[0],0] = bestCentToSplit

這裡是更改其所屬的類別,其中bestClustAss = splitClustAss.copy()是進行k-means後所返回的矩陣,其中第一列為類別,第二列為SSE值,因為當k=2是k-means返回的是類別0,1兩類,因此這裡講類別為1的更改為其質心的長度,而類別為0的返回的是該簇原先的類別。

舉個例子:
例如:目前劃分成了0,1兩個簇,而要求劃分成3個簇,則在演算法進行時,假設對1進行劃分得到的SSE最小,則將1劃分成了2個簇,其返回值為0,1兩個簇,將返回為1的簇改成2,返回為0的簇改成1,因此現在就有0,1,2三個簇了。

  centList[bestCentToSplit] = bestNewCents[0,:].tolist()[0]#replace a centroid with two best centroids 
        centList.append(bestNewCents[1,:].tolist()[0])
        clusterAssment[nonzero(clusterAssment[:,0].A == bestCentToSplit)[0],:]= bestClustAss#reassign new clusters, and SSE

其中bestNewCents是k-means的返回簇中心的值,其有兩個值,分別是第一個簇,和第二個簇的座標(k=2),這裡將第一個座標賦值給 centList[bestCentToSplit] = bestNewCents[0,:].tolist()[0],將另一個座標新增到centList中 centList.append(bestNewCents[1,:].tolist()[0])

執行與結果

>>> from numpy import *
>>> import kMeans
>>> dat = mat(kMeans.loadDataSet('testSet2.txt'))
>>> cent,assment=kMeans.biKmeans(dat,3)
sseSplit, and notSplit:  570.722757425 0.0
the bestCentToSplit is:  0
the len of bestClustAss is:  60
sseSplit, and notSplit:  68.6865481262 38.0629506357
sseSplit, and notSplit:  22.9717718963 532.659806789
the bestCentToSplit is:  0
the len of bestClustAss is:  40

可以看到進行了兩次的劃分,第一次最好的劃分是在0簇,第二次劃分是在1簇。
視覺化如下圖所示:

這裡寫圖片描述

Mini Batch k-Means

在原始的K-means演算法中,每一次的劃分所有的樣本都要參與運算,如果資料量非常大的話,這個時間是非常高的,因此有了一種分批處理的改進演算法。
使用Mini Batch(分批處理)的方法對資料點之間的距離進行計算。
Mini Batch的好處:不必使用所有的資料樣本,而是從不同類別的樣本中抽取一部分樣本來代表各自型別進行計算。n 由於計算樣本量少,所以會相應的減少執行時間n 但另一方面抽樣也必然會帶來準確度的下降。