JP2017151497A

JP2017151497A - 時系列モデルパラメータの推定方法

Info

Publication number: JP2017151497A
Application number: JP2016030599A
Authority: JP
Inventors: 岡本　知樹; Tomoki Okamoto; 知樹岡本; 恒竹田; Hisashi Takeda
Original assignee: Tokyo Electric Power Co Holdings Inc
Current assignee: Tokyo Electric Power Co Holdings Inc
Priority date: 2016-02-22
Filing date: 2016-02-22
Publication date: 2017-08-31

Abstract

【課題】単位時間に対して取得されていないデータを有する時系列データであっても、時系列モデルパラメータの推定ができる時系列モデルパラメータの推定方法を提供する。【解決手段】時間順で取得された複数系列の時系列データにおいて、取得されたデータ毎の当該データと一つ前に取得されたデータとの差分からなる差分時系列データ、及び、これらの時系列データに共通する最小単位時間に基づいて、時系列モデルパラメータを推定する時系列モデルパラメータの推定方法であって、系列数に対応する多変量正規分布の平均及び分散共分散行列、並びに、差分毎に、当該差分を求めたときのデータと一つ前に取得されたデータとの最小単位時間に基づく時間間隔、及び、複数系列のうち、異なる２つの系列毎の差分の重複する最小単位時間に基づく時間間隔により、差分時系列データの尤度関数を決定し、この尤度関数の値を最大化する平均及び分散共分散行列を求める。【選択図】図４

Description

本発明は、時系列モデルパラメータの推定方法に関する。

多次元時系列データ（系列数＝ｎ、時間数＝Ｔ）から時系列モデルパラメータを推定する手法については、欠損が無いケースでは、広く知られている時系列分析手法（Auto Regressive Integrated Moving Average モデル）がある（例えば、特許文献１参照）。この手法は、時系列データの差分（系列数＝ｎ、時間数＝Ｔ−１）が多変量正規分布に従うことにより、標本である差分データから多変量正規分布の母数（平均、分散、共分散）を推定するものである。具体的には、
１．ｎ次元データが多変量正規分布に従うことから尤度を計算
２．時間毎のｎ次元データが独立で同一な分布に従うことから、時間毎に計算された（Ｔ−１）個の尤度の積を計算
３．この積を最大化するような多変量正規分布の母数（平均、分散、共分散）を計算し、推定値とする
の手順を経る。平均、分散、共分散の推定値は、それぞれ標本平均、標本分散、標本共分散に一致し、数式で表せることがわかっている。標本平均がドリフト項、標本分散及び標本共分散が拡散項を表すパラメータの推定値となる。

特開平１１−２０３３６１号公報

多次元時系列データに欠損がないときは上述した方法によりパラメータの推定値を求めることができるが、実際のデータには、系列毎に様々な欠損が発生する場合や、サンプリング間隔が系列毎やサンプリング毎に異なる場合がある。このような場合、各系列で存在するデータの時間が不揃いとなり、ｎ個のデータの組みの間だけの相関を考慮すればよいという仮定が成立しないため、欠損を伴う（Ｔ−１）個のｎ次元データから多変量正規分布の母数の推定に上述した１〜３の従来技術の手順をそのまま適用することができないという課題があった。なお、以降の説明において、サンプリング間隔が異なる場合であって、最小単位時間で計測した場合に、データが存在しない場合も「欠損」と呼ぶこととする。

本発明はこのような課題に鑑みてなされたものであり、欠損を伴う多次元時系列データであっても、時系列モデルパラメータの推定ができる時系列モデルパラメータ推定方法を提供することを目的とする。

前記課題を解決するために、第１の本発明に係る時系列モデルパラメータの推定方法は、時間順で取得された複数系列の時系列データにおいて、取得されたデータ毎に、当該データと一つ前に取得されたデータとの差分からなる差分時系列データ、及び、前記複数系列の時系列データに共通する最小単位時間に基づいて、前記時系列データの時系列モデルパラメータを推定する時系列モデルパラメータの推定方法であって、前記複数系列の系列数に対応する多変量正規分布の平均及び分散共分散行列、並びに、前記時系列データの差分毎に、当該差分を求めたときの前記データと前記一つ前に取得されたデータとの前記最小単位時間に基づく時間間隔、及び、前記複数系列のうち、異なる２つの系列毎の前記差分の重複する前記最小単位時間に基づく時間間隔により、前記差分時系列データの尤度関数を決定し、前記尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求めることを特徴とする。

