CN103261908A - 用于估计动脉/组织/静脉动态系统的感兴趣的量的系统和方法 - Google Patents

Abstract

本发明涉及一种用于通过将软概率方法应用于灌注成像来估计血液动力学参数的系统和方法。这样的方法还可以估计动脉输入函数或互补分布函数，并且因此更普遍的是可以估计任何感兴趣的量。本发明特别地从已知的方法中脱颖而出，这是因为本发明需要引入具有生理或血液动力学性质的严格的软信息的先验，而不会通过任意假设或不期望的假设来约束或强制执行所期望的估计。

Description

用于估计动脉/组织/静脉动态系统的感兴趣的量的系统和方法

本发明涉及一种用于通过将软概率方法应用于灌注加权成像技术来估计血液动力学参数的系统和方法。而且，这样的方法允许估计互补累积密度函数或动脉输入函数，随后且更普遍地，允许估计任何感兴趣的量。本发明不同于已知的方法的地方尤其在于：本发明需要引入生理或血液动力学的严格的软信息的先验而无需通过任意的假设和不期望的假设来约束或强制执行所需要的估计。

本发明特别依赖于灌注加权核磁共振成像(PW-MRI)或灌注计算机断层扫描(CT)。这些技术允许快速获得有关诸如大脑或心脏的器官的血液动力学的有用的信息。此信息对于帮助执业医师在诸如脑急中风的紧急病理治疗时做出诊断和治疗决策是特别关键的。

这些技术依赖于核磁共振或计算机断层扫描设备。此设备发送诸如脑的身体部位多个数字图像序列。为此，所述设备在身体的所述部位上施加高频电磁波的组合，然后测量由某些原子重新发射的信号。这样，该设备允许确定成像体积的每一个点（或体元）的化学组成且随后确定生物组织的种类。

通过专门的处理单元来分析成像序列。该处理单元最终通过合适的人机界面将由灌注加权图像估计出的血液动力学参数传送给执业医师。这样，执业医师可以做出诊断并决定哪一种治疗决策是合适的。

核磁共振或计算机断层扫描灌注加权图像是通过静脉注入造影剂(例如，用于核磁共振成像的钆螯合剂)并通过记录其在图像的每一个体元内随时间变化的团而获得的。为了简洁的目的，我们将省去识别体元的索引x,y,z。例如，与将坐标x,y,z的体元的信号表示为S_x,y,z(t)不同的是，，我们将简单地将其表示为S(t)。应理解的是，下文描述的操作和计算通常针对所感兴趣的每一个体元进行，以便最终获得表示待估计的血液动力学参数的图像或地图。

标准模型允许将所测得的随时间t变化的信号强度S(t)与所述造影剂的浓度C(t)联系起来。

例如，在灌注计算机断层扫描中，每一个体元的信号与浓度成正比：S(t)=k·C(t)+S₀。在通过核磁共振进行的灌注成像中，存在指数关系S(t)=S₀·e^-k·TE·C(t)。在这两种情形中，S₀表示造影剂到达之前的平均信号强度。关于核磁共振成像，k是常数，取决于顺磁磁化率与组织内的造影剂浓度之间的关系，且TE是回波时间。每一个体元的常数k的值是未知的，对于所有感兴趣的体元，其被设定成任意值。因而，得到的是相对的估计值而不是绝对值。尽管这样，该相对信息仍然是相关的，因为主要感兴趣的这些值在空间内的相对变化，尤其在正常组织与病变组织之间的相对变化。

通常来说，我们用S(t)=Ψ(C(t),Θ_S)表示将理论信号S(t)与造影剂的理论浓度C(t)关联在一起的模型，Θ_S是所述模型的自由参数的矢量。例如，对通过核磁共振进行的灌注加权成像或计算机断层扫描，我们有Θ_S=(S₀)。

封闭在每一个体元内的组织的体积内的造影剂的每次质量守恒写为

C_a(t)是供给该组织体积的动脉内的造影剂浓度，已知为动脉输入函数或AIF。BF是组织体积内的血流，且C_v(t)是排空该组织体积的静脉内的造影剂浓度，已知为输出函数或VOF。

假设动脉/组织/静脉动态系统是线性非时变的，我们有

其中h(t)是系统脉冲响应-或组织内的造影剂的通过时间的概率密度函数-且

表示卷积。那么，具有初始条件C(t=0)=0的前述微分方程的形式解可写为：

其中R(t)是互补累积密度函数或残余函数，定义为

其中H是单位阶跃函数。脉冲响应和互补累积密度函数限定了另一个血液动力学参数：组织或MTT内的平均通过时间：

$\mathrm{MTT}=\underset{0}{\overset{+\∞}{\∫}}t.h\left(t\right)\mathrm{dt}=\underset{0}{\overset{+\∞}{\∫}}R\left(t\right)\mathrm{dt}$ （如果 $\underset{t\→+\∞}{\lim }t\·h\left(t\right)=0$ ）

还可以由关系式BV=BF·MTT限定血流BV。

在通过核磁共振进行的灌注加权成像中，诸如BF、MTT或BV的血液动力学参数以及互补累积密度函数一般以如下方式估计。

对每一个体元来说，通过利用关系式

在时间点t_i,i=1,N取样获得的实验灌注信号S_exp(t)被转换成浓度-时间曲线C(t)。对所有体元来说，常数k被设定成非零的任意值（如k·TE=1）。通过例如采用造影剂到达之前常数S₀的平均值来对该常数进行估计。请注意，只有当与造影剂到达时间或BAT相比灌注信号采集开始得足够早时，这才是可能的。由浓度C(t)且假设相关的理论动脉输入函数C_a(t)是已知的，通过数值去卷积来估计乘积BF·R(t)。

为了获得用于对浓度-时间曲线C(t)进行去卷积的理论动脉输入函数C_a(t)，已经提出了好几种方法。

在第一种方法中，执业医师手动选择全身实验动脉输入函数。该函数可以在例如用于脑灌注成像的大脑侧动脉或内颈动脉中测得，或可以从额外的测量中获得，例如光学测量。虽然其允许获得具有高信噪比的信号，该方法还是存在许多缺陷。首先，该方法需要人为干预和/或额外的测量。这在临床紧急情况中是不期望的且这使得过程和最终的结果更难以复制。第二且最重要的是，该全身动脉输入函数并不与每一个体元的局部动脉输入函数相匹配。就延迟（因为与在脉管系统的上游采集的全身动脉输入函数相比，局部动脉输入函数通常是滞后的）和分散（因为造影剂在脉管系统的下游的扩散比在上游的慢）来说，全身动脉输入函数与局部动脉输入函数是不同的。目前已知的是这些现象最终会对血液动力学参数估计产生相当大的影响，原因是由于卷积的对称性，这些缺点直接影响互补累积密度函数的估计。因而，例如并不能最终获得局部动脉输入函数与局部静脉输出函数之间真实的平均通过时间（MTI）的估计，而仅仅是获得全身动脉输入函数与静脉输出函数之间的平均通过时间。为了克服这些不一致，一些作者引入了新的描述性参数，诸如

以量化全身动脉输入函数与局部动脉输入函数之间的延迟，即使它们不属于最初的标准灌注模型（在标准模型中动脉输入函数是每一个体元内的真实局部动脉输入函数）。其他方法趋于使局部动脉输入函数的这些不一致对血液动力学参数的估计的影响最小。然而，这些方法在全身问题中引入了新的未知量且仅仅是逃避这一问题。

根据第二种方法，全身动脉输入函数是经由信号处理技术（诸如数据聚类或独立成分分析（ICA））而从灌注加权图像中自动获得的。虽然该方法允许避免人为干预，但其并不能解决全身动脉输入函数所固有的延迟和分散问题并引入了新的未知量（如，其可能获得静脉输出函数而不是动脉输入函数）。

根据第三种方法，局部动脉输入函数是通过信号处理技术和选择标准而从灌注加权图像自动获得的。例如，寻找待估计其血液动力学参数或互补累积密度函数的当前组织体元的紧邻域中的“最佳”函数。第三种方法的目的是通过至少在一定程度上克服延迟和分散问题而最终获得偏差较小且更准确的估计。然而，先验和后验均不能保证这样获得的局部动脉输入函数是感兴趣的体元的“真实的局部函数”的相关近似值。例如，该“真实”函数可能位于讨论中的邻域之外（如果太小的话），或正相反，该“真实”函数会与另一个动脉输入函数混淆（如果太大的话）。此外，该“最佳”局部动脉输入函数是“正常”动脉输入函数（即，稍早/极早的造影剂到达时间，具有大的幅度等）。但是，目的是准确地区分正常动脉输入函数与病理动脉输入函数（例如缺血性动脉输入函数）。因此，即使最后的结果可能优于使用全身方法获得的结果，但局部动脉输入函数具有不确定性，更不必说血液动力学参数或互补累积密度函数在很大程度上仍不确定。

为了通过由上述方法获得的理论动脉输入函数C_a(t)对实验浓度-时间曲线C(t)进行去卷积，首先随时间对标准卷积模型

进行离散化，例如根据矩形方法的近似值：

$\∀i=1,N,C\left(t_i\right)=\mathrm{BF}\·\underset{0}{\overset{t_i}{\∫}}{C}_a\left(\mathrm{\τ}\right)\·R(t-\mathrm{\τ})\mathrm{d\τ}\≈\mathrm{BF}\·\mathrm{\Δt}\·\underset{k=0}{\overset{i}{\mathrm{\Σ}}}{C}_a\left(t_i\right)\·R(t_i-t_k)$

其中，Δt是取样时间段。这样，如果我们使

$b=|\begin{array}{c}R\left(t_1\right)\\ R\left(t_2\right)\\ ...\\ R\left(t_N\right)\end{array}$

$c=|\begin{array}{c}C\left(t_1\right)\\ C\left(t_2\right)\\ ...\\ C\left(t_N\right)\end{array}d=\mathrm{BF}.b$

那么，我们可获得线性系统Ad=c。

实际中，矩阵A是严重病态的且几乎是奇异矩阵，这使得不获得无意义的解和偏离的估计就不能获得该线性系统的数值解。因此，必须求助于不同的方法，以获得矩阵A的（例如）伪逆矩阵

并且随后，通过获得d的估计

这些用于获得伪逆的方法包括基于A的奇异值截断的方法（截断奇异值分解法或(T)SVD），例如，sSVD法（简单奇异值分解法）、cSVD法（循环奇异值分解法）以及oSVD法（振荡指数奇异值分解法）或频域内的Hunt去卷积。

更一般地，可以使诸如||Ad-c||²+||Γd||²的标准最小化，其中||Γd||²是正则项，其支持特定解并允许通过

获得d的估计。这些方法包括Tikhonov正则化和小波变换法等。

一旦得到

那么可以通过

获得BF的估计

这是因为，通过定义，

然而，例如在基于奇异值分解的方法中，BF通常被估计为

以用于补偿这些方法所固有的对的系统性低估。接着，一旦通过获得了b的估计

且通过例如下述矩形近似法由

获得了

的估计。最后，通常通过

获得BV的估计，即使乘积的估计并不是估计的乘积。

可以发现基于这些动脉输入函数的许多变化，例如，实验动脉输入函数可以先被拟合到参数或半参数理论模型C_a(t,Θ_a)，其中Θ_a是参数的矢量，以便人为增大信噪比。信号也可以被人为过采样，以便使数值去卷积更稳定或克服在造影剂循环信号（第一次通过）与再循环信号（第二次通过）之间存在时间重叠时可能因这种再循环产生的问题。但是这些变化仍基于数值去卷积方法，例如截断奇异值分解。

