ランダムフォレストはクラスタリングに用いられる

## 「ランダムフォレストはクラスタリングに用いられる」

と、[wikipedia](https://ja.wikipedia.org/wiki/%E3%83%A9%E3%83%B3%E3%83%80%E3%83%A0%E3%83%95%E3%82%A9%E3%83%AC%E3%82%B9%E3%83%88)を初めいろんなところに書いてあるが具体的な方法はあまり見かけない.  

### 具体的な方法が書かれている文献

T. Shi, S. Horvath, "Unsupervised Learning With Random Forest Predictors", 2006, 
https://horvath.genetics.ucla.edu/html/RFclustering/RFclustering/RandomForestHorvath.pdf.  

L. Breiman, A. Cutler, "Manual--Setting Up, Using, And Understanding Random Forests V4.0", 2003,
 https://www.stat.berkeley.edu/~breiman/Using_random_forests_v4.0.pdf.  

[Shi-Horvath, 2006] によると手法が提案されたのは [Breiman-Cutler, 2003].  
ただしそちらは Fortran プログラムのマニュアルとして書かれていて, アルゴリズムの理解には前者の方が読みやすい.  

これ以外のランダムフォレストでクラスタリングする具体的な方法は見つからなかった[^clustering].  

[^clustering]: ノードの分割で直接クラスタを生成するような方法を想像していたので「wikipedia その他が言及してるのって本当にこのアルゴリズム ?」という思いがある.  

### アルゴリズム

クラスタリングしたいデータ $X$ に対して以下を行う.  

1. 異なる変数同士が独立で変数毎には $X$ と同分布な人工データ $\widetilde{X}$ を生成,  
2. $X$ と $\widetilde{X}$ をランダムフォレストで分類,  
3. 2 で作ったモデルを元に, $X_i$ と $X_j$ が同じノードに入る割合 $P_{ij}$ で類似度行列 $P$ を作成,  
5. $D = \sqrt{1-P}$ を使って距離行列ベースのクラスタリング手法でクラスタリング.  

$\widetilde{X}$ の生成には,  

* 変数毎に $X$ から無作為抽出する,  
* 変数毎に $X$ の分布を適当に近似した乱数を取る,  

といった方法を利用できる.  
  
ランダムフォレストが利くのは 3 の類似度行列を生成する部分で, 同じノードに入る割合を見るので決定木のアンサンブル系のモデル[^ensemble]なら置き換えられる.  
それ以外だと類似度を測る別の指標が必要になる.  

[^ensemble]: ブースティングの場合は個々の木の役割が違うので「$X_i$ と $X_j$ が同じノードに入る割合」が「$X_i$ と $X_j$ の類似度」を表さないかもしれないが, アルゴリズムの適用はできるので何らかの意味でのクラスタリングは行われる.  

### 性質

距離行列がサイズ [レコード数]^2 の正方行列になるのであまりレコードの多いデータには適用できない.  
反面, 高次元データには強くてレコード数より説明変数次元が大きければ[レコード数]^2 までデータを圧縮できる利点がある (DNA データ等).  
ユークリッドでない距離に基いてクラスタリングが行われるため単純な方法と異なる結果が得られることも期待されるらしい.  

### 実験

ランダムフォレストを使ったクラスタリングの動作確認と他のクラスタリング (k-means) との比較の実験.  

```python
import numpy as np
import matplotlib.pyplot as plt
from sklearn.ensemble import RandomForestClassifier
from sklearn.cluster import AgglomerativeClustering

# 比較用
from sklearn.cluster import KMeans

rng = np.random.default_rng()
```

#### データ生成

クラスタリングするデータを生成する.  
混合ガウス分布のパラメータをランダムで設定してそれに従うデータを同じくランダムに生成する.  
標準正規分布に従う乱数データを平均・標準偏差のパラメータ配列で変換すると速い.  
20 次元データのプロットは無理なので最初の 2 次元を抽出してプロットしてみる.  

```python
# 混合ガウス分布のパラメータ設定
m = 10
d = 20
a = rng.uniform(0.5, 1.5, size=m) # 混合比
a = a/np.sum(a)
mu = rng.normal(size=[m, d]) # 平均
mu = mu/np.linalg.norm(mu, axis=1, keepdims=True)
t = 0.05*rng.normal(size=[m, d, d]) # "共標準偏差" (t@t.T が共分散になる)

# データ生成
n = 10**3
c = rng.choice(m, size=n, p=a)
X_org = mu[c]+(np.random.normal(size=[n, 1, d])@t[c])[:, 0]

# 2 次元に射影して混合ガウス分布のでのクラスで色分けしてプロット
plt.scatter(X_org[:, 0], X_org[:, 1], s=5, c=c, cmap=plt.cm.jet)
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2988004/eefed61c-f8c3-7e06-be3d-685947d6a95e.png)

2 次元に制限しているので分離不可能に見えるが実際は 20 次元なので見た目よりはクラス毎に固まっているはず.  

#### ランダムフォレストモデルの作成

各変数について全データからランダムにサンプリングする方法で人工データを生成してランダムフォレストで 2 値分類する.  
汎化性能はあまり気にしなくていいと思われるが, ノードを分割しすぎると類似度行列がほとんど 0 になってしまう恐れがあるのでパラメータは適度に調整しておく.  

```python
# 人工データ生成
ii = rng.choice(n, size=[n, d])
X_gen = X_org[ii, np.arange(d)]

# RandomForest で2値分類
# オリジナルデータならラベル 1
X = np.vstack([X_org, X_gen])
y = np.hstack([np.ones(n), np.zeros(n)])
rfc = RandomForestClassifier(min_samples_split=0.05)
rfc.fit(X, y)
print(rfc.score(X, y))
```

    0.958

train スコア (accuracy) は 0.958 で高め.  

#### クラスタリング

ランダムフォレストモデルから距離行列を作成して階層型クラスタリングの `sklearn.cluster.AgglomerativeClustering` で距離行列を基にクラスタリングする.  
距離行列に基いた階層型クラスタリングの詳細はリファレンスか解説記事を参照.  

```python
# オリジナルデータの類似度行列作成
# T: X_org の各行データが RandomForest の各木でどの葉に落ちるか表す行列, shape=[n, n_estimator]
# P: 類似度行列, shape=[n, n]
T = rfc.apply(X_org)
P = np.mean(T[:, np.newaxis]==T, axis=-1)

# 類似度行列から距離行列を作ってクラスタリング
# distance matrix に基づくクラスタリングに対応した関数が必要
# AgglomerativeClustering だと affinity='precomputed' が必要で linkage は 'ward' 以外なら何でも使えるらしい
k = m # クラスタ数
D = np.sqrt(1-P)
clst = AgglomerativeClustering(k, affinity='precomputed', linkage='average')
clst.fit(D)

# 同じく射影してクラスタリング結果で色分けしてプロット
plt.scatter(X_org[:, 0], X_org[:, 1], s=5, c=clst.labels_, cmap=plt.cm.jet)
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2988004/9bc52a95-77c6-505c-5206-a745db983878.png)


