K-means聚类详解

K-means 算法是一种基于相似属性将一组数据点划分为不同集群或组的方法。它是一种无监督学习算法，这意味着它不需要标记数据来查找数据集中的模式。

  • K-Means使用步骤
  • 初始中心点怎么确定
  • K值怎么确定
  • 「K均值聚类」的Tips

什么是聚类？

聚类的目标是将项目分成组，使得组中的对象比组外的对象更相似，我们可以在 k-means 中使用我们想要的任何相似度函数来比较两个点。

如何定义集群中的相似性？

比较两个数据点相似性的最常见和最直接的方法是使用距离。与使用勾股定理求对角线长度的方式相同，我们可以计算两个数据点之间的距离。

这种比较相似度的方法称为欧几里得距离。K-means 聚类中另一个常用的相似度函数称为曼哈顿距离。

这两个距离公式对用于聚类任务的各种数据集都很有效。然而，给定的相似度函数可能并非对数据点的每个分布都有效。这就提出了一个问题，“一个好的相似度函数的特性是什么？”。

良好相似函数的特征

可以使用三个基本属性来评估相似度函数是否良好。当然，还有更多的标准可以考虑，但这些是最重要的。在符号下方d(x,y)_d_ ( _x_ ,_和_) 读作“x 和 y 之间的距离”。

对称性： $d(x,y) = d(y,x)$

对称性很重要，否则我们可以说 A 看起来像 B，但 B 看起来一点也不像 A。

正可分离性：$d(x,y) = 0 \quad if \quad and \quad only \quad if \quad x=y$

正可分性的属性很重要，因为否则可能会有两个不同的数据点 A 和 B，我们无法区分。

三角不等式：$d(x,y) \leq d(x,z) + d(z,y)$

三角不等式很重要，否则我们可以说 A 看起来像 B，B 看起来像 C，但 A 看起来不像 C。

常用聚类方法概述

本文重点介绍 K 均值聚类，但这只是现有的众多聚类算法之一。

1.分区聚类

K-means 是分区聚类算法（也称为基于质心的聚类）的一个例子。这意味着它接受用户提供的集群数量_k_，并将数据划分为许多分区。在分区聚类中，每个数据点只能属于一个集群，任何集群都不能为空。

许多分区聚类算法是非确定性的。这意味着，即使您保持输入固定，每次运行算法时也可能最终得到不同的集群。

2.层次聚类

层次聚类下的算法通过自上而下或自下而上构建层次结构来将对象分配给集群。

  • 自顶向下的方法称为分裂聚类。它的工作原理是从一个集群中的所有点开始，然后在每一步拆分最不相似的集群，直到每个数据点都在一个单独的集群中。
  • 自底向上的方法称为凝聚聚类。这种方法迭代地合并集群中两个最相似的点，直到只有一个大集群。

与分区聚类方法不同，层次聚类是确定性的。这意味着集群分配在同一数据集上的运行之间不会有所不同。

自上而下和自下而上的方法都会产生一个基于树的点层次结构，称为 dendrogram。事实证明，此树状图可用于选择要使用的集群数量。您可以在任何深度切割树以生成不相交树的子集，每个树代表一个集群。

比较分区聚类和层次聚类

两种聚类方法各有优缺点。首先我们来看看使用K-means（partitional clustering）的优缺点。

K-means 的优点 K-means 的缺点
易于实施 难以预测聚类的数量（K-Value）
k-Means 可能比具有大量变量的层次聚类更快（如果 K 很小） 初始种子对最终结果有很大影响
k-Means 可能产生比层次聚类更紧密的聚类 数据的顺序对最终结果有影响
重新计算质心时，实例可以更改集群（移动到另一个集群） 对缩放敏感：重新缩放数据集（标准化或标准化）将 完全改变结果。

现在，让我们来看看层次聚类的优点和缺点。

