LoginSignup
35
20

More than 5 years have passed since last update.

明治大学Advent Calendar 2018Day 6
@ZoneTsuyoshi

Toeplitz Inverse Covariance-based Clustering ~時系列クラスタリング~

Last updated at Posted at 2018-12-05

　今回の記事では，time series clustering（時系列クラスタリング）に関する概要を述べ，その後に KDD 2017 paper である TICC (Toeplitz Inverse Covariance-based Clustering) という手法について解説する．
　また，説明の際に Markov Random Filed（マルコフ確率場, MRF），Expectation Maximumzation algorithm（EM algorithm）, Dynamic Programming（動的計画法, DP）, ADMM を用いる．MRF の説明は軽く行い，EM, DP, ADMM に関しては理解していなくても読み進められると思う．

Time Series Clustering

　Time series clustering には，大きく３つのタスク（種類）がある．

  • Whole time series clustering
    • 個々の時系列集合に対するクラスタリングのこと
    • (ex) 患者100人の脳波データを似ているデータごとにクラスタリング
  • Subsequence time series clustering
    • 一本の時系列データのセグメント（時系列データをある間隔ごとに区切ること）に対するクラスタリングのこと
    • (ex) 100秒ごとに区切られた一人の患者の脳波データをクラスタリング
  • Time point clustering
    • 一本の時系列データの各時点ごとのクラスタリング．ただし，どのクラスにも属さない（ノイズと判断される）時点もある
    • (ex) 一人の患者の脳波データの各時点の値に対するクラスタリング

　今回，説明していく TICC は time point clustering，subsequence time series clustering を同時に行える手法である1．

TICC の手法の概要

　TICC は時系列データのセグメントとクラスタリングを同時に行う手法である．各クラスターごとに固有の多層 MRF を有すると仮定し，EM algorithm 的に最適解を探索して行く．

Markov Random Field

　Markov Random Filed (MRF) とは，各特徴量をノードと捉え，特徴量間の関係をエッジとして捉えるグラフィカルモデルの一種である．特徴量$x_1,x_2$間に関係性がなければ，エッジの重みが0，すなわち$x_1,x_2$の間にエッジがないと考える．逆に，$x_1,x_2$間のエッジの重みが非0ならば，特徴量$x_1,x_2$間の関係性を重みで表す2．

mrf_graph.png

　今回考える多層 MRF は次のように図式化される．各レイヤーが一つの時点における MRF を表しており，その時点における特徴量をノード，特徴量間の関係をエッジとして表している．また，時点を跨る特徴量間の関係を異なるレイヤー間のエッジで表している．

スクリーンショット 2018-09-05 22.06.05.png

学習概要

　最適化問題を解く際に，EM algorithm 的に E step, M step を収束するまで繰り返す形式となる．E step においてDP，M step において ADMM を使用していく3．

最適化問題の定式化

最適化問題のためのセットアップ

　長さ $T$ の時系列データ $\{x_t\}_{t=1}^T$ を考える．各時点の観測は $n$ 次元であるとし，これらを窓 $w$ 毎にまとめて横に並べたベクトル

$$X_t=(x_{t-w+1}^T,\cdots,x_t^T)\in\mathbb{R}^{nw}$$

を考える．構成としては，図でイメージするとわかりやすい．

subclust.png

　ここから最適化問題をセットアップしていく．クラスター数を $K$ とし，各クラスターは固有の精度行列 $\Theta_i\in\mathbb{R}^{nw\times nw}$ を有するとする．精度行列とは，ガウス型逆共分散行列のことであるが，MRF のネットワークの繋がりを記述したものである．例えば，$(j,k)$ 要素が0でなければ，ネットワークにおいては，変数 $x_j,x_k$ の間にエッジがあるということである．すなわち，クラスター毎にネットワーク構造が異なり，そのネットワーク構造の違いをもとにクラスタリングしようというのである．ここで，$\Theta_i$ はブロック Toeplitz 行列になっており，

\Theta_i=\left[\begin{array}{cccccc}
A^{(0)}&(A^{(1)})^T&(A^{(2)})^T&\cdots&\cdots&(A^{(w-1)})^T\\
A^{(1)}&A^{(0)}&(A^{(1)})^T&\ddots&&\vdots\\
A^{(2)}&A^{(1)}&\ddots&\ddots&\ddots&\vdots\\
\vdots&\ddots&\ddots&\ddots&(A^{(1)})^T&(A^{(2)})^T\\
\vdots&&\ddots&A^{(1)}&A^{(0)}&(A^{(1)})^T\\
A^{(w-1)}&\cdots&\cdots&A^{(2)}&A^{(1)}&A^{(0)}
\end{array}\right]

