五、代码实现（python）
以下代码来自Peter Harrington《Machine Learing in Action》。
本例代码实现二分k均值算法。
代码如下（保存为kMeans.py)：
# -- coding: utf-8 --
from numpy import *
def loadDataSet(fileName):
# 获取数据集
dataMat = []
fr = open(fileName)
for line in fr.readlines():
curLine = line.strip().split(' ')
fltLine = map(float,curLine)
dataMat.append(fltLine)
return dataMat
def distEclud(vecA, vecB):
# 根据式()计算vecA, vecB两点间的欧氏距离
return sqrt(sum(power(vecA - vecB, 2)))
def randCent(dataSet, k):
# 随机生成k个质心
n = shape(dataSet)[1] # 获取数据集特征数量，即列数
centroids = mat(zeros((k,n))) # 初始化一个k行n列的矩阵，元素为0，用于存储质心
for j in range(n):
minJ = min(dataSet[:,j]) # 获取数据集第j列的最小值
rangeJ = float(max(dataSet[:,j]) - minJ) # 计算数据集第j列中，最大值减最小值的差
# 随机生成k行1列的数组，元素在0到1之间，乘以rangeJ再加上minJ，则可得随机生成的第j列中最小值与最大值之间的一个数
centroids[:,j] = mat(minJ + rangeJ * random.rand(k,1))
return centroids
def kMeans(dataSet, k, distMeas=distEclud, createCent=randCent):
# kMeans函数接受4个输入参数，数据集及簇的数目为必选参数，计算距离默认为欧氏距离，创建初始质心默认为随机生成
m = shape(dataSet)[0] # 获取数据集数量，即行数
clusterAssment = mat(zeros((m,2))) # 初始化一个m行2列的矩阵，元素为0，第一列存储当前最近质心，第二列存储数据点与质心的距离平方
centroids = createCent(dataSet, k) # 创建k个初始质心
clusterChanged = True
while clusterChanged:
clusterChanged = False
for i in range(m):
# 循环m个数据，寻找距离第i个数据最近的质心
minIndex = -1 # 初始化最近质心
minDist = inf # 初始化第i个数据与最近质心的最小距离为无穷大
for j in range(k):
# 循环k个质心，寻找离第i个数据最近的质心
distJI = distMeas(centroids[j,:],dataSet[i,:]) #计算第i行数据与第j个质点的欧氏距离
if distJI < minDist:
minIndex = j # 更新最近质心为第j个
minDist = distJI # 更新第i个数据与最近的质心的最小距离
if clusterAssment[i,0] != minIndex: clusterChanged = True # 若其中clusterAssment存储的质心与此次结果不一样，则需迭代，直至没有质心的改变
clusterAssment[i,:] = minIndex,minDist**2 # 更新clusterAssment数据
for cent in range(k):
ptsInClust = dataSet[nonzero(clusterAssment[:,0].A==cent)[0]] # 获取属于第cent个质心的所有数据
centroids[cent,:] = mean(ptsInClust, axis=0) # 计算属于第cent个质心的所有数据各列的平均值，更新第cent个质心
return centroids, clusterAssment # centroids为当前k个质心，clusterAssment为各个数据所属质心及距离该质心的距离平方
def biKmeans(dataSet, k, distMeas=distEclud):
# biKmeans函数接受3个输入参数，数据集及簇的数目为必选参数，计算距离默认为欧氏距离
m = shape(dataSet)[0] # 获取数据集数量，即行数
clusterAssment = mat(zeros((m,2))) # 初始化一个m行2列的矩阵，元素为0，第一列存储当前最近质心，第二列存储数据点与质心的距离平方
centroid0 = mean(dataSet, axis=0).tolist()[0] # 将所有点作为一个簇，计算数据集各列的平均值，作为初始簇的质心
centList = [centroid0] # centList存储各个质心
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)):
# 循环簇的个数，寻找使SSE下降最快的簇的划分
ptsInCurrCluster = dataSet[nonzero(clusterAssment[:,0].A==i)[0],:] # 获取属于第i个质心的所有数据
centroidMat, splitClustAss = kMeans(ptsInCurrCluster, 2, distMeas) # 将第i个簇二分为2个簇
sseSplit = sum(splitClustAss[:,1]) # 计算第i个簇二分为2个簇后的SSE值
sseNotSplit = sum(clusterAssment[nonzero(clusterAssment[:,0].A!=i)[0],1]) # 计算剩余数据的SSE值
if (sseSplit + sseNotSplit) < lowestSSE:
# 若二分后总体SSE值下降，则更新簇的信息
bestCentToSplit = i
bestNewCents = centroidMat # 第i个簇二分后的质心
bestClustAss = splitClustAss.copy()# 第i个簇二分后的结果
lowestSSE = sseSplit + sseNotSplit # 更新当前SSE值
bestClustAss[nonzero(bestClustAss[:,0].A == 1)[0],0] = len(centList) # 更新簇的分配结果，将二分后第二个簇分配到新簇
bestClustAss[nonzero(bestClustAss[:,0].A == 0)[0],0] = bestCentToSplit # 更新簇的分配结果，将二分后第一个簇分配到被划分簇
centList[bestCentToSplit] = bestNewCents[0,:].tolist()[0] # 更新簇的分配结果，更新未划分簇质心为二分后第一个簇的质心
centList.append(bestNewCents[1,:].tolist()[0]) # 更新簇的分配结果，添加新质心为二分后第二个簇的质心
clusterAssment[nonzero(clusterAssment[:,0].A == bestCentToSplit)[0],:]= bestClustAss # 将被划分簇数据更新为划分后簇
return mat(centList), clusterAssment