CN101002104B - 动脉输入和组织残留函数在灌注磁共振成像中的盲确定 - Google Patents

动脉输入和组织残留函数在灌注磁共振成像中的盲确定 Download PDF

Info

Publication number
CN101002104B
CN101002104B CN2005800225254A CN200580022525A CN101002104B CN 101002104 B CN101002104 B CN 101002104B CN 2005800225254 A CN2005800225254 A CN 2005800225254A CN 200580022525 A CN200580022525 A CN 200580022525A CN 101002104 B CN101002104 B CN 101002104B
Authority
CN
China
Prior art keywords
voxel
input function
interest
area
artery input
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
CN2005800225254A
Other languages
English (en)
Other versions
CN101002104A (zh
Inventor
托芬·塔克斯特
雷内特·格鲁纳
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
UNIFOB Universitetsforskning Bergen (Bergen Univ Res Foundation)
Original Assignee
UNIFOB Universitetsforskning Bergen (Bergen Univ Res Foundation)
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by UNIFOB Universitetsforskning Bergen (Bergen Univ Res Foundation) filed Critical UNIFOB Universitetsforskning Bergen (Bergen Univ Res Foundation)
Publication of CN101002104A publication Critical patent/CN101002104A/zh
Application granted granted Critical
Publication of CN101002104B publication Critical patent/CN101002104B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01RMEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00Arrangements or instruments for measuring magnetic variables
    • G01R33/20Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48NMR imaging systems
    • G01R33/54Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
    • G01R33/56Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution
    • G01R33/5601Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution involving use of a contrast agent for contrast manipulation, e.g. a paramagnetic, super-paramagnetic, ferromagnetic or hyperpolarised contrast agent

Landscapes

  • Physics & Mathematics (AREA)
  • Health & Medical Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Radiology & Medical Imaging (AREA)
  • Engineering & Computer Science (AREA)
  • Signal Processing (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Condensed Matter Physics & Semiconductors (AREA)
  • General Physics & Mathematics (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

本发明提供一种对有脉管的动物对象的感兴趣区域的灌注磁共振成像方法,所述方法包括:向所述对象的脉管系统中给予造影剂药剂;在从造影剂到达所述感兴趣区域之前到至少所述造影剂穿过所述感兴趣区域的第一通过的末端期间内在一系列时间值(t)上确定所述感兴趣区域的体素(i)的磁共振信号强度si(t);由确定的所述信号强度的值si(t)和一动脉输入函数v(t),对每个所述体素确定组织残留函数ri(t)的值;可选择的,由确定的ri(t)的值产生所述感兴趣区域的图象;改进包括由si(t)产生体素特定动脉函数vi(t)并使用所述体素特定动脉函数来确定组织残留函数ri(t)的值。

Description

动脉输入和组织残留函数在灌注磁共振成像中的盲确定
技术领域
本发明涉及一种磁共振(MR)成像,尤其是灌注成像的方法,和成像后数据处理的方法及其程序。
背景技术
MR成像可用于检查体内的身体结构和组织功能及存活性。一种在用于例如由于中风、肿瘤、梗塞等而导致的有脑、心或肝损伤的患者中特别有价值的成像技术,叫做灌注成像。
在这种技术中,一些MR造影剂药剂(例如,象Amersham及Schering所销售的Omniscan
Figure 058225254_0
或Magnevist一样的螯化钆)被给予到患者的脉管系统中,在包括造影剂药剂穿过感兴趣区域的组织的期间内,收集来自感兴趣区域的MR图象。为此,使用快速图象获取顺序,例如自旋回波、梯度回讯(GRASS或FLASH)、平面回波、RARE、混杂、半激励等等。在本领域中这样的顺序和给予MR造影剂药剂用于关注成像是众所公知的(例如见1988年Wehrli等人编辑的“Biomedical Magnetic Resonance Imaging”)。
在临床实践中,通常检查关注图象系列并定性地评估结果。
然而,在许多情形中,希望有对例如区域血流量、区域血容量、区域平均通行时间、到达的区域时间和到顶点的区域时间的定量结果。(例如见Rempp等人的Radiology,193:637-641,1994,以及Vonken等人的MRM,41:343-350,1999)。
使用传统的示踪剂动力学理论,分析动态的MR灌注图象以对每个体素提取组织残留函数r(t)。由r(t),可以确定上面提到的区域参数的定量的数值。然而,局部造影剂信号,c(t),其可以由在时间t的MR信号强度对时间0(即,在造影剂到达关注组织之前)的MR信号强度的比的开方而得出,是r(t)和动脉输入函数v(t)的卷积。
r(t)描述在时间t时仍存在与组织区域中的造影剂,因而是一个依赖组织的生理参数,例如血容量和平均通行时间的函数。v(t)描述了造影剂怎样被输送到组织体素中,这样给出了在器官中的脉管结构的印象。常规上,v(t)假定是空间上不变的,即,在由一体素ci(t)确定ri(t)中,v(t)确定为独立值i而通过反卷积ri(t)*v(t)来确定ri(t)。这样
c i ( t ) = k . ln s i ( t ) s i ( o ) = r i ( t ) * v ( t )
(其中,k是依赖于成像技术中所用的参数的常数,而si(t)是体素i在时间t的信号强度。)
发明内容
现在我们已经惊奇地发现通过由MR图象强度信号si(t)确定vi(t)的体素特定值,临床信息意义更大的ri(t)值和因此上面提到的区域参数的值可被确定。另外,通过将vi(t)表现为时间的函数,可以评估器官的血液供应模式。
因而从本发明提供的一个方面看,在对有脉管的动物对象中的感兴趣区域的灌注磁共振成像方法中,所述方法包括:将造影剂药剂给予到所述对象的脉管系统中;在从造影剂到达所述感兴趣区域之前到至少所述造影剂穿过所述感兴趣区域的第一通过的末端期间内在一系列时间值(t)上确定对所述感兴趣区域的体素(i)的磁共振信号强度si(t);由确定的所述信号强度的值si(t)和动脉输入函数v(t),对每个所述体素确定组织残留函数ri(t)的值;可选择由确定的ri(t)的值产生所述感兴趣区域的图象;改进包括由si(t)产生体素特定动脉函数vi(t)并使用所述体素特定动脉函数来确定组织残留函数ri(t)的值。
从另一方面来看,本发明提供一种灌注磁共振成像数据处理方法,所述方法包括:获得对其中脉管中已被给予磁共振造影剂药剂的有脉管对象中的感兴趣区域的体素(i)的一组磁共振信号强度值si(t),所述信号强度值是对在从造影剂到达所述感兴趣区域之前到至少所述造影剂穿过所述感兴趣区域的第一通过的末端期间的一系列时间值(t);由所述该组信号强度值si(t)和一动脉输入函数vi(t),对每个所述体素确定组织残留因子ri(t)的值;可选择由确定的ri(t)的值产生所述感兴趣区域的图象;改进包括由si(t)产生体素特定动脉函数vi(t)并使用所述体素特定动脉函数来确定组织残留函数ri(t)的值。
从又一方面看,本发明提供一种灌注磁共振成像方法或灌注磁共振成像数据处理方法,其中对于造影剂泄漏效应及或造影剂再循环效应,校正了信号强度值。
从再一方面看,本发明提供一种灌注磁共振成像方法或灌注磁共振成像数据处理方法,其中对每个体素的造影信号被分成其没有泄漏或具有泄漏成分的第一通过和第二通过部分。
从另一方面看,本发明提供一种灌注磁共振成像方法或灌注磁共振成像数据处理方法,其中使用理查德森-露西反卷积算法来估计组织残留函数并而后估计动脉输入函数,然后使用维纳反卷积滤子从对动脉输入函数的估计中找到对组织残留函数的新的估计。
从另一方面看,本发明提供一种灌注磁共振成像方法或灌注磁共振成像数据处理方法,其中在迭代循环中依次使用理查德森-露西反卷积算法和维纳反卷积滤子。
从另一方面看,本发明提供一种灌注磁共振成像方法或灌注磁共振成像数据处理方法,其中使用理查德森-露西反卷积算法来估计组织残留函数而后估计动脉输入函数。
从另一方面看,本发明提供一种灌注磁共振成像方法或灌注磁共振成像数据处理方法,其中在迭代循环中使用理查德森-露西反卷积算法。
从另一方面看,本发明提供一种灌注磁共振成像方法或灌注磁共振成像数据处理方法,其中产生一参数灌注图象。
从另一方面看,本发明提供一种计算机软件产品,其携带的指令在用于本发明的灌注磁共振成像方法的数据处理装置上被执行时,将使该装置构建为对体素(i)接受一组磁共振信号强度值si(t)作为输入,所述体素是在其中脉管中已被给予磁共振造影剂药剂的有脉管对象的感兴趣区域中,所述信号强度值是对在从造影剂到达所述感兴趣区域之前到至少所述造影剂穿过所述感兴趣区域的第一通过的末端期间内的一系列时间值(t),及由此产生体素特定动脉函数vi(t)、体素特定组织残留因子ri(t),及可选的,感兴趣区域的图象。到计算机程序的输入可以是原始信号强度值本身或者优选为归一化的信号强度值(si(t)/si(0))或者归一化信号强度值的对数。
从再一方面看,本发明还提供一种用于本发明的灌注磁共振成像方法中的数据处理装置,该装置构建为对体素(i)接受一组磁共振信号强度值si(t)作为输入,所述体素是在其中脉管中已被给予磁共振造影剂药剂的有脉管对象的感兴趣区域中,所述信号强度值是对在从造影剂到达所述感兴趣区域之前到至少所述造影剂穿过所述感兴趣区域的第一通过的末端期间内的一系列时间值(t),及由此产生体素特定动脉函数vi(t)、体素特定组织残留因子ri(t),及可选的,感兴趣区域的图象。
该数据处理装置优选为这样的装置,其也设置为控制磁共振成像仪。
本发明的方法比起用于获得参数灌注图象(例如rCBV、rCBF、rMTT和到达的区域时间)的现有技术提供了显著的改进,因为可自动产生体素特定的动脉输入函数。相比之下,在现有技术方法中,需要一个手动的,交互式的步骤且所用的动脉输入函数对所有体素是一样的。对于本发明的方法的执行,仅需要依照成像装置和成像次序的特性,以及可选的,造影剂的特性,来设定在算法中使用的参数的值。因而可以预加载适当的参数设定到用来完成参数灌注图象产生的计算机程序中,从而操作者仅需要向计算机指明所用的造影剂的性质。实际上由机器制造商提供的程序可能包含一查询表,从而之后可以自动选择计算所需要的参数的适当的值。虽然在下面会被更详细地解释,在本发明的处理和方法中的数据处理通常涉及对每个体素的归一化的信号强度函数(通常与ln si(t)/si(0)成正比)的傅立叶变换以得出函数ci(ω)的值,将ci(ω)对数复变换以得出函数
Figure S05822525420070108D000041
逆傅立叶变换
Figure S05822525420070108D000042
以得出
Figure S05822525420070108D000043
分成短时间函数
Figure S05822525420070108D000045
和长时间函数
Figure S05822525420070108D000046
以及反转该过程而将
Figure S05822525420070108D000047
变成vi(t)。随后通过使用维纳滤子或其他反卷积方法的反卷积得到r(t)。将
Figure S05822525420070108D000048
分成
Figure S05822525420070108D000049
Figure S05822525420070108D0000410
可能通常涉及部分函数跟随最初过滤或在函数
Figure S05822525420070108D0000411
的加权下在部分重叠的位置上的迭代拟合。
以传统的方式,从确定的rj(t)的值,可以完成产生上面提到的原始参数,及/或感兴趣区域的图象。
一种可选择的策略,其同时也是本发明的一部分,是使用盲迭代反卷积方法来对每个体素同时估计组织残留函数和动脉输入函数体素。一种用于迭代算法的可能性是结合两种反卷积方法;即理查德森-露西反卷积算法以孤立出卷积的函数之一,然后再次使用理查德森-露西函数以得到第二个函数,然后维纳反卷积滤子算法以获得对第一个函数等等的更好的估计。以这种方式,可以消除每种方法的限制。或者,可以在所有的迭代中使用理查德森.露西反卷积。下面要进一步描述这些方法。
在本发明的方法中使用的快速成像技术可便利地成为任何在每个切片之间能用少于2秒钟,优选1.5秒钟的时间延迟获得图象的技术,优选是图象获取时间少于1秒,更优选是少于100毫秒。为此可以使用上述的成像技术。优选使用基于自旋回波(T2成像)或梯度回波(T2 *成像)信号产生的回波平面顺序,特别是梯度回波回波平面成像。
优选使用T2 *成像。
尽管可以使用例如Gd DTPA、Gd DTPA-BMA和Gs HPD03A这样的分布到细胞外流体中的常规MR造影剂在本发明的方法中,在要成像的器官或组织不是脑的情况下,优选使用血液淤积造影剂,即保留在脉管中的那种。这种造影剂在MR造影剂制造商,例如Schering、Amersham、Nycomed等的专利出版物中有很多记载。
脉管内造影剂的组织浓度c(t),依赖组织流速F,组织函数r(t),和动脉输入函数v(t)。该函数关系是在r(t)和v(t)之间的卷积。
c(t)=F.r(t)*v(t)=u(t)*v(t)    (1)
其中u(t)=F.r(t)。残留函数r,描述在即刻注射后的时间t时仍存在于组织区域中的造影剂的部分,因而是组织的生理参数(血容量和平均通行时间)的函数。流速F,可以从反卷积的造影剂浓度函数的最初高度得出,因为在时间t=0时(r(0)=1),在组织中没留有造影剂。在方程式(1)中,忽略了造影剂的脉管外泄漏和再循环。在具有完整的血脑屏障的灌注研究中,当仅仅分析代表造影剂的第一通过的时间点测量时,这是充分的近似。
方程式(1)对每个体素(i)都是成立的,其具有它特有的特定动脉输入函数和组织残留函数
ci(t)=ui(t)*vi(t)    (2)
体素特定动脉输入函数包括在其到被观察组织的通行上的,从注入脉管w(t)的动脉输入函数的任何延迟δ,及离差g。
vi(t)=w(t)*gi(t)*δi(t-t0i)    (3)
其中t0i是造影剂在体素(i)中的到达时间。因而在文献中描述的脉管运送函数由g和δ之间的卷积来表达。
在使用方程式(2)在MR灌注成像之前,ci(t)必须与MR信号强度si(t)的相应改变联系起来。对于回波平面成像实验,在组织的松弛特性和造影剂浓度之间的线性关系是良好的近似。
c i ( t ) = - k TE · ln s i ( t ) s i ( 0 ) - - - ( 4 )
其中TE是回波时间,si(0)是在注射造影剂之前的恒定信号强度。比例常数K,依赖于磁场强度、造影剂类型和成像顺序参数。
复倒频谱是对观测到的信号的特定的频率表达。如果观测到的信号是几个单独的信号的集合,所述单独的信号属于复倒频谱中的单独波段,通过过滤该复倒频谱可对信号进行分解。在本发明中,复倒频谱分析可应用到方程式(2)中的卷积上以便确定体素特定动脉输入函数。
通过对关注信号的同型变换而获得复倒频谱表达。存在几个同型变换,导致不同的倒频谱,象双对数倒频谱(bicepstrum)、光谱根倒频谱、通用倒频谱和复倒频谱。用来计算复倒频谱的同型变换,是复对数函数,log,其具有复指数函数,Ae=Acos(Φ)+iAsin(Φ),作为反函数。可以用不同的方法来计算这个同型变换,包括傅立叶变换方法和提取复数根的方法。在傅立叶方法中,傅立叶变换,F被应用到方程式(2)中的卷积中
Ci(z=ejw)=F[ci(t)]=F[vi(t)]·F[ui(t)]    (5)
G ^ i ( z ) = log ( C i ( z ) ) = log ( [ v i ( t ) ] ) + log ( [ u i ( t ) ] ) - - - ( 6 )
Figure S05822525420070108D000062
的反傅立叶变换(F-1)表达,给出想要的复倒频谱表达, c ^ i ( t ) = v ^ i ( t ) + u ^ i ( t ) . 这样,可以在方程式(2)中的两个卷积的函数之间建立线性的关系。然而,在计算上要求的相位展开规划和由于有限傅立叶表达引起的混淆,在使用这种方法而不是下面要描述的复数根提取中会是一个问题,而复数根提取不会受到混淆的影响,优选用在本发明的方法中。
在方程式(2)中的造影剂信号,
Figure S05822525420070108D000065
tε[1,T]是实数,稳定和有限的。在这种情况中,ci(t)的z变换,Ci(z),可以在包括单位圆的区域内被分析,并可以被表达成为N=Mi+Mo阶的复数多项式,
C i ( z ) = | A i , | Π k = 1 M i ( 1 - a i , k z - 1 ) Π k = 1 M o ( 1 - b i , k z - 1 ) - - - ( 7 )
系数|ai,k|≤1及|bi,k|≤1,从而在方程式(7)中的因数分解项分别表达在单位圆内侧和外侧的Mi零点和Mo零点。方程式(7)的复对数是
C ^ i ( z ) = log | A i | + r . log ( z ) + Σ k = 1 M i log ( 1 - a i , k z - 1 ) + Σ k = 1 M i log ( 1 - b i , k z ) - - - ( 8 )
这样,通过使用数字多项式根算法,可以直接计算方程式(8)的复零点ai,k,和1/bi,k。可以根据下式计算复倒频谱
c ^ i ( t ) = log | A i | t = 0 - &Sigma; k = 1 M i a i , k t t t < 0 - &Sigma; k = 1 Mo b i , k t t t > 0
典型的50个时间点测量的灌注成像系列,在体内估计中给出收敛和稳定的结果。通过仅包括代表造影剂药剂的第一通过的图象帧(例如,22时间点测量),可有利地减少计算的时间。
在复倒频谱中的线性关系 c ^ i ( t ) = v ^ i ( t ) + u ^ i ( t ) , 意味着如果两个导出的函数位于复倒频谱中的不同波段内,那么可以通过简单的理想的带通过滤而分开它们。正确的分开对于确定体素特定动脉输入函数是重要的。在经过理想的低通过滤以后,复倒频谱可假定成 c ^ i ( t ) = v ^ i ( t ) . 接着,基于傅立叶方法,可以应用反同型变换到被过滤的复倒频谱数据上以获得vi(t)。第一步是将v变换成对数傅立叶域,F-1(Vi)=Vi(eiw),然后是普通傅立叶域,Vi(z=ei)。按照这一方法所有估计的动脉输入函数都可以被归一化到空间域中的单位面积上,因为从倒频谱分开中不能唯一地确定任何对动脉输入函数的乘法系数。这一归一化的可能的物理解释是,同样数量的造影剂,在这种情况中是单位面积的造影剂药剂,通过每个图象体素。这样,假定正确再现了体素特定输入函数的形态,仅仅一个比例因子将它从相应的手工选择的动脉输入函数上分开。
根据下式,简单频率域维纳过滤是优选用于反卷积体素特定动脉输入函数和组织残留函数的方法。
U i ( e iw ) = C i ( e iw ) V i * ( e iw ) | V i ( e iw ) | 2 + q i - - - ( 9 )
参数qi是算出的造影剂浓度的信噪比的倒数的开方。这个维纳滤子给出初始信号的后验概率估值的最大值,其谬误噪音被当作独立的高斯过程模式,而信噪比对所有频率都被取作为常数。可以用其他的标准的反卷积方法,例如基于正则化、奇异值分解方法或最大相似性的方法,来代替维纳过滤方法。
当在特定情况中,可以执行本发明的方法而不考虑造影剂泄漏出脉管或者造影剂再循环的可能性,通常希望考虑这些因素,特别是在造影剂是分布到细胞外流体中(而不是停留在血液中的血液淤积药剂)的简单的螯化物时。
除了正常脑组织,几乎所有组织都具有造影剂穿过毛细管膜的实质泄漏。许多脑肿瘤不会具有完整的血脑屏障,而泄漏造影剂到脉管外空间中。泄漏经常对造影剂信号产生明显的贡献,在确定体素特定动脉输入函数中不应被忽视。为了估计泄漏项,必须考虑第一通过造影剂信号的较晚的部分。这一较晚的信号部分通常受到再循环的造影剂的信号的污染。因此,在确定体素特定动脉输入函数中,值得考虑再循环效应。
因而在优选实施例中,每个体素的造影剂信号被分成没有泄漏和有泄漏成分的第一通过和第二通过部分。然后,通过在复倒频谱域中过滤而获得每个体素特定动脉输入函数的初始估计。然而,已经看到,如果所有的初始值大于零,最终的结果不会明显受到初始化的影响。此后,可以使用一维迭代露西-维纳反卷积方法来改善对每个体素特定动脉输入函数的初始计算值以及其相应的组织残留函数和泄漏项。
在这个实施例中,对泄漏和再循环使用标准模型(见Vonken等人的MRM 43:820-827(2000)),先去掉造影剂信号的由于再循环和泄漏引起的附加项。使用造影剂信号在空间或傅立叶域中的简单特性和皮卡尔迭代,可以进行对这些附加项的估计。在去掉附加项以后,分离出第一通过信号。使用复倒频谱分析获得对在每个体素中的动脉输入函数的初始估计。
得到的对动脉输入函数的初始估计可被用作向非参数确定性图象约束反卷积算法的输入(见Kundur等人的IEEE  Signal Processing Magazine 13:43-64(1996))。在这些算法中,使用了动脉输入函数和组织残留函数这两者的非负和已知的有限支持的确定性约束。对体素特定组织残留函数和动脉输入函数的估计被迭代更新到局部最优。对每个迭代,在前的图象约束被加强。
有很多中非参数确定性图象约束反卷积算法(看Kundur等人的IEEESignal Processing Magazine 13:43-64(1996))。最早提出的和最被人熟知的一些是迭代盲反卷积算法、非负和支集约束递归逆过滤算法和模拟退火算法。在Kundur(上面的)中给出了对这些算法和其他一些的很好的回顾。
最近的在非参数确定性图象约束反卷积算法上的发展是迭代链算法(看Tsumuraya等人的J.Opt.Soc.Am.A 13:1532-1536 (1996))。这些算法顺次结合了两个不同的反卷积算法而利用了两者的优点。在马上要描述的优选实施例中,使用一维迭代露西-维纳反卷积算法来提高对每个体素特定动脉输入函数和它相应的组织残留函数的初始估计。选择迭代露西-维纳反卷积算法是有利的,因为没有参数一定要被手工设定并且在每个迭代步骤中的两个反卷积算法倾向于抵消每个算法在单独使用时在反卷积结果上具有的负效果。这一点已经在模拟数据和体内数据上被证实。与用前述的算法所获得的结果相比,结果受到噪音的影响要小得多。
对单个体素i的,有泄漏但没有再循环的造影剂信号模型,是:
c i ( t ) = F i &CenterDot; r i ( t ) * v i ( t ) + &kappa; i &CenterDot; &Integral; 0 t F i &CenterDot; r i ( &tau; ) * v i ( &tau; ) d&tau; + v i ( t ) - - - ( 10 )
其中ci(t)是造影剂信号,Fi是血流量,ri(t)是组织残留函数,vi(t)是动脉输入函数,κi是泄漏常数而vi(t)是噪音信号,其通常假定是不相关的和高斯分布的。
在简单的标准现有技术模型中,造影剂的再循环部分,m∈<0,1>是假定对所有体素共有的参数。同样由于造影剂一次穿过循环系统的再循环引起附加的动脉输入函数,v0(t),也假定为对所有体素是共有的。由于身体和心脉管系统的结构,这些假设被认为是合理的。再循环时间,Ti,对每个体素是特定的。使用这些假设,具有一次造影剂再循环但没有泄漏的体素特定造影剂信号模型是:
ci(t)=Fi·ri(t)*vi(t)+m·v0(t)*[Fi·ri(t)*vi(t+Ti)]+vi(t)    (11)
结合方程式(10)和(11),给出具有泄漏和再循环时间Ti的体素特定造影剂信号模型:
c ij ( t ) = F i &CenterDot; r i ( t ) * v i ( t ) + m &CenterDot; ( F i &CenterDot; r i ( t ) * vo ( t ) * v i ( t + T i ) ) + &kappa; i &CenterDot; &Integral; 0 t F i &CenterDot; r i * v i ( &tau; ) d&tau; - - - ( 12 )
+ &kappa; i &CenterDot; &Integral; 0 t F i &CenterDot; r i ( &tau; ) * v 0 ( &tau; ) * v i ( &tau; + T i ) d&tau; + v i ( t )
方程式(12)的傅立叶域可以被因式分解成便利的形式
C i ( &omega; ) = F i &CenterDot; R i &prime; ( &omega; ) &CenterDot; V i ( &omega; ) * ( 1 + &kappa; i j&omega; ) &CenterDot; ( 1 + m &CenterDot; V 0 ( &omega; ) &CenterDot; e - j T i &omega; ) + N i ( &omega; ) - - - ( 13 )
其中大写字母表示由相同的小写字母给出的相应的时间信号的傅立叶变换。当可以忽略造影剂泄漏和再循环时,方程式(13)简化成复数乘积
Ci(ω)=Fi·Ri(ω)·Vi(ω)+Ni(ω)    (14)
这是前面所用的用于估计在复倒频谱域中的动脉输入函数vi(t)的标准的起始表达。噪音项,Ni(ω),被忽略。然而,通常不能忽略造影剂泄漏和再循环,因而对vi(t)的估计必须基于方程式(13)。因此,在估计动脉输入函数之前,必须从造影剂信号ci(t)中去掉造影剂泄漏和再循环的贡献。否则,函数1+κi/jω和1+m·e-jTω会使对缓慢改变的函数vi(t)的估计变得错误。
表示没有泄漏的第一通过造影剂信号
di(t)=Fi·r(t)*vi(t)    (15)
并忽略噪音项,方程式(13)可以被简化成
C i ( &omega; ) = D i ( &omega; ) * ( 1 + &kappa; i j&omega; ) &CenterDot; ( 1 + m &CenterDot; V 0 ( &omega; ) &CenterDot; e - j ij &omega; ) - - - ( 16 )
这样,对每个造影剂信号,对di(t)的最优估计对应了对全体参数m、全体再循环函数v0(t)和所有体素特定泄漏常数κi及局部再循环时间Ti,i=1,2,N,其中N是在感兴趣区域内的体素的个数的最优估计。
为了简化估计问题,v0(t)被简单选择为高斯或伽马函数,其中对给定获取顺序和所有患者手工设定参数。这是合理的第一级近似,因为对所有人来说脉管系统特性是相似的。在初步的实验中,对所用的获得参数,具有平均值为零和标准差为2的高斯分布为最好。
可以最优化全部三个剩下的参数,但是这里再循环部分m,基于初步实验被定为值0.05,m在区间[0.025-0.200]之间改变。再一次,这是合理的第一级近似,因为对所有人来说脉管系统的渗透性特性是相似的。用小的递增步骤最优化m的值是直截了当的,但可能明显增加总的计算时间。
使用简单的带5个迭代的皮卡尔松弛算法可以最优化在空间或傅立叶域中的参数κi和Ti(看Kreyszig“Advanced engineering Mathematics”,JohnWiley,NY,1999,1156页)。可以假定造影剂信号曲线的晚些的部分代表稳态,其中在脉管内和脉管外空间之间存在着造影剂的平衡。通过忽略任何再循环效果和使用稳态假设(看Vonken(前面的))找出κi的初始值。最优化的准则仅仅是在具有泄漏和再循环的估计的函数和观测到的浓度信号之间的总的距离。
对每个估计的κi,Ti的最优值是通过在预先指定范围内试验所有可能的Ti的整数值来确定。Ti的上限和下限可以依赖于要检查的身体部分来指定。
一旦完成了最优化,就知道了第一通过造影剂信号,di(t)。在用di(t)替换ci(t)以后,使用象前面那样的复倒频谱分析,方程式(14)然后可用作用于确定体素特定动脉输入函数vi(t)的起始点。
在现在描述的优选实施例中,迭代露西-维纳反卷积算法被引入以获得对体素特定动脉输入函数和体素特定组织残留函数这两者的更好的估计。来自复倒频谱过滤的动脉输入函数仅用作初始估计。
露西-维纳反卷积算法是两个反卷积算法的结合。第一部分是理查德森-露西反卷积算法,其工作在对具有附加的泊松分布噪音的图象的空间域中。它是由下式给出的最大期望(EM)算法
r j + 1 ( x ) = { [ c ( x ) r j ( x ) * v ( x ) ] * v ( - x ) } r j ( x ) - - - ( 17 )
其中j是迭代次数。表示体素的标记i,已经被抑制了。这一算法是天体物理学界非常喜欢和广泛使用以对通过大气记录的有斑点的天文图象反卷积。它满足了图象的正性约束。能量守恒以及自动支撑。因为这个模型具有内在的正性约束,因此可能解的数量减少了并可获得生理学上有意义的结果。
此外,可以选择通过将动脉输入函数归一化到相同(单位)面积,强制造影剂的保持,以保证相同量的造影剂通过所有体素。
另外,动脉输入函数象第一通过造影剂信号一样被限制在同样的时间集,即在预给药之后和在用泄漏模型确认为非第一通过的那些之前的时间点。组织残留函数被限制为第一时间点。
因为它是EM算法,其保证在单独使用时会收敛到最优值。给定输入函数,该最优值是最大后验概率或可能性。假定附加噪音是高斯的,导致不同的EM算法,其不保证图象值的正性。
该最优值总是局部的。它到全体最优值的距离由动脉输入函数的初始估计决定。
基于频率域的反卷积方法没有其中导致局部信息在空间域中损失的空间改变收敛速度和振铃的问题。因此,基于频率域的反卷积方法可以恢复在空间域中的局部信息损失。
当知道动脉输入函数v(t)时,频率域维纳反卷积滤子会最小化对造影剂信号c(t)的均方差。维纳滤子由下式给出
R ( &omega; ) = C ( &omega; ) V ( &omega; ) + | | N ( &omega; ) | | 2 | | R ( &omega; ) | | 2 - - - ( 18 )
在统计公式中,当噪音是高斯分布时,该算法使组织残留函数的后验概率最小化。
从方程式(18)可以看到,虽然想要得到对组织残留函数的估计,为了这一目的需要这一函数本身的功率谱。这使得构建最优滤子较困难且对噪音敏感。然而,通过先使用理查德森-露西反卷积方法来获得对组织残留函数的估计,构建维纳滤子可以更容易且对噪音的敏感性更低。
在迭代循环中依次使用这两个算法。首先,使用复倒频谱分析将组织残留函数和动脉输入函数这两者初始化成大于零的常数值。然后使用理查德森-露西算法估计组织残留函数。其次,使用这一估计来找出对具有相同算法的动脉输入函数的新的估计:
r j + 1 ( x ) = { [ c ( x ) r j ( x ) * v j ( x ) ] * v j ( - x ) } r j ( x ) - - - ( 19 )
v j + 1 ( x ) = { [ c ( x ) v j ( x ) * r j + 1 ( x ) ] * r j + 1 ( - x ) } v j ( x ) - - - ( 20 )
然后,使用维纳反卷积滤子来找出对组织残留函数等等的新的估计。最后,通过使用维纳反卷积滤子找出动脉输入函数的新的估计来完成循环。当计算维纳滤子时使用对动脉输入函数或组织残留函数的先前的估计。
R j + 2 ( &omega; ) = C ( &omega; ) V j + 1 ( &omega; ) + | | N ( &omega; ) | | 2 | | R j + 1 ( &omega; ) | | 2 - - - ( 21 )
V j + 2 ( &omega; ) = C ( &omega; ) R j + 1 ( &omega; ) + | | N ( &omega; ) | | 2 | | V j + 1 ( &omega; ) | | 2 - - - ( 22 )
通常,在约6次迭代后结束迭代。
在本发明方法的优选实施例中,考虑了造影剂从脉管中的泄漏以及造影剂再循环的影响。对于除了脑以外的器官和组织,在使用常规ECF-分布造影剂,例如GdDTPA时,由相当程度的造影剂泄漏。此外,即使在脑灌注成像中,因为脑肿瘤可能有漏的毛细脉管,可能发生穿过血脑屏障的造影剂泄漏。在优选实施例中,对每个体素自动估计泄漏和再循环效应。然后从测得的造影剂信号中去掉它们的信号贡献。剩下来的信号是造影剂的第一通过信号。然后这个信号被用于对体素特定动脉输入函数和体素特定组织残留函数的自动的,迭代的估计。使用本发明这一实施例的方法,维纳滤子常数被自动估计使得如果用于初始化的话仅需要限定在复倒频谱域中的滤子;然而有了对体素特定动脉输入函数和组织残留函数的迭代计算,复倒频谱滤子的精确值和形态不是那么重要了,可选择忽略将信号变换成复倒频谱。
在本发明的另一方面中,理查德森-露西反卷积被用于所有盲估计。在分离出第一通过造影剂信号以后,设定对两个卷积的函数的初始估计。这些估计或者是通过倒频谱分析获得,或者是通过设定函数值为大于零,即非负的常数而获得。对于初始化,优选为卷积的造影剂信号的总能量(曲线下的面积)在对组织残留函数和动脉输入函数初始化中守恒,即流速对组织残留函数的面积的标定等于被卷积的造影剂信号的面积被动脉输入函数的面积除。
在优选实施例中,可选择通过将动脉输入函数归一化成相同(单位)面积来强制造影剂的守恒,以保证相同量的造影剂通过所有体素。
在另一优选实施例中,动脉输入函数象第一通过造影剂信号一样被限制在同样的时间集,即在预给药之后和在用泄漏模型确认为非第一通过的那些之前的时间点。组织残留函数被限制到第一时间点。
此后,通过先改善对动脉输入函数的估计开始迭代。这可以通过将动脉输入函数初始化成归一化的卷积的造影剂浓度来实现。该估计被用于改善对组织残留函数等的估计。理查德森-露西反卷积被用在每次迭代中。每个函数被迭代更新到一局部最优。通常,在约6次迭代后结束迭代。
用于完成本发明的计算机软件产品可以被提供在物理载体形式中例如盘、磁带或其他存储装置,或者可以作为从遥远位置处传来的信号,例如通过象因特网这样的网络被提供。
附图说明
下面参考接下来的示例和附图来描述本发明,在附图中:
图1示出模拟分散效果的两个动脉输入函数v(t)对时间的图示;
图2示出没有噪音的动脉输入函数的图示(估计的(实线)和模拟的(点划线)。a)和b)表(1)的情形I,c)和d)表(1)的情形II。3s的延迟到达时间在c)和d)中。);
图3示出具有加入的获得噪音的动脉输入函数的图示(估计的(平均值±两个标准差,实线)和模拟的(点划线)。a)和b)表(1)的情形I,c)和d)表(1)的情形II。3s的延迟到达时间在c)和d)中。);
图4示出没有噪音的被流速标定的组织残留函数的图示(反卷积的(实线)和指数模拟的(点划线)。a)和b)表(1)的情形I,c)和d)表(1)的情形II。动脉输入函数的延迟到达时间,3s,用在c)和d)中。);
图5示出具有加入的获得噪音的被流速标定的组织残留函数的图示(反卷积的(实线)和指数模拟的(点划线)。a)和b)表(1)的情形I,c)和d)表(1)的情形II。动脉输入函数的延迟到达时间,3s,用在c)和d)中。);
图6示出没有噪音的反卷积的造影剂浓度函数的形状和初始高度(流速)的图示(组织残留函数具有指数形状((a)和(b)),三角形形状((c)和(d)),以及盒形。左栏;反卷积的(实线)和初始(点划线)浓度函数(F=60ml/100g/min,MTT=4s)。右栏;对具有3s(o)、4s()、5s(Δ)MTT的增大的输入流速值算得的流速值(平均值±标准差));
图7示出具有减小的信噪比的估计的流速的图示(50个体素的平均值(实线)和标准差(点划线)a)在重复最优化维纳滤子常数时b)具有不变的维纳滤子常数,q=0.01);
图8以切片示出体内反卷积造影剂浓度作为时间的函数的图象(在用a反卷积每个体素以后,在挑选的8个时间点测量中的信号强度的变化),手工选择的动脉输入函数,b)盲估的,体素特定动脉输入函数);
图9示出的图象是脑血流量(上排),到峰值时间(中排),平均通行时间(下排)(手工选择的动脉输入函数(左)盲估的动脉输入函数(右)。在流速图象中颜色更深的值表示具有更高流速的体素,而在到峰值时间和平均通行时间图象中更深的值表示更短的时间);
图10示出组织标签遮罩(label masks)的图象(a)预造影剂灌注图象,b)叠加在同一图象上的组织标签遮罩。白质;黑色遮罩刚好在侧脑室上方。灰质;六个剩下的黑色遮罩。两个最前面的遮罩;皮质灰色物质。就在侧脑室后面的两个遮罩;头尾状核。两个最后面的遮罩;豆状核。)
图11示出参数灌注图象。(黑度随着数量和流速而线性增大。亮度随着时间线性增大。到峰值时间和平均通行时间图象按比例线性变化以示出相关图象特征。a)具有再循环、泄漏和对动脉输入函数的改进的估计的完整的模型。b)具有再循环、泄漏和对动脉输入函数的复倒频谱估计的完整的模型。c)没有再循环和泄漏以及对动脉输入函数的复倒频谱估计的模型。第一行:血量。第二行:血流量。第三行:平均通行时间。第四行:延迟时间);
图12示出组织残留函数的时程。(黑度随着流速而线性增大。a)具有再循环、泄漏和对动脉输入函数的改进的估计的完整的模型。b)具有再循环、泄漏和对动脉输入函数的复倒频谱估计的完整的模型。c)没有再循环和泄漏以及对动脉输入函数的复倒频谱估计的模型。行一:帧18。行二:帧20。行二:帧23。行二:帧26);及
图13示出泄漏和再循环图象(黑度随着泄漏量线性增大;亮度随着泄漏系数的大小及时间而线性增大。所有图象按比例线性变化以示出相关图像特征。大于0.1的泄漏系数设定为0.1。第一行,泄漏系数图象。第二行,稳态造影剂信号图象。第三行,再循环时间图象。a,c)具有再循环、泄漏和对动脉输入函数的改进的估计的完整的模型。b,d)具有泄漏和对动脉输入函数的改进的估计,但没有再循环项的模型。e)具有再循环、泄漏和对动脉输入函数的改进的估计的完整的模型。f)具有再循环和对动脉输入函数的改进的估计,但没有泄漏项的模型)。
图14示出泄漏和再循环的模型。选定的感兴趣区域在来自肿瘤患者的数据中:a)在肿瘤区域中的结果,b)在皮质灰质中的结果,c)对结合的再循环和泄漏函数在肿瘤区域中的收敛特性,d)在灌注时间数列中叠加在第一图象上的感兴趣区域。
图15示出泄漏和再循环模型。选择的感兴趣区域在来自具有纤维肌性发育异常的患者的数据中的健康灰质。
图16示出到峰值时间,估计的泄漏和再循环参数。a)相对的到峰值时间(rTTP),即卷积的造影剂信号的到峰值时间,b)绝对的到峰值时间(TTP),即在反卷积后的到峰值时间(这里是迭代盲反卷积),c)Ti,即第二通过的到峰值时间,d)泄漏参数,κi。较亮的值大于较暗的值。
图17示出计算机模拟的结果。在a)改变动脉输入函数但组织残留函数为固定的(MTT=4s,指数的)时,b)对固定的动脉输入函数(a=3,b=1.8)和固定的组织残留函数(MTT=4s,指数的)改变噪音状况时,c)改变组织残留函数的形状但动脉输入函数(a=3,b=1.8)为固定时,估计的流值。在d中示出使用不同形状的组织残留函数(具有输入平均通行时间4s)再现的平均通行时间。
图18示出体内动脉输入函数。在图象帧之间的时差是TR=1.442秒。更高的值是暗的。A)左前叶肿瘤。B)右中间脑纤维肌性发育异常。
图19示出体内组织残留函数。在图象帧之间的时差是TR=1.442秒。第一时间集是组织残留函数的最大值,r(t)=1。更小的值在颜色上更亮。A)左前叶肿瘤。B)右中间脑纤维肌性发育异常。
图20示出对具有右中间脑动脉纤维肌性发育异常的患者的体内分析。a)选择动脉输入函数作为具有最早和最大信号变化的体素(上面的更暗的曲线),随后伽马拟合(颜色较浅的曲线)。选择了七个体素(最高亮的体素在较小的图象中),b)建议的平均通行时间方法,c)常规的到峰值时间方法,d)建议的到峰值时间方法。更暗的值表示更短的时间。
图21为感兴趣区域。具有右中间脑动脉纤维肌性发育异常的患者。在表12中示出MTT和TTP的平均值。
具体实施方式
示例1
模拟
按照Calamante等人的MRM44:466-473(2000)和
Figure 058225254_3
stergaard等人的MRM36:715-725(1996)的描述产生综合的数据。参数被调整到在典型的体内灌注研究中获得的模拟数据。在每个情况中,产生包含五十个图象的时间数列。由于每个体素的时间进展在所描述的算法中被分别处理,综合数据经常被分成不同空间段,其每个表示或者不同的物理(即噪音)或者不同的生理(即组织特性)条件。在所有估计中,流速被设定为反卷积的造影剂浓度曲线的最大值,而在出现这一最大值之前的占用时间给出到峰值时间。
必须作出对动脉输入函数和组织残留函数的模型。伽马变量函数被用于模拟动脉输入函数。
v i ( t ) = 0 if t < t oi ( t - t oi ) a e - ( t - t oi ) / b if t > t oi
其中toi是药剂出现在关注的体素(i)中的时间,而a和b是模型相关参数。组织残留函数,ri(t),通常以指数衰减函数为模型。当与组织流速Fi相乘时,它可表达为
u i ( t ) = F i &CenterDot; r i ( t ) = F i &CenterDot; e - t / MTT i - - - ( 23 )
其中MTTi表示通过关注体素的相应的平均通行时间。
通过在傅立叶域中的简单的相乘,得到在ui(t)和vi(t)之间的卷积,
ci(t)=F-1[F[vi(t)]F[ui(t)]]    (24)
其中F和F-1表示傅立叶变换和其逆变换。因为这一方法可能受到有限、离散傅立叶变换的振铃磁感伪影的影响,与Calamante等人(前面的)的解析表达进行了比较。
在方程式(24)中的造影剂浓度曲线,c(t),随后被变换成信号-时间曲线。基线信号强度被设定为si(0)=600MR单位。为了获得逼真的噪音分布,从临床灌注测量中抽取一部分设置在目标边界外的数据。在去掉平均噪音值以后,这些具有噪音的五十时间帧测量(64×64象素)被加入到综合数据组中。
在动脉输入函数中的分散和延迟的效果被分别研究以简化对模拟结果的解释。通过将延迟Δt0,Δt0∈[0,5]秒,加到动脉输入函数来模拟药剂注射的延迟的到达时间。通过对动脉输入函数选择不同的参数值来模拟分散。在动脉输入函数下方的面积在每个情况中被归一化到单位面积。这样,动脉输入函数越宽,最大振幅越小。或者,也可以通过将动脉输入函数以传送函数卷积来实现分散(见Calamante等人的(前面的))。在表达延迟和分散结果中使用了两个不同的动脉输入函数(看附图的图1)。对表1的两个动脉输入函数使用了相同的指数的组织残留函数。在当前的算法中,应重建模拟的动脉输入函数的不同形状。另一方面,任何被延迟的到达时间,被加到组织残留函数上。
动脉输入函数vi(t)和组织残留函数(被流速定比)ui(t),在复倒频谱中是可分开的。为了展示这一点,将不同的动脉输入函数与指数衰减的盒形或三角形组织残留函数分别卷积。
 动脉输入函数  组织残留函数
 I     a=3,  b=1.8II    a=3,  b=1.2  F=60ml/100g/min,MTT=4sF=60ml/100g/min,MTT=4s
