CN103091656B - 基于正则化约束多项式拟合磁共振线圈灵敏度的计算方法 - Google Patents

Abstract

本发明涉及一种基于正则化约束多项式拟合磁共振线圈灵敏度的计算方法，其包括以下步骤：A1、对相控阵列线圈Coil_l(l＝1,2,…,L)进行扫描，得到复合图像f_l的求和图像计算灵敏度系数C_l(p)′；A2、对所述灵敏度系数C_l(p)′进行基于Tikhonov正则约束的多项式拟合，生成灵敏度图C_l′；A3、对相控阵列线圈Coil_l进行扫描，得到相应复数图S_l，根据复数图S_l和灵敏度图C_l′，生成整个视野域的复合图像。本发明无需对被检体进行预扫描，对灵敏度图的异常数据进行校正，提高了线圈灵敏度估计的准确性，从而提高图像重建质量，而且该计算方法快速简单，极大地增强了用户的使用效果，提高了工作效率，具有很好的应用价值。

Description

基于正则化约束多项式拟合磁共振线圈灵敏度的计算方法

【技术领域】

本发明涉及MRI快速成像领域，尤其涉及一种基于正则化约束多项式拟合磁共振线圈灵敏度的计算方法。

【背景技术】

关于多通道线圈并行成像的研究开始于1987年，起初敏感度编码被提出来替代所有的梯度相位编码，由于可以获得的有效敏感度编码数是有限的，这个方法的应用受到了一定限制。

1993年，Ra.J.B等人提出更加实际有效的混合编码方案，即同时利用敏感度编码与梯度相位编码。在经历了初始的基本理论研究阶段后，研究者提出了并行成像的实际应用方案，它是利用相控阵线圈空间信息取代梯度编码信息，各个线圈单元同时对K空间进行欠采样，使得每个线圈单元采集的K空间信号数量大大减少，成像速度大大提高。并行成像方法不仅可以提高成像速度，而且还能提高图像的对比度，因此在心脏、动态增强、血管成像等对比成像要求很高的检查中具有广泛的实用价值。由于各个线圈是欠采样采集K空间数据，根据奈奎斯特采样定理，对各个线圈单元采集的数据直接进行图像重建会发生混叠，需要结合线圈的灵敏度矩阵进行重建。当前并行成像的图像重建算法主要有两类，一类是基于图像域的重建，此类代表算法如SENSE（SensitivityEncoding）；另一类是基于K空间域的重建，此类代表算法如SMASH（Simultaneous Acquisition of Spatial Harmonics）。随后出现了将上述方法进一步优化或扩展的研究，如AUTO-SMASH，VD-AUTO-SMASH，GRAPPA，GENERALIZED SMASH，PILS，SPACE-RIP，TSENSE，k-t SENSE，TGRAPPA和k-tGRAPPA等。值得一提的是SENSE和GRAPPA算法已被广泛应用在商业化MRI系统上。

SENSE算法为Pruessmann教授于1999年提出，利用各个线圈采集到的K空间数据，经傅里叶变换（Fast Fourier Transform，FFT）可以重建出一幅FOV（视野）降低含混迭伪影的图像，利用敏感度信息，对重叠图像进行展开和组合，最后得到无重叠完整的图像。SENSE意为灵敏度编码，其简要流程如附图1所示。

SENSE技术使用多个表面接收线圈，将它们以一定结构放置在受检物体周围不同的位置，在施加脉冲序列后各线圈同时采集磁共振信号。每个线圈都有独立的接收通道，如前置放大器、模数转换和快速离散傅里叶变换模块，从而可得到多个磁共振图像。对其中的每个线圈而言，SENSE技术在保持K空间最大值不变的情况下，增加相位编码之间的距离，减少K空间的采样步数。这样K空间采样密度随之也降低，导致傅里叶成像中的视野FOV相应缩小，在相位编码方向出现像素重叠，如附图1中经FFT生成图像的步骤所示。每个接收线圈都生成一幅具有重叠伪影的磁共振图像，称之为中间图像。因此，每个接收线圈有自己独特的空间敏感特征，暗含受检体磁共振信号的位置。在一定的条件下完全可利用该特征进行空间编码，重建原始图像。

