はじめに
こちらの記事ではクラスタリングについて、以下の内容をまとめました。
- K-means法の実行(KMeansクラス)
- エルボー法、シルエット分析(silhouette_score()関数)
- 結果の可視化(PCAクラス)
- 凝集型クラスタリングの実行(AgglomerativeClusteringクラス)
- デンドログラムの描画(dendrogram()関数)
クラスタリング
教師なし学習のひとつで、データの類似性を計算してクラスタというグループに分類すること。各データの所属するクラスタについての正解ラベルが与えられるのではなく、データの特徴を学習してグループ分けをするので、結果の妥当性を人間が判断する必要がある。
クラスタリングは手法によって次の4タイプに分けられる。
クラスタリングの種類
| カテゴリ | 内容 | 主な手法 |
|---|---|---|
| 分割型 | 指定されたk個のクラスタに分割する | K-means法 |
| 階層型 | ツリー構造を作成してクラスタを切り出す | 凝集型 分裂型 |
| 密度ベース型 | 高密度の領域をクラスタとみなし、低密度の領域は外れ値とみなす | DBSCAN |
| モデルベース型 | データが特定の確率分布から生成されたと仮定する | GMM |
※階層的クラスタリングの分裂型を分割型ということもある
k-means法
k-means法は、次の流れでクラスタリングを行う。
- クラスタ数を事前に決める
- クラスタ数の分だけ仮の中心点を用意する
- 各データにとって最も近い中心点を見つけ、同じ中心点に割り当てられたデータを同じクラスタとする
- クラスタに所属する全てのデータと、そのクラスタの中心点との距離の二乗の合計(クラスタ内誤差二乗和,SSE)を計算する
- SSEが最小となるように中心点を動かしながら決めて分類する
k-means法の実装コード
ポルトガルの卸売業者の顧客データを例に使って、この卸売業者の顧客を購買の特徴の類似性から分類する。
実装はこのような流れで行う。
KMeansクラス
k-means法を実行するには、sklearn.clusterモジュールのKMeansクラスクラスの持つfit()メソッドを利用する。
KMeansクラスをインスタンス化する際はクラスタ数を指定する必要がある。クラスタ数の最適な値を求めるには、クラスタ数を変化させながらfit()を実行し、挙動を見て決める。(エルボー法、シルエット分析)
コンストラクタ引数
-
n_clusters(int)：クラスタ数を指定する -
random_state(int)：クラスタ中心点の初期位置を決める乱数のseedを固定する -
n_init：異なる中心点で距離計算を実行する回数(int)
メソッド
| メソッド名 | 内容 |
|---|---|
fit() |
クラスタリングを実行し、その結果を自身の属性labes_にNumpy配列で格納し、自身のオブジェクトを返す。 |
transform() |
すべてのサンプルについて、すべてのクラスタの中心との距離を調べ、2次元Numpy配列に格納して返す。この配列は1行が1つのサンプルについての情報で、各クラスタの中心との距離をすべて列挙したものになるので、行数=サンプル数、列数=クラスター数になる。 |
predict() |
すべてのサンプルについて、どのクラスタに属するかを調べ、1次元Numpy配列にクラスタ番号を格納して返す。この配列は 要素数=サンプル数 となる。 |
fit_transform() |
fit()とtransform()を同時に行う。戻り値はtransform()と同じ。 |
fit_predict() |
fit()とpredict()を同時に行う。戻り値はpredict()と同じ。 |
1. データの前処理
1-1. データの取得
データは以下のGitHubへのリンクからCSVファイルをダウンロードし、同一フォルダ内に配置する。ファイル名にスペースが含まれている場合は、_に変更する。(Wholesale_customers_data.xlsx)
https://github.com/VaderSame/Wholesale-Customers-Dataset/blob/main/Wholesale%20customers%20data.xlsx
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.preprocessing import StandardScaler # データの標準化を行うクラス
from sklearn.cluster import KMeans # k-means法を実行するクラス
%matplotlib inline
df = pd.read_excel("Wholesale_customers_data.xlsx")
df.head()
各列の内容は以下の通り。
| 列名 | 内容 |
|---|---|
| Channel | 業態 (1:サービス業、2:小売業) |
| Region | 地域 (1:リスボン市、2:ポルト市、3:その他) |
| Fresh | 生鮮食品の年間購入額 |
| Milk | 乳製品の年間購入額 |
| Grocery | 食料雑貨品の年間購入額 |
| Frozen | 冷凍食品の年間購入額 |
| Detergents_Paper | 洗剤と紙製品の年間購入額 |
| Delicassen | 惣菜の年間購入額 |
1-2. データの確認
基本統計量
df.describe()
欠損値の有無
df.isnull().sum().any()
# False
1-3. 不要なデータの削除
"Channel"カラムや"Region"カラムは、データの値の種類を数値で表しているだけで、数値の大小に意味を持たないので、分析の対象から外す。
df = df.drop(["Channel","Region"],axis=1)
df.head()
1-4. 特徴量の標準化
k-means法で距離を計算する際に、スケールの大きい特徴量が他の特徴量に影響を与えるので、すべての特徴量の平均値と標準偏差を揃える必要がある。
データの標準化を行うにはsklearn.preprocessingモジュールのStandardScalerクラスを利用する。
StandardScalerクラスのメソッド
| メソッド名 | 内容 |
|---|---|
fit() |
2次元データ(DataFrameなど)を受け取り、その各列について平均と標準偏差を計算し、その結果を内部に保持する |
transform() |
内部に保持された平均(μ)と標準偏差(σ)を使って、標準化(x - μ) / σを行い、その結果を2次元Numpy配列で返す |
fit_transform() |
fit()とtransform()をまとめて実行する |
sc = StandardScaler()
df_scaled = pd.DataFrame(
sc.fit_transform(df), # dfを標準化した2次元Numpy配列
columns=df.columns)
df_scaled.describe()
fit_transform()メソッドにdfを渡すとデータを標準化した2次元Numpy配列を戻すので、それを再びDataFrame化し、df_scaledとする。
2. 最適なクラスタ数の決定
最適なクラスタ数を決定するには、エルボー法を使って大体の範囲を絞り、その範囲の中ので最適な値を選択するのにシルエット分析を使う。
2-1. エルボー法
エルボー法は、k-meansの過程で計算するSSE(クラスタ内誤差二乗和:各データと所属クラスタの中心点との距離の二乗の合計)の値の変化を見て、最適なクラスタ数の値の範囲を絞る手法である。
SSE の値はクラスタ数を増やすほど小さくなるが、ある値を超えると SSE の減少率が低下してくる。この時のクラスタ数 k が最適なクラスタ数に近い値であると考える手法である。
以下はエルボー法を実装するコードである。
SSEは、学習済みのモデルのinertia_属性から参照できる。
# SSEの値を格納するリストを作成
sse_list =[]
# クラスタ数を1ずつ増やしながらモデルの作成と学習を行う
# 各モデルのSSEの値をリストに追加してく
for n in range(2,31):
model = KMeans(n_clusters=n,random_state=0,n_init=10)
model.fit(df_scaled)
sse_list.append(model.inertia_)
# 横軸：クラスタ数, 縦軸：SSEの線グラフをプロット
# 値の減少幅が低下している時のクラスタ数を見つける
plt.plot(np.arange(2,31),np.array(sse_list))
plt.xlabel("cluster")
plt.ylabel("SSE")
plt.show()
実行結果：
SSEの値の減少幅が低下しているのは、クラスタ数が 5 の辺りであるとわかる。
2-2. シルエット分析
シルエット分析とは、$i$番目のデータについて、次の2つの値を用いてクラスタリングの評価を行うことである。
- $a_{i}$ : 同じクラスタ内にあるすべてのデータとの平均距離
- $b_{i}$ : 隣接するクラスタ内にあるすべてのデータとの平均距離
$a_{i}$は、クラスタ内のデータの密集度を表し、$b_{i}$はクラスタ間の乖離度を表す。
$a_{i}$,$b_{i}$について、次のことが言える。
- $a_{i}$が低いほど クラスタ内は密集度が高い
- $b_{i}$が高いほど クラスタ間が乖離度が高い
よって
- $a_{i}$が低く かつ $b_{i}$が高いほど良いモデルである
この$a_{i}$と$b_{i}$の関係性を表す指標がシルエット係数 $s_{i}$ で、$b_{i}$と$a_{i}$の差を、$b_{i}$と$a_{i}$のどちらか大きい方で割った値である。
s_{i}=\frac {b_{i}-a_{i}}{max(a_{i},b_{i})}
$s_{i}$は $-1 <s_{i}< 1$ の範囲の値をとり、$s_{i}$の値とモデルの評価との関係性は以下。
| $s_{i}$ | モデルの評価 |
|---|---|
| 1に近い | 良い |
| -1に近い | 悪い |
シルエット係数を求めるには、sklearn.metricsモジュールのsilhouette_score()関数を使用する。
silhouette_score(X,labels)の引数
-
X: 2次元の説明データ -
labels: 各データの属するクラスタ番号を要素とする1次元Numpy配列(predict()の戻り値)
エルボー法でクラスタ数は5に近い値であるということがわかっているので、その周辺の全ての値でシルエット分析を行う。
from sklearn.metrics import silhouette_score
for n in range(2,8):
model = KMeans(n_clusters=n,random_state=0,n_init=10)
labels = model.fit_predict(df_scaled)
score = silhouette_score(df_scaled,labels)
print(f"クラスタ数:{n}, シルエットスコア:{score}")
実行結果：
クラスタ数:2, シルエットスコア:0.5715166500367507
クラスタ数:3, シルエットスコア:0.5620284080698716
クラスタ数:4, シルエットスコア:0.41041538163747193
クラスタ数:5, シルエットスコア:0.3715836236351362
クラスタ数:6, シルエットスコア:0.37203143575408226
クラスタ数:7, シルエットスコア:0.3177928222685055
クラスタ数2と3のスコアが僅差で高いことがわかった。
3. モデルの作成
クラスタ数3でクラスタリングを実行する。
model = KMeans(n_clusters=3, random_state=0, n_init=10)
model.fit(df_scaled)
model.labels_
array([1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 1, 0, 0, 0, 1, 0, 1, 1, 1, 1, 1,
1, 2, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 0, 1, 1, 1, 0, 0,
0, 0, 0, 2, 0, 0, 1, 1, 1, 0, 1, 1, 0, 0, 1, 1, 1, 2, 1, 0, 1, 2,
1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 0, 1, 1, 2, 2, 1,
1, 1, 1, 1, 2, 1, 0, 1, 1, 1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 1, 0,
・・・(省略 )
])
クラスタ数を3と指定したKMeansクラスのfit()メソッドを実行すると、属性labels_に分割後のクラスラベルが1次元Numpy配列で格納される。
4. クラスタリングの評価
クラスタリング結果からモデルの評価を行う。
4-1. シルエットスコア
前述のシルエット分析で使ったシルエットスコアで評価する。
silhouette_score(df_scaled, model.labels_)
# 0.3425521614309284
4-2. 次元削減 → 散布図
多次元のデータを主成分分析で2次元に変換し、散布図に描画する。
主成分分析にはsklearn.decompositionモジュールのPCAクラスを利用する。
コンストラクタ引数
-
n_components: 変換後の次元数
メソッド
-
fit(data): 主成分分析を実行する -
transform(data): 各データを主成分軸へ射影した座標を要素とする2次元Numpy配列を返す -
fit_transform(data):fit()とtransform()を同時に実行する
from sklearn.decomposition import PCA
# PCAの作成と実行
pca = PCA(n_components=2)
pca_data = pca.fit_transform(df_scaled)
pca_data
n_components=2と指定したので2つの新しい軸へ射影した座標を、データの数だけ並べた2次元Numpy配列を取得できる。
[[ 1.93290546e-01 -3.05099963e-01]
[ 4.34419896e-01 -3.28412623e-01]
[ 8.11143234e-01 8.15095701e-01]
[-7.78647832e-01 6.52753734e-01]
... (省略)
[ 3.46570362e+00 -1.03983801e+00]
[-9.18022733e-01 -3.00465906e-02]
[-1.10513709e+00 -8.61337874e-01]]
# 散布図の描画
plt.scatter(pca_data[:,0], pca_data[:,1], c=model.labels_)
plt.xlabel("PC1")
plt.ylabel("PC2")
plt.show()
4-3. グループ化 → 棒グラフ
クラスタごとにグループ化して平均を棒グラフに表示する。
# クラスタの列を追加
df_scaled["cluster"] = model.labels_
# クラスタの値でグループ化しその平均をプロット
cluster_mean = df_scaled.groupby("cluster").mean()
cluster_mean.plot(kind="bar")
plt.show()
凝集型クラスタリング
階層的クラスタリングは、データを階層構造として構築してまとめていく手法。次のような特徴が挙げられる。
- 樹形図(テンドログラム)でクラスタの結合過程を可視化できる
- クラスタ数を事前に指定する必要がない
階層的クラスタリングには凝集型と分裂型がある。
- 凝集型(ボトムアップ): 小さなクラスタをマージしていく
- 分裂型(トップダウン): 大きな一つのクラスタを分割していく
AgglomerativeClusteringクラス
凝集型クラスタリングを実行するには、sklearn.clusterのAgglomerativeClusteringクラスを利用する。
コンストラクタ引数
-
n_clusters(int): クラスタ数(Noneも指定できる) -
metric(str) : データ間の距離の定義-
"manhattann", "l1": マンハッタン距離 -
"euclidean","l2": ユークリッド距離 -
"cosine": コサイン距離
-
-
linkage(str) : クラスタ間の距離として使う基準(連結方法)-
"ward": クラスタ内の分散を最小化する
(デフォルト,metricの指定が"euclidean"で固定される) -
"complete": 2つのクラスタに属する各データ間の最大距離 -
"average": 2つのクラスタに属する各データ間の平均距離 -
"single": 2つのクラスタに属する各データ間の最小距離
-
-
compute_distances(bool): クラスタ間の距離計算をするかどうか -
distance_threshold(float) : この値(閾値)以下の距離はマージされない(指定した場合はn_clustersをNoneにする必要がある)
# 先ほど使用したデータセットをそのまま利用します
from sklearn.cluster import AgglomerativeClustering
# モデルの作成・実行
ac = AgglomerativeClustering(linkage="ward",
compute_distances=True)
df_scaled["cluster"] = ac.fit_predict(df_scaled)
属性
-
children_: マージされた2つのノードのインデックスの並びを要素とするnumpy配列。(n_samples - 1, 2 )の行列になる。 -
distance_: 全てのマージについて、マージした2つのノード間の距離を計算し、その値を格納する1次元numpy配列。要素数はn_samples - 1となる。(※AgglomerativeClusteringのコンストラクタでcompute_distances=Trueと指定しないとこの属性は作られない)
例えば、次のようなマージをしたとする。
- インデクス番号(0,1)のノードがマージされ、インデクス番号10のノードがクラスタとして作られる
- (2,3),(4,5),(6,7),(8,9)のノードも同様にマージされ、11,12,13,14のノードが作られる
- (10,11)のノードがマージされ、15のノードが作られる
- (12,13)のノードがマージされ、16のノードが作られる
- (16,14)のノードがマージされ、17のノードが作られる
- (15,17)のノードがマージされ、18のノードが作られ、終了する
children_は、このマージされる2つのノードを全て格納した2次元numpy配列で、distances_は、マージされる時の2つのノードの距離を全て格納した1次元numpy配列である。
| merge_index | children_ |
新しいノード | distances_ |
|---|---|---|---|
| 0 | [ 0, 1 ] | 10 | ノード0とノード1の距離 |
| 1 | [ 2, 3 ] | 11 | ノード2とノード3の距離 |
| 2 | [ 4, 5 ] | 12 | ノード4とノード5の距離 |
| 3 | [ 6, 7 ] | 13 | ノード6とノード7の距離 |
| 4 | [ 8, 9 ] | 14 | ノード8とノード9の距離 |
| 5 | [10, 11 ] | 15 | ノード10とノード11の距離 |
| 6 | [12, 13 ] | 16 | ノード12とノード13の距離 |
| 7 | [16, 14 ] | 17 | ノード16とノード14の距離 |
| 8 | [15, 17 ] | 18 | ノード15とノード17の距離 |
ここからわかるのは、
- マージの合計回数
-
children_の要素数 -
distances_の要素数 - n_samples - 1
がすべて同じ値になるということである。
dendrogram()関数
デンドログラムの描画には、scipy.cluster.hierarchyモジュールのdendrogram関数を利用する。
引数
-
第1引数: 連結行列(手作りする必要がある) -
truncate_mode(str) : デンドログラムの表示を簡略化する-
"level":p引数で指定した階層数を表示する -
"lastp":p引数で指定したクラスター数を上位から表示する
-
-
p(int) :truncate_modeで指定する値 -
color_threshold(float) : 色分けの閾値(結合距離の値)
連結行列
連結行列(linkage_matrix)とは、各行が1回のマージについての情報を4列で保持する、( n_samples - 1, 4 )の形状の行列である。
連結行列の各列は次の4つの情報を持つ
- 第0列 : マージされた左ノード =
children_[i][0] - 第1列 : マージされた右ノード =
children_[i][1] - 第2列 : マージされた2つのノードの距離 =
distances_[i] - 第3列 : そのマージでできた新しいノードが包含するノード数
第3列はマージした左ノードと右ノードがそれぞれ持っているノード数の合計を計算するが、マージするノードの値がリーフノードか内部ノードかによって挙動が変わる。
子ノードのインデックス = node_index とすると、
- リーフノード(node_index < n_samples) :
1 + 1 - 内部ノード(node_index >= n_samples) :
linkage_matrix[node_index - n_samples, 3]
マージするノードが内部ノードの場合は、連結行列の第3列の値を参照するが、その行番号は「マージするノードのインデックスから、リーフノードの個数を引いた値」と同じになる。
from scipy.cluster.hierarchy import dendrogram
children = ac.children_
distances = ac.distances_
n_samples = df_scaled.shape[0]
# 連結行列を0で初期化
linkage_matrix = np.zeros((n_samples - 1, 4))
# 各ノードが包含するサンプル数のリストを、リーフノード分だけ1で初期化
counts = np.ones(n_samples)
# 連結行列の作成
for i, (index_l, index_r) in enumerate(children):
count_l = counts[index_l] if index_l < n_samples else linkage_matrix[index_l - n_samples, 3]
count_r = counts[index_r] if index_r < n_samples else linkage_matrix[index_r - n_samples, 3]
linkage_matrix[i, 0] = index_l
linkage_matrix[i, 1] = index_r
linkage_matrix[i, 2] = distances[i]
linkage_matrix[i, 3] = count_l + count_r
counts = np.append(counts, count_l + count_r)
# デンドログラムの描画
fig, ax = plt.subplots(figsize=(15, 6))
dendrogram(linkage_matrix)
ax.set_xlabel("Node Index")
ax.set_ylabel("Distance")
plt.show()
truncate_modeに"level"を渡して、pで階層の数を指定すると、その値より深い階層のノードが省略され、すっきりしたデンドログラムになる。
dendrogram(linkage_matrix, truncate_mode="level", p=5)
color_threshold=20と指定すると、クラスタ間の距離が20以下で結合したノードは同じクラスタとして同じ色で表示する。
dendrogram(linkage_matrix, color_threshold=20)
