JP2004537357A

JP2004537357A - 関節疾患および関節疾患の時間変化を定量的に評価するためのシステムおよび方法

Info

Publication number: JP2004537357A
Application number: JP2003517824A
Authority: JP
Inventors: マージャッタソフィアトッターマンサーラ; タメツ−ペナジョセ; アシュトンエドワード; パーカーケビン
Original assignee: ヴァーチャルスコピックスリミテッドライアビリティカンパニー
Priority date: 2001-07-27
Filing date: 2002-07-26
Publication date: 2004-12-16
Also published as: WO2003012724A1; WO2003012724A9; CA2452547A1; EP1412909A4; EP1412909A1; US20030035773A1

Abstract

ヒトまたは動物の関節では、特定の対象が関節疾患の指標またはバイオマーカーの働きをする。関節の３次元画像（１０２）でバイオマーカーを識別し定量化する（１０４）。複数の３次元画像を経時的に撮影することができ（１０６）、その画像中でバイオマーカーを経時的に追跡することができる（１１２）。統計的セグメント化技法を使用して、最初の画像中のバイオマーカーを識別し、その識別を残りの画像に引き渡す。

Description

【技術分野】
【０００１】
本発明は、関節疾患および関節疾患の時間変化を定量的に評価するためのシステムおよび方法を対象とし、バイオマーカー（biomarker）を使用した該システムおよび方法を対象とする。
【背景技術】
【０００２】
かなりの割合の人々が、変形性関節症、他の変性および外傷後疾患などの関節疾患に苦しんでいる。加えて、不幸にして機能減退の進行にまで至る、膝、肩、肘、手首、足首および他の複合関節、ならびにそれらを支持する靱帯および構造の数多くの損傷がある。このような状況の評価、および新しい治療法による改善を含むこのような状況の時間変化の追跡では、定量的情報を得ることが必要である。過去には、痛みまたは不快感といった主観的な尺度が使用された。Ｘ線フィルム上の画像およびディジタルＸ線画像の測定によってより客観的な尺度を得ることができるが、それらは伝統的に、画像の手動トレースまたはカリパス測定によって評価されている。ＭＲＩおよびＣＴスキャナの３Ｄ画像が使用可能であれば、通常は画像ワークステーションにインタフェースされたマウスまたはトラック・ボールを使用した関心の物体のトレースによって、より詳細な手動評価を得ることができる。ヒトおよび動物の膝の変形性関節症で実施される測定の例には、軟骨の厚さ、軟骨の体積、軟骨および骨の画像強度、ならびに軟骨のＴ２緩和時間が含まれる。
【０００３】
これまでの研究におけるこれらの測定は、ユーザが関心の構造を識別し、境界または領域をトレースするか、またはアクティブ・コンタ（active contour）を初期化する必要がある手動または半手動システムを必要とする（例えば、非特許文献１、非特許文献２、非特許文献３、非特許文献４、非特許文献５参照）。
【０００４】
従来技術は、全般的な異常または全般的な経時変化を評価することはできる。しかし、従来の測定は、微妙な異常または微妙な変化を評価、定量化するのにあまり適してはおらず、複雑なトポロジー（topology）または形状を正確に記述することができない。さらに、未処理画像の手動および半手動測定は、高い空間間変動および同一観察者内変動を受ける。さらに、トレーシングが一連の２Ｄ画像に基づくときに手動および半手動測定は、３Ｄ画像中にでこぼこした不規則な境界を生み出す傾向がある。
【０００５】
【特許文献１】
米国特許第６１６９８１７号明細書（Ｐａｒｋｅｒ他）
【非特許文献１】
Eckenstein F., Gavazzeni H.S., Sittek H., Haubner, M., Losch, A., Milz, S., Englmeier, K-H., Schulte, E., Putz, R, Reiser, M., "Determination of Knee Joint Cartilage Thickness using Three-Dimensional Magnetic Resonance Chondro-Crassometry (3D MR-CCM)," Magnetic Resonance in Medicine 36: 256-265, 1996
【非特許文献２】
Solloway, S., Hutchinson, C.E., Waterton, J.C., Taylor, C., "The Use of Active Shape Models for Making Thickness Measurements of Articular Cartilage from MR Images," Magnetic Resonance in Medicine 37: 943-952, 1997
【非特許文献３】
Stammberger, T., Eckstein, F., Englmeier, K-H., Reiser, M. "Determination of 3D Cartilage Thickness Data from MR Imaging: Computational Method and Reproducibility in the Living," Magnetic Resonance in Medicine 41: 529-536, 1999
【非特許文献４】
Ghosh, S., Ries, M., Lane, N., Majundar, S. "Segmentation of High Resolution Articular Cartilage MR Images," 46th Annual Meeting, Orthopaedic Research Society, March 12-15,2000, Orlando Florida
【非特許文献５】
Dardzinski, B.J., Mosher, T.J., Li, S., Van Slyke, M.A., Smith, M.B.,"Spatial Variation of T2 in Human Articular Cartilage, Radiology 205: 546-550, 1997
【非特許文献６】
Peet, F.G., Sahota, T.S. "Surface Curvature as a Measure of Image Texture" IEEE Transactions on Pattern Analysis and Machine Intelligence 1985 Vol PAMI-7 G: 734-738
【非特許文献７】
Struik, D.J., Lectures on Classical Differential Geometry, 2nd ed., Dover, 1988
【非特許文献８】
Duda, R.O., Hart, P.E., Pattern Classification and Scene Analysis, Wiley & Sons, 1973
【非特許文献９】
Jain, A.K., Fundamentals of Digital Image Processing, Prentice Hall, 1989
【非特許文献１０】
Matheny, A., Goldgof, D. "The Use of Three and Four Dimensional Surface Harmonics for Nonrigid Shape Recovery and Representation," IEEE Transactions on Pattern Analysis and Machine Intelligence 1995, 17: 967-981
【非特許文献１１】
Chen, C.W, Huang, T.S., Arrot, M. "Modeling, Analysis, and Visualization of Left Ventricle Shape and Motion by Hierarchical Decomposition," IEEE Transactions on Pattern Analysis and Machine Intelligence 1994, 342-356
【非特許文献１２】
Pena, J.-T., Totterman, S.M.S., Parker, K.J. "MRI Isotropic Resolution Reconstruction from Two Orthogonal Scans," SPIE Medical imaging, 2001
【発明の開示】
【発明が解決しようとする課題】
【０００６】
当技術分野では、未処理画像の手動および半手動測定ならびに２Ｄ画像の使用に関連した上記の問題の解決が求められていることは明らかである。
【０００７】
したがって本発明の主な目的は、関節および関節疾患のより正確な定量化を提供することにある。本発明の他の目的は、関節疾患の時間変化のより正確な定量化を提供することにある。本発明の他の目的は、上で述べたニーズに対処することにある。
【０００８】
上記の目的および他の目的を達成するために、本発明は、重要な構造または下位構造、その常態および異常の識別、ならびに関節疾患および病態の敏感な指標であるその具体的なトポロジー的、形態的、放射線学的および薬物動態学的特性の識別を対象とする。構造の異常および常態、そのトポロジー的および形態的特性、ならびに放射線学的および薬物動態学的パラメータはバイオマーカーと呼ばれ、バイオマーカーの特定の測定は、関節疾患の定量的評価に役立つ。
【０００９】
発明者らは、以下の新しいバイオマーカーは、ヒトおよび動物の変形性関節症の敏感な指標であることを発見した。
【００１０】
・軟骨下骨板（subchondral bone plate）の形状
・軟骨層（layers of cartilage）およびその相対寸法
・軟骨層内での信号強度分布
・関節軟骨（articulating cartilage）表面間の接触面積
・軟骨形状の表面トポロジー
・骨髄水腫の強度
・骨間分離距離
・半月板（meniscus）の形状
・半月板の表面積
・半月板と軟骨の接触面積
・軟骨構造特性
・軟骨表面特性
・半月板構造特性
・半月板表面特性
・パンヌス（pannus）構造特性
・関節液（joint fluid）特性
・骨増殖体（osteophyte）特性
・骨特性
・溶解性病変（lytic lesion）特性
・プロテーゼ接触（prosthesis contact）特性
・プロテーゼの摩耗
・関節間隔（joint spacing）特性
・脛骨内側軟骨（tibia medial cartilage）の体積
・脛骨外側軟骨（tibia lateral cartilage）の体積
・大腿骨軟骨（femur cartilage）の体積
・膝蓋骨軟骨（patella cartilage）の体積
・脛骨内側軟骨（tibia medial cartilage）の湾曲
・脛骨外側軟骨（tibia lateral cartilage）の湾曲
・大腿骨軟骨（femur cartilage）の湾曲
・膝蓋骨軟骨（patella cartilage）の湾曲
・軟骨屈曲（cartilage bending）エネルギー
・軟骨下骨板（subchondral bone plate）の湾曲
・軟骨下骨板（subchondral bone plate）屈曲エネルギー
・半月板の体積
・骨増殖体の体積
・軟骨Ｔ２病変（cartilage t2 lesion）の体積
・骨髄水腫（bone marrow edema）の体積および数
・滑液（synovial fluid）の体積
・滑膜（synovial）の厚化
・軟骨下骨嚢胞（subchondral bone cyst）の体積
・脛骨の運動学的平行移動（kinematic tibial translation）
・脛骨の運動学的回転
・脛骨の運動学的外反
・椎体間の距離
・骨組みの減退（subsidence of cage）の程度
・角度測定による脊柱前湾症の程度
・椎体間のオフセットの程度
・大腿骨（femoral bone）特性
・膝蓋骨（patella）特性
バイオマーカーを抽出する好ましい技法は、統計学に基づく推論を用いる技法である（特許文献１参照）。特許文献１における開示は、その全体が参照によって本明細書の開示に組み込まれる。形状およびトポロジーを定量化する好ましい方法は、先行文献（湾曲分析については、例えば、非特許文献６、７を参照；形状およびトポロジーの記述子については、例えば、非特許文献８、９を参照；球面調和関数については、例えば、非特許文献１０、１１を参照。）によって明確に記述されている形態的およびトポロジー的な式を用いる方法である。
【００１１】
これらの形態およびトポロジー測定は、過去において、関節バイオマーカーに対しては適用されていなかった。
【００１２】
それぞれの関節バイオマーカーの定量的測度を得ることができる。これは、湾曲、トポロジーおよび形状のうちの１つまたは複数であることができる。
【発明を実施するための最良の形態】
【００１３】
次に、本発明の好ましい一実施形態を、図面を参照して詳細に説明する。
【００１４】
図１に、バイオマーカーおよびバイオマーカーの傾向を経時的に識別するプロセスの概要を示す。ステップ１０２では器官の３次元画像を撮影する。ステップ１０４では、その画像中の少なくとも１つのバイオマーカーを識別する。そのための技法については図２を参照して後に説明する。ステップ１０４ではさらに、そのバイオマーカーの少なくとも１つの定量的測定値を得る。ステップ１０６では、その器官の同じ領域の複数の３次元画像を経時的に撮影する。場合によってはステップ１０６をステップ１０４の前に実施することができる。これらのステップの順序はそのときの都合によって決めることができる。ステップ１０８では、経時的に撮影したこれらの複数の画像中の同じ１つまたは複数のバイオマーカーおよびそれらの定量的測定値を識別する。そのための技法については図３を参照して後に説明する。複数の画像でバイオマーカーを識別することは、ステップ１１０における、４次元、すなわち３次元空間および１つの時間次元に、その器官のモデルを作成することを可能にする。ステップ１１２では、このモデルから、その１つまたは複数のバイオマーカーの展開像（development）を経時的に追跡することができる。
【００１５】
バイオマーカーを抽出する好ましい方法は、特許文献１に明確に記載されている統計学に基づく推論を用いる方法である。この特許文献１の開示内容はその全体が参照によって本明細書の開示に組み込まれている。磁気共鳴画像法などの方法によって得た未処理（生）の画像データから、被写体が４次元（空間および時間）に再構築されて、視覚化される。この再構築および視覚化は、まず最初に一連の画像中の最初の画像を、その画像データの平均値および分散の統計的推定、ならびに十分に類似した画素（ボクセル）の接合によって複数の領域に分割し、次いで、画像の構成要素の既知の動き特性（例えば筋肉および腱のばね定数）を使用してそれぞれの領域の画像間の剛体動きおよび変形動き（rigid and deformational motion）を推定することによって、これらの領域を残りの画像に外挿することによって実施される。この被写体および被写体の複数の領域を、３次元空間と１つの時間次元から成る４次元（４Ｄ）仮想現実環境でレンダリングし、相互作用させることができる。４つの次元は
このセグメント化（segmentation）を、図２を参照して説明する。まず最初に、ステップ２０１で、ＭＲＩなどによって一連の画像を撮影する。このようにして未処理（生）の画像データを得る。次いで、ステップ２０３で、この一連の画像中の最初の画像の未処理（生）データをコンピューティング装置に入力する。次に、ステップ２０５で、それぞれのボクセルに対して、画像データの局所平均値（local mean value）および領域分散（region variance）を推定する。ステップ２０７で、ステップ２０５で推定した平均値および分散の比較によって、ボクセル間の連結性（connectivity）を推定して、複数の領域を形成する。連結性を推定した後に、分離する必要がある領域を決定し、ステップ２０９でそれらの領域を分離する。ステップ２１１のセグメンテーションの緩和（segmentation relaxation）によって、それらの領域の正確さをさらに高めることができる。次いで、結合する必要がある領域を決定し、ステップ２１３でそれらの領域を結合する。ステップ２１５で再びセグメンテーションの緩和を実行する。このようにして、未処理（生）の画像データを、ステップ２１７の最終結果であるセグメント化された画像に変換する。これらのプロセスの詳細は前掲の特許文献１に出ている。
【００１６】
４Ｄモデル（３つの空間次元および１つの時間次元）の生成を、図３を参照して説明する。ステップ３０１で、一連の画像中の最初の画像、およびこの画像から先に説明したようにして得た完全にセグメント化された画像を入力すると、動き追跡／推定アルゴリズムが、セグメント化された画像を１つのフレームから別のフレームに渡すのに必要な情報を提供する。３Ｄ動きの推定においては、理想的には、剛体構成要素と非剛体構成要素の両方の存在を考慮に入れなければならない。本発明によれば、その画像中の選択された特徴ポイントの位置決めの後に、それぞれのボクセルの動きベクトルを推定する。
【００１７】
関節に存在する多くの構造の動きを考慮するため、本発明のアプローチでは、セグメント化された画像中に見られるさまざまな構造の物質特性のアプリオリな（a priori）知識を使用することによって、軟組織の局所変形を考慮に入れる。このような知識は、ステップ３０３で、適当なデータベースの形態で入力される。さらに、剛体構造の動きおよび軟組織の動きに対して異なる戦略を適用することができる。選択した複数のポイントを登録した後、選択したポイントの動きベクトルを内挿することによって、画像中のすべてのボクセルの動きベクトルを計算する。それぞれのボクセルの動きベクトルが推定されてしまうと、一連の画像の次の画像のセグメント化は、前の画像のセグメント化の伝搬（propagation）に過ぎない。この技法を、一連の画像のすべての画像が解析されるまで繰り返す。なお、解析の都合によっては時間の定義および配列の順序を逆にすることができる。
【００１８】
画像解析用の、およびタイム・エボリューション解析（time-evolution analysis）用の有限要素モデル（ＦＥＭ）が知られている。本発明は、同様のアプローチに従い、その解剖学的構造の物理特性をモデル化した質量とばねとからなるメッシュ（mesh）の総エネルギーを最小化することによって、ポイント対応（point correspondence）を回復する。本発明では、このメッシュは、画像中の単一の構造に拘束されず、代わりに、この体積測定（volumetric）の画像全体を自由にモデル化することができる。このとき、トポロジー的特性はセグメント化された最初の画像によって供給され、物理特性は、アプリオリな特性およびセグメント化された最初の画像によって供給される。この動き推定アプローチは、一連の画像中の２つの連続した画像間のＦＥＭベースのポイント対応回復アルゴリズムである。メッシュのそれぞれのノードは、追跡しようとする画像の自動的に選択された特徴ポイントであり、ばねの剛性は、セグメント化された最初の画像およびヒトの解剖学的構造のアプリオリな知識、ならびに筋肉、骨などの一般的な生体力学的特性から計算される。
【００１９】
変形可能な多くのモデルは、ばねによって連結されたポイント質量を駆動するベクトル力フィールド（vector force field）を画像から抽出できると仮定する。このようなモデルの大部分は、半自動特徴抽出アルゴリズムを構築するこのアプローチを使用する。本発明も同様のアプローチを使用し、ｔ＝ｎにサンプリングされた画像は、３つの動的スカラー・フィールドの集合、Φ（ｘ，ｔ）＝｛ｇ_ｎ（ｘ），｜∇ｇ_ｎ（ｘ）｜，∇^２ｇ_ｎ（ｘ）｝、すなわちグレー・スケール画像値、この画像値の勾配（gradient）の絶対値、およびこの画像値のラプラス演算子（Laplacian）の集合であると仮定する：
したがって、Φ（ｘ，ｔ）の変化は、スカラー・フィールド・エネルギーＵ_Φ（ｘ）∝（ΔΦ（ｘ））^２の２次の変化を引き起こす。さらに、この画像の根底にある構造は、これらのスカラー・フィールドと平衡状態にある、ばねによって連結されたポイント質量のメッシュとしてモデル化されると仮定する。平衡状態は外部力フィールド（external force field）の存在を仮定するが、力フィールドの形状は重要ではない。ポイント質量の分布は時間変化するものとする。時刻ｔ＝ｎの後の時間間隔Δｔにおける総エネルギーの変化は下式によって与えられる。
【００２０】
【数１】