表1:在模拟中使用的参数值
模拟是在3、4和5秒这三个不同的平均通行时间下完成的。
反卷积算法对噪音敏感。为了检验维纳滤子反卷积的噪音敏感性,使用蒙特卡洛模拟。应用的在表(1)(情形1)中给出的动脉输入函数和组织残留函数加有综合产生的噪音。产生平均值为零的独立高斯噪音,其被应用而产生信噪比,SNRs,在从10到10000的范围内。通过在五十个单独的时间数列上求平均而计算出结果。
示例2
MRI实验
包括在这一研究中的体内测量是从正在进行的对职业深海潜水员(N=120)的主要研究中选出的。深海潜水员脑受损的机制仍不清楚,然而以前的研究推断出大多数结果可以用脉管阻塞来解释。在这一研究中,选择的对象被认为是健康的,因为并没有确认有明显的病理。灌注获得是在SiemensVision 1.5 T全身医用扫描仪上完成的。选择了标准的,商业的梯度-回波平面成像(GE-EPI)顺序,参数为TR/TE=1442ms/60.70ms,FoV=230×230,矩阵=128×128。造影剂,二亚乙基钆三胺戊乙酸,(Gd-DTPA,即,可从Schering AG买到的Magnevist
Figure 058225254_4
),在12秒后自动注射。总共测量50个时间点。图象数据被传送到离线计算机上用于后处理。除了所建议的方法,使用手工选择的动脉输入函数来完成反卷积。该过程是在C中进行,包括几个数字处方程序。基于XITE包和MATLAB的图象处理工具被用于图象操作,显象和估计。
为了让该方法达到最好的效果,带通滤子和维纳反卷积滤子都要被最优化。在这一研究中,通过检查计算的复倒频谱和已反卷积的造影剂浓度函数,ui(t),来最优化这些参数。
在所有体内研究中,对所有体素选择相同的带通滤子值。假定,变化最慢的函数,动脉输入函数,位于复倒频谱的早些部分,而变化更快一些的组织残留函数位于晚些时候。典型地,在复倒频谱中选择四个或五个最低的时间集以分离动脉输入函数。通常,在分离动脉输入函数时如果一个额外的时间集被包括或排除的话,过滤的结果不太受到影响。
当从估计的体素特定动脉输入函数反卷积造影剂浓度函数时,使用维纳滤子常数,q。这个常数与观测到的造影剂浓度函数的信噪比直接相关。通过对q选择一过小的值,噪音被增大。另一方面,q值过大,不会完成反卷积。在模拟中,在不加入噪音时选择0.009的小q值。这个q值随着噪音值的增大而增大:对SNR∈[10000,10],q∈[0.009,0.15]。在体内数据中对每个体素选择q=0.05这样的值。
示例3
结果
模拟的数据被用于研究用于计算体素特定动脉输入函数和组织残留函数的新的盲的方法的表现。体内记录用于将该方法的表现与其中手工选择单一的动脉输入函数的传统方法相比较。
当没有加入噪音时,图2,和获得噪音被加入时,图3,在两种模拟中精确地估计模拟的动脉输入函数的时程。如上所述,在到达时间上的诱发延迟(在两幅图的右边的图象)没有被包括在估计的动脉输入函数中。
模拟的残留函数的形状在反卷积的造影剂浓度函数中仅仅是部分重现,ui(t)=Fi.ri(t)。指数函数(点划线)的初始振幅被更高斯形状的函数所代替(图4和图5)。类似地,反卷积过程无法重现线性组织残留函数的初始尖锐边缘,或者盒形组织残留函数的恒定平直段(图6)。在组织残留函数中(图(4c,d)-(5c,d)),如预期的那样,可以看到在一些模拟中被加入到动脉输入函数上(图(1c,d)-(2c,d))的延迟的到达时间。
当输入流速值F=30ml/100g/分钟和F=60ml/100g/分钟,使用指数形状的组织残留函数,重复模拟表现出在模拟的动脉输入函数的延迟到达的增加和在算得的到峰值时间值TTP的增加之间的几乎精确的匹配,(表2)。
对造影剂到达中的每个模拟的延迟,计算类似的流速值。虽然独立于延迟,与输入值相比,流速值总是被低估(表2)。当对组织残留函数使用其他形状时,可以看到同样的趋势(图6)。对组织残留函数的每个形状,在输入和算得的组织流速值之间发现线性的关系。这样的事实,即没有发现对组织残留函数的所有选定的模型的独特的关系,显示出反卷积过程依赖于组织残留函数的形状。估计的动脉输入函数不受不同形状的输入组织残留函数的影响,意味着在复倒频谱中的相当准确的分离。
在所有模拟中应用相同的维纳滤子常数。这样在上述模拟中,通过最优化该常数而获得的任何效果被忽略。加入不同水平的高斯噪音,蒙特卡洛模拟显示出输入流速值(F=60ml/100g/分钟)被反卷积过程低估了(图7)。通过重复最优化维纳滤子常数q,在信噪比值的整个范围内获得所建议的方法的稳定性(图7a)。在图7a中,在信噪比的呈现出的范围内,q增大了六倍,而在图7b中使用常数1=0.01。
 TTP[s]   Fest[ml/100g/min](F=30ml/100g/min)    TTP[s]   Fest[ml/100g/min](F=60ml/100g/min)
 Δt0=0s  4.4±0.6     20.2±1.7  4.7±0.7     41.0±5.2
 Δt0=1s  5.2±0.5     20.5±1.7  5.5±0.7     41.4±4.6
 Δt0=2s  6.2±0.6     20.0±1.6  7.0±1.1     40.8±4.8
 Δt0=3s  7.6±0.6     20.0±1.7  8.5±1.8     41.4±4.9
 Δt0=4s  8.6±0.5     20.1±1.7  9.2±1.5     41.0±3.9
 Δt0=5s  9.9±0.8     19.5±1.3  10.5±1.3     40.6±3.7
