説明

血流動態解析装置およびその制御プログラム

【課題】血管の血流動態に関する情報を正確かつ安定に算出する。
【解決手段】造影濃度曲線算出部2bは、時系列の撮影データから組織の時間濃度曲線Cを算出する。動脈の時間濃度曲線算出部2cは、組織の時間濃度曲線Cから動脈の時間濃度曲線aiを算出する。動脈の時間濃度曲線のカーブフィッティング部2dは、ガンマ確率密度関数の定数倍の関数をカーブモデルとして適用し、動脈の時間濃度曲線aのカーブフィッティングを行う。組織の時間濃度曲線のカーブフィッティング部2eは、不完全ガンマ関数の差をカーブモデルとして適用し、かつ、動脈の時間濃度曲線のカーブフィッティングで得られたパラメータの値の全てまたは一部を用いて、組織の時間濃度曲線Cのカーブフィッティングを行う。組織の時間濃度曲線のカーブフィッティング部2eは、カーブフィッティングにより得たパラメータの値を用いて、局所血流動態に関する情報を算出する。

【発明の詳細な説明】
【技術分野】
【0001】
本発明は、血流動態解析装置およびその制御プログラムに関し、特に、血管の血流動態に関する情報を正確かつ安定に算出することができる血流動態解析装置およびその制御プログラムに関する。
【背景技術】
【0002】
従来、X線CT(Computed Tomography)装置やMRI(Magnetic Resonance Imaging)装置では、静脈から造影剤を投与し、時系列画像データを収集して、その画像を元に解析し、組織内の血流動態に関する情報を得る検査が行われている。これはパーヒュージョン検査と呼ばれ、造影剤の撮影断面への集中の度合いが画像の濃度変化として収集できることを利用したものである。
【0003】
例えば、肺循環や造影剤投与のばらつきをなくすために、組織に流入する動脈の時間濃度曲線(Time Intensity Cur:TIC):aを入力関数として、測定した組織の時間濃度曲線Cとのデコンボリューション(Deconvolution:逆畳み込み積分)を行い、組織固有のインパルス応答関数(伝達関数)R(residue function)を求め、得られたインパルス応答関数から、血流動態を定量的に表す指標である血流量(脳の場合、CBF:Cerebral Blood Flow)、平均通過時間:MTT(Mean Transit Time)、および、血液量(脳の場合、CBV:Cerebral Blood Volume)などのパラメータを算出している。この解析方法は、デコンボリューション法と呼ばれている。
【0004】
また、非特許文献1では、動脈の時間濃度曲線をガンマ関数でフィッティングし、組織の時間濃度曲線をガンマ関数でフィッティングし、ガンマ関数でそれぞれフィッティングされた動脈の時間濃度曲線と組織の時間濃度曲線とのデコンボリューションを行い、得られたインパルス応答関数から、血流量、平均通過時間、および血液量を求めるようにしている。
【先行技術文献】
【非特許文献】
【0005】
【非特許文献1】T.Hackander et al.,The Use of gamma functions in the calculation of organ perfusion functions for non-diffusible radioactive tracers,phys.Med.Biol.,1988,vol.33,no.1,pp53-59,UK
【発明の概要】
【発明が解決しようとする課題】
【0006】
従来のデコンボリューション法では、組織の時間濃度曲線のデコンボリューションが数値的に不安定であり、インパルス応答関数の算出精度がある一定以上にはならない。また、インパルス応答関数のピークの面積は正確に求めることができるものの、ピークの幅や高さを独立に正確に求めることができない。つまり、デコンボリューション法で求めたい値は、ピークの高さ(CBFに相当)またはピークの幅(MTTに相当)であり、十分な精度が得られないという課題があった。
【0007】
また、非特許文献1の技術においても、デコンボリューションを行っているため、数値的に不安定であり、また求めるパラメータの数が多いため、正確性に欠ける課題があった。
【0008】
本発明はこのような状況に鑑みてなされたものであり、その目的は、血管の血流動態に関する情報を正確かつ安定に算出することが可能な血流動態解析装置およびその制御プログラムを提供することである。
【課題を解決するための手段】
【0009】
請求項1記載の本発明の特徴は、造影剤が注入された被検体を医用画像撮影装置で撮影した時系列の撮影データから、各画素の画素値の変化を計測して組織の時間濃度曲線を算出する組織の時間濃度曲線算出手段と、組織の時間濃度曲線算出手段により算出された組織の時間濃度曲線から、動脈領域を設定し、動脈領域での時間濃度曲線を算出する動脈の時間濃度曲線算出手段と、動脈の時間濃度曲線算出手段により算出された動脈の時間濃度曲線を、複数のパラメータを持つ動脈カーブモデルに当てはめる動脈の時間濃度曲線のカーブフィッティング手段と、組織の時間濃度曲線算出手段により算出された組織の時間濃度曲線を、複数のパラメータを持つ組織カーブモデルに当てはめる組織の時間濃度曲線のカーブフィッティング手段と、組織の時間濃度曲線のカーブフィッティング手段により算出されたパラメータの値を用いて、血流動態に関する情報を算出する血流情報算出手段とを備えることであり、組織の時間濃度曲線のカーブフィッティング手段は、動脈の時間濃度曲線のカーブフィッティング手段により算出されたパラメータの値の全てまたは一部を、組織カーブモデルのパラメータとして用いることを特徴とする。
【0010】
請求項4記載の本発明の特徴は、造影剤が注入された被検体を医用画像撮影装置で撮影した時系列の撮影データから、各画素の画素値の変化を計測して組織の時間濃度曲線を算出する組織の時間濃度曲線算出ステップと、組織の時間濃度曲線算出ステップで算出された組織の時間濃度曲線から、動脈領域を設定し、動脈領域での時間濃度曲線を算出する動脈の時間濃度曲線算出ステップと、動脈の時間濃度曲線算出ステップで算出された動脈の時間濃度曲線を、複数のパラメータを持つ動脈カーブモデルに当てはめる動脈の時間濃度曲線のカーブフィッティングステップと、組織の時間濃度曲線算出ステップで算出された組織の時間濃度曲線を、複数のパラメータを持つ組織カーブモデルに当てはめる組織の時間濃度曲線のカーブフィッティングステップと、組織の時間濃度曲線のカーブフィッティングステップで算出されたパラメータの値を用いて、血流動態に関する情報を算出する血流情報算出ステップとを血流動態解析装置が備えるコンピュータに実行させることであり、組織の時間濃度曲線のカーブフィッティングステップでは、動脈の時間濃度曲線のカーブフィッティングステップで算出されたパラメータの値の全てまたは一部を、組織カーブモデルのパラメータとして用いることを特徴とする。
【発明の効果】
【0011】
本発明によれば、血管の血流動態に関する情報を正確かつ安定に算出することが可能な血流動態解析装置およびその制御プログラムを提供することができる。
【図面の簡単な説明】
【0012】
【図1】本発明を適用した血流動態解析装置1の構成例を示す図である。
【図2】局所血流解析アプリケーションのアルゴリズムの一例を示す図である。
【図3】組織の血流動態に関する情報の算出処理を説明するフローチャートである。
【図4】図3のステップS5における、第1の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細を説明するフローチャートである。
【図5】第1の実施の形態での組織の時間濃度曲線のカーブフィッティングを説明する図である。
【図6】カーブフィッティングに多項式を用いる例を説明する図である。
【図7】図3のステップS5における、第2の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細を説明するフローチャートである。
【図8】第2の実施の形態での組織の時間濃度曲線のカーブフィッティングを説明する図である。
【図9】図3のステップS5における、第3の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細を説明するフローチャートである。
【図10】第3の実施の形態での組織の時間濃度曲線のカーブフィッティングを説明する図である。
【図11】平均通過時間MTTを算出した結果の一例を示す図である。
【発明を実施するための形態】
【0013】
以下、本発明の実施の形態について図面を参照して詳細に説明する。
【0014】
図1は、本発明を適用した血流動態解析装置1の構成例を示す図である。
【0015】
血流動態解析装置1は、CPU(Central Processing Unit)1a、ROM(Read Only Memory)1b、RAM(Random Access Memory)1c、および入出力インターフェイス1dが、バス1eを介して接続されている。入出力インターフェイス1dには、記憶部1f、通信部1g、表示部1h、入力部1i、およびリムーバブルメディア1jが接続されている。
【0016】
CPU1aは、入力部1iからの入力信号に基づいて血流動態解析装置1を起動するためのブートプログラムをROM1bから読み出して実行し、記憶部1fに格納されている各種オペレーティングシステムを読み出す。またCPU1aは、入力部1iからの入力信号に基づいて各種の制御を行ったり、ROM1bや記憶部1fに記憶されたプログラムおよびデータを読み出してRAM1cにロードしたり、あるいはRAM1cから読み出されたプログラムのコマンドに基づいて、データ演算または加工などの一連の処理を実行する。
【0017】
記憶部1fは、半導体メモリや磁気ディスクなどで構成されており、CPU1aで実行されるプログラムやデータを記憶する。記憶部1fには、CPU1aが実行するプログラムとして、例えば、局所血流動態に関する情報を解析するためのアプリケーションが用意される。以下、局所血流動態に関する情報を解析するこのアプリケーションを局所血流解析アプリケーションという。
【0018】
通信部1gは、LAN(Local Area Network)カードやモデムなどで構成されており、血流動態解析装置1をローカルエリアネットワークやインターネットといった通信媒体に接続することを可能にする。すなわち通信部1gは、通信媒体から受信したデータを、入出力インターフェイス1dおよびバス1eを介してCPU1aに送信し、CPU1aからバス1eおよび入出力インターフェイス1dを介して受信したデータを、通信媒体に送信する。
【0019】
表示部1hは、例えば液晶ディスプレイであり、CPU1aからバス1eおよび入出力インターフェイス1dを介して受信した信号に基づいて、CPU1aの処理結果などを表示する。
【0020】
入力部1iは、血流動態解析装置1の操作者が各種の操作を入力するキーボードやマウスなどの入力デバイスにより構成されており、操作者の操作に基づいて入力信号を生成し、入出力インターフェイス1dおよびバス1eを介してCPU1aに送信する。
【0021】
リムーバブルメディア1jは、例えば光ディスクやフレキシブルディスクなどであり、図示せぬディスクドライブによって読み出されたデータが、入出力インターフェイス1dおよびバス1eを介してCPU1aに送信され、CPU1aからバス1eおよび入出力インターフェイス1dを介して受信したデータが、ディスクドライブによって書き込まれる。
【0022】
本実施の形態の血流動態解析装置1は、X線CT装置やMRI装置などの医用画像撮影装置(いわゆる医用モダリティ)により、ダイナミックス撮影により収集された画像データから血流動態を解析する装置である。このため、血流動態解析装置1は、かかる画像データを入手できる環境にあればよく、医用画像撮影装置と一体に構成されていてもよいし、医用画像撮影装置とは別体で構成されていてもよい。別体で構成される場合には、収集された画像データはリムーバブルメディア1jで、または通信部1gを介して図示せぬ医用画像撮影装置から血流動態解析装置1に送られる。
【0023】
また、血流動態解析装置1は、血流動態解析結果を、通信部1gを介して、PACS(Picture Archiving and Communication System)に転送し、そこに保存させることができる。これによって、特定の画像データの血流動態解析結果の検索や閲覧などを容易に行うことが可能となる。
【0024】
医用画像撮影装置には、X線CT装置やMRI装置のほか、超音波診断装置やPET(Positron Emission Tomography)装置などが含まれる。
【0025】
図2は、局所血流解析アプリケーションのアルゴリズムの一例を示す図である。
【0026】
図2に示すように、局所血流解析アプリケーション2は、時系列データ取得部2a、造影濃度曲線算出部2b、動脈の時間濃度曲線算出部2c、動脈の時間濃度曲線のカーブフィッティング部2d、および組織の時間濃度曲線のカーブフィッティング部2eのアルゴリズムから成り立っている。
【0027】
時系列データ取得部2aは、造影剤を静脈注入(ボーラス注入)した被検者について、医用画像撮影装置が所定時間撮影して生成した時系列の撮影データを取得し、取得した時系列の撮影データを造影濃度曲線算出部2bに出力する。
【0028】
造影濃度曲線算出部2bは、時系列データ取得部2aから出力される時系列の撮影データを入力し、各画素の画素値の変化を計測して造影濃度曲線を算出する。
【0029】
つまり、被検者に造影剤を比較的短時間で静脈注入すると、造影剤は静脈を通り心臓に戻り、肺を循環して再び心臓に戻り、動脈へ流入し、各臓器の組織に到達し、再び静脈を通り心臓へ戻る。ここで、例えば、医用画像撮影装置が脳の領域のCT画像やMRI画像を時系列で撮影していた場合、造影剤濃度の変化に応じた各画素の画素値の変化を計測することができる。造影濃度曲線算出部2bは、この計測した各画素の画素値の変化に基づいて、造影濃度曲線を算出することができる。
【0030】
例えば、CT画像の場合、次式(1)に示すように、各時相の画素値Sから造影剤到着前の画素値Sを引き算することにより、造影剤濃度に比例する量(TDC:Time Density Curve)が算出される。式(1)において、iは、撮影時相を表している。
【0031】
=S−S ・・・(1)
また例えば、MRI画像の場合、次式(2)に示すように、各時相の画素値Sと造影剤到着前の画素値Sの比の対数(ΔR2)から、造影剤濃度に比例する量(TCC:Time Concentration Curve)が算出される。式(2)において、iは、撮影時相を表し、TEは、エコー時間を表している。
【0032】
ΔR2=−log(S/S)/TE ・・・(2)
造影濃度曲線算出部2bは、上記の式(1)または式(2)で算出した造影濃度曲線を、組織の時間濃度曲線Cとし、算出値を動脈の時間濃度曲線算出部2c、および、組織の時間濃度曲線のカーブフィッティング部2eに出力する。
【0033】
この造影濃度曲線算出部2bは、時系列の撮影データから、各画素の画素値の変化を計測して組織の時間濃度曲線を算出する組織の時間濃度曲線算出手段として機能する。
【0034】
動脈の時間濃度曲線算出部2cは、造影濃度曲線算出部2bから出力される組織の時間濃度曲線Cを入力し、その曲線から、CT画像やMRI画像の解析単位となるボクセルでの動脈の時間濃度曲線を算出する。つまり、CT画像やMRI画像の解析単位となるボクセルの大きさは、0.5mm乃至3mm程度の大きさであり、組織の中の主要な動脈であれば、血管にほぼ包含されるボクセルが存在する。
【0035】
動脈の時間濃度曲線算出部2cは、このようなボクセルを自動的あるいはユーザからの指定によって特定することで、動脈領域を設定し、その動脈領域での時間濃度曲線(ここでは、i番目の撮影時刻での濃度値ai)を算出する動脈の時間濃度曲線算出手段として機能する。
【0036】
動脈時間濃度曲線算出部2cは、算出した動脈の時間濃度曲線aiを、動脈の時間濃度曲線のカーブフィッティング部2dに出力する。
【0037】
動脈の時間濃度曲線のカーブフィッティング部2dは、動脈の時間濃度曲線算出部2cから出力される動脈の時間濃度曲線aを入力する。動脈の時間濃度曲線のカーブフィッティング部2dは、次式(3)によるガンマ確率密度関数の定数倍の関数(ガンマ確率密度関数に比例する関数)を動脈の時間濃度曲線のカーブモデル(以下、動脈カーブモデルと称する)として適用し、入力した動脈の時間濃度曲線aのカーブフィッティング(動脈カーブモデルへの当てはめ)を行う動脈の時間濃度曲線のカーブフィッティング手段として機能する。
【数1】