なるほど, 分からん.  

#### 混合行列による評価[^evaluation]

[^evaluation]: 1 次元の指標で評価することも考えたがちょうどいいものが見つからなかったので保留. 今回は「真のクラス」も存在するので通常のクラスタリング用の指標だと不十分になる.  

20 次元データのプロットでは評価しにくいので多クラス混合行列を作って眺めてみる.  
クラス分類ではないので対角成分が大きければいいわけではなく, 各行各列で突出して大きい成分があれば良い結果と見なされることに注意.  

```python
# 真の (混合ガウス分布としての) クラスとクラスタリングによるクラスを one-hot 化してそれぞれの組み合わせになるレコードをカウント
oh_org = c[:, np.newaxis]==np.arange(m)
oh_clst = clst.labels_[:, np.newaxis]==np.arange(m)
Z = np.sum((oh_org[..., np.newaxis])*(oh_clst[:, np.newaxis]), axis=0)
print(Z)
```

    [[  0   0   0   0   0   1 112   0   0   0]
     [ 64   2   0   1   0   0   0   0   0   0]
     [  0   0   0   0   2   1   0   0  79   0]
     [  0   0   1   0   0  95   0   0   0   0]
     [  0   0 117   0   0   0   0   0   1   0]
     [  1   1   0   1   0   0   0   0   0 110]
     [  0  75   0   0   0   0   0   0   0   0]
     [  0   0   0   0  86   0   0   0   0   0]
     [  0   7   0   0   0   0   0 103   0   1]
     [  1   2   2 131   0   0   2   0   0   1]]

各行各列で大きい成分がちょうど 1 つずつあるので良い感じじゃないかな ?  

#### 他のクラスタリングとの比較

k-means `sklearn.cluster.KMeans`  と比較してみる.  

```python
# k-means と比較
clst_km = KMeans(k)
clst_km.fit(X_org)

# 同様にプロット
plt.figure(figsize=[s/1.25 for s in plt.rcParams['figure.figsize']], dpi=plt.rcParams['figure.dpi']/1.25)
plt.scatter(X_org[:, 0], X_org[:, 1], s=5, c=clst_km.labels_, cmap=plt.cm.jet)
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2988004/48897a6d-ced0-b613-24f4-2f9effb92fdc.png)


```python
# 同様に混合行列を計算
oh_org = c[:, np.newaxis]==np.arange(m)
oh_clst_km = clst_km.labels_[:, np.newaxis]==np.arange(m)
Z_km = np.sum((oh_org[..., np.newaxis])*(oh_clst_km[:, np.newaxis]), axis=0)
print(Z_km)
```

    [[113   0   0   0   0   0   0   0   0   0]
     [  0   0  61   0   1   1   0   1   3   0]
     [  0  81   1   0   0   0   0   0   0   0]
     [  0   0   0   1   0   0   0   0   0  95]
     [  0   0   1 117   0   0   0   0   0   0]
     [  0   0   0   0   0   1   1 111   0   0]
     [  0   0   0   0   0   1   0   0  74   0]
     [  0   0   1   0  84   0   0   1   0   0]
     [  1   0   0   0   1   0 109   0   0   0]
     [  0   0   3   0   0 135   0   1   0   0]]

クラスタリングするのに簡単すぎるデータだったかもしれないがランダムフォレストでクラスタリングできることの確認はできたといえそうである.