【００２１】
上式で、α、βおよびγは、個々のフィールドの変化の寄与に対する重み付けであり、ηは、ひずみエネルギーの増大に重み付けを与え、Ｋは、ＦＥＭ剛性行列であり、ΔＸは、ＦＥＭノード変位行列（node displacement matrix）である。この式の解析によれば、画像フィールドまたはメッシュ・ポイント分布が変化するとその系の総エネルギーが増大する。したがって、ｇ_ｎからｇ_ｎ＋１へのポイント対応は、その総エネルギー変化が最小であるメッシュ構成によって与えられる。したがって、ポイント対応は下式によって与えられる。
【００２２】
【数２】

【００２３】
上式で
【００２４】
【数３】

【００２５】
であり、この表記において、ｍｉｎ_ｐｑはｑを最小にするｐの値である。
【００２６】
以上に示した諸式を使用して、画像中のすべてのボクセルの動き（ポイント対応）を推定することができると思われるが、ボクセルの数が一般に百万を超え、また、これらの式が複雑であることから、全面的な最小化（global minimization）は困難である。この問題を単純化するため、ステップ３０５で、画像から選択したポイントを用いて粗いＦＥＭメッシュを構築する。このエネルギー最小化は選択されたポイントのポイント対応を与える。
【００２７】
このようなポイントの選択は単純なものではない。第１に、実用上、ポイントの数は、一般に≡１０^４と非常に少なくなければならない。このとき、選択したポイントが画像動き全体を記述するように注意を払わなければならない。第２に、領域境界の追跡は領域の動きの正確な記述には十分なので、領域境界は、重要な特徴である。第３に、領域境界では勾配の絶対値が大きく、そのラプラス演算子は、ゼロ交差ポイントにあり、したがって領域境界を追跡しやすい特徴にする。したがって、ＦＥＭの構築では、セグメント化された境界ポイントが選択される。
【００２８】
境界ポイントは画像ポイントについての小規模の部分集合だが、それでも実用上、境界ポイントの数は多すぎる。このポイントの数を減らすため、境界ポイントの条件付き無作為抽出をポイント抽出ステップに使用する。その制約条件は、すでに選択したポイントに近すぎるポイントは選択しないということから成る。この制約条件によって、境界に沿ってポイントをより均一に選択することができる。最後に、それぞれの領域の内部のポイントでの動き推定の誤差を低減するため、同じ距離制約条件を使用して画像のさらにいくつかのポイントを無作為に選択する。実験結果によれば、２５６×２５６×３４ボクセルの一般的な体積測定画像の動きを推定し記述するのには、５０００個から１００００個のポイントがあれば十分である。選択したポイントのうち、７５％が境界ポイントとして任意に選択され、残りの２５％は内部のポイントである。もちろん、必要に応じて他の割合を使用することもできる。
【００２９】
追跡するポイントの集合を選択したら、次のステップは、ステップ３０７でそれらのポイントに対してＦＥＭメッシュを構築することある。このメッシュは、それぞれの領域の物質特性および相互作用特性をコード化することによって許される動きの種類を制限する。最初のステップは、すべてのノード・ポイント（nodal point）、その近傍のノード・ポイントを見つけることである。近傍のノード・ポイントを見つけるこの操作は、メッシュのボロノイ図（Voronoi diagram）を構築することに対応することを当業者は理解しよう。その２重の、デローニー三角形分割（Delaunay triangulation）は、与えられたノード構成に対する可能な最良の四面体有限要素を表す。ボロノイ図は、膨張（dilation）アプローチによって構築される。このアプローチの下で、離散的な体積中のそれぞれのノード・ポイントを膨張させる。このような膨張は２つの目的を達成する。第１に、膨張した１つのポイントが他のポイントと接触したときにそれはテストされ、そのためその近傍のポイントを識別することができる。第２に、すべてのボクセルをメッシュのポイントに関連づけることができる。
【００３０】
すべてのポイントｘｉを近傍のポイントｘｊと関連づけたら、この２つのポイントは、ばね定数
【００３１】
【数４】

