JP4769983B2 - 異常検出装置および異常検出方法 - Google Patents

Description

本発明は、異常検出装置および異常検出方法に関し、特に、任意の被検出対象から得られた任意の多次元特徴データから通常とは異なる異常状態を高い精度で検出することができる異常検出装置および異常検出方法に関するものである。

従来、異常状態の監視システムの一種として、セキュリティ分野における映像監視や老人介護のモニタリングシステムなど、カメラによる異常監視システムが多く利用されている。しかし、これらは人が監視を行うシステムであり、手間と費用がかかっていた。そこで、コンピュータによる異常状態の自動検出システムが求められていた。

異常状態の一種である、人や車など任意の対象物の異常動作の認識の研究としては、本発明者の一部らが発表した下記の特許文献１に、立体高次局所自己相関特徴（以下、CHLACデータとも記す）を用いた異常動作認識を行う技術が開示されている。特許文献１の異常動作検出装置は、ビデオカメラから入力される動画像データからフレーム間差分データを生成し、複数のフレーム間差分データからなる３次元データから立体高次局所自己相関によって特徴データを抽出し、過去の特徴データから主成分分析手法により得られた主成分ベクトルに基づく部分空間と、最新の特徴データとの距離を計算する。そして、この距離が所定値よりも大きい場合に異常と判定する。対象物の正常動作を部分空間として学習し、異常動作をそこからの逸脱として検出することにより、例えば画面内に複数人いる場合でも１人の異常動作を検出可能である。
特開２００６−０７９２７２号公報

上記した従来の異常検出方法においては、主成分分析手法により得られた主成分ベクトルに基づく正常領域を示す部分空間と最新の特徴データとの距離によって異常か否かを判定していた。しかし、主成分分析手法により得られた主成分ベクトルに基づく正常領域を示す部分空間の近傍や部分空間内にも本来異常と判定すべき特徴データが存在する可能性があるという問題点があった。

本発明の目的はこのような問題点を解決し、動作に関する特徴データに限らず、任意の被検出対象から得られた任意の多次元特徴データ（特徴ベクトル）から通常とは異なる異常状態を高い精度で検出することができる異常検出装置および異常検出方法を提供する点にある。

本発明の異常検出装置は、入力されたデータから多次元の特徴データを抽出する特徴データ抽出手段と、特徴データ空間において予め求められている正常空間を示す錐状体と前記特徴データ抽出手段により抽出された特徴データとの角度を計算する角度計算手段と、前記角度が所定値よりも大きい場合に異常と判定する異常判定手段とを備えたことを主要な特徴とする。

また、上記した異常検出装置において、前記正常空間を示す錐状体は、複数の基底ベクトルにより各辺が規定された多角錐である点にも特徴がある。
また、上記した異常検出装置において、前記多角錐は、ノンネガティブマトリックスファクトライゼーション（Non-negative Matrix Factorization）法により得られた基底ベクトルによって各辺が規定された、正常部分空間内における正常サンプルの分布範囲を近似する多角錐である点にも特徴がある。

また、上記した異常検出装置において、前記多角錐は、正常サンプルから主成分分析によって求めた主成分ベクトルに基づいて生成された複数の基底ベクトルによって各辺が規定された、正常部分空間内における正常サンプルの分布範囲を近似する多角錐である点にも特徴がある。

また、上記した異常検出装置において、前記正常空間を示す錐状体は、円錐の中心を示す方向ベクトル情報および円錐の広がりを示す角度情報によって規定された円錐である点にも特徴がある。また、上記した異常検出装置において、更に、特徴データのスケールを正規化するスケール正規化手段を備えた点にも特徴がある。

また、上記した異常検出装置において、更に、学習用の入力データに基づき、前記特徴データ抽出手段によって抽出された複数の特徴データから正常空間を示す錐状体を算出する学習手段を備えた点にも特徴がある。

また、上記した異常検出装置において、更に、特徴データを複数のクラスタに分けるクラスタ分割手段を備え、前記学習手段は分割された各クラスタ毎に正常空間を示す錐状体を算出し、前記異常判定手段は、前記各クラスタ毎に判定を行い、前記特徴データが全てのクラスタにおいて異常と判定された場合にのみ当該特徴データを異常と判定する点にも特徴がある。

また、上記した異常検出装置において、前記特徴データ抽出手段は、入力された複数の画像フレームデータからなる動画像データからフレーム間差分データを生成し、複数の前記フレーム間差分データからなる３次元データから立体高次局所自己相関によって動作特徴データを生成する点にも特徴がある。

本発明の異常検出方法は、入力されたデータから特徴データを抽出するステップ、特徴データ空間において予め求められている正常空間を示す錐状体と前記特徴データ抽出手段により抽出された特徴データとの角度を計算するステップ、前記角度が所定値よりも大きい場合に異常と判定するステップを含むことを主要な特徴とする。

本発明によれば、以下のような効果がある。
（１）従来異常であるのに正常と判定されていた事象についても異常の判定が可能となり、異常判定の精度が向上する。
（２）特徴抽出や異常判定のための計算量が少なく、計算量は対象数に依らず一定であり、実時間処理が可能である。

（３）さまざまな監視対象に対応可能であり、特徴ベクトルとしてはCHLACデータを始めとして任意のものを採用可能である。
（４）CHLACを採用して角度で異常判定を行う場合には、対象のスケールにロバストに異常を検出できる。
（５）異常検出に限らず一般の認識課題（顔、人認識など）にも適用可能である。