根据某些比较性研究，在已经就打算模拟典型实时数据的合成数据进行过测试的不同方法中，就偏差（相对于真实值的系统误差）、精确度（估计相对于真实值的标准偏差随输入信号的信噪比的变化）以及相对于随患者、组织种类、病理等的不同而实际出现的各种互补累积密度函数R(t)和动脉输入函数C_a(t)的稳固性而言，基于奇异值分解的去卷积及其变化形式（线性截断奇异值分解（sSVD）、截断循环（cSVD）或平滑截断循环（oSVD））最终得到了BF和MTI的最佳估计。

然而，在与动脉输入函数的选择有关的上述问题之上，此类数值方法存在明显的固有问题。

首先，BF.b的估计并不随时间减小，而是以有时呈负值的程度发生振荡。但R(t)是时间t时保留在体元内的造影剂的量，因此其必然是正的下降函数。Ad hoc法（例如oSVD）允许降低这些异常振荡，但异常振荡仍然存在，因为异常振荡是奇异值分解所固有的。这就是为什么在此类方法中血流BF是通过估计的，而血流BF应该被估计为

因为在标准灌注模型内，

通过取最大值,而不是原始值，只能寄希望于消除这些振荡对BF估计的影响。因此，这些方法不是完全令人满意的。具体地，通过更严格的互补累积密度函数的估计方法所能获得的BF估计的准确性以及随之而来的血液动力学参数的准确性仍然是未知的。因此，这些数值去卷积方法与标准灌注模型相矛盾，因为它们提供了不能满足所述模型的特性的解。

为了摆脱此问题并获得互补累积密度函数的生理上可接受的估计，已经引入了这些互补累积密度函数的参数模型R(t,Θ_R)，Θ_R是所述模型的参数的矢量。这些模型例如通过采用贝叶斯法被拟合到实验信号上。然而，该方法在这一点上似乎不成熟。的确，应当自行选择互补累积密度函数的先验非参数估计，以便确定适于描述它们的参数或半参数模型，因为Monte-Carlo模拟已经显示出：如果这些模型不完全适于准确地描述实际中可能出现的所有种类的互补累积密度函数，那么所得到的血液动力学参数的估计（诸如MTI或BF）会变得异常。因而，选择互补累积密度函数的理论模型是至关重要的且可以仅仅通过将模型应用于实验数据而恰当地做出选择。因此，有必要使用非参数法来估计所述互补累积密度函数，以便在下一步骤中用允许得到互补累积密度函数的生理上可接受的估计的典型的参数或半参数方法替代它们。

但是，如前所述，通过诸如基于奇异值分解的方法获得的估计与根据标准灌注模型的非常定义

相矛盾。这些估计不是生理上和身体上可接受的。因此，无法将参数或半参数理论模型拟合到这些估计上。更不必说无法将若干模型彼此比较且不可以选择最适合于描述实验累积密度函数的那些模型。因此，无法再在建模和对灌注现象的理解方面取得进展，因为诸如奇异值分解方法的数值方法存在固有的缺点。

此外，互补累积密度函数和血液动力学参数的非参数估计中存在的问题看起来不是由经验性动脉输入函数获得的经验性浓度-时间曲线的“简单”去卷积问题。的确，如果实际上在此问题内给定了动脉输入函数且已知具有绝对的确定性和极大的准确性，那么可能发生这种情况。但至多仅会提供已知仅取决于测量噪声且必须被预先转换成浓度-时间曲线的经验性动脉测量信号。换句话说，假设经验性动脉输入函数等于理论动脉输入函数，则忽略实验信号上的测量噪声和由这些信号的浓度-时间曲线的估计或转化产生的不确定性。

标准卷积模型

仅包括通常不能直接测得的理论信号。实际上，所测得的信号通常是理论信号和测量噪声的总和。因此，实验性动脉和组织灌注信号分别写为S_exp(t)=S_th(t)+ξ_t和其中ξ_t和

是对测量噪声建模的零-平均随机过程。此外，我们有和

其中C_th(t)是理论浓度-时间曲线且C_ath(t)是核磁共振灌注加权成像内的理论动脉输入函数，或更通常是S_th(t)=Ψ(C_th(t),Θ_S)，如前所述。因此，应用于实验信号的理论标准灌注模型必须被写为

即，在核磁共振灌注加权成像内为：

$\ln \frac{S_{\exp }\left(t\right)-{\mathrm{\ξ}}_t}{S_0}=\mathrm{BF}\·\ln \frac{S_{\mathrm{aexp}}\left(t\right)-{\mathrm{\ξ}}_t^a}{S_{0a}}\⊗R\left(t\right)$

或等同地为：

$\ln \frac{S_{\exp }\left(t\right)+{\mathrm{\ξ}}_t}{S_0}=\mathrm{BF}\·\ln \frac{S_{\mathrm{aexp}}\left(t\right)+{\mathrm{\ξ}}_t^a}{S_{0a}}\⊗R\left(t\right)$

而不是：

$C_{\exp }\left(t\right)=\mathrm{CBF}\·\left[C_{a_{\exp }}\left(t\right)\right]\⊗R\left(t\right)+{\mathrm{\ξ}}_t$

其中C_exp(t)是实验浓度-时间曲线，且

是理论动脉输入函数，该理论动脉输入函数是有误差的隐含的数学卷积模型，大多数去卷积模型都基于该数学卷积模型。由这一点可以看到，优选将标准灌注模型写为：

