• 検索結果がありません。

画像再構成

ドキュメント内 博士論文 (ページ 38-52)

5.1 filtered back projection

filterd back projection (FBP) 法は,単純逆投影法におけるボケを改善するため補正

filterを掛けて逆投影する方法である.filter処理を実空間で行う重畳積分逆投影法と周波

数空間で行う方法がある.現在主流である重畳積分逆投影法は,実空間上で各方向からの 投影データに補正関数 (filter) を重畳して新しい補正投影データを作り,逆投影を行う.

補正関数として,Ramachandran filter,Shepp&Logan filter,Chesler filterなどが用いられる [48].この再構成法では,高集積部が存在するとストリークアーチファクトを発生させ,

近傍の集積がマイナス値になるなどの問題がある[48].

被者体のある断面位置 (x,y) における放射性核種の濃度は d(x,y)で表わされる.x-y平面において,x軸に平行なy軸方向の投影データP(x) は,x軸に対するy軸の積分で あり,式 (5.1)で表わすことができる.

𝑃(𝑥) = ) 𝑑U𝑥,𝑦Y𝑑𝑦 (5.1)

x-y座標系を角度θ回転させたr-t座標系におけるt方向の投影データをp (r,θ)とす る.f (x,y)をp(r,θ)に対応させる変換式は式 (5.2)で表わされる.

𝑝U𝑟,𝜃Y = ) ) 𝑓U𝑥,𝑦Y𝛿(𝑥 cos 𝜃 + 𝑦 sin 𝜃 − 𝑟)𝑑𝑥𝑑𝑦

ƒ ƒ

(5.2)

δは原点では無限大となりそれ以外ではゼロとなるデルタ関数であり,rは式 (5.3)で表わ される.

𝑟 = 𝑥 𝑐𝑜𝑠 𝜃 + 𝑦 𝑠𝑖𝑛 𝜃 (5.3)

5.2 投影データ収集の概略

層像を再構成することを利用し,断面の各位置における濃度関数を復元することによ り,d(x,y) は式 (5.4)で表わされる.

𝑑U𝑥,𝑦Y = ) 𝑝(𝑟,𝜃)

-(5.4)

図 5.3 重畳積分逆投影法の概略

しかしながら,多数の投影データを逆投影することで直線の合成で作成した濃度に等 しくなるため,逆投影画像は真の画像に対してボケを有する.ボケを修正し,逆投影を行 う方法として重畳積分逆投影法を用いる.重畳積分逆投影法は式(5.5)で表わされる (図 5.3).ここで,⊗は重畳積分を示している.

𝑑U𝑥,𝑦Y = ) 𝑝U𝑟𝜃Y

-⊗ 𝑔(𝑟)𝑑𝜃 (5.5)

投影データと補正フィルタの重畳積分を,周波数空間内で行う方法としてフィルタ補正逆 投影法があり,重畳積分逆投影法と等価な方法である.ボケの補正フィルタには,Ramp フィルタ,Ramachandranフィルタ,Shepp-Loganフィルタがあり,目的に応じて使い分け る.

5.2 maximum likelihood expectation maximization

この方法は,統計学的推定に基づいた解析的手法である.取得した投影データになる ための最も確からしい値を割り当てる方法である[49].期待値最大化法を用い,まず初期 値を仮定し,その結果を用いて次々に計算を繰り返し (iteration) 最も確からしい値を算 出する.ML-EM法は全ての投影方向のデータを同時に利用するため計算に時間がかか る.FBP法と比較して,高いカウント領域でのストリークアーチファクトは出現せず,

低カウント領域でのS/Nが改善される.しかしながら,FBP 法とは異なる方法であるた め,得られる値が必ずしも一致しないため,定量検査においては注意が必要である.

xijを画素jから出たγ線が検出器iで検出される計測値とし,検出器iにおける投影 データをyiとすると,yiは式(5.6) で表される (図5.4).

𝑦 = ˆ 𝑥‡‰

(5.6)