采用多线圈技术进行并行采集时，每个相控阵线圈各自独立地采集欠采样的K空间数据。对于二维磁共振成像，假设同时采集的线圈个数为L，不考虑噪声，则每个线圈采集到的信号的离散表达式为：

$s_l(k_x,k_y)=\underset{x=0}{\overset{N_x}{\mathrm{\Σ}}}\underset{y=0}{\overset{N_y}{\mathrm{\Σ}}}{C}_l(x,y)\mathrm{\ρ}(x,y)e^{-i(k_{x}x+k_{y}y)}---\left(1\right)$

式（1）中，(k_x,k_y)为图像空间位置(x,y)所对应的k空间坐标，s_l(k_x,k_y)为第l个相控阵线圈采集的k空间数据，l＝1,2,...,L；C_l(x,y)为第l个相控阵线圈在图像(x,y)处所对应的敏感度系数，ρ为待重建图像，N_x和N_y分别为图像在x与y方向上的采样点数目即重建图像的大小。

SENSE方法就是在图像空间利用多个相控阵线圈所成的图像结合各自特有空间敏感度信息来产生一幅如同传统成像方式一样的全视野，不含混迭伪影的图像。并行重建主要通过这种欠采样方式减少相位梯度编码的数量以相应地缩短扫描时间，达到快速成像。由此可知，为了求得重建图像ρ(x,y)，需要知道每个相控阵线圈的准确敏感度系数。

然而，一般无法预先确切知道每个相控阵线圈的敏感度；并且，因为信噪比的要求，目前主要由预先扫描的不含欠采样伪影的各线圈图像除以各图像平方和的平方根采集得到的图像作为敏感度估计，若将预先扫描的数组线圈图像记为f_l，则有下式：

$C_l(x,y)=f_l(x,y)/\sqrt{\underset{l}{\mathrm{\Σ}}{|f_l(x,y)}^2|}---\left(2\right)$

在使用单个线圈进行传统磁共振成像中，一幅图像某个位置上像素点的值应等于该点的像素值与相应位置上的空间敏感度值的乘积。而在单线圈成像时，线圈各个位置上具有相近的空间敏感度值，所以，所有位置上的空间敏感度信息值可近似为1。而采用多线圈并行重建技术时，在全采样的情况下，一幅图像某个位置上像素点的值应等于该点的像素值与相应线圈位置上的空间敏感度值的乘积，由于所成图像不发生混叠图像，所以与另外线圈的空间敏感度值无关；但在加速采样的情况下，以加速因子R=2为例，在混叠图像上同样位置下的像素点的值应等于两个位置上像素值的叠加。SENSE算法的核心就是利用这个关系把混叠的图像分离出来，最后重建出一幅全视野的，不含混叠伪影的图像。

以上的重建关系可以用下列方程组来表示：

$\left(\begin{array}{c}{S_1}^A(x,y)\\ {S_2}^A(x,y)\\ M\\ {S_L}^A(x,y)\end{array}\right)=C\left(\begin{array}{c}\mathrm{\ρ}(x,y)\\ \mathrm{\ρ}(x,y+k)\\ M\\ \mathrm{\ρ}(x,y+(R-1)k)\end{array}\right)C=\left(\begin{array}{c}{C}_1(x,y)\mathrm{\Λ}{C}_1(x,y+(R-1)k)\\ C_2(x,y)\mathrm{\Λ}{C}_2(x,y+(R-1)k)\\ M\\ C_L(x,y)\mathrm{\Λ}{C}_L(x,y+(R-1)k)\end{array}\right)---\left(3\right)$