また、第２の本発明に係る時系列モデルパラメータの推定方法は、時間順で取得された複数系列の時系列データにおいて、取得されたデータ毎に、当該データと一つ前に取得されたデータとの差分からなる差分時系列データ、及び、前記複数系列の時系列データに共通する最小単位時間に基づいて、前記時系列データの時系列モデルパラメータを推定する時系列モデルパラメータの推定方法であって、前記複数系列の系列数に対応する多変量正規分布の平均及び分散共分散行列、並びに、前記時系列データの差分毎に、当該差分を求めたときの前記データと前記一つ前に取得されたデータとの前記最小単位時間に基づく時間間隔、及び、前記複数系列のうち、異なる２つの系列毎の前記差分の重複する前記最小単位時間に基づく時間間隔により、前記差分時系列データの尤度関数を決定するステップと、前記複数系列の時系列データにおいて、前記最小単位時間で前記データが取得されたと仮定したときと、取得された前記時系列データとを比較して、前記データが取得されていない欠損部分にデータを補填して、前記尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求めるステップと、を有することを特徴とする。

また、第３の本発明に係る時系列モデルパラメータの推定方法は、時間順で取得された複数系列の時系列データにおいて、取得されたデータ毎に、当該データと一つ前に取得されたデータとの差分からなる差分時系列データ、及び、前記複数系列の時系列データに共通する最小単位時間に基づいて、前記時系列データの時系列モデルパラメータを推定する時系列モデルパラメータの推定方法であって、前記複数系列のうち、一部の系列の差分時系列データを取り出し、前記一部の系列の系列数に対応する多変量正規分布の平均及び分散共分散行列、並びに、前記差分時系列データの前記差分毎に、当該差分を求めたときの前記データと前記一つ前に取得されたデータとの前記最小単位時間に基づく時間間隔、及び、前記複数系列のうち、異なる２つの系列毎の前記差分の重複する前記最小単位時間に基づく時間間隔により、前記一部の系列の尤度関数を決定し、当該尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求める第１のステップと、前記複数の系列のうち、残りの未知の平均、分散共分散行列の要素が１つ以上含まれるように、一部の系列の差分時系列データを前記第１のステップと同様に取り出し、前記第１のステップで取得された前記平均及び前記分散共分散行列を既知として尤度関数を決定し、当該尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求める第２のステップと、を有し、全ての系列の差分時系列データに対する平均及び分散共分散行列が求まるまで、前記第２のステップを繰り返すことを特徴とする。

また、第４の本発明に係る時系列モデルパラメータの推定方法は、時間順で取得された複数系列の時系列データにおいて、取得されたデータ毎に、当該データと一つ前に取得されたデータとの差分からなる差分時系列データ、及び、前記複数系列の時系列データに共通する最小単位時間に基づいて、前記時系列データの時系列モデルパラメータを推定する時系列モデルパラメータの推定方法であって、前記複数系列のうち、一部の系列の差分時系列データを取り出し、前記一部の系列の系列数に対応する多変量正規分布の平均及び分散共分散行列、並びに、前記差分時系列データの前記差分毎に、当該差分を求めたときの前記データと前記一つ前に取得されたデータとの前記最小単位時間に基づく時間間隔、及び、前記複数系列のうち、異なる２つの系列毎の前記差分の重複する前記最小単位時間に基づく時間間隔により、前記一部の系列の尤度関数を決定し、前記一部の系列の時系列データにおいて、前記最小単位時間で前記データが取得されたと仮定したときと、取得された前記時系列データとを比較して、前記データが取得されていない欠損部分にデータを補填して、前記尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求める第１のステップと、前記複数の系列のうち、残りの未知の平均、分散共分散行列の要素が１つ以上含まれるように、一部の系列の差分時系列データを前記第１のステップと同様に取り出し、前記第１のステップで取得された前記平均及び前記分散共分散行列を既知として尤度関数を決定し、前記残りの系列の一部の系列の時系列データにおいて、前記最小単位時間で前記データが取得されたと仮定したときと、取得された前記時系列データとを比較して、前記データが取得されていない欠損部分にデータを補填して、前記尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求める第２のステップと、を有し、全ての系列の差分時系列データに対する平均及び分散共分散行列が求まるまで、前記第２のステップを繰り返すことを特徴とする。

本発明に係る時系列モデルパラメータの推定方法を以上のように構成すると、単位時間に対して取得されていないデータを有する多次元時系列データであっても、時系列モデルパラメータの推定ができる。

欠損があるときの時系列データに対する差分時系列データを説明するための説明図である。３次元の４時刻分の差分時系列データを説明するための説明図である。時系列モデルパラメータ推定システムの構成を示す説明図である。時系列モデルパラメータの推定方法のフローチャートである。

以下、本発明の好ましい実施形態について図面を参照して説明する。

