sktimeのRocketClassifierを調べてみると、ROCKETという変換器を通したデータを分類器に食わせているだけだった

# はじめに

[過去記事](https://qiita.com/takahashi-ry/items/e94e630a2088e8626c51)でsktimeの時系列分類のユーザガイドで、`RocketClassifier`というモデルを使用してたのを目にしました。そもそも時系列分類という分野自体も聞いたことがなかった浅学な筆者。これを機にモデルについて少し調べてみることにしました。

また、本記事作成のために書いたコードは以下のGistにあります。よろしければ、そちらも併せてご覧ください。

https://gist.github.com/takahashi-r/a4a94f99332101fa1f67f8a4137afd7e

# sktimeとは

sktimeについては簡単にですが過去に記事を書いたのでそちらをご覧ください。

https://qiita.com/takahashi-ry/items/e94e630a2088e8626c51

READMEから一部抜粋です。

>🚀 Version 0.13.4 out now! Check out the release notes here.
>
> sktime is a library for time series analysis in Python. It provides a unified interface for multiple time series learning tasks. Currently, this includes time series classification, regression, clustering, annotation and forecasting. It comes with time series algorithms and scikit-learn compatible tools to build, tune and validate time series models.

要約：
sktime は、Python で時系列分析を行うためのライブラリで、分類、回帰、クラスタリング、アノテーション、および予測ができる。つまり何でもできる。

## install 

### pip 

```
pip install sktime
```

### conda

```
conda install -c conda-forge sktime
```

# RocketClassifierとは

まずは[公式doc](https://www.sktime.net/en/stable/api_reference/auto_generated/sktime.transformations.panel.rocket.Rocket.html#sktime.transformations.panel.rocket.Rocket)を見に行ってみました。

>Classifier wrapped for the Rocket transformer using RidgeClassifierCV.
>
>This classifier simply transforms the input data using the Rocket [1] transformer and builds a RidgeClassifierCV estimator using the transformed data.

おっと？早速書いてありますね。ざっくり訳すと、

>>データをRocketという変換機を使用して変換し、変換されたデータを使用して RidgeClassifierCVで推定を行っている

`RidgeClassifierCV` はsklearnの線形分類モデルなのでsktimeとは関係がありません。ずいぶん早いタイトル回収ですが、`RocketClassifier`は新しいモデルではなくて、入力データをRocketで変換して線形分類モデルに食わせたものだということがわかりました。（正確に言うと`StandardScaler`も間に噛んでいるんですが、`RocketClassifier`の実装にはそこまで本質的な話ではないので割愛）

さて、Rocketっていったいどんな変換器なのか見ていきましょう！

# ROCKETってどんな変換器？

Let's [公式doc](https://www.sktime.org/en/stable/api_reference/auto_generated/sktime.transformations.panel.rocket.Rocket.html#rocket)

一発目にRocketの意味が書いてあります。

> RandOm Convolutional KErnel Transform
>> ランダム畳み込みカーネル変換

畳み込みカーネルと言えば、畳み込みニューラルネットワーク（CNN: Convolutional Neural Network）がパッと思いつきますが、それとなにか関係あるのでしょうか？ランダムがいったい何を指すのか気になりますね。
参考論文名が記載されているので、詳しくはそっちを見ていきましょう！

# ROCKETの論文を読んでみた

論文のタイトルは「[ROCKET: Exceptionally fast and accurate time series classification using random convolutional kernels](https://arxiv.org/abs/1910.13051)」

abstractをざっくり要約すると

>> 時系列分類のほとんどの方法は計算が複雑です。つまり、小さなデータセットでもたくさんのトレーニング時間が必要で、大きなデータセットになるとより扱いにくくなるということです。ランダムな畳み込みカーネルを使用する単純な線形分類器が、既存の方法の計算コスト大幅に上回り、精度が同程度であることを示します。

CNNなどでは特徴抽出を行うためのカーネルは学習によって手に入れていましたが、Rocketではランダム生成することによっていい感じのカーネルが手に入るといったことが言われているみたいですね。カーネルを学習しなくていいので、その分の計算コストが浮くのがうれしい㌽みたいです。

## ランダム畳み込みカーネルの作り方

じゃあ、ランダムカーネルとやらはいったいどうやって作ってるんだっていう話をこれからしていきたいと思います。
カーネルの重みってなんだっけ？ダイレーションってなんだっけ？そもそもカーネルってなんだよって人は[付録](#付録-カーネルについて)で簡単にまとめてありますので、そちらを先にご覧ください。

ではソースコードと一緒にカーネルの生成方法を紐解いていきましょう。

カーネルを生成している関数は`_generate_kernels`なので、以降このメソッドから該当箇所を切り抜いていきます。

https://github.com/sktime/sktime/blob/4194aacaa105ef59807ea56c1bb3497b37eb1613/sktime/transformations/panel/rocket/_rocket.py#L116-L189


### weights：重み

重みは標準積分布からサンプリングされます。サンプリング後、重みは平均が0になるように調整を行います。

```
_weights = np.random.normal(0, 1, _num_channel_indices * _length).astype(np.float32)
# loop中
_weights[a3:b3] = _weights[a3:b3] - _weights[a3:b3].mean()
```

※`_num_channel_indices`はもう少し下に説明があります

### lengths：長さ

長さは「7,9,11」から等確率でランダムに選択されます。

```lengths.py
candidate_lengths = np.array((7, 9, 11), dtype=np.int32)
lengths = np.random.choice(candidate_lengths, num_kernels).astype(np.int32) # num_kernelsはカーネルの数でdefaultで10,000
```

### biases：バイアス

バイアスは一様分布からサンプリングします。

```biases.py
biases[i] = np.random.uniform(-1, 1)
```

### dilations：ダイレーション

ダイレーションは指数分布からのサンプリングを行います。

```dilations.py
dilation = 2 ** np.random.uniform(
    0, np.log2((n_timepoints - 1) / (_length - 1))
) # n_timepointsはデータの長さ
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1272355/60ed7474-edf2-6fb0-9868-f0a93f8b5f6d.png)

```dilation_plot.py
def f(n_input, kernel_length): 
    return 2 ** np.random.uniform(0, np.log2((n_input - 1) / (kernel_length - 1)))

n_input = 200
kernel_length = 7
plt.hist([f(n_input, kernel_length) for _ in range(1000)])
plt.title("dilation")
```


### paddings：パディング

パディングは1/2でパディングを追加するかどうかを決めています。パディングが追加される場合系列の両端も畳み込みの対象となります。

```padding.py
padding = ((_length - 1) * dilation) // 2 if np.random.randint(2) == 1 else 0 # _lengthはカーネルの長さ→[7,9,11]のどれか
```

### num_channel_indices

データが多次元の場合、全ての特徴量を畳み込むわけではなく、畳み込みの対象となる特徴量もランダムに指定しています。各カーネル毎に特徴量を **いくつ** 使うのかを決めているのが、この`num_channel_indices`です。

```num_channel_indices.py
num_channel_indices = np.zeros(num_kernels, dtype=np.int32)
for i in range(num_kernels):
    limit = min(n_columns, lengths[i])
    num_channel_indices[i] = 2 ** np.random.uniform(0, np.log2(limit + 1))
```

上記のコードを括りだして適当にサンプリングした結果は以下のようになります。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1272355/60545dd6-85b2-82b2-a365-6bdff57f4f27.png)

カーネルの長さが最大で11なので、`limit`が大きくても11にしかなりません。$\log_2 (11+1) = 3.584...$なので、多次元で畳み込むことのできる特徴量のサイズは12程度になります。

```num_channel_indices_plot.py
num_kernels = 10000
n_columns = 20

candidate_lengths = np.array((7, 9, 11), dtype=np.int32)
lengths = np.random.choice(candidate_lengths, num_kernels).astype(np.int32)
num_channel_indices = np.zeros(num_kernels, dtype=np.int32)
for i in range(num_kernels):
    limit = min(n_columns, lengths[i])
    num_channel_indices[i] = 2 ** np.random.uniform(0, np.log2(limit + 1))

channel_indices = np.zeros(num_channel_indices.sum(), dtype=np.int32)
plt.hist(num_channel_indices)
plt.title("num_channel_indices")
```

### channel_indices

先ほど、各カーネル毎に使用する特徴量の数を決めました。今度は、**どの**特徴量を使用するのかを決めます。

```channel_indices.py
# replace=False は非復元抽出を行うという意味
# 例えば5次元データに対して、_num_channel_indicesが3の場合
# [0, 1, 2, 3, 4] から、3つを等確率で非復元抽出するという操作を行う
channel_indices[a2:b2] = np.random.choice(
    np.arange(0, n_columns), _num_channel_indices, replace=False
)
```

# 終わりに

本記事では省きましたが、カーネルのダイレーションを一様分布からサンプリングした場合や、カーネルサイズを5,000にしてみると精度がどうなるのか、ほかにも比較検討が行われているので、気になるかたはぜひ論文を読んでみてください。

# おまけ

**？**「畳み込みカーネルは特徴的な形状を探すやつだから、ランダム生成でいい形のカーネルが出てくるのがポイントだと思ったんだけど、例えば100次元とかのときってランダム生成でいい感じのカーネルが見つからない気がするんだけどどうなの？」
**僕**「特に言及されていなかった気がします。関連の論文見てないんでよくわかんないですけど、とりあえず実際にやってみます」

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1272355/eed61cca-c45e-6943-ce28-1fb73eb18be9.png)

**僕**「あ、ほんとだ。精度良くないですね」
**？**「多次元拡張はなんかアイディア加えないと難しい印象だなぁ。論文読んでないから詳しくはわかんないけどね」


# 付録 カーネルについて

カーネルを理解するには畳み込みを理解しないといけません。畳み込みとは、データから特徴を抽出する作業のことをいいます。そしてその特徴抽出はカーネルを使って行います。つまり、カーネルはデータから「特徴抽出するためのフィルタ」のような役割をするものです。

画像をご覧ください。特徴抽出はカーネルの重み付きの和で表されます。

```math
2 * 2 + 4 * 0.5 + 5 * 0.3 = 7.5
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1272355/e6f93864-7702-6194-4b44-c8f2d8cf67b5.png)

以降、カーネルを横にスライドさせて次々と重み付きの和を計算していくのが畳み込みです。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1272355/775e9725-d1ee-9130-ca66-27758b843727.png)

さて、ここで本記事で取り扱うカーネルを制御するパラメータのようなものが、画像のどこに当たるのかを説明していきます。

## 重み

画像中のカーネルの数字が重みに当たります。CNNではカーネルの重みも推定すべきパラメータになりますが、Rocketでは重みは確率で取り出すだけなので、重みを学習することはありません。

## 長さ

長さよりもサイズといったほうがピンときやすいかもですね。上記画像で行くと3です。

## バイアス

CNNの場合で考えてみましょう。$w$をカーネル、$X$を画像データ、$b$をバイアスだとすると

```math
\begin{align}
Y_{conv} &= \sum_{c} W[:,:,c] * X[:,:,c] \\
H &= \phi(Y_{conv} + b )
\end{align} 
```

となります。ここで、$\phi$は活性化関数を表しています。プーリングはこのHに対して行っていく操作です。Rocketでは活性化関数を通さずにプーリングを行うことに注意しましょう（活性化関数に恒等関数を指定したという見方もできますが）。

ちなみに画像中では$Y_{conv}$部分しか計算していないので、バイアスは出てきていません。

## パディング

パディングとは畳み込みを行うときにカーネルの中心をどの部分から始めるのかということを決めるパラメータです。例えば、画像中の畳み込みを例に挙げると、カーネルの中心である重み0.5はデータの左から二番目の要素から始まっています。この場合のパディングは0になります。

ではパディング1とはどういう状態になるのでしょうか。パディング1はカーネルの中心が左端に合わせられるときのことを言います。この場合、データからカーネルがはみ出すのですが、はみ出した分は0埋めを行います。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1272355/bbbe18a2-d14c-0d04-0a43-d7bc1b413591.png)

これによって何が起きるのかというと、以下のような畳み込みを行った場合$Y_{conv}$のサイズが元のデータよりも小さくなるのです。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1272355/5d1c423c-5499-f35c-ce0a-c3e519a401a2.png)

