机器学习基础_十七_聚类-k均值

1. 引言

聚类(clustering)是一种无监督学习方法,主要任务是将数据集划分成组,这些组叫作簇(cluster),尽可能使得一个簇内的数据点非常相似且不同簇内的数据点非常不同。与分类算法类似,聚类算法为每个数据点分配(或预测)一个数字,表示这个点属于哪个簇,但不同的是聚类无需预先标注数据而是基于数据内在结构进行划分,聚类多应用于用户画像、图像分割、异常检测等领域。

图片

1.1. 常见聚类算法概述

下面学习几种常见的聚类算法:

[!note] 原图缺失说明
原文此处附有一张 scikit-learn 常见聚类算法比较图(对比 k 均值、层次聚类、DBSCAN 等算法在模拟数据集上的表现),导入时图片丢失。可参考 scikit-learn 官方文档「Comparing different clustering algorithms on synthetic datasets」查看该对比图。

2. k 均值聚类

k 均值是最简单也最常用的聚类算法之一。它试图找到代表数据特定区域的簇中心(cluster center)。算法交替执行以下两个步骤:将每个数据点分配给最近的簇中心,然后将每个簇中心设置为所分配的所有数据点的平均值。如果簇的分配不再发生变化,那么算法结束。

flowchart TD
    A[初始化簇中心] --> B[分配:每个数据点<br>分配给最近的簇中心]
    B --> C[更新:簇中心设为<br>所分配点的平均值]
    C --> D{簇的分配<br>是否变化?}
    D -- 是 --> B
    D -- 否 --> E[算法结束]

下面在一个模拟数据集上对算法进行演示:

import numpy as np
from sklearn.datasets import make_blobs
from sklearn.cluster import KMeans
from sklearn.metrics import pairwise_distances
import matplotlib as mpl
import matplotlib.pyplot as plt
from cycler import cycler
from matplotlib.colors import colorConverter

# 绘制散点图的工具方法
def discrete_scatter(x1, x2, y=None, markers=None, s=10, ax=None, 
                     labels=None, padding=.2, alpha=1, c=None, markeredgewidth=None):
    if ax is None:
        ax = plt.gca()
    if y is None:
        y = np.zeros(len(x1))
    unique_y = np.unique(y)
    if markers is None:
        markers = ['o', '^', 'v', 'D', 's', '*', 'p', 'h', 'H', '8', '<', '>'] * 10
    if len(markers) == 1:
        markers = markers * len(unique_y)
    if labels is None:
        labels = unique_y
    lines = []
    current_cycler = mpl.rcParams['axes.prop_cycle']
    for i, (yy, cycle) in enumerate(zip(unique_y, current_cycler())):
        mask = y == yy
        if c is None:
            color = cycle['color']
        elif len(c) > 1:
            color = c[i]
        else:
            color = c
        if np.mean(colorConverter.to_rgb(color)) < .4:
            markeredgecolor = "grey"
        else:
            markeredgecolor = "black"
        lines.append(ax.plot(x1[mask], x2[mask], markers[i], markersize=s, 
                             label=labels[i], alpha=alpha, c=color, 
                             markeredgewidth=markeredgewidth, 
                             markeredgecolor=markeredgecolor)[0])
    if padding != 0:
        pad1 = x1.std() * padding
        pad2 = x2.std() * padding
        xlim = ax.get_xlim()
        ylim = ax.get_ylim()
        ax.set_xlim(min(x1.min() - pad1, xlim[0]), max(x1.max() + pad1, xlim[1]))
        ax.set_ylim(min(x2.min() - pad2, ylim[0]), max(x2.max() + pad2, ylim[1]))
    return lines