１．時系列データ｛ｙ_i,t｝の定義
まず、系列ｎ、時間数Ｔとするときの時系列データを｛ｙ_i,t｝（ｉ＝１〜ｎ、ｔ＝１〜Ｔ）とすると、この時系列データ｛ｙ_i,t｝の隣接するデータ同士の差分データ（差分時系列データ）である｛ｘ_i,t｝（ｉ＝１〜ｎ、ｔ＝１〜（Ｔ−１））を、ｘ_i,t＝ｙ_i,(t+1)−ｙ_i,tと定義する。ここで、ｔは最小単位時間に基づいて振り出されている。時系列データが、対数を取るなど所定の変換をした結果として「拡散項が正規分布に従うドリフト（定数）項付きランダムウォークモデル」に従う場合には、以下に示す式（１）が成立する。ここで、ｘは差分時系列データであり、Ｎは正規分布であり、μは平均であり、Σは分散共分散行列である。

このような差分時系列データから、全てのｉ，ｊ，ｔ，ｕに対して、平均Ｅ、分散Ｖ及び共分散Ｃｏｖについては、以下に示す式（２）〜（５）が成立する。なお、ｊはｉとは異なる系列（ｉ≠ｊ）を示している。

２．欠損データがある場合について
次に、ｉ番目（ｉ＝１〜ｎ）のデータの一部に欠損があるときの差分時系列データについて説明する。図１は、ｙ_i,t(s)とｙ_i,t(s+1)の間のデータが欠損している場合を示している。ここで、ｓは非欠損データだけを並べた際に、何番目であるかを示し、ｔ（ｓ）はｓの関数でｓ番目の時刻を表している。例えば、この図１においては、ｙ_i,t(s)の次にｙ_i,t(s+1)が取得できたことを示している。このとき、欠損データ間の差分データは取得することができず、ｙ_i,t(s)とｙ_i,t(s+1)の差分データ（ｘ_i,s ^*とする）だけを取得することができる。すなわち、欠損データがあるときの差分時系列データ（ｘ_i,s ^*）は、以下に示す式（６）で表される。

ここで、式（７）に示すように、差分データの和の期待値Ｅは、各々の差分データの期待値の和で表され、また、差分データの和の分散Ｖは、各々の差分データの和とそれらの差分データの共分散で表される。上述した式（２）〜（５）の性質を用いると、以下に示す式（７）、（８）が成立する。

この式（７）、（８）の関係より、上述した式（６）で表される差分時系列データの平均Ｅ及び分散Ｖはそれぞれ、式（９）、（１０）で表される。

また、欠損を伴う２つの時系列データにおいて、差分時系列データの共分散は、以下に示す式（１１）で表される。ここで、２つの時系列データをｉ，ｊとし、非欠損データの順番をｓ，ｕで表すものとする。なお、式（１１）において、時刻ｔのときが重複しており、（ｔ−１），（ｔ＋１）のときは重複していないことを示している。また、時刻の異なるデータの共分散の値は０となる。

この式（１１）より、欠損を伴う２つの時系列データそれぞれの差分時系列データ（ｘ_i,s ^*とｘ_j,u ^*）間の共分散は、以下に示す式（１２）で表される。

したがって、式（９）、（１０）、（１２）より、欠損を伴う２つの時系列データの差分時系列データは、式（１３）に示す正規分布に従うことになる。

３．欠損がある差分時系列データからのモデルパラメータの最尤推定
上述したように、本実施形態においては、差分時系列データが多変量正規分布に従うことから、当該正規分布の確率密度関数は、当該データを既知とした場合のモデルパラメータの関数（尤度関数）とみることができる。したがって、これを最大化するモデルパラメータ（μ_i，σ_i ²，σ_ij）の値を求めることで、最尤推定量が得られる。なお、差分データの系列間は、系列数を次元とする多変量正規分布に従うことから、式（１）に示す分散共分散行列Σ（＝｛σ_i ²，σ_ij｝）が半正定値行列でなければならない。つまり、Σの全固有値が非負となる条件下で上記尤度関数を最大化することになる。

上述した欠損がある差分時系列データについて、図２に示すように、３次元の４時刻分を例にして説明する。ここでは、ｘ、ｙ、ｚの３つの系列があり、上段に示すように、それぞれ４時刻分のデータが発生するが、下段に示すように、一部のデータに欠損が発生している。具体的には、系列ｘでは、ｘ₁が取得できずにｘ₁及びｘ₂の合計であるｘ₁ ^*が取得され、ｘ₃が取得できずにｘ₃及びｘ₄の合計であるｘ₂ ^*が取得されている。また、系列ｙでは、ｙ₁はｙ₁ ^*として取得されているが、ｙ₂が取得できずにｙ₂及びｙ₃の合計であるｙ₂ ^*が取得され、ｙ₄がｙ₃ ^*として取得されている。また、系列ｚでは、ｚ₁はｚ₁ ^*として取得されているが、ｚ₂〜ｚ₄が取得できずにｚ₂〜ｚ₄の合計であるｚ₂ ^*が取得されている。

上述した欠損を含む差分時系列データをＸとし、平均をμとし、分散共分散行列をΣとすると、差分時系列データＸは、平均μ及び分散共分散行列Σの正規分布に従い、以下に示す式（１４）で表される。

