1次元のデータをカーネル密度推定(KDE)でクラスタリングする

次のような1次元のデータをクラスタリングしたくなることがあるかもしれません。

```python
data=[3016,1556,348,1185,341,1189,342,425,343,419,349,1190,340,424,345,422,347,
      423,344,1187,343,1189,354,1166,350,423,345,1188,345,425,342,425,342,424,
      345,421,347,1184,346,1187,347,419,345,425,344,421,346,422,346,424,345,
      8269,3014,1561,344,1183,345,1188,343,432,335,428,344,1185,343,421,346,
      425,342,427,342,1187,343,1185,345,1187,342,424,344,1188,343,419,349,424,
      344,424,346,424,343,1188,343,1189,340,422,346,422,348,422,344,428,341,
      419,350,8267,3013,1562,343,1186,353,1179,342,419,347,424,346,1187,343,
      427,341,420,348,419,353,1186,340,1188,344,1189,336,422,348,1186,344,425,
      344,423,344,421,349,419,349,1189,339,1185,343,428,340,428,340,422,347,
      422,346,426,343]
```

このデータは人間なら次のように分割することになるでしょう。

```json
[348 341 342 343 349 340 345 347 344 343 354 350 345 345 342 342 345 347
 346 347 345 344 346 346 345 344 345 343 335 344 343 346 342 342 343 345
 342 344 343 349 344 346 343 343 340 346 348 344 341 350 343 353 342 347
 346 343 341 348 353 340 344 336 348 344 344 344 349 349 339 343 340 340
 347 346 343]
[425 419 424 422 423 423 425 425 424 421 419 425 421 422 424 432 428 421
 425 427 424 419 424 424 424 422 422 422 428 419 419 424 427 420 419 422
 425 423 421 419 428 428 422 422 426]
[1185 1189 1190 1187 1189 1166 1188 1184 1187 1183 1188 1185 1187 1185
 1187 1188 1188 1189 1186 1179 1187 1186 1188 1189 1186 1189 1185]
[1556 1561 1562]
[3016 3014 3013]
[8269 8267]
```

しかし、同じことを計算機で実現する方法は自明ではありません。

クラスタリングといえばk-meansが定番ですが、これは多次元データ向けのクラスタリング手法になります。ソートしてから次元を増やせば適用できそうですが、おそらくオーバースペックでしょう。

そこで、1次元データに対してカーネル密度推定(KDE)を適用し、密度の極小値で区間分割することでクラスタリングしてみました。

## 実装

カーネル密度推定(KDE)にはscikit-learnを使っています。また、極小値の検出にSciPyを使いました。

```python
import numpy as np
from sklearn.neighbors import KernelDensity
from scipy.signal import argrelextrema

data=[3016,1556,348,1185,341,1189,342,425,343,419,349,1190,340,424,345,422,347,
      423,344,1187,343,1189,354,1166,350,423,345,1188,345,425,342,425,342,424,
      345,421,347,1184,346,1187,347,419,345,425,344,421,346,422,346,424,345,
      8269,3014,1561,344,1183,345,1188,343,432,335,428,344,1185,343,421,346,
      425,342,427,342,1187,343,1185,345,1187,342,424,344,1188,343,419,349,424,
      344,424,346,424,343,1188,343,1189,340,422,346,422,348,422,344,428,341,
      419,350,8267,3013,1562,343,1186,353,1179,342,419,347,424,346,1187,343,
      427,341,420,348,419,353,1186,340,1188,344,1189,336,422,348,1186,344,425,
      344,423,344,421,349,419,349,1189,339,1185,343,428,340,428,340,422,347,
      422,346,426,343]

a=np.array(data)
a2=np.log2(a)
minval, maxval = a2.min(), a2.max()

kde = KernelDensity(kernel='gaussian', bandwidth=0.1).fit(a2.reshape(-1,1))
s = np.linspace(minval, maxval, 1024)
e = kde.score_samples(s.reshape(-1,1))

mi = argrelextrema(e, np.less)[0]

for i in range(0,len(s[mi])):
    print(a[(a2 >= minval) * (a2 <= s[mi][i])])
    minval = s[mi][i]
print(a[(a2 >= minval)])
```

これを実行すると上の結果が得られます。`data`を変えながら挙動を見ていくと「それっぽく」分割してくれます。

`bandwidth`は大きすぎても小さすぎてもうまくいかないので、問題の性質にマッチする値を指定する必要があります。グリッドサーチで良い`bandwidth`を探せば自動化できるのかもしれません。

`a2`で元の値のlogを取ってからKDEにかけていることに注意してください。今回はデータの特性として絶対値の差分より倍率を見たいデータだったので、こうした方が望む結果に近づきそうだと考えました。

差分の方が重要なデータであればlogを取らない方が良いかもしれません。その場合は`bandwidth`と`s`の分割数を調整する必要があるでしょう。

## 参考URL

下記の記事が本稿の元ネタです。

- [machine learning \- How would one use Kernel Density Estimation as a 1D clustering method in scikit learn? \- Stack Overflow](https://stackoverflow.com/questions/35094454/how-would-one-use-kernel-density-estimation-as-a-1d-clustering-method-in-scikit/)

下記の記事ではグリッドサーチで適切な`bandwidth`を探しています。

- [scikit\-learnでカーネル密度推定 \| 分析ノート](https://analytics-note.xyz/statistics/scikit-learn-kde/)