時系列のクラスタリングとアンサンブル
鈴木雄大, 池野心平, サトヤキ・ロイ
こんにちは、NRIの時系列予測チームの鈴木です。
第４弾までの記事で、iterativeモデルとglobalモデルを導入し、それぞれ優れた点があることを紹介しました。
今回は第５弾として、時系列データをクラスタリングによって分類し、アンサンブルに活かす方法をご紹介します。
1. 同じデータセットでも、生データを見ていくと分割して扱いたくなる
Walmartコンペティションの学習データから、50事例ほどランダムにピックアップしてきました。同一のデータセットに含まれるデータですが、明らかにIdによって特徴がバラバラです。
第２弾で見たように、各種iterativeモデルは系列IDによって適合性が大きくブレます（e.g. 短期間のトレンドに依存しすぎて折り返されたときに外す、ハズレ値に弱く極端なスパイクがあると引っ張られる、など）。
第４弾でご紹介したglobalモデルは、これらをすべて１つのモデルで予測しようとしますが、学習の結果Idによって向き不向きが出てしまうのは明らかです。（e.g. 説明変数の寄与が甚大な系列と軽微な系列がある。 大域的に目的変数のラグ１への依存度が高いlgbmモデルが出来上がったとき、ラグ１よりも周期的な要素や、ホワイトノイズ的な要素が強い系列に対する予測が外れる、など）。
「このIdにはモデルAを重視しよう」といった形で、系列ごとにモデルのアンサンブル比率をコントロールすることで、正確かつロバストな予測が得られるようになります。
ここでは系列を大まかなグループに仕分け、各グループ内でモデルの比率を制御する方法を考えます。
図１をよくみると、ある程度似た雰囲気を持つ系列が横に並んでいることがわかります。例えば１段目は細かな変動と２回あるスパイク、２段目は大きな季節性といった具合です。
実はこれは、時系列に特化したある高速なクラスタリング手法で分類した結果からサンプリングしています。次項でこの手法をご紹介します。
(ただしここでは意味解釈を優先して、SKUではなくDept単位での粗いクラスタリングに留めています。)
Canonical time-series characteristics
クラスタリング手法には色々ありますが、ここでは時系列に特化した方法を扱います。
大別すると一般に２つの方法があります。
- 時系列の「形 (shape)」 が似ているものを１つのクラスタにまとめる
- 時系列の「特性(characteristics)」が似ているものを１つのクラスタにまとめる
1. は「Dynamic Time Warping（DTW: 動的時間収縮法)」と呼ばれる手法です。「２つの折れ線を横にスライドさせると重なる」ものを「距離が近い」とみなして、２つの折れ線グラフ同士の距離を算出します。こうして定義された距離を用いてクラスタリングが可能です。しかしDTW距離は定義から計算時間が非常に長いため、今回は採用しません。
2. は各時系列の「特性」（たとえば自己相関の強さ、周期性の強さなど）を複数取り出し、これを座標(z1, z2, ...)とする空間の点とみなす方法です。この空間上で距離を測ることで時系列同士をクラスタリングします。特性の似た特性を持つ点を「近い」とみなすということです今回はこちらの方法を採用します。
時系列データの特性を測る指標は非常にたくさんあるのですが、特に効果的なものをピックアップしている研究(*1)があり、これを活用しましょう。
このクラスタリングの方法を具体的な実行コードで解説します。まず多系列・時系列データを以下の形式のオブジェクトtrとして準備します。
head(tr)
| Id | .index | .value |
| 1 | 2010/2/5 | 881833.4 |
| 1 | 2010/2/12 | 1457182 |
| 1 | 2010/2/19 | 1118257 |
| 1 | 2010/2/26 | 681391.6 |
| 1 | 2010/3/5 | 762652.6 |
| 1 | 2010/3/12 | 803886.9 |
# settings ----
library(tidyverse)
library(magrittr)
library(progress)
library(Rcatch22) # 特性抽出を行ってくれるパッケージ
feat <- tr
feat %<>% nest(df = c(.index, .value))
feat %<>%
mutate(.feat = map(df, function(df) {
catch22_all(df$.value)
}))
feat %<>%
select(Id, .feat) %>%
unnest(.feat) %>%
arrange(names, Id)
feat %>% head %>% print
| Id | names | values |
| 1 | CO_Embed2_Dist_tau_d_expfit_meandiff | 0.093 |
| 10 | CO_Embed2_Dist_tau_d_expfit_meandiff | 0.077 |
| 11 | CO_Embed2_Dist_tau_d_expfit_meandiff | 0.093 |
| 12 | CO_Embed2_Dist_tau_d_expfit_meandiff | 0.554 |
| 13 | CO_Embed2_Dist_tau_d_expfit_meandiff | 0.128 |
| ... | ... | ... |
catch22によって各Idから、特に有用な22の特性指標が一気に計算されます。（この計算は非常に高速です。）
あとはこれを用いてクラスタリングしていきましょう。
feat %<>% pivot_wider(names_from=names, values_from=values) %>% filter(complete.cases(.)) feat %<>% mutate(across(where(is.numeric), scale)) pca <- feat %>% select(-Id) %>% prcomp tsne <- Rtsne::Rtsne( pca$x, pca=FALSE, perplexity = floor((nrow(pca$x) - 1) / 3) ) info.tsne <- tibble(tsne1 = tsne$Y[, 1], tsne2 = tsne$Y[,2]) hc.tsne = hclust(dist(tsne$Y)) info.tsne$hclust = factor(cutree(hc.tsne, 5)) hc.tsne.cent = info.tsne %>% group_by(hclust) %>% select(tsne1, tsne2) %>% summarize_all(mean) p <- ggplot(info.tsne, aes(x = tsne1, y = tsne2, colour = hclust)) + geom_point(alpha = 0.75) + theme_bw() + ggrepel::geom_label_repel(aes(label = hclust), data = hc.tsne.cent) + guides(colour = FALSE) clst <- bind_cols( info.tsne, feat %>% select(Id) ) clst %<>% select(Id, clst = hclust)
このように、各Idを特性指標で点にマップすることで、クラスタリングを行うことができます。
実務上で使うときにはある程度工夫の余地があり、例えば
- 生の時系列ではなく、製品分類などで括りを使い合計値に対してクラスタリングを施す（意味解釈や逐次更新バッチ時の安定性を優先）
- 平滑化処理を行い、細かなノイズ(fluctuation)以外の部分でクラスタリングする
- クラスタリング手法を変更する(PCA+kmeansなど)
Walmartコンペティションデータにおいては、（技巧より実務への応用性を重視するポリシーのため）Dept単位でまとめた系列に対してクラスタリングを施しました。
2010年と2011年12月の大きな山を持ち、それ以外は平坦な１つ目のクラスタ、減少傾向を持つ２つ目のクラスタ、全期間で平坦な３つ目のクラスタなどが摘出されています。
このクラスタを使って、以下のような解釈の基、モデルの配分を決めました。
- prophetは系列の安定性が高い場合に有用
- SARIMAXは人間の目から見て予測が困難なケースで使用（black box forecast）
- スパイクが多く立つ系列はetsのようにシンプルな平滑化か、lgbmの外部変数の活用で処理すべき
- クラスタ４のように予測しやすい系列は、多くのモデルをアンサンブルさせることで予測する
このクラスター別のアンサンブル予測をSubmitすることでスコアが
| モデル名 | private | public |
| クラスターアンサンブル | 2291.69 | 2235.81 |
となり、private / public共に精度上昇によってリーダーボードを越すことができます。
これで当初の目的であった、1位/688チームが達成されます。
次回の内容
多系列・時系列予測モデルを第１回から第５回まで行ってきました。次回はこれを見通しよく整理していきます。
またここまで深層学習は出てきていませんでしたが、実際には時系列予測に関する深層学習はいくつか使いやすいものが実装されています。これも紹介する予定です。
(*1 Lubba, C.H.; Sethi, S.S.; Knaute, P.; Schultz, S.R.; Fulcher, B.D.; Jones, N.S. Catch22: Canonical time-series characteristics. Data
Min. Knowl. Discov. 2019, 33, 1821–185)
