JP5738778B2 - 最適モデル推定装置、方法、及びプログラム - Google Patents

Description

本発明は、最適モデル推定装置、方法、及びプログラムに係り、特に、最適な状態空間モデルを推定する最適モデル推定装置、方法、及びプログラムに関する。

従来より、線形モデルの残差(観測値と予測値の差)とその共分散の関係を線形方程式で書き下し、最小二乗法でシステムノイズと観測ノイズを推定する手法が提案されている（非特許文献１）。

Brian J. Odelson et al., "A New Autocovariance Least-Squares Method for Estimating Noise Covariances", Technical report, TWMCC,2003年.

しかしながら、上記の特許文献１に記載の技術では、線形モデルを対象としているため、非線形モデルのモデリングを行うことができない、という問題が生じる。

本発明は、上記の事情を鑑みてなされたもので、非線形モデルを対象として、最適な状態空間モデルを推定することができる最適モデル推定装置、方法、及びプログラムを提供することを目的とする。

上記の目的を達成するために本発明に係る最適モデル推定装置は、システムノイズを用いて状態ベクトルＸを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いて定められた状態空間モデルについて、前記状態ベクトルＸの推定誤差分散の共分散行列を、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和で正規化した正規化推定誤差分散共分散行列と、前記システムノイズの共分散行列を前記総和で正規化した正規化システムノイズ共分散行列と、前記観測ノイズの分散値を前記総和で正規化した正規化観測ノイズ分散とを用いたシグマ点列の発生、状態ベクトルＸの状態更新、及び観測更新を行うＮｏｒｍａｌｉｚｅｄ ＵＫＦ（Ｕｎｓｃｅｎｔｅｄ Ｋａｌｍａｎ Ｆｉｌｔｅｒ）に従って、前記観測値の時系列データに基づいて前記状態ベクトルＸを逐次推定すると共に、予測値を時系列に算出し、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段と、予め用意された前記正規化システムノイズ共分散行列及び前記正規化観測ノイズ分散の組み合わせの各々を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を特定し、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を用いた前記状態空間モデルを、最適モデルの推定結果とする最適モデル推定手段と、を含んで構成されている。

本発明に係る最適モデル推定方法は、状態推定手段及び最適モデル推定手段を含む最適モデル推定装置における最適モデル推定方法であって、前記最適モデル推定装置は、状態推定手段によって、システムノイズを用いて状態ベクトルＸを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いて定められた状態空間モデルについて、前記状態ベクトルＸの推定誤差分散の共分散行列を、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和で正規化した正規化推定誤差分散共分散行列と、前記システムノイズの共分散行列を前記総和で正規化した正規化システムノイズ共分散行列と、前記観測ノイズの分散値を前記総和で正規化した正規化観測ノイズ分散とを用いたシグマ点列の発生、状態ベクトルＸの状態更新、及び観測更新を行うＮｏｒｍａｌｉｚｅｄ ＵＫＦ（Ｕｎｓｃｅｎｔｅｄ Ｋａｌｍａｎ Ｆｉｌｔｅｒ）に従って、前記観測値の時系列データに基づいて前記状態ベクトルＸを逐次推定すると共に、予測値を時系列に算出し、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出するステップと、最適モデル推定手段によって、予め用意された前記正規化システムノイズ共分散行列及び前記正規化観測ノイズ分散の組み合わせの各々を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を特定し、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を用いた前記状態空間モデルを、最適モデルの推定結果とするステップと、を含んで実行することを特徴とする。

本発明に係るプログラムは、コンピュータを、システムノイズを用いて状態ベクトルＸを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いて定められた状態空間モデルについて、前記状態ベクトルＸの推定誤差分散の共分散行列を、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和で正規化した正規化推定誤差分散共分散行列と、前記システムノイズの共分散行列を前記総和で正規化した正規化システムノイズ共分散行列と、前記観測ノイズの分散値を前記総和で正規化した正規化観測ノイズ分散とを用いたシグマ点列の発生、状態ベクトルＸの状態更新、及び観測更新を行うＮｏｒｍａｌｉｚｅｄ ＵＫＦ（Ｕｎｓｃｅｎｔｅｄ Ｋａｌｍａｎ Ｆｉｌｔｅｒ）に従って、前記観測値の時系列データに基づいて前記状態ベクトルＸを逐次推定すると共に、予測値を時系列に算出し、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段、及び予め用意された前記正規化システムノイズ共分散行列及び前記正規化観測ノイズ分散の組み合わせの各々を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を特定し、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を用いた前記状態空間モデルを、最適モデルの推定結果とする最適モデル推定手段として機能させるためのプログラムである。