と表せる．ただし，$A^{(j)}\in\mathbb{R}^{n\times n}$ である．ブロック Toeplitz 行列とは，このようにブロック対角が全て同じで，ブロック二重対角も（転置を除いて）同じといったように，対角からの外れ具合で同じ値を取るブロック行列のことを指す4．

image.png

MRF をイメージするために再度，上の画像を見て欲しい．この場合は，window size $w=3$ であり，

\Theta_i=\left[\begin{array}{ccc}
A^{(0)}&(A^{(1)})^T&(A^{(2)})^T\\
A^{(1)}&A^{(0)}&(A^{(1)})^T\\
A^{(2)}&A^{(1)}&A^{(0)}
\end{array}\right]

となる．$A^{(0)}$ はそれぞれのレイヤー（時点）におけるエッジを表している．例えば，画像の青いレイヤーに注目した時，どの変数同士が繋がっているかを表している．これが対角に並んでいるということは，紫や橙のレイヤーでも同じ変数同士の繋がりを持っているということを意味する．$A^{(1)}$ は隣のレイヤーの変数との繋がりを表している．例えば，画像の青・紫レイヤーに注目したときに，青レイヤーの変数と紫レイヤーの変数の繋がりを表している．$A^{(2)}$ 以降も同様に，２つ以上離れたレイヤーの変数間の繋がりを表現している．

　さて，最適化問題のセットアップに戻ろう．各クラスター $i$ に関して割り当て集合 $P_i$ を考える．例えば，$X_2,X_3,X_7$ がクラスター $i$ に属している場合，$P_i={2,3,7}$ となる．集合用語で述べると，${1,\cdots,T}$ を $P_1,\cdots,P_K$ に直和分解することになる．

定式化

　では，最適化問題を定式化しよう．最適化問題は次式のようになる：

\mathrm{arg}\min_{\Theta\in\mathcal{T},P}\sum_{i=1}^K\left[\|\lambda\circ\Theta_i\|_1+\sum_{X_t\in P_i}\left(-ll(X_t,\Theta_i)+\beta\mathbf{1}\{X_{t-1} \notin P_i\}\right)\right].

これを TICC 問題と呼ぶことにする．第１項はスパースな MRF を得るための l1 norm である．第２項は，

$$ll(X_t,\Theta_t)=-\frac{1}{2}(X_t-\mu_i)^T\Theta_i(X_t-\mu_i)+\frac{1}{2}\log\det\Theta_i-\frac{n}{2}\log(2\pi)$$

で記述できる対数尤度関数である．正則化項がなければ，この項を最小化する問題となる．第３項は，temporal consistency と呼ばれる項であり，前の時点と別のクラスターに遷移するのに罰則 $\beta$ を課している．すなわち，クラスターが変わりすぎないようにするための項である．

部分問題

　TICC 問題は直接最適化するのが困難である．そこで，幾つかの部分問題に分割し，それらを順番に解くことで $P,\Theta$ の最適解を見つけていく．

クラスター割り当て（E step）

　クラスター割り当てにおいては，部分問題

$$\mathrm{arg}\min_{P}\sum_{i=1}^K\sum_{X_t\in P_i}\left(-ll(X_t,\Theta_i)+\beta\mathbf{1}{X_{t-1} \notin P_i}\right)$$

を解く．この問題は，動的計画法（DP）を用いれば，大域最適解を見つけることができる．動的計画法とは，あるコストを最小にするパスを見つける際に用いられるアルゴリズムである．具体的には，途中時点では最もコストの小さいパスだけを保存し，最後にコスト最小パスをバックワードで見つけるものである．

image.png

画像を見て欲しい． 例えば，$-ll(2,1)$ に行くには $K$ 通りのパスがある．$-ll(1,1)$ から来た場合のコストは $-ll(1,1)$，$-ll(1,2)$ から来た場合のコストは $-ll(1,2)+\beta$ といった具合になる．最終時点に到達したときに最小コストのパスを見つけたいので，途中時点においてもその時点に至るまでのパスの中で最小コストのものだけ保存しておけば十分だろう．

精度行列の更新（M step）

　精度行列の更新においては，次の部分問題（TGL; Toeplitz Graphical Lasso 問題）

\mathrm{arg}\min_{\Theta_i\in\mathcal{T}}-\log\det\Theta_i+\mathrm{tr}(S_i\Theta_i)+\frac{1}{|P_i|}\|\lambda\circ\Theta_i\|_1

を各 $\Theta_i$ について解けば良い（各 $\Theta_i$ は並列に解くことができるので，以後簡単のため $\Theta$ と表記する）．ここで，$S_i$ はクラスター $i$ の標本共分散である．この形の問題に対しては，ADMM のような手法が有効であることが知られている．具体的には，$\Theta=Z$ なる $Z$ を導入してラグランジュ型の問題に落とし込んで行く：

