(58)【調査した分野】(Int.Cl.,DB名)
前記気流速度分布算出工程は、左右の前記肺野領域のそれぞれを体幹軸方向の複数のブロック領域に分割し、左右の前記肺野領域のそれぞれについて体幹軸方向の気流速度の分布を表す特徴量を算出する請求項1〜7の何れか一項に記載の胸部診断支援情報生成方法。
【発明を実施するための形態】
【0020】
以下、図面を参照して本発明の実施の形態を詳細に説明する。ただし、発明の範囲は、図示例に限定されない。
【0021】
−第1の実施の形態−
〔胸部診断支援システム100の構成〕
まず、構成を説明する。
図1に、本実施の形態における胸部診断支援システム100の全体構成を示す。
図1に示すように、胸部診断支援システム100は、撮影装置1と、撮影用コンソール2とが通信ケーブル等により接続され、撮影用コンソール2と、診断用コンソール3とがLAN(Local Area Network)等の通信ネットワークNTを介して接続されて構成されている。胸部診断支援システム100を構成する各装置は、DICOM(Digital Image and Communications in Medicine)規格に準じており、各装置間の通信は、DICOMに則って行われる。
【0022】
〔撮影装置1の構成〕
撮影装置1は、例えば、呼吸運動に伴う肺の膨張及び収縮の形態変化、心臓の拍動等の、周期性(サイクル)を持つ胸部の動態を撮影する装置である。動態撮影は、人体の胸部に対し、X線等の放射線を連続照射して複数の画像を取得(即ち、連続撮影)することにより行う。この連続撮影により得られた一連の画像を動態画像と呼ぶ。また、動態画像を構成する複数の画像のそれぞれをフレーム画像と呼ぶ。
撮影装置1は、
図1に示すように、放射線源11、放射線照射制御装置12、放射線検出部13、読取制御装置14等を備えて構成されている。
【0023】
放射線源11は、被写体Mを挟んで放射線検出部13と対向する位置に配置され、放射線照射制御装置12の制御に従って、被写体Mに対し放射線(X線)を照射する。
放射線照射制御装置12は、撮影用コンソール2に接続されており、撮影用コンソール2から入力された放射線照射条件に基づいて放射線源11を制御して放射線撮影を行う。撮影用コンソール2から入力される放射線照射条件は、例えば、連続照射時のパルスレート、パルス幅、パルス間隔、1撮影あたりの撮影フレーム数、X線管電流の値、X線管電圧の値、フィルタ種等である。パルスレートは、1秒あたりの放射線照射回数であり、後述するフレームレートと一致している。パルス幅は、放射線照射1回当たりの放射線照射時間である。パルス間隔は、連続撮影において、1回の放射線照射開始から次の放射線照射開始までの時間であり、後述するフレーム間隔と一致している。
【0024】
放射線検出部13は、FPD等の半導体イメージセンサにより構成される。FPDは、例えば、ガラス基板等を有しており、基板上の所定位置に、放射線源11から照射されて少なくとも被写体Mを透過した放射線をその強度に応じて検出し、検出した放射線を電気信号に変換して蓄積する複数の画素(検出素子)がマトリックス状(二次元状)に配列されている。各画素は、例えばTFT(Thin Film Transistor)等のスイッチング部により構成されている。
【0025】
読取制御装置14は、撮影用コンソール2に接続されている。読取制御装置14は、撮影用コンソール2から入力された画像読取条件に基づいて放射線検出部13の各画素のスイッチング部を制御して、当該各画素に蓄積された電気信号の読み取りをスイッチングしていき、放射線検出部13に蓄積された電気信号を読み取ることにより、画像データを取得する。この画像データがフレーム画像である。そして、読取制御装置14は、取得したフレーム画像を撮影用コンソール2に出力する。画像読取条件は、例えば、フレームレート、フレーム間隔、画素サイズ、画像サイズ(マトリックスサイズ)等である。フレームレートは、1秒あたりに取得するフレーム画像数であり、パルスレートと一致している。フレーム間隔は、連続撮影において、1回のフレーム画像の取得動作開始から次のフレーム画像の取得動作開始までの時間であり、パルス間隔と一致している。
【0026】
ここで、放射線照射制御装置12と読取制御装置14は互いに接続され、互いに同期信号をやりとりして放射線照射動作と画像の読み取りの動作を同調させるようになっている。この他、後述するオフセット補正に用いるオフセット補正係数を算出するための複数の暗画像を取得するキャリブレーション時は、放射線照射動作と同期せず、放射線が照射されない状態で、リセット〜蓄積〜データ読取〜リセットの一連の画像の読み取り動作を行うが、一連の動態撮影前、一連の動態撮影後のいずれかのタイミングで行うようにしてもよい。
【0027】
〔撮影用コンソール2の構成〕
撮影用コンソール2は、放射線照射条件や画像読取条件を撮影装置1に出力して撮影装置1による放射線撮影及び放射線画像の読み取り動作を制御するとともに、撮影装置1により取得された動態画像を撮影技師によるポジショニングの確認や診断に適した画像であるか否かの確認用に表示する。
撮影用コンソール2は、
図1に示すように、制御部21、記憶部22、操作部23、表示部24、通信部25を備えて構成され、各部はバス26により接続されている。
【0028】
制御部21は、CPU(Central Processing Unit)、RAM(Random Access Memory
)等により構成される。制御部21のCPUは、操作部23の操作に応じて、記憶部22に記憶されているシステムプログラムや各種処理プログラムを読み出してRAM内に展開し、展開されたプログラムに従って後述する撮影制御処理を始めとする各種処理を実行し、撮影用コンソール2各部の動作や、撮影装置1の放射線照射動作及び読み取り動作を集中制御する。
【0029】
記憶部22は、不揮発性の半導体メモリやハードディスク等により構成される。記憶部22は、制御部21で実行される各種プログラムやプログラムにより処理の実行に必要なパラメータ、或いは処理結果等のデータを記憶する。例えば、記憶部22は、
図2に示す撮影制御処理を実行するための撮影制御処理プログラムを記憶している。また、記憶部22は、胸部の動態画像を撮影するための放射線照射条件及び画像読取条件を記憶している。各種プログラムは、読取可能なプログラムコードの形態で格納され、制御部21は、当該プログラムコードに従った動作を逐次実行する。
【0030】
操作部23は、カーソルキー、数字入力キー、及び各種機能キー等を備えたキーボードと、マウス等のポインティングデバイスを備えて構成され、キーボードに対するキー操作やマウス操作により入力された指示信号を制御部21に出力する。また、操作部23は、表示部24の表示画面にタッチパネルを備えても良く、この場合、タッチパネルを介して入力された指示信号を制御部21に出力する。
【0031】
表示部24は、LCD(Liquid Crystal Display)やCRT(Cathode Ray Tube)等のモニタにより構成され、制御部21から入力される表示信号の指示に従って、操作部23からの入力指示やデータ等を表示する。
【0032】
通信部25は、LANアダプタやモデムやTA(Terminal Adapter)等を備え、通信ネットワークNTに接続された各装置との間のデータ送受信を制御する。
【0033】
〔診断用コンソール3の構成〕
診断用コンソール3は、撮影用コンソール2から動態画像を取得し、取得した動態画像に基づいて画像解析を行い、解析結果を表示する動態画像処理装置である。
診断用コンソール3は、
図1に示すように、制御部31、記憶部32、操作部33、表示部34、通信部35を備えて構成され、各部はバス36により接続されている。
【0034】
制御部31は、CPU、RAM等により構成される。制御部31のCPUは、操作部33の操作に応じて、記憶部32に記憶されているシステムプログラムや、各種処理プログラムを読み出してRAM内に展開し、展開されたプログラムに従って、後述する画像解析処理を始めとする各種処理を実行し、診断用コンソール3各部の動作を集中制御する。制御部31は、後述する画像解析処理を実行することにより、抽出手段、気流速度算出手段、気流速度分布算出手段としての機能を実現する。
【0035】
記憶部32は、不揮発性の半導体メモリやハードディスク等により構成される。記憶部32は、制御部31で画像解析処理を実行するための画像解析処理プログラムを始めとする各種プログラムやプログラムにより処理の実行に必要なパラメータ、或いは処理結果等のデータを記憶する。これらの各種プログラムは、読取可能なプログラムコードの形態で格納され、制御部31は、当該プログラムコードに従った動作を逐次実行する。
【0036】
操作部33は、カーソルキー、数字入力キー、及び各種機能キー等を備えたキーボードと、マウス等のポインティングデバイスを備えて構成され、キーボードに対するキー操作やマウス操作により入力された指示信号を制御部31に出力する。また、操作部33は、表示部34の表示画面にタッチパネルを備えても良く、この場合、タッチパネルを介して入力された指示信号を制御部31に出力する。
【0037】
表示部34は、LCDやCRT等のモニタにより構成され、制御部31から入力される表示信号の指示に従って、操作部33からの入力指示やデータ等を表示する。
【0038】
通信部35は、LANアダプタやモデムやTA等を備え、通信ネットワークNTに接続された各装置との間のデータ送受信を制御する。
【0039】
〔胸部診断支援システム100の動作〕
次に、上記胸部診断支援システム100における動作について説明する。
【0040】
(撮影装置1、撮影用コンソール2の動作)
まず、撮影装置1、撮影用コンソール2による撮影動作について説明する。
図2に、撮影用コンソール2の制御部21において実行される撮影制御処理を示す。撮影制御処理は、制御部21と記憶部22に記憶されている撮影制御処理プログラムとの協働により実行される。
【0041】
まず、撮影技師により撮影用コンソール2の操作部23が操作され、撮影対象(被写体M)の患者情報(患者の氏名、身長、体重、年齢、性別等)の入力が行われる(ステップS1)。
【0042】
次いで、放射線照射条件が記憶部22から読み出されて放射線照射制御装置12に設定されるとともに、画像読取条件が記憶部22から読み出されて読取制御装置14に設定される(ステップS2)。
【0043】
次いで、操作部23の操作による放射線照射の指示が待機される(ステップS3)。ここで、撮影実施者は、撮影装置1において被写体Mのポジショニング等の撮影準備を行うとともに、安静呼吸の動態を撮影するために被験者(被写体M)に楽にするように指示し、安静呼吸を促す。撮影準備が整った時点で、操作部23を操作して放射線照射指示を入力する。
【0044】
操作部23により放射線照射指示が入力されると(ステップS3;YES)、放射線照射制御装置12及び読取制御装置14に撮影開始指示が出力され、動態撮影が開始される(ステップS4)。即ち、放射線照射制御装置12に設定されたパルス間隔で放射線源11により放射線が照射され、放射線検出部13によりフレーム画像が取得される。予め定められたフレーム数の撮影が終了すると、制御部21により放射線照射制御装置12及び読取制御装置14に撮影終了の指示が出力され、撮影動作が停止される。撮影されるフレーム数は、少なくとも1呼吸サイクルが撮影できる枚数である。
【0045】
撮影により取得されたフレーム画像は順次撮影用コンソール2に入力され、各フレーム画像に対して補正処理が行われる(ステップS5)。ステップS5の補正処理においては、オフセット補正処理、ゲイン補正処理、欠陥画素補正処理の3つの補正処理が行われる。まず最初に、取得された各フレーム画像に対してオフセット補正処理が行われ、取得された各フレーム画像に重畳された暗電流に起因するオフセット値が除去される。オフセット補正処理では、例えば、取得した各フレーム画像の各画素値(濃度値。以下、信号値という)から、予め記憶されたオフセット補正係数を減算する処理が行われる。ここで、オフセット補正係数は、予め放射線非照射時に取得した複数のフレーム画像を平均化した画像である。次いで、ゲイン補正処理が行われ、各フレーム画像の各画素に対応する各検出素子の個体差や読み出しアンプのゲインムラによって生じる画素毎のばらつきが除去される。ゲイン補正処理では、例えば、オフセット補正後の各フレーム画像に、予め記憶されたゲイン補正係数を乗算する処理が行われる。ここで、ゲイン補正係数は、放射線検出部13に一様に放射線を照射した時に取得した複数のオフセット補正済みフレーム画像を平均化した画像と、このときの放射線照射条件で期待される出力信号値の関係から、補正後の各画素の信号値が一様となるように予め算出され、記憶された係数である。次いで、欠陥画素補正処理が行われ、周囲の画素と比較して感度が非線形な画素や、感度がない欠落画素が除去される。欠陥画素補正処理では、例えば、予め記憶された欠陥画素位置情報マップに従って、欠陥画素位置情報マップに登録された各欠陥画素において、欠陥画素の信号値をその近傍の欠陥でない画素の信号値の平均値で置き換える処理が行われる。ここで、欠陥画素位置情報マップは、放射線検出部13に一様に放射線を照射した時に取得したオフセット補正、ゲイン補正済みのフレーム画像から、予め複数の欠陥画素が認識され、その欠陥画素の位置が登録されたマップである。上記オフセット補正係数及びゲイン補正係数、欠陥画素位置情報マップは、ビニングやダイナミックレンジ等の収集モードに応じて、予め、それぞれ最適な値が記憶されており、それぞれの収集モードにおいて対応する最適な値が読み出されるようになっている。
【0046】
次いで、補正処理後の各フレーム画像と、撮影順を示す番号と対応付けて記憶部22に記憶されるとともに(ステップS6)、表示部24に表示される(ステップS7)。ここで、各フレーム画像を記憶する直前に、各フレーム画像の各画素の信号値を真数から対数に変換する対数変換処理を行ってから記憶しても良い。撮影技師は、表示された動態画像によりポジショニング等を確認し、撮影により診断に適した画像が取得された(撮影OK)か、再撮影が必要(撮影NG)か、を判断する。そして、操作部23を操作して、判断結果を入力する。尚、撮影により取得された各フレーム画像は、全撮影の終了後に纏めて入力するようにしても良い。
【0047】
操作部23の所定の操作により撮影OKを示す判断結果が入力されると(ステップS8;YES)、動態撮影で取得された一連のフレーム画像のそれぞれに、動態画像を識別するための識別IDや、患者情報、放射線照射条件、画像読取条件、撮影順を示す番号、撮影日時等の情報が付帯され(例えば、DICOM形式で画像データのヘッダ領域に書き込まれ)、通信部25を介して診断用コンソール3に送信される(ステップS9)。そして、本処理は終了する。一方、操作部23の所定の操作により撮影NGを示す判断結果が入力されると(ステップS8;NO)、記憶部22に記憶された一連のフレーム画像が削除され(ステップS10)、本処理は終了する。尚、このケースに於いては、再撮影が実行されることとなる。
【0048】
(診断用コンソール3の動作)
次に、診断用コンソール3における動作について説明する。
診断用コンソール3においては、通信部35を介して撮影用コンソール2から動態画像の一連のフレーム画像が受信されると、制御部31と記憶部32に記憶されている画像解析処理プログラムとの協働により
図3に示す画像解析処理が実行される。
【0049】
以下、
図3を参照して画像解析処理の流れについて説明する。
まず、基準となる一のフレーム画像(基準画像という)から肺野領域が抽出される(ステップS21)。
基準画像としては、安静呼気位のフレーム画像とすることが好ましい。安静呼気位では、安静呼吸時において肺野領域の面積が最も小さくなるので、基準画像の各小領域を他のフレーム画像に対応付けたときに、小領域が他のフレーム画像の肺野外の領域に対応付けられることがないためである。
安静呼気位のフレーム画像は、一連のフレーム画像の中から横隔膜の位置が最も高い位置にある画像を抽出することで取得することができる。また、まず各フレーム画像から肺野領域を抽出し、肺野領域の面積が最小のフレーム画像(肺野領域内の画素数が最も少ない画像)を基準画像とすることとしてもよい。
肺野領域の抽出方法は何れの方法であってもよい。例えば、フレーム画像(基準画像)の各画素の信号値(濃度値)のヒストグラムから判別分析によって閾値を求め、この閾値より高信号の領域を肺野領域候補として1次抽出する。次いで、1次抽出された肺野領域候補の境界付近でエッジ検出を行い、境界付近の小領域でエッジが最大となる点を境界に沿って抽出すれば肺野領域の境界を抽出することができる。
【0050】
次いで、時間軸方向のローパスフィルタリング処理が施される(ステップS22)。肺野領域には、換気と血流の双方の動きによる信号値の変化が現れる。呼吸換気に伴う信号値の時間による変化幅は、肺血流に伴う信号値変化幅の約10倍となる。よって、換気の解析をする場合は、血流等に伴う高周波数な信号変化はノイズとなる。この血流等に伴う信号変化を除去するために、信号値の時間変化に対して、カットオフ周波数1.0Hzでローパスフィルタ処理が行われる。
【0051】
次いで、各フレーム画像の肺野領域が複数の小領域に分割され、各フレーム画像の小領域が互いに対応付けられる(ステップS23)。各小領域の画素の位置は、制御部31のRAMに記憶される。
【0052】
ここで、呼吸サイクルは、呼気期と吸気期により構成される。
図4は、1つの呼吸サイクル(深呼吸時)において撮影された複数の時間位相T(T=t0〜t6)のフレーム画像を示す図である。
図4に示すように、呼気期は、横隔膜が上がることによって肺から空気が排出され、肺野の領域が小さくなる。最大呼気位では、横隔膜の位置が最も高い状態となる。吸気期は、横隔膜が下がることにより肺に空気が取り込まれ、
図4に示すように胸郭中の肺野の領域が大きくなる。最大吸気位では、横隔膜の位置が最も下がった状態となる。即ち、肺野領域の同一部分の位置は呼吸運動に応じて時間とともに変化するため、各フレーム画像間で肺野の同一の部分(特に、下部領域(横隔膜付近))を示す画素位置がずれることになる。
しかし、安静呼吸時に撮影された画像においては、上記の位置ずれは小さく、後述する解析結果が狂ってしまうほどの位置ずれはおきない。
図5の画像D1は、安静呼気位(安静呼吸時において横隔膜の位置が最も高くなったタイミング)のフレーム画像である。
図5の画像D2は、安静吸気位(安静呼吸時において横隔膜の位置が最も低くなったタイミング)のフレーム画像である。即ち、
図5の画像D1とD2とは、呼吸1サイクルにおいて最も形状の差の大きいタイミングで撮影された画像である。しかし、
図5の画像D1、D2間においては、最も位置ずれの大きい肺野領域の下部領域においても位置ずれはわずかであることがわかる(画像D2のA11は画像D1のA1と同じ画素位置を示し、画像D2のA2は画像D1のA1と肺野における同一部分を描画した領域を示している)。
【0053】
そこで、ステップS23における具体的な処理としては、まず、基準画像の抽出された肺野領域を複数の小領域(例えば、2mm×2mmの矩形領域)に分割する(
図5参照)。次いで、他のフレーム画像の肺野領域を基準画像の各小領域と同じ画素位置の小領域(放射線検出部13の同じ検出素子から出力される信号値の領域)に分割する。次いで、各フレーム画像間の同じ画素位置の各小領域を互いに対応付ける。各フレーム画像の肺野領域は基準画像と同じ画素位置の領域(複数の小領域からなる領域)となる。この処理では、フレーム画像の小領域への分割及び対応付けを高速に行うことが可能となる。
【0054】
なお、ここでは安静呼吸時の動態画像を使用しているが、撮影された動態画像が深呼吸時の画像である場合は、
図6に示すように、肺野の同一の部分を示す画素位置が大幅にずれる(画像D4のA31は画像D3のA3と同じ画素位置を示し、画像D4のA4は画像D3のA3と肺野における同一部分を描画した領域を示している)。そのため、安静呼吸時と同様に各フレーム画像における、基準画像の各小領域と同じ画素位置の領域をその小領域に対応する領域としてしまうと、後述する解析結果が診断に利用し得ないものとなる。そこで、このような場合には、各フレーム画像間の対応点を抽出する対応点抽出処理(ローカルマッチング処理)及び非線形歪変換処理(ワーピング処理)を行って、肺野領域の同一部分が描画された領域を各フレーム画像間で対応付ける。また、安静呼吸時に於いて、処理速度をあまり重視せず、解析精度をより高めたい場合にも、これらの処理を行うこととしてもよい。
【0055】
ローカルマッチング処理では、まず、撮影順が1番(最初)のフレーム画像から肺野領域を抽出し、この抽出した肺野領域を、例えば、2mm角の矩形からなる小領域に分割する。
次いで、撮影順が1番のフレーム画像をP1、これと隣接するフレーム画像(撮影順が隣接するフレーム画像(即ち、時間的に隣接するフレーム画像。以下同様。))をP2とし、P2に、P1の各小領域の探索領域を設定する。ここで、P2の探索領域は、P1における各小領域における中心点の座標を(x,y)とすると、同一の中心点(x,y)をもち、P1の小領域よりも縦横の幅が大きくなるように設定する(例えば、1.5倍)。そして、P1の各領域毎に、P2の探索範囲で最もマッチング度合いが高くなる領域を求めることで、P1の各小領域に対するP2上での対応位置を算出する。マッチング度合いとしては、最小二乗法や相互相関係数を指標に用いる。そして、P2の肺野領域をP1の各小領域の対応位置で分割する。
次いで、P2を、新たにP1とみなし、撮影順がP2の次のフレーム画像を新たなP2とみなして、P1の各小領域におけるP2の対応位置を算出する。以上の処理を繰り返すことで、各フレーム画像の各小領域が隣接するフレーム画像のどの位置に対応するかが求まる。求めた処理結果は、制御部31のRAMに記憶される。
【0056】
次いで、ワーピング処理が行われる。具体的には、撮影順が1番のフレーム画像をP1、これと撮影順が隣接する(時間的に隣接する)フレーム画像をP2とし、上記ローカルマッチング処理で算出された隣接するフレーム画像間の各小領域の対応位置に基づいて、各小領域毎にP1からP2へのシフトベクトルを算出する。次いで、算出されたシフトベクトルを多項式でフィッティングして、この多項式を用いて各小領域における各画素のシフトベクトルを算出する。そして、算出された各画素のシフトベクトルに基づいて、ワーピング処理を行い、P2の各小領域内の各画素の位置をP1のフレーム画像の対応する画素の位置にシフトする。次いで、ワーピング処理されたP2を、新たにP1とみなし、撮影順がP2の次のフレーム画像を新たなP2とみなして、上記処理を行う。以上の処理を撮影順の早い隣接フレーム画像間から順次繰り返すことで、全てのフレーム画像の各小領域の位置を撮影順が1のフレーム画像(この処理での基準画像)に略一致させることが可能となる。各フレーム画像間における小領域の位置の対応関係は、制御部31のRAMに記憶される。
【0057】
次いで、一連のフレーム画像のステップS23で対応付けられた各小領域毎に、フレーム間差分処理が実行され、フレーム間差分値が算出される(ステップS24)。
フレーム間差分値は、そのフレーム画像が撮影されたタイミングにおける信号変化量を示す値である。呼吸により息を吸ったり吐いたりすれば、その息の流れに応じて肺の密度が変化し、これによってX線透過量(つまり、画素の出力信号値)が変化する。よって、信号変化量は、そのタイミングにおける気流速度を示す値(気流速度を表す特徴量)とみなすことができる。
具体的には、まず、各フレーム画像の各小領域内の画素の信号値(平均信号値)が算出される。次いで、撮影順が隣接するフレーム画像間で各小領域の信号値の差分を算出するフレーム間差分処理が行われる。ここでは、各小領域毎に、フレーム番号NとN+1(Nは1、2、3・・・)のフレーム画像について、N+1−Nの差分値が算出される。
【0058】
次いで、吸気期のフレーム画像群と呼気期のフレーム画像群への分類が行われる(ステップS25)。各小領域においてフレーム間差分値の符号が正の期間が吸気期、負の期間が呼気期である。
【0059】
次いで、各小領域において、呼気期のフレーム間差分値の最大値(絶対値の最大値)が呼気最大気流速度として取得され、各小領域の値が呼気最大気流速度からなる呼気最大気流速度画像が作成される。また、各小領域において、吸気期のフレーム間差分値の最大値(絶対値の最大値)が吸気最大気流速度として取得され、各小領域の値が吸気最大気流速度からなる呼気最大気流速度画像が作成される(ステップS26)。
【0060】
次いで、呼気最大気流速度画像、吸気最大気流速度画像のそれぞれにおいて、肺野領域が体幹軸方向(頭−足方向、上下方向)の複数のブロック領域に分割され、各ブロック領域で最大気流速度の代表値(ここでは、平均値)が算出される(ステップS27)。各ブロック領域の最大気流速度の代表値は、各ブロック領域の気流速度を表す特徴量である。
肺野領域は、例えば、
図7Aに示すように、左肺野、右肺野のそれぞれが上、中、下の3つのブロック領域L1〜L3、R1〜R3に分割される。
また、ここでは各ブロック領域における最大気流速度の代表値として平均値を用いているが、平均値に限定せず、最小値、最大値、中央値、積分値を用いることとしてもよい。
【0061】
そして、呼気最大気流速度画像、吸気最大気流速度画像のそれぞれにおいて、体幹軸方向に隣接するブロック領域間で最大気流速度平均値の差分値(下のブロック領域の最大気流速度平均値−上のブロック領域の最大気流速度平均値)が算出され、算出された差分値に基づいて、肺野における体幹軸方向における最大気流速度分布を表す特徴量(最大気流速度分布特徴量と呼ぶ)が算出される(ステップS28)。
最大気流速度分布特徴量は、呼気と吸気のそれぞれについて、左右の肺野毎に、下記の[数1]により算出される。また、左右の肺野毎だけでなく、肺野全体を対象として下記の[数1]により最大気流速度分布特徴量を算出することとしてもよい。
【数1】
特徴量の算出後、解析結果(最大気流速度分布特徴量を含む)が表示部34に表示される(ステップS29)。
【0062】
正常な肺野では、呼気と吸気共に、上肺野に向かうにつれて徐々に気流速度が低下する。一方、肺気腫等の換気疾患のある肺野では、疾患の重症度が大きくなるにつれて上記傾向が失われていく。即ち、正常な肺野では、体幹軸方向に隣接するブロック領域間の気流速度の差分値(下のブロック領域の最大気流速度平均値−上のブロック領域の最大気流速度平均値)の値がほぼすべて0より大きくなり、算出された最大気流速度分布特徴量の値は約0.9以上となる。異常のある肺野では、上記隣接ブロック領域間の気流速度の差分値の値に負の値が混じるので、算出された最大気流速度分布特徴量の値は0.9より小さくなる。異常が重症であるほど算出された特徴量の値は小さくなる。
このように、肺野の体幹軸方向における最大気流速度分布を表す最大気流速度分布特徴量によって、吸気気流速度と呼気気流速度の間に差がみられない換気疾患であっても医師が容易に把握することが可能となり、診断の判断材料とすることができる。
【0063】
図8A〜
図8Cに、ステップS29における解析結果の表示例を示す。解析結果としては、呼気と吸気それぞれの最大気流速度分布特徴量のほか、最大気流速度画像を表示することが好ましい。最大気流速度画像を表示することで、実際に各肺野の最大気流速度の分布がどのようになっているのかを医師が容易に(直感的に)把握することができるからである。
最大気流速度画像は、
図8Aに示すように、フレーム画像(基準画像)上の各小領域を最大気流速度に応じた彩度(又は輝度)で色分けして示した画像としてもよいし、
図8B に示すように、フレーム画像(基準画像)上の各ブロック領域を最大気流速度の平均値に応じた彩度(又は輝度)で色分けして示した画像としてもよい。最大気流速度画像の近傍には、最大気流速度と彩度(又は輝度)との対応関係を表示することが好ましい。また、
図8Cに示すように、最大気流速度画像においては、各ブロック領域の最大気流速度の平均値を数値で表示することとしてもよい。
図8Cにおいて、a〜fはそれぞれブロック領域[1]〜[6]の最大気流速度の平均値を表す数値を示している。
なお、呼気最大気流速度画像と吸気最大気流速度画像の双方を表示することが好ましいが、何れか結果の悪い方のみを表示することとしてもよい。
【0064】
また、ステップS29においては、各小領域毎の気流速度比(吸気の最大気流速度/呼気の最大気流速度)のヒストグラム、及び肺野全体における気流速度比の傾向を示す指標値(平均値、標準偏差、半値幅等)を算出(作成)して併せて表示することが好ましい(
図11〜
図13の34d参照)。最大気流速度分布特徴量は、吸気気流速度と呼気気流速度の間に差がみられない換気疾患を判別可能とするものであるが、気流速度比のヒストグラムや指標値を表示することにより、吸気気流速度と呼気気流速度との間の差が大きい換気疾患についての診断支援情報も提供することができるようになるからである。また、ヒストグラム等と併せて、フレーム画像(基準画像)上の各小領域を気流速度比に応じた彩度(又は輝度)で色分けして示した画像を併せて表示することとしてもよい(
図11〜
図13の34a、34b参照)。このようにすれば、肺野内の気流速度比が大きい局所的な異常個所を医師が容易に把握することが可能となるからである。
【0065】
また、ステップS29においては、併せて、オリジナル動画像と換気動画像のいずれか又は双方を表示することとしてもよい。これらを表示することにより、解析結果と併せて実際の肺の動きを医師が観察可能となるからである。オリジナル動画像は、撮影装置1において撮影された動態画像のフレーム画像を順次切り替え表示(パラパラ捲り表示)したものである。換気動画像は、ステップS22でローパスフィルタにかけたフレーム画像を順次切り替え表示したものである。
【0066】
なお、上記第1の実施の形態においては、右肺野、左肺野のそれぞれを体幹軸方向に対して3分割することでブロック領域に分割することとして説明したが、
図7Bに示すように、3分割以上することとしてもよい。ブロック領域を細かくすることで、最大気流速度特徴量の精度を向上することができる。また、正常な気流速度の傾向として、内側と外側でも気流速度の傾向に違いがあるので、内側から外側への気流速度の変化を考慮に入れて、
図7Cに示すように、臨床に即したブロック領域で分割することとしてもよい。これにより、最大気流速度特徴量の算出精度を向上させることができる。
【0067】
−第2の実施の形態−
次に、第2の実施の形態について説明する。
第2の実施の形態においては、診断用コンソール3の記憶部32に正常肺野の呼気、吸気それぞれの最大気流速度分布のテンプレート(正常気流速度分布テンプレート)が記憶されている点、及び診断用コンソール3の制御部31において実行される画像解析処理の内容が第1の実施の形態と異なる。その他、胸部診断支援システム100の全体構成、各装置の構成、及び撮影装置1と撮影用コンソール2の動作については、第1の実施の形態において説明したものと同様であるので説明を援用する。
【0068】
以下、第2の実施の形態における画像解析処理(画像解析処理Bとする)について説明する。
図9に、第2の実施の形態における画像解析処理Bのフローチャートを示す。画像解析処理Bは、制御部31と記憶部32に記憶されている画像解析処理Bプログラムとの協働により実行される。
【0069】
まず、ステップS31〜ステップS37の処理が実行されることにより、一連のフレーム画像から呼気期及び吸気期それぞれの最大気流速度画像が作成され、それぞれ体幹軸方向の複数のブロック領域に分割され、各ブロック領域の最大気流速度の代表値(ここでは、平均値)が算出される。ステップS31〜ステップS37の処理は、ステップS21〜S27の処理と同様であるので説明を援用する。
【0070】
次いで、記憶部32から呼気、吸気それぞれの、正常気流速度分布テンプレートが読み出される(ステップS38)。
ここで、換気機能が正常な肺野の気流速度には「下肺野→上肺野」にかけてなだらかに低下していくという傾向が見られる。そこで、正常気流速度分布テンプレートは、この傾向を満たすように、各ブロック領域に最大気流速度の平均値を割り当てた画像を用いることができる。例えば、
図7Aにおける領域L1(R1)、L2(R2)、L3(R3)に対し、L1(R1)=3、L2(R2)=2、L3(R3)=1の値を割り当てることができる。
或いは、複数人の予め正常ということがわかっている複数の肺野の動態画像に対し、上述のステップS31〜S37と同様の処理を行って、各ブロック領域毎の最大気流速度の平均値を算出し、これを平均化したものを正常気流速度分布テンプレートとしてもよい。
【0071】
なお、最大気流速度画像と正常気流速度分布テンプレートのブロック領域は互いに対応している。また、正常気流速度分布テンプレートは、下記の条件に応じて記憶部32に複数記憶しておき、操作部33からの入力や患者情報に応じた正常気流速度分布テンプレートを読み出して使用することとしてもよい。
条件:性別・年齢、胴囲・胸囲、身長・体重、肺活量、手術歴(ペースメーカー有無、肺切除有無等)、撮影時呼吸法(安静呼吸、深呼吸(努力呼吸))、過去の病歴等
【0072】
次いで、呼気の右肺野、呼気の左肺野、吸気の右肺野、吸気の左肺野のそれぞれについて、それぞれの肺野内の最大気流速度画像と正常気流速度分布テンプレートの値を用いて、最大気流速度画像と正常気流速度分布テンプレートとの相互相関値(相互相関係数)が算出される(ステップS39)。各肺野の相互相関値は、各肺野における最大気流速度分布と正常な気流速度分布との相関を示しており、各肺野における最大気流速度の分布を表す特徴量(最大気流速度分布特徴量)となる。
【0073】
ステップS39においては、呼気右肺野の相互相関値、呼気左肺野の相互相関値、吸気右肺野の相互相関値、吸気左肺野の相互相関値の4つの相互相関値とその積が算出される。
例えば、呼気右肺野の相互相関値は、呼気最大気流速度画像の右肺野内の各ブロック領域の最大気流速度の平均値と、呼気の正常気流速度分布テンプレートの右肺野内の各ブロック領域の平均値とを用いて、下記の [数2]によって求めることができる(例えば、公知文献1:高木幹雄、下田陽久著、「新編画像解析ハンドブック」、東京大学出版会、20004年参照)。
【数2】
ここで、j=1、2・・・は下肺野から順に付与されるものとする。
【0074】
即ち、呼気右肺野の最大気流速度の相互相関値は、呼気の最大気流速度画像の右肺野内の各ブロック領域の最大気流速度の平均値と、呼気の正常気流速度分布テンプレートの右肺野内の各ブロック領域の平均値との共分散を算出し、算出された共分散を、呼気の最大気流速度画像の右肺野内の各ブロック領域の呼気最大気流速度の平均値の標準偏差と呼気正常気流速度分布テンプレートの右肺野内の各ブロック領域の平均値の標準偏差の積で除算することにより求めることができる。
呼気左肺野の相互相関値、吸気右肺野の相互相関値、吸気左肺野の相互相関値についても、呼気右肺野の部分をそれぞれ呼気左肺野、吸気右肺野、吸気左肺野に置き換え、用いる最大気流速度画像及び正常気流速度分布テンプレートを対応するものに変更するだけで、算出方法は同じである。
なお、特徴量としては、各肺野毎に相互相関値を用いるのではなく、上述の公知文献1に記載のSSDA法による類似度を求めることとしてもよい。
【0075】
正常な肺野では、呼気と吸気共に、上肺野に向かうにつれて徐々に気流速度が低下する。一方、肺気腫等の換気疾患のある肺野では、疾患の重症度が大きくなるにつれて上記傾向が失われていく。よって、各肺野毎に、最大気流速度画像と正常気流速度分布テンプレートとの相互相関値を求めることにより、正常の肺野では相互相関値が1に近くなり、異常が重症であるほど1から遠ざかる値となる。各肺野の積は、異常が重症な肺野ほど1から遠ざかる値となる。このように、肺野の体幹軸方向における最大気流速度分布を表す最大気流速度分布特徴量を表示することで、吸気気流速度と呼気気流速度の間に差がみられない換気疾患であっても医師が容易に把握することが可能となる。
【0076】
特徴量の算出後、解析結果(最大気流速度分布特徴量を含む)が表示部34に表示される(ステップS40)。
ステップS40における解析結果は、第1の実施の形態の解析結果と同様に、
図8A 〜
図8Cに示す表示とすることができる。即ち、最大気流速度分布特徴量のほか、最大気流速度画像を表示することが好ましい。また、各小領域毎の気流速度比(吸気の最大気流速度/呼気の最大気流速度)のヒストグラム、及び肺野全体における気流速度比の傾向を示す指標値(平均値、標準偏差、半値幅等)を算出(作成)して併せて表示することが好ましい(
図11〜
図13の34d参照)。また、ヒストグラム等と併せて、フレーム画像(基準画像)上の各小領域を気流速度比に応じた彩度(又は輝度)で色分けして示した画像を併せて表示することとしてもよい(
図11〜
図13の34a、34b参照)。また、オリジナル動画像と換気動画像のいずれか又は双方を表示することとしてもよい。
【0077】
更に、ステップS40においては、解析した正常気流速度分布テンプレートを併せて表示することとしてもよい。例えば、肺野を正常気流速度分布テンプレートの値に応じた彩度(又は輝度)で色分けして示した画像として表示してもよいし、正常気流速度分布テンプレートの各ブロック領域の値を数値で表示することとしてもよい。
なお、呼気と吸気の双方の正常気流速度分布テンプレートを表示することが好ましいが、何れか結果の悪い方のみを表示することとしてもよい。
【0078】
また、ステップS40においては、最大気流速度画像と正常気流速度分布テンプレートの対応するブロック領域毎の「差の絶対値」又は「差の2乗値」を算出して表示することとしてもよい。これにより、どのブロック領域が異常であるかを医師が認識することができる。「差の絶対値」又は「差の2乗値」についても、呼気と吸気の双方を表示することが好ましいが、何れか結果の悪い方のみを表示することとしてもよい。また、差分を取る前に、最大気流速度画像と正常気流速度分布テンプレート画像のレベルを合わせるために、それぞれを平均値、最大値、最小値、又は中央値で除算することで正規化してもよい。
【0079】
なお、ブロック領域への分割については、上述のように、各肺野を上、中、下の3分割にする他、第1の実施の形態と同様、
図7Bに示すように、体幹軸方向に3分割以上することとしてもよい。また、
図7Cに示すように、臨床に即したブロック領域で分割することとしてもよい。更に、
図10A、
図10Bに示すように、肺野領域を体幹軸方向及びこれと略直交する方向の複数のブロック領域に分割することとしてもよい。これにより、肺野の体幹軸方向の気流速度の分布に加えて肺野の内側から外側にかけての気流速度の分布を考慮して特徴量を算出することができる。
【0080】
ここで、上記画像解析処理Bにより算出される最大気流速度分布特徴量(相互相関値)の有効性についての検証結果について説明する。
表1、表2に、
図11に示す正常な肺野について画像解析処理Bで算出した相互相関値(上段)と、
図12に示す肺気腫の肺野(GOLD3)について画像解析処理Bで算出した相互相関値(下段)と、を示す。
図12に示すように、この肺気腫の肺野は、気流速度比を用いた解析では正常のような結果を示すものである。
表1は、各肺野を
図7Aに示すように上、中、下の3つのブロック領域に分割した場合を示す。表2は、各肺野を
図7Bに示すように20のブロック領域に分割した場合を示す。また、いずれの場合も、換気機能が正常な肺野の気流速度の傾向に基づいて作成したテンプレートを用いて相互相関値を算出している。なお、相互相関値が0未満である場合は、全て0とする。
【表1】
【表2】
【0081】
表1に示すように、左右の各肺野を3分割して正常気流速度分布テンプレートとの相互相関値を算出した場合、上段の正常な肺野の場合は、呼気右肺野0.99、呼気左肺野0.98、吸気右肺野0.99、吸気左肺野0.99であり、全ての相互相関値が1に近い値となった。下段の肺気腫の肺野の場合は、呼気右肺野0.92、呼気左肺野0.98、吸気右肺野0.98、吸気左肺野0.87となった。また、正常な肺野、肺気腫の肺野のそれぞれについて、各相互相関値の積をとると、正常な肺野は0.96、肺気腫の肺野は0.78となり、正常な肺野と肺気腫の肺野の値の差が際立つ結果となった。
表2に示すように、左右の各肺野を20分割して正常気流速度分布テンプレートとの相互相関値を算出した場合、上段の正常な肺野の場合は、呼気右肺野0.95、呼気左肺野0.96、吸気右肺野0.96、吸気左肺野0.98であり、全ての相互相関値が1に近い値となった。下段の肺気腫の肺野の場合は、呼気右肺野0.58、呼気左肺野0.55、吸気右肺野0.65、吸気左肺野0.58となり、0.6付近の値となった。また、正常な肺野、肺気腫の肺野のそれぞれについて、各相互相関値の積をとると、正常な肺野は0.87、肺気腫の肺野は0.12となり、正常な肺野と肺気腫の肺野の値の差が際立つ結果となった。
本願発明者等の検討によれば、各相互相関値の積が左右の各肺野を3分割程度で0.9以上、20分割程度で0.8以上あれば正常な肺野の目安となり、換気重症度が進む程積が0に近くなる。上記表1、表2においては、各相互相関値の積でみると、双方とも正常な肺野は目安以上、肺気腫の肺野は目安未満の値が得られており、診断に有効であるといえる。従って、最大気流速度分布特徴量である相互相関値を生成して出力することによって、解析対象の肺野が正常か否かを判断するための判断材料を提供することが可能となる。なお、左右の肺野を3分割した表1の結果と、20分割した表2の結果とを比べると、20分割したほうが正常と肺気腫の肺野の相互相関値が離れた値となっており、異常であることが明確に表された結果となった。
【0082】
また、最大気流速度分布特徴量の有効性を評価するため、予めスパイロ検査結果によりCOPDの気流制限の度合いが分類されたサンプル(正常:29名、GOLD1又は2:20名、GOLD3又は4:31名(GOLD3:23名、GOLD4:8名))について上記画像解析処理Bにより最大気流速度分布特徴量(相互相関値)を算出し、最大気流速度分布特徴量の値が各グループ間で統計的有意差が認められるかをWilcoxon検定を用いて評価した。Wilcoxon検定はデータの有効性を評価するための公知の統計手法である(公知文献2:市原清志著「バイオサイエンスの統計学」、南江堂、1990年)。Wilcoxon検定では、値が0.05未満で有意差があるとされる。
【0083】
表3に、左右の各肺野を
図7Aに示すように上、中、下の3つのブロック領域に分割し、換気機能が正常な肺野の気流速度の傾向に基づいて作成した正常気流速度分布テンプレートを用いて算出した相互相関値についてのWilcoxon検定結果を示す。この場合、正常な肺野の相互相関値とGOLD3又はGOLD4の肺野の相互相関値との間に有意差が見られた。
【表3】
表4に、各肺野を
図7Bに示すように20のブロック領域に分割し、換気機能が正常な肺野の気流速度の傾向に基づいて作成した正常気流速度分布テンプレートを用いて算出した相互相関値についてのWilcoxon検定結果を示す。この場合、正常とGOLD3、正常とGOLD4、GOLD1又は2とGOLD4、GOLD3とGOLD4との間に有意差が見られた。
【表4】
【0084】
以上説明したように、胸部診断支援システム100によれば、診断用コンソール3の制御部31は、通信部35により撮影用コンソール2から胸部の動態を撮影することにより生成された複数のフレーム画像が入力されると、複数のフレーム画像のうちの一のフレーム画像から肺野領域を抽出し、抽出された肺野領域を複数の小領域に分割して複数のフレーム画像間にわたり小領域を対応付け、小領域毎に気流速度を表す特徴量を算出する。そして、肺野領域を体幹軸方向の複数のブロック領域に分割し、当該分割された各ブロック領域に含まれる複数の小領域の気流速度を表す特徴量に基づいて、各ブロック領域の気流速度を表す特徴量を算出し、当該算出した特徴量に基づいて、肺野領域の体幹軸方向の気流速度の分布を表す特徴量を算出する。
【0085】
正常な肺野であれば下肺野から上肺野にいくに従って気流速度が低下する。一方、肺気腫等の換気疾患のある肺野では、疾患の重症度が大きくなるにつれて上記傾向が失われていく。従って、肺野領域の体幹軸方向の気流速度の分布を表す特徴量を算出することで、吸気気流速度と呼気気流速度の間に差がみられない換気疾患であっても医師が容易に把握することができるような診断支援情報を提供することが可能となる。
【0086】
また、左右の肺野領域のそれぞれを体幹軸方向の複数のブロック領域に分割し、左右の肺野領域のそれぞれについて体幹軸方向の気流速度の分布を表す特徴量を生成することで、左右それぞれの肺野領域における体幹軸方向の気流速度の分布の傾向を医師が容易に把握できるような診断支援情報を提供することが可能となる。
【0087】
また、肺野領域を体幹軸方向の3以上のブロック領域に分割して肺野領域の体幹軸方向の気流速度の分布を表す特徴量を算出することで、体幹軸方向の気流速度の分布が肺野の上部に向かって単調減少しているか否かを精度よく表す診断支援情報を提供することが可能となる。
【0088】
また、肺野領域を前記体幹軸方向及びこれと略直交する方向の複数のブロック領域に分割することで、肺野の上下方向の気流速度の分布に加えて肺野の内側から外側にかけての気流速度の分布をも考慮した診断支援情報を提供することが可能となる。
【0089】
肺野領域の体幹軸方向の気流速度の分布を表す特徴量としては、例えば、肺野領域の気流速度と、正常な肺野の気流速度分布を示すテンプレートとの相互相関値を算出することができる。これにより、肺野領域の気流速度分布が正常との相関が強いか否かを数値で医師に提供することが可能となる。
【0090】
なお、上述した本実施の形態における記述は、本発明に係る好適な胸部診断支援システムの一例であり、これに限定されるものではない。
【0091】
例えば、上記の説明では、本発明に係るプログラムのコンピュータ読み取り可能な媒体としてハードディスクや半導体の不揮発性メモリ等を使用した例を開示したが、この例に限定されない。その他のコンピュータ読み取り可能な媒体として、CD-ROM等の可搬
型記録媒体を適用することが可能である。また、本発明に係るプログラムのデータを通信回線を介して提供する媒体として、キャリアウエーブ(搬送波)も適用される。
【0092】
その他、胸部診断支援システム100を構成する各装置の細部構成及び細部動作に関しても、本発明の趣旨を逸脱することのない範囲で適宜変更可能である。