层次聚类的优点 层次聚类的缺点
输出比 k-means 返回的非结构化集群集信息更多的层次结构。 一旦一个点被分配给一个 集群，它就不能再四处移动了。
易于实施 时间复杂度：不适合大数据集
初始种子对最终结果有很大影响
数据的顺序对最终结果有影响

我已经概述了两个聚类算法系列。还有其他类型的聚类算法我没有介绍，例如基于密度的聚类。有关聚类算法的更完整概述，请访问此资源。

K-means Clustering 如何在视觉上工作？

此可视化显示了使用 3 个集群运行的 k-means 算法。首先，初始化三个质心。这些是每个集群的初始中心，由蓝色、红色和绿色大点表示。

接下来，计算给定数据点与三个质心中的每一个之间的距离。这是针对所有数据点完成的。每个点都分配给它最接近的质心的集群。当我在可视化中单击重新分配点时会发生这种情况。

然后，通过对相应集群中每个数据点的坐标求平均值来重新计算质心。当我单击Update Centroids时会发生这种情况。

这个重新分配点和更新质心的过程一直持续到质心不再移动。

您可以使用此可视化自己直观地了解 k-means！在这里试试。

K-Means伪代码

K-Means聚类步骤是一个循环迭代的算法，非常简单易懂：

  1. 假定我们要对N个样本观测做聚类，要求聚为K类，首先选择K个点作为初始中心点；
  2. 接下来，按照距离初始中心点最小的原则，把所有观测分到各中心点所在的类中；
  3. 每类中有若干个观测，计算K个类中所有样本点的均值，作为第二次迭代的K个中心点；
  4. 然后根据这个中心重复第2、3步，直到收敛（中心点不再改变或达到指定的迭代次数），聚类过程结束。

⭐️使用 sklearn 中的 K-Means

https://github.com/ljpzzz/machinelearning/blob/master/classic-machine-learning/kmeans_cluster.ipynb

KMeans类的主要参数有：

  1. n_clusters: 即我们的k值，一般需要多试一些值以获得较好的聚类效果。k值好坏的评估标准在下面会讲。
  2. max_iter： 最大的迭代次数，一般如果是凸数据集的话可以不管这个值，如果数据集不是凸的，可能很难收敛，此时可以指定最大的迭代次数让算法可以及时退出循环。
  3. n_init：用不同的初始化质心运行算法的次数。由于K-Means是结果受初始值影响的局部最优的迭代算法，因此需要多跑几次以选择一个较好的聚类效果，默认是10，一般不需要改。如果你的k值较大，则可以适当增大这个值。
  4. init： 即初始值选择的方式，可以为完全随机选择’random’,优化过的’k-means++’或者自己指定初始化的k个质心。一般建议使用默认的’k-means++’。
  5. algorithm：有“auto”, “full” or “elkan”三种选择。”full”就是我们传统的K-Means算法， “elkan”是我们原理篇讲的elkan K-Means算法。默认的”auto”则会根据数据值是否是稀疏的，来决定如何选择”full”和“elkan”。一般数据是稠密的，那么就是 “elkan”，否则就是”full”。一般来说建议直接用默认的”auto”
  6. n_jobs：表示任务使用CPU数量
  7. random_state：表示随机数生成器的种子。
  8. verbose：0表示不输出日志信息；1表示每隔一段时间打印一次日志信息。如果大于1，打印次数频繁。
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets.samples_generator import make_blobs

# X为样本特征，Y为样本簇类别， 共1000个样本，每个样本2个特征，共4个簇，簇中心在[-1,-1], [0,0],[1,1], [2,2]， 簇方差分别为[0.4, 0.2, 0.2]
X, y = make_blobs(n_samples=1000, n_features=2, centers=[[-1,-1], [0,0], [1,1], [2,2]], cluster_std=[0.4, 0.2, 0.2, 0.2], 
                  random_state =9)
plt.scatter(X[:, 0], X[:, 1], marker='o')
plt.show()