もちろんパディングを2以上にすると出力サイズのほうが入力サイズよりも大きくなります。

## ストライド

ストライドはカーネルをずらしていくと幅を指定するパラメータです。画像中の順々にずれていく畳み込みはストライドが1であるといいます。ストライドを2にするとセルを二つずらすことになります。ちなみにRocketではストライドは1固定です。

## ダイレーション

ダイレーションなのかディレーションなのかどっちかわからなかったので、とりあえずダイレーションとさせてください。

ダイレーションは畳み込みの対象となるデータのセルの幅を拡張するものです。なにを言ってるのかわからないと思うので、まずは画像をご覧ください。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1272355/d439d0ac-537e-9df5-11aa-5924440dc59c.png)

これはダイレーション2の場合の畳み込み計算の一部を表したものです。もともと畳み込みの対象となるデータは「2,4,5」だった（ダイレーション1）のが、「2,5,3」に代わっていることがわかります。

CNNのダイレーションについては下記記事の「ディレーション」の部分を参照ください。
https://dreamer-uma.com/pytorch-cnn/


# 参考
- [sktime.org](https://www.sktime.org/en/stable/index.html)
- [ROCKET: Exceptionally fast and accurate time series classification using random convolutional kernels](https://arxiv.org/abs/1910.13051)