# 绘制在模拟数据集上k均值聚类过程图形工具方法
def plot_kmeans_algorithm():
    X, y = make_blobs(random_state=1)
    with mpl.rc_context(rc={'axes.prop_cycle': cycler('color', ["#9696f1", 
                                                              "#f80404", 
                                                              "#07aa07"])}):
        fig, axes = plt.subplots(3, 3, figsize=(10, 8), subplot_kw={'xticks': (), 'yticks': ()})
        axes = axes.ravel()
        axes[0].set_title("Input data")
        discrete_scatter(X[:, 0], X[:, 1], ax=axes[0], markers=['o'], c='w')
        axes[1].set_title("Initialization")
        init = X[:3, :]
        discrete_scatter(X[:, 0], X[:, 1], ax=axes[1], markers=['o'], c='w')
        discrete_scatter(init[:, 0], init[:, 1], [0, 1, 2], ax=axes[1], 
                         markers=['^'], markeredgewidth=2)
        axes[2].set_title("Assign Points (1)")
        km = KMeans(n_clusters=3, init=init, max_iter=1, n_init=1).fit(X)
        centers = km.cluster_centers_
        labels = np.argmin(pairwise_distances(init, X), axis=0)
        discrete_scatter(X[:, 0], X[:, 1], labels, markers=['o'], 
                         ax=axes[2])
        discrete_scatter(init[:, 0], init[:, 1], [0, 1, 2], 
                         ax=axes[2], markers=['^'], markeredgewidth=2)
        axes[3].set_title("Recompute Centers (1)")
        discrete_scatter(X[:, 0], X[:, 1], labels, markers=['o'], 
                         ax=axes[3])
        discrete_scatter(centers[:, 0], centers[:, 1], [0, 1, 2], 
                         ax=axes[3], markers=['^'], markeredgewidth=2)
        axes[4].set_title("Reassign Points (2)")
        km = KMeans(n_clusters=3, init=init, max_iter=1, n_init=1).fit(X)
        labels = km.labels_
        discrete_scatter(X[:, 0], X[:, 1], labels, markers=['o'], 
                         ax=axes[4])
        discrete_scatter(centers[:, 0], centers[:, 1], [0, 1, 2], 
                         ax=axes[4], markers=['^'], markeredgewidth=2)
        km = KMeans(n_clusters=3, init=init, max_iter=2, n_init=1).fit(X)
        axes[5].set_title("Recompute Centers (2)")
        centers = km.cluster_centers_
        discrete_scatter(X[:, 0], X[:, 1], labels, markers=['o'], 
                         ax=axes[5])
        discrete_scatter(centers[:, 0], centers[:, 1], [0, 1, 2], 
                         ax=axes[5], markers=['^'], markeredgewidth=2)
        axes[6].set_title("Reassign Points (3)")
        labels = km.labels_
        discrete_scatter(X[:, 0], X[:, 1], labels, markers=['o'], 
                         ax=axes[6])
        markers = discrete_scatter(centers[:, 0], centers[:, 1], [0, 1, 2], 
                                   ax=axes[6], markers=['^'], 
                                   markeredgewidth=2)
        axes[7].set_title("Recompute Centers (3)")
        km = KMeans(n_clusters=3, init=init, max_iter=3, n_init=1).fit(X)
        centers = km.cluster_centers_
        discrete_scatter(X[:, 0], X[:, 1], labels, markers=['o'], 
                         ax=axes[7])
        discrete_scatter(centers[:, 0], centers[:, 1], [0, 1, 2], 
                         ax=axes[7], markers=['^'], markeredgewidth=2)
        axes[8].set_axis_off()
        axes[8].legend(markers, ["Cluster 0", "Cluster 1", "Cluster 2"], loc='best')

# 调用方法进行聚类并绘图
plot_kmeans_algorithm()
plt.show()

输出图形结果:

图片