ここで、差分時系列データＸは、取得された全てのデータのベクトルである。そして、平均μは、式（９）を用いて説明したように、各系列の平均（例えば、系列ｘにおけるμ_x）に各時刻の長さ（個数）をかけた値である。具体的には、図２に示すように、データｘ₁ ^*は、ｘ₁₁とｘ₁₂に対応しているため、その時刻の長さ（個数）は「２」となり、その平均は２μ_xとなる。他のデータについても同様である。

次に、分散共分散行列Σは、式（１３）を用いて説明したように、分散は、平均と同様に、各系列の分散（例えば、系列ｘにおけるσ_x ²）に各時刻の長さ（個数）をかけた値である。具体的には、図２に示すように、ｘ₁ ^*の分散は２σ_x ²となる。また、共分散は、重複する時刻がないときは「０」となり、重複する時刻があるときは、共分散の値に重複する時刻の長さ（個数）をかけたものである。具体的には、図２に示すように、ｘ₁ ^*とｙ₃ ^*とは重複する時刻がないため「０」となり、ｘ₁ ^*とｙ₂ ^*とは２番目の時刻のデータが重複するため、系列ｘと系列ｙとの共分散σ_xyにその重複する時刻の長さ（個数）である「１」をけたσ_xyとなる。

以上より求められた平均μ及び分散共分散行列Σより、尤度関数Ｌを求めると、以下に示す式（１５）で表され、また、対数尤度ＬＬは式（１６）で表され、この対数尤度ＬＬを最大化する平均μ及び分散共分散行列Σを求めることで時系列モデルパラメータの値を決定することができる。なお、式（１５）及び式（１６）の第１項の「７」は、図２に示すように、取得されたデータが、ｘ₁ ^*〜ｚ₂ ^*の数が７個であることが理由である。

なお、系列数３を次元とする多変量正規分布の分散共分散行列は、以下に示す式（１７）で表され、この分散共分散行列の固有値が非負となる条件下での尤度関数Ｌの最大化を行うことが必要である。

４．欠損が全くない場合
ここまで説明したように、時系列データに欠損がある場合でも、得られた差分データから時系列モデルパラメータを推定することができる。そこで、この方法が、欠損が全くない場合でも適用可能なことについて説明する。

欠損を考慮した場合の差分時系列データｘ_i,t(s)において、欠損がないときはｔ（ｓ）＝ｓとなる。したがって、欠損がないときの差分時系列データをｘ_i,s（ｉ＝１〜Ｎ，ｓ＝１〜Ｓ（＝Ｔ−１））とおくと、差分時系列データＸ、平均μは以下に示す式（１８）、（１９）で表される。なお、ｍは推定された平均を示す。

また、推定されるべき分散共分散行列を以下の式（２０）に示すＶとすると、全体の分散共分散行列Σは式（２１）で表される。なお、Ｏは０行列である。

以上より、対数尤度ＬＬは以下に示す式（２２）で表される。

この式（２２）は、多変量正規分布に従うＮ次元ベクトルがＳ個あって、それらが独立である場合の対数尤度ＬＬにほかならない。この場合の平均ｍ及び推定された分散共分散行列Ｖにかかる最尤推定量については、解析解が知られており、以下の式（２３）、（２４）で示される通り、標本平均、標本分散共分散行列に一致する。また、この標本分散共分散行列が半正定値行列となることも知られている。

５．一次元の場合
ここまでの説明では複数の系列の時系列データについて説明してきたが、ここでは１次元の場合について説明する。上述した式（６）より、１次元の差分データｘ_is ^*は、以下に示す式（２５）のように正規分布に従う。なお、１次元であるため、ｉ＝１となる。また、差分時系列データＸ、平均μ及び分散共分散行列Σも示す。

そして、この式（２５）より、尤度関数Ｌ及び対数尤度ＬＬは式（２６）、（２７）で表される。

ここで、式（２７）に示される対数尤度ＬＬを最大化する平均μ_i及びσ_i ²は、解析的に求めることができ、最尤推定量は、式（２８）、（２９）となる。なお、もし欠損が全くない場合には、以下の式（２８），（２９）は標本平均、標本分散に一致する。

６．導出した尤度関数を最大化するモデルパラメータの計算手法（一括使用による解法）
以上のようにして導出された尤度関数に基づいて、この尤度関数を最大化するモデルパラメータ（平均、分散、共分散）を求める方法として、まず、一括使用による解法を説明する。この一括使用による解法では、尤度関数から数理最適化手法（例えば、準ニュートン法等）を用いてモデルパラメータを求める方法である。但し、制約条件（半正定値条件）及び目的関数（対数尤度最大化）がともに非線形となるため、モデルパラメータを直接変えながらモデルパラメータの最尤推定量を探索するのではなく、欠損部分に入れるデータを変化させながら、都度見かけ上欠損を無くし、その上で、式（２３），（２４）に示した解析解に従って計算し、しかる後、欠損を考慮した尤度関数の値を計算することを繰り返しながら探索する。なお、この一括使用による解法では、式（２４）が半正定値行列であることから、半正定値条件を充足しながら、対数尤度を最大化するモデルパラメータを探索することができる。

