特許第5773398号(P5773398)IP Force 特許公報掲載プロジェクト 2022.1.31 β版

知財求人 - 知財ポータルサイト「IP Force」

▶ オレア メディカルの特許一覧

特許5773398動脈/組織/静脈の動態系の関心対象量を推定するためのシステムおよび方法
(19)【発行国】日本国特許庁(JP)
(12)【公報種別】特許公報(B2)
(11)【特許番号】5773398
(24)【登録日】2015年7月10日
(45)【発行日】2015年9月2日
(54)【発明の名称】動脈/組織/静脈の動態系の関心対象量を推定するためのシステムおよび方法
(51)【国際特許分類】
   A61B 5/055 20060101AFI20150813BHJP
   A61B 6/03 20060101ALI20150813BHJP
   A61B 5/05 20060101ALI20150813BHJP
【FI】
   A61B5/05 382
   A61B6/03 375
   A61B5/05 383
   A61B5/05ZDM
【請求項の数】14
【全頁数】45
(21)【出願番号】特願2013-533262(P2013-533262)
(86)(22)【出願日】2011年10月11日
(65)【公表番号】特表2013-539704(P2013-539704A)
(43)【公表日】2013年10月28日
(86)【国際出願番号】FR2011052374
(87)【国際公開番号】WO2012049421
(87)【国際公開日】20120419
【審査請求日】2013年6月10日
(31)【優先権主張番号】1058251
(32)【優先日】2010年10月11日
(33)【優先権主張国】FR
(73)【特許権者】
【識別番号】510139128
【氏名又は名称】オレア メディカル
(74)【代理人】
【識別番号】110000970
【氏名又は名称】特許業務法人 楓国際特許事務所
(72)【発明者】
【氏名】パウトー、ファブリース
【審査官】 田中 洋介
(56)【参考文献】
【文献】 A.BENAVOLI et al.,Hard-Constrained versus Soft-Constrained Parameter Estimation,IEEE TRANSACTIONS ON AEROSPACE AND ELECTRONIC SYSTEM,2006年10月,Vol.42 No.4,pp.1224-1239
【文献】 L.Willats et al.,Modelling the Bolus Dispersion from DSC-MRI data,Proceedings of ISMRM,2007年,p.1445
【文献】 VOLKER J. SCHMID ET AL,"Quantitative Analysis of Dynamic Contrast-Enhanced MR Images Based on Bayesian P-Splines",IEEE TRANSACTIONS ON MEDICAL IMAGING,2009年 6月,Vol.28,No.6,p.789-798
【文献】 FRANCESCA ZANDERIGO ET AL,"Nonlinear Stochastic Regularization to Characterize Tissue Residue Function in Bolus-Tracking MRI: Assessment and Comparison With SVD,Block-Circulant SVD,and Tikhonov",IEEE TRANSACTIONS ON BIOMEDICAL ENGINEERING,2009年 5月,Vol.56,No.5,p.1287-1297
(58)【調査した分野】(Int.Cl.,DB名)
A61B 5/055
A61B 5/00
A61B 6/00−6/14
A61B 8/00− 8/15
G06T 1/00−7/60
JSTPlus(JDreamIII)
JMEDPlus(JDreamIII)
(57)【特許請求の範囲】
【請求項1】
器官のボクセルと呼ばれる単位容積において複数の動脈/組織/静脈の動態系の中から関心対象量(14)を推定(57)するために、灌流強調イメージング解析システムの処理ユニット(4)によって実行され、実験による灌流データ(15)から前記関心対象量を推定するステップを含んでいる方法であって、前記ステップが、
前記ボクセルにおける動脈/組織/静脈の動態系の関心対象量の推定に関係するパラメータに鑑みた実験による灌流データの直接確率分布の割り当て(54)を行うこと、および
前記動態系のインパルス応答に関する完全不確実情報の導入(51b)を行ない、最大エントロピーの原理を適用することで、前記量の結合事前確率分布の割り当て(53)を行なうことによって、前記関心対象量の事後周辺確率分布(56)をベイズ法に従って評価することを含むことを特徴とする方法。
【請求項2】
器官のボクセルと呼ばれる単位容積において複数の動脈/組織/静脈の動態系の中から関心対象量(14)を推定(57)するために、灌流強調イメージング解析システムの処理ユニット(4)によって実行され、実験による灌流データ(15)から前記関心対象量を推定するステップを含んでおり、前記動態系は、線形であり、時不変であり、関係C(t)=BF・C(t)*R(t)によって形式的に決定され、ここでC(t)はボクセルを循環する造影剤の濃度であり、C(t)は前記ボクセルに流れ込む動脈における前記造影剤の濃度であり、BFは前記ボクセルにおける血流量であり、*は畳み込み積を表わし、R(t)は前記ボクセルにおける移動時間の相補累積分布関数である方法であって、前記ステップが、
前記ボクセルにおける動脈/組織/静脈の動態系の関心対象量の推定に関係するパラメータに鑑みた実験による灌流データの直接確率分布の割り当て(54)を行うこと、および
前記ボクセルにおける相補累積分布関数R(t)に関する完全不確実情報の導入(51b)を行ない、最大エントロピーの原理を適用することで、前記量の結合事前確率分布の割り当て(53)を行なうことによって、前記関心対象量の事後周辺確率分布(56)をベイズ法に従って評価することを含むことを特徴とする方法。
【請求項3】
前記量の結合事前確率分布の割り当てが、ボクセルに流れ込む動脈における前記造影剤の濃度−時間曲線に関する完全不確実情報の導入(51c)によってさらに達成されることを特徴とする請求項2に記載の方法。
【請求項4】
該当の複数のボクセルについての逐次の反復によって実行されることを特徴とする請求項1〜3のいずれか一項に記載の方法。
【請求項5】
推定した関心対象量に関する信頼区間によって表わされる補足の情報を計算するステップ(58)を含むことを特徴とする請求項1〜4のいずれか一項に記載の方法。
【請求項6】
推定した関心対象量に関する信頼区間への賭け率によって表わされる補足の情報を計算するステップ(59)を含むことを特徴とする請求項1〜5のいずれか一項に記載の方法。
【請求項7】
該当のボクセルにおける動脈/組織/静脈の動態系の関心対象量の推定の問題に関係するパラメータに鑑みた実験による灌流データの直接確率分布の割り当て(54)、および
前記量の結合事前確率分布の割り当て(53)の産物の適切さの指標によって表わされる補足の情報を計算するステップ(60)を含むことを特徴とする請求項4に記載の方法。
【請求項8】
推定した関心対象量(14)を、該関心対象量をユーザ(6)へと提供することができるヒューマン・マシン・インターフェイス(5)へともたらすステップを含むことを特徴とする請求項4〜7のいずれか一項に記載の方法。
【請求項9】
前記推定した関心対象量に関する補足の情報を、該情報をユーザ(6)へと提供することができるヒューマン・マシン・インターフェイス(5)へともたらすステップを含むことを特徴とする請求項5〜8のいずれか一項に記載の方法。
【請求項10】
実験による灌流データが、実験による灌流信号の値のベクトルまたは該ベクトルの濃度−時間曲線の値のベクトルへの変換を含むことを特徴とする請求項1〜9のいずれか一項に記載の方法。
【請求項11】
記憶手段と、外界との通信手段と、処理手段とを備える処理ユニット(4)であって、
通信手段が、外界から実験による灌流データ(15)を受け取ることができ、
処理手段が、器官のボクセルと呼ばれる単位容積において複数の動脈/組織/静脈の動態系の中から関心対象量(14)を推定するための請求項1〜10のいずれか一項に記載の方法を実行するように構成されていることを特徴とする処理ユニット。
【請求項12】
通信手段が、推定された関心対象量(14)を、該関心対象量をユーザ(6)へと提供することができるヒューマン・マシン・インターフェイス(5)に適したフォーマットでもたらす(57)ことを特徴とする請求項11に記載の処理ユニット。
【請求項13】
通信手段が、推定された関心対象量(14)に関する補足の情報を、該情報をユーザ(6)へと提供することができるヒューマン・マシン・インターフェイス(5)に適したフォーマットでもたらす(58、59)ことを特徴とする請求項12に記載の処理ユニット。
【請求項14】
請求項11〜13のいずれか一項に記載の処理ユニット(4)と、前記処理ユニット(4)によって実行される請求項1〜10のいずれか一項に記載の方法に従って、推定された量(14)をユーザ(6)へと提供することができるヒューマン・マシン・インターフェイス(5)とを備える灌流強調イメージング解析システム。
【発明の詳細な説明】
【技術分野】
【0001】
本発明は、灌流強調イメージング(perfusion−weighted imaging)へと確定的でない確率論的手法(soft probabilistic methods)を適用することによって血行動態パラメータを推定するためのシステムおよび方法に関する。さらに、そのような方法は、相補累積分布関数または動脈入力関数の推定を可能にし、その後に、より一般的に、任意の関心対象量の推定を可能にする。本発明は、とりわけ、必要とされる推定を恣意的な望ましくない仮定によって制約または行使することがなく、不確実な事前の生理学的または血行動態的な情報の導入を必要とする点で、公知の方法から相違する。
【背景技術】
【0002】
本発明は、詳細には、灌流強調磁気共鳴イメージング(PW−MRI)または灌流コンピュータ断層撮影(CT)に依存する。それらの技術は、脳または心臓などの器官の血行動態についての有用な情報を迅速に取得することを可能にする。この情報は、急性脳梗塞などの病変の救急処置において診断および治療上の決定を行おうとする医師を助けるためにとくに重要である。
【0003】
それらの技術は、核磁気共鳴装置またはコンピュータ断層撮影装置に依存する。この装置は、脳などの体の一部分の複数のデジタル画像系列をもたらす。この目的のため、前記装置は、高周波数の電磁波の組み合わせを体の前記部分へと照射し、特定の原子によって再放出される信号を測定する。このやり方で、この装置は、画像化された塊の各点(すなわち、ボクセル)の化学組成を明らかにし、したがって生物学的組織の種類を明らかにすることを可能にする。
【0004】
画像系列は、専用処理ユニットによって解析される。この処理ユニットは、最終的には、適切なヒューマン・マシン・インターフェイスによって、灌流強調画像からの血行動態パラメータの推定を医師へともたらす。このやり方で、医師は、診断を行い、どの治療決定が適切であるかを判断することができる。
【0005】
核磁気共鳴またはコンピュータ断層撮影灌流強調画像は、造影剤(例えば、磁気共鳴イメージングにおけるガドリニウムキレート)を静脈に注入し、画像の各ボクセルにおいて長時間にわたってこの造影剤のボーラスを記録することによって取得される。簡潔にするため、ボクセルを識別するための添え字x,y,zを省略する。例えば、座標x,y,zのボクセルの信号を、Sx,y,z(t)と表わす代わりに、単にS(t)と表わす。後述される操作および演算が、推定すべき血行動態パラメータを表わす画像またはマップを最終的に取得するために、着目される各ボクセルについて一般的に実行されることを、理解すべきである。
【0006】
標準モデルが、或る時間tにわたって測定された信号強度S(t)を、前記造影剤の濃度C(t)に関連付けることを可能にする。
【0007】
例えば、灌流コンピュータ断層撮影では、各ボクセルの信号が濃度に正比例し、すなわちS(t)=k・C(t)+Sである。核磁気共鳴による灌流イメージングにおいては、指数関数の関係S(t)=S・e−k・TE・C(t)が存在する。どちらの場合も、Sが、造影剤の到達前の平均の信号強度を表わしている。核磁気共鳴イメージングに関しては、kが、組織における常磁性磁化率と造影剤の濃度との間の関係に依存する定数であり、TEが、エコー時間である。各ボクセルの定数kの値は、未知であるため、着目されるすべてのボクセルについて任意の値に設定される。したがって、絶対的な推定ではなく、相対的な推定が得られる。しかしながら、この相対的な情報は、空間におけるそれらの値の相対的な変化(詳細には、正常な組織と病んだ組織との間での違い)が主たる関心の対象であるため、依然として有意義である。
【0008】
一般に、理論上の信号S(t)を造影剤の理論上の濃度C(t)へと結び付けるモデルが、S(t)=Ψ(C(t),Θ)と称され、Θが、このモデルの自由パラメータのベクトルである。例えば、磁気共鳴またはコンピュータ断層撮影による灌流強調イメージングにおいては、Θ=(S)である。
【0009】
各々の時点における各ボクセルに囲まれた組織領域の造影剤の質量の保存が、
【0010】
【数1】
【0011】
と表わされる。C(t)が、この組織領域へと送り込まれる動脈の造影剤の濃度であり、動脈入力関数またはAIFとして知られている。BFが、この組織領域における血流量であり、C(t)が、この組織領域から流れ出る静脈の造影剤の濃度であり、静脈出力関数またはVOFとして知られている。
【0012】
動脈/組織/静脈の動態系が、線形であり、時間変化しないと仮定すると、
【0013】
【数2】
【0014】
であり、ここでh(t)は、系のインパルス応答または組織における造影剤の移動時間の確率密度関数であり、「x」を○で囲んだ符号(以下、「*」で表すことがある)は、畳み込み積を表わす。したがって、初期条件C(t=0)=0における上述の微分方程式の形式解は、C(t)=BF・C(t)*R(t)と表わされ、ここでR(t)は、
【0015】
【数3】
【0016】
と定義される相補累積分布関数または剰余関数であり、ここでHは、ヘビサイド単位ステップ一般化関数である。インパルス応答および相補累積分布関数から、別の血行動態パラメータ、すなわち組織における平均通過時間またはMTTが定められる。
【0017】
【数4】
【0018】
また、血液容量BVを関係BV=BF・MTTによって定めることもできる。
【0019】
核磁気共鳴による灌流強調イメージングにおいて、BF、MTT、またはBVなどの血行動態パラメータ、ならびに相補累積分布関数は、現在のところ以下のように推定されている。
【0020】
各々のボクセルについて、時点t,i=1,Nにおいてサンプリングされた実験による灌流信号Sexp(t)が、関係
【0021】
【数5】
【0022】
を使用することによって濃度−時間曲線C(t)へと変換される。定数kが、すべてのボクセルについて、ゼロでない任意の値(例えば、k・TE=1)に設定される。定数Sが、例えば造影剤の到着前の平均を取ることによって推定される。これが、灌流信号の取得が造影剤のボーラス到着時間(すなわち、BAT)と比べて充分に早期に開始される場合にかぎって可能であることに、注意すべきである。濃度C(t)から、関連の理論上の動脈入力関数C(t)が既知であると仮定し、積BF・R(t)が、数値逆畳み込みによって推定される。
【0023】
濃度−時間曲線C(t)の逆畳み込みのための理論上の動脈入力関数C(t)を得るために、いくつかの手法が提案されている。
【0024】
第1の手法においては、医師が、実験による大域的な動脈入力関数を手作業で選択する。これを、脳の灌流イメージングにおいて、例えば対側性シルビウス動脈(contralateral sylvian artery)または内頸動脈において測定することができ、あるいは例えば光学的な測定など、追加の測定から得ることができる。高い信号対雑音比の信号を得ることができるにせよ、この手法は、多くの欠点を有している。第1に、人間の介入および/または追加の測定が必要である。これは、臨床の緊急の状況において望ましくなく、手順および最終的な結果の再現を、より困難にする。第2に、とりわけ、この大域的動脈入力関数は、各ボクセルの局所的動脈入力関数に一致しない。局所的動脈入力関数からの相違が、遅延(局所的動脈入力関数は、一般に血管系の上流において取得される大域的動脈入力関数と比べて遅れるため)および分散(造影剤の伝播が、血管系の上流よりも下流において遅いため)に関して存在する。ここで、畳み込み積の対称性により、これらの欠点が相補累積分布関数の推定に直接影響するため、これらの現象が、最終的に、血行動態パラメータの推定にかなりの影響を有することが既知である。したがって、例えば、局所的動脈入力関数と局所的静脈出力関数との間の真の平均通過時間(MTT)の推定を最終的に得ることができず、大域的動脈入力関数と静脈出力関数との間の平均通過時間しか得ることができない。この乖離を克服するために、一部の者は、大域的動脈入力関数と局所的動脈入力関数との間の遅延を定量化する
【0025】
【数6】
【0026】
というパラメータなどの新たな記述パラメータを、たとえそれらが(動脈入力関数が各ボクセルの真の局所的動脈入力関数である)元の標準的な灌流モデルに属さないにせよ、導入している。他の方法は、局所的動脈入力関数のこれらの乖離について、血行動態パラメータの推定への影響を最小限にしようと試みている。しかしながら、それらは大域的問題に新たな未知の要素を導入し、それを回避しているにすぎない。
【0027】
第2の手法によれば、大域的動脈入力関数が、データクラスタリングまたは独立成分分析(ICA)などの信号処理技術によって灌流強調画像から自動的に取得される。この手法は、人間の介在を回避できるようにするかもしれないが、大域的動脈入力関数につきものの遅延および分散の問題を解決しておらず、新たな未知の要素を導入している(例えば、動脈入力関数の代わりに静脈出力関数を得ることができる)。
【0028】
第3の手法によれば、局所的動脈入力関数が、信号処理技術および選択基準によって灌流強調画像から自動的に取得される。例えば、<<最良の>>関数が、血行動態パラメータまたは相補累積分布関数を推定すべき現在の組織ボクセルのすぐ近傍において捜索される。この第3の手法の目的は、遅延および分散の問題を少なくとも或る程度は克服することによって、より偏りが少なく、より正確な推定を、最終的に得ることである。しかしながら、このようにして得られた局所的動脈入力関数が対象のボクセルの<<真の局所的関数>>の適切な近似であることが、演繹的にも帰納的にも保証されない。例えば、この<<真の>>関数が、該当の近傍の外側に位置する(該当の近傍が小さすぎる場合)可能性があり、反対に、別の動脈入力関数と混同される(該当の近傍が大きすぎる場合)可能性がある。さらに、この<<最良の>>局所的動脈入力関数は、<<正常な>>動脈入力関数(すなわち、短い/高い早さの造影剤到着時間を有する、大きな振幅を有する、など、)の中から求められる。しかしながら、目的は、例えば虚血性の動脈入力関数など、病気の動脈入力関数から、正常な動脈入力関数を正確に区別することである。結果として、たとえ最終的な結果が大域的な手法よりも良好になりうる場合でも、局所的な動脈入力関数についての不確かさが依然として大きく、血行動態パラメータまたは相補累積分布関数についての不確かさはなおさらである。
【0029】
上述の方法から得られるとおりの理論上の動脈入力関数C(t)によって実験による濃度−時間曲線C(t)の逆畳み込みを行なうために、標準的な畳み込みモデルC(t)=BF・C(t)*R(t)が、最初に例えば矩形法の近似に従って時間について離散化される。
【0030】
【数7】
【0031】
ここで、Δtはサンプリング周期である。このやり方で、下記のように線形システムAd=cに帰着する。
【0032】
【数8】
【0033】
実際には、行列Aが、深刻に悪い状態であり、ほとんど特異であるため、無意味な解および異常な推定を得ることなく線形システムを数値的に解くことはできない。したがって、例えば行列Aの疑似逆行列A〜−1を得、次いで
【0034】
【数9】
【0035】
によってdの推定d^を得るために、種々の方法に頼らなければならない。逆行列を得るためのそれらの方法として、sSVD法(単純特異値分解)、cSVD法(円形特異値分解)、およびoSVD法(振動特異値分解)などのAの特異値の打ち切りにもとづく方法(打ち切り特異値分解または(T)SVD)や、周波数域におけるHunt逆畳み込みが挙げられる。
【0036】
より一般的には、
【0037】
【数10】
【0038】
などの基準を最小にすることができ、ここで
【0039】
【数11】
【0040】
は、特定の解を有利にし、
【0041】
【数12】
【0042】
によってdの推定を得ることを可能にする正則化項である。そのような方法として、Tikhonov正則化およびウェーブレット変換にもとづく方法などが挙げられる。
【0043】
ひとたびd^が得られると、定義により
【0044】
【数13】
【0045】
であるため、BF^=d^(t)=d^(0)によってBFの推定BF^を得ることができる。しかしながら、BFは、多くの場合、例えば特異値分解にもとづく方法において、それらの方法に固有のd^(0)の系統的な過小評価を補償するために、
【0046】
【数14】
【0047】
として推定される。次いで、
【0048】
【数15】
【0049】
によってbの推定b^が得られ、例えば矩形法の近似に従うことにより、
【0050】
【数16】
【0051】
によって
【0052】
【数17】
【0053】
の推定が得られる。最後に、典型的には、たとえ積の推定量が推定量の積でなくても、
【0054】
【数18】
【0055】
によってBVの推定が得られる。
【0056】
これらの動脈入力関数にもとづく方法について、多数の変種を見つけることができ、例えば信号対雑音比を人工的に高めるために、実験による動脈入力関数を、最初にパラメトリックまたはセミパラメトリック理論モデルC(t,Θ)(ここで、Θはパラメータのベクトルである)へと適合させることができる。同様に、数値逆畳み込みをより安定にし、あるいは造影剤の循環信号(第1の通過)と再循環信号(第2の通過)との間に時間的な重なり合いが存在する場合に、再循環に起因して生じうる問題を克服するために、信号を人工的にオーバーサンプリングすることができる。しかしながら、それらの方法も、依然として打ち切り特異値分解などの数値逆畳み込み法にもとづいている。
【0057】
典型的な実データを模擬するように意図された合成データについて試験した種々の方法についてのいくつかの比較研究によれば、特異値分解にもとづく逆畳み込みおよびその変種(線形打ち切り特異値分解(sSVD)、打ち切り円形特異値分解(cSVD)、または平滑化打ち切り円形特異値分解(oSVD))が、実務において患者、組織の種類、病状、などに応じて生じうる種々の相補累積分布関数R(t)および動脈入力関数C(t)に対して、偏り(真の値に対する系統誤差)、精度(入力信号の信号対雑音比の関数としての真の値に対する推定の標準偏差)、およびロバストさに関してBFおよびMTTパラメータの最良の推定を最終的にもたらす。
【0058】
しかしながら、上述した動脈入力関数の選択に関する問題に加えて、この種の数値法は、激しい固有の問題を抱えている。
【0059】
第一に、BF・bの推定d^が、時間とともに減少するのではなく、ときには負の値をとりうるような程度にまで振動する。しかしながら、時刻tにおいてボクセルに残る造影剤の量R(t)は、当然ながら正であり、減少関数である。oSVDなどの特別な方法は、それらの異常な振動を軽減することができるが、それらは特異値分解に固有のものであるがゆえに依然として残る。これが、この種の方法において、血流量BFが
【0060】
【数19】
【0061】
によって推定される一方で、標準的な灌流モデルにおいては
【0062】
【数20】
【0063】
であるがゆえに、
【0064】
【数21】
【0065】
として推定されるべき理由である。始点における値の代わりに最大値をとることによって、それらの振動のBFの推定に対する影響だけを消去することが期待される。したがって、これらの方法は、完璧に満足できるものにはなり得ない。詳細には、相補累積分布関数、したがって血行動態パラメータのより厳密な推定方法によって到達できるBFの推定の精度は、未知のままである。したがって、これらの数値逆畳み込み法は、モデルの特性を満たさない解をもたらすため、標準的な灌流モデルと矛盾する。
【0066】
この問題を克服し、相補累積分布関数の生理学的に容認できる推定を得るために、それらの相補累積分布関数のパラメトリックモデルR(t,Θ)が導入されており、ここでΘは、このモデルのパラメータのベクトルである。これらのモデルが、例えばベイズ法を使用することによって実験による信号に合わせられる。しかしながら、この手法は、現時点においては時期尚早であるように見受けられる。実際、相補累積分布関数を記述するために適したパラメータまたはセミパラメトリックモデルを決定するために、相補累積分布関数の事前の非パラメータ的推定が思いのままでなければならない。なぜならば、それらのモデルが、実務において生じうるあらゆる種類の相補累積分布関数を正しく表わすことに完璧に適していないと、得られるMTTまたはBFなどの血行動態パラメータの推定が異常になることが、モンテカルロシミュレーションによって示されているからである。したがって、相補累積分布関数の理論上のモデルの選択が重要であり、モデルを実験データへと適用することによってのみ正しく行なうことが可能である。したがって、前記相補累積分布関数を推定するための非パラメータ法が、おそらくは次の段階における古典的なパラメータまたはセミパラメトリックモデルによる置き換えによって相補累積分布関数の生理学的に容認できる推定を得ることを可能にするために、必要である。
【0067】
しかしながら、上述のように、特異値分解にもとづく方法などの方法によって得られる推定は、標準的な灌流モデルのもとでの
【0068】
【数22】
【0069】
の定義そのものに矛盾する。それらは、生理学的にも、物理学的にも、容認することができない。したがって、パラメータまたはセミパラメトリック理論モデルをそれらの推定へと当てはめることが不可能である。さらに、いくつかのモデルを互いに比較し、実験による相補累積分布関数を表わすために最も適したモデルを選択することも、不可能である。したがって、特異値分解法などの数値法につきものの欠陥ゆえに、灌流現象のモデル化および理解における進歩はもはや不可能である。
【0070】
さらに、相補累積分布関数および血行動態パラメータの非パラメータ的推定の根底にある問題は、実験に基づく動脈入力関数による、実験に基づく濃度−時間曲線の<<単純な>>逆畳み込みの問題ではないように思われる。実際、絶対的な確かさおよび無限の精度で知られた動脈入力関数が問題において実際に与えられるならば、そのようであるかもしれない。しかしながら、最良でも、最大でも測定雑音までしか知ることができず、前もって濃度−時間曲線へと変換しなければならない実験的に測定された動脈信号が与えられるにすぎない。換言すると、実験的に測定された動脈入力関数が理論上の動脈入力関数に等しいと仮定することによって、実験による信号における測定雑音ならびにそれらの信号からの濃度−時間曲線の推定または変換に起因する不確かさは無視される。
【0071】
標準的な畳み込みモデルC(t)=BF・C(t)*R(t)は、直接的に測定することが一般に不可能である理論上の信号だけを含む。実際、測定される信号は、典型的には、理論上の信号と測定雑音との合計である。したがって、実験による組織および動脈の灌流信号は、それぞれSexp(t)=Sth(t)+ξおよびSaexp(t)=Sath(t)+ξと表わされ、ここでξおよびξは、測定雑音をモデル化する平均ゼロの確率過程である。さらに
【0072】
【数23】
【0073】
および
【0074】
【数24】
【0075】
であり、ここでCth(t)は、磁気共鳴灌流強調イメージングにおける理論上の濃度−時間曲線であり、Cath(t)は、理論上の動脈入力関数であり、あるいは、より一般的には、上述のようにSth(t)=Ψ(Cth(t),Θ)である。したがって、実験による信号へと適用される理論上の標準的な灌流モデルは、Cth(t)=BF・Cath(t)*R(t)と表わされなければならず、すなわち核磁気共鳴灌流強調イメージングにおいては以下のように表わされなければならず、
【0076】
【数25】
【0077】
あるいは同等には、下記のように表わされなければならず、
【0078】
【数26】
【0079】
下記のようにではない。
【0080】
【数27】
【0081】
ここで、Cexp(t)は、実験による濃度−時間曲線であり、Caexp(t)は、実験による動脈入力関数であり、逆畳み込み法の大部分は、この誤りのある間接的な数学的畳み込みモデルにもとづいている。この時点で、必ずしも正でないランダム変数の対数をとることを避けるために、標準的な灌流モデルを
【0082】
【数28】
【0083】
と表わすことが好ましいことを、理解できるであろう。
【0084】
逆畳み込みされる信号の測定の不確かさが、逆畳み込みプロセスの最終的な結果に大きく影響することが知られており、それらの不確かさに起因する入力信号の非常に小さい変動が、最終的な結果のかなりの変動につながりかねない。まさにこれらの問題を克服し、このような不安定さを軽減するために、Tikhonov正則化などの逆畳み込み法または特異値分解にもとづく方法が導入されている。さらに、現在のところ完全に無視されている測定の雑音および不確かさの動脈入力関数への影響が、さらにより顕著であり、動脈の測定雑音ξおよびS0aへの不確かさが、今や畳み込み行列Aにおいて生じ、結果として増幅され、伝播する。動脈入力関数への誤差および不確かさを無視することは、血行動態パラメータの推定に重大な誤差を引き起こし、それらの推定に精度の錯覚も引き起こす。いくつかの方法は、現在まで避けられたままであるこの問題を最小限にするために、動脈入力関数における測定雑音の消去を目的としている。推定の誤差を支配し、定量化できるよう、動脈入力関数の不確かさを血行動態パラメータおよび相補累積分布関数の推定に伝播させることができる方法を思いのままにできることが、好ましいと考えられる。
【0085】
さらに、標準的な灌流モデルをCth(t)=BF・Cath(t)*R(t)と表わすことを、避けることが好ましいと考えられる。実際、標準的な灌流モデルが、最初に動脈/組織/静脈の動態系のインパルス応答h(t)を定め、そこから相補累積分布関数R(t)が、
【0086】
【数29】
【0087】
と計算される。したがって、標準的な灌流モデルをインパルス応答h(t)の関数として
【0088】
【数30】
【0089】
と表わすことが、例えば、矩形法による近似
【0090】
【数31】
【0091】
のように相補累積分布関数の推定R^を最終的に得るために、h^によって測定時点t,j=1,Nにおけるインパルス応答h(t)を推定すべくこのモデルを当てはめるために、好都合である。詳細には、インパルス応答h(t)を相補累積分布関数の推定R^から推定しない方がよい。数値的な観点から、最初に導関数(h(t))を推定し、第2の段階において不定積分(R(t))を推定することが好ましく、逆ではないことを、容易に理解できるであろう。しかしながら、相補累積分布関数R(t)がインパルス応答h(t)の代わりに直接推定されることを、確認することができる。
【0092】
さらに、逆畳み込みの問題は、不良設定であり、無限に多くの可能解の余地がある。さらに、血行動態パラメータ、インパルス応答、相補累積分布関数、または動脈入力関数の推定に関して灌流強調イメージングにおいて直面される不良設定問題は、測定雑音の伝播の問題および実験による信号の濃度−時間曲線への変換の問題によって3重にされる逆畳み込みの問題である。
【0093】
おそらくは1つの解だけしか認めない良設定問題へと帰着させるために、先験的に可能な無限の解の中から求められる解への事前の情報、制約を加える必要がある。これが、特定の種類の事前情報を程度の差こそあるが明示的に、あるいは程度の差こそあるが直接的に注入する多数の逆畳み込みおよび推定の方法を、技術水準において見つけることができる理由である。
【0094】
例として、上述のとおりの古典的なTikhonov正則化を挙げることが可能である。その最も一般的なバージョンにおいては、行列ΓがΓ=αIに等しい。これは、大きなユークリッドノルムの解を不利にし、スカラーαが、この制約におけるこの不利の重みの大きさを表わしている。
【0095】
他方で、灌流強調イメージングにおいて古典的である打ち切り特異値分解にもとづく方法は、正則化パラメータの役割を果たす所与のしきい値を下回る畳み込み行列の特異値をキャンセルすることからなる。それらの小さな特異値は、信号の高周波成分に関係しており、したがって特異値の打ち切りは、低域通過フィルタとして機能する。しかしながら、他方では、相補累積分布関数は、t=0での不連続性ゆえに、無限の帯域幅の信号である(相補累積分布関数が
【0096】
【数32】
【0097】
について定義されていることを思い出されたい)。他方で、小さな特異値の打ち切りは、追加の生理学的情報に直接対応するだけでなく、とりわけ追加の代数的情報、すなわち所与の線形部分空間の帰属関係に対応する。高すぎる周波数の振動が残ることが確認される。詳細には、動脈入力関数による相補累積分布関数の逆畳み込みによって推定された理論上の信号は、測定雑音に追従する傾向を有する(過剰適合)。さらに、打ち切りのしきい値が、明示的な生理学的制約の重みを意味せず(解を何らかの線形部分空間に属するように強いることを意味するにすぎない)、その効果が複雑で直接的でなく、さまざまな信号の種類に応じて最も適切な値を決定するための一般的基準を与えることが難しい。反対に、基準に対する特異値の切り捨ての作用が、きわめて間接的であるため、例えばoSVD法における相補累積分布関数の粗さなど、そのような基準の適切な最適化の可能性が保証されない。
【0098】
さらに、現在の血行動態パラメータまたは相補累積分布関数の推定方法は、それらの推定の精度(すなわち、推定量の標準偏差)の評価を可能にしておらず、それらの推定およびそれらの推定の精度の評価における信頼度の評価も可能にしていない。詳細には、逆畳み込みおよび推定の方法による標準的な灌流モデルの適合度の定量化が困難である。したがって、打ち切り特異値分解にもとづく方法が、実験による信号の過剰適合の傾向を有することに、注目することができる。推定による信号が滑らかでなく、むしろ反対に、測定雑音に従う傾向を有する。結果として、古典的なχの統計値(すなわち、二乗誤差の合計)などのモデルの適合度の指標に関して、低い値が得られる可能性がある。そのような低い値は、高品質の適合を示しているが、過剰適合ゆえに事実ではない。したがって、それらの統計値は、この場合には不適当であって誤解につながりかねず、過剰適合の検出および考慮を望む場合には避けるべきである。したがって、技術水準による逆畳み込みおよび血行動態パラメータの推定の方法に適合度の<<良好な>>指標を導入することは、自明ではない。したがって、そのような適合度の指標が一意に定められ、過剰適合を考慮に入れ、標準的な灌流モデルの適合度だけをデータへと定量化することができる方法を思いのままにすることが望まれる。
【0099】
要約すると、血行動態パラメータ、インパルス応答、または相補累積分布関数を推定するための現在の方法は、多数の方法論的な欠陥および制御できない多数の近似を免れない。したがって、結果の解釈が難しくなり、灌流強調イメージングが充分に活用されていない。
【発明の概要】
【0100】
本発明の目的は、公知の方法の使用から生じるすべての欠点に回答することにある。本発明の主たる目的は、必ずしも満足されず、とりわけ実験的に検証することが不可能であり、あるいは生理学的または血行動態的な種類ではないその場しのぎの制約を導入することなく、生理学的または血行動態的な種類の特定の制約を満たす解の中からインパルス応答hの推定h^を探し出すことを可能にする新規な方法を提供することにある。さらに、それらの方法は、特別な方法に頼ることを必要とせずに、そのような制約の重みを自動的かつ独自のやり方で決定することを可能にする。
【0101】
したがって、本発明は、血行動態パラメータ、相補累積分布関数、または動脈入力関数の推定の問題が最終的に良設定であることを主張することを可能にする。本発明は、その独特な解決策(おそらくは複数)を提供する。
【0102】
本発明によってもたらされる主要な利点の中でも、以下が可能であることを挙げることができる(ただし、これらに限られるわけではない)。
【0103】
・血行動態パラメータ、インパルス応答、相補累積分布関数、または動脈入力関数についての定性的または半定性的で不確実な生理学的情報を定量的なやり方で変換することができる。
【0104】
・動脈入力関数および実験による信号における不確かさおよび誤差を明示的に考慮に入れ、関心対象量の推定に伝播させることができる。
【0105】
・それらのパラメータ、インパルス応答、相補累積分布関数、または動脈入力関数の推定を、偏り(系統誤差)、精度(統計誤差)、線形性、および実務において生じうる種々の状況に対するロバストさに関して、改善することができる。
【0106】
・血行動態パラメータの推定、インパルス応答、相補累積分布関数、または動脈入力関数について、それらの推定において有することができる信頼を改善および明確化するために、信頼区間を得ることができ、それらの信頼区間についての賭け率さえ得ることができる。
【0107】
・インパルス応答、相補累積分布関数、または動脈入力関数について、生理学的および血行動態的な観点から容認することができ、標準的な灌流モデルに準拠している非パラメータ的な推定を得ることができ、したがって次の段階において、これらのインパルス応答、相補累積分布関数、または動脈入力関数について、パラメトリックまたはセミパラメトリックな理論上のモデルの適合、比較、および最終的には選択を行なうことができる。
【0108】
・動脈入力関数が対象の各ボクセルの組織に関して実際にどの程度まで容認できる動脈入力関数であるかを定量化することができる。
【0109】
・大域的灌流モデルの適合度の客観的かつ定量的な指標をもたらすことができ、最終的には比較および最も適した大域的灌流モデルの選択を可能にすることができる。
【0110】
この目的のため、本発明は、器官のボクセルと呼ばれる単位容積において複数の動脈/組織/静脈の動態系の中から関心対象量を推定するための方法を具現化する。そのような方法は、灌流強調イメージング解析システムの処理ユニットによって実行されるように意図され、実験による灌流データから前記関心対象量を推定するステップを含んでいる。必ずしも満たされず、とりわけ実験的に検証することが不可能であり、あるいは生理学的または血行動態的な種類のものではない特別な制約を導入することなく、特定の生理学的または血行動態的な制約を満たす解の中から前記関心対象量を探し求めるために、本発明による推定のステップは、
・該当のボクセルにおける動脈/組織/静脈の動態系の関心対象量の推定に関係するパラメータに鑑みた灌流データの直接確率分布を割り当て、
・前記動態系のインパルス応答に関する完全不確実情報を導入し、最大エントロピーの原理を適用することで、前記量の結合事前確率分布を割り当てることによって前記関心対象量の事後周辺分布をベイズ法に従って評価することを含む。
【0111】
さらに本発明は、器官のボクセルと呼ばれる単位容積において複数の動脈/組織/静脈の動態系の中から関心対象量を推定するための方法であって、前記動態系が、線形であり、時不変であり、関係C(t)=BF・C(t)*R(t)によって形式的に決定され、ここでC(t)はボクセルを循環する造影剤の濃度であり、C(t)は前記ボクセルに流れ込む動脈における前記造影剤の濃度であり、BFは前記ボクセルにおける血流であり、*は畳み込み積を表わし、R(t)は前記ボクセルにおける移動時間の相補累積分布関数である、方法を具現化する。そのような方法は、灌流強調イメージング解析システムの処理ユニットによって実行されるように意図され、実験による灌流データから前記関心対象量を推定するステップを含んでいる。本発明によれば、上述したように、前記推定のステップが、
・該当のボクセルにおける動脈/組織/静脈の動態系の関心対象量の推定に関係するパラメータに鑑みた灌流データの直接確率分布を割り当て、
・前記ボクセルにおける相補累積分布関数R(t)に関する完全不確実情報を導入し、最大エントロピーの原理を適用することで、前記量の結合事前確率分布を割り当てることによって前記関心対象量の事後周辺分布をベイズ法に従って評価することを含む。
【0112】
好ましい実施の形態によれば、両方の場合において、本発明は、ボクセルに流れ込む動脈における前記造影剤の濃度−時間曲線に関する完全不確実情報の導入による前記量の結合事前確率分布の割り当てをさらに企画する。
【0113】
好都合には、本発明による方法を、該当の複数のボクセルについての逐次の反復によって実行することができる。
【0114】
関心対象量の推定の精度を評価するために、本発明による方法は、推定した関心対象量に関する信頼区間によって表わされる補足の情報を計算するステップを含むことができる。
【0115】
好ましい実施の形態によれば、そのような方法が、推定した関心対象量に関する信頼区間への賭け率によって表わされる補足の情報を計算するステップをさらに含むことができる。
【0116】
さらに、本発明による方法は、
・該当のボクセルにおける動脈/組織/静脈の動態系の関心対象量の推定の問題に関係するパラメータに鑑みた実験による灌流データの直接確率分布の割り当て、および
・前記量の結合事前確率分布の割り当ての産物の適切さの指標によって表される補足の情報を計算するステップを含むことができる。
【0117】
推定を提供するために、本発明による方法は、推定した関心対象量を、該関心対象量をユーザへと提供することができるヒューマン・マシン・インターフェイスへともたらすステップを含むことができる。
【0118】
補足の情報が計算される場合には、そのようなプロセスが、前記関心対象量に関する補足の情報を、該情報をユーザへと提供することができるヒューマン・マシン・インターフェイスへともたらすステップをさらに含むことができる。
【0119】
好ましい実施の態様によれば、本発明において、実験による灌流データが、実験による灌流信号の値のベクトルまたは該ベクトルの濃度−時間曲線の値のベクトルへの変換を含む。
【0120】
さらに本発明は、第2の目的に従い、記憶手段と、外界との通信手段と、処理手段とを備える処理ユニットに関する。通信手段が、外界から実験による灌流データを受け取ることができる。処理手段が、器官のボクセルと呼ばれる単位容積において複数の動脈/組織/静脈の動態系の中から関心対象量を推定する本発明による方法を実行するように構成されている。
【0121】
本発明による方法によって生成された推定をユーザに提供するために、そのような処理ユニットの通信手段は、推定された関心対象量を、該関心対象量をユーザへと提供することができるヒューマン・マシン・インターフェイスに適したフォーマットでもたらすことができる。
【0122】
この実施の態様によれば、さらに通信手段が、推定された関心対象量に関する補足の情報を、該情報をユーザへと提供することができるヒューマン・マシン・インターフェイスに適したフォーマットでもたらすことができる。
【0123】
第3の目的によれば、さらに本発明は、本発明による処理ユニットと、前記処理ユニットによって実行される本発明による方法に従って推定された量をユーザへと提供することができるヒューマン・マシン・インターフェイスとを備える灌流強調イメージング解析システムに関する。
【0124】
他の特徴および利点は、以下の説明を検討し、添付の図面を点検することによって、さらに明らかになる。
【図面の簡単な説明】
【0125】
図1図1は、灌流強調画像を解析するためのシステムの実施の形態を示している。
図2図2は、灌流強調画像を解析するためのシステムの実施の形態を示している。
図3図3は、造影剤の注入前における、核磁気共鳴イメージング装置によって得られた人間の脳のスライスの灌流画像を示している。
図4図4は、脳組織内での造影剤の循環の最中における、核磁気共鳴イメージング装置によって得られた人間の脳のスライスの灌流画像を示している。
図5図5は、人間の脳のボクセルに関する核磁気共鳴による灌流強調信号S(t)を示している。
図6図6は、人間の脳のボクセルを循環する造影剤の典型的な濃度−時間曲線C(t)を示している。
図7図7は、典型的な動脈入力関数C(t)を示している。
図8図8は、本発明による方法を示している。
図9図9は、本発明に従って推定された脳血液量のマップを示している。
図10図10は、本発明に従って推定された脳虚血の場合における脳血流のマップを示している
図11図11は、本発明に従って推定された平均通過時間のマップを示している。
図12図12は、脳血流が信頼区間に属する確率のマップを示している。
【発明を実施するための形態】
【0126】
図1が、灌流強調画像の解析のためのシステムを示している。核磁気共鳴イメージング装置またはコンピュータ断層撮影装置1が、コンソール2によって制御される。ユーザが、装置1を制御するためのパラメータ11を選択することができる。装置1によって生成される情報10から、人間または動物の体の一部分の複数の一連のデジタル画像12が得られる。好ましい例として、人間の脳の観察からのデジタル画像を使用して、先行技術からの解決策および本発明を例証する。他の器官も考えられる。
【0127】
一連の画像12を、随意によりサーバ3に保存し、患者の医療記録13を構成してもよい。そのような記録13は、灌流強調または拡散強調画像など、さまざまな種類の画像を含むことができる。一連の画像12を、専用の処理ユニット4を使用して解析することができる。そのような処理ユニットは、画像を集めるべく外界と通信するための手段を備える。さらに、この通信のための手段は、処理ユニットが専用のヒューマン・マシン・インターフェイス5によって医師6または研究者に灌流強調画像12からの血行動態パラメータ14の推定を最終的にもたらすことを可能にする。分析システムのユーザ6は、診断の承認または却下、自身が適切であると考える治療行為の決定、さらなる検査の実行、などを行なうことができる。随意により、ユーザは、設定16を通じて処理ユニット4の動作を設定することができる。例えば、表示のしきい値を定めることができ、あるいは閲覧を希望する推定パラメータを選択することができる。
【0128】
図2が、前処理ユニット7が一連の画像12を解析して各ボクセルの灌流データ15を取り出す解析システムの実施の形態を示している。このようにして、血行動態パラメータ14の推定を担当する処理ユニット4が、この動作から解放され、外界との通信のための手段によって受信される灌流データ15からの推定方法を実行する。
【0129】
図3が、人間の脳の厚さ5mmのスライスの画像12の典型例を示している。この画像は、核磁気共鳴によって得られている。この技術を使用して、各々のスライスについて、1.5×1.5×5mmの寸法を有する128×128のボクセルの行列を得ることができる。双一次補間を使用して、画像20などの458×458ピクセルのフラット画像を生成することができる。
【0130】
図4が、図3に関して示した画像と同様の画像20を示している。しかしながら、この画像は、造影剤の注入後に得られている。この画像は、灌流強調脳画像の典型例である。図3に示した同じ画像と対照的に、動脈がはっきりと現れている。公知の技術に従い、血行動態パラメータを推定するために、異常のある半球の対側にある半球内の1つ以上の動脈入力関数21を選択することが可能である。
【0131】
図5bが、図2に関連して説明した前処理ユニット7によってもたらされるデータ15などの核磁気共鳴による灌流強調信号S(t)の例を示している。したがって、灌流強調信号は、造影剤の注入後のボクセルの時間tにつれての変化を表わしている。例えば、図5bは、そのような信号を50秒の期間にわたって示している。縦座標は、信号強度を任意単位にて示している。そのような信号を得るために、図1による処理ユニット4(あるいは、図2による前処理ユニット7)が、例えば図5aに示されるような時点t,t,…,t,…,tにおける核磁気共鳴による一連のn枚の灌流強調画像I1,I2,…,Ii,…,Inを解析する。所与のボクセル(例えば、ボクセルV)について、造影剤の注入後のボクセルの時間tにつれての変化を表わす灌流強調信号S(t)が決定される。
【0132】
図6が、図5bに示したとおりの灌流強調信号から導出された濃度−時間曲線を示している。すでに述べたように、灌流強調信号と関連の濃度−時間曲線との間には、関係が存在する。すなわち、核磁気共鳴による灌流強調イメージングにおいては、指数関数の関係S(t)=S・e−k・TE・C(t)が存在し、ここでSは、造影剤の到着前の平均の信号強度であり、TEは、エコー時間であり、kは、常磁性磁化率と組織における造影剤の濃度との間の関係に依存する定数である。
【0133】
図6が、ボクセル内の造影剤の濃度の時間変化の観察を可能にする。ボクセルにおける造影剤の1回目の通過において高い振幅のピークが存在し、その後に造影剤の再循環(2回目の通過)の現象に関係したより低い振幅のピークが存在することを、見て取ることができる。
【0134】
図7が、図4に関連して示したボクセル21などの動脈ボクセルにおける造影剤の循環を表わす典型的な動脈入力関数C(t)を示している。図7は、詳細には、造影剤の1回目の通過の後の再循環の現象がきわめて弱いことを示している。
【0135】
図8が、器官内のボクセルの複数の動脈/組織/静脈の動態系において関心対象量を推定するための本発明による方法の例を示している。そのような方法を、図1および2に関連して説明した相応に構成された灌流強調イメージング解析システムの処理ユニットによって実行することができる。
【0136】
本発明による方法は、主として、血行動態パラメータ、サンプリング時点における理論上のインパルス応答の値、または相補累積分布関数の値などの推定すべき種々の関心対象量について、1つ以上の事後周辺分布を割り当てるステップ56を含む。さらに、前記推定を計算するステップ57も含む。
【0137】
そのような事後周辺分布を割り当てるために、処理ユニットの設定プロセス50を行なうことが必要である。処理ユニット自身が、好ましくは、1つ以上の設定セッティングによって、この設定を実行することができる。設定は、1つ以上の事後周辺分布からなるライブラリの構築ももたらすことができ、このライブラリは、あらかじめ確立されて前記ユニットのプログラムメモリに保存される。本発明によれば、前記ライブラリを、外部の処理ユニットによって使用または提供されるときに充実したものにすることができる。この外部の処理ユニットは、前記設定セッティングから前記設定を実行することができ、前記ライブラリを出力するために処理ユニットと協働することができる。
【0138】
したがって、本発明による方法は、割り当て56の前に実行される設定ステップを含むことができ、その中でも、以下が必要かつ充分である。
【0139】
・対象ボクセルの動脈/組織/静脈の動態系の関心対象量の推定の問題に関係するすべてのパラメータを鑑みた実験データへの直接確率分布の割り当て54。
【0140】
・それらすべてのパラメータについての結合事前確率分布の割り当て53。
【0141】
本発明は、さまざまな応用の事例において1つ以上の関心対象量の推定を可能にする。
【0142】
・理論上の動脈入力関数が、絶対的な確かさおよび無限の精度で知られると仮定される。
【0143】
・局所的動脈入力関数は、やはり推定されるべき未知の遅延にて大域的動脈入力関数から相違する。
【0144】
・動脈入力信号だけが、おそらくは遅延まで測定される。
【0145】
・動脈入力関数が与えられず、動脈入力信号が測定されない。これが最も現実的な事例である。
【0146】
設定ステップは、該当の応用の事例に依存することができる。
【0147】
図8が、動脈入力信号がおそらくは遅延τまで測定されるときの第1の応用の例における本発明による方法を示している。
【0148】
実験の灌流モデルMを、以下のように表わすことができる。
【0149】
【数33】
【0150】
b=[R(t),…,R(t)]は、未知の相補累積分布関数の値のベクトルである。すでに説明したように、このベクトルを、インパルス応答h=[h(t)=0,…,h(t)]の値のベクトルから表現することができる。例えば、左矩形近似法によれば、
【0151】
【数34】
【0152】
であり、右矩形近似法によれば、
【0153】
【数35】
【0154】
であり、
c=[c(t),…,c(t)]は、ボクセルにおける造影剤の未知の理論上の濃度の値のベクトルであり、
s=[S(t),…,S(t)]は、灌流強調イメージングによる信号強度の実験による測定値の実または複素ベクトルであり、
Θは、S(t)をC(t)に結び付けるパラメータのベクトルであり、
ξ=[ξ(t),…,ξ(t)]は、前記測定信号における測定雑音のベクトルであり、
=[S(t),…,S(t)]は、実験による動脈入力信号S(t)の測定値のベクトルであり、
ΘSaは、S(t)をC(t)に結び付けるパラメータのベクトルであり、
ξ=[ξ(t),…,ξ(t)]は、前記動脈信号における測定雑音のベクトルであり、
a=[C(t),…,C(t)]は、未知の理論上の動脈入力関数C(t)の値のベクトルであり、
Aは、矩形近似法、台形規則、あるいはSimpson法、Boole法、Gauss−Legendre法などのより高次の方法などによって畳み込み積分を数値的に近似することによって得られる未知の理論上の動脈入力関数aの値のベクトルに関する畳み込み行列である。Aは、cSVDおよびoSVDなどのような特異値の打ち切りにもとづく逆畳み込み法において使用されるような循環畳み込み行列であってもよい。
【0155】
モデルMの血行動態パラメータのベクトルをΘと称し、例えばここではΘ=(BF,τ)である。
【0156】
図1、2、および8に関連し、処理ユニット4によって実行される方法は、実験による信号sの測定値のベクトルについての情報Iおよび実験による動脈入力信号sの測定値のベクトルについての情報ISaをそれぞれ導入することからなる2つの初期設定ステップ51dおよび51eを含むことができる。
【0157】
例えば、実数値の測定雑音ベクトルの組の最初の2つのモーメント
【0158】
【数36】
【0159】
および
【0160】
【数37】
【0161】
に専念すると、最大エントロピー(Lebesgue基準測度のもとでの微分Shannonエントロピー)の原理により、ベクトルξおよびξが、標準偏差σおよびσSaをそれぞれ有する相互に独立な白色の定常ガウス確率過程と見なされる必要がある。
【0162】
より一般的には、測定雑音ベクトルの組(ξ,ξ)を特徴付けるパラメータを(E,ESa)と称する。例えば、(E,ESa)=(σ,σSa)である。
【0163】
本発明による方法は、ベクトルaおよびb、ベクトルΘ、パラメータのベクトルΘおよびΘSa、ならびにベクトルの組(E,ESa)に鑑みて測定ベクトルの組(s,s)に結合直接確率分布を割り当てる設定ステップ54を含む。したがって、この結合直接確率分布は、p(s,s|a,b,Θ,Θ,ΘSa,E,ESa,I,ISa,M)と表わされる。
【0164】
例えば、
【0165】
【数38】
【0166】
である(ξおよびξを、標準偏差σおよびσSaをそれぞれ有する相互に独立な白色の定常ガウス確率過程と考えた場合)。あるいは、本発明は、所望であれば、前記結合直接確率分布を別なやり方で、未知の理論上の相補累積分布関数の値のベクトルbにて直接にではなく、インパルス応答の値のベクトルhにて表現することを可能にできる。この目的のために、bをhに関して表現するだけでよい。したがって、結合直接確率分布54が、例えば下記のように表わされる。
【0167】
【数39】
【0168】
一般性を失うことなく、上述した理由でベクトルhが優先的に注目されることを後に考慮すべきである。
【0169】
同じやり方で、複素信号の測定値のベクトルの組についての結合直接確率分布を、それらの実部および虚部の直接確率分布を乗算することによって割り当てることもできる。
【0170】
さらに、本発明は、例えば核磁気共鳴による灌流強調イメージングにおいて
【0171】
【数40】
【0172】
を使用することによって灌流信号sおよびsの値のベクトルを濃度cexpの値のベクトルへと正確に変換できるときに変種を提供する。例えば、造影剤の到着前の平均の灌流信号強度Sを正確かつ他のパラメータとは別個独立に推定できるとき、そのようにすることが可能である。この実施の形態によれば、他のすべてのパラメータに鑑みたベクトルの組(cexp,cexpa)についてのガウス確率分布を下記のように割り当てることが依然として正当であることを示すことができ、
【0173】
【数41】
【0174】
σおよびσは、今やcexpおよびcexpaのそれぞれの測定雑音の標準偏差である。
【0175】
より一般的には、実験による灌流信号sの値のベクトルならびにその濃度−時間曲線cexpの値のベクトルへの変換を指して、用語<<実験灌流データ>>が使用される。その後に、一般性を失うことなく、実験灌流データをsおよびsと称する。
【0176】
さらに、この方法は、モデルMの血行動態パラメータΘのついての情報IΘ、インパルス応答hについての情報I、および動脈入力関数aについての情報Iをそれぞれ導入することからなる3つの設定ステップ51a、51b、および51cを含む。設定ステップ51a〜51eにおいて導入される情報は、本発明によれば、処理ユニットを設定し、すなわち処理ユニットが事後周辺分布の割り当て56を行なうことができ、次いで関心対象量の推定57を行なうことができるようにするための設定パラメータを構成する。それらの情報または設定パラメータから、本発明による方法は、p(a,h,Θ,Θ,ΘSa,E,ESa|IΘ,I,I,I,ISa,M)と表わすことができる結合事前確率分布を割り当てるステップ53を含むことができる。
【0177】
最大エントロピーの原理を適用し、この分布を、下記のように典型的に因数分解することができる。
【0178】
【数42】
【0179】
したがって、そのような方法は、事前確率分布p(Θ|IΘ,M)を割り当てることからなるステップ52aを含む。
【0180】
例えば、無情報事前分布を割り当てることができる。例えば、情報IΘが、BFおよびτがそれぞれ区間[BFmin,BFmax]および[τmin,τmax]に属することを知ることだけからなる場合、事前確率分布p(Θ|IΘ,M)を、下記のように表わすことができる。
【0181】
【数43】
【0182】
他端において、過去の実験から得られ、例えば陽電子放出断層撮影(PET)または動脈スピンラベリング(ASL)などの定量的な灌流強調イメージング技術から得られる相対頻度サンプリング分布または事後周辺確率分布などの情報確率分布を割り当てることも可能である。
【0183】
本発明による方法は、情報I(あるいは、相補累積分布関数だけが関心の対象である場合には、I)からの事前確率分布
【0184】
【数44】
【0185】
を割り当てることからなる設定ステップ52bをさらに含む。この情報は、確実および不確実情報で作られる。
【0186】
本発明によれば、確実情報が、不確実情報から区別される。確実情報は、確かであると考えられ、すなわち確率が1に等しいあらゆるブール命題(Boolean proposition)に相当する。例えば、<<この曲線は滑らかである>>または<<この信号はこのモデルに従う>>などのブール命題は、確実情報を構成する。対照的に、不確実情報は、ブール命題を示すことからなるが、そのようなブール命題を特定の確率で示すにすぎないあらゆるブール命題に関する。最後に、これは、<<この曲線はほぼ滑らかである>>または<<この信号はこのモデルにほぼ従う>>などのブール命題を導入することを意味する。以後、確実でないあらゆる不確実情報を<<完全不確実情報>>と呼ぶ。
【0187】
確実情報と不確実情報との間の相違を強調するために、最も単純かつ最も極端な事例を検討する。ガウス分布N(x,σ)に従って独立かつ同一に分布したサンプルx,…,xから実際の量xを推定したいと仮定する。実数aを設定する。事前の確実な定量的情報<<x=a>>を考える。この確実情報は、ディラック事前確率分布p(x|a)=δ(x−a)へと変換され、すなわちx=aであればp(x|a)=1であり、x≠aであればp(x|a)=0である。
【0188】
次に、対応する不確実な定量的情報、すなわち<<xがaにほぼ近い>>、つまりは
【0189】
【数45】
【0190】
を考える。
【0191】
最大エントロピーの原理を適用し、この不確実な定量的情報は、ガウス事前確率分布へと変換される。
【0192】
【数46】
【0193】
確実情報および対応する不確実情報が、異なる事前確率分布へと変換されることを理解できるであろう。詳細には、不確実情報についての事前分布において、事前確率分布p(ε)を割り当てることによって<<ほぼ>>を定量的なやり方で先験的に変換することを可能にする新たな(ハイパー)パラメータεが現れる。詳細には、確実情報は、εを0へと向けることによる不確実情報の限定の事例と理解される。しかしながら、本発明は、ε>0である完全不確実情報だけを取り扱う。
【0194】
確実および対応する完全不確実な事前確率分布は常に異なっており、通常はそれ自身が異なる関心対象量の推定をもたらす。上述の極端な例では、実際に、確実な場合にはx=aが既知であるためxを推定する必要がない。他方で、不確実な場合には、ガウス分布の数学的期待値を推定するという古典的な問題に帰着する。
【0195】
次に、第2の例として、<<yが(セミ)パラメトリックモデルM(x,Θ)に従う>>という古典的な確実な定量的情報を考える(ここで、Θは、Mのパラメータのベクトルであり、yは実数である)。本質的に、この定量的情報は、ディラック事前確率分布へと変換される。
【0196】
【数47】
【0197】
ここで、zがz〜N(y,σ)などの雑音の多い実験データである場合、尤度関数または直接確率分布
【0198】
【数48】
【0199】
が得られ、ここから例えばベイズ規則を適用することによって関心対象量Θおよびσを推定することができる。
【0200】
【数49】
【0201】
次に、<<yが(セミ)パラメトリックモデルM(x,Θ)にほぼ従う>>という対応する不確実な定量的情報を検討する。この不確実情報は、以下のように変換される。
【0202】
【数50】
【0203】
最大エントロピーの原理を適用し、この不確実情報は、下記の事前確率分布へと変換される。
【0204】
【数51】
【0205】
上述のようにz〜N(y,σ)である場合、尤度関数または直接確率分布
【0206】
【数52】
【0207】
が得られ、ここから例えばベイズ規則を適用することによって関心対象量を推定することができる。
【0208】
【数53】
【0209】
この第2の例から、確実情報および対応する完全不確実情報が、第1の例のように事前確率分布のレベルにおいてのみならず、尤度関数のレベルにおいても相違しうることを、見て取ることができる。やはり、一般に、関心対象量の異なる推定が最終的に得られる。この例は、本発明に従って(セミ)パラメトリックモデルおよび確定的な血行動態パラメータ推定法を<<不確実に>>するやり方を大まかに示している。
【0210】
定性的情報を説明するための第3の例を検討する。S(t)が実際の信号であり、s=[S(t),…,S(t)]がその値のベクトルであるとする。<<S(t)が滑らかである>>という確実な定性的情報を検討する。この定性的情報が、
【0211】
【数54】
【0212】
という定量的情報へと正規に(canonically)変換され、ここで
【0213】
【数55】
【0214】
は、関数のLノルムを表わす。したがって、S(t)は線形であり、S(t)=at+bなどの(a,b)が存在する。実際、直線よりも滑らかなものは存在しない。この確実情報が、事前確率分布p(S(t)|t,a,b)=δ[S(t)−at−b]へと最終的に変換される。したがって、実務においては、組(a,b)の推定に帰着する。
【0215】
ここで、<<S(t)がほぼ滑らかである>>という対応する完全不確実な定性的情報を考える。これが、
【0216】
【数56】
【0217】
という完全不確実な定量的情報へと変換される。関数空間を扱っているため、この不確実情報を変換する事前確率プロセスを得るために最大エントロピーの原理を直接適用することができないように思われる。しかしながら、
【0218】
【数57】
【0219】
という定量的情報の離散時間バージョンを考えることができ、ここでds/dtは、sからds/dt=Dsなどの二次導関数の有限差分数値近似を使用することによって得られる。
【0220】
したがって、以下のようになる。
【0221】
【数58】
【0222】
次いで、最大エントロピーの原理の適用により、ベクトルsについての多変量ガウス事前確率分布がもたらされる。
【0223】
【数59】
【0224】
やはり、確実情報と完全不確実情報との間の相違を評価することができる。確実情報の場合には、sが直線y=at+b上に位置するように制約される一方で、完全不確実情報の場合には、推定しなければならないε>0の値に応じて直線から離れるように多少なりとも移動する。第1の場合が、2つのパラメータ(a,b)の推定に帰着する一方で、第2の場合は、N+1個のパラメータの推定に帰着する(ベクトルsについてのN個のパラメータおよびεハイパーパラメータ)。
【0225】
本発明による不確実および完全不確実情報の考え方を明確に定義したので、インパルス応答h(t)についての確実および不確実情報に戻る。
【0226】
確実な定量的情報は、h(t)=h(0)=0を含み、例えば右矩形法による
【0227】
【数60】
【0228】
を含む。両方の情報により、推定すべき値/パラメータの数をN−2に減らすことができる。理解できるとおり、賢明な選択は、値h’=[h(t),…,h(tN−1)]を保ち、h(t)を
【0229】
【数61】
【0230】
によって表現することからなる。さらに、確実な定量的情報
【0231】
【数62】
【0232】
も考慮することができる。したがって、それらの確実な定量的情報を純粋に定性的な完全不確実情報と組み合わせることによってベクトルh’の結合事前確率分布を割り当てることが残る。例えば、上述のように、<<h(t)がほぼ滑らかである>>という完全不確実な定性的情報Iを考える。すでに述べたように、この定性的情報Iを、最初に
【0233】
【数63】
【0234】
という完全不確実な定量的情報へと変換することができる。
【0235】
サンプリング時点t=1,Nにおける時間離散化の後で、hの二次導関数の二次数値近似が、例えば下記によって与えられる。
【0236】
【数64】
【0237】
例えば下記の四次の式など、さらに高次の数値近似も使用することが可能である。
【0238】
【数65】
【0239】
それらの数値近似を、行列表示にて
【0240】
【数66】
【0241】
と表わすことができ、ここでDは、例えば二次の近似の場合に次元Nの正方行列
【0242】
【数67】
【0243】
であり、h’は上記にて定めたとおりである。
【0244】
したがって、h’の二次導関数のユークリッドノルムの平方を、
【0245】
【数68】
【0246】
と表わすことができる。
【0247】
したがって、インパルス応答h(t)が多少なりとも滑らかになると定性的に仮定することを、
【0248】
【数69】
【0249】
が多少なりとも大きくなると仮定することによって定量的に変換することができる。したがって、同じユークリッドノルム
【0250】
【数70】
【0251】
を有するすべての連続確率分布の中から境界点h(t)およびh(t)が別々に取り扱われる事前確率分布p(h’|h(t),h(t),I)が探し求められる。ハイパー象限[0,+∞]N−2を台とするそれらのすべての分布の中から最も不確か、したがって最も正当であるがゆえに最大の微分Shannonエントロピー(Lebesgue基準測度のもとで)を有する分布を選択することからなる最大エントロピーの原理の適用により、一定のベクトルの数学的期待値M=(μ,…,μおよび共分散行列ε−1/σを有する[0,+∞]N−2についての条件付き多変量打ち切りガウス分布(あるいは、相補累積分布関数の値のベクトルについて[0,1]N−2で打ち切られた多変量ガウス分布)が得られる。
【0252】
【数71】
【0253】
ここでC(μ,ε,σ)は、正規化定数
【0254】
【数72】
【0255】
である。
【0256】
したがって、2つの新たなハイパーパラメータ、すなわち我々の灌流モデルの正則化パラメータの役割を果たすスカラー数学的期待値μおよび逆分数分散εが、我々の大域的灌流モデルに現れる。εは、実験データsおよびsによってもたらされる確実な定量的情報に対する事前の不確実な定性的情報Iの不確実さを定量化する。
【0257】
さらに、定義により、下記のとおりである。
【0258】
【数73】
【0259】
しかしながら、遅延τがサンプリング時点t,i=1,Nにほとんど等しくないという事実を考慮するために、
【0260】
【数74】
【0261】
とすることができ、ここで
【0262】
【数75】
【0263】
かつC1,1(μ1,1,ε1,1,σ)は正規化定数
【0264】
【数76】
【0265】
であり、
【0266】
【数77】
【0267】
だけが示される。
【0268】
最後に、情報Iについてハイパー象限[0,+∞]を台とするhの事前確率分布を、以下のように表わすことができる。
【0269】
【数78】
【0270】
一般に、例えばE=(μ,ε)またはE=(μ,ε,μ1,1,ε1,1)など、hの事前分布のハイパーパラメータのベクトルがE(または、E)と称される。
【0271】
の事前確率分布を割り当てることが残る。例えば、無情報非正則Bayes−Laplace−Lhoste−Jeffreys分布p(μ,ε|I,M)∝ε−1またはp(μ,ε,μ1,1,ε1,1|I,M)∝ε−1ε1,1−1を割り当てることができる。
【0272】
同じやり方で、
【0273】
【数79】
【0274】
という完全不確実な定量的情報Iを導入でき、より一般的には、
【0275】
【数80】
【0276】
という完全不確実情報I,k=1,2,…を導入することができる。
【0277】
最後に、最大エントロピーの原理を適用することによってhの新たな事前多変量打ち切りガウス確率分布が得られる。
【0278】
上述のような不確実情報I,I,…が唯一の完全不確実情報であり、我々の問題に正当に導入することができる唯一の制約であり、それらが検証できる任意の仮説ではなく、あるいは実験によるものでなく、むしろ対照的に、それらは基本的な生理学的特性の表現にすぎず、それなしでは血行動態パラメータ、インパルス応答、または相補累積分布関数の推定の問題が実際に絶対的に無意味になりかねず、それらが我々の問題にとって論理的に必要かつ充分であることに、注目すべきである。他のあらゆる情報は、対照的に、潜在的に実験によって検証されることができる単純仮説であると考えられる。
【0279】
本発明によれば、それでもやはり、単純作業仮説である他の完全不確実情報を導入することができる。例えば、インパルス応答が所与の関数形式にほぼ従うにすぎない旨を、確実情報によってこの形式に正確に従うように制約するのではなく、示すことができる。この種の半定量的な完全不確実情報の追加は、提案された関数形式でインパルス応答をどの程度まで表わすことができるかを判断できるようにする。したがって、例えば<<h(t)がパラメトリックまたはセミパラメトリックモデルM:f(t,Θ)にほぼ従う>>という上述のような完全不確実情報Iを導入したいと仮定する(ここで、Θはこのモデルのパラメータのベクトルである)。すでに述べたように、そのようなパラメトリックモデルは、例えば下記のように与えられる2パラメータΓモデルなど、すでに提案されている。
【0280】
【数81】
【0281】
この場合に、パラメータMTTを、上述のようにインパルス応答の第1のモーメント(または、相補累積分布関数の積分)を数値的に推定する必要なく、直接推定できることに注目すべきである。
【0282】
hが所与の関数形式f(t,Θ)にほぼ従う旨を示すことは、剰余のベクトルのユークリッドノルム
【0283】
【数82】
【0284】
が多少なりとも大きい旨を定量的に示すことを意味する。やはり最大エントロピーの原理を適用し、同じやり方で、hの事前確率分布がハイパー象限[0,+∞]で打ち切られた多変量分布であることが分かる。
【0285】
【数83】
【0286】
ここでC(ε,σ,Θ)は、正規化定数
【0287】
【数84】
【0288】
である。
【0289】
同じやり方で、<<h(t)が所与の値h=[h(t),…,h(t)]のベクトルにほぼ従う>>などの不確実情報Iを導入することができる。
【0290】
最大エントロピーの原理を適用し、事前確率分布
【0291】
【数85】
【0292】
を得る。
【0293】
ここでC(ε,σ,h)は、正規化定数
【0294】
【数86】
【0295】
である。
【0296】
さらに、本発明は、インパルス応答h(t)(または、相補累積分布関数)についてのいくつかの情報およびそれらの対応する事前確率分布を組み合わせることを可能にする。したがって、p(h|E,I,M),…,p(h|E,I,M)がそれぞれの正則化パラメータE,…,E(例えば、E=(ε,μ)、E=(ε,Θ)、など)で情報I,…,Iを変換するn個の事前確率分布である場合、それらn個の情報を考慮するhの事前分布を、E=(E,…,E)およびI=I∧…∧Iとすることによって
【0297】
【数87】
【0298】
と表わすことができる。
【0299】
情報Iを局所的動脈入力関数へとエンコードするために、本発明による方法は、事前確率分布
【0300】
【数88】
【0301】
を割り当てることからなるステップ52cをさらに含む。そのような分布は、動脈入力関数についての確実および/または不確実情報を導入して組み合わせることによって、インパルス応答hに関する分布と同じやり方で割り当てられる。例えば、動脈入力関数がほぼ滑らかであり、正であり、単峰であり、双峰であり、原点においてゼロであり、ゼロに漸近であり、あるいは所与の範囲を有する旨を指定することができ、もしくは所与のパラメータまたはセミパラメトリックモデルC(t,Θ)にほぼ従うことだけを示すことができ、ここでΘは、例えば11パラメータの<<3ガンマ>>モデルなどの、パラメータのベクトルである。
【0302】
【数89】
【0303】
すでに述べたように、ガウス事前確率分布
【0304】
【数90】
【0305】
を得る。
【0306】
ここでC(ε,σSa,Θ)は、正規化定数
【0307】
【数91】
【0308】
である。
【0309】
同様に、ベクトルaが局所的動脈入力関数などの所与の値のベクトルa=[a(t),…,a(t)]にほぼ近いことだけを示す(そうでない場合、動脈入力関数が与えられる事例に帰着する)ことができ、事前確率分布
【0310】
【数92】
【0311】
を得る、ここでC(ε,σ,a)は、正規化定数
【0312】
【数93】
【0313】
である。
【0314】
本発明による方法は、典型的には以下のように因数分解される事前確率分布p(E,ESa,Θ,ΘSa|I,ISa,M)を割り当てることからなる設定ステップ52dをさらに含む。
【0315】
【数94】
【0316】
例えば、無情報事前確率分布が得られる。
【0317】
【数95】
【0318】
ステップ52b〜52dにおいて導入された種々のハイパーパラメータを考慮し、結合事前確率分布53を、p(a,E,h,E,Θ,Θ,ΘSa,E,ESa|IΘ,I,I,I,ISa,M)と書き直すことができる。
【0319】
以下の表現を簡単にするために、処理ユニットの設定パラメータとして入力された情報一式をI=(IΘ,I,I,I,ISa,M)とする。
【0320】
上述のとおりの直接確率分布54および結合事前確率分布53に鑑み、ベイズ規則を適用することによってすべてのパラメータの結合事後確率分布55が得られる。
【0321】
【数96】
【0322】
初期設定が行なわれると、今や本発明は、ベクトルΞ=(a,E,h,E,Θ,Θ,ΘSa,E,ESa)のすべての要素中からθと称するべき関心対象量の推定を可能にする。例えば、θ=BFまたはθ=τ、あるいはθ=h(t)などである。
【0323】
したがって、本発明による方法は、θの事後周辺分布を評価することからなるステップ56、すなわち
【0324】
【数97】
【0325】
を含む。
【0326】
この事後周辺分布から、θの推定θ^を入手57することができる。例えば、この分布の数学的期待値
【0327】
【数98】
【0328】
をとることによって、二次損失関数
【0329】
【数99】
【0330】
のもとでのベイズ推定量が得られる。同じやり方で、
【0331】
【数100】
【0332】
によって最大事後推定量θ^、すなわちMAPを得ることができる。
【0333】
例として、各々のサンプリング時点tにおけるインパルス応答の値h(t),i=1,Nの事後周辺確率分布を、他のすべての時点を周辺化することによって得ることができ、
【0334】
【数101】
【0335】
その後にh^(t)またはh^(t)などのそれらの値の推定を得ることができる。
【0336】
同じやり方で、必ずしもサンプリング時点tに等しくない任意の時点xにおけるインパルス応答h(x)の値の推定h^(x)(または、R^(x))を入手57することも可能である。この目的のために、上述のような相補累積分布関数R(t)の値の表現にh(x)を導入し、その周辺確率分布を計算すれば充分である。そのような追加の時点x,…,xを推定の問題に導入することは、考慮される時点の数が多いほど
【0337】
【数102】
【0338】
などの積分の数値近似ならびにdh/dtまたはdh/dtの数値近似がより良好になるため、さらに望ましい。その後に、得られる推定もより良好である。
【0339】
したがって、
【0340】
【数103】
【0341】
であるため、例えばここでは右矩形法などの数値積分法を適用することによって、<<平均>>推定
【0342】
【数104】
【0343】
または最確推定
【0344】
【数105】
【0345】
などといったこのパラメータの推定を得ることができる。
【0346】
さらに、本発明による方法は、パラメータθの推定の精度の推定を得、さらにはそれらの推定についての信頼区間を得るステップ58を含むことができる。本発明によれば、そのような方法が、前記信頼区間についての賭け率を得るステップ59をさらに備えることができる。例えば、推定θ^の精度を、θの事後周辺確率分布の共分散行列によって定量化することができる。
【0347】
【数106】
【0348】
したがって、例えばθに対する<<1シグマ>>における信頼(ハイパー)区間
【0349】
【数107】
【0350】
ならびにθがこの区間に属する確率
【0351】
【数108】
【0352】
または同等のオッズp/(1−p)が得られる。
【0353】
また、パラメータBV=BF・MTT、相補累積分布関数の値のベクトルb、静脈出力関数の値のベクトルv=Ah、または理論上の濃度−時間曲線の値のベクトルc=BF・Abについての推定、信頼区間、およびそれら信頼区間への賭け率を得ることも、いくつかのランダム変数の結合確率分布に鑑みてそれらの任意の関数の確率分布関数を計算できるがゆえに可能である。例えば、数学的期待値の線形性によって、相補累積分布関数R(t)の値の推定が、例えば右矩形法に従うことによってインパルス応答
【0354】
【数109】
【0355】
および
【0356】
【数110】
【0357】
のそれらからすぐに得られる。同じやり方で、cおよびΘについての結合確率分布p(c,Θ|s,s,I)から、理論上の灌流信号Sth(t)の値の確率分布sth=[Sth(t),…,Sth(t)]を、sth=Ψ(c,Θ)であるがゆえに導出することもできる。この分布から、上述した方法と同様の方法に従って、推定sth^を入手57することができ、信頼ハイパー区間58およびそれらハイパー区間についての賭け率59も入手することができる。
【0358】
本発明によれば、信頼区間またはそれらの信頼区間の賭け率が、処理ユニットの設定のセッティング62を可能にすることができる。したがって、設定パラメータを変更し、より高品質の推定をもたらすことができる。
【0359】
さらに、本発明による方法は、理論上および実験による灌流強調信号の間の剰余r(t)=S(t)−sth^(t),i=1,Nを計算するステップ60を含むことができる。したがって、本発明は、それらのベクトルの間の種々の統計量または距離D(s,sth^)の計算を可能にし、その最も古典的な1つは、誤差の平方和(SSE)
【0360】
【数111】
【0361】
である。それらの種々の統計量は、灌流モデルの実験データへの適合度を定量化できるようにする。このようにして、関心対象の各ボクセルにおけるモデルの<<誤差マップ>>が得られる。上述の理由(すなわち、過剰適合の問題)で、灌流モデルの実験データsおよびsへの適合度の定量化を、関心対象の各ボクセルの灌流モデルを鑑みて実験データ(s,s)の確率を計算することによってきわめて好都合に達成することができる。
【0362】
【数112】
【0363】
この場合には、誤差マップが1−p(s,s|I)にもとづくと考えられる。
【0364】
さらに本発明によれば、剰余r(t)が実際に独立であるかどうか、等しく分布しているかどうか、ガウスであるかどうか、などをチェックするために、Q−Qプロットまたはポアンカレの帰還写像などの種々の統計テストまたは種々のグラフィカル診断技術の適用61を行なうこともできる。したがって、本発明によれば、反復プロセスおよび試行錯誤によって、灌流現象のモデル化、理解、および処理を進歩させるとともに、血行動態パラメータ、インパルス応答、相補累積分布関数、動脈入力関数、または静脈出力関数の良好な推定を精密に得るために、設定プロセス50(詳細には、理論上の灌流モデル)の修正および洗練62を行なうことができる。
【0365】
次に、図8に関連して、本発明による方法を、理論上の動脈入力関数が遅延τまで絶対的な確かさおよび無限の精度で与えられると仮定される第2の応用の例について説明する。
【0366】
したがって、実験による灌流モデルMを
【0367】
【数113】
【0368】
と表わすことができ、ここですべての量は、今や既知であると仮定される理論上の動脈入力関数の値のベクトルであるa=[C(t),…,C(t)]を除き、すでに定めたとおりである。
【0369】
したがって、結合直接確率分布が、p(s|a,h,Θ,Θ,E,I)と表わされ、結合事前確率分布が、p(h,E,Θ,Θ,E|I)と表わされ、I=(IΘ,I,I,M)である。これらの分布が、上述のように割り当てられる。ベイズ規則は、下記のようになる。
【0370】
【数114】
【0371】
続いて、すでに説明した同じやり方で、任意のパラメータ
【0372】
【数115】
【0373】
の推定、信頼区間、およびそれらの信頼区間の賭け率が得られる。
【0374】
次に、図8に関連して、本発明による方法を、理論上の動脈入力関数が与えられず、動脈入力信号も測定されない第3の応用の例について説明する。
【0375】
したがって、実験による灌流モデルMを
【0376】
【数116】
【0377】
と表わすことができ、ここですべての量は、すでに定めたとおりである。
【0378】
結合直接確率分布が、やはりp(s|a,h,Θ,Θ,E,I)と表わされ、結合事前確率分布は、今やp(a,E,h,E,Θ,Θ,E|I)と表わされ、I=(IΘ,I,I,I,M)である。ベイズ規則は、下記のようになる。
【0379】
【数117】
【0380】
続いて、すでに説明した同じやり方で、任意のパラメータ
【0381】
【数118】
【0382】
の推定、信頼区間、およびそれらの信頼区間の賭け率が得られる。
【0383】
該当の応用の例がどれであっても、公知の方法と対照的に、本発明による推定方法が、関心対象量について先験的に有する定性的、定量的、または半定量的な情報を、実験による測定によってもたらされるそれらの量についての事後情報を一意に決定するために、ベイズ確率理論へと変換することだけからなる点で、正確な方法であることに、注目すべきである。本発明は、問題を解決するために論理的に必要かつ充分な種々の不確実な定性的制約をエンコードしている最も不確かな確率分布だけを導入するため、検証されることができ、あるいは(実際には可能であっても)検証されることができない(詳細には、インパルス応答または相補累積分布関数についての)恣意的な仮説が、不要である。
【0384】
したがって、動脈入力関数が与えられず、せいぜい測定されるにすぎない場合でも、それらの方法は、提案された動脈入力信号が関心対象の各ボクセルの<<真の>>局所的動脈入力関数に効果的に対応することができるか否かをテストすることを可能にし、実際に、動脈入力信号が適しておらず、真の局所的動脈入力関数に対応しない場合に、問題の解であると同時に種々の事前制約を満足する(または、少なくともその確率が先験的に無視できない)パラメータ(血行動態パラメータ、相補累積分布関数、など)の組が、存在しない可能性がある。この場合、確率論は、提案された動脈入力信号を雑音と解釈し、標準偏差σ/εが、より適した動脈入力信号から典型的に得られる標準偏差と比べてはるかに大きくなる。
【0385】
したがって、本発明は、技術水準による大域的または局所的動脈入力関数の種々の選択および推定方法をテストするというきわめて興味深い応用を有する。それらの大域的または局所的動脈入力関数から得られた推定(詳細には、それらの正則化パラメータEの推定)が、或るボクセルから他のボクセルへと過剰に頻繁に異常であることが明らかになる場合、それらの関数の新たなより適切な(局所的)選択法を導入するか、あるいは動脈入力関数が与えられること、または動脈入力信号が事前に測定されることを必要としない方法に頼るか、のいずれかが必要であると結論付けられ、これがまさに上述した第3の方法の目的である。
【0386】
応用の例として、図1または2において説明したような灌流強調解析イメージングシステムによる本発明の主たる実行のステップに言及することができる。
【0387】
・処理ユニット4(または、前処理ユニット7)によって、患者の記録が開かれ、あるいは一連の画像が取り込まれ、関心対象の一連の画像が選択され、詳細には図5aによって示されるように各ボクセルについての灌流信号S(t)を得るための時間につれての灌流強調画像I1〜Inが選択される。
【0388】
・関心対象のスライスまたは領域のユーザ6による特定を可能にするためのヒューマン・マシン・インターフェイス5による画像のプレビュー。
【0389】
・本発明による推定方法の実行を可能にするための設定パラメータ(入力情報)からの処理ユニット4の設定。
【0390】
・推定すべき1つ以上の関心対象量の選択。
【0391】
・人間の脳などの器官についてのBF^、τ^、またはMTT^などの血行動態パラメータ14の処理ユニット4による推定。
【0392】
・動脈入力関数が与えられないときのE、E、ΘまたはE、ESaまたはΘSaなどの他のパラメータの随意による推定。
【0393】
・随意による、インパルス応答の値のベクトルの推定h^、相補累積分布関数の推定b^、動脈入力関数の値のベクトルの推定a^、静脈出力関数の値のベクトルの推定v^、理論上の濃度の値のベクトルの推定c^、理論上の信号の推定sth^、または剰余の推定r^。
【0394】
・の前記推定による関心対象量14が、内容を医師へと示すために、例えば各ピクセルの輝度または色が計算された値に応じて決まるマップとして最終的に表示されるよう、ヒューマン・マシン・インターフェイス5へと引き渡される。
【0395】
・ユーザによって選択された関心対象のいくつかのボクセルについて、インパルス応答、相補累積分布関数、動脈入力関数、静脈出力関数、理論上の濃度、理論上の信号、または剰余が、随意により表示される。
【0396】
・血行動態パラメータなどの関心対象の1つ以上のパラメータについての信頼マップまたは賭け率マップが随意により表示される。
【0397】
・ユーザによって選択された関心対象のいくつかのボクセルに関して、いくつかのインパルス応答、相補累積分布関数、動脈入力関数、などについての信頼区間またはそれら区間への賭け率が随意により表示される。
【0398】
・実験データと大域的非パラメータ灌流モデルとの間の1つ以上の距離の誤差マップが随意により表示され、詳細には大域的灌流モデルに鑑みた実験データの確率が表示される。
【0399】
・1つ以上の血行動態パラメータ、インパルス応答(または、相補累積分布関数)、または局所的動脈入力関数の異常な分布を特徴とする関心対象の前記病気の領域の補助された選択。
【0400】
・おそらくは病変の領域に関係しており、医師による治療行為(例えば、血栓を分解するための経静脈血栓溶解)の決定の対象となりうる異常な灌流の組織領域部分の処理ユニットによる推定。
【0401】
・病変の容積と異常な灌流の領域との間の比など、医師による診断および治療の決定(例えば、血栓を分解するための経静脈血栓溶解)の精緻化を可能にするための何らかの量の処理ユニットによる推定。
【0402】
すでに述べたように、処理ユニット4の設定プロセス50を、ユニットそのものによって実行することができる(設定プロセス50の実行)。あるいは、前記設定は、推定を望む関心対象量に応じた結合事後確率分布のライブラリの保存および選択するこを含むことができる。このライブラリの構築を、処理ユニット4と協働することができる専用のユニットによって達成することができる。
【0403】
あるいは、前記設定の精緻化のために、いくつかの関心対象量についての信頼区間およびそれら信頼区間の賭け率の推定に続いて、繰り返しを行なうことができる。実験データと大域的非パラメータ灌流モデルとの間の距離の提供、詳細には大域的灌流モデルに鑑みた実験データの確率の表示は、設定の更新をさらに含むことができる。
【0404】
したがって、本発明は、パラメータの推定を、各ボクセルの輝度または色が計算された値に例えば線形なやり方で依存する<<パラメータマップ>>の形態で表示することを目的とする。同様に、本発明は、おそらくは、それらの推定の標準偏差を<<信頼度マップ>>の形態で表示するとともに、対応する信頼区間についての賭け率を<<賭け率マップ>>の形態で表示することを目的とする。インパルス応答、相補累積分布関数、動脈入力関数、静脈出力関数、濃度−時間曲線、灌流信号、または剰余の値のベクトルの推定に関して、本発明は、それらをユーザによって選択された各ボクセルについて時系列の形態で表示することを目的とする。最後に、本発明は、実験による信号と非パラメータ灌流モデルとの間の距離またはこのモデルに鑑みた実験データの確率を<<誤差マップ>>の形態で表示することを目的とする。
【0405】
図9〜12が、本発明に従って推定された血行動態パラメータ14あるいはそれらに関する標準偏差または確率などのいくつかの関心対象量について、マップの形態での表示の態様を説明することを可能にしている。
【0406】
したがって、核磁気共鳴イメージングによって分析された人間の脳について、図9が、脳血液容量CBVの推定を眺めることを可能にする。そのようなマップ(458×458ピクセル)が、おそらくは虚血性の領域80の強調を可能にする。実際、適切なインターフェイス5によって、右後大脳動脈の領域において反対側の半球と比べてパラメータCBVが明らかに増加していることを見て取ることができる。図9に示されるマップを読み取ることによって、虚血の後の血管拡張が明らかになる。
【0407】
図10が、虚血性脳梗塞の場合における脳血流の推定に関するマップ(458×458ピクセル)を示すことを可能にする。マップを分析することによって、虚血の後の右後大脳動脈の領域における反対側の半球と比べたパラメータCBF(脳血流量)の減少を見て取ることができる。そのようなマップが、おそらくは虚血性の領域80の強調を可能にする。
【0408】
図11が、平均通過時間MTTの推定に関するマップ(458×458ピクセル)を示すことを可能にする。マップを分析することによって、虚血の後の右後大脳動脈の領域81における反対側の半球と比べたMTTの明らかな増加を見て取ることができる。
【0409】
図12が、脳血流量が信頼区間[CBF^−σ^CBF,CBF^+σ^CBF]に属する確率の推定に関するマップ(458×458ピクセル)を示すことを可能にする。マップを分析することによって、異常なモデルの適合を除いて、確率が0.68のあたりに中心を有することを見て取ることができる。これは、CBFの事後確率分布がガウスの法則に従う場合に期待することができる値である。
【0410】
上述のマップのおかげで、本発明は、技術水準の公知の技術を用いた場合には入手することができなかったさまざまな有用な情報をユーザに提供することを可能にする。この提供は、図1または2に示したとおりの処理ユニット4を、推定したパラメータ14を例えば図9〜13によって示されるとおりのマップの形態でユーザ6へと示すことができるヒューマン・マシン・インターフェイス5に適したフォーマットで、外界との通信のための手段がそのパラメータ14をもたらすことができるように、構成することによって可能にされる。
【0411】
本発明によれば、もたらされる情報がより多数かつより公正である。医師にとって利用可能な情報が、診断における判断および治療上の意志決定における医師の自信を高めると考えられる。
【0412】
本発明によるシステムの性能を向上させるために、処理ユニットが、血行動態パラメータ、相補累積分布関数、または動脈入力関数の推定が必要とされる画像のボクセルについての演算を並列化するための手段を備えることができる。これは、グラフィカル・プロセッサ・ユニット(GPU)などのハードウェア技術、コンピュータクラスタ、または並列モンテカルロ法などのソフトウェア技術、などを使用することによって達成することができる。あるいは、本発明による処理ユニットが、離れて位置する演算手段に頼ることができる。このやり方で、演算の時間をさらに大幅に減らすことができる。
図1
図2
図3
図4
図5a-5b】
図6
図7
図8
図9-10】
図11
図12