上图是输入数据(Input data)与 k 均值算法的三个步骤。簇中心用三角形表示,数据点用圆形表示。颜色表示簇成员。在上面的代码中,实例化 KMeans 算法时,用 n_clusters=3 参数指定要寻找三个簇,所以通过声明三个随机数据点为簇中心来将算法初始化(上图中"Initialization")。然后开始迭代算法。首先,每个数据点被分配给距离最近的簇中心(上图中"Assign Points (1)“/“分配数据点(1)”)。接下来,将簇中心修改为所分配点的平均值(上图中"Recompute Centers (1)”/“重新计算中心(1)”)。然后将这一过程再重复两次。在第三次迭代之后,为簇中心分配的数据点保持不变,因此算法结束。

给定新的数据点,k 均值会将其分配给最近的簇中心。下面的代码示例展示了上图中学到的簇中心的边界:

import numpy as np
from sklearn.datasets import make_blobs
from sklearn.cluster import KMeans
import matplotlib as mpl
import matplotlib.pyplot as plt
from matplotlib.colors import ListedColormap, colorConverter

cm3 = ListedColormap(["#9696f1","#f80404","#07aa07"])

# 绘制散点图的工具方法
def discrete_scatter(x1, x2, y=None, markers=None, s=10, ax=None, 
                     labels=None, padding=.2, alpha=1, c=None, markeredgewidth=None):
    if ax is None:
        ax = plt.gca()
    if y is None:
        y = np.zeros(len(x1))
    unique_y = np.unique(y)
    if markers is None:
        markers = ['o', '^', 'v', 'D', 's', '*', 'p', 'h', 'H', '8', '<', '>'] * 10
    if len(markers) == 1:
        markers = markers * len(unique_y)
    if labels is None:
        labels = unique_y
    lines = []
    current_cycler = mpl.rcParams['axes.prop_cycle']
    for i, (yy, cycle) in enumerate(zip(unique_y, current_cycler())):
        mask = y == yy
        if c is None:
            color = cycle['color']
        elif len(c) > 1:
            color = c[i]
        else:
            color = c
        if np.mean(colorConverter.to_rgb(color)) < .4:
            markeredgecolor = "grey"
        else:
            markeredgecolor = "black"
        lines.append(ax.plot(x1[mask], x2[mask], markers[i], markersize=s, 
                             label=labels[i], alpha=alpha, c=color, 
                             markeredgewidth=markeredgewidth, 
                             markeredgecolor=markeredgecolor)[0])
    if padding != 0:
        pad1 = x1.std() * padding
        pad2 = x2.std() * padding
        xlim = ax.get_xlim()
        ylim = ax.get_ylim()
        ax.set_xlim(min(x1.min() - pad1, xlim[0]), max(x1.max() + pad1, xlim[1]))
        ax.set_ylim(min(x2.min() - pad2, ylim[0]), max(x2.max() + pad2, ylim[1]))
    return lines

# 绘制决策边界工具方法
def plot_2d_classification(classifier, X, fill=False, ax=None, eps=None, 
                           alpha=1, cm=cm3):
    if eps is None:
        eps = X.std() / 2.
    if ax is None:
        ax = plt.gca()
    x_min, x_max = X[:, 0].min() - eps, X[:, 0].max() + eps
    y_min, y_max = X[:, 1].min() - eps, X[:, 1].max() + eps
    xx = np.linspace(x_min, x_max, 1000)
    yy = np.linspace(y_min, y_max, 1000)
    X1, X2 = np.meshgrid(xx, yy)
    X_grid = np.c_[X1.ravel(), X2.ravel()]
    decision_values = classifier.predict(X_grid)
    ax.imshow(decision_values.reshape(X1.shape), extent=(x_min, x_max, 
                                                         y_min, y_max), 
             aspect='auto', origin='lower', alpha=alpha, cmap=cm)
    ax.set_xlim(x_min, x_max)
    ax.set_ylim(y_min, y_max)
    ax.set_xticks(())
    ax.set_yticks(())

