JP4945220B2

JP4945220B2 - 血管狭窄部位探索装置、プログラム、及び画像診断装置

Info

Publication number: JP4945220B2
Application number: JP2006313138A
Authority: JP
Inventors: 義強楊
Original assignee: Toshiba Corp; Toshiba Medical Systems Corp; Toshiba Medical Systems Engineering Co Ltd
Current assignee: Toshiba Corp; Canon Medical Systems Corp
Priority date: 2006-11-20
Filing date: 2006-11-20
Publication date: 2012-06-06
Anticipated expiration: 2026-11-20
Also published as: JP2008125736A

Description

本発明は、画像内の血管から狭窄部位を探索する技術に関する。

近年の医療診断では、血管の狭窄程度を定量的に把握することが重要な診断項目の一つとなっている。例えば、頸部の血管狭窄が著しいときには脳梗塞が起こるおそれがあることが知られている。また、冠状動脈ツリーの場合には、狭窄は心筋層の血流を著しく妨げるおそれがある。

そのため、従来、医師や技師等がデジタルアンギオ画像のような血管造影像から手作業により狭窄部位を探索する方法や、血管の輪郭を抽出して輪郭形状より狭窄部位を探索する方法（例えば、「特許文献１参照。」）により、狭窄部位が探索されていた。

血管の輪郭抽出方法としては、血管の正常部分と狭窄部分の夫々において血管走行方向に対して垂直な線分を操作者に引かせることにより、その垂直線分の濃度プロファイルを抽出し、この濃度プロファイルに基づき血管辺縁を決定し、狭窄率を算出する画像処理がある。この画像処理により得られた輪郭形状より血管上の各血管断面の径と血管断面径の平均値を求め、各点における血管断面径と平均値との比を求め、その比が所定値以下であれば狭窄部位としている。

このような従来の狭窄部位探索方法によっては、血管の輪郭を抽出する画像処理の精度が狭窄部位の探索精度を左右する。そのため、画像処理を短時間に抑えようとすると血管輪郭の抽出精度が落ち、狭窄部位を精度良く探索できないおそれがある。また、血管輪郭の抽出精度を上げようとすると、画像処理に時間がかかり、狭窄部位探索終了までに長時間を要してしまうおそれがある。

特開平０７−２７１９７６号公報

本発明は、以上の問題点に鑑みてなされたものであって、その目的は、画像から血管の狭窄部位を探索する技術に関し、特に探索処理時間の短縮を図りつつ探索精度を向上させ、効率のよい狭窄部位探索を行う技術を提供することにある。

上記課題を解決するために、本発明は、血流速度を画素値に含む機能画像と部位の形態を表す形態画像とを記憶する画像記憶手段と、前記機能画像に基づき血管各所の血流速度を算出する血流速度算出手段と、前記血流速度算出手段により算出された血流速度から閾値を作成する閾値作成手段と、前記閾値を記憶する閾値記憶手段と、前記血流速度と前記閾値とを比較する速度比較手段を含み、前記閾値以上の値を示す血流速度に対応する血管位置を第１の狭窄部位候補として探索する第一の探索手段と、前記血流速度の血管位置に対する変化特性に基づき、血管に沿って並ぶ血流速度列を位置毎に各周波数成分に分解する分解手段を含み、各位置の各周波数成分の傾向を基に、高周波数成分の大きさがウエーブレット分解レベルの増加に従い減少する血管位置を第二の狭窄部位候補として探索する第二の探索手段と、前記第一の狭窄部位候補及び前記第二の狭窄部位候補の双方を満たす血管位置を狭窄部位として、その位置情報を前記形態画像上に付加する付加手段と、を備えること、を特徴とする。

本発明によっては、機能画像に基づき血管各所の血流速度を算出しておき、この血流速度の遅速に基づき第一の狭窄部位候補を探索し、さらに前記血管速度の血流変化に基づき第二の狭窄部位候補を探索し、第一の狭窄部位候補及び第二の狭窄部位候補の双方を満たす血管位置を狭窄部位とした。

診断装置を用いた観察によって直接得られる血流速度を利用して狭窄部位を探索するため、血管輪郭抽出の精度に依らず精度の高い狭窄部位探索を行うことが可能となり、そのため処理時間が短縮する。

以下、本発明に係る血管狭窄部位抽出の好適な各実施形態について、図面を参照しながら詳細に説明する。

図１は、本発明に係る血管狭窄部位抽出を実施する血管狭窄部位探索装置の構成を示すブロック図である。図１に示す血管狭窄部位探索装置１は、画像内から血管の狭窄部位を探索してその位置を表示する装置である。狭窄部位は、血管幅が急速に狭くなっている。そのため、狭窄部位付近では、血流速度が速く、かつ狭窄部位前後の血流速度と比べてその血流速度は特異的であり血流速度が速くなっている。血管狭窄部位探索装置１は、血流速度の遅速と空間的推移に基づき狭窄部位を探索する。

この血管狭窄部位探索装置１は、画像記憶部１０と、血管位置特定部２０と、遅速探索部３０と、血流変化探索部４０と、狭窄部位付加部５０と、狭窄率算出部６０と、表示部７０とを備える。これら各構成は、コンピュータにより構成される。又はこれら各構成は、それぞれが専用の回路で構成される。