７．導出した尤度関数を最大化するモデルパラメータの計算手法（逐次使用による解法）
ここでは、ｎ個の系列の中から一部の系列を抜き出して、モデルパラメータを計算し、ついで、残りの未知のパラメータが１つ以上含まれるように一部の系列を前回同様に抜き出し、既に計算されているモデルパラメータを固定値として未知のパラメータを計算し、上記の処理を繰り返して全てのモデルパラメータを計算するという方法である。なお、この方法の場合、全てのモデルパラメータが得られた後、そのモデルパラメータにより半正定値条件を満足するか否かを確認し、半正定値条件を満たしていない場合は、モデルパラメータから定まる分散共分散行列が半正定値条件を満たすように変形し補正解を求めることが必要である。

例えば、まず、ｎ個の系列を個々に１次元系列ととらえ、上述した式（２８），（２９）により、平均と分散の最尤推定量を計算し、ついで、ｎ個の系列から２系列を取り出す全ての組み合わせを考え、それぞれの組み合わせを２次元系列ととらえ、１次元系列として得た平均と分散の推定値を固定した対数尤度から、共分散の最尤推定量を計算し、最後に、半正定値条件を確認し、この条件を満たしていない場合には適宜変形するという方法がある。なお、このような解析的な計算方法の代わりに、上述した数理最適化手法による計算方法を用いることも可能である。

なお、半正定値条件を満たさない場合の変形方法としては、上述の方法により得られた分散共分散行列において、条件を満たさない推定分散共分散行列に対し、条件を満たす分散共分散行列を別に考え、両者の差（例えば、何らかの行列ノルム）が最小となるものを見いだす方法がある。あるいは、上述の方法により得られた分散共分散行列において、条件を満たさない推定分散共分散行列を固有値分解し、対角行列の対角成分の値を非負となるよう調整し、分解式に当てはめて分散共分散行列を修正するという方法がある。

８．時系列モデルパラメータ推定システム
次に、これまで説明した時系列モデルパラメータの推定を行う時系列モデルパラメータ推定システム１００について図３及び図４を用いて説明する。この時系列モデルパラメータ推定システム１００は、図３に示すように、ＣＰＵやＲＡＭ、ＲＯＭ等を有し、プログラムを実行することにより時系列モデルパラメータの推定を行う処理部１１０と、所定の情報を入力する入力部１２０と、推定結果等を記憶するためのメモリーやハードディスク等からなる記憶部１３０と、推定結果等を出力する出力部１４０と、を有して構成されている。

図４を用いて、処理部１１０で実行される処理について説明する。処理部１１０は、時系列モデルパラメータ推定処理が実行されると、まず、差分時系列データ（欠損を含む場合がある）を読み込む（ステップＳ２００）。なお、予め時系列データ（例えば、取得時刻と取得されたデータの組み合わせ）を記憶部１３０に記憶しておいてこれを読み込むように構成してもよいし、時系列データを読み込んで、この時系列データから差分時系列データを生成するように構成してもよい。そして、一括使用か逐次使用かを判断する（ステップＳ２０２）。どちらの方法を選択するかを予め記憶部１３０に設定しておいてもよいし、入力部１２０からの入力により選択するように構成してもよい。

ステップＳ２０２において、一括使用が選択された場合には、上述したように、欠損が無い状態での時系列データから求められる差分時系列データを生成し（ステップＳ２０４）、この欠損のない時系列データから求められる差分時系列データから平均、分散、共分散を計算し（ステップＳ２０６）、その結果から尤度（対数尤度）を計算する（ステップＳ２０８）。このステップＳ２０４〜Ｓ２０８の処理を、欠損データを代えて繰り返し実行し、尤度が最大になったか否かを判断し（ステップＳ２１０）、最大になったと判断したときは、そのときの平均、分散、共分散を時系列モデルパラメータの推定値として出力する（ステップＳ２２４）。出力方法として、ディスプレイやプリンタ等からなる出力部１４０に出力（表示）してもよいし、記憶部１３０に記憶させてもよい。