$\left\{\begin{array}{c}{S}_{\exp }\left(t\right)=S_0{e}^{-C_{\mathrm{th}}\left(t\right)}+{\mathrm{\ξ}}_t\\ S_{\mathrm{aexp}}\left(t\right)={S_{0a}e}^{-C_{\mathrm{ath}}\left(t\right)}+{\mathrm{\ξ}}_t^a\\ C_{\mathrm{th}}\left(t\right)=\mathrm{BF}.C_{\mathrm{ath}}\left(t\right)\⊗R\left(t\right)\end{array}\right.$

以避免对非必要的正随机变量求对数。

已知待去卷积的信号的测量的不确定性对去卷积过程的最终结果产生了相当大的影响：由于这些不确定性，输入信号的极小量变化可能会使最后的结果产生相当大的变化。为了准确地克服这些问题且减少这些不稳定性，已经引入了去卷积法，例如Tikhonov正则法或基于奇异值分解的方法。何况测量噪声和不确定性对动脉输入函数的影响（目前完全被忽略）甚至更可观：此时，在卷积矩阵A中出现了动脉测量噪声

和对S_0a的不确定性，且其随后被放大并传输。忽略这些误差和对动脉输入函数的不确定性会造成与血液动力学参数估计相关的严重误差以及对这些估计的准确性的错觉。某些方法的目的是消除动脉输入函数上的测量噪声以使到现在都仍在逃避的该问题最小化。优选的是一种自行允许传输与动脉输入函数相关、与血液动力学参数的估计和互补累积密度函数相关的不确定性以便控制并量化估计误差的方法。

此外，还优选避免将灌注模型写为

的确，标准灌注模型首先界定了动脉/组织/静脉动态系统的脉冲响应h(t)，互补累积密度函数R(t)由此通过

计算。因此，方便的是将标准灌注模型写为脉冲响应h(t)的函数，即为以将此模型拟合，从而由

估计测量时间点t_j,j=1,N的脉冲响应h(t)，以便最终获得互补累积密度函数的估计

如

（通过矩形法获得的近似值）。具体地，最好不通过互补累积密度函数的估计

来估计脉冲响应h(t)。从数值点的观点来讲，容易理解的是优选先估计导数（h(t)），以便在第二步骤中估计反导数（R(t)），而不是相反的顺序。尽管这样，仍可核查直接估出的是互补累积密度函数R(t)而不是脉冲响应h(t)。

此外，去卷积问题是不适定问题且容许无穷多个解。何况在用于估计血液动力学参数、脉冲响应、互补累积密度函数或动脉输入函数的灌注加权成像内面临的不适定问题是与测量噪声传输的问题和实验信号转换成浓度-时间曲线的问题三重叠加在一起的去卷积问题。

为了将其归结为可能容许唯一解的适定问题，必须增加在先信息以及在先前可能的无穷解中寻求的对解的约束。这就是为什么可以在现有技术中找到许多去卷积和估计方法，每一种方法多少都明确地或直接地引入特定种类的在先信息。

举个例子，可以提及先前所述的典型的Tikhonov正则法。在其最通用的版本中，矩阵Γ等于Γ=αI_N。这消除了此约束条件的较大欧几里得范数的解和用于量化此罚则的加权的标量α。

另一方面，基于截断奇异值分解的方法是灌注加权成像内的典型方法，这些方法在于删除低于充当正则化参数的给定阈值的、卷积矩阵的奇异值。那些小的奇异值与信号的高频分量有关，使得奇异值的截断点用作低通滤波器。但一方面，由于在t=0处的中断（互补累积密度函数被界定在]-∞,+∞[内），互补累积密度函数是无限带宽信号。另一方面，小的奇异值的截断并不是直接的且仅对应于额外的生理信息，但最重要的是对应于额外的代数信息，即，给定的线性子空间的成员。核查是否保持着过高的频率振荡。具体地，通过互补累积密度函数与动脉输入函数的去卷积而估计的理论信号具有与测量噪声相符的趋势（过拟合）。此外，截断阈值并未达到明确的生理约束条件的加权（而仅仅是迫使该解归属于一些线性子空间），其作用是复杂的而不是直接的，且难以给出一般的标准来确定其最合适的值，这取决于不同种类的信号。相反，并不能保证可以合适地优化这样的标准（例如，oSVD法内的互补累积密度函数的粗略度），因为奇异值的截断对该标准所起的作用是非常间接的。

此外，现有的血液动力学参数或互补累积密度函数估计方法并不允许估计那些估计的精确度（即，估计量的标准偏差），也不允许估计可能存在于那些估计和存在于那些估计的精确度的估计中的置信度。具体地，难以通过去卷积和估计方法来量化标准灌注模型的拟合优度。因此，应注意到，基于截断奇异值分解的方法具有过拟合实验信号的趋势。所估计的信号并不是平滑的，而是正相反，具有与测量噪声相符的趋势。因此，可以获得用于测量模型拟合优度的低值，例如典型的χ²统计（即，平方误差的和）。这种低的值表明高质量的拟合，而不是过拟合所造成的情形。由此，如果想要检测并考虑过拟合的话，这些统计在此情形中是不合适的且具误导性，应该被避免。因此，引入拟合优度的“良好”测量以便用于根据现有技术的去卷积和血液动力学参数估计方法并非显而易见。因此，期望的是一种自行的方法，对于该方法，这样的拟合优度测量是唯一界定的，且该方法将过拟合考虑在内并仅仅量化标准灌注模型对数据的拟合优度。

总之，当前用于估计血液动力学参数、脉冲响应或互补累积密度函数的方法存在许多不受控制的方法学错误和近似值。因此，对结果的解读很困难，且灌注加权成像并未被充分利用。

本发明的目的是提供对采用已知方法所产生的所有缺点的解决方案。本发明的主要目的是提供一种的新的方法，其允许在满足生理或血液动力学种类的某些约束条件的解中寻求脉冲响应h的估计

而不需引入并不一定需要被满足且最重要的是不可能通过实验证实或属于非生理或非血液动力学种类的约束条件ad hoc。此外，这些方法允许以自动的且唯一的方式确定这样的约束条件的加权，而不需要求助于ad hoc法。

因而，本发明可以宣称最终解决了血液动力学参数、互补累积密度函数或动脉输入函数的问题。本发明提供了其独特的解决方法（可能多种）。

在本发明提供的主要优势中，非限定性的，我们可以引述的可能的优势是：

-以量化的方式转化定性或半定性的关于血液动力学参数、脉冲响应、互补累积密度函数或动脉输入函数的软生理信息；

-明确地将关于动脉输入函数以及实验信号的不确定性和误差考虑在内，并在感兴趣的量的估计上传播这些不确定性和误差；

-就相对于实际中可能出现的不同情况的偏差（系统误差）、精确度（统计误差）、线性度和稳固度而言，改进那些参数、脉冲响应、互补累积密度函数或动脉输入函数的估计；

-获得关于血液动力学参数估计、脉冲响应、互补累积密度函数或动脉输入函数的置信区间-且甚至获得所述置信区间的可能性，以便改进且澄清那些估计的置信度；

-获得可被生理和血液动力学观点接受且符合标准灌注模型的脉冲响应、互补累积密度函数或动脉输入函数的非参数估计，以便允许在下一步骤中拟合、比较且最终选择用于所述脉冲响应、所述互补累积密度函数或所述动脉输入函数的参数或半参数理论模型；

-量化到何种程度，动脉输入函数实际上是可接受的用于感兴趣的每一个体元内的组织的动脉输入函数；

-提供全身灌注模型的拟合优度的客观且量化的测量，这允许最终比较并选择最合适的全身灌注模型。

为了此目的，本发明具体体现为一种用于估计器官的单元体积-称为体元-内的多个动脉/组织/静脉动态系统中的感兴趣的量的方法。这样的方法的目的是通过灌注加权成像分析系统的处理单元来实施且包括从实验灌注数据估计所述感兴趣的量的步骤。为了在实现某些生理或血液动力学约束条件的解中寻求所述感兴趣的量而不引入不一定需要被满足且最重要的是不可能从实验中得到证实或是非生理或非血液动力学种类的约束条件ad hoc，根据本发明的用于估计的步骤包括：根据贝叶斯方法，通过以下过程评估所述感兴趣的量的后验边缘概率分布：

将直接概率分布分配给实验灌注数据，并且给定估计所述体元内的所述动脉/组织/静脉动态系统的所述感兴趣的量的过程中涉及的参数；

通过引入与所述动态系统的脉冲响应相关的严格的软信息，将联合先验概率分布分配给所述量。

本发明还具体体现为一种用于估计器官的单元体积-称为体元-内的多个动脉/组织/静脉动态系统中的感兴趣的量的方法，所述动态系统是线性时不变的且通过关系式

确定，其中C(t)是在体元内循环的造影剂的浓度，C_a(t)是供给所述体元的动脉内的造影剂的浓度，BF是所述体元内的血流，

代表卷积，且R(t)是所述体元内的通过时间的互补累积密度函数。这样的方法目的在于通过灌注加权成像分析系统的处理单元来实施且包括从实验灌注数据估计所述感兴趣的量的步骤。根据本发明且如前所述，所述估计步骤包括：根据贝叶斯方法，通过以下过程评估所述感兴趣的量的后验边缘概率分布：

将直接概率分布分配给实验灌注数据，并且给定估计当前考虑的体元内的所述动脉/组织/静脉动态系统的所述感兴趣的量的过程中涉及的参数；

通过引入与所述体元内的互补累积密度函数R(t)相关的严格的软信息，将联合先验概率分布分配（53）给所述量。

根据优选的实施方案，在这两种情况下，本发明还计划通过引入将与供给所述体元的所述动脉内的所述造影剂的浓度-时间曲线有关的严格的软信息来将所述联合先验概率分布分配给所述量。

有利地，根据本发明的方法通过逐次迭代正在考虑的多个体元来实施。

为了估计感兴趣的量的估计的精确度，根据本发明的方法包括用于计算由与所估计的感兴趣的量有关的置信区间表示的补充信息的步骤。

根据优选的实施方案，这样的方法可以包括用于计算由与所估计的感兴趣的量有关的置信区间的可能性（betting odd）表示的补充信息的步骤。

还计划，根据本发明的方法可以包括用于计算由以下的乘积的充分性量度表示的补充信息的步骤：

将直接概率分布分配给实验灌注数据，并且给定估计正在考虑中的所述体元内的所述动脉/组织/静脉动态系统的所述感兴趣的量的问题中涉及的参数；

将联合先验概率分布分配给所述量。

为了提供任何估计，根据本发明的方法可以包括用于将所估计的感兴趣的量传送至能够将其呈现给用户的人机界面的步骤。

如果计算出补充信息，那么这样的方法还可以包括用于将与所估计的感兴趣的量有关的任何补充信息传送至能够将其呈现给用户的人机界面的步骤。

根据优选的实施方式，本发明提出：所述实验灌注数据包含在实验灌注信号的值的矢量中或包含在从所述矢量到浓度-时间曲线的值的矢量的转换中。

根据本发明的第二个目的，本发明还涉及处理单元，其包括用于存储的装置、用于与外界通信的装置以及用于处理的装置。用于通信的装置能够从外界接收实验灌注数据。用于处理的装置适于实施根据本发明的用于估计器官的单元体积-称为体元-内的多个动脉/组织/静脉动态系统中的感兴趣的量的方法。

为了向用户提供由根据本发明的方法产生的估计，用于通信的装置以适合于人机界面的格式传送所估计的感兴趣的量，所述人机界面能够将所述感兴趣的量呈现给用户。

根据此实施模式，用于通信的装置以适合于人机界面的格式传送与所估计的感兴趣的量有关的补充信息，所述人机界面能够将所述补充信息呈现给用户。

根据第三个目的，本发明还涉及一种灌注加权成像系统，其包括根据本发明的处理单元和能够根据本发明且通过所述处理单元实施的方法将所估计的量呈现给用户的人机界面。

在阅读下面的描述并浏览附图时，本发明的其他特征和优势将更清楚地显现，附图包括：

-图1和2显示了用于分析灌注加权图像的系统的两个实施方式；

-图3和4分别显示了在注入造影剂之前和在所述造影剂在所述脑组织内循环期间、通过核磁共振成像装置获得的人脑切片的灌注图像；

-图5a和5b显示了通过核磁共振获得的与人脑的体元有关的灌注加权信号S(t)；

-图6显示了在人脑的体元内循环的造影剂的典型的浓度-时间曲线；

-图7显示了典型的动脉输入函数C_a(t)；

-图8显示了根据本发明的方法；

-图9显示了根据本发明估计得到的脑血容量图；

-图10显示了在脑缺血的情形中，根据本发明估计得到的脑血流图；

-图11显示了根据本发明估计得到的平均通过时间图；

-图12显示了脑血流属于置信区间的概率图。

图1阐释了一种用于分析灌注加权图像的系统。核磁共振成像或计算机断层扫描装置1受控制台2的控制。用户可以选择用于控制该设备1的参数11。从由装置1产生的信息10获得了人体或动物体的一部分的多个数字图像序列12。作为优选实施例，我们将阐释现有技术的解决方法和使用由人脑观察到的数字图像的本发明。也可以考虑其他器官。

图像序列12可以被随意地存储在服务器3中且构成了患者的医疗记录13。这样的记录13可以包括各种类型的图像，例如灌注加权或扩散加权图像。使用专门的处理单元4分析图像序列12。所述处理单元包括用于与外界通信的装置以收集图像。所述用于通信的装置还允许处理单元最后通过专门的人机界面5将从灌注加权图像12获得的对血液动力学参数14的估计提供给执业医师或研究者6。分析系统用户6可以证实或反驳诊断、决定他认为合适的治疗行为、进行进一步的调研等等。任选地，用户可以通过设置16来配置处理单元4的操作。例如，其可以界定显示阈值或选择他想浏览的所估计的参数。

图2阐释了分析系统的实施方式，其中预处理单元7分析图像序列12以检索每一个体元的灌注数据15。因而，负责估计血液动力学参数14的处理单元不再担当这一功能，且通过用于与外界通信的装置接收到的灌注数据15来执行估计方法。

图3显示了5mm厚的人脑切片的图像12的典型示例。该图像通过核磁共振获得。对每一个切片来说，使用该技术可以获得128×128体元的矩阵，其尺寸是1.5×1.5×5mm。使用双线性插值，可以产生458×458像素的扁平图像，例如图像20。

图4显示了图像20，其类似于图3所示的图像。然而，该图像20是在注入造影剂后获得的。该图像20是灌注加权脑图像的典型示例。所显示出的动脉无疑与图3中描绘的相同图像相反。根据已知的技术，可以选择病变半球对侧的半球中的一个或若干个动脉输入函数21以便估计血液动力学参数。

图5b阐释了通过核磁共振获得的灌注加权信号S(t)的示例，其作为通过图2所示的预处理单元7传送的数据15。因此，灌注加权信号表示在注入造影剂后体元随时间的演变。例如，图5b描绘了在50秒时间段内的这样的信号。纵坐标描绘了任意单位的信号强度。为了获得这样的信号，根据图1的处理单元4（或可选择地，根据图2的预处理单元7）分析通过核磁共振在时间点t₁、t₂、……t_i、……t_n处获得的n个灌注加权图像序列I1、I2、……、Ii、……，正如图5a描绘的。确定给定的体元（例如体元V）的灌注加权信号S(t)，其表示在注入造影剂后体元随时间t的演变。

图6显示了由灌注加权信号得到的浓度-时间曲线，正如图5b描绘的。正如之前描述过的，灌注加权信号与相关的浓度-时间曲线之间存在关系。因此，在通过核磁共振的灌注加权成像中，存在指数关系S(t)=S₀.e^-k.TE.C(t)，其中S₀是造影剂到达之前的信号的平均强度，TE是回波时间，且k是常数，该常数取决于顺磁磁化率与组织内的造影剂浓度之间的关系。

图6允许查看体元内的造影剂浓度随时间的演变。已观察到在体元内造影剂第一次通过时存在高的幅度峰值，然后是与造影剂再循环（第二次通过）现象有关的较低的幅度峰值。

图7阐释了典型的动脉输入函数C_a(t)，其表示动脉体元（例如，图4呈现的体元21）内的造影剂的循环。图7特别显示了：在造影剂第一次通过后，再循环现象非常弱。

图8显示了根据本发明的用于估计器官内的体元的多个动脉/组织/静脉动态系统中的感兴趣的量的方法的实施例。这样的方法可以通过结合图1和2所描述的且相应调节后的灌注加权成像分析系统的处理单元来实施。

根据本发明的方法主要包括用于将一个或多个边缘后验分布分配给待估计的多个感兴趣的量的步骤56，所述感兴趣的量为例如血液动力学参数、取样时间点处的理论脉冲响应值或互补累积密度函数值。该方法还包括用于计算所述估计的步骤57。

为了分配这样的后验边缘分布，需要配置50处理单元。处理单元自身优选地可以通过一个或多个配置设置值进行该配置。该配置还可以产生包括一个或多个边缘后验分布的库的构型，所述库是预先建立的且存储在处理单元的程序存储器中。本发明规定了所述库可以被富集，因为其是被能够通过所述配置设置值执行所述配置且能够与处理单元配合以便输出所述库的外部处理单元使用和传送的。

因而，根据本发明的方法可以包括在分配56之前执行的配置步骤，其中，下面的过程是充分必要条件：

-将直接概率分布分配54给实验数据，并且给定估计正在讨论中的体元内的动脉/组织/静脉动态系统的感兴趣的量的问题所涉及的所有参数；

-将联合先验概率分布分配53给所有那些参数。

本发明允许估计不同应用情形中的一个或多个感兴趣的量：

-假设理论动脉输入函数是已知的且具有绝对确定性和无限准确性；

-局部动脉输入函数以未知的延迟区别于全身动脉输入函数，该延迟也有待估计；

-动脉输入信号仅被测量，可能长达延迟；

-未给出动脉输入函数且未测量动脉输入信号–这是最实际的情形。

配置步骤可以取决于讨论中的应用情形。

图8阐释了根据本发明的、用于动脉输入信号被测量（可能长达延迟τ）的第一应用实施例的方法。

实验灌注模型M可以被写为：

$M:\left\{\begin{array}{c}s=\mathrm{\Ψ}(c,{\mathrm{\Θ}}_S)+\mathrm{\ξ}\\ s_a=\mathrm{\Ψ}(a,{\mathrm{\Θ}}_{S_a})+{\mathrm{\ξ}}_a\\ c=\mathrm{BF}\·\mathrm{Ab}(t-\mathrm{\τ})=\mathrm{BF}\·B(t-\mathrm{\τ})a\end{array}\right.$

b=[R(t₁),...,R(t_N)]^T是未知的互补累积密度函数的值的矢量。正如前面解释过的，此矢量可以由脉冲响应的值的矢量h=[h(t₁)=0,...,h(t_N)]^T来表示。例如，根据左矩形近似法

或根据右矩形近似法 $b={[\mathrm{1,1}-\mathrm{\Δt}\·h\left(t_2\right),...,1-\mathrm{\Δt}\·\underset{i=2}{\overset{N}{\mathrm{\Σ}}}h\left(t_N\right)]}^T;$

c=[C(t₁),...,C(t_N)]^T是体元内的造影剂的未知理论浓度的值的矢量；

s=[S(t₁),...,S(t_N)]^T是由灌注加权成像获得的信号的强度的实验测量值的实数矢量或复数矢量；

Θ_S是使S(t)与C(t)相关的参数的矢量；

ξ=[ξ(t₁),...,ξ(t_N)]^T是所测得的信号的测量噪声的矢量；

s_a=[S_a(t₁),...,S_a(t_N)]^T是实验动脉输入信号S_a(t)的测量结果的矢量；

是使S_a(t)与S_a(t)相关的参数的矢量；

ξ_a=[ξ_a(t₁),...,ξ_a(t_N)]^T是所述动脉信号的测量噪声的矢量；

a=[C_a(t₁),...,C_a(t_N)]^T是未知的理论动脉输入函数C_a(t)的值的矢量；

A是与通过数值逼近卷积积分获得的未知的理论动脉输入函数a的值的矢量相关的卷积矩阵，该卷积积分是通过矩形近似法、梯形法或通过较高阶的方法（诸如Simpson法、Boole法、Gauss-Legendre法等）得到的。A还可以是循环卷积矩阵，例如，基于奇异值的截断、cSVD以及oSVD的去卷积法中使用的那些矩阵。

我们将模型M的血液动力学参的矢量表示为Θ，例如此处Θ=(BF,τ)。

结合图1、2和8，通过处理单元4实施的方法可以包括两个初始配置步骤51d和51e，其分别引入关于实验信号s的测量值的矢量的信息I_S和关于实验动脉信号s_a的测量值的矢量的信息

。

如果我们例如紧跟着实值测量噪声矢量偶的刚开始两个时刻E(ξ,ξ_a)≡(0,0)和

那么最大熵原理（根据Lebesgue参照测量的微分Shannon熵）要求矢量ξ和ξ_a必须被视为是分别具有标准偏差σ_S和的相互独立的、白的、稳定的高斯随机过程。

更普遍地，我们将表征该测量噪声矢量偶(ξ,ξ_a)的参数表示为 $(E_S,E_{S_a}).$ 例如 $(E_S,E_{S_a})=({\mathrm{\σ}}_S,{\mathrm{\σ}}_{S_a}).$

根据本发明的方法包括用于将联合直接概率分布分配给测量矢量偶(s,s_a)的配置步骤54，并且给定矢量a和b、矢量Θ、参数Θ_S和的矢量以及矢量偶

那么所述联合直接概率分布写作 $p(s,s_a|a,b,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},E_S,E_{S_a},I_S,I_{S_a},M).$

例如，我们有

$p(s,s_a|a,b,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},{\mathrm{\σ}}_S,{\mathrm{\σ}}_{S_a},I_S,I_{S_a},M)\∝$