まず、本発明の実施例としては、一例としてCHLACデータを用いた異常動作の検出について開示するが、本発明は、CHLACデータに限らず、任意の被検出対象から得られた任意の多次元特徴データ（特徴ベクトル）から、異常動作に限らず、通常とは異なる異常状態、即ち発生確率の小さな事象を高い精度で検出することができるものである。

次に、「異常」の定義であるが、異常の事象を全て列挙することができないように、異常そのものを定義することはできない。そこで本明細書においては、異常を「正常ではないもの」として定義する。正常は、例えば特徴の統計的な分布を考えた場合に分布の集中する領域であるものとすれば、統計的な分布から学習可能である。そして、その分布から大きく逸脱する動作を異常とする。

例えば、通路や道路などに設置された防犯カメラでは、歩く動作のような一般的な動作は正常動作として学習され、認識されるが、挙動不審のような動作は、歩く動作のような周期運動でもなく、発現する確率の非常に低い事象であるため、異常として認識されることになる。

本発明においては、特徴ベクトル空間内において正常な特徴ベクトルが１つあるいは複数の（超立体の）多角錐形状あるいは円錐形状の領域を合わせた正常空間領域内に分布していることを見出し、特徴ベクトルが、上記した１つあるいは複数の多角錐あるいは円錐形状領域の内の１つの内部に存在しているか否かを判定し、いずれかの多角錐あるいは円錐形状の内部に存在していれば正常、そうでなければ異常と判定する。以下に本発明の実施例の一例としてCHLACデータを用いた異常動作の検出について説明する。

図１は、本発明の実施例１（多角錐）の判定方法を示す説明図である。図１は説明を簡単にするために、CHLAC特徴データ空間を３次元（実際には例えば２５１次元）としたものであり、正常動作空間は正常動作のCHLAC特徴データを含む１個の多角錐の内部に存在している。

図１（ｂ）は判定すべきサンプル特徴ベクトルｘおよび多角錐の各辺を表す基底ベクトルと交差する平面で切った場合の断面における配置例を示しており、図１（ｃ）は原点およびサンプル特徴ベクトルｘから一番近い多角錐の表面へ降ろした垂線を含む平面で切った場合の断面を示している。

特徴ベクトルｘから一番近い多角錐の表面が多角錐の辺（基底ベクトル）の間にある平面である場合（図１（ｂ）におけるＡの範囲）には、角度θはこの平面と特徴ベクトルｘとのなす角度となる。また、特徴ベクトルｘから一番近い多角錐の表面が辺（基底ベクトル）である場合（図１（ｂ）におけるＢの範囲）には、角度θはこの辺（基底ベクトル）と特徴ベクトルｘとのなす角度となる。

従って、角度θが同じ値となる特徴ベクトルの位置を結んだ線は図１（ｂ）に示す線Ｌのような位置となる。異常動作の特徴データは、正常動作空間を示す多角錐の表面との角度θが大きくなるので、この角度θが所定の閾値以上か否かによって異常を判定する。

図５は、本発明による異常検出装置の構成を示すブロック図である。ビデオカメラ１０は対象となる人や装置の動画像フレームデータをリアルタイムで出力する。ビデオカメラ１０はモノクロでもよいしカラーカメラであってもよい。コンピュータ１１は例えば動画像を取り込むためのビデオキャプチャー回路を備えた周知のパソコン（ＰＣ）であってもよい。本発明は、パソコンなどの周知の任意のコンピュータ１１に後述するプログラムを作成してインストールすることにより実現される。

モニタ装置１２はコンピュータ１１の周知の出力装置であり、例えば異常動作が検出されたことをオペレータに表示するために使用される。キーボード１３およびマウス１４は、オペレータが入力に使用する周知の入力装置である。なお、実施例においては、例えばビデオカメラ１０から入力された動画像データを実時間で処理してもよいし、一旦動画像ファイルに保存してから順次読み出して処理してもよい。

図６は、本発明の異常検出処理の検出時の処理内容を示すフローチャートである。Ｓ１０においては、ビデオカメラ１０からのフレームデータの入力が完了するまで待ち、フレームデータを入力（メモリに読み込む）する。このときの画像データは例えば２５６階調のグレイスケールデータである。

Ｓ１１においては、入力されたフレームデータから後述する方法によって特徴ベクトルを生成する。特徴ベクトルは例えば後述する２５１次元のCHLAC特徴データであってもよい。Ｓ１２においては、下記の数式１に示すように、特徴ベクトルについてベクトルノルムによる正規化を行う。即ち、特徴ベクトルｘをベクトルノルム（ベクトルの長さ、｜ｘ｜）で除算する。

Ｓ１３においては、未処理のクラスタを選択する。クラスタの数および各クラスタの基底ベクトルｘb、および閾値は予め後述する学習処理によって求められている。

Ｓ１４においては、以下に示す方法によって多角錐の表面と入力特徴ベクトルとの角度を求める。なお、上記した基底ベクトルを求めておけば処理が高速化されるが、以下の方法は基底ベクトルとして学習時に得られている全ての正常な特徴ベクトルをそのまま用いてもよい。基底ベクトルの集合をＸb、判定すべき正規化された入力ベクトルをｙとすると、凸多角錐（Ｓ）および多角錐への角度θは次の数式２で表される。

cosθをsinθに変換することによって数式２は数式３に示すように変形可能であり、数式２の最後の式を最小化するαを求めれば、Ｘ_bαが多角錐表面への射影ベクトルとなり、θが求まる。但し、入力ベクトルと多角錐内ベクトルの成す最小角度θがπ／２より大きいとき、つまり∀λ＞π／２のときは、最後の式の解はα＝０となる。そこで、得られた解がα＝０の場合は最小角度をθ＝π／２（上限）とする。

