臓器動態の定量化方法、装置、臓器位置の予測方法、装置、放射線照射方法、装置及び臓器異常検出装置
【課題】金属マーカーでは不可能な、臓器内部領域全体にわたって任意の点での局所的な動態追跡を可能とする。
【解決手段】微小に変動・変形した臓器（１０）の複数のＣＴ画像を入力データとし、臓器内に分布した血管や気管１２を抽出して細線化し、細線化した画像１４を用いて分岐点１６や連結点の座標を抽出し、これらを特徴点として複数のＣＴ画像間の個々の点の移動を３次元空間中で追跡することにより、臓器（１０）の連動を計測する。
【解決手段】微小に変動・変形した臓器（１０）の複数のＣＴ画像を入力データとし、臓器内に分布した血管や気管１２を抽出して細線化し、細線化した画像１４を用いて分岐点１６や連結点の座標を抽出し、これらを特徴点として複数のＣＴ画像間の個々の点の移動を３次元空間中で追跡することにより、臓器（１０）の連動を計測する。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、放射線診断・治療計画及び治療において、呼吸等に伴って変動・変形する肺等の臓器に対する臓器内各部の対応付け（位置合わせ）や、腫瘍等の関心体積の変形や位置を、呼吸波形の位相との相関として把握するために、高速３次元ＣＴ（いわゆる４次元ＣＴ）を用いて、ある時間間隔で取得された各時刻の３次元デジタル画像を用いて、各画素位置の動態を定量化するための臓器動態の定量化方法、装置、これを利用した臓器位置の予測方法、装置、放射線照射方法、装置、及び、臓器異常検出装置に関する。
【背景技術】
【０００２】
放射線治療の照射対象となる腫瘍の位置は、呼吸や日々の臓器の変動・変形により動く。放射線治療の際には、決められた位置に照射目標となる臓器（ターゲット）が来るように、一般にＸ線透過像（２次元画像）を用いて位置合わせを行なっている。
【０００３】
呼吸に伴って動く臓器に対しては、例えば呼吸同期法を用いて、呼吸波形のある領域のみ照射を行なうことにより、呼吸移動により生じる余分な照射を減らしている。
【０００４】
しかし、この呼吸波形は、体表面上に設定した、ある点の１次元的な位置を示しているものであり、又、その照射中にも臓器は３次元的に動いている。従って、より高精度な照射を行なうためには、臓器の各位置の動態（運動の方向及び速度）を把握しなければならないが、これを２次元画像から解析することは困難である。
【０００５】
特に、粒子線治療において、肺のように腫瘍と正常組織の密度が大きく異なる場合、照射ビーム（以下単にビームと称する）に対して垂直な方向だけでなく、ビーム軸方向の臓器の変動がビームの飛程に大きく影響し、その結果、腫瘍への過小照射、正常組織への過大照射につながる危険性がある。そのため、臓器の動態を３次元的に把握する技術が求められているが、現在臨床では、このような技術は確立されていない。
【０００６】
近年、臓器の３次元的な位置を定量化する研究が、いくつか報告されている。それらの多くは、Ｘ線ＣＴ画像を用いて抽出された対象臓器の表面上に、動きを追うための目印となる参照点を設定し、２つのＣＴ画像上の参照点同士が一致するように、臓器を弾性体モデルとして、その変形を計算している。
【０００７】
しかしながら、この場合、特徴点を人為的に設定しなければならず、又、臓器内部の各位置の対応は必ずしも現実の動きを反映していない。又、具体的な動きの追跡のために金属マーカーが一般的に用いられているが、患者の体内へ埋め込むため、負担が大きく、部位や個数が制限され、臓器の動態を細かく捉えるのが困難である。
【０００８】
一方、非特許文献１には、有限要素法により肺の３次元的な変形を計算することが記載されている。
【０００９】
又、非特許文献２では、肺の３次元ＣＴ画像の画素値の勾配を比較することにより、内部の変形を定量化しており、腫瘍部分の関心体積の変形についても論じている。
【００１０】
更に、非特許文献３には、呼吸運動に伴う肺の変形での臓器内の対応付けのために、血管分岐点の動きを追跡することが記載されている。
【００１１】
【非特許文献１】Ｔ．Ｚhangら「Ｔechnical note：Ａ novel boundary condition using contact elements for finite element based deformable image registration」Ｍedical Ｐhysics ３１（２００４）２４１２−２４１５
【非特許文献２】Ｄ．Ｓarrutら「Ｎonrigid registration method to assess reproducibility of breath-holding with ＡＢＣ in lung cancer」Ｉnternational Ｊournal of Ｒadiation Ｏncology，Ｂiology，Ｐhysics，６１（２００５）５９４−６０７
【非特許文献３】四方秀則ら「肺の非剛体レジストレーションのための血管分岐点位置決めアルゴリズム」 電子情報通信学会論文誌Ｄ−II Ｖol.Ｊ８４-Ｄ-II Ｎo.０７ pp．１−１１
【発明の開示】
【発明が解決しようとする課題】
【００１２】
しかしながら、非特許文献１に記載の技術は、臓器の外形と内部の変形が必ずしも一致していないという問題点を有する。
【００１３】
又、非特許文献２に記載の技術は、内部の解剖学的特徴を追っているわけではなく、関心体積変形の正確さについては問題がある。
【００１４】
又、非特許文献３に記載の技術は、単に、血管分岐点の位置決めを行なっているだけである。
【００１５】
本発明は、前記従来の問題点を解消するべくなされたもので、臓器の動態を正確に定量化することを課題とする。
【課題を解決するための手段】
【００１６】
本発明は、微少に変動・変形した臓器の複数のＣＴ画像を入力データとし、臓器内に分布した血管や気管を抽出して細線化し、細線化した画像を用いて分岐点や連結点の座標を抽出し、これらを特徴点として複数のＣＴ画像間の個々の点の移動を３次元空間中で追跡することにより、臓器の運動を計測するようにして、前記課題を解決したものである。
【００１７】
又、前記特徴点の追跡を、複数のＣＴ画像及び特徴点座標を用いたポイントパターンマッチングにより行なうようにしたものである。
【００１８】
又、複数の特徴点の運動を補間することにより、任意点の運動が分かるようにしたものである。
【００１９】
本発明は、又、前記の方法により得られた動的な情報と、臓器の周期的な運動の情報との相関を取ることにより、臓器の位置を予測することを特徴とする臓器位置の予測方法を提供するものである。
【００２０】
本発明は、又、前記の方法により予測された臓器の位置に、放射線を照射することを特徴とする放射線照射方法を提供するものである。
【００２１】
本発明は、又、微少に変動・変形した臓器の複数のＣＴ画像を入力する手段と、臓器内に分布した血管や気管を抽出して細線化する手段と、細線化した画像を用いて分岐点や連結点の座標を抽出する手段と、これらを特徴点として複数のＣＴ画像間の個々の点の移動を３次元空間中で追跡することにより、臓器の運動を計測する手段と、を備えることにより、前記課題を解決したものである。
【００２２】
又、前記特徴点の追跡を、複数のＣＴ画像及び特徴点座標を用いたポイントパターンマッチングにより行う手段を備えたものである。
【００２３】
又、複数の特徴点の運動を補間して、任意点の運動を求める手段を備えたものである。
【００２４】
本発明は、又、前記の装置により得られた動的な情報と、臓器の周期的な運動の情報との相関を取ることにより、臓器の位置を予測する手段を備えたことを特徴とする臓器位置予測装置を提供するものである。
【００２５】
本発明は、又、前記の装置により予測された臓器の位置に、放射線を照射する手段を備えたことを特徴とする放射線照射装置を提供するものである。
【００２６】
本発明は、又、前記の装置により得られた動的な情報に基づいて、臓器の局所的に異常な変形を捉える手段を備えたことを特徴とする臓器異常検出装置を提供するものである。
【発明の効果】
【００２７】
本発明においては、肺の場合、気管支や血管のように連結する解剖学的特徴に着目し、それらのトポロジー（分岐や連結）は、臓器の運動や変形でも保持されていることを利用して、臓器の動態を定量的に計測する。即ち、動態の前後で対象のトポロジーが保存されることから、それらの特徴点座標を抽出し、個々の点の移動を３次元空間中で追跡することで、局所的にも高い精度を有する臓器の運動（方向、速度）を計測する。
【００２８】
例えば、気管支に沿って存在する血管のＣＴ値が、肺野のＣＴ値に比較して大きいことを利用し、図１（Ａ）に示す如く、ＣＴ画像から血管や気管１２を抽出する。
【００２９】
次いで、抽出された血管や気管部を、細線化アルゴリズムを利用して、太さ１画素の細線化画像１４に変換し、血管や気管の特徴点（図１では分岐点１６）座標を特定し、これを動態追跡の特徴点とする。これにより、従来動きの追跡に用いられる金属マーカーでは不可能な、肺内部領域全体にわたって、任意の点での局所的な動態追跡が可能となる。
【００３０】
次いで、ある時間間隔で取得された２つのＣＴ画像に対して、これらの特徴点の変位を追跡するために、２つのＣＴ画像及び特徴点座標を用いて、ポイントパターンマッチングを行なう。例えば確率的弛緩法を用いてこの操作を行なうことにより、特徴点の変位が、その近傍の特徴点の変位と互いに矛盾が無いように、各特徴点の変位量（方向と長さ）を求めることができる。即ち、抽出した特徴点の追跡から、図１（Ｂ）に示す如く、各点の運動が分かる。
【００３１】
更に、それらの点から３次元空間での補間を行なうことで、任意点の運動が分かる。従って、臓器（図では肺１０）内の全ての画素に対してその変位量が求められる。
【００３２】
以上の操作を連続的に取得されたＣＴ画像に対して逐次行なうことにより、連続的に動態定量化を行なうことができる。従って、動く対象に対する４次元的な線量計算が可能となり、呼吸動態臓器に対する、より高精度な放射線治療を行うことが可能となる。
【００３３】
放射線治療のための治療計画（線量計算、評価による放射線の照射方向・回数等の照射方法を決定する操作）等では、腫瘍の部分を関心体積として設定し、その関心体積やその周りの正常組織への放射線量計算を行なうが、従来は静的な３次元ＣＴ画像に対して行われていたものを、本発明を用いれば、４次元ＣＴを用いて取得されたＣＴ画像を用いて、動的な対象に対してこれを行なうことが可能となり、関心体積の変形が自動的に計算でき、各時刻で取得されたＣＴ画像に対して設定する必要が無くなる。これにより、４次元的な治療計画を、より現実的なものとすることができる。
【００３４】
又、近年開発され、医療現場へ普及しつつある４次元ＣＴから逐次得られる時系列の３次元デジタル画像を用いて、本発明の操作を逐次的に行なうことにより、臓器の３次元的な動きを、ある時間領域にわたって得ることができる。
【００３５】
従って、本発明により得られた動的な情報を、呼吸波形のような参照情報と相関をとることにより、任意の時点での臓器の動き（３次元の方向とその速度）を予測することができる。
【００３６】
又、本発明は、動きの特徴点を得る方法と、その特徴点の動きを追跡する方法に分かれているため、例えば血管造影法を用いて肝臓のＸ線ＣＴ画像を取得する等のように、特徴点を得る方法を工夫すれば、様々な部位に応用することが可能である。
【００３７】
このようにして、臓器の３次元的な動きを把握することにより、従来位置の不確定さのために設定せざるを得なかった照射マージンを、より縮小して、正常組織への不必要な放射線被曝を抑えることができ、晩期障害等の放射線障害を減らすことができる。
【００３８】
臓器の解剖学的な変性は、臓器の運動に歪みや局所的に異常な変形を引き起こすので、本発明を用いて、正常時あるいは正常部分との動態を定量的に比較して、局所的な運動の歪みや異常性を検出することで、医学的診断に応用できる。例えば、ＣＴ値では判別し難い腫瘍の存在や肺正常組織への放射線障害である繊維化の兆候を知る手懸りとなる。
【発明を実施するための最良の形態】
【００３９】
以下図面を参照して、本発明の実施形態を詳細に説明する。
【００４０】
本発明の第１実施形態を実施するための装置は、例えばコンピュータを用いて、図２に示す如く構成され、肺領域抽出手段２０と、細線化手段２２と、分岐点抽出手段２４と、特徴点追跡手段２６と、を含んでいる。
【００４１】
前記細線化手段２２では、あるＣＴ値を閾値として２値化を行ない、それによって得られた１画素値の図形を削っていく。しかし、肺の周りの肉や骨の部分が２値化によって残ってしまうと、それも細線化処理の対象となり、その部分は肺内の気管支や血管に比べて体積（画素数）がかなり大きいため、計算時間が大幅に長くなる。
【００４２】
そこで、計算時間を短くし、又、必要な肺領域のみを計算させるために、前記肺領域抽出手段２０で、まず肺領域を抽出する。
【００４３】
具体的には、図３に示したような方法で行なうことができる。処理は、図３（Ａ）に例示する、それぞれの２次元ＣＴスライス画像毎に行なう。
【００４４】
まず、あるスライス画像について、ある設定した矩形領域でスムージング処理を行う。これを全ての２次元領域について行う。この処理によって、図３（Ｂ）に示す如く、肺内のＣＴ値が平均的に小さい領域と、肺の周りのＣＴ値の高い領域に分けられる。次に、この画像について、図３（Ｃ）に示す如く、ある適当なしきい値により二値化を行い、ＣＴ値の小さい領域のみ抽出する。その後、左右の肺領域を指定することによって体外部など不要領域を削除して、図３（Ｄ）に示す如く、肺領域のみを抽出する。以上の操作をすべてのスライスについて行うことにより、肺領域の３次元画像を得る。ここで得られた領域のみを後の細線化処理の対象とする。
【００４５】
肺領域の抽出結果の例を図４に示す。ここでは、平均化の矩形領域サイズを１５、ＣＴしきい値を−５００とした。
【００４６】
前記細線化手段２２は、厚さ・太さのある図形を、トポロジーを変化させずに太さ１の中心線に対応する線図形に変換する操作（細線化）を行う。この細線化によって、肺内の血管（気管支）の分岐点の座標を画素単位で特定することができ、その結果得られた分岐点を動きの特徴点として追跡することにより、肺の動きを定量化することができる。トポロジーを変えないということは、連続した図形の数や、図形の穴や空洞の数といった図形の性質が変化しないことを意味する。ここでは、ユークリッド距離変換を伴う逐次型アルゴリズムによる細線化法を用いている。
【００４７】
具体的には、まず（１）３次元ＣＴ画像について、適当なＣＴ値を閾値として２値化を行なう。ここでは、閾値より大きい画素を１画素として細線化の対象図形とし、閾値以下の画素を０画素（背景）とする。
【００４８】
次に、（２）前記肺領域抽出手段２０で抽出した肺領域の画像をかけて、不要な領域を消去する（画素値を０とする）。
【００４９】
その後、（３）画像全体にユークリッド二乗距離変換処理を行なう。これによって、図形上の全画素に、それぞれの画素から背景までの最大二乗距離が入る。ここで、二乗距離を用いるのは整数演算で済ませるためであり、その後の処理において距離の大小関係のみを把握すれば良いからである。
【００５０】
次に、（４）距離値の小さい画素から消去判定を行なっていく。この消去判定ではトポロジーが変化しないものを消去可能としている。まず、最小距離値（最も外側）の画素の中で消去可能性を調べ、そこで消去不可能な画素を一時保存画素とする。一時保存は、現時点ではトポロジー保存の観点から消去不可能であるが、消去処理が進むに連れて消去可能となる可能性があるために用いる。消去可能な画素について、端点（隣接する２６近傍の正数画素数が１）の場合は、永久保存（１画素として確定）する。その他の消去可能な画素については、その隣接する２６近傍の正数画素数によって３個毎にまとめてグループに分け、正数画素数の小さい画像画素（即ち、より孤立して突出した部分）から消去判定を行なう。これにより、境界面上の無意味な凹凸によって不要な枝が生じるのを抑えることができる。
【００５１】
次いで、（５）ここで消去可能な画素を一時保存し、消去可能なものの中で端点を永久保存し、その他を消去する。
【００５２】
次に、（６）消去した画素の６近傍画素（画素（ボクセル）同士が面で接する画素）について、ユークリッド二乗距離値を持つ画素を消去判定画素のグループに含める。
【００５３】
以後、上記（４）以降の操作を消去可能な画素が無くなるまで繰り返す。その結果得られた正確数値を全て１にすることにより、細線化画像が得られる。
【００５４】
この細線化法においてユークリッド距離を用いる理由は、より一般的且つ単純な方法と考えられる３つの直交軸方向に関する境界画素（６近傍画素距離）に対する消去判定に比べ、図形の回転によって細線化結果に及ぼす影響が小さくなるためである。肺の血管（気管支）のような、様々な方向を向いた図形に対する細線化においては、ユークリッド距離を用いた方が、より正確に中心線を得ることができると考えられる。
【００５５】
細線化の例を図５に示す。ここでは、５１２（ピクセル）×５１２（ピクセル）×１８１（スライス）の３次元ＣＴ画像について、細線化を行なった。２値化の閾値ＣＴ値は−４００とした。左は、ＣＴ画像中のある３０スライス分のＭＩＰ（Maximum Intensity Projection）像で、右は同じ部分の細線化画像である。右図のコントラストは、スライス深さに依存している。
【００５６】
前記分岐点抽出手段２４の分岐点抽出に際して、肺の動きを定量化するために追跡すべき特徴点として、本実施形態では、肺内の血管（気管支）の分岐点を用いる。その座標を特定するために、太さ１となっている細線化画像１４中の線図形の分岐点１６を取り出し、これを特徴点とする。分岐点抽出では、細線化画像中の各画素について、「２６近傍の１画素数が３」、且つ、「その画素を中心とした３×３×３領域において中心画素値を０としたときの成分指数が３」である画素を分岐点とした。ここで、成分指数とは、３次元画像中の２６近傍で繋がった１画素の集合の数である。
【００５７】
画像によっては、抽出した分岐点数が非常に多くなる場合がある。この場合、後述の特徴点追跡における計算時間の短縮や計算機メモリの制限、追跡エラーの除去のために、ある近傍領域に複数の分岐点が存在するものは削除する。
【００５８】
上の細線化画像に対して分岐点抽出を行なった結果を図６に示す。ここでは、分岐点のうちで各軸方向に±２画素近傍に複数存在するものを削除した。画像全体では２４００点程度の分岐点が得られており、これは、後で行なう特徴点追跡処理において十分な数である。
【００５９】
前記特徴点追跡手段２６の特徴点追跡では、ある微小時間間隔で取得された２つの３次元ＣＴ画像（ｆ０、ｆ１）間の物体の変形を定量化するために、それぞれの画像から得られた特徴点を対応点候補としてポイントパターンマッチングを行ない、画像ｆ０に対する画像ｆ１の対応点を求める。その結果得られる、対応点同士を結んだベクトルが、変形を表わす変位ベクトル１８となる。ここでは、確率的弛緩法によるポイントパターンマッチングを行なう。これは、Barnardらによって報告（Stephen T.Barnard and William B.Thompson,“Disparity Analysis of Image,”IEEE TRANSACTION OF PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VoL.PAMI-2, No.4, July 1980, PP333-340）された方法を、３次元画像に適用、修正したものである。
【００６０】
具体的には、図７に示す如く、画像ｆ０、ｆ１から抽出された対応点候補を、それぞれａｉ（ｉ＝１，２，…，ｎa）、ｂｋ（ｋ＝１，２，…，ｎb）とする。ａｉに対して、ａｉを中心として（２ｒ＋１）3画素の大きさの立方体内に存在する画像ｆ１上の対応点候補ｂｋを対応点対候補とする。変位ベクトル１８は、ａｉを始点とし、ｂｋを終点とするベクトルとして定義されるが、この段階で、ａｉは０個以上のｂｋと関係付けられるので、ａｉに対して一般的に複数の変位ベクトル候補が存在する。変位ベクトル候補の数が０の場合は、探索範囲内にｂｋが存在しない場合である。
【００６１】
ａｉに対して複数の変位ベクトル候補が存在するということは、それだけ曖昧さが存在していることを意味する。この曖昧さを、ａｉに対して得られている変位ベクトル候補とａｉの近傍のａｊ（ｊ≠ｉ）に対して得られている変位ベクトル候補が互いに矛盾しないように候補を絞っていく処理を弛緩法で行なう。ここで、変位ベクトル候補をラベルλｋで表わす。
【００６２】
【数１】
【００６３】
各ラベルの妥当性を、ラベル確率として定義する。結果として、高いラベル確率を持つ変位ベクトルが、求める動きとして判断される。
【００６４】
ａｉがＬ個の対応点対候補ｂｋｌ（ｌ＝１，２，…，Ｌ）を持つとき、ａｉはＬ＋１個のラベルを持つ。このラベル集合をΛｉとする。
【００６５】
【数２】
【００６６】
ここで、λｍは、対応点対候補が存在しないラベルを表わす。
【００６７】
（ｘｉ，ｙｉ，ｚｉ）及び（ｘｋｌ，ｙｋｌ，ｚｋｌ）を中心とする1辺（２ｄ＋１）画素の立方体を考え、ａｉとｂｋの類似度を次式で定義する。
【００６８】
【数３】
【００６９】
但し、Ｃは適当に設定される定数である。Ｓｉ（λｋｌ）の値域は［０，１］であるから、これを用いてａｉのλｋｌに対する初期ラベル確率Ｐｉ(0)（λｋｌ）を、次のように定義する。
【００７０】
【数４】
【００７１】
ａｉのラベルλｋｌのラベル確率Ｐｉ（λｋｌ）は、λｋｌに対して高いラベル確率を持つａｊがａｉの近傍に存在するならば、即ち、近傍の特徴点が同じ変位ベクトルを持つ確率が大きいならば、ａｉがλｋｌを持つことに関する矛盾は少ないと考えられる。従って、その場合、より高い確率を持つように更新する。そこで、まずａｉのラベルλｋｌの局所的無矛盾度Ｑｉ（λｋｌ）を次式で定義する。
【００７２】
【数５】
【００７３】
ｌ’については、ａｉのラベルλｋｌ＝（Δｘｋｌ，Δｙｋｌ，Δｚｋｌ）に対して、ある程度の位置ずれＴを許容して、次式を満たすａｊのラベルλｋｌ’であるｌ’について考える。
【００７４】
【数６】
【００７５】
ここでは、Ｔ＝１としている。従って、（５）式は、ａｉの近傍にあるａｊのうちで、ほぼ同じ方向のベクトルを持つラベル確率の和を表わしている。ここで近傍点ａｊとしては、次式を満たすように選択される。
【００７６】
【数７】
【００７７】
ここで、Ｒは、近傍点ａｊを探索する範囲を表わす。
【００７８】
ラベル確率の更新は、Ｑｉ（λｋｌ）を用いて、次のようになされる。
【００７９】
【数８】
ここで、
【数９】
であり、Ａ、Ｂは適当に設定される正の定数である。
【００８０】
（８）式、（９）式、（１０）式によるラベル確率の更新を、各ａｉ（ｉ＝１，２，…，ｎａ）に対して並列的（独立）に行なう処理を、１回の反復処理とし、ラベル確率の変化が無くなる（収束する）まで反復する。反復処理は１０回程度のオーダーである。
【実施例】
【００８１】
特徴点追跡の試験的適用として、先に細線化処理で用いたＣＴ画像と、それを少し計算（ある一点を基準に５％膨張させ、ＣＴ値の揺らぎとして、平均値０、標準偏差３０のランダムな値を各画素に加算）により変形させた画像を用いて、それぞれ肺領域抽出、細線化、分岐点抽出を行ない、得られた結果を用いて、特徴点追跡を行なった。ここで、確率的弛緩法によるポイントパターンマッチングにおける各パラメータの値は、ｒ＝１５、ｄ＝２、Ｃ＝０．０００１、Ｒ＝１５、Ａ＝０．３、Ｂ＝３．０である。
【００８２】
初期ＣＴ画像上の各対応候補点における計算されたラベル確率のうちで、最大のものを採用した。その結果、対応点対候補を持つ点数２４８６について、変位ベクトルを持つ点は１３４３点であり、残りの１１４３点は対応点無しであった。
【００８３】
得られた変位ベクトルを図８に示す。ベクトル数としては、ボクセル全体の変位の定量化のために行なわれるグリッド補間には十分な数が得られている。
【００８４】
計算の結果得られた変位ベクトルの妥当性を評価するために、そのベクトル成分の理想値からのズレを求めた。計算された変位ベクトルの成分はほぼ理想値に近いが、理想値からズレているものがある。
【００８５】
ズレの原因としては、（１）細線化に行なう際に生じる中心線の座標の揺らぎがある。図形形状が複雑であったり、真の中心線がボクセルとボクセルの間に来てしまったり、又、消去判定を行なうボクセルの順序等のアルゴリズムの影響等も含めて、１ボクセル程度の揺らぎが出てしまう。即ち、抽出された特徴点の座標に１程度の揺らぎが存在する。変位ベクトル成分の理想値についても、ボクセルの座標から計算される理想値なので、特徴点の座標を特定する段階で、既に揺らぎが含まれている。
【００８６】
別の原因としては、（２）特徴点追跡における誤認識が考えられる。ポイントパターンマッチングにおいて、間違った点を対応点としてしまうことが有り得る。只、ズレの分布を見る限り、ほとんどのベクトルのズレが±１程度に収まっているので、ほぼ正確に計算されていると思われる。ｘ、ｙ、ｚ方向のボクセルサイズは、今回のデータの場合、それぞれ０．６２５、０．６２５、１ｍｍであるので、本手法における動きの定量化の揺らぎも、それと同程度と結論付けられる。
【００８７】
変位ベクトルは、空間中に分散した限られた数のボクセルに対する動きを示しているが、局所的（変位ベクトル同士の間隔より小さい領域）に変位ベクトルと異なった変位は、ほとんど起こらないと考えられるので、グリッド補間によって、空間中の全てのボクセルに対して変位量が求められる。
【００８８】
動きの定量化を用いた応用例として、動きを考慮した線量分布評価等が挙げられる。その場合、本手法を用いれば、ある時刻のＣＴ画像が設定したターゲットボリュームについて、変形後のものを自動的に計算することができる。その例を図９に示す。ここでは、簡単のために直方体の領域を与えた。実線Ａで示す直方体は、設定したボリューム、破線Ｂで示す直方体は、変位ベクトルからグリッド補間することにより計算されたボリュームであり、ほぼ正確に変形した直方体が得られる。
【００８９】
なお、前記説明においては、肺が対象臓器とされていたが、対象臓器はこれに限定されず、血管造影等の手法を用いて特徴点の抽出方法を工夫することにより、肺以外の肝臓等の臓器にも同様に適用できる。
【００９０】
次に、図１０を参照して、本発明を利用して臓器位置を予測し、放射線を照射するようにした本発明の第２実施形態を説明する。
【００９１】
本実施形態においても、図１０の左方に示すように、ＣＴガントリ３２内のベッド３４上の患者３０に対して、Ｘ線管球３６からＸ線３７を照射し、Ｘ線検出器３８により得られた画像に対して、第１実施形態と同様に、４次元ＣＴ（４Ｄ−ＣＴ）４０を用いて、臓器の動態画像（４次元ＣＴ画像）を得る。このとき、例えば患者３０の体表面に固定した呼吸センサ５０から、運動情報取得手段５２で呼吸波形を同時に計測することにより、両者の時間的な相関を得ておく。
【００９２】
そして、処理部６０で、４次元ＣＴ画像から第１実施形態の手法による動態定量化を行い、各位置の動態情報を得る。これにより、上記で計測した呼吸波形と各位置の運動との相関を得て、任意の呼吸位相での臓器の位置を予測することができる。
【００９３】
次に、上で得られた相関データを用いて、動きを考慮した４次元の治療計画（腫瘍や正常組織への線量評価・照射方法等の条件の最適化）を行う。
【００９４】
具体的には、図１０の右方に示したように、放射線治療装置７０のベッド７２上の患者３０の体表面に固定した呼吸センサ５０から、運動情報取得手段５２により呼吸波形をリアルタイムで計測する。
【００９５】
そして、放射線を照射すべき呼吸位相を指定し、前記相関データから臓器位置をリアルタイムで予測して、治療計画に合わせた必要なタイミングで、照射制御装置７４を介して放射線照射装置７６から治療用放射線を照射する。
【００９６】
次に、図１１（手順を示す流れ図）及び図１２（概念図）を参照して、本発明を利用して臓器異常を検出するようにした本発明の第３実施形態を説明する。
【００９７】
本実施形態においては、図１１に手順を示す如く、第１実施形態の手法により得られた各画素（ボクセル）ごとの変位ベクトルを、図１２に示す如く、そのベクトルの成分によって領域を区分する。即ち、地図上の各地点の標高によって等高線を得るように、ベクトルマップの成分の大きさによって領域を分けることにより、異なる変位の領域を抽出することができる。これによって、異常領域を検出することができる。
【図面の簡単な説明】
【００９８】
【図１】本発明による動態定量化を示す概念図
【図２】本発明の第１実施形態を実施するための装置の基本的な構成を示すブロック図
【図３】第１実施形態における肺領域の抽出手順を示す図
【図４】同じく抽出結果を示す図
【図５】同じく細線化の手順を示す図
【図６】同じく分岐点抽出結果を示す図
【図７】同じく特徴点追跡の手順を示す図
【図８】同じく計算された変位ベクトルを示す図
【図９】同じく変形後のボリュームの計算例を示す図
【図１０】本発明の第２実施形態の手順を示す流れ図
【図１１】本発明の第３実施形態の手順を示す流れ図
【図１２】本発明の第３実施形態を示す概念図
【符号の説明】
【００９９】
１０…肺
１２…血管や気管
１４…細線化画像
１６…分岐点
１８…変位ベクトル
２０…肺領域抽出手段
２２…細線化手段
２４…分岐点抽出手段
２６…特徴点追跡手段
３０…患者
３２…ＣＴガントリ
３６…Ｘ線管球
３８…Ｘ線検出器
４０…４次元ＣＴ（４Ｄ−ＣＴ）
５０…呼吸センサ
５２…運動情報取得手段
６０…処理部
７０…放射線治療装置
７４…照射制御装置
７６…放射線照射装置
【技術分野】
【０００１】
本発明は、放射線診断・治療計画及び治療において、呼吸等に伴って変動・変形する肺等の臓器に対する臓器内各部の対応付け（位置合わせ）や、腫瘍等の関心体積の変形や位置を、呼吸波形の位相との相関として把握するために、高速３次元ＣＴ（いわゆる４次元ＣＴ）を用いて、ある時間間隔で取得された各時刻の３次元デジタル画像を用いて、各画素位置の動態を定量化するための臓器動態の定量化方法、装置、これを利用した臓器位置の予測方法、装置、放射線照射方法、装置、及び、臓器異常検出装置に関する。
【背景技術】
【０００２】
放射線治療の照射対象となる腫瘍の位置は、呼吸や日々の臓器の変動・変形により動く。放射線治療の際には、決められた位置に照射目標となる臓器（ターゲット）が来るように、一般にＸ線透過像（２次元画像）を用いて位置合わせを行なっている。
【０００３】
呼吸に伴って動く臓器に対しては、例えば呼吸同期法を用いて、呼吸波形のある領域のみ照射を行なうことにより、呼吸移動により生じる余分な照射を減らしている。
【０００４】
しかし、この呼吸波形は、体表面上に設定した、ある点の１次元的な位置を示しているものであり、又、その照射中にも臓器は３次元的に動いている。従って、より高精度な照射を行なうためには、臓器の各位置の動態（運動の方向及び速度）を把握しなければならないが、これを２次元画像から解析することは困難である。
【０００５】
特に、粒子線治療において、肺のように腫瘍と正常組織の密度が大きく異なる場合、照射ビーム（以下単にビームと称する）に対して垂直な方向だけでなく、ビーム軸方向の臓器の変動がビームの飛程に大きく影響し、その結果、腫瘍への過小照射、正常組織への過大照射につながる危険性がある。そのため、臓器の動態を３次元的に把握する技術が求められているが、現在臨床では、このような技術は確立されていない。
【０００６】
近年、臓器の３次元的な位置を定量化する研究が、いくつか報告されている。それらの多くは、Ｘ線ＣＴ画像を用いて抽出された対象臓器の表面上に、動きを追うための目印となる参照点を設定し、２つのＣＴ画像上の参照点同士が一致するように、臓器を弾性体モデルとして、その変形を計算している。
【０００７】
しかしながら、この場合、特徴点を人為的に設定しなければならず、又、臓器内部の各位置の対応は必ずしも現実の動きを反映していない。又、具体的な動きの追跡のために金属マーカーが一般的に用いられているが、患者の体内へ埋め込むため、負担が大きく、部位や個数が制限され、臓器の動態を細かく捉えるのが困難である。
【０００８】
一方、非特許文献１には、有限要素法により肺の３次元的な変形を計算することが記載されている。
【０００９】
又、非特許文献２では、肺の３次元ＣＴ画像の画素値の勾配を比較することにより、内部の変形を定量化しており、腫瘍部分の関心体積の変形についても論じている。
【００１０】
更に、非特許文献３には、呼吸運動に伴う肺の変形での臓器内の対応付けのために、血管分岐点の動きを追跡することが記載されている。
【００１１】
【非特許文献１】Ｔ．Ｚhangら「Ｔechnical note：Ａ novel boundary condition using contact elements for finite element based deformable image registration」Ｍedical Ｐhysics ３１（２００４）２４１２−２４１５
【非特許文献２】Ｄ．Ｓarrutら「Ｎonrigid registration method to assess reproducibility of breath-holding with ＡＢＣ in lung cancer」Ｉnternational Ｊournal of Ｒadiation Ｏncology，Ｂiology，Ｐhysics，６１（２００５）５９４−６０７
【非特許文献３】四方秀則ら「肺の非剛体レジストレーションのための血管分岐点位置決めアルゴリズム」 電子情報通信学会論文誌Ｄ−II Ｖol.Ｊ８４-Ｄ-II Ｎo.０７ pp．１−１１
【発明の開示】
【発明が解決しようとする課題】
【００１２】
しかしながら、非特許文献１に記載の技術は、臓器の外形と内部の変形が必ずしも一致していないという問題点を有する。
【００１３】
又、非特許文献２に記載の技術は、内部の解剖学的特徴を追っているわけではなく、関心体積変形の正確さについては問題がある。
【００１４】
又、非特許文献３に記載の技術は、単に、血管分岐点の位置決めを行なっているだけである。
【００１５】
本発明は、前記従来の問題点を解消するべくなされたもので、臓器の動態を正確に定量化することを課題とする。
【課題を解決するための手段】
【００１６】
本発明は、微少に変動・変形した臓器の複数のＣＴ画像を入力データとし、臓器内に分布した血管や気管を抽出して細線化し、細線化した画像を用いて分岐点や連結点の座標を抽出し、これらを特徴点として複数のＣＴ画像間の個々の点の移動を３次元空間中で追跡することにより、臓器の運動を計測するようにして、前記課題を解決したものである。
【００１７】
又、前記特徴点の追跡を、複数のＣＴ画像及び特徴点座標を用いたポイントパターンマッチングにより行なうようにしたものである。
【００１８】
又、複数の特徴点の運動を補間することにより、任意点の運動が分かるようにしたものである。
【００１９】
本発明は、又、前記の方法により得られた動的な情報と、臓器の周期的な運動の情報との相関を取ることにより、臓器の位置を予測することを特徴とする臓器位置の予測方法を提供するものである。
【００２０】
本発明は、又、前記の方法により予測された臓器の位置に、放射線を照射することを特徴とする放射線照射方法を提供するものである。
【００２１】
本発明は、又、微少に変動・変形した臓器の複数のＣＴ画像を入力する手段と、臓器内に分布した血管や気管を抽出して細線化する手段と、細線化した画像を用いて分岐点や連結点の座標を抽出する手段と、これらを特徴点として複数のＣＴ画像間の個々の点の移動を３次元空間中で追跡することにより、臓器の運動を計測する手段と、を備えることにより、前記課題を解決したものである。
【００２２】
又、前記特徴点の追跡を、複数のＣＴ画像及び特徴点座標を用いたポイントパターンマッチングにより行う手段を備えたものである。
【００２３】
又、複数の特徴点の運動を補間して、任意点の運動を求める手段を備えたものである。
【００２４】
本発明は、又、前記の装置により得られた動的な情報と、臓器の周期的な運動の情報との相関を取ることにより、臓器の位置を予測する手段を備えたことを特徴とする臓器位置予測装置を提供するものである。
【００２５】
本発明は、又、前記の装置により予測された臓器の位置に、放射線を照射する手段を備えたことを特徴とする放射線照射装置を提供するものである。
【００２６】
本発明は、又、前記の装置により得られた動的な情報に基づいて、臓器の局所的に異常な変形を捉える手段を備えたことを特徴とする臓器異常検出装置を提供するものである。
【発明の効果】
【００２７】
本発明においては、肺の場合、気管支や血管のように連結する解剖学的特徴に着目し、それらのトポロジー（分岐や連結）は、臓器の運動や変形でも保持されていることを利用して、臓器の動態を定量的に計測する。即ち、動態の前後で対象のトポロジーが保存されることから、それらの特徴点座標を抽出し、個々の点の移動を３次元空間中で追跡することで、局所的にも高い精度を有する臓器の運動（方向、速度）を計測する。
【００２８】
例えば、気管支に沿って存在する血管のＣＴ値が、肺野のＣＴ値に比較して大きいことを利用し、図１（Ａ）に示す如く、ＣＴ画像から血管や気管１２を抽出する。
【００２９】
次いで、抽出された血管や気管部を、細線化アルゴリズムを利用して、太さ１画素の細線化画像１４に変換し、血管や気管の特徴点（図１では分岐点１６）座標を特定し、これを動態追跡の特徴点とする。これにより、従来動きの追跡に用いられる金属マーカーでは不可能な、肺内部領域全体にわたって、任意の点での局所的な動態追跡が可能となる。
【００３０】
次いで、ある時間間隔で取得された２つのＣＴ画像に対して、これらの特徴点の変位を追跡するために、２つのＣＴ画像及び特徴点座標を用いて、ポイントパターンマッチングを行なう。例えば確率的弛緩法を用いてこの操作を行なうことにより、特徴点の変位が、その近傍の特徴点の変位と互いに矛盾が無いように、各特徴点の変位量（方向と長さ）を求めることができる。即ち、抽出した特徴点の追跡から、図１（Ｂ）に示す如く、各点の運動が分かる。
【００３１】
更に、それらの点から３次元空間での補間を行なうことで、任意点の運動が分かる。従って、臓器（図では肺１０）内の全ての画素に対してその変位量が求められる。
【００３２】
以上の操作を連続的に取得されたＣＴ画像に対して逐次行なうことにより、連続的に動態定量化を行なうことができる。従って、動く対象に対する４次元的な線量計算が可能となり、呼吸動態臓器に対する、より高精度な放射線治療を行うことが可能となる。
【００３３】
放射線治療のための治療計画（線量計算、評価による放射線の照射方向・回数等の照射方法を決定する操作）等では、腫瘍の部分を関心体積として設定し、その関心体積やその周りの正常組織への放射線量計算を行なうが、従来は静的な３次元ＣＴ画像に対して行われていたものを、本発明を用いれば、４次元ＣＴを用いて取得されたＣＴ画像を用いて、動的な対象に対してこれを行なうことが可能となり、関心体積の変形が自動的に計算でき、各時刻で取得されたＣＴ画像に対して設定する必要が無くなる。これにより、４次元的な治療計画を、より現実的なものとすることができる。
【００３４】
又、近年開発され、医療現場へ普及しつつある４次元ＣＴから逐次得られる時系列の３次元デジタル画像を用いて、本発明の操作を逐次的に行なうことにより、臓器の３次元的な動きを、ある時間領域にわたって得ることができる。
【００３５】
従って、本発明により得られた動的な情報を、呼吸波形のような参照情報と相関をとることにより、任意の時点での臓器の動き（３次元の方向とその速度）を予測することができる。
【００３６】
又、本発明は、動きの特徴点を得る方法と、その特徴点の動きを追跡する方法に分かれているため、例えば血管造影法を用いて肝臓のＸ線ＣＴ画像を取得する等のように、特徴点を得る方法を工夫すれば、様々な部位に応用することが可能である。
【００３７】
このようにして、臓器の３次元的な動きを把握することにより、従来位置の不確定さのために設定せざるを得なかった照射マージンを、より縮小して、正常組織への不必要な放射線被曝を抑えることができ、晩期障害等の放射線障害を減らすことができる。
【００３８】
臓器の解剖学的な変性は、臓器の運動に歪みや局所的に異常な変形を引き起こすので、本発明を用いて、正常時あるいは正常部分との動態を定量的に比較して、局所的な運動の歪みや異常性を検出することで、医学的診断に応用できる。例えば、ＣＴ値では判別し難い腫瘍の存在や肺正常組織への放射線障害である繊維化の兆候を知る手懸りとなる。
【発明を実施するための最良の形態】
【００３９】
以下図面を参照して、本発明の実施形態を詳細に説明する。
【００４０】
本発明の第１実施形態を実施するための装置は、例えばコンピュータを用いて、図２に示す如く構成され、肺領域抽出手段２０と、細線化手段２２と、分岐点抽出手段２４と、特徴点追跡手段２６と、を含んでいる。
【００４１】
前記細線化手段２２では、あるＣＴ値を閾値として２値化を行ない、それによって得られた１画素値の図形を削っていく。しかし、肺の周りの肉や骨の部分が２値化によって残ってしまうと、それも細線化処理の対象となり、その部分は肺内の気管支や血管に比べて体積（画素数）がかなり大きいため、計算時間が大幅に長くなる。
【００４２】
そこで、計算時間を短くし、又、必要な肺領域のみを計算させるために、前記肺領域抽出手段２０で、まず肺領域を抽出する。
【００４３】
具体的には、図３に示したような方法で行なうことができる。処理は、図３（Ａ）に例示する、それぞれの２次元ＣＴスライス画像毎に行なう。
【００４４】
まず、あるスライス画像について、ある設定した矩形領域でスムージング処理を行う。これを全ての２次元領域について行う。この処理によって、図３（Ｂ）に示す如く、肺内のＣＴ値が平均的に小さい領域と、肺の周りのＣＴ値の高い領域に分けられる。次に、この画像について、図３（Ｃ）に示す如く、ある適当なしきい値により二値化を行い、ＣＴ値の小さい領域のみ抽出する。その後、左右の肺領域を指定することによって体外部など不要領域を削除して、図３（Ｄ）に示す如く、肺領域のみを抽出する。以上の操作をすべてのスライスについて行うことにより、肺領域の３次元画像を得る。ここで得られた領域のみを後の細線化処理の対象とする。
【００４５】
肺領域の抽出結果の例を図４に示す。ここでは、平均化の矩形領域サイズを１５、ＣＴしきい値を−５００とした。
【００４６】
前記細線化手段２２は、厚さ・太さのある図形を、トポロジーを変化させずに太さ１の中心線に対応する線図形に変換する操作（細線化）を行う。この細線化によって、肺内の血管（気管支）の分岐点の座標を画素単位で特定することができ、その結果得られた分岐点を動きの特徴点として追跡することにより、肺の動きを定量化することができる。トポロジーを変えないということは、連続した図形の数や、図形の穴や空洞の数といった図形の性質が変化しないことを意味する。ここでは、ユークリッド距離変換を伴う逐次型アルゴリズムによる細線化法を用いている。
【００４７】
具体的には、まず（１）３次元ＣＴ画像について、適当なＣＴ値を閾値として２値化を行なう。ここでは、閾値より大きい画素を１画素として細線化の対象図形とし、閾値以下の画素を０画素（背景）とする。
【００４８】
次に、（２）前記肺領域抽出手段２０で抽出した肺領域の画像をかけて、不要な領域を消去する（画素値を０とする）。
【００４９】
その後、（３）画像全体にユークリッド二乗距離変換処理を行なう。これによって、図形上の全画素に、それぞれの画素から背景までの最大二乗距離が入る。ここで、二乗距離を用いるのは整数演算で済ませるためであり、その後の処理において距離の大小関係のみを把握すれば良いからである。
【００５０】
次に、（４）距離値の小さい画素から消去判定を行なっていく。この消去判定ではトポロジーが変化しないものを消去可能としている。まず、最小距離値（最も外側）の画素の中で消去可能性を調べ、そこで消去不可能な画素を一時保存画素とする。一時保存は、現時点ではトポロジー保存の観点から消去不可能であるが、消去処理が進むに連れて消去可能となる可能性があるために用いる。消去可能な画素について、端点（隣接する２６近傍の正数画素数が１）の場合は、永久保存（１画素として確定）する。その他の消去可能な画素については、その隣接する２６近傍の正数画素数によって３個毎にまとめてグループに分け、正数画素数の小さい画像画素（即ち、より孤立して突出した部分）から消去判定を行なう。これにより、境界面上の無意味な凹凸によって不要な枝が生じるのを抑えることができる。
【００５１】
次いで、（５）ここで消去可能な画素を一時保存し、消去可能なものの中で端点を永久保存し、その他を消去する。
【００５２】
次に、（６）消去した画素の６近傍画素（画素（ボクセル）同士が面で接する画素）について、ユークリッド二乗距離値を持つ画素を消去判定画素のグループに含める。
【００５３】
以後、上記（４）以降の操作を消去可能な画素が無くなるまで繰り返す。その結果得られた正確数値を全て１にすることにより、細線化画像が得られる。
【００５４】
この細線化法においてユークリッド距離を用いる理由は、より一般的且つ単純な方法と考えられる３つの直交軸方向に関する境界画素（６近傍画素距離）に対する消去判定に比べ、図形の回転によって細線化結果に及ぼす影響が小さくなるためである。肺の血管（気管支）のような、様々な方向を向いた図形に対する細線化においては、ユークリッド距離を用いた方が、より正確に中心線を得ることができると考えられる。
【００５５】
細線化の例を図５に示す。ここでは、５１２（ピクセル）×５１２（ピクセル）×１８１（スライス）の３次元ＣＴ画像について、細線化を行なった。２値化の閾値ＣＴ値は−４００とした。左は、ＣＴ画像中のある３０スライス分のＭＩＰ（Maximum Intensity Projection）像で、右は同じ部分の細線化画像である。右図のコントラストは、スライス深さに依存している。
【００５６】
前記分岐点抽出手段２４の分岐点抽出に際して、肺の動きを定量化するために追跡すべき特徴点として、本実施形態では、肺内の血管（気管支）の分岐点を用いる。その座標を特定するために、太さ１となっている細線化画像１４中の線図形の分岐点１６を取り出し、これを特徴点とする。分岐点抽出では、細線化画像中の各画素について、「２６近傍の１画素数が３」、且つ、「その画素を中心とした３×３×３領域において中心画素値を０としたときの成分指数が３」である画素を分岐点とした。ここで、成分指数とは、３次元画像中の２６近傍で繋がった１画素の集合の数である。
【００５７】
画像によっては、抽出した分岐点数が非常に多くなる場合がある。この場合、後述の特徴点追跡における計算時間の短縮や計算機メモリの制限、追跡エラーの除去のために、ある近傍領域に複数の分岐点が存在するものは削除する。
【００５８】
上の細線化画像に対して分岐点抽出を行なった結果を図６に示す。ここでは、分岐点のうちで各軸方向に±２画素近傍に複数存在するものを削除した。画像全体では２４００点程度の分岐点が得られており、これは、後で行なう特徴点追跡処理において十分な数である。
【００５９】
前記特徴点追跡手段２６の特徴点追跡では、ある微小時間間隔で取得された２つの３次元ＣＴ画像（ｆ０、ｆ１）間の物体の変形を定量化するために、それぞれの画像から得られた特徴点を対応点候補としてポイントパターンマッチングを行ない、画像ｆ０に対する画像ｆ１の対応点を求める。その結果得られる、対応点同士を結んだベクトルが、変形を表わす変位ベクトル１８となる。ここでは、確率的弛緩法によるポイントパターンマッチングを行なう。これは、Barnardらによって報告（Stephen T.Barnard and William B.Thompson,“Disparity Analysis of Image,”IEEE TRANSACTION OF PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VoL.PAMI-2, No.4, July 1980, PP333-340）された方法を、３次元画像に適用、修正したものである。
【００６０】
具体的には、図７に示す如く、画像ｆ０、ｆ１から抽出された対応点候補を、それぞれａｉ（ｉ＝１，２，…，ｎa）、ｂｋ（ｋ＝１，２，…，ｎb）とする。ａｉに対して、ａｉを中心として（２ｒ＋１）3画素の大きさの立方体内に存在する画像ｆ１上の対応点候補ｂｋを対応点対候補とする。変位ベクトル１８は、ａｉを始点とし、ｂｋを終点とするベクトルとして定義されるが、この段階で、ａｉは０個以上のｂｋと関係付けられるので、ａｉに対して一般的に複数の変位ベクトル候補が存在する。変位ベクトル候補の数が０の場合は、探索範囲内にｂｋが存在しない場合である。
【００６１】
ａｉに対して複数の変位ベクトル候補が存在するということは、それだけ曖昧さが存在していることを意味する。この曖昧さを、ａｉに対して得られている変位ベクトル候補とａｉの近傍のａｊ（ｊ≠ｉ）に対して得られている変位ベクトル候補が互いに矛盾しないように候補を絞っていく処理を弛緩法で行なう。ここで、変位ベクトル候補をラベルλｋで表わす。
【００６２】
【数１】
【００６３】
各ラベルの妥当性を、ラベル確率として定義する。結果として、高いラベル確率を持つ変位ベクトルが、求める動きとして判断される。
【００６４】
ａｉがＬ個の対応点対候補ｂｋｌ（ｌ＝１，２，…，Ｌ）を持つとき、ａｉはＬ＋１個のラベルを持つ。このラベル集合をΛｉとする。
【００６５】
【数２】
【００６６】
ここで、λｍは、対応点対候補が存在しないラベルを表わす。
【００６７】
（ｘｉ，ｙｉ，ｚｉ）及び（ｘｋｌ，ｙｋｌ，ｚｋｌ）を中心とする1辺（２ｄ＋１）画素の立方体を考え、ａｉとｂｋの類似度を次式で定義する。
【００６８】
【数３】
【００６９】
但し、Ｃは適当に設定される定数である。Ｓｉ（λｋｌ）の値域は［０，１］であるから、これを用いてａｉのλｋｌに対する初期ラベル確率Ｐｉ(0)（λｋｌ）を、次のように定義する。
【００７０】
【数４】
【００７１】
ａｉのラベルλｋｌのラベル確率Ｐｉ（λｋｌ）は、λｋｌに対して高いラベル確率を持つａｊがａｉの近傍に存在するならば、即ち、近傍の特徴点が同じ変位ベクトルを持つ確率が大きいならば、ａｉがλｋｌを持つことに関する矛盾は少ないと考えられる。従って、その場合、より高い確率を持つように更新する。そこで、まずａｉのラベルλｋｌの局所的無矛盾度Ｑｉ（λｋｌ）を次式で定義する。
【００７２】
【数５】
【００７３】
ｌ’については、ａｉのラベルλｋｌ＝（Δｘｋｌ，Δｙｋｌ，Δｚｋｌ）に対して、ある程度の位置ずれＴを許容して、次式を満たすａｊのラベルλｋｌ’であるｌ’について考える。
【００７４】
【数６】
【００７５】
ここでは、Ｔ＝１としている。従って、（５）式は、ａｉの近傍にあるａｊのうちで、ほぼ同じ方向のベクトルを持つラベル確率の和を表わしている。ここで近傍点ａｊとしては、次式を満たすように選択される。
【００７６】
【数７】
【００７７】
ここで、Ｒは、近傍点ａｊを探索する範囲を表わす。
【００７８】
ラベル確率の更新は、Ｑｉ（λｋｌ）を用いて、次のようになされる。
【００７９】
【数８】
ここで、
【数９】
であり、Ａ、Ｂは適当に設定される正の定数である。
【００８０】
（８）式、（９）式、（１０）式によるラベル確率の更新を、各ａｉ（ｉ＝１，２，…，ｎａ）に対して並列的（独立）に行なう処理を、１回の反復処理とし、ラベル確率の変化が無くなる（収束する）まで反復する。反復処理は１０回程度のオーダーである。
【実施例】
【００８１】
特徴点追跡の試験的適用として、先に細線化処理で用いたＣＴ画像と、それを少し計算（ある一点を基準に５％膨張させ、ＣＴ値の揺らぎとして、平均値０、標準偏差３０のランダムな値を各画素に加算）により変形させた画像を用いて、それぞれ肺領域抽出、細線化、分岐点抽出を行ない、得られた結果を用いて、特徴点追跡を行なった。ここで、確率的弛緩法によるポイントパターンマッチングにおける各パラメータの値は、ｒ＝１５、ｄ＝２、Ｃ＝０．０００１、Ｒ＝１５、Ａ＝０．３、Ｂ＝３．０である。
【００８２】
初期ＣＴ画像上の各対応候補点における計算されたラベル確率のうちで、最大のものを採用した。その結果、対応点対候補を持つ点数２４８６について、変位ベクトルを持つ点は１３４３点であり、残りの１１４３点は対応点無しであった。
【００８３】
得られた変位ベクトルを図８に示す。ベクトル数としては、ボクセル全体の変位の定量化のために行なわれるグリッド補間には十分な数が得られている。
【００８４】
計算の結果得られた変位ベクトルの妥当性を評価するために、そのベクトル成分の理想値からのズレを求めた。計算された変位ベクトルの成分はほぼ理想値に近いが、理想値からズレているものがある。
【００８５】
ズレの原因としては、（１）細線化に行なう際に生じる中心線の座標の揺らぎがある。図形形状が複雑であったり、真の中心線がボクセルとボクセルの間に来てしまったり、又、消去判定を行なうボクセルの順序等のアルゴリズムの影響等も含めて、１ボクセル程度の揺らぎが出てしまう。即ち、抽出された特徴点の座標に１程度の揺らぎが存在する。変位ベクトル成分の理想値についても、ボクセルの座標から計算される理想値なので、特徴点の座標を特定する段階で、既に揺らぎが含まれている。
【００８６】
別の原因としては、（２）特徴点追跡における誤認識が考えられる。ポイントパターンマッチングにおいて、間違った点を対応点としてしまうことが有り得る。只、ズレの分布を見る限り、ほとんどのベクトルのズレが±１程度に収まっているので、ほぼ正確に計算されていると思われる。ｘ、ｙ、ｚ方向のボクセルサイズは、今回のデータの場合、それぞれ０．６２５、０．６２５、１ｍｍであるので、本手法における動きの定量化の揺らぎも、それと同程度と結論付けられる。
【００８７】
変位ベクトルは、空間中に分散した限られた数のボクセルに対する動きを示しているが、局所的（変位ベクトル同士の間隔より小さい領域）に変位ベクトルと異なった変位は、ほとんど起こらないと考えられるので、グリッド補間によって、空間中の全てのボクセルに対して変位量が求められる。
【００８８】
動きの定量化を用いた応用例として、動きを考慮した線量分布評価等が挙げられる。その場合、本手法を用いれば、ある時刻のＣＴ画像が設定したターゲットボリュームについて、変形後のものを自動的に計算することができる。その例を図９に示す。ここでは、簡単のために直方体の領域を与えた。実線Ａで示す直方体は、設定したボリューム、破線Ｂで示す直方体は、変位ベクトルからグリッド補間することにより計算されたボリュームであり、ほぼ正確に変形した直方体が得られる。
【００８９】
なお、前記説明においては、肺が対象臓器とされていたが、対象臓器はこれに限定されず、血管造影等の手法を用いて特徴点の抽出方法を工夫することにより、肺以外の肝臓等の臓器にも同様に適用できる。
【００９０】
次に、図１０を参照して、本発明を利用して臓器位置を予測し、放射線を照射するようにした本発明の第２実施形態を説明する。
【００９１】
本実施形態においても、図１０の左方に示すように、ＣＴガントリ３２内のベッド３４上の患者３０に対して、Ｘ線管球３６からＸ線３７を照射し、Ｘ線検出器３８により得られた画像に対して、第１実施形態と同様に、４次元ＣＴ（４Ｄ−ＣＴ）４０を用いて、臓器の動態画像（４次元ＣＴ画像）を得る。このとき、例えば患者３０の体表面に固定した呼吸センサ５０から、運動情報取得手段５２で呼吸波形を同時に計測することにより、両者の時間的な相関を得ておく。
【００９２】
そして、処理部６０で、４次元ＣＴ画像から第１実施形態の手法による動態定量化を行い、各位置の動態情報を得る。これにより、上記で計測した呼吸波形と各位置の運動との相関を得て、任意の呼吸位相での臓器の位置を予測することができる。
【００９３】
次に、上で得られた相関データを用いて、動きを考慮した４次元の治療計画（腫瘍や正常組織への線量評価・照射方法等の条件の最適化）を行う。
【００９４】
具体的には、図１０の右方に示したように、放射線治療装置７０のベッド７２上の患者３０の体表面に固定した呼吸センサ５０から、運動情報取得手段５２により呼吸波形をリアルタイムで計測する。
【００９５】
そして、放射線を照射すべき呼吸位相を指定し、前記相関データから臓器位置をリアルタイムで予測して、治療計画に合わせた必要なタイミングで、照射制御装置７４を介して放射線照射装置７６から治療用放射線を照射する。
【００９６】
次に、図１１（手順を示す流れ図）及び図１２（概念図）を参照して、本発明を利用して臓器異常を検出するようにした本発明の第３実施形態を説明する。
【００９７】
本実施形態においては、図１１に手順を示す如く、第１実施形態の手法により得られた各画素（ボクセル）ごとの変位ベクトルを、図１２に示す如く、そのベクトルの成分によって領域を区分する。即ち、地図上の各地点の標高によって等高線を得るように、ベクトルマップの成分の大きさによって領域を分けることにより、異なる変位の領域を抽出することができる。これによって、異常領域を検出することができる。
【図面の簡単な説明】
【００９８】
【図１】本発明による動態定量化を示す概念図
【図２】本発明の第１実施形態を実施するための装置の基本的な構成を示すブロック図
【図３】第１実施形態における肺領域の抽出手順を示す図
【図４】同じく抽出結果を示す図
【図５】同じく細線化の手順を示す図
【図６】同じく分岐点抽出結果を示す図
【図７】同じく特徴点追跡の手順を示す図
【図８】同じく計算された変位ベクトルを示す図
【図９】同じく変形後のボリュームの計算例を示す図
【図１０】本発明の第２実施形態の手順を示す流れ図
【図１１】本発明の第３実施形態の手順を示す流れ図
【図１２】本発明の第３実施形態を示す概念図
【符号の説明】
【００９９】
１０…肺
１２…血管や気管
１４…細線化画像
１６…分岐点
１８…変位ベクトル
２０…肺領域抽出手段
２２…細線化手段
２４…分岐点抽出手段
２６…特徴点追跡手段
３０…患者
３２…ＣＴガントリ
３６…Ｘ線管球
３８…Ｘ線検出器
４０…４次元ＣＴ（４Ｄ−ＣＴ）
５０…呼吸センサ
５２…運動情報取得手段
６０…処理部
７０…放射線治療装置
７４…照射制御装置
７６…放射線照射装置
【特許請求の範囲】
【請求項１】
微少に変動・変形した臓器の複数のＣＴ画像を入力データとし、
臓器内に分布した血管や気管を抽出して細線化し、
細線化した画像を用いて分岐点や連結点の座標を抽出し、
これらを特徴点として複数のＣＴ画像間の個々の点の移動を３次元空間中で追跡することにより、臓器の運動を計測することを特徴とする臓器動態の定量化方法。
【請求項２】
前記特徴点の追跡を、複数のＣＴ画像及び特徴点座標を用いたポイントパターンマッチングにより行うことを特徴とする請求項１に記載の臓器動態の定量化方法。
【請求項３】
複数の特徴点の運動を補間して、任意点の運動を求めることを特徴とする請求項１又は２に記載の臓器動態の定量化方法。
【請求項４】
請求項１乃至３のいずれかに記載の方法により得られた動的な情報と、臓器の周期的な運動の情報との相関を取ることにより、臓器の位置を予測することを特徴とする臓器位置の予測方法。
【請求項５】
請求項４に記載の方法により予測された臓器の位置に、放射線を照射することを特徴とする放射線照射方法。
【請求項６】
微少に変動・変形した臓器の複数のＣＴ画像を入力する手段と、
臓器内に分布した血管や気管を抽出して細線化する手段と、
細線化した画像を用いて分岐点や連結点の座標を抽出する手段と、
これらを特徴点として複数のＣＴ画像間の個々の点の移動を３次元空間中で追跡することにより、臓器の運動を計測する手段と、
を備えたことを特徴とする臓器動態の定量化装置。
【請求項７】
前記特徴点の追跡を、複数のＣＴ画像及び特徴点座標を用いたポイントパターンマッチングにより行う手段を備えたことを特徴とする請求項６に記載の臓器動態の定量化装置。
【請求項８】
複数の特徴点の運動を補間して、任意点の運動を求める手段を備えたことを特徴とする請求項６又は７に記載の臓器動態の定量化装置。
【請求項９】
請求項６乃至８のいずれかに記載の装置により得られた動的な情報と、臓器の周期的な運動の情報との相関を取ることにより、臓器の位置を予測する手段を備えたことを特徴とする臓器位置の予測装置。
【請求項１０】
請求項９に記載の装置により予測された臓器の位置に、放射線を照射する手段を備えたことを特徴とする放射線照射装置。
【請求項１１】
請求項６乃至８のいずれかに記載の装置により得られた動的な情報に基づいて、臓器の局所的に異常な変形を捉える手段を備えたことを特徴とする臓器異常検出装置。
【請求項１】
微少に変動・変形した臓器の複数のＣＴ画像を入力データとし、
臓器内に分布した血管や気管を抽出して細線化し、
細線化した画像を用いて分岐点や連結点の座標を抽出し、
これらを特徴点として複数のＣＴ画像間の個々の点の移動を３次元空間中で追跡することにより、臓器の運動を計測することを特徴とする臓器動態の定量化方法。
【請求項２】
前記特徴点の追跡を、複数のＣＴ画像及び特徴点座標を用いたポイントパターンマッチングにより行うことを特徴とする請求項１に記載の臓器動態の定量化方法。
【請求項３】
複数の特徴点の運動を補間して、任意点の運動を求めることを特徴とする請求項１又は２に記載の臓器動態の定量化方法。
【請求項４】
請求項１乃至３のいずれかに記載の方法により得られた動的な情報と、臓器の周期的な運動の情報との相関を取ることにより、臓器の位置を予測することを特徴とする臓器位置の予測方法。
【請求項５】
請求項４に記載の方法により予測された臓器の位置に、放射線を照射することを特徴とする放射線照射方法。
【請求項６】
微少に変動・変形した臓器の複数のＣＴ画像を入力する手段と、
臓器内に分布した血管や気管を抽出して細線化する手段と、
細線化した画像を用いて分岐点や連結点の座標を抽出する手段と、
これらを特徴点として複数のＣＴ画像間の個々の点の移動を３次元空間中で追跡することにより、臓器の運動を計測する手段と、
を備えたことを特徴とする臓器動態の定量化装置。
【請求項７】
前記特徴点の追跡を、複数のＣＴ画像及び特徴点座標を用いたポイントパターンマッチングにより行う手段を備えたことを特徴とする請求項６に記載の臓器動態の定量化装置。
【請求項８】
複数の特徴点の運動を補間して、任意点の運動を求める手段を備えたことを特徴とする請求項６又は７に記載の臓器動態の定量化装置。
【請求項９】
請求項６乃至８のいずれかに記載の装置により得られた動的な情報と、臓器の周期的な運動の情報との相関を取ることにより、臓器の位置を予測する手段を備えたことを特徴とする臓器位置の予測装置。
【請求項１０】
請求項９に記載の装置により予測された臓器の位置に、放射線を照射する手段を備えたことを特徴とする放射線照射装置。
【請求項１１】
請求項６乃至８のいずれかに記載の装置により得られた動的な情報に基づいて、臓器の局所的に異常な変形を捉える手段を備えたことを特徴とする臓器異常検出装置。
【図１】
【図２】
【図７】
【図１０】
【図１１】
【図１２】
【図３】
【図４】
【図５】
【図６】
【図８】
【図９】
【図２】
【図７】
【図１０】
【図１１】
【図１２】
【図３】
【図４】
【図５】
【図６】
【図８】
【図９】
【公開番号】特開２００６−３１４６４３（Ｐ２００６−３１４６４３Ａ）
【公開日】平成１８年１１月２４日（２００６．１１．２４）
【国際特許分類】
【出願番号】特願２００５−１４１８８０（Ｐ２００５−１４１８８０）
【出願日】平成１７年５月１３日（２００５．５．１３）
【出願人】（３０１０３２９４２）独立行政法人放射線医学総合研究所 (149)
【出願人】（５０４１４５３６４）国立大学法人群馬大学 (352)
【Ｆターム（参考）】
【公開日】平成１８年１１月２４日（２００６．１１．２４）
【国際特許分類】
【出願日】平成１７年５月１３日（２００５．５．１３）
【出願人】（３０１０３２９４２）独立行政法人放射線医学総合研究所 (149)
【出願人】（５０４１４５３６４）国立大学法人群馬大学 (352)
【Ｆターム（参考）】
[ Back to top ]