其中，L为线圈个数，S_L ^A为单个相控阵线圈所成的混叠图像，ρ为重建图像，C为空间敏感度信息，R为加速因子，k为步长，k＝N_y/R，N_y为重建图像y方向的大小。

这个重建过程也可以用一个线性公式来进行表达：

S^A＝Cρ                                         (4)

通常情况下，在这个线性系统中要求加速因子小于相控阵线圈的个数，即R<L，这样公式(4)中的方程组即为超定方程组，该方程组的解通常为在最小二乘意义下求得的解，那么所求的重建图像为

ρ^LS＝(C^HC)^-1C^HS^A                                   (5)

其中，C^H为矩阵C的赫尔米特（Hermitian Matrix，也称埃尔米特矩阵）转置。

《磁共振医学》1999年第42卷第5期952—962页提出了灵敏度编码的最初成像算法，利用相控阵列线圈Coil_l自校准生成灵敏度图C_l′，即利用复合图像f_l的平方图像f_l ²的求和平方根图像计算出灵敏度图C_l′；

《中国医疗器械杂志》2005年第29卷第3期196—198页提出了灵敏度编码的改进算法，首先利用相控阵列线圈Coil_l自校准生成灵敏度系数C_l(p)′，然后利用最小二乘法对灵敏度系数C_l(p)′进行曲线拟合生成灵敏度图C_l′；

中国专利ZL03138669.5，其名称为线圈灵敏图生成方法和平行成像方法以及磁共振成像装置，提出通过对相控阵列线圈Coil_l（l＝1,2,...,L）进行整个视野域校准扫描得到复合图像f_l的求和图像计算灵敏度因数C_l(p)′。此外，利用加权最小二乘法对绝对值求和图像Am的像素值Am(p)的平方进行曲线拟合，生成灵敏度图C_l′。

常规的灵敏度编码算法，如《磁共振医学》1999年第42卷第5期952—962页提出的灵敏度编码成像算法，只利用相控阵列线圈Coil_l自校准生成灵敏度图C_l′，即利用复合图像f_l的平方图像f_l ²的求和平方根图像计算出灵敏度图C_l′，尽管提供了良好的信噪比，但该方法的均匀性太差，灵敏度图的表面也不平滑，且其需要对被检体进行预扫描，致使扫描时间过长。

现有技术方案中，尚未提出如何消除灵敏度异常值对重建图像造成的伪影。如专利号为ZL03138669.5的中国专利和《中国医疗器械杂志》2005年第29卷第3期196—198页提出的灵敏度编码改进算法。提出使用最小二乘法对灵敏度系数进行二维曲线拟合，但该方法没有解决最小二乘法拟合过程中对噪声比较敏感，拟合中使用的正规方程组病态问题，致使灵敏度数据发生异常，从而对最终重建图像产生很大伪影的问题。

或者说，磁共振成像具有较高的软组织对比度与空间分辨率，并能根据需要灵活选择成像参数与成像层面，已经广泛应用于临床。然而，由于数据采集时间较长而导致成像速度较慢，是磁共振成像的主要不足之处。随着梯度场强已经接近极限，扫描速度也受着严格的制约，现有技术提出了多通道采集与并行成像算法，使得成像可以不再依赖梯度性能的提高就可以大大加快数据采集。多通道并行成像技术主要是利用相控阵线圈中单个相控阵线圈的空间敏感度差异来编码空间信息，降低成像所必需的梯度编码数，从而获得了更快的扫描速度。

但是，在数据采集重建过程中，运动常常会使线圈数据发生异常，一般无法预先确切知道每个相控阵线圈的敏感度，从而对最终重建图像质量产生很大影响。

因此，现有技术需要改进。

【发明内容】

有鉴于此，有必要提出一种新型的基于正则化约束多项式拟合磁共振线圈灵敏度的计算方法，即基于正则化约束的采用多项式拟合的提升磁共振线圈灵敏度的计算方法。