本発明によれば、状態推定手段によって、システムノイズを用いて状態ベクトルＸを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いて定められた状態空間モデルについて、前記状態ベクトルＸの推定誤差分散の共分散行列を、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和で正規化した正規化推定誤差分散共分散行列と、前記システムノイズの共分散行列を前記総和で正規化した正規化システムノイズ共分散行列と、前記観測ノイズの分散値を前記総和で正規化した正規化観測ノイズ分散とを用いたシグマ点列の発生、状態ベクトルＸの状態更新、及び観測更新を行うＮｏｒｍａｌｉｚｅｄ ＵＫＦ（Ｕｎｓｃｅｎｔｅｄ Ｋａｌｍａｎ Ｆｉｌｔｅｒ）に従って、前記観測値の時系列データに基づいて前記状態ベクトルＸを逐次推定すると共に、予測値を時系列に算出し、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する。そして、最適モデル推定手段によって、予め用意された前記正規化システムノイズ共分散行列及び前記正規化観測ノイズ分散の組み合わせの各々を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を特定し、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を用いた前記状態空間モデルを、最適モデルの推定結果とする。

このように、各々正規化した推定誤差分散共分散行列、システムノイズの共分散行列、及び観測ノイズの分散値を用いたＮｏｒｍａｌｉｚｅｄ ＵＫＦに従って状態ベクトルの逐次推定を行って、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を特定することにより、非線形モデルを対象として、最適な状態空間モデルを推定することができる。

以上説明したように、本発明の最適モデル推定装置、方法、及びプログラムによれば、各々正規化した推定誤差分散共分散行列、システムノイズの共分散行列、及び観測ノイズの分散値を用いたＮｏｒｍａｌｉｚｅｄ ＵＫＦに従って状態ベクトルの逐次推定を行って、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を特定することにより、非線形モデルを対象として、最適な状態空間モデルを推定することができる、という効果が得られる。

本発明の実施の形態に係る最適モデル推定装置の構成を示す概略図である。 観測値の時系列データを説明するための図である。 出力データを示す図である。 本発明の実施の形態に係る最適モデル推定装置における最適モデル推定処理ルーチンの内容を示すフローチャートである。 （Ａ）テスト１に対する各モデルの評価結果を示す図、及び（Ｂ）テスト２に対する各モデルの評価結果を示す図である。

以下、図面を参照して本発明の実施の形態を詳細に説明する。

＜発明の概要＞
経済時系列データに代表される現実の時系列は時間と共に大きく変動し、一般的に解析が困難である。この様な時系列は非定常と呼ばれ、想定するモデルも非線形になる。本発明は、非線形非定常時系列モデルに対する一つの状態推定アルゴリズムであるUnscented Kalman Filter(UKF)において、非線形時系列モデル中に存在する様々なノイズの分散値が未知である場合でも、その総和を1に抑えることで安定的な推定を可能にするＮｏｒｍａｌｉｚｅｄ ＵＫＦ（Ｕｎｓｃｅｎｔｅｄ Ｋａｌｍａｎ Ｆｉｌｔｅｒ）に関するものである。

＜システム構成＞
本発明の実施の形態に係る最適モデル推定装置１００は、観測値の時系列データが入力され、状態空間の最適モデルを推定する。この最適モデル推定装置１００は、ＣＰＵと、ＲＡＭと、後述する最適モデル推定処理ルーチンを実行するためのプログラムを記憶したＲＯＭとを備えたコンピュータで構成され、機能的には次に示すように構成されている。図１に示すように、最適モデル推定装置１００は、入力部１０と、演算部２０と、出力部３０とを備えている。