図 5.4 画素jを通り検出器iで検出される過程

xijの期待値は,画素jにおける画像の値λi,画素jから出たγ線が検出される確率Cij

とすると式(5.7) で表すことができる.

𝐸 (𝑥𝑖𝑗) = 𝐶𝑖𝑗𝜆𝑖 (5.7)

xijの統計的変動がポアソン分布に従う場合,画像の値λに対してxが測定される条件 付き確率P(x|λ) は式(5.8) で表すことができる (尤度関数).

𝑃(𝑥|𝜆) = Ž Ž 𝑒𝑥𝑝U−𝐶‡‰𝜆Y (𝐶‡‰𝜆)•‘

𝑥‡‰!

(5.8)

ML-EM法による画像再構成は式 (5.9) で表される.

𝜆= 𝜆”AH

∑ 𝐶 ‡‰ˆ 𝐶‡‰𝑦

∑ 𝐶B ‡‰𝜆B”AH

(5.9)

ここで,kは繰り返し回数である.jは再構成画像の座標であり,1から画素の最後mま で通し番号をとる例えばマトリクスサイズが64×64ならj = 1~4096である.iは検出器上 の画素番号であり,角度方向のデータも含めてn個の一連のデータである.投影方向数 nvが72ならi = 1~4608 (= 64×72) となる.

λjはある画素jに対する逆投影 (back-projection) のRI濃度 (カウント) であり,λj’は 画素jを通り検出器iに入る光子の投影 (forward-projection) のカウントである.yiは検出 器iに入射した投影カウントであり,Cijは画素jから出た光子 (γ線) が検出器iに到達す る割合 (確率) である.yiは実際に測定される投影データである.

式 (5.9) は,forward-projection,back-projection,総確率での規格化,再構成値の更新 の計算要素に分解できる.

forward-projectionの計算項は,式 (5.9)の第二項の分母であり,式 (5.10)である.

ˆ 𝐶‡‰—

˜

‰—™H

𝜆B (5.10)

back-projectionの項は,式 (5.9)の第二項の分子であり,式 (5.11)である.

𝑦𝐶‡‰ (5.11)

総確率での規格化する項は,式 (5.9)の第一項の分母であり,式 (5.12)である.

1

‡™H𝐶‡‰ (5.12)

式 (5.9)は画素jに着目したカウントの式であることから,再構成を行うにはすべて の画素についてこの演算を行う.

計算ステップを以下に示す.

(1) 検出確率Cijを計算する (最初に一回だけ計算する).

(2) 初期画像を仮定する.

(3) 検出器iから画素jに対するBack-projectionを計算する.

(4) 画素jを通り検出器iに入る光子のForward-projectionを計算する.

(5) ステップ (3) (4) をすべての投影角度について計算し,両者の比を取って加算する.

(6) 全確率で規格化する.

(7) 以上の計算値を初期画像の画素jに掛け算して更新画像を作成する.

(8) 次の画素 (j + 1) に移ってステップ (3) に戻る.すべての画素の計算が終わったら更 新画像を初期画像としてステップ (3) に戻る.

さらに,検出確率に減弱,散乱,コリメータ開口径に対する補正項を織り込むことで 定量性を向上させることができる.計算ループを繰り返すことによってλはRI分布画像 に近づいていく.計算ループの打ち切りに関しては明確なルールが存在しないため,実験 的 (経験的) に行っているのが現状である.また初期値としては「正の値であること」と いう制限はあるが,一様分布を仮定することで対処する.

5.3 ordered subset expectation maximization

ordered subset expectation maximization (OSEM) 法は,ML-EM法の投影データをいく つかのグループ (subset) に分けることにより計算を高速化したものであ[50].あるsubset の投影データを計算し初期値し,次のsubsetで計算をし,次々に再構成を繰り返す

(iteration).これにより,データ収束のために時間が短縮される.Subsetとiterationの理論

