(57)【要約】【課題】システムから得られた時系列の多変量データを用いて、次元間の非線形な因果関係を推定することを可能とする技術を提供する。【解決手段】因果推定装置において、時系列多次元数値ベクトルのデータを入力する入力部と、入力された時系列多次元数値ベクトルのデータを用いて、ある時刻のデータを過去の時刻のデータから予測する非線形回帰モデルを学習する回帰モデル学習部と、前記非線形回帰モデルを用いて、前記時系列多次元数値ベクトルのデータにおける次元jの次元iに対する因果の強さを計算する因果推定部と、前記因果推定部により計算された因果の強さを出力する出力部とを備える。【選択図】図１
2019144779.tif 000002
【特許請求の範囲】
【請求項１】
時系列多次元数値ベクトルのデータを入力する入力部と、
入力された時系列多次元数値ベクトルのデータを用いて、ある時刻のデータを過去の時刻のデータから予測する非線形回帰モデルを学習する回帰モデル学習部と、
前記非線形回帰モデルを用いて、前記時系列多次元数値ベクトルのデータにおける次元jの次元iに対する因果の強さを計算する因果推定部と、
前記因果推定部により計算された因果の強さを出力する出力部と
を備えることを特徴とする因果推定装置。
【請求項２】
前記因果推定部は、
前記非線形回帰モデルにおける時刻t-pの次元jの誤差項の変動が時刻tの次元iに与える影響を用いて前記因果の強さを計算する、又は、
時刻t-pにおける次元jに微小量Δが与えられた際の時刻tにおける次元iの前記非線形回帰モデルによる予測値と、当該微小量が与えられなかった場合の時刻tにおける次元iの前記非線形回帰モデルによる予測値との誤差を用いて前記因果の強さを計算する、又は、
前記非線形回帰モデルによる次元iの予測値における時刻t-pの次元jの値を含む項を用いて前記因果の強さを計算する
ことを特徴とする請求項１に記載の因果推定装置。
【請求項３】
前記回帰モデル学習部は、スパース項を考慮したスパースモデリングにより前記非線形回帰モデルの学習を行う
ことを特徴とする請求項１又は２に記載の因果推定装置。
【請求項４】
前記回帰モデル学習部は、ニューラルネットワークを用いて前記非線形回帰モデルの学習を行う
ことを特徴とする請求項１ないし３のうちいずれか１項に記載の因果推定装置。
【請求項５】
前記回帰モデル学習部は、前記非線形回帰モデルの計算時に、当該非線形回帰モデルの各パラメータの重要度を計算し、
前記因果推定部は、前記重要度を用いて前記因果の強さを計算する
ことを特徴とする請求項１又は２に記載の因果推定装置。
【請求項６】
因果推定装置が実行する因果推定方法であって、
時系列多次元数値ベクトルのデータを入力する入力ステップと、
入力された時系列多次元数値ベクトルのデータを用いて、ある時刻のデータを過去の時刻のデータから予測する非線形回帰モデルを学習する回帰モデル学習ステップと、
前記非線形回帰モデルを用いて、前記時系列多次元数値ベクトルのデータにおける次元jの次元iに対する因果の強さを計算する因果推定ステップと、
前記因果推定ステップにより計算された因果の強さを出力する出力ステップと
を備えることを特徴とする因果推定方法。
【請求項７】
コンピュータを、請求項１ないし５のうちいずれか１項に記載の因果推定装置における各部として機能させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、システムから収集した時系列数値データを分析し、データ間の因果関係を推定する技術に関するものである。なお、本明細書で使用する「因果」とは、あくまでデータ上に表れる関係性に基づいた因果であり、例えばデータAが変動した後にはデータBにも変動が見られるといったような事実から推定されるものである。データ上の因果がその背後にある「真の因果」を表しているとは限らないが、システムの挙動把握や異常時の原因推定においては、データ上の因果も十分に有用であると考えられるため、本発明において推定の対象とする。
【背景技術】
【０００２】
システムから時系列の多変量データが得られる場合に、得られたデータに基づいてデータ間の因果関係を推定することが、システムの挙動の理解やシステムに異常が発生した場合の原因解明に向けて重要である（非特許文献１、非特許文献２）。
【０００３】
特に、対象データが時系列である場合、過去のデータを用いて将来のデータを予測するようなベクトル自己回帰（VAR）に基づいた、グレンジャー因果性（非特許文献３）やインパルス応答関数（非特許文献４）による因果推定は、入力データが多次元に渡る場合でも、少ない時間で計算が可能である。特に、後者のインパルス応答関数に関しては、因果の強さを定量的に評価できるといった利点がある。 
【先行技術文献】
【非特許文献】
【０００４】
【非特許文献１】Kobayashi, Satoru, Kensuke Fukuda, and Hiroshi Esaki. "Causation mining in network logs." ACM SIGCOMM CoNEXT 2016 Student Workshop. 2016.
【非特許文献２】Gonzalez, Jose Manuel Navarro, Javier Andion Jimenez, and Juan Carlos Duenas Lopez. "Root Cause Analysis of Network Failures Using Machine Learning and Summarization Techniques." IEEE Communications Magazine 55.9 (2017): 126-131.
【非特許文献３】Barnett, Lionel, Adam B. Barrett, and Anil K. Seth. "Granger causality and transfer entropy are equivalent for Gaussian variables." Physical review letters 103.23 (2009): 238701.
【非特許文献４】Pesaran, H. Hashem, and Yongcheol Shin. "Generalized impulse response analysis in linear multivariate models." Economics letters 58.1 (1998): 17-29.
【非特許文献５】Koop, Gary, M. Hashem Pesaran, and Simon M. Potter. "Impulse response analysis in nonlinear multivariate models." Journal of econometrics 74.1 (1996): 119-147.
【非特許文献６】Shimizu, Shohei, et al. "A linear non-Gaussian acyclic model for causal discovery." Journal of Machine Learning Research 7.Oct (2006): 2003-2030.
【発明の概要】
【発明が解決しようとする課題】
【０００５】
VARを用いた一般的なインパルス応答関数分析は、線形な回帰に基づいている。しかし、システムから得られるデータは、線形な関係だけではなく、非線形な関係も多く含まれると考えられる。特に、syslogの出現有無のようなものがデータに含まれている場合、あるsyslogとあるsyslogが同時に出現した場合（AND）、もしくはどちらか一方が出現した場合（OR）、別のsyslogが出現するといった、非線形な因果関係が考えられる。
【０００６】
非線形なインパルス応答関数についての理論的な議論は（非特許文献５）で与えられているが、実用において、システムデータにおける複雑な関係性を十分に表すことができ、かつインパルス応答関数についても理論的に導出可能であるような非線形回帰を実現する具体的な方法については提案されていない。
【０００７】
インパルス応答関数以外で、多変量データ間の因果関係について推定する手法としては、PCアルゴリズム（非特許文献１）やLiNGAM（非特許文献６）等が提案されているが、PCアルゴリズムについては、密な因果関係を持つ場合には計算量が非常に多くなってしまい、かつ因果の強さについては推定できず、またLiNGAMについては、線形な関係性を前提としている。そのため、多次元に渡るデータにおける非線形な因果推定をどのように実現するかが課題になっている。
【０００８】
本発明は上記の点に鑑みてなされたものであり、システムから得られた時系列の多変量データを用いて、次元間の非線形な因果関係を推定することを可能とする技術を提供することを目的とする。
【課題を解決するための手段】
【０００９】
開示の技術によれば、時系列多次元数値ベクトルのデータを入力する入力部と、
入力された時系列多次元数値ベクトルのデータを用いて、ある時刻のデータを過去の時刻のデータから予測する非線形回帰モデルを学習する回帰モデル学習部と、
前記非線形回帰モデルを用いて、前記時系列多次元数値ベクトルのデータにおける次元jの次元iに対する因果の強さを計算する因果推定部と、
前記因果推定部により計算された因果の強さを出力する出力部と
を備えることを特徴とする因果推定装置が提供される。
【発明の効果】
【００１０】
開示の技術によれば、システムから得られた時系列の多変量データを用いて、次元間の非線形な因果関係を推定することを可能とする技術が提供される。
【図面の簡単な説明】
【００１１】
【図１】本発明の実施の形態における因果推定装置１００の構成図である。
【図２】因果推定装置１００のハードウェア構成図である。
【図３】実施例１における処理手順を示すフローチャートである。
【図４】実施例３と実施例４の組み合わせで因果の強さを算出する例を示す図である。
【図５】シミュレーションにおいて生成したデータのイメージである。
【図６】シミュレーションにおいてN=100とした場合の精度評価の結果を示す図である。
【図７】シミュレーションにおいてN=500とした場合の精度評価の結果を示す図である。
【発明を実施するための形態】
【００１２】
以下、図面を参照して本発明の実施の形態（本実施の形態）を説明する。以下で説明する実施の形態は一例に過ぎず、本発明が適用される実施の形態は、以下の実施の形態に限られるわけではない。
【００１３】
（システム構成）
図１に、本実施の形態における因果推定装置１００の構成例を示す。図１に示すように、本実施の形態における因果推定装置１００は、入力部１０１、記憶部１０２、因果推定部１０３、回帰モデル学習部１０４、及び出力部１０５を有する。
【００１４】
入力部１０１は、時系列の多次元数値ベクトルデータ、各種パラメータ等の外部の情報を因果推定装置１００に入力させる。記憶部１０２は、入力部１０１より入力されたデータ、モデル及びパラメータ等を保持する。因果推定部１０３は、次元間の因果の強さの計算を行う。回帰モデル学習部１０４は、非線形回帰モデルの学習を行う。出力部１０５は、因果推定部１０３により計算された次元間の因果の強さの出力を行う。回帰モデル学習部１０４と因果推定部１０３における処理については、後述の実施例１〜６において詳細に説明する。
【００１５】
（ハードウェア構成例）
上述した因果推定装置１００は、例えば、コンピュータに、本実施の形態で説明する処理内容を記述したプログラムを実行させることにより実現可能である。
【００１６】
すなわち、因果推定装置１００は、コンピュータに内蔵されるＣＰＵやメモリ等のハードウェア資源を用いて、因果推定装置１００で実施される処理に対応するプログラムを実行することによって実現することが可能である。上記プログラムは、コンピュータが読み取り可能な記録媒体（可搬メモリ等）に記録して、保存したり、配布したりすることが可能である。また、上記プログラムをインターネットや電子メール等、ネットワークを通して提供することも可能である。
【００１７】
図２は、本実施の形態における上記コンピュータのハードウェア構成例を示す図である。図２のコンピュータは、それぞれバスＢで相互に接続されているドライブ装置１５０、補助記憶装置１５２、メモリ装置１５３、ＣＰＵ１５４、インタフェース装置１５５、表示装置１５６、及び入力装置１５７等を有する。
【００１８】
当該コンピュータでの処理を実現するプログラムは、例えば、ＣＤ−ＲＯＭ又はメモリカード等の記録媒体１５１によって提供される。プログラムを記憶した記録媒体１５１がドライブ装置１５０にセットされると、プログラムが記録媒体１５１からドライブ装置１５０を介して補助記憶装置１５２にインストールされる。但し、プログラムのインストールは必ずしも記録媒体１５１より行う必要はなく、ネットワークを介して他のコンピュータよりダウンロードするようにしてもよい。補助記憶装置１５２は、インストールされたプログラムを格納すると共に、必要なファイルやデータ等を格納する。
【００１９】
メモリ装置１５３は、プログラムの起動指示があった場合に、補助記憶装置１５２からプログラムを読み出して格納する。ＣＰＵ１５４は、メモリ装置１５３に格納されたプログラムに従って、モデル学習装置１００に係る機能を実現する。インターフェース装置１５５は、ネットワークに接続するためのインターフェースとして用いられる。表示装置１５６はプログラムによるＧＵＩ（Graphical User Interface）等を表示する。入力装置１５７はキーボード及びマウス、ボタン、又はタッチパネル等で構成され、様々な操作指示を入力させるために用いられる。なお、表示装置１５６を備えないこととしてもよい。
【００２０】
以下、因果推定装置１００の動作例を実施例１〜６として説明する。以下に示す実施例１が基本的な動作例であり、実施例２〜６では、実施例１と異なる点を主に説明している。
【００２１】
（実施例１）
実施例１では、入力された時系列の多次元数値ベクトルデータを用いて非線形回帰モデルx_t = c+f(x_t-τ,x_t-τ+1,…,x_t-1)+ε_tを推定し，モデルのインパルス応答関数によって次元間の因果推定を行う例を説明する。
【００２２】
実施例１における因果推定装置１００の動作を図３のフローチャートに沿って説明する。
【００２３】
Ｓ１０１）入力部１０１より、システムより収集された時系列多次元数値ベクトルデータセットX={x_1、…、x_T}が入力される。収集されるデータの例としては、各時刻における、各インターフェース上におけるトラヒック量や、CPU・メモリ負荷、テンプレート化されたsyslogのID出現回数などある。
【００２４】
Ｓ１０２）回帰モデル学習部１０４は、入力されたXを用いて、非線形回帰モデルx_t = c+f(x_t-τ,x_t-τ+1,…,x_t-1)+ε_t（ただしcは定数項，fは任意の非線形関数，ε_tは時刻tにおける誤差項）を学習する。ここで、モデル式z=f(y)としては累乗モデルz=a*y^b、指数モデルz=a*b^y等任意のモデルが考えられる。学習方法としても、最小二乗法を用いた回帰（Bohme, J. "Estimation of source parameters by maximum likelihood and nonlinear regression." Acoustics, Speech, and Signal Processing, IEEE International Conference on ICASSP'84.. Vol. 9. IEEE, 1984.）など、任意の方法を用いて良い。モデル及び学習方法の選択については、あらかじめ定めておき記憶部１０２に保存してもよいし、入力部１０１から入力することで選択してもよい。
【００２５】
Ｓ１０３）因果推定部１０３は、学習したモデルを元に、非線形回帰モデルのインパルス応答関数を計算する。ここで、インパルス応答関数とは、時刻t-pのデータの次元jに対して与えられたショックが、時刻tのデータの次元iに対してどの程度影響を及ぼすかを表すものであり、x_{t,i}のε_{t-p,j}での偏微分∂x_{t,i}／∂ε_{t-p,j}（p時刻前の次元jの誤差項の変動が次元iに与える影響を示す）で定義される。一般的なインパルス応答関数に関する議論は非特許文献５で与えられているが、ここでは簡単のため、モデル式fが任意のyで可微分かつ誤差項ε_tが次元間で独立である場合について説明する。任意のpに対するインパルス応答関数は、以下のように再帰的に計算できる。
【００２６】
まず、データセットx_t-τ,x_t-τ+1,…,x_t-1が与えられた元での、次元jのショックに対する時刻p後の次元iのインパルス応答関数をIRF_{i,j}(p，x_t-τ,x_t-τ+1,…,x_t-1)として定義する。これは、後述する通りp>0においてはインパルス応答関数がデータx_t-τ,x_t-τ+1,…,x_t-1に依存するためである。p=0である場合のインパルス応答関数は、定義より
【００２７】
【数１】
2019144779.tif 000003
のように定数で与えられる。次に、p=1の場合には、微分の連鎖律及び上式により、
【００２８】
【数２】
2019144779.tif 000004
となる。ここで、f_i(・)はf(・)において次元iの値を与える関数である。p=2については、
【００２９】
【数３】
2019144779.tif 000005
となることから、IRF_{i,j}(p, x_t-τ,x_t-τ+1,…,x_t-1)については
【００３０】
【数４】
2019144779.tif 000006
のように一般化できる。上式はx_t-τ,x_t-τ+1,…,x_t-1に依存するため、非特許文献５における議論と同様、期待値を取ることにより次元jのショックに対する時刻p後の次元iのインパルス応答関数IRF_{i,j}(p)を
【００３１】
【数５】
2019144779.tif 000007
として求めることができる。ここで，E[・]は・の期待値を表している。期待値の計算方法としては、x_tの事前分布に基づいて数値積分を行う方法や、収集されたデータセットX上における
【００３２】
【数６】
2019144779.tif 000008
の平均を取る方法などがある。IRFの計算においては、回帰モデルの微分
【００３３】
【数７】
2019144779.tif 000009
が必要であるが、これは、各モデルに対応した微分式を予め記憶部１０２に保存する方法や、モデルを入力部１０１より入力する場合、微分式についても入力させる方法、数値的に計算する方法などがある。
【００３４】
因果推定部１０３は、計算したインパルス応答関数IRF_{i,j}(0),…, IRF_{i,j}(p_max)を元に、次元jから次元iに向けた因果の強さを計算する。p_maxについては、予め定めた値を記憶部１０２に保存する方法や、入力部１０１より与える方法が考えられる。計算方法としては、単純にIRF_{i,j}(0),…, IRF_{i,j}(p_max)の中の一つを用いる方法や、総和を取る方法、重み付け平均を取る方法、絶対値が最大となる値を採用する方法など、様々な方法がある。
【００３５】
Ｓ１０４）因果推定部１０３は、全ての次元間について因果の強さを計算し、出力として、次元数をNとした場合に、i行j列の要素が次元jから次元iに向けた因果の強さとなるようなN×Nの行列を出力部１０５より出力する。
【００３６】
（実施例２）
実施例２において、因果推定装置１００の動作の全体の流れは実施例１で説明した図３に示す流れと同じであるが、Ｓ１０３における因果の強さの計算方法が実施例１と異なる。
【００３７】
実施例２では、因果推定部１０３は、実施例１のようなインパルス応答関数ではなく、次元jに微小量を与えた際の次元iの予測値の変化を元に因果の強さを計算する。時刻t-pにおける次元jに微小量Δが与えられた際の時刻tの次元iの予測値と、微小量が与えられなかった場合の予測値との誤差をDIFF_{i,j}(p, x_t-τ,x_t-τ+1,…,x_t-1)とおくと、これは
【００３８】
【数８】
2019144779.tif 000010
として与えられる。インパルス応答関数と同様、x_t-τ,x_t-τ+1,…,x_t-1に依存するため，因果の強さを求める際には期待値を取り、
【００３９】
【数９】
2019144779.tif 000011
として、DIFF_{i,j}(1)，... ，DIFF_{i,j}(p_max)を用いて実施例１と同様に因果の強さを決定し、全ての次元間についての因果の強さが出力部１０５より出力される。
【００４０】
（実施例３）
実施例３において、因果推定装置１００の動作の全体の流れは実施例１で説明した図３に示す流れと同じであるが、Ｓ１０３における因果の強さの計算方法が実施例１と異なる。
【００４１】
実施例３における因果推定部１０３は、インパルス応答関数ではなく次に説明する方法で因果の強さを求める。
【００４２】
因果推定部１０３は、f_i(x_t-τ,x_t-τ+1,…,x_t-1)においてx_{t-p,j}を含む項{a_1 g_1(x_t-τ,x_t-τ+1,…,x_t-1), …,a_M g_M(x_t-τ,x_t-τ+1,…,x_t-1)}（a_mは定数，g_mは関数）のみを抽出し、定数a_mや関数g_mにおけるx_{t-p,j}の次数を用いて、次元jから次元iへの因果の強さを決定する。
【００４３】
例えば、fが累乗モデルであり、f_i(x_t-τ,x_t-τ+1,…,x_t-1)においてx_{t-p,j}を含む項がa*x_{t-p,j}^b*g(x_t-τ,x_t-τ+1,…,x_t-1)として与えられるとする。ここで、gはx_{t-p,i}以外の変数に関する関数である。この時、次元jの値が時刻p後の次元iに与える影響を定数a, bや関数gを用いて表現する。
【００４４】
例えば、単純に係数aを影響の強さとしても良いし、a*bのような積の形で与えても良い。また、g(x_t-τ,x_t-τ+1,…,x_t-1)は変数x_t-τ,x_t-τ+1,…,x_t-1に依存する関数であるが、実施例１及び実施例２と同様に，期待値を取りaやbに掛け合わせる方法を用いてもよい。そのような計算をp=1,…,p_maxについて行い、それらの値を用いて実施例１と同様に因果の強さを決定する。
【００４５】
（実施例４）
実施例４において、因果推定装置１００の動作の全体の流れは実施例１で説明した図３に示す流れと同じであるが、Ｓ１０２における学習方法が実施例１と異なる。なお、実施例４は、実施例２、３にも適用可能である。
【００４６】
実施例４では、回帰モデル学習部１０４が非線形回帰を行う際に、スパースなモデリングを行うために、学習時にスパース項L(x_t-τ,x_t-τ+1,…,x_t-1)を考慮した学習を行う。これは、非線形回帰の過学習によって誤ったパラメータ推定をおこなってしまい、本来は存在しない因果について存在すると誤って推定したり、存在する因果を見逃してしまう事を防ぐためである。
【００４７】
回帰モデル学習部１０４が実行するスパース項を考慮した学習方法としては、最小二乗法を用いた回帰において、目的関数に対してL2ノルム項λL_2(x_t-τ,x_t-τ+1,…,x_t-1)=λΣ_{i=1}^τ||x_{t-i}||^2をペナルティ項として加えた最小化を行う方法（λは事前に与えた定数か、又は入力部１０１より入力）や、L1ノルム項λL_1(x_t-τ,x_t-τ+1,…,x_t-1)=λΣ_{i=1}^τ||x_{t-i}||^1を加えた最小化について近接勾配法（Beck, Amir, and Marc Teboulle. "A fast iterative shrinkage-thresholding algorithm for linear inverse problems." SIAM journal on imaging sciences 2.1 (2009): 183-202.）を用いて解く方法などがある。
【００４８】
（実施例５）
実施例５において、因果推定装置１００の動作の全体の流れは実施例１で説明した図３に示す流れと同じであるが、Ｓ１０２における学習方法等が実施例１と異なる。実施例５は、実施例２〜４にも適用可能である。
【００４９】
実施例５では、回帰モデル学習部１０４は、非線形回帰をニューラルネットワークを用いて行う。ニューラルネットワークは、多様な非線形回帰を単純なモデル化で実現できるといった利点や、実施例１において必要となる微分項についても連鎖則を用いて容易に計算できるといった利点がある。
【００５０】
非線形回帰x_t = c+f(x_t-τ,x_t-τ+1,…,x_t-1)+ε_tをニューラルネットワークを用いて行う場合、入力層がx_t-τ,x_t-τ+1,…,x_t-1を入力とするτ×N次元ノード、出力層がx_tを出力とするN次元ノードとなるようなニューラルネットワークを設計し、ニューラルネットワークのパラメータについては、データセットXを用いて学習を行うことで取得し、記憶部１０２に保存する。
【００５１】
ニューラルネットワークにおける中間層の数や次元数、活性化関数、学習時のパラメータ（バッチサイズ，学習エポック数等）については、予め定めておき記憶部１０２に保存しておく方法や、入力部１０１より与えて指定する方法がある。
【００５２】
なお、実施例１のＳ１０３において必要となるx_{t-p,j}によるx_{t,i}=f_i(x_t-τ,x_t-τ+1,…,x_t-1)の微分
【００５３】
【数１０】
2019144779.tif 000012
については、バックプロパゲーション法（Goh, A. T. C. "Back-propagation neural networks for modeling complex systems." Artificial Intelligence in Engineering 9.3 (1995): 143-151.）により計算が可能である。このようにバックプロパゲーション法により微分を計算することで、実施例１におけるＳ１０３と同様にして、インパルス応答関数を計算し、次元jから次元iに向けた因果の強さを計算する。
【００５４】
また、次元数が多いなどで微分の計算量が大きくなってしまうといった場合には、微分の計算の代わりに、実施例３のように係数のみを用いて因果の強さを計算する事としてもよい。例えば、p時刻前の次元jが次元iに与える因果の強さについて、入力層のx_{t-p,j}と、出力層のx_{t,i}とを結ぶリンク上の重みの積を全ての経路分について足し合わせるといった方法を用いることができる。
【００５５】
図４には、一例として、入力層、中間層、出力層を有する３層のニューラルネットワークにおいて、p=1時刻前の次元j=1が次元i=3に与える因果の強さが、入力層のx_{t-1,1}と、出力層のx_{t,3}とを結ぶリンク上の重みの積を全ての経路分について足し合わせた値w^1_11 * w^1_12 * w^1_13 + w^2_13 * w^2_23 * w^2_33で計算されることが示されている。
【００５６】
（実施例６）
実施例６は、実施例１、実施例２、実施例３のぞれぞれと、Ｓ１０２における回帰モデルの計算方法、Ｓ１０３における因果の強さを計算する方法が異なり、その他の処理は同様である。
【００５７】
実施例６では、因果推定部１０３が、実施例１、実施例２、実施例３において因果の強さを計算する際に、非線形回帰モデルの各パラメータの重要度を考慮した計算を行う。ここで重要度とは、各パラメータの非線形回帰に対する寄与の強さを表すようなものであり、寄与が強いパラメータほど因果の強さを推定する上で重要になると仮定している。パラメータの重要度としては、例えばモデルのデータに対するフィッシャー情報量（Jauffret, Claude. "Observability and Fisher information matrix in nonlinear regression." IEEE Transactions on Aerospace and Electronic Systems 43.2 (2007).）を用いる方法などある。
【００５８】
実施例６において、回帰モデル学習部１０４は、回帰モデルの計算時に、各パラメータθ_1,…,θ_Kの重要度F_1,…,F_Kについても計算し、記憶部１０２へと保存する。
【００５９】
因果推定部１０３は、各パラメータの重要度についても考慮して実施例１、実施例２、実施例３における因果の強さの計算を行う。その方法としては、非線形回帰モデルを用いて実施例１〜３で説明した方法で因果の強さの計算を行う際に、例えば単純に非線形回帰モデルのパラメータの値θ_kに対して重要度F_kをかけた値θ_k*F_kを新たなパラメータθ'_kとみなす方法や、重要度に対して閾値を与え、F_kが閾値未満の場合、θ_k=0とみなしてしまうなどといった方法がある。
【００６０】
（効果について）
実施例を用いて説明した本発明に係る技術により、システムから得られた時系列の多変量データを用いて次元間の非線形な因果関係を定量的に評価することが可能となる。
【００６１】
ここでは、効果を示すために、一例として、実施例１、４、５を組み合わせる事により、ニューラルネットワークを用いてスパースな学習を行った非線形回帰モデルにおいて、インパルス応答関数を用いた因果推定を行った結果を示す。
【００６２】
ここでは、ラグτ=1で以下の通り因果関係を与えたN個のsyslog idの出現x_i, i=1,…Nに関するデータをシミュレーションにより生成し、因果推定装置１００による因果推定を行った。
【００６３】
データに対して与えられた因果のイメージを図５に示す。この例では、各時刻tにおいて、idがi=1,…,N/2であるsyslog（例：図５のsyslog id=1, 2）については、後述する通り一時刻前の出現有無に依存するベルヌーイ分布で出現確率が決まり、idがi=N/2+1,…,Nのsyslog（例：図５のsyslog id=51, 52）については、一時刻前のi=1,…,N/2のsyslogの出現に依存し出現有無が確定する。以下、各時刻のsyslog出現のルールについてより具体的に説明する。
【００６４】
全てのi(i<N/2)についてx_{i,t}=1の場合、確率q_contのベルヌーイ分布で、x_{i,t}=0の場合、確率q_iのベルヌーイ分布でq_{i,t+1}が決まる。ここではq_cont=0.7、q_iはi%2=1の場合0.5、i%2=0の場合0.01とした。全てのi(i%2=1∧i<N/2)について、x_{i,t}=1かつx_{i+1,t}=1であればx_{i+N/2,t+1}=1であり、x_{i,t}=1又はx_{i+1,t}=1であればx_{i+N/2+1,t+1}=1である。すなわち、i→i+N/2, i→i+N/2+1, i+1→i+N/2, i+1→i+N/2+1 (i<N/2)の因果関係があることになる。図５の例では、1→51, 1→52, 2→51, 2→52の因果関係が示されている。x_t, t=1,…,Tを観測データXとし、因果推定装置１００により上記の因果関係を推定した。
【００６５】
因果推定の評価については、データ取得期間Tを1000,10000,100000と変えた場合の評価を実施した。1000,10000,100000は、毎分のデータ取得とした場合、約16時間、1週間、2ヶ月強分のデータ量に相当する。
【００６６】
k番目のデータx_kからl番目のデータx_lへの因果有無を，実施例５を用いてニューラルネットワークにより学習した非線形回帰モデル（τ＝１）において、実施例１を用いて計算したIRF_{l,k}(1)に閾値を与えて因果の有無を判断し、閾値を変化させた時のPR-AUCを比較した。ここで、PR-AUCとは、閾値に依存して決まるprecision（因果があると判断した組のうち、実際に因果があった組の割合）及びrecall（実際に因果がある組のうち、因果があると判断できた組の割合）について、閾値を変化させた場合の変化を、横軸にrecall、縦軸にprecisionを取った際にプロットされるPR曲線について、その曲線の下側の面積を表したものであり，PR-AUCが高いほど推定精度が高いとされる。
【００６７】
比較対象は回帰モデルとして線形VARを用いた際のIRF（従来技術における非特許文献４）である。また、ニューラルネットワークのモデルとしては、活性化関数はsigmoid、重み減衰をλとして与え、実施例４におけるL2ノルム項を追加した学習により、パラメータがスパースになるような学習を行った。比較対象のモデルについては、中間層を1層とし、次元数を入力次元のrh倍としたモデル（DNN）と、入力層と出力層のみのモデル（2-layer NN。L2ノルム付き線形VARに相当する。）を比較した。
【００６８】
データの次元数N=100の場合の結果を図６に示し、N=500の場合の結果を図７に示す。横軸はデータ取得期間であり、T=1000の場合、T=10000の場合、及びT=100000の場合の３パターンについて評価している。縦軸はPR-AUCであり、PR-AUCが高いほど因果推定の精度が高いと言える。なお、N=100の場合は、L2ノルム項の係数λは10^-4、N=500の場合には10^-5としている。図６、７に示すように、線形VARや2層のニューラルネットワークに比べて、中間層を与えた非線形ニューラルネットワークによる因果推定は少ないデータ取得期間で高精度に因果推定を行うことができており、非線形な回帰によって、非線形な因果関係の推定が高精度にできていることが確認できる。
【００６９】
（実施例のまとめ）
以上、説明したように、実施例１においては、システム監視データが、N次元の時系列多次元数値ベクトルとして表現される場合に、時刻tにおけるデータをx_t =(x_{t,1},…,x_{t,N})とする。因果推定装置１００は、収集したデータセットX={x_1,…,x_T}を用いて、時刻tのデータを時刻t-τ〜時刻t-1のデータにより表現する非線形回帰モデルx_t = c+f(x_t-τ,x_t-τ+1,…,x_t-1)+ε_t（ただしcは定数項，fは任意の非線形関数，ε_tは時刻tにおける誤差項）を学習し、監視データにおける次元jの次元iに対する因果の強さを、p時刻前の次元jの誤差項の変動が次元iに与える影響∂x_{t,i} / ∂ε_{t-p,j}（p=1,…,p_max）を用いて計算する。
【００７０】
実施例２では、実施例１において、因果推定装置１００は、監視データにおける次元jの次元iに対する因果の強さを、偏微分を用いて計算する代わりに、時刻p前の次元jに微小量Δを与えた際の次元iの予測値の変化量x'_{t,i} - x_{ t,i }（x'_{t,i}は時刻p前の次元jに微小量Δを与えた際の次元iの予測値，x_{t,i}は与えなかった場合の予測値，p=1,…,p_max）を用いて計算する。
【００７１】
実施例３では、実施例１において、因果推定装置１００は、監視データにおける次元jの次元iに対する因果の強さを、偏微分を用いて計算する代わりに、f_i(x_t-τ,x_t-τ+1,…,x_t-1)においてx_{t-p,j}を含む項{a_1 g_1(x_t-τ,x_t-τ+1,…,x_t-1), …,a_M g_M(x_t-τ,x_t-τ+1,…,x_t-1)}（a_mは定数，g_mは関数）のみに着目し、定数a_mや関数g_mを用いて計算する。
【００７２】
実施例４では、実施例１において、因果推定装置１００は、非線形回帰を行う際に，スパースなモデリングをおこなうために、学習時にスパース項L(x_t-τ,x_t-τ+1,…,x_t-1)を考慮して学習を実施する。
【００７３】
実施例５では、実施例１において、因果推定装置１００は、非線形回帰をニューラルネットワークを用いて行う。
【００７４】
実施例６では、実施例１、実施例２、あるいは実施例３において、因果推定装置１００は、学習した非線形回帰モデルにおける各パラメータθ_1,…,θ_Kの重要度F_1,…,F_Kを定義し、因果の強さの計算時にパラメータの重要度についても考慮した計算を行う。
【００７５】
上述したとおり、本発明の実施の形態により、時系列多次元数値ベクトルのデータを入力する入力部と、入力された時系列多次元数値ベクトルのデータを用いて、ある時刻のデータを過去の時刻のデータから予測する非線形回帰モデルを学習する回帰モデル学習部と、前記非線形回帰モデルを用いて、前記時系列多次元数値ベクトルのデータにおける次元jの次元iに対する因果の強さを計算する因果推定部と、前記因果推定部により計算された因果の強さを出力する出力部とを備えることを特徴とする因果推定装置が提供される。
【００７６】
前記因果推定部は、例えば、前記非線形回帰モデルにおける時刻t-pの次元jの誤差項の変動が時刻tの次元iに与える影響を用いて前記因果の強さを計算する、又は、時刻t-pにおける次元jに微小量Δが与えられた際の時刻tにおける次元iの前記非線形回帰モデルによる予測値と、当該微小量が与えられなかった場合の時刻tにおける次元iの前記非線形回帰モデルによる予測値との誤差を用いて前記因果の強さを計算する、又は、前記非線形回帰モデルによる次元iの予測値における時刻t-pの次元jの値を含む項を用いて前記因果の強さを計算する。
【００７７】
前記回帰モデル学習部は、スパース項を考慮したスパースモデリングにより前記非線形回帰モデルの学習を行うこととしてもよい。
【００７８】
前記回帰モデル学習部は、ニューラルネットワークを用いて前記非線形回帰モデルの学習を行うこととしてもよい。
【００７９】
前記回帰モデル学習部は、前記非線形回帰モデルの計算時に、当該非線形回帰モデルの各パラメータの重要度を計算し、前記因果推定部は、前記重要度を用いて前記因果の強さを計算することとしてもよい。
【００８０】
また、本発明の実施の形態により、因果推定装置が実行する因果推定方法であって、時系列多次元数値ベクトルのデータを入力する入力ステップと、入力された時系列多次元数値ベクトルのデータを用いて、ある時刻のデータを過去の時刻のデータから予測する非線形回帰モデルを学習する回帰モデル学習ステップと、前記非線形回帰モデルを用いて、前記時系列多次元数値ベクトルのデータにおける次元jの次元iに対する因果の強さを計算する因果推定ステップと、前記因果推定ステップにより計算された因果の強さを出力する出力ステップとを備えることを特徴とする因果推定方法が提供される。
【００８１】
また、本発明の実施の形態により、コンピュータを、上記の因果推定装置における各部として機能させるためのプログラムが提供される。
【００８２】
以上、本実施の形態について説明したが、本発明はかかる特定の実施形態に限定されるものではなく、特許請求の範囲に記載された本発明の要旨の範囲内において、種々の変形・変更が可能である。
【符号の説明】
【００８３】
１００  因果推定装置
１０１  入力部
１０２  記憶部
１０３  因果推定部
１０４  回帰モデル学習部
１０５  出力部
１５０  ドライブ装置
１５１  記録媒体
１５２  補助記憶装置
１５３  メモリ装置
１５４  ＣＰＵ
１５５  インターフェース装置
１５６  表示装置
１５７  入力装置
【図１】
2019144779.tif 000013
【図２】
2019144779.tif 000014
【図３】
2019144779.tif 000015
【図４】
2019144779.tif 000016
【図５】
2019144779.tif 000017
【図６】
2019144779.tif 000018
【図７】
2019144779.tif 000019