JP2014041547A - 時系列データ解析装置、方法、及びプログラム - Google Patents

時系列データ解析装置、方法、及びプログラム Download PDF

Info

Publication number
JP2014041547A
JP2014041547A JP2012184449A JP2012184449A JP2014041547A JP 2014041547 A JP2014041547 A JP 2014041547A JP 2012184449 A JP2012184449 A JP 2012184449A JP 2012184449 A JP2012184449 A JP 2012184449A JP 2014041547 A JP2014041547 A JP 2014041547A
Authority
JP
Japan
Prior art keywords
state
noise
normalized
covariance matrix
observation
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
JP2012184449A
Other languages
English (en)
Inventor
Shinya Murata
眞哉 村田
Noriko Takaya
典子 高屋
Masashi Uchiyama
匡 内山
Kunio Kashino
邦夫 柏野
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nippon Telegraph and Telephone Corp
Original Assignee
Nippon Telegraph and Telephone Corp
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nippon Telegraph and Telephone Corp filed Critical Nippon Telegraph and Telephone Corp
Priority to JP2012184449A priority Critical patent/JP2014041547A/ja
Publication of JP2014041547A publication Critical patent/JP2014041547A/ja
Pending legal-status Critical Current

Links

Images

Abstract

【課題】非線形モデルを対象として、ノイズの統計量が未知であっても、通常のUFKと同程度の精度で、ロバストな解析結果を得る。
【解決手段】状態空間モデリング部31が、システムノイズを用いた状態更新式、及び観測ノイズを用いた観測方程式によって定められた状態空間モデルを定義し、状態推定部32が、状態空間モデル及び観測値の時系列データに基づいてUFKを行う際に、システムノイズの共分散行列Qの対角成分及び観測ノイズの分散値rの総和tr(Q)+rで正規化された正規化共分散行列^P'、及び総和tr(Q)+rを用いて変更された重みに基づいてシグマ点列を発生させ、^P'、tr(Q)+rで正規化された正規化共分散行列 ̄P’、正規化共分散行列^P’、正規化システムノイズQ’、及び正規化観測ノイズr’を用いて状態更新及び観測更新を行って、状態ベクトルX及び共分散行列Pを逐次推定する。
【選択図】図1

Description