本发明的技术方案如下：一种基于正则化约束多项式拟合磁共振线圈灵敏度的计算方法，其包括以下步骤：A1、对相控阵列线圈Coil_l（l＝1,2,...,L）进行扫描，得到复合图像f_l的求和图像计算灵敏度系数C_l(p)′；A2、对所述灵敏度系数C_l(p)′进行基于Tikhonov正则约束的多项式拟合，生成灵敏度图C_l′；A3、对相控阵列线圈Coil_l进行扫描，得到相应的复数图像S_l，根据复数图像S_l和灵敏度图C_l′，生成整个视野域的复合图像。

优选的，所述计算方法具体执行以下步骤：B1、开始并行成像处理；B2、进行扫描，从而得到与相控阵列线圈Coil_l相应的复数图像f_l(l=1,2,…,L)；B3、生成每个复数图像f_l的平方图像f_l ²，并且生成平方图像f_l ²的求和平方根图像As，B4、计算相控阵列线圈Coil_l中的像素p的灵敏度系数C_l(p)′，其中，f_l(p)为复数图像中的像素p的像素值，As(p)为绝对值求和图像As中的像素p的像素值；B5、基于多项式拟合及Tikhonov正则约束根据灵敏度系数C_l(p)′生成灵敏度图C_l′；B6、对相控阵列线圈Coil_l进行扫描，得到相应的复数图像S_l（l=1,2,…,L）；B7、根据灵敏度图C_l′和复数图像S_l生成整个视野域的复数图像。

优选的，所述计算方法中，所述Tikhonov正则约束采用以下约束计算公式：arg minλ||X||₂使得||AX-b||₂≤ε；其中，λ为正则约束项参数，ε为噪声因子，||·||₂分别表示矩阵范数为2。

优选的，所述计算方法中，采用岭迹法确定所述正则约束项参数λ。

优选的，所述计算方法中，采用非线性的最速下降法、共轭梯度法或者高斯牛顿法计算所述约束计算公式。

优选的，所述计算方法中，所述矩阵范数为1。

优选的，所述计算方法中，所述多项式拟合采用以下步骤：C1、由已知数据画出函数的散点图，确定拟合多项式的次数n；C2、列表计算(j＝0,1,...,2n)和(j＝0,1,...,2n)；C3、列出法方程组，计算a₀,a₁,...a_n；C4、得到拟合多项式

优选的，所述计算方法中，所述对相控阵列线圈Coil_l进行扫描，是对其进行具有较小视野域的非完整视野域扫描。

优选的，所述计算方法中，步骤A3之后，还执行以下步骤A4：输出所述复合图像到本地显示终端或远程显示终端。

优选的，所述计算方法中，所述对相控阵列线圈Coil_l（l＝1,2,...,L）进行扫描，其中，自动校准信号行数为32，加速因子为4。

基于以上方案，本发明提出了一种基于Tikhonov正则化约束的多项式拟合线圈灵敏度计算方法，该方法无需对被检体进行预扫描，在对线圈灵敏度进行多项式拟合过程中使用Tikhonov正则化约束，对灵敏度图的异常数据进行校正，提高了线圈灵敏度估计的准确性，从而提高图像重建质量，而且该计算方法快速简单，极大地增强了用户的使用效果，提高了工作效率，具有很好的应用价值。

【附图说明】

图1是现有技术的SENSE算法示意图；

图2是一个实施例的并行成像算法的流程图；

图3a是又一个实施例的真实脑部数据重建结果图的参考图像；

图3b是图3a的实施例采用本发明的重建方法得到的图像；

图3c是图3a的实施例采用SENSE传统估计灵敏度方法得到的图像。

【具体实施方式】

下面结合附图，对本发明的具体实施方式进行详细描述。

