CN101305908A - 并行采集图像重建的信噪比计算方法 - Google Patents

Abstract

本发明公开了一种并行采集图像重建的信噪比计算方法，其包括：确定并行采集图像重建的重建表达式，使图像重建满足线性运算ρ＝Mr，其中ρ为重建图像，M为重建表达式，r为原始数据；根据重建表达式确定图像重建的加权系数；根据加权系数和原始数据计算并行采集图像的信噪比(SNR)。并行采集图像的SNR，不仅与线圈的几何形态和位置有关，还受重建方法及采样模式的影响。在本发明中，重建图像SNR的计算考虑到所有采集通道上所有相位编码线上读方向对应位置的原始数据的贡献，因此能够更准确地反映并行采集图像重建的SNR损失情况，尤其是能够准确地反映重建中结合的参考线的多少所导致的SNR的变化情况。使得对于并行采集图像重建的SNR损失的评估更加准确。

Description

并行采集图像重建的信噪比计算方法

技术领域

本发明涉及磁共振成像(MRI，Magnetic Resonance Imaging)技术，特别是涉及一种并行采集(PA，Parallel Acquisition)图像重建的信噪比(SNR，Signal to Noise Ratio)计算方法。

背景技术

在MRI技术中，成像的速度是非常重要的参数。早期的检查常常需要花费若干个小时，随后由于在场强、梯度硬件以及脉冲序列方面技术的改进，成像的速度有了一个较大提升。但是，场梯度快速变换以及高密度连续射频(RF，Radio Frequency)脉冲会带来人体生理极限无法承受的特殊吸收率(SAR，Specific Absorption Rate)以及器官组织的致热量，因此，成像速度的提升遇到了瓶颈。

随后，研究人员发现，借助应用复杂的电脑图像重建算法以及相配合的线圈阵列，磁共振的成像速度可以被大大的提升，该种技术通常被称为并行成像(parallel imaging)技术。并行采集图像重建是一种用于并行采集的图像重建(reconstruction)技术，它利用相控阵线圈的空间灵敏度(sensitivity)差异进行空间编码，并用相控阵线圈同时采集，获得比常规磁共振成像快2-6倍甚至更高的成像速度。采用并行成像技术，对磁共振成像系统提出新的要求，如需要多个接收通道、多元阵列线圈及线圈灵敏度校准、用特殊的数据处理和图像重建方法等等。

并行成像可以提高图像采集速度，成像速度的提高是通过降低K空间填充率实现的。如果K空间填充率低于采样定理的限制会导致直接利用傅立叶重建后的图像出现伪影。一般的磁共振成像是通过采集物体的频域信息通过傅立叶变换得到图像的。根据采样定理，在图像域物体的重复周期反比于频域的采样间隔。如果图像的空间重复周期小于图像本身的大小，重建后的图像就会重叠在一起，这种现象在信号处理中称为混叠。

对于多线圈采集，虽然每个线圈采集到的K空间信息是不够的，但可以利用不同线圈采集到信号的差异经过处理得到一幅完整的图像。并行成像去混叠伪影的重建算法大致上可以分为两类：空间谐波并行采集(SMASH，SiMultaneous Acquisition of Spatial Harmonics)和灵敏度编码并行采集技术(SENSE，SENSitivity Endoding parallel acquisitiontechnique)。其中，SMASH方法是一种利用各个通道线圈的灵敏度函数组成空间谐波并进行辅助编码的方法。一般线圈的灵敏度函数在空间是缓慢变化的，并可近似认为是高斯分布函数，于是可以利用各通道线圈的灵敏度函数的线性组合构成一定频率的空间谐波。并利用空间谐波函数凑出没有实际采集的相位编码线。Griswold等人对SMASH做了进一步的改进，提出了广义自校准行采集(GRAPPA，Generalized Autocalibrating Partially ParallelAcquisitions)算法。与SMASH拟和最终图像满采样时k空间欠采的数据不同，GRAPPA拟和每个通道欠采的k空间数据。GRAPPA提高了重建图像的精度，但需要更多的处理时间。目前已经提出了在GRAPPA算法中考虑参考线进行重建的方案，即在拟合欠采的K空间数据时，将低频满采样的数据替换对应的拟合出来的数据，以提高重建图像的SNR。