入力部１０は、入力された観測値の時系列データを受け付ける。観測値は、スカラー値であり、すなわち、時系列データは、１次元時系列である。例えば、図２に示すように、エポックｔ_nと観測値（ｔ_n）との組（ｎ＝１、・・・、Ｎ）からなる時系列データである。

また、入力部１０は、オペレータから、後述する状態空間モデルの定義の入力を受け付ける。

演算部２０は、状態空間モデリング部２１、状態推定部２２、最適モデル推定部２３、及び推定結果データ生成部２４を備えている。

状態空間モデリング部２１は、オペレータによる入力部１０からの状態空間モデルの定義の入力を受け付けて、以下に説明するように、時系列データを状態空間で表現した非線形ガウス型モデルを定義する。

まず、以下の（１）式、（２）式に示す非線形ガウス型モデリングを行う。

ただし、Ｘ_ｔは時刻ｔにおける状態ベクトルで、通常、非観測量が取られる。状態ベクトルＸ_tを時系列データy_tを用いて推定することを目的とする。wはシステムノイズベクトルであり、正規分布に従う白色ノイズベクトルである。μは観測ノイズであり、正規分布に従う白色ノイズである。また、Ｎ（ｍ、ｖ）は平均ｍ、分散ｖの正規分布を表す。したがって、Ｑはシステムノイズの共分散行列に相当し、ｒは観測ノイズの分散値に相当する。上記（２）式に示されるように、y_tは時刻tにおける観測値で、状態X_tの線形変換Hにより生成される。上記（２）式は、観測ノイズを用いて状態ベクトルＸtと観測値ｙtとの関係を示している。f()が非線形箇所であり、f()の型はモデルにより異なるが、循環変動+簡単なトレンドモデルでは以下の（３）式のように状態が定義され、状態更新式が組まれる。

ここで、x_t, ・・・,x_t−M+1,a₁,・・・,a_Mまでは自己回帰過程で表現された循環変動成分に関する変量、T_tがトレンド成分の変量である。また、状態更新式のｆ（Ｘ_t-1）、ｗは、以下の（４）式で表わされる。

ここでｗはシステムノイズベクトルであるが、通常ω₁とω_2×M+1以外は0である。この循環変動+簡単なトレンドモデルでは、上記（２）式の観測方程式のＨは、1行（2×M+1）列のベクトルで、1番目、（2×M+1）番目の要素が1であり、その他は0である。

このように、状態空間モデリング部２１は、オペレータからの入力を受け付けて、上記（４）式に示すようなｆ（）の型を定義すると共に、上記（２）式のＨを定義する。

また、状態空間モデリング部２１は、後述する正規化されたシステムノイズの共分散行列の対角成分及び観測ノイズの分散値を0〜1まで振り、正規化されたシステムノイズの共分散行列及び観測ノイズの分散値の組み合わせの各々を用意する。例えば、共分散行列Ｑの対角成分の各値及び分散値ｒの値の総和が１となるように正規化されたシステムノイズの共分散行列及び観測ノイズの分散値の組み合わせを複数用意する。

状態推定部２２は、入力された時系列データと、状態空間モデリング部２１でモデリングされた状態空間モデルとに基づき、状態ベクトルの逐次推定を行う。モデルは非線形であるので、Unscented Kalman Filter(UKF)をベースにしたNormalized UKFにより状態ベクトルX_tの推定を行う。通常のUKFアルゴリズムは、以下の（５）式〜（１５）式で表されるシグマ点列の発生、状態更新、及び観測更新（フィルタ）からなる。

ただし、λは、予め定められた定数であり、Ｌは、状態ベクトルＸの次元数である。￣ｙは観測値の予測値であり、￣ｙ＝Ｈ￣Ｘである。ｖは予測誤差であり、Ｐは状態ベクトルＸの推定誤差分散の共分散行列である。