# 绘制kmeans边界图形工具方法
def plot_kmeans_boundaries():
    X, y = make_blobs(random_state=1)
    init = X[:3, :]
    km = KMeans(n_clusters=3, init=init, max_iter=2, n_init=1).fit(X)
    discrete_scatter(X[:, 0], X[:, 1], km.labels_, markers=['o'])
    discrete_scatter(km.cluster_centers_[:, 0], km.cluster_centers_[:, 1], 
                     [0, 1, 2], markers=['^'], markeredgewidth=2)
    plot_2d_classification(km, X, cm=cm3, alpha=.4)

# 调用方法进行聚类并绘图
plot_kmeans_boundaries()
plt.show()

输出 k 均值算法找到的簇中心和簇边界:

图片

2.1. 使用 scikit-learn 实现 k 均值聚类

在 scikit-learn 库中使用 k 均值聚类算法。将 k 均值聚类算法应用于上图中的模拟数据,首先将 KMeans 类实例化,设置要寻找的簇个数 3(如果不指定 n_clusters,它的默认值是 8)。然后对数据调用 fit 方法:

from sklearn.datasets import make_blobs
from sklearn.cluster import KMeans

# 生成模拟的二维数据
X, y = make_blobs(random_state=1)

# 构建聚类模型
kmeans = KMeans(n_clusters=3)
kmeans.fit(X)

# 算法会为X中的每个训练数据点分配一个簇标签。标签存储在kmeans.labels_属性中
print("Cluster memberships:\n{}".format(kmeans.labels_))

输出:

Cluster memberships:[0 1 1 1 2 2 2 1 0 0 1 1 2 0 2 2 2 0 1 1 2 1 2 0 1 2 2 0 0 2 0 0 2 0 1 2 1
1 1 2 2 1 0 1 1 2 0 0 0 0 1 2 2 2 0 2 1 1 0 0 1 2 2 1 1 2 0 2 0 1 1 1 2 0
0 1 2 2 0 1 0 1 1 2 0 0 0 0 1 0 2 0 0 1 1 2 2 0 2 0]

2.2. 用 predict 分配新数据点

因为要找的是 3 个簇,所以簇的编号是 0 到 2。可以用 predict 方法为新数据点分配簇标签。预测时会将最近的簇中心分配给每个新数据点,但现有模型不会改变。对训练集运行 predict 会返回与 labels_ 相同的结果:

from sklearn.datasets import make_blobs
from sklearn.cluster import KMeans

# 生成模拟的二维数据
X, y = make_blobs(random_state=1)

# 构建聚类模型
kmeans = KMeans(n_clusters=3)
kmeans.fit(X)

# 算法会为X中的每个训练数据点分配一个簇标签。标签存储在kmeans.labels_属性中
print("Cluster memberships:\n{}".format(kmeans.labels_))

# 对训练数据点应用predict方法,结果与kmeans.labels_相同
print(kmeans.predict(X))

输出结果:

Cluster memberships:[0 2 2 2 1 1 1 2 0 0 2 2 1 0 1 1 1 0 2 2 1 2 1 0 2 1 1 0 0 1 0 0 1 0 2 1 2
2 2 1 1 2 0 2 2 1 0 0 0 0 2 1 1 1 0 1 2 2 0 0 2 1 1 2 2 1 0 1 0 2 2 2 1 0
0 2 1 1 0 2 0 2 2 1 0 0 0 0 2 0 1 0 0 2 2 1 1 0 1 0][0 2 2 2 1 1 1 2 0 0 2 2 1 0 1 1 1 0 2 2 1 2 1 0 2 1 1 0 0 1 0 0 1 0 2 1 2
2 2 1 1 2 0 2 2 1 0 0 0 0 2 1 1 1 0 1 2 2 0 0 2 1 1 2 2 1 0 1 0 2 2 2 1 0
0 2 1 1 0 2 0 2 2 1 0 0 0 0 2 0 1 0 0 2 2 1 1 0 1 0]