一方、ステップＳ２０２において、逐次使用が選択された場合には、全系列データから一部の系列データを抜き出し、その系列に対して、上述したように、時系列モデルパラメータ（平均、分散共分散行列）を計算する（ステップＳ２１２）。ついで、残りの未知のパラメータが１つ以上含まれるように一部の系列を前回同様に抜き出し（ステップＳ２１４）、既に計算されている系列の時系列モデルパラメータを既知の値として未知の時系列モデルパラメータを計算し（ステップＳ２１６）、全ての系列の計算が終了するまでステップＳ２１４〜Ｓ２１６を繰り返す（ステップＳ２１８）。そして、全ての系列が計算されたと判断すると、上述した処理により算出された最終的な時系列モデルパラメータが半正定値条件を満たすか否かを判断し（ステップＳ２２０）、満たさないときは補正解を算出して補正し（ステップＳ２２２）、結果を平均、分散、共分散を時系列モデルパラメータの推定値として出力する（ステップＳ２２４）。

それでは、具体的な実施例として、下記の表１に示すように、ｘ，ｙ，ｚの３つの系統における７時刻分の時系列データに対し、上述した一括使用による解法及び逐次使用による解法に基づいて時系列モデルパラメータを求める場合について説明する。なお、以下の表１において、「？」はそのデータが欠損していることを示している。また、各系列における「差」の値は、当該時刻におけるデータと、一つ前に取得されたデータとの差分を示している。

また、この表１に示す時系列データ（差分時系列データ）に対して、全ての差分を、ベクトルにした差分時系列データＸ、並びに、この差分時系列データＸから求められる平均μ及び分散共分散行列Σは以下のように表される。なお、Ｄは差分時系列データＸのデータ数を示している（ｘ_tに対して４個の差分データが取得され、ｙ_tに対して６個の差分データが取得され、ｚ_tに対して５個の差分データが取得されているため、合計でＤ＝１５となる）。

また、以降の説明において、対数尤度ＬＬを、データ数Ｄと差分時系列データＸをパラメータに有する、平均μ及び分散共分散行列Σの関数として以下のように定義し、この対数尤度ＬＬが最大になる時系列モデルパラメータを求めるものとする。

但し、この対数尤度ＬＬにおいて、系列x_t〜z_tの3系列（3次元）の多変量正規分布の分散共分散行列をＶとしたとき、この分散共分散行列Ｖは、次式のように表され、その固有値が非負となる（すなわち、半正定値行列となる）条件下での対数尤度ＬＬの最大化を行う。

（第１の実施例）
まず、第１の実施例として、上述した差分時系列データＸに対し、一括解法（一括使用による解法）で、時系列モデルパラメータを推定する場合について説明する。

対数尤度ＬＬの最大化にあたって、上述した分散共分散行列Ｖの半正定値条件を充足させるために、欠損した時系列データに対して何らかの値を与えて欠損無しの状態を生成し、この場合に理論的に得られる平均μ及び分散共分散行列Ｖの推定値を用いて対数尤度ＬＬを求める。ここで、欠損無しの場合の理論値は必ず半正定値条件を満たすことから、欠損を補填する値を変えて計算を行うことで、対数尤度ＬＬを最大とする平均μ及び分散共分散行列Ｖを見いだすことにすれば、半正定値条件の下で対数尤度ＬＬの最大化を行うこととなる。

まず、表１の値より、ｘ、ｙ、ｚ各系列の平均μ_x、μ_y、μ_z、並びに、分散共分散行列Ｖの値を、式（２３）、（２４）を用いて求める。

上述したように、欠損データであるｘ_３，ｘ₄，ｚ₄を変化させながら、対数尤度ＬＬが最大となる平均μ及び分散共分散行列Ｖを求める。

（第２の実施例）
次に、第２の実施例として、逐次解法（逐次使用による解法）で、時系列モデルパラメータを推定する場合について説明する。ここでは、ｘ_t、ｙ_t、ｚ_tからなる系列数３の時系列データに対し、それぞれ１次元データとしてモデルパラメータを算出し、次に、３つの系列から２つの系列を選択し、１次元の結果を用いることで、２次元データとしてモデルパラメータを算出し、それらの結果から、最終的なモデルパラメータの推定値を求める場合について説明する。

［１］ｘ_tだけの一次元データとしたときのモデルパラメータ
以下の表３は、表１から系列ｘ_tだけの時系列データ及びその差分を抽出したものである。

この表２から、式（２８）（２９）により系列ｘ_tの平均μ_x及び分散σ² _xを求める。これをｍ_x及びｓ_x ²（但しｓ_x≧０）とする。

［２］ｙ_tだけの一次元データとしたときのモデルパラメータ
表１から系列ｙ_tだけの時系列データ及びその差分を抽出し、［１］と同様の手法により平均μ_y及び分散σ_y ²を求める。これをｍ_y及びｓ_y ²（但しｓ_y≧０）とする。

［３］ｚ_tだけの一次元データとしたときのモデルパラメータ
表１から系列ｚ_tだけの時系列データ及びその差分を抽出し、［１］と同様の手法により平均μ_z及び分散σ_z ²を求める。これをｍ_z及びｓ_z ²（但しｓ_z≧０）とする。

