[k-means++]仕組みの解説とpythonでの自作コード例

def init(k,X_):
    #=======
    #k: クラスタ数
    #X_ : データ点の座標を格納したnp.array
    #=======
    
    clusters = [] #初期重心を管理するリスト

    #1. データから一点を選択
    centroid_id = np.random.choice([i for i in range(len(X_))])
    clusters.append(X_[centroid_id])
    
    X_ = np.delete(X_, centroid_id, axis=0) #選択した点を入力データのリストから除外
    
    # 4. k 個のクラスタ中心を得られるまで計算を繰り返す。
    while len(clusters) < k:
        dists = []

        #2. 各データ点 と各クラスタ中心との距離を計算し、最も近いものを取り出す
        for x_i in range(len(X_)):
            dist = float('INF')
            for cent_j in range(len(clusters)):
                dist = min(dist, distance(X_[x_i], clusters[cent_j]) ) #今まで見た最短距離と今見てる距離との小さい方を選択
            dists.append(dist)
        
        len_dist = len(dists)
        
        #3. 確率分布にしたがって、データ点を1点選択してクラスタ中心とする。
        prob = dists / np.sum(dists) #確率分布
        
        new_c = np.random.choice([i for i in range(len_dist)], p=prob) #確率分布から点を一点取り出す
        clusters.append(X_[new_c])
        
        X_ = np.delete(X_, new_c, axis=0) #取り出した要素を消去
        #距離のリストを初期化
        dists = []
  
    return clusters

import numpy as np
import matplotlib.pyplot as plt
from sklearn.cluster import KMeans
from sklearn.datasets import make_blobs
from sklearn.cluster import kmeans_plusplus
import copy

np.random.seed(16)

X_, _ = make_blobs(1000, 2,centers=3, random_state=10, cluster_std=1, center_box=(-20,20))
plt.scatter(X_[:,0], X_[:,1])

#自作コード
centers = init(k, X_)

plt.scatter(X_[:,0], X_[:,1]) #データ点の散布図
for i in range(k):
  c = centers[i]
  #print(type(c), c)
  plt.scatter(c[0], c[1], label=f'{i}', marker="X", s=100) #初期重心の位置

plt.legend()

#sklearn
centers_sk, _ = kmeans_plusplus(X, n_clusters=k, n_local_trials=1)

plt.scatter(X[:,0], X[:,1])
for i in range(k):
  c = centers_sk[i]
  print(type(c), c)
  plt.scatter(c[0], c[1], label=f'{i}', marker="X", s=100)

plt.legend()