説明

確率モデル更新システム、確率モデル更新装置、確率モデル更新方法およびプログラム

【課題】3以上の次元を有する計測情報を高精度に確率モデル化し、かつ、ストリーム処理によって確率モデルを逐次更新することを可能にする。
【解決手段】受信部11は、2人以上のユーザの行動を検出する2以上のセンサから、ユーザの行動を示すセンサ情報を受信し、センサ情報記億部12は、センサ情報を記憶する。モデル化部13は、時間を変数とし、ユーザとセンサを変化しないと仮定し、センサ情報記億部12が記憶するセンサ情報を確率モデル化する。モデル記憶部14は、確率モデルを記憶する。モデル更新部15は、モデル化部13が確率モデル化に用いなかったセンサ情報に基づいて、モデル記憶部14が記憶する確率モデルを更新する。欠損値予測部16および異常値検出部17は、センサ情報と確率モデルとを比較して、それぞれ欠損値を予測し、異常値を検出する。

【発明の詳細な説明】
【技術分野】
【0001】
本発明は、取得したデータに基づいて確率モデルを更新する確率モデル更新システム、確率モデル更新装置、確率モデル更新方法およびプログラムに関する。
【背景技術】
【0002】
近年、安価なセンサの出現や、ネットワーク環境の普及によって、多種多様なセンサを容易に用いることができるようになった。それに伴い、顧客の購買行動の分析や地域の交通状況を分析などについて、マルチセンサの情報に基づいてモデル化する技術が盛んに研究されている。業務のマネジメントの観点では、主観的に捉えきれない変化を検出する方法として、センサデータを使ったマネジメント支援方法が検討されている。
【0003】
しかし、非定型業務においては、ミスやトラブルの種類は多種多様であるため、あらかじめミスやトラブルの種類を想定して検出することは難しい。そのため、業務の定常状態をモデル化しそのモデルとの差を検知することが考えられる。
【0004】
非特許文献1では、各センサ情報の分布を個別に仮定することができ、かつ多様な分布(ガウス分布、ベルヌーイ分布、ポアソン分布、指数分布など)を扱える。そのため、従来の手法に比べてマルチセンサ環境のデータに対して高精度のモデル化が可能であり、欠損値予測や異常値検出の性能が高い。
【先行技術文献】
【非特許文献】
【0005】
【非特許文献1】K.Hayashi, T.Takenouchi, T.Shibata, Y.Kamiya, D.Kato, K.Kunieda, K.Yamada and K.Ikeda, "Exponential Family Tensor Factorization for Missing-Values Prediction and Anomaly Detection", IEEE International Conference on Data Mining 2010.
【発明の概要】
【発明が解決しようとする課題】
【0006】
しかしながら、非特許文献1の技術は、バッチ処理方式であり、計算量が膨大である。たとえば、1年などの長期間かつ大規模なデータを用いてモデル化を行うには、大規模コンピュータ・クラスタを用いて数日かかる。したがって、業務状態のリアルタイムモニタリングなど、センサデータのストリーム処理が必要な場面では使用できない。
【0007】
本発明は、上述のような事情に鑑みてなされたもので、3以上の次元を有する計測情報を高精度に確率モデル化し、かつ、ストリーム処理によって確率モデルを逐次更新することを可能にする確率モデル更新システム、確率モデル更新装置、確率モデル更新方法およびプログラムを提供することを目的とする。
【課題を解決するための手段】
【0008】
本発明の第1の観点に係る確率モデル更新装置は、
2以上の対象のそれぞれについて2以上の特性を計測する計測装置から、3以上の次元を有する計測情報を受信する受信手段と、
前記計測情報を記憶する記憶手段と、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化手段と、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新手段と、
を備えることを特徴とする。
【0009】
本発明の第2の観点に係る確率モデル更新システムは、
2以上の対象のそれぞれについて2以上の特性を計測する計測装置と、前記計測装置と接続する確率モデル更新装置とで構成される確率モデル更新システムであって、
前記計測装置は、それぞれ、
前記対象のそれぞれについて2以上の特性を計測する計測手段と、
前記計測手段が計測した前記特性から3以上の次元を有する計測情報を生成する生成手段と、
前記計測情報を送信する送信手段と、を備え、
前記確率モデル更新装置は、
前記計測装置から、前記計測情報を受信する受信手段と、
前記計測情報を記憶する記憶手段と、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化手段と、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新手段と、
を備えることを特徴とする。
【0010】
本発明の第3の観点に係る確率モデル更新方法は、
2以上の対象のそれぞれについて2以上の特性を計測する計測装置がそれぞれ実行する
前記対象のそれぞれについて2以上の特性を計測する計測ステップと、
前記計測ステップで計測した前記特性から3以上の次元を有する計測情報を生成する生成ステップと、
前記計測情報を送信する送信ステップと、
確率モデル更新装置が実行する
前記計測装置から、前記計測情報を受信する受信ステップと、
前記計測情報を記憶する記憶ステップと、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化ステップと、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新ステップと、
を備えることを特徴とする。
【0011】
本発明の第4の観点に係るプログラムは、コンピュータを、
2以上の対象のそれぞれについて2以上の特性を計測する計測装置から、3以上の次元を有する計測情報を受信する受信手段、
前記計測情報を記憶する記憶手段、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化手段、および、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新手段、
として機能させることを特徴とする。
【発明の効果】
【0012】
本発明によれば、3以上の次元を有する計測情報を高精度に確率モデル化し、かつ、ストリーム処理によって確率モデルを逐次更新することが可能になる。
【図面の簡単な説明】
【0013】
【図1】本発明の実施の形態に係る確率モデル更新システムの構成例を示すブロック図である。
【図2】実施の形態に係る確率モデル更新装置の機能構成例を示す図である。
【図3】指数型分布における関数ψとその導関数をまとめた表を示す図である。
【図4】M次のデータテンソルのための推定アルゴリズムの疑似コードを示す図である。
【図5】ETFのオンラインアルゴリズムの疑似コードを示す図である。
【図6】距離ベースの外れ値の概念を示す図である。
【図7】実施の形態に係る確率モデル更新装置のモデル化処理およびモデル更新処理の動作の一例を示すフローチャートである。
【図8】本発明の実施の形態に係る確率モデル更新装置のハードウェア構成例を示すブロック図である。
【発明を実施するための形態】
【0014】
本実施の形態では、本発明をユーザの行動を検出するマルチセンサシステムに適用する。本実施の形態では、対象はユーザであり、計測装置はセンサであり、計測情報はセンサ情報である。また、特性はユーザの各種行動であり、計測情報の次元は、「ユーザ」、「特性」、「時間」である。
【0015】
以下、本発明を実施するための形態について図を参照して詳細に説明する。なお図中、同一または同等の部分には同一の符号を付す。
【0016】
図1は、本発明の実施の形態に係る確率モデル更新システムの構成例を示すブロック図である。確率モデル更新システム100は、ネットワーク上の確率モデル更新装置1と複数のセンサ2とで構成される。
【0017】
センサ2は、ユーザの行動を検出して、確率モデル更新装置1に送信する。センサ2は、たとえば、ユーザの位置を検出するセンサや、ユーザの動きを検出するセンサなどに加え、ユーザが使用する端末なども含む。端末の場合は、ユーザが端末を操作した情報をユーザの行動を示すセンサ情報として確率モデル更新装置1に送信する。確率モデル更新装置1は、それぞれのセンサ2からセンサ情報を受信し、これに基づいてユーザの行動の確率モデル化を行う。また、確率モデル更新装置1は、定期的に確率モデルを更新する。
【0018】
図2は、実施の形態に係る確率モデル更新装置の機能構成例を示す図である。確率モデル更新装置1は、受信部11、センサ情報記憶部12、モデル化部13、モデル記憶部14、モデル更新部15、欠損値予測部16および異常値検出部17を備える。
【0019】
受信部11は、複数のセンサ2からセンサ情報を受信する。
【0020】
センサ情報記憶部12は、受信部11が受信したセンサ情報を記憶する。
【0021】
モデル化部13は、センサ情報記憶部12が記憶した所定の期間のセンサ情報に基づいて、確率モデルを生成するモデル化処理を行う。モデル化処理の詳細は、後述する。
【0022】
モデル記憶部14は、モデル化部13のモデル化処理によって生成された確率モデルを記憶する。
【0023】
モデル更新部15は、モデル記憶部14が記憶する確率モデルを定期的に更新するモデル更新処理を行う。モデル更新処理の詳細は、後述する。
【0024】
欠損値予測部16は、受信部11が受信したセンサ情報とモデル記憶部14が記憶する確率モデルとに基づいて欠損値の予測を行う。欠損値の予測方法については、後述する。
【0025】
異常値検出部17は、センサ情報記憶部12が記憶するセンサ情報とモデル記憶部14が記憶する確率モデルとに基づいて異常値の検出を行う。異常値の検出方法については、後述する。
【0026】
ここで、モデル化部13およびモデル更新部15がそれぞれ行う、モデル化処理およびモデル更新処理について説明する。モデル化処理およびモデル更新処理は、指数分布族のテンソル因子分解(ETF:Exponential family Tensor Factorization)のオンラインアルゴリズムによって行われる。テンソルとは、多次元の配列である。M次元の配列をM次のテンソルと呼ぶ。特別な場合として、1次、2次のテンソルはそれぞれベクトルおよび行列である。
【0027】
モデル化処理およびモデル更新処理の説明を行うにあたって、まず、代表的なテンソル因子分解のアルゴリズムと、指数分布族のテンソル因子分解のバッチアルゴリズムとを説明する。以下では、理解を容易にするため3次のテンソルで考える。しかし、結果はより高次のテンソルに容易に拡張できる。以下、テンソルは、アンダーラインを引いた文字で表記する。
【0028】
まず、代表的なテンソル因子分解としてTucker分解を説明する。をD×D×D次の観測値テンソルとする。観測入力D≡Dを有する。Tucker分解は、を分解したコアテンソルと、3つの因子行列U(m)(m=1,2,3)を生成する方法を提供する。の(i,j,k)番要素を以下の式で表す。
【数1】