［４］ｘ_t，ｙ_tの二次元データとしたときのモデルパラメータ
以下の表３は、表１から系列ｘ_t，ｙ_tの時系列データ及びその差分を抽出したものである。

系列ｘ_t，ｙ_tにおいて、差分データをベクトル化したものをＸ_xyとし、そのときの平均をμ_xy、分散共分散行列をΣ_xyとすると、以下のように表される。ここで、μ_x、μ_y、σ_x ²、σ_y ²は、一次元のデータとして求められた値を使用する。また、系列ｘ_t，ｙ_tにおいては、５個の差分データが取得されているため、Ｄ＝１０となる。

そして、ＬＬ＝ＬＬ（μ_xy，Σ_xy；１０，Ｘ_xy）が最大となるσ_xyを求める。なお、系列
ｘ_t及び系列ｙ_tの２系列（２次元）多変量正規分布の分散共分散行列Ｖ_xyは、以下のように表される。

ここで、分散共分散行列Ｖ_xyが半正定値行列、すなわち、σ_xyが次式の範囲であることが必要である。

以上より、（μ_x、μ_y、σ_x ²、σ_y ²）を（ｍ_x、ｍ_y、ｓ_x ²、ｓ_y ²）に固定し、σ_xyを前式の範囲で変化させながら、対数尤度ＬＬが最大となるσ_xyを求める。これをs_xyとする。

［５］ｘ_t，ｚ_tの二次元データとしたときのモデルパラメータ
表１から系列ｘ_t，ｚ_tの時系列データ及びその差分を抽出し、［４］と同様の手法により共分散σ_xzを求める。これをｓ_xzとする。

［６］ｙ_t，ｚ_tの二次元データとしたときのモデルパラメータ
表１から系列ｙ_t，ｚ_tの時系列データ及びその差分を抽出し、［４］と同様の手法により共分散σ_yzを求める。これをｓ_yzとする。

［７］推定値の補正
上述した（１）〜（６）から、本実施例における時系列モデルパラメータの推定値は以下のようになる。

ここで、分散共分散行列Ｖを固有値分解すると以下のようになる。

この分散共分散行列Ｖにおいて、右辺の２番目の行列には、対角成分に固有値が並んでいるが、負となる成分が１つ以上ある場合には、この分散共分散行列Ｖは半正定値行列となっていない。このため、半正定値行列となるように、負の成分を補正する。以下の例は、正の微小数ε（＝０．０００１）に補正した場合を示している（補正後の分散共分散行列をＶ^*とする）。

このが最終的な推定値となる。

以上のように、本実施形態に係る時系列モデルパラメータの推定方法によれば、多次元時系列データの一部に欠損があったとしても、差分データのパターンに応じた同時分布を生成できるようになったため、尤度さらにはこれを最大化する母数の値を計算することが可能になった。これにより、欠損を伴う多次元時系列データであっても、欠損していない残りの全てのデータを使用した時系列モデルパラメータの推定が可能となる。なお、欠損を伴う多次元時系列データだけでなく、最小単位時間は共通するが、この最小単位時間の整数倍（系列毎に倍数が異なる場合や、１つの系列内でも倍数が一定でない場合、及びその組み合わせも含む）でデータが取得される時系列データに対しても、単位時間に対して取得されていないデータを上記欠損と考えれば、本発明を適用することができる。

１００時系列モデルパラメータ推定システム
１１０処理部
１２０入力部
１３０記憶部
１４０表示部

Claims