コンピュータにより構成される血管狭窄部位探索装置１は、例えば演算制御部（ＣＰＵ：Central Processing Unitやグラフィックチップ）、主記憶部（ＲＡＭ：Random Access Memory）、外部記憶部（ＨＤＤ：Hard Disk Drive）を共通線で電気的に接続して相互にデータ入出力可能としたコンピュータと、このコンピュータにデータ入出力可能に接続されたモニタ及び入力インターフェースで構成される。

外部記憶部に血管狭窄部位抽出を実現させるプログラムを記憶し、そのプログラムを適宜主記憶部に展開して演算制御部が実行する。モニタは、ＬＣＤ（Liquid Crystal Display）ディスプレイあるいはＣＲＴ(Cathode Ray Tube)ディスプレイで構成される。入力インターフェースは、キーボードやホイールトラッキング機能を有するマウス等である。

この血管狭窄部位探索装置１は、ワークステーションで実施してもよいし、コンピュータを搭載した画像診断装置で実施してもよい。画像診断装置としては、例えばＰＥＴ−ＣＴ装置、Ｘ線ＣＴ装置、Ｘ線診断装置、ＭＲ装置、核医学診断装置（ＰＥＴ装置又はＳＰＥＣＴ装置）等の機能画像及び形態画像が撮影可能な単一の機器又は混合のシステムが挙げられる。

画像記憶部１０は、外部記憶部を含み構成される。画像記憶部１０には、医用画像が記憶されている。この医用画像は、被検体の同一箇所を撮影して得られた被検体を基準にした座標系において同一領域を有する機能画像及び形態画像の２種類が記憶されている。

医用画像は、被検体を画像診断装置で撮影することで得られた被検体内の画像である。この医用画像は、３次元の実空間で表されるボリュームデータ（ボクセルデータ）である。ボリュームデータは、複数の２次元で表される画像データを読み込み、補間処理を行うことにより生成される。この医用画像は、血管狭窄部位探索装置１に３次元画像生成部を備えるようにして予め生成しておいてもよいし、他のワークステーション等のコンピュータやボリュームデータの生成が可能な画像診断装置において生成し、ネットワーク又はＣＤ−ＲＯＭ等の可搬記憶媒体を介して画像記憶部１０に記憶させるようにしてもよい。

形態画像は、ＰＥＴ−ＣＴ装置のＣＴ撮影手段、Ｘ線ＣＴ装置、又はＸ線診断装置で撮影されたＣＴ画像、若しくはＭＲＩ装置で撮影されたＭＲＩ画像である。形態画像を構成するボクセルは、座標上に３次元に配列され、その画素値によりボクセルが対応する位置にある部位の形態を表している。

機能画像は、ＰＥＴ−ＣＴ装置のＰＥＴ撮影手段、又は核医学診断装置で撮影された核医学画像、若しくはＭＲＩ装置で撮影されたＭＲＩ画像である。特に本実施形態の機能画像は、画像領域内の血流速度を視覚的に表現した画像である。各ボクセルには、血流速度に応じた信号強度を示す値が画素値として付与される。

この機能画像は、例えばＭＲ装置を用いてＰＳ（ＰｈａｓｅＳｈｉｆｔ）法によりＭＲ信号を取得することにより生成される。フローエンコードパルスと呼ばれる双極性の傾斜磁場を含むパルスシーケンスによってＭＲ信号を収集して複素差分すると、血流速度に比例した位相シフトのシフト量を示す信号強度のＭＲ信号が得られる。この信号強度を画素値とすることで、血流速度を視覚的に表現した機能画像が生成される。

血管位置特定部２０は、演算制御部を含み構成される。この血管位置特定部２０は、機能画像上で血管の位置を特定する。血管狭窄部位探索装置１は、機能画像を用いて血流速度の遅速と位置上の変化とから狭窄部位を探索するため、この血管位置特定部２０により血流速度情報を取得する領域を機能画像上で特定する。血管位置特定部２０は、血管の形態が表れている形態画像から血管位置を特定し、血管の中心線を抽出する。複雑に分岐する血管中心線を単線に単離し、単離した各血管中心線を被検体の座標系に従って機能画像上に当てはめる。

遅速探索部３０は、演算制御部を含み構成される。この遅速探索部３０は、機能画像から血管各所の血流速度を算出し、血流速度の遅速から第１の狭窄部位候補を探索する。血流速度の算出では、血管中心線上の各ボクセルの画素値をパラメータとする。所定以上の値を有する血流速度情報を選別し、選別した血流速度情報に対応する血管位置を第一の狭窄部位候補とする。選別した血流速度情報に対応する血管位置は、その血流速度情報を取得するために用いた画素値を有する血管位置である。

血流変化探索部４０は、演算制御部を含み構成される。この血流変化探索部４０は、血流変化情報を取得し、血流変化情報から第２の狭窄部位候補を探索する。血流変化情報は、血流速度情報の位置に対する変化特性であり、血流速度の空間的推移を示している。血流変化情報をウェーブレット変換して位置毎の周波数成分に分解し、その高周波数成分の振幅がウェーブレット分解レベルの増加に従い減少する血管位置を第二の狭窄部位候補とする。

