获取案例链接、直播课件、数据集请访问数据酷客社区版块。
K-Means 是一种最经典和常用的聚类方法。它通过多轮迭代的方式不断更新不同类样本的中心，计算样本到每个中心的距离，然后更新样本所属的类。最终能够把样本划分到 K 个类中。本案例中，我们首先使用 Python 实现 K-Means 算法，基于一份随机数据集，使用动画演示聚类过程和优化目标的变化。然后将 K-Means 应用于图像分割问题。最后我们还将使用一份中文新闻数据集，用 K-Means 算法进行自动新闻主题聚类，并使用柱状图和词云图对聚类结果进行可视化分析。
1 实现 K-Means 算法
K-Means 算法的基本运行流程为：
- 1 随机选择 $k$ 个点作为初始中心
- 2 Repeat：
- 2.1 将每个样本指派到最近的中心，形成 $k$ 个类
- 2.2 重新计算每个类的中心为该类样本均值
- 3 直到中心不发生变化
在整个算法中，2.1 步骤运算量最大，因为该步骤需要计算每一个样本到 $k$ 个中心的距离。假设我们的数据集表示为 Pandas 的 DataFrame 格式 X。要实现每个样本到中心的距离计算，可以有多种实现方式。最简单的是通过 for 循环遍历 X 的每一行，计算每一行到中心的距离。为了提高运算效率，我们可以借助 DataFrame 的 apply 函数加快运算。Numpy 是 Python 中矩矩阵运算的高效工具，我们也可以使用矩阵运算的方式来实现样本到中心的距离计算。下面，我们分别来看下这三种实现方法。
1.1 使用 iterrows 遍历的方式实现
假设我们使用欧式距离计算样本到中心的距离。对于样本
使用最简单的方式来实现，先用一个函数 point_dist 计算一个样本到中心的距离。这里我们使用 Numpy 的线性代数模块 linalg 中的 norm 方法。
import numpy as np
def point_dist(x,c): #定义距离计算函数
return np.linalg.norm(x-c)然后使用 iterrows 方法遍历样本计算样本到中心的距离，定义 k_means_iterrows 方法实现 K-Means 算法。
def k_means1(X,k):
centers = X.sample(k).values #从数据集随机选择 K 个样本作为初始化的类中心，k 行 d 列
X_labels = np.zeros(len(X)) #样本的类别
error = 10e10
while(error > 1e-6):
for i,x in X.iterrows():#指派样本类标签
X_labels[i] = np.argmin([point_dist(x,centers[i,:]) for i in range(k)])
centers_pre = centers
centers = X.groupby(X_labels).mean().values #更新样本均值，即类中心
error = np.linalg.norm(centers_pre - centers)#计算error
return X_labels, centers用一个简单的随机数据集来测试时间性能。Sklearn 中的 datasets 模块的 make_blobs 函数能够自动生成一些供测试聚类算法的随机数据集。它能够根据输入的参数生成数据集和对应的类标签。常用的参数如下表：
from sklearn import datasets
import pandas as pd
X, y = datasets.make_blobs(n_samples=5000, n_features=8, cluster_std = 0.5,centers=3,random_state=99)
X_df = pd.DataFrame(X)在该数据集上用我们实现的 k_means1 方法运行 K-Means 聚类。使用 iPython 提供的魔法命令 %time 记录运行时间。
%time labels,centers = k_means1(X_df,3) # for 循环CPU times: user 5.16 s, sys: 39.3 ms, total: 5.19 s
Wall time: 5.14 s
1.2 使用 apply 遍历的方式实现
提高运算效率，可以使用 DataFrame 的 apply 函数，它可以对数据框中的每一行执行一个复杂的函数。在我们的例子中，是计算每一行与每一个中心的距离。
def k_means2(X,k):
#初始化 K 个中心，从原始数据中选择样本
centers = X.sample(k).values
X_labels = np.zeros(len(X)) #样本的类别
error = 10e10
while(error > 1e-6):
#********#
X_labels = X.apply(lambda r : np.argmin([point_dist(r,centers[i,:]) for i in range(k)]),axis=1)
centers_pre = centers
centers = X.groupby(X_labels).mean().values #更新样本均值，即类中心
error = np.linalg.norm(centers_pre - centers)#计算error
return X_labels, centers
%time labels,centers = k_means2(X_df,3) # apply 运算CPU times: user 3.24 s, sys: 19.9 ms, total: 3.26 s Wall time: 3.23 s
1.3 矩阵运算的方式实现
数据集表示成
已经聚类中心，计算样本到中心距离，并将样本划分到距离最小的类的流程如下图所示。
使用 Numpy 实现上述计算流程的代码为：
for i in range(k):
D[:,i] = np.sqrt(np.sum(np.square(X - C[i,:]),axis=1))
labels = np.argmin(D,axis=1)得到样本的类标签后，聚类中心的更新流程为：1）根据类标签对样本进行分组；2）将聚类中心更新为每一组样本的均值。Python 实现的代码为：
C = X.groupby(labels).mean().values现在我们更新 K-Means 算法的实现，函数名为 k_means。
import pandas as pd
import numpy as np
def k_means(X,k):
C = X.sample(k).values #从数据集随机选择 K 个样本作为初始化的类中心，k 行 d 列
X_labels = np.zeros(len(X)) #记录样本的类别
error = 10e10 #停止迭代的阈值
while(error > 1e-6):
D = np.zeros((len(X),k)) #样本到每一个中心的距离，n 行 k 列
for i in range(k):
D[:,i] = np.sqrt(np.sum(np.square(X - C[i,:]),axis=1))
labels = np.argmin(D,axis=1)
C_pre = C
temp_C = X.groupby(labels).mean() #更新样本均值，即类中心
C = np.zeros((k,X.shape[1]))
for i in temp_C.index:
C[i,:] = temp_C.loc[i,:].values
if C.shape == C_pre.shape:
error = np.linalg.norm(C_pre - C)#计算error
else:
print(C.shape, C_pre.shape)
return labels, C%time labels,centers = k_means(X_df,3)# 矩阵运算CPU times: user 150 ms, sys: 0 ns, total: 150 ms Wall time: 149 ms
1.4 聚类结果可视化
下面我们使用一份随机生成的二维数据集，使用我们上一小节实现的 k_means 完成聚类，然后使用不同颜色标注不同类的样本以及类中心。
color_dict ={0:"#E4007F",1:"#007979",2:"blue",3:"orange"}#洋红，深绿，蓝色，橘色from sklearn import datasets
import matplotlib.pyplot as plt
import pandas as pd
%matplotlib inline
X, y = datasets.make_blobs(n_samples=1000, n_features=2, cluster_std = 1.5,centers=4,random_state=999)
X_df = pd.DataFrame(X,columns=["x1","x2"])
labels,centers= k_means(X_df,4)
fig, ax = plt.subplots(figsize=(8, 8)) #设置图片大小
for i in range(len(centers)):
ax.scatter(X_df[labels == i]["x1"],X_df[labels == i]["x2"],color=color_dict[i],s=50,alpha=0.4)
ax.scatter(centers[int(i),0],centers[int(i),1],color="r",s=100,marker="+")
plt.xlabel("$x_1$")
plt.ylabel("$x_2$")Text(0, 0.5, '
1.5 使用动画展示 K-Means 聚类过程
要动态展示 K-Means 聚类过程，我们需要在每一步迭代中记录每一个类的中心，以及每一个类的样本集合。创建 k_means_steps，在完成聚类的同时，将每一步迭代的每类样本和中心返回。
def k_means_steps(X,k):
#初始化 K 个中心，从原始数据中选择样本
#********#
samples_list = [] #记录每一个中间迭代中每一类样本
centers_list = [] #记录每一个中间迭代中每一类样本中心
#********#
C = X.sample(k).values
labels = np.zeros(len(X)) #样本的类别
error = 10e10
while(error > 1e-6):
D = np.zeros((len(X),k)) #样本到每一个中心的距离
for i in range(k):
D[:,i] = np.sqrt(np.sum(np.square(X - C[i,:]),axis=1))
labels = np.argmin(D,axis=1)
C_pre = C
C = X.groupby(labels).mean().values #更新样本均值，即类中心
#********# 记录当前迭代地每一类的样本集合和中心
samples,centers2 = [],[]
for i in range(k):
samples.append(X[labels == i])
centers2.append(C[i,:])
samples_list.append(samples)
centers_list.append(centers2)
#********#
if C.shape == C_pre.shape:
error = np.linalg.norm(C_pre - C)#计算error
else:
print(C.shape, C_pre.shape)
return labels, C,samples_list,centers_list #********# 返回最终的聚类结果，聚类中心，每一步的聚类结果和聚类中心我们可以借助 matplotlib.animation 动画模块来实现下面的 init_draw 函数是动画最开始时绘制的内容，包含数据。update_draw 则是每次更新的内容。
labels,centers,samples_list,centers_list= k_means_steps(X_df,4)
fig, ax = plt.subplots(figsize=(8, 8))
samples_obj = []
centers_obj = []
def init_draw(): # 展现样本数据
ax.set_title("K-Means聚类过程:")
for i in range(len(centers)):
samples_obj.append(ax.scatter(samples_list[0][i]["x1"],samples_list[0][i]["x2"],color=color_dict[i],s=50,alpha=0.6))
centers_obj.append(ax.scatter(centers_list[0][i][0],centers_list[0][i][1],color="r",s=100,marker="+"))
plt.xlabel("$x_1$")
plt.ylabel("$x_2$")
def update_draw(t): # 实现动画中每一帧的绘制函数，i为第几帧
ax.set_title("K-Means聚类过程:" + str(t))
samples,centers = samples_list[t],centers_list[t]
for i in range(len(centers)):
samples_obj[i].set_offsets(samples[i])
centers_obj[i].set_offsets(centers[i])
plt.close()
#演示决策面动态变化
import matplotlib.animation as animation
from IPython.display import HTML
animator = animation.FuncAnimation(fig, update_draw, frames= range(1,len(centers_list)), init_func=init_draw,interval=2000)
HTML(animator.to_jshtml())
1.6 失真度量 J 的变化
我们来查看随着迭代的进行，K-Means聚类模型的优化目标，即失真度量
import pandas as pd
import numpy as np
def k_means_inertia(X,k):
#初始化 K 个中心，从原始数据中选择样本
C = X.sample(k).values
labels = np.zeros(len(X)) #样本的类别
inertia_list = [] #*****记录优化目标****#
error = 10e10
while(error > 1e-6):
D = np.zeros((len(X),k)) #样本到每一个中心的距离
for i in range(k):
D[:,i] = np.sqrt(np.sum(np.square(X - C[i,:]),axis=1))
labels = np.argmin(D,axis=1)
inertia_list.append(np.square(np.min(D,axis=1)).sum()) #****记录当前步骤的失真度量****#
C_pre = C
temp_C = X.groupby(labels).mean() #更新样本均值，即类中心
C = np.zeros((k,X.shape[1]))
for i in range(len(temp_C)):
C[i,:] = temp_C.loc[i,:].values
if C.shape == C_pre.shape:
error = np.linalg.norm(C_pre - C)#计算error
return labels,C,inertia_list在随机数据上进行聚类，并将失真度量的变化以折线图的形式绘制出来。
X, y = datasets.make_blobs(n_samples=1000, n_features=2, cluster_std=1,centers=3,random_state=99)
X_df = pd.DataFrame(X,columns=["x1","x2"])
labels,centers,inertia_list = k_means_inertia(X_df,3)
fig, ax = plt.subplots(figsize=(9, 6)) #设置图片大小
t = range(len(inertia_list))
plt.plot(t,inertia_list,c="#E4007F",marker="o",linestyle='dashed')
plt.xlabel("t")
plt.ylabel("inertia")
plt.xticks(t)
plt.title("K-Means算法优化目标的变化")Text(0.5, 1.0, 'K-Means算法优化目标的变化')
2 使用 K-Means 算法进行图像分割
我们先加载一张测试图片，将图片打印出来。使用 PIL.Image.open 方法来打开图片，然后使用 matplotlib 中的 imshow 方法将图片可视化。
from PIL import Image
fig, ax = plt.subplots(figsize=(6, 5)) #设置图片大小
path = './input/timg.jpg'
img = Image.open(path)
plt.imshow(img)
plt.box(False) #去掉边框
plt.axis("off")#不显示坐标轴(-0.5, 599.5, 514.5, -0.5)
将一张图片转换成表格形式。每一行为一个像素，三列分别为像素的 R，B，G取值。获取图片的每一个像素
Image 类的
getpixel 方法。
import pandas as pd
def image_dataframe(image): #将图片转换成DataFrame,每个像素对应每一行，每一行包括三列
rbg_values = []
for i in range(image.size[0]):
for j in range(image.size[1]):
x,y,z= image.getpixel((i,j)) # 获取图片的每一个像素 (i,j)(i,j) 的 RBG 值
rbg_values.append([x,y,z])
return pd.DataFrame(rbg_values,columns=["R","B","G"]),img.size[0],img.size[1]
img_df,m,n = image_dataframe(img)打印转换的数据框如下：
img_df.head()
print(m,n,m*n,len(img_df))600 515 309000 309000
使用我们在上一节实现的 K-Means 算法对像素进行聚类。
labels, _ = k_means(img_df,2)将生成的灰度图可视化，对图像可视化使用 plt.imshow 方法。
fig, ax = plt.subplots(figsize=(6, 5)) #设置图片大小
labels = labels.reshape((m,n))
pic_new = Image.new("L",(m,n))
#根据类别向图片中添加灰度值
for i in range(m):
for j in range(n):
pic_new.putpixel((i,j),int(256/(labels[i][j] + 1)))
plt.imshow(pic_new)
plt.box(False) #去掉边框
plt.axis("off")#不显示坐标轴(-0.5, 599.5, 514.5, -0.5)
实现一个函数 img_from_labels ，将像素聚类类别标签，转换成一张灰度图。
def img_from_labels(labels,m,n):
labels = labels.reshape((m,n))
pic_new = Image.new("L",(m,n))
#根据类别向图片中添加灰度值
for i in range(m):
for j in range(n):
pic_new.putpixel((i,j),int(256/(labels[i][j] + 1)))
return pic_new调整聚类数量
Matplotlib 将生成的灰度图绘制出来。
fig, ax = plt.subplots(figsize=(18, 10)) #设置图片大小
img = Image.open(path) #显示原图
plt.subplot(2,3,1)
plt.title("原图")
plt.imshow(img)
plt.box(False) #去掉边框
plt.axis("off")#不显示坐标轴
for i in range(2,7):
plt.subplot(2,3,i)
plt.title("k=" + str(i))
labels, _ = k_means(img_df,i)
pic_new = img_from_labels(labels,m,n)
plt.imshow(pic_new)
plt.box(False) #去掉边框
plt.axis("off")#不显示坐标轴
3 使用 K-Means 进行中文新闻聚类
首先，我们使用 Pandas 的 read_csv 方法将中文新闻语料读取进来， encoding 参数和 sep 参数分别设置为 "utf8" 和 "t" 。
import pandas as pd
news = pd.read_csv("./input/chinese_news_cutted_train_utf8.csv",sep="t",encoding="utf8")
news.head()
news["分类"].value_counts()汽车 544
科技 511
旅游 510
健康 492
文化 491
房地产 480
财经 469
新闻 458
体育 437
教育 437
娱乐 370
女人 322
Name: 分类, dtype: int64
3.1 将新闻表示成向量格式
使用 sklearn.feature_extraction.text 模块的 TfidfVectorizer，将词列表格式的新闻转换成向量形式。向量中的每一个维度代表字典中的一个词，维度的取值代表词在对应文档中的 TF-IDF 取值。
#加载停用词
stop_words = []
file = open("./input/stopwords.txt")
for line in file:
stop_words.append(line.strip())
file.close()
from sklearn.feature_extraction.text import TfidfVectorizer
vectorizer = TfidfVectorizer(stop_words=stop_words,min_df=0.01,max_df=0.5,max_features=500)
news_vectors = vectorizer.fit_transform(news["分词文章"])3.2 使用 k_means 方法对新闻进行聚类
在真实的文本应用中，文本表示成向量后，由于词典的词通常会很大，我们得到的是稀疏矩阵。例如上一步 news_vectors 就是一个稀疏存储格式。
type(news_vectors)scipy.sparse.csr.csr_matrix
在本案例中，为了简便我们直接使用 todense 方法将稀疏向量转换成稠密矩阵。(!注意，实际中很少这么操作！)
news_df = pd.DataFrame(news_vectors.todense())现在我们可以直接使用第一节实现的 k_means 方法对新闻进行聚类，这里我们将聚类个数设置成 12 。
labels, _ = k_means(news_df,12)在 news 中新建 labels 列保存将得到的样本聚类结果。
news["labels"] = labels3.3 使用柱状图查看新闻聚类的主题分布
在新闻数据集中，由于 分类 一列已经标注了新闻的主题类别，我们可以聚类得到的每一组新闻中的主题分布，从而定性地观察聚类效果。将聚类结果按照聚类结果进行分组，可以使用 DataFrame 的 groupby 方法。然后使用 value_counts 方法统计每一组的不同主题的新闻数量分布。
fig, ax = plt.subplots(figsize=(18, 12)) #设置图片大小
cluster_topics = news.groupby(["labels"])["分类"].value_counts()
for i in range(12):
plt.subplot(3,4,i+1)
plt.title("cluster " + str(i))
topic_dist = pd.Series(data=0,index = news["分类"].unique())
topic_dist += cluster_topics[i]
topic_dist.plot(kind="bar")
plt.xticks(rotation=45)
plt.tight_layout()
从上图可以看到，部分分组对应于特定主题，然而在另一上部分分组中，不同主题的新闻分布较为分散。
3.4 使用词云图查看新闻聚类结果
首先，使用 groupby 方法，根据聚类结果对新闻进行分组。
cluster_contents = news[["labels","分词文章"]].groupby(["labels"])将 wordcloud 模块导入词云类 WordCloud，新建一个词云对象 cloud 。为了显示效果，构造词云对象时需要设置停用词词表。由于我们显示的是中文，为了避免乱码，还需要制定显示字体。其中字体文件为 ./input/simfang.ttf 。
from wordcloud import WordCloud
stop_words_set = set(stop_words)
cloud = WordCloud(font_path = './input/simfang.ttf',background_color='white', max_words=100,stopwords=stop_words_set,width=600,height=400)
def series_to_word_list(ts):
results = []
for words in ts.values:
results.extend(words.split(" "))
return " ".join(results)遍历每一个聚类的新闻样本，得到新闻内容的内容，然后显示成词云图。
fig, ax = plt.subplots(figsize=(18, 12)) #设置图片大小
for cluster, group in cluster_contents:
plt.subplot(3,4,cluster+1)
plt.title("cluster " + str(cluster))
wc = cloud.generate_from_text(series_to_word_list(group["分词文章"]))
plt.imshow(wc)
plt.axis("off")
plt.tight_layout()
可以看到，我们可以直接根据词云图对得到的聚类进行解释。
4 总结
本案例中我们使用首先使用三种方式实现了 K-Means 算法并对不同实现的时间性能进行了对比，结果发现向量化实现能够大大提高运行效率。然后我们使用 K-Means 算法进行图像分割，展示了不同 K 取值下生成的灰度图的变化。最后，我们使用 K-Means 在一份中文新闻数据集进行了主题聚类。
本案例使用的主要 Python 工具如下：
→ 机器学习十讲视频回放
数据科学人工智能：视频回放|机器学习第一讲：机器学习介绍zhuanlan.zhihu.com
→ 机器学习十讲直播案例
数据科学人工智能：直播案例 | 使用KNN对新闻主题进行自动分类zhuanlan.zhihu.com
1万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