与SMASH在频域处理的方案不同，SENSE方法通过在图像域求解线性方程组的方法来消除欠采样引起的伪影。由于图像空间周期引起的重叠效应。SENSE方法采集数据时在k空间的中心区域满采样，在外围区域欠采样。于是原始的K空间数据被分为两部分：均匀欠采样的数据和低频满采样的数据。均匀欠采样的数据用来生成混叠的图像，而低频满采样的数据用来生成模糊的组织图像进而得到实时的线圈灵敏度分布和加权矩阵，最后将利用均匀欠采样的数据生成的混叠图像与利用低频满采样数据得到的加权矩阵合成，得到高分辨率的无混叠图像。这里，用于得到线圈灵敏度分布和加权矩阵的低频满采样数据称为参考数据，K空间中低频满采样的相位编码线称为参考线。在本申请人的中国专利“并行采集图像重建的方法和装置”中，描述了一种在SENSE算法中考虑参考线进行重建的方案，即按照混合采样模式生成包含K空间欠采样数据和低频满采样数据相结合的重建数据，并在图像重建时考虑混合采样模式。

在本申请人的申请号为200410082376.8的中国专利“磁共振成像快速广义自校准并行采集图像重建算法”中描述了一种广义的SNR损失评估算法，这种算法既可应用于图像空间，又可应用于K空间并行成像方法中。然而上述SNR损失评估算法并没有考虑参考线对于SNR的贡献。

发明内容

本发明的一个目的在于提出一种并行采集图像重建的SNR计算方法，在计算并行采集图像重建的SNR时考虑参考线的贡献。

为了实现上述目的，本发明的技术方案包括：

一种并行采集图像重建的信噪比计算方法，包括：确定并行采集图像重建的重建表达式，使图像重建满足线性运算ρ＝Mr，其中ρ为重建图像，M为重建表达式，r为原始数据；根据所述的重建表达式确定图像重建的加权系数；根据所述加权系数和所述原始数据计算并行采集图像的信噪比。

所述方法还包括：根据并行采集图像的信噪比与标准采集图像的信噪比的比值作为相对信噪比评估并行采集图像重建的信噪比损失。

其中，所述重建表达式为重建矩阵。

其中，所述原始数据为对K空间均匀欠采样获得的数据。

其中，所述原始数据为对K空间均匀欠采样和低频满采样相结合的混合采样获得的数据。

其中，对于基于K空间的重建方式，所述确定并行采集图像重建的重建表达式包括：确定并行采集图像重建的重建矩阵：M＝C·F·G，其中G表示利用拟合系数对均匀欠采样的K空间位置进行填充同时保留低频满采样的数据的k空间数据填充矩阵；C表示利用线圈灵敏度的复数共轭进行等效于平方和合并SOS的通道合并运算矩阵；F表示傅立叶变换矩阵。

其中，对于基于图像空间的重建方式，所述确定并行采集图像的重建表达式包括：确定并行成像采样矩阵：E＝U·F^-1·P，其中E表示并行成像采样矩阵，U表示欠采样矩阵，F表示傅立叶变换矩阵，P表示线圈灵敏度分布矩阵；对所述并行成像采样矩阵求逆得到并行采集图像重建的重建矩阵：M＝inv(E′E)·E′，其中M表示重建矩阵，E表示所述并行成像采样矩阵。

其中，所述欠采样矩阵U表示均匀欠采样矩阵或者对K空间进行均匀欠采样和低频满采样的混合模式欠采样矩阵。