数式３の最後の式は周知の凸２次計画問題であるので、結局、凸２次計画問題を解くことによってθを求めることができる。コンピュータで凸２次計画問題を解くためには、例えばサイバネットシステム株式会社が販売している科学計算用ソフトウェアＭＡＴＬＡＢ（登録商標）システムのOptimization Toolboxのquadprog関数を使用することができる。なお、凸２次計画問題のアルゴリズムについては例えば下記の非特許文献１、２に記載されている。
Coleman, T.F. and Y. Li, OA Reflective Newton Method for Minimizing a Quadratic Function Subject to Bounds on some of the Variables,O SIAM Journal on Optimization, Vol. 6, Number 4, pp.1040-1058, 1996. Gill, P. E. and W. Murray, and M.H. Wright, Practical Optimization, Academic Press, London, UK,1981.

なお、角度θは以下の方法によっても求められる。数式３の下から２行目の式は周知のＮＮＬＳ（Non-negative Least Squares）問題であるので、ＮＮＬＳ問題を解くことによってθを求めることができる。

コンピュータでＮＮＬＳ問題を解くためには、例えばサイバネットシステム株式会社が販売している科学計算用ソフトウェアＭＡＴＬＡＢ（登録商標）システムのOptimization Toolboxのnnls関数を使用することができる。

なお、ＮＮＬＳ問題のアルゴリズムについては例えば下記の非特許文献３、４に記載されている。
Bro, R., de Jong, S., 1997. A fast non-negativity constrained linear least squares algorithm. Journal of Chemometrics 11, 393-401. Lawson, C., Hanson, R., 1995. Solving Least Squares Problems. SIAM, Philadelphia, PA.

Ｓ１５においては、算出された角度θが異常判定の閾値未満か否かが判定され、判定結果が否定の場合にはＳ１６に移行するが、肯定の場合にはＳ１７に移行する。閾値は学習時に決定する。Ｓ１６においては、全クラスタについて処理が完了したか否かが判定され、判定結果が否定の場合にはＳ１３に移行するが、肯定の場合にはＳ１８に移行する。

Ｓ１７においては、「正常」を出力し、またＳ１８においては、「異常」を出力する。Ｓ１９においては、処理の終了が指示されたか否かが判定され、判定結果が否定の場合にはＳ１０に移行して処理を繰り返す。以上のような検出時の処理によって入力されたサンプルデータが異常か否かが判定される。

図７は、本発明の実施例１の異常検出処理の学習時の処理内容を示すフローチャートである。この処理は図６の検出処理に先立って実行しておく。Ｓ３０においては、動作学習用のサンプルフレームデータを入力する。Ｓ３１においては、Ｓ１１と同様に特徴ベクトルを生成する。Ｓ３２においては、Ｓ１２と同様にベクトルノルムによる正規化を行う。Ｓ３３においては、特徴ベクトルをMean Shift（ミーンシフト）法によりクラスタリングする。

図１４は、Ｓ３３のクラスタリング処理の内容を示すフローチャートである。発明者らはvon Mises-Fisher分布を用いたMean Shift法を発明した。以下これを説明する。Ｓ１２０においては、未処理のサンプル点を１つ選択して参照点ξの初期値とする。Ｓ１２１においては、数式４に示す更新式に基づいて参照点ξの更新処理を行う。κはvon Mises-Fisher分布における集中度を表すパラメータであり、実験により設定する。

Ｓ１２２においては、ξに変化が無くなったか否か、即ち更新後のξと更新前のξとの差の絶対値が所定値以下か否かが判定され、判定結果が否定の場合にはＳ１２１に移行するが、肯定の場合にはＳ１２３に移行する。Ｓ１２３においては、全てのサンプル点について処理が完了したか否かが判定され、判定結果が否定の場合にはＳ１２０に移行するが、肯定の場合にはＳ１２４に移行する。以上の処理により、各サンプル点と対応する参照点ξが求まる。Ｓ１２４においては、近い位置にある点を統合する。即ち、参照点ξ相互の差が所定値以下である参照点を持つサンプルの集合を１つのクラスタと見なす。

図７に戻って、Ｓ３４においては、未処理のクラスタを１つ選択する。Ｓ３５においては、多角錐の基底ベクトルを例えば下記の方法によって求める。正常な特徴ベクトルの分布（最も外側の境界面）を示す多角錐は各頂点を表す特徴ベクトル＝基底ベクトルによって表現でき、多角錐の内部にある特徴ベクトルを取り除いても多角錐の形状は不変である。そこで、周知のリーブワンアウト法（Leave One Out Method）を用いて多角錐の内部に存在する（多角錐の頂点とならない）特徴ベクトルを取り除き、基底ベクトルを得る。

図２はリーブワンアウト法を示す説明図である。図２に示すようにサンプルの集合から非独立サンプルｘ_tを取り除いても、多角錐の形状は変わらない。従って以下に示す方法によって非独立サンプルｘ_tを取り除き、独立サンプルｘ_bのみを残す。

図１６は、Ｓ３５のリーブワンアウト法による基底ベクトル抽出処理の内容を示すフローチャートである。Ｎ個の特徴ベクトルからなるサンプル集合Ｘ＝ｘ₁,...,ｘ_Nがあったとする。Ｓ１５０においては、変数ｉを１にセットし、Ｓ１５１においては、凸２次計画問題またはＮＮＬＳ問題を解くことによってｘ_iを除いた集合Ｘ~によって形成される多角錐とベクトルｘ_iとの最小角度ψ_iを算出する。Ｓ１５２においては、変数ｉに１を加算し、Ｓ１５３においては、変数ｉがＮより大きいか否かが判定され、判定結果が否定の場合にはＳ１５１に移行するが、肯定の場合にはＳ１５４に移行する。