Tikhonov正则化方法的核心思想是采用加和的方式，在最小均方误差代价函数的基础上直接集成一个关于原始图像x′的先验概率函数；图像重建问题常采用正则化方法，根据对于解的先验知识来构造附加约束条件或改变求解策略，使问题的解变得确定和稳定，例如，通常以图像的平滑性作为解的约束条件，实践证明，这是一种对于不适定问题求解的有效方法。本发明提出了一种线圈灵敏度多项式拟合过程中应用Tikhonov正则约束的并行重建算法，一个实施例如下：一种基于正则化约束多项式拟合磁共振线圈灵敏度的计算方法，其包括以下步骤：

A1、对相控阵列线圈Coil_l（l＝1,2,...,L）进行扫描，得到复合图像f_l的求和图像然后计算灵敏度系数C_l(p)′；例如，步骤A1包括以下步骤：A11、开始并行成像处理；A12、进行扫描，从而得到与相控阵列线圈Coil_l相应的复数图像f_l(l=1,2,…,L)；A13、生成每个复数图像f_l的平方图像f_l ²，并且生成平方图像f_l ²的求和平方根图像As，A14、计算相控阵列线圈Coil_l中的像素p的灵敏度系数C_l(p)′，其中，f_l(p)为复数图像中的像素p的像素值，As(p)为绝对值求和图像As中的像素p的像素值。

A2、通过对所述灵敏度系数C_l(p)′进行基于Tikhonov正则约束的多项式拟合，生成灵敏度图C_l′；也就是说，基于多项式拟合及Tikhonov正则约束根据灵敏度系数C_l(p)′生成灵敏度图C_l′。优选的，所述Tikhonov正则约束采用以下约束计算公式：arg minλ||X||₂使得||AX-b||₂≤ε；其中，λ为正则约束项参数，ε为噪声因子，||·||₂分别表示矩阵范数为2，或者，所述矩阵范数为1，将基于Tikhonov正则化的||·||₂约束替换成||·||₁稀疏约束，即约束计算公式为arg minλ||X||₁使得||AX-b||₁≤ε。优选的，采用岭迹法确定所述正则约束项参数λ。优选的，采用非线性的最速下降法、共轭梯度法或者高斯牛顿法计算所述约束计算公式。

优选的，所述多项式拟合采用以下步骤：C1、由已知数据画出函数的散点图，确定拟合多项式的次数n；C2、列表计算和C3、列出法方程组，计算a₀,a₁,...a_n；C4、得到拟合多项式其中，多项式拟合方法说明如下。

假设给定数据点(x_i,y_i)(i=0,1,…,m)，Φ为所有次数不超过n(n≤m)的多项式构成的函数类，现求一使得

$I=\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{[p_n\left(x_i\right)-y_i]}^2=\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{(\underset{k=0}{\overset{n}{\mathrm{\&Sigma;}}}{a}_k{x}_i^k-y_i)}^2=\min ---\left(6\right)$

当拟合函数为多项式时，称为多项式拟合，满足式（6）的p_n(x)称为最小二乘拟合多项式。特别地，当n=1时，称为线性拟合或直线拟合。显然，为a₀,a₁,...a_n的多元函数，因此上述问题即为求I＝I(a₀,a₁,...a_n)的极值问题。由多元函数求极值的必要条件，得

$\frac{\&PartialD;I}{{\&PartialD;a}_j}=2\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}(\underset{k=0}{\overset{n}{\mathrm{\&Sigma;}}}{a}_k{x}_i^k-y_i)x_i^j=0,j=\mathrm{0,1},...,n---\left(7\right)$

即，

$\underset{k=0}{\overset{n}{\mathrm{\&Sigma;}}}\left(\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^{j+k}\right)a_k=\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^j{y}_i,j=\mathrm{0,1},...,n---\left(8\right)$

其中，式（8）是关于a₀,a₁,...a_n的线性方程组，用矩阵表示为