其中，对于直角坐标系下的采样，所述原始数据r包括不同采集通道、不同相位编码线在读出方向上同一位置的所有数据。

其中，对于非直角坐标系下的采样，所述原始数据包括所有采集通道的所有k空间数据。

其中，对于直角坐标系下的采样，所述根据重建表达式确定图像重建的加权系数包括：若读出方向为水平方向，相位编码方向为竖直方向，对于重建图像中的水平方向位置为jRO，竖直方向位置为iPE的像素点(jRO，iPE)，确定其加权系数为： $w_{\mathrm{jRO}}^{\mathrm{iPE}}=(w_1,w_2,\·\·\·,w_{\mathrm{nCh}\×\mathrm{nPE}}),$ 相应的重建数学表达式为： ${\mathrm{Ima}}_{\mathrm{jRO}}^{\mathrm{iPE}}=w_{\mathrm{jRO}}^{\mathrm{iPE}}\·r_{\mathrm{jRO}},$ 其中，nCh表示采集通道数，nPE表示相位编码行数，w_jRO ^iPE表示计算该像素点的加权系数行向量，(w₁，w₂，…，w_nCh×nPE)为重建矩阵M_jRO中的第iPE行的所有数值，r_iRO为K空间读出方向上位置同为jRO的不同通道，不同相位编码位置的原始数据组成的列向量，重建矩阵M_jRO为对r_jRO进行重建的重建矩阵，Ima_jRO ^iPE为重建图像中位置为(jRO，iPE)的像素点。

其中，所述根据加权系数和原始数据计算并行采集图像的信噪比包括： ${\mathrm{SNR}}_{\mathrm{jRO}}^{\mathrm{iPE}}=\frac{|w_{\mathrm{jRO}}^{\mathrm{iPE}}\·r_{\mathrm{jRO}}|}{\sqrt{w_{\mathrm{jRO}}^{\mathrm{iPE}}\·\mathrm{\ψ}\·{\left(w_{\mathrm{jRO}}^{\mathrm{iPE}}\right)}^{\′}}\sqrt{\mathrm{mPE}}},$ 其中SNR_jRO ^iPE表示点(jRO，iPE)并行采集图像的信噪比，w_jRO ^iPE表示图像重建的加权系数，mPE为实际采集的相位编码线数，ψ表示噪声相关矩阵。

由上述技术方案可以看出，在本发明的SNR计算方法中，将图像重建表达为线性运算ρ＝Mr，根据重建矩阵M确定图像重建的加权系数，然后据图像重建的加权系数和原始数据计算并行采集图像的SNR。因此，可以实现在SNR计算中，不单考虑线圈的几何形态，还充分考虑了原始数据和采样模式，由于图像中每一像素点的信号的计算考虑到所有采集通道上所有相位编码线上对应位置的原始数据的贡献，能够更准确地反映并行采集图像重建的SNR损失情况，尤其是能够准确地反映重建中结合的参考线的多少所导致的SNR的变化情况。使得对于并行采集图像重建的SNR损失的评估更加准确。

附图说明

下面将通过参照附图详细描述本发明的优选实施例，使本领域的普通技术人员更清楚本发明的上述及其它特征和优点，相同的标号表示相同的部件，附图中：

图1是根据本发明实施例一的并行采集图像重建的SNR计算方法的流程图；

图2是根据本发明实施例一的并行采集图像重建的SNR计算方法的流程图；

图3是采用本发明的方法计算的SNR随着参考线的变化而变化的示意图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

根据本发明，在计算并行采集图像重建的SNR时，确定并行采集图像重建的重建表达式，根据重建表达式确定图像重建的加权系数，然后根据图像重建的加权系数和原始数据计算并行采集图像的SNR。

下面通过两个具体实施例对本发明进行详细阐述。

图1是根据本发明实施例一的并行采集图像重建的SNR计算方法的流程图。本实施例应用于基于K空间的并行采集图像重建，例如GRAPPA重建方法。从图1可以看出，本实施例主要包括以下步骤：