狭窄部位付加部５０は、演算制御部を含み構成される。この狭窄部位付加部５０は、第一の狭窄部位候補及び第二の狭窄部位候補の双方を満たす血管位置を狭窄部位とし、形態画像上にその位置情報を付加する。位置情報は、狭窄部位領域に相当するボクセルを白色等の単色で塗り潰して識別可能とした印である。第１の狭窄部位候補は、血流速度が他と比べて速い位置が選別されたものであるが、この第１の狭窄部位候補には血管が細い末梢血管も含まれる。従って第２の狭窄部位候補と合わせることで、精度良く狭窄部位を抽出ことができる。

狭窄率算出部６０は、演算制御部を含み構成される。この狭窄率算出部６０は、狭窄部位の狭窄率を算出する。狭窄部位の血管断面径と狭窄部位の近傍の正常血管の血管断面径の比を算出する。

表示部７０は、モニタを含み構成される。この表示部７０は、狭窄部位が位置情報によって示された形態画像と、算出された狭窄率とを表示する。

図２に、この血管狭窄部位探索装置１の概略動作を示す。図２は、血管狭窄部位探索装置１の血管狭窄部位抽出動作を示すフローチャートである。

まず、血管位置特定部２０は、画像記憶部１０から形態画像を読み出し（Ｓ０１）、血管中心線を抽出する（Ｓ０２）。血管中心線を抽出すると、画像記憶部１０から機能画像を読み出し（Ｓ０３）、血管中心線を機能画像に当てはめる（Ｓ０４）。

機能画像上に血管中心線が当てはめられることにより、機能画像上で血管位置が特定されると、遅速探索部３０は、機能画像上の血管中心線各所の血流速度情報を取得する（Ｓ０５）。血流速度情報を取得すると、各血流速度情報と所定の閾値とを比較し（Ｓ０６）、閾値以上の値を示す血流速度情報を有する血管位置を第一の狭窄部位候補として記憶する（Ｓ０７）。

第一の狭窄部位候補を記憶すると、血流変化探索部４０が、血流速度情報の位置に対する変化特性を示す血流変化情報を位置毎の周波数成分に分解する（Ｓ０８）。各分解レベルの高周波数成分の傾向を位置毎に解析して血流速度が前後の血管位置に比べて特異な血管位置を検索し（Ｓ０９）、該当する血管位置を第二の狭窄部位候補として記憶する（Ｓ１０）。

第二の狭窄部位候補が記憶されると、狭窄部位付加部５０が、第一の狭窄部位候補及び第二の狭窄部位候補の双方を満たす血管位置を狭窄部位とし、形態画像上にその位置情報を付加する（Ｓ１１）。

狭窄部位が位置情報により示されると、狭窄率算出部６０が、その位置情報で示される狭窄部位の血管断面径と狭窄部位の近傍にある正常血管の血管断面径とを算出して、その比をとって狭窄率を算出する（Ｓ１２）。

狭窄率が算出されると、表示部７０が、狭窄部位が位置情報によって示された形態画像と、算出された狭窄率とを表示する（Ｓ１３）。

このような血管狭窄部位探索装置１における血管位置特定についてさらに詳細に説明する。

図３は、血管位置特定部２０の血管位置特定の動作を示すフローチャートである。血管位置の特定では、形態画像から血管像を抽出し、抽出された血管像から血管中心線を抽出し、中心線のネットワークを単線に単離し、単離された各血管中心線を機能画像上にマッチングする。

まず、Ｓ２１において、血管位置特定部２０は、血管の始点座標と終点座標を取得する。

血管位置特定部２０は、位置指定用のカーソルをモニタ上に描画する。具体的には、グラフィックメモリにカーソルを書き込む。操作者が入力インターフェースを用いてカーソルを移動させ、狭窄部位を挟み込むように始点及び終点をマークすると、血管位置特定部２０は、始点座標情報及び終点座標情報を記憶する。

血管の始点座標情報と終点座標情報を取得すると、血管像抽出部２１は、血管の仮の中心線を決定する（Ｓ２２）。

図４は、血管の軸方向断面を示す模式図である。血管の仮の中心線の決定では、ｉ番目の仮の中心点は、ｉ−１番目、ｉ−２番目・・・の移動平均方向Ｌ１にあるｉ−１番目の隣のボクセルを含む平面Ｐ１上に定められる。この平面Ｐ１は、ｉ番目ｉ−１番目、ｉ−２番目・・・の移動平均方向Ｌ１と直交する平面であり、血管をわぎりしたものである。

血管像抽出部２１は、この平面Ｐ１の濃度を示す濃度プロファイルＧ１を生成する。図５に濃度プロファイルＧ１を示す。濃度プロファイルＧ１は、ｉ番目ｉ−１番目、ｉ−２番目・・・の移動平均方向Ｌ１と直交する平面Ｐ１の濃度分布を示している。濃度は、平面Ｇ１上の各ボクセルの画素値である。