Ｓ１５４においては、変数ｉを１にセットし、Ｓ１５５においては、角度ψ_iが所定の閾値ｓ（正の小さな値）より大きいか否かが判定され、判定結果が否定の場合にはＳ１５７に移行するが、肯定の場合にはＳ１５６に移行する。Ｓ１５６においては、ベクトルｘ_iを基底ベクトルｘ_bとして登録する。なお、ここでψ_iがあまりに大きい値（閾値ｔ以上）のときには外れ値として基底ベクトルには登録しない。結局は、ｓ＜ψ_i＜ｔのものだけを登録する。
Ｓ１５７においては、変数ｉに１を加算し、Ｓ１５８においては、変数ｉがＮより大きいか否かが判定され、判定結果が否定の場合にはＳ１５５に移行するが、肯定の場合には処理を終了する。リーブワンアウト法を適用することによって基底ベクトルの数は元の特徴ベクトル数と較べると非常に少なくなる。

図７に戻って、Ｓ３６においては、全クラスタについて処理が完了したか否かが判定され、判定結果が否定の場合にはＳ３４に移行するが、肯定の場合には処理を終了する。以上のような学習処理によって、正常範囲を示す１個あるいは複数個の多角錐の基底ベクトルが求められる。

次に、処理の細部について説明する。図１２は、Ｓ１１の特徴ベクトル生成処理の内容を示すフローチャートである。Ｓ１００においては、差分データを生成し、２値化する。Ｓ１００においては動画像データに対して「動き」の情報を検出し、背景など静止しているものを除去する目的で差分データを生成する。差分の取り方としては、隣接するフレーム間の同じ位置の画素の輝度の変化を抽出するフレーム間差分方式を採用するが、フレーム内における輝度の変化部分を抽出するエッジ差分、あるいは両方を採用してもよい。

更に、「動き」に無関係な色情報やノイズを除去するために自動閾値選定による２値化を行う。２値化の方法としては、例えば一定閾値、下記非特許文献５に開示されている判別最小二乗自動閾値法、閾値０及びノイズ処理方式（濃淡画像において差が０以外を全て動き有り＝１とし、公知のノイズ除去法によってノイズを除去する方法）を採用可能である。以上の前処理により入力動画データは画素値に「動いた（１）」「動かない（０）」の論理値をもつフレームデータの列となる。
大津展之、"判別および最小２乗規準に基づく自動しきい値選定法"電子通信学会論文誌Ｄ、J63-D-4、P349-356、1980年。

Ｓ１０１においては、１フレーム分の画素データに関する相関パターンカウント処理を行い、フレームCHLACデータを生成する。詳細については後述するが、２５１次元の特徴ベクトルデータが生成される。時系列の二値差分データからの動作特徴の抽出には、立体高次局所自己相関(CHLAC) 特徴を用いる。Ｎ次自己相関関数は次の数式５のようになる。

ここで、f は時系列画素値（差分値）であり、参照点（注目画素）ｒおよび参照点からみたＮ個の変位ａ_i (i = 1，．．．,N) は差分フレーム内の二次元座標に時間も成分として持つ三次元のベクトルである。更に、時間方向の積分範囲はどの程度の時間方向の相関を取るかのパラメータである。但し、Ｓ１０１におけるフレームCHLACデータは１フレーム毎のデータであり、これを時間方向に所定期間だけ積分（加算）したものがCHLAC特徴データとなる。

変位方向、次数のとり方により高次自己相関関数は無数に考えられるが、これを局所領域に限定したものが高次局所自己相関関数である。立体高次局所自己相関特徴では変位方向を参照点ｒを中心とする３×３×３画素の局所領域内、即ち参照点ｒの２６近傍に限定している。特徴量の計算は１組の変位方向に対して式１の積分値が１つの特徴量になる。従って変位方向の組み合わせ（＝マスクパターン）の数だけ特徴量が生成される。

特徴量の数、つまり特徴ベクトルの次元はマスクパターンの種類に相当する。２値画像の場合、画素値１を何回乗算しても１であるので、二乗以上の項は乗数のみが異なる１乗の項と重複するものとして削除する。また数式５の積分操作（平行移動：スキャン）で重複するパターンは１つの代表パターンを残して他を削除する。式１右辺の式は参照点（f(r)：局所領域の中心）を必ず含むので、代表パターンとしては中心点を含み、パターン全体が３×３×３画素の局所領域内に収まるものを選択する。

この結果、中心点を含むマスクパターンの種類は、選択画素数が１個のもの：１個、２個のもの：２６個、３個のもの：２６×２５／２＝３２５個の計３５２個あるが、式１の積分操作（平行移動：スキャン）で重複するパターンを除くと、マスクパターンの種類は２５１種類となる。即ち、１つの３次元データに対する立体高次局所自己相関特徴ベクトルは２５１次元となる。

なお、画素の値が多値の濃淡画像の場合には、例えば画素値をａとすると、相関値はａ（０次）≠ａ×ａ（１次）≠ａ×ａ×ａ（２次）となり、選択画素が同じでも乗数の異なるものを重複削除できない。従って、多値の場合には、２値の場合より選択画素数が１の場合に２個、選択画素数が２の場合に２６個増加し、マスクパターンの種類は計２７９種類となる。