表2:估计的到峰值时间、TTP和流速、Fest(平均值±标准差)。模拟参数如表1中的情形I。加入了获得噪音。
通过与其中在反卷积中使用片内手工选择的动脉输入函数的常规的方法所获得的结果相比较,来评估建议的在体内测量的盲估计方法的表现。为了比较,反卷积前此动脉输入函数被归一化到单位面积。在两种情况中都计算描述血流量F,平均通行时间MTT和到峰值时间MTP的参数图。
反卷积的造影剂浓度的时程显示出,当使用新方法时,与传统方法相比,在更少的时间帧内可见到造影剂药剂的第一通过,(图8)。这暗示出在两种方法之间在平均通行时间和到峰值时间值上的差异。
体素方面比较显示出,比起传统方法,60%的目标体素在平均通行时间上有大于5%的下降(相对于与传统方法平均下降为19.6%)。在20%的目标体素中看到在平均通行时间上的轻微增加,平均值为19%。当使用新方法时,在40%的目标图象体素中看到了在流速值上超过5%的增大(相对于传统方法平均增加21.2%),而35%的剩下的体素显示出减小的值(相对于传统方法平均下降19.2%)。
与传统方法相比,在新方法中约所有目标图象体素的37%在到峰值时间上有超过5%的增大(平均增加50%),而在两种方法之间,所有其他目标体素的差异小于5%。参数图象是有噪音的(图9),但当使用新方法时(右栏),比起传统方法(左栏),不同的组织结构明显地更容易被看出。
对于定量分析,选择代表正常前白质(WM)和正常皮质灰质(GM)的感兴趣区域(表3)。获得的灰质向白质的流比是相似的,使用盲估方法时,(GM/WM=2.4±0.7),使用手工选择的动脉输入函数方法时,(GM/WM=2.3±0.7)。当使用盲估动脉输入函数时,平均通行时间值更小。当使用盲估方法时,(21.7±0.8s),比起传统方法,(0.8±0.8s),在灰质和白质区域之间的在到峰值时间上的差异更大,使得更容易从视觉上区别组织结构。
               盲估            常规
 GM     WM  GM     WM
 CBF  2.3±0.6     0.9±0.2  2.2±0.6     0.9±0.3
 TTP  29.4±0.5     31.1±0.7  28.9±0.5     29.7±0.9
 MTT  4.2±0.6     5.0±1.1  4.6±0.5     5.2±0.8