$\left[\begin{array}{cccc}m+1& \underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i& ...& \underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^n\\ \underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i& \underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^2& ...& \underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^{n+1}\\ .& .& & .\\ .& .& & .\\ .& .& & .\\ \underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^n& \underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^{n+1}& ...& \underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^{2n}\end{array}\right]\left[\begin{array}{c}{a}_0\\ a_1\\ .\\ .\\ .\\ a_n\end{array}\right]=\begin{array}{c}\left[\begin{array}{c}\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{y}_i\\ \underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i{y}_i\\ .\\ .\\ .\\ \underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^n{y}_i\end{array}\right]---\left(9\right)\end{array}$

式（8）或式（9）称为正规方程组或法方程组。从式（4）中解出a_k(k=0,1,…，n)，从而可得多项式如下：

$p_n\left(x\right)=\underset{k=0}{\overset{n}{\mathrm{\&Sigma;}}}{a}_k{x}^k---\left(10\right)$

可证明，式（10）中的p_n(x)满足式（6），即p_n(x)为所求的拟合多项式。

因此，把称为最小二乘拟合多项式p_n(x)的平方误差，记作：

${\left|\right|r\left|\right|}_2^2=\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{[p_n\left(x_i\right)-y_i]}^2$

由式(7)可得

${\left|\right|r\left|\right|}_2^2=\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{y}_i^2-\underset{k=0}{\overset{n}{\mathrm{\&Sigma;}}}{a}_k\left(\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^k{y}_i\right)---\left(11\right)$

本发明的多项式拟合的一般方法包括以下几步：

(1)由已知数据画出函数粗略的图形——散点图，确定拟合多项式的次数n；

(2)列表计算(j＝0,1,...,2n)和(j＝0,1,...,2n);

(3)写出正规方程组，求出a₀,a₁,...a_n；

(4)写出拟合多项式

并且，Tikhonov正则化约束说明如下。

从式（9）可以看出，a₀,a₁,...a_n的求解可以归结为一个Ax=b逆问题。对于低秩的SENSE敏感度拟合过程中正规方程组A为病态矩阵，所以使得拟合中的微小扰动很容易使系统解产生巨大改变。拟合多项式往往遭到噪声的破坏而发生异常，此时用最小二乘估计得到解是基于残差最小下的解。

但是由于灵敏度图破坏数据和异常值的存在往往会造成大残差的情况，而最小二乘解通常对具有大残差的数据点十分敏感，这些异常值能够破坏最小二乘解。可见，基于最小二乘意义下求得的解并不具有鲁棒性（Robustness，也称健壮性）。

针对上述问题，提出了一种具有鲁棒性的Tikhonov正则化约束的重建算法。通过鲁棒性约束对数据进行校正，使校正后的数据具有较高的可靠性，从而降低异常值对方程解的敏感性，即

arg minλ||X||₂subject to||AX-b||₂≤ε

(12)

其中λ为正则约束项参数，ε为噪声因子，||·||₂分别表示矩阵的2范数。在本发明中，为了便于操作，同时考虑准确性，例如，采用岭迹法（Ridge Trace）确定正则约束项参数λ。求解方程式（12）有很多种非线性迭代方法，比如最速下降法，共轭梯度法，高斯牛顿法等。对于本发明各各实施例，在任一个迭代搜索方法中，最关键的是在每一步迭代中找到搜索方向，共轭梯度法是一种线性搜索方法，是目前求解线性方程中最重要的一个方法，它可以很快的线性收敛到局部最优解。因此，本发明各实施例在求解非线性最优化问题时，优选的，采用共轭梯度法来求解。

A3、通过对相控阵列线圈Coil_l进行扫描，得到相应的复数图像S_l，根据复数图像S_l和灵敏度图C_l′，生成整个视野域的无伪影复合图像。例如，步骤A3包括以下步骤：A31、对相控阵列线圈Coil_l进行扫描，得到相应的复数图像S_l（l=1,2,…,L）；A32、根据灵敏度图C_l′和复数图像S_l生成整个视野域的复数图像。优选的，步骤A3之后，还执行以下步骤A4：输出所述复合图像到本地显示终端或远程显示终端。这样，就可以直观地、迅速地获知计算结果。