【0038】
ここで言うガンマ確率密度関数とは、時間βに1回の割合で通過できる何重かの障壁αの全部を通過する所要時間の確率密度を表わしている。動脈関数は、肺の障壁を通過して分散された関数と考えることができるため、このガンマ確率密度関数を、動脈カーブモデルとして適用し、動脈カーブモデルへの当てはめ(近似)を行うようにする。
【0039】
動脈の時間濃度曲線a(t)は、上記の式(3)で表されるガンマ確率密度関数の定数倍の関数でのカーブフィッティングによって、次式(4)により表わされる。式(4)において、tは、時刻を表し、t0,α,β,Kは、いずれもパラメータを表している。αは、障壁の数であり、βは、障壁αの全部を通過する所要時間であり、Kは、線形パラメータである。
【数2】

【0040】
動脈の時間濃度曲線のカーブフィッティング部2dは、上記の式(4)における4つのパラメータt0,α,β,Kの全てを変数とし、次式(5)に示す、動脈カーブモデルと動脈の濃度値aとの2乗誤差を最小にするt0,α,β,Kを算出する。実際には、t0には、造影剤が対象臓器に到達する前の特定の時刻を用い、α,β,Kを算出する方法を用いる場合が多い。
【数3】

【0041】
動脈の時間濃度曲線のカーブフィッティング部2dは、算出したパラメータα,β,Kの値を、組織の時間濃度曲線のカーブフィッティング部2eに出力する。
【0042】
組織の時間濃度曲線のカーブフィッティング部2eは、造影濃度曲線算出部2bから出力される組織の時間濃度曲線Cを入力するとともに、動脈の時間濃度曲線のカーブフィッティング部2dから出力されたパラメータt0,α,β,Kの値を入力する。
【0043】
組織の時間濃度曲線のカーブフィッティング部2eは、予め定めた関数を組織の時間濃度曲線のカーブモデル(以下、組織カーブモデルと称する)として適用し、入力した組織の時間濃度曲線Cのカーブフィッティングを行う組織の時間濃度曲線のカーブフィッティング手段として機能する。
【0044】
ここで、カーブフィッティングの際、動脈の時間濃度曲線のカーブフィッティング部2dから出力されたパラメータt0,α,β,Kの値の全てまたは一部を、組織カーブモデルのパラメータとして用いる。なお、組織の時間濃度曲線のカーブフィッティング処理の詳細は、後述する。
【0045】
組織の時間濃度曲線のカーブフィッティング部2eは、カーブフィッティングにより得たパラメータの値を用いて、局所血流動態に関する情報である、平均通過時間MTT(Mean Transit Time)、局所脳血流量CBF(Cerebral Blood Flow)、および局所脳血液量CBV(Cerebral Blood Volume)を算出する血流情報算出手段としても機能する。
【0046】
次に、図3のフローチャートを参照して、本発明の実施の形態における組織の血流動態に関する情報の算出処理について説明する。この処理では、特に、脳組織を例に挙げ、その血流動態に関する情報を算出する処理について説明する。
【0047】
ステップS1において、時系列データ取得部2aは、医用画像撮影装置が所定時間撮影して生成した時系列の撮影データを取得する。
【0048】
ステップS2において、造影濃度曲線算出部2bは、ステップS1の処理で取得された時系列の撮影データから、各画素の画素値の変化を計測して造影濃度曲線を算出する。この造影濃度曲線を、組織の時間濃度曲線Cとする。
【0049】
ステップS3において、動脈の時間濃度曲線算出部2cは、ステップS2の処理で算出された組織の時間濃度曲線Cから、CT画像やMRI画像の解析単位となる動脈領域を設定し、その動脈領域での時間濃度曲線aiを算出する。
【0050】
ステップS4において、動脈の時間濃度曲線のカーブフィッティング部2dは、上記の式(3)で表わしたガンマ確率密度関数の定数倍の関数を、動脈カーブモデルとして適用し、ステップS3の処理で算出された動脈の時間濃度曲線aのカーブフィッティングを行う。動脈の時間濃度曲線aのカーブフィッティングによって、t0,α,β,Kのパラメータの値が算出される。
【0051】
ステップS5において、組織の時間濃度曲線のカーブフィッティング部2eは、予め定めた関数を、組織カーブモデルとして適用し、ステップS2の処理で算出された組織の時間濃度曲線Cのカーブフィッティングを行う。組織の時間濃度曲線Cのカーブフィッティング処理の詳細は、図4、図7、および図9のフローチャートを参照して後述するが、ステップS4の処理で算出されたt0,α,β,Kのパラメータの値の全てまたは一部が、組織カーブモデルのパラメータとして用いられる。
【0052】
そして、組織の時間濃度曲線のカーブフィッティング部2eは、組織の時間濃度曲線Cのカーブフィッティングにより得られたパラメータの値を用いて、血流動態に関する情報である、平均通過時間MTT、局所脳血液量CBV、および局所脳血流量CBFを算出する。
【0053】
以上のように、血流動態解析装置1は、動脈の時間濃度曲線aをガンマ確率密度関数の定数倍の関数を用いてカーブフィッティングし、そのカーブフィッティングにより得られたパラメータの値の全てまたは一部を、組織の時間濃度曲線Cのカーブフィッティングのパラメータに用いることによって、安定かつ正確に、平均通過時間MTT、局所脳血液量CBV、および局所脳血流量CBFを算出することが可能となる。
【0054】
また、血流動態解析装置1は、算出した平均通過時間MTT、局所脳血液量CBV、および局所脳血流量CBFに基づいて、対応する画像を表示部1hに表示させることができ、操作者は、血流動態の状態を容易に把握することが可能となる。
【0055】
次に、図4のフローチャートを参照して、図3のステップS5における、第1の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細について説明する。
【0056】
第1の実施の形態では、図5に示すように、組織の時間濃度曲線のカーブモデルに不完全ガンマ関数の差を用いる点に特徴を有する。図5において、横軸は時間を表わし、縦軸は濃度を表わしており、a(t)は、ガンマ確率密度関数の定数倍の関数でカーブフィッティングした動脈の時間濃度曲線を表わし、c(t)は、不完全ガンマ関数の差でカーブフィッティングした組織の時間濃度曲線を表わし、f(t)は、組織のインパルス応答関数を表わしている。
【0057】
ステップS11において、組織の時間濃度曲線のカーブフィッティング部2eは、次式(6)で表わされる不完全ガンマ関数の差、または、式(6)の変形である次式(7)を組織カーブモデルとして適用する。なお、実際の組織の時間濃度曲線Cのカーブフィッティングには、式(7)を微分した次式(8)を適用することもできる。
【数4】

