2章:ホテリングの $T^2$ 法による異常検知

はじめに

この記事は,異常検知と変化検知を学ぶアドベントカレンダー2018 – Takala’s Memory2日目の記事です.MLPシリーズの 「異常検知と変化検知」 の2章まとめとなります.

素人がoutput駆動の勉強用に書いています.不備や勘違いもあるかもしれません.差し支えなければ間違え等については,ご指摘・ご教授頂けると助かります.よろしくお願いします.

また,本記事において使用している図表はすべて,上記参考書を基にして作成したものです.

コンテンツ

  • 訓練データのモデリング
  • 異常度の定義
  • ホテリングの$T^2$法を用いた異常判定

訓練データのモデリング

異常検知の基本的な枠組みは,

  • 訓練データから正常時のモデルを作成する
  • そのモデルに基づき,新たな観測データに対する異常度を定義する
  • 異常度が閾値異常であれば異常と判定する

である.今回扱う,古典的な外れ値検出手法(異常検知手法)であるホテリングの$T^2$法は,正常時のモデルというものを正規分布で表現するところからスタートする.

具体的な手順を見ていこう.

まず,$M$次元の$N$個の観測値からなるデータ:

$$
\mathcal{D} = \left\{ \boldsymbol{x}^{(n)} | n=1,2,..,N \right\}
$$

を考える.このとき,ホテリングの$T^2$法では,このデータ$\mathcal{D}$が,異常標本を含まないか,含んでいたとしても圧倒的に少数であるという前提のもと,各標本が独立に,次の確率密度関数に従うと仮定する.
$$
\mathcal{N}(\boldsymbol{x} | \boldsymbol{\mu}, \boldsymbol{\Sigma})
=\frac{\boldsymbol{\Sigma}^{-1/2}}{(2\pi)^{M/2}}
\exp \left\{
-\frac{1}{2} \left( \boldsymbol{x} – \boldsymbol{\mu} \right)^\mathsf{T}
\boldsymbol{\Sigma}^{-1}(\boldsymbol{x} – \boldsymbol{\mu})
\right\}
$$

関数の形からなんとなくわかると思うが,これは多次元の正規分布(ガウス分布)を表す確率密度関数である.1次元の場合の正規分布のパラメータは平均・分散であったが,多次元の場合には平均ベクトル$\boldsymbol{\mu}$および,分散共分散行列$\boldsymbol{\Sigma}$(*1)が確率密度関数のパラメータとなる.(このある確率密度関数でデータを表現するという行為がまさしくモデリングである.)

これらのパラメータを最尤推定法によって算出する.つまり,モデルのパラメータを決定するということだ.

すぐわかると思うが,正規分布のパラメータを最尤推定法に求める作業は非常に簡単で,結局のところデータ$\mathcal{D}$の平均(ベクトル)・分散(共分散行列)を求めることに他ならない(*2).

結局,パラメータの最尤推定量(データの平均ベクトル・分散共分散行列)を,$\hat{\boldsymbol{\mu}}, \hat{\boldsymbol{\Sigma}}$とすると,データ$\mathcal{D}$を表現する確率密度関数は以下のように表現される.

$$
p(\boldsymbol{x} | \mathcal{D}) = \mathcal{N}(\boldsymbol{x} | \hat{\boldsymbol{\mu}}, \hat{\boldsymbol{\Sigma}})
$$

ここまでで,データを確率密度関数で表現し,正常時のモデルというものを作成した.次のステップで異常度を定義する.

異常度の定義

さきほどの確率密度関数$p(\boldsymbol{x} | \mathcal{D})$を用いて異常度$a(\boldsymbol{x}’)$を

$$
a(\boldsymbol{x}’ ) = (\boldsymbol{x}’ – \hat{\boldsymbol{\mu}})^\mathsf{T}
\hat{\boldsymbol{\Sigma}}^{-1}(\boldsymbol{x}’ – \hat{\boldsymbol{\mu}})
$$
と定義する.これは,観測データ$\boldsymbol{x}’$が,どれだけ標本平均$\hat{\boldsymbol{\mu}}$から離れているかを表すもので,マハラノビス距離(の2乗)と呼ばれる.

一般的な多次元におけるユークリッド距離(の2乗)と異なるところは,$\hat{\boldsymbol{\Sigma}}^{-1}$の部分である.$\hat{\boldsymbol{\Sigma}}$が分散共分散行列であることを考えると,各軸を標準偏差で重みづけをしていると直感的に考えことができる.ある軸について,標準偏差が大きいほど(ばらつきが大きいほど),その軸の重みは大きくなる(*3).

さてここまでで,ある観測データ集合(訓練データ)$\mathcal{D}$に対し,新たな観測値$\boldsymbol{x}’$の異常度$a(\boldsymbol{x}’)$を定義した.残された問題は,異常度$a(\boldsymbol{x}’)$の閾値をどのように設定するかということになる.

異常判定

ホテリングの$T^2$法

いきなりであるが,さきほどのマハラノビス距離として定義された異常度$a$について以下の定理が成り立つことが知られている.(この証明は難しいらしい,異常検知という枠組みでは道具として使えればOKだそうな)

$M$次元正規分布$\mathcal{N}(\boldsymbol{x} | \boldsymbol{\mu}, \boldsymbol{\Sigma})$からの,$N$個の独立標本$\left\{ \boldsymbol{x}^{(i)}| i=1,2,…,N \right\}$に基づき,標本平均$\hat{\boldsymbol{\mu}}$および,標本分散共分散行列$\hat{\boldsymbol{\Sigma}}$を求める.このとき,新たな観測値$\boldsymbol{x}’$について,以下が成立する.
1. $\boldsymbol{x}’ – \hat{\boldsymbol{\mu}}$は,平均$\boldsymbol{0}$で分散共分散行列が$\frac{N+1}{N}\boldsymbol{\Sigma}$の$M$次元正規分布に従う.

