整形外科は四肢と躯幹の骨・関節・軟部組織を対象とし、そのいずれもが3次元的な形態を有している。治療対象はその形態に変形が生じた場合が多く、3次元的変形を経時的に追跡して計測する手法が必要になる。

　しかし従来整形外科領域の画像計測では撮影対象物の2次元投影像に過ぎない単純レントゲン写真が用いられることが多かったため、撮影時に肢位が変化した場合には同一の像が得られず、3次元的変形の正確な追跡には不向きであった。一方CTやMRIなどの3次元画像は3次元的位置情報を有しているため正確な3次元的計測に使用可能なはずである。しかし2次元画像では計測対象部位が同一平面上にあるか否かにかかわらず1枚の平面上に表示されるため観察者にとって把握が容易であるのに対して、3次元画像は通常多数枚の2次元断層像の集合として表示されるため、対象部位が同一断層内に入っていない場合は把握が困難になる。

　任意の測定対象部位を経時的に追跡するためには、初回に撮影された画像上の測定対象部位が、異なる時期に撮影された画像上でも同定可能である必要がある。そこで本研究では、コンピュータを用いた画像マッチングという手法を用いることで初回画像上の任意点が他の画像上でどの点に対応するかを求める手法を開発した。また経時的画像を用いて追跡精度の評価を行った。

　本研究で用いるマッチングとは2つの画像を一致させる手法であり、rigid matchingとdeformable matchingに分けられる。前者は平行移動・回転移動の組み合わせで各画素間の位置関係は変えないのに対し、後者は各画素間の位置関係の変化を許している。後者はさらに対象画像に対して仮定する変形規則によって分類されるが、本研究では計算量及び生体の特性を考慮し、弾性体と仮定するelastic matchingを用いた。これは2つの対象画像のうち一方を弾性体と仮定し、画像を弾性変形させて一致させるものであり、これにより変形を伴う画像間でのマッチングが可能となる。Fig.1にマッチングを用いた変形追跡の概略を示す。

Fig.1 初回画像上の点のマッチングによる追跡。初回画像をマッチングにより2回目画像に一致させることで、初回画像上の測定対象点●の座標を2回目画像上での対応点◆の座標に変換する。その点の座標を2回目画像に挿入することにより、初回画像上の測定対象点を2回目画像上で同定することが可能になる。

　各マッチングは画像の一致度を表す式であるエネルギー関数を求めてそれを最適化(最小化)することで行う。各マッチングの方法により関数が異なるため以下では分けて解説する。なお画像の一致度としては2つの画像の対応する点同志の画素値の差の二乗平均(以下E_d)を用いた。

　またrigid、elasticの両マッチングにおいてエネルギー関数の最適化に最急降下法を用いた。これは関数値Eがその偏微分のベクトルであるgradient(ここでは∇E)の反対方向に進むと最も減少率が高いことを利用した方法である。

　rigid matchingではエネルギー関数はE_dと等しく、初回画像を1つの剛体と考えて回転移動と平行移動を行うことで最適化した。

　elastic matchingでは初回画像を弾性体と仮定し、全画素をx、y、z方向に個々に移動させてエネルギー関数(以下でE_eと表す)を最適化した。この仮定は画像の断裂や極端な体積変化など立体としての連続性を無視した変形を防ぐことを目的とし、弾性変形時の歪みエネルギーを制御項としてRで表した。E_eを画像の一致度E_dと制御項Rの和とすることで画像一致の際の変形制御が可能になる。

　E_dは、初回及び2回目画像上の点xにおける画素値をv0(x)、v1(x)とし,初回画像上の画素の変位ベクトルをu、画素数をNとすると次式で表せる。

　Rは、対象とする画像をボアソン比vの線形弾性体と仮定し、横弾性係数をG垂直歪みを、剪断歪みをとおくと下式で表せる。

　なおポアソン比は形態変化に伴って体積が変化する物質ではl/4〜1/3であることから、今回の実験では1/3と仮定した。

　以上の式からE_eは下式で表される。

　対象画像上で画像内の構造が重ならない場合にもマッチングを行うために、原画像から低解像度の画像を作成してそこからマッチングを開始し、解像度を上げていくmultiresolution techniqueという手法を用いた。これは局所解(エネルギー関数が極小にはなるが最小ではない解)の回避にも役立つ。今回は解像度がl/16の画像から開始し、1/8、1/4、1/2、原画像の順でマッチングを行うことで、最大16ボクセル未満の変形までは対応可能とした。