濃度プロファイルＧ１を生成すると、血管位置特定部２０は、濃度プロファイルＧ１上を、濃度の極小値ｄ_ｍｉｎを求めてサーチする。形態画像の場合、血管壁と血管内部とでは、血管壁のほうが画素値が高く得られ、血管内部は血管壁の輝度値よりも低い傾向にあるためである。尚、形態画像がサブストラクションの場合は、極大値を求めてサーチする。

極小値ｄ_ｍｉｎが得られると、血管位置特定部２０は、この極小値ｄ_ｍｉｎを画素値に有するボクセルをｉ番目の仮の中心Ｃｋとする。血管位置特定部２０は、この仮の中心Ｃｋの決定プロセスを始点座標情報に対応する始点から終点座標に対応する終点まで繰り返す。

仮の中心線が決定すると、血管位置特定部２０は、血管縁を確定する（Ｓ２３）。

血管抽出部２１は、仮の中心Ｃｋ毎の濃度プロファイルＧ１から、それぞれ極小値ｄ_ｍｉｎを持つ点とその隣り合う位置の濃度値の差が極大になる濃度値ｄ_１を有する点、換言すると一次微分の最大点を決める。この一次微分の最大点は、仮の血管内壁縁となる。このとき、仮の中心の濃度値をｄ_ｍｉｎ、一次微分の最大点の濃度をｄ_１とする。仮の中心の濃度値をｄ_ｍｉｎ、一次微分の最大点の濃度をｄ_１から、バックグラウンド領域を決定するための収束判定幅εを以下式１から求める。

式１：ε＝（ｄ_１−ｄ_ｍｉｎ）×β
βは、形態画像の条件により決定される収束判定用のパラメータであり、予め適宜な値に決められている。

この後、血管位置特定部２０は、濃度プロファイル上で隣り合う位置の濃度差を順次演算し、この濃度差が収束判定幅εに収まる点の集まりをバックグラウンド領域とする。このバックグラウンド領域は、血管内を示す領域である。

次いで、バックグラウンド領域の濃度値ｄ_ｂを平均値演算により求め、血管壁か否かを弁別するための閾値ｄ_ｔｈを式２から求める。

式２：ｄ_ｔｈ＝ｄ_ｂ−（ｄ_ｂ−ｄ_ｍｉｎ）×α
αは、閾値を決定するパラメータで、予め適宜に決められている。

閾値ｄ_ｔｈが定まると、この閾値ｄ_ｔｈに最も近い値を画素値に有するボクセルを濃度プロファイルからサーチし、サーチしたボクセルを血管縁Ｅとして確定する。

血管縁Ｅが確定すると、血管位置特定部２０は、抽出された血管縁Ｅを基に血管中心線を抽出する（Ｓ２４）。

図６は、血管位置特定部２０による血管中心線の抽出を示す模式図である。血管位置特定部２０は、確定した血管像の各断面領域において、血管縁にある任意のボクセルＢｓから垂線Ｌ_２を下ろし、垂線Ｌ_２と血管縁との交点のボクセルＢｅを求める。血管位置特定部２０は、ボクセルＢｓ付近の血管縁を構成する座標データから、ボクセルＢｓの接線ベクトルを求め、この接線ベクトルに垂直なベクトルを垂線Ｌとする。

この垂線Ｌに沿ってサーチ点Ｓを移動させ、サーチ点Ｓの座標データが血管縁を構成する一点の座標データと一致すれば、当該座標データの点をボクセルＢｅとする。ボクセルＢｅが求められると、血管位置特定部２０は、ボクセルＢｓとボクセルＢｅの座標データからボクセルＢｓとボクセルＢｅの中点となるボクセルＢｃを演算する。

このボクセルＢｃの演算を血管縁にある各ボクセルをボクセルＢｓとして演算する。この結果、ボクセルＢｃが複数存在すれば、それら全てのボクセルＢｃの重心となるボクセルＢｇを演算し、この重心のボクセルＢｇを一断面における血管の中心とする。

血管中心線が抽出されると、血管位置特定部２０は、血管中心線を補正し、血管中心線を確定する（Ｓ２５）。

血管中心線の各点は、連結性を有するはずである。血管位置特定部２０は、血管中心線の連結性を２６近傍検索により調べ、非連結部分に血管中心となる点を追加する。

血管位置特定部２０は、血管中心線上の隣り合うｉ番目のボクセルＢｇ_ｉ（ｘ，ｙ，ｚ）とｉ＋１番目のボクセルＢｇ_ｉ＋１（ｘ＋１，ｙ＋１，ｚ＋１）の距離を座標データから算出する。この距離が１ボクセル分でなければ、即ちボクセルＢｇ_ｉ（ｘ，ｙ，ｚ）とボクセルＢｇ_ｉ＋１（ｘ＋１，ｙ＋１，ｚ＋１）とが隣接していなければ、ボクセルＢｇ_ｉ（ｘ，ｙ，ｚ）とボクセルＢｇ_ｉ＋１（ｘ＋１，ｙ＋１，ｚ＋１）との中点に新たな血管中心点を追加する。