本発明は、時系列データ解析装置、方法、及びプログラムに関する。
時系列データの解析において、データの生成プロセスをシステムモデルと観測モデルとで記述する状態空間モデリングがよく用いられる。例えば、線形モデルの残差(観測値と予測値との差)とその共分散との関係を線形式で書き下し、最小二乗法でシステムノイズと観測ノイズとを推定する手法が提案されている(例えば、非特許文献1参照)。
しかし、実際の問題ではモデルが非線形になることが多い。非線形モデルに対する有効な状態推定アルゴリズムにUnscented Kalman Filter(UKF)が存在する。
BrianJ.Odelsonetal., "A New Autocovariance Least-Squares Method for Estimating Noise Covariances", Technical report, TWMCC, 2003.
しかしながら、従来手法であるUFKでは、システムノイズ及び観測ノイズの分散値の未知性及び初期値任意性により、解析結果がこれらノイズの設定値に大きく依存してしまい、ロバストな解析結果を得られない場合がある、という問題があった。
本発明は、上記の事情を鑑みてなされたもので、非線形モデルを対象として、ノイズの統計量が未知であっても、通常のUFKと同程度の精度で、ロバストな解析結果を得ることができる時系列データ解析装置、方法、及びプログラムを提供することを目的とする。
上記目的を達成するために、本発明の時系列データ解析装置は、システムノイズを用いて状態ベクトルXを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルと観測値との関係を示す観測方程式によって定められた状態空間モデル、並びに前記観測値の時系列データに基づいて、シグマ点列の発生、状態ベクトルの状態更新、及び観測更新を行うUFK(UnscentedKalmanFilter)によって、前記状態ベクトルX及び該状態ベクトルの推定誤差分散の共分散行列Pを逐次推定する際に、前記状態更新後の前記共分散行列Pの予測値 ̄Pをシステムノイズの共分散行列Qの対角成分及び観測ノイズの分散値rの総和で正規化した正規化共分散行列 ̄P’、前記観測更新後の前記共分散行列Pの推定値^Pを前記総和で正規化した正規化共分散行列^P’、前記システムノイズの共分散行列Qを前記総和で正規化した正規化システムノイズQ’、及び前記観測ノイズの分散値rを前記総和で正規化した正規化観測ノイズr’を用いると共に、前記総和を用いて変更された重みに基づいて、前記シグマ点列を発生させる状態推定手段と、前記状態推定手段の推定結果に基づいて、前記正規化システムノイズQ’と前記正規化観測ノイズr’との比率をパラメトライズして得られる尤度が最大となるときのパラメータを、前記状態空間モデルのパラメータとして同定するモデル同定手段と、を含んで構成されている。
本発明の時系列データ解析装置によれば、状態空間モデルとして、システムノイズを用いて状態ベクトルXを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルと観測値との関係を示す観測方程式によって定められたモデルを定義する。そして、状態推定手段が、この状態空間モデルと観測値の時系列データとに基づいて、シグマ点列の発生、状態ベクトルの状態更新、及び観測更新を行うUFK(UnscentedKalmanFilter、無香カルマンフィルタ)によって、状態ベクトルX及び状態ベクトルの推定誤差分散の共分散行列Pを逐次推定する。この際、状態推定手段は、状態更新後の共分散行列Pの予測値 ̄Pをシステムノイズの共分散行列Qの対角成分及び観測ノイズの分散値rの総和で正規化した正規化共分散行列 ̄P’、観測更新後の共分散行列Pの推定値^Pを総和で正規化した正規化共分散行列^P’、システムノイズの共分散行列Qを総和で正規化した正規化システムノイズQ’、及び観測ノイズの分散値rを総和で正規化した正規化観測ノイズr’を用いる。これにより、UFKアルゴリズム内の全ての共分散行列Pの推定値^Pが正規化される。さらに、状態推定手段は、システムノイズの共分散行列Qの対角成分及び観測ノイズの分散値rの総和を用いて変更された重みに基づいて、シグマ点列を発生させる。これにより、推定精度をテイラー展開の意味で2次まで保証し、UKFと同等の精度を得ることができる。
また、モデル同定手段が、状態推定手段の推定結果に基づいて、正規化システムノイズQ’と正規化観測ノイズr’との比率を、例えば0〜1で振るなどしてパラメトライズして得られる尤度が最大となるときのパラメータを、状態空間モデルのパラメータとして同定する。
このように、システムノイズ及び観測ノイズを、システムノイズの共分散行列Qの対角成分及び観測ノイズの分散値rの総和で正規化し、さらに、総和を用いて変更された重みに基づいてシグマ点列を発生させることにより、ノイズの統計量が未知であっても、通常のUFKと同程度の精度で、ロバストな解析結果を得ることができる。
また、前記状態推定手段は、前記総和が収束するまで前記UFKを繰り返し行い、前記総和が収束した際の前記状態ベクトルXの推定値^X及び前記推定値^Pの各々を平滑化した平滑化推定値~X及び~Pを、最終的な推定結果とすることができる。
また、本発明の時系列データ解析方法は、状態推定手段と、モデル同定手段とを含む時系列データ解析装置における時系列データ解析方法であって、前記状態推定手段は、システムノイズを用いて状態ベクトルXを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルと観測値との関係を示す観測方程式によって定められた状態空間モデル、並びに前記観測値の時系列データに基づいて、シグマ点列の発生、状態ベクトルの状態更新、及び観測更新を行うUFK(UnscentedKalmanFilter)によって、前記状態ベクトルX及び該状態ベクトルの推定誤差分散の共分散行列Pを逐次推定する際に、前記状態更新後の前記共分散行列Pの予測値 ̄Pをシステムノイズの共分散行列Qの対角成分及び観測ノイズの分散値rの総和で正規化した正規化共分散行列 ̄P’、前記観測更新後の前記共分散行列Pの推定値^Pを前記総和で正規化した正規化共分散行列^P’、前記システムノイズの共分散行列Qを前記総和で正規化した正規化システムノイズQ’、及び前記観測ノイズの分散値rを前記総和で正規化した正規化観測ノイズr’を用いると共に、前記総和を用いて変更された重みに基づいて、前記シグマ点列を発生させ、前記モデル同定手段は、前記状態推定手段の推定結果に基づいて、前記正規化システムノイズQ’と前記正規化観測ノイズr’との比率をパラメトライズして得られる尤度が最大となるときのパラメータを、前記状態空間モデルのパラメータとして同定する方法である。
また、本発明の時系列データ解析プログラムは、コンピュータを、上記の時系列データ解析装置を構成する各手段として機能させるためのプログラムである。
以上説明したように、本発明の時系列データ解析装置、方法、及びプログラムによれば、システムノイズ及び観測ノイズを、システムノイズの共分散行列Qの対角成分及び観測ノイズの分散値rの総和で正規化し、さらに、総和を用いて変更された重みに基づいてシグマ点列を発生させることにより、ノイズの統計量が未知であっても、通常のUFKと同程度の精度で、ロバストな解析結果を得ることができる、という効果が得られる。
本実施の形態に係る時系列データ解析装置の構成を示す概略図である。 観測値の時系列データのデータ構造の一例を示す図である。 本実施の形態における時系列データ解析処理ルーチンの内容を示すフローチャートである。 実験1における人工時系列データを示すグラフである。 実験1における1期先の予測結果を示すグラフである。 実験1におけるモデルの同定結果を示す表である。 実験2における1期先の予測結果を示すグラフである。 実験2におけるAR係数の推定結果を示すグラフである。 実験2におけるモデルの同定結果を示す表である。
以下、図面を参照して本発明の実施の形態を詳細に説明する。
<本実施の形態の概要>
本実施の形態に係る時系列データ解析装置では、まず、通常のUKFアルゴリズムに存在する状態の共分散行列の状態更新後の予測値、及び観測更新後の推定値の各々を、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和で割る。すなわち、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和を1に正規化する。この操作によりアルゴリズムに存在する全ての共分散行列の推定値は正規化されたものになる。次に、正規化された共分散行列推定値から発生されたシグマ点列を使い時間更新及び観測更新を行って状態の期待値及び共分散行列の推定を行うが、ここで各シグマ点列の重みを適切に変更することで、推定精度をテイラー展開の意味で2次まで保証する。これはUKFと同等の精度である。本実施の形態では説明を簡潔にするため観測モデルは線形を仮定するが、非線形観測モデルに対しても同様なスケーリング処理を行うことでUKFと同等な精度の状態推定を行うことができる。最終的なモデルの同定には対数尤度を使用する。
上記概要に従った本実施の形態における手法を、以下では、Noises Normalized UKF&URTS Smoother(NN−UKF&URTSS)と呼ぶ。
<システム構成>
本実施の形態に係る時系列データ解析装置10は、CPUと、RAMと、後述する時系列データ解析処理ルーチンを実行するためのプログラムを記憶したROMとを備えたコンピュータで構成されている。このコンピュータは、機能的には、図1に示すように、入力部20と、演算部30と、出力部40とを含んだ構成で表すことができる。
入力部20は、入力された観測値の時系列データを受け付ける。観測値は、スカラー値であり、すなわち、時系列データは、1次元時系列である。観測値は、例えば、図2に示すように、時刻tと観測値yとの組(t=1,・・・,N)からなる時系列データである。
また、入力部20は、オペレータから、後述する状態空間モデルを定義するための設定値の入力を受け付ける。
演算部30は、状態空間モデリング部31と、状態推定部32と、モデル同定部33とを含んだ構成で表すことができる。
状態空間モデリング部31は、オペレータにより入力部20から入力された設定値を受け付けて、下記(1)式及び(2)式に示す状態空間モデルを定義する。
Figure 2014041547
ただし、Xは時刻tにおける状態ベクトルで、通常、非観測量が取られる。状態ベクトルXを時系列データyを用いて推定することを目的とする。μはシステムノイズベクトル、wは観測ノイズベクトルであり、共に正規分布に従う白色ノイズベクトルである。また、N(m、v)は平均m、分散vの正規分布を表す。したがって、Qはシステムノイズの共分散行列に相当し、rは観測ノイズの分散値に相当する。上記(2)式に示されるように、yは時刻tにおける観測値で、状態Xの線形変換Hにより生成される。上記(2)式は、観測ノイズを用いて状態ベクトルXと観測値yとの関係を示している。f(・)が非線形箇所であり、これによりシステムモデルは非線形になる。f(・)の型はモデルにより異なるが、簡単な循環変動モデルでは、下記(3)式に示すように状態が定義され、状態更新式が組まれる。
Figure 2014041547
ここで、x,・・・,xt−M+1は自己回帰過程で表現された循環変動成分、a,・・・,aは各循環変動成分にかかる係数である。以下、この係数をAR係数と呼ぶ。この状態の定義を使用したシステムモデルにおける状態更新式のf(Xt−1)及びwは、下記(4)式で表わされる。
Figure 2014041547
ここでwは観測ノイズベクトルであるが、通常w以外の要素は非常に小さな値を設定する(≒0)。この循環変動モデルでは、観測方程式(2)式のHは1行2×M列のベクトルで、1番目の要素が1でその他は0になる。
このように、状態空間モデリング部31は、オペレータからの入力を受け付けて、上記(4)式に示すようなf(・)の型を定義すると共に、上記(2)式のHを定義する。
状態推定部32は、入力された時系列データと、状態空間モデリング部31でモデリングされた状態空間モデルとに基づき、状態ベクトルの逐次推定を行う。モデルは非線形であるので、上記概要で説明したNN−UKF&URTSSによって、状態ベクトルXの逐次状態推定(オンライン状態推定)を行う。NN−UKFのアルゴリズムは、下記(5)式〜(17)式で表されるシグマ点列の発生、状態更新(1期先予測)、及び観測更新(フィルタ)からなる。
Figure 2014041547
ただし、κは予め定められた定数であり、nは状態ベクトルXの次元数である。vは予測誤差であり、Pは状態ベクトルXの推定誤差分散の共分散行列である。また、 ̄A(数式中では、記号A上に“ ̄”)は状態更新により予測されたAの予測値を示し、^A(数式中では、記号A上に“^”)は観測更新により推定されたAの推定値を示す。なお、 ̄P'、^P'、Q’、r’は下記(18)式〜(21)式に示すように、共分散行列の予測値 ̄P、共分散行列の推定値^P、システムノイズの共分散行列Q、及び観測ノイズの分散値rを、tr(Q)+rで割って正規化したものである。
Figure 2014041547
ここで、tr(Q)+rはシステムノイズμの共分散行列Qの対角成分及び観測ノイズwの分散値rの総和であり、アルゴリズム中のαとの関係はtr(Q)+r=1/αである。このスケーリング操作により、アルゴリズムに存在するシステムノイズの共分散行列の対角成分及び観測ノイズの分散値は正規化され、全て1以下になる。また、(6)式及び(9)式は、通常のUKFのシグマ点列の重みであるが、これを(7)式及び(10)式のように変更することで、状態の推定精度をテイラー展開の意味で2次まで保証している。多くの場合、tr(Q)+rの値は未知であるので最初はtr(Q)+r=1、すなわちα=1に設定してNN−UKFを実行する。そして、NN−UKFの実行毎に、下記(22)式及び(23)式に基づいて、tr(Q)+rを推定する。
Figure 2014041547
状態推定部32は、NN−UKF実行時に使用した(tr(Q)+r)nowの値と、NN−UKF実行終了時に推定した(tr(Q)+r)newの値との差が所定の閾値以下であれば、tr(Q)+rの値が収束したとみなし、この実行によって推定された状態の結果を採用する手法を取る。最後の観測値yt=Nまでの処理が完了すると、(13)式を使用したNoises Normalized Unscented RTS Smoother(NN−URTSS)により、下記(24)式〜(26)式に示す平滑化を行う。
Figure 2014041547
ここで、St−1は平滑化ゲイン、~Xt−1及び~Pt−1(数式中では、記号XまたはP上に“~”)はそれぞれ状態の期待値及び共分散行列の平滑化推定値であり、この推定値が最終的にアウトプットされる。
モデル同定部33は、状態推定部32により正規化されたシステムノイズの共分散行列Q’及び観測ノイズの分散値r’を利用して、Q’とr’との比率γをパラメトライズして振って、γ毎に下記(27)式〜(30)式に示す対数尤度LL(γ)を求める。そして、対数尤度LL(γ)が最大になったγから得られるQ’及びr’を、状態推定モデルのパラメータとして推定することにより、モデルの同定を行う。
Figure 2014041547
ここで、M(γ)がシステムノイズμの共分散行列の対角成分及び観測ノイズwの分散値の総和の推定値となる。
出力部40は、状態推定部32により推定された状態ベクトルの平滑化推定値~Xt−1及びその共分散行列の平滑化推定値~Pt−1、並びにモデル同定部33により同定されたモデルを示すパラメータの推定値を出力する。
<時系列データ解析装置の作用>
次に、本実施の形態に係る時系列データ解析装置10の作用について説明する。まず、オペレータにより、状態空間モデルを定義するための設定値が時系列データ解析装置10に入力されると、時系列データ解析装置10によって、入力された設定値が、メモリ(図示省略)へ格納される。また、時刻t〜tの観測値からなる時系列データが、時系列データ解析装置10に入力されると、時系列データ解析装置10によって、入力された時系列データが、メモリへ格納される。そして、時系列データ解析装置10によって、図3に示す時系列データ解析処理ルーチンが実行される。
まず、ステップ100で、メモリに格納された状態空間モデルを定義するための設定値及び観測値の時系列データを取得する。そして、ステップ102で、状態空間モデリング部31が、取得した設定値を用いて、上記(1)式及び(2)式に示す状態空間モデルを定義する。
次に、ステップ104で、状態推定部32が、NN−UKF&URTSSを実行するための初期設定として、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和(tr(Q)+r)、並びに時系列データの時刻を示すtをそれぞれ1にセットする。
次に、ステップ106で、状態推定部32が、NN−UKFのアルゴリズムに従って、上記ステップ100で取得した観測値の時系列データ、及び上記ステップ104で設定された初期設定を用いて、シグマ点列の発生、状態更新、及び観測更新を逐次行って、状態ベクトル及びその共分散行列を推定する。
次に、ステップ108で、状態推定部32が、上記(22)式及び(23)式によりtr(Q)+rを推定し、上記ステップ106で使用した(tr(Q)+r)nowの値と、本ステップで推定した(tr(Q)+r)newの値との差が所定の閾値以下か否かを判定することにより、tr(Q)+rの値が収束したか否かを判定する。tr(Q)+rの値が収束していない場合には、ステップ106へ戻り、推定された状態ベクトル^X、推定され正規化された共分散行列^P'、及び推定されたtr(Q)+rを用いて、NN−UKFによる状態ベクトル及びその共分散行列の推定を繰り返す。
一方、tr(Q)+rの値が収束した場合には、ステップ110へ移行し、tがNになったか否かを判定することにより、最後の観測値yまで処理が完了したか否かを判定する。tがNに達していない場合には、ステップ112へ移行して、tを1インクリメントして、ステップ106へ戻る。t=Nとなった場合には、ステップ114へ移行する。
ステップ114では、状態推定部32が、上記(13)式を使用したNN−URTSSにより、上記(24)式〜(26)式に示す平滑化を行い、状態ベクトルの平滑化推定値~Xt−1及びその共分散行列の平滑化推定値~Pt−1を得る。
次に、ステップ116で、モデル同定部33が、Q’とr’との比率γをパラメトライズして振って、γ毎に上記(27)式〜(30)式に示す対数尤度LL(γ)を求め、対数尤度LL(γ)が最大になったγから得られるQ’及びr’を、状態推定モデルのパラメータとして推定することにより、モデルの同定を行う。
次に、ステップ118で、出力部40が、上記ステップ114で得られた状態ベクトルの平滑化推定値~Xt−1及びその共分散行列の平滑化推定値~Pt−1、並びに上記ステップ106で同定されたモデルを示すパラメータの推定値を出力し、時系列データ解析処理ルーチンを終了する。
<実験結果>
ここで、下記(31)式(3次のAR過程、AR係数は1.5,-0.9,0.3、観測ノイズの分散値は1)で生成された人工時系列データを使用したNN−UKF&URTSSの実験結果について説明する。
Figure 2014041547
まず、実験1ではAR係数を既知として循環変動成分のみを推定した。図4に上記(31)式に基づいて生成された人工時系列データ(AR(3))、図5に実験1による1期先の予測結果、図6にシステムノイズと観測ノイズの値をλ(=システムノイズ/(システムノイズ+観測ノイズ))とした各モデルの同定結果を示す。図4及び図5から1期先の予測が精度良く行われていることがわかり、図6から観測ノイズの分散値(図中、ノイズ推定値)の真値1を推定することに成功していることがわかる。
次に、実験2ではAR係数も未知とし、非線形システムモデル((1)式〜(4)式)を使用して状態推定を行った。図7に実験2による1期先の予測結果、図8に3つのAR係数の逐次推定(オンライン推定)結果を示す。図8からAR係数のそれぞれの真値(1.5,-0.9,0.3)を適切に推定していることがわかる。また、図9はモデルの同定結果であり、本実施の形態の手法であるNN−UKF&URTSSの有効性を示している。
以上説明したように、本実施の形態に係る時系列データ解析装置によれば、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和を1で正規化して、UKFアルゴリズムに存在する状態ベクトルの誤差分散の共分散行列の状態更新後の予測値、及び観測更新後の推定値の各々を正規化し、正規化された共分散行列の推定値から発生されたシグマ点列を使って時間更新及び観測更新を行うことにより、ノイズの統計量が未知であっても、ロバストな推定結果を提供することが可能になる。また、システムノイズの共分散行列の対角成分及び観測ノイズの分散値の総和に基づいて適切に変更された重みに基づいてシグマ点列を発生させることにより、推定精度をテイラー展開の意味で2次まで保証することができ、従来のUKFと同等の精度を得ることができる。
なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
例えば、状態空間モデルの定義を予め設定しておき、状態空間モデルを定義するための設定値の入力を不要としてもよい。
また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能である。
10 時系列データ解析装置
20 入力部
30 演算部
31 状態空間モデリング部
32 状態推定部
33 モデル同定部
40 出力部