# 训练 & 预测
from sklearn.cluster import KMeans
km = KMeans(n_clusters=2, random_state=9, n_init=10, max_iter=300, init='k-means++', verbose=1)
y_pred = km.fit_predict(X)
plt.scatter(X[:, 0], X[:, 1], c=y_pred)
plt.show()

# 观察在不同的k值下Calinski-Harabasz分数，分数值𝑠越大则聚类效果越好
from sklearn import metrics
s = metrics.calinski_harabaz_score(X, y_pred)

如何在 Python 中从头开始编写 K-means？

我们的 k-means 实现将分为五个辅助方法和一个运行算法的主循环。让我们一一介绍这些功能。

  • 成对距离
  • 初始化中心
  • 更新分配
  • 更新中心
  • 计算损失
  • 主循环
  • 完整的实现

计算成对距离

1
2
3
4
5
6
7
8
9
10
11
12
13
14
def pairwise_dist(self, x, y):
        """
        Args:
            x: N x D numpy array
            y: M x D numpy array
        Return:
            dist: N x M array, where dist2[i, j] is the euclidean distance between 
            x[i, :] and y[j, :]
        """
        xSumSquare = np.sum(np.square(x),axis=1);
        ySumSquare = np.sum(np.square(y),axis=1);
        mul = np.dot(x, y.T);
        dists = np.sqrt(abs(xSumSquare[:, np.newaxis] + ySumSquare-2*mul))
        return dists

pairwise_dist 函数与前面描述的相似性函数等效。这是我们比较两点相似性的指标。在这里，我使用的是欧几里得距离。我使用的公式可能看起来与欧几里德距离函数的常规公式不同。这是因为我们正在执行矩阵操作，而不是使用两个单个向量。在这里深入阅读。

初始化中心

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
def _init_centers(self, points, K, **kwargs):
        """
        Args:
            points: NxD numpy array, where N is # points and D is the dimensionality
            K: number of clusters
            kwargs: any additional arguments you want
        Return:
            centers: K x D numpy array, the centers. 
        """
        row, col = points.shape
        retArr = np.empty([K, col])
        for number in range(K):
            randIndex = np.random.randint(row)
            retArr[number] = points[randIndex]
        
        return retArr

该函数接收点数组并随机选择其中的 K 个作为初始质心。该函数仅返回 K 个选定点。

更新分配

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
def _update_assignment(self, centers, points):
        """
        Args:
            centers: KxD numpy array, where K is the number of clusters, and D is the dimension
            points: NxD numpy array, the observations
        Return:
            cluster_idx: numpy array of length N, the cluster assignment for each point

        Hint: You could call pairwise_dist() function.
        """
        row, col = points.shape
        cluster_idx = np.empty([row])
        distances = self.pairwise_dist(points, centers)
        cluster_idx = np.argmin(distances, axis=1)

        return cluster_idx

更新分配函数负责选择每个点应该属于哪个集群。首先，我使用 pairwise_dist 函数计算每个点和每个质心之间的距离。然后，我得到每一行的最小距离的索引。最小距离的索引也是给定数据点的聚类分配索引，因为我们希望将每个点分配给最近的质心。

更新中心

1
2
3
4
5
6
7
8
9
10
11
12
13
14
def _update_centers(self, old_centers, cluster_idx, points):
        """
        Args:
            old_centers: old centers KxD numpy array, where K is the number of clusters, and D is the dimension
            cluster_idx: numpy array of length N, the cluster assignment for each point
            points: NxD numpy array, the observations
        Return:
            centers: new centers, K x D numpy array, where K is the number of clusters, and D is the dimension.
        """
        K, D = old_centers.shape
        new_centers = np.empty(old_centers.shape)
        for i in range(K):
            new_centers[i] = np.mean(points[cluster_idx == i], axis = 0)
        return new_centers

更新中心功能负责对属于给定集群的所有点进行平均。该平均值是相应聚类的新质心。该函数返回新中心的数组。