血管中心線が確定すると、血管位置特定部２０は、血管中心線のネットワークを単線に単離する（Ｓ２６）。

図７は、血管中心線の単離を示す模式図である。血管位置特定部２０は、確定した血管中心線の座標データを、被検体の頭から足側、被検体の腹側から背中側、被検体の左側の順番に辿り、血管中心線の座標データが終端に行き着くと、その辿った線を単離する。即ち、血管中心線として確定した座標データの集まりから分離する。血管位置特定部２０は、血管中心線として確定した座標データが残っている限り、この単離を繰り返す。

血管中心線が単線に単離されると、血管位置特定部２０は、この血管中心線を機能画像に当てはめて、機能画像上において血管中心線を確定する（Ｓ２７）。

形態画像と機能画像は、被検体内の同一箇所を撮影したものである。血管位置特定部２０は、形態画像の空間解像度と機能画像の空間解像度から、空間解像度の比を求め、その比に合わせて形態画像における血管中心線の座標データを変倍する。

血管狭窄部位探索装置１は、このＳ２１〜Ｓ２７の動作を行うことにより、機能画像上において血管中心線を単線毎に抽出する。

尚、血管中心線の特定は、その他の公知の技術を用いて行ってもよく、例えばＳｈａｐｅＤｒｉｖｅｎＦｌｏｗ法（参考文献：Vessel Segmentation Using A Shape Driven Flow. In Proceedings of MICCAI, St.Malo, France, (2004)）により血管の輪郭を抽出してから血管中心線を求めてもよく、また血管中心線の補正をＬｏｃａｌＳＳ−ｃｏｄｉｎｇ法（参考文献：Efficient Skeletonization of Volumetric Objects, IEEE TRANSACTIONS ON VISUALIZATION AND COMPUTER GRAPHICS, Vol. 5, No.3, pp.196 - 209 , July-Sept (1999)）を用いて行ってもよい。

血管中心線が機能画像上において特定されると、血管狭窄部位探索装置１は、血管中心線の各所における血流速度から狭窄部位の候補を選出する。この血管狭窄部位探索装置１における第１の狭窄部位探索についてさらに詳細に説明する。図８は、遅速探索部３０の第１の狭窄部位候補を探索する動作を示すフローチャートである。

まず、遅速探索部３０は、機能画像の血管各所から血流速度情報を取得する（Ｓ３１）。

血管速度情報を取得する位置は、血管中心線上の各点である。血管中心線上の血流が、血管の曲がりくねりや分岐による影響の最も少ないからである。ここで、ＰＳＭＲＡ画像等の機能画像は、血流速度に比例して位相シフトしたプロトンスピンのシフト量を多値階調化したものである。従って、遅速探索部３０は、フローエンコード量に血管中心線上の各点に対応するボクセルの画素値を乗じることで血管各所の血流速度情報を取得する。すなわち、遅速探索部３０は、以下の式３を血管中心線上の各点について演算して血流速度を得る。

式３：Ｖ（ｘ，ｙ，ｚ）＝Ｉ（ｘ，ｙ，ｚ）×ＦＥＰ／１８０００

Ｖ（ｘ，ｙ，ｚ）は、血管中心線上の座標（ｘ，ｙ，ｚ）で示される位置の血流速度情報である。Ｉ（ｘ，ｙ，ｚ）は、血管中心線上の座標（ｘ，ｙ，ｚ）で示される位置の画素値である。ＦＥＰは、フローエンコード量（ｃｍ／ｓｅｃ）である。「１８０００」は、スケール値である。尚、フローエンコード量は、位相シフトをｃｍ／ｓｅｃに換算して用いる。スケール値は、階調やフローエンコード量に伴って適宜用意される。

血管各所の血流速度情報を取得すると、遅速探索部３０は、血流速度情報の値との比較対象となる閾値を生成する（Ｓ３２）。

遅速探索部３０は、算出された血管各所における血流速度情報の平均値を算出し、この平均値を閾値とする。狭窄部位では、他の血管各所に対して血液の流れが速くなっているためである。

閾値を生成すると、遅速探索部３０は、血流速度情報と閾値とを比較する（Ｓ３３）。比較により、閾値を上回る血流速度情報が選出されると、この選出された血流速度情報を有する血管位置を第１の狭窄部位として記憶する（Ｓ３４）。

血管狭窄部位探索装置１は、このＳ３１〜Ｓ３４の動作を行うことにより、血流速度から狭窄部位の候補を探索する。

第１の狭窄部位候補が探索されると、血管狭窄部位探索装置１は、血管中心線に沿って並ぶ血流速度の血流変化から第２の狭窄部位候補を選出する。この血管狭窄部位探索装置１における第２の狭窄部位探索についてさらに詳細に説明する。

図９は、血流変化探索部４０の第２の狭窄部位候補を探索する動作を示すフローチャートである。

まず、血流変化探索部４０は、血流変化情報を読み出し、読み出した血流変化情報をウェーブレット変換する（Ｓ４１）。

図１０は、血流変化情報を示し、図１１は、血流変化情報を位置毎に周波数分解したものである。

