自分用のクラスタリングメモ
クラスタリングの実装を覚えたくて自分用のメモにブログとして残しました。
コードはこちら。
データセットは Student Performance Data Set を使用しています。
なお文中のdf_dummyは、pd.get_dummiesを使ってOneHotEncodingしたものになります。
実装メインなので、理論は各リンク先の記事を参照ください。
k-means法＋PCAによる次元削減
from sklearn.cluster import KMeans
k = 3
SEED = 42
kmeans_model = KMeans(n_clusters=k, random_state=SEED).fit(df_dummy.select_dtypes(exclude='object'))
labels = kmeans_model.labels_
colors = [plt.get_cmap("tab10")(x) for x in labels]
PCAを使うと、累積寄与率（この主成分軸一つで、データの何割を説明することができているか）の高い順に主成分得点を得ることができます。
display(pd.DataFrame(pca.explained_variance_ratio_, index=[f"PC{i}" for i in range(59)]).T)
クラスタリングの確認のため、PCAによる第1主成分軸をX軸、第2主成分軸をY軸に散布図を表示してみます。
from sklearn.decomposition import PCA
pca = PCA(copy=True, n_components=None, whiten=False)
feature = pca.fit_transform(df_dummy.select_dtypes(exclude='object'))
# グラフ化
plt.figure(figsize=(6, 6))
plt.scatter(feature[:, 0], feature[:, 1], alpha=0.8, color=colors)
plt.show()
最適なクラスター数の調べ方
k-means法によるクラスタリングはクラスタ数を指定する必要があるが、クラスタ数をいくつにするかは明確な基準がありません。
完璧な方法ではないが、「エルボー法」と「シルエット分析」により適切なクラスタ数の検討が可能です。
エルボー法
distortions = []
num = 20
for i in range(1,num+1):
km = KMeans(n_clusters=i,
init='k-means++',
n_init=10,
max_iter=300,
random_state=0)
km.fit(df_dummy.select_dtypes(exclude='object'))
distortions.append(km.inertia_) # km.fitするとkm.inertia_が得られる
plt.plot(range(1,num+1),distortions,marker='o')
plt.xlabel('Number of clusters')
plt.ylabel('Distortion')
plt.show()
Y軸のSSE値が"ヒジ"のように鋭角に曲がった点が最適なクラスター数らしいです。
クラスタ数5が適切と判断し、クラスタ数5のk-means法のクラスタリングを見ると以下のようになります。
シルエット法
クラスタ内のサンプルがどの程度密にグルーピングされているか、グラフで可視化しました。
クラスタ数を2～20で変化させてグラフ化してみました。
from sklearn.metrics import silhouette_samples
from matplotlib import cm
SEED = 42
num = 20
fig = plt.figure(figsize=(10,20))
for k in range(2, num+1):
kmeans_num = KMeans(n_clusters=k, random_state=SEED).fit_predict(df_dummy.select_dtypes(exclude='object'))
cluster_labels = np.unique(kmeans_num) # y_kmの要素の中で重複を無くす
n_clusters=cluster_labels.shape[0] # 配列の長さを返す。つまりここでは n_clustersで指定した3となる
# シルエット係数を計算
silhouette_vals = silhouette_samples(df_dummy.select_dtypes(exclude='object'),
kmeans_num,metric='euclidean') # サンプルデータ, クラスター番号、ユークリッド距離でシルエット係数計算
y_ax_lower, y_ax_upper= 0,0
yticks = []
# グラフ
plt.subplot((num//3+1), 3, k-1)
for i,c in enumerate(cluster_labels):
c_silhouette_vals = silhouette_vals[kmeans_num==c] # cluster_labelsには 0,1,2が入っている（enumerateなのでiにも0,1,2が入ってる（たまたま））
c_silhouette_vals.sort()
y_ax_upper += len(c_silhouette_vals) # サンプルの個数をクラスターごとに足し上げてy軸の最大値を決定
color = cm.jet(float(i)/n_clusters) # 色の値を作る
plt.barh(range(y_ax_lower,y_ax_upper), # 水平の棒グラフのを描画（底辺の範囲を指定）
c_silhouette_vals, # 棒の幅（1サンプルを表す）
height=1.0, # 棒の高さ
edgecolor='none', # 棒の端の色
color=color) # 棒の色
yticks.append((y_ax_lower+y_ax_upper)/2) # クラスタラベルの表示位置を追加
y_ax_lower += len(c_silhouette_vals) # 底辺の値に棒の幅を追加
silhouette_avg = np.mean(silhouette_vals) # シルエット係数の平均値
plt.axvline(silhouette_avg,color="red",linestyle="--") # 係数の平均値に破線を引く
plt.title(f"cluster {k}")
plt.yticks(yticks,cluster_labels + 1) # クラスタレベルを表示
plt.ylabel('Cluster')
plt.xlabel('silhouette coefficient')
plt.tight_layout() # グラフ間の隙間を調整して座標部分が重なりを解消できる
全てのクラスタで、ある程度のサンプルが赤の点線(全クラスタの平均シルエット係数)を超えていればよいクラスタリングらしいです。
また、各クラスタの縦の長さがだいたい等しければ、よいクラスタリングらしいです。
クラスタ数10が適切と判断し、クラスタ数10のk-means法のクラスタリングを見ると以下のようになります。
その他の次元削減
その他の次元削減手法として、t-SNEとUMAPを試してみました。
t-SNE で次元削減
perplexityという引数の違いで結果が大きく変わるので、perplexityを2～10で変化させてみました。
from sklearn.manifold import TSNE
fig = plt.figure(figsize=(10,10))
num = 10
for i in range(2, num):
plt.subplot((num//3+1), 3, i-1)
plt.title(f"perplexity: {i}")
model_tsne = TSNE(n_components=2, perplexity=i)
feature = model_tsne.fit_transform(df_dummy.select_dtypes(exclude='object'))
plt.scatter(feature[:, 0], feature[:, 1], alpha=0.8, color=colors)
plt.tight_layout()
UMAPでクラスタリング
n_neighborsという引数の違いで結果が大きく変わるので、n_neighborsを2～10で変化させてみました。
import umap.umap_ as umap
from scipy.sparse.csgraph import connected_components
fig = plt.figure(figsize=(10,10))
num = 10
for i in range(2, num):
plt.subplot((num//3+1), 3, i-1)
plt.title(f"n_neighbors: {i}")
model_umap = umap.UMAP(n_components=2,n_neighbors=i)
feature = model_umap.fit_transform(df_dummy.select_dtypes(exclude='object'))
plt.scatter(feature[:, 0], feature[:, 1], alpha=0.8, color=colors)
plt.tight_layout()
余談ですが、UMAPをインストールするときはコマンドに注意しましょう。
pip install umap # 不正解
pip install umap-learn # 正解
以上になります、最後までお読みいただきありがとうございました。
Discussion