\mathrm{minimize\quad}-\log\det\Theta+\mathrm{tr}(S\Theta)+\|\lambda\circ Z\|_1,\\
\mathrm{subject\ to\quad}\Theta=Z,Z\in\mathcal{T}.

ここで，$\mathcal{T}$ は $nw\times nw$ Toeplitz 行列全体を示す．こうすることで，拡張ラグランジュ法が適用可能になる．拡張ラグランジュ法は，次のように制約条件の二乗ノルムを最適化問題に敢えて加えることで，最適化が容易になるテクニックだと思ってもらえれば良い（詳細は最後に述べる文献参照）．

\mathcal{L}_p(\Theta,Z,U)=-\log\det\Theta+\mathrm{tr}(S\Theta)+\|\lambda\circ Z\|_1+\frac{\rho}{2}\|\Theta-Z+U\|_F^2

ここで，$Z$ は consensus 変数，$\rho>0$ は ADMM 罰則パラメータ，$U\in\mathbb{R}^{nw\times nw}$ は scaled 双対変数と呼ばれている．

ADMM はこのように拡張ラグランジュ法を使った手法の一つと思ってもらえればよく，この定式化によって，TGL 部分問題は次の３ステップを収束するまで繰り返せば良い．

\mathrm{(a)}\quad\Theta^{k+1}:=\mathrm{arg}\min_\Theta\mathcal{L}_p(\Theta,Z^k,U^k)\\
\mathrm{(b)}\quad Z^{k+1}:=\mathrm{arg}\min_{Z\in\mathcal{T}}\mathcal{L}_p(\Theta^{k+1},Z,U^k)\\
\mathrm{(c)}\quad U^{k+1}:=U^k+(\Theta^{k+1}-Z^{k+1})

$\Theta, Z$ の更新も閉形式（cos, exp などのよく知られた関数で記述できる形式）で書くことができ，陽に各ステップ計算できる．また，TGL 部分問題は凸問題であることから，ADMM は大域解への収束が保証されている．

実験結果

Macro F1-score

　F1-score とは，precision と recall の調和平均で定義される．クラスター $i$ の予測に関して，precision とは，クラスター $i$ に属すると予測したもののうち実際に属していた割合である．一方，recall とは，クラスター $i$ に属するもののうちクラスター $i$ と予測された割合である．一般に，precision と recall はトレードオフの関係であり，F1-score を用いることがある．この F1-score を各クラスターに対して計算し，その平均を取ったものが Macro F1-score である．Macro F1-score は高い方が良く，0以上1以下である．
　擬似生成したデータにおいて，各手法で Macro F1-score を比べた結果が次の表である．Temporal Sequence にある「1,2,1」とは，正解系列が1,2,1のクラスター遷移をしているということである．ここで，１つのクラスターに滞在している時点数は$100\times$（クラスター数）である．

image.png

ここで，GMM は混合ガウス法，EEV は正則化 GMM，DTW は動的時間伸縮法，GAK,Euclidean はそれぞれ global alignment kernel，Euclidean distance metric，Neural Gas は自己組織化マップ，K-means は Kmeans 法である．また，1クラスターに滞在するサンプル数を変化させた時の Macro F1-score の遷移は次のようになったと述べられている．

image.png

TICC が他の手法を上回っていることがわかる．

Case Study

Case Study として，車の走行状態のクラスタリングが示されている．7変数（ブレーキ，二次元加速度，方向，速度，RPM，ガソリン）を用いて，クラスター数5つとして TICC を適用したところ表のように，各クラスターの解釈が容易であった．また，GPS 位置に対して，クラスターを色で当てはめたところ，ブレーキを踏んでいる状態や曲がっている状態等のクラスタリングに成功している．

image.png

image.png

終わりに

　今回は TICC という時系列クラスタリング手法について述べてきた．時系列クラスタリングは，確立した手法が少なく，未だ未だ発展可能性を秘めている分野である．今回の手法は，https://github.com/davidhallac/TICC に公開されているので，ぜひ試してみてほしい．

参考文献

  1. 正確には，一本の時系列データの segment と subsequence time series clustering を同時に行う手法である． ↩

  2. MRFの写真は， Kaikai Liu, Ting An, Hao Cai, Lirida Naviner, Jean-Francois Naviner, Herv´e Petit. A General Cost-effective Design Structure forProbabilistic-Based Noise-Tolerant Logic Functionsin Nanometer CMOS Technology. 2013. より拝借した． ↩

  3. このように E step, M step を順次繰り返して最適化する手法を EM algorithm と呼ぶと思ってくれれば良い． ↩

  4. Toeplitz 行列と言った場合は，対角成分，二重対角成分が全て同じ行列のことを言う．空間的・時間的距離に応じて値を一緒にしてパラメータ数を減らしたい場合などに用いられる． ↩

35
20
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
35
20