第11回 k-means②：スケール対応とクラスタ数の選択¶
本日の到達目標¶
- 初期値依存性の問題と k-means++ による対処法を説明できる
- ミニバッチ k-means のアルゴリズムと通常 k-means との精度/速度トレードオフを説明できる
- シルエット分析でクラスタ数を選択する手順を説明できる
- 発展：WCSS が \(k\) について単調に減ることを示し、だからこそ WCSS 最小化では \(k\) を決められない理由を説明できる
1. 初期値依存性の問題¶
1.1 ランダム初期化の落とし穴¶
Lloyd's アルゴリズム（第10回）は初期重心の選び方によって最終結果が異なる。
最悪ケースの例：
局所最適解に陥ることがある → 大域最適は保証されない
1.2 k-means++ による初期化¶
アルゴリズム（Arthur & Vassilvitskii, 2007）：
- 一様ランダムに1点を選んで最初の重心 \(\mu_1\) とする
- 各データ点 \(x_i\) について、最も近い既存重心との距離 \(D(x_i)\) を計算
- 確率 \(\propto D(x_i)^2\) でデータ点を選び、次の重心にする
- \(k\) 個の重心が決まるまで繰り返す
直観：既存の重心から「遠いデータ点」ほど次の重心に選ばれやすい
→ 初期重心が互いに離れた良い配置になる
理論的保証：k-means++ の期待WCSSは最適解の \(O(\log k)\) 倍以内
導出①：なぜ \(D\) ではなく \(D^2\) に比例させるのか
一様ランダムに初期重心を選ぶと何が起きるかから考えます。点が密集した塊が1つ、離れた小さな塊が1つあるとき、一様に選べば大きい塊から2つとも選ばれる確率が高い。すると密集塊に重心が2つ入り、小さい塊には1つも入りません。第10回の導出②より反復は改善方向にしか動かないので、この配置から抜け出せません。
k-means++ は、既存の重心からの最近距離 \(D(\mathbf{x})\) を使い、確率
で選びます。ここで指数が2であることには理由があります。目的関数 WCSS は
でした。いま重心が決まっている途中段階では、点 \(\mathbf{x}\) が WCSS に与えている寄与はちょうど \(D(\mathbf{x})^2\) です。つまり分子の \(D(\mathbf{x})^2\) は、その点が現在いくら「損」を出しているかそのものです。
いま最も損を出している点ほど、次の重心に選ばれやすい。目的関数と同じ尺度で重み付けしているので、\(D\) の1乗でも3乗でもなく2乗なのです。第10回で「WCSS は二乗和である」と決めたことが、ここに効いています。
この授業で扱う範囲外：期待 WCSS が最適解の \(O(\log k)\) 倍以内に収まるという保証の証明。結果だけを使います。
使った数学：確率の重み付け、目的関数との対応。
【カードに記入・3分】：k-means++ の確率 \(\propto D^2\) で選ぶことで、どのような状況を回避できますか？
2. ミニバッチ k-means¶
2.1 標準 k-means の計算量の問題¶
\(n = 10^7\)、\(d = 10\)（PCA後の次元）、\(k = 5\) の場合： - 1反復 = \(O(nkd) = O(5 \times 10^8)\) の距離計算 - 収束まで50反復とすると \(2.5 \times 10^{10}\) の演算 → 数分〜数十分
2.2 ミニバッチ k-means のアルゴリズム¶
基本アイデア：全データを毎回使わず、ランダムにサンプリングした ミニバッチ（例：1000件）で更新する。
アルゴリズム： 1. ミニバッチ \(B \subset \{1, \ldots, n\}\) をランダムサンプリング（\(|B| \ll n\)） 2. \(B\) の各点を最近傍重心に割り当て 3. 重心を「指数移動平均的」に更新： $\(\boldsymbol{\mu}_j \leftarrow \boldsymbol{\mu}_j + \eta (x_i - \boldsymbol{\mu}_j)\)$ （\(x_i\) はクラスタ \(j\) に割り当てられた \(B\) 内の点） 4. 収束まで繰り返す
計算量：\(O(T \cdot b \cdot k \cdot d)\) — \(b = |B|\) がミニバッチサイズ
2.3 通常 k-means vs ミニバッチ k-means¶
| 観点 | 通常 k-means | ミニバッチ k-means |
|---|---|---|
| 1反復のコスト | \(O(nkd)\) | \(O(bkd)\)（\(b \ll n\)） |
| WCSS精度 | 高い | やや高い（通常数%の劣化） |
| 収束速度 | 遅い | 速い |
| 適用場面 | \(n\) が中規模 | \(n\) が大規模（\(10^6\)以上） |
列指向ストレージとの接続：ミニバッチでは毎回全行ではなく一部の行を読む。行指向だと各行全列を読む必要があるが、列指向ならミニバッチで必要な列だけを部分読みできる（行選択と列選択の組み合わせ）。
3. クラスタ数 k の選択¶
3.1 エルボー法（Elbow Method）¶
手順：\(k = 1, 2, \ldots, K_{\max}\) に対して WCSS を計算し、グラフを描く
問題点：実際のデータでは「肘」が曖昧なことが多い → より客観的な指標が必要
導出②：WCSS 最小化では \(k\) を決められない（\(k\) について単調に減るから）
エルボー法が「最小値」ではなく「折れ曲がり」を探すのには理由があります。
\(k\) 個に分けた最適な分割の WCSS を \(W(k)\) とします。\(k+1\) 個に分けるとき、\(k\) 個のときの分割のどれか1つを2つに割るという分け方が必ず作れます。点を2つの側に分けたそれぞれの平方和の合計は、第10回の導出①より、1つの重心でまとめたときの平方和以下です（重心が2つに増えて、それぞれが自分の側の平均になるため）。
つまり \(k+1\) 個での達成可能な値の中に \(W(k)\) 以下のものが必ずあり、\(W(k+1)\) はその最小値なので
WCSS は \(k\) について単調非増加です。極端な場合、\(k = n\) では各点が自分だけのクラスタになり重心が自分自身と一致するので
したがって「WCSS を最小にする \(k\)」を選ぶと、必ず \(k = n\) になります。 これはクラスタリングとして無意味です。WCSS はそもそも \(k\) を選ぶ基準になりえません。
だから見るのは値そのものではなく減り方の変化（肘）になり、より客観的にしたいならシルエット係数のように「密集度と分離度の釣り合い」を測る別の量が要ります。シルエット係数が \(k\) について単調でないのは、\(b(i)\)（隣接クラスタとの距離）が \(k\) を増やすと小さくなる方向に働くからです。
この「選べる範囲が広がれば最小値は下がらない」という論法は、第6回の \(R^2\)、第13回の PLS 成分数でも同じ形で出てきます。
使った数学：最小値の定義（達成可能な値の集合が広がれば最小値は下がらない）。
3.2 シルエット分析（Silhouette Analysis）¶
シルエット係数（各データ点 \(i\) に対して）：
- \(a(i)\)：\(i\) と同じクラスタ内の他の点との平均距離（小さいほど密）
- \(b(i)\)：\(i\) に最も近い隣接クラスタ内の点との平均距離（大きいほど分離）
| \(s(i)\) の値 | 解釈 |
|---|---|
| ≈ +1 | クラスタ内に密、隣接クラスタと分離が良い（理想） |
| ≈ 0 | 2つのクラスタの境界付近にある |
| ≈ −1 | 誤ったクラスタに割り当てられている可能性 |
導出③：シルエット係数が必ず \(-1\) と \(+1\) の間に入る理由
表の3行は暗記する必要がありません。定義から場合分けするだけで出ます。\(a(i) \geq 0\)、\(b(i) \geq 0\) です。
\(b(i) \geq a(i)\) のとき（隣のクラスタの方が遠い＝良い状態）、分母は \(b(i)\) なので
\(0 \leq a(i)/b(i) \leq 1\) なので \(0 \leq s(i) \leq 1\)。\(a(i) = 0\)（自分のクラスタ内で完全に一致）のとき \(s(i) = 1\) です。
\(a(i) > b(i)\) のとき（隣のクラスタの方が近い＝おかしい状態）、分母は \(a(i)\) なので
\(b(i)/a(i)\) は \(0\) 以上 \(1\) 未満なので \(-1 \leq s(i) < 0\)。
合わせて \(-1 \leq s(i) \leq 1\)、\(a(i) = b(i)\) のとき \(s(i) = 0\) です。
分母が \(\max(a,b)\) である理由がここで分かります。\(a\) でも \(b\) でもなく大きい方で割ることで、どちらの場合でも比が1を超えず、目盛りが \([-1, 1]\) に収まるのです。大きさの異なるクラスタ間でも同じ物差しで比べられるように作られています。
使った数学：場合分けと不等式。
平均シルエット係数（全データ点の平均）：
→ \(k\) を変えて \(\bar{s}\) が最大になる \(k\) を選ぶ
3.3 シルエットプロット¶
各クラスタ内の点のシルエット係数を並べて可視化：
クラスタ1 ████████████░░░░░░░░ ← 係数が揃って高い → 良い分離
クラスタ2 ██████████████████ ← 係数が高い
クラスタ3 ████░░░░░░░░░░ ← 係数にばらつき → このkは不適切かも
─────────────── s(i)
-1 0 +1
3.4 その他の指標（参考）¶
| 指標 | 概要 |
|---|---|
| Gap統計量 | ランダムデータの WCSS との比較 |
| Calinski-Harabasz指数 | クラスタ間分散/クラスタ内分散の比 |
| Davies-Bouldin指数 | クラスタの「似ていなさ」の指標（小さいほど良い） |
4. 再現実験のための乱数シード管理¶
k-means はランダム性（初期化・ミニバッチのサンプリング）を含む → 再現性の確保が必要
管理すべき要素：
- 乱数シード（random_state = 42 等）
- データのスナップショット版（第3回）
- k-means++ の初期化方法（同一）
DAGへの組み込み：
実験設定ファイル:
random_seed: 42
init: k-means++
data: snapshot_20240101.parquet
実験A (k=3): → WCSS_A, シルエット_A
実験B (k=4): → WCSS_B, シルエット_B
実験C (k=5): → WCSS_C, シルエット_C
↓
比較表 + k 選択
→ 設定ファイルを保存することで同じ実験が誰でも再現できる
5. 基盤リンクの整理¶
| 基盤概念 | 本回での役割 | 具体的な使い方 |
|---|---|---|
| 列指向 | ミニバッチで必要列だけ読む | 大規模時のI/Oをさらに削減 |
| 整合性 | 乱数シード・データ版を固定 | k の比較実験の再現性 |
| DAG | 反復実験（k選択）をワークフローとして管理 | 各 k の実験を並列化できる |
| 近似 | ミニバッチやサンプリングで近似最適化 | 厳密な最適解の代わりに高速化 |
6. 本日のまとめ¶
- k-means++：\(D^2\) 比例の確率で初期重心を選ぶことで局所最適を回避
- ミニバッチ k-means：サブサンプルで更新することで大規模データに対応（精度は若干劣化）
- シルエット分析：\(s(i) = (b(i) - a(i)) / \max(a(i), b(i))\) で各点の「クラスタへの帰属度」を評価
- 再現実験には乱数シードとデータスナップショットの固定が必須
次回予告¶
第12回では PLS（部分最小二乗法）①：共分散最大化により予測と次元削減を同時に行う原理を学ぶ。
共線性がある説明変数群（時間窓特徴量等）に対して重回帰が不安定になる状況への代替手法。
事後課題¶
以下の問いに答えること（A4 半枚）：
- k-means のエルボー法では k を客観的に決めることが難しい場合がある。シルエット分析がより客観的な理由を、それぞれの数式を参照しながら説明してください。
- ミニバッチ k-means の「重心を指数移動平均的に更新する」方式が、標準的な「クラスタ内平均で更新する」方式と収束先が近くなる理由を直観的に説明してください。