# YADING: Fast Clustering of Large-Scale Time Series Data
Navigation: [[時系列クラスタリング]] | [[密度ベースクラスタリング]] | [[Rui Ding]] | [[Dongmei Zhang]]
> [!abstract] 概要
> ビッグデータの時代において、高速でスケーラブルな分析技術は、データ分析におけるリアルタイムで対話的な体験を可能にする基盤技術として重要性を増している。時系列は多様な応用分野で広く利用できる。時系列の個数が多く(例: 数百万)、各時系列の次元も高い(例: 数千)ため、大規模な時系列のクラスタリングは難しく、対話的な探索を支えるリアルタイム処理はいっそう難しい。本論文では、大規模な時系列を高速かつ高品質に自動でクラスタリングする、新しいエンドツーエンドの時系列クラスタリングアルゴリズム YADING を提案する。具体的には、YADING は三つの段階からなる。入力データセットのサンプリング、サンプリングしたデータセット上でのクラスタリング、残りの入力データのサンプリング側で得たクラスタへの割り当てである。特にサンプル数の下限と上限を理論的に証明し、YADING の高い性能を保証するとともに、入力データセットとサンプルの分布の一貫性を担保する。また、類似度尺度に L1 ノルムを、クラスタリング手法に複数密度アプローチを選ぶ。理論的な限界により、この選択は位相の揺れとランダムノイズに起因する時系列の変動に対して YADING が頑健であることを保証する。評価では、典型的な規模(各 1,000 次元の時系列 100,000 本)のデータセットで、YADING は最先端のサンプリングベースのクラスタリングアルゴリズム DENCLUE 2.0 の約 40 倍、DBSCAN と CLARANS の約 1,000 倍高速であることを示した。YADING は Microsoft のプロダクトチームでもサービス性能の分析に使われており、本論文ではその 2 事例を紹介する。
## 論文情報
- **タイトル**: YADING: Fast Clustering of Large-Scale Time Series Data
- **著者**: [[Rui Ding]]、[[Qiang Wang]]、[[Yingnong Dang]]、[[Qiang Fu]]、[[Haidong Zhang]]、[[Dongmei Zhang]]([[Microsoft Research]]、北京)
- **媒体**: PVLDB Vol. 8, No. 5(第 41 回 VLDB、2015 年 8 月 31 日〜9 月 4 日、ハワイ)
- **発表年**: 2015
## 概要
YADING は、大規模時系列を対話的な速さで自動クラスタリングするエンドツーエンドの手法である。データ削減(ランダムサンプリングと PAA による次元削減)、サンプル上の複数密度クラスタリング(L1 距離)、残りの系列の割り当ての 3 段階からなる。サンプル数の上下限が入力規模に依存しないこと、密度半径を kdis 曲線の変曲点から自動推定できること、割り当てを三角不等式による枝刈りで削減することが技術の要点である。全体の計算量は入力の系列数と次元にほぼ線形である。
## 問題設定
- 対象は $N$ 本・各 $D$ 次元の時系列 $\mathcal{T}_{N\times D}$ である。監視の現場ではデータセットが約 1 万〜10 万本・約 100〜1,000 次元に及び、電子商取引では数百万本になる。
- クラスタリングは、サーバ群の共通の性能プロファイルの把握や、時系列をカテゴリ属性へ変換して他のカテゴリ属性と併せて分析するために必要になる。例として、ネットワーク使用のスパイクが共通するサーバ群が同じルータを共有するなら、そのルータの問題を疑える。
- 対話的な探索では高速性が要る。既存研究の多くは精度向上のための高コストな類似度尺度に注力しており、スケーラビリティと性能に応える研究は少ない。
- 入力データの分布に仮定を置かない。形状が不規則で密度の異なるクラスタと、どのクラスタにも属さないノイズが混在する。
## 提案手法
全体像は 3 段階である。計算量は $O(Nsd + ND + s^2 \log s + sD \log D)$ で、$s$ はサンプル数($N$ に依存しない)、$d$ は削減後の次元である。$s^2\log s$ と $sD\log D$ の項は $N$ が大きいとき無視でき、実質は $N$ と $D$ に線形になる。
### データ削減
- **ランダムサンプリング**: 入力分布を仮定しないためランダムサンプリングを選ぶ。サンプル数 $s$ の決め方に二つの補題を与える。
- 補題 1(下限): 母集団比 $p_i$ のクラスタから信頼水準 $1-\alpha$ でサンプル内に $m$ 個以上を得るには $s \ge \frac{m + z_\alpha(z_\alpha/2 + \sqrt{m + z_\alpha^2/4})}{p_i}$ が要る。$p_i>1\%$、$m=5$、信頼水準 95% なら $s_l \ge 1{,}030$ で、これを下回ると 1% 未満のクラスタが見落とされる確率が無視できない。
- 補題 2(上限): 許容誤差 $\epsilon$、信頼水準 $1-\alpha$ のとき、全クラスタで母集団比の偏差 $|p_i-p_i'|<\epsilon$ を保つには $s\ge z_{\alpha/2}^2/(4\epsilon^2)$ で足りる。$\epsilon=0.01$、95% なら約 9,600 である。$p_i(1-p_i)\le 1/4$ で全比率を一括して抑えるため緩い上限である。
- どちらも $N$ に依存しない。実験に基づき、実運用では $s=2{,}000$ を採る。
- **PAA による次元削減**: 系列を $d$ 個のフレームに分け、各フレームの平均で表す。フレーム数 $d$ は標本化定理の考えから周波数の上限を推定して決める。各系列の自己相関曲線の最初の極小の位置 $y'$($g_i(y')<0$ のとき半周期に対応)から代表周波数を求め、その 80 パーセンタイルを周波数上限の近似にする。自己相関は FFT で求める。計算量は $O(sD\log D + ND)$ である。
### 複数密度クラスタリング
- **L1 距離を採る理由**: 計算が軽く、衝撃性ノイズに頑健である。位相のずれが小さければ L1 距離も小さい。位相差が大きい系列も、小さなずれをもつ系列の連鎖でつながれば密度ベースで同じクラスタに入る。
- **密度半径**: 各点の $k$ 近傍までの距離 kdis を降順に並べた曲線を kdis 曲線とし、最頻の kdis 値を密度半径と定義する。補題 4・5 により、kdis 曲線の変曲点の Y 値が密度半径に等しい。複数の密度があれば変曲点も複数現れる。
- **推定アルゴリズム**: 左右の傾きの差が最小の点を変曲点として選び、その左右の部分曲線で再帰する。有意な変曲点がなくなれば止める。kdis 曲線の作成が $O(ds^2)$、推定が $O(s\log s)$ である。
- **クラスタリング**: 密度半径を小さいものから順に DBSCAN($k=4$)を適用し、抽出した点を除いて次の半径へ進む。残りをノイズとする。計算量は $O(s^2\log s)$ である。
### 割り当てと枝刈り
- 未ラベルの系列は、最近傍のラベル付き系列との距離がそのクラスタの密度半径未満ならそのクラスタに入れ、そうでなければノイズとする。素朴には $O(Nsd)$ で全体で最大の項になる。
- 三角不等式による枝刈り: 未ラベル点 $a$ がラベル点 $b$ から密度半径 $\varepsilon$ より遠い距離 $dis$ にあれば、$b$ の近傍のうち $b$ から $dis-\varepsilon$ 以内の点も $a$ から $\varepsilon$ より遠い。この計算を省く。
- 整列近傍グラフ(Sorted Neighbor Graph, SNG)を、コア点ごとに全サンプル点への距離を昇順で保持する構造として持つ。二分探索で省ける近傍を飛ばす。実際に 2〜4 倍の高速化になる。最悪計算量は $O(Nsd)$ のままである。
原論文の擬似コードは、フレーム数の推定が Table 1、密度半径の推定が Table 2、複数密度の DBSCAN が Table 3、割り当てが Table 4 である。本ページでは文章で要約し、擬似コードは転載しない。
### 位相の揺れへの頑健性(補題 3)
初期位相 $a$ の系列と、位相が $\delta_i\in[0,\Delta]$ 一様にずれた $n$ 本の系列が同じクラスタになる確率は、$\varepsilon$ を距離しきい値、$m$ を系列長、$M$ を定数として $P(E_n)\ge 1-n\left(1-\frac{\varepsilon}{mMk\Delta}\right)^n$ と下から抑えられる。データ数 $n$ が増えるほど、許容できる位相の揺れ $\Delta$ が大きくなる。ランダムノイズについての解析は著者のプロジェクトサイトに置かれている。
## 新規性
- 時系列クラスタリングの各段(サンプリング、次元削減、密度推定、割り当て)にパラメータの自動推定を組み込み、人手のパラメータ設定を不要にしたエンドツーエンドの手法である。
- サンプル数の上下限を、入力規模に依存しない形で理論的に導いた。既存の DENCLUE 2.0 は分布を保つ根拠なしにサンプル数を選ぶとしている。
- 密度半径の定義と、kdis 曲線の変曲点による複数密度の自動推定。従来の密度推定は手動か低速であった。
- 三角不等式に基づく割り当ての枝刈りと、その支えとなるデータ構造 SNG。
- L1 距離と複数密度クラスタリングの組み合わせについて、位相の揺れへの頑健性の理論限界を与えた。
## 実験設定
- **計算機**: 2.4 GHz Intel Xeon E5-2665、128 GB RAM。比較の公平のため YADING と比較手法はすべて単一スレッド実装である。
- **データ**: 5 種の確率モデル(AR(1)、強制振動、ドリフト、ピーク、ランダムウォーク)から生成した合成データ。テンプレート A(15 グループ、グループの大きさ 0.1〜30%)を RQ1・RQ2 とノイズの評価に、テンプレート B(8 グループ、うち 2 つは強制振動の位相 β を変えて位相の揺れを作る)を位相の評価に使う。実データは UCR の StarLightCurves(3 クラスタ、9,236 本、1,024 次元。当時の UCR で最大規模)を使う。
- **指標**: 実行時間、NMI(正規化相互情報量)、正しく特定できたクラスタ数 NCICluster。
- **比較手法**: エンドツーエンドで DENCLUE 2.0(YADING と同じサンプル数を使用)、DBSCAN(MinPts=4、$\varepsilon$ は最も有意な密度半径)、CLARANS(クラスタ数は正解を与える)。複数密度の推定は DECODE と比較する。
- **研究設問**: RQ1 効率、RQ2 サンプル数と精度の関係、RQ3 時系列の変動への頑健性。
## 実験結果
- **データセット**: 使用した実データは Table 5(StarLightCurves の 3 クラスタ・9,236 本・1,024 次元)にまとめられている。
- **RQ1(効率)**: 系列数 $N$ を $10^4$〜$10^5$、次元 $D$ を 64〜1,024 に変えた 15 個の合成データで、YADING は DENCLUE 2.0 の約 1 桁、DBSCAN と CLARANS の約 3 桁高速で、$N$ と $D$ にほぼ線形である。100 万本・2,000 次元でも 91.7 秒であった。
- **StarLightCurves**: Table 6 のとおり YADING は 3.1 秒で、NMI は DBSCAN・CLARANS と同程度である。3 手法とも NMI が高くないのは、ラベル 1 と 3 の系列が近く、まとめて 1 つのクラスタにするためである。DENCLUE 2.0 は高次元向けに調整しても品質が低い。
![[wiki/sources/_attachments/2015__VLDB__YADING-Fast-Clustering-of-Large-Scale-Time-Series-Data/table6-starlight.png]]
- **複数密度の推定**: Table 7 のとおり、ラベルではなく kdis 曲線の急な変化の数を正解とすると、全 15 データセットで密度は計 58 個である。YADING は 55 個、DECODE は 40 個を特定し、DECODE の 40 個のうち 37 個は YADING も特定した。推定時間は DECODE が平均 500 秒超、YADING が 100 ミリ秒未満である。
![[wiki/sources/_attachments/2015__VLDB__YADING-Fast-Clustering-of-Large-Scale-Time-Series-Data/table7-density-estimation.png]]
- **RQ2(サンプル数)**: 100 回繰り返した平均で、Table 8(NMI)と Table 9(NCICluster)のとおり NMI はサンプル数 200 でも大きくは下がらないが、NCICluster は $s=200$ で 7.0、$s=500$ で 8.9 と大きく落ち、全体は 12 である。$s=200$・$500$ のデータセット 1 では 15 グループ中 8 個・6 個を見落とし、1% 未満のグループがすべて含まれる。$s=1{,}000$ では最小の 4 グループを見落とすが、0.5% と 0.8% のグループは残る。見落としたグループはサンプル内で 5 個未満であり、下限 1,030 の理論と整合する。$s\ge 2{,}000$ で精度はほぼ飽和する。
| サンプル数 | 200 | 500 | 1K | 2K | 5K | 全数 |
|---|---|---|---|---|---|---|
| NMI(DS1 10K) | 0.915 | 0.932 | 0.948 | 0.952 | 0.956 | 0.965 |
| NMI(DS2 100K) | 0.857 | 0.919 | 0.939 | 0.943 | 0.946 | 0.960 |
| NCICluster(DS1) | 7.0 | 8.9 | 10.6 | 11.3 | 11.6 | 12 |
| NCICluster(DS2) | 6.3 | 8.7 | 10.3 | 11.3 | 11.3 | 12 |
- **RQ3(頑健性)**: ランダムノイズを含む 15 データセットの平均 NMI は Table 10 のとおりで、YADING が最も高く分散も小さい。位相の揺れを含むテンプレート B の 2 データセットで、YADING の NMI は 0.982 と 0.988 であり、位相の揺れをもつ 2 組の系列を正しく分けた。
| 手法 | YADING | DENCLUE 2.0 | DBSCAN | CLARANS |
|---|---|---|---|---|
| 平均 NMI | 0.925±0.027 | 0.523±0.057 | 0.820±0.071 | 0.804±0.018 |
### 図
- 実サービスの時系列を PCA で 2 次元に射影すると、エネルギーの 91% を保ったまま、密度の異なる不規則な形のグループとノイズが見える(Figure 1)。この観察が複数密度クラスタリングの動機になる。
![[wiki/sources/_attachments/2015__VLDB__YADING-Fast-Clustering-of-Large-Scale-Time-Series-Data/fig01-pca-projection.png]]
- 4 近傍距離の曲線(4-dis)には 3 つの変曲点があり、対応する値が 3 つの密度半径の近似になる(Figure 2)。
![[wiki/sources/_attachments/2015__VLDB__YADING-Fast-Clustering-of-Large-Scale-Time-Series-Data/fig02-4dis-curve.png]]
- 割り当ての枝刈りは、ラベル点 $b$ から遠い未ラベル点 $a$ に対し、$b$ の近くの点との距離計算を省く(Figure 3)。
![[wiki/sources/_attachments/2015__VLDB__YADING-Fast-Clustering-of-Large-Scale-Time-Series-Data/fig03-assignment-pruning.png]]
- 系列数と次元に対する実行時間の比較(Figure 4)。縦軸は対数である。
![[wiki/sources/_attachments/2015__VLDB__YADING-Fast-Clustering-of-Large-Scale-Time-Series-Data/fig04-performance.png]]
- 位相の揺れとノイズを含む 1 グループの系列(灰色)と基準(黒)。YADING はこれらを同じクラスタにまとめた(Figure 5)。
![[wiki/sources/_attachments/2015__VLDB__YADING-Fast-Clustering-of-Large-Scale-Time-Series-Data/fig05-phase-perturbation.png]]
## 実運用の適用事例
YADING を核に、時系列をその場でグループ化する対話的な多次元分析ツールを作り、Microsoft のオンラインサービスのチームが使った。
- **チーム A(性能カウンタ、Figure 6)**: 46,000 台のサーバの CPU 使用率と記憶域使用率(各 1,008 点)を、約 4.5 秒で 21 グループにした。CPU 使用率のグループは、極めて低い(グループ 2、25.63%)、中程度、高い(異なるパターン)に分かれた。低使用率のグループを記憶域使用率などの属性で切り替えると、その属性の時系列がその場で再クラスタリングされる。
![[wiki/sources/_attachments/2015__VLDB__YADING-Fast-Clustering-of-Large-Scale-Time-Series-Data/fig06-cpu-groups.png]]
- **チーム B(エラー分布、Figure 7)**: サーバごとに約 100 種のエラー事象の 24 時間の発生回数を曲線にしたエラー分布(約 1 万台)を 1 秒未満で 16 グループにした。健全なグループ 1 と、問題があるらしいグループ 2 を得た。グループ 2 を配置のトポロジで切ると全サーバが特定のサーバクラスタ X に属し、さらに、そのサーバ群に入れたソフトウェアのバージョンの問題を突き止められた。過去に約 1 日かけた同種の調査が数分で済んだ。過去の性能問題のエラー分布をクラスタリングして説明と対処を付ける知識ベースの構想も挙がった。
![[wiki/sources/_attachments/2015__VLDB__YADING-Fast-Clustering-of-Large-Scale-Time-Series-Data/fig07-error-distribution-groups.png]]
## 考察
- 現実装は単一スレッドである。割り当て段は並列化しやすく、局所性鋭敏ハッシュ(LSH)による加速は今後の課題とされる。
- 入力系列は等長を前提とする。可変長には整数因子のダウンサンプリング等の前処理が要る。
- L1 距離以外にも、距離の公理(非負性・対称性・三角不等式)を満たす尺度は枠組みに載る。$L_p$ ノルムは超球の体積係数を替えるだけで使え、ARIMA 由来のモデルパラメータ空間への写像も可能である。DTW と Pearson 相関は距離の公理を満たさないため、簡単には載らない。
- ストリームには、L1 距離の逐次更新と距離行列で多くの段が対応できる。自己相関を逐次更新できない PAA の次元削減だけが難しく、将来課題である。
- iSAX による索引は、粒度の人手指定と非可逆圧縮による密度推定の偏りの懸念から採らなかった。$R^*$ 木は高次元で素朴な総当たりより遅かった。
## 強み / 弱点・課題
- 強み: 入力規模に依存しないサンプル数の理論保証、パラメータの自動推定、大規模データでの速さ、実運用での対話的分析の実績。
- 強み: 位相の揺れとノイズへの頑健性に理論限界を与えた。
- 弱点: 実データの定量評価は StarLightCurves 1 件のみで、実運用の 2 事例は定性的な報告である(機密のため詳細は伏せられている)。
- 弱点: 割り当て段は最悪 $O(Nsd)$ のままで、枝刈りの効果はデータ分布に依存する(2〜4 倍)。
- 弱点: 比較手法は単一スレッド実装で、DENCLUE 2.0 は高次元向けに著者が調整した。CLARANS にはクラスタ数の正解を与えるなど、設定の有利不利が混在する。
- 弱点: サンプル数の補題は正規近似などの仮定を置く(緩い上限であることは論文も認める)。位相の揺れの補題は解析関数の仮定を置く。
- 課題: 可変長系列、DTW や相関のような非距離尺度、ストリームの次元削減。
## 関連
- 概念: [[時系列クラスタリング]] / [[密度ベースクラスタリング]] / [[時系列次元削減]] / [[サンプリング手法]] / [[次元の呪い]]
- ソース: [[@1996__KDD__A Density-Based Algorithm for Discovering Clusters in Large Spatial Databases with Noise]] / [[@2016__SIGMOD Record__k-Shape - Efficient and Accurate Clustering of Time Series]] / [[@2005__Pattern-Recognit__Clustering of time series data - a survey]]
- エンティティ: [[Rui Ding]] / [[Qiang Wang]] / [[Yingnong Dang]] / [[Qiang Fu]] / [[Haidong Zhang]] / [[Dongmei Zhang]] / [[Microsoft Research]]
## 出典
- `.raw/papers/2015__VLDB__YADING-Fast-Clustering-of-Large-Scale-Time-Series-Data.pdf`(PVLDB Vol. 8, No. 5)