算法实现
一、K均值算法
K均值算法是典型的原型聚类算法,此类算法假设聚类结构能够通过一组原型刻画。通常情况下,算法先对原型进行初始化,然后对原型进行迭代更新求解。
K均值算法主要过程概括为:
1>初始化均值向量{
μ
1
,
μ
2
,
.
.
,
μ
k
μ_1,μ_2,..,μ_k
μ1,μ2,..,μk},
k
k
k值由用户给出,为聚类簇数。
2>对样本进行归类:对给定数据集
X
X
X中每个样本
x
x
x利用距离度量
d
=
∣
∣
x
−
μ
∣
∣
2
d=||x-μ||_2
d=∣∣x−μ∣∣2计算每个样本到均值向量的距离,把样本
x
x
x归为其距离最近的均值向量所对应的簇。比如:对第
i
i
i个样本
x
i
x_i
xi计算其到每个均值向量
μ
1
,
μ
2
,
.
.
.
,
μ
k
μ_1,μ_2,...,μ_k
μ1,μ2,...,μk的距离,如果他到第
j
j
j簇
C
j
C_j
Cj,就把
x
i
x_i
xi划分到
C
j
C_j
Cj中。
3>更新均值向量:现在我们每个均值向量
μ
j
μ_j
μj对应的簇
C
j
C_j
Cj下都有我们第二步通过归类所得到的样本。我们对每个簇下的样本进行求平均得到:
μ
j
′
=
1
∣
C
j
∣
∑
x
∈
C
i
x
μ_j^{'}=\frac{1}{|C_j|}\sum_{x∈C_i}x
μj′=∣Cj∣1∑x∈Cix。这样我们就可以得到把
μ
j
′
μ_j^{'}
μj′当做新的均值向量
μ
j
μ_j
μj。
4>跳至第二步,进行迭代,直至满足迭代终止条件为止。一般迭代终止条件为第二步中所有样本的类所属簇不在改变,或者由用户设定迭代次数。
二、K均值算法的实现
首先,根据我们上面的算法可以知道主要是对样本进行归类、更新均值向量,这样我们就需要编写两个函数。然后还需要将以上两个函数合并,给出迭代终止条件。
1、导入需要使用的模块
import numpy as np
from scipy.io import loadmat
import pandas as pd
import matplotlib.pyplot as plt
2、初始化均值向量
这里我们还需要编写一个能够初始化均值向量的函数。聚类的最终结果一定程度上取决于我们初始均值向量的选取,我们这里采用随机方法,产生初始均值向量。
def initial_centriod(X,k):
m,d = X.shape
initial_centriod=np.zeros((k,d))
index=np.random.randint(0,m,k)
for i in range(k):
initial_centriod[i,:]=X[index[i],:]
return initial_centriod
代码第三行创建一个
k
∗
d
k*d
k∗d零矩阵用于存放均值向量,其中
k
k
k是我们设定的聚类簇数,
d
d
d为样本的维度或者说是我们研究的特征空间的维度,矩阵每一行就是用来存放初始均值向量,因为我们要产生
k
k
k个均值向量,所以矩阵为k行。
第四行我们产生一个大小为
k
k
k,范围为
1
—
m
1—m
1—m的数组,
m
m
m为样本的大小。最后对零向量每一行进行遍历,使其每一行为
m
m
m个样本中随机产生的一个。例如,我们把
k
k
k设置为3,一共有20个样本,那么index
就有可能是数组
[
1
,
3
,
4
]
[1,3,4]
[1,3,4],到底下对零矩阵每一行进行遍历时,初始均值向量的第一行就会是第1个样本,第二行为第3个样本,第三行为第4个样本。也就是说我们把第1,3,4个样本作为初始均值向量。
3、对每个样本进行归类
def find_closest_centriod(X,centriod):
m=X.shape[0]
k = centriod.shape[0]
index_centriod=np.zeros(m)
for i in range(m):
minidist = 1000000000000
for j in range(k):
dist = (centriod[j]-X[i]).T @ (centriod[j]-X[i])
if dist < minidist:
minidist = dist
index_centriod[i]=j
return index_centriod
这里代码第4行的index_centriod
它的维度为训练集大小
m
m
m,用来存放着
m
m
m个样本的所属簇的编号。然后我们队每一个样本进行遍历,在对每一个均值向量进行遍历,找出每一个样本最近的均值向量,第8行即是在计算当前遍历样本与当前遍历均值向量的距离平方大小。找出
x
i
x_i
xi的最近均值向量
μ
j
μ_j
μj后,即把index_centriod
的第i个位置设置为均值向量
μ
j
μ_j
μj对应的簇的编号j。
4、计算每个簇下样本的平均值,得到新的均值向量
def compute_centroid(X,index,k):
m,d = X.shape
centriod_mat = np.zeros((k,d))
for j in range(k):
l=[]
for i in range(m):
if index[i]==j:
l.append(X[i])
centriod_mat[j,:]=(sum(l))/len(l)
return centriod_mat
这里我们遍历每个簇,对于每个簇创建一个列表,用于每个簇对应的样本,遍历每个上面分类得到的关于每个样本所属类别的数组,如果数组第 i i i个位置对应的数为第 j j j个簇对应的编号,那么就相当于我们的第i个样本 x i x_i xi属于第 j j j簇,就将其存放在列表中,最后对 j j j簇对应的列表中的所有样本求平均,即为新的均值向量 μ j μ_j μj,将其存放在已经创建的矩阵的第 j j j行,这个矩阵的第 i i i行存放着第 i i i簇的均值向量。
5、将两个函数整合在一起,即主函数
def run_k_means(X,centriod,max_iter):
old_centriod=centriod
for i in range(max_iter):
old_index=find_closest_centriod(X,old_centriod)
k = old_centriod.shape[0]
new_centriod=compute_centroid(X,old_index,k)
new_index=find_closest_centriod(X,new_centriod)
if list(new_index)==list(old_index):
break
else:
old_centriod=new_centriod
return new_index,new_centriod
参数centriod
为初始均值向量,max_iter
为用户设置的迭代最大次数。如果我们发现迭代前后两次对样本分类结果一样,就立即停止迭代,否则执行到最大迭代次数max_iter
。
对数据进行聚类
一、处理数据并查看
先查看数据的维度:
data = loadmat('F:\\MachineLearning\data\ex7data2.mat')
X=data['X']
print(X.shape)
得到
(300,2)
即数据为二维,可以对其进行可视化:
#查看数据分布
data2 = pd.DataFrame(X,columns=['x1','x2'])
fig1,ax1=plt.subplots(figsize=(9,6))
ax1.scatter(data2['x1'],data2['x2'],s=20,color='blue')
plt.xlabel('x1')
plt.ylabel('x2')
plt.show()
得到数据分布图:
很明显,我们可以把数据分为三类。
二、对数据进行聚类
1、首先,我们产生初始均值向量,用上面的initial_centriod
函数:
centriod = initial_centriod(X,3)
2、接着,运行主函数进行聚类,并处理数据:
index,centriod=run_k_means(X,centriod,10)
data2['index']=index
data_index1=data2[data2['index']==0]
data_index2=data2[data2['index']==1]
data_index3=data2[data2['index']==2]
,这里我们设置最大迭代次数为10,然后得到了聚类原型,也就是最后一次迭代的均值向量centriod
,以及每个样本聚类的结果index
。然后我们给数据添加一列,即为index
,给出每个样本对应的聚类编号。然后分别把属于簇1、2、3的样本分别存放data_index1
、data_index2
、data_index3
。
3、查看聚类结果
fig2,ax2 = plt.subplots(figsize=(9,6))
ax2.scatter(data_index1['x1'],data_index1['x2'],s=20,color='blue',label='cluster1')
ax2.scatter(data_index2['x1'],data_index2['x2'],s=20,color='red',label='cluster2')
ax2.scatter(data_index3['x1'],data_index3['x2'],s=20,color='green',label='cluster3')
plt.legend()
plt.show()
这里我们将簇1用蓝色表示,簇2用红色表示,簇3用绿色表示:
发现聚类结果良好。
当然也会发现例如一下这样的情况:
这样的聚类结果显然是不能令人满意的,这主要是因为我们的初始均值向量的选择会影响聚类结果,然而我们的均值向量是随机产生的,会有一定可能性(可能性较小)选择了“不好”的初始均值向量。如果遇到这种情况我们可以多次运行程序,当然也可以根据数据分布自行选择。
利用聚类算法进行图像压缩
一、查看原图,并产看数据
首先,查看我们需要处理的原图:
from IPython.display import Image,display
display(Image(filename='F:/MachineLearning/data/bird_small.png'))
得到原图:
接着,处理数据:
先查看数据维度
image_data=loadmat('F:/MachineLearning/data/bird_small.mat')
A=image_data['A']
print(A.shape)
(128, 128, 3)
即数据为三维数组。前两个维度即128和128表示我们的平面图由128128个像素点,第三个维度3表示我们每个像素点的红、蓝、绿的强度。所以我们有128128个样本,样本的维度为3。而我们利用聚类算法进行压缩的原理主要是将我们128*128个样本处理成k个簇,再把每个样本代表的簇的原型代替原先样本,利用这些点组成新的图。
二、利用聚类算法进行压缩图片
A=A/255
X=np.reshape(A,(A.shape[0]*A.shape[1],A.shape[2]))
initial_image_centriod=initial_centriod(X,16)
index,centriod=run_k_means(X,initial_image_centriod,10)
X_recovered = centriod[index.astype(int),:]
X_recovered=np.reshape(X_recovered,(A.shape[0],A.shape[1],A.shape[2]))
plt.imshow(X_recovered)
plt.show()
第一行:归一化数据,让imshow函数可以处理。第二行:我们把样本展开,形成128*128行,3列的样本矩阵,可以让程序处理。第三行:初始化聚类均值向量,这里我们聚类簇数设置为16,簇数决定了压缩图片的还原原图的质量,一般这个值设置的越大,质量越好,同时图片内存也越大。第四行:运行程序。第五行:将原始数据点用聚类原型代替。第六行:再次改变维度使变为图片维度。第七行:形成压缩图片。最后,显示图片。得到图:
这是我们聚类簇数为16的图。如果我们将第三行聚类簇数设置为32、8和4,会看到如下三图:
会发现参数为32好于参数16,参数为16明显好于参数为8,参数为8明显好于参数为4。这不难理解,聚类簇数相当压缩图片所使用的色彩数,色彩越少,图片质量肯定越差。
三、用sklearn中K-Means库处理此问题
from skimage import io
pic=io.imread('F:/MachineLearning/data/bird_small.png')
pic=pic/255
print(pic.shape)
io.imshow(pic)
plt.show
这里我们直接从原图得到数据,第一行:导入需要的模块。第二行:io.imread
将我们的原图处理为数组形式。第三行:归一化。第四行:查看数据维度,为(128, 128, 3)
。第五、六行再利用得到数据重构原图。我们会得到:
这实际上就是原图。
from sklearn.cluster import KMeans
#系数解释:n_cluster:聚类数。n_init:使用初始点次数。n_jobs:是否使用多线程,-1使用。
Cluster = KMeans(n_clusters=4,n_init=100,n_jobs=-1)
Cluster.fit(X)
centriod_pic=Cluster.cluster_centers_
print(centriod_pic)
C=Cluster.predict(data_pic)
pic_recovered=centriod_pic[C].reshape((128,128,3))
上面是调用K-Means库的情况。第一行:导入模块。第三行:设置参数。第四行:拟合数据。第五行:读取聚类原型矩阵。第七行:读取聚类结果C,即每个样本属于簇组成的数组。第八行:再把数据处理成可以图片数组。
这里我们可以得到centriod_pic
,聚类原型矩阵:
[[0.48994379 0.40168338 0.32354671]
[0.1287319 0.13038871 0.12085462]
[0.91127734 0.8588342 0.74082312]
[0.79460449 0.63997646 0.41903748]]
最后查看原图和压缩图片的差别: