(57)【特許請求の範囲】
【請求項１】
造影剤が注入された被検体を医用画像撮影装置で撮影した時系列の撮影データから、各画素の画素値の変化を計測して組織の時間濃度曲線を算出する組織の時間濃度曲線算出手段と、
前記組織の時間濃度曲線算出手段により算出された前記組織の時間濃度曲線から、動脈領域を設定し、前記動脈領域での時間濃度曲線を算出する動脈の時間濃度曲線算出手段と、
前記動脈の時間濃度曲線算出手段により算出された前記動脈の時間濃度曲線を、複数のパラメータを持ち、ガンマ確率密度関数の定数倍の関数である動脈カーブモデルに当てはめる動脈の時間濃度曲線のカーブフィッティング手段と、
前記組織の時間濃度曲線算出手段により算出された前記組織の時間濃度曲線を、複数のパラメータを持ち、不完全ガンマ関数の差を含む関数、組織のインパルス応答に矩形関数、或いは、ガンマ確率密度関数の定数倍の関数を用いた、いずれかの組織カーブモデルに当てはめる組織の時間濃度曲線のカーブフィッティング手段と、
前記組織の時間濃度曲線のカーブフィッティング手段により算出されたパラメータの値を用いて、血流動態に関する情報を算出する血流情報算出手段と
を備え、
前記組織の時間濃度曲線のカーブフィッティング手段は、前記動脈の時間濃度曲線のカーブフィッティング手段により算出されたパラメータの値の全てまたは一部を、前記組織カーブモデルのパラメータとして用いる
ことを特徴とする血流動態解析装置。
【請求項２】
造影剤が注入された被検体を医用画像撮影装置で撮影した時系列の撮影データから、各画素の画素値の変化を計測して組織の時間濃度曲線を算出する組織の時間濃度曲線算出ステップと、
前記組織の時間濃度曲線算出ステップで算出された前記組織の時間濃度曲線から、動脈領域を設定し、前記動脈領域での時間濃度曲線を算出する動脈の時間濃度曲線算出ステップと、
前記動脈の時間濃度曲線算出ステップで算出された前記動脈の時間濃度曲線を、複数のパラメータを持ち、ガンマ確率密度関数の定数倍の関数である動脈カーブモデルに当てはめる動脈の時間濃度曲線のカーブフィッティングステップと、
前記組織の時間濃度曲線算出ステップで算出された前記組織の時間濃度曲線を、複数のパラメータを持ち、不完全ガンマ関数の差を含む関数、組織のインパルス応答に矩形関数、或いは、ガンマ確率密度関数の定数倍の関数を用いた、いずれかの組織カーブモデルに当てはめる組織の時間濃度曲線のカーブフィッティングステップと、
前記組織の時間濃度曲線のカーブフィッティングステップで算出されたパラメータの値を用いて、血流動態に関する情報を算出する血流情報算出ステップと
を血流動態解析装置が備えるコンピュータに実行させ、
前記組織の時間濃度曲線のカーブフィッティングステップでは、前記動脈の時間濃度曲線のカーブフィッティングステップで算出されたパラメータの値の全てまたは一部を、前記組織カーブモデルのパラメータとして用いることを特徴とする制御プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、血流動態解析装置およびその制御プログラムに関し、特に、血管の血流動態に関する情報を正確かつ安定に算出することができる血流動態解析装置およびその制御プログラムに関する。
【背景技術】
【０００２】
従来、Ｘ線CT（Computed Tomography）装置やMRI（Magnetic Resonance Imaging）装置では、静脈から造影剤を投与し、時系列画像データを収集して、その画像を元に解析し、組織内の血流動態に関する情報を得る検査が行われている。これはパーヒュージョン検査と呼ばれ、造影剤の撮影断面への集中の度合いが画像の濃度変化として収集できることを利用したものである。
【０００３】
例えば、肺循環や造影剤投与のばらつきをなくすために、組織に流入する動脈の時間濃度曲線（Time Intensity Cur：TIC）：ａｉを入力関数として、測定した組織の時間濃度曲線Ｃｉとのデコンボリューション（Deconvolution：逆畳み込み積分）を行い、組織固有のインパルス応答関数（伝達関数）Ｒｉ（residue function）を求め、得られたインパルス応答関数から、血流動態を定量的に表す指標である血流量（脳の場合、CBF：Cerebral Blood Flow）、平均通過時間：MTT（Mean Transit Time）、および、血液量（脳の場合、CBV：Cerebral Blood Volume）などのパラメータを算出している。この解析方法は、デコンボリューション法と呼ばれている。
【０００４】
また、非特許文献１では、動脈の時間濃度曲線をガンマ関数でフィッティングし、組織の時間濃度曲線をガンマ関数でフィッティングし、ガンマ関数でそれぞれフィッティングされた動脈の時間濃度曲線と組織の時間濃度曲線とのデコンボリューションを行い、得られたインパルス応答関数から、血流量、平均通過時間、および血液量を求めるようにしている。
【先行技術文献】
【非特許文献】
【０００５】
【非特許文献１】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
【発明の概要】
【発明が解決しようとする課題】
【０００６】
従来のデコンボリューション法では、組織の時間濃度曲線のデコンボリューションが数値的に不安定であり、インパルス応答関数の算出精度がある一定以上にはならない。また、インパルス応答関数のピークの面積は正確に求めることができるものの、ピークの幅や高さを独立に正確に求めることができない。つまり、デコンボリューション法で求めたい値は、ピークの高さ（CBFに相当）またはピークの幅（MTTに相当）であり、十分な精度が得られないという課題があった。
【０００７】
また、非特許文献１の技術においても、デコンボリューションを行っているため、数値的に不安定であり、また求めるパラメータの数が多いため、正確性に欠ける課題があった。
【０００８】
本発明はこのような状況に鑑みてなされたものであり、その目的は、血管の血流動態に関する情報を正確かつ安定に算出することが可能な血流動態解析装置およびその制御プログラムを提供することである。
【課題を解決するための手段】
【０００９】
請求項１記載の本発明の特徴は、血流動態解析装置において、造影剤が注入された被検体を医用画像撮影装置で撮影した時系列の撮影データから、各画素の画素値の変化を計測して組織の時間濃度曲線を算出する組織の時間濃度曲線算出手段と、組織の時間濃度曲線算出手段により算出された組織の時間濃度曲線から、動脈領域を設定し、動脈領域での時間濃度曲線を算出する動脈の時間濃度曲線算出手段と、動脈の時間濃度曲線算出手段により算出された動脈の時間濃度曲線を、複数のパラメータを持ち、ガンマ確率密度関数の定数倍の関数である動脈カーブモデルに当てはめる動脈の時間濃度曲線のカーブフィッティング手段と、組織の時間濃度曲線算出手段により算出された組織の時間濃度曲線を、複数のパラメータを持ち、不完全ガンマ関数の差を含む関数、組織のインパルス応答に矩形関数、或いは、ガンマ確率密度関数の定数倍の関数を用いた、いずれかの組織カーブモデルに当てはめる組織の時間濃度曲線のカーブフィッティング手段と、組織の時間濃度曲線のカーブフィッティング手段により算出されたパラメータの値を用いて、血流動態に関する情報を算出する血流情報算出手段とを備えることであり、組織の時間濃度曲線のカーブフィッティング手段は、動脈の時間濃度曲線のカーブフィッティング手段により算出されたパラメータの値の全てまたは一部を、組織カーブモデルのパラメータとして用いることを特徴とする。
【００１０】
請求項２記載の本発明の特徴は、制御プログラムにおいて、造影剤が注入された被検体を医用画像撮影装置で撮影した時系列の撮影データから、各画素の画素値の変化を計測して組織の時間濃度曲線を算出する組織の時間濃度曲線算出ステップと、組織の時間濃度曲線算出ステップで算出された組織の時間濃度曲線から、動脈領域を設定し、動脈領域での時間濃度曲線を算出する動脈の時間濃度曲線算出ステップと、動脈の時間濃度曲線算出ステップで算出された動脈の時間濃度曲線を、複数のパラメータを持ち、ガンマ確率密度関数の定数倍の関数である動脈カーブモデルに当てはめる動脈の時間濃度曲線のカーブフィッティングステップと、組織の時間濃度曲線算出ステップで算出された組織の時間濃度曲線を、複数のパラメータを持ち、不完全ガンマ関数の差を含む関数、組織のインパルス応答に矩形関数、或いは、ガンマ確率密度関数の定数倍の関数を用いた、いずれかの組織カーブモデルに当てはめる組織の時間濃度曲線のカーブフィッティングステップと、組織の時間濃度曲線のカーブフィッティングステップで算出されたパラメータの値を用いて、血流動態に関する情報を算出する血流情報算出ステップとを血流動態解析装置が備えるコンピュータに実行させることであり、組織の時間濃度曲線のカーブフィッティングステップでは、動脈の時間濃度曲線のカーブフィッティングステップで算出されたパラメータの値の全てまたは一部を、組織カーブモデルのパラメータとして用いることを特徴とする。
【発明の効果】
【００１１】
本発明によれば、血管の血流動態に関する情報を正確かつ安定に算出することが可能な血流動態解析装置およびその制御プログラムを提供することができる。
【図面の簡単な説明】
【００１２】
【図１】本発明を適用した血流動態解析装置１の構成例を示す図である。
【図２】局所血流解析アプリケーションのアルゴリズムの一例を示す図である。
【図３】組織の血流動態に関する情報の算出処理を説明するフローチャートである。
【図４】図３のステップＳ５における、第１の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細を説明するフローチャートである。
【図５】第１の実施の形態での組織の時間濃度曲線のカーブフィッティングを説明する図である。
【図６】カーブフィッティングに多項式を用いる例を説明する図である。
【図７】図３のステップＳ５における、第２の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細を説明するフローチャートである。
【図８】第２の実施の形態での組織の時間濃度曲線のカーブフィッティングを説明する図である。
【図９】図３のステップＳ５における、第３の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細を説明するフローチャートである。
【図１０】第３の実施の形態での組織の時間濃度曲線のカーブフィッティングを説明する図である。
【図１１】平均通過時間MTTを算出した結果の一例を示す図である。
【発明を実施するための形態】
【００１３】
以下、本発明の実施の形態について図面を参照して詳細に説明する。
【００１４】
図１は、本発明を適用した血流動態解析装置１の構成例を示す図である。
【００１５】
血流動態解析装置１は、CPU(Central Processing Unit)１ａ、ROM（Read Only Memory）１ｂ、RAM（Random Access Memory）１ｃ、および入出力インターフェイス１ｄが、バス１ｅを介して接続されている。入出力インターフェイス１ｄには、記憶部１ｆ、通信部１ｇ、表示部１ｈ、入力部１ｉ、およびリムーバブルメディア１ｊが接続されている。
【００１６】
CPU１ａは、入力部１ｉからの入力信号に基づいて血流動態解析装置１を起動するためのブートプログラムをROM１ｂから読み出して実行し、記憶部１ｆに格納されている各種オペレーティングシステムを読み出す。またCPU１ａは、入力部１ｉからの入力信号に基づいて各種の制御を行ったり、ROM１ｂや記憶部１ｆに記憶されたプログラムおよびデータを読み出してRAM１ｃにロードしたり、あるいはRAM１ｃから読み出されたプログラムのコマンドに基づいて、データ演算または加工などの一連の処理を実行する。
【００１７】
記憶部１ｆは、半導体メモリや磁気ディスクなどで構成されており、CPU１ａで実行されるプログラムやデータを記憶する。記憶部１ｆには、CPU１ａが実行するプログラムとして、例えば、局所血流動態に関する情報を解析するためのアプリケーションが用意される。以下、局所血流動態に関する情報を解析するこのアプリケーションを局所血流解析アプリケーションという。
【００１８】
通信部１ｇは、LAN（Local Area Network）カードやモデムなどで構成されており、血流動態解析装置１をローカルエリアネットワークやインターネットといった通信媒体に接続することを可能にする。すなわち通信部１ｇは、通信媒体から受信したデータを、入出力インターフェイス１ｄおよびバス１ｅを介してCPU１ａに送信し、CPU１ａからバス１ｅおよび入出力インターフェイス１ｄを介して受信したデータを、通信媒体に送信する。
【００１９】
表示部１ｈは、例えば液晶ディスプレイであり、CPU１ａからバス１ｅおよび入出力インターフェイス１ｄを介して受信した信号に基づいて、CPU１ａの処理結果などを表示する。
【００２０】
入力部１ｉは、血流動態解析装置１の操作者が各種の操作を入力するキーボードやマウスなどの入力デバイスにより構成されており、操作者の操作に基づいて入力信号を生成し、入出力インターフェイス１ｄおよびバス１ｅを介してCPU１ａに送信する。
【００２１】
リムーバブルメディア１ｊは、例えば光ディスクやフレキシブルディスクなどであり、図示せぬディスクドライブによって読み出されたデータが、入出力インターフェイス１ｄおよびバス１ｅを介してCPU１ａに送信され、CPU１ａからバス１ｅおよび入出力インターフェイス１ｄを介して受信したデータが、ディスクドライブによって書き込まれる。
【００２２】
本実施の形態の血流動態解析装置１は、Ｘ線CT装置やMRI装置などの医用画像撮影装置（いわゆる医用モダリティ）により、ダイナミックス撮影により収集された画像データから血流動態を解析する装置である。このため、血流動態解析装置１は、かかる画像データを入手できる環境にあればよく、医用画像撮影装置と一体に構成されていてもよいし、医用画像撮影装置とは別体で構成されていてもよい。別体で構成される場合には、収集された画像データはリムーバブルメディア１ｊで、または通信部１ｇを介して図示せぬ医用画像撮影装置から血流動態解析装置１に送られる。
【００２３】
また、血流動態解析装置１は、血流動態解析結果を、通信部１ｇを介して、PACS(Picture Archiving and Communication System)に転送し、そこに保存させることができる。これによって、特定の画像データの血流動態解析結果の検索や閲覧などを容易に行うことが可能となる。
【００２４】
医用画像撮影装置には、Ｘ線CT装置やMRI装置のほか、超音波診断装置やPET（Positron Emission Tomography）装置などが含まれる。
【００２５】
図２は、局所血流解析アプリケーションのアルゴリズムの一例を示す図である。
【００２６】
図２に示すように、局所血流解析アプリケーション２は、時系列データ取得部２ａ、造影濃度曲線算出部２ｂ、動脈の時間濃度曲線算出部２ｃ、動脈の時間濃度曲線のカーブフィッティング部２ｄ、および組織の時間濃度曲線のカーブフィッティング部２ｅのアルゴリズムから成り立っている。
【００２７】
時系列データ取得部２ａは、造影剤を静脈注入（ボーラス注入）した被検者について、医用画像撮影装置が所定時間撮影して生成した時系列の撮影データを取得し、取得した時系列の撮影データを造影濃度曲線算出部２ｂに出力する。
【００２８】
造影濃度曲線算出部２ｂは、時系列データ取得部２ａから出力される時系列の撮影データを入力し、各画素の画素値の変化を計測して造影濃度曲線を算出する。
【００２９】
つまり、被検者に造影剤を比較的短時間で静脈注入すると、造影剤は静脈を通り心臓に戻り、肺を循環して再び心臓に戻り、動脈へ流入し、各臓器の組織に到達し、再び静脈を通り心臓へ戻る。ここで、例えば、医用画像撮影装置が脳の領域のCT画像やMRI画像を時系列で撮影していた場合、造影剤濃度の変化に応じた各画素の画素値の変化を計測することができる。造影濃度曲線算出部２ｂは、この計測した各画素の画素値の変化に基づいて、造影濃度曲線を算出することができる。
【００３０】
例えば、CT画像の場合、次式（１）に示すように、各時相の画素値Ｓｉから造影剤到着前の画素値Ｓ０を引き算することにより、造影剤濃度に比例する量（TDC：Time Density Curve）が算出される。式（１）において、ｉは、撮影時相を表している。
【００３１】
Ｃｉ＝Ｓｉ−Ｓ０                                        ・・・（１）
また例えば、MRI画像の場合、次式（２）に示すように、各時相の画素値Ｓｉと造影剤到着前の画素値Ｓ０の比の対数（ΔR2＊）から、造影剤濃度に比例する量（TCC：Time Concentration Curve）が算出される。式（２）において、ｉは、撮影時相を表し、TEは、エコー時間を表している。
【００３２】
ΔＲ２＊＝−log（Ｓｉ／Ｓ０）／TE                        ・・・（２）
造影濃度曲線算出部２ｂは、上記の式（１）または式（２）で算出した造影濃度曲線を、組織の時間濃度曲線Ｃｉとし、算出値を動脈の時間濃度曲線算出部２ｃ、および、組織の時間濃度曲線のカーブフィッティング部２ｅに出力する。
【００３３】
この造影濃度曲線算出部２ｂは、時系列の撮影データから、各画素の画素値の変化を計測して組織の時間濃度曲線を算出する組織の時間濃度曲線算出手段として機能する。
【００３４】
動脈の時間濃度曲線算出部２ｃは、造影濃度曲線算出部２ｂから出力される組織の時間濃度曲線Ｃｉを入力し、その曲線から、CT画像やMRI画像の解析単位となるボクセルでの動脈の時間濃度曲線を算出する。つまり、CT画像やMRI画像の解析単位となるボクセルの大きさは、０．５ｍｍ乃至３ｍｍ程度の大きさであり、組織の中の主要な動脈であれば、血管にほぼ包含されるボクセルが存在する。
【００３５】
動脈の時間濃度曲線算出部２ｃは、このようなボクセルを自動的あるいはユーザからの指定によって特定することで、動脈領域を設定し、その動脈領域での時間濃度曲線（ここでは、ｉ番目の撮影時刻での濃度値ａi）を算出する動脈の時間濃度曲線算出手段として機能する。
【００３６】
動脈時間濃度曲線算出部２ｃは、算出した動脈の時間濃度曲線ａiを、動脈の時間濃度曲線のカーブフィッティング部２ｄに出力する。
【００３７】
動脈の時間濃度曲線のカーブフィッティング部２ｄは、動脈の時間濃度曲線算出部２ｃから出力される動脈の時間濃度曲線ａｉを入力する。動脈の時間濃度曲線のカーブフィッティング部２ｄは、次式（３）によるガンマ確率密度関数の定数倍の関数（ガンマ確率密度関数に比例する関数）を動脈の時間濃度曲線のカーブモデル（以下、動脈カーブモデルと称する）として適用し、入力した動脈の時間濃度曲線ａｉのカーブフィッティング（動脈カーブモデルへの当てはめ）を行う動脈の時間濃度曲線のカーブフィッティング手段として機能する。
【数１】
0005591512.tif 000002
【００３８】
ここで言うガンマ確率密度関数とは、時間βに１回の割合で通過できる何重かの障壁αの全部を通過する所要時間の確率密度を表わしている。動脈関数は、肺の障壁を通過して分散された関数と考えることができるため、このガンマ確率密度関数を、動脈カーブモデルとして適用し、動脈カーブモデルへの当てはめ（近似）を行うようにする。
【００３９】
動脈の時間濃度曲線ａ（t）は、上記の式（３）で表されるガンマ確率密度関数の定数倍の関数でのカーブフィッティングによって、次式（４）により表わされる。式（４）において、ｔは、時刻を表し、t0，α，β，Ｋは、いずれもパラメータを表している。αは、障壁の数であり、βは、障壁αの全部を通過する所要時間であり、Ｋは、線形パラメータである。
【数２】
0005591512.tif 000003
【００４０】
動脈の時間濃度曲線のカーブフィッティング部２ｄは、上記の式（４）における４つのパラメータt0，α，β，Ｋの全てを変数とし、次式（５）に示す、動脈カーブモデルと動脈の濃度値ａｉとの２乗誤差を最小にするt0，α，β，Ｋを算出する。実際には、t0には、造影剤が対象臓器に到達する前の特定の時刻を用い、α，β，Ｋを算出する方法を用いる場合が多い。
【数３】
0005591512.tif 000004
【００４１】
動脈の時間濃度曲線のカーブフィッティング部２ｄは、算出したパラメータα，β，Ｋの値を、組織の時間濃度曲線のカーブフィッティング部２ｅに出力する。
【００４２】
組織の時間濃度曲線のカーブフィッティング部２ｅは、造影濃度曲線算出部２ｂから出力される組織の時間濃度曲線Ｃｉを入力するとともに、動脈の時間濃度曲線のカーブフィッティング部２ｄから出力されたパラメータt0，α，β，Ｋの値を入力する。
【００４３】
組織の時間濃度曲線のカーブフィッティング部２ｅは、予め定めた関数を組織の時間濃度曲線のカーブモデル（以下、組織カーブモデルと称する）として適用し、入力した組織の時間濃度曲線Ｃｉのカーブフィッティングを行う組織の時間濃度曲線のカーブフィッティング手段として機能する。
【００４４】
ここで、カーブフィッティングの際、動脈の時間濃度曲線のカーブフィッティング部２ｄから出力されたパラメータt0，α，β，Ｋの値の全てまたは一部を、組織カーブモデルのパラメータとして用いる。なお、組織の時間濃度曲線のカーブフィッティング処理の詳細は、後述する。
【００４５】
組織の時間濃度曲線のカーブフィッティング部２ｅは、カーブフィッティングにより得たパラメータの値を用いて、局所血流動態に関する情報である、平均通過時間MTT（Mean Transit Time）、局所脳血流量CBF（Cerebral Blood Flow）、および局所脳血液量CBV（Cerebral Blood Volume）を算出する血流情報算出手段としても機能する。
【００４６】
次に、図３のフローチャートを参照して、本発明の実施の形態における組織の血流動態に関する情報の算出処理について説明する。この処理では、特に、脳組織を例に挙げ、その血流動態に関する情報を算出する処理について説明する。
【００４７】
ステップＳ１において、時系列データ取得部２ａは、医用画像撮影装置が所定時間撮影して生成した時系列の撮影データを取得する。
【００４８】
ステップＳ２において、造影濃度曲線算出部２ｂは、ステップＳ１の処理で取得された時系列の撮影データから、各画素の画素値の変化を計測して造影濃度曲線を算出する。この造影濃度曲線を、組織の時間濃度曲線Ｃｉとする。
【００４９】
ステップＳ３において、動脈の時間濃度曲線算出部２ｃは、ステップＳ２の処理で算出された組織の時間濃度曲線Ｃｉから、CT画像やMRI画像の解析単位となる動脈領域を設定し、その動脈領域での時間濃度曲線ａiを算出する。
【００５０】
ステップＳ４において、動脈の時間濃度曲線のカーブフィッティング部２ｄは、上記の式（３）で表わしたガンマ確率密度関数の定数倍の関数を、動脈カーブモデルとして適用し、ステップＳ３の処理で算出された動脈の時間濃度曲線ａｉのカーブフィッティングを行う。動脈の時間濃度曲線ａｉのカーブフィッティングによって、t0，α，β，Ｋのパラメータの値が算出される。
【００５１】
ステップＳ５において、組織の時間濃度曲線のカーブフィッティング部２ｅは、予め定めた関数を、組織カーブモデルとして適用し、ステップＳ２の処理で算出された組織の時間濃度曲線Ｃｉのカーブフィッティングを行う。組織の時間濃度曲線Ｃｉのカーブフィッティング処理の詳細は、図４、図７、および図９のフローチャートを参照して後述するが、ステップＳ４の処理で算出されたt0，α，β，Ｋのパラメータの値の全てまたは一部が、組織カーブモデルのパラメータとして用いられる。
【００５２】
そして、組織の時間濃度曲線のカーブフィッティング部２ｅは、組織の時間濃度曲線Ｃｉのカーブフィッティングにより得られたパラメータの値を用いて、血流動態に関する情報である、平均通過時間MTT、局所脳血液量CBV、および局所脳血流量CBFを算出する。
【００５３】
以上のように、血流動態解析装置１は、動脈の時間濃度曲線ａｉをガンマ確率密度関数の定数倍の関数を用いてカーブフィッティングし、そのカーブフィッティングにより得られたパラメータの値の全てまたは一部を、組織の時間濃度曲線Ｃｉのカーブフィッティングのパラメータに用いることによって、安定かつ正確に、平均通過時間MTT、局所脳血液量CBV、および局所脳血流量CBFを算出することが可能となる。
【００５４】
また、血流動態解析装置１は、算出した平均通過時間MTT、局所脳血液量CBV、および局所脳血流量CBFに基づいて、対応する画像を表示部１ｈに表示させることができ、操作者は、血流動態の状態を容易に把握することが可能となる。
【００５５】
次に、図４のフローチャートを参照して、図３のステップＳ５における、第１の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細について説明する。
【００５６】
第１の実施の形態では、図５に示すように、組織の時間濃度曲線のカーブモデルに不完全ガンマ関数の差を用いる点に特徴を有する。図５において、横軸は時間を表わし、縦軸は濃度を表わしており、ａ（ｔ）は、ガンマ確率密度関数の定数倍の関数でカーブフィッティングした動脈の時間濃度曲線を表わし、ｃ（ｔ）は、不完全ガンマ関数の差でカーブフィッティングした組織の時間濃度曲線を表わし、ｆ（ｔ）は、組織のインパルス応答関数を表わしている。
【００５７】
ステップＳ１１において、組織の時間濃度曲線のカーブフィッティング部２ｅは、次式（６）で表わされる不完全ガンマ関数の差、または、式（６）の変形である次式（７）を組織カーブモデルとして適用する。なお、実際の組織の時間濃度曲線Ｃｉのカーブフィッティングには、式（７）を微分した次式（８）を適用することもできる。
【数４】
0005591512.tif 000005
【００５８】
ステップＳ１２において、組織の時間濃度曲線のカーブフィッティング部２ｅは、ステップＳ１１の処理で適用したカーブモデルを用いて、組織の時間濃度曲線Ｃｉのカーブフィッティングを行う。なお、組織カーブモデルの複数のパラメータのうち、α，β，Ｋのパラメータは、図３のステップＳ４の処理で算出された値を用い、組織の時間濃度曲線Ｃｉのカーブフィッティングにより新たに決定するパラメータは、ｔｄ、α１、Ｋ４である。ｔｄは、理想濃度曲線からの遅れ時間であり、α１は、障壁の数であり、Ｋ４は、線形パラメータである。
【００５９】
また、カーブフィッティングの高速算出方法として、多項式ｙをカーブフィッティングに用いることができる。例えば、組織の時間濃度曲線のピーク時刻をｔｐとすると、次式（９）で表わすことができる。
【数５】
0005591512.tif 000006
【００６０】
また、上記の式（９）を展開すると、次式（１０）を得ることができる。
【数６】
0005591512.tif 000007
【００６１】
そして、上記の式（１０）から、ｔｐ−ｔｄ−ｔ０とα１の間の関係を多項式で近似することができる。
【００６２】
図６は、横軸ｔｐ−ｔｄ−ｔ０と縦軸α１の関係の一例を示すグラフである。図６の例の場合、ｔｐ−ｔｄ−ｔ０とα１の間の関係を、次式（１１）で近似できることを示している。
【数７】
0005591512.tif 000008
【００６３】
すなわち、上記の式（１１）を用いて、組織の時間濃度曲線Ｃｉのピーク時刻ｔｐから、極めて容易に、α１を求めることができることを示している。Ｋ４は、線形パラメータのため、最適化法の適用は不要である。
【００６４】
このように、ｔｐ−ｔｄ−ｔ０とα１の間の関係を近似した多項式を用いて、高速に組織の時間濃度曲線Ｃｉをカーブフィッティングすることができる。
【００６５】
また、カーブフィッティングの他の算出方法として、上記の式（１１）で求めたパラメータα１の値を用いて、次式（１２）をカーブフィッティングに適用することもできる。
【数８】
0005591512.tif 000009
【００６６】
上記の式（１２）で決定しなければならないパラメータは、ｔｄのみになり、１変数の最適化であるため、高速に値を決定することができる。
【００６７】
ステップＳ１３において、組織の時間濃度曲線のカーブフィッティング部２ｅは、ステップＳ１２の処理で得られたパラメータｔｄ、α１、Ｋ４の値、および、図３のステップＳ４の処理で得られたパラメータt0，α，β，Ｋの値から、次式（１３）に示すようにして、平均通過時間MTTを算出する。
【００６８】
ＭＴＴ＝α１・β                        ・・・（１３）
ステップＳ１４において、組織の時間濃度曲線のカーブフィッティング部２ｅは、次式（１４）に示すように、組織の時間濃度曲線Ｃｉのピーク部のカーブ下の面積（AUC：Area Under the Curve）ＡＣｃを、動脈の時間濃度曲線ａiのピーク部のカーブ下の面積ＡＣａで除算して局所脳血液量CBVを算出する。
【００６９】
ＣＢＶ＝ＡＣｃ／ＡＣａ                  ・・・（１４）
ステップＳ１５において、組織の時間濃度曲線のカーブフィッティング部２ｅは、次式（１５）に示すように、局所脳血液量CBVを、平均通過時間MTTで除算して局所脳血流量CBFを算出する。
【００７０】
ＣＢＦ＝ＣＢＶ／ＭＴＴ                    ・・・（１５）
なお、以上の処理において、ステップＳ１３とステップＳ１４の処理は順不同である。
【００７１】
以上のように、第１の実施の形態によれば、ガンマ確率密度関数を組織カーブモデルとして適用し、動脈の時間濃度曲線ａｉのカーブフィッティングより得られたパラメータを用いて、組織の時間濃度曲線Ｃｉのカーブフィッティング（近似）を行うようにすることで、安定かつ正確に、組織の局所脳血流に関する情報である、局所脳血液量CBV、局所脳血流量CBF、および平均通過時間MTTを算出することが可能となる。本実施例では、組織の入力関数（カーブフィッティングされた動脈の時間濃度曲線）と出力関数（インパルス応答関数）が双方ともガンマ確率密度関数で表わせると仮定しており、式（６）や式（７）は、造影剤が組織に入り、そして出る間に組織にとどまっている造影剤の濃度を表わすこととなり、物理的意味が明確であり、より正確に解析を行うことができる。
【００７２】
次に、図７のフローチャートを参照して、図３のステップＳ５における、第２の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細について説明する。
【００７３】
第２の実施の形態では、図８に示すように、組織のインパルス応答ｆ（ｔ）に矩形（BOX）関数を用いる点に特徴を有する。図８において、横軸は時間を表わし、縦軸は濃度を表わしており、ａ（ｔ）は、ガンマ確率密度関数の定数倍の関数でカーブフィッティングした動脈の時間濃度曲線を表わし、ｃ（ｔ）は、組織の時間濃度曲線のカーブモデルでカーブフィッティングした組織の時間濃度曲線を表わしている。
【００７４】
ここで、カーブフィッティングした動脈の時間濃度曲線を次式（１６）で表わし、矩形関数（インパルス応答）を次式（１７）で表わすと、それらのコンボリューション（畳み込み積分）は、次式（１８）で表わすことができる。
【数９】
0005591512.tif 000010
【００７５】
ステップＳ２１において、組織の時間濃度曲線のカーブフィッティング部２ｅは、上記の式（１８）を組織カーブモデルとして適用する。Γ（ａ）は、ガンマ関数であり、Ｐ（ａ，ｘ）は、不完全ガンマ関数である。従って、この組織カーブモデルは、不完全ガンマ関数を含む関数である。
【００７６】
ステップＳ２２において、組織の時間濃度曲線のカーブフィッティング部２ｅは、ステップＳ２１の処理で適用したカーブモデルを用いて、組織の時間濃度曲線Ｃｉのカーブフィッティングを行う。なお、組織カーブモデルの複数のパラメータのうち、α，β，Ｋのパラメータは、図３のステップＳ４の処理で算出された値を用い、組織の時間濃度曲線Ｃｉのカーブフィッティングにより新たに決定するパラメータは、ｕ、ｗ、ｔｄである。
【００７７】
例えば、矩形関数の幅が大きくないと仮定した場合、上記の式（１８）は、次式（１９）に示すように変形することができる。
【数１０】
0005591512.tif 000011
【００７８】
上記の式（１９）において、ａ（ｔ−ｔｄ）を動脈の時間濃度曲線ａｉに、Ｆ（ｔ）を組織の時間濃度曲線Ｃｉに置き換えると、次式（２０）に示す連立方程式を得ることができる。パラメータｕ、ｗは、線形演算で解くことができ、パラメータｔｄは、黄金比分割法などの最適化法を用いて解くことができる。
【数１１】
0005591512.tif 000012
【００７９】
ステップＳ２３において、組織の時間濃度曲線のカーブフィッティング部２ｅは、ステップＳ２２の処理で得られたパラメータｕ、ｗの値から、次式（２１）に示すようにして、平均通過時間MTT、局所脳血液量CBV、および局所脳血流量CBFを算出する。
【００８０】
ＭＴＴ＝ｗ
ＣＢＶ＝ｕ
ＣＢＦ＝ｕ／ｗ                          ・・・（２１）
以上のように、第２の実施の形態によれば、カーブフィッティングされた動脈の時間濃度曲線と矩形関数をコンボリューションした関数を、組織カーブモデルとして適用し、動脈の時間濃度曲線ａｉのカーブフィッティングより得られたパラメータを用いて、組織の時間濃度曲線Ｃｉのカーブフィッティングを行うようにすることで、安定かつ正確に、組織の局所脳血流動態に関する情報である、局所脳血液量CBV、局所脳血流量CBF、および平均通過時間MTTを算出することが可能となる。
【００８１】
次に、図９のフローチャートを参照して、図３のステップＳ５における、第３の実施の形態での組織の時間濃度曲線のカーブフィッティング処理の詳細について説明する。
【００８２】
第３の実施の形態では、図１０に示すように、組織のインパルス応答ｆ（ｔ）にガンマ確率密度関数の定数倍の関数を用いる点に特徴を有する。図１０において、横軸は時間を表わし、縦軸は濃度を表わしており、ａ（ｔ）は、ガンマ確率密度関数の定数倍の関数でカーブフィッティングした動脈の時間濃度曲線を表わし、ｃ（ｔ）は、組織の時間濃度曲線のカーブモデルでカーブフィッティングした組織の時間濃度曲線を表わしている。
【００８３】
ここで、カーブフィッティングした動脈の時間濃度曲線を次式（２２）で表わし、ガンマ確率密度関数の定数倍の関数（インパルス応答）を次式（２３）で表わすと、それらのコンボリューションは、次式（２４）で表わすことができる。
【数１２】
0005591512.tif 000013
【００８４】
ステップＳ３１において、組織の時間濃度曲線のカーブフィッティング部２ｅは、上記の式（２４）を組織カーブモデルとして適用する。この組織カーブモデルは、動脈の時間濃度曲線ａｉのカーブフィッティングに適用する動脈カーブモデルと同じである。
【００８５】
ステップＳ３２において、組織の時間濃度曲線のカーブフィッティング部２ｅは、ステップＳ２１の処理で適用したカーブモデルを用いて、組織の時間濃度曲線Ｃｉのカーブフィッティングを行う。なお、組織カーブモデルの複数のパラメータのうち、α，β，Ｋのパラメータは、図３のステップＳ４の処理で算出された値を用い、組織の時間濃度曲線Ｃｉのカーブフィッティングにより新たに決定するパラメータは、ｔｄ、α１、Ｋ２である。
【００８６】
ステップＳ３３において、組織の時間濃度曲線のカーブフィッティング部２ｅは、ステップＳ３２の処理で得られたパラメータα１の値、および、図３のステップＳ４の処理で得られたパラメータα，β，Ｋの値から、次式（２５）に示すようにして、平均通過時間MTT、局所脳血液量CBV、および局所脳血流量CBFを算出する。
【数１３】
0005591512.tif 000014
【００８７】
なお、局所脳血流量CBFは、次式（２６）の関係を用いて算出するようにしてもよい。
【００８８】
ＣＢＦ＝ＣＢＶ／ＭＴＴ                    ・・・（２６）
以上のように、第３の実施の形態によれば、カーブフィッティングされた動脈の時間濃度曲線とガンマ確率密度関数をコンボリューションした関数を、組織カーブモデルとして適用し、動脈の時間濃度曲線ａｉのカーブフィッティングより得られたパラメータを用いて、組織の時間濃度曲線Ｃｉのカーブフィッティングを行うようにすることで、安定かつ正確に、組織の局所血流動態に関する情報である、局所脳血液量CBV、局所脳血流量CBF、および平均通過時間MTTを算出することが可能となる。
【００８９】
図１１は、従来のデコンボリューション法と本発明の方法とで数値実験を行い、平均通過時間MTTを算出した結果の一例を示している。横軸は、組織の時間濃度曲線の数値モデルに与えた通過時間を表し、縦軸は、平均通過時間MTTの算出値を表している。
【００９０】
従来のデコンボリューション法では、平均通過時間MTTの値が小さい領域において、実際の平均通過時間MTTより大きな値が得られており、組織の時間濃度曲線の数値モデルに与えた通過時間と算出結果の線形性があまり良いとは言えない。つまり、時間濃度曲線のデコンボリューションは、数値的に不安定であり、測定ノイズが拡大されて結果に重畳されることが知られている。
【００９１】
そこで、結果の安定化を図るため、様々な技術が使われているが、結果の正確性を重視するか、あるいは、安定性を重視するかはトレードオフの関係にあり、デコンボリューションの理論では、両者を同時に向上させることができない。
【００９２】
一方、本実施の形態によれば、数値的なデコンボリューションを行うことなく、動脈の時間濃度曲線ａｉをカーブフィッティングすることによって得たパラメータを利用して、組織の時間濃度曲線Ｃｉのカーブフィッティングを行い、その結果得たパラメータの値を用いて、平均通過時間MTTおよび局所脳血液量CBVを算出することができ、パラメータ算出の安定性の向上を期待することができる。
【００９３】
従って、算出された平均通過時間MTTは、理論値に近い値が得られており、デコンボリューション法より線形性が改善されていることがわかる。また、CBF＝CBV／MTTの関係から、局所脳血流量CBFについても、正確かつ安定した値を得ることができる。
【００９４】
さらに、動脈の時間濃度曲線ａｉのカーブフィッティングより得られたパラメータを、組織の時間濃度曲線Ｃｉのカーブフィッティングに用いることにより、動脈の時間濃度曲線の依存症を小さくすることができる。従って、造影剤の注入時間を長くすることが可能になり、患者への負担を軽減することができる。
【００９５】
以上においては、対象組織として、脳組織を例に挙げ説明したが、本実施の形態では、脳組織に限られるものではなく、肝臓、心臓、あるいは肺などの他の組織を対象とすることも勿論可能である。
【００９６】
なお、肝臓や肺など、呼吸による動きがある組織、あるいは、心臓では、組織の動き補正を行う必要がある。また動き補正だけでなく、ヘマトクリット値（一定量の血液中に含まれる赤血球の割合）の補正、脳の比重による補正、その他の補正を適宜行う必要があるが、本実施の形態では、特に言及しない。
【００９７】
この発明は、上記実施の形態そのままに限定されるものではなく、実施段階ではその要旨を逸脱しない範囲で構成要素を変形して具体化したり、上記実施の形態に開示されている複数の構成要素を適宜組み合わせたりすることにより種々の発明を形成できる。例えば、実施の形態に示される全構成要素から幾つかの構成要素を削除してもよい。さらに、異なる実施の形態に亘る構成要素を適宜組み合わせても良い。
【符号の説明】
【００９８】
１  血流動態解析装置
２  局所血流解析アプリケーション
２ａ  局所血流解析アプリケーション
２ｂ  造影濃度曲線算出部
２ｃ  動脈の時間濃度曲線算出部
２ｄ  動脈の時間濃度曲線のカーブフィッティング部
２ｅ  組織の時間濃度曲線のカーブフィッティング部
【図１】
0005591512.tif 000015
【図２】
0005591512.tif 000016
【図３】
0005591512.tif 000017
【図４】
0005591512.tif 000018
【図５】
0005591512.tif 000019
【図６】
0005591512.tif 000020
【図７】
0005591512.tif 000021
【図８】
0005591512.tif 000022
【図９】
0005591512.tif 000023
【図１０】
0005591512.tif 000024
【図１１】
0005591512.tif 000025