時間順で取得された複数系列の時系列データにおいて、取得されたデータ毎に、当該データと一つ前に取得されたデータとの差分からなる差分時系列データ、及び、前記複数系列の時系列データに共通する最小単位時間に基づいて、前記時系列データの時系列モデルパラメータを推定する時系列モデルパラメータの推定方法であって、
前記複数系列の系列数に対応する多変量正規分布の平均及び分散共分散行列、並びに、前記時系列データの差分毎に、当該差分を求めたときの前記データと前記一つ前に取得されたデータとの前記最小単位時間に基づく時間間隔、及び、前記複数系列のうち、異なる２つの系列毎の前記差分の重複する前記最小単位時間に基づく時間間隔により、前記差分時系列データの尤度関数を決定し、
前記尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求めることを特徴とする時系列モデルパラメータの推定方法。
時間順で取得された複数系列の時系列データにおいて、取得されたデータ毎に、当該データと一つ前に取得されたデータとの差分からなる差分時系列データ、及び、前記複数系列の時系列データに共通する最小単位時間に基づいて、前記時系列データの時系列モデルパラメータを推定する時系列モデルパラメータの推定方法であって、
前記複数系列の系列数に対応する多変量正規分布の平均及び分散共分散行列、並びに、前記時系列データの差分毎に、当該差分を求めたときの前記データと前記一つ前に取得されたデータとの前記最小単位時間に基づく時間間隔、及び、前記複数系列のうち、異なる２つの系列毎の前記差分の重複する前記最小単位時間に基づく時間間隔により、前記差分時系列データの尤度関数を決定するステップと、
前記複数系列の時系列データにおいて、前記最小単位時間で前記データが取得されたと仮定したときと、取得された前記時系列データとを比較して、前記データが取得されていない部分にデータを補填して、前記尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求めるステップと、を有することを特徴とする時系列モデルパラメータの推定方法。
時間順で取得された複数系列の時系列データにおいて、取得されたデータ毎に、当該データと一つ前に取得されたデータとの差分からなる差分時系列データ、及び、前記複数系列の時系列データに共通する最小単位時間に基づいて、前記時系列データの時系列モデルパラメータを推定する時系列モデルパラメータの推定方法であって、
前記複数系列のうち、一部の系列の差分時系列データを取り出し、前記一部の系列の系列数に対応する多変量正規分布の平均及び分散共分散行列、並びに、前記差分時系列データの前記差分毎に、当該差分を求めたときの前記データと前記一つ前に取得されたデータとの前記最小単位時間に基づく時間間隔、及び、前記複数系列のうち、異なる２つの系列毎の前記差分の重複する前記最小単位時間に基づく時間間隔により、前記一部の系列の尤度関数を決定し、当該尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求める第１のステップと、
前記複数の系列のうち、残りの未知の平均、分散共分散行列の要素が１つ以上含まれるように、一部の系列の差分時系列データを前記第１のステップと同様に取り出し、前記第１のステップで取得された前記平均及び前記分散共分散行列を既知として尤度関数を決定し、当該尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求める第２のステップと、を有し、
全ての系列の差分時系列データに対する平均及び分散共分散行列が求まるまで、前記第２のステップを繰り返すことを特徴とする時系列モデルパラメータの推定方法。
時間順で取得された複数系列の時系列データにおいて、取得されたデータ毎に、当該データと一つ前に取得されたデータとの差分からなる差分時系列データ、及び、前記複数系列の時系列データに共通する最小単位時間に基づいて、前記時系列データの時系列モデルパラメータを推定する時系列モデルパラメータの推定方法であって、
前記複数系列のうち、一部の系列の差分時系列データを取り出し、前記一部の系列の系列数に対応する多変量正規分布の平均及び分散共分散行列、並びに、前記差分時系列データの前記差分毎に、当該差分を求めたときの前記データと前記一つ前に取得されたデータとの前記最小単位時間に基づく時間間隔、及び、前記複数系列のうち、異なる２つの系列毎の前記差分の重複する前記最小単位時間に基づく時間間隔により、前記一部の系列の尤度関数を決定し、前記一部の系列の時系列データにおいて、前記最小単位時間で前記データが取得されたと仮定したときと、取得された前記時系列データとを比較して、前記データが取得されていない部分にデータを補填して、前記尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求める第１のステップと、
前記複数の系列のうち、残りの未知の平均、分散共分散行列の要素が１つ以上含まれるように、一部の系列の差分時系列データを前記第１のステップと同様に取り出し、前記第１のステップで取得された前記平均及び前記分散共分散行列を既知として尤度関数を決定し、前記残りの系列の一部の系列の時系列データにおいて、前記最小単位時間で前記データが取得されたと仮定したときと、取得された前記時系列データとを比較して、前記データが取得されていない部分にデータを補填して、前記尤度関数の値を最大化する前記差分時系列データの平均及び分散共分散行列を求める第２のステップと、を有し、
全ての系列の差分時系列データに対する平均及び分散共分散行列が求まるまで、前記第２のステップを繰り返すことを特徴とする時系列モデルパラメータの推定方法。
前記全ての系列の差分時系列データに対する前記差分時系列データの平均及び分散共分散行列が半正定値条件を満たしていないときは、半正定値条件を満たすように、前記平均及び前記分散共分散行列を補正することを特徴する請求項３または４に記載の時系列モデルパラメータの推定方法。
前記差分時系列データの平均は、前記系列毎の平均と前記差分を求めたときの前記データと前記一つ前に取得されたデータとの前記最小単位時間に基づく時間間隔との積であり、
前記差分時系列データの分散は、前記系列数に対応する多変量正規分布の分散と当該差分を求めたときの前記データと前記一つ前に取得されたデータとの前記最小単位時間に基づく時間間隔との積であり、
前記差分時系列データの共分散のうち、前記異なる２つの系列の同一の時刻の値は、前記系列数に対応する多変量正規分布の共分散と前記異なる２つの系列毎の前記差分の重複する前記最小単位時間に基づく時間間隔との積であり、
前記差分時系列データの共分散のうち、前記異なる２つの系列の異なる時刻の値は０であることを特徴とする請求項１〜５のいずれか一項に記載の時系列モデルパラメータの推定方法。