【0058】
ステップS12において、組織の時間濃度曲線のカーブフィッティング部2eは、ステップS11の処理で適用したカーブモデルを用いて、組織の時間濃度曲線Cのカーブフィッティングを行う。なお、組織カーブモデルの複数のパラメータのうち、α,β,Kのパラメータは、図3のステップS4の処理で算出された値を用い、組織の時間濃度曲線Cのカーブフィッティングにより新たに決定するパラメータは、t、α、Kである。tは、理想濃度曲線からの遅れ時間であり、αは、障壁の数であり、Kは、線形パラメータである。
【0059】
また、カーブフィッティングの高速算出方法として、多項式yをカーブフィッティングに用いることができる。例えば、組織の時間濃度曲線のピーク時刻をtとすると、次式(9)で表わすことができる。
【数5】

【0060】
また、上記の式(9)を展開すると、次式(10)を得ることができる。
【数6】

【0061】
そして、上記の式(10)から、t−t−tとαの間の関係を多項式で近似することができる。
【0062】
図6は、横軸t−t−tと縦軸αの関係の一例を示すグラフである。図6の例の場合、t−t−tとαの間の関係を、次式(11)で近似できることを示している。
【数7】

【0063】
すなわち、上記の式(11)を用いて、組織の時間濃度曲線Cのピーク時刻tから、極めて容易に、αを求めることができることを示している。Kは、線形パラメータのため、最適化法の適用は不要である。
【0064】
このように、t−t−tとαの間の関係を近似した多項式を用いて、高速に組織の時間濃度曲線Cをカーブフィッティングすることができる。
【0065】
また、カーブフィッティングの他の算出方法として、上記の式(11)で求めたパラメータαの値を用いて、次式(12)をカーブフィッティングに適用することもできる。
【数8】