可以看到,聚类算法与分类算法有些相似,每个元素都有一个标签。但并不存在真实的标签,因此标签本身并没有先验意义。以人脸图像聚类为例:聚类的结果可能是,算法找到的第 3 个簇仅包含某个人的面孔,但只有在查看图片之后才能知道这一点,而且数字 3 是任意的。算法给的唯一信息就是所有标签为 3 的人脸都是相似的。

对于上面的模拟数据集上运行的聚类算法,意味着不应该为其中一组的标签是 0、另一组的标签是 1 这一事实赋予任何意义,不需要在意它。再次运行该算法可能会得到不同的簇编号,原因在于初始化的随机性质。

3. k 均值与分解方法的比较

k 均值是一种聚类算法,但在 k 均值和分解方法(比如 PCA 和 NMF)之间存在一些相似之处。PCA 试图找到数据中方差最大的方向,而 NMF 试图找到累加的分量,这通常对应于数据的"极值"或"部分"。两种方法都试图将数据点表示为一些分量之和。与之相反,k 均值则尝试利用簇中心来表示每个数据点。可以将其看作仅用一个分量来表示每个数据点,该分量由簇中心给出。这种观点将 k 均值看作是一种分解方法,其中每个点用单一分量来表示,被称为k 均值聚类

方法 表示数据的方式 特点
PCA 方差最大的方向(主成分) 寻找数据中方差最大的方向
NMF 累加的非负分量 通常对应于数据的"极值"或"部分"
k 均值 最近的簇中心(单一分量) 每个数据点仅用一个簇中心来表示

下面比较 PCA、NMF 和 k 均值,分别显示提取的分量,以及利用 100 个分量对测试集中人脸的重建。对于 k 均值,重建就是在训练集中找到的最近的簇中心:

import numpy as np
import matplotlib.pyplot as plt
from sklearn.model_selection import train_test_split
from sklearn.decomposition import PCA
from sklearn.datasets import fetch_lfw_people
from sklearn.decomposition import NMF
from sklearn.cluster import KMeans

people = fetch_lfw_people(min_faces_per_person=20, resize=0.7)
image_shape = people.images[0].shape

mask = np.zeros(people.target.shape, dtype=np.bool_)
for target in np.unique(people.target):
    mask[np.where(people.target == target)[0][:50]] = 1

X_people = people.data[mask]
y_people = people.target[mask]

# 将灰度值缩放到0到1之间,而不是在0到255之间
# 以得到更好的数据稳定性
X_people = X_people / 255.

X_train, X_test, y_train, y_test = train_test_split(X_people, y_people, stratify=y_people, random_state=0)

nmf = NMF(n_components=100, random_state=0)
nmf.fit(X_train)

pca = PCA(n_components=100, random_state=0)
pca.fit(X_train)

kmeans = KMeans(n_clusters=100, random_state=0)
kmeans.fit(X_train)

X_reconstructed_pca = pca.inverse_transform(pca.transform(X_test))
X_reconstructed_kmeans = kmeans.cluster_centers_[kmeans.predict(X_test)]
X_reconstructed_nmf = np.dot(nmf.transform(X_test), nmf.components_)

# 绘图
fig, axes = plt.subplots(3, 5, figsize=(8, 8), subplot_kw={'xticks': (), 'yticks': ()})
fig.suptitle("Extracted Components")
for ax, comp_kmeans, comp_pca, comp_nmf in zip(axes.T, kmeans.cluster_centers_, pca.components_, nmf.components_):
    ax[0].imshow(comp_kmeans.reshape(image_shape))
    ax[1].imshow(comp_pca.reshape(image_shape), cmap='viridis')
    ax[2].imshow(comp_nmf.reshape(image_shape))

axes[0, 0].set_ylabel("kmeans")
axes[1, 0].set_ylabel("pca")
axes[2, 0].set_ylabel("nmf")