本実施の形態で用いるNormalized UKFでは、上記の通常のアルゴリズムに存在する￣P,＾Pを、全ノイズの分散値の総和、つまりシステムノイズｗの共分散行列Ｑの対角成分の総和と観測ノイズμの分散値rとの和traceＱ+rで割る。例えば、上記（１５）式の￣Pは以下の（１６）式で定義される￣P’に置き換えられる。

この操作によって、上記（１１）式のシステムノイズ共分散行列Ｑや上記（１２）式の観測ノイズの分散値rは、以下の（１７）式、（１８）式に示すように正規化される。つまり、システムノイズ共分散行列Ｑの要素及び観測ノイズの分散値rは全て１以下になる。また、システムノイズ共分散行列Ｑの要素及び観測ノイズの分散値の総和が１となる。

この一連の正規化の操作によってアルゴリズムの形は適宜変更されるが、上記（８）式、（９）式に含まれるαの設定値を小さくとることで通常のUKFのアルゴリズムと一致させることができる。これは以下のように示される。

ここでM・α²が1になるようにαを設定することで、正規化数Mのシグマ点列発生への影響が消える。つまり以下の（２３）式、（２４）式のようになる。

上記（２４）式は状態が正規分布に従うと仮定した時の通常のシグマ点列発生である。この場合L+κ=3に設定することが多い。

このように、状態推定部２２は、入力された観測値の時系列データに基づいて、通常のＵＫＦのアルゴリズムにおいて、正規化された状態ベクトルＸの推定誤差分散の共分散行列￣P’,＾P’、正規化されたシステムノイズ共分散行列Ｑ’、観測ノイズの分散値r’に置き換えたNormalized UKFのアルゴリズム（シグマ点列の発生、状態更新、観測更新）に従って、状態ベクトルＸの逐次推定を行うと共に、予測値を時系列に算出し、時系列に算出された予測値と対応する観測値との予測誤差を各々算出する。

また、状態推定部２２は、状態ベクトルＸの逐次推定を、状態空間モデリング部２１で用意した、正規化されたシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の組み合わせの各々を用いた場合についてそれぞれ行う。

最適モデル推定部２３は、正規化されたシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の組み合わせの各々について、状態推定部２２によって逐次行われたシグマ点列の発生、状態更新、及び観測更新の結果を利用して、以下の（２５）式〜（２９）式に基づき、モデルの性能を意味する赤池情報量ＡＩＣが最小になったところで、最適モデルを出力する。つまり、モデルの同定と同時にノイズ分散値の最適推定も行う。

ただし、ｖは上記（１３）式で示した予測誤差であり、ＡＩＣは赤池情報量である。ここで、正規化定数Mは以下の（３０）式で推定される。

最小のＡＩＣを記録したときの正規化されたシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の組み合わせが、ノイズ分散値の最適な組み合わせであると特定され、特定されたシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の組み合わせを適用した状態空間モデルが、最適な推定モデルの推定結果とされる。

推定結果データ生成部２４は、最適モデル推定部２３から最適モデルの推定結果を受け取ると、図３に示すように、最適モデルについて時刻ｔ＝Ｎにおける観測更新後の状態ＸNを含む推定モデルパラメータ及び推定システムノイズ分散値を示す推定結果データを生成する。パラメータ推定結果には、状態ｘ_t，・・・，ｘ_t-M+1、パラメータａ_t，・・・，ａ_t-M+1，Ｔ_t，システムノイズ分散Ｑ’、観測ノイズ分散ｒ’が含まれている。

出力部３０により、推定結果データ生成部２４によって生成された推定結果を示すデータがユーザに出力される。

＜最適モデル推定装置の作用＞
次に、本実施の形態に係る最適モデル推定装置１００の作用について説明する。まず、オペレータにより、状態空間モデルの定義が最適モデル推定装置１００に入力されると、最適モデル推定装置１００によって、入力された状態空間モデルの定義が、メモリ（図示省略）へ格納される。また、時刻ｔ₁〜ｔ_Nの観測値からなる時系列データが、最適モデル推定装置１００に入力されると、最適モデル推定装置１００によって、入力された時系列データが、メモリへ格納される。そして、最適モデル推定装置１００によって、図４に示す最適モデル推定処理ルーチンが実行される。