【0066】
上記の式(12)で決定しなければならないパラメータは、tのみになり、1変数の最適化であるため、高速に値を決定することができる。
【0067】
ステップS13において、組織の時間濃度曲線のカーブフィッティング部2eは、ステップS12の処理で得られたパラメータt、α、Kの値、および、図3のステップS4の処理で得られたパラメータt0,α,β,Kの値から、次式(13)に示すようにして、平均通過時間MTTを算出する。
【0068】
MTT=α・β ・・・(13)
ステップS14において、組織の時間濃度曲線のカーブフィッティング部2eは、次式(14)に示すように、組織の時間濃度曲線Cのピーク部のカーブ下の面積(AUC:Area Under the Curve)ACを、動脈の時間濃度曲線aiのピーク部のカーブ下の面積ACで除算して局所脳血液量CBVを算出する。
【0069】
CBV=ACc/ACa ・・・(14)
ステップS15において、組織の時間濃度曲線のカーブフィッティング部2eは、次式(15)に示すように、局所脳血液量CBVを、平均通過時間MTTで除算して局所脳血流量CBFを算出する。
【0070】
CBF=CBV/MTT ・・・(15)
なお、以上の処理において、ステップS13とステップS14の処理は順不同である。
【0071】
以上のように、第1の実施の形態によれば、ガンマ確率密度関数を組織カーブモデルとして適用し、動脈の時間濃度曲線aのカーブフィッティングより得られたパラメータを用いて、組織の時間濃度曲線Cのカーブフィッティング(近似)を行うようにすることで、安定かつ正確に、組織の局所脳血流に関する情報である、局所脳血液量CBV、局所脳血流量CBF、および平均通過時間MTTを算出することが可能となる。本実施例では、組織の入力関数(カーブフィッティングされた動脈の時間濃度曲線)と出力関数(インパルス応答関数)が双方ともガンマ確率密度関数で表わせると仮定しており、式(6)や式(7)は、造影剤が組織に入り、そして出る間に組織にとどまっている造影剤の濃度を表わすこととなり、物理的意味が明確であり、より正確に解析を行うことができる。
【0072】
次に、図7のフローチャートを参照して、図3のステップS5における、第2の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細について説明する。
【0073】
第2の実施の形態では、図8に示すように、組織のインパルス応答f(t)に矩形(BOX)関数を用いる点に特徴を有する。図8において、横軸は時間を表わし、縦軸は濃度を表わしており、a(t)は、ガンマ確率密度関数の定数倍の関数でカーブフィッティングした動脈の時間濃度曲線を表わし、c(t)は、組織の時間濃度曲線のカーブモデルでカーブフィッティングした組織の時間濃度曲線を表わしている。
【0074】
ここで、カーブフィッティングした動脈の時間濃度曲線を次式(16)で表わし、矩形関数(インパルス応答)を次式(17)で表わすと、それらのコンボリューション(畳み込み積分)は、次式(18)で表わすことができる。
【数9】