【0029】
ここで、zqrsの(q,r,s)番要素、uiq(m)は因子行列U(m)の(i,q)番要素、εijkは観測ノイズである。行列U(m)は、m次のの上の相関関係の構造を表す。Tucker分解は、二乗誤差の和(εijk)を最小にすることによってパラメータおよび{U(m)}を推定する。高次特異値分解(HOSVD:Higher-order singular value decomposition)は、Tucker分解の解法の1つである。以下では、テンソルの要素を行列の形式に並べた展開(行列形式)を考える。結果の行列もまた、m次元のに従った構造を保持する。詳細は、T. G. Kolda and B. W. Bader, “Tensor decompositions and applications," Sandia National Laboratories, Albuquerque, NM and Livermore, CA, Tech. Rep., 2007.を参照。HOSVDは、因子行列U(m)を展開テンソルX(m)の左特異ベクトルとして推定する。
【0030】
後の便宜のため、ベクトルと行列を使って、数1を以下の数2の形に書く。数式中のアルファベットまたはギリシャ文字の上部に矢印が付いている記号は、ベクトルを表す。以下、これらの上部に矢印が付いている文字を○→と表記する。同様に、数式中の上部に〜が付いている文字を○〜と表記し、上部に^が付いている文字を○^と表記する。それぞれ、○にはアルファベットまたはギリシャ文字が入る。x→は、要素が任意に並べ替えられたの要素で与えられるD次元ベクトルである。これは、行列のベクトル化になぞらえている。
【数2】