まず、ステップＳ１０１において、入力された観測値の時系列データを取得する。そして、ステップＳ１０２において、入力された状態空間モデルの定義に基づいて、上記（１）式、（２）式の非線形ガウス型モデルを設定すると共に、正規化されたシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の組み合わせを複数生成する。

そして、ステップＳ１０３において、上記ステップＳ１０２で生成された、正規化されたシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の組み合わせの各々について、当該組み合わせを設定したＮｏｒｍａｌｉｚｅｄ ＵＫＦのアルゴリズムに従って、上記ステップＳ１０１で取得した観測値の時系列データを用いて、シグマ点列の発生、状態更新、及び観測更新を逐次行って、状態ベクトルＸの逐次推定を行う。観測値の時系列データを用いた状態ベクトルＸの逐次推定を、正規化されたシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の全ての組み合わせについて行う。

次のステップＳ１０４では、上記ステップＳ１０３で正規化されたシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の全ての組み合わせについて行われたＮｏｒｍａｌｉｚｅｄ ＵＫＦのアルゴリズムに従ったシグマ点列の発生、状態更新、及び観測更新の結果に基づいて、上記（２９）式に従って、最小のＡＩＣを記録したシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の組み合わせを特定し、特定されたシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の組み合わせが適用された状態空間モデルを、最適モデルの推定結果とする。

そして、ステップＳ１０５において、最適モデルについて上記ステップＳ１０３で求められた時刻ｔ＝Ｎにおける観測更新後の状態ベクトルＸ、及び上記ステップＳ１０５で特定された最適なシステムノイズの共分散行列Ｑ’及び観測ノイズの分散値ｒ’の組み合わせを含む推定結果データを生成し、出力部３０により出力し、最適モデル推定処理ルーチンを終了する。

＜実験結果＞
次に、人工データを利用した実験結果について説明する。テスト１では、時系列モデルが自己回帰項のみ(2次AR、AR係数は1.5,−0.9)で、この項にシステムノイズが0.07、観測ノイズが0.03のる問題とした。テスト2は、時系列モデルが自己回帰項＋トレンド項のモデルでそれぞれの項にシステムノイズが0.07、0.03のり、観測ノイズがのらない問題とした。それぞれ人工データを発生し、本発明の手法でモデルの同定がどれほどできるか評価した。

テスト１に対する評価結果を図５（Ａ）に示す、テスト２に対する評価結果を図５（Ｂ）に示す。テスト１、２の双方とも、システムノイズと観測ノイズの値をλ（＝システムノイズ／（システムノイズ＋観測ノイズ））を用いて振った結果、AICが最小になったモデルが正解モデルになっており(システムノイズと観測ノイズがλの値に対応しており、自己回帰係数も正しく推定されている)、本発明の手法の有効性を示している。

以上説明したように、本実施の形態に係る最適モデル推定装置によれば、各々正規化した、推定誤差分散共分散行列、システムノイズの共分散行列、及び観測ノイズの分散値を用いたＮｏｒｍａｌｉｚｅｄ ＵＫＦに従って状態ベクトルＸの逐次推定を行って、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を特定することにより、非線形モデルを対象として、最適な状態空間モデルを推定することができる。

また、本実施の形態で用いる状態推定アルゴリズムは、非線形ガウス型でモデル化された時系列の状態を高精度で推定するUnscented Kalman Filter(UKF)をベースにして、モデルに存在する全ノイズ分散値の総和を1に正規化する手法を用いる。通常、ノイズの分散値は未知であるので、解析者の経験や勘により設定する場合がほとんどであり、状態の推定結果にも大きく影響を与える。本実施の形態では、モデル中に存在する全てのノイズ分散値の総和を1に抑えることで、安定的な解析、状態推定を可能にする。

また、赤池情報量等に代表されるモデル選択基準を使うことで、ノイズの統計量の推定までを含んだ最適なモデルの同定も可能になる。

なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。

例えば、状態空間モデルの定義を予め設定しておき、状態空間モデルの定義の入力を不要としてもよい。

また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能である。