Claims (4)

  1. システムノイズを用いて状態ベクトルXを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルと観測値との関係を示す観測方程式によって定められた状態空間モデル、並びに前記観測値の時系列データに基づいて、シグマ点列の発生、状態ベクトルの状態更新、及び観測更新を行うUFK(Unscented Kalman Filter)によって、前記状態ベクトルX及び該状態ベクトルの推定誤差分散の共分散行列Pを逐次推定する際に、前記状態更新後の前記共分散行列Pの予測値 ̄Pをシステムノイズの共分散行列Qの対角成分及び観測ノイズの分散値rの総和で正規化した正規化共分散行列 ̄P’、前記観測更新後の前記共分散行列Pの推定値^Pを前記総和で正規化した正規化共分散行列^P’、前記システムノイズの共分散行列Qを前記総和で正規化した正規化システムノイズQ’、及び前記観測ノイズの分散値rを前記総和で正規化した正規化観測ノイズr’を用いると共に、前記総和を用いて変更された重みに基づいて、前記シグマ点列を発生させる状態推定手段と、
    前記状態推定手段の推定結果に基づいて、前記正規化システムノイズQ’と前記正規化観測ノイズr’との比率をパラメトライズして得られる尤度が最大となるときのパラメータを、前記状態空間モデルのパラメータとして同定するモデル同定手段と、
    を含む時系列データ解析装置。
  2. 前記状態推定手段は、前記総和が収束するまで前記UFKを繰り返し行い、前記総和が収束した際の前記状態ベクトルXの推定値^X及び前記推定値^Pの各々を平滑化した平滑化推定値~X及び~Pを、最終的な推定結果とする請求項1記載の時系列データ解析装置。
  3. 状態推定手段と、モデル同定手段とを含む時系列データ解析装置における時系列データ解析方法であって、
    前記状態推定手段が、システムノイズを用いて状態ベクトルXを非線形に時間更新するための状態更新式、及び観測ノイズを用いて状態ベクトルと観測値との関係を示す観測方程式によって定められた状態空間モデル、並びに前記観測値の時系列データに基づいて、シグマ点列の発生、状態ベクトルの状態更新、及び観測更新を行うUFK(Unscented Kalman Filter)によって、前記状態ベクトルX及び該状態ベクトルの推定誤差分散の共分散行列Pを逐次推定する際に、前記状態更新後の前記共分散行列Pの予測値 ̄Pをシステムノイズの共分散行列Qの対角成分及び観測ノイズの分散値rの総和で正規化した正規化共分散行列 ̄P’、前記観測更新後の前記共分散行列Pの推定値^Pを前記総和で正規化した正規化共分散行列^P’、前記システムノイズの共分散行列Qを前記総和で正規化した正規化システムノイズQ’、及び前記観測ノイズの分散値rを前記総和で正規化した正規化観測ノイズr’を用いると共に、前記総和を用いて変更された重みに基づいて、前記シグマ点列を発生させ、
    前記モデル同定手段が、前記状態推定手段の推定結果に基づいて、前記正規化システムノイズQ’と前記正規化観測ノイズr’との比率をパラメトライズして得られる尤度が最大となるときのパラメータを、前記状態空間モデルのパラメータとして同定する
    時系列データ解析方法。
  4. コンピュータを、請求項1または請求項2記載の時系列データ解析装置を構成する各手段として機能させるための時系列データ解析プログラム。