表3:灌注参数
示例4
成像研究
从其中男性深海潜水员经过用常规MR成像和灌注MR成像的扩展MR协议检查的数据库中选择对象。从数据库中,选择没有明显病理的十个最年轻的男性潜水员。这些男性的平均年龄是39.4(范围在33-43岁)。
用SiemensVision 1.5 Tesla全身医用扫描仪记录该MR图象。使用动态磁化造影剂成像,所有的患者经过灌注MR检查。使用的是梯度-回波回波-平面顺序,参数为TR/TE 1442ms/60.70ms,翻转(flip)角90度,视场230×230象素而矩阵大小为128×128象素。切片厚度,5mm,具有0.3mm的片间间隔。扫描时间:对11个切片是1分钟21秒,每个切片50次获取。切片取向平行于在前叶基部和脑尾部边界之间划出的直线。造影剂:1mmol/mlgadobutrol。造影剂的量:体重的0.2ml/kg。用速率为5ml/s的电动注射器通过口径18的在肘脉中的i.v.导管完成所有检查。造影剂注射延迟是12秒,其后直接接着20ml的0.9mg/ml氯化钠的注射。
灌注图象被传送的离线计算机上用于后处理。图象处理工具是来自挪威,奥斯陆大学,信息学系的XITE公共领域图象处理包。
第一通过药剂信号和第二通过药剂信号在时间上不会交叠。因此,通过简单地使用第一通过药剂部分在体素特定动脉输入函数确定中,可以避免再循环效应。由于灌注记录是在健康男性上完成的,假定血脑屏障是完整的。因此,在确定体素特定动脉输入函数中不需要包括泄漏项。对所有反卷积,通过简单地在频率域中用同一经验发现的信噪比的负二次方进行维纳过滤而完成反卷积步骤。
对动脉输入函数的手工选择:当确定作为解剖感兴趣区域时,在相同切片中手工选择动脉输入函数体素。从前脑动脉周围选择六个具有最早、最快和显著的信号下降的体素。
解剖区域的选择:从最能确定深灰质核的切片选择所有感兴趣区域。通过从每个头尾状核和豆状核,从左和右前脑皮层及左和右前白质中选择约20个体素来手工提取这些感兴趣区域(看图10a,b)、排除在感兴趣区域内的具有与动脉输入函数相似的信号的体素。它们一定包含了动脉或静脉。
跟在反卷积之后,感兴趣区域遮罩被用于对每个感兴趣区域计算rTTP、rMTT、rCBF和rCBV的平均值和标准差。
信号和统计分析:本发明的盲反卷积方法将每个动脉输入函数归一化到单位面积上。因此,对不同的组织,仅比较rCBF和rCBV的比率是有意义的。然而,请注意,通过合并一种方法以自动地确定大动脉的动脉输入函数可以容易地获得正确的归一化。
理想地,反卷积不应该影响体素特定大脑的脑量。然而,所有反卷积方法都降低信噪比并可能不同地影响不同组织的残留函数。结果是,在灰质和白质的脑血量之间的比值可能在反卷积之后改变。对本发明盲反卷积方法、常规反卷积方法和没有反卷积,计算该比值。
使用两侧Wilcoxon符号秩实验来确定在两种反卷积方法之间的有统计意义的重要差异。对rTTP、rMTT、及rCBF和rCBV的速率进行这样的实验。在每个尾部的小于0.0125的p值被认为是有意义的。当计算显著水平以避免大多数平分时使用两位数。
区域的到峰值的时间:与传统方法相比,用新方法,对灰质和白质两者,到峰值的时间略微变长(1.8%和6.0%,表4)。对新方法来说,在灰质和白质在峰值时间上的相应的差异明显变大(36.8%,表4)。
区域的平均通行时间:对两种方法,对灰质的平均通行时间没有明显的不同。相反,与传统方法相比,新的方法的在灰质和白质的平均通行时间之间的时间差异明显更长(21.3%,表5)。
区域的脑血流量:对两种方法,在灰质和白质之间的流速比没有明显不同(表6)。
区域的脑血量:与传统方法相比,新的方法的在灰质和白质的血量比略微低些。相应的,与没有反卷积的从观测到的数据的量的计算相比,传统方法的在灰质和白质的血量比要低一些(表7)。
                      本发明方法                            传统方法