步骤101：确定基于K空间的并行采集图像重建的重建矩阵，使图像重建满足线性运算ρ＝Mr，其中ρ为重建图像，M为重建矩阵，r为原始数据。

在直角坐标系采样时，对于直角坐标系下的采样，原始数据r包括不同采集通道、不同相位编码线在读出方向上同一位置的所有数据。一般来说，读出方向为水平方向，相位编码方向为竖直方向。在本实施例中，ρ是[nPE，1]的向量，表示重建图像中第jRO列；r为[nPE×nCh，1]的向量，表示从所有采集通道1至nCh上读出的第jRO列上的所有原始数据。

原始数据r可以是对K空间均匀欠采样获得的数据，也可以是对K空间均匀欠采样和低频满采样相结合的混合采样获得的数据。

在本实施例中，按照公式(1)确定并行采集图像重建的重建矩阵：

M＝C·F·G    (1)

其中，G表示利用拟合系数对均匀欠采样的K空间位置进行填充同时保留低频满采样的数据的k空间数据填充矩阵；C表示利用线圈灵敏度的复数共轭进行等效于平方和合并SOS的通道合并运算矩阵；F表示傅立叶变换矩阵。在本实施例中，F和G均为[nPE×nCh，nPE×nCh]的矩阵，C为[nPE，nPE×nCh]的矩阵。

步骤102：根据重建矩阵确定图像重建的加权系数。

设读出方向为水平方向，相位编码方向为竖直方向，对于重建图像中的水平方向位置为jRO，竖直方向位置为iPE的像素点(jRO，iPE)，确定其加权系数为：

$w_{\mathrm{jRO}}^{\mathrm{iPE}}=(w_1,w_2,\·\·\·,w_{\mathrm{nCh}\×\mathrm{nPE}})---\left(2\right)$

相应的重建数学表达式为：

${\mathrm{Ima}}_{\mathrm{jRO}}^{\mathrm{iPE}}=w_{\mathrm{jRO}}^{\mathrm{iPE}}\·r_{\mathrm{jRO}}---\left(3\right)$

其中，nCh表示采集通道数，nPE表示相位编码行数，w_jRO ^iPE表示计算该像素点的加权系数行向量，(w₁，w₂，…，w_nCh×nPE)为重建矩阵M_jRO中的第iPE行的所有数值，r_jRO为K空间读出方向上位置同为jRO的不同通道，不同相位编码位置的原始数据组成的列向量，重建矩阵M_jRO为对r_jRO进行重建的重建矩阵，Ima_jROiPE为重建图像中位置为(jRO，iPE)的像素点。

步骤103：根据确定的图像重建的加权系数和原始数据计算并行采集图像的SNR。

按照公式(4)计算并行采集图像重建的SNR：

${\mathrm{SNR}}_{\mathrm{jRO}}^{\mathrm{iPE}}=\frac{|w_{\mathrm{jRO}}^{\mathrm{iPE}}\·r_{\mathrm{jRO}}|}{\sqrt{w_{\mathrm{jRO}}^{\mathrm{iPE}}\·\mathrm{\ψ}\·{\left(w_{\mathrm{jRO}}^{\mathrm{iPE}}\right)}^{\′}}\sqrt{\mathrm{mPE}}}---\left(4\right)$

其中，SNR_jRO ^iPE表示点(jRO，iPE)并行采集图像的SNR，w_jRO ^iPE表示图像重建的加权系数，mPE为实际采集的相位编码线数，ψ表示噪声相关矩阵。

步骤104：将并行采集图像的SNR与标准采集图像的SNR的比值作为相对SNR评估并行采集图像重建的SNR损失。

利用公式(5)计算并行采集图像的SNR与标准采集图像的SNR的比值：

${\mathrm{SNR}}^{\mathrm{rel}}=\frac{{\mathrm{SNR}}_{\mathrm{jRO}}^{\mathrm{iPE}}}{{\mathrm{SNR}}^{\mathrm{opt}}}---\left(5\right)$