${\left({\mathrm{\σ\σ}}_a\right)}^{-N}\exp \{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\frac{{[S\left(t_i\right)-\mathrm{\Ψ}(\mathrm{BF}.\mathrm{Ab}(t_i-\mathrm{\τ}),{\mathrm{\Θ}}_S)]}^2}{{2\mathrm{\σ}}^2}+\frac{{[S_a\left(t_i\right)-\mathrm{\Ψ}(a\left(t_i\right),{\mathrm{\Θ}}_{S_a})]}^2}{{{2\mathrm{\σ}}_a}^2}\}$

如果将ξ和ξ_a视为是分别具有标准偏差σ_S和

的相互独立的、白的、稳定的高斯随机过程，本发明可选择地允许以不同的方式表达所述联合直接概率分布（如有需要），但不是直接以未知的理论互补累积密度函数的值的矢量b来表达，而是以脉冲响应的值的矢量h来表达。基于此目的，就h而言，应该恰当表达b。那么，联合直接概率分布54可被写作，例如 $p(s,s_a|a,h,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},{\mathrm{\σ}}_S,{\mathrm{\σ}}_{S_a},I_S,I_{S_a},M)\∝$

${\left({\mathrm{\σ\σ}}_a\right)}^{-N}\exp \{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\frac{{[S\left(t_i\right)-\mathrm{\Ψ}(\mathrm{BF}.\mathrm{Ab}(t_i-\mathrm{\τ}),{\mathrm{\Θ}}_S)]}^2}{{2\mathrm{\σ}}^2}+\frac{{[S_a\left(t_i\right)-\mathrm{\Ψ}(a\left(t_i\right),{\mathrm{\Θ}}_{S_a})]}^2}{{{2\mathrm{\σ}}_a}^2}\}$

在不损失通用性的前提下，我们随后应考虑到，基于前面给出的原因，将优先对矢量h感兴趣。

这样，可以进一步将复信号的测量值的矢量偶的实部和虚部的直接概率分布相乘，由此将联合直接概率分布分配给复信号的测量值的矢量偶。

本发明还提供了这样一种变型：灌注信号s和s_a的值的矢量可以被准确转化成浓度的值c_exp的矢量，例如在使用核磁共振的灌注加权成像中，

当例如可以独立于其他参数准确地估计造影剂到达之前的平均灌注信号强度S₀时，这么做是可能的。根据此实施方案，可以显示出，在给定例如如下所示的所有其他参数的情况下，将高斯概率分布分配给矢量偶(c_exp,c_expa)仍是恰当的：

$p(c_{\exp },c_{\mathrm{expa}}|a,h,\mathrm{\Θ},\mathrm{\σ},{\mathrm{\σ}}_a,I_S,I_{S_a},M)\∝$

${\left({\mathrm{\σ\σ}}_a\right)}^{-N}\exp \{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\frac{{[c_{\exp }\left(t_i\right)-\mathrm{BF}.\mathrm{Ab}\left(t_i\right)]}^2}{{2\mathrm{\σ}}^2}+\frac{{[c_{\mathrm{expa}}\left(t_i\right)-a\left(t_i\right)]}^2}{{{2\mathrm{\σ}}_a}^2}\}$

σ和σ_a此时分别是关于c_exp和c_expa的测量噪声的标准偏差。

更普遍地，我们将使用术语“实验灌注数据”来表示实验灌注信号s的值的矢量及其到浓度-时间曲线的值的矢量c_exp的转换。此后，在不损失通用性的情况下，我们将s和s_a表示为实验灌注数据。

此外，该方法包括三个配置步骤51a、51b以及51c，其分别引入关于模型M的血液动力学参数Θ的信息I_Θ的片段、关于脉冲响应h的信息I_h的片段以及关于动脉输入函数a的信息I_a的片段。根据本发明，在步骤51a至51e中引入的信息片段构成了用于配置（即:使处理单元能够分配56先验边缘分布且随后估计57感兴趣的量）的配置参数。从这些信息片段-或配置参数-，根据本发明的方法可以包括用于分配联合先验概率分布的步骤53，该联合先验概率分布可写作：

$p(a,h,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},E_S,E_{S_a}|I_{\mathrm{\Θ}},I_h,I_a,I_S,I_{S_a},M)$

应用最大熵原理，此分布通常可以被因式分解为：

$p(a,h,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},E_S,E_{S_a}|I_{\mathrm{\Θ}},I_h,I_a,I_S,I_{S_a},M)$

$=p\left(\mathrm{\Θ}\right|I_{\mathrm{\Θ}},M)\·p(h,E_S|I_h,I_S,M)\·p(a,E_{S_a}|I_a,I_{S_a},M)\·p({\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a}|I_S,I_{S_a},M)$

因此该方法包括步骤52a，其分配先验概率分布p(Θ|I_Θ,M)。

例如，可以分配非信息型的先验分布。例如，如果信息I_Θ仅能告知BF和τ分别属于区间[BF_min,BF_max]和[τ_min,τ_max]，那么先验概率分布p(Θ|I_Θ,M)可以被表达为：

$p\left(\mathrm{\Θ}\right|I_{\mathrm{\Θ}},M)=p(\mathrm{BF},\mathrm{\τ}|I_{\mathrm{\Θ}},M)=$

${\mathrm{\χ}}_{[{\mathrm{BF}}_{\min },{\mathrm{BF}}_{\max }]}\left(\mathrm{BF}\right)\·{\mathrm{\χ}}_{[{\mathrm{\τ}}_{\min },{\mathrm{\τ}}_{\max }]}\left(\mathrm{\τ}\right)/[\log \left({\mathrm{BF}}_{\max }\right)-\log \left({\mathrm{BF}}_{\min }\right)]/({\mathrm{\τ}}_{\max }-{\mathrm{\τ}}_{\min })/\mathrm{BF}$

另一方面，还可以分配信息型概率分布，例如，在过去的实验中通过例如量化灌注加权成像技术（正电子成像术（PET）或动脉自旋标记（ASL））获得的相对频率取样分布或边缘后验概率分布。

根据本发明的方法还包括配置步骤52b，其从信息I_h的片段（或I_b，如果仅仅对互补累积密度函数感兴趣）分配先验概率分布p(h,Ε_S|I_h,I_S,M)=p(h|Ε_S,I_h,M)·p(Ε_S|I_S,M)。此信息片段由硬信息片段和软信息片段构成。

根据本发明，硬信息区别于软信息。硬信息片段对应于任何被认为是确定性的布尔命题-即，其概率等于1。例如，诸如“该曲线是平滑的”或“该信号符合该模型”的布尔命题构成了硬信息片段。相比之下，软信息片段关注于任何指示且仅指示特定概率的布尔命题。最后，这相当于引入诸如“该曲线大致是平滑的”或“该信号大致符合该模型”的布尔命题。此后，我们要求“严格的软信息”，非硬信息的任何软信息。

为了强调硬信息与软信息之间的差异，我们考虑了最简单的且最极端的情形。假设想要从根据高斯分布N(x,σ²)的独立且恒等分布的样本x₁,...,x_n估计实量x。确定一实数a。考虑先前的硬定量信息片段”x=a”。该硬信息转换成Dirac先验概率分布p(x|a)=δ(x-a)，即，如果x=a，那么其是p(x|a)=1；且如果x≠a，其是p(x|a)=0。

现在，我们考虑对应的软定量信息“x大致接近a”，即，“||x-a||²＝(x-a)²大致是较大的”，其中||||代表欧几里得范数。

应用最大熵原理，该软信息片段转换成高斯先验概率分布：

$p\left(x\right|a,\mathrm{\ϵ})=\frac{1}{\sqrt{2\mathrm{\π}}\mathrm{\ϵ}}\exp [-\frac{{(x-a)}^2}{{2\mathrm{\ϵ}}^2}]$ 且ε>0

可以看出，硬信息片段和对应的软信息片段转换成不同的先验概率分布。具体地，在软信息片段的先验分布中，新的（超级）参数ε似乎允许以量化的方式、通过分配先验概率分布p(ε)来转换先验的“大致”。特别地，我们发现硬信息片段是软信息片段的极限情况（使ε趋向0）。然而，本发明仅涉及严格的软信息，对此，ε>0。