对象   年龄(岁) 灰质 白质 时差 灰质 白质 时差
    1     43   13.6±0.5   15.3±2.9     1.7   1 2.9±0.7   14.0±1.1     1.1
    2     40   13.7±0.7   16.2±4.2     2.6   13.7±0.6   15.9±5.7     2.2
    3     40   16.0±0.9   17.6±1.1     1.7   15.9±0.8   17.6±2.3     1.8
    4     38   11.3±0.9   13.3±1.1     2.0   11.2±0.8   12.8±0.1     1.7
    5     33   16.3±0.8   20.7±4.8     4.4   1 6.0±0.8   1 8.6±2.7     2.7
    6     40   14.6±1.1   18.4±2.9     3.8   1 5.0±0.8   17.6±1.7     2.6
    7     41   12.0±0.7   14.6±1.4     2.5   11.8±0.8   14.2±2.0     2.4
    8     42   11.7±1.3   13.9±2.6     2.2   10.8±1.4   11.5±2.0     0.7
    9     41   19.2±1.3   22.4±2.8     3.3   18.7±1.1   21.5±2.9     2.7
    10     36   16.4±1.1   19.7±3.2     3.3   16.3±0.8   18.5±2.9     2.2
    均值     39.4   14.5±2.5  17.2±3.1     2.8±0.9   14.2±2.6   16.2±3.1     2.0±0.7
表4.区域的到峰值时间。单位,秒。显著性水平:灰质,p=0.024。白质,p=0.002。时差,p=0.002
                   本发明方法                        传统方法
对象   年龄(岁) 灰质 白质 时差 灰质 白质 时差
    1     43     5.0±0.7     6.6±0.7     1.6   5.4±0.6   6.5±1.2     1.1
    2     40     5.6±0.8     7.5±0.7     1.9   5.7±0.7   7.2±1.7     1.5
    3     40     6.1±0.9     7.6±0.9     1.5   5.9±0.8   7.3±1.4     1.4
    4     38     4.8±0.7     6.2±1.3     1.3   5.4±0.7   6.3±1.1     0.9
    5     33     6.3±0.9     7.9±0.9   1.6   6.4±0.7   7.8±1.4     1.5
    6     40     6.7±1.0     7.9±0.8   1.2   6.5±0.8   7.4±1.1     0.9
    7     41     5.1±0.9     6.2±0.9   1.1   5.0±0.7   6.2±1.4     1.2
    8     42     6.6±1.1     7.6±0.6   1.0   6.7±1.0   7.4±1.3     0.6
    9     41     6.9±1.0     7.8±0.8   0.9   6.5±0.7   7.3±1.4     0.8
    10     36     6.4±0.8     7.4±0.9   1.1   6.4±0.6   7.4±1.3     1.0
    均值     39.4     5.9±0.8     7.3±0.7   1.3±0.3   6.0±0.6   7.1±0.5     1.1±0.3
表5.区域的平均通行时间。单位,秒。显著性水平:灰质,p=0.41。白质,p=0.022。时差,p=0.001
                     本发明方法                           传统方法
对象   年龄(岁) 灰质 白质 时差 灰质 白质 时差
    1     43     2.3±0.4     0.9±0.2     2.5   2.3±0.5   0.9±0.2     2.6
    2     40     2.0±0.5     0.6±0.2     3.3   2.3±0.5   0.7±0.2     3.3
    3     40     2.2±0.5     0.9±0.3     2.5   2.5±0.6   1.0±0.3     2.5
    4     38     2.8±0.6     1.4±0.4     2.0   3.2±0.7   1.5±0.4     2.1
    5     33     2.3±0.4     1.1±0.2     2.1   2.5±0.5   1.1±0.3     2.4
    6     40     2.3±0.4     0.8±0.2     2.7   3.0±0.6   1.1±0.2     2.9
    7     41     2.4±0.4     1.3±0.3     1.9   2.6±0.6   1.2±0.3     2.2
    8     42     1.5±0.4     0.8±0.2     2.0   1.3±0.5   0.9±0.5     1.5
    9     41     1.8±0.5     0.7±0.2     2.5   2.2±0.5   0.9±0.2     2.6
    10     36     2.2±0.3     0.9±0.3     2.5   2.5±0.4   0.9±0.2     2.8
    均值     39.4     2.2±0.4     0.9±0.3     2.4±0.4   2.4±0.5   1.0±0.2     2.5±0.5