【００３２】
を有するばねによって連結されていると仮定する。ただしｌおよびｍは物質を識別する。このばね定数は、連結されたこれらのポイントの物質の相互作用特性によって定義される。物質の相互作用特性は、その物質の既知の特性に従ってユーザがあらかじめ定義する。連結されたポイントが同じ領域に属する場合、ばね定数は
【００３３】
【数５】

【００３４】
となり、その領域の物質の弾性から導かれる。連結されたポイントが互いに別の領域に属する場合、ばね定数は、境界での物質と物質の間の平均相互作用力から導き出される。画像化している物体がヒトの肩である場合、ばね定数は以下のような表から導き出すことができる。
【００３５】
【表１】

【００３６】
理論的には、相互作用は、隣接する任意の２つの領域間で定義されなければならない。しかし実際には、画像中の主要な解剖学的構成要素間の相互作用だけを定義し、残りを任意定数として残すことは許容される近似である。このような近似において、もたらされる誤差は、先に記載したいくつかの仮定でもたらされる他の誤差に比べればそれほど重大なものではない。
【００３７】
骨のおおよそのサイズおよび画像強度（輝度）は通常アプリオリに分かっているため、ばね定数は自動的に割り当てることができる。アプリオリな予想と一致するセグメント化された画像領域では、骨に対して比較的に堅い弾性定数が割り当てられる。軟組織には比較的に軟かい弾性定数が割り当てられる。
【００３８】
メッシュを構成した後、ステップ３０９で、一連の画像中の次の画像を入力し、ステップ３１１で、一連の画像中のこれらの２つの連続する画像間のエネルギーを最小化する。エネルギーＵを最小化する問題は２つの別個の問題、すなわち、剛体構造の動きに関連したエネルギーを最小化する問題と、変形可能な組織の動きに関連したエネルギーを最小化する問題とに分離することができる。これらのエネルギーはともに同じエネルギー関数を使用するが、異なる戦略に依存する。
【００３９】
剛体（構造の）動きの推定は、メッシュ変形エネルギー（ΔＸ^ＴＫΔＸ）／２に対する剛体（構造の）動きの寄与がゼロに極めて近いという事実に依存する。セグメント化およびその解剖学的構造のアプリオリな知識は、いずれのポイントが剛体に属するかを指示する。このようなポイントが個々のすべての剛体領域に対して選択される場合、剛体（構造の）動きエネルギーの最小化は、それぞれの剛体領域Ｒ_ｉに対して、その領域自体のエネルギー
【００４０】
【数６】

