最適モデル推定装置、方法、及びプログラム
【課題】誤差の分散値を推定せずに、高精度にモデリングを行うことができるようにする。
【解決手段】状態空間モデリング部２２によって、観測値の時系列データに基づいて、カルマンフィルタのモデルパラメータを推定する。正規化カルマンフィルタデザイン部２２において、状態ベクトルＸの推定誤差分散の共分散行列を予測誤差分散で正規化した正規化推定誤差分散共分散行列と、正規化誤差λとを用いた正規化カルマンフィルタが定められる。状態推定部２３によって、定められた正規化カルマンフィルタに従って、予測値を時系列に算出すると共に、対応する観測値との予測誤差を各々算出する。最適モデル推定部２４によって、予め用意された正規化誤差λの各値を用いたときに状態推定部２３によって各々算出された予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた正規化カルマンフィルタを、最適モデルの推定結果とする。
【解決手段】状態空間モデリング部２２によって、観測値の時系列データに基づいて、カルマンフィルタのモデルパラメータを推定する。正規化カルマンフィルタデザイン部２２において、状態ベクトルＸの推定誤差分散の共分散行列を予測誤差分散で正規化した正規化推定誤差分散共分散行列と、正規化誤差λとを用いた正規化カルマンフィルタが定められる。状態推定部２３によって、定められた正規化カルマンフィルタに従って、予測値を時系列に算出すると共に、対応する観測値との予測誤差を各々算出する。最適モデル推定部２４によって、予め用意された正規化誤差λの各値を用いたときに状態推定部２３によって各々算出された予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた正規化カルマンフィルタを、最適モデルの推定結果とする。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、最適モデル推定装置、方法、及びプログラムに係り、特に、カルマンフィルタを用いた最適モデルを推定する最適モデル推定装置、方法、及びプログラムに関する。
【背景技術】
【０００２】
時系列データを状態空間でモデリングする際、モデル誤差の分散値と観測誤差の分散値とが未知となっている場合がほとんどである。状態空間では、カルマンフィルタによる状態の逐次推定及び推定された状態に基づく予測が可能であるが、誤差の分散値に関しては推測による設定が多い。
【０００３】
例えば、状態の推定誤差の状態空間モデルを考え、予測誤差の自己相関関数を基に誤差の分散値を推定する手法が知られている（非特許文献１）。
【先行技術文献】
【非特許文献】
【０００４】
【非特許文献１】“Ａ ｎｅｗ ａｕｔｏｃｏｖａｒｉａｎｃｅ ｌｅａｓｔ−ｓｑｕａｒｅｓ ｍｅｔｈｏｄ ｆｏｒ ｅｓｔｉｍａｔｉｎｇ ｎｏｉｓｅ ｃｏｖａｒｉａｎｃｅ”，Ｏｄｅｌｓｏｎ ｅｔ ａｌ.Ｊｏｕｒｎａｌ ｏｆ Ａｕｔｏｍａｔｉｃａ， ｖｏｌ.４２，ｐｐ.３０３−３０８，２００６.
【発明の概要】
【発明が解決しようとする課題】
【０００５】
しかしながら、上記の特許文献１に記載の技術では、誤差の分散値を推定しているため、推定した誤差の分散値が間違っている場合、誤差の分散値の間違った設定によりモデリングの精度が下がり、予測精度も低下してしまう、という問題が生じる。
【０００６】
本発明は、上記の事情を鑑みてなされたもので、誤差の分散値を推定せずに、高精度にモデリングを行うことができる最適モデル推定装置、方法、及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【０００７】
上記の目的を達成するために本発明に係る最適モデル推定装置は、観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定する状態空間モデリング手段と、前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段と、予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とする最適モデル推定手段と、を含んで構成されている。
【０００８】
本発明に係る最適モデル推定方法は、状態空間モデリング手段、状態推定手段、及び最適モデル推定手段を含む最適モデル推定装置における最適モデル推定方法であって、前記最適モデル推定装置は、前記状態空間モデリング手段によって、観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定するステップと、前記状態推定手段によって、前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出するステップと、前記最適モデル推定手段によって、予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とするステップと、を含んで実行することを特徴とする。
【０００９】
本発明に係るプログラムは、コンピュータを、観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定する状態空間モデリング手段、前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段、及び予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とする最適モデル推定手段として機能させるためのプログラムである。
【００１０】
本発明によれば、状態空間モデリング手段によって、観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定する。状態空間モデリング手段によって推定された前記モデルパラメータに基づいて、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタが定められる。
【００１１】
そして、状態推定手段によって、定められた正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する。最適モデル推定手段によって、予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とする。
【００１２】
このように、モデル誤差分散及び観測誤差分散の何れか一方を予測誤差分散で正規化した正規化誤差λを用いた正規化カルマンフィルタに従って、観測値の時系列データに基づいて予測値を算出し、予測誤差に基づいて最適な正規化誤差λを特定することにより、誤差の分散値を推定せずに、高精度にモデリングを行うことができる。
【発明の効果】
【００１３】
以上説明したように、本発明の最適モデル推定装置、方法、及びプログラムによれば、モデル誤差分散及び観測誤差分散の何れか一方を予測誤差分散で正規化した正規化誤差λを用いた正規化カルマンフィルタに従って、観測値の時系列データに基づいて予測値を算出し、予測誤差に基づいて最適な正規化誤差λを特定することにより、誤差の分散値を推定せずに、高精度にモデリングを行うことができる、という効果が得られる。
【図面の簡単な説明】
【００１４】
【図１】本発明の実施の形態に係る最適モデル推定装置の構成を示す概略図である。
【図２】観測値の時系列データを説明するための図である。
【図３】（Ａ）入力される観測値を示す図、（Ｂ）パラメータ推定結果を示す図、及び（Ｃ）出力される予測値を示す図である。
【図４】本発明の実施の形態に係る最適モデル推定装置における最適モデル推定処理ルーチンの内容を示すフローチャートである。
【発明を実施するための形態】
【００１５】
以下、図面を参照して本発明の実施の形態を詳細に説明する。
【００１６】
＜システム構成＞
本発明の実施の形態に係る最適モデル推定装置１００は、観測値の時系列データが入力され、状態空間の最適モデルを推定する。この最適モデル推定装置１００は、ＣＰＵと、ＲＡＭと、後述する最適モデル推定処理ルーチンを実行するためのプログラムを記憶したＲＯＭとを備えたコンピュータで構成され、機能的には次に示すように構成されている。図１に示すように、最適モデル推定装置１００は、入力部１０と、演算部２０と、出力部３０とを備えている。
【００１７】
入力部１０は、入力された観測値の時系列データを受け付ける。観測値は、スカラー値であり、すなわち、時系列データは、１次元時系列である。例えば、図２に示すように、エポックｔnと観測値（ｔn）との組（ｎ＝１、・・・、Ｎ）からなる時系列データである。
【００１８】
演算部２０は、状態空間モデリング部２１、正規化カルマンフィルタデザイン部２２、状態推定部２３、最適モデル推定部２４、推定結果データ生成部２５、及び予測値算出部２６を備えている。
【００１９】
状態空間モデリング部２１は、入力部１０により受け付けた観測値の時系列データに基づいて、以下に説明するように、時系列データを状態空間で表現したガウス型モデルのモデルパラメータを推定する。
【００２０】
まず、一般に、カルマンフィルタを表わすガウス型モデルは以下の（１）式、（２）式で表現される。
【００２１】
【数１】
【００２２】
ただし、Ｘｔは時刻ｔにおける状態ベクトルで、通常、非観測量（モデルの状態及び未知パラメータ）が取られる。ｗは正規白色ノイズのベクトルであり、モデル誤差が含まれる。上記（１）式に示されるように、状態ベクトルＸｔは線形（非線形）作用ｆとノイズベクトルｗにより時間更新される。
【００２３】
また、ｙｔは時刻ｔにおける観測値であり、状態Ｘｔの線形（非線形）作用ｇにより生成され、正規白色観測誤差μがのる。上記（２）式は、観測誤差を用いて状態ベクトルＸtと観測値ｙtとの関係を示している。
【００２４】
上記（１）式、（２）式で示される状態空間表現は、ＡＲ（自己回帰モデル）に代表される様々な時系列モデルを含む、汎用性が高いモデルである。以降、時系列モデルとして、以下の（３）式、（４）式で表わされるＡＲ＋ドリフトモデルを用いた場合を例として説明する。
【００２５】
【数２】
【００２６】
ただし、ωがモデル誤差に相当し、μが観測誤差に相当する。ここで、Ｎ（ｍ、ｖ）は平均ｍ、分散ｖの正規分布を表す。したがって、ｑはモデル誤差分散に相当し、ｒは観測誤差分散に相当する。このモデルの状態空間表現は、以下の（５）式〜（８）式で示される。
【００２７】
【数３】
【００２８】
Ｈは、状態を観測値に変換する入出力行列を表す。
【００２９】
状態空間モデリング部２１は、観測値の時系列データに基づいて、上記（５）式〜（８）式で表わされるガウス型モデルのモデルパラメータ（ｘt-1，・・・，ｘt-M、ａ1，・・・，ａM，ｄt，ｑ，ｒ）を推定する。
【００３０】
正規化カルマンフィルタデザイン部２２は、状態空間モデリング部２１によって推定されたガウス型モデルのモデルパラメータに基づいて、オペレータによる入力部１０からの入力を受け付けて、以下に説明するように、正規化拡張カルマンフィルタをデザインする。
【００３１】
状態空間モデリング部２１により状態空間表現された時系列データは、以下のアルゴリズムで表現される。まず、状態ベクトルを時間更新するアルゴリズムが、以下の（９）式〜（１１）式で表現される。
【００３２】
【数４】
【００３３】
上記（１１）式で表わされる^Ｆは、ｆの微分をオペレータが計算して入力することにより、求められる。
【００３４】
また、得られた観測値を用いる観測更新のアルゴリズムが、以下の（１２）式〜（１５）式で表される。
【００３５】
【数５】
【００３６】
Ｘｔ｜ｓはｓ期におけるｔ期の状態推定値を意味する。＾Ｘは状態Ｘの期待値である。Ｐ’は状態Ｘの正規化推定誤差分散の共分散行列であり、以下の（１６）式に示すように、推定誤差分散の共分散行列Ｐを予測誤差ｑ＋ｒで正規化したものである。
【００３７】
【数６】
【００３８】
この正規化プロセスにより、通常のカルマンフィルタのアルゴリズムからモデル誤差分散ｑと観測誤差分散ｒとが消え（Ｑ’＝ｄｉａｇ（λ、０、…、０））、代わりに以下の（１７）で定義するλが導入される。
【００３９】
【数７】
【００４０】
λは予測誤差分散ｑ＋ｒにおけるモデル誤差分散ｑの割合であり、０≦λ≦１の値をとるパラメータである。このλの値を設定することで、上記（９）式〜（１５）式で表現された正規化拡張カルマンフィルタのデザインが完了する。以下、λを、正規化誤差と称する。
【００４１】
このように、上記（１６）式に示す状態Ｘの正規化推定誤差分散の共分散行列と、上記（１７）式に示す正規化誤差とを用いた正規化拡張カルマンフィルタのアルゴリズムがデザインされる。
【００４２】
状態推定部２３は、予め用意されたλの各値（例えば、０から１まで０．１ずつ加算した各値）及び次数Ｍの各値（例えば、１から１０までの自然数）の組み合わせについて、正規化カルマンフィルタデザイン部２２によってデザインされ、かつ、当該λの値及び次数Ｍの値が設定された、上記（９）式〜（１５）式で表現された正規化拡張カルマンフィルタのアルゴリズムに従って、観測値の時系列データを用いて、時間更新及び観測更新を逐次行う。
【００４３】
最適モデル推定部２４は、λの値及び次数Ｍの値の各組み合わせについて、状態推定部２３によって逐次行われた時間更新及び観測更新の結果を利用して、以下の（１８）式に基づき、時系列データの最適モデルを推定する。
【００４４】
【数８】
【００４５】
ただし、ｖｔは上記（１５）式で示した予測誤差であり、Σｖｔは予測誤差ｖｔの分散である。ＡＩＣは赤池情報量であり、通常、ＡＲを想定した場合は次数Ｍ（何期前のデータまで使用するか）の関数である。しかしながら、本実施の形態で用いる正規化拡張カルマンフィルタのアルゴリズムに対しては、ＡＩＣは、Ｍとλの関数になる。ここで、上記（１８）で用いる予測誤差ｑ＋ｒには、以下の（１９）式に従って推定された値を適用すればよい。
【００４６】
【数９】
【００４７】
ただし、Ｐt’[1,1]は、行列Ｐt’の要素（１，１）である。
【００４８】
上記のように、Ｍとλをパラメトライズし、最小のＡＩＣ（Ｍ、λ）を記録したときのＭとλの値の組み合わせが、観測値の時系列データに対して最適な組み合わせであると特定され、特定されたＭとλの値の組み合わせを適用した正規化拡張カルマンフィルタが、最適な推定モデルの推定結果とされる。
【００４９】
推定結果データ生成部２５は、最適モデル推定部２４から最適モデルの推定結果を受け取ると、最適モデルについて時刻ｔ＝Ｎにおける観測更新後の状態ＸN、λ、Ｍを含む、パラメータ推定結果を示す推定結果データを生成する。図３（Ｂ）に示すように、パラメータ推定結果には、状態ｘ1，・・・，ｘM、パラメータａ1，・・・，ａM，ｄt，λ、モデル誤差分散ｑ、観測誤差分散ｒが含まれている。モデル誤差分散ｑｅｓｔｉｍａｔｅｄと観測誤差分散ｒｅｓｔｉｍａｔｅｄは、以下の（２０）式、（２１）式により推定される。
【００５０】
【数１０】
【００５１】
上記（２０）式、（２１）式で表される、モデル誤差の分散ｑｅｓｔｉｍａｔｅｄと観測誤差の分散ｒｅｓｔｉｍａｔｅｄは、時系列データのダイナミクスを数理的に表現している。
【００５２】
予測値算出部２６は、最適モデル推定部２４から最適モデルの推定結果（時刻ｔ＝Ｎにおける観測更新後の状態ＸN|N、λ、Ｍ）を受け取ると、上記（９）式に従って、＾ＸＮ＋１｜Ｎを算出し、時刻ｔ＝Ｎ＋１における予測値Ｈ＾ＸＮ＋１｜Ｎを算出する。また、入力部１０により当該時刻ｔ＝Ｎ＋１の観測値を受け取ると、上記（１２）式、（１３）式に従って観測更新を行うと共に、上記（９）式、（１０）式に従って時間更新を行って、次の時刻ｔ＝Ｎ＋２における予測値Ｈ＾ＸＮ＋２｜Ｎ＋１を算出する。このように、図３（Ｃ）に示すように、エポックｔnと予測値（ｔn）との組（ｎ＝Ｎ＋１、Ｎ＋２、・・・）が順次算出される。
【００５３】
出力部３０により、推定結果データ生成部２５によって生成された推定結果を示すデータがユーザに出力される。また、出力部３０により、予測値算出部２６によって算出された予測値がユーザに出力される。
【００５４】
＜最適モデル推定装置の作用＞
次に、本実施の形態に係る最適モデル推定装置１００の作用について説明する。まず、時刻ｔ1〜ｔNの観測値からなる時系列データが、最適モデル推定装置１００に入力されると、最適モデル推定装置１００によって、入力された時系列データが、メモリ（図示省略）へ格納される。そして、最適モデル推定装置１００によって、図４に示す最適モデル推定処理ルーチンが実行される。
【００５５】
まず、ステップＳ１０１において、入力された観測値の時系列データを取得する。そして、ステップＳ１０２において、上記ステップＳ１０１で取得した観測値の時系列データの状態空間表現となるように、上記（５）式〜（８）式で示すガウス型モデルのモデルパラメータを推定する。
【００５６】
ステップＳ１０３では、上記ステップＳ１０２で推定されたガウス型モデルのモデルパラメータを出力部３０によりオペレータに対して出力し、オペレータの手計算により求められた上記（１１）式の値の入力を入力部１０により受け付けることにより、上記（９）式〜（１３）式で表わされる正規化拡張カルマンフィルタをデザインする。
【００５７】
そして、ステップＳ１０４において、予め用意されたλの値及び次数Ｍの値の各組み合わせについて、当該値の組み合わせを設定した正規化拡張カルマンフィルタのアルゴリズムに従って、上記ステップＳ１０１で取得した観測値の時系列データを用いて、時間更新及び観測更新を逐次行う。観測値の時系列データを用いた、時間更新及び観測更新を、λの値及び次数Ｍの値の全ての組み合わせについて行う。
【００５８】
次のステップＳ１０５では、上記ステップＳ１０４でλの値及び次数Ｍの値の全ての組み合わせについて行われた正規化拡張カルマンフィルタのアルゴリズムに従った時間更新及び観測更新の結果に基づいて、上記（１８）式に従って、最小のＡＩＣ（Ｍ，λ）を記録したλの値と次数Ｍの値の組み合わせを特定し、特定されたλの値と次数Ｍの値の組み合わせが適用された、正規化拡張カルマンフィルタの状態空間モデルを、最適モデルの推定結果とする。
【００５９】
そして、ステップＳ１０６において、最適モデルについて上記ステップＳ１０４で求められた時刻ｔ＝Ｎにおける観測更新後の状態ＸN|N、及び上記ステップＳ１０５で特定された最適なλの値と次数Ｍの値の組み合わせを示す推定結果データを生成し、出力部３０により出力する。ステップＳ１０７では、最適モデルについて上記ステップＳ１０４で求められた時刻ｔ＝Ｎにおける観測更新後の状態ＸN|N、及び上記ステップＳ１０５で特定された最適なλの値と次数Ｍの値の組み合わせに基づいて、次時刻ｔ＝Ｎ＋１の予測値を算出して、出力部３０により出力する。
【００６０】
次のステップＳ１０８では、入力部１０により入力された次時刻の観測値を取得する。そして、ステップＳ１０９において、上記（１２）式、（１３）式に従って、上記ステップＳ１０８で取得した観測値に基づいて観測更新を行って、上記（９）式、（１０）式に従って時間更新を行うことにより、当該時刻における予測値を算出して、出力部３０により出力する。
【００６１】
ステップＳ１１０では、予測を終了するか否かを判定する。例えば、予め定めた予測終了時刻に到達していない場合には、予測を終了しないと判定して、上記ステップＳ１０８へ戻り、次時刻の観測値を取得する。一方、予め定めた予測終了時刻に到達した場合には、予測を終了すると判定し、最適モデル推定処理ルーチンを終了する。
【００６２】
以上説明したように、本実施の形態に係る最適モデル推定装置によれば、モデル誤差分散を予測誤差分散で正規化した正規化誤差λを用いた正規化拡張カルマンフィルタのアルゴリズムに従って、観測値の時系列データに基づいて予測値を算出し、予測誤差に基づいて最適な正規化誤差λを特定することにより、誤差の分散値を推定せずに、高精度にモデリングを行うことができる。
【００６３】
また、通常は未知であり、かつ推定が困難であるモデル誤差と観測誤差に対してロバストな時系列モデリングを行うことが可能となる。その結果、誤差分散の設定ミスによるモデリング及び予測精度の低下を回避することができ、精度の高い予測が実現される。
【００６４】
また、求めた正規化誤差λからモデル誤差、観測誤差それぞれの分散値を推定することができるため、状態空間モデルの全未知パラメータが推定可能ということになり、対象時系列のダイナミクスを数理的に把握することができる。
【００６５】
カルマンフィルタのアルゴリズム全体を正規化することにより、未知量であったモデル誤差と観測誤差の分散を考慮する必要がなくなり、代わりに予測誤差分散に対するモデル誤差分散の割合を示す正規化誤差λがパラメータとしてアルゴリズムに導入される。赤池情報量に基づくモデル推定の段階で正規化誤差λをパラメイトライズすることで、正規化誤差λの最適値を推定することが可能となる。正規化誤差λの最適値を使用することで、モデル誤差と観測誤差の推定分散値も算出することができる。
【００６６】
なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
【００６７】
例えば、カルマンフィルタにより１時刻先を予測する場合を例に説明したが、２時刻先、３時刻先といった長期の予測を行うようにしてもよい。
【００６８】
また、モデル誤差分散ｑを予測誤差分散ｑ＋ｒで正規化した値を、正規化誤差λとする場合を例に説明したが、これに限定されるものではなく、観測誤差分散ｒを予測誤差分散ｑ＋ｒで正規化した値を、正規化誤差λとして用いて、正規化拡張カルマンフィルタをデザインするようにしてもよい。
【００６９】
また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能である。
【符号の説明】
【００７０】
１０ 入力部
２０ 演算部
２１ 状態空間モデリング部
２２ 正規化カルマンフィルタデザイン部
２３ 状態推定部
２４ 最適モデル推定部
２５ 推定結果データ生成部
２６ 予測値算出部
３０ 出力部
１００ 最適モデル推定装置
【技術分野】
【０００１】
本発明は、最適モデル推定装置、方法、及びプログラムに係り、特に、カルマンフィルタを用いた最適モデルを推定する最適モデル推定装置、方法、及びプログラムに関する。
【背景技術】
【０００２】
時系列データを状態空間でモデリングする際、モデル誤差の分散値と観測誤差の分散値とが未知となっている場合がほとんどである。状態空間では、カルマンフィルタによる状態の逐次推定及び推定された状態に基づく予測が可能であるが、誤差の分散値に関しては推測による設定が多い。
【０００３】
例えば、状態の推定誤差の状態空間モデルを考え、予測誤差の自己相関関数を基に誤差の分散値を推定する手法が知られている（非特許文献１）。
【先行技術文献】
【非特許文献】
【０００４】
【非特許文献１】“Ａ ｎｅｗ ａｕｔｏｃｏｖａｒｉａｎｃｅ ｌｅａｓｔ−ｓｑｕａｒｅｓ ｍｅｔｈｏｄ ｆｏｒ ｅｓｔｉｍａｔｉｎｇ ｎｏｉｓｅ ｃｏｖａｒｉａｎｃｅ”，Ｏｄｅｌｓｏｎ ｅｔ ａｌ.Ｊｏｕｒｎａｌ ｏｆ Ａｕｔｏｍａｔｉｃａ， ｖｏｌ.４２，ｐｐ.３０３−３０８，２００６.
【発明の概要】
【発明が解決しようとする課題】
【０００５】
しかしながら、上記の特許文献１に記載の技術では、誤差の分散値を推定しているため、推定した誤差の分散値が間違っている場合、誤差の分散値の間違った設定によりモデリングの精度が下がり、予測精度も低下してしまう、という問題が生じる。
【０００６】
本発明は、上記の事情を鑑みてなされたもので、誤差の分散値を推定せずに、高精度にモデリングを行うことができる最適モデル推定装置、方法、及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【０００７】
上記の目的を達成するために本発明に係る最適モデル推定装置は、観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定する状態空間モデリング手段と、前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段と、予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とする最適モデル推定手段と、を含んで構成されている。
【０００８】
本発明に係る最適モデル推定方法は、状態空間モデリング手段、状態推定手段、及び最適モデル推定手段を含む最適モデル推定装置における最適モデル推定方法であって、前記最適モデル推定装置は、前記状態空間モデリング手段によって、観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定するステップと、前記状態推定手段によって、前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出するステップと、前記最適モデル推定手段によって、予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とするステップと、を含んで実行することを特徴とする。
【０００９】
本発明に係るプログラムは、コンピュータを、観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定する状態空間モデリング手段、前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段、及び予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とする最適モデル推定手段として機能させるためのプログラムである。
【００１０】
本発明によれば、状態空間モデリング手段によって、観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定する。状態空間モデリング手段によって推定された前記モデルパラメータに基づいて、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタが定められる。
【００１１】
そして、状態推定手段によって、定められた正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する。最適モデル推定手段によって、予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とする。
【００１２】
このように、モデル誤差分散及び観測誤差分散の何れか一方を予測誤差分散で正規化した正規化誤差λを用いた正規化カルマンフィルタに従って、観測値の時系列データに基づいて予測値を算出し、予測誤差に基づいて最適な正規化誤差λを特定することにより、誤差の分散値を推定せずに、高精度にモデリングを行うことができる。
【発明の効果】
【００１３】
以上説明したように、本発明の最適モデル推定装置、方法、及びプログラムによれば、モデル誤差分散及び観測誤差分散の何れか一方を予測誤差分散で正規化した正規化誤差λを用いた正規化カルマンフィルタに従って、観測値の時系列データに基づいて予測値を算出し、予測誤差に基づいて最適な正規化誤差λを特定することにより、誤差の分散値を推定せずに、高精度にモデリングを行うことができる、という効果が得られる。
【図面の簡単な説明】
【００１４】
【図１】本発明の実施の形態に係る最適モデル推定装置の構成を示す概略図である。
【図２】観測値の時系列データを説明するための図である。
【図３】（Ａ）入力される観測値を示す図、（Ｂ）パラメータ推定結果を示す図、及び（Ｃ）出力される予測値を示す図である。
【図４】本発明の実施の形態に係る最適モデル推定装置における最適モデル推定処理ルーチンの内容を示すフローチャートである。
【発明を実施するための形態】
【００１５】
以下、図面を参照して本発明の実施の形態を詳細に説明する。
【００１６】
＜システム構成＞
本発明の実施の形態に係る最適モデル推定装置１００は、観測値の時系列データが入力され、状態空間の最適モデルを推定する。この最適モデル推定装置１００は、ＣＰＵと、ＲＡＭと、後述する最適モデル推定処理ルーチンを実行するためのプログラムを記憶したＲＯＭとを備えたコンピュータで構成され、機能的には次に示すように構成されている。図１に示すように、最適モデル推定装置１００は、入力部１０と、演算部２０と、出力部３０とを備えている。
【００１７】
入力部１０は、入力された観測値の時系列データを受け付ける。観測値は、スカラー値であり、すなわち、時系列データは、１次元時系列である。例えば、図２に示すように、エポックｔnと観測値（ｔn）との組（ｎ＝１、・・・、Ｎ）からなる時系列データである。
【００１８】
演算部２０は、状態空間モデリング部２１、正規化カルマンフィルタデザイン部２２、状態推定部２３、最適モデル推定部２４、推定結果データ生成部２５、及び予測値算出部２６を備えている。
【００１９】
状態空間モデリング部２１は、入力部１０により受け付けた観測値の時系列データに基づいて、以下に説明するように、時系列データを状態空間で表現したガウス型モデルのモデルパラメータを推定する。
【００２０】
まず、一般に、カルマンフィルタを表わすガウス型モデルは以下の（１）式、（２）式で表現される。
【００２１】
【数１】
【００２２】
ただし、Ｘｔは時刻ｔにおける状態ベクトルで、通常、非観測量（モデルの状態及び未知パラメータ）が取られる。ｗは正規白色ノイズのベクトルであり、モデル誤差が含まれる。上記（１）式に示されるように、状態ベクトルＸｔは線形（非線形）作用ｆとノイズベクトルｗにより時間更新される。
【００２３】
また、ｙｔは時刻ｔにおける観測値であり、状態Ｘｔの線形（非線形）作用ｇにより生成され、正規白色観測誤差μがのる。上記（２）式は、観測誤差を用いて状態ベクトルＸtと観測値ｙtとの関係を示している。
【００２４】
上記（１）式、（２）式で示される状態空間表現は、ＡＲ（自己回帰モデル）に代表される様々な時系列モデルを含む、汎用性が高いモデルである。以降、時系列モデルとして、以下の（３）式、（４）式で表わされるＡＲ＋ドリフトモデルを用いた場合を例として説明する。
【００２５】
【数２】
【００２６】
ただし、ωがモデル誤差に相当し、μが観測誤差に相当する。ここで、Ｎ（ｍ、ｖ）は平均ｍ、分散ｖの正規分布を表す。したがって、ｑはモデル誤差分散に相当し、ｒは観測誤差分散に相当する。このモデルの状態空間表現は、以下の（５）式〜（８）式で示される。
【００２７】
【数３】
【００２８】
Ｈは、状態を観測値に変換する入出力行列を表す。
【００２９】
状態空間モデリング部２１は、観測値の時系列データに基づいて、上記（５）式〜（８）式で表わされるガウス型モデルのモデルパラメータ（ｘt-1，・・・，ｘt-M、ａ1，・・・，ａM，ｄt，ｑ，ｒ）を推定する。
【００３０】
正規化カルマンフィルタデザイン部２２は、状態空間モデリング部２１によって推定されたガウス型モデルのモデルパラメータに基づいて、オペレータによる入力部１０からの入力を受け付けて、以下に説明するように、正規化拡張カルマンフィルタをデザインする。
【００３１】
状態空間モデリング部２１により状態空間表現された時系列データは、以下のアルゴリズムで表現される。まず、状態ベクトルを時間更新するアルゴリズムが、以下の（９）式〜（１１）式で表現される。
【００３２】
【数４】
【００３３】
上記（１１）式で表わされる^Ｆは、ｆの微分をオペレータが計算して入力することにより、求められる。
【００３４】
また、得られた観測値を用いる観測更新のアルゴリズムが、以下の（１２）式〜（１５）式で表される。
【００３５】
【数５】
【００３６】
Ｘｔ｜ｓはｓ期におけるｔ期の状態推定値を意味する。＾Ｘは状態Ｘの期待値である。Ｐ’は状態Ｘの正規化推定誤差分散の共分散行列であり、以下の（１６）式に示すように、推定誤差分散の共分散行列Ｐを予測誤差ｑ＋ｒで正規化したものである。
【００３７】
【数６】
【００３８】
この正規化プロセスにより、通常のカルマンフィルタのアルゴリズムからモデル誤差分散ｑと観測誤差分散ｒとが消え（Ｑ’＝ｄｉａｇ（λ、０、…、０））、代わりに以下の（１７）で定義するλが導入される。
【００３９】
【数７】
【００４０】
λは予測誤差分散ｑ＋ｒにおけるモデル誤差分散ｑの割合であり、０≦λ≦１の値をとるパラメータである。このλの値を設定することで、上記（９）式〜（１５）式で表現された正規化拡張カルマンフィルタのデザインが完了する。以下、λを、正規化誤差と称する。
【００４１】
このように、上記（１６）式に示す状態Ｘの正規化推定誤差分散の共分散行列と、上記（１７）式に示す正規化誤差とを用いた正規化拡張カルマンフィルタのアルゴリズムがデザインされる。
【００４２】
状態推定部２３は、予め用意されたλの各値（例えば、０から１まで０．１ずつ加算した各値）及び次数Ｍの各値（例えば、１から１０までの自然数）の組み合わせについて、正規化カルマンフィルタデザイン部２２によってデザインされ、かつ、当該λの値及び次数Ｍの値が設定された、上記（９）式〜（１５）式で表現された正規化拡張カルマンフィルタのアルゴリズムに従って、観測値の時系列データを用いて、時間更新及び観測更新を逐次行う。
【００４３】
最適モデル推定部２４は、λの値及び次数Ｍの値の各組み合わせについて、状態推定部２３によって逐次行われた時間更新及び観測更新の結果を利用して、以下の（１８）式に基づき、時系列データの最適モデルを推定する。
【００４４】
【数８】
【００４５】
ただし、ｖｔは上記（１５）式で示した予測誤差であり、Σｖｔは予測誤差ｖｔの分散である。ＡＩＣは赤池情報量であり、通常、ＡＲを想定した場合は次数Ｍ（何期前のデータまで使用するか）の関数である。しかしながら、本実施の形態で用いる正規化拡張カルマンフィルタのアルゴリズムに対しては、ＡＩＣは、Ｍとλの関数になる。ここで、上記（１８）で用いる予測誤差ｑ＋ｒには、以下の（１９）式に従って推定された値を適用すればよい。
【００４６】
【数９】
【００４７】
ただし、Ｐt’[1,1]は、行列Ｐt’の要素（１，１）である。
【００４８】
上記のように、Ｍとλをパラメトライズし、最小のＡＩＣ（Ｍ、λ）を記録したときのＭとλの値の組み合わせが、観測値の時系列データに対して最適な組み合わせであると特定され、特定されたＭとλの値の組み合わせを適用した正規化拡張カルマンフィルタが、最適な推定モデルの推定結果とされる。
【００４９】
推定結果データ生成部２５は、最適モデル推定部２４から最適モデルの推定結果を受け取ると、最適モデルについて時刻ｔ＝Ｎにおける観測更新後の状態ＸN、λ、Ｍを含む、パラメータ推定結果を示す推定結果データを生成する。図３（Ｂ）に示すように、パラメータ推定結果には、状態ｘ1，・・・，ｘM、パラメータａ1，・・・，ａM，ｄt，λ、モデル誤差分散ｑ、観測誤差分散ｒが含まれている。モデル誤差分散ｑｅｓｔｉｍａｔｅｄと観測誤差分散ｒｅｓｔｉｍａｔｅｄは、以下の（２０）式、（２１）式により推定される。
【００５０】
【数１０】
【００５１】
上記（２０）式、（２１）式で表される、モデル誤差の分散ｑｅｓｔｉｍａｔｅｄと観測誤差の分散ｒｅｓｔｉｍａｔｅｄは、時系列データのダイナミクスを数理的に表現している。
【００５２】
予測値算出部２６は、最適モデル推定部２４から最適モデルの推定結果（時刻ｔ＝Ｎにおける観測更新後の状態ＸN|N、λ、Ｍ）を受け取ると、上記（９）式に従って、＾ＸＮ＋１｜Ｎを算出し、時刻ｔ＝Ｎ＋１における予測値Ｈ＾ＸＮ＋１｜Ｎを算出する。また、入力部１０により当該時刻ｔ＝Ｎ＋１の観測値を受け取ると、上記（１２）式、（１３）式に従って観測更新を行うと共に、上記（９）式、（１０）式に従って時間更新を行って、次の時刻ｔ＝Ｎ＋２における予測値Ｈ＾ＸＮ＋２｜Ｎ＋１を算出する。このように、図３（Ｃ）に示すように、エポックｔnと予測値（ｔn）との組（ｎ＝Ｎ＋１、Ｎ＋２、・・・）が順次算出される。
【００５３】
出力部３０により、推定結果データ生成部２５によって生成された推定結果を示すデータがユーザに出力される。また、出力部３０により、予測値算出部２６によって算出された予測値がユーザに出力される。
【００５４】
＜最適モデル推定装置の作用＞
次に、本実施の形態に係る最適モデル推定装置１００の作用について説明する。まず、時刻ｔ1〜ｔNの観測値からなる時系列データが、最適モデル推定装置１００に入力されると、最適モデル推定装置１００によって、入力された時系列データが、メモリ（図示省略）へ格納される。そして、最適モデル推定装置１００によって、図４に示す最適モデル推定処理ルーチンが実行される。
【００５５】
まず、ステップＳ１０１において、入力された観測値の時系列データを取得する。そして、ステップＳ１０２において、上記ステップＳ１０１で取得した観測値の時系列データの状態空間表現となるように、上記（５）式〜（８）式で示すガウス型モデルのモデルパラメータを推定する。
【００５６】
ステップＳ１０３では、上記ステップＳ１０２で推定されたガウス型モデルのモデルパラメータを出力部３０によりオペレータに対して出力し、オペレータの手計算により求められた上記（１１）式の値の入力を入力部１０により受け付けることにより、上記（９）式〜（１３）式で表わされる正規化拡張カルマンフィルタをデザインする。
【００５７】
そして、ステップＳ１０４において、予め用意されたλの値及び次数Ｍの値の各組み合わせについて、当該値の組み合わせを設定した正規化拡張カルマンフィルタのアルゴリズムに従って、上記ステップＳ１０１で取得した観測値の時系列データを用いて、時間更新及び観測更新を逐次行う。観測値の時系列データを用いた、時間更新及び観測更新を、λの値及び次数Ｍの値の全ての組み合わせについて行う。
【００５８】
次のステップＳ１０５では、上記ステップＳ１０４でλの値及び次数Ｍの値の全ての組み合わせについて行われた正規化拡張カルマンフィルタのアルゴリズムに従った時間更新及び観測更新の結果に基づいて、上記（１８）式に従って、最小のＡＩＣ（Ｍ，λ）を記録したλの値と次数Ｍの値の組み合わせを特定し、特定されたλの値と次数Ｍの値の組み合わせが適用された、正規化拡張カルマンフィルタの状態空間モデルを、最適モデルの推定結果とする。
【００５９】
そして、ステップＳ１０６において、最適モデルについて上記ステップＳ１０４で求められた時刻ｔ＝Ｎにおける観測更新後の状態ＸN|N、及び上記ステップＳ１０５で特定された最適なλの値と次数Ｍの値の組み合わせを示す推定結果データを生成し、出力部３０により出力する。ステップＳ１０７では、最適モデルについて上記ステップＳ１０４で求められた時刻ｔ＝Ｎにおける観測更新後の状態ＸN|N、及び上記ステップＳ１０５で特定された最適なλの値と次数Ｍの値の組み合わせに基づいて、次時刻ｔ＝Ｎ＋１の予測値を算出して、出力部３０により出力する。
【００６０】
次のステップＳ１０８では、入力部１０により入力された次時刻の観測値を取得する。そして、ステップＳ１０９において、上記（１２）式、（１３）式に従って、上記ステップＳ１０８で取得した観測値に基づいて観測更新を行って、上記（９）式、（１０）式に従って時間更新を行うことにより、当該時刻における予測値を算出して、出力部３０により出力する。
【００６１】
ステップＳ１１０では、予測を終了するか否かを判定する。例えば、予め定めた予測終了時刻に到達していない場合には、予測を終了しないと判定して、上記ステップＳ１０８へ戻り、次時刻の観測値を取得する。一方、予め定めた予測終了時刻に到達した場合には、予測を終了すると判定し、最適モデル推定処理ルーチンを終了する。
【００６２】
以上説明したように、本実施の形態に係る最適モデル推定装置によれば、モデル誤差分散を予測誤差分散で正規化した正規化誤差λを用いた正規化拡張カルマンフィルタのアルゴリズムに従って、観測値の時系列データに基づいて予測値を算出し、予測誤差に基づいて最適な正規化誤差λを特定することにより、誤差の分散値を推定せずに、高精度にモデリングを行うことができる。
【００６３】
また、通常は未知であり、かつ推定が困難であるモデル誤差と観測誤差に対してロバストな時系列モデリングを行うことが可能となる。その結果、誤差分散の設定ミスによるモデリング及び予測精度の低下を回避することができ、精度の高い予測が実現される。
【００６４】
また、求めた正規化誤差λからモデル誤差、観測誤差それぞれの分散値を推定することができるため、状態空間モデルの全未知パラメータが推定可能ということになり、対象時系列のダイナミクスを数理的に把握することができる。
【００６５】
カルマンフィルタのアルゴリズム全体を正規化することにより、未知量であったモデル誤差と観測誤差の分散を考慮する必要がなくなり、代わりに予測誤差分散に対するモデル誤差分散の割合を示す正規化誤差λがパラメータとしてアルゴリズムに導入される。赤池情報量に基づくモデル推定の段階で正規化誤差λをパラメイトライズすることで、正規化誤差λの最適値を推定することが可能となる。正規化誤差λの最適値を使用することで、モデル誤差と観測誤差の推定分散値も算出することができる。
【００６６】
なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
【００６７】
例えば、カルマンフィルタにより１時刻先を予測する場合を例に説明したが、２時刻先、３時刻先といった長期の予測を行うようにしてもよい。
【００６８】
また、モデル誤差分散ｑを予測誤差分散ｑ＋ｒで正規化した値を、正規化誤差λとする場合を例に説明したが、これに限定されるものではなく、観測誤差分散ｒを予測誤差分散ｑ＋ｒで正規化した値を、正規化誤差λとして用いて、正規化拡張カルマンフィルタをデザインするようにしてもよい。
【００６９】
また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能である。
【符号の説明】
【００７０】
１０ 入力部
２０ 演算部
２１ 状態空間モデリング部
２２ 正規化カルマンフィルタデザイン部
２３ 状態推定部
２４ 最適モデル推定部
２５ 推定結果データ生成部
２６ 予測値算出部
３０ 出力部
１００ 最適モデル推定装置
【特許請求の範囲】
【請求項１】
観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定する状態空間モデリング手段と、
前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段と、
予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とする最適モデル推定手段と、
を含む最適モデル推定装置。
【請求項２】
状態空間モデリング手段、状態推定手段、及び最適モデル推定手段を含む最適モデル推定装置における最適モデル推定方法であって、
前記最適モデル推定装置は、
前記状態空間モデリング手段によって、観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定するステップと、
前記状態推定手段によって、前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出するステップと、
前記最適モデル推定手段によって、予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とするステップと、
を含んで実行することを特徴とする最適モデル推定方法。
【請求項３】
コンピュータを、
観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定する状態空間モデリング手段、
前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段、及び
予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とする最適モデル推定手段
として機能させるためのプログラム。
【請求項１】
観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定する状態空間モデリング手段と、
前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段と、
予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とする最適モデル推定手段と、
を含む最適モデル推定装置。
【請求項２】
状態空間モデリング手段、状態推定手段、及び最適モデル推定手段を含む最適モデル推定装置における最適モデル推定方法であって、
前記最適モデル推定装置は、
前記状態空間モデリング手段によって、観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定するステップと、
前記状態推定手段によって、前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出するステップと、
前記最適モデル推定手段によって、予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とするステップと、
を含んで実行することを特徴とする最適モデル推定方法。
【請求項３】
コンピュータを、
観測値の時系列データに基づいて、モデル誤差を用いて状態ベクトルＸを時間更新するための状態更新式、及び観測誤差を用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いた状態空間モデルで表わされるカルマンフィルタのモデルパラメータを推定する状態空間モデリング手段、
前記状態空間モデリング手段によって推定された前記モデルパラメータに基づいて定められた、前記状態ベクトルＸの推定誤差分散の共分散行列を、モデル誤差分散及び観測誤差分散の和である予測誤差分散で正規化した正規化推定誤差分散共分散行列と、前記モデル誤差分散及び前記観測誤差分散の何れか一方を前記予測誤差分散で正規化した正規化誤差λとを用いた状態ベクトルＸの時間更新及び観測更新を行う正規化カルマンフィルタに従って、前記時系列データの各観測値を用いて予測値を時系列に算出すると共に、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段、及び
予め用意された前記正規化誤差λの各値を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化誤差λを特定し、最適な正規化誤差λを用いた前記正規化カルマンフィルタを、最適モデルの推定結果とする最適モデル推定手段
として機能させるためのプログラム。
【図１】
【図２】
【図３】
【図４】
【図２】
【図３】
【図４】
【公開番号】特開２０１３−６１７６８（Ｐ２０１３−６１７６８Ａ）
【公開日】平成２５年４月４日（２０１３．４．４）
【国際特許分類】
【出願番号】特願２０１１−１９９４７５（Ｐ２０１１−１９９４７５）
【出願日】平成２３年９月１３日（２０１１．９．１３）
【出願人】（０００００４２２６）日本電信電話株式会社 (13,992)
【Ｆターム（参考）】
【公開日】平成２５年４月４日（２０１３．４．４）
【国際特許分類】
【出願日】平成２３年９月１３日（２０１１．９．１３）
【出願人】（０００００４２２６）日本電信電話株式会社 (13,992)
【Ｆターム（参考）】
[ Back to top ]