DEIM Forum 2017 G3-4

尤度に基づくストリームデータの異常検出手法の感度に関する実験
山岸

祐己y

和巳y

斉藤

y 静岡県立大学経営情報イノベーション研究科 〒 4228526 静岡県静岡市駿河区谷田 521
E-mail: yyamagissygmail. om, yyk-saitou-shizuoka-ken.a .jp
あらましストリームデータにおける異常検出において，様々なバースト検知手法が考案されているが，一般的には

1-カテゴリー情報を扱うことを前提としているため，複数カテゴリー情報の分析には向いていない．特に代表的な手法は Kleinberg のバースト検知であるが，バースト検知そのものについても，想定されるバースト期間に応じて感度のパラメータを慎重に設定する必要がある．よって我々は，複数カテゴリー情報を持ったストリームデータにおける，多様な区間長の異常検出手法を提案する．提案手法は，複数カテゴリー情報の出現確率分布の変化に関する尤度最大化と尤度比検定に基づいて設計されており，ラフなパラメータ設定においても多様な区間長の異常（バースト）検出を可能としている．
キーワード 異常検出，バースト検知，ストリームデータ，時系列データ
1.

はじめに

ストリームデータにおける異常検出手法は，Kleinberg のバー
スト検知 [1℄ をはじめとして幅広く研究がなされており [2℄ [3℄，

定した状態で，提案手法の方が様々な区間長の異常に対応できることを人工データによる比較実験で示す．
2.

問題設定

近年では異常（バースト）の数理モデル化の研究も進んでいJ -カテゴリーの情報を持つストリームデータについて，それる [4℄．Kleinberg のバースト検知は，パラメータを多様に調整

らの確率分布における変化の異常と，異常の持続区間を推定すできることを考慮すると汎用性が高いように思えるが，バースることを考える．n 番目の観測ステップ tn でのカテゴリー情報

ト状態を検知する感度のパラメータ ( ) を決定した時点で想定を sn とすると，観測されたストリームデータ D はされるバースト（異常）区間の長さもある程度固定されてしまうため，想定よりも長期的な異常，すなわち潜在的な異常を検出できない恐れがある．更に，これらのバースト検知や異常検出手法 [1℄ [2℄ [3℄ は，一般的に 1-カテゴリーの情報を扱うために設計されているので，複数カテゴリーの情報をもつストリー
ムデータに対してはカテゴリー毎に独立して適応する必要があり，複数カテゴリー間の関係を考察することが難しい．例えば，

D = f( s1 ; t1 );    ; ( sN ; tN )g;(1)のように表せる．ここで， sn =2 f1;    ; J g である．便宜上， sn

はカテゴリー j 2 f1;    ; J g の J -次元ベクトルダミー変数として

8>>>< 1 if sn = j;
sn j = >>>: 0 otherwise:;(2)ある 2-カテゴリーがほぼ同時期にバーストを示したとしても，のように変換する．カテゴリー j が出現する確率 p j が多項分

そのバーストの度合いは各カテゴリーによって相対的に決めら布に従っていると仮定すると，ストリームデータ D に対するれているため，カテゴリー間の直接的な比較をすることはでき

対数尤度は，出現確率のパラメータベクトル p = f p1 ;    ; pJ g にない．また，これらの手法は観測データの時刻間隔に依存する

よって

ところが大きいため，観測時刻間隔がもともと一定になっているデータや，観測時刻間隔があまり変化しないデータには向いL(D; p) =ていない．
よって我々は，上記の異常検出手法が扱うようなデータ形式に加え，カテゴリーが複数含まれていたり，観測時刻間隔が一定だったりするような多種データを扱うことを前提として，多様な区間長の異常（バースト）を検出する手法を提案する．提

N X
J
X

sn j log p j ;;n=1 j=1(3)のように計算できる．式 (3) の最尤推定量はp j =

PN

n=1 sn; j

N

;

(4)

のように与えられる．

案手法は，代表的なバースト検知の研究 [1℄ [5℄ と同様の回顧的

ここからは，このパラメータベクトル p が，確率分布の変化(Retrospe tive) な観点から，データの確率分布に対する尤度最に従って K 箇所で分割された階段関数の形をとることを考える．

大化と尤度比検定に基づいて考案したものである．Kleinberg のすなわち，k 2 f1;    ; K g 番目の分割点の時刻 T k (t1 < T k < tN )手法におけるパラメータは，バーストのスケーリングとバース
ト検知の感度の 2 つを要するのに対し，提案手法におけるパラにおいてパラメータベクトルが pk から pk+1 に切り替わることを仮定する． K 個の分割点を持つ集合 CK = fT 1 ;    ; T K g とし，メータは，異常検出の感度を決定するものただ 1 つであるため，

便宜上 T k 1 < T k , T 0 = t1 , T K +1 = tN とする．更に，CK による D

汎用性が高い設計となっている．更に，感度のパラメータを固の分割を

Dk = fn; T k 1 < tn <= T k g;

(5)

3. 1

貪

欲

法

まず，貪欲法 (A1) の手順について説明する．このアルゴリズ

すなわちN = f1; 2;    ; N g = f1g [ D1 [    [ DK +1 ;ムは，バックトラッキングをしないデータの 2 分割の繰り返し(6)とし，jDk j は (T k 1 ; T k ℄ における観測ステップ数とする．ここで，任意の k 2= f1;    ; K + 1g について jDk j , 0 ，つまり少なくとも一つの tn 2 Dk が存在することを条件とする．この時点で，この分割点検出問題は，部分集合 CK  T の探索問題というこである．つまり，既に選択された (k1) 個の分割点

Ck 1 を固定したまま k 番目の分割点 T k を Ck 1 に新たに追加することを繰り返す．一般的に，2(LR(Ck )LR(Ck 1 )) は，N が十分に大2

きいとき  分布に従うことが知られているため，我々はこのアルゴリズムの終了条件として 2 検定を採用する．この 2 検とになる．ここで，T は観測ステップの集合 T = ft1 ; t2 ;    ; tN g

定の危険率は事前に設定する必要がある．貪欲法アルゴリズム

である．

の手順は以下となる．

分割点 CK によって与えられる D の対数尤度は，パラメータ

ベクトル集合 PK +1 = f p1 ;    ; pK +1 g を用いてL(D; PK +1 ; CK ) =

K +1
X
k=1

L(Dk ; pk );

A1-1.
A1-2.
(7)

推定量は
p k j =T k = arg maxtn 2T fLR(Ck 1 [ ftn g)g を探索する．

Ck = Ck 1 [ fT k g のように更新する．
LR(Ck 1 )) が，設定された危険率と自2
由度 J 1 における  の棄却限界値よりも小さければ，CK を
A1-3.
A1-4.

のように計算できる．つまり，式 (7) の k と j についての最尤

k = 1, C0 = ; のように初期化する．

もし 2(LR(Ck )出力して終了する．

P

n2Dk sn j
;

jDk j

;

;

(8)

K +1 X X
JX

sn j log p k j ;;k=1 n2Dk j=1;

k = k + 1 とし，A1-2 に戻る．

ここで，A1-3 での Ck の各分割点は，T i 1 < T i for i = 2;    ; kとなる．これらを式 (7) に代入すると，L(D; P K +1 ; CK ) =A1-5.

を満たすように再インデックスする．明らかに，このアルゴリ(9)が導ける．従って，この分割点検出問題は，式 (9) を最大化す
る CK の探索問題に帰着できる．しかし，式 (9) では分割点集合

CK の導入によってどれだけ対数尤度が改善したかという直接的な評価をすることができない．この問題において，分割を考慮しない，すなわちパラメータベクトルの変化が無いことを仮定したときの対数尤度からの改善度合いを評価することは重要となるため，尤度比最大化問題として目的関数を構築する．もズムの計算量は O(NK ) と高速であるため，大規模な N にたい
しても実用的な計算時間で結果を得ることが可能である．しかし，先程も説明したように，このアルゴリズムはバックトラッ
キングを行わないため，プアーな局所解に陥ってしまうことが危惧される．
3. 2

局所探索法

次に，局所探索法 (A2) について説明する．このアルゴリズ
ムは，A1 で得られた解 CK から始まり，分割点の改善を 1 つず

し，パラメータベクトルに一切の変化がない，すなわち C0 = ;つ試みるものである．つまり，分割点 T k を一度取り去り，残っとするならば，式 (9) はることを k = 1 から K まで繰り返す．ここで，  n  は集合差をL(D; P 1 ; C0 ) =JXX

sn j log p 1 j ;;n2N j=1;(10)となる．ここで，p 1 j =;P

n2N sn j

jNj

;

;

た CK n fT k g を固定して，よりよい尤度を得られる T k0 を探索す表す．もし，全ての k (k = 1;    ; K ) に対して分割点の置換が行われない，すなわち，全ての k に対して T k0 = T k ならば，これ以上の改善は望めないとして処理を終了する．局所探索法のア(11)

ルゴリズムは以下となる．

である．よって，K 個の分割点を持つ場合と分割点を一切持た

A2-1.

ない場合との対数尤度比はA2-2.

LR(CK ) = L(D; P K +1 ; CK ) L(D; P 1 ; C0 );(12)のように与えられる．最終的に，この分割点検出問題は LR(CK )を最大化する CK の探索問題に帰着できる．
3.

分割点検出法

A2-3.

k = 1, h = 0 のように初期化する．
T k0 = arg maxtn 2T fLR(CK n fT k g [ ftn g)g を探索する．
もし T 0 = T ならば h = h + 1 とし，さもなければ h = 0kkとして CK = CK n fT k g [ fT k0 g のように更新する．
A2-4.
A2-5.

もし h = K ならば CK を出力して終了する．

もし k = K ならば k = 1, さもなければ k = k + 1 とし，A2-2 に戻る．

式 (12) を網羅的に解くと最適解が保証されるが，計算量が明らかに，このアルゴリズムの計算量は改善が終わらない限り

O( N ) となってしまうため，ある程度大きい N に対して K >=3増え続けてしまうが，ある程度大規模な問題に対しても，せいとなってしまうと，実用的な計算時間で解くことができない．

ぜい貪欲法アルゴリズムの計算量 O(NK ) の数倍程度で終了すしたがって，我々は任意の K について解くための高速な解法をることを我々は既に実験によって示している [6℄．

K

提案する．以下では，まず貪欲法 (A1) と局所探索法 (A2) を説

3. 3

明し，更にそれらを組み合わせた提案解法について説明する．

もし，計算量を最低限に抑えることを目的として，単純に貪

提案解法

表1

欲法アルゴリズムと局所探索法アルゴリズムを組み合わせると，

CK を得る．
A2 で CK を改善する．

C2.

p 1

p 2

p 3

p1 (

)

1 =3

1 =3

1 =3

(

)

2 =4

1 =4

1 =4

(

)

6 =8

1 =8

1 =8

(

)

1 =3

1 =3

1 =3

D1
p2 D2
p3 D3
p4 D4

A1 で

C1.

3-カテゴリーの出現確率設定;;;となる．確かに，これだけでも十分な近似解が期待できるが，分割点数 K が貪欲法アルゴリズムによって決定されてしまうた

上手く機能しないという特徴が見られる．感度のパラメータはめ，問題に対して不適切な分割点数のまま局所改善を行ってし固定しているため，パラメータを変化させてさらに検証する余まう恐れが大いにある．したがって我々は，不必要な分割点は地はあるが，異常区間の長さが予測できない状況を考慮すると，極力追加せず，且つ必要な分割点は極力追加することを目的と提案手法のほうが多様な区間長の異常に適応できる可能性がした，アルゴリズムの反復的な組み合わせを提案する．提案解高いと言える．更に，各手法に要した計算時間の比較 (Intel(R)法の手順は以下となる．

Xeon(R) X5690 3.47GHz) を図 5,6 に示す．図より，提案手法

は個々のサンプルの問題に依存して計算時間がブレてしまって
P1.

A1-1 から処理を開始する．

P2.

A1-4 の処理後に k > 2 ならば， k を

P3.

C CK として出力する．
CK を A2 で改善し，改善した CK を Ck として出力する．

P4.

A1-5 から処理を再開させ，ステップ I2 へ戻る．

=

いるが，全体的な計算時間については，全ての場合において提案手法の方が短いことがわかる．

この手順では，分割点が追加される度に局所探索法アルゴリズ
ムを行うため，更なる計算量の増加が予想されるが，ある程度大規模な問題に対しても，せいぜい貪欲法アルゴリズムの計算
量 O(NK ) の数倍から十数倍程度で終了することを我々は既に実験によって示している [6℄．結局のところ，提案解法において事前に設定が必要となるのは，貪欲法アルゴリズムにおける 2
検定の危険率のみであり，これが分割点数を大きく左右するため，本手法における感度のパラメータということになる．
4.

実

験

4. 1

人工データによる比較実験

我々は，ストリームデータの異常検出において，尤度の最大

sum of absolute error of probability

2000
1800
1600
1400
1200
1000
800
600
400
200

化が有効な方法であると仮定し，その最適化問題についての効0100

率的な解法を提案した．ここでは，隠れマルコフモデルを用い

500

1000

5000

10000

time steps of anomaly periods

た Kleinberg( [1℄) の手法を技術水準とし，提案手法との比較実
験を行う．

図1

真の出現確率との絶対誤差の比較 (Kleinberg の手法)

実験で用いるのは 3-カテゴリーをランダムに発生させた人工ストリームデータで，出現確率を階段関数の形で変化させたものである．真の分割点数は K = 3 であり，異常（バースト）区

間は D2 と D3 とした．異常区間の両端部分の観測ステップ数

を jD1 j = 10000, jD4 j = 10000 で固定したまま，異常区間の観測

ステップ数は 5 パターン (100, 500, 1000, 5000, 10000) に変化さ
せ，それぞれ 100 サンプルずつ生成した．Kleinberg の手法のスケーリングパラメータ（今回は s = 1:5 ）に則って，出現確率の変化は表 1 のように設定した．なお，Kleinberg の手法の感= 1:0, 提案手法のアルゴリズムにおける 2
検定の危険率（感度のパラメータ）は p = 0:0001 で固定した．

度のパラメータは各手法による異常検出結果を基にした推定出現確率と，真の出現確率との絶対誤差の比較を図 1 から 4 に示す．図より，提案手法は極端に短い異常区間 jD2 j = jD3 j = 100 においては上手

参考までに，各異常区間における両手法の代表的な結果を
図 7 から 16 に示す
4. 2

現実データによる実験

我々は，国内の大規模化粧品レビューサイト  osme（注 1）から，各アイテムにつけられた今までのレビューの点数情報を取得
し，最も多くのレビューを有する 2 アイテム (Oshima Tsubaki
Camellia Hair Care Oil (Oshima Tsubaki) と Conditioner Essential (Albion)) に対して提案手法を適応した．人工データによる実験のときと同様，2 の危険率は p = 0:0001 だが， osme のレビュー点数の範囲は 0 から 7 ( j = 1;    ; 8 として考える) であるため，自由度は J

1= 7 となる．図 17 と 18 に 2 アイテムの得点確率分布の分割結果を示す．両図より，提案手法は適切な分割数と様々な分割区間長によって，レビュー得点の確率分

く機能していないが，それ以外の場合においての性能は良好である．逆に，Kleinberg の手法は，異常区間が長くなってくると（注 1）：http://www. osme.net/

300

sum of absolute error of probability

sum of absolute error of probability

2000
1800
1600
1400
1200
1000
800
600
400
200
0

250

200

150

100

50

0
100

500

1000

5000

10000

100

time steps of anomaly periods
図2

真の出現確率との絶対誤差の比較 (提案手法)

図4

5000

10000

真の出現確率との絶対誤差の比較の詳細 (提案手法)

0.6

250

processing time (sec.)sum of absolute error of probability

1000

0.7

300

200

150

100

50

0.5

0.4

0.3

0.2

0.1

0

0
100

500

1000

5000

20200

10000

21000

22000

30000

40000

time steps of stream data N

time steps of anomaly periods
図3500

time steps of anomaly periods

図5

真の出現確率との絶対誤差の比較の詳細 (Kleinberg の手法)

計算時間の比較 (Kleinberg の手法)

布における変化を表現できているように見える．更に，分割点

験では，両手法の異常検出に関する感度のパラメータを固定しが集中している期間や短い区間は得点分布の信頼性が有意に低た状態で，提案手法の方が様々な区間長の異常に対応できるこく，長い区間は得点分布の信頼性が有意に高いということも自とを示した．今回の比較実験では，計算時間においても提案手然に考えられる．

法が僅かに優れている結果となった．現実データを使った実験

5.

おわりに

では，レビュー点数の分布変化の視覚化に成功し，得点分布の信頼性指標における有用性を示した．

代表的な異常（バースト）検出手法特有の問題を補うため，
カテゴリーが複数含まれていたり，観測時刻間隔が一定だった
りするような多種データを扱うことを前提として，多様な区間長の異常（バースト）を検出する手法を提案した．我々は，謝辞 本研究は，JSPS 特別研究員奨励費 15K00311 の支援を受けて行ったものである．

データの確率分布に対する尤度最大化と尤度比検定に基づいた
文問題設定を行い，その問題を高速且つ高精度に解くための解法を提案した．人工データを使った Kleinberg の手法との比較実[1℄

献Jon M. Kleinberg. Bursty and hierar hi al stru ture in streams. In

Pro eedings of the Eighth ACM SIGKDD International Conferen e

0.8

0.6

0.7

0.5

0.6

probability

processing time (sec.)

0.7

0.4

0.3

0.5

0.4

0.2

0.3

0.1

0.2

0

1
2
3
proposed

0.1
20200

21000

22000

30000

40000

0.5

1

time steps of stream data N
図6

1.5

time step n

計算時間の比較 (提案手法)

図8

異常区間の観測ステップ

=

2
×104

100 における出現確率の絶対誤差の中央値に最も近い結果 (提案手法)

0.8
0.7

0.8

0.7

1
2
3
Kleinberg

0.6

0.5

probability

probability

0.6

1
2
3
Kleinberg

0.4
0.3

0.5

0.4

0.2
0.3
0.1
0.2
0
0.5

1

1.5

time step n
図7

異常区間の観測ステップ

=

2
0.1

×104

0.5

100 における出現確率の絶対誤差の中央値に最も近い結果 (Kleinberg の手法)

1

1.5

time step n
図9異常区間の観測ステップ

=

2
×104

500 における出現確率の絶対誤差の中央値に最も近い結果 (Kleinberg の手法)

on Knowledge Dis overy and Data Mining, July 23-26, 2002, Edmonton, Alberta, Canada, pp. 91101, 2002.
[2℄

ÆYunyue Zhu and Dennis E. Shasha. E

ient elasti

timelines. In Pro eedings of the 23rd Annual International ACM SI-

burst dete tion

GIR Conferen e on Resear h and Development in Information Retrieval, pp. 4956, 2000.

in data streams. In Pro eedings of the Ninth ACM SIGKDD Interna-

tional Conferen e on Knowledge Dis overy and Data Mining, Washington, DC, USA, August 24 - 27, 2003, pp. 336345, 2003.
[3℄

Aaron Sun, Daniel Dajun Zeng, and Hsin hun Chen. Burst dete tionfrom multiple data streams: A network-based approa h. IEEE Trans.

Systems, Man, and Cyberneti s, Part C, Vol. 40, No. 3, pp. 258267,
2010.
[4℄

Hua-Wei Shen, Dashun Wang, Chaoming Song, and Albert-László
Barabási.

Modeling and predi ting popularity dynami s via rein-for ed poisson pro esses. In Pro eedings of the Twenty-Eighth AAAI

Conferen e on Arti ial Intelligen e, July 27 -31, 2014, Québe City,
Québe , Canada., pp. 291297, 2014.
[5℄

Russell C. Swan and James Allan. Automati

generation of overview[6℄

Yuki Yamagishi,

Seiya Okubo,

Kazumi Saito,

Kouzou Ohara,

Masahiro Kimura, and Hiroshi Motoda. A method to divide stream
data of s ores over review sites. In PRICAI 2014: Trends in Arti ialIntelligen e - 13th Pa i Rim International Conferen e on Arti ialIntelligen e, Gold Coast, QLD, Australia, De ember 1-5, 2014. Proeedings, pp. 913919, 2014.

0.8

0.7

0.8

1
2
3
proposed

0.7

0.6

probability

probability

0.6

0.5

0.4

0.4

0.3

0.2

0.2

0.1

0.1
0.5

1

1.5

time step n
異常区間の観測ステップ

=

2
×10

0.5

500 における出現確率の絶対誤差の

図 12

2

異常区間の観測ステップ

=

×10

4

1000 における出現確率の絶対誤差

0.8

1
2
3
Kleinberg

0.7

1
2
3
Kleinberg

0.6

probability

0.6

0.5

0.4

0.5

0.4

0.3

0.3

0.2

0.2

0.1

0.1
0.5

1

1.5

time step n
図 11

1.5

の中央値に最も近い結果 (提案手法)

0.8

0.7

1

time step n4

中央値に最も近い結果 (提案手法)

probability

0.5

0.3

図 10

1
2
3
proposed

異常区間の観測ステップ

=

2
×10

0.5

1000 における出現確率の絶対誤差の中央値に最も近い結果 (Kleinberg の手法)

1

1.5

2

time step n

4

図 13

異常区間の観測ステップ

=

2.5

3
×10

4

5000 における出現確率の絶対誤差の中央値に最も近い結果 (Kleinberg の手法)

0.8

0.7

0.8

1
2
3
proposed

0.7

0.6

probability

probability

0.6

0.5

0.4

0.5

0.4

0.3

0.3

0.2

0.2

0.1
0.5

1

1.5

2

2.5

time step n
図 14

1
2
3
proposed

異常区間の観測ステップ

=

0.1

3
×10

0.5

1

4

5000 における出現確率の絶対誤差

1.5

2

2.5

3

3.5

time step n
図 16の中央値に最も近い結果 (提案手法)

異常区間の観測ステップ

=

4

×104

10000 における出現確率の絶対誤差の中央値に最も近い結果 (提案手法)

0.8

0.7

1
2
3
Kleinberg

0.6

0
1
2
3
4
5
6
7

0.5

0.4
0.5

probability

probability

0.6

0.4

0.3

0.3

0.2

0.2

0.1

0.1
0.5

1

1.5

2

2.5

3

time step n
図 15

異常区間の観測ステップ

=

3.5
×10

4

0
0.5

4

10000 における出現確率の絶対誤差の中央値に最も近い結果 (Kleinberg の手法)

1

1.5

time step n
図 17

2
×10

4

Oshima Tsubaki Camellia Hair Care Oil (Oshima Tsubaki) の結果

0.35

0

1

2

3

4

5

6

7

0.3

probability

0.25

0.2

0.15

0.1

0.05

0
0.5

1

1.5

time step n
図 18

Conditioner Essential (Albion) の結果

2×104