１０ 入力部
２０ 演算部
２１ 状態空間モデリング部
２２ 状態推定部
２３ 最適モデル推定部
２４ 推定結果データ生成部
１００ 最適モデル推定装置

Claims (3)

1. システムノイズを用いて状態ベクトルＸを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いて定められた状態空間モデルについて、前記状態ベクトルＸの推定誤差分散の共分散行列を、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和で正規化した正規化推定誤差分散共分散行列と、前記システムノイズの共分散行列を前記総和で正規化した正規化システムノイズ共分散行列と、前記観測ノイズの分散値を前記総和で正規化した正規化観測ノイズ分散とを用いたシグマ点列の発生、状態ベクトルＸの状態更新、及び観測更新を行うＮｏｒｍａｌｉｚｅｄ ＵＫＦ（Ｕｎｓｃｅｎｔｅｄ Ｋａｌｍａｎ Ｆｉｌｔｅｒ）に従って、前記観測値の時系列データに基づいて前記状態ベクトルＸを逐次推定すると共に、予測値を時系列に算出し、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段と、
   予め用意された前記正規化システムノイズ共分散行列及び前記正規化観測ノイズ分散の組み合わせの各々を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を特定し、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を用いた前記状態空間モデルを、最適モデルの推定結果とする最適モデル推定手段と、
   を含む最適モデル推定装置。
2. 状態推定手段及び最適モデル推定手段を含む最適モデル推定装置における最適モデル推定方法であって、
   前記最適モデル推定装置は、
   状態推定手段によって、システムノイズを用いて状態ベクトルＸを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いて定められた状態空間モデルについて、前記状態ベクトルＸの推定誤差分散の共分散行列を、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和で正規化した正規化推定誤差分散共分散行列と、前記システムノイズの共分散行列を前記総和で正規化した正規化システムノイズ共分散行列と、前記観測ノイズの分散値を前記総和で正規化した正規化観測ノイズ分散とを用いたシグマ点列の発生、状態ベクトルＸの状態更新、及び観測更新を行うＮｏｒｍａｌｉｚｅｄ ＵＫＦ（Ｕｎｓｃｅｎｔｅｄ Ｋａｌｍａｎ Ｆｉｌｔｅｒ）に従って、前記観測値の時系列データに基づいて前記状態ベクトルＸを逐次推定すると共に、予測値を時系列に算出し、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出するステップと、
   最適モデル推定手段によって、予め用意された前記正規化システムノイズ共分散行列及び前記正規化観測ノイズ分散の組み合わせの各々を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を特定し、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を用いた前記状態空間モデルを、最適モデルの推定結果とするステップと、
   を含んで実行することを特徴とする最適モデル推定方法。
3. コンピュータを、
   システムノイズを用いて状態ベクトルＸを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルＸと観測値との関係を示す観測方程式を用いて定められた状態空間モデルについて、前記状態ベクトルＸの推定誤差分散の共分散行列を、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和で正規化した正規化推定誤差分散共分散行列と、前記システムノイズの共分散行列を前記総和で正規化した正規化システムノイズ共分散行列と、前記観測ノイズの分散値を前記総和で正規化した正規化観測ノイズ分散とを用いたシグマ点列の発生、状態ベクトルＸの状態更新、及び観測更新を行うＮｏｒｍａｌｉｚｅｄ ＵＫＦ（Ｕｎｓｃｅｎｔｅｄ Ｋａｌｍａｎ Ｆｉｌｔｅｒ）に従って、前記観測値の時系列データに基づいて前記状態ベクトルＸを逐次推定すると共に、予測値を時系列に算出し、時系列に算出された前記予測値と対応する観測値との予測誤差を各々算出する状態推定手段、及び
   予め用意された前記正規化システムノイズ共分散行列及び前記正規化観測ノイズ分散の組み合わせの各々を用いたときに前記状態推定手段によって各々算出された前記予測誤差に基づいて、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を特定し、最適な正規化システムノイズ共分散行列及び正規化観測ノイズ分散を用いた前記状態空間モデルを、最適モデルの推定結果とする最適モデル推定手段
   として機能させるためのプログラム。