表6.区域的脑血流量。对比率的显著性水平,p=0.097。
    对象  年龄(岁)  本发明方法   传统方法   无反卷积
    1     43     2.0     2.3     2.4
    2     40     2.5     2.5     3.3
    3     40     2.0     2.1     2.2
    4     38     1.6     1.9     1.7
    5     33     1.8     1.9     1.9
    6     40     2.3     2.6     2.6
    7     41     1.6     1.8     2.2
    8     42     1.7     1.4     2.3
    9     41     2.3     2.4     2.8
    10     36     2.3     2.5     2.6
    均值     39.4     2.0±0.3     2.1±0.4     2.4±0.5
表7.区域脑血量比。显著性水平;新方法和传统方法,p=0.042;传统方法和无反卷积,p=0.092。
示例5
模拟和成像
如上所述产生综合数据。参数被调整成在典型的,体内磁共振T2*灌注检查期间获得的模拟数据。完成这些模拟以调查本发明的方法在存在再循环时表现如何。具有8s延迟和30s时间常数的分散的指数函数与上面给出的两个动脉函数(和组织残留函数)中的任一卷积。
使用Siemens Vision 1.5 T全身医用扫描仪在一个在脑左前叶上具有肿瘤的患者身上完成体内成像。灌注检查是临床估计,常规程序的一部分。
选择标准的,商用的梯度-回波平面影像(GE-EPT)顺序,参数TR/TE=1442ms/60.70ms,FoV=230×230,矩阵=128×128。在12秒后,自动注射造影剂,二亚乙基钆三胺戊乙酸,(Gd-DTPA,Schering)。记录了11个平行于头颅底部的截面并对于每个体素测量了50个时间点。图象数据被传送到离线计算机上用于后处理。该过程在C中进行,包括几个数字处方程序。基于XITE包和MATLAB的图象处理工具被用于图象操作,显象和估计。
在将所要需要的数据进行了对数变换之后,每个图象数列的3到14帧被用于计算每个体素的基线组织信号。这个基线常数信号被从对数变换的体素信号中减去以给出体素特定造影信号。
用Linux下的脚本完成该计算,其中小的处理步骤被储存到带有奔腾II主板的个人电脑上的盘上。在这样的构架内对每个部分的计算时间约为20分钟。在新式计算机上的简单最优化将以50-100的倍数缩短该计算时间。
为了在小的局部区域内估计组织特定参数,在屏幕上交互式地划出小的组织遮罩。对每个组织种类标记20个以上的体素。对灰质和白质的遮罩在左侧顶骨区域。肿瘤的区域1位于外周,而区域2位于具有最大造影剂信号的部分的中心。
仅使用对体素特定动脉输入函数的复倒频谱估计,包含泄漏和再循环项减小了组织特定血液动力学量和时间参数的变化(表8和9,模型1和2)。用迭代链反卷积对体素特定动脉输入函数估计和组织残留函数估计的精炼,进一步减小了组织特定血液动力学量和时间参数的变化(表8和9,模型2和5)。在变化的减小中,泄漏项比再循环项更重要(表8和9,模型2、3和5)。
在组织特定血液动力学量和时间参数上的变化的减小,随完整的模型及迭代链估计,使得不同的组织成分更容易在参数图象中视觉上分开(图11)。特别是,具有高血量、高血流量、混合的到峰值的时间以及相当长的平均通行时间,肿瘤变得更明显了。同样,组织残留函数的时程也变得平稳很多并且持续时间明显更短(图12)。
对灰质、白质和两个肿瘤区域,估计的泄漏系数都较低并且是相似的(表10)。忽略再循环对估计的泄漏系数有很小的影响(表10,图13a,b)。这样,在第一通过药剂已经通过以后(图13,c,d),在肿瘤区域中的明显的造影剂信号的持续性,不得不用由于高血流量引起的比在正常脑皮层中更大的泄漏来解释(看表8)。在稳态泄漏图象中在脑膜薄膜和脉络丛中的高的造影剂信号(图13,c,d),符合这样的事实,即在这些组织中的毛细管可透过造影剂。
从心脏和背部的循环路线对白质来说要长于灰质。这是因为脑动脉在进入白质之前通过灰质进入脑。因此,对白质来说,再循环时间应该长于灰质。
用迭代链估计,整个模型的组织特定再循环时间与这一预测一致(表10和图13e)。相比之下,当忽略泄漏项时,对灰质、白质和两个肿瘤区域的再循环时间都变长并且更加相似(表10和图13f)。
模型     灰质      白质   肿瘤区域1   肿瘤区域2
    脑血量(相对单位)
    12     3.36±0.103.03±0.09     1.00±0.051.00±0.05     5.28±0.224.35±0.09     8.22±0.395.94±0.18
    345     3.03±0.092.97±0.082.94±0.08     1.00±0.041.00±0.041.00±0.04     4.35±0.124.41±0.104.23±0.10     6.05±0.145.75±0.155.67±0.15
    脑血流量(相对单位)
    12345     4.21±0.163.47±0.123.35±0.083.12±0.073.10±0.07     1.00±0.061.00±0.061.00±0.051.00±0.051.00±0.05     7.29±0.435.08±0.174.67±0.124.40±0.114.29±0.10     8.11±0.355.83±0.235.40±0.135.47±0.115.34±0.10
表8.血液动力学量参数。平均值±标准差。体素数字,灰质,54;白质,59,肿瘤区域1,30;肿瘤区域2,22。所有的值都是相对于白质的平均值。模型1,仅仅是复倒频谱估计,既没有包括再循环也没有包括泄漏。模型2,仅仅是复倒频谱估计,包括了再循环和泄漏。模型3,迭代估计,包括了再循环,没包括泄漏。模型4,迭代估计,没包括再循环,包括了泄漏。模型5,迭代估计,包括再循环和泄漏。
 模型      灰质     白质   肿瘤区域1   肿瘤区域2
    到峰值的时间(帧)
    12345     4.26±0.084.56±0.086.19±0.056.52±0.086.54±0.07     5.87±0.245.81±0.157.52±0.137.34±0.137.36±0.13     3.30±0.133.27±0.083.93±0.115.30±0.085.40±0.09     3.77±0.433.73±0.112.73±0.115.73±0.095.82±0.08
    平均通行(帧)
    12345     4.04±0.113.56±0.084.02±0.073.90±0.063.91±0.06     5.40±0.194.34±0.144.81±0.194.35±0.124.36±0.12     3.71±0.133.51±0.104.22±0.094.14±0.064.20±0.06     5.01±0.114.15±0.125.05±0.164.34±0.114.38±0.09
表9.血液动力学参数。平均值±标准差。模型3,迭代估计,包括再循环,不包括泄漏。模型4,迭代估计,不包括再循环,包括泄漏。模型5,迭代估计,包括再循环和泄漏。在两个连续帧之间的时距:1440ms。
模型     灰质     白质   肿瘤区域1   肿瘤区域2
    泄漏系数
    56     0.019±0.0010.019±0.001   0.019±0.0020.021±0.002   0.015±0.0020.015±0.002   0.022±0.0010.022±0.001
    再循环时   间(帧)
    57     15.8±0.418.1±0.3   18.0±0.618.3±0.5   15.5±0.218.3±0.3   15.0±0.217.5±0.1
表10.泄漏和再循环。皮卡尔迭代。平均值±标准差。模型5,包括再循环和泄漏。模型6,包括再循环,不包括泄漏。模型7,不包括再循环,包括泄漏。在两个连续帧之间的时距:1440ms。
示例6
模拟
为了估计理查德森.露西迭代盲反卷积方法的目的,模拟第一通过信号的综合数据。在少数体内测量中估计泄漏和再循环模型的影响。
对于第一通过数据,使用在傅立叶域中的简单相乘,完成在综合产生的动脉输入函数和组织残留函数之间的卷积,c(t)=u(t)*v(t)。用时间t的伽马变量函数,作动脉输入函数的模型
v ( t ) = c 0 ( t - t 0 ) &alpha; e - ( t - t 0 ) / &beta; - - - ( 23 )
其中,t0是造影剂到达时间,c0、α和β参数确定伽马函数的形状。对所有t≥t0的值,动脉输入函数的表达都是正确的,在其他情况,则为零。选择值α=3和β[1.2,2.0]以研究具有不同形状动脉输入函数。对所有动脉输入函数,c0被用于将动脉输入函数归一化到单位面积上。用线性或指数形状来模拟被脑血流量CBF定比的组织残留函数,
u exp ( t ) = CBF &CenterDot; e - 1 MTT , t &GreaterEqual; 0 - - - ( 24 )
u lin ( t ) = CBF &CenterDot; ( 1 - 1 MTT ) , t &le; 2 MTT - - - ( 25 )
Hbox(t)=CBF    t≤MTT    (26)
对所有其他的时间点,为零(
Figure 058225254_5
stergaard等人)。对固定的平均通行时间MTT=4秒,研究了较大范围的流速值,CBF[10,60]ml/100g/分钟。
卷积的数据c(t),通过下式,被变换成MR信号强度,s(t),
s ( t ) = s 0 ( t ) e - k &CenterDot; TE &rho; &kappa; H &CenterDot; c ( t ) + N &CenterDot; &eta; ( t ) - - - ( 27 )
假定在脉管造影剂浓度和观测到的磁化造影剂之间为线性关系。TE是该顺序的回波时间而K是假定对所有组织都共有的实验特定常数(这里设定为一),ρ是脑组织的密度,KH是用于在毛细管和大脉管之间的血细胞比容差异的修正因子(McCallum Optics Communications 75:101-105(1990))。基线信号强度,s0,被设定为通常在对正常灰质的灌注检查中得到的值(500MR单位)。最后,加入来自在真实体内检查中外部对象边界的图象噪音η(t),t=1,…50,使得64×64图象体素具有单独的噪音变化。在加入图象噪音之前,减去所有噪音体素的平均值。用不同比例因子(N=0,1,2.5,4)加入噪音以获得不同信噪比。这分别给了预造影剂噪音标准差,σs=0,15,30,50。
示例7
MRI实验
从具有不同脑脉管异常的两个患者上获得的数据被用于检验所建议的这套方法;一个患者在脑的左前叶中有肿瘤(anaplastic coligoastrocytoma),一个患者遭受纤维肌性发育异常。
体内灌注检查都是临床估计、常规程序的一部分。没有在患者身上作实验研究,而获得的数据的使用与医院保险单一致。选择标准的,卖方提供的梯度-回波回波平面成像顺序(GE-EPI,Siemens Visionl.5T,埃尔兰根,德国),选择参数为TR/TE 1442ms/60.70ms,翻转角90度,视场230×230,象素矩阵大小为128×128象素。使用MR兼容电动注射器(Medirad,匹兹堡)在12秒后自动注射造影剂Gadovist(Schering,德国),其剂量按照体重(0.2mmol/kg)。与用于造影剂注射的一样,用相同的注射速度(5ml/s)完成随后的盐水冲洗。对总计五十个时间点重复测量十一个图象切片,其与真正的轴获取轻微地倾斜以避免在额窦中的磁感伪影。
示例8
在MATLAB(Mathworks,Natic,MA,美国)中产生综合数据。随后使用NordicNeuroLab AS(NordicNeuralab AS,卑尔根,挪威)提供的另外市售可得的nICE软件中的专门的扩展的工具处理和显示这些数据和获得的图象数据。作者使用C/C++完成了建议的这套方法,包括几个数字处方程序(Press等人,Numerical Recipes in c,剑桥大学出版社(1992))。在2.2GHz奔腾III英特尔处理器上使用当前工具的计算时间,对具有50个时间点测量的体内图象切片(128×128体素)图象切片,是15秒。
在对数变换该模拟的或获得的灌注时间数列以后,平均约十个预造影剂图象以计算每个体素的基线组织信号。预造影剂图象范围和对第一通过的限制对每个时间数列来说是被交互选择的,但是在数列中对所有体素是共有的。从被对数变换的体素信号中减去基线以给出体素特定造影剂信号。不进行平滑。
为了研究泄漏和再循环模型的影响,从患者数据中选择灰质、白质和病理区域中的小的(6×6体素)感兴趣区域,对其求均值和接续不同的后续处理步骤。计算在不同区域中对m、κi、Ti的估计。血量作为在减去泄漏和再循环后得到的在第一通过信号下面的面积来计算。这个血量与作为在所有获得的时间点下的面积计算出的、作为选定的代表第一通过的时间点下的面积计算出的、和作为伽马拟合的第一通过函数下的面积计算出的血量相比较。
在反卷积以后,对所有时间数列,在体素方面计算脑血流量、脑血量(定义为反卷积的造影剂浓度的面积)、平均通行时间和到峰值时间的血液动力学图。分别将血流量和血量确定为反卷积造影剂信号的最大值和面积。随后由中间体积定理,MTT=CBV/CBF,计算平均通行时间,而到峰值时间被设定为其中出现估计的动脉输入函数的最大值的那个时间点。在体内数据中的动脉输入函数的归一化,影响血流量和血量的绝对值时,平均通行时间值与当使用动脉输入函数的真实面积时获得的那些是可比拟的。
对模拟的数据,计算出的血液动力学参数被和真实的输入值相比较。模拟使得可以检验用于再现(i)不同形状的动脉输入函数、(ii)不同流动条件、(iii)不同形状的组织残留函数(即不同的脉管结构)、(iv)不同信噪比的迭代算法。
对于体内数据,体素特定动脉输入函数和反卷积的组织残留函数被显现为时间的函数。另外,通过在计算出的血液动力学图中交互选择组织区域(14×14体素)来完成对感兴趣区域的分析。这些结果与在相同区域中使用以平截奇异值分解(固定的平截阈值为0.2)(
Figure 058225254_6
stergaard等人的MGM 50:164-174(2003))的标准方法获得的值相比较。在标准方法中,其在下面指的是整体方法,动脉造影剂信号在所有体素中都被用作动脉输入函数。通过找出在造影剂信号强度上具有最早和最大变化的切片内图象体素来半自动地选择整体动脉输入函数。随后一伽马变量函数被拟合到这个动脉信号中以避免再循环的影响,图20。因为相容性的原因,这个整体动脉输入函数在反卷积前被归一化到单位面积。
示例9
结果
新的方法对泄漏和再循环给出了良好的体素方面的近似,否则其污染第一通过信号(选出的例子图14,图15)。在几次迭代内获得了估计参数m、κi、Ti的收敛(图14c),因而对所有患者和所有图象体素,迭代的次数保持为固定(总计5次迭代)。它足以保持参数m对一个患者的所有体素是不变的,但对每个患者该参数必须被最优化。例如,在肿瘤患者中,值为m=0.05,而在具有纤维肌性发育异常的患者中它增大到m=0.15。用固定的高斯模型取代在人体内造影剂的分散v0是粗糙的,但在这一研究中对所有患者都已足够。
当如这个研究建议的那样确定造影剂第一通过时,与使用伽马拟合的确定相比,血量较少了。(表11)
模型     灰质     白质   肿瘤区域1   肿瘤区域2
    1234     130±595±792±582±5     98±369±463±455±4     384±10259±14250±12211±11     322±10228±13222±11183±11
