運動中のオブジェクトの大域的な形状を追跡するシステムおよび方法
ここに記載されているのは、運動中のオブジェクトの大域的な形状を追跡するシステムおよび方法である。ここではこの大域的な形状の初期輪郭に沿って１つ以上の制御点が定められる。運動中の１つ以上の制御点の各々は、オブジェクトの運動に伴って追跡される。運動中の制御点の位置の不確かさは、多数の手法を使用して表される。大域的な形状を制約するための不確かさは、予備の形状モデルを使用して利用される。択一的な１実施形態では、制御点毎に複数の外観モデルが形成され、各モデルによって形成される動きベクトルが組み合わされてオブジェクトの形状が追跡される。オブジェクトの形状の運動は、ディスプレイおよびカラーベクトルを使用して視覚的に追跡することができる。
【発明の詳細な説明】
【技術分野】
【０００１】
関連する出願との相互参照
本出願は、２００３年３月７日付け提出の米国暫定特許出願番号第60/452,669号、２００３年５月２７日付け提出の米国暫定特許出願番号第60/473,425号、および２００３年１０月３日付け提出の米国暫定特許出願番号第60/508,367号に優先権を主張するものであり、ここではこれらの全体を参考のために組み込む。
【０００２】
発明の分野
本発明は、運動中の形状を追跡する方法に関し、また殊に異分散ノイズ（heteroscedastic noise）がある場合に線形制約（linear constraint）を有する形状を追跡する方法に関する。
【０００３】
発明の背景
ほとんどの視覚的追跡アプリケーションに対して測定データは、不確定であり、また欠落していることもある。すなわち、画像はノイズおよび歪みを伴って記録されるのであり、その一方で閉塞（occlusion）により、関心領域の一部が観測できないようになることがあるのである。不確かさは、大域的に均一になり得るが、ほとんどの実世界のシナリオにおいて、これは本質的に異分散である。すなわち非等方であり不均質である。よい例が心エコー図（超音波心臓データ）である。超音波は、反射アーチファクトの影響、例えば薄膜から得られる鏡面反射体の影響を受けやすい。「見る方向」が１つであるため、鏡面構造体の垂直な面により、強いエコーが形成されるが、傾斜したまた「軸からずれた」表面により、弱いエコーしか形成されないか、またはまったくエコーが形成されないことがある（音響的な「ドロップアウト」）。心エコー図に対し、このようなドロップアウトは、超音波ビームに対して組織表面が平行である心臓の領域において発生し得る。
【０００４】
利用し易く、比較的コストがかからず、また非侵襲性であることに起因して心臓超音波画像は、心機能を評価するために広く使用されている。殊に心室の運動の分析は、虚血および梗塞の程度を評価するのに有効である。心内膜壁の検出またはセグメンテーションは、左心室の収縮および弾性の定量化に対する第１ステップである。いくつかの既存の方法の例に含まれるのは、ピクセルベースのセグメンテーション／クラスタリングアプローチ（例えばカラーキネシス（Color Kinesis））、オプティカルフローの種々の変形形態、変形可能なテンプレート（deformable template）およびマルコフ確率過程／フィールド、およびアクティブ輪郭（active contour）／スネーク（snake）である。いくつかの方法は２次元、３次元または４次元（３Ｄ＋時間）空間で使用される。
【０００５】
しかしながら既存のセグメンテーションおよび検出法は、心内膜壁の正確な局所的運動を再生しようとはしておらず、また多くの場合、壁に沿った動きの成分が無視される。このような単純化された処理は、目下の輪郭の法線に沿ってのみサーチする輪郭追跡器でも使用されている。これは局所的な壁の異常を検出するのには不適切である。それは異常な左心室の局所的な動きは、輪郭の法線から外れることが多く、また（音波検査者の手の動きまたは患者の呼吸による動きに起因する）平行移動または回転運動のような大域的な動きも、輪郭において法線を外れた局所的な動きを発生させることはいうまでない。望ましいのは、心内膜壁の大域的な形状と、その局所的な動きとを追跡して、局所的な壁の動きの異常を検出することである。この情報は、虚血および梗塞のさらなる診断に使用可能である。
【０００６】
一般的に共分散は、背後にある異分散ノイズを反映する画像フィーチャ（image feature）またはフロー推定（flow estimate）に割り当てることができる。全体的にノイズレベルが低くデータがきれいである場合、異分散的な性質は無視することができ、また大域的な不確かさを、局所的な推定の代わりに使用することができる。しかしながら極めてノイズの多い入力に対して、殊に空間的に変化する構造的なノイズを有する入力に対して、局所共分散行列において符号化される情報は、背後にある画像構造またはオブジェクトを高い信頼性でロバストに推定するのに重要になる。
【０００７】
追跡のフレームワークにおいてモデル制約を課すことはふつうに行われている。例としては、小球体またははパラメタ化された楕円などの簡単なモデル、および判別テンプレート（discriminative template）のような複雑なモデルが含まれる。ほとんどの実践的なケースにおいて部分空間モデルは形状追跡に適している。それは、主要な形状変化を捕捉するモードの数は限られており、また通常、この形状を表すのに使用されるフィーチャコンポーネントの元々の数よりも格段に小さいからである。さらにＰＣＡ（Principal Component Analysis）ベースの固有形状部分空間によって、任意の複雑な形状変化を捕捉することができる。ここでこれは元々の空間において、極めて単純なパラメトリックモデルを有しているとしても、極めて非線形である。
【０００８】
測定ベクトルが異分散ノイズによって影響を受ける場合、制約的な部分空間への正射影は、そろっていないだけではなく、情報のロスという意味で極めて多くのダメージを受ける。ノイズが等方的であると共に均一である特殊な場合にのみ、これをそろえることができるのである。
【０００９】
しかしながら部分空間制約追跡への既存の取り組みのほとんどは、測定値における異分散ノイズを考慮していない。「点分散モデル（Point Distribution Model）」または「アクティブ形状モデル（Active Shape Model）」において、ＰＣＡベースの部分空間形状モデルは、ランドマーク点（landmark point）対応関係を有するトレーニング形状（training shape）に基づいて導出される。結果的に得られる固有形状の部分空間によって、トレーニングデータセットにおける最も重要な変化が捕捉される。検出時にモデルは摂動されて、候補の位置におけるテスト画像との一致のための合成画像が形成される。しかしながら測定ノイズは、この処理においてモデル化されていなかった。
【００１０】
異分散ノイズ特性が利用可能である場合であっても、これらは通例、部分空間モデルフィッティング中に廃棄されていた。例えば、測定においてフルの共分散行列が捕捉される公知の１アプローチにおいて、かなりその場限りの特別な閾値判定が適用されて、測定値の平均値が、モデルの共分散によって定められる超楕円体の制約に限定される。この操作は、測定ノイズには依存しない。
【００１１】
別の公知のアプローチは、２ステップアプローチを適用して、カルマンフィルタリングフレームワークにおいて形状空間制約を課している。この形状空間は、線形変換されたアフィン部分空間または固有部分空間である。しかしながら形状空間への投影は正射影であり、測定の異分散のノイズを考慮していない。したがってこのアプローチは投影中に情報をロスすることになる。
【００１２】
別の公知のアプローチではガウス分布を使用して関心対象のオブジェクト（この場合、顔である）の外観を適応的にモデリングしている。ここでこれはＥＭアルゴリズムを使用して学習される。本発明と同様に、局所的な不確かさは共分散行列に捕捉される。異なるのは、本発明では殊に、非等方性の不確かさが存在する場合に投影にわたる交わりの決定的な選択および部分空間モデルの制約を検討することである。
【００１３】
別の公知のアプローチでは、オプティカルフロー推定中に暗黙的に部分空間制約使用しており、ここではフローの不確かさも利用される。別の応用に対する別のフレームワークにおいてではあるが、本発明によって認識されるのは、「部分空間投影処理においてフローベクトルの信頼性が高ければ高いほど、大きな影響がある」ことである。
【００１４】
別の公知のアプローチでは、基本行列および楕円のフィッティングに対して異分散回帰を適用している。このフィッティングは、元々の空間において、パラメータ化されたモデルで行われる。本発明では、形状の変化のパラメータ化（極めて複雑であり、高度に非線形である）は回避される。その代わりに本発明では、例えばＰＣＡを介して、部分空間線形確率モデルを形成し、またフィッティングしたデータの平均および共分散の両方について閉じた形の解決手段を得る。
【００１５】
またＭ推定器またはRANSACによるロバストモデルマッチングも適用されており、モデルの異常値であるデータ成分の影響が制限または除去される。ここでも（空間的または時間的に）局所的に変化する不確かさは、これらのフレームワークにおいて使用されない。
【００１６】
関連する別のアプローチに含まれるのは、データインピュテーション（data imputation）、すなわち、欠落したデータを信頼性のある値で「埋める」ことである。この領域における取り組みは、音声認識、医用画像分析、および社会科学に対して広い応用を有する統計学において定着している。しかしながらデータインピュテーションの定式化の問題は、ふつう０−１の利用可能性を前提としている。すなわちデータ成分は欠落しているか、利用可能であるかのいずれかを前提としているのである。部分空間モデルの制約と、形状ダイナミックおよび測定データノイズの異分散の性質についての情報および形状ダイナミックスについての情報とを融合するための統一的なフレームワークが必要である。
【００１７】
発明の要約
本発明は、運動中のオブジェクトの大域的な形状を追跡するシステムおよび方法に関する。ここでは大域的な形状の初期輪郭に沿って１つ以上の制御点が定められる。オブジェクトが運動に伴って上記の１つ以上の制御点の各々が追跡される。運動中の制御点の位置の不確かさが推定される。この不確かさを表す１形態は共分散行列である。部分空間形状制約モデルを使用する場合、この不確かさは、非正射影射影および／または情報フュージョンを使用して利用される。後続の各輪郭が表示される。
【００１８】
本発明はまたオブジェクトの形状の運動を視覚的に追跡するためのシステムに関する。１つ以上の第１のカラーベクトルが生成されて、この形状の輪郭に沿った制御点の収縮が表される。１つ以上の第２のカラーベクトルが生成されて、この形状の輪郭に沿った制御点の拡張が表される。上記の第１および第２のカラーベクトルは周期的に表示され、この形状の動きが表示される。
【００１９】
本発明はまた運動中のオブジェクトの大域的な形状を追跡する方法に関する。１つ以上の制御点を大域的な形状に沿って定める。オブジェクトの運動に伴って１つ以上の制御点の各々を追跡する。各制御点に対して複数の外観モデルを形成する。この形状を追跡する各モデルによって形成される動きベクトルを組み合わる。
【００２０】
図面の簡単な説明
添付の図面に関連して以下に本発明の有利な実施形態をより詳細に説明する。ここでは類似の参照番号は類似の要素を示す。
【００２１】
図１は、本発明にしたがって形状追跡法を実現する例示的なシステムのブロック図であり、
図２は、音響的なドロップアウトのエリアと、推定された局所的な壁運動の不確かさとを示す心エコー図画像であり、
図３は、制御点の位置の不確かさが極めて非等方であり不均質である心内膜輪郭を示す左心室の心エコー図であり、
図４は、相異なるα値で強く適合化したＰＣＡ（ＳＡ（strongly adapted）−ＰＣＡ）モデルおよびインクリメント式ＰＣＡモデルの例を示しており、
図５ａ〜５ｃは、左心室の心エコー図の画像を示しており、ここでは心内膜壁が初期化され（ａ）、ＩＰＣＡ（ｂ）およびＳＡ−ＰＣＡ（ｃ）を使用して追跡され、
図６ａおよび６ｂは左心室の心エコー図の画像を示しており、ここでは心内膜壁の動きは本発明にしたがって追跡され、
図７ａ〜７ｄはテスト輪郭を示しており、これは本発明にしたがい、先端の４室側から見た左心室の心内膜壁の心エコー図画像を表しており、
図８ａ〜８ｄはテスト輪郭を示しており、これは本発明にしたがい、胸骨傍短軸（parasternal short axis）側から見た左心室の心内膜壁の心エコー図画像を表しており、
図９は、本発明による画像視覚化方法を例示する心エコー図画像を示しており、
図１０は、本発明の１実施形態による、マルチモデル制御点ベースの追跡器を示している。
【００２２】
詳細な説明
本発明は、異分散ノイズがある場合に線形制約によって形状を追跡する方法に関する。このような方法が利用される１例では、心筋層壁の大域的な形状と、その局所的な動きとが追跡されて、心臓における局所的な壁運動の異常が検出される。この方法は、心臓の心内膜壁または心外膜壁を追跡するためにも使用可能である。本発明が、頭部の運動、顔の特徴、手の運動または別の身体運動などの人間の特徴部分の運動を識別するなどのように形状追跡が有効な別の応用にも使用でき、しかしながらこれらに限定されるないことを当業者は理解されたい
本発明はまた２次元、３次元および４次元（３Ｄ＋時間）における心臓、肺または徐々に成長する腫瘍などの解剖学的構造の医学的な分析にも使用可能である。
【００２３】
本発明を説明するため、左心室の心内膜壁を追跡する例を述べる。図１は、本発明による左心室心内膜壁形状追跡方法を使用する心エコーシステムの例示的なアーキテクチャを示している。超音波トランスデューサのような医用センサ１０２が使用されて患者に検査が行われる。センサ１０２は、所定の健康診断に見合った医療用測定値を得るために使用される。例えば心臓に問題を経験した患者は心エコー検査を受けて、これが所定の心疾患の診断に役立てられる。超音波システムにより、さまざまな視点から心臓の２次元、３次元および４次元（３Ｄ＋時間）の画像が供給される。
【００２４】
センサ１０２から得られた情報は、ワークステーションまたはパーソナルコンピュータとすることの可能なプロセッサ１０４に伝達される。プロセッサ１０４は、センサデータを画像に変換し、これがディスプレイ１０８に伝達される。またディスプレイ１０８は、この画像に関連する別の図形情報または情報のテーブルも伝達することができる。本発明によれば、プロセッサ１０４には、心内膜壁の初期輪郭を表すデータも供給されている。このデータは、医師または音波検査者などのユーザによって手動で供給されるか、またはプロセッサ１０４によって自動的に供給される。輪郭は個別の点の列からなっており、その運動はプロセッサ１０４によって追跡され、ディスプレイ１０８に示される。個別の点がどのように追跡されるかについての詳細を以下詳述する。
【００２５】
医用センサ１０２からのデータに加えてプロセッサ１０４は別のデータ入力を受け取ることができる。例えばこのプロセッサは、プロセッサ１０４に付随したデータベース１０６からデータを受け取ることができる。このようなデータに含まれ得るのは、心内膜壁に対する輪郭形状の候補を表す部分空間モデルである。これらの部分空間モデルは、複数の患者の代表である左心室の画像か、または統計的な情報に基づいてコンピュータで生成した輪郭形状のモデルとすることが可能である。プロセッサ１０４は、ベイズカーネルマッチングまたはオプティカルフローベース法などの既知のアプローチを使用して、輪郭形状の個々の点を追跡する。追跡中に累積される誤差は、マルチテンプレートアダプティブマッチングフレームワークを使用して修正される。追跡の不確かさは各点において共分散行列の形で表され、これは引き続いて、非正射影の投影を使用することにより、部分空間形状制約によって十分に活用される。
【００２６】
図２は典型的な心エコー図画像を示している。音響的なドロップアウトを有する左心室の心内膜壁の部分は、実線の楕円２０８によって示されている。局所的な壁運動の推定は、点線の楕円２０２，２０４によって示されている。音響的なドロップアウトがあるため、心内膜壁はこの画像においてつねに最も強いエッジになるわけではない。本発明によれば、部分空間モデルの制約と、形状ダイナミックおよび測定ノイズの異分散的な性質についての情報および形状ダイナミクスについての情報との融合に対して、統合化されたフレームワークが使用される。この部分空間モデルは、所定の部分空間分布、例えばガウス分布、または単純な部分空間の制約、例えば固有空間モデルの形を取り得る。
【００２７】
本発明では心内膜壁を表す輪郭上の個々の制御点が追跡される。この追跡は、ベイズカーネルマッチングアプローチまたはフローベースのアプローチを使用することによって行うことができる。ベイズカーネルマッチングアプローチの例は、共同発明者Dorin Comaniciuによって著されたBayesian Kernel Tracking, Annual Conf. of the German Society for Pattern Recognition （DAGM'02）, Zurich, Switzerland, 438-445, 2002なる題名の論文に記載されており、その全体を参考のために組み込む。個々の点を追跡するオプティカルフローベースのアプローチの例は、"Density Estimation-Based Information Fusion for Multiple Motion Computation"なる名称の同時係属出願の明細書第10/681,02に記載されており、その全体を参考のために組み込む。本発明によれば、追跡中の誤差の累積は、心運動の周期的な性質を利用し、マルチテンプレートアダプティブマッチングフレームワークを使用することによって修正される。
【００２８】
追跡の不確かさは各点において共分散行列によって表され、これは続いて、非正射影の投影を使用することにより、部分空間形状制約によって十分に利用される。
【００２９】
追跡フレームワークは、イメージ列にわたる２つのステップの繰り返し処理である。制御点を有する初期輪郭は第１のイメージに（手動または自動で）引かれ、つぎに後続の画像毎に順番に各制御点がまず独立に追跡される。ここでは非等方性の不確かさも記録される。第２のステップとして新しい輪郭が、非正射影の投影を使用して、実現可能な部分空間に射影される。この実現可能な部分空間は、トレーニング輪郭に基づいて学習され、また追跡器に利用可能な初期輪郭を使用して都度のケースに適合される。また初期化の信頼度も考慮される（すなわち手動の初期化は、完全自動化の初期化よりも信頼性が高い）。
【００３０】
上に示したように、追跡処理中に複数のテンプレートが使用される。複数のテンプレートを使用することにより、結果的に目下のケースに対する形状の統計がより正確に表現されることになる。ベイズカーネルマッチングアプローチでは、第１のテンプレートは、初期化されたフレームから選択される。後続のテンプレートは、それが既存のテンプレートとは十分に異なっている共に、位置決めのために十分な情報を与える場合に追加される。ここでこれは、それ自体とのカーネルマッチングによって測定される。より多くの情報を与えるパッチは、それ自体とマッチングした場合に、より高い信頼性を有するパッチである。
【００３１】
１つ以上のテンプレートを使用することの決定は、心臓の運動が周期的であり、したがって１周期内の異なる外観のパターンが、後続の周期においても再度現れるという観察に基づく。複数のテンプレートとマッチングするため、マッチングはテンプレート毎に行われ、最良のマッチングを有するものが選択される。または計算を節約するため、前にマッチしたテンプレートの近くのテンプレートだけをマッチングすることができる。ここでの周期的な運動が使用される。
【００３２】
位置決めのマッチングについての不確かさは、最適な位置の近くで尤度マップ（likelihood map）によって計算され、ここでこれはカーネルマッチング処理を使用して同様に計算される。この尤度面（likelihood surface）はつぎに共分散行列を推定するために使用される。ここでこれは例えば、ヘッセの行列の重み付けされた推定の逆行列または２次元のガウス分布とのフィッティングによって行われる。
【００３３】
フローベースのアプローチの場合、通例のオプティカルフローの実現では隣りのフレームだけを使用しているため、長いシーケンスにわたる追跡処理は、誤差の累積およびドリフトに影響されやすくなっている。ここでは複数のテンプレートがフロー計算に使用されており、カーネルマッチングに対するのと同じである。新しいテンプレートは、フローの不確かさが大きい場合にはいつでもどこでも追加されるの対して、局所的な傾斜は既存のテンプレートから区別可能である。
【００３４】
制御点毎に新たな位置が得られた後、つぎのステップでは、人間の心臓の「正しい」形状変化を捉える統計的なモデルによって、全体的な形状を制約する。点分散モデルとしても知られているＰＣＡベースの形状モデルまたはアクティブ形状モデルが使用される。
【００３５】
心臓輪郭にわたる不確かさは、均一でもなく、すなわちエリアによっては、例えば信号のドロップアウトに起因して他のエリアよりも悪く、また等方性でもない、例えば、エッジに沿った位置決めの方が傾斜方向に沿った位置決めよりも悪いのである。図３は、このような輪郭にわたる非等方性かつ不均一なノイズの例を示している。図からわかるように個々の点、例えば点３０２，３０４，３０６および３０８ははじめ、心内膜壁の輪郭に沿っていると識別されている。点毎に確かさの測定が行われ、また各点のまわりに楕円体が形成される。この惰円体は、この所定の点が正しい位置にある確かさのレベルを表している。楕円体が大きければ大きいほど、この所定の点の位置についての不確かさのレベルが高くなる。図３からわかるように点３０２は、それを取り囲む比較的小さい楕円体を有し、輪郭におけるこの点の位置についての確かさの程度が高いをことを示している。点３０６および３０８は、それらの周りに比較的大きな楕円体を有しており、これらの点の位置についての不確かさが高いことを示している。
【００３６】
新しい輪郭的にわたる一般的な不確かさのモデルを考慮することによって、最適な射影はもはや正射影ではなくなる。部分空間における最適な解は実際に、形状モデル部分空間における交差分布における最尤度の形状である。さらにこの部分空間に分布モデルがある場合、このような追加情報を無視する理由はない。以下では、非正射影の投影の仕方について詳しい分析を示す。また部分空間モデル分布が利用可能な場合、このようなモデルの情報と、不確かな入力との融合のさせ方を示す。
【００３７】
本発明は、部分空間モデルベースフュージョンアプローチを使用する。同じｎ次元変数ｘのノイズを伴った２つの測定値が与えられる場合、各々は多次元のガウス分布ｐ１およびｐ２によって特徴付けられ、ｘの最尤度推定は、２つの重心に対するマハラノビス距離の最小和を有する点である。すなわち、
ｘ* ＝ argmin ｄｍ
であり、ここで
ｄｍ ＝ (ｘ−ｘ１)ＴＣ１−１(ｘ−ｘ１)＋(ｘ−ｘ２)ＴＣ２−１(ｘ−ｘ２） （１）
である。
【００３８】
ｘについての導関数をとると、
ｘ* ＝ Ｃ(Ｃ１−１ｘ＋Ｃ２−１ｘ２)
Ｃ ＝ (Ｃ１−１＋Ｃ２−１) （２）
を得る。
【００３９】
これは、ｘのＢＬＵＥ（best linear unbiased estimate）として知られている。
【００４０】
上記の２つのガウシアンのうちの１つが次元ｐの部分空間にある、例えばＣ２が特異値を有するとする。Ｃ２ ＝ ＵΛＵＴの特異値分解（singular value decomposition）により、（ただしＵ ＝ [ｕ１，ｕ２，…，ｕｎ]，ここでｕｉは正規直交であり、Λ ＝ diag｛λ１，λ２，…，λｐ，０，…，０｝である）、式（１）におけるｘ２に対するマハラノビス距離をつぎの標準的な形に書き換える。すなわち、
【００４１】
【数１】
である。
【００４２】
λｉが０なると、Ｕ０Ｔｘ ＝ ０，ただしＵ０ ＝ [ｕｐ＋１，ｕｐ＋２，…，ｕｎ]でない限り、ｄｍ，２は無限大になる。（ここで一般性を失うことなく仮定したのは、部分空間が元々の空間の原点を通ることである。ｘ２が部分空間にあるため、Ｕ０Ｔｘ２ ＝ ０である。）
したがって式（３）においてゼロでないλｉに相応する項だけ残しておけばよい。すなわち、
【００４３】
【数２】
であり、ここで
Ｃ＋はＣの疑似逆行列（pseudoinverse）、Ｕｐ ＝ [ｕ１，ｕ２，…，ｕｐ]である。
【００４４】
さらにＵ０Ｔｘ ＝ ０であるため、ｘを別の形で表して、この制約を反映することができる。すなわち、１×ｐベクトルｙに対して、
ｘ ＝ ＵＵＴｘ ＝ Ｕ([Ｕｐ｜Ｕ０]Ｔｘ) ＝ Ｕ[ｙ｜０] ＝ Ｕｐｙ （５）
である。式（１）はいまや、特異値にかかわらずつぎの一般的な形をとる。すなわち、
ｄｍ ＝ (Ｕｐｙ−ｘ１)ＴＣ１−１(Ｕｐｙ−ｘ１)＋(Ｕｐｙ−ｘ２)ＴＣ２＋(Ｕｐｙ−ｘ２) （６）
である。
【００４５】
ｙについての導関数をとると、
ｙ* ＝ Ｃｙ＊ＵｐＴ(Ｃ１−１ｘ１＋Ｃ２＋ｘ２) （７）
Ｃｙ＊ ＝ [ＵｐＴ(Ｃ１−１＋Ｃ２＋)Ｕｐ]−１ （８）
ｘ* ＝ Ｕｐｙ* ＝ Ｃｘ＊(Ｃ１−１ｘ１＋Ｃ２＋ｘ２) （９）
Ｃｘ＊ ＝ ＵｐＣｙ＊ＵｐＴ （１０）
である。
【００４６】
式（７）は情報空間への融合を示しており、これに続いて部分空間への変換が行われる。式（９）は、元々の空間に戻す座標変換である。Ｃｘ＊およびＣｙ＊は、ｘ*およびｙ*に対する相応の共分散行列であることは示すことができる。
【００４７】
この解は、正しい逆行列の代わりに疑似逆行列を使用することによって式（２）単純に一般化したものではないことに注意されたい。これはｘ*が部分空間にあることを制約しない。
【００４８】
択一的には式（７）および（８）を
ｙ* ＝ (ＵｐＴＣ１−１Ｕｐ＋Λｐ−１)−１(ＵｐＴC１−１ｘ１＋Λｐ−１ｙ２) （１１）
と書き表すことができる。
【００４９】
ここでｙ２はＵｐによって張られた部分空間におけるｘ２の変換された座標であり、Λｐ ＝ diag{λ１，λ2，…，λｐ}である。
【００５０】
式（１１）は、２つの分布の部分空間におけるＢＬＵＥフュージョンとして示すこともでき、その１つはＮ(ｙ２，Λｐ)であり、別の１つは部分空間におけるＮ(ｘ１，Ｃ１)の「交点」（射影ではない！）Ｎ((ＵｐＴＣ１−１Ｕｐ)−１ＵｐＴＣ１−１ｘ１，(ＵｐＴC１−１Ｕｐ)−１)である。
【００５１】
トレーニングサンプルの大きなプールから学習した統計的形状モデルを使用して、所定の心臓から輪郭を導くことは、時として問題となることがある。本発明によれば、トレーニングサンプルは、目下の心臓の形状モデルを得るために使用され、一般的な心臓の形状モデルを得るために使用されるのではない。したがって、目下のケースに対して知られていることに向かって一般的なモデルを適合させようとする強いモチベーションがある。本発明によれば、患者の心臓の心内膜壁の（手動または自動検出による）初期輪郭は検証されて、既存のＰＣＡモデルに適合化される。
【００５２】
心内膜壁の実際の輪郭を決定してその運動を追跡する際に２つのアプローチが考慮される。すなわち、初期輪郭が決定的（deterministic）である１アプローチと、初期輪郭が不確かである別の１アプローチとが考慮されるのである（これは、初期輪郭が自動検出アルゴリズムからのものである場合に発生する。これはまた不確かさも発生させる）。
【００５３】
初期輪郭が１点である（決定的に確定している）と仮定される場合、ＳＡ−ＰＣＡ（strongly-adpted-PCA）モデルが使用されてこの点の運動が追跡される。ここで仮定されるのは、旧いＰＣＡモデル（現在のケースを除く）と、現在のケースに対して初期化された輪郭とが、一緒になって現在のケースの変化を表すが、相対的なエネルギー（すなわち代表的な出力）はそれぞれａおよび１−ａ、ただし０＜ａ＜１となることである。言い換えると、現在のケースの形状変化の１部分は、一般的なモデルによって表されるのに対して、残りは初期輪郭の方向に捉えられるとを仮定するのである。
【００５４】
ＰＣＡモデルは、その平均、固有値行列および固有ベクトル行列で表され、またそれぞれｘｍ，ΛおよびＵで表される。元々の完全な共分散行列Ｃが格納される場合（これは、元々の次元数（dimensionality）が許容されないほどに高すぎないケースである）、適合された平均および共分散行列は単に、寄与する２つのソースの重み付けされた和である。すなわち、
ｘｍ，ｎｅｗ ＝ αｘｍ＋(１−α)ｘ （１２）
Ｃｎｅｗ ＝ α(Ｃ＋(ｘｍ−ｘｍ，ｎｅｗ）(ｘｍ−ｘｍ，ｎｅｗ）Ｔ)
＋(1−α) (ｘｍ−ｘｍ，ｎｅｗ）(ｘｍ−ｘｍ，ｎｅｗ）Ｔ
＝ αＣ＋α(１−α) (ｘｍ−ｘｍ，ｎｅｗ）Ｔ （１３）
である。
【００５５】
固有分析をＣｎｅｗに実行して、新しい部分空間モデルを得ることができる。
【００５６】
Ｃが格納されておらず、｛ｘｍ, Λ, Ｕ｝だけが部分空間において利用可能な場合、直接的な代数的操作によって以下のようにして、適合化された固有分析結果｛ｘｍ，ｎｅｗ, Λｎｅｗ, Ｕｎｅｗ｝に到達することができる。すなわち、
初期輪郭ｘは、ｘＳ ＝ ＵＴｘｄ（ただしｘｄ ＝ ｘ−ｘｍ）の部分空間成分を、またｘｒ ＝ (ｘ−ｘｍ)−ＵｘＳｉの残差ベクトルを有する。ｘｒを正規化してノルム１を有するようにしたものをｘｒｕとする（またはｘｒのノルムがゼロの場合ゼロとする）。
【００５７】
結合されたエネルギーを表す適合された固有ベクトル行列はつぎの形を有する。すなわち、
Ｕｎｅｗ ＝ [Ｕ, ｘｒｕ]Ｒ （１４）
である。
【００５８】
ＲおよびΛｎｅｗは、以下の固有分析問題の解になる。すなわち、
【００５９】
【数３】
であり、ここでｅｒ ＝ ｘｒｕＴ(ｘ−ｘｍ)は残りのエネルギーである。
【００６０】
上記の式は、新しいデータと旧いデータとの間の調整可能なエネルギー比を有し、ＩＰＣＡよりも一般的である。αを点の総数に対するモデルにおける点の数の比と設定した場合、これらはＩＰＣＡに等しくなる。通例、これは１に極めて近くなる。それはトレーニングセットにおける輪郭の数はふつう大きいからである。αを比較的小さな値（例えば０．５）に設定した場合、ＰＣＡモデルは目下のケースに向かって強力に適合化される。このことがその名前の由来である。図４は、相異なるα値を有するＩＰＣＡおよびＳＡ−ＰＣＡの簡単な２次元の図を示している。
【００６１】
点４０４は目下のケースを表している。各ｘ４０２は、特定のモデルに相応するトレーニング点を表している。楕円４０６は、元々のモデル分布を表している。楕円４０８は、強く適合化されたＰＣＡモデルに相応するインクリメント式ＰＣＡモデルを示しており、ここではα ＝ ０．９９である。楕円４１０は、α ＝ ０．５の強く適合化されたＰＣＡモデルを示している。楕円４１２は、α ＝ ０．１の強く適合化されたＰＣＡモデルを示している。各楕円は、相応する分布の等しく９０％で確からしい輪郭を示している。
【００６２】
実際に目下の心臓から得られる輪郭は、一般的なトレーニングセットにおける輪郭よりも、同じ心臓の初期輪郭に似ていることのほうがはるかに多い。殊に目下の心臓がトレーニングセットで良好に表されない異常な形を有する場合はそうである。ここのシステムでは、αを０．５にセットした。これによって初期輪郭からの強い影響が可能である。すなわち、モデルエネルギーの５０％は、初期輪郭から得られる（９９個の例でトレーニングされたモデルによってＩＰＣＡが適合される場合、これはわずかに１％に減少する）。
【００６３】
図５ａ〜５ｃはＩＰＣＡとＳＡ−ＰＣＡとの比較を示している。胸骨傍の短軸側から見たものは、（凹形部分のある）異常な形を有し、これに対してトレーニングセットは円形が大多数を占めている。初期輪郭（図５ａ）においてではあるが極めてわずかな重み付け（＜０．０１％）で記録されているインクリメント式ＰＣＡモデルは、目下の形状の凹形の性質を捉えておらず、この輪郭を典型的な円形（図５ｂ）に制約している。この結果、実際にインクリメント式ステップのない旧いPCAモデルを使用して得られたのと同じになっている。α ＝ ０．５の適応形ＰＣＡモデルでは、真の境界により一層良好にフィットしている（図５ｃ）。
【００６４】
フュージョンアプローチによる非正射影の投影と、ＳＡ−ＰＣＡモデル適合化との間の微妙ではあるが重要な相互作用はつぎのようになっている。モデル平均および共分散を有するフュージョンは、部分空間にはあるがモデル分散とはあまりに離れている輪郭を濾波するために必要である。しかしながらこのような（正射影よりも）強い制約は、必然的に低い尤度または異常な心臓輪郭を変化させる。ＳＡ−ＰＣＡモデルは、与えられた初期輪郭において供給される追加情報を使用することにより、目下のケースの方向に向かって強力にモデルをシフトすることを介してこのジレンマの間を取り持つのである。
【００６５】
ダイナミック処理によって定められる予想およびノイズ測定値からの情報をカルマンフィルタによって融合する。形状追跡に適用する場合、付加的な大域的制約が必要であり、この制約によって全体的な形状を有利な範囲に安定させる。本発明では、統合化されたフュージョンフレームワークを利用して、部分空間モデルの制約をカルマンフィルタに組み込む。
【００６６】
カルマンフィルタに対して、測定値の更新の式はつぎの形である。すなわち、
ｘｋ＋１｜ｋ＋１ ＝ ｘｋ＋１｜ｋ＋Ｋ(ｚｋ＋１−Ｈｚｋ＋１｜ｋ) （１６）
であり、ここで
Ｋ ＝ Ｐｋ＋１｜ｋＨＴ(ＨＰｋ＋１｜ｋＨＴ＋Ｒ)−１ （１７）
Ｐｋ＋１｜ｋ＋１ ＝ (Ｉ−ＫＨ)Ｐｋ＋１｜ｋ （１８）
Ｐｋ＋１｜ｋ ＝ ＳＰｋ｜ｋＳＴ＋Ｑ （１９）
である。
【００６７】
ここでＰは推定誤差共分散であり、ｘｉ｜ｊは時間jにおいて与えられた状態の時間iにおける状態推定である。この測定モデルはｚｋ ＝ Ｈｘｋ＋ｒｋであり、ただしｒｋは共分散Ｒを有する測定値ノイズを表す。このシステム／処理モデルはｘｋ＋１ ＝ Ｓｘｋ＋ｑｋであり、ここでｑｋは共分散Ｑを有するシステムノイズを表す。
【００６８】
上記のアプローチを使用することにより、部分空間制約および異常分散ノイズを有するカルマンフィルタの更新の式は、
ｘｋ＋１｜ｋ＋１ ＝ Ｐｋ＋１｜ｋ＋１((Ｐｋ｜ｋ＋Ｑ)＋ｘｋ｜ｋ＋Ｒ−１ｚ＋Ｃ２＋ｘ２) （２０）
Ｐｋ＋１｜ｋ＋１ ＝ Ｕｐ[ＵｐＴ(((Ｐｋ｜ｋ＋Ｑ)＋＋Ｒ−１＋Ｃ２＋)Ｕｐ]−１ＵｐＴ （２１）
によって与えられる。ここで仮定するのは、上記のシステムノイズ共分散Ｑは部分空間において制約されることである。情報空間において利用可能なすべての知識を組み合わせる解の対称性を観察する。これらの式によって、システムダイナミック、部分空間制約およびノイズ情報の統合化されたフュージョンが得られる。これらは、追跡システムに影響を及ぼすさまざまな不確かさの完全な表現を表す。
【００６９】
本発明の択一的な１実施形態によれば、制御点ベースのオブジェクトの表現と、フレームにわたる統合モデルへのロバストな融合とに基づく追跡技法が使用される。制御点ベースのオブジェクト表現の集合は、異なる時点に取得されて維持される。制御点によって示唆される推定の運動は、ロバストに結合されてオブジェクトのつぎの位置が決定される。このオブジェクトの視覚的な追跡は時間にわたっていくつかのモデルを維持することによって達成される。結果的に得られるのは確率密度関数のパラメトリックでない表現であり、これはオブジェクトの外形を特徴付ける。追跡は、オプティカルフローによるその不確かさと、運動の推定とを各モデルから独立に得ることによって行われる。各制御点に対する最終的な推定は、可変帯域幅密度フュージョンＶＢＤＦ（Variable-Bandwidth Density-based Fusion）のようなロバストな融合手法を使用して計算される。ＶＢＤＦにより、変位密度関数（displacement density function）の最も重要なモードの位置が計算され、その間にその不確かさが考慮される。ＶＤＢＦプロシージャにより、複数のデータソースと、変位の推定における異常値とが管理される。閉塞は、大きな残差に対し、推定の不確かさを介して処理されるのがふつうである。モデルは、フローが独立して計算される複数の領域に分割される。残りの位置合わせ誤差は推定の共分散行列のスケールを計算するために使用される。したがって信頼のできない変の影響が低減される。
【００７０】
閉塞および外形の変化によるオブジェクト追跡の試みは本発明において、マルチモード制御点ベースのアプローチを介して処理される。２次元外形モデルに対して複数のモデルを維持することは、これを単峰形分布（unimodal distribution）に制約せず、またＶＤＢＦメカニズムにより、ロバストに複数の推定が統合されて、制御点毎に最も優位な運動が決定される。追跡中の変化をモデリングするため、時間にわたるオブジェクトの外形の複数のサンプルが維持される。図１０には、本発明による例示的なマルチモデルの制御点ベースの追跡器が示されている。最も上の３つのフレーム１００２，１００４および１００６は、モデルセットの目下のサンプルを示しており、重なり合った制御点の関連する集合を有する。制御点ベースのアプローチは、大域的な表現よりもロバストであり、したがって照明の変化および中断に対してさほど敏感でない。本発明の別の利点は、部分的な閉塞が制御点レベルでマッチングの尤度を分析することによって処理できることである。
【００７１】
各制御点は独立して処理され、その位置および共分散行列は現在の画像においてすべてのモデルテンプレートについて推定される。例えば、制御点の１つはグレイな矩形１０１４によって示され、各モデルについてのその位置および不確かさはＩｎｅｗ１００８によって示されている。ＶＢＤＦロバストフュージョンプロシージャが適用されて、関連した不確かさと共に最も優勢な運動（モード）が決定される。これはフレーム１０１０に示したとおりである。各制御点の推定位置における変化が、閉塞または外形の変化によることに注意されたい。
【００７２】
目下のフレームにおける制御点の位置はさらに、大域的なパラメトリックな運動モデルによって制約される。類似性変換モデルおよびそのパラメタは、各制御点における信頼度（confidence）を使用して推定される。したがって信頼性の高い制御点は、大域的な運動の推定により大きく貢献するのである。残差、基準の外観が比較的低い場合、現在のフレームがモデルセットに加えられる。閾値を選択して、オブジェクトが重要な閉塞を有する箇所に画像が加えられないようにする。モデルにおけるテンプレートの数は可変でも固定でもよい。テンプレートの数が固定の場合、所定のテンプレートが廃棄される（例えば、最も旧いテンプレートを廃棄する）方式を提供することが可能である。
【００７３】
ＶＢＤＦ推定器は、アダプティブカーネル帯域幅（adaptive kernel bandwidth）を有するノンパラメトリック密度推定（nonparametric density estimation）に基づいている。ＶＢＤＦ推定器は、密度関数の最重要モードの位置として定義される。モード計算は、マルチスケール最適化フレームにおける可変帯域幅平均シフト法（variable-bandwidth mean shift technique）に基づいている。
【００７４】
ｘｉ ∈Ｒｄ，ｉ＝１…ｎを利用可能なｄ次元推定値とし、ここで各々は共分散行列Ｃｉによって与えられる付随した不確かさを有する。密度関数の最重要モードは、マルチスケール方式で反復的に決定される。帯域幅行列Ｈｉ ＝ Ｃｉ＋α２Ｉは、各点ｘｉに付随しており、ここでＩは識別行列（identity matrix）であり、パラメタαによって分析のスケールが決定される。位置ｘにおけるサンプル点密度推定器は、
【００７５】
【数４】
によって定められ、ここでＤはｘとｘｉとの間のマハラノビス距離を表し、
【００７６】
【数５】
である。
【００７７】
ｘにおける可変帯域幅平均シフトベクトルは、
【００７８】
【数６】
によって与えられ、ここでＨｈは、データ依存の重み付けωｊ(ｘ)によって重み付けされた帯域幅マトリクスの調和平均を表し、
【００７９】
【数７】
である。
【００８０】
目下の位置ｘにおいて計算されたデータ依存の重み付けは、
【００８１】
【数８】
と表され、ここでこれは
【００８２】
【数９】
を満たすことに注意されたい。
【００８３】
点ｘ＋ｍ(ｘ)に対応する密度はつねに、ｘに対応する密度よりも大きいかまたは等しいことを示すことができる。したがって平均シフトを使用して目下の位置を繰り返して更新することによって、基礎密度（underlying density）の停留点（stationary point）に収束するヒルクライミングプロシージャが得られる。
【００８４】
ＶＢＤＦ推定器は、いくつかのスケールで適応形平均シフトプロシージャを繰り返して適用することによって最重要モードを探し出す。点ｘｉの広がりに対してパラメタαを大きく選択することにより、大きなスケールから開始される。この場合、密度面（density surface）は単峰形（unimodal）であり、したがって決定されるモードは、大域的に最も密度の高い領域に対応することになる。パラメタαの値を低減しまた前のスケールにおいて決定されたモードから平均シフトの繰り返しを開始して、このプロシージャを繰り返す。最終ステップに対して各点に関連する帯域幅行列は、共分散行列に等しく、すなわちＨｉ ＝ Ｃｉである。
【００８５】
ＶＢＤＦ推定器は、複数のソースモデルを処理することの可能な情報フュージョンに対する強力なツールである。局所的な近傍における複数の点は複数の動きを示すため、これは動き推定に対して重要である。最重要モードは、最も関係がある動きに対応する。
【００８６】
ｎ個のモデルＭ０，Ｍ１…Ｍｎがあるとする。位置がｘｉｊ，ｉ＝１…ｎ，ｊ＝１…ｃによって表されるｃ個のコンポーネントを画像毎に維持する。新たな画像が利用可能になると、コンポーネント毎およびモデル毎に位置および不確かさが推定される。このステップはいくつかの手法、例えば、画像の相関、時空間エネルギー（spatio-temporal energy）の正則化（regularization）または空間傾斜（spatial gradient）に基づく手法を使用して行うことができる。ＶＢＤＦ法を使用した場合、得られる結果はコンポーネント毎の動き推定ｘｉｊおよびその不確かさＣｉｊである。したがってｘｉｊは、モデルｉについてコンポーネントｊの局所的推定を表す。共分散行列のスケールも、適合する残差から推定される。各制御点が閉塞（occlude）しており、したがって閉塞が制御点レベルで処理される場合には、このことによって共分散行列のサイズが増大する。
【００８７】
ＶＢＤＦロバストフュージョン法は、目下のフレームにおいてコンポーネントｊに対して最も関係する位置ｘjを決定するために適用される。スケールにわたるモードトラッキングによって、
【００８８】
【数１０】
が得られ、ここで重み付けωｉは（２６）のように定義され、
【００８９】
【数１１】
である。ここでこれらが
【００９０】
【数１２】
を満たすことに注意されたい。
【００９１】
点ｘ＋ｍ(ｘ)に対応する密度はつねに、ｘに対応する密度よりも大きいかまたは等しいことを示すことができる。したがって平均シフトを使用して目下の位置を繰り返して更新することによって、基礎密度の停留点に収束するヒルクライミングプロシージャが得られる。
【００９２】
ＶＢＤＦ推定器は、いくつかのスケールで適応形平均シフトプロシージャを繰り返して適用することによって最重要モードを探し出す。点ｘｉの広がりに対してパラメタαを大きく選択することにより、大きなスケールから開始される。この場合、密度面は単峰形であり、したがって決定されるモードは、大域的に最も密度の高い領域に対応することになる。パラメタαの値を低減しまた前のスケールにおいて決定されたモードから平均シフトの繰り返しを開始して、このプロシージャを繰り返す。最終ステップに対して各点に関連する帯域幅行列は、共分散行列に等しく、すなわちＨｉ ＝ Ｃｉである。
【００９３】
ＶＢＤＦ推定器は、複数のソースモデルを処理することの可能な情報フュージョンに対する強力なツールである。局所的な近傍における複数の点は複数の動きを示すため、これは動き推定に対して重要である。最重要モードは、最も関係がある動きに対応する。
【００９４】
以下では、複数の制御点モデルをトラッキングする例を本発明にしたがって説明する。ここではｎ個のモデルＭ０，Ｍ１…Ｍｎがあるとする。位置がｘｉｊ，ｉ＝１…ｎ，ｊ＝１…ｃによって表されるｃ個の制御点を画像毎に維持する。新たな画像が利用可能になると、制御毎およびモデル毎に位置および不確かさが推定される。このステップはいくつかの手法、例えば、画像の相関、時空間エネルギーの正則化または空間傾斜に基づく手法を使用して行うことができる。ＶＢＤＦ法を使用した場合、得られる結果は制御点毎の動き推定ｘｉｊおよびその不確かさＣｉｊである。したがってｘｉｊは、モデルｉについてコンポーネントｊの局所的推定を表す。共分散行列のスケールも、適合する残差から推定される。各制御点が閉塞（occlude）しており、したがって閉塞が制御点レベルで処理される場合には、このことによって共分散行列のサイズが増大する。
【００９５】
ＶＢＤＦロバストフュージョン法は、目下のフレームにおいてコンポーネントｊに対して最も関係する位置ｘjを決定するために適用される。スケールにわたるモードトラッキングによって、
【００９６】
【数１３】
が得られ、ここで重み付けωｉは（２８）のように定義され、
【００９７】
【数１４】
である。
【００９８】
各制御点の計算に続いて、推定の共分散行列によって与えられる重み付けによって、重み付けされた矩形のフィッティングが行われる。画像パッチは、４つのパラメタによって定義される相似変換によって結びつけられる。ダイナミックに変わる制御点位置ｘの相似変換は、式
【００９９】
【数１５】
によって特徴付けられ、ここでｔｘ，ｔｙは平行移動のパラメタであり、ａ，ｂは２Ｄのの回転とスケーリングのパラメタである。
【０１００】
最小化された判定条件は、基準位置ｘ0jと、推定されたｘj（目下のフレームにおけるｊ番目の制御点位置）との間のマハラノビス距離の総和である。
【０１０１】
【数１６】
【０１０２】
最小化は、標準的な重み付け最小自乗法を介して行われる。共分散行列は制御点毎に使用されるため、不確かさの大きな点の影響が低減されることに注意されたい。
【０１０３】
トラッキングされた制御点に矩形をフィッティングした後、この矩形内でダイナミックに変化するコンポーネントの候補が均一的に再サンプルされる。ここで仮定するのは、この矩形に対する各制御点の相対位置があまり変化しないことである。ある制御点のオプティカルフローによって計算されるトラック位置と再サンプル位置との距離が、許容される閾値よりも大きい場合、このトラック位置は、異常値とみなされ、再サンプルされた点に置き換えられる。十分な個数の制御点が小さな残差を有する場合、目下の画像がモデルの集合に加えられる。複数のモデルと目下のフレームとの間の残差中央値が、あらかじめ定められた閾値Ｔｈと比較される。
【０１０４】
フレームｉにおいてコンポーネントｊが位置ｘijを有するモデルＭ０，Ｍ１…Ｍｎの集合が与えられる場合、ここのオブジェクトトラッキングアルゴリズムは、以下のステップによってまとめることができる。すなわち、
【０１０５】
【数１７】
【０１０６】
ここで提案したマルチテンプレートフレームワークは、形状トラッキングのコンテキストにおいて直接適用可能である。トラッキングされた点が、スプラインによってモデル化された形状の制御点を表す場合、複数位置推定のロバストフュージョンを使用することによって、形状の位置推定の信頼性が増大する。形状空間（shape space）が、学習された部分空間制約によって制約される場合、補正も少なくなる。輪郭が利用可能な場合、トラッキングに使用されるモデルテンプレートは、形状間の距離に基づいてモデル集合からオンラインで選択することができる。
【０１０７】
図６ａおよび６ｂは左心室の心エコー図の画像を示しており、ここでは心内膜壁の運動が本発明にしたがってトラッキングされている。画像６０２，６１０に対するフレーム１は、壁部の初期輪郭を示しており、これはこの画像においてドットで示されている。後続のフレーム６０４，６０６，６０８，６１２，６１４，６１６は、この壁部の運動が時間についてどのようにトラッキングされるかを示している。各ドットに対する測定は本発明にしたがって行われる。
【０１０８】
例
この例において、心エコー図画像において手動でトレースされた心臓輪郭がトレーニングセットとして使用される。先端の４室（開いた輪郭）および胸骨傍短軸（閉じた輪郭）の両側から見た図が、図７ａ〜７ｄおよび図８ａ〜８ｄに示したように検査される。輪郭毎にランドマーク点が付される。図７ｂおよび図８ｂは、一緒に引かれたトレーニング輪郭のセットを示している。各輪郭がベクトルであり、このベクトルがそのコンポーネントとして、順序づけされたランドマーク点の座標を有する箇所に固有分析アプローチを適用する（開いた輪郭に対して３４個のコンポーネントまた閉じた輪郭に対して３６個のコンポーネント）。この場合、列がトレーニングベクトルである行列にＰＣＡを実行する。つぎに固有値を使用して、部分空間におけるモデル共分散行列として対角行列を形成する。測定値共分散行列は、さまざまなシナリオをテストするために調整される。
【０１０９】
図７ａおよび８ａはテスト輪郭を示しており、ここでは実線の曲線７０４，８０２は、グラウンドトルースであり、破線の曲線７０２，８０４は、ノイズを有する測定値である。比較的不確定な点を一層調整することによって、目下のノイズを有する輪郭に最も近い、固有輪郭部分空間における輪郭を探すことが望ましい。図７ｃおよび図８ｃは、図７ａおよび図８ａの輪郭に垂直な投影を行った場合の結果を示している（等方性の共分散は、ランドマーク点のまわりの小さな円７１２，８１２によって示されている）。図７ｄおよび８ｄは結果を示している。この結果は、はるかにグラウンドトルースに近いことがわかる。トレーニングデータがテストデータとはかなり異なるため、これは完全な一致にはならない。このため、テスト輪郭におけるいくつかの小さな形状の変形は、部分空間において実現できないことがある。
【０１１０】
本発明によれば、輪郭は視覚化ツールを使用して診断に役立てることができる。医師が心エコー図から心臓疾患をより簡単に診断できるようにするため、ＬＶ輪郭は、リアルタイムにカラーベクトルで表示される。ベクトルの長さが変位の大きさを表すのに対し、ベクトルの方向はＬＶ壁輪郭における点の運動方向を表す。例えば、オレンジ色は、運動が収縮である場合を示すために使用され、暗い青は、運動が拡張である場合を示すために使用される。図９は表示結果の一例を示している。ＬＶ輪郭における点の運動が追跡されている。この運動は時間および空間領域の両方においてガウシアンによって平滑化される。医師がよりわかりやすくするためにベクトルの長さは３倍される。
【０１１１】
この視覚化法により、医師は心内膜の各セグメントの運動の大きさおよび方向を用意に見ることができる。この視覚化方法と大域的な運動補償とを組み合わせることによって、医師はＬＶのすべてのセグメントにおいて収縮方向および大きさを簡単に見ることができる。虚血性の領域または他の以上な領域を直接かつ人間の目で識別することができる。
【０１１２】
健康診断中に機能の感度を決定する方法に対する実施形態を説明したが、上に教示したことを考え合わせれば、当業者は変更および変形をなし得ることに注意されたい。したがって開示した本発明の特定の実施形態において種々の変更を行うことができ、またこれらの変更が添付の請求項によって定められる本発明の範囲内にあることを理解されたい。ここまで特許法によって要求される詳細で本発明を説明したが、請求され特許法によって保護されるものは添付の請求項に記載されている。
【図面の簡単な説明】
【０１１３】
【図１】本発明にしたがって形状追跡法を実現する例示的なシステムのブロック図である。
【図２】音響的なドロップアウトのエリアと、推定された局所的な壁運動の不確かさとを示す心エコー図画像である。
【図３】制御点の位置の不確かさが極めて非等方であり不均質である心内膜輪郭を示す左心室の心エコー図である。
【図４】相異なるα値で強く適合化したＰＣＡ（ＳＡ−ＰＣＡ）モデルおよびインクリメント式ＰＣＡモデルの例を示す図である。
【図５】心内膜壁が初期化され（ａ）、ＩＰＣＡ（ｂ）およびＳＡ−ＰＣＡ（ｃ）を使用して追跡される左心室の心エコー図の画像である。
【図６】本発明にしたがって心内膜壁の動きを追跡した左心室の心エコー図の画像である。
【図７】本発明にしたがって先端の４室側から見た左心室の心内膜壁の心エコー図画像を表すテスト輪郭である。
【図８】本発明にしたがって胸骨傍短軸側から見た左心室の心内膜壁の心エコー図画像を表すテスト輪郭である。
【図９】本発明による画像視覚化方法を例示する心エコー図画像である。
【図１０】本発明の１実施形態による、マルチモデル制御点ベースの追跡器を示す図である。
【技術分野】
【０００１】
関連する出願との相互参照
本出願は、２００３年３月７日付け提出の米国暫定特許出願番号第60/452,669号、２００３年５月２７日付け提出の米国暫定特許出願番号第60/473,425号、および２００３年１０月３日付け提出の米国暫定特許出願番号第60/508,367号に優先権を主張するものであり、ここではこれらの全体を参考のために組み込む。
【０００２】
発明の分野
本発明は、運動中の形状を追跡する方法に関し、また殊に異分散ノイズ（heteroscedastic noise）がある場合に線形制約（linear constraint）を有する形状を追跡する方法に関する。
【０００３】
発明の背景
ほとんどの視覚的追跡アプリケーションに対して測定データは、不確定であり、また欠落していることもある。すなわち、画像はノイズおよび歪みを伴って記録されるのであり、その一方で閉塞（occlusion）により、関心領域の一部が観測できないようになることがあるのである。不確かさは、大域的に均一になり得るが、ほとんどの実世界のシナリオにおいて、これは本質的に異分散である。すなわち非等方であり不均質である。よい例が心エコー図（超音波心臓データ）である。超音波は、反射アーチファクトの影響、例えば薄膜から得られる鏡面反射体の影響を受けやすい。「見る方向」が１つであるため、鏡面構造体の垂直な面により、強いエコーが形成されるが、傾斜したまた「軸からずれた」表面により、弱いエコーしか形成されないか、またはまったくエコーが形成されないことがある（音響的な「ドロップアウト」）。心エコー図に対し、このようなドロップアウトは、超音波ビームに対して組織表面が平行である心臓の領域において発生し得る。
【０００４】
利用し易く、比較的コストがかからず、また非侵襲性であることに起因して心臓超音波画像は、心機能を評価するために広く使用されている。殊に心室の運動の分析は、虚血および梗塞の程度を評価するのに有効である。心内膜壁の検出またはセグメンテーションは、左心室の収縮および弾性の定量化に対する第１ステップである。いくつかの既存の方法の例に含まれるのは、ピクセルベースのセグメンテーション／クラスタリングアプローチ（例えばカラーキネシス（Color Kinesis））、オプティカルフローの種々の変形形態、変形可能なテンプレート（deformable template）およびマルコフ確率過程／フィールド、およびアクティブ輪郭（active contour）／スネーク（snake）である。いくつかの方法は２次元、３次元または４次元（３Ｄ＋時間）空間で使用される。
【０００５】
しかしながら既存のセグメンテーションおよび検出法は、心内膜壁の正確な局所的運動を再生しようとはしておらず、また多くの場合、壁に沿った動きの成分が無視される。このような単純化された処理は、目下の輪郭の法線に沿ってのみサーチする輪郭追跡器でも使用されている。これは局所的な壁の異常を検出するのには不適切である。それは異常な左心室の局所的な動きは、輪郭の法線から外れることが多く、また（音波検査者の手の動きまたは患者の呼吸による動きに起因する）平行移動または回転運動のような大域的な動きも、輪郭において法線を外れた局所的な動きを発生させることはいうまでない。望ましいのは、心内膜壁の大域的な形状と、その局所的な動きとを追跡して、局所的な壁の動きの異常を検出することである。この情報は、虚血および梗塞のさらなる診断に使用可能である。
【０００６】
一般的に共分散は、背後にある異分散ノイズを反映する画像フィーチャ（image feature）またはフロー推定（flow estimate）に割り当てることができる。全体的にノイズレベルが低くデータがきれいである場合、異分散的な性質は無視することができ、また大域的な不確かさを、局所的な推定の代わりに使用することができる。しかしながら極めてノイズの多い入力に対して、殊に空間的に変化する構造的なノイズを有する入力に対して、局所共分散行列において符号化される情報は、背後にある画像構造またはオブジェクトを高い信頼性でロバストに推定するのに重要になる。
【０００７】
追跡のフレームワークにおいてモデル制約を課すことはふつうに行われている。例としては、小球体またははパラメタ化された楕円などの簡単なモデル、および判別テンプレート（discriminative template）のような複雑なモデルが含まれる。ほとんどの実践的なケースにおいて部分空間モデルは形状追跡に適している。それは、主要な形状変化を捕捉するモードの数は限られており、また通常、この形状を表すのに使用されるフィーチャコンポーネントの元々の数よりも格段に小さいからである。さらにＰＣＡ（Principal Component Analysis）ベースの固有形状部分空間によって、任意の複雑な形状変化を捕捉することができる。ここでこれは元々の空間において、極めて単純なパラメトリックモデルを有しているとしても、極めて非線形である。
【０００８】
測定ベクトルが異分散ノイズによって影響を受ける場合、制約的な部分空間への正射影は、そろっていないだけではなく、情報のロスという意味で極めて多くのダメージを受ける。ノイズが等方的であると共に均一である特殊な場合にのみ、これをそろえることができるのである。
【０００９】
しかしながら部分空間制約追跡への既存の取り組みのほとんどは、測定値における異分散ノイズを考慮していない。「点分散モデル（Point Distribution Model）」または「アクティブ形状モデル（Active Shape Model）」において、ＰＣＡベースの部分空間形状モデルは、ランドマーク点（landmark point）対応関係を有するトレーニング形状（training shape）に基づいて導出される。結果的に得られる固有形状の部分空間によって、トレーニングデータセットにおける最も重要な変化が捕捉される。検出時にモデルは摂動されて、候補の位置におけるテスト画像との一致のための合成画像が形成される。しかしながら測定ノイズは、この処理においてモデル化されていなかった。
【００１０】
異分散ノイズ特性が利用可能である場合であっても、これらは通例、部分空間モデルフィッティング中に廃棄されていた。例えば、測定においてフルの共分散行列が捕捉される公知の１アプローチにおいて、かなりその場限りの特別な閾値判定が適用されて、測定値の平均値が、モデルの共分散によって定められる超楕円体の制約に限定される。この操作は、測定ノイズには依存しない。
【００１１】
別の公知のアプローチは、２ステップアプローチを適用して、カルマンフィルタリングフレームワークにおいて形状空間制約を課している。この形状空間は、線形変換されたアフィン部分空間または固有部分空間である。しかしながら形状空間への投影は正射影であり、測定の異分散のノイズを考慮していない。したがってこのアプローチは投影中に情報をロスすることになる。
【００１２】
別の公知のアプローチではガウス分布を使用して関心対象のオブジェクト（この場合、顔である）の外観を適応的にモデリングしている。ここでこれはＥＭアルゴリズムを使用して学習される。本発明と同様に、局所的な不確かさは共分散行列に捕捉される。異なるのは、本発明では殊に、非等方性の不確かさが存在する場合に投影にわたる交わりの決定的な選択および部分空間モデルの制約を検討することである。
【００１３】
別の公知のアプローチでは、オプティカルフロー推定中に暗黙的に部分空間制約使用しており、ここではフローの不確かさも利用される。別の応用に対する別のフレームワークにおいてではあるが、本発明によって認識されるのは、「部分空間投影処理においてフローベクトルの信頼性が高ければ高いほど、大きな影響がある」ことである。
【００１４】
別の公知のアプローチでは、基本行列および楕円のフィッティングに対して異分散回帰を適用している。このフィッティングは、元々の空間において、パラメータ化されたモデルで行われる。本発明では、形状の変化のパラメータ化（極めて複雑であり、高度に非線形である）は回避される。その代わりに本発明では、例えばＰＣＡを介して、部分空間線形確率モデルを形成し、またフィッティングしたデータの平均および共分散の両方について閉じた形の解決手段を得る。
【００１５】
またＭ推定器またはRANSACによるロバストモデルマッチングも適用されており、モデルの異常値であるデータ成分の影響が制限または除去される。ここでも（空間的または時間的に）局所的に変化する不確かさは、これらのフレームワークにおいて使用されない。
【００１６】
関連する別のアプローチに含まれるのは、データインピュテーション（data imputation）、すなわち、欠落したデータを信頼性のある値で「埋める」ことである。この領域における取り組みは、音声認識、医用画像分析、および社会科学に対して広い応用を有する統計学において定着している。しかしながらデータインピュテーションの定式化の問題は、ふつう０−１の利用可能性を前提としている。すなわちデータ成分は欠落しているか、利用可能であるかのいずれかを前提としているのである。部分空間モデルの制約と、形状ダイナミックおよび測定データノイズの異分散の性質についての情報および形状ダイナミックスについての情報とを融合するための統一的なフレームワークが必要である。
【００１７】
発明の要約
本発明は、運動中のオブジェクトの大域的な形状を追跡するシステムおよび方法に関する。ここでは大域的な形状の初期輪郭に沿って１つ以上の制御点が定められる。オブジェクトが運動に伴って上記の１つ以上の制御点の各々が追跡される。運動中の制御点の位置の不確かさが推定される。この不確かさを表す１形態は共分散行列である。部分空間形状制約モデルを使用する場合、この不確かさは、非正射影射影および／または情報フュージョンを使用して利用される。後続の各輪郭が表示される。
【００１８】
本発明はまたオブジェクトの形状の運動を視覚的に追跡するためのシステムに関する。１つ以上の第１のカラーベクトルが生成されて、この形状の輪郭に沿った制御点の収縮が表される。１つ以上の第２のカラーベクトルが生成されて、この形状の輪郭に沿った制御点の拡張が表される。上記の第１および第２のカラーベクトルは周期的に表示され、この形状の動きが表示される。
【００１９】
本発明はまた運動中のオブジェクトの大域的な形状を追跡する方法に関する。１つ以上の制御点を大域的な形状に沿って定める。オブジェクトの運動に伴って１つ以上の制御点の各々を追跡する。各制御点に対して複数の外観モデルを形成する。この形状を追跡する各モデルによって形成される動きベクトルを組み合わる。
【００２０】
図面の簡単な説明
添付の図面に関連して以下に本発明の有利な実施形態をより詳細に説明する。ここでは類似の参照番号は類似の要素を示す。
【００２１】
図１は、本発明にしたがって形状追跡法を実現する例示的なシステムのブロック図であり、
図２は、音響的なドロップアウトのエリアと、推定された局所的な壁運動の不確かさとを示す心エコー図画像であり、
図３は、制御点の位置の不確かさが極めて非等方であり不均質である心内膜輪郭を示す左心室の心エコー図であり、
図４は、相異なるα値で強く適合化したＰＣＡ（ＳＡ（strongly adapted）−ＰＣＡ）モデルおよびインクリメント式ＰＣＡモデルの例を示しており、
図５ａ〜５ｃは、左心室の心エコー図の画像を示しており、ここでは心内膜壁が初期化され（ａ）、ＩＰＣＡ（ｂ）およびＳＡ−ＰＣＡ（ｃ）を使用して追跡され、
図６ａおよび６ｂは左心室の心エコー図の画像を示しており、ここでは心内膜壁の動きは本発明にしたがって追跡され、
図７ａ〜７ｄはテスト輪郭を示しており、これは本発明にしたがい、先端の４室側から見た左心室の心内膜壁の心エコー図画像を表しており、
図８ａ〜８ｄはテスト輪郭を示しており、これは本発明にしたがい、胸骨傍短軸（parasternal short axis）側から見た左心室の心内膜壁の心エコー図画像を表しており、
図９は、本発明による画像視覚化方法を例示する心エコー図画像を示しており、
図１０は、本発明の１実施形態による、マルチモデル制御点ベースの追跡器を示している。
【００２２】
詳細な説明
本発明は、異分散ノイズがある場合に線形制約によって形状を追跡する方法に関する。このような方法が利用される１例では、心筋層壁の大域的な形状と、その局所的な動きとが追跡されて、心臓における局所的な壁運動の異常が検出される。この方法は、心臓の心内膜壁または心外膜壁を追跡するためにも使用可能である。本発明が、頭部の運動、顔の特徴、手の運動または別の身体運動などの人間の特徴部分の運動を識別するなどのように形状追跡が有効な別の応用にも使用でき、しかしながらこれらに限定されるないことを当業者は理解されたい
本発明はまた２次元、３次元および４次元（３Ｄ＋時間）における心臓、肺または徐々に成長する腫瘍などの解剖学的構造の医学的な分析にも使用可能である。
【００２３】
本発明を説明するため、左心室の心内膜壁を追跡する例を述べる。図１は、本発明による左心室心内膜壁形状追跡方法を使用する心エコーシステムの例示的なアーキテクチャを示している。超音波トランスデューサのような医用センサ１０２が使用されて患者に検査が行われる。センサ１０２は、所定の健康診断に見合った医療用測定値を得るために使用される。例えば心臓に問題を経験した患者は心エコー検査を受けて、これが所定の心疾患の診断に役立てられる。超音波システムにより、さまざまな視点から心臓の２次元、３次元および４次元（３Ｄ＋時間）の画像が供給される。
【００２４】
センサ１０２から得られた情報は、ワークステーションまたはパーソナルコンピュータとすることの可能なプロセッサ１０４に伝達される。プロセッサ１０４は、センサデータを画像に変換し、これがディスプレイ１０８に伝達される。またディスプレイ１０８は、この画像に関連する別の図形情報または情報のテーブルも伝達することができる。本発明によれば、プロセッサ１０４には、心内膜壁の初期輪郭を表すデータも供給されている。このデータは、医師または音波検査者などのユーザによって手動で供給されるか、またはプロセッサ１０４によって自動的に供給される。輪郭は個別の点の列からなっており、その運動はプロセッサ１０４によって追跡され、ディスプレイ１０８に示される。個別の点がどのように追跡されるかについての詳細を以下詳述する。
【００２５】
医用センサ１０２からのデータに加えてプロセッサ１０４は別のデータ入力を受け取ることができる。例えばこのプロセッサは、プロセッサ１０４に付随したデータベース１０６からデータを受け取ることができる。このようなデータに含まれ得るのは、心内膜壁に対する輪郭形状の候補を表す部分空間モデルである。これらの部分空間モデルは、複数の患者の代表である左心室の画像か、または統計的な情報に基づいてコンピュータで生成した輪郭形状のモデルとすることが可能である。プロセッサ１０４は、ベイズカーネルマッチングまたはオプティカルフローベース法などの既知のアプローチを使用して、輪郭形状の個々の点を追跡する。追跡中に累積される誤差は、マルチテンプレートアダプティブマッチングフレームワークを使用して修正される。追跡の不確かさは各点において共分散行列の形で表され、これは引き続いて、非正射影の投影を使用することにより、部分空間形状制約によって十分に活用される。
【００２６】
図２は典型的な心エコー図画像を示している。音響的なドロップアウトを有する左心室の心内膜壁の部分は、実線の楕円２０８によって示されている。局所的な壁運動の推定は、点線の楕円２０２，２０４によって示されている。音響的なドロップアウトがあるため、心内膜壁はこの画像においてつねに最も強いエッジになるわけではない。本発明によれば、部分空間モデルの制約と、形状ダイナミックおよび測定ノイズの異分散的な性質についての情報および形状ダイナミクスについての情報との融合に対して、統合化されたフレームワークが使用される。この部分空間モデルは、所定の部分空間分布、例えばガウス分布、または単純な部分空間の制約、例えば固有空間モデルの形を取り得る。
【００２７】
本発明では心内膜壁を表す輪郭上の個々の制御点が追跡される。この追跡は、ベイズカーネルマッチングアプローチまたはフローベースのアプローチを使用することによって行うことができる。ベイズカーネルマッチングアプローチの例は、共同発明者Dorin Comaniciuによって著されたBayesian Kernel Tracking, Annual Conf. of the German Society for Pattern Recognition （DAGM'02）, Zurich, Switzerland, 438-445, 2002なる題名の論文に記載されており、その全体を参考のために組み込む。個々の点を追跡するオプティカルフローベースのアプローチの例は、"Density Estimation-Based Information Fusion for Multiple Motion Computation"なる名称の同時係属出願の明細書第10/681,02に記載されており、その全体を参考のために組み込む。本発明によれば、追跡中の誤差の累積は、心運動の周期的な性質を利用し、マルチテンプレートアダプティブマッチングフレームワークを使用することによって修正される。
【００２８】
追跡の不確かさは各点において共分散行列によって表され、これは続いて、非正射影の投影を使用することにより、部分空間形状制約によって十分に利用される。
【００２９】
追跡フレームワークは、イメージ列にわたる２つのステップの繰り返し処理である。制御点を有する初期輪郭は第１のイメージに（手動または自動で）引かれ、つぎに後続の画像毎に順番に各制御点がまず独立に追跡される。ここでは非等方性の不確かさも記録される。第２のステップとして新しい輪郭が、非正射影の投影を使用して、実現可能な部分空間に射影される。この実現可能な部分空間は、トレーニング輪郭に基づいて学習され、また追跡器に利用可能な初期輪郭を使用して都度のケースに適合される。また初期化の信頼度も考慮される（すなわち手動の初期化は、完全自動化の初期化よりも信頼性が高い）。
【００３０】
上に示したように、追跡処理中に複数のテンプレートが使用される。複数のテンプレートを使用することにより、結果的に目下のケースに対する形状の統計がより正確に表現されることになる。ベイズカーネルマッチングアプローチでは、第１のテンプレートは、初期化されたフレームから選択される。後続のテンプレートは、それが既存のテンプレートとは十分に異なっている共に、位置決めのために十分な情報を与える場合に追加される。ここでこれは、それ自体とのカーネルマッチングによって測定される。より多くの情報を与えるパッチは、それ自体とマッチングした場合に、より高い信頼性を有するパッチである。
【００３１】
１つ以上のテンプレートを使用することの決定は、心臓の運動が周期的であり、したがって１周期内の異なる外観のパターンが、後続の周期においても再度現れるという観察に基づく。複数のテンプレートとマッチングするため、マッチングはテンプレート毎に行われ、最良のマッチングを有するものが選択される。または計算を節約するため、前にマッチしたテンプレートの近くのテンプレートだけをマッチングすることができる。ここでの周期的な運動が使用される。
【００３２】
位置決めのマッチングについての不確かさは、最適な位置の近くで尤度マップ（likelihood map）によって計算され、ここでこれはカーネルマッチング処理を使用して同様に計算される。この尤度面（likelihood surface）はつぎに共分散行列を推定するために使用される。ここでこれは例えば、ヘッセの行列の重み付けされた推定の逆行列または２次元のガウス分布とのフィッティングによって行われる。
【００３３】
フローベースのアプローチの場合、通例のオプティカルフローの実現では隣りのフレームだけを使用しているため、長いシーケンスにわたる追跡処理は、誤差の累積およびドリフトに影響されやすくなっている。ここでは複数のテンプレートがフロー計算に使用されており、カーネルマッチングに対するのと同じである。新しいテンプレートは、フローの不確かさが大きい場合にはいつでもどこでも追加されるの対して、局所的な傾斜は既存のテンプレートから区別可能である。
【００３４】
制御点毎に新たな位置が得られた後、つぎのステップでは、人間の心臓の「正しい」形状変化を捉える統計的なモデルによって、全体的な形状を制約する。点分散モデルとしても知られているＰＣＡベースの形状モデルまたはアクティブ形状モデルが使用される。
【００３５】
心臓輪郭にわたる不確かさは、均一でもなく、すなわちエリアによっては、例えば信号のドロップアウトに起因して他のエリアよりも悪く、また等方性でもない、例えば、エッジに沿った位置決めの方が傾斜方向に沿った位置決めよりも悪いのである。図３は、このような輪郭にわたる非等方性かつ不均一なノイズの例を示している。図からわかるように個々の点、例えば点３０２，３０４，３０６および３０８ははじめ、心内膜壁の輪郭に沿っていると識別されている。点毎に確かさの測定が行われ、また各点のまわりに楕円体が形成される。この惰円体は、この所定の点が正しい位置にある確かさのレベルを表している。楕円体が大きければ大きいほど、この所定の点の位置についての不確かさのレベルが高くなる。図３からわかるように点３０２は、それを取り囲む比較的小さい楕円体を有し、輪郭におけるこの点の位置についての確かさの程度が高いをことを示している。点３０６および３０８は、それらの周りに比較的大きな楕円体を有しており、これらの点の位置についての不確かさが高いことを示している。
【００３６】
新しい輪郭的にわたる一般的な不確かさのモデルを考慮することによって、最適な射影はもはや正射影ではなくなる。部分空間における最適な解は実際に、形状モデル部分空間における交差分布における最尤度の形状である。さらにこの部分空間に分布モデルがある場合、このような追加情報を無視する理由はない。以下では、非正射影の投影の仕方について詳しい分析を示す。また部分空間モデル分布が利用可能な場合、このようなモデルの情報と、不確かな入力との融合のさせ方を示す。
【００３７】
本発明は、部分空間モデルベースフュージョンアプローチを使用する。同じｎ次元変数ｘのノイズを伴った２つの測定値が与えられる場合、各々は多次元のガウス分布ｐ１およびｐ２によって特徴付けられ、ｘの最尤度推定は、２つの重心に対するマハラノビス距離の最小和を有する点である。すなわち、
ｘ* ＝ argmin ｄｍ
であり、ここで
ｄｍ ＝ (ｘ−ｘ１)ＴＣ１−１(ｘ−ｘ１)＋(ｘ−ｘ２)ＴＣ２−１(ｘ−ｘ２） （１）
である。
【００３８】
ｘについての導関数をとると、
ｘ* ＝ Ｃ(Ｃ１−１ｘ＋Ｃ２−１ｘ２)
Ｃ ＝ (Ｃ１−１＋Ｃ２−１) （２）
を得る。
【００３９】
これは、ｘのＢＬＵＥ（best linear unbiased estimate）として知られている。
【００４０】
上記の２つのガウシアンのうちの１つが次元ｐの部分空間にある、例えばＣ２が特異値を有するとする。Ｃ２ ＝ ＵΛＵＴの特異値分解（singular value decomposition）により、（ただしＵ ＝ [ｕ１，ｕ２，…，ｕｎ]，ここでｕｉは正規直交であり、Λ ＝ diag｛λ１，λ２，…，λｐ，０，…，０｝である）、式（１）におけるｘ２に対するマハラノビス距離をつぎの標準的な形に書き換える。すなわち、
【００４１】
【数１】
である。
【００４２】
λｉが０なると、Ｕ０Ｔｘ ＝ ０，ただしＵ０ ＝ [ｕｐ＋１，ｕｐ＋２，…，ｕｎ]でない限り、ｄｍ，２は無限大になる。（ここで一般性を失うことなく仮定したのは、部分空間が元々の空間の原点を通ることである。ｘ２が部分空間にあるため、Ｕ０Ｔｘ２ ＝ ０である。）
したがって式（３）においてゼロでないλｉに相応する項だけ残しておけばよい。すなわち、
【００４３】
【数２】
であり、ここで
Ｃ＋はＣの疑似逆行列（pseudoinverse）、Ｕｐ ＝ [ｕ１，ｕ２，…，ｕｐ]である。
【００４４】
さらにＵ０Ｔｘ ＝ ０であるため、ｘを別の形で表して、この制約を反映することができる。すなわち、１×ｐベクトルｙに対して、
ｘ ＝ ＵＵＴｘ ＝ Ｕ([Ｕｐ｜Ｕ０]Ｔｘ) ＝ Ｕ[ｙ｜０] ＝ Ｕｐｙ （５）
である。式（１）はいまや、特異値にかかわらずつぎの一般的な形をとる。すなわち、
ｄｍ ＝ (Ｕｐｙ−ｘ１)ＴＣ１−１(Ｕｐｙ−ｘ１)＋(Ｕｐｙ−ｘ２)ＴＣ２＋(Ｕｐｙ−ｘ２) （６）
である。
【００４５】
ｙについての導関数をとると、
ｙ* ＝ Ｃｙ＊ＵｐＴ(Ｃ１−１ｘ１＋Ｃ２＋ｘ２) （７）
Ｃｙ＊ ＝ [ＵｐＴ(Ｃ１−１＋Ｃ２＋)Ｕｐ]−１ （８）
ｘ* ＝ Ｕｐｙ* ＝ Ｃｘ＊(Ｃ１−１ｘ１＋Ｃ２＋ｘ２) （９）
Ｃｘ＊ ＝ ＵｐＣｙ＊ＵｐＴ （１０）
である。
【００４６】
式（７）は情報空間への融合を示しており、これに続いて部分空間への変換が行われる。式（９）は、元々の空間に戻す座標変換である。Ｃｘ＊およびＣｙ＊は、ｘ*およびｙ*に対する相応の共分散行列であることは示すことができる。
【００４７】
この解は、正しい逆行列の代わりに疑似逆行列を使用することによって式（２）単純に一般化したものではないことに注意されたい。これはｘ*が部分空間にあることを制約しない。
【００４８】
択一的には式（７）および（８）を
ｙ* ＝ (ＵｐＴＣ１−１Ｕｐ＋Λｐ−１)−１(ＵｐＴC１−１ｘ１＋Λｐ−１ｙ２) （１１）
と書き表すことができる。
【００４９】
ここでｙ２はＵｐによって張られた部分空間におけるｘ２の変換された座標であり、Λｐ ＝ diag{λ１，λ2，…，λｐ}である。
【００５０】
式（１１）は、２つの分布の部分空間におけるＢＬＵＥフュージョンとして示すこともでき、その１つはＮ(ｙ２，Λｐ)であり、別の１つは部分空間におけるＮ(ｘ１，Ｃ１)の「交点」（射影ではない！）Ｎ((ＵｐＴＣ１−１Ｕｐ)−１ＵｐＴＣ１−１ｘ１，(ＵｐＴC１−１Ｕｐ)−１)である。
【００５１】
トレーニングサンプルの大きなプールから学習した統計的形状モデルを使用して、所定の心臓から輪郭を導くことは、時として問題となることがある。本発明によれば、トレーニングサンプルは、目下の心臓の形状モデルを得るために使用され、一般的な心臓の形状モデルを得るために使用されるのではない。したがって、目下のケースに対して知られていることに向かって一般的なモデルを適合させようとする強いモチベーションがある。本発明によれば、患者の心臓の心内膜壁の（手動または自動検出による）初期輪郭は検証されて、既存のＰＣＡモデルに適合化される。
【００５２】
心内膜壁の実際の輪郭を決定してその運動を追跡する際に２つのアプローチが考慮される。すなわち、初期輪郭が決定的（deterministic）である１アプローチと、初期輪郭が不確かである別の１アプローチとが考慮されるのである（これは、初期輪郭が自動検出アルゴリズムからのものである場合に発生する。これはまた不確かさも発生させる）。
【００５３】
初期輪郭が１点である（決定的に確定している）と仮定される場合、ＳＡ−ＰＣＡ（strongly-adpted-PCA）モデルが使用されてこの点の運動が追跡される。ここで仮定されるのは、旧いＰＣＡモデル（現在のケースを除く）と、現在のケースに対して初期化された輪郭とが、一緒になって現在のケースの変化を表すが、相対的なエネルギー（すなわち代表的な出力）はそれぞれａおよび１−ａ、ただし０＜ａ＜１となることである。言い換えると、現在のケースの形状変化の１部分は、一般的なモデルによって表されるのに対して、残りは初期輪郭の方向に捉えられるとを仮定するのである。
【００５４】
ＰＣＡモデルは、その平均、固有値行列および固有ベクトル行列で表され、またそれぞれｘｍ，ΛおよびＵで表される。元々の完全な共分散行列Ｃが格納される場合（これは、元々の次元数（dimensionality）が許容されないほどに高すぎないケースである）、適合された平均および共分散行列は単に、寄与する２つのソースの重み付けされた和である。すなわち、
ｘｍ，ｎｅｗ ＝ αｘｍ＋(１−α)ｘ （１２）
Ｃｎｅｗ ＝ α(Ｃ＋(ｘｍ−ｘｍ，ｎｅｗ）(ｘｍ−ｘｍ，ｎｅｗ）Ｔ)
＋(1−α) (ｘｍ−ｘｍ，ｎｅｗ）(ｘｍ−ｘｍ，ｎｅｗ）Ｔ
＝ αＣ＋α(１−α) (ｘｍ−ｘｍ，ｎｅｗ）Ｔ （１３）
である。
【００５５】
固有分析をＣｎｅｗに実行して、新しい部分空間モデルを得ることができる。
【００５６】
Ｃが格納されておらず、｛ｘｍ, Λ, Ｕ｝だけが部分空間において利用可能な場合、直接的な代数的操作によって以下のようにして、適合化された固有分析結果｛ｘｍ，ｎｅｗ, Λｎｅｗ, Ｕｎｅｗ｝に到達することができる。すなわち、
初期輪郭ｘは、ｘＳ ＝ ＵＴｘｄ（ただしｘｄ ＝ ｘ−ｘｍ）の部分空間成分を、またｘｒ ＝ (ｘ−ｘｍ)−ＵｘＳｉの残差ベクトルを有する。ｘｒを正規化してノルム１を有するようにしたものをｘｒｕとする（またはｘｒのノルムがゼロの場合ゼロとする）。
【００５７】
結合されたエネルギーを表す適合された固有ベクトル行列はつぎの形を有する。すなわち、
Ｕｎｅｗ ＝ [Ｕ, ｘｒｕ]Ｒ （１４）
である。
【００５８】
ＲおよびΛｎｅｗは、以下の固有分析問題の解になる。すなわち、
【００５９】
【数３】
であり、ここでｅｒ ＝ ｘｒｕＴ(ｘ−ｘｍ)は残りのエネルギーである。
【００６０】
上記の式は、新しいデータと旧いデータとの間の調整可能なエネルギー比を有し、ＩＰＣＡよりも一般的である。αを点の総数に対するモデルにおける点の数の比と設定した場合、これらはＩＰＣＡに等しくなる。通例、これは１に極めて近くなる。それはトレーニングセットにおける輪郭の数はふつう大きいからである。αを比較的小さな値（例えば０．５）に設定した場合、ＰＣＡモデルは目下のケースに向かって強力に適合化される。このことがその名前の由来である。図４は、相異なるα値を有するＩＰＣＡおよびＳＡ−ＰＣＡの簡単な２次元の図を示している。
【００６１】
点４０４は目下のケースを表している。各ｘ４０２は、特定のモデルに相応するトレーニング点を表している。楕円４０６は、元々のモデル分布を表している。楕円４０８は、強く適合化されたＰＣＡモデルに相応するインクリメント式ＰＣＡモデルを示しており、ここではα ＝ ０．９９である。楕円４１０は、α ＝ ０．５の強く適合化されたＰＣＡモデルを示している。楕円４１２は、α ＝ ０．１の強く適合化されたＰＣＡモデルを示している。各楕円は、相応する分布の等しく９０％で確からしい輪郭を示している。
【００６２】
実際に目下の心臓から得られる輪郭は、一般的なトレーニングセットにおける輪郭よりも、同じ心臓の初期輪郭に似ていることのほうがはるかに多い。殊に目下の心臓がトレーニングセットで良好に表されない異常な形を有する場合はそうである。ここのシステムでは、αを０．５にセットした。これによって初期輪郭からの強い影響が可能である。すなわち、モデルエネルギーの５０％は、初期輪郭から得られる（９９個の例でトレーニングされたモデルによってＩＰＣＡが適合される場合、これはわずかに１％に減少する）。
【００６３】
図５ａ〜５ｃはＩＰＣＡとＳＡ−ＰＣＡとの比較を示している。胸骨傍の短軸側から見たものは、（凹形部分のある）異常な形を有し、これに対してトレーニングセットは円形が大多数を占めている。初期輪郭（図５ａ）においてではあるが極めてわずかな重み付け（＜０．０１％）で記録されているインクリメント式ＰＣＡモデルは、目下の形状の凹形の性質を捉えておらず、この輪郭を典型的な円形（図５ｂ）に制約している。この結果、実際にインクリメント式ステップのない旧いPCAモデルを使用して得られたのと同じになっている。α ＝ ０．５の適応形ＰＣＡモデルでは、真の境界により一層良好にフィットしている（図５ｃ）。
【００６４】
フュージョンアプローチによる非正射影の投影と、ＳＡ−ＰＣＡモデル適合化との間の微妙ではあるが重要な相互作用はつぎのようになっている。モデル平均および共分散を有するフュージョンは、部分空間にはあるがモデル分散とはあまりに離れている輪郭を濾波するために必要である。しかしながらこのような（正射影よりも）強い制約は、必然的に低い尤度または異常な心臓輪郭を変化させる。ＳＡ−ＰＣＡモデルは、与えられた初期輪郭において供給される追加情報を使用することにより、目下のケースの方向に向かって強力にモデルをシフトすることを介してこのジレンマの間を取り持つのである。
【００６５】
ダイナミック処理によって定められる予想およびノイズ測定値からの情報をカルマンフィルタによって融合する。形状追跡に適用する場合、付加的な大域的制約が必要であり、この制約によって全体的な形状を有利な範囲に安定させる。本発明では、統合化されたフュージョンフレームワークを利用して、部分空間モデルの制約をカルマンフィルタに組み込む。
【００６６】
カルマンフィルタに対して、測定値の更新の式はつぎの形である。すなわち、
ｘｋ＋１｜ｋ＋１ ＝ ｘｋ＋１｜ｋ＋Ｋ(ｚｋ＋１−Ｈｚｋ＋１｜ｋ) （１６）
であり、ここで
Ｋ ＝ Ｐｋ＋１｜ｋＨＴ(ＨＰｋ＋１｜ｋＨＴ＋Ｒ)−１ （１７）
Ｐｋ＋１｜ｋ＋１ ＝ (Ｉ−ＫＨ)Ｐｋ＋１｜ｋ （１８）
Ｐｋ＋１｜ｋ ＝ ＳＰｋ｜ｋＳＴ＋Ｑ （１９）
である。
【００６７】
ここでＰは推定誤差共分散であり、ｘｉ｜ｊは時間jにおいて与えられた状態の時間iにおける状態推定である。この測定モデルはｚｋ ＝ Ｈｘｋ＋ｒｋであり、ただしｒｋは共分散Ｒを有する測定値ノイズを表す。このシステム／処理モデルはｘｋ＋１ ＝ Ｓｘｋ＋ｑｋであり、ここでｑｋは共分散Ｑを有するシステムノイズを表す。
【００６８】
上記のアプローチを使用することにより、部分空間制約および異常分散ノイズを有するカルマンフィルタの更新の式は、
ｘｋ＋１｜ｋ＋１ ＝ Ｐｋ＋１｜ｋ＋１((Ｐｋ｜ｋ＋Ｑ)＋ｘｋ｜ｋ＋Ｒ−１ｚ＋Ｃ２＋ｘ２) （２０）
Ｐｋ＋１｜ｋ＋１ ＝ Ｕｐ[ＵｐＴ(((Ｐｋ｜ｋ＋Ｑ)＋＋Ｒ−１＋Ｃ２＋)Ｕｐ]−１ＵｐＴ （２１）
によって与えられる。ここで仮定するのは、上記のシステムノイズ共分散Ｑは部分空間において制約されることである。情報空間において利用可能なすべての知識を組み合わせる解の対称性を観察する。これらの式によって、システムダイナミック、部分空間制約およびノイズ情報の統合化されたフュージョンが得られる。これらは、追跡システムに影響を及ぼすさまざまな不確かさの完全な表現を表す。
【００６９】
本発明の択一的な１実施形態によれば、制御点ベースのオブジェクトの表現と、フレームにわたる統合モデルへのロバストな融合とに基づく追跡技法が使用される。制御点ベースのオブジェクト表現の集合は、異なる時点に取得されて維持される。制御点によって示唆される推定の運動は、ロバストに結合されてオブジェクトのつぎの位置が決定される。このオブジェクトの視覚的な追跡は時間にわたっていくつかのモデルを維持することによって達成される。結果的に得られるのは確率密度関数のパラメトリックでない表現であり、これはオブジェクトの外形を特徴付ける。追跡は、オプティカルフローによるその不確かさと、運動の推定とを各モデルから独立に得ることによって行われる。各制御点に対する最終的な推定は、可変帯域幅密度フュージョンＶＢＤＦ（Variable-Bandwidth Density-based Fusion）のようなロバストな融合手法を使用して計算される。ＶＢＤＦにより、変位密度関数（displacement density function）の最も重要なモードの位置が計算され、その間にその不確かさが考慮される。ＶＤＢＦプロシージャにより、複数のデータソースと、変位の推定における異常値とが管理される。閉塞は、大きな残差に対し、推定の不確かさを介して処理されるのがふつうである。モデルは、フローが独立して計算される複数の領域に分割される。残りの位置合わせ誤差は推定の共分散行列のスケールを計算するために使用される。したがって信頼のできない変の影響が低減される。
【００７０】
閉塞および外形の変化によるオブジェクト追跡の試みは本発明において、マルチモード制御点ベースのアプローチを介して処理される。２次元外形モデルに対して複数のモデルを維持することは、これを単峰形分布（unimodal distribution）に制約せず、またＶＤＢＦメカニズムにより、ロバストに複数の推定が統合されて、制御点毎に最も優位な運動が決定される。追跡中の変化をモデリングするため、時間にわたるオブジェクトの外形の複数のサンプルが維持される。図１０には、本発明による例示的なマルチモデルの制御点ベースの追跡器が示されている。最も上の３つのフレーム１００２，１００４および１００６は、モデルセットの目下のサンプルを示しており、重なり合った制御点の関連する集合を有する。制御点ベースのアプローチは、大域的な表現よりもロバストであり、したがって照明の変化および中断に対してさほど敏感でない。本発明の別の利点は、部分的な閉塞が制御点レベルでマッチングの尤度を分析することによって処理できることである。
【００７１】
各制御点は独立して処理され、その位置および共分散行列は現在の画像においてすべてのモデルテンプレートについて推定される。例えば、制御点の１つはグレイな矩形１０１４によって示され、各モデルについてのその位置および不確かさはＩｎｅｗ１００８によって示されている。ＶＢＤＦロバストフュージョンプロシージャが適用されて、関連した不確かさと共に最も優勢な運動（モード）が決定される。これはフレーム１０１０に示したとおりである。各制御点の推定位置における変化が、閉塞または外形の変化によることに注意されたい。
【００７２】
目下のフレームにおける制御点の位置はさらに、大域的なパラメトリックな運動モデルによって制約される。類似性変換モデルおよびそのパラメタは、各制御点における信頼度（confidence）を使用して推定される。したがって信頼性の高い制御点は、大域的な運動の推定により大きく貢献するのである。残差、基準の外観が比較的低い場合、現在のフレームがモデルセットに加えられる。閾値を選択して、オブジェクトが重要な閉塞を有する箇所に画像が加えられないようにする。モデルにおけるテンプレートの数は可変でも固定でもよい。テンプレートの数が固定の場合、所定のテンプレートが廃棄される（例えば、最も旧いテンプレートを廃棄する）方式を提供することが可能である。
【００７３】
ＶＢＤＦ推定器は、アダプティブカーネル帯域幅（adaptive kernel bandwidth）を有するノンパラメトリック密度推定（nonparametric density estimation）に基づいている。ＶＢＤＦ推定器は、密度関数の最重要モードの位置として定義される。モード計算は、マルチスケール最適化フレームにおける可変帯域幅平均シフト法（variable-bandwidth mean shift technique）に基づいている。
【００７４】
ｘｉ ∈Ｒｄ，ｉ＝１…ｎを利用可能なｄ次元推定値とし、ここで各々は共分散行列Ｃｉによって与えられる付随した不確かさを有する。密度関数の最重要モードは、マルチスケール方式で反復的に決定される。帯域幅行列Ｈｉ ＝ Ｃｉ＋α２Ｉは、各点ｘｉに付随しており、ここでＩは識別行列（identity matrix）であり、パラメタαによって分析のスケールが決定される。位置ｘにおけるサンプル点密度推定器は、
【００７５】
【数４】
によって定められ、ここでＤはｘとｘｉとの間のマハラノビス距離を表し、
【００７６】
【数５】
である。
【００７７】
ｘにおける可変帯域幅平均シフトベクトルは、
【００７８】
【数６】
によって与えられ、ここでＨｈは、データ依存の重み付けωｊ(ｘ)によって重み付けされた帯域幅マトリクスの調和平均を表し、
【００７９】
【数７】
である。
【００８０】
目下の位置ｘにおいて計算されたデータ依存の重み付けは、
【００８１】
【数８】
と表され、ここでこれは
【００８２】
【数９】
を満たすことに注意されたい。
【００８３】
点ｘ＋ｍ(ｘ)に対応する密度はつねに、ｘに対応する密度よりも大きいかまたは等しいことを示すことができる。したがって平均シフトを使用して目下の位置を繰り返して更新することによって、基礎密度（underlying density）の停留点（stationary point）に収束するヒルクライミングプロシージャが得られる。
【００８４】
ＶＢＤＦ推定器は、いくつかのスケールで適応形平均シフトプロシージャを繰り返して適用することによって最重要モードを探し出す。点ｘｉの広がりに対してパラメタαを大きく選択することにより、大きなスケールから開始される。この場合、密度面（density surface）は単峰形（unimodal）であり、したがって決定されるモードは、大域的に最も密度の高い領域に対応することになる。パラメタαの値を低減しまた前のスケールにおいて決定されたモードから平均シフトの繰り返しを開始して、このプロシージャを繰り返す。最終ステップに対して各点に関連する帯域幅行列は、共分散行列に等しく、すなわちＨｉ ＝ Ｃｉである。
【００８５】
ＶＢＤＦ推定器は、複数のソースモデルを処理することの可能な情報フュージョンに対する強力なツールである。局所的な近傍における複数の点は複数の動きを示すため、これは動き推定に対して重要である。最重要モードは、最も関係がある動きに対応する。
【００８６】
ｎ個のモデルＭ０，Ｍ１…Ｍｎがあるとする。位置がｘｉｊ，ｉ＝１…ｎ，ｊ＝１…ｃによって表されるｃ個のコンポーネントを画像毎に維持する。新たな画像が利用可能になると、コンポーネント毎およびモデル毎に位置および不確かさが推定される。このステップはいくつかの手法、例えば、画像の相関、時空間エネルギー（spatio-temporal energy）の正則化（regularization）または空間傾斜（spatial gradient）に基づく手法を使用して行うことができる。ＶＢＤＦ法を使用した場合、得られる結果はコンポーネント毎の動き推定ｘｉｊおよびその不確かさＣｉｊである。したがってｘｉｊは、モデルｉについてコンポーネントｊの局所的推定を表す。共分散行列のスケールも、適合する残差から推定される。各制御点が閉塞（occlude）しており、したがって閉塞が制御点レベルで処理される場合には、このことによって共分散行列のサイズが増大する。
【００８７】
ＶＢＤＦロバストフュージョン法は、目下のフレームにおいてコンポーネントｊに対して最も関係する位置ｘjを決定するために適用される。スケールにわたるモードトラッキングによって、
【００８８】
【数１０】
が得られ、ここで重み付けωｉは（２６）のように定義され、
【００８９】
【数１１】
である。ここでこれらが
【００９０】
【数１２】
を満たすことに注意されたい。
【００９１】
点ｘ＋ｍ(ｘ)に対応する密度はつねに、ｘに対応する密度よりも大きいかまたは等しいことを示すことができる。したがって平均シフトを使用して目下の位置を繰り返して更新することによって、基礎密度の停留点に収束するヒルクライミングプロシージャが得られる。
【００９２】
ＶＢＤＦ推定器は、いくつかのスケールで適応形平均シフトプロシージャを繰り返して適用することによって最重要モードを探し出す。点ｘｉの広がりに対してパラメタαを大きく選択することにより、大きなスケールから開始される。この場合、密度面は単峰形であり、したがって決定されるモードは、大域的に最も密度の高い領域に対応することになる。パラメタαの値を低減しまた前のスケールにおいて決定されたモードから平均シフトの繰り返しを開始して、このプロシージャを繰り返す。最終ステップに対して各点に関連する帯域幅行列は、共分散行列に等しく、すなわちＨｉ ＝ Ｃｉである。
【００９３】
ＶＢＤＦ推定器は、複数のソースモデルを処理することの可能な情報フュージョンに対する強力なツールである。局所的な近傍における複数の点は複数の動きを示すため、これは動き推定に対して重要である。最重要モードは、最も関係がある動きに対応する。
【００９４】
以下では、複数の制御点モデルをトラッキングする例を本発明にしたがって説明する。ここではｎ個のモデルＭ０，Ｍ１…Ｍｎがあるとする。位置がｘｉｊ，ｉ＝１…ｎ，ｊ＝１…ｃによって表されるｃ個の制御点を画像毎に維持する。新たな画像が利用可能になると、制御毎およびモデル毎に位置および不確かさが推定される。このステップはいくつかの手法、例えば、画像の相関、時空間エネルギーの正則化または空間傾斜に基づく手法を使用して行うことができる。ＶＢＤＦ法を使用した場合、得られる結果は制御点毎の動き推定ｘｉｊおよびその不確かさＣｉｊである。したがってｘｉｊは、モデルｉについてコンポーネントｊの局所的推定を表す。共分散行列のスケールも、適合する残差から推定される。各制御点が閉塞（occlude）しており、したがって閉塞が制御点レベルで処理される場合には、このことによって共分散行列のサイズが増大する。
【００９５】
ＶＢＤＦロバストフュージョン法は、目下のフレームにおいてコンポーネントｊに対して最も関係する位置ｘjを決定するために適用される。スケールにわたるモードトラッキングによって、
【００９６】
【数１３】
が得られ、ここで重み付けωｉは（２８）のように定義され、
【００９７】
【数１４】
である。
【００９８】
各制御点の計算に続いて、推定の共分散行列によって与えられる重み付けによって、重み付けされた矩形のフィッティングが行われる。画像パッチは、４つのパラメタによって定義される相似変換によって結びつけられる。ダイナミックに変わる制御点位置ｘの相似変換は、式
【００９９】
【数１５】
によって特徴付けられ、ここでｔｘ，ｔｙは平行移動のパラメタであり、ａ，ｂは２Ｄのの回転とスケーリングのパラメタである。
【０１００】
最小化された判定条件は、基準位置ｘ0jと、推定されたｘj（目下のフレームにおけるｊ番目の制御点位置）との間のマハラノビス距離の総和である。
【０１０１】
【数１６】
【０１０２】
最小化は、標準的な重み付け最小自乗法を介して行われる。共分散行列は制御点毎に使用されるため、不確かさの大きな点の影響が低減されることに注意されたい。
【０１０３】
トラッキングされた制御点に矩形をフィッティングした後、この矩形内でダイナミックに変化するコンポーネントの候補が均一的に再サンプルされる。ここで仮定するのは、この矩形に対する各制御点の相対位置があまり変化しないことである。ある制御点のオプティカルフローによって計算されるトラック位置と再サンプル位置との距離が、許容される閾値よりも大きい場合、このトラック位置は、異常値とみなされ、再サンプルされた点に置き換えられる。十分な個数の制御点が小さな残差を有する場合、目下の画像がモデルの集合に加えられる。複数のモデルと目下のフレームとの間の残差中央値が、あらかじめ定められた閾値Ｔｈと比較される。
【０１０４】
フレームｉにおいてコンポーネントｊが位置ｘijを有するモデルＭ０，Ｍ１…Ｍｎの集合が与えられる場合、ここのオブジェクトトラッキングアルゴリズムは、以下のステップによってまとめることができる。すなわち、
【０１０５】
【数１７】
【０１０６】
ここで提案したマルチテンプレートフレームワークは、形状トラッキングのコンテキストにおいて直接適用可能である。トラッキングされた点が、スプラインによってモデル化された形状の制御点を表す場合、複数位置推定のロバストフュージョンを使用することによって、形状の位置推定の信頼性が増大する。形状空間（shape space）が、学習された部分空間制約によって制約される場合、補正も少なくなる。輪郭が利用可能な場合、トラッキングに使用されるモデルテンプレートは、形状間の距離に基づいてモデル集合からオンラインで選択することができる。
【０１０７】
図６ａおよび６ｂは左心室の心エコー図の画像を示しており、ここでは心内膜壁の運動が本発明にしたがってトラッキングされている。画像６０２，６１０に対するフレーム１は、壁部の初期輪郭を示しており、これはこの画像においてドットで示されている。後続のフレーム６０４，６０６，６０８，６１２，６１４，６１６は、この壁部の運動が時間についてどのようにトラッキングされるかを示している。各ドットに対する測定は本発明にしたがって行われる。
【０１０８】
例
この例において、心エコー図画像において手動でトレースされた心臓輪郭がトレーニングセットとして使用される。先端の４室（開いた輪郭）および胸骨傍短軸（閉じた輪郭）の両側から見た図が、図７ａ〜７ｄおよび図８ａ〜８ｄに示したように検査される。輪郭毎にランドマーク点が付される。図７ｂおよび図８ｂは、一緒に引かれたトレーニング輪郭のセットを示している。各輪郭がベクトルであり、このベクトルがそのコンポーネントとして、順序づけされたランドマーク点の座標を有する箇所に固有分析アプローチを適用する（開いた輪郭に対して３４個のコンポーネントまた閉じた輪郭に対して３６個のコンポーネント）。この場合、列がトレーニングベクトルである行列にＰＣＡを実行する。つぎに固有値を使用して、部分空間におけるモデル共分散行列として対角行列を形成する。測定値共分散行列は、さまざまなシナリオをテストするために調整される。
【０１０９】
図７ａおよび８ａはテスト輪郭を示しており、ここでは実線の曲線７０４，８０２は、グラウンドトルースであり、破線の曲線７０２，８０４は、ノイズを有する測定値である。比較的不確定な点を一層調整することによって、目下のノイズを有する輪郭に最も近い、固有輪郭部分空間における輪郭を探すことが望ましい。図７ｃおよび図８ｃは、図７ａおよび図８ａの輪郭に垂直な投影を行った場合の結果を示している（等方性の共分散は、ランドマーク点のまわりの小さな円７１２，８１２によって示されている）。図７ｄおよび８ｄは結果を示している。この結果は、はるかにグラウンドトルースに近いことがわかる。トレーニングデータがテストデータとはかなり異なるため、これは完全な一致にはならない。このため、テスト輪郭におけるいくつかの小さな形状の変形は、部分空間において実現できないことがある。
【０１１０】
本発明によれば、輪郭は視覚化ツールを使用して診断に役立てることができる。医師が心エコー図から心臓疾患をより簡単に診断できるようにするため、ＬＶ輪郭は、リアルタイムにカラーベクトルで表示される。ベクトルの長さが変位の大きさを表すのに対し、ベクトルの方向はＬＶ壁輪郭における点の運動方向を表す。例えば、オレンジ色は、運動が収縮である場合を示すために使用され、暗い青は、運動が拡張である場合を示すために使用される。図９は表示結果の一例を示している。ＬＶ輪郭における点の運動が追跡されている。この運動は時間および空間領域の両方においてガウシアンによって平滑化される。医師がよりわかりやすくするためにベクトルの長さは３倍される。
【０１１１】
この視覚化法により、医師は心内膜の各セグメントの運動の大きさおよび方向を用意に見ることができる。この視覚化方法と大域的な運動補償とを組み合わせることによって、医師はＬＶのすべてのセグメントにおいて収縮方向および大きさを簡単に見ることができる。虚血性の領域または他の以上な領域を直接かつ人間の目で識別することができる。
【０１１２】
健康診断中に機能の感度を決定する方法に対する実施形態を説明したが、上に教示したことを考え合わせれば、当業者は変更および変形をなし得ることに注意されたい。したがって開示した本発明の特定の実施形態において種々の変更を行うことができ、またこれらの変更が添付の請求項によって定められる本発明の範囲内にあることを理解されたい。ここまで特許法によって要求される詳細で本発明を説明したが、請求され特許法によって保護されるものは添付の請求項に記載されている。
【図面の簡単な説明】
【０１１３】
【図１】本発明にしたがって形状追跡法を実現する例示的なシステムのブロック図である。
【図２】音響的なドロップアウトのエリアと、推定された局所的な壁運動の不確かさとを示す心エコー図画像である。
【図３】制御点の位置の不確かさが極めて非等方であり不均質である心内膜輪郭を示す左心室の心エコー図である。
【図４】相異なるα値で強く適合化したＰＣＡ（ＳＡ−ＰＣＡ）モデルおよびインクリメント式ＰＣＡモデルの例を示す図である。
【図５】心内膜壁が初期化され（ａ）、ＩＰＣＡ（ｂ）およびＳＡ−ＰＣＡ（ｃ）を使用して追跡される左心室の心エコー図の画像である。
【図６】本発明にしたがって心内膜壁の動きを追跡した左心室の心エコー図の画像である。
【図７】本発明にしたがって先端の４室側から見た左心室の心内膜壁の心エコー図画像を表すテスト輪郭である。
【図８】本発明にしたがって胸骨傍短軸側から見た左心室の心内膜壁の心エコー図画像を表すテスト輪郭である。
【図９】本発明による画像視覚化方法を例示する心エコー図画像である。
【図１０】本発明の１実施形態による、マルチモデル制御点ベースの追跡器を示す図である。
【特許請求の範囲】
【請求項１】
運動中のオブジェクトの大域的な形状を追跡する方法において、
当該の大域的な形状に沿って１つ以上の制御点を定めるステップと、
前記オブジェクトの運動に伴って当該の１つ以上の制御点の各々を追跡するステップと、
運動中の制御点の位置の不確かさを表すステップと、
当該の不確かさを利用し、予備の形状モデルを使用して前記大域的な形状を制約するステップとを含むことを特徴とする、
運動中のオブジェクトの大域的な形状を追跡する方法。
【請求項２】
前記の１つ以上の制御点の各々を追跡するステップは、ベイズカーネルマッチングアプローチを使用する、
請求項１に記載の方法。
【請求項３】
前記の１つ以上の制御点の各々を追跡するステップは、オプティカルフローベースのアプローチを使用する、
請求項１に記載の方法。
【請求項４】
マルチテンプレートアダプティブマッチングフレームワークを使用して、前記の１つ以上の制御点の追跡中に累積誤差を修正するステップをさらに含む、
請求項１に記載の方法。
【請求項５】
前記の予備の形状モデルは、部分空間形状制約モデルである、
請求項１に記載の方法。
【請求項６】
ＳＡ−ＰＣＡ（strong adapted Principal Component Analysis）モデルを使用して、前記の形状を制約する、
請求項１に記載の方法。
【請求項７】
前記の不確かさを利用するステップはさらに共分散行列を使用し、非正射影投影または情報フュージョンを利用し、部分空間形状制約モデルを使用することを含む、
請求項１に記載の方法。
【請求項８】
運動中のオブジェクトの大域的な形状を追跡するシステムにおいて、
当該の大域的な形状の初期輪郭に沿って１つ以上の制御点を定める手段と、
前記オブジェクトの運動に伴って当該の１つ以上の制御点の各々を追跡する手段と、
運動中の制御点の位置の不確かさを表す手段と、
当該の不確かさを利用し、予備の形状モデルを使用して前記大域的な形状を制約する手段とを有することを特徴とする、
運動中のオブジェクトの大域的な形状を追跡するシステム。
【請求項９】
前記の追跡手段により、ベイズカーネルマッチングアプローチが使用される、
請求項８に記載のシステム。
【請求項１０】
前記の追跡手段により、オプティカルフローベースのアプローチが使用される、
請求項８に記載のシステム。
【請求項１１】
マルチテンプレートアダプティブマッチングフレームワークを使用して、前記の１つ以上の制御点の追跡中に累積誤差を修正する手段がさらに含まれている、
請求項８に記載のシステム。
【請求項１２】
ＳＡ−ＰＣＡ（strong adapted Principal Component Analysis）モデルが使用されて、前記の輪郭の形状が制約される、
請求項８に記載のシステム。
【請求項１３】
前記の不確かさを利用するステップはさらに共分散行列を使用し、非正射影投影または情報フュージョンを利用し、部分空間形状制約モデルを使用することを含む、
請求項８に記載のシステム。
【請求項１４】
オブジェクトの形状の運動を視覚的に追跡する方法において、
第１のカラーの１つ以上のベクトルを生成して、前記の形状の輪郭に沿った制御点の収縮を表すステップと、
第２のカラーの１つ以上のベクトルを生成して、前記の形状の輪郭に沿った制御点の拡張を表すステップと、
前記の第１および第２のカラーのベクトルを周期的に表示して、前記形状の運動を表示するステップとを有することを特徴とする、
オブジェクトの形状の運動を視覚的に追跡する方法。
【請求項１５】
オブジェクトの形状の運動を視覚的に追跡するシステムにおいて、
第１のカラーの１つ以上のベクトルを生成して、前記の形状の輪郭に沿った制御点の収縮を表す手段と、
第２のカラーの１つ以上のベクトルを生成して、前記の形状の輪郭に沿った制御点の拡張を表す手段と、
前記の第１および第２のカラーのベクトルを周期的に表示して、前記形状の運動を表示する手段とを有することを特徴とする、
オブジェクトの形状の運動を視覚的に追跡するシステム。
【請求項１６】
運動中のオブジェクトの大域的な形状を追跡する方法において、
当該の大域的な形状に沿って１つ以上の制御点を定めるステップと、
前記オブジェクトの運動に伴って当該の１つ以上の制御点の各々を追跡するステップと、
制御点毎に複数の外観モデルを形成するステップと、
各モデルによって形成される動きベクトルを組み合わせて形状を追跡するステップとを有することを特徴とする、
運動中のオブジェクトの大域的な形状を追跡する方法。
【請求項１７】
ロバスト情報フュージョンを使用して、各制御点の位置の推定を計算するステップをさらに含む、
請求項１６に記載の方法。
【請求項１８】
前記の複数の外観モデルを、前記の１つ以上の制御点の位置の不確かさに基づいて形成する、
請求項１６に記載の方法。
【請求項１９】
前記のロバスト情報フュージョンは、可変帯域幅密度ベースフュージョンＶＢＤＦ（Variable-Bandwidth Density-based Fusion）である、
請求項１８に記載の方法。
【請求項２０】
ロバスト情報フュージョンは、
平均ベクトルおよび共分散行列によって表される多数の測定値を入力するステップと、
可変帯域幅密度推定値を計算することによって前記測定値を一緒に融合するステップと、
可変帯域幅密度推定値の最重要ノードを検出するステップとを有することを特徴とする、
請求項１９に記載の方法。
【請求項１】
運動中のオブジェクトの大域的な形状を追跡する方法において、
当該の大域的な形状に沿って１つ以上の制御点を定めるステップと、
前記オブジェクトの運動に伴って当該の１つ以上の制御点の各々を追跡するステップと、
運動中の制御点の位置の不確かさを表すステップと、
当該の不確かさを利用し、予備の形状モデルを使用して前記大域的な形状を制約するステップとを含むことを特徴とする、
運動中のオブジェクトの大域的な形状を追跡する方法。
【請求項２】
前記の１つ以上の制御点の各々を追跡するステップは、ベイズカーネルマッチングアプローチを使用する、
請求項１に記載の方法。
【請求項３】
前記の１つ以上の制御点の各々を追跡するステップは、オプティカルフローベースのアプローチを使用する、
請求項１に記載の方法。
【請求項４】
マルチテンプレートアダプティブマッチングフレームワークを使用して、前記の１つ以上の制御点の追跡中に累積誤差を修正するステップをさらに含む、
請求項１に記載の方法。
【請求項５】
前記の予備の形状モデルは、部分空間形状制約モデルである、
請求項１に記載の方法。
【請求項６】
ＳＡ−ＰＣＡ（strong adapted Principal Component Analysis）モデルを使用して、前記の形状を制約する、
請求項１に記載の方法。
【請求項７】
前記の不確かさを利用するステップはさらに共分散行列を使用し、非正射影投影または情報フュージョンを利用し、部分空間形状制約モデルを使用することを含む、
請求項１に記載の方法。
【請求項８】
運動中のオブジェクトの大域的な形状を追跡するシステムにおいて、
当該の大域的な形状の初期輪郭に沿って１つ以上の制御点を定める手段と、
前記オブジェクトの運動に伴って当該の１つ以上の制御点の各々を追跡する手段と、
運動中の制御点の位置の不確かさを表す手段と、
当該の不確かさを利用し、予備の形状モデルを使用して前記大域的な形状を制約する手段とを有することを特徴とする、
運動中のオブジェクトの大域的な形状を追跡するシステム。
【請求項９】
前記の追跡手段により、ベイズカーネルマッチングアプローチが使用される、
請求項８に記載のシステム。
【請求項１０】
前記の追跡手段により、オプティカルフローベースのアプローチが使用される、
請求項８に記載のシステム。
【請求項１１】
マルチテンプレートアダプティブマッチングフレームワークを使用して、前記の１つ以上の制御点の追跡中に累積誤差を修正する手段がさらに含まれている、
請求項８に記載のシステム。
【請求項１２】
ＳＡ−ＰＣＡ（strong adapted Principal Component Analysis）モデルが使用されて、前記の輪郭の形状が制約される、
請求項８に記載のシステム。
【請求項１３】
前記の不確かさを利用するステップはさらに共分散行列を使用し、非正射影投影または情報フュージョンを利用し、部分空間形状制約モデルを使用することを含む、
請求項８に記載のシステム。
【請求項１４】
オブジェクトの形状の運動を視覚的に追跡する方法において、
第１のカラーの１つ以上のベクトルを生成して、前記の形状の輪郭に沿った制御点の収縮を表すステップと、
第２のカラーの１つ以上のベクトルを生成して、前記の形状の輪郭に沿った制御点の拡張を表すステップと、
前記の第１および第２のカラーのベクトルを周期的に表示して、前記形状の運動を表示するステップとを有することを特徴とする、
オブジェクトの形状の運動を視覚的に追跡する方法。
【請求項１５】
オブジェクトの形状の運動を視覚的に追跡するシステムにおいて、
第１のカラーの１つ以上のベクトルを生成して、前記の形状の輪郭に沿った制御点の収縮を表す手段と、
第２のカラーの１つ以上のベクトルを生成して、前記の形状の輪郭に沿った制御点の拡張を表す手段と、
前記の第１および第２のカラーのベクトルを周期的に表示して、前記形状の運動を表示する手段とを有することを特徴とする、
オブジェクトの形状の運動を視覚的に追跡するシステム。
【請求項１６】
運動中のオブジェクトの大域的な形状を追跡する方法において、
当該の大域的な形状に沿って１つ以上の制御点を定めるステップと、
前記オブジェクトの運動に伴って当該の１つ以上の制御点の各々を追跡するステップと、
制御点毎に複数の外観モデルを形成するステップと、
各モデルによって形成される動きベクトルを組み合わせて形状を追跡するステップとを有することを特徴とする、
運動中のオブジェクトの大域的な形状を追跡する方法。
【請求項１７】
ロバスト情報フュージョンを使用して、各制御点の位置の推定を計算するステップをさらに含む、
請求項１６に記載の方法。
【請求項１８】
前記の複数の外観モデルを、前記の１つ以上の制御点の位置の不確かさに基づいて形成する、
請求項１６に記載の方法。
【請求項１９】
前記のロバスト情報フュージョンは、可変帯域幅密度ベースフュージョンＶＢＤＦ（Variable-Bandwidth Density-based Fusion）である、
請求項１８に記載の方法。
【請求項２０】
ロバスト情報フュージョンは、
平均ベクトルおよび共分散行列によって表される多数の測定値を入力するステップと、
可変帯域幅密度推定値を計算することによって前記測定値を一緒に融合するステップと、
可変帯域幅密度推定値の最重要ノードを検出するステップとを有することを特徴とする、
請求項１９に記載の方法。
【図１】
【図２】
【図３】
【図４】
【図１０】
【図２】
【図３】
【図４】
【図１０】
【公表番号】特表２００６−５２４５３４（Ｐ２００６−５２４５３４Ａ）
【公表日】平成１８年１１月２日（２００６．１１．２）
【国際特許分類】
【出願番号】特願２００６−５０６９５９（Ｐ２００６−５０６９５９）
【出願日】平成１６年３月８日（２００４．３．８）
【国際出願番号】ＰＣＴ／ＵＳ２００４／００７０６８
【国際公開番号】ＷＯ２００４／０８１８７５
【国際公開日】平成１６年９月２３日（２００４．９．２３）
【出願人】（５９３０７８００６）シーメンス コーポレイト リサーチ インコーポレイテツド (47)
【氏名又は名称原語表記】Ｓｉｅｍｅｎｓ Ｃｏｒｐｏｒａｔｅ Ｒｅｓｅａｒｃｈ，Ｉｎｃ．
【住所又は居所原語表記】７５５ Ｃｏｌｌｅｇｅ Ｒｏａｄ Ｅａｓｔ，Ｐｒｉｎｃｅｔｏｎ， ＮＪ ０８５４０，Ｕｎｉｔｅｄ Ｓｔａｔｅｓ ｏｆ Ａｍｅｒｉｃａ
【出願人】（５９３０６３１０５）シーメンス メディカル ソリューションズ ユーエスエー インコーポレイテッド (156)
【氏名又は名称原語表記】Ｓｉｅｍｅｎｓ Ｍｅｄｉｃａｌ Ｓｏｌｕｔｉｏｎｓ ＵＳＡ，Ｉｎｃ．
【住所又は居所原語表記】５１ Ｖａｌｌｅｙ Ｓｔｒｅａｍ Ｐａｒｋｗａｙ，Ｍａｌｖｅｒｎ，ＰＡ １９３５５−１４０６，Ｕ．Ｓ．Ａ．
【Ｆターム（参考）】
【公表日】平成１８年１１月２日（２００６．１１．２）
【国際特許分類】
【出願日】平成１６年３月８日（２００４．３．８）
【国際出願番号】ＰＣＴ／ＵＳ２００４／００７０６８
【国際公開番号】ＷＯ２００４／０８１８７５
【国際公開日】平成１６年９月２３日（２００４．９．２３）
【出願人】（５９３０７８００６）シーメンス コーポレイト リサーチ インコーポレイテツド (47)
【氏名又は名称原語表記】Ｓｉｅｍｅｎｓ Ｃｏｒｐｏｒａｔｅ Ｒｅｓｅａｒｃｈ，Ｉｎｃ．
【住所又は居所原語表記】７５５ Ｃｏｌｌｅｇｅ Ｒｏａｄ Ｅａｓｔ，Ｐｒｉｎｃｅｔｏｎ， ＮＪ ０８５４０，Ｕｎｉｔｅｄ Ｓｔａｔｅｓ ｏｆ Ａｍｅｒｉｃａ
【出願人】（５９３０６３１０５）シーメンス メディカル ソリューションズ ユーエスエー インコーポレイテッド (156)
【氏名又は名称原語表記】Ｓｉｅｍｅｎｓ Ｍｅｄｉｃａｌ Ｓｏｌｕｔｉｏｎｓ ＵＳＡ，Ｉｎｃ．
【住所又は居所原語表記】５１ Ｖａｌｌｅｙ Ｓｔｒｅａｍ Ｐａｒｋｗａｙ，Ｍａｌｖｅｒｎ，ＰＡ １９３５５−１４０６，Ｕ．Ｓ．Ａ．
【Ｆターム（参考）】
[ Back to top ]