例如，上述各实施例中，所述对相控阵列线圈Coil_l进行扫描，是对其进行具有较小视野域的非完整视野域扫描；又如，对相控阵列线圈Coil_l进行扫描，是相控阵列线圈Coil_l中的每个线圈，对其扫描范围进行扫描，其视野域不是完整的整个视野域，仅是其中的一部分，各线圈的视野域叠加形成整个视野域。例如，上述各实施例中，对相控阵列线圈Coil_l（l＝1,2,...,L）进行扫描，其中，自动校准信号行数为32，加速因子为4。或者，自动校准信号行数为16，加速因子为3。或者，自动校准信号行数为24，加速因子为6。

优选的，一个实施例的流程图如图2所示，示出了线圈灵敏度多项式拟合过程中应用Tikhonov正则约束的并行重建算法的详细流程图，为整个发明的技术路线，具体说明如下。

B1、开始并行成像处理；

B2、进行扫描，从而得到与相控阵列线圈Coil_l相应的复数图像f_l(l=1,2,…,L)；

B3、生成每个复数图像f_l的平方图像f_l ²，并且生成平方图像f_l ²的求和平方根图像As:;

B4、计算相控阵列线圈Coil_l中的像素p的灵敏度系数C_l(p)′，其中，f_l(p)为复数图像中的像素p的像素值，As(p)为绝对值求和图像As中的像素p的像素值；

B5、基于多项式拟合及Tikhonov约束根据灵敏度系数C_l(p)′生成灵敏度图C_l′；

B6、进行具有较小视野域（FOV）的扫描，从而得到相控阵列线圈Coil_l相应的复数图像S_l（l=1,2,…,L）；

B7、根据灵敏度图C_l′和复数图像S_l生成整个视野域的复数图像；

B8、结束。

为了体现本算法的优越性给出了自动校准信号(Auto-calibration Signal，ACS)的ACS行数=32，加速因子R=4时，应用本发明方法重建结果。通常可以理解为，加速因子R越大，图像的信噪比越差；ACS行数越多，图像的质量越好。图3a为模板参考图像；图3b给出了本发明的重建图像；图3c传统估计灵敏度方法。通过对比各图，同样可以看出，本发明的计算方法无论在伪影的消除还是边缘的细节分辨都明有明显优势。

本发明的关键点是一种相控阵列线圈灵敏度的计算方法，该方法包括：通过对相控阵列线圈Coil_l（l＝1,2,...,L）进行扫描得到的复合图像f_l的求和图像计算灵敏度因数C_l(p)′，然后通过多项式拟合线圈灵敏度过程中应用具有鲁棒性的Tikhonov正则约束生成灵敏度图C_l′；这样，本发明提出了一种基于Tikhonov正则约束及多项式拟合的线圈灵敏度生成方法，该算法无需对被检体进行预扫描，应用多项式拟合及Tikhonov正则化对灵敏度图的异常数据进行拟合校正，即应用Tikhonov正则化约束的多项式拟合方法对灵敏度图的异常数据进行校正，使校正后的数据具有较高的可靠性，提高了接受线圈空间敏感性值的精确度，提高线圈灵敏度估计的准确性，从而有效的消除破坏数据在重建图像中产生的伪影，提高重建图像的分辨率，从而提高重建图像质量。

需要补充说明的是，本发明的计算方法，其直接目的不是获得诊断结果或健康状况，而只是对扫描图像进行重建和处理，以获取作为中间结果的信息，以及计算该信息的方法，从而有效地消除破坏数据在重建图像中产生的伪影，提高重建图像的分辨率。根据目前的医学知识和本发明所公开的内容，从所获得的信息本身不能够直接得出疾病的诊断结果。