【００４１】
を最小化する剛体（構造の）動きの回転Ｒ_ｉおよび平行移動Ｔ_ｉを見つけることによって達成される。上式で、
【００４２】
【数７】

【００４３】
であり、
【００４４】
【数８】

【００４５】
は、剛体領域Ｒ_ｉに属する各ポイントの最適変位行列である。この剛体（構造の）動きの最小化問題のそれぞれの剛体領域に対する自由度は、回転行列で３、平行移動行列で３、の６でしかない。したがって、一連の画像中の任意の２つの画像間の差が十分に小さい場合には、６次元最速降下法によって１２の成分（回転９、平行移動３）を見つけることができる。
【００４６】
剛体（構造の）動きパラメータを見つけたら、系の総エネルギーＵの最小化によって変形（可能な組織の）動きを推定する。この最小化は、剛体エネルギーの最小化ほどには単純化することができず、追加の考慮事項なしで、変形可能な３Ｄ物体の自由度はメッシュ全体のノード・ポイントの数の３倍である。この問題のこの性質によって、メッシュのそれぞれのノードに対して単純な勾配降下技法（gradient descent technique）を使用することができる。位置および動きエネルギーから、系のラグランジアン（動きエネルギーから位置エネルギーを差し引いたものとして物理学で定義される動きポテンシャル）を使用して、局所駆動力がエネルギー・フィールドの勾配である系のすべてのノードに対してオイラー−ラグランジュの方程式を導き出すことができる。メッシュのすべてのノードに対して局所エネルギーは下式によって与えられる。
【００４７】
【数９】