【0075】
ステップS21において、組織の時間濃度曲線のカーブフィッティング部2eは、上記の式(18)を組織カーブモデルとして適用する。Γ(a)は、ガンマ関数であり、P(a,x)は、不完全ガンマ関数である。従って、この組織カーブモデルは、不完全ガンマ関数を含む関数である。
【0076】
ステップS22において、組織の時間濃度曲線のカーブフィッティング部2eは、ステップS21の処理で適用したカーブモデルを用いて、組織の時間濃度曲線Cのカーブフィッティングを行う。なお、組織カーブモデルの複数のパラメータのうち、α,β,Kのパラメータは、図3のステップS4の処理で算出された値を用い、組織の時間濃度曲線Cのカーブフィッティングにより新たに決定するパラメータは、u、w、tである。
【0077】
例えば、矩形関数の幅が大きくないと仮定した場合、上記の式(18)は、次式(19)に示すように変形することができる。
【数10】

【0078】
上記の式(19)において、a(t−t)を動脈の時間濃度曲線aに、F(t)を組織の時間濃度曲線Cに置き換えると、次式(20)に示す連立方程式を得ることができる。パラメータu、wは、線形演算で解くことができ、パラメータtは、黄金比分割法などの最適化法を用いて解くことができる。
【数11】

【0079】
ステップS23において、組織の時間濃度曲線のカーブフィッティング部2eは、ステップS22の処理で得られたパラメータu、wの値から、次式(21)に示すようにして、平均通過時間MTT、局所脳血液量CBV、および局所脳血流量CBFを算出する。
【0080】
MTT=w
CBV=u
CBF=u/w ・・・(21)
以上のように、第2の実施の形態によれば、カーブフィッティングされた動脈の時間濃度曲線と矩形関数をコンボリューションした関数を、組織カーブモデルとして適用し、動脈の時間濃度曲線aのカーブフィッティングより得られたパラメータを用いて、組織の時間濃度曲線Cのカーブフィッティングを行うようにすることで、安定かつ正確に、組織の局所脳血流動態に関する情報である、局所脳血液量CBV、局所脳血流量CBF、および平均通過時間MTTを算出することが可能となる。
【0081】
次に、図9のフローチャートを参照して、図3のステップS5における、第3の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細について説明する。
【0082】
第3の実施の形態では、図10に示すように、組織のインパルス応答f(t)にガンマ確率密度関数の定数倍の関数を用いる点に特徴を有する。図10において、横軸は時間を表わし、縦軸は濃度を表わしており、a(t)は、ガンマ確率密度関数の定数倍の関数でカーブフィッティングした動脈の時間濃度曲線を表わし、c(t)は、組織の時間濃度曲線のカーブモデルでカーブフィッティングした組織の時間濃度曲線を表わしている。
【0083】
ここで、カーブフィッティングした動脈の時間濃度曲線を次式(22)で表わし、ガンマ確率密度関数の定数倍の関数(インパルス応答)を次式(23)で表わすと、それらのコンボリューションは、次式(24)で表わすことができる。
【数12】