図１０に示すように、血流変化情報は、血管中心線上の各位置を横軸とし、血流速度情報を縦軸とした血流速度の位置に対する変化特性を示す血流変化曲線ｆ（ｘ，ｙ，ｚ）である。血流変化探索部４０は、単線の血管中心線上の血流速度情報をボクセル毎に読み出し、血流速度情報を血管位置の並びに従って並べることにより、血流変化情報を取得する。この血流変化曲線ｆ（ｘ，ｙ，ｚ）は、図中の位置Ａに示される位置で血流速度が特異的に変化しており、図中の範囲Ｂでは血流速度が安定的に変化している。この位置Ａが狭窄部位の第２の候補であり、血流変化探索部４０は、この位置Ａを周波数解析により探索する。

血流速度が特異的に変化している位置Ａを探索すべく、血流変化探索部４０は、血流変化曲線ｆ（ｘ，ｙ，ｚ）をウェーブレット変換し、図１１に示すような各分解レベルのウェーブレット詳細化係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_ｓ＝２ ^ｊを取得する。周波数分解としては、血流変化情報の横軸情報を保持でき、かつ各周波数成分のパワーを読み取ることのできる周波数領域及び空間領域を変数としたウェーブレット変換等の周波数分解が望ましい。

ウェーブレット変換は、低周波数領域をその領域において高周波数領域と低周波数領域へ分解する工程をツリー状に繰り返すオクターブ分割として捉えることができる。例えば、ウェーブレット変換によりデータ列を分解レベルｊまで分解することで、血流変化曲線ｆ（ｘ，ｙ，ｚ）から、血管中心線上の位置（ｘ，ｙ，ｚ）毎の各分解レベルのウェーブレット詳細化係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_ｓ＝２ ^ｊが取得される。ｓは、スケールである。

このウェーブレット詳細化係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_ｓ＝２ ^ｊは、分解レベルｊに対応したスケールの周波数成分の各位置における大きさを示す。この周波数成分の大きさは、その周波数スケールにおける血流速度の空間的変位量である。ウェーブレット変換では、分解レベルが進む毎に、対象のデータ列を高周波数成分と低周波成分とに分解し、さらに次の分解レベルで低周波成分を高周波成分と低周波成分に分解していく。従って、ウェーブレット詳細化係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_ｓ＝２ ^ｊは、その分解レベルが進むつれ、より低周波の成分を示している。

周波数分解して位置（ｘ，ｙ，ｚ）毎の各周波数成分を得ると、Ｓ４２〜Ｓ４４において、血流変化探索部４０は、同一血管位置の各周波数成分を各位置において比較して、血流速度の周波数成分が低周波成分になればなるほど小さくなる血管位置を探索する。

血流速度が前後の位置と比べて特異的な血管位置（図中Ａ）では、血流速度が急峻になるため高周波成分が大きい（図１１参照）。一方、この特異的な血管位置の低周波成分は、相対的に小さくなる。即ち、血流速度が前後の位置と比べて特異的な血管位置（図中Ａ）では、血流速度の周波数成分が低周波成分になればなるほど小さくなる。

一般的に、リプシッツ（Ｌｉｐｓｃｈｉｔｚ）定数をＬ＝Ｓ＾αとしたとき、ウェーブレット詳細化係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_ｓ＝２ ^ｊとリプシッツ定数とは、以下の式４に従うことが分かっている（参考文献：S.Mallat and W.L.Hwang,“Singularity detection and processing with wavelets,”IEEE Trans. Inform. Theory,vol.38,pp.617-643,Mar. 1992.）
式４：｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_ｓ＝２ ^ｊ≦Ｓ^α

流速が前後位置と比較して特異的な部分は、αが負数となるから、上記式にも示されるとおり、各分解レベルのウェーブレット詳細化係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_ｓ＝２ ^ｊが分解レベルｊの増加によって小さくなっていく。

血流変化探索部４０は、まず、ウェーブレット変換により得られたウェーブレット詳細化係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_ｓ＝２ ^ｊのデータ列に表れている複数の振幅範囲のそれぞれについて、ウェーブレット最大係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_{ｍａｘｉｍｕｍ}を算出する（Ｓ４２）。ウェーブレット最大係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_{ｍａｘｉｍｕｍ}は、一つの振幅範囲における最大振幅値である。血流速度の変位量は、狭窄部位で最大となるからである。

次に、このウェーブレット最大係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_{ｍａｘｉｍｕｍ}を有する血管位置（ｘ，ｙ，ｚ）における各分解レベルのウェーブレット最大係数｜Ｗｓｆ（ｘ，ｙ，ｚ）｜_{ｍａｘｉｍｕｍ}を比較する（Ｓ４３）。比較の結果、分解レベルの増加に従って小さくなっていれば（Ｓ４３，Ｙｅｓ）、Ｓ４４において、この血管位置（ｘ，ｙ，ｚ）を第２の狭窄部位候補として記憶する。