JP2012184449A 2012-08-23 2012-08-23 時系列データ解析装置、方法、及びプログラム Pending JP2014041547A (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP2012184449A JP2014041547A (ja) 2012-08-23 2012-08-23 時系列データ解析装置、方法、及びプログラム

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2012184449A JP2014041547A (ja) 2012-08-23 2012-08-23 時系列データ解析装置、方法、及びプログラム

Publications (1)

Publication Number Publication Date
JP2014041547A true JP2014041547A (ja) 2014-03-06

Family

ID=50393743

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2012184449A Pending JP2014041547A (ja) 2012-08-23 2012-08-23 時系列データ解析装置、方法、及びプログラム

Country Status (1)

Country Link
JP (1) JP2014041547A (ja)

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106610980A (zh) * 2015-10-22 2017-05-03 日本电气株式会社 用于对时空序列数据进行分类/预测的设备和方法
CN108255786A (zh) * 2017-11-28 2018-07-06 中南大学 一种称重结果的干扰补偿计算方法及系统
JP2018116441A (ja) * 2017-01-17 2018-07-26 富士通株式会社 処理装置、調整パラメータの予測モデル推定方法、及び調整パラメータの予測モデル推定プログラム
CN112162563A (zh) * 2020-09-15 2021-01-01 郑州轻工业大学 基于自适应弱敏无迹Kalman滤波的直升机状态估计方法
CN112468116A (zh) * 2020-12-01 2021-03-09 哈尔滨工程大学 一种基于Gibbs采样器的自适应平滑方法
CN117879540A (zh) * 2024-03-12 2024-04-12 西南应用磁学研究所(中国电子科技集团公司第九研究所) 基于改进卡尔曼滤波的磁罗盘传感器自适应信号滤波方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2003161778A (ja) * 2001-11-26 2003-06-06 Mitsubishi Electric Corp 目標追尾装置及び方法
JP2004070452A (ja) * 2002-08-02 2004-03-04 Nippon Hoso Kyokai <Nhk> 状態推定装置、その方法及びそのプログラム、並びに、現在状態推定装置及び未来状態推定装置
US20060251304A1 (en) * 2005-03-21 2006-11-09 Charles Florin System and method for kalman filtering in vascular segmentation

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2003161778A (ja) * 2001-11-26 2003-06-06 Mitsubishi Electric Corp 目標追尾装置及び方法
JP2004070452A (ja) * 2002-08-02 2004-03-04 Nippon Hoso Kyokai <Nhk> 状態推定装置、その方法及びそのプログラム、並びに、現在状態推定装置及び未来状態推定装置
US20060251304A1 (en) * 2005-03-21 2006-11-09 Charles Florin System and method for kalman filtering in vascular segmentation

Cited By (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106610980A (zh) * 2015-10-22 2017-05-03 日本电气株式会社 用于对时空序列数据进行分类/预测的设备和方法
CN106610980B (zh) * 2015-10-22 2022-03-18 日本电气株式会社 用于对时空序列数据进行分类/预测的设备和方法
JP2018116441A (ja) * 2017-01-17 2018-07-26 富士通株式会社 処理装置、調整パラメータの予測モデル推定方法、及び調整パラメータの予測モデル推定プログラム
CN108255786A (zh) * 2017-11-28 2018-07-06 中南大学 一种称重结果的干扰补偿计算方法及系统
CN108255786B (zh) * 2017-11-28 2021-07-16 中南大学 一种称重结果的干扰补偿计算方法及系统
CN112162563A (zh) * 2020-09-15 2021-01-01 郑州轻工业大学 基于自适应弱敏无迹Kalman滤波的直升机状态估计方法
CN112162563B (zh) * 2020-09-15 2023-01-31 郑州轻工业大学 基于自适应弱敏无迹Kalman滤波的直升机状态估计方法
CN112468116A (zh) * 2020-12-01 2021-03-09 哈尔滨工程大学 一种基于Gibbs采样器的自适应平滑方法
CN112468116B (zh) * 2020-12-01 2023-06-16 哈尔滨工程大学 一种基于Gibbs采样器的自适应平滑方法
CN117879540A (zh) * 2024-03-12 2024-04-12 西南应用磁学研究所(中国电子科技集团公司第九研究所) 基于改进卡尔曼滤波的磁罗盘传感器自适应信号滤波方法

Similar Documents

Publication Publication Date Title
Särkkä et al. Gaussian filtering and smoothing for continuous-discrete dynamic systems
JP2014041547A (ja) 時系列データ解析装置、方法、及びプログラム
JP6768155B2 (ja) 状態推定装置
US10620597B2 (en) Gray box model estimation for process controller
Dyvak et al. Improving the computational implementation of the parametric identification method for interval discrete dynamic models
Piga et al. LPV system identification under noise corrupted scheduling and output signal observations
JP2014041565A (ja) 時系列データ解析装置、方法、及びプログラム
JP5261740B2 (ja) コンクリート構造物の中性化深さ予測装置および中性化深さをコンピュータに計算させるためのプログラム
Qi et al. Estimation of distribution function for control valve stiction estimation
Bansal et al. A new stochastic simulation algorithm for updating robust reliability of linear structural dynamic systems subjected to future Gaussian excitations
Figueroa et al. An approach for identification of uncertain Wiener systems
JP5738778B2 (ja) 最適モデル推定装置、方法、及びプログラム
JP2013061768A (ja) 最適モデル推定装置、方法、及びプログラム
Schöbi et al. PC-Kriging: A new meta-modelling method and its application to quantile estimation
Evangelidis et al. Quantitative verification of Kalman filters
Larsson et al. Estimation of continuous-time stochastic system parameters
JP2016018323A (ja) パラメータ推定方法、装置、及びプログラム
Lindström Estimating parameters in diffusion processes using an approximate maximum likelihood approach
Pandurangan et al. The use of polynomial chaos for parameter identification from measurements in nonlinear dynamical systems
JP5710513B2 (ja) 信号抽出装置、方法、及びプログラム
Lwin Parameter estimation in first-order autoregressive model for statistical process monitoring in the presence of data autocorrelation
JP6954346B2 (ja) パラメータ推定装置、パラメータ推定方法、及びプログラム
Ardekani et al. Bayesian damage identification of simply supported beams from elastostatic data
Bell Yes, the CAPM is absurd: OLS is misunderstood and incorrectly modeled mathematically
Renard et al. Bayesian Kriging with lognormal data and uncertain variogram parameters

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20140815

A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20150515

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20150609

A02 Decision of refusal

Free format text: JAPANESE INTERMEDIATE CODE: A02

Effective date: 20151201