【0084】
ステップS31において、組織の時間濃度曲線のカーブフィッティング部2eは、上記の式(24)を組織カーブモデルとして適用する。この組織カーブモデルは、動脈の時間濃度曲線aのカーブフィッティングに適用する動脈カーブモデルと同じである。
【0085】
ステップS32において、組織の時間濃度曲線のカーブフィッティング部2eは、ステップS21の処理で適用したカーブモデルを用いて、組織の時間濃度曲線Cのカーブフィッティングを行う。なお、組織カーブモデルの複数のパラメータのうち、α,β,Kのパラメータは、図3のステップS4の処理で算出された値を用い、組織の時間濃度曲線Cのカーブフィッティングにより新たに決定するパラメータは、t、α、Kである。
【0086】
ステップS33において、組織の時間濃度曲線のカーブフィッティング部2eは、ステップS32の処理で得られたパラメータαの値、および、図3のステップS4の処理で得られたパラメータα,β,Kの値から、次式(25)に示すようにして、平均通過時間MTT、局所脳血液量CBV、および局所脳血流量CBFを算出する。
【数13】

【0087】
なお、局所脳血流量CBFは、次式(26)の関係を用いて算出するようにしてもよい。
【0088】
CBF=CBV/MTT ・・・(26)
以上のように、第3の実施の形態によれば、カーブフィッティングされた動脈の時間濃度曲線とガンマ確率密度関数をコンボリューションした関数を、組織カーブモデルとして適用し、動脈の時間濃度曲線aのカーブフィッティングより得られたパラメータを用いて、組織の時間濃度曲線Cのカーブフィッティングを行うようにすることで、安定かつ正確に、組織の局所血流動態に関する情報である、局所脳血液量CBV、局所脳血流量CBF、および平均通過時間MTTを算出することが可能となる。
【0089】
図11は、従来のデコンボリューション法と本発明の方法とで数値実験を行い、平均通過時間MTTを算出した結果の一例を示している。横軸は、組織の時間濃度曲線の数値モデルに与えた通過時間を表し、縦軸は、平均通過時間MTTの算出値を表している。
【0090】
従来のデコンボリューション法では、平均通過時間MTTの値が小さい領域において、実際の平均通過時間MTTより大きな値が得られており、組織の時間濃度曲線の数値モデルに与えた通過時間と算出結果の線形性があまり良いとは言えない。つまり、時間濃度曲線のデコンボリューションは、数値的に不安定であり、測定ノイズが拡大されて結果に重畳されることが知られている。
【0091】
そこで、結果の安定化を図るため、様々な技術が使われているが、結果の正確性を重視するか、あるいは、安定性を重視するかはトレードオフの関係にあり、デコンボリューションの理論では、両者を同時に向上させることができない。
【0092】
一方、本実施の形態によれば、数値的なデコンボリューションを行うことなく、動脈の時間濃度曲線aをカーブフィッティングすることによって得たパラメータを利用して、組織の時間濃度曲線Cのカーブフィッティングを行い、その結果得たパラメータの値を用いて、平均通過時間MTTおよび局所脳血液量CBVを算出することができ、パラメータ算出の安定性の向上を期待することができる。
【0093】
従って、算出された平均通過時間MTTは、理論値に近い値が得られており、デコンボリューション法より線形性が改善されていることがわかる。また、CBF=CBV/MTTの関係から、局所脳血流量CBFについても、正確かつ安定した値を得ることができる。
【0094】
さらに、動脈の時間濃度曲線aのカーブフィッティングより得られたパラメータを、組織の時間濃度曲線Cのカーブフィッティングに用いることにより、動脈の時間濃度曲線の依存症を小さくすることができる。従って、造影剤の注入時間を長くすることが可能になり、患者への負担を軽減することができる。
【0095】
以上においては、対象組織として、脳組織を例に挙げ説明したが、本実施の形態では、脳組織に限られるものではなく、肝臓、心臓、あるいは肺などの他の組織を対象とすることも勿論可能である。
【0096】
なお、肝臓や肺など、呼吸による動きがある組織、あるいは、心臓では、組織の動き補正を行う必要がある。また動き補正だけでなく、ヘマトクリット値(一定量の血液中に含まれる赤血球の割合)の補正、脳の比重による補正、その他の補正を適宜行う必要があるが、本実施の形態では、特に言及しない。
【0097】
この発明は、上記実施の形態そのままに限定されるものではなく、実施段階ではその要旨を逸脱しない範囲で構成要素を変形して具体化したり、上記実施の形態に開示されている複数の構成要素を適宜組み合わせたりすることにより種々の発明を形成できる。例えば、実施の形態に示される全構成要素から幾つかの構成要素を削除してもよい。さらに、異なる実施の形態に亘る構成要素を適宜組み合わせても良い。
【符号の説明】
【0098】
1 血流動態解析装置
2 局所血流解析アプリケーション
2a 局所血流解析アプリケーション
2b 造影濃度曲線算出部
2c 動脈の時間濃度曲線算出部
2d 動脈の時間濃度曲線のカーブフィッティング部
2e 組織の時間濃度曲線のカーブフィッティング部