以上所述实施例仅表达了本发明的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对本发明专利范围的限制；并且，上面列出的各个技术特征，其相互组合所能够形成各个实施方案，应被视为属于本发明说明书记载的范围。对于本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进，这些都属于本发明的保护范围。因此，本发明专利的保护范围应以所附权利要求为准。

Claims (10)

1.一种基于正则化约束多项式拟合磁共振线圈灵敏度的计算方法，其特征在于，包括以下步骤：

A1、对相控阵列线圈Coil_l(l＝1,2,…，L)进行扫描，得到复合图像f_l的求和图像计算灵敏度系数C_l(p)′；

A2、对所述灵敏度系数C_l(p)′进行基于Tikhonov正则约束的多项式拟合，生成灵敏度图C_l′；

A3、对相控阵列线圈Coil_l进行扫描，得到相应的复数图像S_l，根据复数图像S_l和灵敏度图C_l′，生成整个视野域的复合图像。

2.根据权利要求1所述计算方法，其特征在于，具体执行以下步骤：

B1、开始并行成像处理；

B2、进行扫描，从而得到与相控阵列线圈Coil_l相应的复数图像f_l(l＝1,2,…,L)；

B3、生成每个复数图像f_l的平方图像f_l ²，并且生成平方图像f_l ²的求和平方根图像As， $\mathrm{As}=\sqrt{\underset{l=1}{\overset{L}{\mathrm{\&Sigma;}}}{\left|f_l\right|}^2};$

B4、计算相控阵列线圈Coil_l中的像素p的灵敏度系数C_l(p)′，其中，f_l(p)为复数图像中的像素p的像素值，As(p)为绝对值求和图像As中的像素p的像素值；

B5、基于多项式拟合及Tikhonov正则约束根据灵敏度系数C_l(p)′生成灵敏度图C_l′；

B6、对相控阵列线圈Coil_l进行扫描，得到相应的复数图像S_l(l＝1,2,…,L)；

B7、根据灵敏度图C_l′和复数图像S_l生成整个视野域的复合图像。

3.根据权利要求1或2所述计算方法，其特征在于，所述Tikhonov正则约束采用以下约束计算公式：

arg minλ||X||₂使得||AX-b||₂≤ε

其中，λ为正则约束项参数，ε为噪声因子，||·||₂分别表示矩阵范数为2。

4.根据权利要求3所述计算方法，其特征在于，采用岭迹法确定所述正则约束项参数λ。

5.根据权利要求3所述计算方法，其特征在于，采用非线性的最速下降法、共轭梯度法或者高斯牛顿法计算所述约束计算公式。

6.根据权利要求3所述计算方法，其特征在于，所述矩阵范数为1，将基于Tikhonov正则化的||·||₂约束替换成||·||₁稀疏约束，即约束计算公式为argminλ||X||₁使得||AX-b||₁≤ε。

7.根据权利要求3所述计算方法，其特征在于，所述多项式拟合采用以下步骤：

C1、由已知数据画出函数的散点图，确定拟合多项式的次数n；

C2、列表计算 $\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^j,(j=\mathrm{0,1,2}n)$ 和 $\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}{x}_i^j{y}_i,(j=\mathrm{0,1,}\&CenterDot;\&CenterDot;\&CenterDot;,2n);$

C3、列出法方程组，计算a₀,a₁,…a_n；

C4、得到拟合多项式

8.根据权利要求1所述计算方法，其特征在于，所述对相控阵列线圈Coil_l进行扫描，是对其进行具有较小视野域的非完整视野域扫描。

9.根据权利要求1所述计算方法，其特征在于，步骤A3之后，还执行以下步骤A4：输出所述复合图像到本地显示终端或远程显示终端。

10.根据权利要求1所述计算方法，其特征在于，所述对相控阵列线圈Coil_l(l＝1,2,…,L)进行扫描，其中，自动校准信号行数为32，加速因子为4。