(57)【要約】【課題】監視対象の状態異常を速やかに検知することができる異常検知システムおよびプログラムを提供する。【解決手段】監視対象１に対する観測データを収集して時系列観測データとして保存する観測データ更新部３１と、観測データを、訓練用データおよび検証用データのいずれかに分類する観測データ分類部３３と、訓練用データに基づき、監視対象の線形状態空間モデルのモデルパラメータを同定するシステム同定部３５と、モデルパラメータ、および検証用データを入力として、監視対象の状態変数の確率分布の推定値を計算する推定部３６と、推定値に基づき、監視対象の異常度を計算する異常度計算部３７とを備え、観測データ更新部３１は、観測データを収集すると、新しく収集した観測データを時系列観測データに追加するとともに、時系列観測データのデータ数がしきい値よりも大きい場合に、最も前に収集した観測データを破棄する。【選択図】図１
2019191836.tif 000002
【特許請求の範囲】
【請求項１】
線形状態空間モデルによりモデル化が可能な監視対象の異常を検知する異常検知システムであって、
前記監視対象に対する観測データを収集して時系列観測データとしてデータベースに保存する観測データ更新部と、
前記時系列観測データの少なくとも一部の観測データを、訓練用データおよび検証用データのいずれかに分類する観測データ分類部と、
前記訓練用データに基づき、前記監視対象の線形状態空間モデルのモデルパラメータを同定するシステム同定部と、
前記システム同定部により同定されるモデルパラメータ、および前記検証用データを入力として、前記監視対象の状態変数の確率分布の推定値を計算する推定部と、
前記推定部が計算した前記推定値に基づき、前記監視対象の異常度を計算する異常度計算部とを備え、
前記観測データ更新部は、監視対象の観測データを収集すると、新しく収集した観測データを前記時系列観測データに追加するとともに、前記時系列観測データのデータ数がしきい値よりも大きい場合に、前記時系列観測データにおける最も前に収集した観測データを破棄することを特徴とする、異常検知システム。
【請求項２】
前記推定部は、カルマンフィルタに基づいて前記確率分布の推定値を計算する、請求項１に記載の異常検知システム。
【請求項３】
コンピュータを請求項１または２に記載の異常検知システムとして機能させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、線形状態空間モデルによってモデル化が可能な監視対象に関する異常検知システムおよびプログラムに関するものである。
【背景技術】
【０００２】
近年、電気設備および機械設備には、状態監視機能が備えられ、いわゆるＣＢＭ(Condition Based Maintenance：状態基準保全)実現への取組みが進んでいる。
【０００３】
一般的な異常検知システムの原理については、非特許文献１において統一的に記述されており、特に時系列データの異常検知についての高度なアプローチのひとつとして、状態空間モデルによる異常検知方法が論じられている。
【０００４】
状態空間モデルによる異常検知方法とは、監視対象を状態空間モデルで表現し、当該監視対象に対する観測データから監視対象の内部状態を推定するとともに、内部状態の推定値に基づく異常度を算出するものである。
【０００５】
図５は、非特許文献１において、「状態空間モデルによる異常検知」として記載された異常検知手法に基づく異常検知システムの構成を表したものである。
【０００６】
符号１は監視対象を表す。符号２は監視対象１に対する観測データを一定周期で収集記録する観測データ収集部を表す。符号２１は観測データ収集部により生成される時系列観測データＤが保存されるデータベースを表す。符号２２と符号２３とは、それぞれ時系列観測データＤを分割して生成される訓練用データＤｔｒと検証用データＤｖｒが保存されるデータベースを表す。
【０００７】
符号３は、従来の異常検知システムを表す。異常検知システム３は、システム同定部３５、カルマンフィルタ推定部３６および異常度計算部３７を備える。
【０００８】
監視対象１の時刻ｔにおける内部状態（状態変数）はｚ（ｔ）というｍ次元ベクトルで表される。観測可能な状態量（以下「観測データ」と呼ぶ）はｘ（ｔ）というＭ次元ベクトルで表される。このとき、監視対象に対応する線形状態空間モデルは、監視対象に対する入力が無い場合に、下記の式で表すことができる。
【０００９】
【数１】
2019191836.tif 000003
【００１０】
2019191836.tif 000004
【００１１】
図６は、非特許文献１に掲載された異常検知アルゴリズムをフローチャートの形で示したものである。
【００１２】
この異常検知アルゴリズムによれば、十分なサンプリング数（データ数）から成る監視対象１の時系列観測データＤが最初に用意され、観測データ収集部２は、時系列観測データＤを訓練用データＤｔｒと検証用データＤｖｒとに分割する（Ｓ６１）。
【００１３】
ここで、訓練用データＤｔｒのデータ数をＴ１、検証用データＤｖｒのデータ数をＴ２とし、時系列観測データＤのデータ数をＴ１＋Ｔ２とする。
【００１４】
次にシステム同定部３５は、上記の訓練用データＤｔｒを用いて、監視対象１のモデルパラメータＡ、Ｃ、Ｑ、Ｒ、Ｑ０、ｚ０を、部分空間同定法などのシステム同定法を用いて求める（ステップＳ６２）。この処理の詳細については後述する。ここで、Ｑは内部状態の遷移のばらつきを表すｍ×ｍの共分散行列である。Ｒは観測データのばらつきを表すＭ×Ｍの共分散行列である。Ｑ０は時刻０でのＱである。ｚ０は時刻０でのｚ（ｔ）である。
【００１５】
次に、システム同定部３５は、式Ａμ０＝ｚ０により、μ０を得る（ステップＳ６３）。
【００１６】
カルマンフィルタ推定部３６は、検証用データＤｖｒに含まれる観測値ｘ（ｔ）（ｔ＝１，．．．，Ｔ２）について、非特許文献１で定義されるカルマンフィルタによる状態推定演算を実行する。まずは、ｔ＝１とする（ステップＳ６４）。カルマンフィルタ推定部３６は、ｔがＴ２以下であるか判定する（ステップＳ６５）。ここでは、ｔがＴ２以下であるため、次のステップＳ６６に進む。
【００１７】
カルマンフィルタ推定部３６は、下記の式でｔ＝１として、Ｋ１、μ１、Ｖ１およびＱ１を求める（ステップＳ６６）。
【００１８】
【数２】
2019191836.tif 000005
【００１９】
カルマンフィルタ推定部３６は、式Ａμ１＝ｚ１から、ｚ１を求める。
【００２０】
更に異常度計算部３７は、数３により定義した異常度ａ（ｘ（ｔ））を算出することによって、監視対象１の異常状態を定量的に評価する（ステップＳ６７）。
【００２１】
【数３】
2019191836.tif 000006
なお上述した式で、記号「≡」は「定義」を表している。
【００２２】
さらに、カルマンフィルタ推定部３６は、ｔを１増加させる（ステップＳ６８）。
【００２３】
ここまで説明したステップＳ６５からＳ６８を繰り返すことで、カルマンフィルタ推定部３６は、各時刻ｔ＝１，．．．，Ｔ２における状態変数ｚ（ｔ）の確率分布（平均と分散（μｔ、Ｖｔ）)の推定値を計算することができる。また異常度計算部３７は、カルマンフィルタ推定部３６が計算した推定値に基づき、監視対象１の異常度を計算することができる。
【００２４】
図７は、非特許文献１で説明される部分空間同定法によるシステム同定（図６のステップＳ６２）をより詳細なフローチャートの形で示したものである。このシステム同定は、異常検知システム３において、システム同定部３５によって実装されるものである。
【００２５】
ここで、長さＴのＭ次元時系列観測データＤが与えられているとする。また窓幅ｗと部分時系列の数Ｎについて、ＮがＭｗより十分大きくなるように選ばれているとする。そして時系列のおおよその中間点を時刻ｔに設定する。
【００２６】
システム同定部３５は、次式で示す行列Ｘ（ｔ）、ＸｐおよびＸｆを求める（ステップＳ７１）。
【００２７】
【数４】
2019191836.tif 000007
【００２８】
システム同定部３５は、次式で定義する行列Ｓｐｆ、Ｓｐｐ、ＳｆｐおよびＳｆｆを求める（ステップＳ７２）。
【００２９】
【数５】
2019191836.tif 000008
【００３０】
システム同定部３５は、次式で定義する行列Ｗを求める。そしてＷに対して特異値が大きい順にｍ個の、長さが１に規格化された左右特異ベクトルを求める（ステップＳ７３）。
【００３１】
【数６】
2019191836.tif 000009
【００３２】
2019191836.tif 000010
【００３３】
【数７】
2019191836.tif 000011
【００３４】
システム同定部３５は、下記の式を用いて、αｉおよびβｉから状態変数系列の推定値を求める（ステップＳ７５）。
【００３５】
【数８】
2019191836.tif 000012
【００３６】
システム同定部３５は、下記の式で定義される、観測データからなるデータＤ’と状態変数の推定値Ｚ’とを求める。なお、Ｔ’＝Ｔ−ｗとする。
【００３７】
【数９】
2019191836.tif 000013
【００３８】
システム同定部３５は、下記の式で定義するベクトルＸ、Ｚ、Ｚ＋およびＺ−を求める（ステップＳ７７）。
【００３９】
【数１０】
2019191836.tif 000014
【００４０】
システム同定部３５は、モデルパラメータＡ、Ｑ、ＣおよびＲの推定値を下記の式から求める。
【００４１】
【数１１】
2019191836.tif 000015
【００４２】
ここまで説明したように、図７に示す部分空間同定法によるシステム同定によれば、線形状態空間モデルのパラメータＡ，Ｃ，Ｑ，Ｒ、Ｑ０、ｚ０の推定値を、時系列観測データｘ（ｔ）から求めることができる。
【先行技術文献】
【非特許文献】
【００４３】
【非特許文献１】井出剛、「入門機械学習による異常検知」、コロナ社、２０１５年
【発明の概要】
【発明が解決しようとする課題】
【００４４】
上述した従来の異常検知アルゴリズムは、予め用意された時系列観測データを、訓練用データと検証用データとに分割し、監視対象の状態空間モデルに係るモデルパラメータの推定を、オフラインで一括処理する必要があった。この一括処理は例えば、時系列観測データが一定量蓄積された後に行われていた。
【００４５】
したがって、例えば経年劣化によって監視対象の特性に変動が生じた場合に、言い換えれば状態空間モデルのモデルパラメータに変化が生じた場合に、監視対象のシステム同定をやり直す必要があった。
【００４６】
そこで本発明は、監視対象の状態異常を速やかに検知することができる異常検知システムおよびプログラムを提供することを目的とするものである。
【課題を解決するための手段】
【００４７】
本発明は、非特許文献１において示された状態空間モデルによる異常検知アルゴリズムを発展させ、同アルゴリズムのオンライン化を図ったものである。
【００４８】
第１の観点による異常検知システムは、線形状態空間モデルによりモデル化が可能な監視対象の異常を検知する異常検知システムであって、前記監視対象に対する観測データを収集して時系列観測データとしてデータベースに保存する観測データ更新部と、前記時系列観測データの少なくとも一部の観測データを、訓練用データおよび検証用データのいずれかに分類する観測データ分類部と、前記訓練用データに基づき、前記監視対象の線形状態空間モデルのモデルパラメータを同定するシステム同定部と、前記システム同定部により同定されるモデルパラメータ、および前記検証用データを入力として、前記監視対象の状態変数の確率分布の推定値を計算する推定部と、前記推定部が計算した前記推定値に基づき、前記監視対象の異常度を計算する異常度計算部とを備え、前記観測データ更新部は、監視対象の観測データを収集すると、新しく収集した観測データを前記時系列観測データに追加するとともに、前記時系列観測データのデータ数がしきい値よりも大きい場合に、前記時系列観測データにおける最も前に収集した観測データを破棄することを特徴とするものである。
【００４９】
また、第２の観点による異常検知システムにおいては、前記複数の管理単位それぞれで有する、電力供給機器の供給余力に基づいて、前記管理単位それぞれの達成困難度を求めることが好ましい。
【００５０】
上述したように本発明の解決手段を異常検知システムとして説明してきたが、本発明はこれらに実質的に相当する装置、プログラムまたは当該プログラムを記録した記憶媒体としても実現し得るものであり、本発明の範囲にはこれらも包含されるものと理解されたい。
【００５１】
例えば、第３の観点によるプログラムは、第１または第２の観点による異常検知システムにおいて、線形状態空間モデルによりモデル化が可能な監視対象の異常を検知するためのプログラムである。
【発明の効果】
【００５２】
本発明によれば、監視対象の特性に変動が生じた場合に、監視対象のモデルパラメータの再学習処理のため、異常検知システムの動作を停止しなければならないという制約をなくすことができる。
【００５３】
そして本発明の異常検知システムは、観測データを収集した後、監視対象の内部状態の推定値に基づく異常度を速やかに算出することができる。
【図面の簡単な説明】
【００５４】
【図１】本発明に係る異常検知システムの実施形態の構成を示す図である。
【図２】本実施形態の観測データのサンプリング処理を示す図である。
【図３】本発明に係る異常検知システムの実施形態のアルゴリズムを表す図である。
【図４Ａ】本実施形態の異常検知システムの実施例を表す図である。
【図４Ｂ】或る時点での各ポインタの指示位置と、新しい観測データがサンプリングされた後の各ポインタの指示位置とを表す図である。
【図５】従来の異常検知システムの構成を示す図である。
【図６】従来の異常検知システムのアルゴリズムを表す図である。
【図７】部分空間同定法の学習アルゴリズムを示す図である。
【発明を実施するための形態】
【００５５】
以下、本発明の実施形態について、図面を参照して説明する。図１は、異常検知システムの実施形態の構成を表したものである。符号１は、監視対象を表す。符号３は、本実施形態の異常検知システムを表す。
【００５６】
本実施形態の異常検知システム３は、観測データ更新部３１と、観測データが保存されるデータベース３２と、観測データを訓練用データおよび検証用データのいずれかに分類する観測データ分類部３３と、訓練用データおよび検証用データがそれぞれ保存されるデータベース３４ａおよび３４ｂと、システム同定部３５と、カルマンフィルタ推定部３６と異常度計算部３７とを備える。
【００５７】
図２は、監視対象１の観測量をサンプリングして、観測データｘ（ｔ）として保存する過程を示す図である。任意の絶対的な時刻ｋにおける監視対象１の観測データをξ（ｋ）で表すと、上記の線形状態空間モデルの観測データｘ（ｔ）は、次式で定義される。
【００５８】
【数１２】
2019191836.tif 000016
【００５９】
以下、異常検知システム３が行う異常検出のアルゴリズムを、図３を参照して説明する。監視対象１の新たな観測データξ（ｋ＋Ｔ１＋Ｔ２）が得られたとき（ステップＳ３１）、観測データｘ（ｔ）の時系列観測データを、次式により更新する（ステップＳ３２）。
【００６０】
【数１３】
2019191836.tif 000017
【００６１】
上述した処理により、時系列観測データＤ０がＤ１に更新される。図２に示すように、更新後の時系列観測データＤ１では、更新前の時系列観測データの中で最も古い観測データ、すなわち最も前に収集された観測データであったξ（ｋ）が破棄されている。観測データ更新部３１は、このデータ破棄処理を、時系列観測データのデータ数がしきい値（この実施形態ではＴ１＋Ｔ２）よりも大きい場合に行う。更新前の時系列観測データＤ０と更新後の時系列観測データＤ１とのデータ数は同じである。
【００６２】
次に、観測データ分類部３３は、時系列観測データＤ１を訓練用データＤｔｒと検証用データＤｖｒに分ける（ステップＳ３３）。なお、時系列観測データＤ１のすべてを訓練用データＤｔｒと検証用データＤｖｒのいずれかに分けることができ、また時系列観測データＤ１の一部を訓練用データＤｔｒと検証用データＤｖｒのいずれかに分けることもできる。
【００６３】
その後、システム同定部３５は、上述し図７に示したフローチャートで詳細に示すシステム同定を行い、新たなモデルパラメータＡ、Ｃ、Ｑ、Ｒ、Ｑ０、ｚ０を得る（ステップＳ３４）。
【００６４】
次に、システム同定部３５は、式Ａμ０＝ｚ０により、μ０を得る（ステップＳ３５）。
【００６５】
カルマンフィルタ推定部３６は、検証用データＤｖｒに含まれる観測データｘ（ｔ）（ｔ＝１，．．．，Ｔ２）に対して、カルマンフィルタによる状態推定演算を実行する。まずは、ｔ＝１とする（ステップＳ３６）。カルマンフィルタ推定部３６は、ｔがＴ２以下であるか判定する（ステップＳ３７）。ここでは、ｔがＴ２以下であるため、次のステップＳ３８に進む。
【００６６】
カルマンフィルタ推定部３６は、下記の式でｔ＝１として、Ｋ１、μ１、Ｖ１およびＱ１を求める（ステップＳ３８）。
【００６７】
【数１４】
2019191836.tif 000018
【００６８】
カルマンフィルタ推定部３６は、式Ａμ１＝ｚ１から、ｚ１を求める。
【００６９】
更に異常度計算部３７は、下記の式により定義した異常度ａ（ｘ（ｔ））を算出することにより、監視対象１の異常状態を定量的に評価する（ステップＳ３９）。
【００７０】
【数１５】
2019191836.tif 000019
【００７１】
異常度ａ（ｘ（ｔ））が所定のしきい値を越える場合には、例えば異常度計算部３７はオペレータに警告を行うことができる。このように、異常検知システム３は、速やかに状態異常を検知してオペレータに警告することができる。
【００７２】
さらに、カルマンフィルタ推定部３６は、ｔを１増加させる（ステップＳ４０）。
【００７３】
ここまで説明したステップＳ３７からＳ４０を繰り返すことで、カルマンフィルタ推定部３６は、各時刻ｔ＝１，．．．，Ｔ２における状態変数ｚ（ｔ）の確率分布（平均と分散（μｔ、Ｖｔ）)の推定値を計算することができる。また異常度計算部３７は、カルマンフィルタ推定部３６が計算した前記推定値に基づき、前記監視対象の異常度を計算することができる。
【００７４】
ｔ＝Ｔ２＋１であるときに、カルマンフィルタ推定部３６は、ステップＳ３７で、ｔがＴ２以下でないと判定する。その後処理はステップＳ３１に戻る。
【００７５】
なお、本実施形態では、カルマンフィルタ推定部３６がカルマンフィルタに基づいて状態推定演算を行っているが、他の同定法で状態推定演算を行うこともできる。
【００７６】
その後、監視対象に対する本実施形態では一定のサンプリング周期が経過し、監視対象１の新たな観測データξ（ｋ＋Ｔ１＋Ｔ２＋１）が得られたときには（ステップＳ３１）、前記と同様に、観測データｘ（ｔ）の時系列データを、次式により更新する（ステップＳ３２）。
【００７７】
【数１６】
2019191836.tif 000020
【００７８】
これにより、図２に示すように時系列観測データＤ１がＤ２に更新される。その後、上述したシステム同定と異常度計算（ステップＳ３８）とを含むステップＳ３３からＳ４０を実行することにより、異常検知システム３は新たな時系列観測データＤ２に基づく異常度を得ることができる。
【００７９】
このように、監視対象の新たな観測データが得られるごとに、時系列観測データを更新し、上述したシステム同定と異常度計算を繰り返すことによって、最新の観測データに基づいた異常度をリアルタイムで得ることができる。
【００８０】
前記の如く、本実施形態の異常検知システムを適用することによって、異常検知アルゴリズムにおいて、状態空間モデルのシステム同定をオンライン化することができ、システム同定のための期間を別途設ける必要が無くなる。例えば具体的には、時系列観測データを車両外のシステムに移して車両外システムで異常度を求める必要がなくなる。
【００８１】
更に、監視対象の直近の観測データに基づいて、状態空間モデルのシステム同定が行われるため、監視対象１の特性変動にも追従可能となる。そのため、本実施形態の異常検知システム３によれば、状態異常の検知精度を向上させることができる。
【実施例】
【００８２】
図４Ａは、本発明の異常検知システムの実施例における構成を示したものである。
【００８３】
符号４は異常検知システムを表す。この実施例の異常検知システムの基本的な構成は図１に示した異常検知システム３と同じであるが、図１において観測データが保存されるデータベース３２と訓練用データが保存されるデータベース３４ａと検証用データが保存されるデータベース３４ｂとは、具体的なメモリ装置４２として実装されている。また、図４Ａには、図１に示す観測データ分類部３３は設けられておらず、後述するように、観測データ更新部３１が観測データを訓練用データおよび検証用データのいずれかに分類する。
【００８４】
メモリ装置４２は、観測データ１点分（ある時刻に取得した観測データ）のエリアごとにアドレッシングされている。符号４２ａは訓練用データＤｔｒの開始位置を指示するポインタを表す。符号４２ｂは検証用データＤｖｒの開始位置を指示するポインタを表す。符号４２ｃは、次に取得する新しい観測データを書き込む位置を指示するポインタを表す。
【００８５】
ここで、観測データ用のメモリエリアは有限であるため、論理的なリングバッファを構成するように、前記の各ポインタのアドレッシングが管理されているものとする。
【００８６】
図４Ｂは、或る時点での各ポインタの指示位置と、新しい観測データがサンプリングされて観測データ更新部４１によって時系列観測データが更新された後の各ポインタの指示位置とを表している。例えば図４Ｂにおいて左側に示す、或る時点での時系列観測データは、図２の時系列観測データＤ０に対応する。また、図４Ｂにおいて右側に示す、時系列観測データが更新された後の時系列観測データは、図２の時系列観測データＤ１に対応する。
【００８７】
図４Ｂに示すポインタ４２ａ’は、ポインタ４２ａが指すデータの次のデータを指している。ポインタ４２ｂ’および４２ｃ’についても同様である。新しい観測データがサンプリングされた際に、観測データ更新部４１は、ポインタをこのようにずらして、最も古い観測データ、すなわち最も前に収集された観測データが参照されないようにする。そして、更に新しく収集する観測データが最も古い観測データが記憶されている位置に書き込まれるようにする。
【００８８】
上記のように、各ポインタのアドレッシングを用いて時系列観測データの一部が更新される。これにより、訓練用データと検証用データのそれぞれの内容を一サンプリング時間ずつシフトすることが可能となり、各データの更新およびシステム同定部で使用される行列演算を効率よく実行することができる。
【００８９】
以上、異常検知システム３，４について説明したが、異常検知システム３，４として機能させるためにコンピュータを好適に用いることが可能である。そのようなコンピュータは、異常検知システム３，４の各機能を実現する処理内容を記述したプログラムを該コンピュータの記憶部に格納しておき、該コンピュータのＣＰＵによってこのプログラムを読み出して実行させることで実現することができる。
【００９０】
また、このプログラムは、コンピュータ読取り可能媒体に記録されていてもよい。コンピュータ読取り可能媒体を用いれば、コンピュータにインストールすることが可能である。ここで、プログラムが記録されたコンピュータ読取り可能媒体は、非一過性の記録媒体であってもよい。非一過性の記録媒体は、特に限定されるものではないが、例えば、ＣＤ−ＲＯＭまたはＤＶＤ−ＲＯＭなどの記録媒体であってもよい。
【００９１】
上述の実施形態は代表的な例として説明した。しかしながら、本発明の趣旨および範囲内で、多くの変更および置換ができることは当業者に明らかである。したがって、本発明は、上述の実施形態によって制限するものと解するべきではなく、特許請求の範囲から逸脱することなく、種々の変形または変更が可能である。例えば、実施形態の構成図に記載の複数の構成ブロックを１つに組み合わせたり、あるいは１つの構成ブロックを分割したりすることが可能である。
【産業上の利用可能性】
【００９２】
本発明は、線形状態空間モデルによるモデリングが可能なシステムであれば電子装置のみに適用されるものではなく、適切なセンサを備えることにより機械系（例えば駆動系）の異常診断システムへも適用可能である。
【符号の説明】
【００９３】
１…監視対象
３…異常検知システム
３１…観測データ更新部
３２…データベース
３３…観測データ分類部
３４ａ…データベース
３４ｂ…データベース
３５…システム同定部
３６…カルマンフィルタ推定部
３７…異常度計算部
４…異常検知システム
４１…観測データ更新部
４２…メモリ装置
４２ａ，４２ａ’…ポインタ
４２ｂ，４２ｂ’…ポインタ
４２ｃ，４２ｃ’…ポインタ
４５…システム同定部
４６…カルマンフィルタ推定部
４７…異常度計算部
【図１】
2019191836.tif 000021
【図２】
2019191836.tif 000022
【図３】
2019191836.tif 000023
【図４Ａ】
2019191836.tif 000024
【図４Ｂ】
2019191836.tif 000025
【図５】
2019191836.tif 000026
【図６】
2019191836.tif 000027
【図７】
2019191836.tif 000028