其中SNR^opt表示最大SNR，SNR^rel表示相对SNR，即相对于测量时间的SNR效率。

由于步骤103计算出了重建图像中每个象素点的SNR，由此生成SNR损失图像，该SNR损失图像准确地反映了并行采集图像重建所带来的SNR，利用SNR损失图像，可以对并行采集图像重建的SNR进行更准确的估计。

图2是根据本发明实施例二的并行采集图像重建的SNR计算方法的流程图。本实施例应用于基于图像空间的并行采集图像重建，例如SENSE重建方法。从图2可以看出，本实施例主要包括以下步骤：

步骤201：确定并行采集图像重建的重建矩阵，使图像重建满足线性运算ρ＝Mr，其中ρ为重建图像，M为重建矩阵，r为原始数据。

在本实施例中，按照公式(6)确定并行成像采样矩阵：

E＝U·F^-1·P    (6)

其中，E表示并行成像采样矩阵，U表示欠采样矩阵，F表示傅立叶变换矩阵，P表示线圈灵敏度分布矩阵。欠采样矩阵U即可以是均匀欠采样矩阵，也可以是对K空间进行均匀欠采样和低频满采样的混合模式欠采样矩阵。

按照公式(7)对并行成像采样矩阵求逆得到并行采集图像重建的重建矩阵：

M＝inv(E′·E)·E′(7)

其中，M表示重建矩阵，E表示所述并行成像采样矩阵。

步骤202至204与步骤102至104完全相同，在此不予赘述。

以上两个实施例中都公开的是直角坐标系下的采样情况，对于非直角坐标系下的采样，原始数据r包括所有采集通道的所有k空间数据。

由以上两个实施例可以看出，在本发明的SNR计算方法中，将图像重建表达为线性运算ρ＝Mr，根据重建矩阵M确定图像重建的加权系数，然后据图像重建的加权系数和原始数据计算并行采集图像的SNR。并行采集图像的SNR，不仅与线圈的几何形态和位置有关，还受重建方法及采样模式的影响。在本发明中，重建图像SNR的计算考虑到所有采集通道上所有相位编码线上读方向对应位置的原始数据的贡献，因此能够更准确地反映并行采集图像重建的SNR损失情况，尤其是能够准确地反映重建中结合的参考线的多少所导致的SNR的变化情况。使得对于并行采集图像重建的SNR损失的评估更加准确。

图3是采用本发明的方法计算的SNR随着参考线的变化而变化的示意图。图中的横坐标是重建图像中像素点的位置，纵坐标为SNR效率，图中每条曲线代表不同的参考线所导致的重建后图像的SNR的不同，从图中可以看出，在重建中考虑越多的参考线，重建图像的SNR效率越高。

需要说明的是，本发明可以应用于任意加权系数和采样方案已知的任何阵列线圈结合算法，不限定在上述两个实施例中的GRAPPA方法和SENSE方法。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (12)

1、一种并行采集图像重建的信噪比计算方法，其特征在于，包括：

确定并行采集图像重建的重建表达式，使图像重建满足线性运算ρ＝Mr，其中ρ为重建图像，M为重建表达式，r为原始数据；

根据所述的重建表达式确定图像重建的加权系数；

根据所述加权系数和所述原始数据计算并行采集图像的信噪比。

2、根据权利要求1所述的方法，其特征在于，还包括：

根据并行采集图像的信噪比与标准采集图像的信噪比的比值作为相对信噪比评估并行采集图像重建的信噪比损失。

3、根据权利要求1所述的方法，其特征在于，所述重建表达式为重建矩阵。

4、根据权利要求1至3中任一项所述的方法，其特征在于，所述原始数据为对K空间均匀欠采样获得的数据。

5、根据权利要求1至3中任一项所述的方法，其特征在于，所述原始数据为对K空间均匀欠采样和低频满采样相结合的混合采样获得的数据。

