(57)【要約】【課題】固有値に無限大値が含まれる共分散行列をコンピュータ上で演算する場合であっても、演算量を増大させることなく計算の精度を向上させることのできる技術を提供する。【解決手段】共分散行列を行列Ａと対角行列Ｄと行列Ａの転置行列の積の形で表して演算を行う処理装置１において、入力部８は、対角行列Ｄを、無限大値と対角行列Ｄ１との積と、対角行列Ｄ２との和で表したとき、対角行列Ｄ１及び対角行列Ｄ２を入力させる。第１の格納部４及び第２の格納部５はそれぞれ対角行列Ｄ１及び対角行列Ｄ２を格納する。第３の格納部６は行列Ａを格納する。演算部７は、対角行列Ｄ１、対角行列Ｄ２及び行列Ａを用いて、共分散行列を下三角行列Ｌと対角行列Ｅと下三角行列Ｌの転置行列の積の形で表し、対角行列Ｅを無限大値と対角行列Ｅ１の積と、対角行列Ｅ２との和で表したときの下三角行列Ｌ、対角行列Ｅ１及び対角行列Ｅ２を計算する。出力部９は得られた下三角行列Ｌを出力する。【選択図】図２
2008176683.tif 000002
【特許請求の範囲】
【請求項１】
共分散行列を行列Ａと対角行列Ｄと行列Ａの転置行列の積の形で表して演算を行う処理装置であって、
対角行列Ｄを、無限大値と対角行列Ｄ１との積と、対角行列Ｄ２との和で表したとき、対角行列Ｄ１及び対角行列Ｄ２を入力させる入力手段と、
対角行列Ｄ１及び対角行列Ｄ２を格納する第１の格納手段と、
行列Ａを格納する第２の格納手段と、
前記第１の格納手段及び第２の格納手段に格納された対角行列Ｄ１、対角行列Ｄ２及び行列Ａを用いて、前記共分散行列を下三角行列Ｌと対角行列Ｅと下三角行列Ｌの転置行列の積の形で表し、対角行列Ｅを無限大値と対角行列Ｅ１の積と、対角行列Ｅ２との和で表したときの下三角行列Ｌ、対角行列Ｅ１及び対角行列Ｅ２を計算する演算手段と、
得られた下三角行列Ｌを出力する出力手段と
を備えることを特徴とする処理装置。
【請求項２】
前記第１の格納手段は、対角行列Ｄ１を格納する手段と対角行列Ｄ２を格納する手段とを含む
を備えることを特徴とする請求項１記載の処理装置。
【請求項３】
前記第１の格納手段は、対角行列Ｄ１及び対角行列Ｄ２の対角要素を格納する要素格納手段と、該要素格納手段に格納されている対角要素がそれぞれ対角行列Ｄ１あるいは対角行列Ｄ２のいずれを表すかを示す情報を格納する格納手段とを含む
ことを特徴とする請求項１記載の処理装置。
【請求項４】
共分散行列を行列Ａと対角行列Ｄと行列Ａの転置行列の積の形で表して演算を行うための共分散行列演算プログラムであって、
対角行列Ｄを、無限大値と対角行列Ｄ１との積と、対角行列Ｄ２との和で表したとき、対角行列Ｄ１及び対角行列Ｄ２を入力させ、
対角行列Ｄ１、対角行列Ｄ２及び行列Ａを格納手段に格納し、
前記格納手段に格納された対角行列Ｄ１、対角行列Ｄ２及び行列Ａを用いて、前記共分散行列を下三角行列Ｌと対角行列Ｅと下三角行列Ｌの転置行列の積の形で表し、対角行列Ｅを無限大値と対角行列Ｅ１の積と、対角行列Ｅ２との和で表したときの下三角行列Ｌ、対角行列Ｅ１及び対角行列Ｅ２を計算し、
得られた下三角行列Ｌを出力する
処理をコンピュータに実行させることを特徴とする共分散行列演算プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、時系列データに基づいて状態変数値を推定する場合に用いられる共分散行列を演算する技術に関する。
【背景技術】
【０００２】
従来から、時系列データに基づいて現在の状態変数値を推定する問題に対しては、カルマンフィルタが広く用いられている。カルマンフィルタとは、解析時までの時系列データを用いて時刻更新アルゴリズムと観測更新アルゴリズムを繰り返して行い、最適推定を行う手法である。カルマンフィルタアルゴリズムの実行に際しては、状態変数値推定の共分散行列の演算が行われる（例えば、非特許文献１）。
【０００３】
ここで、例えば、カルマンフィルタを用いてランダムに位置移動する運動体の位置座標を推定するときに、上記の時間更新において、推定対象とされている位置移動する運動体の位置に対する状態遷移推定の誤差が大きい場合がある。このような場合には、カルマンフィルタの共分散行列の演算においては、入力値として与える状態遷移雑音の共分散行列の固有値の１つあるいは複数が、特異的すなわち無限大発散することとして演算を行う必要が生じる。この他にも、初期時刻における運動体の位置が確定していない等の場合においても、本来は、共分散行列の固有値に特異値が含まれることとして行列演算を行う必要がある。
【０００４】
固有値に特異値が含まれる場合は、コンピュータでは固有値に無限大値を含む行列の演算ができないため、カルマンフィルタに入力する共分散行列をコンピュータで演算するための各種の手法が従来から提案されている。
【０００５】
かかる手法の一つに、無限大値を用いる代わりに、数値的に安定が保たれる程度、具体的には固有値のうち最大値と最小値の比が１００対１程度に収まる範囲で十分に大きい値を代入して演算する方法がある。
【０００６】
他の手法としては、２つの行列（Ｓ及びＴとする）を用いて、共分散行列をＳとＴの逆行列Ｔ−１との積に置き換えて演算する方法がある（例えば、非特許文献２）。
更に他の手法としては、特異性を持つ共分散行列に対しムーア・ペンローズ（Moore-Penrose）型の一般逆行列を用いて演算する方法がある（例えば、特許文献１）。
【０００７】
なお、公知技術として、カルマンフィルタの計算を行列に対する直交・三角行列化の計算に帰着させる方法が用いられており、特に、特異性を持たない対称行列（Σとする）に対して修正コレスキー分解を行うことにより、Σ＝ＡＤＡＴを満たす行列Ａ及び対角行列Ｄを求める方法（例えば、非特許文献３）や、グラム・シュミット法により特異性を持たない行列に対して直交・三角行列化を行う方法（例えば、非特許文献４）が用いられている。
【特許文献１】特許第２８４６９５８号公報
【非特許文献１】片山徹著、「応用カルマンフィルタ」、朝倉書店、２０００年１月、ｐ．１６４−１７２
【非特許文献２】C. C. Paige、「Covariance Matrix Representation in Linear Filtering」、Linear Algebra and Its Role in Systems Theory、B. N. Datta、Ed. Providence、American Mathematical Society、１９８５、ｐｐ．３０９−３２１
【非特許文献３】伊理正夫著、「岩波講座  応用数学１［基礎１］線形代数Ｉ」、岩波書店、１９９７年９月、ｐ．６９−７２
【非特許文献４】森正武、杉原正顕、室田一雄著、「岩波講座  応用数学１２［方法２］線形計算」、岩波書店、１９９７年１１月、ｐ．６８−７１
【発明の開示】
【発明が解決しようとする課題】
【０００８】
共分散行列の行列要素に特異値が含まれる場合における上記の各種手法のうち、無限大値に替えて十分に大きい値を代入する方法によれば、計算処理上の数値的不安定の発生や算出精度が低下するという問題がある。また、共分散行列を上記のとおりＳＴ−１に置き換えて演算する方法によれば、算出精度の低下を防ぐことができる一方で、計算量が増大し、演算時間が長期化するという問題がある。また、ムーア・ペンローズ型の一般逆行列を使用する方法に関しては、一般逆行列は、その元の行列の固有値を有限値から無限大へ極限をとったものとは必ずしも一致せず、このため、数学的に正しい方法とは言えない。
【０００９】
本発明は、固有値に無限大値が含まれる共分散行列をコンピュータ上で演算する場合であっても、演算量を増大させることなく計算の精度を向上させることのできる技術を提供することを目的とする。
【課題を解決するための手段】
【００１０】
上記課題を解決するために、本発明に係る処理装置は、共分散行列を行列Ａと対角行列Ｄと行列Ａの転置行列の積の形で表して演算を行う処理装置であって、対角行列Ｄを、無限大値と対角行列Ｄ１との積と、対角行列Ｄ２との和で表したとき、対角行列Ｄ１及び対角行列Ｄ２を入力させる入力手段と、対角行列Ｄ１及び対角行列Ｄ２を格納する第１の格納手段と、行列Ａを格納する第２の格納手段と、前記第１の格納手段及び第２の格納手段に格納された対角行列Ｄ１、対角行列Ｄ２及び行列Ａを用いて、前記共分散行列を下三角行列Ｌと対角行列Ｅと下三角行列Ｌの転置行列の積の形で表し、対角行列Ｅを無限大値と対角行列Ｅ１の積と、対角行列Ｅ２との和で表したときの下三角行列Ｌ、対角行列Ｅ１及び対角行列Ｅ２を計算する演算手段と、得られた下三角行列Ｌを出力する出力手段とを備えた構成とする。
【００１１】
上記の処理装置においては、共分散行列を行列Ａと対角行列Ｄと行列Ａの転置行列の積（ＡＤＡＴ）で表し、行列式中の対角行列Ｄに関し、Ｄ＝∞×Ｄ１＋Ｄ２の形で表す。共分散行列に含まれる無限大値をとる固有値は、上記対角行列Ｄ１の外に出されるため、ＡＤＡＴ＝ＬＥＬＴを満たす下三角行列Ｌ及び対角行列Ｅを求める計算を、演算量を増大させることなくコンピュータを用いて高精度で求めることが可能となる。
【００１２】
前記第１の格納手段は、対角行列Ｄ１を格納する手段と対角行列Ｄ２を格納する手段とを含む構成としてもよい。あるいは、前記第１の格納手段は、対角行列Ｄ１及び対角行列Ｄ２の対角要素を格納する要素格納手段と、該要素格納手段に格納されている対角要素がそれぞれ対角行列Ｄ１あるいは対角行列Ｄ２のいずれを表すかを示す情報を格納する格納手段とを含む構成としてもよい。
【００１３】
本発明は、上記の処理装置に限らない。例えば、上記の処理方法やその方法を実行するプログラム等であっても、本発明に含まれる。
【発明の効果】
【００１４】
本発明によれば、共分散行列の固有値に無限大値が含まれる場合であっても、コンピュータにおいて、演算量を増大させることなく高精度の演算を実行することが可能となる。
【発明を実施するための最良の形態】
【００１５】
以下、本発明の好適な実施の形態について、図面を参照しながら詳細に説明する。
図１は、カルマンフィルタの更新処理についての概念図である。まず、図１を参照して、カルマンフィルタアルゴリズムに用いる各種の行列について説明する。
【００１６】
図１に示す行列のうち、行列Ｆｉは時刻ｉにおける状態遷移行列、Ｇｉは時刻ｉにおける駆動行列、Ｈｉは時刻ｉの観測行列である。そして、行列Ｐｉは時刻ｉにおける推定誤差共分散行列であり、Ｑｉ及びＲｉは、それぞれ時刻ｉにおける非負定値対称行列及び正定値対称行列である。各行列の次数は、Ｒｉがｍ１行×ｍ１列、Ｈｉがｍ１行×ｍ２列、Ｐｉがｍ２行×ｍ２列、Ｆｉがｍ２行×ｍ２列、Ｇｉがｍ２行×ｍ３列、Ｑｉがｍ３行×ｍ３列とし、以降の説明においては、ｍ＝ｍ１＋ｍ２、ｎ＝ｍ１＋ｍ２＋ｍ３とする。
【００１７】
カルマンフィルタを適用するシステムに推定誤差共分散行列の初期値
【数１】
2008176683.tif 000003
【００１８】
及び時系列データが入力されると、以後はカルマンフィルタアルゴリズムにより、時刻ｉにおける推定誤差共分散行列Ｐｉと与えられた各行列（Ｆｉ、Ｇｉ、Ｈｉ、Ｑｉ及びＲｉ）とから観測更新及び時間更新の計算を行うことにより、時刻（ｉ＋１）における予測誤差共分散行列Ｐｉ＋１が得られる。Ｐｉ＋１により状態推定値が求められ、更にＰｉ＋１と与えられる時刻（ｉ＋１）の行列群に対して観測更新及び時間更新により、時刻（ｉ＋２）における予測誤差共分散行列Ｐｉ＋２を求めることができる。以後同様に、更新を繰り返して時刻ｔにおける予測誤差共分散行列Ｐｔを求め、状態推定値を求めていく。
【００１９】
先に従来の共分散行列の演算処理における問題点について説明したとおり、初期値として入力される共分散行列Ｐ０や時系列データ（Ｆｉ、Ｇｉ、Ｈｉ、Ｑｉ及びＲｉからなる行列群）の値によっては、数値的に安定した演算結果を得ることができない。そこで、数値的に安定した結果を得るために、公知の技術として平方根フィルタが用いられている。平方根フィルタアルゴリズムによれば、上記の行列のうち、行列Ｑｉ、Ｒｉ及びＰｉについてはそれぞれ平方根形を用いて計算を行う。以下、平方根フィルタに入力するための共分散行列の演算方法について説明する。
【００２０】
図２は、本実施形態に係る共分散行列の演算を行う情報処理装置のブロック図である。情報処理装置１は、入力部８、第１の格納部４、第２の格納部５、第３の格納部６、演算部７及び出力部９を含んで構成される。
【００２１】
入力部８は、情報処理装置１に入力される時系列データ及び情報処理装置における共分散行列の演算結果を入力させる。具体的には、共分散行列に修正コレスキー分解を行うことにより得られる形、すなわち行列Ａと対角行列Ｄと行列Ａの転置行列ＡＴとの積の形（ＡＤＡＴ）で表した場合の、行列Ａ及び対角行列Ｄを入力させる。行列Ａは、上記の時系列で与えられる行列群及びＰｉを要素に含む。
【００２２】
ここで、対角行列Ｄの対角要素のうち、無限大値からなる要素については無限大値と第１の対角行列Ｄ（１）積で表し、他の要素については第２の対角行列Ｄ（２）の対角要素とし、対角行列Ｄは第１の対角行列Ｄ（１）と第２の対角行列Ｄ（２）との和で表されることとする。すなわち、第１及び第２の対角行列Ｄ（１）、Ｄ（２）は以下の（１）式を満たす。
Ｄ＝∞×Ｄ（１）＋Ｄ（２）      （１）
第１の格納部４及び第２の格納部５は、置換処理により得られた行列Ｄ（１）、Ｄ（２）の対角要素を格納する。第３の格納部６は、上記行列Ａの行列要素を格納する。
【００２３】
演算部７は、それぞれの格納部に格納されている値を用いて、カルマンフィルタの計算を行い、ＡＤＡＴ＝ＬＥＬＴを満たす下三角行列Ｌ及び対角行列Ｅを求める。詳細については後述する。出力部９は、演算結果を外部に出力するとともに、演算部７の演算の結果得られる共分散行列Ｐｉ＋１を出力する。カルマンフィルタにおいては、次のステップの計算に用いるため、得られたＰｉ＋１は入力部８に与えられる。
【００２４】
図３は、第１の格納部４及び第２の格納部５の構造の例を示す図である。
図３（ａ）は、第１及び第２の格納部のそれぞれに、第１及び第２の対角行列の対角要素の値を保持する場合のメモリの状態を示す。第１及び第２の格納部４、５は、上記の例においてはそれぞれ行列Ｄの次数と等しい個数の領域を備えており、それぞれに第１及び第２の対角行列の対角要素を格納する。
【００２５】
図３（ｂ）は、第１の格納部４には、対角要素の値を格納し、第２の格納部５には、第１の格納部４に格納されている要素は第１あるいは第２の対角行列のいずれの要素であるかを識別するための値を格納する場合のメモリの状態を示す。図３（ｂ）に示す例では、第１の格納部４に第１の対角行列の要素が格納される場合は、第２の格納部５の対応する領域に例えば「１」を格納し、第１の格納部４に第２の対角行列の要素が格納される場合は、第２の格納部５の対応する領域に例えば「０」を格納することとしている。例えば、行列Ｄの１行には無限大値が含まれるため、第１の格納部４にはＤ１（１）の値が格納されており、第２の格納部５には、第１の格納部４に保持されている値はＤ（１）の値であることを示す情報である「１」が格納されている。行列Ｄの２行には無限大値は含まれないため、第１の格納部４にはＤ２（２）の値が格納されており、第２の格納部５には、第１の格納部４に保持されている値はＤ２（２）の値であることを示す情報「０」が格納されている。
【００２６】
なお、図３（ｂ）に示すような対角要素の保持方法は、２つの対角行列Ｄ（１）、Ｄ（２）の要素ｄｉ（１）及びｄｉ（２）の内積は０であることから実施可能である。
対角行列Ｄの要素のうち、無限大値をとる対角要素については上記第１の対角行列Ｄ（１）に含め、残りの対角要素については第２の対角行列Ｄ（２）に含めることで、共分散行列に無限大値が含まれる場合であっても、コンピュータで行列演算をすることが可能となる。また、従来の方法、例えば無限大値の代わりに数値的安定性を保つことのできる程度に大きい値を代入して以後の計算を行う方法や、共分散行列を行列ＳとＴの逆行列であるＴ−１との積ＳＴ−１に置き換えて計算する方法と比較して、演算量を増大させることなく精度の高い値を得ることができる。
【００２７】
ところで、上記の修正コレスキー分解により得られる行列については、ｍ行×ｎ列の行列Ａは
【数２】
2008176683.tif 000004
（ｍ＝ｍ１＋ｍ２、ｎ＝ｍ１＋ｍ２＋ｍ３）と表される。（２）式で表されるｍ行×ｎ列行列Ａと正規直交行列Θとにより、ｍ行×ｎ列の行列ＡΘは、以下のように表現される。
【数３】
2008176683.tif 000005
【００２８】
なお、（３）式中の
【数４】
2008176683.tif 000006
はｍ１行×ｍ１列行列であり、
【数５】
2008176683.tif 000007
はｍ２行×ｍ１列行列である。
【００２９】
上記の（３）式に示すｍ行×ｎ列の行列ＡΘのうち、共分散行列を含む部分を取り出したｍ行×ｍ列の行列を行列Ｌとおくと、行列Ｌは、下の（４）式のように表される。
【数６】
2008176683.tif 000008
【００３０】
上記の修正コレスキー分解により得られたＡＤＡＴに対し、直交・三角行列化の計算を行ってＡＤＡＴ＝ＬＥＬＴ（Ｅはｍ次対角行列）を満たすｍ次下三角行列Ｌを求めることで、カルマンフィルタの誤差推定共分散行列
【数７】
2008176683.tif 000009
【００３１】
が得られる。直交・三角行列化の計算には、例えばグラム・シュミット法を用いる。
【００３２】
以下、図面を参照して、直交・三角行列化の具体的な方法について説明する。
図４は、グラム・シュミット法による直交・三角行列化のアルゴリズムを示す図であり、図５から図７は、その方法を示したフローチャートである。
【００３３】
入力される行列は、上記のｍ行×ｎ列行列Ａ及びｎ次対角行列Ｄである。対角行列Ｄについては、上記（１）式に示すとおり、第１の対角行列Ｄ（１）と第２の対角行列Ｄ（２）との和に置換されている。
【００３４】
行列Ａ及び対角行列Ｄの入力に対して出力される行列は、ｍ次下三角行列Ｌ、ｍ行×ｎ列重み付き直交行列Ｂ及びｍ次対角行列Ｅである。行列Ｄが２つの対角行列の和に置換されているので、行列Ｅについても、無限大値を含む対角要素については第１の対角行列Ｅ（１）に、無限大値以外の対角要素については第２の対角行列Ｅ（２）に含まれ、Ｅ＝∞×Ｅ（１）＋Ｅ（２）と表される。
【００３５】
図４から図７において用いている変数のうち、「ｉ」及び「ｊ」は、演算すべき行列要素を識別するために用いるための変数であり、「ｋ」は、第１及び第２の対角行列Ｅ（１）、Ｅ（２）（あるいはＤ（１）、Ｄ（２））を識別するための変数である。
【００３６】
行列Ｂは、
Ａ＝ＬＢ          （５）
で定義される。このとき、ＡＤＡＴ＝（ＬＢ）Ｄ（ＬＢ）Ｔ＝ＬＢＤＢＴＬＴより、
ＢＤＢＴ＝Ｅ          （６）
の関係が導かれる。
【００３７】
また、ＡＤＡＴ＝ＬＥＬＴであり、両辺はそれぞれ
【数８】
2008176683.tif 000010
及び
【数９】
2008176683.tif 000011
と変形できることから、（７）式が導かれる。
【数１０】
2008176683.tif 000012
【００３８】
そして、（６）式の両辺については、それぞれ以下のように変形することができる。
【数１１】
2008176683.tif 000013
【数１２】
2008176683.tif 000014
これより、
【数１３】
2008176683.tif 000015
【００３９】
が導かれる。また、
ΘΘＴ＝Ｉ（行列Ｉは単位行列）          （９）
である。
【００４０】
ｍ行×ｎ列の直交行列Ｂについては、以下の式のように表すこととする。
【数１４】
2008176683.tif 000016
【数１５】
2008176683.tif 000017
【００４１】
は行列Ｂのｉ行におけるｎ個の行列要素（ｂｉ１、ｂｉ２、…、ｂｉｎ）を表す。
【００４２】
以上の前提条件の下で、グラム・シュミット法による直交・三角行列化の計算を行う。
まず、図５のステップＳ１で、変数ｉを初期化してｉ＝１とし、ステップＳ２で、ｉの値が行列Ｌ及びＥの次数であるｍ以下の条件を満たすか否かを判定する。ｉがｍ以下の値をとる場合は、ステップＳ３に進む。
【００４３】
ステップＳ３で、変数ｋを初期化してを代入する。上記のとおり、変数ｋは行列要素に無限大値を含むか否かに応じて分離した、第１及び第２の対角行列を識別するために用いられている。そして、ステップＳ４で、変数ｋが２以下の場合はステップＳ５に進み、第ｋの対角行列Ｅ（ｋ）のｉ行ｉ列の行列要素ｅｉ（ｋ）を初期化して０を設定し、ステップＳ６でｋに１を加算してステップＳ４に戻る。ステップＳ４でｋが２を超えると判定すると、ステップＳ７に進む。
【００４４】
ステップＳ７で、変数ｊを初期化して１を設定し、ステップＳ８で、ｊの値が行列Ｌ及びＥの次数であるｍ以下の条件を満たすか否かを判定する。ｊがｍ以下であるときは、ステップＳ９に進み、行列Ｌのｉ行ｊ列の要素ｌｉｊを初期化して０を設定し、ステップＳ１０でｊに１を加算してステップＳ８に戻る。ステップＳ８で、ｊがｍを超えると判定すると、行列Ｌのｉ行目の行列要素については全て初期化を終えたとして、図６のステップＳ１１に進む。
【００４５】
ステップＳ１１で、上記の処理で初期化により０が設定されている行列Ｌのｉ行の行列要素のうち、ｉ列の要素ｌｉｉについては１を設定する。ステップＳ１２で、行列Ｂのｉ行の要素を設定する。ここで、行列Ｌのｉ行の各要素については、先にステップＳ９及びステップＳ１１において設定した値（ｌｉｉ＝１を、ｌｉｊ＝０（ｊ≠０））を設定している。設定した値を用いると上記（５）は、
Ａ＝ＬＢ＝Ｂ
となることから、ステップＳ１２においては、行列Ｂのｉ行の要素は、
【数１６】
2008176683.tif 000018
【００４６】
と設定される。そして、ステップＳ１３でｉに１を加算してステップＳ２に戻る。ステップＳ２において、ｉがｍよりも大きい場合は、全ての行列要素について計算が完了したと判断して、初期値の設定を終えて図７のステップＳ１４に進む。
【００４７】
以上の処理は、求めるべき行列Ｅ、Ｌ及びＢの初期化を行うものであり、図４のアルゴリズムにおいては１行目から６行目に示す処理に対応している。
ステップＳ１４で、変数ｊを初期化して１を設定し、ステップＳ１５で、ｊの値がｍ以下であるとの条件を満たすか否かを判定する。ｊがｍ以下の値をとる場合は、ステップＳ１６に進み、ｋを初期化して１を設定して更にステップＳ１７に進む。
【００４８】
ステップＳ１７で、ｋが２以下であるか否かを判定する。ｋが２以下の場合はステップＳ１８に進み、（６）式より、第ｋの対角行列Ｅ（ｋ）のｊ行の対角要素を下の（１１）式のように設定する。
【数１７】
2008176683.tif 000019
【００４９】
そして、ステップＳ１９に進み、ステップＳ１８で設定したｅｊ（ｋ）が０であるか否かを判定する。「０」が設定されている場合は、ｍ次下三角行列Ｌ及びｍ×ｎ重み直交行列Ｂの行列要素については先に設定した初期値を用いることとして、ここでは特に設定し直さず、ステップＳ２０でｋに１を加算してステップＳ１７に戻る。ｅｊ（ｋ）に０以外の値が設定されている場合は、ステップＳ２１に進む。
【００５０】
ステップＳ２１で、変数ｉについて初期化してｊ＋１を設定し、ステップＳ２２で、ｉの値がｍ以下であるとの条件を満たすか否かを判定する。ｉがｍ以下の値をとる場合は、ステップＳ２３に進む。
【００５１】
ステップＳ２３で、行列Ｌ及びＢの要素について、下の（１２）式のように設定する。
【数１８】
2008176683.tif 000020
【００５２】
そして、ステップＳ２４でｉに１を加算し、ステップＳ２２に戻る。
ステップＳ２２で、ｉがｍを超えた場合は、ステップＳ２５に進み、ｊに１を加算してステップＳ１５に戻る。ステップＳ１５でｊがｍを超えた場合は、処理を終了する。
【００５３】
以上説明したように、図４及び図５から図７に示す方法によれば、変数ｊについて対角行列Ｅ（１）及びＥ（２）の対角要素を（１１）式に示すとおりに設定し、設定された値が０以外の場合には、対応するｉ行ｊ列の下三角行列Ｌの要素、及びｉ列の直交行列Ｂの要素について（１２）式に示すとおりに設定し、ｊの値がｍを超えると、処理を終了する。
【００５４】
ステップＳ１５以降の処理は、グラム・シュミット法を適用しており、図４のアルゴリズムにおいては、７行目以降の処理に対応している。
なお、上記の実施形態においては行列Ｌを下三角行列としているが、求めるべき三角行列を上三角行列Ｕとしても、Ｕ＝ＬＴであることから、同様の方法により共分散行列の演算が可能である。
【００５５】
また、上記の実施形態においては、カルマンフィルタに適用するための共分散行列の演算方法について例に挙げて説明しているが、これに限るものではなく、一般の共分散行列の計算に対しても適用が可能である。
【００５６】
更には、直交・三角行列化の方法としてはグラム・シュミット法について説明しているが、他の方法、例えば、ハウスホルダ法やギブンズ法を用いることも可能である。
【００５７】
（付記１）
共分散行列を行列Ａと対角行列Ｄと行列Ａの転置行列の積の形で表して演算を行う処理装置であって、
対角行列Ｄを、無限大値と対角行列Ｄ１との積と、対角行列Ｄ２との和で表したとき、対角行列Ｄ１及び対角行列Ｄ２を入力させる入力手段と、
対角行列Ｄ１及び対角行列Ｄ２を格納する第１の格納手段と、
行列Ａを格納する第２の格納手段と、
前記第１の格納手段及び第２の格納手段に格納された対角行列Ｄ１、対角行列Ｄ２及び行列Ａを用いて、前記共分散行列を下三角行列Ｌと対角行列Ｅと下三角行列Ｌの転置行列の積の形で表し、対角行列Ｅを無限大値と対角行列Ｅ１の積と、対角行列Ｅ２との和で表したときの下三角行列Ｌ、対角行列Ｅ１及び対角行列Ｅ２を計算する演算手段と、
得られた下三角行列Ｌを出力する出力手段と
を備えることを特徴とする処理装置。
（付記２）
前記第１の格納手段は、対角行列Ｄ１を格納する手段と対角行列Ｄ２を格納する手段とを含む
ことを特徴とする付記１記載の処理装置。
（付記３）
前記第１の格納手段は、対角行列Ｄ１及び対角行列Ｄ２の対角要素を格納する要素格納手段と、該要素格納手段に格納されている対角要素がそれぞれ対角行列Ｄ１あるいは対角
行列Ｄ２のいずれを表すかを示す情報を格納する格納手段とを含む
ことを特徴とする付記１記載の処理装置。
（付記４）
前記演算手段は、グラム・シュミット法により計算を行い、該グラム・シュミット法による計算においては、
対角行列Ｅ１及び対角行列Ｅ２の対角要素として０を設定し、下三角行列Ｌの対角要素として１を、下三角行列Ｌの対角要素以外の要素については０を設定することにより初期化を行い、
対角行列Ｄ１に基づいて、対角行列Ｅ１のｊ行の対角要素を設定し、
設定された対角行列Ｅ１のｊ行の対角要素の値が０でないときは、下三角行列Ｌの（ｊ＋１）行ｊ列からｍ行ｊ列までの要素を、対角行列Ｅ１のｊ行の対角要素の値及び対角行列Ｄ１に基づいて順次設定してゆき、
対角行列Ｄ２に基づいて、対角行列Ｅ２のｊ行の要素を設定し、
設定された対角行列Ｅ２のｊ行の対角要素の値が０でないときは、下三角行列Ｌの（ｊ＋１）行ｊ列からｍ行ｊ列までの要素を、対角行列Ｅ２のｊ列の対角要素の値及び対角行列Ｄ２に基づいて順次設定してゆく
処理を実行する
ことを特徴とする付記２または３記載の処理装置。
（付記５）
前記行列Ａ、対角行列Ｄ１及び対角行列Ｄ２は、カルマンフィルタの共分散行列の演算に使用される
ことを特徴とする付記１から３のいずれか１つに記載の処理装置。
（付記６）
共分散行列を行列Ａと対角行列Ｄと行列Ａの転置行列の積の形で表して演算を行うための共分散行列演算プログラムであって、
対角行列Ｄを、無限大値と対角行列Ｄ１との積と、対角行列Ｄ２との和で表したとき、対角行列Ｄ１及び対角行列Ｄ２を入力させ、
対角行列Ｄ１、対角行列Ｄ２及び行列Ａを格納手段に格納し、
前記格納手段に格納された対角行列Ｄ１、対角行列Ｄ２及び行列Ａを用いて、前記共分散行列を下三角行列Ｌと対角行列Ｅと下三角行列Ｌの転置行列の積の形で表し、対角行列Ｅを無限大値と対角行列Ｅ１の積と、対角行列Ｅ２との和で表したときの下三角行列Ｌ、対角行列Ｅ１及び対角行列Ｅ２を計算し、
得られた下三角行列Ｌを出力する
処理をコンピュータに実行させることを特徴とする共分散行列演算プログラム。
【図面の簡単な説明】
【００５８】
【図１】カルマンフィルタの更新処理についての概念図である。
【図２】本実施形態に係る共分散行列の演算を行う情報処理装置のブロック図である。
【図３】第１の格納部及び第２の格納部の構造を示す図である。
【図４】グラム・シュミット法による直交・三角行列化のアルゴリズムを示す図である。
【図５】グラム・シュミット法による直交・三角行列化方法を示すフローチャート（その１）である。
【図６】グラム・シュミット法による直交・三角行列化方法を示すフローチャート（その２）である。
【図７】グラム・シュミット法による直交・三角行列化方法を示すフローチャート（その３）である。
【図８】グラム・シュミット法による直交・三角行列化のアルゴリズムの変形例方法を示す図である。
【図９】グラム・シュミット法による直交・三角行列化方法の変形例を示すフローチャートである。
【符号の説明】
【００５９】
１      情報処理装置
４      第１の格納部
５      第２の格納部
６      第３の格納部
７      演算部
８      入力部
９      出力部
【図１】
2008176683.tif 000021
【図２】
2008176683.tif 000022
【図３】
2008176683.tif 000023
【図４】
2008176683.tif 000024
【図５】
2008176683.tif 000025
【図６】
2008176683.tif 000026
【図７】
2008176683.tif 000027
【図８】
2008176683.tif 000028
【図９】
2008176683.tif 000029