【００４８】
上式で、Ｇ_ｍはボロノイ図の近傍を表す。
【００４９】
したがって、すべてのノードに対して自由度３の問題があり、その最小化は、局所ノードのエネルギーを繰り返し低減する単純な勾配降下技法を使用して実行される。局所ノードの勾配降下式は、
【００５０】
【数１０】

【００５１】
である。
【００５２】
上式で、メッシュ・エネルギーの勾配は解析によって計算可能であり、フィールドのエネルギーの勾配は、２つの異なる分解能で、画像から数値的に推定され、ｘ（ｎ＋ｌ）は次ノードの位置であり、ｖは勾配の寄与に対する重み係数である。
【００５３】
この最小化のすべてのステップで、それぞれのノードに対するプロセスは、近傍のノードの前の変位を考慮する。このプロセスは、総エネルギーが局所最小に達するまで繰り返される。小さな変形では局所最小は全体最小に近いか、またはそれに等しい。このように見つけられた変位ベクトルは、そのノード・ポイントの推定される動きを表す。
【００５４】
以上に説明した最小化プロセスから、サンプリングされた変位フィールドΔＸが得られたら、その変位フィールドを使用して、一連の１つの画像から次の画像へセグメント化を追跡するのに必要な密集した（dense）動きフィールドを推定する（ステップ３１３）。この密集した動きは、メッシュのすべての近傍のノードの寄与に重みを付けることによって推定される。一定速度モデルを仮定する。ボクセルｘの時刻ｔの推定される速度はｖ（ｘ，ｔ）＝Δｘ（ｔ）／Δｔである。密集した動きフィールドは下式によって推定される。
【００５５】
【数１１】

【００５６】
上式で、
【００５７】
【数１２】

【００５８】
であり、ｋ^ｌ，ｍは、ボクセルｘおよびｘ_ｊに関連した物質ｌとｍの間のばね定数または剛性、Δｔは、一連の画像中の連続する画像間の時間間隔、｜ｘ−ｘ_ｊ｜はボクセル間の単純なユークリッド距離であり、内挿は、ボクセルｘに最も近いノードの近傍のノードを使用して実行される。この内挿は、その物質特性
【００５９】
【数１３】

【００６０】
によってすべての近傍のノードの寄与に重み付けする。したがって推定されるボクセルの動きは、均質なあらゆる領域に対して同様であり、これは領域境界であっても同じである。
【００６１】
次いで、ステップ３１５で、一連の画像の次の画像をセグメント化データで満たす。これは、１つの画像で決定された領域を次の画像に繰り越すことを意味する。そうするためには、この次の画像のすべてのボクセルに対して速度を推定する。これは、下式によって与えられる、推定される動きの逆写像（reverse mapping）によって達成される。
【００６２】
【数１４】