计算损失

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
def _get_loss(self, centers, cluster_idx, points):
        """
        Args:
            centers: KxD numpy array, where K is the number of clusters, and D is the dimension
            cluster_idx: numpy array of length N, the cluster assignment for each point
            points: NxD numpy array, the observations
        Return:
            loss: a single float number, which is the objective function of KMeans. 
        """
        dists = self.pairwise_dist(points, centers)
        loss = 0.0
        N, D = points.shape
        for i in range(N):
            loss = loss + np.square(dists[i][cluster_idx[i]])
        
        return loss

损失函数是我们评估聚类算法性能的指标。我们的损失只是每个点与其聚类质心之间的平方距离之和。在我们的实现中，我们首先调用成对距离来获得每个点和每个中心之间的距离矩阵。我们使用 cluster_idx 为每个点选择与集群对应的适当距离2。

主循环

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
def __call__(self, points, K, max_iters=100, abs_tol=1e-16, rel_tol=1e-16, verbose=False, **kwargs):
        """
        Args:
            points: NxD numpy array, where N is # points and D is the dimensionality
            K: number of clusters
            max_iters: maximum number of iterations (Hint: You could change it when debugging)
            abs_tol: convergence criteria w.r.t absolute change of loss
            rel_tol: convergence criteria w.r.t relative change of loss
            verbose: boolean to set whether method should print loss (Hint: helpful for debugging)
            kwargs: any additional arguments you want
        Return:
            cluster assignments: Nx1 int numpy array
            cluster centers: K x D numpy array, the centers
            loss: final loss value of the objective function of KMeans
        """
        centers = self._init_centers(points, K, **kwargs)
        for it in range(max_iters):
            cluster_idx = self._update_assignment(centers, points)
            centers = self._update_centers(centers, cluster_idx, points)
            loss = self._get_loss(centers, cluster_idx, points)
            K = centers.shape[0]
            if it:
                diff = np.abs(prev_loss - loss)
                if diff < abs_tol and diff / prev_loss < rel_tol:
                    break
            prev_loss = loss
            if verbose:
                print('iter %d, loss: %.4f' % (it, loss))
        return cluster_idx, centers, loss

现在，在主循环中，我们可以组合所有实用函数并实现伪代码。首先，中心用_init_centers随机初始化。然后，对于指定的迭代次数，我们重复update_assignment和update_centers步骤。每次迭代后，我们计算总损失并将其与之前的损失进行比较。如果差异小于我们的阈值，则算法执行完成。

完整的 K-Means 类实现

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
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
from __future__ import absolute_import
from __future__ import print_function
from __future__ import division

import sys
import matplotlib
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import axes3d
from tqdm import tqdm
# Load image
import imageio


# Set random seed so output is all same
np.random.seed(1)