立体高次局所自己相関特徴の性質としては、変位方向を局所領域に限定しているためにデータに対する加法性があり、画面全体を通して積分しているためにデータ内の位置不変性がある。

Ｓ１０２においては、フレームCHLACデータをフレーム対応に保存する。Ｓ１０３においては、現在のCHLACデータにＳ１０１で求めた最新のフレームCHLACデータを加算すると共に、所定期間以上古いフレームに対応するフレームCHLACデータを現在のCHLACデータから減算して新たなCHLACデータを生成し、保存する。

図１９は、本発明による動画像の実時間処理の内容を示す説明図である。動画像のデータはフレームのシーケンスとなっている。そこで時間方向に一定幅の時間窓を設定し、窓内のフレーム集合を１つの３次元データとする。そして、新たなフレームが入力される度に時間窓を移動させ、古いフレームを削除することで有限な３次元データが得られる。この時間窓の長さは、認識しようとする動作の１周期と等しいか、より長く設定することが望ましい。

なお、実際には画像フレームデータは差分を取るために１フレームのみ保存され、フレームと対応するフレームCHLACデータが時間窓分だけ保存される。即ち、図１９において時刻ｔに新たなフレームが入力された時点では、直前の時間窓（t-1，t-n-1）に対応するフレームCHLACデータはすでに算出されている。但し、フレームCHLACデータを計算するためには直近の３つの差分フレームが必要であるが、（t-1）フレームは端なのでフレームCHLACデータは（t-2）フレームに対応するものまで算出されている。

そこで、新たに入力されたｔフレームを用いて（t-1）フレームに対応するフレームCHLACデータを生成してCHLACデータに加算する。また、最も古い（t-n-1）フレームと対応するフレームCHLACデータをCHLACデータから減算する。このような処理により、時間窓と対応するCHLAC特徴データが更新される。

図１３は、Ｓ１０１の立体高次局所自己相関特徴抽出処理の内容を示すフローチャートである。Ｓ１１０においては、２５１個の相関パターンカウンタをクリアする。Ｓ１１１においては、未処理の画素を１つ選択する（フレーム内で注目画素を順にスキャンする）。Ｓ１１２においては、未処理のマスクパターンを１つ選択する。

図１７は、３次元画素空間における自己相関処理座標を示す説明図である。図１７においては、ｔ−１フレーム、ｔフレーム、ｔ＋１フレームの時間的に隣接する３つの差分フレームのｘｙ平面を並べて図示してある。本発明においては、注目する参照画素を中心とする３×３×３（＝２７）画素の立方体の内部の画素について相関を取る。

マスクパターンは、相関を取る画素の組合せを示す情報であり、マスクパターンによって選択された画素のデータは相関値の計算に使用されるが、マスクパターンによって選択されなかった画素は無視される。前記したように、マスクパターンでは注目画素（中心の画素）は必ず選択される。また、２値画像で０次〜２次までの相関値を考えた場合、３×３×３画素の立方体において重複を排除した後のパターン数は２５１個となる。

図１８は、自己相関マスクパターンの例を示す説明図である。図１８（１）は注目画素のみの最も簡単な０次のマスクパターンである。（２）はハッチングを施した２つの画素が選択されている１次マスクパターン例、（３）、（４）はハッチングを施した３つの画素が選択されている３次マスクパターン例であり、この他に多数のパターンがある。

図１３に戻って、Ｓ１１３においては、前記した数式５を用いて相関値を計算する。数式５のf(r)f(r+a₁)…f(r+a_N)はマスクパターンと対応する座標の差分２値化３次元データの画素値を掛け合わせる（＝相関値、０または１）ことに相当する。また、数式５の積分操作は注目画素をフレーム内で移動（スキャン）させて相関値をマスクパターンと対応するカウンタによって足し合わせる（１をカウントする）ことに相当する。

Ｓ１１４においては、相関値は１か否かが判定され、判定結果が肯定の場合にはＳ１１５に移行するが、否定の場合にはＳ１１６に移行する。Ｓ１１５においては、マスクパターンと対応する相関パターンカウンタを＋１する。Ｓ１１６においては、全てのパターンについて処理が完了したか否かが判定され、判定結果が肯定の場合にはＳ１１７に移行するが、否定の場合にはＳ１１２に移行する。

Ｓ１１７においては、全ての画素について処理が完了したか否かが判定され、判定結果が肯定の場合にはＳ１１８に移行するが、否定の場合にはＳ１１１に移行する。Ｓ１１８においては、パターンカウンタ値の集合を２５１次元フレームCHLACデータとして出力する。

図３は、本発明の実施例２（近似的多角錐１）の判定方法を示す説明図である。図３は説明を簡単にするために、CHLAC特徴データ空間を３次元（実際には例えば２５１次元）としている。実施例１においては、正常な特徴ベクトルの分布範囲を包含する多角錐を求めたが、実施例２においては、学習時に正常な特徴ベクトルの分布範囲を近似する多角錐を求め、この近似的多角錐を用いて、実施例１と同様に異常か否かを判定するものである。近似的多角錐の求め方について以下に実施例を説明する。