【0031】
ここで、○の中に×の記号は行列のクロネッカ積を表す。また、Wは、ベクトル化されたノイズであるε→のもとで、ベクトル化されたコアテンソルであるz→の線形マッピングを表す。この形式では、Tucker分解は、標準線形モデルとして見られる。標準線形モデルとの主な違いは、Wがクロネッカ積を通して{U(m)}で表現されることである。
【0032】
行列表記を用いて、数1は数3に書き直せる。
【数3】

【0033】
(1)とE(1)は、テンソルとノイズに関係する展開テンソルである。X(2)およびX(3)も数3と同じように書き直せる。後述する{U(m)}に関する期待対数尤度の勾配を導くためにこの表現を用いる。
【0034】
続いて、指数型分布族のテンソル因子分解のバッチアルゴリズムを説明する。Tucker分解では、数1において誤差εijkの二乗和を最小にすることによって、パラメータU(1)、U(2)、U(3)およびを推定する。この推定は、確率の観点では、球状ガウスノイズεの想定の下で、最尤推定解として考えることができる。しかしながら、この仮定はデータが異質な分布を有する場合は適当ではない。この問題を扱うために、モデルを数4に示すように一般化する。
【数4】

【0035】
ここで、θ→≡Wz→は、自然母数と呼ばれるD次元のベクトルであり、hはx→の想定された指数型分布族を特定する指標である。分布関数Exponは指数型分布族である。
【数5】

【0036】
関数F(x)は、基底因数、ψ(θ)は対数分配関数すなわち数6である。
【数6】