class KMeans(object):

    def __init__(self):  # No need to implement
        pass

    def pairwise_dist(self, x, y):  # [5 pts]
        """
        Args:
            x: N x D numpy array
            y: M x D numpy array
        Return:
                dist: N x M array, where dist2[i, j] is the euclidean distance between 
                x[i, :] and y[j, :]
                """
        xSumSquare = np.sum(np.square(x),axis=1);
        ySumSquare = np.sum(np.square(y),axis=1);
        mul = np.dot(x, y.T);
        dists = np.sqrt(abs(xSumSquare[:, np.newaxis] + ySumSquare-2*mul))
        return dists

    def _init_centers(self, points, K, **kwargs):  # [5 pts]
        """
        Args:
            points: NxD numpy array, where N is # points and D is the dimensionality
            K: number of clusters
            kwargs: any additional arguments you want
        Return:
            centers: K x D numpy array, the centers. 
        """
        row, col = points.shape
        retArr = np.empty([K, col])
        for number in range(K):
            randIndex = np.random.randint(row)
            retArr[number] = points[randIndex]
        
        return retArr

    def _update_assignment(self, centers, points):  # [10 pts]
        """
        Args:
            centers: KxD numpy array, where K is the number of clusters, and D is the dimension
            points: NxD numpy array, the observations
        Return:
            cluster_idx: numpy array of length N, the cluster assignment for each point

        Hint: You could call pairwise_dist() function.
        """
        row, col = points.shape
        cluster_idx = np.empty([row])
        distances = self.pairwise_dist(points, centers)
        cluster_idx = np.argmin(distances, axis=1)

        return cluster_idx

    def _update_centers(self, old_centers, cluster_idx, points):  # [10 pts]
        """
        Args:
            old_centers: old centers KxD numpy array, where K is the number of clusters, and D is the dimension
            cluster_idx: numpy array of length N, the cluster assignment for each point
            points: NxD numpy array, the observations
        Return:
            centers: new centers, K x D numpy array, where K is the number of clusters, and D is the dimension.
        """
        K, D = old_centers.shape
        new_centers = np.empty(old_centers.shape)
        for i in range(K):
            new_centers[i] = np.mean(points[cluster_idx == i], axis = 0)
        return new_centers

    def _get_loss(self, centers, cluster_idx, points):  # [5 pts]
        """
        Args:
            centers: KxD numpy array, where K is the number of clusters, and D is the dimension
            cluster_idx: numpy array of length N, the cluster assignment for each point
            points: NxD numpy array, the observations
        Return:
            loss: a single float number, which is the objective function of KMeans. 
        """
        dists = self.pairwise_dist(points, centers)
        loss = 0.0
        N, D = points.shape
        for i in range(N):
            loss = loss + np.square(dists[i][cluster_idx[i]])
        
        return loss

    def __call__(self, points, K, max_iters=100, abs_tol=1e-16, rel_tol=1e-16, verbose=False, **kwargs):
        """
        Args:
            points: NxD numpy array, where N is # points and D is the dimensionality
            K: number of clusters
            max_iters: maximum number of iterations (Hint: You could change it when debugging)
            abs_tol: convergence criteria w.r.t absolute change of loss
            rel_tol: convergence criteria w.r.t relative change of loss
            verbose: boolean to set whether method should print loss (Hint: helpful for debugging)
            kwargs: any additional arguments you want
        Return:
            cluster assignments: Nx1 int numpy array
            cluster centers: K x D numpy array, the centers
            loss: final loss value of the objective function of KMeans
        """
        centers = self._init_centers(points, K, **kwargs)
        for it in range(max_iters):
            cluster_idx = self._update_assignment(centers, points)
            centers = self._update_centers(centers, cluster_idx, points)
            loss = self._get_loss(centers, cluster_idx, points)
            K = centers.shape[0]
            if it:
                diff = np.abs(prev_loss - loss)
                if diff < abs_tol and diff / prev_loss < rel_tol:
                    break
            prev_loss = loss
            if verbose:
                print('iter %d, loss: %.4f' % (it, loss))
        return cluster_idx, centers, loss

这是我们在 Python 中完整的 K-means 类实现。我鼓励您复制此代码（或自己实现！）并在您自己的机器上运行它。这是巩固您对 K-means 算法理解的最佳方式。scikit-learn中的KMeans聚类实现

初始中心点怎么确定

在k-means算法步骤中，有两个地方降低了SSE：

  1. 把样本点分到最近邻的簇中，这样会降低SSE的值；
  2. 重新优化聚类中心点，进一步的减小了SSE。

这样的重复迭代、不断优化，会找到局部最优解（局部最小的SSE），如果想要找到全局最优解需要找到合理的初始聚类中心。

那合理的初始中心怎么选？

方法有很多，譬如先随便选个点作为第1个初始中心C1，接下来计算所有样本点与C1的距离，距离最大的被选为下一个中心C2，直到选完K个中心。这个算法叫做K-Means++，可以理解为 K-Means的改进版，它可以能有效地解决初始中心的选取问题，但无法解决离群点问题。

K值怎么确定