【００６３】
上式で、Ｈは、次の画像の同じボクセル空間Ｓ（ｘ）に含まれるポイントの数である。この写像は、時刻ｔ＋Δｔに、その空間すべてを満たすわけではないが、写像された近傍ボクセル間の単純な内挿を使用してその空間を埋めることができる。次の画像のすべてのボクセルに対して速度が推定されれば、その画像のセグメント化は、単純に、Ｌ（ｘ，ｔ）およびＬ（ｘ，ｔ＋Δｔ）は、時刻ｔおよびｔ＋Δｔのボクセルｘのセグメント化ラベルであるとすると、Ｌ（ｘ，ｔ＋Δｔ）＝Ｌ（ｘ−ｖ（ｘ，ｔ＋Δｔ）Δｔ，ｔ）である。
【００６４】
このように展開したセグメント化を、ステップ３１７で、ステップ２１１および２１５で実施したようなリラクゼーション・ラベリングによって調整し、画像中のメッシュ・ノードに微調整を実施する。次いで、ステップ３１９で、一連の画像の最後の画像がセグメント化されたと判定されない限り、ステップ３０９で次の画像が入力される。ステップ３１９で一連の画像の最後の画像がセグメント化されたと判定された場合には、ステップ３２１でこの操作は終了となる。
【００６５】
以上に説明した操作は、図４のブロック図に示したようなシステムに実装することができる。システム４００は、画像データ、物質特性のデータベースなどを入力するための入力装置４０２を含む。入力装置４０２を通して入力された情報はワークステーション４０４で受け取られる。ワークステーション４０４は、ハードディスク・ドライブなどの記憶装置４０６と、先に開示した処理を実行して４Ｄデータを提供する処理ユニット４０８と、この４Ｄデータを例えば表面レンダリングによって表示用に準備するグラフィックス・レンダリング・エンジン４１０とを有する。出力装置４１２は、レンダリング・エンジン４１０によってレンダリングされた画像を表示するためのモニタ、またはビデオ・レコーダなどの画像を記録するための別の記憶装置、あるいはその両方を含むことができる。ワークステーション３０４およびグラフィックス・レンダリング・エンジン４１０の例は、ＳｉｌｉｃｏｎＧｒａｐｈｉｃｓＩｎｄｉｇｏワークステーションおよびＩｒｉｘＥｘｐｌｏｒｅｒ３Ｄグラフィックス・エンジンである。
【００６６】
識別されたバイオマーカーの形状およびトポロジーは、解析幾何学において周知の任意の適当な技法によって定量化することができる。形状およびトポロジーを定量化する好ましい方法は、上述の文献によって明確に記述されている形態的およびトポロジー的な式を用いる方法である。
【００６７】
新しいバイオマーカーの定量的測定の一例として、成人の膝を、１．５テスラＭＲＩシステム（1.5 Tesla MRI system）を用いて、面内（in-plane）分解能０．３ｍｍ、スライス厚２．０ｍｍでスキャンした。大腿骨、脛骨およびフィビア（fibia）の軟骨を、特許文献１（前掲）の統計的推論技法を使用してセグメント化した。これらの軟骨構造の特性は、形態的およびトポロジー的測定を適用することによって、得られた。このような測定の１つが局所表面湾曲の推定である。局所表面湾曲を測定する技法は解析幾何学においてよく知られている。例えば、Ｓ（ｘ，ｙ，ｚ）を、外向き法線および平均曲率を有する構造の表面とした場合、局所量は、非特許文献７（前掲）の８３ページに出ている２次方程式の根から決定することができる。この測定は、従来の２Ｄ画像断面の手動トレーシングを使用しては実施できない方法で、軟骨の定量的かつ再現可能な非常に微妙な特性評価を与える。
【００６８】
図５に、表面湾曲の定量的高次測度のグレー・スケール・グラフを、軟骨の表面の位置の関数として示す。この図は、上部大腿骨から膝に向かって軟骨の内面のほうを見下ろしたものである。暗から明へのわずかな相違は、高次の測定である局所湾曲の定量的測定を指示している。
【００６９】
これらのデータは、その後に間隔をもって個々の画像がスキャンされるたびに、繰り返し経時的に解析される。経時的傾向の表現方法には好ましい２つのタイプがある。１つのタイプでは、繰り返される高次測定値が、動画を形成するような早いシーケンスでされる連続する測定を上書きする態様で、図５のように表示される。その相補的な表現では、高次測度が時間の関数として与えられる傾向プロットが描かれる。例えば、特定の軟骨局所領域について、局所湾曲の平均値および標準偏差（または範囲）を、時間の関数としてプロットすることができる。
【００７０】
これらの測定の正確さ、および小さな下位構造の微妙な変化に対する測定の感度は、画像化システム（imaging system）の分解能に大幅に依存する。残念ながら、大部分のＣＴ、ＭＲＩおよび超音波システムの面外軸ないし「ｚ」軸の分解能はよいとは言えない。これらのシステムの面内分解能は一般に、１ミリメートル弱、離れた物体を解像できる程度だが、面外（スライス厚）は一般に１．５ｍｍ以上にセットされる。微妙な変化および小さな欠陥を高次構造測定を使用して評価するためには、３本の直交軸すべてで分解能が１ミリメートルよりもよいことが好ましい。これは、高分解能ボクセル・データ・セットを生み出す直交方向すなわち面外方向の高分解能スキャンの融合によって達成することができる（非特許文献１２。この文献はその全体が参照によって本明細書の開示に組み込まれる）。構造内の微妙な欠陥の評価に加えて、この高分解能ボクセル・データ・セットは、薄い構造、湾曲した構造、または蛇行した構造のより正確な測定を可能にする。
【００７１】
療法に対する人または動物の反応の追跡、あるいは病気の進行の監視では、バイオマーカーの経時的な傾向を精確に監視することが望ましい。繰り返されるスキャンを独立に精査し、関心を持たせるバイオマーカーを手動または半手動でトレースし、または測定しなければならないため、従来のやり方でこれを実行することは難しい。そのそれぞれの時間間隔は、測定を繰り返す新たな手間のかかるプロセスを表すことになる。統計的推論を用いてバイオマーカーを識別し、そのバイオマーカーをスキャンから別のスキャンへ経時的に追跡する、前掲特許文献１に記載されているものなどの４Ｄアプローチをとることは非常に有利である。すなわち、それは、関心を持たせるバイオマーカーの最初のセグメント化の結果を、後の時間間隔で撮影されたスキャンのデータ・セットに引き継ぐことである。１つのスキャンから次のスキャンへバイオマーカーの境界を追跡するために、探索が実施される。このアプローチの正確さ、精度および再現性は、１つのスキャンからそれ以降のスキャンへの境界情報の自動追跡のない、または引継ぎのない、手動または半手動画像測定を実行する際の正確さ、精度および再現性よりも優れている。
【００７２】
先に一覧した新しいバイオマーカーの定量的評価は、関節の状態、特に関節疾患の進行の客観的な測定を提供する。これらのバイオマーカーの正確な測定を経時的に得ることも非常に有用であり、これは特に、新しい療法に対する反応の程度を判断し、または加齢に伴う傾向を評価するのに役立つ。（軟骨の単純な厚さまたは体積などの）従来のバイオマーカーの手動および半手動トレーシングは、高い固有の変動性（バラツキ）を有し、そのため、連続するスキャンをトレースするときには、この変動性のために微妙な傾向が隠されることがある。このことは、従来の方法では、全般的な変化しか確認できず、時には非常に長い期間が経たないと分からない全般的な変化しか確認できないことを意味する。本発明の発明者は、統計的検定を使用してバイオマーカーを抽出し、このバイオマーカーを、１つの時間間隔から次の時間間隔へ境界を自動的に引き渡す態様で、非常に正確かつ再現可能なセグメント化が提供する、４Ｄ物体として経時的に取り扱うことによって、そこから経時的な傾向を検出できることを発見した。このように、それ自体が微妙な病状を捕らえる選択されたバイオマーカーと、経時的な正確さおよび信頼性を増大させる４Ｄアプローチとの組合せは、以前には得られなかった感度を生み出す。
【００７３】
本発明の好ましい一実施形態を説明してきたが、本明細書の開示を検討した当業者なら、本発明の範囲内で他の実施形態を実現できることを容易に理解しよう。例えば、適当な任意の画像化技術を使用することができる。したがって、本発明は、請求項によってのみ限定されるものと解釈すべきである。
【図面の簡単な説明】
【００７４】
【図１】本発明のプロセスの好ましい実施形態の概要を示す流れ図である。
【図２】図１のプロセスで使用されるセグメント化プロセスの流れ図である。
【図３】経時的に撮影した複数の画像中でセグメント化された１つの画像を追跡するプロセスを示す図である。
【図４】図１〜３のプロセスを実装することができるシステムのブロック図である。
【図５】好ましい実施形態に従って形成されたバイオマーカーの画像を示す図である。