これにより、血管狭窄部位抽出装置１は、血流速度から第１の狭窄部位候補を記憶し、かつ血流変化から第２の狭窄部位候補を記憶する。第１の狭窄部位候補と第２の狭窄部位候補が探索されると、狭窄部位付加部５０により、狭窄部位が確定され、形態画像上にその位置情報が付加される。この狭窄部位の確定動作をさらに詳細に説明する。

図１２は、この狭窄部位の位置情報付加の動作を示すフローチャートである。

まず、狭窄部付加部５０は、第１の狭窄部位候補と第２の狭窄部位候補となった座標データを比較する（Ｓ５１）。比較の結果、第１の狭窄部位候補と第２の狭窄部位候補に共通の座標データがあれば（Ｓ５１，Ｙｅｓ）、当該座標データを狭窄部位の座標データとして記憶する（Ｓ５２）。

次に、形態画像を読み出し（Ｓ５３）、記憶した狭窄部位の座標データに相当するボクセルとその周辺範囲のボクセルに単色の画素値等の位置情報を付する（Ｓ５４）。

血管狭窄部位探索装置１は、このＳ５１〜Ｓ５４の動作を行うことにより、狭窄部位の位置情報を形態画像上に付与する。狭窄部位の位置情報が付与されると、血管狭窄部位探索装置１は、狭窄部位の狭窄率を算出する。狭窄部位の狭窄率算出についてさらに詳細に説明する。

図１３は、狭窄率算出部６０による狭窄率算出動作を示すフローチャートである。また、図１４は、狭窄率算出を説明するための血管断面の模式図である。

狭窄率算出部６０は、狭窄部位の周辺にある正常血管の位置を選択する（Ｓ６１）。狭窄率算出部６０は、狭窄部位の座標Ｐｓと血管中心線の位置情報を読み出す。狭窄部位の座標Ｐｓを起点として、血管中心線に沿って所定ボクセル分離れた位置を正常血管の座標Ｐｎとする。

狭窄部位の座標Ｐｓを読み出し、また正常血管の座標Ｐｎが特定されると、狭窄率算出部６０は、狭窄部位の血管断面径Ｄｓ及び正常血管の血管断面径Ｄｎを算出する（Ｓ６２）。

狭窄率算出部６０は、血管位置特定部２０により特定された血管縁の位置データを読み出し、狭窄部位の座標Ｐｓ及びその前後の血管中心線上の座標の移動平均線を算出し、その移動平均線と直交する断面を求める。この断面に含まれる血管縁を構成する任意の起点ボクセルから断面に沿った垂線を下ろし、垂線と血管縁との交点の終点ボクセルを求める。

狭窄率算出部６０は、起点ボクセル付近の血管縁線を構成する座標データから、起点ボクセルの接線を求め、この接線に垂直なベクトルを垂線とする。起点ボクセルと終点ボクセルが特定されると、この起点ボクセル及び終点ボクセルの座標から起点ボクセルと終点ボクセルの距離を算出する。

狭窄部位算出部６０は、任意の起点ボクセルを血管縁に沿って複数設定し、各起点ボクセル毎に距離を求め、平均値をとる。この平均値を狭窄部位の血管断面径Ｄｓとする。正常血管の血管断面径Ｄｎについても同様に算出する。

狭窄部位及び正常血管の血管断面径を求めると、狭窄率算出部６０は、これら血管断面径の比をとり、狭窄率を取得する（Ｓ６３）。すなわち、狭窄率＝（狭窄部位の血管断面径Ｄｓ）／（正常血管の血管断面径Ｄｎ）×１００を演算する。狭窄率は、百分率で表すようにしてもよい。

血管狭窄部位探索装置１は、狭窄部位の探索及び狭窄率の算出を行うと、狭窄部位の位置情報が付与された形態画像及び狭窄率をモニタに表示する。図１５は、表示部７０により表示された表示画面を示す図である。

モニタには、形態画像が表示される。形態画像には、狭窄部位Ｓｐが白色等の単色で塗り潰されて表示される。また、その狭窄部位Ｓｐの近傍に狭窄率Ｓｒを示す文字列が表示される。

このように本実施形態に係る血管狭窄部位探索装置１は、血流速度の遅速から第一の狭窄部位候補を探索し、血管速度の血管位置に対する変化特性に基づき第二の狭窄部位候補を探索し、第一の狭窄部位候補及び前記第二の狭窄部位候補の双方を満たす血管位置を狭窄部位とした。

血管の輪郭を抽出するという画像処理を直接の狭窄部位探索手段とせず、診断装置による観察によって直接得られる血流速度や血流変化を直接の狭窄部位探索手段としているため、血管輪郭抽出の精度に依らず精度の高い狭窄部位探索を行うことが可能となる。

尚、本実施形態では、医用画像を３次元画像として説明したが、２次元画像であってもよく、２次元画像であっても本発明をそのまま適用して狭窄部位を探索することができる。

