(57)【要約】  少なくとも第一と第二デジタル画像又は同じ対象の一連の断面画像のような生物医学画像を登録する方法において、第一画像又は画像セットの中で、一定数のランドマーク、いわゆる特徴が、一定数のピクセル又はボクセルを選択することで個別化される。特徴として選択された各ピクセル又はボクセルの位置は、オプティカルフローベクトルを決定することによって第一から第二画像又は画像セット内で追跡される。第一と第二画像又は画像セットは、逆オプティカルフローベクトルを第二画像又は画像セットのピクセル又はボクセルに適用することで登録される。本発明は、第一画像又は第一の一連の断面画像の各ピクセル又はボクセルの周囲近傍のピクセル又はボクセルを定義することに存する自動追跡可能なランドマーク選択ステップを提供する。各目標ピクセル又はボクセルは、この目標ピクセル又はボクセルと、ウインドウの各ピクセル又はボクセル又はその一部の外観を記述するパラメータの関数として、そして数字で表したマトリックス又は前記ウインドウのピクセル又はボクセルを表現する前記数字で表したマトリックスの変形の一つ以上の特徴的なパラメータの関数として計算される一つ以上の特徴的なパラメータを決定する。有効に追跡可能なランドマークと同時に起こるピクセル又はボクセルが、前記目標ピクセル又はボクセルの前記特徴的なパラメータの関数として決定される。
2009512527.tif 000002
【特許請求の範囲】
【請求項１】
ａ）ピクセル又はボクセルからなる２次元又は３次元アレイで構成された、少なくとも第一と第二デジタル画像又は同じ対象の一連の断面画像を提供するステップと、
ｂ）第一画像又は画像セットの中に、ランドマークとして設定された一定数のピクセル又はボクセルを選択し、追跡すべき前記特徴のリストを生成することで、一定数のいわゆる特徴であるランドマークを定義するステップと、
ｃ）第一画像から第二画像又は、第一から第二の画像位置間のオプティカルフローベクトルによって決定された後の時点で得られた画像又は画像セット又は、特徴として選択された各ピクセル又はボクセルの後の時点で得られた画像又は画像セットから、特徴として選択された各ピクセル又はボクセルの位置を追跡するステップと、
ｄ）第一画像と第二画像或いは、逆オプティカルフローベクトルを前記第二画像又は前記画像セットのピクセル又はボクセルの位置に適用することによって後の時点で得られた画像又は画像セット、を登録するステップと、
からなる、対象の移動によって生じる画像アーチファクトを低減する、生物医学画像の登録方法であって、
Ｂ１）第一画像又は第一の一連の断面画像の各ピクセル又はボクセル周りに隣接した、限られた数で構成されるピクセル又はボクセルを定義するステップと、
Ｂ２）各目標ピクセル又はボクセルに対して、外観を記載する数値パラメータ、すなわち前記ピクセル又はボクセル窓の前記目標ピクセル又はボクセルとピクセル又はボクセルの部分とからなるいわゆる数値外観パラメータの関数として計算され、或いは、前記窓のピクセル又はボクセルを表現する数値パラメータのマトリックス又は数値パラメータの前記マトリックスの変形で構成される一つ以上の特徴パラメータの関数として計算される、一つ以上の特徴のあるパラメータを定義するステップと、
Ｂ３）前記目標ピクセル又はボクセルの特徴的なパラメータの関数としての有効に追跡可能なランドマークに存するピクセル又はボクセルを定義するステップと、
からなる自動追跡可能なランドマーク選択ステップが実行されることを特徴とする登録方法。
【請求項２】
もし前記目標ピクセル又はボクセルが追跡すべき有効なランドマークであった場合に、定義するための関数が、前記目標ピクセル又はボクセルの特徴パラメータの各一つ又は組み合わせの値と閾値とを比較することであることを特徴とする請求項１記載の方法。
【請求項３】
前記各特徴パラメータに対する閾値を求め、前記各特徴パラメータを対応する閾値と比較して、単一の比較結果の関数としての全体的評価パラメータから、有効に追跡可能なランドマークと有効でない追跡可能なランドマークの質を求めることを特徴とする請求項２に記載の方法。
【請求項４】
前記各特徴パラメータに対して閾値の関数としての全体的閾値が求められ、前記全体的評価パラメータと前記全体的閾値とを比較することによって、有効に追跡可能なランドマークと有効でない追跡可能なランドマークの質を求めることを特徴とする請求項３に記載の方法。
【請求項５】
有効に追跡可能なランドマークか否かに関して比較した各結果を数値で表現するために可変有効性が定義され、
部分又は各特性と対応する閾値との比較の前記可変有効性の値の関数として、全体的評価パラメータが求められ、そして
全体的な評価パラメータの関数として、有効に追跡すべきランドマークと有効でない追跡すべきランドマークの質が求められる
ことを特徴とする請求項３又は４に記載の方法。
【請求項６】
包括的評価パラメータから有効に追跡可能なランドマークと有効でない追跡可能なランドマークを求めるための関数が、前記包括的評価パラメータと包括的閾値との比較であることを特徴とする請求項５に記載の方法。
【請求項７】
特徴パラメータセットが一定数のサブセットに細分化されて、各サブセットに対し二次特徴パラメータがそのサブセットの特徴パラメータの関数として求められ、そして、
各二次特徴パラメータに対し閾値が定義されて、各二次特徴パラメータに対し有効に追跡可能なランドマークと有効でない追跡可能なランドマークの質が、前記二次特徴パラメータと対応する閾値との関数によって求められる
ことを特徴とする請求項１ないし６のいずれか１項に記載の方法。
【請求項８】
各二次特徴パラメータと対応する閾値が、特徴パラメータの対応するサブセットの閾値の関数として求められることを特徴とする請求項７に記載の方法。
【請求項９】
前記サブセットの特徴パラメータが、少なくとも共通の特徴を有し、
特に前記目標ピクセル又はボクセル又は前記ピクセル又はボクセル窓のピクセル又はボクセルのそれぞれ又は一部の外観を記載する数値的外観パラメータからなり、そして
前記窓のピクセル又はボクセルを表現する数値行列の一つ以上の特徴パラメータの同じ関数によって、又は前記数値行列の同じ変換によって求められる
ことを特徴とする請求項７又は８に記載の方法。
【請求項１０】
前記包括的評価パラメータが二次特徴パラメータの関数として定義され、そして、
有効に追跡可能なランドマークと有効でない追跡可能なランドマークの質が、前記包括的評価パラメータの関数と、前記二次特徴パラメータに対する閾値の関数として存する包括的閾値と、によって求められることを特徴とする請求項７ないし１０のいずれか１項に記載の方法。
【請求項１１】
有効に追跡可能なランドマークと有効でない追跡可能なランドマークの質を求めるために、各特徴パラメータ及び／又は各二次特徴パラメータ及び／又は各包括的閾値が加重され、そして、その加重が特徴パラメータ又は二次特徴パラメータの関連性に従って求められることを特徴とする請求項２ないし１０のいずれか１項に記載の方法。
【請求項１２】
前記包括的評価パラメータ及び／又は前記二次特徴パラメータ及び／又は包括的閾値を求めるための関数が、一つ以上の特徴パラメータの係数の線形結合か、二乗値の合計か、前記特徴パラメータ又は前記第二パラメータ又は前記閾値が処理される非線形関数に存することを特徴とする請求項１ないし１１のいずれか１項に記載の方法。
【請求項１３】
前記包括的評価パラメータ及び／又は前記二次特徴パラメータ及び／又は前記包括的閾値を求めるための関数が、各特徴パラメータ又は第二パラメータ又は閾値を加重するためにも用いられることを特徴とする請求項１２に記載の方法。
【請求項１４】
前記目標ピクセル又はボクセルの有効に追跡可能なランドマークと有効でない追跡可能なランドマークの質が、分類アルゴリズムを用いて画像の各ピクセル又はボクセルの特徴パラメータを処理することによって求められることを特徴とする請求項１に記載の方法。
【請求項１５】
前記各目標ピクセルのIDが、ベクトル処理されるためにコード化され、前記ベクトルの成分が、前記ピクセル又はボクセルの一つ以上の特徴パラメータに存することを特徴とする請求項１４に記載の方法。
【請求項１６】
特徴パラメータセットが一定数のサブセットに細分化され、
前記特徴パラメータの代わりに、二次特徴パラメータが前記分類アルゴリズムによって処理され、そして、前記二次特徴パラメータが、特徴パラメータのサブセットの関数として求められることを特徴とする請求項１４又は１５に記載の方法。
【請求項１７】
各ピクセル又はボクセルに対して前記分類アルゴリズムによって処理された入力ベクトル成分が、各特徴パラメータと対応する閾値との比較結果、又は各二次特徴パラメータと対応する閾値との比較結果を数値的に表現する可変有効性に存することを特徴とする請求項１４又は１５に記載の方法。
【請求項１８】
前記分類アルゴリズムが、クラスタ化アルゴリズム又は予測アルゴリズムであることを特徴とする請求項１４ないし１７のいずれか１項に記載の方法。
【請求項１９】
前記分類アルゴリズムが、人工の神経ネットワークであることを特徴とする請求項１９記載の方法。
【請求項２０】
一定数の有効なランドマークが、有効に追跡可能なランドマークとして識別された既知の場合の一定数の画像を提供するステップと、
上記ステップＢ１）とＢ２）に存する自動追跡可能ランドマーク選択ステップを適用することによって、一定数の画像の中で有効に追跡可能であると識別された前記ランドマークに対応して前記ピクセル又はボクセルに対する特徴パラメータセットを求めるステップと、
有効に追跡可能であると識別された各ランドマークに一意的に関連付けられ、成分として前記有効に追跡可能なランドマークと一致する前記ピクセル又はボクセルの特徴パラメータを備えたベクトルを生成するステップと、
変数の所定の数値によって有効に追跡可能な質を表現して、前記数値を、有効に追跡可能なランドマークと一致する各ピクセル又はボクセルをコード化するベクトルに関連付けるステップと、
有効に追跡可能なランドマークと一致する各ピクセル又はボクセルをコード化する各ベクトルが、分類アルゴリズムのためのトレーニングデータベースの登録を形成するステップと、
前記データベースによって分類アルゴリズムをトレーニングするステップと、
鍛えられた分類アルゴリズムの入力に前記目標ピクセル又はボクセルの前記特徴パラメータを備えたベクトルを提供することによって、目標ピクセル又はボクセルの有効に追跡可能なランドマークの質を求めるステップと、
からなることを特徴とする請求項１４ないし１９のいずれか１項に記載の方法。
【請求項２１】
前記トレーニングデータベースが、有効に追跡可能ではないと知られたランドマークと一致する画像のピクセル又はボクセルのベクトル形式に組織化され、
有効に追跡可能でないランドマークの質を、追跡可能なランドマークの質を表現する前記変数の値とは異なる変数の所定の数値によって表現する特徴パラメータが備えられ、
一方で有効に追跡可能でないランドマークの質を表現する前記数値が、既知の有効に追跡可能でないランドマークと一致した各ピクセル又はボクセルをコード化するベクトルに加えられたことを特徴とする請求項２０に記載の方法。
【請求項２２】
画像のピクセル又はボクセルが、いわゆるLucas & Kanade自動特徴追跡法又はアルゴリズムによって、事前に又は平行して処理されることを特徴とする請求項１ないし２２のいずれか１項に記載の方法。
【請求項２３】
前記Lucas & Kanadeアルゴリズムが、
ａ）ボクセルの３次元アレイによって形成された第一３次元画像ボリュームＩを提供するステップと、
ｂ）目標ボクセルを取り囲み、例えば前記目標ボクセルを中心とした3×３×３アレイのような調整可能なサイズを有する隣接ボクセルを各目標ボクセルに対して定義するステップと、
ｃ）デカルト座標系を定義し、前記デカルト座標系の各軸に対し、隣接ボクセルに対して第一画像ボリュームの各ボクセルの強度勾配によって形成されるいわゆる勾配ボリュームを計算するステップと、
ここで傾斜ボリュームは次式で定義され、
ア）ｘ方向の勾配ボリューム：
【数１】
2009512527.tif 000024
イ）ｙ方向の勾配ボリューム：
【数２】
2009512527.tif 000025
ウ）ｚ方向の勾配ボリューム：
【数３】
2009512527.tif 000026
ｄ）各目標ボクセルに対し、次式で定義されるいわゆる勾配行列を計算するステップと、
【数４】
2009512527.tif 000027
ｅ）前記３次元画像ボリュームの各目標ボクセルの各勾配行列に対して、以下の等式を適用することで最小固有値？ｍを計算するステップと、
【数５】
2009512527.tif 000028
ここで
【数６】
2009512527.tif 000029
に対して、前記勾配行列の最小固有値を計算する
【数７】
2009512527.tif 000030
ｆ）前記最小固有値？ｍの閾値を定義するステップと、
ｇ）以下の基準を満足する、対応する勾配行列が最小固有値？ｍを有する画像Ｉの追跡すべき特徴としての各ボクセルを選択するステップと、
ｉ）？ｍが前記閾値より大きく、
ii）？ｍが固有値の全ての最小値？ｍの最大値の割合より小さくなく、
iii）もし特徴として選択された別のボクセルが、選択された特徴として維持されている目標ボクセル回りの定義された隣接ボクセルに存在する場合、勾配行列がより大きな固有値の最小値？ｍを有するボクセルの一つのみが特徴として選択され、他の一つは追跡可能な特徴のリストから廃棄され、
iv）特徴として選択されたボクセルの回りの３次元隣接ボクセルの平均信号強度が、調整
可能な平均信号強度閾値より大きい、
以上の各ステップを備えることを特徴とする請求項２２に記載の方法。
【請求項２４】
特徴として選択されたボクセル回りの隣接ボクセルの平均信号強度に対する閾値が約１０に設定されることを特徴とする請求項２３に記載の方法。
【請求項２５】
前記特徴パラメータが、選択窓のピクセル又はボクセルの強度値に存し、また目標ピクセル又はボクセルの強度値に存することを特徴とする請求項１ないし２４のいずれか１項に記載の方法。
【請求項２６】
前記特徴パラメータが、選択された窓のピクセル又はボクセルの画像データを備えた数値行列の特異値に代替的に又は付加的に存することを特徴とする請求項１ないし２５のいずれか１項に記載の方法。
【請求項２７】
前記特徴パラメータが、選択された窓のピクセル又はボクセルを表現する数値行列の勾配行列の固有値に代替的に又は付加的に存することを特徴とする請求項１ないし２６のいずれか１項に記載の方法。
【請求項２８】
前記特徴パラメータが、選択された窓のピクセル又はボクセルを表現する数値行列のヘッシアン行列の固有値に代替的に又は付加的に存することを特徴とする請求項１ないし２７のいずれか１項に記載の方法。
【請求項２９】
前記特徴パラメータが、一つ以上のウエーブレットに基づく関数を用いて前記数値行列を処理することで得られた前記窓のピクセル又はボクセルを表現する数値行列のウエーブレット変換係数の一つ以上又はその組み合わせに代替的に又は付加的に存することを特徴とする請求項１ないし２８のいずれか１項に記載の方法。
【請求項３０】
前記特徴パラメータが、いわゆる行列の同時生起変換に代替的に又は付加的に存することを特徴とする請求項１ないし２９のいずれか１項に記載の方法。
【請求項３１】
前記特徴パラメータが、前記窓のピクセル又はボクセルを表現する数値行列の自己相関変換の一つ以上の係数に代替的に又は付加的に存することを特徴とする請求項１ないし３０のいずれか１項に記載の方法。
【請求項３２】
前記特徴パラメータが、
前記窓のピクセル又はボクセルを表現する数値行列の固有値又は特異値及び／又は、
前記窓のピクセル又はボクセルを表現する数値行列の勾配行列又はヘッシアン行列の固有値及び／又は、
ウェーブレット変換係数の一つ以上の係数及び／又は、
自己相関変換の一つ以上の係数及び／又は、
前記窓のピクセル又はボクセルを表現する数値行列の同時生起行列の特異値
の組み合わせに存することを特徴とする請求項１ないし３１のいずれか１項に記載の方法。
【請求項３３】
同一対象の少なくとも二つ以上の画像が得られるか、又は異なる時間に得られた少なくとも二つ以上の画像が提供され、
第一画像ボリュームＩにおいて有効に追跡可能なランドマークと一致して区別された各ボクセルに対して、前記特徴が、第二の遅い時間に得られた同一対象の追加の画像Ｊのボクセルアレイにおけるその位置に対して追跡され、
前記追跡が、いわゆる「Lucas＆Kanade特徴追跡」のピラミッド形の実施によって行われることを特徴とする請求項１ないし３２のいずれか１項に記載の方法。
【請求項３４】
前記いわゆる「Lucas＆Kanade特徴追跡」のピラミッド形の実施が、
画像Ｉにおける特徴として選択されたボクセルの座標をそれぞれ有する、画像Iのボクセルの３次元アレイに対応する画像ボリュームIにおけるポイントとしてのポイントｕと、画像Ｊのボクセルの３次元アレイに対応する画像ボリュームＪにおける対応するポイントｖとを定義するステップと、
元のボリュームを表現する下部構造（Ｌ＝０）と共に、レベルＬ  ０，．．．，ｍ  におけるボリュームを表現するＩLとＪLとを用いて画像ボリュームＩとＪの二つの理想的なピラミッドを構築するステップと、
最上階Ｌmにおいて初期化されるいわゆる包括的ピラミッド形の運動ベクトルｇを定義するステップと、
【数８】
2009512527.tif 000031
ｉ）uLとしてポイントｕの位置を定義し、以下の式で前記位置を算出し、ここでＬは
ピラミッドの実際の階又はレベルの値を有し、
【数９】
2009512527.tif 000032
ここでｐは実際のボリューム位置であり、
ii）ボリューム勾配に対する以下の式に従って、デカルト座標ｘ、ｙ、ｚの各方向にお
けるボリューム勾配を計算し、
【数１０】
2009512527.tif 000033
iii）以下の等式に従って、勾配行列を計算するために勾配ボリュームを使用し、
【数１１】
2009512527.tif 000034
ここで
【数１２】
2009512527.tif 000035
ここで数値？は、エリアサイズ又は追跡される特徴を表現する目標ボクセルに影響を与える隣接ボクセルを定義し、
iv）各レベルLに対し、次式で定義される反復ベクトルを初期化し、
【数１３】
2009512527.tif 000036
ｖ）ｋ＝１から最大カウントＫまで、又は最小置換
【数１４】
2009512527.tif 000037
まで、以下のボリューム差を計算し、
【数１５】
2009512527.tif 000038
vi）以下の等式に従って、不整合なベクトルを計算し、
【数１６】
2009512527.tif 000039
vii）以下の式に従って、オプティカルフローベクトルを求め、
【数１７】
2009512527.tif 000040
ここで、Ｇ-1はステップiii）で求めた行列Ｇの逆勾配行列であり、
viii）上で定義されたオプティカルフローベクトルを用いて以下のように反復運動ベクト
ルを計算し、
【数１８】
2009512527.tif 000041
そしてこれをレベルＬに対する最後のオプチカルフローと等しくし、
【数１９】
2009512527.tif 000042
ix）最後のレベルＬ＝０が次式で定義される最終オプチカルフローベクトルに達するま
でステップｉ）からviii）を各レベルに対し繰り返し、
【数２０】
2009512527.tif 000043
ｘ）以下の式でボリュームにおける追跡すべき特徴を表現するポイントｕに対応するボ
リュームＪにおけるポイントｖの座標を求め、
ｖ＝ｕ＋ｄ
xi）画像Ｉにおいて選択された特徴に対応する各ボクセルに対して上記方法を繰り返
す、
上記各ステップを備えることを特徴とする請求項３３記載の方法。
【請求項３５】
画像Ｉにおける選択された特徴に対応するボクセルのポイントＵに対応する特徴を表現するボクセルに対応すると認定された画像Ｊにおける各ポイントｖに反転オプティカルフローベクトルを適用するか、又は
画像Ｉにおいて選択された特徴であると認定されたボクセルにそれぞれ対応するポイントｕにオプティカルフローベクトルを適用することによって、
二つの画像ＩとＪを登録することを特徴とする請求項３４に記載の方法。
【請求項３６】
ａ）ＭＲＩ又は超音波で得られた同一対象の少なくとも第一と第二ディジタル画像又はディジタル処理された画像又は断面画像セットを提供するステップと、
ここで前記画像はピクセル又はボクセルの二次元又は三次元アレイによって形成され、
ここで同一組織又は組織領域又は同一の身体構造上の領域の操作が、前記同一の組織又は組織領域又は前記身体構造上の領域に造影剤が存在する状態で、二度目又はそれ以降の時間で実施され、
ｂ）ランドマークとして設定されたピクセル又はボクセルの一定数を選択し、追跡すべき前記特徴のリストを生成することによって、第一画像又は画像セット内に一定数のランドマーク、いわゆる特徴を定義するステップと、
ｃ）特徴として選択された各ピクセル又はボクセルに対する第一から第二画像又は画像のセットからオプティカルフローベクトルを求めることによって、第一から第二画像又は画像のセットから特徴として選択された各ピクセル又はボクセルの位置を追跡するステップと、
ｄ）第二画像又は画像セットのピクセル又はボクセルに反転オプティカルフローを適用することによって、第一と第二画像又は画像セットを登録するステップと、
からなる造影剤で強調した画像診断法、特に造影剤で強調したＭＲＩ又は超音波画像診断法との組合わせで提供されることを特徴とする請求項３３ないし３５のいずれか１項に記載の方法。
【請求項３７】
前記特徴選択ステップの後で、前記特徴追跡ステップを実行する前に、
Ｂ１）第一画像又は第一の断面画像セットの各ピクセル又はボクセルの回りに隣接するピクセル又はボクセルを定義するステップと、
ここで前記隣接するピクセル又はボクセルは、制限された数のピクセル又はボクセルを備え、単一画像の場合に二次元隣接部又はボクセルが選択され、断面画像セットの場合に三次元隣接部又はボクセルが選択され、
Ｂ２）各ピクセル又はボクセルに対し、前記隣接部又はボクセルの全ピクセル又はボクセルの平均信号強度を求めるステップと、
Ｂ３）平均信号強度閾値を定義するステップと、
Ｂ４）ステップＢ２で各隣接部又はボクセルに対して求められた平均信号強度値と、前記平均信号強度値とを、前記平均信号強度閾値に対して比較するステップと、
Ｂ５）前記隣接部又はボクセルの前記平均信号強度値が前記ステップＢ４での閾値より大きい場合、前記ピクセル又はボクセルを追跡すべき特徴として定義することにより、追跡すべき特徴のリストに加えるステップと、
からなる追加の自動特徴選択ステップを実行することを特徴とする請求項３３ないし３６のいずれか１項に記載の方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、
ａ）ピクセル又はボクセルからなる２次元又は３次元アレイで構成された、少なくとも第一と第二デジタル画像又は同じ対象の一連の断面画像を提供するステップと、
ｂ）第一画像又は画像セットの中に、ランドマークとして設定された一定数のピクセル又はボクセルを選択し、追跡すべき前記特徴のリストを生成することで、一定数のいわゆる特徴であるランドマークを定義するステップと、
ｃ）第一画像から第二画像又は、第一から第二の画像位置間のオプティカルフローベクトルによって決定された後の時点で得られた画像又は画像セット又は、特徴として選択された各ピクセル又はボクセルの後の時点で得られた画像又は画像セットから、特徴として選択された各ピクセル又はボクセルの位置を追跡するステップと、
ｄ）第一と第二画像又は、逆オプティカルフローベクトルを第二画像又は画像セットのピクセル又はボクセルの位置に適用することによって後の時点で得られた画像又は画像セットを登録するステップと、
からなる対象の移動によって生じる画像アーチファクトを低減する生物医学画像の登録方法に関する。
【背景技術】
【０００２】
この種類の時間系列の一部としての画像登録方法は、例えば米国特許６、５５３、１５２によって知られている。
【０００３】
米国特許６、５５３、１５２は、視覚的ランドマーク識別と、人間オペレータによる二つの画像の対応する画素のマーキングとによって、二つの画像に対する画像登録方法を提供している。
【０００４】
この方法は、ランドマークの選択を完全にオペレーターの技能にゆだねる。そして、顕著な又は一意の認識可能な構造が画像内に存在していなければ、実際に実行することは非常に困難である。例えば磁気共鳴映像（MRI）や超音波やレントゲン写真画像といった診断画像を考えると、画像の構造的特徴によってランドマークの視覚的識別をすることは、極めて困難であり、大きなエラーが生ずる可能性がある。
【０００５】
種々の異なるアルゴリズムが輪郭を識別するために試みられており、例えばFischer,   B.;Modersitzki,   J.   Curvature   Based   Registration   with  Applications  to  MRMammography;  Springer-Verlag  Berlin  Heidelberg:  ICCS 2002,  LNCS 2331,  2002;  pp 202- 206に開示されている。また、２回以上測定されたか又は異なるモダリティで調査された組織／領域の放射線画像内のランドマークの識別に関しては、例えばShi  J,  Tomasi  C.  Good  Features  to  Track.  1994  IEEE  Conference  on  Computer Vision and Pattern Recognition  (CVPR'94)  1994;  593- 600や、Tommasini  T,  Fusiello A,  Trucco E,  Roberto V. Making Good Features Track Better Proc. IEEE Int. Conf. on Computer Vision and Pattern Recognition  (CVPR  '98) 1998;  178-  183や、Ruiter,  N.  V.;  Muller,  T.  0.; Stotzka,  R.;  Kaiser,  W.  A.  Elastic  registration of x ray mammograms and three-dimensional magnetic resonance imaging data. J Digit Imaging 2001, 14,  52-55  が開示されている。
【０００６】
異なる時間に撮影された同一対象の二枚の平面画像を登録して、時間に対する対称のう動きを計算することは、二次元登録と呼ばれる。同様のことを画像の三次元セット、例えばMRI又はコンピュータ断層撮影（CT）スキャナの断面画像セットに対して行うアルゴリズムは、三次元と称される。ランドマークや対象の移動は、これら画像セットのあらゆる方向に生じる場合がある。
【０００７】
比較される画像の異なるモダリティの数に依存して、アルゴリズムは単一様式の、そして複数様式の登録アルゴリズムに分割される。
【０００８】
二枚の連続したMRI画像を比較することは、単一様式の登録である。
【０００９】
通常、単一様式登録アルゴリズムは、両画像内の明確なランドマークや複雑な形状を識別しようとし、ほとんどの場合、二次元上で動作する。二番目のステップにおいて、両画像の対応するランドマークの位置が比較され、それらの運動ベクトルが計算される。この情報は、動きアーチファクトを除去するために、二番目の画像の全ての画素を新しい位置にシフトバックするために用いられる。
【００１０】
「剛体レジストレーション」は、単一ユニットとしての二次元又は三次元画像／ボリュームを一定の方向に移動させる。「非剛体レジストレーション」は、一定のエリア／ボリュームの単一又は複数の画素を異なる方向へ移動させる。
【００１１】
例えば非常に弾力性のある乳房組織に関し、動きアーチファクトを除去するために非剛体レジストレーションが必要である。
【００１２】
骨のような一定の構造をとらない特別な柔らかさゆえに、乳房に直接隣接した部分は、異なる方向に移動可能である。従って、どの登録アルゴリズムも、他の既存のアプローチも、乳房組織全体に広がるランドマークの識別が必要である。
【００１３】
アルゴリズムの一つのグループは、限定された数の発見すべきランドマークを定義する。これらは、その数価によって分類される。乳房のMRIイメージに関し、MRIにおいて脂肪と空気の間に高いコントラストがあるため、これらのアプローチはしばしば乳房の外側境界や乳房構造において多数のランドマークを見つける（Huwer, S.; Rahmel, J.; Wangenheim, A. v. Data-Driven Registration for Lacal Deformations. Pattern Recognition Letters 1996, 17, 951-957）。
【００１４】
圧縮パネルを用いるか、少なくとも両乳房を腹臥位に横たえて固定することで、通常圧縮装置に直接接する皮膚のすべての部分は良く固定される（図１）。対照的に、心臓の鼓動や呼吸によって、両乳房の内部領域は僅かに動き続ける。外側の皮膚と空気の境界と比較して、脂肪と柔組織間のコントラストは非常に低いので、乳房の中心にある有効なランドマークの数は、これらのアルゴリズムにおいて低いままである。
【００１５】
ランドマーク選択における一般的な問題は、選択したランドマークが非常に少ないことが、不十分な又は不正確な登録につながるという事実に存する。しかし、人間の観察者もプログラムも、例えば乳房中央の異なる組織タイプのより低いコントラストによって、より再割付けが不確実であるランドマークを含めざるを得ないので、追加のランドマークを選択することは、必ずしも正確な登録を保証しない。ランドマークの数を増やすことは、常に登録の計算の複雑さを著しく増加させる。
【００１６】
いわゆる非剛体レジストレーションに特に適した特徴選択と特徴追跡アルゴリズムは、いわゆるLucas & Kanadeアルゴリズムと、Lucas BD,  Kanade T.  An Iterative Stereo Vision, IJCAI81; 674-679 1981 及び付録１として下に同封したBouguet J-Y. Pyramidal Implementation of the Lucas Kanade Feature Tracker, TechReport as Part of Open Computer Vision Library   Documentation, Intel Corporationに詳細が開示された、特にピラミッド形の実施である。Lucas & Kanade技術とそのピラミッド形の特徴追跡実施は、画像内の信頼できるランドマークの自動認識や、異なる時間に撮影された二枚の画像間の追跡に、特に適している。選択された特徴の置換又は移動は、シフトとして又はいわゆるオプティカルフローベクトルとして求められる。この技術は、例えばロボット工学におけるコンピュータ・ビジョンに広く適用され、デジタル画像処理といくつかの基礎コースで教えられるコンピュータビジョンにおける当業者の一般的な知識の一部である。
【発明の開示】
【発明が解決しようとする課題】
【００１７】
しかしLucas & Kanade技術を実用的に特に乳房のデジタルの又はデジタル処理された診断画像に適用することを考慮すると、前記技術は興味の範囲外も含めた非常に多くの特徴を識別し、例えば乳房まわりの空気のノイズの影響を受ける。
【００１８】
以下の例は、上記問題をより詳細に明らかにする。
【００１９】
コントラストを強調した（CE）磁気共鳴乳房X線撮影（MRM）は、浸潤性腫瘍に対する鋭敏な感度によって有用で重要な検査である。それは乳房柔組織の密度と独立な高い陰性適中率を示す。最近のマルチセンターな研究は、画像解析技術に従って、感度と特異性が異なることを示す。１．５Tスキャナを用いて、９６、９３、又は８６％の感度と３０、５０、又は７０％の特異性がそれぞれ計算された。従って、MRMの主な欠点は、偽の確定診断が高い割合で残存していることである。
【００２０】
通常MRイメージングは、患者が腹臥位の状態で、製造業者の二重乳房コイル（米国の場合、一重乳房コイル装置も多い）を用いて１．０又は１．５テスラの撮像装置で実施される。このプロトコルはあらゆる方向において３０から１８０秒の間の収集時間で、少なくともダイナミックコントラスト強調されたＴ１加重したシーケンスを含む。gadopentetate dimeglumine又は他の常磁性MR造影剤のボーラス注入の前に、少なくとも一回測定が実施され、その後に数回の測定が実施される。
【００２１】
乳房組織の領域への実際の造影剤の摂取を計算するために、造影剤の使用前の各ボクセルの信号強度は、造影剤の使用後の信号強度から差し引かれなければならない。MR画像のスライス厚さによる部分容積効果と同様に、呼吸と心臓の鼓動によって生じる両乳房の最低限の動きによって、異なる時間に撮影された同じ画像位置のボクセルは、乳房組織の正確に同じボリュームを示さない。従って（腫瘍の場合を除いて）引き算した画像は完全に黒くはならない。小さな移動の効果は、乳房の外側境界で最も良く現れる。脂肪組織の信号強度は大きく、周囲の空気の信号強度はほとんどゼロであることから、非常に小さな移動が脂肪を先の画像の空気の位置に置換する。脂肪の大きな信号強度から空気のゼロの信号強度を引き算することで、CM用途における大きな信号上昇をシミュレートする。結果として、引き算画像が、画像化された乳房の少なくとも部分まわりの強い造影剤の摂取を表す太く白い線を示す（図２−乳房の白い外側境界が表現されている）
【００２２】
デカルト座標系によって定義された三次元又は容積測定用の画像収集を考慮すると、ｘとｙ方向のあらゆる移動に加えて、（胸に向かっての）ｚ方向の幾らかの移動が常に見られる。これは、一連の検査中の患者の弛緩によって生じる。
【００２３】
アーチファクトを低減するために、乳房は通常何らかの圧縮方法（製造業者に依存する）によって乳房コイルの内側に固定される。この事実にもかかわらず、最小限の移動は常に発生する。従って、乳房を固定しているにもかかわらず、画像には常にｘ、ｙ、ｘ方向に僅かな動きアーチファクトが見られる。この動きアーチファクトは、明るい画素として引き算画像の中で視覚化される。動きアーチファクトがなければ、腫瘍性組織が占める領域を除き、引き算画像は絶対的に黒くなる。
【００２４】
周囲の空気のノイズ入りの信号内で、一定数のランドマークが乳房の外側で検出される。各特徴とみなされるものは追跡しなければならず、コンピュータに不必要な負荷を発生させるという事実の他にも、ノイズは変量効果であり、ノイズの中で見つけられたすべての特徴又はランドマークは、二番目の又はそれに続く画像において対応するランドマークを有さない。しかし実は、ノイズの中から発見されたランドマークに対応する二番目の、又はそれに続く画像又は画像セットにおけるランドマークを調査するアルゴリズムは、一定量のランドマークを第二画像又は画像セット内に再割付けする。これが、画像の不正確な登録につながり、結果に影響を与える。
【００２５】
現在知られている方法によれば、「特徴選択」は以下のように実行できる。
Ｂ１）第一画像又は第一画像セットの各ピクセル又はボクセル周りに隣接するピクセル又はボクセルを定義する。ここで、前記隣接するピクセル又はボクセルは、限られた数のピクセル又はボクセルからなり、単一画像の場合には二次元隣接部又はボクセルが選択され、断面画像セットの場合には、三次元隣接部又はボクセルが選択される。
Ｂ２）各ピクセル又はボクセルに対して、前記隣接部又はボクセルのすべてのピクセル又はボクセルの平均信号強度を定義する。
Ｂ３）平均信号強度閾値を定義する。
Ｂ４）ステップＢ２で各隣接部又はボクセルに対して定義した平均信号強度を、前記平均信号強度閾値と比較する。
Ｂ５）前記隣接部又はボクセルの平均信号強度がステップＢ４の閾値より大きい場合に、そのピクセル又はボクセルは追跡すべき特徴として定義され、追跡すべき特徴のリストに追加される。
【００２６】
平均信号強度閾値は、調整可能であり、変更可能である。経験的又は実験的基準は、例えばテストセットやサンプル画像に対しこの方法を実行し、異なる閾値に対して得られた結果を比較することによって、このような閾値を求めるために適用可能である。
【００２７】
自動特徴選択が、例えば二次元画像又は二次元の場合の一般化としての三次元画像に対していわゆるLukas＆Kanadeアルゴリズムを応用することによって、現在用いられている。良く知られたLukas＆Kanadeアルゴリズムによる自動特徴選択法は、Lucas BD,  Kanade T. An Iterative Image Registration Technique with an Application to Stereo Vision, IJCAI81; 674-679 に非常に詳細に記載されている。
【００２８】
三次元画像や三次元データボリュームを表す断面画像セットに用いられる自動特徴選択は、Lukas＆Kanadeの二次元特徴検出アルゴリズムの実質的な拡張であって、以下により詳細に説明する。
【００２９】
本明細書と請求項において、「特徴」という用語は「有効に追跡可能なランドマーク」という用語と等価である。
【００３０】
追跡すべきランドマーク又は特徴を自動的に選択するための上記方法を用いて、対象の画像の登録が実行可能である。登録のおかげで、動きの差異、又は動きアーチファクトを補償することができる。しかし、登録された画像シーケンスは、十分な外観を有さない。上記方法で画像を処理した後、画像はシャープさが減ったように見え、又は他の効果が見える。
【００３１】
さらに、有効なランドマーク又は特徴の自動選択のための方法は、手作業で、即ち人間オペレータの技能に基づいて実行されるのを避けて、客観的な基準に従って、画像内で追跡されるべき有効なランドマークを選択することを可能にする有用なツールではある。しかし、有効なランドマークの識別は、画像の各ピクセルを目標ピクセルとして考慮し、目標ピクセルの周囲の一定数のピクセルを備えた各目標ピクセルまわりの窓を定義することによって求められた勾配行列の固有値の関数であるところの特定の基準にのみ委ねられている。
【００３２】
この方法ステップは、一種の境界検出であって、境界の特徴を有するとして識別された目標ピクセルの中で、幾つかのみが登録プロセスにおいて有効に追跡可能なランドマークとして選択される。
【００３３】
本発明の目的は、特に二次元又は三次元の生物医学的診断画像に対して画像を登録する方法を提供することに存する。この方法は、二次元及び三次元で動きを追跡可能にすることによって、また特にデジタル画像又はデジタル処理された画像内の対象領域の外側のノイズによるアーチファクトを低減することによって、また隣接する組織タイプの信号強度における高低から独立な、乳房全体にわたりランドマークを広げることによって、既知の登録方法の欠点を効率的に克服している。
【００３４】
本発明は、特に画像内の追跡すべき有効なランドマークを選択するステップを向上させることを目的としている。
【課題を解決するための手段】
【００３５】
本発明による画像を登録するための方法は、
ａ）組織又は組織領域又は解剖した領域のＭＲＩ，超音波又はレントゲン写真の走査により得られた第一のデジタル又はデジタル処理された画像又は断面画像セットを提供し、
ここで前記画像はピクセル又はボクセルの二次元又は三次元アレイによって形成され、
同じ組織又は組織領域又は同じ解剖した領域のＭＲＩ，超音波又はレントゲン写真の走査により得られた少なくとも第二デジタル又はデジタル処理された画像又は同じ解剖した領域の断面画像セットを提供するステップと、
ｂ）ランドマーク又は特徴として設定された一定数のピクセル又はボクセルを選択し、追跡すべき前記特徴のリストを生成することによって、第一画像又は画像セットの中に、一定数のランドマーク、いわゆる特徴を定義するステップと、
ｃ）特徴として選択された各ピクセル又はボクセルに対する第一から第二画像又は画像セットのオプティカルフローベクトルを求めることによって、第一から第二画像又は画像セット間で、特徴としての各ピクセル又はボクセルの位置を追跡するステップと、
ｄ）逆オプティカルフローベクトルを第二画像又は画像セットのピクセル又はボクセルに適用することによって、第一と第二画像又は画像セットを登録するステップと、からなる。
【００３６】
ここにおいて、
Ｂ１）第一画像又は第一の一連の断面画像の各ピクセル又はボクセル周りに隣接した、限られた数で構成されるピクセル又はボクセルを定義するステップと、
Ｂ２）各目標ピクセル又はボクセルに対し、
外観を記載する数値パラメータ、すなわち前記ピクセル又はボクセル窓の前記目標ピクセル又はボクセルとピクセル又はボクセルの部分とからなるいわゆる数値外観パラメータの関数として計算され、又は、
前記窓のピクセル又はボクセルを表現する数値パラメータのマトリックス又は数値パラメータの前記マトリックスの変形で構成される一つ以上の特徴パラメータの関数として計算される、一つ以上の特徴のあるパラメータを定義するステップと、
Ｂ３）前記目標ピクセル又はボクセルの特徴的なパラメータの関数としての有効に追跡可能なランドマーク又は特徴に存するピクセル又はボクセルを定義するステップと、
からなる自動追跡可能なランドマーク選択ステップが実行されることを特徴とする。
【００３７】
本発明の第一実施例によれば、目標ピクセル又はボクセルが追跡すべき有効なランドマーク又は特徴であるかを求めるための関数は、それぞれの値又は目標ピクセル又はボクセルの特徴パラメータと閾値との比較に存する。
【００３８】
閾値は、実験によって求められ、テストセット又はサンプル画像に対して方法を実行し、異なる閾値に対して得られた結果を比較して求められる。
【００３９】
目標ピクセル又はボクセルが追跡すべき有効なランドマーク又は特徴であるかを求めるための代替的な方法は、分類アルゴリズムの入力ベクトルであるところのベクトル形式の目標ピクセル又はボクセルの特徴パラメータを組織化することに存する。
【００４０】
コンポーネントとして目標ピクセル又はボクセルの一定数の特徴パラメータを備えたベクトルを処理するために、前記分類アルゴリズムは、前記タスクを実行するために具体的に訓練される必要がある。
【００４１】
従って、上記代替方法による本発明の方法は、
一定数の有効なランドマークが、有効に追跡可能なランドマーク又は特徴として識別された既知の場合の一定数の画像を提供するステップと、
上記ステップＢ１）とＢ２）に存する自動追跡可能ランドマーク選択ステップを適用することによって、一定数の画像の中で有効に追跡可能であると識別された前記ランドマークに対応して前記ピクセル又はボクセルに対する特徴パラメータセットを求めるステップと、
有効に追跡可能であると識別された各ランドマークに一意的に関連付けられ、成分として前記有効に追跡可能なランドマークと一致する前記ピクセル又はボクセルの特徴パラメータを備えたベクトルを生成するステップと、
変数の所定の数値によって有効に追跡可能な質を表現し、
前記数値を、有効に追跡可能なランドマークと一致する各ピクセル又はボクセルをコード化するベクトルに関連付けるステップと、
有効に追跡可能なランドマークと一致する各ピクセル又はボクセルをコード化する各ベクトルが、分類アルゴリズムのためのトレーニングデータベースの登録を形成するステップと、
前記データベースによって分類アルゴリズムをトレーニングするステップと、
鍛えられた分類アルゴリズムの入力に前記目標ピクセル又はボクセルの前記特徴パラメータを備えたベクトルを提供することによって、目標ピクセル又はボクセルの有効に追跡可能なランドマークの質を求めるステップと、を備える。
【００４２】
改善した本発明によれば、前記トレーニングデータベースは、有効に追跡可能ではないと知られたランドマークと一致する画像のピクセル又はボクセルのベクトル形式に組織化され、有効に追跡可能でないランドマークの質を、追跡可能なランドマークの質を表現する前記変数の値とは異なる変数の所定の数値によって表現する特徴パラメータも備え、一方で有効に追跡可能でないランドマークの質を表現する前記数値を、既知の有効に追跡可能でないランドマークと一致した各ピクセル又はボクセルをコード化するベクトルに加える。
【００４３】
分類アルゴリズムを用いた有効に追跡可能なランドマークを求めるためのこの種の自動的方法は、例えば適応度や他の同様な統計指標のような分類の信頼性を示す統計的尺度をも提供するという更なる長所を有する。
【００４４】
分類アルゴリズムとして、あらゆる種類のこれらアルゴリズムを用いることができる。分類アルゴリズムの特別な場合としては、いわゆるクラスタ化アルゴリズム又は人工の神経ネットワークがある。
【００４５】
目標ターゲットピクセル又はボクセルや、前記ピクセル又はボクセル窓のピクセル又はボクセルのそれぞれ、又は一部の外観を表現する数値外観パラメータの関数として、又は前記窓のピクセル又はボクセルを表現する数値行列又は前記数値行列の転換の一つ以上の特徴パラメータの関数として計算される特徴数値パラメータに関し、多くの異なる関数が多くの異なる組み合わせの中で使用可能である。
【００４６】
特徴又は有効に追跡可能なランドマークに対応する目標ピクセル又はボクセルの特徴パラメータを求めるための第一の関数は、前記隣接部又はボクセルのすべてのピクセル又はボクセルの平均強度値に存する。
【００４７】
目標ピクセル又はボクセルのもう一つの特徴パラメータは、ステップＢ１）で定義された前記窓のピクセル又はボクセル、即ち目標ピクセル又はボクセルに新設するピクセル又はボクセルの一定数を表現する数値行列の特徴パラメータに存する。前記特徴パラメータは、前記窓のピクセル又はボクセルの画像データを備えた数値行列の特異値である。
【００４８】
代替的に、又は上述の特徴パラメータとの組合わせにおいて、目標ピクセル又はボクセルのもう一つの特徴パラメータが、前記窓のピクセル又はボクセルを表現する数値行列の勾配行列の固有値に存する。
【００４９】
目標ピクセル又はボクセルのもう一つの特徴パラメータが、前記窓のピクセル又はボクセルを表現する数値行列のヘッシアン行列の固有値に存する。
【００５０】
代替的に又は上記開示されたものとの組合わせで提供された目標ピクセル又はボクセルの特徴パラメータのためのもう一つの選択が、前記窓のピクセル又はボクセルに対応するデータの数値行列のウェーブレット変換係数の一つ以上に存する。
【００５１】
この場合、いくつかのウェーブレットを基礎にした関数が、代替的に又は組み合わせとして用いるために選択可能である。
【００５２】
ウェーブレット変換のより詳細な記載は、付録１で与えられている。この付録は、インターネット上で利用可能なように公開されており、"Wavelet for Kids, A tutorial introduction"というタイトルで、Duke大学のBrandi Vidakovic とPeter Muellerが著している。この文献では、ウェーブレット理論が要約され、議論され、画像処理への幾つかの応用が開示されている。画像処理でのウェーブレットの実行を開示した章から明らかに、ウェーブレット分解はパラメータを得ることを可能にする。分解の各レベルに対して、ウェーブレット変換は、平均を表す一つの行列と、いわゆる詳細を表す3つの行列とを生成する。上記行列の一つ以上から、例えば行列の要素の平均を取る以外にも、行列の特異値を取ることなどのよって、幾つかのパラメータを抽出することが可能である。これらのパラメータのすべて又は幾つかは、目標ピクセル又はボクセルの特徴パラメータとして用いることができ、また、選択された窓に含まれる周辺ピクセル又はボクセルと関係する各目標ピクセル又はボクセルを表すベクトルの成分を形成するために用いることができる。
【００５３】
代替的に又は上記開示されたものとの組合わせで提供された目標ピクセル又はボクセルの特徴パラメータのためのもう一つの選択が、前記窓のピクセル又はボクセルを表す数値行列の自己相関変換係数の一つ以上に存する。
【００５４】
画像処理における自己相関は、一般的に画像テクスチャを特徴づけるための道具として用いられ、二つの画像の数学的評価から成る。この二つの画像は、異なる時刻に撮影されるか、又は第一画像から場所を移動した結果としての第二画像として生成される。自己相関は、これら二つの画像間の関係を測定する。この数学的評価は、目標ピクセル又はボクセルのコード化において、考慮されるパラメータの数を減らす可能性を提供する。
【００５５】
目標ピクセル又はボクセルのもう一つの特徴パラメータは、前記窓のピクセル又はボクセルを表す前記数値行列の同時生起行列（又はその特異値）に存する。
【００５６】
前記同時生起行列は、固定された空間的関係によって分離された一対のピクセルに対する中間調の二次元ヒストグラムである。この行列は、一対のピクセルの共有の確率分布を近似する。一般的にコントラストやエントロピーのような、テクスチャ測定を計算するために用いられる。
【００５７】
目標ピクセル又はボクセルの特徴パラメータの上記例によれば、各目標ピクセル又はボクセルは、前記窓のピクセル又はボクセルを表現する数値行列の固有値又は特異値及び／又は、前記窓のピクセル又はボクセルを表現する数値行列の勾配行列又はヘッシアン行列の固有値及び／又は、ウェーブレット変換係数の一つ以上の係数及び／又は、自己相関変換の一つ以上の係数及び／又は、同時生起行列の一つ以上の入力又は特異値及び／又は前記窓のピクセル又はボクセルと前記目標ピクセル又はボクセルの平均強度値、の二つ以上の組み合わせ又はサブコンビネーションによって記載される。
【００５８】
既に開示したように、前記特徴パラメータは、前記特徴パラメータのそれぞれに対し、又はその組合わせに対し閾値を定義することによって、画像における目標ピクセル又はボクセルが、有効に追跡可能なランドマークとして考慮可能かどうかを判断するために用いることができる。
【００５９】
もう一つの基準は、前記特徴パラメータをベクトル形式に組織化することに存し、そこで各パラメータは前記ベクトルの成分であり、前記ベクトルの係数を求め、この値を前記係数に対する閾値と比較することで、組織化する。
【００６０】
代替的に、前記ベクトルは既に上で開示した分類アルゴリズムによって処理することもできる。
【００６１】
各ピクセル又はボクセルと、選択されたピクセル又はボクセル窓によって定義された選択された隣接部又はボクセルによって記載された数値間の関係は、単に各ピクセル又はボクセルの外観を表す値に存するオリジナルの数値行列の前記異なる変換の前記特異値及び／又は固有値によって、及び／又は、ウェーブレット変換の一つ以上の係数及び／又は自己相関変換の一つ以上の係数及び／又は同時生起行列の一つ以上の入力又は特異値
によって要約される。そして、各ベクトルにおいて、この関係は強調表示に特に適した、又は選択された窓内の画像のピクセル又はボクセル間の関係の特定の種類に感受性のある数値によって定義される。
【００６２】
従って、ピクセル又はボクセルに対する有効に追跡可能なランドマークの質を求めるためのこれら特徴パラメータの一部又は全てを考慮することは、このピクセル又はボクセルと残りの画像との関係に関しより良くより完全に考慮する助けになる。
【００６３】
すでに述べたように、上記方法は磁気共鳴画像（ＭＲＩ）又はＸ線ＣＴにおける同一の又は異なる種類のスキャナで撮影された二枚の診断画像又は3次元画像ボリュームを登録するのに用いることができる。
【００６４】
上述の登録方法は、例えばＭＲＩや超音波画像や特にコントラストを強調した磁気共鳴マンモグラフィー（ＭＲＭ）といった造影剤強調診断画像との組み合わせによっても提供される。この方法において、乳房組織領域の実際の造影剤摂取を計算するために、造影剤使用前の各ボクセルの信号強度が、造影剤使用後の信号強度から差し引かれる必要がある。このように、画像化された組織に造影剤の存在しない時に撮影された第一画像と、画像化された組織に造影剤が存在しかん流した第二の時間又は後の時間とを登録することは、登録が実行されなかった場合に患者の動きにより引き算画像に現れるアーチファクトを非常に成功裏に抑制するのに有用である。前記二枚の画像の引き算は、造影剤が溜まった組織領域に対応するピクセル又はボクセルの信号強度を高めることに繋がる。造影剤の溜まっていない他の領域は、低い信号強度のピクセル又はボクセルであって、特に暗い灰色か黒に見える。
【発明を実施するための最良の形態】
【００６５】
画像の登録は、特に医療診断の分野において、幾つかの種類の調査の実行を可能にする重要なプロセスである。
【００６６】
いわゆる画像化された領域のかん流挙動の測定において、登録が重要であることが以下の例により明らかになる。この場合、異なる時刻に撮影された身体上の領域の一連の画像が得られる。最初の画像は、造影剤が身体上の領域にない状態で撮影される。続く画像は、造影剤が身体上の領域に達した状態で撮影される。造影剤は、身体上の領域の周辺組織に対して血管やリンパ管の流れの応答を強調し、この流れがより良く見えるようにする。腫瘍性組織や他の病巣のような病変によって影響された組織は、強化された血管新生によって特徴付けられるので、画像化された組織のかん流挙動の測定は、この病理組織を個性化し、画像においてより良く見えるようにレンダリングするための手段である。画像化された身体上の領域の幾つかの領域における強化された血管新生は、この領域において造影剤が他の領域より集中するという結果をもたらす。また、造影剤は投与した時から血管やリンパ管のシステムによって組織に拡散されるまでに時間を要するので、強化された血管新生は、この強化された血管新生が存在する領域において、造影剤をより迅速に供給し、同時に、画像化された身体上の領域の通常に血管新生された領域に比較して、造影剤の最大ピークに達した後に、造影剤の集中のより迅速な減少を提供する。
【００６７】
比較的長い時間間隔に広がった画像の時間系列が得られるので、画像化された身体上の領域の動きは、一枚の画像とそれに続く画像間で生じることができる。この動きは、患者の動きによって、また例えば心臓の鼓動や呼吸といった生理上の理由によって決定される。
【００６８】
かん流挙動は、造影剤の画像化された身体上の領域への投与後に得られた時間系列の各画像から前記画像化された身体上の領域に造影剤の存在しない時点での画像を引き算することによって求められる。図２に示すように、造影剤が存在する領域が明るいままであるのに対して、この引き算が実質的に造影剤のかん流の生じない領域のピクセルを黒かほとんど黒色にするという事実のほかに、（上で定義した）各引き算画像に対して平均強度を求め、グラフに一連の画像の各引き算画像の平均強度をプロットすることで、造影剤が検出された領域のかん流挙動の尺度である、いわゆるかん流曲線を引くことができる。
【００６９】
診断に関するフォローアップ又は治療が施された際の画像化された領域の病状の進展の
フォローアップを実行するために、画像化が異なる時間間隔で実行された場合に、登録も重要である。ここで、画像を撮影するための装置に対して患者はいつでも位置決めしなおす必要があるので、先のかん流調査の場合よりも、一つの画像とそれに続く画像間の時間間隔はより長くなり、一つの画像とその他の画像間の動きは、よりドラマチックになる。
【００７０】
他の場合において、異なる画像化技術で撮影された身体上の領域の画像は登録される必要がある。相関画像又は組み合わせ画像がこの対象に関するより完全で詳細な情報を、又はより良いより明確な組合わせ画像を提供できるように、異なる画像技術はしばしば身体上の領域の画像化された対象の異なる様相や質を強調する又は暴露するので、これは重要である。
【００７１】
強固な構造が画像化される場合には、登録は一連の画像の一枚の他に対する単純な平行移動又は回転平行移動とみなすことができる。画像化された対象が柔組織からなるか、柔組織に埋め込まれている場合、画像化された対象や画像化された領域に歪みが生じるので、登録にはこの効果も考慮することが必要である。
【００７２】
診断画像の時間系列を登録する一般的なケースは、かん流挙動測定法を適用することで実行される乳房中央を調査する場合である。
【００７３】
図１に患者の乳房を固定するための一般的な装置を示す。通常、MR造影は１．０又は１．５テスラの撮像装置に製造業者の二重乳房コイル（米国の場合、単一乳房コイルである場合も多い）を用いて行われ、その間患者はうつ向けになっている。このプロトコルはあらゆる方向において３０から１８０秒の間の収集時間で、少なくともダイナミックコントラスト強調されたＴ１加重したシーケンスを含む（個々のプロトコルは異なっていてもよいし、追加的に脂肪抑圧を用いることもできる）。0.1 mmol/kgから0.2 mmol/kgの gadopentetate dimeglumine又は他の常磁性MR造影剤のボーラス注入の前に、少なくとも一回測定が実施され、その後に数回の測定が実施される。
【００７４】
少なくともMRI画像の最初のボリュームが造影剤が肘正中静脈に注射されて乳房組織に到達する前の時間ｔ０に撮影される。少なくとも画像のもう一つのボリュームが造影剤注射後の時間ｔ１に直接撮影される。通常、時間遅れのある連続したボリュームが撮影される（例えば、造影剤投与後、６回以下の同一ボリュームが連続する）。
【００７５】
図２に、断面画像の３次元ボリュームから取り出された２枚の画像を示す。一枚は造影剤投与前の時刻ｔ０に撮影されたものであり、もう一枚は造影剤投与後の時刻ｔ１に撮影されたものである。引き算画像、即ち第二画像から第一画像をピクセル毎に引き算した結果の像も、同様に示す。ここで登録は行われておらず、例えば矢印で示された領域のように、造影剤の摂取されていない領域のピクセルは、輝いたままであることは明らかである。
【００７６】
図３に、従来の登録方法とも共通の、登録プロセスの主なステップを示す。
【００７７】
立方体で表現された画像ボリュームが異なる時刻ｔ０とｔ１に画像化される。ここでは三角と円で表現された二つの対象が、ランドマークとして選択される。動きのせいで、画像ボリュームｔ０での二つの対象の位置は、ｔ１での画像の同じ対象の位置とは一致しない。
【００７８】
従って、第一ステップは、画像の中で追跡すべきランドマークを個別化することである。三角と円で表現された追跡すべきランドマークがひとたび個別化されると、いわゆる動きベクトルの測定を可能にする追跡ステップが実行される。
【００７９】
これらベクトルによって、画像化された身体構造上の領域の動きを大域的運動として、また画像化された身体構造上の領域の一定の領域の局部的運動として表現する、画像中に拡がるベクトル場を生成することが可能である。逆動きベクトル場を画像に適用すると、第一画像の時刻ｔ０と第二画像の時刻ｔ１との間で生じる動きの効果は、補償される。乳房のかん流試験の例におけるこの登録の効果は、図１９で見ることができる。ここで、左側の連続画像は、時刻ｔ０とｔ１で撮影された乳房画像が登録されていない場合に関し、下側の画像は、図２の一つと対応する引き算画像である。右側の連続画像は、時刻ｔ０とｔ１に撮影された２枚の画像を示すが、時刻ｔ１に撮影された画像は、登録方法で処理済みのものである。引き算画像に対応する最後の画像において、左側の引き算画像ではピクセルが輝いていた多くの領域が消失し、又は対応するピクセルが完全に黒か灰色の様相を呈していることが明らかに見ることができる。
【００８０】
本発明の目的は、以下のステップでもクリティカルな、画像の登録プロセスの第一ステップを改良することである。画像内の信頼できるランドマークを定義することは、より良い登録結果を得るためにクリティカルである。したがって、有効に追跡可能なランドマークとなる可能性の非常に高い画像のピクセル又はボクセルの特徴的性質を求めるためのより精密で完全な方法が本発明の目的である。
【００８１】
図４を参照して、単純化のために二次元の場合について、デジタル又はデジタル処理された画像が、n, m = 1, ... , NであるピクセルアレイP(n, m)によって形成される。各ピクセルは、一定の外観を有するモニター又は画面又は印刷された紙上の小さなドットである。ドットの外観は、数値データによって表現可能である。いわゆるグレースケール画像において、各ドットは前記グレースケールの一定のグレーレベルに対応した一定の強度を有する。カラー画像の場合には、ピクセルの外観を数値データで完全に表現するために、より多くのパラメータが通常必要である。カラー画像におけるピクセルの外観を一義的に表現するために、幾つかのシステムが定義されている。そのうち一つのシステムは、いわゆる既知のHSV（色調、彩度、明暗度）又はいわゆるRGB（赤、緑、青）システムである。もしカラー画像を処理した場合、これらの値が対応する数値データに置き換える必要があることが当業者には自明であるため、本単純化した例においては強度I(n, m)のみが示される。
【００８２】
従って、図５に示されたn, m = 1, ... , 10の数値データのアレイI(n, m)は、図１のピクセルアレイの概要図としてのn, m = 1, ... , 10のピクセルP(n, m)で構成される画像に対応し、そして、数値データのアレイは行列である。
【００８３】
各ピクセルP(n, m)は、例えばグレースケールデジタル画像のピクセルのグレーレベルのように、数値的にその外観を表現した数値データI(n, m)に一義的に関連する。
【００８４】
隣接する断面（二次元）画像の一定数によって形成されるボリューム画像として、三次元画像を定義することができる。この場合、画像はいわゆるボクセルで形成され、そしてアレイは画像データの数値行列も三次元であるように、三次元である。実際には、断面画像セットは定義された厚みを有し、従ってこれら画像は、純粋な二次元画像よりもよりスライスされている。この観点から、常にボクセルと称することがより正確である。それにもかかわらず、モニター画面上の断面画像の視覚的表現は、純粋に二次元画像であって、スライスの厚みは考慮されず、この場合、ピクセルが画像ドットの適切な定義である。
【００８５】
図５、６、７に示すように、画像の目標ピクセルは、その成分がその目標ピクセル周りのピクセルに関する情報を含んだベクトルによってコード化可能である。通常、前記周辺領域は目標ピクセルを中心とした（2? + 1）(  (2? + 1)ピクセルの次元の窓から成る。ここで、?は(? = 1, 2, ..., N)でユーザーによって任意に定義される整数であり、ここでは「ステップ」と称す。
【００８６】
周囲のピクセルで形成された窓の上記定義は、勾配ｘの周辺ピクセルの定義と等価であり、ここでｘは整数であって、この整数はステップにおける目標ピクセルから隣接部までの距離を示す。目標ピクセルとして前記窓の中心にあるピクセルを考慮すれば、勾配１の周辺ピクセルを備えた前記窓は、目標ピクセルに直接隣接したピクセルのシェルであり、勾配２の周辺ピクセルは、目標ピクセルの周囲の二番目に近いピクセルのシェルを備え、画像を形成するピクセルアレイの各方向において、勾配１は目標ピクセルから１距離ステップに対応し、勾配２は、目標ピクセルから２距離ステップに対応する。この定義は、ボクセルアレイで形成された３次元画像に対して準用する。
【００８７】
前記窓の最小寸法は、３×３寸法のピクセルアレイに存し、その中心ピクセルが目標ピクセルであり、上で勾配として定義されたステップ？は、１に等しい。例えば５×５又は７×７ピクセル窓といった、より大きな窓を選択することもできて、それぞれのステップ又は勾配は、？＝２と？＝３である。簡単にするために、本例においては、目標ピクセルが中心にある３×３ピクセルアレイに対応した窓を選択する。
【００８８】
画像のピクセルを表現する数値データ行列を参照して、図５にこの窓を示す。窓Wは九つのピクセルを備え、そのうち一つが目標ピクセルである。図５に示された窓は、目標ピクセルP(2, 2)を中心とし、ピクセルP(1, 1), P(1, 2), P(1, 3), P(2, 1), P(2, 2), P(2, 3), P(3, 1), P(3, 2), P(3, 3)を備えており、各ピクセルは対応する数値、即ち前記ピクセルの信号強度に関連する数値パラメータI(1, 1), I(1, 2), I(1, 3), I(2, 1), I(2, 2), I(2, 3), I(3, 1), I(3, 2), I(3, 3)で表現される。
【００８９】
前記目標ピクセルと画像内の特定の周辺領域との関係が考慮されるように、目標ピクセルP(2, 2)は、窓Wの隣接画素に関する情報も用いてコード化される。前記目標ピクセルP(2, 2)と上で定義した周辺ピクセルとの関係を表すベクトルの成分として、前記ピクセルの強度は、目標ピクセルの強度とともに測定される。
【００９０】
前記ベクトルは、図６に図式的に示される。明確に表すために、目標ピクセルを表現する各ベクトルは、九つの成分を有する。
【００９１】
図７に、目標ピクセルをコード化するための前記ベクトルの一般的表現を示す。この場合、目標ピクセルはその強度によって表現され、強度I(n, m)を有するピクセルP(n, m)として定義される。
【００９２】
図８に図７のベクトルを参照する窓Wと、いわゆる勾配行列におけるその変換を示す。
【００９３】
当業者には自明のことであろうが、二方向が定義可能で、対象の窓における各ピクセルに対し、この二方向において勾配が評価可能であるように、ピクセルを表現する数値アレイと、この場合窓W内のピクセルの強度データアレイI(m-1, n-1), I(m-1, n), I(m-1, n+1), I(m, n-1), I(m, n), I(m, n+1), I(m+1, n-1), I(m, n), I(m+1, n+1) とは、二次元オブジェクトである。
【００９４】
勾配行列は、固有値?p(p = 1, 2)で、この場合?1?  ?2で表現可能なように、対角化可能である。
【００９５】
選択された窓Wのピクセルの強度Iのオリジナル行列は、さらに加工可能であり、前記オリジナル行列に対し、いわゆるヘッシアン行列が計算可能である。また、この場合ヘッシアン行列は、その固有値?pによって表現可能である。
【００９６】
選択された窓の強度値Iの3×３行列を考慮した場合、通常、前記行列は対角化可能でなく、固有値は上述のように有意義なものではない。これを考慮して、より一般的な場合において、図８に示すように選択された窓のピクセルの強度Iのオリジナル行列は、いわゆる特異値siで表現可能である。
【００９７】
選択された窓に対応する強度行列の特異値又は代替的に、勾配行列又はヘッシアン行列の固有値と組合わせて用いることで、対応した目標ピクセルを一義的にコード化するためのベクトルを生成することが可能である。前記ベクトルの成分は、選択された窓のピクセルと、選択された窓のピクセルに対応する強度値の前記行列を処理することで得られる勾配行列又はヘッシアン行列の固有値とに対応した強度値の行列の前記特異値に存する。
【００９８】
上で開示したように、図示された例は、選択された窓のピクセルに対応する強度値のオリジナル行列の変換の特別な場合に限定されるが、その他の変換も代替的に又は組み合わせとして適用可能である。例えば、ウェーブレット分解が実行可能であり、目標ピクセルのコード化ベクトルの成分として、それらのすべて又は少なくとも一部に対し、平均値と詳細値が利用可能である。
【００９９】
代替的に、又は組合わせにおいて、選択された窓のピクセルに対応する強度値のオリジナル行列の自己相関変換が利用可能であり、目標ピクセルのコード化ベクトルの成分として、それらのすべて又は少なくとも一部に対し、得られたパラメータが利用可能である。
【０１００】
ウェーブレット変換も前処理ステップとして利用可能である。得られた数値データの行列は、対応する勾配行列及び／又は対応するヘッシアン行列を用いて処理可能であり、その特異値又は固有値は、代替的に又は組合わせにおいて求めることができる。
【０１０１】
前記窓のピクセル又はボクセルを表現する行列のいわゆる同時生起変換に存する、ピクセルの窓を表現するオリジナル行列の別な変換が利用可能である。この場合、この同時生起行列の特異値又は固有値が、目標ピクセルを表すパラメータとして利用可能である。
【０１０２】
一定数の隣接部からなる選択された窓における目標ピクセルを特徴付けるための幾つかのパラメータの上記記載は、単に単純化のために二次元の場合に限定されているが、上述の選択は、三次元の場合にも適用可能であることを再び強調しておく必要がある。
【０１０３】
上記結果として、画像の各目標ピクセル又はボクセルは、上記パラメータの一つ以上によって表現可能である。これらのパラメータは、暗示的に目標ピクセル又はボクセルと選択された窓のピクセル又はボクセルとの関係を示し、それぞれ選択された窓に対応する画像の一定の特徴を現すことができる。従って、これらパラメータの二つ以上、又はすべてを、ピクセル又はボクセルが画像内の有効に追跡可能なランドマークと一致するかどうかを求めるための数値としてみなすことで、非常に信頼性のあるランドマークの識別を提供することになる。
【０１０４】
目標ピクセル又はボクセルと関連付けられた少なくとも上述のパラメータの一部を利用して有効に追跡可能なランドマークを選択するための第一の方法によれば、例えば閾値の数値の定義を実験データに基づかせることで、閾値を求めることができる。
【０１０５】
いったん閾値が上述の各目標ピクセル又はボクセルに対するベクトルの成分である各パラメータに対して定義されると、前記パラメータのそれぞれ又はパラメータのグループの組み合わせと、対応する閾値とを比較することによって、目標ピクセル又はボクセルが有効に追跡可能なランドマークと対応するかどうかを求めるための基準を設定する関数を定義することができる。
【０１０６】
非常に多数の関数が異なる基準に対して定義可能であり、異なる可能性が考慮される。
【０１０７】
最も単純だが非常に時間を要する解は、ピクセル又はボクセルをコード化するベクトルの成分である各パラメータの値に対して閾値を定義することであり、各パラメータを関連する閾値と比較して、パラメータが閾値より高いか低いかを調べることである。各パラメータに対する比較結果は、パラメータが関連する閾値より低いか高いか、又は比較結果が有効に追跡可能なランドマークかそうでないかを、例えば０か１かの数値によって表現可能である。有効に追跡可能なランドマークに対して、特徴パラメータの幾つかは負の比較結果、即ち、前記特徴パラメータがランドマークが有効に追跡可能でないことを示す比較結果を有することもあるので、特徴パラメータと関連する閾値との各個別の比較結果の数値からなる成分を有するベクトルは、各特徴パラメータとその閾値との比較結果を要約する数値からなる包括的評価パラメータを得るために、さらに処理可能である。
【０１０８】
このことは、特徴パラメータの一つ又はサブセットに対する比較結果がランドマークが有効に追跡可能ではないと示しているにもかかわらず、有効に追跡可能なランドマークを識別することを可能にする。
【０１０９】
上記特徴パラメータと関連する閾値との比較結果の包括的評価は、各比較結果の関数として、又は特徴パラメータの関数として包括的評価パラメータを算出することで求めることができる。
【０１１０】
例えば包括的評価パラメータは、特徴パラメータのセットの平均値として算出することができる。この包括的評価パラメータは、包括的閾値と比較することができる。この場合、包括的評価パラメータに対する包括的閾値は、各特徴パラメータと関連する各単一閾値の関数として、例えば閾値の平均値として、求めることができる。
【０１１１】
包括的評価パラメータを求める他の方法は、特徴パラメータにより形成された成分からなるベクトルを定義し、次に前記ベクトルの係数、即ちその長さを計算することである。包括的評価パラメータである前記ベクトルの係数と比較される包括的閾値は、各単一特徴パラメータと関連する閾値によって形成された成分からなる閾値ベクトルの係数として算出することができる。
【０１１２】
また別な方法は、対応する特徴パラメータの対応する閾値との比較結果が有効に追跡可能なランドマークを示すか否かに依存する、例えば０又は１といった２値を前提とする、比較評価変数によって、各特徴パラメータの対応する閾値との比較結果をパラメータ化することである。
【０１１３】
このように、特徴パラメータと対応する閾値との比較結果が数値的に表現可能であり、そして、包括的結果は、その成分が比較評価パラメータ値であるベクトル形式で表現可能である。包括評価パラメータが二つの離散値０と１のうち一つを有することで定義される上記例に関し、比較評価ベクトルの成分は、０と１の数値のセットから成る。
【０１１４】
上記オプションと同様に、包括的評価パラメータは、この比較評価ベクトルの係数として、又はその成分の平均値として計算可能である。
【０１１５】
上記方法の別の特徴として、特徴パラメータの一つ以上のサブセットが定義可能であり、それは前記サブセットの特徴パラメータの関数としての二次特徴パラメータを計算するために用いられる。
【０１１６】
あらゆる種類の関数が、二次特徴パラメータを求め、又は算出するために利用可能である。
【０１１７】
二次特徴パラメータを求めるための特別な方法によれば、特徴パラメータの各サブセットは、目標ピクセル又はボクセル又は前記ピクセル又はボクセル窓のピクセル又はボクセルのそれぞれ又は一部を表す数値的外観の同じ種類の関数に関する特徴パラメータと、
窓のピクセル又はボクセルを表現する数値行列又は同一の数値行列の同一の変換についての一つ以上の特徴パラメータの同じ種類の関数に関する特徴パラメータと、を備える。
【０１１８】
これによって、特徴パラメータのセット全体が、各サブセットの特徴パラメータの関数として求められた二次特徴パラメータによって要約される特徴パラメータのサブセットに分割される。
【０１１９】
実際の例によれば、目標ピクセル又はボクセルを備えた一定の窓内のピクセル又はボクセルの信号強度と、前記目標ピクセル又はボクセルの信号強度の値は、特徴パラメータのサブセットを形成することとして定義される。これらの値を要約する二次特徴パラメータは、目標ピクセル又はボクセル又は窓に含まれるピクセル又はボクセルの信号強度の平均信号強度として計算可能である。
【０１２０】
同一の平均値関数が、窓のピクセル又はボクセルと目標ピクセル又はボクセルを表す数値画像データの行列の特異値から、及び／又は勾配行列の固有値から、又はヘッシアン行列の固有値から、又は自己相関行列の特異値から、二次特徴パラメータを計算するために利用可能である。
【０１２１】
上記ステップは、一定数のパラメータを有する一次特徴パラメータセットを、より少ない数のパラメータを有する二次特徴パラメータセットに減じることに存する。
【０１２２】
上述のように一次特徴パラメータとして定義可能な特徴パラメータのセットの評価と同様に、二次パラメータは、ピクセル又はボクセルが有効に追跡可能なランドマークに対応するか、有効に追跡可能でないランドマークに対応するかを決定するために用いることができる。前記二次パラメータは、一次特徴パラメータに適用したものと同じステップによって処理可能である。従って、二次パラメータをその二次閾値と比較するステップが実行可能であるのと同様に、二次閾値は定義可能である。またこの場合に、包括的閾値と比較される包括的比較評価パラメータも求めることができる。
【０１２３】
また上述の基準に従って特徴パラメータをグルーピングすることによって得られた前記二次特徴パラメータの場合において、単一二次特徴パラメータの組み合わせと対応する閾値の組み合わせは、例えば加重和又はその他の加重した組み合わせによって作成され、そして、そのように定義された包括的二次特徴パラメータと包括的二次特徴閾値のレベルとで比較が実行可能である。
【０１２４】
原理的に、一次特徴パラメータとその対応する閾値との比較を実行する上述のそれぞれの方法は、二次特徴パラメータのセットに対しても適用できる。特に、例えば二次特徴パラメータが閾値より大きいか小さいかによって０と１を取るような、各特徴パラメータとその関連する閾値との間の比較に関する異なる離散値を含む比較結果ベクトルに対して提供することにメリットがある。
【０１２５】
次に、包括的結果を和として、特に比較結果ベクトルの成分の加重和として求めることができる。
【０１２６】
パラメータ数が減ることによって、計算と評価がより簡単になることが長所である。さらに、二次特徴パラメータを求めるために選択された関数の種類が、正しくない又は統計的許容誤差の外にある一つ以上の前記パラメータの値の範囲を制限することを可能にする。
【０１２７】
比較のために包括的評価パラメータを求める場合において、ピクセル又はボクセルが有効に追跡可能なランドマークかどうかを求める場合の妥当性に相応じて、各一次又は二次特徴パラメータを加重する加重値を定義することが可能である。図２０は、上述のごとく定義された前記ピクセル又はボクセルの特徴を考慮することによって、ピクセル又はボクセルが有効に追跡可能なランドマークと対応するか否かを評価するための上述のステップを表現した概要図である。
【０１２８】
図２０で示されたケースは、一次特徴としてここで定義された目標ピクセル又はボクセルの特徴パラメータの選択されたサブセットの関数としての、各目標ピクセル又はボクセルに対して二次特徴パラメータを求めた場合を指す。
【０１２９】
以下で説明することではあるが、前記処理ステップのメカニズムも、第一特徴パラメータに適用することができる。
【０１３０】
図４ないし８に関連して、先に説明したように、第一特徴パラメータのセットが求められる。目標ピクセルは、I(m, n)で表される。隣接するピクセルは、３×３ピクセル寸法を有する窓に含まれるものとして選択される。一次特性は、選択された窓内の隣接部の強度Iと、目標ピクセルの強度と、目標ピクセルの一つを含む選択された窓に含まれるピクセルの信号強度の行列の特異値ｓｉと、目標ピクセルの一つを含む選択された窓に含まれるピクセルの信号強度の行列の勾配行列の固有値と、目標ピクセルの一つを含む選択された窓に含まれるピクセルの信号強度の行列のヘッシアン行列の固有値と、で構成される。
【０１３１】
上述のように求められたあらゆる別の特徴パラメータは、図示された特徴パラメータに付加したり、置き換えたりすることができるので、図２０の例で選択された特徴パラメータの種類に、本発明は限定されるものではないことを強調する必要がある。
【０１３２】
一次特徴パラメータの前記種類のそれぞれは、特徴パラメータのサブセットとして定義される。これは、中括弧によって図２０に示される。二次特徴パラメータは、関数１として示された関数によって、一次特徴パラメータの各サブセットから求められる。対応する一次特徴パラメータのサブセットから二次特徴パラメータを求めるための関数は、各サブセット間で同一又は異なるものであってよい。
【０１３３】
関数の例は、二乗関数又は係数関数等の和の平均値関数又は加重関数に存することができる。
【０１３４】
他の種類の関数は、サブセットの一次特徴パラメータの少なくとも幾つか、又は一つだけ除いた全部を除去するために、閾値との比較を提供する。各サブセットに対し、二次特徴パラメータはそのように求めることができ、図２０において、この二次特徴パラメータは、Ｐ１，Ｐ２，Ｐ３，Ｐ４で示される。ここで、目標ピクセルが有効に追跡可能なランドマークと認められるか否かを求めるために、上述のすべてのステップが実行可能である。各二次特徴パラメータに対して個々の閾値を求めることができ、そして各二次特徴パラメータと対応する閾値との個々の比較が実行可能である。ピクセルが有効に追跡可能なランドマークか否かの決定は、上述の包括的比較評価プロセスによるこれらの比較結果に基づいて行われる。
【０１３５】
本実施例において、包括的二次特徴パラメータが図２０の関数２によって示される二次特徴パラメータの関数として求められる。関数２は、関数１又は関数１に対して示された例のうち一つと等しくすることができる。
【０１３６】
上述のすでに開示した方法のうち一つに従って求められた包括的閾値は、例えば各二次特徴パラメータに対する個々の閾値を求め、次に、関数２と同じ又は異なる関数によって、包括的二次特徴パラメータをもとめるために、この個々の閾値のセットから包括的閾値を求めることによって求められる。
【０１３７】
包括的二次特徴パラメータと包括的閾値との比較は、比較関数によって実行され、そして、比較結果はピクセルI(m, n)が有効に追跡可能なランドマークと一致するか否かの評価を与える。
【０１３８】
代替案として、全ての一次特徴パラメータから包括的評価パラメータを求めるために関数を選ぶことができ、そして、閾値から求めるための別の関数又は同じ関数が、各一次特徴パラメータに対して包括的閾値を定義することを認めなければならない。比較関数は、次に、この一次特徴パラメータに対する包括的閾値に適用され、また包括的評価特徴パラメータに適用される。
【０１３９】
好ましい実施例において、関数１と関数２として示された異なる関数が、線形結合であって、好ましくは前記一次又は二次特徴パラメータの加重線形結合であって、例えば加重和である。同様に、第二特徴パラメータに対する閾値と包括的閾値とは、各一次パラメータに対する閾値と、各二次パラメータに対する閾値それぞれの加重和として計算可能である。
【０１４０】
特別に長所を有する実施例は、以下のステップからなる。各一次特徴パラメータとそれらの閾値とを比較するステップと、各成分が各一次パラメータと対応する閾値との比較結果からなる比較結果ベクトルを定義するステップ。もし一次パラメータが閾値より小さいとき、値０とし、もし一次パラメータが対応する閾値より大きいとき、値を１とするように比較結果を定義するステップ。比較結果ベクトルを値０と１の連続で形成するステップ。比較結果ベクトルの成分を組合わせるステップ。特に比較結果ベクトルの成分の加重和を求めるステップと、そして、この結果を包括的閾値と比較するステップ。包括的閾値は、実験又は各一次特徴パラメータに対する閾値の正規化された加重和によって求めることができる。比較結果ベクトルの成分の加重和と包括的閾値との比較結果に関連して、ピクセル又はボクセル又はピクセル又はボクセルのグループが有効に追跡可能なランドマークであるかどうかを決定するステップ。
【０１４１】
代替的に、上述のように二次パラメータを形成するために一次パラメータがグループ化された場合に、同じステップが二次特徴パラメータに対して適用可能である。この場合、比較結果ベクトルは、その成分が各二次特徴パラメータと対応する閾値との比較結果として定義される。また、ここで数値、特に０と１が、二次パラメータが対応する閾値より高いか低いかという条件にそれぞれ対応する比較結果を表現するために選択可能である。次に、ピクセル又はボクセル又はピクセル又はボクセルのグループが追跡すべき良いランドマークであるか否かということに関連するピクセル又はボクセル又はピクセル又はボクセルのグループの有効性が、比較結果ベクトルの成分の加重和と包括的閾値との比較によって求められる。
【０１４２】
特徴パラメータは、前記特徴パラメータの種類又はその数学的意味に依存した異なる方法によっても評価可能である。
【０１４３】
従って、例えば一次特徴パラメータのリストに含まれるヘッシアン行列の固有値を考慮すると、前記ヘッシアン行列の固有値の値から、ピクセル又はボクセルの重要性を評価することができる。
【０１４４】
ウェーブレット変換係数を考慮することによっても、画像のピクセル又はボクセルの重要性を評価することができる。この技術は、文献「Wavelet-based Salient Points: Applications to Using Color and Texture Features, E. Loupias, N. Sebe, Visual' 00, pp. 223-232, Lyon, France」と「Image Retrieval using Wavelet-based Salient Points, Q. Tian, N. Sebe, M.S. Lew, E. Loupias, T.S. Huang, Journal of Electronic Imaging, Vol. 10, No.4, pp. 835-849, October, 2001」とにより詳細に開示されている。
【０１４５】
従って、上述の開示内容に加えて、特徴がピクセル又はボクセルのものかを評価するための閾値が、追跡すべき特徴又はランドマーク即ち画像のピクセル又はボクセルの一つに対応し、幾つかの特徴パラメータは、ヘッシアン行列の固有値とウェーブレット変換係数に関する代替的な方法によって評価することができる。
【０１４６】
図２１に、有効に追跡可能なランドマーク又は有効に追跡可能なランドマークでないものの質を評価する別の実施例を示す。ここで、いわゆる分類アルゴリズム、例えばクラスタ化アルゴリズム又は予測アルゴリズム、例えば人工の神経回路網又は統計アルゴリズムが用いられる。
【０１４７】
この実施例を実行可能にするために、第一に、有効に追跡可能なランドマークか否かが知られたランドマークに対応したピクセルのデータベースを構築する必要がある。ここで、i=1, ..., n の各ランドマークLiに対し、成分が前記ランドマークLiが特定された画像のピクセル又はボクセルの前記方法ステップによって求められた一次又は二次特徴パラメータであるベクトルが構築される。
【０１４８】
ベクトルの成分は、更にそれぞれ有効に追跡可能なランドマークか否かの質を示す二つの異なる値を推測可能な有効な変数を備える。
【０１４９】
前記データベースは、図２１に概略的に示される。ここで、Liはランドマーク、Viは対応する特性のベクトルであり、「有効性」は、ランドマークLi又はそれに一致するピクセルが有効に追跡可能か否かを表現するパラメータである。データベースのレコードは、通常の手段で分類アルゴリズムを鍛えるために用いられる。ひとたびアルゴリズムが鍛えられると、ピクセルI(m, n)の特徴パラメータは、それによって処理され、そして分類アルゴリズムは出力として有効な変数の値と、ピクセルI(m, n)が有効に追跡可能なランドマーク又は有効でない追跡可能なランドマークと一致するか否かの予測又は決定を与える。
【０１５０】
本実施例及び図２０の一つ前の実施例に従う評価プロセスの組み合わせも適用可能であることを強調しておくことは重要である。この場合、一次特徴パラメータを用いる代わりに、目標ピクセル又はボクセルと、鍛えられたデータベースのランドマークのための入力データベクトルを構築するために二次特徴パラメータのみが用いられる。
【０１５１】
上述の異なる評価ステップに対して、特に有効に追跡可能であると識別された二つのランドマークが互いに一定の距離にある場合に、別なランドマークを確認する又は廃棄するステップが実行可能であることも考慮されるべきである。これは、前記二つのランドマークの一次特徴パラメータの少なくとも一つ、又は一部又は全部、又は二次特徴パラメータの少なくとも一つ、又は一部又は全部を比較することによって、そして、選択関数を二つのランドマーク間の選択基準を表す比較結果に適用することによって、実行可能である。
【０１５２】
以下において、自動ランドマーク識別の例を図１１ないし１４を用いて示す。これらの図式的な例は、理解を容易にするために、二次元画像を用いている。三次元画像への方法の一般化は、ここでは簡潔に記述された数学的形式主義によって、概要が与えられる。
【０１５３】
三次元と、三次元ボリュームデータを表現する断面画像セットにおける自動特徴選択作業は、Lucas & Kanadeの二次元特徴検出アルゴリズムの実質的拡張である。
【０１５４】
時間の僅かに異なる間に記録された同一患者の二枚の画像ボリュームを、IとJと定義する。
【０１５５】
第一の画像ボリュームIにおいて、二つのボリューム間で追跡すべき特徴とみなされた一定数の選択されたボクセルが識別される必要がある。
【０１５６】
Iをオリジナルの第一画像ボリュームとし、Jを次のボリュームとする。ここで、対応する特徴は見出されていなければならない。
【０１５７】
第一ステップは、第一画像ボリュームにおける最小固有値の計算で構成される。これは、以下のように実行される。
【０１５８】
第一画像ボリュームは、三つの勾配ボリュームIx, Iy, Iz（画像ボリュームを表現するデカルト座標系のx-, y-, z方向）に変換され、x、y、ｚは各シングルボクセルの座標である。
【０１５９】
勾配ボリュームは、一般的にx方向に(2?x + 1)、ｙ方向に(2?y + 1)、ｚ方向に(2?z + 1)（ここで?x, ?y, ?z = 1, 2, 3, 4, ..., nボクセル）のサイズのいわゆる窓に含まれる隣接するボクセルに対する各ボクセルの強度勾配から算出される。
【０１６０】
目標ボクセルを中心とする隣接ボクセルの３×３×３アレイのサイズを考えると、３つの勾配ボリュームは、以下のように定義される。
ａ）ｘ方向の勾配ボリューム
【数１】
2009512527.tif 000003
ｂ）ｙ方向の勾配ボリューム
【数２】
2009512527.tif 000004
ｃ）ｚ方向の勾配ボリューム：
【数３】
2009512527.tif 000005
次のステップとして、各ボクセルに対し、三つすべての先に算出した勾配ボリュームを用いて、勾配行列Ｇが設定される。勾配行列は、以下のように生成される。
【数４】
2009512527.tif 000006
各ボクセルの各勾配行列に対して、以下のように勾配行列から最小固有値?minが算出される。
定義：Graphics Gems IV, Paul S. Heckbert, 1994, Academic Press, S 193-98
【数５】
2009512527.tif 000007
ここで
【数６】
2009512527.tif 000008
なぜなら、Ｇは中心モーメントの行列だからである。
勾配行列Ｇの固有値は、次式で算出される。
【数７】
2009512527.tif 000009
次に、固有値?mの最小値に対する閾値が定義される。
【０１６１】
次に、ボクセルが追跡可能な特徴を表現しているか否かを検討するために、以下の基準が適用される。
１．もし?mが閾値より大きい場合、ボクセルの実際の位置が、追跡すべき特徴のリストに    保存される。
２．もし？ｍが固有値の全ての最小値?mの最大値のパーセンテージより小さい場合、追跡    すべき特徴のリストから落とされる。
３．もし実際の特徴周りの定義された領域（調整可能、例えば３×３×３）において、他    の特徴が存在する場合、固有値のより小さな最小値?mを有する特徴が、追跡すべき特    徴のリストから落とされる。
４．もし特徴位置周りの3次元隣接ボクセルの平均信号強度(例えば? = 1)が、ゼロに近    い非常に低い値（例えば１０）より小さいか等しい場合、その特徴は、周囲の空気の    中で乳房の外側にあるとみなされ、追跡すべき特徴のリストから廃棄される。
５．最後に、特徴の定義可能な最大数のみがキープされる。もしそれ以上の特徴が追跡す    べき特徴のリストに存在する場合、より小さな最小値?mを有する特徴が、そのリスト    から落とされる。
【０１６２】
これらのステップは、二次元の場合について図１１ないし１４で視覚的に説明されている。この場合、ピクセルのアレイが小さな四角のアレイとして示されている。第一の基準を満足する固有値の最小値を有する選択された特徴の位置に対応した異なるピクセルが示されている。これらの特徴と位置は、図３において対応するピクセルに異なる外観を与えることによって強調表示されている。図４において、5×５隣接部が、その5×５隣接部の内側を刻む円によって強調表示されている。
【０１６３】
図１１ないし１８の例において、グリッドはボクセルアレイを表現している。グリッドの各要素は、画像を形成するボクセルである。
【０１６４】
図１１の例において、４種類のボクセルがランドマークとして区別され、対応する勾配行列の最小固有値に従って分類されている。１０２で示されたボクセルは、最大の固有値を有し、２０２で示されたボクセルは、大きな固有値を有し、３０２で示されたボクセルは、小さな固有値を有し、４０２で示されたボクセルは、最小の固有値を有する。
【０１６５】
図１２において、各ボクセルに対して、直径５ピクセル／ボクセルのユークリッド距離内の隣接部が特徴として選択された中央ピクセルの周りに定義される。前記隣接部選択は、図４に円で示されている。
【０１６６】
図１２において、このイベントの第一のケースが、５０２’と記されたより太い円で強調表示されている。ここで、特徴として選択された二つのピクセルが、ピクセルのユークリッド距離として上で定義された同一の隣接部内に存在する。これらの二つのピクセルは２０２’と３０２’で示され、上記定義に従って、特徴２０２’として選択されたピクセルは、第二のピクセル３０２’より大きな固有値を備えた勾配行列を有する。よって、この第二のピクセル３０２’は廃棄される。
【０１６７】
図１３には次のステップを示し、図１１で特徴として選択されたピクセル３０２’はもはや存在しない。図１３において、隣接部を示す円を備えたピクセル２０２’’が設けられている。図１２と同様の状況がここでも存在し、ピクセル２０２’’の隣接部に、ピクセル２０２’’より小さな固有値を有する勾配行列を持ったピクセル３０２’’が存在する。
【０１６８】
従って、この場合においても、特徴として選択されたピクセル３０２’’は追跡すべき特徴のリストから廃棄される。
【０１６９】
最後に、図１４のピクセルアレイにおいて、元の特徴として選択されたピクセルの幾つかのみが存在している。ここで、単一の隣接部内において、一つの特徴以上は見られなくなっている。従って、追跡可能な特徴は、画像／ボリューム中に散開している。
【０１７０】
有効に追跡可能なピクセルを定義するための数学的形式化と基準と同様に、図１１ないし１４に示された本例は、目標ピクセルと選択された窓のピクセルとの強度信号に存し、また前記窓の画像データの勾配行列の固有値に存する特徴パラメータのセットに限られる。
【０１７１】
本発明の次のステップによれば、ピクセルの強度はさらに調査され、そして、所定値より低い強度を有する特徴として選択されたピクセルも、追跡すべき特徴に対応するピクセルのリストから廃棄される。ほとんどゼロに近い非常に低い閾値が、周囲空気内で見つけられた特性に用いられる。このステップは、図示しないが、当業者には自明である。
【０１７２】
ひとたび2次元又は3次元画像の特徴、即ち有効に追跡可能なランドマークとみなされたピクセル又はボクセルが求められると、前記特徴の追跡が実行される。本発明の好ましい実施例によれば、これはいわゆる「Lucas Kanade特徴追跡のピラミッド形の実施」を用いて実行される。詳細な記載は、Jean Yves Bouguetが記した記事「Pyramidal Implementation of the Lucas Kanade Feature Tracker, Description of the Algorithm」Intel Corporation Microprocessor Research Labs.発行、により与えられる。
【０１７３】
上述の出版物によって二次元の例について開示された理論を、以下で3次元の場合に更に適用する。上記方法とアルゴリズムとによって個別化された追跡すべき特徴を表わす各ボクセルに対し、以下のステップを実行する必要がある。
【０１７４】
ボリュームIにおけるポイントとしてｕを、ボリュームＪにおける対応するポイントとしてｖを定義する。第一ステップにおいて、レベルＬ 0, ..., mにおけるボリュームを表現するＩＬとＪＬとを備え、元のボリュームを表わす基底部（L = 0）を備えたボリュームＩとＪに対応する二つのピラミッドが設けられる必要がある。各以下のフロアのボリュームは、一つの平均値に対する直接の隣接部における一定量のピクセルを組合わせることによって寸法が削減される。追加のフロアの数は、ボリュームのサイズに依存する。
【０１７５】
４．２．以下のステップにおいて、いわゆる包括的ピラミッド形動きベクトルｇが最上階で初期化される。
【数８】
2009512527.tif 000010
全てのレベルＬに対して、結果として生じる動きベクトルdLが次に計算される。
これは、以下のステップに従って実行される。
ａ）各ボリュームILに対して、ポイントuLの位置が次式によって計算される。
【数９】
2009512527.tif 000011
ここでｐは実際のボリューム位置である。
ｂ）ボリュームのピラミッド形の表現の各レベルＬにおいて、デカルト座標ｘ、ｙ、ｚの各方向において勾配が計算される。デカルト座標ｘに対し、ボリューム勾配が以下の式によって計算される。この等式は、他の二つのデカルト座標ｙ、ｚに対しても適用される。
【数１０】
2009512527.tif 000012
ｃ）次に、ボリューム勾配が以下の式による勾配行列を計算するために用いられる。
【数１１】
2009512527.tif 000013
ここで値？は、エリアサイズ又はボクセルに影響を与える隣接部を定義する。行列の要素は、上述の定義によるものと同様にして得られる。即ち、
【数１２】
2009512527.tif 000014
ｄ）各レベルＬに対して、反復ベクトル？が初期化される。
【数１３】
2009512527.tif 000015
ｅ）次に、ｋ＝１から最大カウントＫまで、又は最小置換
【数１４】
2009512527.tif 000016
まで、以下のボリューム差を計算する。
【数１５】
2009512527.tif 000017
以下の等式に従って、不整合なベクトルを計算する。
【数１６】
2009512527.tif 000018
次に、次式によってオプティカルフローベクトルが定義される。
【数１７】
2009512527.tif 000019
上で定義したオプティカルフローベクトルを用いて、反復動きベクトルが計算可能である。
【数１８】
2009512527.tif 000020
そしてこれをレベルＬに対する最後のオプティカルフローと等しく設定する。
【数１９】
2009512527.tif 000021
上述のステップを各レベルに対し、最後のレベルＬ＝０になるまで繰り返す。それによって、次のレベルＬ−１に対して、包括的オプティカルフローが上記計算値から算出される。
【数２０】
2009512527.tif 000022
ｆ）最終的なオプティカルフローベクトルは、次式となる。
【数２１】
2009512527.tif 000023
ｇ）次に、ボリュームＩのポイントｕに対応するボリュームＪのポイントｖの座標は、以下のように計算できる。
v = u + d
【０１７６】
図１５ないし１８は、明らかに実際のものの簡単表現である上述の追跡方法の効果の図による概略的な発想を与えようとするものである。
【０１７７】
図１５ないし１８において、異なる時間Ｔ０とＴ１に撮影された二枚の画像ＩとＪが、図１１ないし１４の例のように、ピクセルの二次元アレイとして表現されている。
【０１７８】
追跡すべき特徴に対応し、上述の方法によって選択されたピクセル２０が、画像Ｉの二次元ピクセルアレイにおいて示されている。画像Ｊにおいて同一の位置のピクセル２０’’が、画像化された対象この場合乳房の動きが無い場合におけるピクセル２０に対応している。ピクセル２０’’の一定の隣接ピクセルが、ピクセル２０’’を中心とした７×７ピクセルアレイに対応した灰色の影つきの四角で示されており、隣接部は番号１０が付されている。
【０１７９】
画像Ｉの特徴２０の新しい位置が２０’で示されている。これは、ピクセル２０の位置における乳房組織が、ピクセル２０’の位置に時間とともに移動したことを意味している。明らかに、前記ピクセル２０’は、追跡が実行される方法を説明するためのみに示されている。実際の場合において、このピクセルは未知であり、さもなければ追跡は不必要である。
【０１８０】
既に上で述べた数学的形式による選択された特徴に対応するピクセルを追跡する一般的な方法について、以下に図で与えれ、簡単のために二次元の場合に限定された実施例によって記載する。しかし、追跡の一般的数学的記述を踏まえた当業者には、それでも3次元の場合について追跡する方法を理解し、評価することができる。
【０１８１】
追跡の初めに、画像I内のピクセル２０が、同じ位置、即ち画像Jのピクセル２０’’にあるものと仮定する。このピクセル２０’’を、以下では「追跡されたピクセル」と呼び、一方で特徴として選択された画像Iのピクセル２０を「特徴ピクセル」と呼ぶこととする。まず追跡されたピクセル２０’’の隣接部に対応した一定の領域の反転勾配行列にアルゴリズムを適用し、本例においては画像Jのピクセル２０’’を中心とした前記７×７ピクセルアレイ（この値は調整可能）が計算される。この隣接部が、いわゆる追跡領域に相当する。Lucas & Kanadeの特徴追跡アルゴリズムのピラミッド形の実施を適用することで、これは3次元画像ボリュームに対しても適応されて、詳細は後述するが、画像J上の追跡された特徴の新しい位置を見つけることができる。
【０１８２】
画像J内の追跡された特徴２０’の新しい位置は、オプティカルフローベクトルを決定する。これは、第一の反復ステップにおいて計算されて定義される。図１５と１６の比較から明らかなように、追跡領域は、追跡領域の中心が同じ組織領域２０’の新しい位置により近くなるようにシフトされる。
【０１８３】
追跡は、画像I内の特徴を表わすピクセル２０と画像J内の追跡領域の中心にある新しい識別された位置との間の不整合ベクトルの計算を繰り返すことによって、さらに反復的して適用される。これは、新しいオプティカルフローベクトルを定義することになり、図１７に示すように、ピクセル隣接部又は追跡領域の新たなシフトを生じさせる。
【０１８４】
画像I内の特徴ピクセル２０と追跡領域の中心間の不整合ベクトルを求めるステップは、不整合ベクトルが調整可能な最小値に達するまで続けられる。
【０１８５】
この状況を図１８に示す。図１５ないし１７の図的表現に関連して、ピクセル２０’は追跡領域がシフトするという事実を強調するために位置を変えておらず、追跡領域は最終的に画像Iのピクセル２０としての乳房の同じ位置を表現する正しい新しい位置２０’を中心とする位置に動いている。
【０１８６】
上記例は、ピクセルによって表現されたただ一つの特徴に対して実行されている。上記ステップは、画像内で選択された各特徴に対して実行する必要がある。
【０１８７】
すべての特徴がその新しい位置にまで追跡された後に、画像／ボリューム全体にわたって、ただし乳房組織を表現する領域内でのみ、異なる位置にあるオプティカルフローベクトルの一定量が存在する。他の登録アルゴリズムと比較した最大の長所は、画像／ボリューム全体にわたってオプティカルフローベクトルが散会することである。次に、逆オプティカルフローベクトルを画像／ボリュームJ内の各特徴と、その周囲の隣接部とに適用することによって、モーフィングが実行される。これは、新しい仮想画像／ボリュームJ'をもたらす。
【０１８８】
図１９に、実際の診断画像に対する前記方法の効果を示す。左側の第一の画像セットは、造影剤投与無しでの時刻ｔ０における乳房検査からの代表的なスライスを示す。左側の第二画像は、造影剤投与後の時刻ｔ１における同じスライス位置を撮影したものである。左側の第三画像は、時刻ｔ１のピクセルの信号強度から時刻ｔ０のものをピクセル毎に引き算した結果を視覚化したものであり、いわゆる引き算画像である。時刻ｔ１とｔ０間で乳房の相対的移動が無い場合、造影剤強調が生じない中で、結果の画像では黒色を示すゼロ信号強度のピクセルが多数を占める。時刻ｔ１とｔ０間での移動は、画像／ボリュームIとJの同じ位置のピクセル単位でアーチファクトが発生するという効果を有し、これは組織の同じ領域には対応していない。本例においては、引き算画像は小さな矢印で示されたコントラスト強調された腫瘍を含む乳房を表わしている。乳房のほとんどの種類の悪性および良性腫瘍は、増加した血管新生によって生じる追加的な造影剤摂取を示す。両乳房の中央の白く残った構造は、移動と組み合わされた正常な乳房柔組織の弱い造影剤摂取の結果である。
【０１８９】
図１９の右側は、前記登録アルゴリズムが適用された対応する画像を示す。
【０１９０】
時刻ｔ０における画像はアルゴリズムによって修正されておらず、従って同じものを留めている。しかし、三次元特徴追跡後は、新たな仮想画像ボリュームが算出される。時刻ｔ１’での右側の第二画像は、登録適用後の同じスライス位置の対応する仮想画像を示す。右側の最後の画像は、時刻ｔ１’から時刻ｔ０でのピクセルの信号強度を引く、ピクセル毎の引き算結果を視覚化したものであり、新たな（仮想）引き算画像である。一つの乳房は、小さな矢印で示されたコントラスト強調された腫瘍を表わしている。正常な乳房柔組織を表わす乳房内の領域は、幾分抑制された造影剤摂取を示し、動きアーチファクトの減少により、僅かに少ない輝きを見せている。両乳房上部の白い構造の大部分は、動きアーチファクトの補償によって消失している。
【０１９１】
一枚は造影剤投与前、一枚は投与後の二枚の連続した画像のみが検討されているので、図１９の例は単純化されている。
【０１９２】
実際には、一連の断面画像（例えば５０枚以上）が、造影剤の投与前後の一定時間毎に撮影され、引き算画像が連続した画像の各二枚間で算出される。
【０１９３】
１組の中に一枚以上の断面画像を確保するのに加えて、２組以上（例えば５組以上）の画像を次々に確保することもできる。これは、複数の引き算画像／ボリュームをもたらす。もちろん、開示したアルゴリズムは各組のすべての組み合わせに適用できるので、造影剤投与の前か後かに撮影されたにかかわらず、その間に生じた動きを等しく削除することができる。
【０１９４】
本発明に係る上述の方法は、コンピュータで実行されるプログラムとしてコード化されたアルゴリズムによって実行可能である。
【０１９５】
このプログラムは、一枚以上のフレキシブルディスクや、一枚以上の光磁気ディスク、一枚以上のＣＤ−ＲＯＭ又はＤＶＤ−ＲＯＭその他携帯用記憶装置や、安定したハードディスク等に保存できる。代替的に、ペンドライブ等の固体電子装置にも保存できる。
【図面の簡単な説明】
【０１９６】
【図１】胸部ＭＲＩで使用される乳房固定装置の概要図である。
【図２】引き算プロセスと登録ステップがない場合に生じる動きアーチファクトを示す。上部左右には、同じＭＲスキャナ位置の乳房の二枚の断面画像を示す。第二画像は第一画像の約２．５分後の造影剤投与後に測定した図である。両画像は、断面画像の全体的な3次元セットの代表的な断面図である。
【図３】本発明に係る登録ステップを示す概要図である。
【図４】１０×１０画素アレイで形成されたデジタルの又はデジタル処理された画像の概要図である。
【図５】図４の画像に対応する画像データアレイの概要図であり、そこにおいて各画素の視覚的外観、例えばグレースケールデジタル画像の画素の強度は数値で表現される。
【図６】ピクセルの外観、特に選択された窓のピクセルの強度（図５に示す）を表す数値によって目標ピクセルをコード化するために用いられた最先端技術に従うベクトルの第一実施例を示す概要図である。
【図７】図６に従うピクセルコード化方法の概要図を示し、そこにおいてベクトルの成分はより一般化された方法で表現されている。
【図８】座標ｍ、ｎの目標ピクセルを中心とする一般的な３×３窓と、一般的な勾配行列と、数値的な偏導関数の計算の一つの方法を示す図であり、例えばSobel演算子やFrei-Chen演算子（Digital Image Processing, W. K. Pratt, John Wiley & Sons, Inc. New York, NY, USA, 1991）のような、より複雑な演算子を用いても導関数を計算することができる。
【図９】一般的なヘッシアン行列と数値的二次偏導関数の一つの計算方法を示す。
【図１０】本発明に係る目標ピクセルを表現するベクトルの概要図である。
【図１１】有効に追跡可能なランドマークの自動測定法の一例の概要図である。
【図１２】有効に追跡可能なランドマークの自動測定法の一例の概要図である。
【図１３】有効に追跡可能なランドマークの自動測定法の一例の概要図である。
【図１４】有効に追跡可能なランドマークの自動測定法の一例の概要図である。
【図１５】追跡がどのように実行されるかを示す概要図である。
【図１６】追跡がどのように実行されるかを示す概要図である。
【図１７】追跡がどのように実行されるかを示す概要図である。
【図１８】追跡がどのように実行されるかを示す概要図である。
【図１９】造影剤投与前後の異なる時間に撮影された二枚の画像と、特徴追跡と画像登録を実行した後の、両方の画像の引き算で得られた画像とを、図の左側に示した図2の例と比較して示すことによる本発明の結果を表した図である。図２から図１９までのすべての画像は、登録以外の他の処理は一切しないで同一のMR試験で撮影されたものである。
【図２０】対応する特徴パラメータによって、ピクセルの有効に又は非有効に追跡可能な特徴の質を測定するためのステップからなる本発明にかかる方法の第一実施例を示す図である。
【図２１】対応する特徴パラメータによって、ピクセルの有効に又は非有効に追跡可能な特徴の質を測定するためのステップからなる本発明にかかる方法の第二の実施例を示す図である。
【図１】
2009512527.tif 000044
【図３】
2009512527.tif 000045
【図４】
2009512527.tif 000046
【図５−７】
2009512527.tif 000047
【図８】
2009512527.tif 000048
【図９−１０】
2009512527.tif 000049
【図１１】
2009512527.tif 000050
【図１２】
2009512527.tif 000051
【図１３】
2009512527.tif 000052
【図１４】
2009512527.tif 000053
【図１５】
2009512527.tif 000054
【図１６】
2009512527.tif 000055
【図１７】
2009512527.tif 000056
【図１８】
2009512527.tif 000057
【図２０】
2009512527.tif 000058
【図２１】
2009512527.tif 000059
【図２】
2009512527.tif 000060
【図１９】
2009512527.tif 000061
【国際調査報告】
2009512527.tif 000062
2009512527.tif 000063
2009512527.tif 000064