Claims

患者の関節を評価するための方法であって、
（ａ）前記関節の少なくとも１つの３次元画像を撮影すること、
（ｂ）前記少なくとも１つの３次元画像中の少なくとも１つのバイオマーカーを識別すること、
（ｃ）前記少なくとも１つのバイオマーカーの少なくとも１つの定量的測定値を得ること、および
（ｄ）前記少なくとも１つのバイオマーカーの識別と前記少なくとも１つの定量的測定値とを記憶媒体に記憶すること
を備えることを特徴とする方法。
前記ステップ（ｄ）は、前記少なくとも１つの３次元画像を前記記憶媒体に記憶することを備えることを特徴とする請求項１に記載の方法。
前記ステップ（ｂ）は、前記少なくとも１つのバイオマーカーを識別するために、前記少なくとも１つの３次元画像を統計的にセグメント化することを備えることを特徴とする請求項１に記載の方法。
前記少なくとも１つの３次元画像は、経時的に撮影した前記関節の複数の３次元画像を備えることを特徴とする請求項１に記載の方法。
前記ステップ（ｂ）は、前記少なくとも１つのバイオマーカーを識別するために、前記複数の３次元画像から選択された３次元画像を統計的にセグメント化することを備えることを特徴とする請求項４に記載の方法。
前記ステップ（ｂ）はさらに、前記選択された３次元画像中で識別された前記少なくとも１つのバイオマーカーに従って、前記複数の３次元画像中で前記少なくとも１つのバイオマーカーを識別するために、動き追跡および推定を備えることを特徴とする請求項５に記載の方法。
前記複数の３次元画像および前記複数の３次元画像中で識別された前記少なくとも１つのバイオマーカーを使用して、前記関節および前記少なくとも１つのバイオマーカーのモデルを、３次元空間および１つの時間次元で形成することを特徴とする請求項６に記載の方法。
前記バイオマーカーを前記モデル中で経時的に追跡することを特徴とする請求項７に記載の方法。
前記少なくとも１つの３次元画像の３つの次元の分解能はすべて１ｍｍよりも細かいことを特徴とする請求項１に記載の方法。
前記少なくとも１つの定量的測定値は高次定量的測定値を備えることを特徴とする請求項９に記載の方法。
前記高次定量的測定値は、湾曲、トポロジーおよび形状のうちの少なくとも１つを備えることを特徴とする請求項１０に記載の方法。
前記少なくとも１つのバイオマーカーは、
軟骨下骨板の形状、
軟骨層およびそれらの相対寸法、
軟骨層内での信号強度分布、
関節軟骨表面間の接触面積、
軟骨形状の表面トポロジー、
骨髄水腫の強度、
骨間分離距離、
半月板の形状、
半月板の表面積、
半月板と軟骨の接触面積、
軟骨構造特性、
軟骨表面特性、
半月板構造特性、
半月板表面特性、
パンヌス構造特性、
関節液特性、
骨増殖体特性、
骨特性、
溶解性病変特性、
プロテーゼ接触特性、
プロテーゼの摩耗、
関節間隔特性、
脛骨内側軟骨の体積、
脛骨外側軟骨の体積、
大腿骨軟骨の体積、
膝蓋骨軟骨の体積、
脛骨内側軟骨の湾曲、
脛骨外側軟骨の湾曲、
大腿骨軟骨の湾曲、
膝蓋骨軟骨の湾曲、
軟骨屈曲エネルギー、
軟骨下骨板の湾曲、
軟骨下骨板屈曲エネルギー、
半月板の体積、
骨増殖体の体積、
軟骨Ｔ２病変の体積、
骨髄水腫の体積および数、
滑液の体積、
滑膜の厚化、
軟骨下骨嚢胞の体積、
脛骨の運動学的平行移動、
脛骨の運動学的回転、
脛骨の運動学的外反、
椎体間の距離、
骨組みの減退の程度、
角度測定による脊柱前湾症の程度、
椎体間のオフセットの程度、
大腿骨特性、および
膝蓋骨特性
からなるグループから選択されることを特徴とする請求項１に記載の方法。
ステップ（ａ）は磁気共鳴画像法によって実施されることを特徴とする請求項１に記載の方法。
患者の関節を評価するためのシステムであって、
（ａ）前記関節の少なくとも１つの３次元画像を受け取る入力装置と、
（ｂ）前記入力装置と通信して、前記関節の前記少なくとも１つの３次元画像を受け取り、前記少なくとも１つの３次元画像中の少なくとも１つのバイオマーカーを識別し、前記少なくとも１つのバイオマーカーの少なくとも１つの定量的測定値を得るプロセッサと、
（ｃ）前記プロセッサと通信して、前記少なくとも１つの３次元画像と、前記少なくとも１つのバイオマーカーの識別と、前記少なくとも１つの定量的測定値とを格納する記憶装置と、
（ｄ）前記少なくとも１つの３次元画像と、前記少なくとも１つのバイオマーカーの前記識別と、前記少なくとも１つの定量的測定値とを表示する出力装置と
を備えることを特徴とするシステム。
前記記憶装置はさらに、前記少なくとも１つの３次元画像を記憶することを特徴とする請求項１４に記載のシステム。
前記プロセッサは、前記少なくとも１つの３次元画像を統計的にセグメント化することによって前記少なくとも１つのバイオマーカーを識別することを特徴とする請求項１４に記載のシステム。
前記少なくとも１つの３次元画像は、経時的に撮影した前記関節の複数の３次元画像を備えることを特徴とする請求項１４に記載のシステム。
前記プロセッサは、前記複数の３次元画像から選択された３次元画像を統計的にセグメント化することによって前記少なくとも１つのバイオマーカーを識別することを特徴とする請求項１７に記載のシステム。
前記プロセッサは、前記選択された３次元画像中で識別された前記少なくとも１つのバイオマーカーに従って、前記複数の３次元画像中で前記少なくとも１つのバイオマーカーを識別するために、動き追跡および推定を使用することを特徴とする請求項１８に記載のシステム。
前記複数の３次元画像および前記複数の３次元画像中で識別された前記少なくとも１つのバイオマーカーを使用して、前記関節および前記少なくとも１つのバイオマーカーのモデルを、３次元空間および１つの時間次元で形成することを特徴とする請求項１９に記載のシステム。
前記少なくとも１つの３次元画像の３つの次元の分解能はすべて１ｍｍよりも細かいことを特徴とする請求項１４に記載のシステム。
前記少なくとも１つの定量的測定値は高次定量的測定値を備えることを特徴とする請求項１４に記載のシステム。
前記高次定量的測定値は、湾曲、トポロジーおよび形状のうちの少なくとも１つを備えることを特徴とする請求項２２に記載のシステム。
前記少なくとも１つのバイオマーカーは、
軟骨下骨板の形状、
軟骨層およびそれらの相対寸法、
軟骨層内での信号強度分布、
関節軟骨表面間の接触面積、
軟骨形状の表面トポロジー、
骨髄水腫の強度、
骨間分離距離、
半月板の形状、
半月板の表面積、
半月板と軟骨の接触面積、
軟骨構造特性、
軟骨表面特性、
半月板構造特性、
半月板表面特性、
パンヌス構造特性、
関節液特性、
骨増殖体特性、
骨特性、
溶解性病変特性、
プロテーゼ接触特性、
プロテーゼの摩耗、
関節間隔特性、
脛骨内側軟骨の体積、
脛骨外側軟骨の体積、
大腿骨軟骨の体積、
膝蓋骨軟骨の体積、
脛骨内側軟骨の湾曲、
脛骨外側軟骨の湾曲、
大腿骨軟骨の湾曲、
膝蓋骨軟骨の湾曲、
軟骨屈曲エネルギー、
軟骨下骨板の湾曲、
軟骨下骨板屈曲エネルギー、
半月板の体積、
骨増殖体の体積、
軟骨Ｔ２病変の体積、
骨髄水腫の体積および数、
滑液の体積、
滑膜の厚化、
軟骨下骨嚢胞の体積、
脛骨の運動性平行移動、
脛骨の運動性回転、
脛骨の運動性バルカス、
椎体間の距離、
ケージの沈み込みの程度、
角度測定による脊柱前湾症の程度、
椎体間のオフセットの程度、
大腿骨特性、および
膝蓋骨特性
からなるグループから選択されることを特徴とする請求項１４に記載のシステム。