表11:  相对脑血量计算。平均值±标准差。每个感兴趣区域(ROI)包含6×6体素。模型1)在所有测得的时间帧下的面积。模型2)在固定数量的时间帧(手工去掉交互选择的但对所有图象体素是共同的帧)下的面积。模型3),在伽马拟合的第一通过(非线性Levenberg-Marquant最小平方拟合)下的面积。模型4)在如本文件中所建议那样去掉了再循环和泄漏项之后剩下的第一通过信号下的面积。由于归一化的动脉输入函数,单位ml/100g以一常数定比。
所建议的方法还显示出,与使用其他策略计算的血量相比,在感兴趣区域内的值的变化减小了(表11)。
对灰质和白质,估计的泄漏系数是同样的低,κi=1.6,但对两个肿瘤区域则增大了,κi=5.7(图16)。与灰质相比,再循环时间Ti对白质是更类似的(Ti=20.2秒)。
迭代反卷积
计算机模拟
在没有噪音时,在一些模拟的设定中用迭代盲反卷积精确地再现了输入流动值(图17a)。然而当对动脉输入函数使用不同参数设定时(图17a)可确定在估计的流动值中的较小的变化。
参数a=3和b=1.8的动脉输入函数,以及MTT=4秒的指数组织残留函数,可近似重现在没有噪音情况下的真实的输入流动值(图17b)。平均流动值也在存在更大噪音水平时被重现出来(示出方程(27)中N=0、1、2.5和4时的结果)。当噪音水平较高时(N=4),平均流动值被高估并且在数值中的差异增大(图17b)。
迭代盲反卷积方法对合理噪音环境下(N=1,方程(27))的指数和线性形状的组织残留函数给出的类似的流动估计。对盒形的组织残留函数,流动值一致是更小,(图17c)。对所有的模拟流动值重现了输入平均通行时间值(4秒),(图17d)。
体内分析
在两个患者上估计的动脉输入函数与从完整的放射检查诊断出的患者的血液供应模式相关紧密(图18)。在肿瘤患者中,造影剂到达部分肿瘤区域与到达脑供血动脉一样早。而后,造影剂到达灰质比到达白质要略早一点(图18,帧A ii和iii)。在发育异常的患者中,造影剂到达左侧半球要比到达受发育异常影响的右侧区域早。就像通过比较在左右半球中从造影剂到达到造影剂冲出的帧数能看到的那样,造影剂药剂通过中间右侧区域也较慢。这与根据患者的脉管造影术检查在右半球中确定的侧枝循环符合得很好。再一次,造影剂到达正常灰质比到达正常白质要早(图18,帧B ii和iii)。
当将组织残留函数显现为时间的函数时,可看到详细的造影剂冲洗动力学(图19)。在肿瘤患者中,在肿瘤的中心部分看到的造影剂冲洗较慢的区域(图19参考下文,黑箭头)。同样在发育异常患者中的右半球的区域中看到冲洗与对侧半球相比要更慢。
对于在右半球中具有纤维肌性发育异常的患者,发现参数(a=4.13,b=1.53)的伽马函数来表示具有最早和最大信号变化的体素(图20,7个体素)。这个伽马函数在整体方法中被用作动脉输入函数。使用建议的方法或传统方法获得的在这个患者中的到峰值时间图,示出由右中间脑动脉通常供应的区域到峰值的时间较晚(图20)。在建议的方法中,到峰值的时间上的差异更为明显。
在有纤维肌性发育异常的患者中的感兴趣区域分析在建议的方法和传统的方法之间给出可比拟的血液动力学参数(表12,图20)。用建议的方法所获得的平均通行时间要略微短于用传统方法获得的那些,但如果对目前的信噪比选择小于0.2的平截因子的话可能会相持平。该差异也可以用不同的确定第一通过的策略来解释。对建议的方法,平均通行时间和到峰值时间值的变化要更小。对建议的方法,可以看到在受影响的右半球中的平均通行时间(ROI2、4、6)要长于在健康的白质中的(ROI3、5),其又长于在健康的皮层灰质中的(ROI1)。在传统方法中对应的值不是那么易于辨识,除了最下面的区域(ROI6)。与建议的方法相反,到峰值时间和平均通行时间值这两者在ROI6中都要长于在附近位置的ROI4中的。建议的方法显示出,与ROI4中的相比,在ROI6中的动脉输入函数更晚且更分散(图18,灰色箭头)。
 ROI M1:    MTT(s)     Ml:TTP(s)     M2:MTT(s)     M2:TTP(s)
    123456     4.4±0.34.8±0.54.6±0.55.6±0.64.7±0.45.4±0.7     24.3±0.625.6±1.425.3±0.826.2±1.024.2±0.727.4±1.0     5.8±0.65.4±0.95.2±0.85.1±1.15.5±0.57.2±2.0     24.3±0.526.3±1.425.4±1.027.4±1.524.8±0.928.7±1.0
表12:感兴趣区域(ROI)分析
具有右中间脑动脉纤维肌性发育异常的患者。平均值±标准差。M1:用建议的泄漏和再循环模型在确定第一通过后的迭代反卷积。M2:在通过仅包括50个图象帧的前27个确定第一通过后的平截奇异值分解。

Claims (12)

1.一种对有脉管的动物对象的感兴趣区域的灌注磁共振成像方法,所述方法包括:
在从造影剂到达所述对象的脉管系统中所述感兴趣区域之前到至少所述造影剂穿过所述感兴趣区域的第一通过的末端期间内在一系列时间值t上确定所述感兴趣区域的体素i的磁共振信号强度si(t);
由所确定的信号强度值si(t)和一体素特定动脉输入函数vi(t),对每个所述体素,通过使用理查德森-露西反卷积算法估计组织残留函数ri(t)及体素特定动脉输入函数vi(t),以及通过使用维纳反卷积滤子由所估计的体素特定动脉输入函数vi(t)重新估计组织残留函数ri(t),确定组织残留函数ri(t)的值;以及
由所确定的ri(t)的值产生所述感兴趣区域的图象。
2.一种灌注磁共振成像数据处理方法,所述方法包括:
获得对其中脉管中已被给予磁共振造影剂药剂的有脉管对象中的感兴趣区域的体素i的一组磁共振信号强度值si(t),所述信号强度值是关于在从造影剂到达所述感兴趣区域之前到至少所述造影剂穿过所述感兴趣区域的第一通过的末端期间的一系列时间值t的信号强度值;
由所述组信号强度值si(t)和一体素特定动脉输入函数vi(t),对每个所述体素,通过使用理查德森-露西反卷积算法估计组织残留函数ri(t)及体素特定动脉输入函数vi(t),以及使用通过使用维纳反卷积滤子由所估计的体素特定动脉输入函数vi(t)重新估计组织残留函数ri(t),确定组织残留函数ri(t)的值;以及
由所确定的ri(t)的值产生所述感兴趣区域的图象。
3.如权利要求1或2所述的方法,其中每个体素的造影信号被分成没有泄漏和有泄漏成分的第一通过和第二通过部分。
4.如权利要求1或2所述的方法,其中在一迭代循环中依次使用理查德森-露西反卷积算法和维纳反卷积滤子。
5.如权利要求4所述的方法,其中通过使用维纳反卷积滤子找到对体素特定动脉输入函数的新的估计而结束循环。
6.如权利要求1或2所述的方法,其中使用理查德森-露西反卷积算法来估计组织残留函数,而后估计体素特定动脉输入函数。
7.如权利要求6所述的方法,其中在一迭代循环中使用理查德森-露西反卷积算法。
8.如权利要求7所述的方法,其中通过找到对体素特定动脉输入函数的新的估计而结束循环。
9.如权利要求5所述的方法,其中在六次迭代后结束迭代循环。
10.如权利要求8所述的方法,其中在六次迭代后结束迭代循环。
11.如权利要求1或2所述的方法,其中产生一参数灌注图象。
12.一种灌注磁共振成像装置,包括
磁共振信号强度确定部件,用于在从造影剂到达所述对象的脉管系统中所述感兴趣区域之前到至少所述造影剂穿过所述感兴趣区域的第一通过的末端期间内在一系列时间值t上确定所述感兴趣区域的体素i的磁共振信号强度si(t);
组织残留函数确定部件,用于由所确定的信号强度值si(t)和一体素特定动脉输入函数vi(t),对每个所述体素,通过使用理查德森-露西反卷积算法估计组织残留函数ri(t)及体素特定动脉输入函数vi(t),以及通过使用维纳反卷积滤子由所估计的体素特定动脉输入函数vi(t)重新估计组织残留函数ri(t),确定组织残留函数ri(t)的值;以及
图象产生部件,用于由所确定的ri(t)的值产生所述感兴趣区域的图象。
CN2005800225254A 2004-05-04 2005-05-04 动脉输入和组织残留函数在灌注磁共振成像中的盲确定 Expired - Fee Related CN101002104B (zh)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
US56764704P 2004-05-04 2004-05-04
US60/567,647 2004-05-04
PCT/GB2005/001680 WO2005106517A2 (en) 2004-05-04 2005-05-04 Blind determination of the arterial input and tissue residue functions in perfusion mri