硬先验概率分布和对应的严格的软先验概率分布总是不同的，他们通常产生本身不同的、感兴趣的量的估计。在我们的极端实施例中，实际上不需要在硬的情形中估计x，因为我们已经知道x=a。另一方面，在软的情形中，其实质上是估计高斯分布的数学预期值的典型问题。

我们现在考虑作为第二实施例的典型的硬定量信息片段“y符合（半）参数模型M(x,Θ)”，Θ是M的参数的矢量，且y是实数。通过定义，该定量的信息片段转换成Dirac先验概率分布：

p(y|x,Θ,M)=δ[y-M(x,Θ)]

现在，如果z是噪声实验数据，诸如z～N(y,σ²)，那么我们得到似然函数或直接概率分布

$p\left(z\right|x,\mathrm{\Θ},\mathrm{\σ},M)=-\frac{1}{\sqrt{2\mathrm{\π}}\mathrm{\σ}}\exp \{-\frac{{[z-M(x,\mathrm{\Θ})]}^2}{2{\mathrm{\σ}}^2}\}$

由此，可以通过应用例如贝叶斯规则来估计感兴趣的量Θ和σ：

p(Θ,σ|z,x,M)∝p(Θ,σ|M)p(z|x,Θ,σ,M)

现在考虑对应的软定量信息片段“y大致符合（半）参数模型M(x,Θ)”。该软信息转换成：

“||y-M(x,Θ)||²=[y-M(x,Θ)]²大致是较大的”。

应用最大熵原理，该信息转换成先验概率分布：

$p\left(y\right|x,\mathrm{\ϵ},\mathrm{\Θ},M)=\frac{1}{\sqrt{2\mathrm{\π}}\mathrm{\ϵ}}\exp \{-\frac{{[y-M(x,\mathrm{\Θ})]}^2}{2{\mathrm{\ϵ}}^2}\}$ 且ε>0。

如果如之前的那样，z～N(y,σ²)，那么我们得到似然函数或直接概率分布：

$p\left(z\right|x,\mathrm{\ϵ},\mathrm{\Θ},\mathrm{\σ},M)=\∫p(z,y|x,\mathrm{\ϵ},\mathrm{\Θ},\mathrm{\σ},M)\mathrm{dy}=$

$\∫p\left(z\right|y,\mathrm{\σ},M)p\left(y\right|x,\mathrm{\ϵ},\mathrm{\Θ},M)\mathrm{dy}=$

$\∫\frac{1}{\sqrt{2\mathrm{\π}}\mathrm{\σ}}\exp [-\frac{{(z-y)}^2}{2{\mathrm{\σ}}^2}]\frac{1}{\sqrt{2\mathrm{\π}}\mathrm{\ϵ}}\exp \{-\frac{{[y-M(x,\mathrm{\Θ})]}^2}{2{\mathrm{\ϵ}}^2}\}\mathrm{dy}$

由此，可以通过应用例如贝叶斯规则来估计感兴趣的量Θ和σ：

p(ε,Θ,σ|z,x,M)∝p(ε,Θ,σ|M)p(z|x,ε,Θ,σ,M)

由该第二实施例可以看出，硬信息片段和对应的严格的软信息的片段可能不仅在先验概率分布的水平上（例如在第一个实施例中）是不同的，而且还可能在似然函数的水平上不同。最终获得通常不同的感兴趣的量的估计。该实施例显示了：根据本发明，如何以一般的方式“软化”（半）参数模型和硬的血液动力学参数估计方法。

现在考虑第三实施例以阐释定量信息片段。使S(t)为实时信号且使s=[S(t₁),...,S(t_N)]为其值的矢量。考虑硬定量信息片段“S(t)是平滑的”。该定量信息规范地转换成定量信息其中||||₂代表函数L²范数。这遵循：S(t)是线性，的且存在(a,b)，例如S(t)=at+b。事实上，没有什么比直线更平滑。该硬信息最终转换成先验概率分布p(S(t)|t,a,b)=δ[S(t)-at-b]。实际上，其实质上是估计偶(a,b)。

现在考虑对应的严格的软定量信息片段“S(t)大致是平滑的”。此信息转换成严格的软定量信息

大致是较大的”。由于处于函数空间，为了获得转换此软信息的先前随机过程，似乎不能直接应用最大熵原理。然而，可以考虑定量信息

大致是较大的”的离散-时间版本，其中

是通过使用诸如

的二阶导数的有限差值数值近似方案从s获得的。

因此，我们得到

${\left|\right|\frac{d^{2}s}{{\mathrm{dt}}^2}\left|\right|}^2={\left(\mathrm{Ds}\right)}^T\left(\mathrm{Ds}\right)=s^T\left(D^{T}D\right)s\≡s^T\mathrm{\Σs}$

接着，应用最大熵原理提供矢量s的多变量高斯先验概率分布：

$p\left(s\right|\mathrm{\ϵ})\∝\exp (-\frac{s^T\mathrm{\Σs}}{2{\mathrm{\ϵ}}^2})$

再次，可以评估硬信息与严格的软信息之间的不同。在硬信息的情况下，s_i被约束在直线y=at+b上，而在严格的软信息的情况中，他们大致远离该直线，这取决于必须被估计的ε>0的值。第一种情况实质上是估计两个参数(a,b)，而第二种情况实质上是估计N+1个参数（矢量s的N个参数和ε超级参数）。

根据本发明的软信息和严格的软信息的概念被明确界定，现在回溯关于脉冲响应h(t)的硬信息和软信息。

硬定量信息片段包括h(t₁)=h(0)=0，例如，根据右矩形法，

两种信息片段都允许将待估计的值/参数的数量减小至N-2。正如我们将看到的，明智的选择在于，保持值h'=[h(t₂),...,h(t_N-1)]^T并将h(t_N)表达为除此之外，还可以考虑硬定量信息片段

h(t_i)≥0因此，仍通过将那些硬定量信息片段与纯定量的严格的软信息相结合来分配矢量h'的联合先验概率分布。考虑例如-如前所述-严格的软定量信息

“h(t)大致是平滑的”。正如前面解释过的，此定量信息

可以先转换成严格的软定量信息

大致是较大的”。

在取样时间点t_i=1,N处进行时间离散化后，二阶导数h的二阶数值逼近通过下式给出：

$\∀i=2,N-1,\frac{d^{2}h\left(t_i\right)}{{\mathrm{dt}}^2}=\frac{h\left(t_{i-1}\right)+h\left(t_{i+1}\right)-2h\left(t_i\right)}{\mathrm{\Δ}{t}^2}+O\left(\mathrm{\Δ}{t}^2\right)$

还可以使用高阶数值逼近，如四阶公式：

$\∀i=3,N-2,$

$\frac{d^{2}h\left(t_i\right)}{{\mathrm{dt}}^2}=\frac{-\frac{1}{12}h\left(t_{i-2}\right)+\frac{4}{3}h\left(t_{i-1}\right)-\frac{5}{2}h\left(t_i\right)+\frac{4}{3}h\left(t_{i+1}\right)-\frac{1}{12}h\left(t_{i+2}\right)}{\mathrm{\Δ}{t}^2}+O\left(\mathrm{\Δ}{t}^4\right)$

那些数值逼近可以以矩阵符号被写为

其中D是例如在二阶近似的情形下的尺寸N的方阵D，且h'如上界定，其中，

因此，h'的二阶导数的欧几里得范数的平方可以被写为：

${\left|\right|\frac{d^2{h}^{\′}}{{\mathrm{dt}}^2}\left|\right|}^2={\left(\mathrm{Dh}\right)}^T\left(\mathrm{Dh}\right)=h^T\left(D^{T}D\right)h\≡h^T{W}_{1}h.$

由此，可以将“定性地假设脉冲响应h(t)大致是平滑的”定量地转换为“假设

大致是较大的”。因此，在具有相同的欧几里得范数