使用 K 均值算法，性能可能会因您使用的集群数量而有很大差异。要知道，K设置得越大，样本划分得就越细，每个簇的聚合程度就越高，误差平方和SSE自然就越小。所以不能单纯像选择初始点那样，用不同的K来做尝试，选择SSE最小的聚类结果对应的K值，因为这样选出来的肯定是你尝试的那些K值中最大的那个。

“手肘法”（Elbow Method）

为了使用肘部方法，您只需多次运行 K-means 算法，每次迭代将聚类数增加一个。记录每次迭代的损失，然后制作 num cluster vs loss 的折线图。

下面是肘部方法的简单实现：

1
2
from sklearn import metrics
metrics.calinski_harabaz_score(X, y_pred) 
1
2
3
4
5
6
7
8
9
10
from sklearn.cluster import KMeans

loss = []

for i in range(1, 10):
    kmeans = KMeans(n_clusters=i, max_iter=100).fit(p_list)
    loss.append(kmeans.inertia_ / point_number / K)
    
plt.plot(range(1, 10), loss)
plt.show()

运行此方法将输出一个类似于您在下面看到的图：

现在，为了选择正确数量的簇，我们进行目视检查。损耗曲线开始弯曲的点称为 肘点。肘点代表误差和聚类数量之间的合理权衡。在此示例中，肘点位于x = 3 处。这意味着最佳聚类数为 3。

轮廓系数

数据的平均轮廓是评估集群自然数的另一个有用标准。数据实例的轮廓是衡量它与集群内数据匹配程度以及与相邻集群数据匹配程度的度量。

轮廓值是衡量一个对象与其自己的集群（内聚）相比其他集群（分离）的相似程度。轮廓范围从 -1 到 +1，其中高值表示对象与其自己的集群匹配良好，而与相邻集群匹配不佳。如果大多数对象具有较高的值，则集群配置是合适的。如果许多点具有低值或负值，则聚类配置可能具有过多或过少的聚类。

如果您想实现用于聚类分析的轮廓系数，我建议使用 scikit-learn。访问此资源以获取有关实施的完整指南。

正确使用「K均值聚类」的Tips

  1. 输入数据一般需要做缩放，如标准化。原因很简单，K均值是建立在距离度量上的，因此不同变量间如果维度差别过大，可能会造成少数变量“施加了过高的影响而造成垄断”。
  2. 输出结果非固定，多次运行结果可能不同。首先要意识到K-means中是有随机性的，从初始化到收敛结果往往不同。一种看法是强行固定随机性，比如设定sklearn中的random state为固定值。另一种看法是，如果你的K均值结果总在大幅度变化，比如不同簇中的数据量在多次运行中变化很大，那么K均值不适合你的数据，不要试图稳定结果 [2]。
  3. 运行时间往往可以得到优化，选择最优的工具库。基本上现在的K均值实现都是K-means++，速度都不错。但当数据量过大时，依然可以使用其他方法，如MiniBatchKMeans [3]。上百万个数据点往往可以在数秒钟内完成聚类，推荐Sklearn的实现。
  4. 高维数据上的有效性有限。建立在距离度量上的算法一般都有类似的问题，那就是在高维空间中距离的意义有了变化，且并非所有维度都有意义。这种情况下，K均值的结果往往不好，而通过划分子空间的算法（sub-spacing method）效果可能更好。
  5. 运行效率与性能之间的取舍。

小结

因此不难看出，K-Means优点在于原理简单，运行速度快且容易实现，能够处理的数据量大。

当然，也有一些缺点：

  • 在高维上可能不是最佳选项
  • K值、初始点的选取不好确定；
  • 得到的结果只是局部最优；
  • 受离群值影响大。

一个比较粗浅的结论是，在数据量不大时，可以优先尝试其他算法。当数据量过大时，可以试试HDBSCAN。仅当数据量巨大，且无法降维或者降低数量时，再尝试使用K均值。

参考

K-means for Beginners: How to Build from Scratch in Python

#文本聚类
K-means聚类详解
http://example.com/2021/09/06/2021-09-06-K-means聚类详解/
作者
NSX
发布于
2021年9月6日
许可协议