【特許請求の範囲】
【請求項1】
造影剤が注入された被検体を医用画像撮影装置で撮影した時系列の撮影データから、各画素の画素値の変化を計測して組織の時間濃度曲線を算出する組織の時間濃度曲線算出手段と、
前記組織の時間濃度曲線算出手段により算出された前記組織の時間濃度曲線から、動脈領域を設定し、前記動脈領域での時間濃度曲線を算出する動脈の時間濃度曲線算出手段と、
前記動脈の時間濃度曲線算出手段により算出された前記動脈の時間濃度曲線を、複数のパラメータを持つ動脈カーブモデルに当てはめる動脈の時間濃度曲線のカーブフィッティング手段と、
前記組織の時間濃度曲線算出手段により算出された前記組織の時間濃度曲線を、複数のパラメータを持つ組織カーブモデルに当てはめる組織の時間濃度曲線のカーブフィッティング手段と、
前記組織の時間濃度曲線のカーブフィッティング手段により算出されたパラメータの値を用いて、血流動態に関する情報を算出する血流情報算出手段と
を備え、
前記組織の時間濃度曲線のカーブフィッティング手段は、前記動脈の時間濃度曲線のカーブフィッティング手段により算出されたパラメータの値の全てまたは一部を、前記組織カーブモデルのパラメータとして用いる
ことを特徴とする血流動態解析装置。
【請求項2】
前記動脈カーブモデルには、ガンマ確率密度関数の定数倍の関数が用いられる
ことを特徴とする請求項1に記載の血流動態解析装置。
【請求項3】
前記組織カーブモデルには、不完全ガンマ関数の差を含む関数、または、ガンマ確率密度関数の定数倍の関数が用いられる
ことを特徴とする請求項1または2に記載の血流動態解析装置。
【請求項4】
造影剤が注入された被検体を医用画像撮影装置で撮影した時系列の撮影データから、各画素の画素値の変化を計測して組織の時間濃度曲線を算出する組織の時間濃度曲線算出ステップと、
前記組織の時間濃度曲線算出ステップで算出された前記組織の時間濃度曲線から、動脈領域を設定し、前記動脈領域での時間濃度曲線を算出する動脈の時間濃度曲線算出ステップと、
前記動脈の時間濃度曲線算出ステップで算出された前記動脈の時間濃度曲線を、複数のパラメータを持つ動脈カーブモデルに当てはめる動脈の時間濃度曲線のカーブフィッティングステップと、
前記組織の時間濃度曲線算出ステップで算出された前記組織の時間濃度曲線を、複数のパラメータを持つ組織カーブモデルに当てはめる組織の時間濃度曲線のカーブフィッティングステップと、
前記組織の時間濃度曲線のカーブフィッティングステップで算出されたパラメータの値を用いて、血流動態に関する情報を算出する血流情報算出ステップと
を血流動態解析装置が備えるコンピュータに実行させ、
前記組織の時間濃度曲線のカーブフィッティングステップでは、前記動脈の時間濃度曲線のカーブフィッティングステップで算出されたパラメータの値の全てまたは一部を、前記組織カーブモデルのパラメータとして用いる
ことを特徴とする制御プログラム。

【図1】
image rotate

【図2】
image rotate

【図3】
image rotate

【図4】
image rotate

【図5】
image rotate

【図6】
image rotate

【図7】
image rotate

【図8】
image rotate

【図9】
image rotate

【図10】
image rotate

【図11】
image rotate


【公開番号】特開2011−83437(P2011−83437A)
【公開日】平成23年4月28日(2011.4.28)
【国際特許分類】
【出願番号】特願2009−238678(P2009−238678)
【出願日】平成21年10月15日(2009.10.15)
【出願人】(000003078)株式会社東芝 (54,554)
【出願人】(594164542)東芝メディカルシステムズ株式会社 (4,066)
【Fターム(参考)】