【0037】
すべてのExponを等方ガウス分布とする場合は、対数尤度は従来のTucker分解の損失関数と等価である。数6の両辺を微分することによって、導関数ψ’が、Expon(x|θ)の上で自然母数θから条件付き期待値xへの写像であることがわかる。また、2次導関数ψ”を用いて分散を計算することができる。ψは凸関数であり、Exponは対数凹関数である。
【0038】
図3は、指数型分布における関数ψとその導関数をまとめた表を示す図である。
【0039】
指数型分布族のテンソル因子分解で取り扱うモデルは一般化した線形モデル(P. McCullagh and J. A. Nelder, Generalized Linear Models, Second Edition, 1989.を参照)に強く関係している。たとえば、二値データに対してExpon(x|θ)としてベルヌーイ分布を選べば、ψ’はシグモイド関数になる。その場合、モデルはロジスティック回帰と等価であり、z→、x→、およびWはそれぞれ、入力、出力および回帰係数に等しい。
【0040】
ここでのモデルの鍵となる想定は、数4における指標hによって制御される属性の不均質性である。PCAなどの既存の分解モデルとは異なり、このモデルでは指標h(d=1,...,D)を変えることによって、指数型分布族(数5)から任意に分布関数を選択することができる。これによって、より柔軟なデータのモデリングが可能である。
【0041】
以下では、z→をx→に関係する隠れた変数として扱い、一般化された線形の隠れた変数モデル(P. Huber, E. Ronchetti, and M. P. V. Feser, “Estimation of generalized linear latent variable models," Journal of the Royal Statistical Society. Series B, Statistical Methodology, pp. 893?908, 2004.を参照)の概念に従って、z→に標準ガウス事前分布N(0,1)を想定する。また、各因子行列{U(m)}に球状ガウス事前分布の精度αを追加する。それは、二乗正規化に等価である。最終的に、同時対数尤度Lは次の数7のように書ける。
【数7】

【0042】
第1項および第2項は、自然母数θ→=Wz→で尤度(数4)に対応する。第3項と第4項はそれぞれ、z→および{U(m)}の事前分布である。このモデルを指数型分布族のテンソル因子分解という。
【0043】
パラメータの推定には、以下のEMアルゴリズムを用いる。
【0044】
まず、ガウス分布q(z→)≡N(z→|z→,Σ)にラプラス近似を適用して、事後分布p(z→|x→)を近似する。ラプラス近似では、z→は最大事後確率(MAP)すなわちz→の事後確率分布のモード(最頻値)であり、Σは、z→におけるヘシアンの負逆数である。z→を探索するために、数8に示す勾配とヘシアンによる、勾配法を用いる。
【数8】

【0045】
ここで、ψ’→≡(ψ’(θ→),...,ψ’hD(θ→))、およびΨ”→≡diag(ψ”→)である。負のヘシアンは正定行列なので、任意のWに対して、大域で最大値を見いだすことができる。以上の処理をEステップと呼ぶ。
【0046】
次に、U(m)の周辺MAP推定をとる。以下に示す数9の期待対数尤度の近似を考える。
【数9】

【0047】
そして、数3を用いての勾配を計算する。ここに、展開されたコアテンソルZ(1)と自然母数Θを導入する。
【数10】

【0048】
(1)の疑似逆行列をU(1)と書けば、期待対数尤度(数9)の勾配は、以下の数11のように書ける。
【0049】
【数11】

【0050】
Ψ’(1)は、ψ’→の展開テンソルであり、数12である。
【0051】
【数12】

【0052】
ここに、微分と積分の演算を交換できると仮定する。先に指摘したように、勾配におけるE[Ψ’(1)(Θ(1)]の期待値の計算は、関数ψの非線形性のため一般に扱いにくい。
【0053】
期待値を近似するために、GPの積分と微分に関して有用な性質を適用する。共分散関数がガウスカーネルであるGPによって、ψ(θ)の近似をψ〜(θ)と書く。第1に、θがガウス分布のランダム変数なら、任意のnについて期待値E[ψ〜(θ)θ]は、数13を用いて解析的に解ける。
【数13】

【0054】
第2に、ψのn次導関数dψ/dθは、他のGPとしてのψ〜を用いて直接的に近似できる。
【0055】
前節でz→の事後分布をガウス分布として近似したので、θ→=Wz→の事後分布は、Wを定数とすると、ガウス分布である。そこで、カーネル(数13の1つ目の式)の期待値の結果とGP(数14)の導関数を結合することによって、期待値E[Ψ’(1)(Θ(1)]を計算できる。
【数14】

【0056】
同様の方法で、U(2)およびU(3)についても、勾配を計算することができる。{U(m)|m=1,2,3}の更新のために交互最適化を用いる。すなわち、{U(n)|n≠m}を固定してU(m)に関してLを最大化することを、指標mを変化させて繰り返す。{U(m)}に関する(局所)最適解を得るために疑似ニュートン法を用いる。以上の処理をMステップと呼ぶ。
【0057】
EM反復で収束したのち、観測入力Dの下で自然母数θ→の推定値を得る。そして、予測分布E[x→|D]の平均によって欠損値x→を推定できる。モデルパラメータであるz→の事後分布による周辺化は、ベイズ予測分布の計算に必要だが、扱いにくい。しかしながらこの枠組みでは、ψ’(θ0d→)は与えられたθ0d→の条件付き平均であり、したがって、z→の事後分布でψ’〜(θ0d→)を周辺化することにより予測分布の平均を近似できる。
【数15】

【0058】
再び、期待値E[ψ’〜(θ0d→)]はGPの事後分布を用いて解析的に解ける。また、予測分布の変数または他の高次モーメントを導くことができる。ここでは、その説明を省略する。
【0059】
前述の指数型分布族のテンソル因子分解の説明において、ガウス分布の分散を1と想定したので、ガウス分布でデータを正規化する。すなわち、パラメータの推定の前にデータサンプルの標準偏差で各部を正規化する。EMアルゴリズムを始める前に、HOSVDを用いて{U(m)}を初期化する。M次のデータテンソルのための推定アルゴリズムの疑似コードを図4に示す。図4の疑似コードに現れるC、s、yおよびγは、カーネル行列、入力、出力、およびGPのハイパーパラメータである。
【0060】
MステップにおけるGP近似のために、トレーニング入力sを設定する。x→の各要素は数4で独立に分布すると仮定しているので、期待値は、独立の期待値の積、すなわち数16のように、因数分解できる。
【数16】

【0061】
従って、ψ(θ→)の各次に入力を共通にできる。また、事後分布が稠密であるように、および/または、関数ψおよびψ’が大きい値をとるように、トレーニング入力によって領域をカバーすることが重要である。従って、Nをランダムに選択して、重みψ(θ→)の事後分布の平均θ→=WZ0→の各次元から1つの次元の入力θ(n=1,...,N)をとる。
【0062】
一方、入力がまばらな領域では、GPの特性のためGPの平均はゼロに近い。勾配に基づく最適化を適用するためには、この特性は問題である。なぜなら、コスト関数L(θ)は、θ→∞で発散するからである。この問題を回避するために、GPの入力が領域[min(s),max(s)]からはずれる場合は、GPの平均の代わりにバリア関数を用いる。E[ψ(θ)]およびE[ψ’(θ)]をバリア関数として用いるために、それらにゼロ次デルタ法を適用する。
【0063】
Eステップにおける主たる複雑性は、ラプラス近似における共分散(ヘシアンの逆数)の計算である。それは、K=dim(z→)としてKのオーダーである。MステップのGP近似では、式(12)で定義されるカーネル行列Cの逆行列は、NおよびN’をGPにおける観測の通常の入力とその微分の数として、(N+N’)のオーダー(の回数の計算)が必要である。EM反復の間にこの最も費用のかかる計算を必要とするのは1回だけである。さらに、最大化処理の期待対数尤度の勾配における行列乗算のために、(D(N+N’))+Kのオーダーの回数の計算を必要とする。
【0064】
前述のTucker分解の説明で述べたように、因子行列U(m)はm次の観測テンソルの低次元特徴量と見られる。異常な値を含むテンソルからパラメータを推定すれば、因子行列の(異常値に)対応する部分は、通常の部分に比べて外れ値として捕捉される。因子行列を外れ値検出の入力として用いれば、全データ集合に影響する高いインパクトの異常を発見できる。それは、独立モデルアプローチより強力な検出を達成できる。さらに、異質なテンソルにおける異常値は、pTuckerよりもETFの特性空間で弁別される異常として現れる。なぜなら、論理的に、ETFは指数型分布族の適切な想定の下で、データの標準的な部分を抽出するからである。
【0065】
外れ値を検出する最も単純な方法は、入力間の距離を用いることである。入力間の距離として、"an object O in a dataset T is a DB(p,D) outlier if at leaset fraction p of the objects in T lies greater than destance D frm O." (E. M. Knorr, R. T. Ng, and V. Tucakov, "Distance-based outliers: algorithms and applications," The VLDB Journal, vol. 8, no. 3-4, pp. 237-253, 2000.)による定義を採用することができる。
【0066】
ここで、モデル化部13およびモデル更新部15がそれぞれ行うモデル化処理およびモデル更新処理を実行するためのETFのオンラインアルゴリズムついて説明する。ETFのオンラインアルゴリズムは、バッチアルゴリズムをオンライン処理に拡張したものである。オンライン処理では、データテンソルを複数のスライスに分割し、パラメータを逐次推定する。この拡張によってリアルタイムにデータ処理が可能になり、パラメータ推定の精度に必要な計算費用を削減する。
【0067】
データテンソルのl番の次元の量が相当大きく、バッチアルゴリズムでは計算不可能であることを想定する。そのようなテンソルに対して、まず、データテンソルをl番の次元に沿ってスライスし、をD1×...Dl−1×Dl〜×Dl+1×...DLのテンソル〜と、D1×...Dl−1×(Dl−Dl〜)×Dl+1×...DLのテンソル^に分割する。〜をバッチアルゴリズムで処理できるように、充分小さくDl〜(<Dl)を選ぶ。
【0068】
本実施の形態では、3次のデータテンソルを想定しており、l番の次元は、「時間」である。それ以外の次元は、「特性」と「ユーザ」である。「時間」は恒常的に変化するが、「特性」と「ユーザ」は変化しないと仮定する。なお、次元はこれに限らない。たとえば、ユーザの「位置」や行動の「対象」を次元としてもよい。また、変数とする次元(l番の次元)は「時間」に限らず、行動の「対象」や「ユーザ」としてもよい。
【0069】
パラメータ、{U|n≠l}および{ul,i|i=1,...,Dl〜}を、分割されたテンソル〜についてバッチアルゴリズムで推定したのちに、残りのU^について推定する。周辺尤度は以下の数17および数18のように分解できる。
【数17】


【数18】

【0070】
のi番の行ベクトルul,iは、観測x(l)にのみ依存するので、i=1,...,Dlについて、ul,iを逐次推定できる。ul,iと同様に事後共分散Σを逐次アップデートする。なぜなら、数19に示すように、z→が固定されていても事後確率分布のヘシアンは、ul,iに伴って変化するからである。
【数19】

【0071】
l,iの最大化とΣのアップデートは、バッチアルゴリズムのMステップおよびEステップにそれぞれ対応する。このアルゴリズムを図5に示す擬似コードに要約する。
【0072】
図5の擬似コードにおいて、initializeからfor i=1,...,Dの前までがモデル化処理に相当し、for i=1,...,Dとend forの間の部分がモデル更新処理に相当する。
【0073】
オンラインアルゴリズムの重要な利点は、バッチアルゴリズムより大幅に計算費用が低いことである。交互アップデートを必要としないので、収束する速さは相当に速い。他の利点は、時系列データのオンライン処理である。このアルゴリズムによれば、x(l)がi番の時刻ごとに観測されるとき、ul,iをリアルタイムに推定できる。
【0074】
コンピュータによる処理を比較することは難しいが、オンラインアルゴリズムはバッチアルゴリズムより著しく速い。オンラインアルゴリズムでは1つの変数Uしかないので、バッチあるゴリズムで用いる座標傾斜法を必要としない。そのことは、収束のための計算費用を劇的に削減する。
K>D\lの場合は、逆行列の補助定理
Σ=(Σi−1−1+WΨ”W−1
=Σi−1−Σi−1(Ψ”+WΣi−1−1Σi−1
を用いて、効率的に事後共分散Σiをアップデートできる。それによって、コンピュータ処理は、KのオーダーからD\lのオーダーに減少する。
【0075】
アルゴリズムが収束すると、観測入力x→の下におけるパラメータz→および{U}の推定値を得る。これらのパラメータを主に2つの目的、すなわち、欠損値の予測と異常値の検出に使用する。自然母数θ→が与えられたとして、ベイズ予測分布E[x→|x→]の平均によって、欠損値、すなわち観測値の指標の集合Iに含まれないdについてのx→を予測する。因子行列Uは、データのl番の次元における異常値を検出するのに用いられる。
【0076】
以下に、欠損値予測部16が行う、欠損値予測処理について説明する。
【0077】
事後分布によるz→の周辺化は、ベイズ予測分布の計算で必要になるが、扱いにくい。その代わり、GPの予測分布の平均m’の一階導関数により非線形関数ψ’(θ)を近似し、数20で与えられる予測分布の平均を近似する。
【数20】

【0078】
第2行の変換には、数21の関係を用いている。
【数21】


【数22】


以下の定理により、数22の方程式において、p=1およびq=0の場合の期待値Eq[m’h(d)(θ→)]の解析的形式が与えられる。また、以下の定理を用いて、ETFの予測分布の分散(p=2)またはそれ以上の高次のモーメント(p≧3)が得られる。
【0079】
定理:m(p)を、共分散関数がガウスカーネルであるGPの予測平均関数のp次導関数とする。任意の正の整数p,q≧0について、平均μで分散がσのガウス分布であるp(x)の期待値Ep(x*)[x(p)(x)]は、p,q,μおよびσの関数として陽に表される。
【0080】
以下に、異常値検出部17が行う異常値検出処理について説明する。
【0081】
因子行列Uは、観測テンソルのl番の次元の低次元特徴量である。異常な値を含むテンソルからパラメータを推定すれば、因子行列の(異常値に)対応する部分は、通常の部分に比べてはずれ値として捕捉される。因子行列をはずれ値検出の入力として用いることによって、観測ノイズの影響なしに本質的な異常値を探索できる。
【0082】
外れ値を発見するにはいくつかの方法がある。ここでは、前述のKnorr 他が提案する距離ベースの外れ値を採用する。
【0083】
定義:少なくともTのオブジェクトの部分pがOからの距離rより大きい位置にあれば、
データセットTにおけるオブジェクトOはDB(p,r)外れ値である。
【0084】
図6は、距離ベースの外れ値の概念を示す図である。p=0.995、およびD=1000、すなわち図6の点の数が1000であるとして、点Oを中心とする半径rの多次元球に多くても5つした他の点が含まれない場合、点Oは外れ値として検出される。Tucker分解では、{U}およびの間の尺度は、誤りを引き起こす。すなわち、αUおよび(1/α)は同じΘになる。従って、距離ベースの外れ値に適用する前にUを列方向に正規化する必要がある。たとえば、ユークリッド距離をとる場合、正規化はコサイン距離に類似する。
【0085】
ここでは、たとえば多重ネットワークにおける異常ノードのような異常値のみに着目する。言い換えれば、データテンソルの個別の要素の異常値ではなく、モードの区別できる次元の異常値に着目する。
【0086】
図7は、実施の形態に係る確率モデル更新装置のモデル化処理およびモデル更新処理の動作の一例を示すフローチャートである。確率モデル更新装置1の受信部11は、各センサ2からセンサ情報を受信する(ステップS11)。センサ情報記憶部12は、受信部11が受信したセンサ情報を記憶する(ステップS12)。モデル化部13は、所定の時間(たとえば、1ヶ月)が経過したか否かを判定する(ステップS13)。所定の時間が経過していない場合(ステップS13;NO)、ステップS11に戻り、ステップS11〜ステップS13を繰り返す。所定の時間が経過した場合(ステップS13;YES)、モデル化部13は、センサ情報記憶部12が記憶するセンサ情報に基づいて、ユーザの行動を確率モデル化するモデル化処理を実行する(ステップS14)。モデル記憶部14は、モデル化部13がモデル化処理によって生成された確率モデルを記憶する(ステップS15)。
【0087】
次に、モデル更新部15は、時刻iがDになったか否かを判定する(ステップS16)。時刻iがDになっていない場合(ステップS16;NO)、モデル更新部15は、モデル更新処理を実行し(ステップS17)、ステップS16およびステップS17を繰り返す。時刻iがDになった場合(ステップS16;YES)、処理を終了する。
【0088】
なお、上記のフローチャートでは、ステップS11〜ステップS15のモデル化処理は、1度だけ実行しているが、定期的に(たとえば1年ごと)に実行してもよいし、センサやユーザに変更があるごとに実行してもよい。
【0089】
以上説明したように、3以上の次元を有する計測情報を高精度に確率モデル化し、かつ、ストリーム処理によって確率モデルを逐次更新することが可能になる。たとえば、ユーザの行動の確率モデルを逐次更新し、これに基づいて異常値を検出することにより、業務におけるミスやトラブルの検出精度の向上が期待できる。
【0090】
上述の実施の形態では、本発明をユーザの行動を検出するマルチセンサシステムに適用した例を示した。しかし本発明はこれに限らない。たとえば、気象を対象とする場合、計測装置は、「風速」、「気温」、「湿度」といった特性を計測し、「特性」、「時間」、「計測地点」といった次元の気象情報を計測情報として生成する。あるいは、交通車両を対象とする場合、計測装置は、「速さ」、「台数」、「向き」といった特性を計測し、「特性」、「計測地点」、「時間」といった次元の交通情報を計測情報として生成する。
【0091】
図8は、本発明の実施の形態に係る確率モデル更新装置のハードウェア構成例を示すブロック図である。確率モデル更新装置1は、図8に示すように、制御部31、主記憶部32、外部記憶部33、操作部34、表示部35および送受信部36を備える。主記憶部32、外部記憶部33、操作部34、表示部35および送受信部36はいずれも内部バス30を介して制御部31に接続されている。
【0092】
制御部31はCPU(Central Processing Unit)等から構成され、外部記憶部33に記憶されている制御プログラム30に従って、モデル化部13、モデル更新部15、欠損値予測部16および異常値検出部17の各処理を実行する。
【0093】
主記憶部32はRAM(Random-Access Memory)等から構成され、外部記憶部33に記憶されている制御プログラム30をロードし、制御部31の作業領域として用いられる。
【0094】
外部記憶部33は、フラッシュメモリ、ハードディスク、DVD−RAM(Digital Versatile Disc Random-Access Memory)、DVD−RW(Digital Versatile Disc ReWritable)等の不揮発性メモリから構成され、確率モデル更新装置1の処理を制御部31に行わせるためのプログラムをあらかじめ記憶し、また、制御部31の指示に従って、このプログラムが記憶するデータを制御部31に供給し、制御部31から供給されたデータを記憶する。センサ情報記憶部12およびモデル記憶部14は、外部記憶部33に構成される。
【0095】
操作部34はキーボードやテンキーなどのポインティングデバイス等と、キーボードやテンキー等を内部バス30に接続するインタフェース装置から構成されている。さらに操作部34は音声入力装置を備えることとしてもよい。モデル化処理を実行するまでの所定の時間などの入力を行う場合は、操作部34を介して、指示が制御部31に供給される。
【0096】
表示部35は、CRT(Cathode Ray Tube)またはLCD(Liquid Crystal Display)などから構成され、制御部31から送られてきた情報を表示する。モデル化処理を実行するまでの所定の時間などの入力を行う場合は、表示部35は、操作画面を表示する。
【0097】
送受信部36は、通信ネットワークに接続する網終端装置または無線通信装置、およびそれらと接続するシリアルインタフェースまたはLAN(Local Area Network)インタフェースから構成されている。受信部11は、送受信部36を介して通信ネットワークに接続し、センサ2からセンサ情報を受信する。
【0098】
図2に示す確率モデル更新装置1の受信部11、センサ情報記憶部12、モデル化部13、モデル記憶部14、モデル更新部15、欠損値予測部16および異常値検出部17の処理は、制御プログラム30が、制御部31、主記憶部32、外部記憶部33、操作部34、表示部35および送受信部36などを資源として用いて処理することによって実行する。
【0099】
その他、前記のハードウェア構成やフローチャートは一例であり、任意に変更および修正が可能である。
【0100】
制御部31、主記憶部32、外部記憶部33、操作部34、内部バス30などから構成される時刻同期処理を行う中心となる部分は、専用のシステムによらず、通常のコンピュータシステムを用いて実現可能である。たとえば、前記の動作を実行するためのコンピュータプログラムを、コンピュータが読み取り可能な記録媒体(フレキシブルディスク、CD−ROM、DVD−ROM等)に格納して配布し、当該コンピュータプログラムをコンピュータにインストールすることにより、前記の処理を実行する確率モデル更新装置1を構成してもよい。また、インターネット等の通信ネットワーク上のサーバ装置が有する記憶装置に当該コンピュータプログラムを格納しておき、通常のコンピュータシステムがダウンロード等することで確率モデル更新装置1を構成してもよい。
【0101】
また、確率モデル更新装置1の機能を、OS(オペレーティングシステム)とアプリケーションプログラムの分担、またはOSとアプリケーションプログラムとの協働により実現する場合などには、アプリケーションプログラム部分のみを記録媒体や記憶装置に格納してもよい。
【0102】
また、搬送波にコンピュータプログラムを重畳し、通信ネットワークを介して配信することも可能である。たとえば、通信ネットワーク上の掲示板(BBS, Bulletin Board System)に前記コンピュータプログラムを掲示し、ネットワークを介して前記コンピュータプログラムを配信してもよい。そして、このコンピュータプログラムを起動し、OSの制御下で、他のアプリケーションプログラムと同様に実行することにより、前記の処理を実行できるように構成してもよい。
【0103】
上記の実施形態の一部又は全部は、以下の付記のようにも記載されうるが、以下には限られない。
【0104】
(付記1)
2以上の対象のそれぞれについて2以上の特性を計測する計測装置から、3以上の次元を有する計測情報を受信する受信手段と、
前記計測情報を記憶する記憶手段と、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化手段と、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新手段と、
を備えることを特徴とする確率モデル更新装置。
【0105】
(付記2)
前記計測情報は、前記対象の次元と、前記特性の次元と、時間の次元とを含む前記3以上の次元を有し、
前記確率モデルは、前記対象の次元と、前記特性の次元と、前記時間の次元とを含む前記3以上の次元で示され、
前記モデル化手段は、前記時間の次元を変数とし、前記対象および前記特性を含む他の次元は変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定し、
前記モデル更新手段は、前記パラメータの決定に用いた前記計測情報に含まれない時間の前記計測情報に基づいて、前記確率モデルを更新することを特徴とする付記1に記載の確率モデル更新装置。
【0106】
(付記3)
前記計測情報と前記確率モデルとを比較し、異常値を検出する異常値検出手段をさらに備えることを特徴とする付記1または2に記載の確率モデル更新装置。
【0107】
(付記4)
前記計測情報と前記確率モデルとを比較し、欠損値を予測する欠損値予測手段をさらに備えることを特徴とする付記1ないし3のいずれかに記載の確率モデル更新装置。
【0108】
(付記5)
2以上の対象のそれぞれについて2以上の特性を計測する計測装置と、前記計測装置と接続する確率モデル更新装置とで構成される確率モデル更新システムであって、
前記計測装置は、それぞれ、
前記対象のそれぞれについて2以上の特性を計測する計測手段と、
前記計測手段が計測した前記特性から3以上の次元を有する計測情報を生成する生成手段と、
前記計測情報を送信する送信手段と、を備え、
前記確率モデル更新装置は、
前記計測装置から、前記計測情報を受信する受信手段と、
前記計測情報を記憶する記憶手段と、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化手段と、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新手段と、
を備えることを特徴とする確率モデル更新システム。
【0109】
(付記6)
2以上の対象のそれぞれについて2以上の特性を計測する計測装置がそれぞれ実行する
前記対象のそれぞれについて2以上の特性を計測する計測ステップと、
前記計測ステップで計測した前記特性から3以上の次元を有する計測情報を生成する生成ステップと、
前記計測情報を送信する送信ステップと、
確率モデル更新装置が実行する
前記計測装置から、前記計測情報を受信する受信ステップと、
前記計測情報を記憶する記憶ステップと、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化ステップと、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新ステップと、
を備えることを特徴とする確率モデル更新方法。
【0110】
(付記7)
前記計測情報は、前記対象の次元と、前記特性の次元と、時間の次元とを含む前記3以上の次元を有し、
前記確率モデルは、前記対象の次元と、前記特性の次元と、前記時間の次元とを含む前記3以上の次元で示され、
前記モデル化ステップでは、前記時間の次元を変数とし、前記対象および前記特性を含む他の次元は変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定し、
前記モデル更新ステップでは、前記パラメータの決定に用いた前記計測情報に含まれない時間の前記計測情報に基づいて、前記確率モデルを更新することを特徴とする付記6に記載の確率モデル更新方法。
【0111】
(付記8)
前記確率モデル更新装置が実行する
前記計測情報と前記確率モデルとを比較し、異常値を検出する異常値検出ステップをさらに備えることを特徴とする付記6または7に記載の確率モデル更新方法。
【0112】
(付記9)
前記確率モデル更新装置が実行する
前記計測情報と前記確率モデルとを比較し、欠損値を予測する欠損値予測ステップをさらに備えることを特徴とする付記6ないし8のいずれかに記載の確率モデル更新方法。
【0113】
(付記10)
コンピュータを、
2以上の対象のそれぞれについて2以上の特性を計測する計測装置から、3以上の次元を有する計測情報を受信する受信手段、
前記計測情報を記憶する記憶手段、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化手段、および、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新手段、
として機能させることを特徴とするプログラム。
【符号の説明】
【0114】
1 確率モデル更新装置
2 センサ
11 受信部
12 センサ情報記憶部
13 モデル化部
14 モデル記憶部
15 モデル更新部
16 欠損値予測部
17 異常値検出部
31 制御部
32 主記憶部
33 外部記憶部
34 操作部
35 表示部
36 送受信部
100 確率モデル更新システム

【特許請求の範囲】
【請求項1】
2以上の対象のそれぞれについて2以上の特性を計測する計測装置から、3以上の次元を有する計測情報を受信する受信手段と、
前記計測情報を記憶する記憶手段と、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化手段と、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新手段と、
を備えることを特徴とする確率モデル更新装置。
【請求項2】
前記計測情報は、前記対象の次元と、前記特性の次元と、時間の次元とを含む前記3以上の次元を有し、
前記確率モデルは、前記対象の次元と、前記特性の次元と、前記時間の次元とを含む前記3以上の次元で示され、
前記モデル化手段は、前記時間の次元を変数とし、前記対象および前記特性を含む他の次元は変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定し、
前記モデル更新手段は、前記パラメータの決定に用いた前記計測情報に含まれない時間の前記計測情報に基づいて、前記確率モデルを更新することを特徴とする請求項1に記載の確率モデル更新装置。
【請求項3】
前記計測情報と前記確率モデルとを比較し、異常値を検出する異常値検出手段をさらに備えることを特徴とする請求項1または2に記載の確率モデル更新装置。
【請求項4】
前記計測情報と前記確率モデルとを比較し、欠損値を予測する欠損値予測手段をさらに備えることを特徴とする請求項1ないし3のいずれか1項に記載の確率モデル更新装置。
【請求項5】
2以上の対象のそれぞれについて2以上の特性を計測する計測装置と、前記計測装置と接続する確率モデル更新装置とで構成される確率モデル更新システムであって、
前記計測装置は、それぞれ、
前記対象のそれぞれについて2以上の特性を計測する計測手段と、
前記計測手段が計測した前記特性から3以上の次元を有する計測情報を生成する生成手段と、
前記計測情報を送信する送信手段と、を備え、
前記確率モデル更新装置は、
前記計測装置から、前記計測情報を受信する受信手段と、
前記計測情報を記憶する記憶手段と、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化手段と、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新手段と、
を備えることを特徴とする確率モデル更新システム。
【請求項6】
2以上の対象のそれぞれについて2以上の特性を計測する計測装置がそれぞれ実行する
前記対象のそれぞれについて2以上の特性を計測する計測ステップと、
前記計測ステップで計測した前記特性から3以上の次元を有する計測情報を生成する生成ステップと、
前記計測情報を送信する送信ステップと、
確率モデル更新装置が実行する
前記計測装置から、前記計測情報を受信する受信ステップと、
前記計測情報を記憶する記憶ステップと、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化ステップと、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新ステップと、
を備えることを特徴とする確率モデル更新方法。
【請求項7】
前記計測情報は、前記対象の次元と、前記特性の次元と、時間の次元とを含む前記3以上の次元を有し、
前記確率モデルは、前記対象の次元と、前記特性の次元と、前記時間の次元とを含む前記3以上の次元で示され、
前記モデル化ステップでは、前記時間の次元を変数とし、前記対象および前記特性を含む他の次元は変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定し、
前記モデル更新ステップでは、前記パラメータの決定に用いた前記計測情報に含まれない時間の前記計測情報に基づいて、前記確率モデルを更新することを特徴とする請求項6に記載の確率モデル更新方法。
【請求項8】
前記確率モデル更新装置が実行する
前記計測情報と前記確率モデルとを比較し、異常値を検出する異常値検出ステップをさらに備えることを特徴とする請求項6または7に記載の確率モデル更新方法。
【請求項9】
前記確率モデル更新装置が実行する
前記計測情報と前記確率モデルとを比較し、欠損値を予測する欠損値予測ステップをさらに備えることを特徴とする請求項6ないし8のいずれか1項に記載の確率モデル更新方法。
【請求項10】
コンピュータを、
2以上の対象のそれぞれについて2以上の特性を計測する計測装置から、3以上の次元を有する計測情報を受信する受信手段、
前記計測情報を記憶する記憶手段、
前記計測情報の前記3以上の次元で示される確率モデルについて、前記3以上の次元のうち1の次元を変数とし、変数とした次元以外の次元を変化しないと仮定し、前記計測情報を指数型分布族に当てはめて、前記確率モデルのパラメータを決定するモデル化手段、および、
前記変数の値が、前記パラメータの決定に用いた前記計測情報に含まれない値である前記計測情報に基づいて、前記確率モデルを更新するモデル更新手段、
として機能させることを特徴とするプログラム。

【図1】
image rotate

【図2】
image rotate

【図3】
image rotate

【図4】
image rotate

【図5】
image rotate

【図6】
image rotate

【図7】
image rotate

【図8】
image rotate


【公開番号】特開2013−37471(P2013−37471A)
【公開日】平成25年2月21日(2013.2.21)
【国際特許分類】
【出願番号】特願2011−171836(P2011−171836)
【出願日】平成23年8月5日(2011.8.5)
【出願人】(000004237)日本電気株式会社 (19,353)
【出願人】(504143441)国立大学法人 奈良先端科学技術大学院大学 (226)
【Fターム(参考)】