K-means聚类——简单代码实现、与EM算法的关联_用kmeans聚类的初始化实现em算法

算法实现步骤

参考链接：https://blog.csdn.net/weixin_41090915/article/details/79389636

K-means思想：

（1）随机选择k个类簇的中心
（2）计算每一个样本点到所有类簇中心的距离，选择最小距离作为该样本的类簇
（3）重新计算所有类簇的中心坐标，直到达到某种停止条件（迭代次数/簇中心收敛/最小平方误差）

---

EM算法思想：

（1）假设我们想估计知道A和B两个参数，在开始状态下二者都是未知的;
（2）但如果知道了A的信息就可以得到B的信息，反过来知道了B也就得到了A;
（3）首先赋予A某种初值，以此得到B的估计值，然后从B的当前值出发，重新估计A的取值，这个过程一直持续到收敛为止。

---

两者关联：

从上面两者的步骤可以看到，k-means其实就是EM的体现，两者一开始对于参数都是未知的，然后随机选择一个初始变量，得到对象结果后，再反过来推算整体，反复迭代以上两步，直到得到最优结果。

总体思想还是一个迭代优化过程，有目标函数，也有参数变量，只是k-means多了个隐含变量，确定其他参数估计隐含变量，再确定隐含变量估计其他参数，直至目标函数最优。

---

代码实现：

计算步骤：
1.迭代执行以下操作：遍历每个样本点，计算样本点i到各个质心的距离，找到最小距离，将该质心所在簇编号分配给该样本点；2.遍历完所有样本点后，重新计算每个簇的质心。直到所有样本点的簇分配都不再发生变化时迭代停止；
3.最后返回质心和样本点的簇分配矩阵。

import numpy as np
import math
import matplotlib.pyplot as plt

#1.加载数据
def loadDataSet(fileName, splitChar='\t'):
    dataSet = []
    with open(fileName) as fr:
        for line in fr.readlines():
            curline = line.strip().split(splitChar)
            fltline = list(map(float, curline))
            dataSet.append(fltline)
    return np.array(dataSet)


#2.计算欧氏距离：两个向量每个元素差的平方，求和后再开方
def dist_eclud(vecA, vecB):
    vec_square = []
    for element in vecA - vecB:   #元素想减
        element = element ** 2    #求平方
        vec_square.append(element)   #放入list
    return sum(vec_square) ** 0.5    #list求和后，再开方


# 3.构建k个随机质心
def rand_cent(data_set, k):
    n = data_set.shape[1]    
    centroids = np.zeros((k, n))   #构造质心的shape
    #限定质心在数据集范围之内
    for j in range(n):
        min_j = float(min(data_set[:,j]))    #数据集中的最小值
        range_j = float(max(data_set[:,j])) - min_j     #
        centroids[:,j] = (min_j + range_j * np.random.rand(k, 1))[:,0]
    return centroids


#4.算法的主函数
def Kmeans(data_set, k):    
    m = data_set.shape[0]
    cluster_assment = np.zeros((m, 2))     
    centroids = rand_cent(data_set, k)     #获取k个随机质心
    cluster_changed = True

    while cluster_changed:       
        cluster_changed = False
         #遍历m个原始数据
        for i in range(m):        
            min_dist = np.inf; min_index = -1
            #遍历k个质心  
            for j in range(k):
                #计算m个原始数据与每一个质心的距离
                dist_ji = dist_eclud(centroids[j,:], data_set[i,:])
                #找到最小距离，将该质心所在簇编号分配给该样本点
                if dist_ji < min_dist:              
                    min_dist = dist_ji
                    min_index = j   

            #重新计算每个簇的质心。直到所有样本点的簇分配都不再发生变化时迭代停止。
            if cluster_assment[i,0] != min_index:    
                cluster_changed = True      
            cluster_assment[i,:] = min_index, min_dist**2
            
        #返回质心和样本点的簇分配矩阵
        for cent in range(k):   
            pts_inclust = data_set[np.nonzero(list(map(lambda x:x==cent, cluster_assment[:,0])))]
            centroids[cent,:] = np.mean(pts_inclust, axis=0)

    return centroids, cluster_assment            


#加载数据
data_set = loadDataSet('img/1.txt', splitChar=',')
#计算出最终的质心
my_centroids, my_cluster_assment = Kmeans(data_set, 4)
print("最终4个质心为：")
print(my_centroids)

#原始数据的坐标点
point_x = data_set[:,0]
point_y = data_set[:,1]
#最终质心的坐标点 
cent_x = my_centroids[:,0]
cent_y = my_centroids[:,1]
#画图
fig, ax = plt.subplots(figsize=(10,5))
#画原始数据点，形状为圆圈
ax.scatter(point_x, point_y, s=30, c="r", marker="o", label="sample point")
#画最终数据点，形状为倒三角
ax.scatter(cent_x, cent_y, s=100, c="black", marker="v", label="centroids")
#在右上角显示小图标
ax.legend()
ax.set_xlabel("factor1")   #横轴名称
ax.set_ylabel("factor2")   #纵轴名称
plt.show()

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95

---

结果如下：
可以得到新的4个质心

---

把原始数据和最终质心在图上画出来：

可以看到，质心的位置基本分布在每一团数据点的中心，实现了聚类。

---

---

代码实现：K均值聚类（sklearn调用）

from sklearn.cluster import KMeans
import numpy as np
import pandas as pd
from numpy import *
import matplotlib.pyplot as plt
import time
from sklearn.externals import joblib

def loadDataSet(fileName, splitChar='\t'):
    dataSet = []
    with open(fileName) as fr:
        for line in fr.readlines():
            curline = line.strip().split(splitChar)
            fltline = list(map(float, curline))
            dataSet.append(fltline)
    return dataSet

def createDataSet():
    dataSet = [[0.0, 2.0],
               [0.0, 0.0],
               [1.5, 0.0],
               [5.0, 0.0],
               [5.0, 2.0]]
    return dataSet


# 加载数据
data_set = loadDataSet('img/1.txt', splitChar=',')   #1.加载一个txt数据集
# data_set = createDataSet()   #2.自己创建一个数据集
data_set = mat(data_set)


# n_clusters: 指定K的值
# max_iter: 对于单次初始值计算的最大迭代次数
# n_init: 重新选择初始值的次数
# init: kmeans++表示该初始化策略选择的初始均值向量之间都距离比较远，它的效果较好；
# random表示从数据中随机选择K个样本最为初始均值向量。
# n_jobs: 进程个数，为-1的时候是指默认跑满CPU
# tol：算法收敛的阈值,默认0.001。
# verbose：0表示不输出日志信息；1表示每隔一段时间打印一次日志信息。如果大于1，打印次数频繁。
# cluster_centers_：聚类后的质心

clf=KMeans(n_clusters=3, max_iter=300, n_init=10, init='k-means++',n_jobs=-1,tol=0.001)
clf=clf.fit(data_set)
print(clf.cluster_centers_)


# 保存/加载模型，.m or .pkl文件
joblib.dump(clf, "img/train_model.m")
# model = joblib.load("img/train_model.m")

data=np.array([1,5]).reshape(1, -1)
pred=clf.predict(data)
print(pred)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54