脳部SPECT画像再構成に対してOSEM法が用いられている.OSEM法は,ML-EM 法の原理を利用し,投影データをいくつかのグループ (Subset) に分類し,そのSubset内 の投影データごとに繰り返し (Iteration) 修正することにより演算速度を高速化した方法 である [45-47].セット内の投影データは,Subset間の影響が均等化されるように,一定 の離れた角度の投影データが選択され,また規則的に計算順番が決められる.

図 5.5 投影データが18,Subset数を3とした時の計算順番例

Subset数は投影データ数を割り切れる数になり,例えば,図5.5の投影データが18

では,Subset数を3にすると1つのSubsetには 6の投影データが属する.演算の順序は 投影データS1より画像を修正し,その次にS1から最大角度離れたデータS2で修正し,

順にS6まで演算すると1回の近似が終了する.その反復は式 (5.13)で表すことができ る.

𝜆BšH= 𝜆B+ 𝜆B

‡›œ𝐶‡‰ˆ 𝐶‡‰

‡∈œ Ÿ 𝑦

∑ 𝐶˜ ‡‰ 𝜆B¡ (5.13) ここで,𝑆˜は1つのsubsetに属する投影データの集合を表す.

計算手順を以下に示す.

(1) 検出確率𝐶‡‰を計算する.

(2) 初期画像を仮定する.

(3) 検出器iから画素jに対するBack-projectionを計算する.

(4) 画素jを通り検出器iに入る光子の投影Forward-projectionを計算する.

(5) ステップ(3),(4)をSubsetに属する角度のみで計算し,両者の比を取って加算する.

(6) 全確率で規格化する (ただし,nはSubset内データ数).

(7) 以上の計算値を初期画像の画素jに掛け算して更新画像を作成する.

(8) 次の画素 (j+1) に移ってステップ (3) に戻る.すべての画素の計算が終わったら更新 画像を初期画像として,次のSubsetに移りステップ (3) に戻る.すべてのSubsetの計 算が終わったら,同様に更新画像を初期画像として,ステップ (3) に戻る.

5.4 3-dimensional stereotaxic region of interest template

3DSRTは,statistical parametric mapping (SPM) を用いて解剖学的標準化を行った

SPECT画像上に,脳表のみの血管支配領域を内部にも拡充した普遍的ROI群 (三次元定

位ROIテンプレート) を構築している [51,52].SPMは,ロンドンのハマースミス病院

のFristonらによって開発され,脳核医学において,統計学的に脳部SPECT画像の解剖学

的標準化を行う[48].3方向の平行移動と回転により位置あわせ,線形,非線形変換によ る標準脳への解剖学的標準化,フィルタ処理による平滑化後に,統計解析を行う.脳血流 定量解析において,群間比較を容易に行えるため,非常に有用である.このROI群を設 定することで,解析過程に由来する同一データでのバラつきや解剖学的多様性に由来する 被験者間,観察者間,施設間のバラつきを排した客観的定量評価が視覚的にも数値解析的 にも可能となる.本研究で使用する3DSRTは,一側につき12区域 (脳梁辺縁,中心前,

中心,頭頂,角回,側頭,後大脳,脳梁周囲,レンズ核,視床,海馬,小脳半球) から 構成される左右計24領域の三次元ROIテンプレートとなっている.

5.6 3DSRTによるROIテンプレート

(Takeuchi R,Yonekura YMatsuda Het al: Eur J Nucl Med29 (3)331-3412002.)

5.5 減弱補正

体表面に近い線源と深部ではγ線の内部吸収がことなるため,減弱補正が必要とな る.この補正には,従来,画像再構成前に投影データを補正する前補正法のSorenson 法,再構成後の画像に対して補正する後補正法のChang法が用いられる[38,53].これら の方法は,人体を均一な吸収体と仮定し一定の減弱係数を設定する.簡便であるが,頭部 は脳実質と高吸収体の頭蓋骨が存在する不均一吸収体であるため,正確な補正が困難であ る.精度の高い減弱補正を行うには外部線源やX線CTを用いて被写体内部の各部分の線 減弱係数 µ (cm-1) の分布を測定する方法 (CTAC法) がある[38,54].体内の透過率を測 定するCTデータを用いた減弱補正法では,CTの高空間分解能により精度の高い減弱 mapを作成することが可能である[38,54].