＜近似的多角錐の求め方１＞
図８は、近似的多角錐の求め方１を用いた本発明の異常検出処理の学習時の処理内容を示すフローチャートである。図７に示した実施例１のフローチャートとの相違点はＳ３５とＳ４５が相違しているのみであり、他のステップの処理は同一であるので、Ｓ４５の処理について説明する。近似的多角錐の求め方１においては、正常な特徴ベクトルの分布範囲を近似し、かつ基底ベクトルの数が少ない多角錐を求める。角数は任意である。但し、近似的多角錐１では特徴ベクトルの要素が全て非負であることを前提としている。例えばCHLAC特徴データなどではこの非負の条件を満たしており、この条件を満たさないデータに対しては分布が非負領域になるように回転させるなどの前処理により対処できる可能性がある。なお、本発明における近似的多角錐１以外の手法ではこのような前提を置いておらず、任意の特徴ベクトルに適用可能である。

具体的には、例えば下記の数式６の式を下記の非特許文献６に開示されているような公知のノンネガティブマトリックスファクトライゼーション（Non-negative Matrix Factorization）法によって解き、近似的多角錐の基底ベクトルを求める。なお、Ｘは多角錐を構成している特徴ベクトル群、Ｗは近似的多角錐の基底ベクトル、Ｈは基底ベクトルの結合係数である。
Lee, D., & Seung, H. (1999). Learning the parts of objects by non-negative matrix factorization. Nature, 401, 788-791.

具体的には下記の数式７に示すような処理を行う。即ち、（１）ＷおよびＨの初期値をランダムに設定する。（２）数式７に示すような式に従ってＷ、Ｈを更新する。（３）もし、（３）に示す判定式が収束した（値の変化が更新の前後で所定値以下であった）場合には処理を終了し、Ｗが近似的多角錐の基底ベクトルとなるが、そうでなければ（３）に示す更新式によってＷ、Ｈを正規化し、（２）に移行する。

学習時には以上の処理によって近似的多角錐を求め、検出時には近似的多角錐を用いて実施例１と同様に図６に示す処理によって異常の判定を行う。

図３（ｂ）、（ｃ）は、本発明の実施例３（近似的多角錐２）の判定方法を示す説明図である。図３では説明を簡単にするために、CHLAC特徴データ空間を３次元（実際には例えば２５１次元）としている。実施例３においては、前述した実施例２と同様に近似的多角錐を求めるが、その求め方が異なり、得られた近似的多角錐は正常な特徴ベクトルの分布範囲を近似するものとなる。近似的角錐の求め方について以下に実施例を説明する。

＜近似的多角錐の求め方２＞
図３（ｂ）は、実施例３の判定方法を示す説明図である。実施例３（近似的多角錐の求め方２）においては公知の主成分分析によって多角錐の基底ベクトルを算出する。主成分分析の結果得られる第１固有ベクトルは多角錐の中心方向を示し、第２固有ベクトル以降は第１固有ベクトルと直角な多角錐の広がり方向を示している。そこで、第２固有ベクトル以降の固有ベクトルの軸上における特徴ベクトルの分布から近似的多角錐の基底ベクトルを算出する。

図９は、近似的多角錐の求め方２を用いた本発明の異常検出処理の学習時の処理内容を示すフローチャートである。図７に示した実施例１のフローチャートととの相違点はＳ３５がＳ５５、Ｓ５６に置換されているのみであるので、Ｓ５５、Ｓ５６の処理について説明する。

Ｓ５５においては、主成分分析によって主成分ベクトルを出力する。主成分分析手法自体は周知であるので概略を説明する。まず、多角錐を含む部分空間を構成するために、CHLAC特徴データ群から主成分分析により主成分ベクトルを求める。Ｍ次元のCHLAC特徴ベクトルｘを以下のように表す。

なお、Ｍ＝２５１である。また、主成分ベクトルを列に並べた行列Ｕ（固有ベクトル）を以下の数式９のように表す。

なお、Ｍ＝２５１である。主成分ベクトルを列に並べた行列Ｕは、以下のように求める。共分散行列Ｒを次の数式１０に示す。

行列Ｕはこの共分散行列Ｒを用いて、次の数式１１の固有値問題より求まる。

固有値の対角行列Λを次数式１２で表す。

第Ｋ固有値までの累積寄与率β_Kは、以下の数式１３ように表される。

ここで、累積寄与率β_Kが所定値（例えばβ_K＝0.99）となる次元までの固有ベクトルu₁,...,u_Kを用いて多角錐の基底ベクトルを算出する。なお、累積寄与率β_Kの最適値は監視対象や検出精度にも依存するものと考えられるので実験等により決定する。

Ｓ５６においては、求めた主成分ベクトルから近似的多角錐の基底ベクトルを算出する。以下に主成分ベクトルから近似的多角錐の基底ベクトルを算出する方法を２つ説明する。

＜近似的多角錐の求め方２Ａ＞
図３（ｃ）は、近似的多角錐の求め方２の第１の算出方法を示す説明図である。第１の算出方法は、第ｉ固有ベクトル（ｉ＞１）に射影した特徴ベクトルの分布の中心から標準偏差σのｋ倍の距離の点を基底ベクトルとする。即ち、第（２ｉ−３）番目の基底ベクトル＝ｅ₁＋ｋσｅ_i、第（２ｉ−２）番目の基底ベクトル＝ｅ₁−ｋσｅ_iとする。なお、ｅ_iは第ｉ固有ベクトルでり、ｋは例えば実験によって最適値が決定される任意の実数であるが、例えばｋ＝３であってもよい。

この方法によれば、主成分分析による次元の２倍の数の基底ベクトルが得られる。また、標準偏差σのｋ倍の距離の点を基底ベクトルとすることにより、ノイズやはずれ値サンプル等による影響が軽減される。