Publications (2)

Publication Number Publication Date
CN101002104A CN101002104A (zh) 2007-07-18
CN101002104B true CN101002104B (zh) 2012-01-18

Family

ID=35033284

Family Applications (1)

Application Number Title Priority Date Filing Date
CN2005800225254A Expired - Fee Related CN101002104B (zh) 2004-05-04 2005-05-04 动脉输入和组织残留函数在灌注磁共振成像中的盲确定

Country Status (7)

Country Link
EP (1) EP1747477B1 (zh)
JP (1) JP5171251B2 (zh)
CN (1) CN101002104B (zh)
CA (1) CA2565668A1 (zh)
DK (1) DK1747477T3 (zh)
NO (1) NO20065557L (zh)
WO (1) WO2005106517A2 (zh)

Families Citing this family (27)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7620227B2 (en) * 2005-12-29 2009-11-17 General Electric Co. Computer-aided detection system utilizing temporal analysis as a precursor to spatial analysis
JP5536974B2 (ja) * 2006-11-08 2014-07-02 株式会社東芝 X線診断装置及び画像処理装置
GB0708136D0 (en) * 2007-04-26 2007-06-06 Cancer Res Inst Royal Measure of contrast agent
US20100172842A1 (en) * 2007-05-16 2010-07-08 Yeda Research And Development Co., Ltd. At The Weizmann Institute Of Science Assessment of blood-brain barrier disruption
JP5524860B2 (ja) * 2007-12-28 2014-06-18 ブラッコ・シュイス・ソシエテ・アノニム 医療画像用途における固定化された造影剤の定量分析
US10130342B2 (en) 2007-12-28 2018-11-20 Bracco Suisse Sa Initialization of fitting parameters for perfusion assessment based on bolus administration
JP5454882B2 (ja) * 2009-08-31 2014-03-26 ジーイー・メディカル・システムズ・グローバル・テクノロジー・カンパニー・エルエルシー 血流動態解析装置、磁気共鳴イメージング装置、およびプログラム
US20130039553A1 (en) * 2010-03-09 2013-02-14 Aarhus Universitet Method for obtaining a blood flow parameter
CN101912262B (zh) * 2010-07-22 2012-07-04 中国科学院深圳先进技术研究院 磁共振成像定量参数计算方法及系统
JP6133794B2 (ja) * 2011-03-18 2017-05-24 コーニンクレッカ フィリップス エヌ ヴェKoninklijke Philips N.V. かん流比較及び定量化のためのデータの動的正規化
WO2012154714A1 (en) 2011-05-11 2012-11-15 Tyk America, Inc. Degasser snorkel with serpentine flow path cooling
US9644246B2 (en) 2011-05-11 2017-05-09 Tyk America, Inc. Degasser snorkel with serpentine flow path cooling
WO2014065340A1 (ja) 2012-10-23 2014-05-01 株式会社東芝 画像処理装置、医用画像診断装置及び画像処理方法
JP6100603B2 (ja) * 2013-05-13 2017-03-22 東芝メディカルシステムズ株式会社 医用画像撮影解析装置
US20160125597A1 (en) 2013-06-03 2016-05-05 Samsung Life Public Welfare Foundation Novel magnetic resonance image technique for imaging and evaluating collateral circulation
US10660597B2 (en) * 2013-09-19 2020-05-26 Aarhus Universitet Method for estimating perfusion indices
JP6523326B2 (ja) 2014-04-07 2019-05-29 ブラッコ・シュイス・ソシエテ・アノニムBracco Suisse SA 非基本波分析による音響レベルのIn−Situ推定
FR3027115B1 (fr) * 2014-10-13 2019-05-10 Olea Medical Systeme et procede pour estimer une quantite d'interet d'un systeme dynamique artere/tissu/veine
CN106474634A (zh) * 2015-11-17 2017-03-08 南京中硼联康医疗科技有限公司 基于医学影像数据的几何模型建立方法
JP7049993B2 (ja) * 2015-11-30 2022-04-07 ジーイー・ヘルスケア・アクスイェ・セルスカプ Mri造影剤の併用を含む製剤
WO2017097738A1 (en) 2015-12-10 2017-06-15 Bracco Suisse Sa Detection of immobilized contrast agent with dynamic thresholding
CN107243093B (zh) * 2017-06-07 2020-05-29 上海联影医疗科技有限公司 一种灌注处理的方法及装置
CN108814601B (zh) * 2018-05-04 2021-12-07 浙江工业大学 基于动态对比增强mri的生理参数定量统计优化方法
CN109858513B (zh) * 2018-12-21 2021-01-29 中国科学院自动化研究所 基于脑灰质白质形态特征多重降维的脑认知能力测量方法
CN110236544B (zh) * 2019-05-29 2023-05-02 中国科学院重庆绿色智能技术研究院 基于相关系数的中风灌注成像病变区域检测系统及方法
KR102272741B1 (ko) * 2019-07-11 2021-07-02 가톨릭대학교 산학협력단 4차원 자기공명 혈관조영술의 영상정보 후처리를 통한 3차원 감산 동맥조영술과 3차원 감산 정맥조영술 및 4차원 컬러 혈관조영술의 동시 구현 방법과 의료영상 시스템
CN116342605B (zh) * 2023-05-30 2023-08-11 杭州脉流科技有限公司 Ct灌注影像参数估计方法、装置、设备和存储介质

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1455872A (zh) * 2000-09-12 2003-11-12 安盛药业有限公司 用核自旋极化的mr造影剂对样品进行磁共振研究的方法

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP1385429B1 (en) * 2000-10-25 2008-09-03 The John P. Robarts Research Institute Method and apparatus for calculating blood flow parameters
JP5143986B2 (ja) * 2001-09-19 2013-02-13 ジーイー・メディカル・システムズ・グローバル・テクノロジー・カンパニー・エルエルシー 脳血流計測装置
JP4363833B2 (ja) * 2001-10-16 2009-11-11 株式会社東芝 局所血流動態に関するインデックスを演算する方法及び装置
JP4041012B2 (ja) * 2002-06-03 2008-01-30 ジーイー・メディカル・システムズ・グローバル・テクノロジー・カンパニー・エルエルシー Cbf定量解析方法および装置
JP4013128B2 (ja) * 2002-09-12 2007-11-28 株式会社日立メディコ 血流動態解析装置、方法、及び画像診断装置

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1455872A (zh) * 2000-09-12 2003-11-12 安盛药业有限公司 用核自旋极化的mr造影剂对样品进行磁共振研究的方法

Non-Patent Citations (10)

* Cited by examiner, † Cited by third party
Title
C.KIEFER ET AL.Measurement of cerebral blood flow in MRI using anonparametric deconvolution technique with energycompensation.PROC.INTL.SOC.MAG.RESON.MED.9.2001,1582.
C.KIEFER ET AL.Measurement of cerebral blood flow in MRI using anonparametric deconvolution technique with energycompensation.PROC.INTL.SOC.MAG.RESON.MED.9.2001,1582. *
D.Y.RIABKOV ET AL.Estimation of Kinetic Parameters Without Input Functions:Analysis of Three Methods for Multichannel BlindIdentification.IEEE TRANSACTIONS ON BIOMEDICAL ENGINEERING49.2002,491318-1327.
D.Y.RIABKOV ET AL.Estimation of Kinetic Parameters Without Input Functions:Analysis of Three Methods for Multichannel BlindIdentification.IEEE TRANSACTIONS ON BIOMEDICAL ENGINEERING49.2002,491318-1327. *
F.CALAMANTE ET AL.Delay and Dispersion Effects in Dynamic SusceptibilityContrast MRI: Simulations Using Singular Valuedecomposition.MAGNETIC RESONANCE IN MEDICINE44.2000,44466-473.
F.CALAMANTE ET AL.Delay and Dispersion Effects in Dynamic SusceptibilityContrast MRI: Simulations Using Singular Valuedecomposition.MAGNETIC RESONANCE IN MEDICINE44.2000,44466-473. *
OSTERGAARD L ET AL.High Resolution Measurement of Cerebral Blood FlowUsing Intravascular Tracer Bolus Passages. Part II:Experimental Comparison and Preliminary Results.MAGNETIC RESONANCE IN MEDICINE, ACADEMIC PRESS, DULUTH, MN, US36.1996,36726-736.
OSTERGAARD L ET AL.High Resolution Measurement of Cerebral Blood FlowUsing Intravascular Tracer Bolus Passages. Part II:Experimental Comparison and Preliminary Results.MAGNETIC RESONANCE IN MEDICINE, ACADEMIC PRESS, DULUTH, MN, US36.1996,36726-736. *
吕衍春等.磁共振灌注成像的原理及其在肿瘤诊断中的应用.国外医学临床放射学分册 6.1999,(6),341-345.
吕衍春等.磁共振灌注成像的原理及其在肿瘤诊断中的应用.国外医学临床放射学分册 6.1999,(6),341-345. *

Also Published As

Publication number Publication date
JP5171251B2 (ja) 2013-03-27
NO20065557L (no) 2007-02-05
WO2005106517A2 (en) 2005-11-10
CN101002104A (zh) 2007-07-18
DK1747477T3 (da) 2013-06-17
JP2007536048A (ja) 2007-12-13
EP1747477A2 (en) 2007-01-31
CA2565668A1 (en) 2005-11-10
WO2005106517A3 (en) 2005-12-15
EP1747477B1 (en) 2013-03-13

Similar Documents

Publication Publication Date Title
CN101002104B (zh) 动脉输入和组织残留函数在灌注磁共振成像中的盲确定
US8326400B2 (en) Method of MR imaging
EP2038841B1 (en) Method of perfusion imaging
Chappell et al. Introduction to perfusion quantification using arterial spin labelling
US8099149B2 (en) MRI method for quantification of cerebral perfusion
US11375918B2 (en) White matter fibrography by synthetic magnetic resonance imaging
Lüdemann et al. Brain tumor perfusion: Comparison of dynamic contrast enhanced magnetic resonance imaging using T1, T2, and T2* contrast, pulsed arterial spin labeling, and H215O positron emission tomography
WO2016167047A1 (ja) 磁気共鳴イメージング装置及び画像作成方法
Li et al. Deepvolume: Brain structure and spatial connection-aware network for brain mri super-resolution
Wielopolski et al. Breath-hold MR cholangiopancreatography with three-dimensional, segmented, echo-planar imaging and volume rendering
Haji-Valizadeh et al. Rapid reconstruction of four-dimensional MR angiography of the thoracic aorta using a convolutional neural network
Attenberger et al. Retrospective respiratory triggering renal perfusion MRI
US20240062061A1 (en) Methods for training a cnn and for processing an inputted perfusion sequence using said cnn
Grüner et al. Magnetic resonance brain perfusion imaging with voxel‐specific arterial input functions
Kamesh Iyer et al. Quantitative susceptibility mapping using plug-and-play alternating direction method of multipliers
Risse et al. MR Perfusion in the Lung
Meurée Arterial spin labelling: quality control and super-resolution
Taxt et al. Single voxel vascular transport functions of arteries, capillaries and veins; and the associated arterial input function in dynamic susceptibility contrast magnetic resonance brain perfusion imaging
Friedli et al. Quantitative 3D Vuse Pulmonary MRA
Nissan Development and Applications in Breast and Pancreatic MRI: Tools for Investigating Structural and Physiological Features of Neoplastic Tissue
Sourbron Image analysis in DCE-MRI
da Silva et al. FOD-Swin-Net: angular super resolution of fiber orientation distribution using a transformer-based deep model
Li Acceleration of Subtractive Non-contrast-enhanced Magnetic Resonance Angiography
Roy et al. Intra-bin correction and inter-bin compensation of respiratory motion in free-running five-dimensional whole-heart magnetic resonance imaging
Liu Plug-and-Play Methods for Magnetic Resonance Image Reconstruction

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
C17 Cessation of patent right
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20120118

Termination date: 20140504