2. $\hat{\boldsymbol{\Sigma}}^{-1}$は,$\boldsymbol{x}’ – \hat{\boldsymbol{\mu}}$と統計的に独立である.

3. $T^2 = \frac{N-M}{(N+1)M} a(\boldsymbol{x}’)$によって定義される統計量$T^2$は,自由度$(M,N-M)$のF分布に従う.

4. $N \gg M$の場合は,$\boldsymbol{x}’$は,近似的に自由度$M$,スケール因子$1$のカイ2乗分布に従う.

この統計量$T^2$は,ホテリング統計量としばしば呼ばれる.そしてこの定理は ホテリングの$T^2$法と呼ばれるらしい.

ここで,注目すべきは,定理の4である.多くの場面において系の変数の次元($M$)は,サンプルの数$N$よりかなり小さいと考えられる.つまり上の定理の4が(ほぼ)常に成立するということになる.したがって,異常度$a$は,データの物理的単位や数値によらず,普遍的に,自由度$M$,スケール因子$1$のカイ2乗分布に従うということになる.(*4)

ちなみに,カイ2乗分布の確率密度関数は以下:
$$
\chi^2(u | k,s) = \frac{1}{2s\Gamma(k/2)} \left(\frac{u}{2s} \right)^{(k/2)-1}\exp \left(-\frac{u}{2s} \right)
$$
で与えられる.ここで$s$がスケール因子,$k$が自由度である.$\Gamma$は以下に示すガンマ関数のことである.
$$
\Gamma(z) = \int^{\infty}_{0} dt \ t^{z-1} \mathrm{e}^{-t}
$$
式の形は非常に見にくい.実務上は覚える必要はないだろう.

異常判定モデルの構築

新たな統計量$T^2$というものを定義し,それがカイ2乗分布に従うといことがわかった.あとは,統計の分野の検定の手続きと同様のことをすればよい.(カイ2乗検定)

具体的に,カイ2乗分布に基づく異常判定モデルの構築手順を示していく.

まず,誤報率$\alpha$を定める.これは全観測量の端何%を異常と判断するかという基準である.一般的にはは,$0.05,0.01$であったりすることが多い(*5).この$\alpha$に基づき,カイ2乗分布から方程式
$$
1-\alpha = \int^{a_\mathrm{th}}_{0} dx \chi^2(x|M,1)
$$
により,閾値$a_{\mathrm{th}}$を求めておく.(仰々しい方程式になっているが,実際にはこれを解く必要はなく,自由度と$\alpha$を定めたときの$a_\mathrm{th}$というのは一覧表のようなもので容易に知ることができるはず.)

正常標本が圧倒的だと信じられるデータ(*6)から標本平均および標本共分散行列:
$$
\hat{\boldsymbol{\mu}} = \frac{1}{N} \sum^{N}_{n=1} \boldsymbol{x}^{(n)}
\
\hat{\boldsymbol{\Sigma}} = \frac{1}{N} \sum^{N}_{n=1}
\left(\boldsymbol{x}^{(n)}-\hat{\boldsymbol{\mu}} \right)^{\mathsf{T}}
\left(\boldsymbol{x}^{(n)}-\hat{\boldsymbol{\mu}} \right)
$$
を計算しておく.

新たな観測値$\boldsymbol{x}’$を得るとき,異常度としてのマハラノビス距離
$$
a(\boldsymbol{x}’)= (\boldsymbol{x}’-\hat{\boldsymbol{\mu}})^{\mathsf{T}}
\hat{\boldsymbol{\Sigma}}^{-1}
(\boldsymbol{x}’-\hat{\boldsymbol{\mu}})
$$
を計算し,$a(\boldsymbol{x}’) > a_{\mathrm{th}}$ならば異常と判定する.

この一連の流れがホテリングの$T^2$法による異常検知となる.

少し,狐に鼻をつままれたような感覚になるのは,マハラノビ距離がカイ2乗分布に従うという点だろう.この点については,正規分布とカイ2乗分布の関係として参考書に簡単に議論されている.本記事ではとりあげない.

コメント

  • つまり普通の統計検定とプロセスは同じだよね?
  • 毎回思うが,カイ2乗分布の確率密度関数複雑すぎる...
  • マハラノビス距離の定義の中の,分散共分散行列の意味解釈については釈然としないところがあるので,今後の課題

注釈

(*1) : 「分散共分散行列」,「共分散行列」ともに同じ概念である.

(*2) : 最尤推定法については,初歩的な統計学の教科書には絶対に掲載されている.また,参考書(「異常検知を変化検知」)にも,最尤推定法の式展開が掲載されているが,ベクトル表記がされているため,慣れない人には理解がしにくいかもしれない.

(*3) : この説明には自信がない...

(*4) : これはすごいことのような気がするが,統計を少しでもかじったことのある人なら検定の手続きとしてinputされているので,当然のことのようにも感じてしまう.個人的には「へーそう」ぐらいの感想しか持ちえないのが正直なところ.

(*5) : 統計学の文脈では有意水準と言ったりする.

(*6) : 当たり前だが,正常時のモデルに基づいて異常判定をしたいので,訓練データは正常時の観測値の方が(量として)支配的である必要がある.

1