fig, axes = plt.subplots(4, 5, subplot_kw={'xticks': (), 'yticks': ()}, figsize=(8, 8))
fig.suptitle("Reconstructions")
for ax, orig, rec_kmeans, rec_pca, rec_nmf in zip(axes.T, X_test, X_reconstructed_kmeans, X_reconstructed_pca, X_reconstructed_nmf):
    ax[0].imshow(orig.reshape(image_shape))
    ax[1].imshow(rec_kmeans.reshape(image_shape))
    ax[2].imshow(rec_pca.reshape(image_shape))
    ax[3].imshow(rec_nmf.reshape(image_shape))

axes[0, 0].set_ylabel("original")
axes[1, 0].set_ylabel("kmeans")
axes[2, 0].set_ylabel("pca")
axes[3, 0].set_ylabel("nmf")

plt.show()

输出以下两个图形:

图片

上图是对比 k 均值的簇中心与 PCA 和 NMF 找到的分量:

图片

上图是利用 100 个分量(或簇中心)的 k 均值、PCA 和 NMF 的图像重建的对比——k 均值的每张图像中仅使用了一个簇中心。

利用 k 均值做矢量量化的一个好处在于,可以用比输入维度更多的簇对数据进行编码。以 two_moons 数据为例,利用 PCA 主成分分析或 NMF 非负矩阵分解,对这个数据无能为力,因为它只有两个维度。使用 PCA 或 NMF 将其降到一维,会完全破坏数据的结构。但通过使用更多的簇中心,可以用 k 均值找到一种更具表现力的表示,看如下代码:

import matplotlib.pyplot as plt
from sklearn.cluster import KMeans
from sklearn.datasets import make_moons

X, y = make_moons(n_samples=200, noise=0.05, random_state=0)

kmeans = KMeans(n_clusters=10, random_state=0)
kmeans.fit(X)

y_pred = kmeans.predict(X)

plt.scatter(X[:, 0], X[:, 1], c=y_pred, s=60, cmap='Paired')
plt.scatter(kmeans.cluster_centers_[:, 0], kmeans.cluster_centers_[:, 1], s=60, marker='^', c=range(kmeans.n_clusters), linewidth=2, cmap='Paired')
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")
print("Cluster memberships:\n{}".format(y_pred))
plt.show()

输出结果和图形:

[!warning] 原图缺失说明
此处对应的 two_moons 数据集聚类结果图(k=10 时用 10 个簇中心对数据编码的散点图)在导入时图片链接损坏(URL 异常重复),已移除。运行上方代码即可得到该图。

4. 结论

k 均值聚类是一种简单而强大的无监督学习算法,特别适用于:

  1. 矢量量化(vector quantization):将相似的数据点自动分组到簇中
  2. 数据分组:通过矢量量化减少数据维度
  3. 数据压缩:将簇中心作为数据的代表性特征

与 PCA 和 NMF 等分解方法相比,k 均值具有以下优势:

  • 可以使用比原始维度更多的簇来表示数据
  • 对于非线性结构的数据(如 two_moons 数据集)表现更好
  • 计算效率高,适合大规模数据集

然而,k 均值也有一些局限性:

  • 需要预先指定簇的数量
  • 对初始中心点选择敏感
  • 对非凸形状的簇效果不佳

在实际应用中,可以结合多种方法,根据具体问题和数据特性选择最合适的算法。

[!success] 核心要点

  • 聚类定义:无监督学习,将数据划分成簇,簇内相似、簇间不同,无需预先标注
  • k 均值两步迭代:分配(数据点归到最近簇中心)→ 更新(簇中心取均值),分配不再变化即收敛
  • 三行实现KMeans(n_clusters=3).fit(X)labels_ 查看簇标签、predict(X) 分配新数据点
  • 分解视角:k 均值用单一簇中心表示每个点,与 PCA(最大方差方向)、NMF(累加分量)同属"分量表示"
  • 矢量量化:可用比维度更多的簇编码数据,two_moons 这类非线性结构也能表达
  • 局限:需指定 k、对初始中心敏感、对非凸簇效果不佳