的所有连续概率分布中寻找先验概率分布p(h'|h(t₁),h(t_N),I)（边界点h(t₁)和h(t_N)被单独处理）。应用最大熵原理，在支持超象限[0，+∞]^N-2的所有那些分布中选择具有最高微分Shannon熵（根据Lebesgue参照测量）的分布（因为其是最不确定的并且因此是最诚实的），从而获得[0，+∞]^N-2上的条件多变量截断高斯分布，以及常数矢量数学期望值Μ=(μ₁,...,μ₁)^T和协方差矩阵

（或互补累积密度函数的值的矢量在[0,1]^N-2上截断的多变量高斯分布）：

$p\left(h^{\′}\right|h\left(t_1\right),h\left(t_N\right),{\mathrm{\μ}}_1,{\mathrm{\ϵ}}_1,E_s,I_h^1,M)=C_1({\mathrm{\μ}}_1,{\mathrm{\ϵ}}_1,{\mathrm{\σ}}_S)\exp \{-\frac{{\mathrm{\ϵ}}_1^2{(h-M)}^T{W}_1(h-M)}{2{\mathrm{\σ}}_S^2}\}$

其中C₁(μ₁,ε₁,σ_S)是归一化常数

$C_1({\mathrm{\μ}}_1,{\mathrm{\ϵ}}_1,{\mathrm{\σ}}_S)={\left[\underset{h^{\′}\∈{[0,+\∞]}^{N-2}}{\∫}\exp \right\{-\frac{{\mathrm{\ϵ}}_1^2{(h-M)}^T{W}_1(h-M)}{2{\mathrm{\σ}}_S^2}\left\}d{h}^{\′}\right]}^{-1}$

因此，两个新的超参数出现在所获得的全身灌注模型中，标量数学期望值μ₁和逆函数变量ε₁在灌注模型中起到正则化参数的作用。ε₁量化了先前的软定量信息

相对于由实验数据s和s_a提供的硬定量信息的柔软度。

此外，通过定义：

$p(h\left(t_1\right),h\left(t_N\right)|E_S,I_h^1,M)=p\left(h\left(t_1\right)\right|E_S,I_h^1,M)p\left(h\left(t_N\right)\right|E_S,I_h^1,M)=$

$\mathrm{\δ}\left[h\left(t_1\right)\right]\·\mathrm{\δ}[h\left(t_N\right)-1/\mathrm{\Δt}\·\underset{i=2}{\overset{N-1}{\mathrm{\Σ}}}h\left(t_i\right)].$

然而，为了将“延迟τ几乎从不等于取样时间点t_i,i=1,N”这一事实考虑在内，使

$p\left(h\left(t_1\right)\right|{\mathrm{\μ}}_{\mathrm{1,1}},{\mathrm{\ϵ}}_{\mathrm{1,1}},E_S,I_h^1,M)=C_{\mathrm{1,1}}({\mathrm{\μ}}_{\mathrm{1,1}},{\mathrm{\ϵ}}_{\mathrm{1,1}},{\mathrm{\σ}}_S)\exp \{-\frac{{\mathrm{\ϵ}}_{\mathrm{1,1}}^2{[h\left(t_1\right)-{\mathrm{\μ}}_{\mathrm{1,1}}]}^2}{2{\mathrm{\σ}}_S^2}\}$

其中

和C_1,1(μ_1,1,ε_1,1,σ_S)是归一化常数

$C_{\mathrm{1,1}}({\mathrm{\μ}}_{\mathrm{1,1}},{\mathrm{\ϵ}}_{\mathrm{1,1}},{\mathrm{\σ}}_S)={\left[\underset{h\left(t_1\right)\∈\left[\mathrm{0,1}\right]}{\∫}\exp \right\{-\frac{{\mathrm{\ϵ}}_{\mathrm{1,1}}^2{[h\left(t_1\right)-{\mathrm{\μ}}_{\mathrm{1,1}}]}^2}{2{\mathrm{\σ}}_S^2}\left\}\mathrm{dh}\left(t_1\right)\right]}^{-1},$

从而表明且只表明

最后，支持用于信息的超象限[0，+∞]^N-2的h的先验概率分布可以被写作：

$p\left(h\right|{\mathrm{\μ}}_1,{\mathrm{\ϵ}}_1,E_S,I_h^1,M)=$

$\mathrm{\δ}\left[h\left(t_1\right)\right]\mathrm{\δ}[h\left(t_N\right)-1/\mathrm{\Δt}+\underset{i=2}{\overset{N-1}{\mathrm{\Σ}}}h\left(t_i\right)]C_1({\mathrm{\μ}}_1,{\mathrm{\ϵ}}_1,{\mathrm{\σ}}_S)\exp \{-\frac{{{\mathrm{\ϵ}}_1}^2{(h-M)}^T{W}_1(h-M)}{2{{\mathrm{\σ}}_S}^2}\}$

一般来说，应该将h的先前分布的超参数的矢量表示为Ε_h（或Ε_b），如Ε_h=(μ₁,ε₁)或Ε_h=(μ₁,ε₁,μ_1,1,ε_1,1)。

仍分配Ε_h的先验概率分布。可以分配例如非信息型的不合适的Bayes-Laplace-Lhoste-Jeffreys分布或 $p({\mathrm{\μ}}_1,{\mathrm{\ϵ}}_1,{\mathrm{\μ}}_{\mathrm{1,1}},{\mathrm{\ϵ}}_{\mathrm{1,1}}|I_h^1,M)\∝{\mathrm{\ϵ}}_1^{-1}{\mathrm{\ϵ}}_{\mathrm{1,1}}^{-1}.$

可以以相同的方式引入严格的软定量信息

大致是较大的”且更通常地，引入严格的软信息

k=1,2,...

大致是较大的”。

最后，通过应用最大熵原理，获得了h的新的先验多变量截断高斯概率分布。

值得注意的是，如前所述的软信息片段

……是仅有的严格的软信息片段，在我们的问题中可以被恰当地引入的唯一的约束条件是：他们不是可以被实验证实或否定的任意假设，正相反，他们仅仅表达了基本的生理特性，如果没有这一点，估计血液动力学参数、脉冲响应或互补累积密度函数的问题实际上将是毫无意义的：他们对于我们的问题来说是逻辑上充分必要的。反之，任何其他信息都将是可以通过实验证实的简单的假设。

尽管这样，根据本发明，可以引入其他严格的软信息，其为简单的作业假设。例如，可以指明且仅指明：脉冲响应大致符合给定的函数形式而不约束它（通过硬信息）以便精确地符合此形式。增加此种类的半数量且严格的软信息允许确定脉冲响应可以在何种程度上被所提议的函数形式描述。因此，假设如想要引入如前所述的严格的软信息

是大约符合参数或半参数模型M_h:f(t,Θ_h)>>，Θ_h是所述模型的参数的矢量。如前所述，已经提出这样的参数模型，如通过下述给出的两参数Γ模型：

$M_h:\left\{\begin{array}{c}f(t,{\mathrm{\Θ}}_h)=t^{\frac{\mathrm{MTT}}{\mathrm{\β}}-1}{e}^{-t/\mathrm{\β}}\text{/}{\mathrm{\β}}^{\frac{\mathrm{MTT}}{\mathrm{\β}}}/\mathrm{\Γ}\left(\frac{\mathrm{MTT}}{\mathrm{\β}}\right)\mathrm{MTT}>0,\mathrm{\β}>0\\ {\mathrm{\Θ}}_h=(\mathrm{MTT},\mathrm{\β})\\ \mathrm{\Γ}\left(\right):\mathrm{fonction\; Gamma}{d}^{\′}\mathrm{Euler}\end{array}\right.$

请注意，在此情形中，可以直接估计参数MTI而无需如前所述从数值上估计脉冲响应的第一分量（或互补累积密度函数的积分）。

要求h大致符合给定的函数形式f(t,Θ_h)，这相当于定量地要求残余

的矢量的欧几里得范数大致是较大的。再次应用最大熵原理，以相同的方式，发现h的先验概率分布是在超象限[0，+∞]^N-2上截断的多变量高斯分布

$p\left(h\right|E_h,E_S,{\mathrm{\Θ}}_h,I_h^3,M_h,M)=C_3({\mathrm{\ϵ}}_3,{\mathrm{\σ}}_S,{\mathrm{\Θ}}_h)\exp \{-\frac{{\mathrm{\ϵ}}_3^2}{2{\mathrm{\σ}}_S^2}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{[h\left(t_i\right)-f(t_i,{\mathrm{\Θ}}_h)]}^2\},$

其中C₃(ε₃,σ_S,Θ_h)是归一化参数

$C_3({\mathrm{\ϵ}}_3,{\mathrm{\σ}}_S,{\mathrm{\Θ}}_h)={\left[\underset{h\∈{[0,+\∞]}^N}{\∫}\exp \right\{-\frac{{\mathrm{\ϵ}}_3}{2{{\mathrm{\σ}}_S}^2}{\left|\right|h-f(t,{\mathrm{\Θ}}_h)\left|\right|}^2\left\}\mathrm{dh}\right]}^{-1}$

以同样的方式，可以引入诸如“h(t)大约符合各值的给定矢量”的软信息片段

应用最大熵原理，获得先验概率分布：

$p\left(h\right|E_h,E_S,I_h^4,M)\∝C_4({\mathrm{\ϵ}}_4,{\mathrm{\σ}}_S,\stackrel{\‾}{h})\exp \{-\frac{{{\mathrm{\ϵ}}_4}^2}{2{{\mathrm{\σ}}_S}^2}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{[h\left(t_i\right)-\stackrel{\‾}{h}\left(t_i\right)]}^2\},$

其中

是归一化常数

$C_4({\mathrm{\ϵ}}_4,{\mathrm{\σ}}_S,\stackrel{\‾}{h})={\left[\underset{h\∈{[0,+\∞]}^N}{\∫}\exp \right\{-\frac{{{\mathrm{\ϵ}}_4}^2}{2{{\mathrm{\σ}}_S}^2}{\left|\right|h-\stackrel{\‾}{h}\left|\right|}^2\left\}\mathrm{dh}\right]}^{-1}$

本发明还允许将关于脉冲响应h(t)（或互补累积密度函数）的若干软信息片段与它们对应的先验概率分布相结合。因此，如果

是用于转换信息片段

的n先验概率分布且具有各自的正则化参数

（如

等），那么通过使

和

将n信息考虑在内的h的先验分布可以写作 $p\left(h\right|E_h,I_h,M)=\underset{k=1}{\overset{n}{\mathrm{\Π}}}p\left(h\right|E_h^k,I_h^k,M).$

为了编码关于局部动脉输入函数的信息I_a，根据本发明的方法还包括步骤52c，其分配先验概率分布 $p(a,E_{S_a}|I_a,I_{S_a},M)=p\left(a\right|E_{S_a},I_a,M)\·p\left(E_{S_a}\right|I_{S_a},M).$ 这样的分布通过引入和组合关于动脉输入函数的硬和/或软信息片段、以与脉冲响应h有关的方式相同的方式被分配。例如，可以规定：动脉输入函数大致是平滑的、正的、单峰的、双峰的、在原点是零的、渐进零的、具有给定的区域；或可以指明且仅指明：其大致符合给定的参数或半参数模型C_a(t,Θ_a)，其中Θ_a是参数的矢量，该模型为例如十一参数“三-反差系数”模型：

$M_a:\left\{\begin{array}{c}{C}_a(t,{\mathrm{\Θ}}_a)=\\ \frac{a{(t-t_0)}^{{\mathrm{\α}}_0-1}{e}^{-(t-t_0)/{\mathrm{\β}}_0}}{{\mathrm{\β}}_0^{{\mathrm{\α}}_0}\mathrm{\Γ}\left({\mathrm{\α}}_0\right)}+\frac{b{(t-t_1)}^{{\mathrm{\α}}_1-1}{e}^{-(t-t_1)/{\mathrm{\β}}_1}}{{\mathrm{\β}}_1^{{\mathrm{\α}}_1}\mathrm{\Γ}\left({\mathrm{\α}}_1\right)}+\frac{(1-a-b){(t-t_2)}^{{\mathrm{\α}}_2-1}{e}^{-(t-t_2)/{\mathrm{\β}}_2}}{{\mathrm{\β}}_2^{{\mathrm{\α}}_2}\mathrm{\Γ}\left({\mathrm{\α}}_2\right)}\\ {\mathrm{\Θ}}_a=(a,b,{\mathrm{\α}}_0,{\mathrm{\β}}_0,t_0,{\mathrm{\α}}_1,{\mathrm{\β}}_1,t_1,{\mathrm{\α}}_2,{\mathrm{\β}}_2,t_2)\\ \mathrm{\Γ}\left(\right):\mathrm{fonctionGamma}{d}^{\′}\mathrm{Euler}\end{array}\right.$

如前所述，获得了高斯先验概率分布

$p\left(a\right|E_a^1,E_{S_a},{\mathrm{\Θ}}_a,I_a^1,M_a,M)=C({\mathrm{\ϵ}}_a^1,{\mathrm{\σ}}_{S_a},{\mathrm{\Θ}}_a)\exp \{-\frac{{\left({\mathrm{\ϵ}}_a^1\right)}^2}{2{\mathrm{\σ}}_{S_a}^2}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{[a\left(t_i\right)-C_a(t_i,{\mathrm{\Θ}}_a)]}^2\}$

其中是归一化常数

$C({\mathrm{\ϵ}}_a^1,{\mathrm{\σ}}_{S_a},{\mathrm{\Θ}}_a)={\left[\underset{a\∈{[0,+\∞]}^N}{\∫}\exp \right\{-\frac{{\left({\mathrm{\ϵ}}_a^1\right)}^2{\left|\right|a-C_a(t,{\mathrm{\Θ}}_a)\left|\right|}^2}{2{{\mathrm{\σ}}_{S_a}}^2}\left\}\mathrm{da}\right]}^{-1}$

此外，可以指出且仅指出（否则，参考给定了动脉输入函数的情形）矢量a大致接近各值的给定矢量

例如局部动脉输入函数：获得了先验概率分布

$p\left(a\right|E_a^2,E_{S_a},I_a^2,M)=C_a({\mathrm{\ϵ}}_a^2,{\mathrm{\σ}}_{S_a},\stackrel{\‾}{a})\exp \{-\frac{{\left({\mathrm{\ϵ}}_a^2\right)}^2}{2{{\mathrm{\σ}}_{S_a}}^2}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{[a\left(t_i\right)-\stackrel{\‾}{a}\left(t_i\right)]}^2\},$

其中

是归一化常数

$C_a({\mathrm{\ϵ}}_a^2,{\mathrm{\σ}}_{S_a},\stackrel{\‾}{a})={\left[\underset{a\∈{[0,+\∞]}^N}{\∫}\exp \right\{-\frac{{\left({\mathrm{\ϵ}}_a^2\right)}^2{\left|\right|a-\stackrel{\‾}{a}\left|\right|}^2}{2{{\mathrm{\σ}}_{S_a}}^2}\left\}\mathrm{da}\right]}^{-1}$

根据本发明的方法还包括配置步骤52d，其分配先验概率分布 $p(E_S,E_{S_a},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a}|I_S,I_{S_a},M),$ 该概率分布通常被因式分解为：

$p(E_S,E_{S_a},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a}|I_S,I_{S_a},M)$

$=p\left(E_S\right|I_S,M)\·p\left(E_{S_a}\right|I_{S_a},M)\·p\left({\mathrm{\Θ}}_S\right|I_S,M)\·p\left({\mathrm{\Θ}}_{S_a}\right|I_{S_a},M)$

例如，我们得到如下非信息型先验概率分布：

$p(E_S,E_{S_a},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a}|I_S,I_{S_a},M)$