　以上のアルゴリズムのフローチャートをFig.2に示す。

Fig.2 マッチングのアルゴリズム。左半分はmultiresolution techniqueを、右半分は各解像度におけるマッチングを表している。

　ソフトウエアは以上のアルゴリズムに沿ってC++言語で作成し、作成・実行共にパソコン上で行った(Fig.3)。対象画像を2つ選んだ後のマッチングは自動的に行われるようにし、初回画像上の追跡対象部位をマウスでクリックすることで初回画像と2回目画像での座標が表示されるようにした。

Fig.3 作成したソフトウェア。すべての操作はマウスのみで可能である。図は膝関節ACLを撮影したMRIのマッチング結果であり、表示断面は前額、冠状、矢状断面のいずれも選択可能である。

　1ケ月間隔で2回撮影された2種類の画像を用いた変形追跡の精度測定及び経時的変化の追跡を行った。追跡精度は、任意に作成した変換式を用いて非線形的に変形させた3次元画像を、元の画像に対するマッチングにより復元し、その際の復元精度の測定により求めた。(Fig.4)

Fig.4 復元時の誤差の測定法。元の画像上の点P₀が非線形変形、elastic matchingによりそれぞれ点P₁、P₂と写るとき点P₀とP₂の位置の差が復元誤差となる。

1.膝関節前十字靭帯(ACL)の長さ測定

　対象は33歳男性の膝関節MRI(FOV200×200mm,解像度256×256、スライス厚2mm×48枚)とした。

　精度実験では長さと角度について測定した(Tabiel)。長さの復元誤差は平均0.4mm(標準偏差＝1.3mm)で画素サイズ(最大長2.3mm)より小さかった。

Table 1 膝関節MRIを用いた精度実験の結果

　経時的変化の測定では、まず従来のように肉眼的に初回及び2回目の画像上でACL骨付着部の座標を求め、次いで従来法により求められた初回画像上のACL骨付着部が2回目画像上で対応する点を今回開発したマッチングを用いる手法により求めた。

　結果をTable2に示す。計算には1時間57分を要した。従来法ではマッチングを用いた手法に比べて脛骨側付着部の測定値に1画素以上の差があり、初回と2回目で異なる点の測定を行っていた可能性があった。ACLの長さの初回と2回目画像での差は両手法で画素サイズ以下であり経時的変化はなかったことになる。

Table 2 ACL骨付着部の座標の測定結果(mm)

　Fig.5でACL骨付着部を含むマッチング後の初回画像の断面A1と、ここで得られた骨付着部の座標を含む2回目画像の断面B1はほぼ同一の像となっており、ACLの同一と思われる断面も得られている。

Fig.5膝関節ACL骨付着部及びACL断面の追跡結果。初回画像(A₀)を2回目画像にマッチングさせることで骨付着部の2回目画像内での座標が求められ(A₁)、それをもとにしてACLの同一断面が得られている(B₁)。

2.脛骨遠位部骨梁構造追跡

　対象は34歳男性の脛骨遠位部pQCT画像(FOV51.2×51.2mm,解像度512×512.スライス厚1mm×16枚)とした。

　精度実験では骨梁構造には肉眼で認識可能な構造がないため、画像上に無作為に選んだ48点につき測定した(Table3)。平均復元誤差は長さ0.48mmで画素サイズ(最大長1.01mm)より小さかった。

Table 3 pQCTを用いた基礎実験の結果

　経時的画像の追跡には2時間57分を要した。Fig.6は経時的画像を用いた追跡例であり、初回画像はマッチングにより2回目画像に骨全体の位置および骨梁構造がほぼ一致している。

Fig.6 経時的pQCT画像のマッチング例。初回画像と、1ヵ月後画像との骨梁構造はマッチングにより骨梁構造までほぼ一致している。

　上記の2種類の画像に対する精度実験の結果ではマッチング精度は共にボクセルサイズよりも小さく、画像を用いた測定法としては十分であった。マッチングに要した時間は2〜3時間と、操作が全自動であることを考慮すれば実用可能な範囲内であった。

　経時的に撮影された骨関節3次元画像上の任意点の追跡を、画像のマッチングによって可能にし、従来は定性的であった経時的変形の測定を定量化する手法を構築した。

　また経時的画像においてその追跡精度を測定し、十分な精度が実現されていることが確認された。