＜近似的多角錐の求め方２Ｂ＞
上述した近似的多角錐の求め方２Ａと同様に主成分分析を行うが、基底ベクトルは、例えば第ｉ固有ベクトル（ｉ＞１）に射影した特徴ベクトルの分布から手動で正常範囲を設定する。即ち第（２ｉ−３）番目の基底ベクトル＝ｅ₁＋Ｌ₁ｅ_i、第（２ｉ−２）番目の基底ベクトル＝ｅ₁−Ｌ₂ｅ_iとする。なお、ｅ_iは第ｉ固有ベクトルである。また、Ｌ₁、Ｌ₂は手動で設定したパラメータであり、例えば、各軸に射影されたサンプルの最大値、最小値であってもよい。この方法でも、主成分分析による次元の２倍の数の基底ベクトルが得られる。

学習時には以上の処理によって近似的多角錐を求め、検出時には近似的多角錐を用いて実施例１と同様に図７に示す処理によって異常の判定を行う。

図４は、本発明の実施例４（円錐）の判定方法を示す説明図である。図４は説明を簡単にするために、CHLAC特徴データ空間を３次元（実際には例えば２５１次元）としている。実施例４においては、正常な特徴ベクトルの分布範囲を包含する円錐（方向ベクトルγおよび広がり角φ）を求め、この円錐とサンプル特徴データｘとの角度θを求めて異常か否かを判定するものである。以下に実施例を説明する。

図１０は、本発明の異常検出処理の実施例４の検出時の処理内容を示すフローチャートである。図６に示した実施例１のフローチャートととの相違点はＳ１４がＳ６４からＳ６６に置換されているのみであるので、相違点の処理について説明する。

Ｓ６４においては、正規化特徴ベクトルｙが原点を中心とする球面上で等方的な分布を成すように特徴ベクトルの各要素のスケールを調整する。Ｕ、Λを後述の学習処理における主成分分析により得られた固有ベクトルおよび固有値対角行列とするとスケーリング処理後のベクトル~ｚは数式１４によって表される。

Ｓ６５においてはＳ１２と同様の処理によって、ベクトル~ｚはベクトルノルムにより正規化されてｚとなる（||ｚ||＝１）。Ｓ６６においては、数式１５に示すような、特徴ベクトルｚと方向ベクトルγとの角度を算出し、φを減算したものを錐体（円錐）への角度θとする。

図１１は、本発明の異常検出処理の実施例４の学習時の処理内容を示すフローチャートである。図７に示した実施例１のフローチャートとの相違点はＳ８５〜Ｓ８８のみであるので、相違点の処理について説明する。Ｓ８５においては、主成分分析を行い、固有ベクトルＵおよび固有値対角行列Λを求める。

Ｓ８６においては、Ｓ６４と同様の処理によって特徴ベクトルの各要素のスケールを調整する。Ｓ８７においては、Ｓ６５と同様の処理によってベクトルノルムによる正規化を行う。Ｓ８８においては、1-class SVM（1-class Support Vector Machine）法により、方向ベクトルγおよび広がり角度φを求める。

図１５は、Ｓ８８の方向ベクトル算出処理の内容を示すフローチャートである。1-class SVM法自体は周知であるので概略を説明する。図４（ｃ）において、方向ベクトルγおよび広がり角度φで規定される（超）球面上の円は（超）球と（超）平面の交線として表現できる。従って、サンプルの射影点（ノルム１に正規化されている）と原点とを分離し、かつ原点からの距離が最大の（超）平面ｂを求めれば、この超平面ｂへの法線ベクトルおよびこの（超）平面の原点からの距離から、方向ベクトルγ、角度の範囲φが求まる。即ち、これは1-class SVMの問題と等価になり、以下の数式１６に示す1-class SVMの問題を解けばよいことになる。数式１６は1-class SVMの問題において、最終的に解くべき凸２次計画問題である。

なお、νはユーザが指定するパラメータであり、例えばν＝０．０１であってもよい。上記問題は前述の通り、サイバネットシステム社のMATLABを用いても解くことができるが、例えばここでは具体的には以下のようにも解ける。Ｓ１３０においては、変数ｔ＝１とし、またδを数式１７の条件を満足する範囲でランダムに初期化する。

Ｓ１３１においては、数式１８に示すように、予備変数Ｏ、ＳＶ、ρを計算する。

Ｓ１３２においては、数式１９の条件を満たすｉ∈ＳＶを選択する。

Ｓ１３３においては、ｉはあったか否かが判定され、判定結果が否定の場合にはＳ１３４に移行するが、肯定の場合にはＳ１３６に移行する。Ｓ１３４においては、数式１９の条件を満たすｉ∈1,..,Nを求める。Ｓ１３５においては、ｉはあったか否かが判定され、判定結果が否定の場合には処理を終了するが、肯定の場合にはＳ１３６に移行する。Ｓ１３６においては、数式２０に基づいてｊを求める。

Ｓ１３７においては、数式２１に基づいてδを更新する。

Ｓ１３８においては、δは数式２２に示す条件を満足しているか否かが判定され、判定結果が否定の場合にはＳ１３９に移行するが、肯定の場合にはＳ１４０に移行する。

Ｓ１３９においては、所定の範囲［０，１／νＮ］内に入るようにδ_jを射影する。即ち、もしδ＞1/νNであれば、δ＝1/νNとし、もしδ＜０ならば、δ＝０として、最も近い範囲内の点（端点）に移動する。そして、δ_iを再計算する。Ｓ１４０においては、ｔ＝ｔ＋１としてＳ１３１に移行する。以上の処理によって求められたδ、ρを用いて法線ベクトルγおよびφは以下の数式２３によって求まる。