$=p\left({\mathrm{\σ}}_S\right|I_S,M)\·p\left({\mathrm{\σ}}_{S_a}\right|I_{S_a},M)\·p\left(S_0\right|I_S,M)\·p\left(S_{0_a}\right|I_{S_a},M)$

$={({\mathrm{\σ}}_S\·{\mathrm{\σ}}_{S_a})}^{-1}/({S_0}^{\max }-{S_0}^{\min })/({S_{0_a}}^{\max }-{S_{0_a}}^{\min })$

将在步骤52b至52d中引入的不同的超参数考虑在内，联合先验概率分布53可以改写为 $p(a,E_a,h,E_h,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},E_S,E_{S_a}|I_{\mathrm{\Θ}},I_h,I_a,I_S,I_{S_a},M).$

为了简化下面的表达，我们用

表示作为处理单元的配置参数而输入的信息片段组。

正如之前描述过的，给定直接概率分布54和联合先验概率分布53，通过应用贝叶斯规则得到所有参数的联合后验概率分布55：

$p(a,E_a,h,E_h,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},E_S,E_{S_a}|s,s_a,I)=$

$\frac{p(a,E_a,h,E_h,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},E_S,E_{S_a}|I)\·p(s,s_a|a,h,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},E_S,E_{S_a},I)}{p(s,s_a|I)}\∝$

$p(a,E_a,h,E_h,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},E_S,E_{S_a}|I)\·p(s,s_a|a,h,\mathrm{\Θ},{\mathrm{\Θ}}_S,{\mathrm{\Θ}}_{S_a},E_S,E_{S_a},I)$

完成初始配置后，本发明允许估计感兴趣的量，我们将在矢量

的所有要素中指出θ。例如，θ=BF或，θ=τ或θ=h(t_i)等。

因此，根据本发明的方法包括步骤56，其评估θ的边缘后验分布，即，

$p\left(\mathrm{\θ}\right|s,s_a,I)=\underset{\mathrm{\Ξ}\backslash \mathrm{\θ}}{\∫}p\left(\mathrm{\Ξ}\right|s,s_a,I)d(\mathrm{\Ξ}\backslash \mathrm{\θ}).$

可以由该边缘后验分布获得57θ的估计

例如，通过采用该分布

的数学期望值，根据二次损失函数

获得贝叶斯估计，其中||||是欧几里得范数。可以以同样的方式通过

获得最大后验估计

举例来说，可以获得在每一个取样时间点t_i的脉冲响应的值h(t_i),i=1,N的边缘后验分布，这通过使所有其他时间点以下式边缘化而实现：

$\∀i=1,N,p\left(h\left(t_i\right)\right|s,s_a,I)=\underset{\underset{j\≠i}{h\left(t_j\right)}}{\∫}p\left(h\right|s,s_a,I)\mathrm{dh}\left(t_1\right)...\mathrm{dh}\left(t_{j\≠i}\right)...\mathrm{dh}\left(t_N\right)$

且随后获得那些诸如

或

的值的估计。

以同样的方式，还可以获得57在任何时间点x（不一定等于取样时间点t_i）的脉冲响应h(x)的值的估计

（或

）。基于此目的，足以在如前所述的互补累积密度函数R(t)的值的表达式中引入h(x)并计算其边缘概率分布。甚至期望在该估计问题中引入这样的额外的时间点x₁,...,x_L，因为所考虑的时间点的数目越大，诸如

$\underset{0}{\overset{t}{\∫}}h\left(\mathrm{\τ}\right)\mathrm{d\τ}\underset{0}{\overset{+\∞}{\∫}}h\left(\mathrm{\τ}\right)\mathrm{d\τ}=1\underset{0}{\overset{t}{\∫}}{C}_a\left(\mathrm{\τ}\right)\·[H(t-\mathrm{\τ})-\underset{0}{\overset{t-\mathrm{\τ}}{\∫}}h\left(\mathrm{\υ}\right)\mathrm{d\υ}]\mathrm{d\τ}$ 的积分的数值逼近以及

或

的数值逼近越好。相应地，所得到的估计也更好。

于是，由于通过应用数值积分法（例如，此处为正矩形法），可以获得此参数的估计，例如“平均”估计

或最可能的估计 ${\widehat{\mathrm{MTT}}}^P=\mathrm{\Δt}\·\underset{i=2}{\overset{N}{\mathrm{\Σ}}}{t}_i\·{\hat{h}}^P\left(t_i\right).$

此外，根据本发明的方法可以包括步骤58，其用于获得参数θ的估计的精确度的估计以及甚至获得这些估计的置信区间。本发明使得这样的方法可以进一步包括步骤59，其用于获得所述置信区间的可能性。例如，估计的精确度可以通过θ的边缘后验概率分布的协方差矩阵而量化：

$\widehat{\mathrm{\σ}}{\widehat{\mathrm{\θ}}}^{Q^2}=\underset{\mathrm{\θ}}{\∫}(\mathrm{\θ}-{\widehat{\mathrm{\θ}}}^Q){(\mathrm{\θ}-{\widehat{\mathrm{\θ}}}^Q)}^{T}p\left(\mathrm{\θ}\right|s,s_a,I)\mathrm{d\θ}$

接着，我们得到θ在“单标准差”处的置信（超）区间

和θ属于区间的概率或（等同地）几率

还可以获得参数BV=BF·MTT、互补累积密度函数b的值的矢量、静脉输出函数v=Ah的值的矢量或理论浓度-时间曲线的值的矢量c=BF.Ab的估计、置信区间和置信区间的可能性，因为在给定若干随机变量的联合概率分布的前提下可以计算它们的任意函数的概率密度函数。例如，通过线性化数学期望值，可以遵循例如右矩形法立即从那些脉冲响应 ${\hat{b}}^Q={(\mathrm{1,1}-\mathrm{\Δt}\·{\hat{h}}^Q\left(t_2\right),...,1-\mathrm{\Δt}\·\underset{i=2}{\overset{N}{\mathrm{\Σ}}}{\hat{h}}^Q\left(t_N\right))}^T$ 和 ${\hat{b}}^P={(\mathrm{1,1}-\mathrm{\Δt}\·{\hat{h}}^P\left(t_2\right),...,1-\mathrm{\Δt}\·\underset{i=2}{\overset{N}{\mathrm{\Σ}}}{\hat{h}}^P\left(t_N\right))}^T$ 获得互补累积密度函数R(t)的值的估计。以同样的方式，也可以从c和Θ_S的联合概率分布p(c,Θ_S|s,s_a,I)推出理论灌注信号的值的概率分布s_th=[S_th(t₁),...,S_th(t_N)]^T，因为s_th=Ψ(c,Θ_S)。可以通过此分布、以类似于上述的方法获得57估计

以及置信高级区间58和与这些超区间相关的可能性59。

本发明提出：置信区间或所述置信区间的可能性可以允许设置62处理单元的配置。因此，可以修改配置参数并提供较高质量的估计。

本发明还提供了一种包括步骤60的方法，其用于计算理论和实验灌注加权信号之间的残余

本发明允许计算这些矢量之间的各种统计信息或距离

最典型的是平方误差的总和-SSE-

这些多种统计允许量化灌注模型对实验数据的拟合优度。这样，可以获得感兴趣的每一个体元的模型的“误差图”。基于上述原因（即，过拟合问题），可以通过假定感兴趣的每一个体元的灌注模型并计算实验数据(s,s_a)的概率而特别有利地量化灌注模型对实验数据s和s_a的拟合优度：

$p(s,s_a|I)=\underset{\mathrm{\Ξ}}{\∫}p(s,s_a|\mathrm{\Ξ},I)p\left(\mathrm{\Ξ}\right|I)\mathrm{d\Ξ}.$

在此情形中，误差图将基于1-p(s,s_a|I)。

本发明还提出：可以应用61各种统计学检验或各种图形诊断技术（诸如Q-Q散点图或Poincaré回归映射图）以检查残余r(t_i)实际上是否是独立的、恒等分布的高斯分布等。本发明因此提出：通过迭代过程以及试错法，可以校正和改善62配置过程50，尤其是理论灌注模型，以便在建模、理解和处理灌注现象方面取得进展以及获得血液动力学参数、脉冲响应、互补累积密度函数、动脉输入函数或静脉输出函数的更好、更细致的估计。

现在结合图8描述根据本发明的用于第二应用实施例的方法，其中假设理论动脉输入函数是以绝对确定性和无限准确性给出的（高达延迟τ）。

因此，实验灌注模型M可以写为：

