from numpy import *
import xlrd
import matplotlib.pyplot as plt

# 計算歐氏距離
def distEclud(vector1, vector2):
    '''
    :param vector1: 第j個均值向量
    :param vector2: 第i個樣本
    :return: 距離值
    '''
    return sqrt(sum(power(vector2 - vector1, 2)))

#構建聚簇中心
def randCent(dataSet, k):
    n = shape(dataSet)[1]
    centroids = mat(zeros((k,n)))
    for j in range(n):
        minJ = min(dataSet[:,j])
        maxJ = max(dataSet[:,j])
        rangeJ = float(maxJ - minJ)
        centroids[:,j] = minJ + rangeJ * random.rand(k, 1)
    return centroids

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迴圈會不停的對簇進行劃分，直到得到想要的簇數目為止
    while(len(centList) < k):
        lowestSSE = inf    #python中inf表示正無窮
        #遍歷簇列表centList中的每一個簇
        for i in range(len(centList)):
            #-----嘗試劃分每一簇----
            #對每個簇，將該簇中的所有點看成一個小的資料集ptsInCurrCluster
            ptsInCurrCluster = dataSet[nonzero(clusterAssment[:,0].A==i)[0],:]
            #將ptsInCurrCluster輸入到函式KMeans()中進行處理(k=2)。K-均值演算法會生成兩個質心(簇)，同時給出每個簇的誤差值
            centroidMat, splitClustAss = kMeans(ptsInCurrCluster, 2, distMeas)

            #劃分後的樣本誤差之和
            sseSplit = sum(splitClustAss[:, 1])
            #剩餘資料集的誤差之和
            sseNotSplit = sum(clusterAssment[nonzero(clusterAssment[:, 0].A != i)[0], 1])
            print("sseSplit, and notSplit: ", sseSplit, sseNotSplit)
            #將劃分後的樣本誤差之和+剩餘資料集的誤差之和作為本次劃分的誤差
            if(sseSplit + sseNotSplit) < lowestSSE:
                #決定要劃分某一個簇
                bestCentToSplit = i
                #劃分後的質心
                bestNewCents = centroidMat
                #劃分後的簇分配結果矩陣，包含兩列：第一列記錄簇索引，第二列儲存誤差
                bestClustAss = splitClustAss.copy()
                # 本次劃分的誤差
                lowestSSE = sseSplit + sseNotSplit
            #----更新簇的分配結果----
        #將劃分簇中所有點的簇分配結果進行修改，當使用KMeans()函式並且指定簇數為2時，會得到兩個編號分別為0和1的結果簇
        #需要將這些簇編號修改為劃分簇及新加簇的編號，該過程可以通過兩個陣列過濾器來完成。
        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中
        centList[bestCentToSplit] = bestNewCents[0, :].tolist()[0]
        centList.append(bestNewCents[1, :].tolist()[0])
        # 更新SSE的值(sum of squared errors)
        clusterAssment[nonzero(clusterAssment[:, 0].A == bestCentToSplit)[0],:] = bestClustAss
    return mat(centList),clusterAssment

#k-means 聚類演算法
def kMeans(dataSet, k, distMeans =distEclud, createCent = randCent):
    m = shape(dataSet)[0]
    clusterAssment = mat(zeros((m,2)))    #用於存放該樣本屬於哪類及質心距離
    centroids = createCent(dataSet, k)
    clusterChanged = True
    while clusterChanged:
        clusterChanged = False;
        for i in range(m):
            minDist = inf; minIndex = -1;
            for j in range(k):
                distJI = distMeans(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)
        for cent in range(k):
            ptsInClust = dataSet[nonzero(clusterAssment[:,0].A == cent)[0]]   # 去第一列等於cent的所有列
            centroids[cent,:] = mean(ptsInClust, axis = 0)
    return centroids, clusterAssment

# show your cluster only available with 2-D data
def showCluster(dataSet, k, centroids, clusterAssment):
    numSamples, dim = dataSet.shape
    if dim != 2:
        print ("Sorry! I can not draw because the dimension of your data is not 2!")
        return 1

    mark = ['or', 'ob', 'og', 'ok', '^r', '+r', 'sr', 'dr', '<r', 'pr']
    if k > len(mark):
        print ("Sorry! Your k is too large! please contact Zouxy")
        return 1

        # draw all samples
    for i in range(numSamples):
        markIndex = int(clusterAssment[i, 0])
        plt.plot(dataSet[i, 0], dataSet[i, 1], mark[markIndex])

    mark = ['Dr', 'Db', 'Dg', 'Dk', '^b', '+b', 'sb', 'db', '<b', 'pb']
    # draw the centroids
    for i in range(k):
        plt.plot(centroids[i, 0], centroids[i, 1], mark[i], markersize=12)
    plt.show()

def main():
    ## step 1: load data
    print ("step 1: load data...")
    dataSet = []
    data = xlrd.open_workbook('C:/Users/Microstrong/Desktop/watermelon4.0.xlsx')
    table = data.sheets()[0]
    for line in range(0,table.nrows):
        lineArr = table.row_values(line)
        dataSet.append([float(lineArr[0]), float(lineArr[1])])

    ## step 2: clustering...
    print ("step 2: clustering...")
    dataSet = mat(dataSet)
    k = 3
    centroids, clusterAssment = biKmeans(dataSet, k)

    ## step 3: show the result
    print ("step 3: show the result...")
    showCluster(dataSet, k, centroids, clusterAssment)

if __name__ == '__main__':
 main()