Chang法は,均一吸収体補正法の中でも広く臨床で用いられている.この方法は,補

正を行なう前にあらかじめ補正行列を作成し,実測した投影データから再構成した画像 に,後から補正行列を掛け合わせることで補正する.

5.7 Chang法の原理

図5.7に示すように,一様な線源弱係数 (µ) を持つ断面を想定し,この断面に正確に 輪郭を設定する.この断面内に放射能強度がI0となる線源を仮定する.この線源から輪 郭までの距離は線源と投影方向により ( l1,l2,l3) と変化し,投影方向 (Projection 1) で 検出される線源強度は距離l1と線減弱係数µにより減弱され,Io exp (-µl1)となる.

Projection 2およびProjection 3も距離が異なるのみで同様である.投影方向がN方向のと

き,∑¢‡™H𝐼-exp (−𝜇𝑙)となるため,補正行列C(xiyi)は,

𝐶(𝑥,𝑦) = N

¢‡™H𝐼-exp (−𝜇𝑙) (5.14)

となる.この補正マトリックスの各点は断層画面上の対応する点からの各投影方向への

図5.7では,点線源の場合を考慮したが,実際には線源は一定の大きさを持つことか ら,投影方向毎にこの位置をずらして,全ての線源位置での補正係数を作成できる.

5.6 散乱線補正

散乱線は,γ線と人体構成物質との相互作用であるコンプトン効果による.複数回散 乱したγ線はエネルギーが減少してenergy windowにより除外されるが,一時散乱線は

energy windowに混入するため,脳SPECT画像の画質,定量性向上には減弱補正と同時に

散乱線補正が重要である.散乱線の補正には,energy windowの設定に基づくものとし て,main windowの低エネルギー側にsub windowを設定するdual energy window (DEW) 法と両側に設定するTriple energy window (TEW) 法があり,減弱係数分布から散乱成分を 推定する方法は,投影データの重畳積分からもとめるtransmission dependent convolution subtraction (TDCS) 法がある.[38,54]

TEW法は,多くの施設で採用されている.図5.8に示すように,main window (W1) および,2つのsub window (W2,W3) をW1を挟み込むように設定する.main window内 の散乱線成分はwindow幅と2つのsub windowのwindow幅の違いを考慮し,2つのsub windowからの投影𝑝C(𝑟,𝜃),𝑝©(𝑟𝜃)から推定する.図5.8のように台形の面積により散 乱線成分を計算し,main windowの投影𝑝HU𝑟,𝜃Y から減算する.

𝑝U𝑟,𝜃Y = 𝑝HU𝑟,𝜃Y − ª𝑝C(𝑟,𝜃)

𝑊C +𝑝©(𝑟,𝜃) 𝑊© « ×𝑊H

2

(5.15)

5.8 散乱線補正 (TEW)

5.7 分解能補正

核医学では,コリメータにより,患者から放出されるγ線を一定方向から入射するも のだけに限定している.コリメータの孔径は一定の大きさをもつため,距離が離れると空 間分解能が劣化する.コリメータ開口に伴う補正には,周波数と距離の関係を用いる方法 と逐次近似画像再構成に組み込んだ方法がある[38,54].

Collimator 開口径補正付き逐次近似法はML-EM基本式内の検出確率𝐶‡‰の部分に

collimator broad correction (CBC) の項を組み込んで行う.

𝜆BšH= 𝜆

”‰‡™H𝐶‡‰¬ 𝑦𝐶‡‰

˜™H𝐶‡‰ 𝜆B

‡™H

(5.16)

ここで,𝜆šHは修正される1pixel分の放射能,𝜆は仮定された全pixelの放射能,𝑦 は投

ドキュメント内 博士論文 (ページ 38-52)

関連したドキュメント