$M:\left\{\begin{array}{c}s=\mathrm{\Ψ}(c,{\mathrm{\Θ}}_S)+\mathrm{\ξ}\\ c=\mathrm{BF}\·\mathrm{Ab}(t-\mathrm{\τ})=\mathrm{BF}\·B(t-\mathrm{\τ})a\end{array}\right.$

其中所有量都以如上方式界定，不同的是此时假设理论动脉输入函数的值的矢量a=[C_a(t₁),...,C_a(t_N)]^T是已知的。

因此，联合直接概率分布写为p(s|a,h,Θ,Θ_S,Ε_S,I)，而联合先验概率分布写为p(h,Ε_h,Θ,Θ_S,Ε_S|I)，且I=(I_Θ,I_h,I_S,M)。这些分布被按照如上所述来分配。贝叶斯规则变成

p(h,Ε_h,Θ,Θ_S,Ε_S|s,a,I)∝p(h,Ε_h,Θ,Θ_S,Ε_S|a,I)·p(s|a,h,Θ,Θ_S,Ε_S,I)

随后以与前述相同的方式获得任何参数θ∈Ξ=(h,Ε_h,Θ,Θ_S,Ε_S)的估计、置信区间和这些置信区间的可能性。

现在结合图8描述根据本发明的用于第三应用实施例的方法，其中未给出理论动脉输入函数且未测量动脉输入信号。

于是，实验灌注模型M可以写为：

$\left\{\begin{array}{c}s=\mathrm{\Ψ}(c,{\mathrm{\Θ}}_S)+\mathrm{\ξ}\\ c=\mathrm{BF}\·\mathrm{Ab}\left(t\right)=\mathrm{BF}\·B\left(t\right)a\end{array}\right.$

其中所有量都以如上方式界定。

联合直接概率分布仍写为p(s|a,h,Θ,Θ_S,Ε_S,I)，而联合先验概率分布此时写为p(a,Ε_a,h,Ε_h,Θ,Θ_S,Ε_S|I)，且I=(I_Θ,I_h,I_a,I_S,M)。贝叶斯规则变成p(a,Ε_a,h,Ε_h,Θ,Θ_S,Ε_S|s,I)∝p(a,Ε_a,h,Ε_h,Θ,Θ_S,Ε_S|I)·p(s|a,h,Θ,Θ_S,Ε_S,I)

接着以与前述相同的方式获得任何参数θ∈Ξ=(a,Ε_a,h,Ε_h,Θ,Θ_S,Ε_S)的估计、置信区间和这些置信区间的可能性。

值得注意的是，无论讨论的是什么应用实施例，与已知的方法不同，根据本发明的方法是精确方法，因为本发明的方法包括且仅包括将我们具有感兴趣的量的先验的定性、定量或半定量信息转换成贝叶斯概率理论，以便明确地确定与由实验测量提供的这些量有关的后验信息。没有可以被证实或否定（如果确实其可以被证实的话）的任意假设（尤其是对脉冲响应或对互补累积密度函数而言）是必要的，因为本发明仅仅引入了最不确定的概率分布函数，其对逻辑上充分必要的各种软定性约束条件进行编码，以解决问题。

本发明遵循：假定并未给出而至多是测量了动脉输入函数，这些方法允许测试所提出的动脉输入函数是否能够有效地对应于感兴趣的每一个体元的“真正的”局部动脉输入函数：事实上，如果动脉输入信号是不合适的且并不对应于真正的局部动脉输入函数，那么可能不存在作为对问题的解且同时还满足各种先验约束条件（或至少是具有不可忽略的概率的先验）的参数组（血液动力学参数、互补累积密度函数等）。在此情形中，概率理论将所提出的理论动脉输入信号解释为噪声：与通常从更合适的动脉输入信号获得的相比，标准偏差σ_a/ε_a将大得多。

因此，本发明尤其感兴趣的应用是测试根据现有技术的全身或局部动脉输入函数的不同选择和估计方法。如果能够证实：估计，尤其是从那些全身或局部动脉输入函数获得的那些正则化参数Ε_a的估计在一个体元和另一个体元之间偏差得过大的话，将得出如下结论：需要引入新的、更合适的（局部的）对这些函数的选择方法，或，需要求助于并不要求给出动脉输入函数或不需要提前测量动脉输入信号的方法，这恰恰是上述第三种方法的目的。

提及本发明通过灌注加权分析成像系统执行的主要实施步骤作为一种应用实施例，如图1或2所描述的：

-打开患者记录或打开由输入处理单元4（预处理单元7）输入的图像序列，以便选择感兴趣的图像序列-尤其是在一定时间段内的灌注加权图像I1到In，据此获得每一个体元的灌注信号S(t)，如通过图5a阐释的；

-通过人际界面5预览图像，以便允许用户6识别感兴趣的切片或区域；

-由配置参数（输入信息片段）配置处理单元4，以便允许实施根据本发明的估计方法；

-选择待估计的一个或多个感兴趣的量；

-通过处理单元4估计诸如人脑的器官的血液动力学参数14，诸如

或

-当未给出动脉输入函数时，诸如Ε_h、Ε_S、Θ_S或Ε_a、

或

的其他参数的可选估计，

-脉冲响应的值的矢量

互补累积密度函数

动脉输入函数的值的矢量

静脉输出函数的值的矢量

理论浓度的值的矢量

理论信号

或残余

的可选估计；

-向人机界面5传送所述感兴趣的量14的估计，使得人机界面5最后例如以图的形式显示这些估计，其中每一个像素的强度或颜色取决于计算值，以便将它们的内容呈现给执业医师；

-选择性地显示用户所选择的感兴趣的一些体元的脉冲响应、互补累积密度函数、动脉输入函数、静脉输出函数、理论浓度、理论信号或残余的估计；

-选择性地显示一个或若干个感兴趣的参数诸如血液动力学参数的置信图或可能性图；

-选择性地显示用户所选择的感兴趣的一些体元的一些脉冲响应、互补累积密度函数、动脉输入函数等的置信区间或这些区间的可能性；

-选择性地显示实验数据与全身非参数灌注模型之间的一个或若干个距离的误差图，尤其是在给定全身灌注模型的前提下显示实验数据的概率；

-辅助选择感兴趣的病理区域，其特征在于，脉冲响应（或互补累积密度函数）或局部动脉输入函数的一个或若干个血液动力学参数的异常分布；

-由处理单元估计可能与损伤区域有关的异常灌注组织区域的体积，据此，执业医师可以决定治疗行为（例如，静脉血栓溶解，以便分解血凝块）；

-由处理单元估计一些量，例如损伤区域与异常分布区域的体积比，据此，执业医师可以完善其诊断和治疗决策（例如，静脉血栓溶解，以便分解血凝块）；

如前所述，处理单元4的配置参数50可以通过单元本身（执行过程50）来执行。可选择地，所述配置可以为：存储和选择包含联合后验概率分布的库，这取决于想要估计的感兴趣的量。此库的配置可以通过能够与处理单元4共同配合的专门的单元来实现。

可选择地，在估计一些感兴趣的量的置信区间、置信区间的可能性以便完善所述配置后，可能会出现迭代。通过规定实验数据与全身灌注模型之间的距离，尤其是通过在给定全身灌注模型的前提下显示灌注数据的概率，可以引起配置的升级。

因此，本发明的目的在于以“参数地图”的形式显示参数估计，其中每一像素的强度或颜色取决于所计算的值，例如，以线性方式。而且，本发明旨在以“置信图”的形式显示这些估计的标准偏差以及以“可能性图“的形式显示关于对应的置信区间的可能性。关于脉冲响应、互补累积密度函数、动脉输入函数、静脉输出函数、浓度-时间曲线、灌注信号或残余值的矢量的估计，本发明旨在以由用户选择的每一个体元的时间序列的形式对其进行显示。最后，本发明旨在显示实验信号与非参数灌注模型或实验数据的概率之间的距离，假定此模型为“误差图”的形式。

图9到12允许阐释以图的形式呈现一些感兴趣的量的显示模式，诸如根据本发明所估计的血液动力学参数14或甚至标准偏差或与它们相关的概率。

因此，对通过核磁共振成像分析的人脑来说，图9允许显示脑血容量CBV的估计。这样的图（458×458像素）允许突出可能的缺血区域80。事实上，可以通过合适的界面6观察到：与对侧脑半球相比，右部大脑后动脉区域内的参数CBV明显增大了。通过图9阐释的图自己揭露了缺血后的血管舒张。

图10允许阐释与缺血中风情形中的脑血流的估计有关的图（458×458像素）。通过分析该图，可以观察到：与缺血后的对侧脑半球相比，左部大脑后动脉区域内的参数CBF（脑血容量）减小了。这样的图允许突出可能的缺血区域80。

图11允许阐释与平均通过时间MTI的估计有关的图（458×458像素）。可以通过分析该图观察到：与局部缺血后的对侧脑半球相比，右部大脑后动脉区域内的参数MTI明显增大了。

图12允许描述与“脑血流属于置信区间

的概率的估计有关的图（458×458像素）。通过分析该图，可以观察到：除了异常模型拟合外，概率以约0.68为中心。这是有资格预计CBF的后验概率分布是否符合高斯法则的值。

本发明能够通过上述图为用户提供多种有用的信息，这些信息是无法通过使用现有已知技术获得的。可以通过调整如图1或2所示的处理单元4来实现这一点，因为用于与外界通信的装置能够以适于人机界面5的格式传送所估计的参数，该人机界面5能够以如图9至13所阐释的图的形式向用户6呈现所估计的参数。

通过本发明可提供的信息片段量更大且更公正。执业医师可获得的信息能增加执业医师在确定诊断和治疗决定方面的置信度。

为了改善根据本发明的系统的性能，本发明规定：处理单元可以包括用于并行计算图像像素的装置，这需要估计血液动力学参数、互补累积密度函数或动脉输入函数。这可以通过使用诸如图形处理器单元（GPU）的硬件技术、诸如并行Monte-Carlo法的计算机集群或软件技术等实现。可选择地，根据本发明的处理单元可以依赖于远程计算机装置。这样，可以显著减少计算次数。

Claims (14)

1.一种用于估计（57）器官的单元体积-称为体元-内的多个动脉/组织/静脉动态系统中的感兴趣的量的方法，所述方法通过灌注加权成像分析系统的处理单元（4）来实施并且包括从实验灌注数据（15）估计所述感兴趣的量的步骤，其特征在于，所述步骤包括：根据贝叶斯方法，通过以下过程，评估所述感兴趣的量的后验边缘概率分布（56）：

将直接概率分布分配（54）给所述实验灌注数据，并且给定在估计所述体元内的所述动脉/组织/静脉动态系统的所述感兴趣的量的过程中涉及到的参数；

通过引入（51b）与所述动态系统的脉冲响应相关的严格的软信息，将联合先验概率分布分配（53）给所述量。

2.一种用于估计（57）器官的单元体积-称为体元-内的多个动脉/组织/静脉动态系统中的感兴趣的量的方法，所述动态系统是线性非时变的且形式上通过关系式

确定，其中C(t)是在体元内循环的造影剂的浓度，C_a(t)是供给所述体元的动脉内的造影剂的浓度，BF是所述体元内的血流，

代表卷积，且R(t)是所述体元内的通过时间的互补累积密度函数；所述方法通过灌注加权成像分析系统的处理单元（4）来实施且包括从实验灌注数据（15）估计所述感兴趣的量的步骤，其特征在于，所述步骤包括：根据贝叶斯方法，通过以下过程，评估所述感兴趣的量的后验边缘概率分布（56）：

将直接概率分布分配（54）给所述实验灌注数据，并且给定在估计所述体元内的所述动脉/组织/静脉动态系统的所述感兴趣的量的过程中涉及的参数；

通过引入（51b）与所述体元内的互补累积密度函数R(t)相关的严格的软信息，将联合先验概率分布分配（53）给所述量。

3.根据前述权利要求中任一项所述的方法，其特征在于，将所述联合先验概率分布分配给所述量的过程进一步通过引入（51b）与供给所述体元的所述动脉内的所述造影剂的浓度-时间曲线有关的严格的软信息来实现。

4.根据前述权利要求中任一项所述的方法，其特征在于，所述方法通过逐次迭代考虑中的多个体元来实施。

5.根据前述权利要求中任一项所述的方法，其特征在于，所述方法包括用于计算（58）由与所估计的感兴趣的量有关的置信区间表示的补充信息的步骤。

6.根据前述权利要求中任一项所述的方法，其特征在于，所述方法包括用于计算（59）由与所估计的感兴趣的量有关的置信区间的可能性表示的补充信息的步骤。

7.根据权利要求4所述的方法，其特征在于所述方法包括计算（60）由以下结果的充分性量度表示的补充信息的步骤：

将直接概率分布分配（54）给实验灌注数据，并且给定在估计考虑中的体元内的所述动脉/组织/静脉动态系统的所述感兴趣的量的问题中涉及的参数；

将联合先验概率分布分配（53）给所述量。

8.根据权利要求4到7中任一项所述的方法，其特征在于，所述方法包括用于将所估计的感兴趣的量传送（14）给能够将其呈现给用户（6）的人机界面（5）的步骤。

9.根据权利要求5到8中任一项所述的方法，其特征在于，所述方法包括用于将与所估计的感兴趣的量有关的任何补充信息传送给能够将其呈现给用户（6）的人机界面（5）的步骤。

10.根据前述权利要求中任一项所述的方法，其特征在于，所述实验灌注数据包含在实验灌注信号的值的矢量中或包含在从所述矢量到浓度-时间曲线的值的矢量的转换中。

11.一种处理单元（4），其包括用于存储的装置、用于与外界通信的装置以及用于处理的装置，其特征在于：

所述用于通信的装置能够从外界接收实验灌注数据（15）；

所述用于处理的装置适于实施根据权利要求1到10中任一项所述的用于估计器官的单元体积-称为体元-内的多个动脉/组织/静脉动态系统中的感兴趣的量（14）的方法。

12.根据前述权利要求所述的处理单元，其特征在于，所述用于通信的装置以适于人机界面（5）的格式传送（57）所估计的感兴趣的量，所述人机界面（5）能够将所述所估计的感兴趣的量呈现给用户（6）。

13.根据前述权利要求所述的处理单元，其特征在于，所述用于通信的装置以适于人机界面（5）的格式传送（58、59）与所估计的感兴趣的量（14）有关的补充信息，所述人机界面（5）能够将所述补充信息呈现给用户（6）。

14.一种灌注加权成像系统，其包括根据权利要求11到13中的任一项所述的处理单元（4）和人机界面（5），所述人机界面（5）能够根据符合权利要求1到10中任一项所述且通过所述处理单元（4）实施的方法将所估计的量（14）呈现给用户（6）。

Images