KMeans 向量化实现更新集群质心。 Numpy 亲

Question

我在 python 中实现 kmeans。在一次迭代中，我计算了每 150 个点的中心标签：

label = 
array([0, 1, 2, 3, 4, 5, 6, 7, 3, 1, 5, 7, 1, 2, 5, 5, 5, 0, 5, 4, 0, 4,
       6, 7, 7, 1, 7, 0, 0, 3, 3, 0, 5, 5, 1, 1, 0, 4, 3, 7, 0, 1, 3, 7,
       5, 1, 4, 3, 0, 7, 5, 5, 5, 5, 5, 5, 5, 3, 5, 5, 3, 5, 5, 5, 5, 5,
       5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5,
       5, 5, 5, 5, 5, 3, 5, 5, 5, 5, 1, 5, 5, 5, 5, 5, 5, 5, 3, 5, 5, 5,
       5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5,
       5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5], dtype=int64)

和最初的 8 个中心：

centers =
array([[5.1, 3.5, 1.4, 0.2],
       [4.9, 3. , 1.4, 0.2],
       [4.7, 3.2, 1.3, 0.2],
       [4.6, 3.1, 1.5, 0.2],
       [5. , 3.6, 1.4, 0.2],
       [5.4, 3.9, 1.7, 0.4],
       [4.6, 3.4, 1.4, 0.3],
       [5. , 3.4, 1.5, 0.2]])

X 是虹膜数据 X.shape=(150, 4)：

X =
array(  [5.1, 3.8, 1.5, 0.3],
   [5.4, 3.4, 1.7, 0.2],
   [5.1, 3.7, 1.5, 0.4],
   [4.6, 3.6, 1. , 0.2],
   [5.1, 3.3, 1.7, 0.5],
   [4.8, 3.4, 1.9, 0.2],
   [5. , 3. , 1.6, 0.2],
   [5. , 3.4, 1.6, 0.4],
   [5.2, 3.5, 1.5, 0.2],
   [5.2, 3.4, 1.4, 0.2],
   [4.7, 3.2, 1.6, 0.2],
   [4.8, 3.1, 1.6, 0.2],
   [5.4, 3.4, 1.5, 0.4],
            ...

现在我想根据当前的中心标签来更新中心。这意味着迭代标签中的唯一值。然后提取X中的所有对应点，根据提取的所有点计算中心。最后更新中心。例如，在第一次迭代中，提取 X 中标签为 0 的所有元素。然后计算中心（每个维度的平均值）。然后将中心[0] 更新为新中心。以此类推标签 1、2...

这是原始 kmeans 算法的一次迭代。我的问题是如何以 numpy 矢量化方式编写此步骤而不是循环。

Answer 1

更新中心

您可以使用boolean array indexing 和computation along an axis 仅显式地迭代集群而不是每个数据点。

K = 8
for k in range(K):
    centers[k] = X[label==k].mean(axis=0)

更新标签

这也可以通过遍历所有集群来完成：

distances = np.empty(shape=(X.shape[0], K))
for k in range(K):
    distances[:, k] = np.sqrt(np.sum((X - centers[k])**2, axis=1))
labels = distances.argmin(axis=1)

但也可以通过利用矩阵乘法是成对点积来完成，而无需显式循环。

squared_distances = np.sum(centers**2, axis=1) + (np.sum(X**2, axis=1) - 2*centers @ X.T).T
squared_distances[np.isclose(squared_distances, 0)] = 0  # self-distance can become slightly negative with this method (floating point precision problem)
distances = np.sqrt(squared_distances)
labels = distances.argmin(axis=1)