sphereclusterで、球面クラスタリングを扱えるようになる

## sphereclusterに至った経緯

コサイン距離でクラスタリングしようと思ったら、球面クラスタリングに辿り着いた。

大きな次元のノイジーなデータを扱う際に、距離尺度としてコサイン距離の方のクラスタリングを行いたかった。階層的クラスタリングの距離尺度をコサイン距離にするのが手軽ではあるが、非階層的アルゴリズムに比べると計算量が多く、対象のデータ量が数十万以上になると実用的ではない。

ということで、非階層でコサイン距離を使ったクラスタリングを探していたら、球面クラスタリングの実装である[spherecluster](https://github.com/clara-labs/spherecluster)を発見した。
これは、距離としてコサイン距離を用いているわけではないが、入力ベクトルの方向を使用してクラスタを作成しているため、得たい結果は同様のものである。
という訳で、これに目を通した&使ったので、その中での色々が以下に続く。

ちなみに、[scikit-learnのk-meansで独自の距離関数を定義する](https://qiita.com/Kensuke-Mitsuzawa/items/67c17f6bda8a1097740c)のような記事もあったが、確認時点ではKMeansクラスの仕様が変わっていたため、こちらの手段は利用できなかった。

## spherical k-means clustering (skmeans)

k-meansを超球上のデータに拡張したもの。

### 概要

次の手順でクラスタリングする。

1. 入力ベクトルデータを、それぞれ正規化する(ノルムを１にする)
2. 初期値決定法(k-means++など)で、k個の重心ベクトルを決定する
3. 重心ベクトルを正規化する
4. 入力ベクトルそれぞれに対し一番近い重心ベクトルを求め、近い重心ベクトルごとに対応するクラスタに割り当てる
5. クラスタごとに、所属するベクトルから重心ベクトルを計算する
6. 3-5を収束条件を満たすまで繰り返す

k-meansと違い、1.と3.の手順が加わる
この手順により、計算で使用するベクトルは常に単位方向ベクトルとなる。

初期化法については、理論的なところは深く追っていないため、あまり良くわからない。
k-means++自体は、入力データの分布確率に沿った初期値の置き方をする。そのため、これを使っておけば感覚的には問題にはならないと思っている。

### k-meansとの比較

一様分布から生成したデータに対して、どのように違いが出るかを比較する。

```python:k-meansとskmeansの比較
from sklearn.cluster import KMeans
from spherecluster import SphericalKMeans
import numpy as np
import matplotlib.pyplot as plt

np.random.seed(seed=529)
X = np.random.uniform(low=-1, high=1, size=(300, 2))

km = KMeans(n_clusters=5, init="k-means++").fit(X)
skm = SphericalKMeans(n_clusters=5, init="k-means++").fit(X)

colors = ["r", "g", "y", "m", "b"]
ax = plt.subplot(1, 2, 1, xlim=[-1.1, 1.1], ylim=[-1.1, 1.1])
for color, label in zip(colors, np.unique(km.labels_)):
    polt_data = X[km.labels_ == label]
    plt.plot(plot_data[:, 0], plot_data[:, 1], "%s+" % color)
ax.set_title("k-means")

ax = plt.subplot(1, 2, 2, xlim=[-1.1, 1.1], ylim=[-1.1, 1.1])
for color, label in zip(colors, np.unique(skm.labels_)):
    polt_data = X[skm.labels_ == label]
    plt.plot(plot_data[, 0], plot_data[:, 1], "%s+" % color)
ax.set_title("skmeans")

plt.show()
```

### 結果

少し分かりにくいが、k-meansでは紫クラスタと赤、あるいは緑クラスタは、同方向でも別々のクラスタになっているのに対して、skmeansでは綺麗に原点からの方向ごとにクラスタが出来ている。

<img width="527" alt="スクリーンショット 2017-12-07 10.36.32.png" src="https://qiita-image-store.s3.amazonaws.com/0/221698/b6f4a6ba-914f-6e70-6cc3-12c64b342154.png">


## Mixture of von Mises Fisher clustering model (movMF)

混合ガウスモデルを仮定したクラスタリングの超球への拡張。

普通のガウス分布は、平均と分散で記述される。
それに対して、単位円周上の確率を表すフォンミーゼス分布は平均方向と集中度で記述される。
平均方向は平均、集中度は分散に対応するようなパラメータである。(ただし、集中度は大きいほど、ばらつきが小さくなる)

なお、フォンミーゼス分布を多次元に拡張すると、フォンミーゼスフィッシャー分布と呼ぶらしい。

[von Mises-Fisher分布 - 機械学習の「朱鷺の杜Wiki」](http://ibisforest.org/index.php?von%20Mises-Fisher%E5%88%86%E5%B8%83)
[フォンミーゼスフィッシャー分布 | 高校数学の美しい物語](https://mathtrain.jp/vonmises)

### 概要

混合ガウスモデルと同様の考え方で、複数のフォンモーゼスフィッシャー分布から入力データを再現できると仮定して、分布のパラメータをEMアルゴリズムで推定して、クラスタを作成しているようだ。

入力データを再現するような確率分布を推定し、データに対してそれぞれの確率分布(=クラスタ)に所属する確率を計算するため、ハードクラスタリングだけでなくソフトクラスタリングも可能。

### skmeansとの比較

3つの正規分布に従うデータ80個ずつを混ぜて入力とし、クラスタの傾向を比べてみる。
入力データによっては、正規化できずにエラーが起きることがある。

```python:skmeansとvmfの比較
from sklearn.cluster import Means
from spherecluster import SphericalKMeans
from spherecluster import VonMisesFisherMixture
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D


np.random.seed(seed=531)
X = np.concaatenate((
    (np.random.normal(loc=0.5, scale=0.2, size=(80, 3)),
     np.random.normal(loc=-0.3, scale=0.3, size=(80, 3)),
     np.random.normal(loc=0, scale=0.1, size=(80, 3)) + np.array([0.5, -0.5, 0] * 80).reshape(80, 3),),
    axis=0
)
skm = SphericalKMeans(n_clusters=3, init="k-means++").fit(X)
vmf = VonMisesFisherMixture(n_clusters=3, posterior_type="hard").fit(X)

colors = ["r", "g", "b"]
fig = plt.figure(figsize=(10, 5))
ax = fig.add_subplot(1, 3, 1, aspect="equal", projection="3d", xlim=[-1.1, 1.1], ylim=[-1.1, 1.1], zlim=[-1.1, 1.1])
ax.plot(X[:, 0], X[:, 1], X[:, 2], "k+")

ax = fig.add_subplot(1, 3, 2, aspect="equal", projection="3d", xlim=[-1.1, 1.1], ylim=[-1.1, 1.1], zlim=[-1.1, 1.1])
for color, label in zip(colors, np.unique(skm.labels_)):
    plot_data = X[skm.labels_ == label]
    ax.plot(plot_data[:, 0], plot_data[:, 1], plot_data[:, 2], "%s+" % color)
ax.set_title("skmeans")

ax = fig.add_subplot(1, 3, 3, aspect="equal", projection="3d", xlim=[-1.1, 1.1], ylim=[-1.1, 1.1], zlim=[-1.1, 1.1])
for color, label in zip(colors, np.unique(vmf.labels_)):
    plot_data = X[vmf.labels_ == label]
    ax.plot(plot_data[:, 0], plot_data[:, 1], plot_data[:, 2], "%s+" % color)
ax.set_title("vmf")

plt.show()

```


### 結果

パッと見で違いは分かりにくいが、それぞれ元の分布を分割するようにクラスタリングできている。

<img width="584" alt="スクリーンショット 2017-12-07 13.01.25.png" src="https://qiita-image-store.s3.amazonaws.com/0/221698/54a1cb0f-62d1-a46f-77a1-9be50ed445d2.png">

ちなみに、混同行列チックに比較すると以下のようになる。
入力データは3つの分布から生成した80個ずつで構成されているので、それらを綺麗に分けているvmfの方がクラスタリング手法としては優秀だと分かる。

||||vmf||
|:-:|:-:|:-:|:-:|:-:|
|||クラスタA'|クラスタB'|クラスタC'|
||クラスタA|81|1|0|
|skmeans|クラスタB|0|8|82|
||クラスタC|0|68|0|

## まとめと所感

sphereclusterを使って、超球上でクラスタリングを行うspherical k-meansとMixture of von Mises Fisher clustering modelを試した。

spherical k-meansはk-meansに、Mixture of von Mises Fisher clustering modelは混合ガウスクラスタリングモデルに対応する。

適切なクラスタ数を与えられる場合には、Mixture of von Mises Fisher clustering modelが良さそう。しかし、計算量の問題で若干難がある。
そのため、実用上は高速なspherical k-meansの方が何かと便利そう。


## 参考
- [球面クラスタリング (spherical clustering) - 朱鷺の杜Wiki](http://ibisforest.org/index.php?%E7%90%83%E9%9D%A2%E3%82%AF%E3%83%A9%E3%82%B9%E3%82%BF%E3%83%AA%E3%83%B3%E3%82%B0)
- [text miningや方向情報の特徴抽出に便利な，超球上に分布する点をクラスタリングする方法](https://qiita.com/TomHortons/items/5819f14ba65e732b5bc6)
- [角度データのための統計モデル](http://stat.w3.kanazawa-u.ac.jp/sympo16/kato.pdf)
- [球面k-means++法について
](https://www.jstage.jst.go.jp/article/fss/31/0/31_81/_pdf)