血管狭窄部位探索装置の構成を示すブロック図である。血管狭窄部位抽出動作を示すフローチャートである。血管位置特定の動作を示すフローチャートである。血管の軸方向断面を示す模式図である。濃度プロファイルを示す模式図である。血管中心線の抽出を示す模式図である。血管中心線の単離を示す模式図である。第１の狭窄部位候補を探索する動作を示すフローチャートである。第２の狭窄部位候補を探索する動作を示すフローチャートである。血流変化曲線を示す図である。血流変化曲線をウェーブレット変換して得た各分解レベルのウェーブレット詳細化係数Ｗｓｆ（ｘ，ｙ，ｚ）_ｓ＝２ ^ｊを示す図である。狭窄部位の位置情報付加の動作を示すフローチャートである。狭窄率算出動作を示すフローチャートである。血管断面の模式図である。表示画面を示す図である。

符号の説明

１血管狭窄部位探索装置
１０画像記憶部
２０血管位置特定部
３０遅速探索部
４０血流変化探索部
５０狭窄部位付加部
６０狭窄率算出部
７０表示部

Claims

血流速度を画素値に含む機能画像と部位の形態を表す形態画像とを記憶する画像記憶手段と、
前記機能画像に基づき血管各所の血流速度を算出する血流速度算出手段と、
前記血流速度算出手段により算出された血流速度から閾値を作成する閾値作成手段と、前記閾値を記憶する閾値記憶手段と、
前記血流速度と前記閾値とを比較する速度比較手段を含み、
前記閾値以上の値を示す血流速度に対応する血管位置を第１の狭窄部位候補として探索する第一の探索手段と、
前記血流速度の血管位置に対する変化特性に基づき、
血管に沿って並ぶ血流速度列を位置毎に各周波数成分に分解する分解手段を含み、
各位置の各周波数成分の傾向を基に、高周波数成分の大きさがウエーブレット分解レベルの増加に従い減少する血管位置を第二の狭窄部位候補として探索する第二の探索手段と、
前記第一の狭窄部位候補及び前記第二の狭窄部位候補の双方を満たす血管位置を狭窄部位として、その位置情報を前記形態画像上に付加する付加手段と、
を備えること、
を特徴とする血管狭窄部位探索装置。
前記血管各所を特定する血管位置特定手段をさらに備えること、
を特徴とする請求項１記載の血管狭窄部位探索装置。
前記血管位置特定手段は、
前記形態画像に含まれる血管像の中心線を抽出する中心線抽出手段を含み、
前記第一の探索手段及び前記第二の探索手段は、
前記中心線上の各点から狭窄部位候補を探索すること、
を特徴とする請求項２記載の血管狭窄部位探索装置。
前記形態画像から前記狭窄部位の狭窄率を算出する狭窄率算出手段をさらに備えること、
を特徴とする請求項１記載の血管狭窄部位探索装置。
前記形態画像は、ＣＴ画像若しくはＭＲＩ画像であること、
を特徴とする請求項１記載の血管狭窄部位探索装置。
前記機能画像は、ＭＲＩ画像若しくは核医学画像であること、
を特徴とする請求項１記載の血管狭窄部位探索装置。
血流速度を画素値に含む機能画像を撮影する第１の撮影手段と、
形態画像を撮影する第２の撮影手段と、
前記機能画像と前記形態画像とを記憶する画像記憶手段と、
前記機能画像に基づき血管各所の血流速度を算出する血流速度算出手段と、
前記血流速度算出手段により算出された血流速度から閾値を作成する閾値作成手段と、前記閾値を記憶する閾値記憶手段と、
前記血流速度と前記閾値とを比較する速度比較手段を含み、
前記閾値以上の値を示す血流速度に対応する血管位置を第１の狭窄部位候補として探索する第一の探索手段と、
前記血流速度の血管位置に対する変化特性に基づき、
血管に沿って並ぶ血流速度列を位置毎に各周波数成分に分解する分解手段を含み、
各位置の各周波数成分の傾向を基に、高周波数成分の大きさがウエーブレット分解レベルの増加に従い減少する血管位置を第二の狭窄部位候補として探索する第二の探索手段と、
前記第一の狭窄部位候補及び前記第二の狭窄部位候補の双方を満たす血管位置を狭窄部位として、その位置情報を前記形態画像上に付加する付加手段と、
を備えること、
を特徴とする画像診断装置。
血流速度を画素値に含む機能画像と形態画像とを記憶するコンピュータを、
前記機能画像に基づき血管各所の血流速度を算出する血流速度算出手段と、
前記血流速度算出手段により算出された血流速度から閾値を作成する閾値作成手段と、前記閾値を記憶する閾値記憶手段と、
前記血流速度と前記閾値とを比較する速度比較手段を含み、
前記閾値以上の値を示す血流速度に対応する血管位置を第１の狭窄部位候補として探索する第一の探索手段と、
前記血流速度の血管位置に対する変化特性に基づき、
血管に沿って並ぶ血流速度列を位置毎に各周波数成分に分解する分解手段を含み、
各位置の各周波数成分の傾向を基に、高周波数成分の大きさがウエーブレット分解レベルの増加に従い減少する血管位置を第二の狭窄部位候補として探索する第二の探索手段と、
前記第一の狭窄部位候補及び前記第二の狭窄部位候補の双方を満たす血管位置を狭窄部位として、その位置情報を前記形態画像上に付加する付加手段と、
して機能させること、
を特徴とする血管狭窄部位探索プログラム。