6、根据权利要求1至3中任一项所述的方法，其特征在于，对于基于K空间的重建方式，所述确定并行采集图像重建的重建表达式包括：

确定并行采集图像重建的重建矩阵：M＝C·F·G，其中G表示利用拟合系数对均匀欠采样的K空间位置进行填充同时保留低频满采样的数据的k空间数据填充矩阵；C表示利用线圈灵敏度的复数共轭进行等效于平方和合并SOS的通道合并运算矩阵；F表示傅立叶变换矩阵。

7、根据权利要求1至3中任一项所述的方法，其特征在于，对于基于图像空间的重建方式，所述确定并行采集图像的重建表达式包括：

确定并行成像采样矩阵：E＝U·F^-1·P，其中E表示并行成像采样矩阵，U表示欠采样矩阵，F表示傅立叶变换矩阵，P表示线圈灵敏度分布矩阵；

对所述并行成像采样矩阵求逆得到并行采集图像重建的重建矩阵：M＝inv(E′·E)·E′，其中M表示重建矩阵，E表示所述并行成像采样矩阵。

8、根据权利要求7所述的方法，其特征在于，所述欠采样矩阵U表示均匀欠采样矩阵或者对K空间进行均匀欠采样和低频满采样的混合模式欠采样矩阵。

9、根据权利要求1至3中任一项所述的方法，其特征在于，对于直角坐标系下的采样，所述原始数据r包括不同采集通道、不同相位编码线在读出方向上同一位置的所有数据。

10、根据权利要求1至3中任一项所述的方法，其特征在于，对于非直角坐标系下的采样，所述原始数据包括所有采集通道的所有k空间数据。

11、根据权利要求9中任一项所述的方法，其特征在于，对于直角坐标系下的采样，所述根据重建表达式确定图像重建的加权系数包括：

若读出方向为水平方向，相位编码方向为竖直方向，对于重建图像中的水平方向位置为jRO，竖直方向位置为iPE的像素点(jRO，iPE)，确定其加权系数为：

$w_{\mathrm{jRO}}^{\mathrm{iPE}}=(w_1,w_2,\·\·\·,w_{\mathrm{nCh}\×\mathrm{nPE}}),$

相应的重建数学表达式为：

$\mathrm{Im}{a}_{\mathrm{jRO}}^{\mathrm{iPE}}=w_{\mathrm{jRO}}^{\mathrm{iPE}}\·r_{\mathrm{jRO}},$

其中，nCh表示采集通道数，nPE表示相位编码行数，w_jRO ^iPE表示计算该像素点的加权系数行向量，(w₁，w₂，…，w_nCh×nPE)为重建矩阵M_jRO中的第iPE行，r_jRO为K空间读出方向上位置同为jRO的不同通道，不同相位编码位置的原始数据组成的列向量，矩阵M_jRO为对r_jRO进行重建的重建矩阵，Ima_jRO ^iPE为重建图像中位置为(jRO，iPE)的像素点。

12、根据权利要求11所述的方法，其特征在于，所述根据加权系数和原始数据计算并行采集图像的信噪比包括：

${\mathrm{SNR}}_{\mathrm{jRO}}^{\mathrm{iPE}}=\frac{|w_{\mathrm{jRO}}^{\mathrm{iPE}}\·r_{\mathrm{jRO}}|}{\sqrt{w_{\mathrm{jRO}}^{\mathrm{iPE}}\·\mathrm{\ψ}\·{\left(w_{\mathrm{jRO}}^{\mathrm{iPE}}\right)}^{\′}}\sqrt{\mathrm{mPE}}},$ 其中SNR_jRO ^iPE表示点(jRO，iPE)并行采集图像的信噪比，w_jRO ^iPE表示图像重建的加权系数，mPE为实际采集的相位编码线数，ψ表示噪声相关矩阵。

Images