図１１に戻って、Ｓ８９においては、全クラスタについて処理が完了したか否かが判定され、判定結果が否定の場合にはＳ８４に移行するが、肯定の場合には処理を終了する。実施例４においては、検出時の計算負荷が軽く、高速に処理可能である。

以上、異常動作の検出を行う実施例について説明したが、本発明には以下のような変形例も考えられる。実施例においては、学習フェーズによって予め正常動作領域を示す多角錐あるいは円錐の情報を生成しておく例を開示したが、多角錐あるいは円錐の情報を更新しながら異常動作の検出を行うようにすれば、環境の変化（正常状態の変化）にも追随することが可能である。

更に、例えば１分、１時間あるいは１日などのフレーム間隔よりも長い所定の周期で、一部のデータを使用して多角錐あるいは円錐の情報を更新し、次の更新までは固定した多角錐あるいは円錐の情報を使用して異常動作の検出を行ってもよい。このようにすれば学習における処理量が減少する。

本発明の実施例１（多角錐）の判定方法を示す説明図である。 リーブワンアウト法を示す説明図である。 本発明の実施例２（近似的多角錐１）の判定方法を示す説明図である。 本発明の実施例４（円錐）の判定方法を示す説明図である。 本発明による異常検出装置の構成を示すブロック図である。 本発明の異常検出処理の検出時の処理内容を示すフローチャートである。 本発明の実施例１の異常検出処理の学習時の処理内容を示すフローチャートである。 近似的多角錐の求め方１を用いた本発明の異常検出処理の学習時の処理内容を示すフローチャートである。 近似的多角錐の求め方２を用いた本発明の異常検出処理の学習時の処理内容を示すフローチャートである。 本発明の異常検出処理の実施例４の検出時の処理内容を示すフローチャートである。 本発明の異常検出処理の実施例４の学習時の処理内容を示すフローチャートである。 Ｓ１１の特徴ベクトル生成処理の内容を示すフローチャートである。 Ｓ１０１の立体高次局所自己相関特徴抽出処理の内容を示すフローチャートである。 Ｓ３３のクラスタリング処理の内容を示すフローチャートである。 Ｓ８８の方向ベクトル算出処理の内容を示すフローチャートである。 Ｓ３５のリーブワンアウト法による基底ベクトル抽出処理の内容を示すフローチャートである。 ３次元画素空間における自己相関処理座標を示す説明図である。 自己相関マスクパターンの例を示す説明図である。 本発明による動画像の実時間処理の内容を示す説明図である。

符号の説明

ｘ…サンプルベクトル
θ…異常判定指標である角度
１０…ビデオカメラ
１１…コンピュータ
１２…モニタ装置
１３…キーボード
１４…マウス

Claims (10)

1. 入力されたデータから多次元の特徴データを抽出する特徴データ抽出手段と、
   特徴データ空間において、予め求められている正常空間を示す錐状体における前記特徴データから一番近い表面と、前記特徴データ抽出手段により抽出された特徴データとの角度を計算する角度計算手段と、
   前記角度が所定値よりも大きい場合に異常と判定する異常判定手段と、
   を備えたことを特徴とする異常検出装置。
2. 前記正常空間を示す錐状体は、複数の基底ベクトルにより各辺が規定された多角錐であることを特徴とする請求項１に記載された異常検出装置。
3. 前記多角錐は、ノンネガティブマトリックスファクトライゼーション（Non-negative Matrix Factorization）法により得られた基底ベクトルによって各辺が規定された多角錐であることを特徴とする請求項２に記載された異常検出装置。
4. 前記多角錐は、正常サンプルから主成分分析によって求めた主成分ベクトルに基づいて生成された複数の基底ベクトルによって各辺が規定された、正常部分空間内における正常サンプルの分布範囲を近似する多角錐であることを特徴とする請求項２に記載された異常検出装置。
5. 前記正常空間を示す錐状体は、円錐の中心を示す方向ベクトル情報および円錐の広がりを示す角度情報によって規定された円錐であることを特徴とする請求項１に記載された異常検出装置。
6. 更に、特徴データのスケールを正規化するスケール正規化手段を備えたことを特徴とする請求項５に記載の異常検出装置。
7. 更に、学習用の入力データに基づき、前記特徴データ抽出手段によって抽出された複数の特徴データから正常空間を示す錐状体を算出する学習手段を備えたことを特徴とする請求項１に記載された異常検出装置。
8. 更に、特徴データを複数のクラスタに分けるクラスタ分割手段を備え、
   前記学習手段は分割された各クラスタ毎に正常空間を示す錐状体を算出し、
   前記異常判定手段は、前記各クラスタ毎に判定を行い、前記特徴データが全てのクラスタにおいて異常と判定された場合にのみ当該特徴データを異常と判定することを特徴とする請求項１に記載の異常検出装置。
9. 前記特徴データ抽出手段は、入力された複数の画像フレームデータからなる動画像データからフレーム間差分データを生成し、複数の前記フレーム間差分データからなる３次元データから立体高次局所自己相関によって動作特徴データを生成することを特徴とする請求項１に記載の異常検出装置。
10. 入力されたデータから特徴データを抽出するステップ、
    特徴データ空間において、予め求められている正常空間を示す錐状体における前記特徴データから一番近い表面と、前記特徴データ抽出手段により抽出された特徴データとの角度を計算するステップ
    前記角度が所定値よりも大きい場合に異常と判定するステップ
    を含むことを特徴とする異常検出方法。

Images