CN106526511A - 基于k空间中心鬼影定位的SPEED磁共振成像方法 - Google Patents

Abstract

本发明公开了一种基于k空间中心鬼影定位的SPEED磁共振成像方法,本发明主要包括九大步骤：k空间数据采集、填零傅立叶重建、差分变换、建立低分辨重叠鬼影图、鬼影阶数的确立、双层鬼影模型求解、鬼影的分离、多个鬼影映射图的配准求和、逆滤波重建。采用本发明方法可以将k空间中心部分数据用于SPEED成像时的鬼影阶数的定位，避免了常规SPEED成像方法的最小平方误差求解鬼影阶数的步骤，从而将所需采集的三组k空间欠采样数据减少到两组，进一步缩短了SPEED成像方法的数据采集时间。

Description

基于k空间中心鬼影定位的SPEED磁共振成像方法

技术领域

本发明属于磁共振的图像成像领域，涉及一种基于k空间中心鬼影定位的SPEED磁共振成像方法。

背景技术

磁共振成像(Magnetic Resonance Imaging,MRI)对机体没有不良影响，且具有优良的软组织分辨能力，已在临床疾病检测中得到广泛应用。但是，MRI在临床应用中还常受到数据采集时间过长的限制。研究人员已通过提高MRI硬件性能、采用高效的k空间(频率空间)数据采集轨迹、研制快速序列和并行数据采集等方式来提高MRI数据采集的速度，但是临床实际应用中，已有的研究成果还是不能完全满足对快速成像的需求，例如脑功能成像和心脏动态成像等应用中。

MRI的数据采集是在k空间中进行的，在k空间的相位编码(Phase Encoding,PE)方向上减少数据采集点数，可以有效地缩短数据采集时间。SPEED(Skipped Phase Encodingand Edge Deghosting)成像技术就是一种通过在PE方向减少数据点数来缩短数据采集时间的MRI快速成像方法(QS Xiang,Accelerating MRI by skipped phase encoding andedge deghosting(SPEED)，Magnetic Resonance in Medicine,53:1112-1117,2005)。SPEED通过在k空间进行简单规则的欠采样来采集数据，然后基于解析法来重建图像，其解析求解过程非常快速。SPEED成像方法不但易于实现，而且也易于和现有的采集方式结合，是一种很有应用潜力的成像方式。

目前已申请的关于SPEED快速成像方面的MRI专利有：基于小波域稀疏表示的SPEED快速磁共振成像方法(授权号：ZL 2013102071971.1)，提出基于小波域的数据稀疏特性来提高SPEED快速成像方法的成像质量。基于离散余弦变换的SPEED快速磁共振成像方法(授权号：ZL 201310719667.2)，提出基于离散余弦变换来提高SPEED快速成像方法的成像质量。目前还未能查询到任何基于k空间中心鬼影定位的SPEED快速成像方法的授权发明专利或申请。

国内外已发表的关于SPEED成像方面的文章有：2016年，金朝阳、叶海慧、杜一平和向清三提出了基于离散余弦变换和离散小波变换的SPEED快速成像方法(Jin Z,Ye H,DuYP,Xiang QS.Improving image quality for skipped phase encoding and edgedeghosting(SPEED)by exploiting several sparsifying transforms.MagneticResonance in Medicine.75:2031-2045，2016)，该方法利用离散余弦变换和离散小波变换对数据进行稀疏表示，相比于常规SPEED方法采用基于离散差分变换的数据稀疏表示，获得了更好的成像质量。2013年，金朝阳和向清三提出了通用G-SPEED(General-SPEED)采样方法(Jin Z,Xiang QS.Accelerated MRI by SPEED with generalized samplingschemes.Magnetic Resonance in Medicine.70:1674-1681，2013)，突破了传统SPEED方法的采样间隔周期N必须是质数(例如：N＝5、7、11)的限制，通过秩判据的方式，使得N不但可为质数，也可为合数(例如：N＝2、4、6、8、9)。2009年，常征等人提出EMA-SPEED(EfficientMultiple Acquisition by SPEED)算法(Chang Z,Xiang QS,Ji J,and Yin FF.Efficientmultiple acquisitions by skipped phase encoding and edge deghosting(SPEED)using shared spatial information.Magnetic Resonance in Medicine.61:229-233,2009)，通过共享多个采集间的相似空间信息进一步缩短了SPEED的数据采集时间，从而可获得比单次采集更高的加速比。2007年，基于MRA数据本身就非常稀疏的特性，常征和向清三将SPEED的双层模型简化到单层模型(Chang Z and Xiang QS.Simplified skippedphase encoding and edge deghosting(SPEED)for imaging sparse objects withapplications to MRA.Med Phys.34:3173-3182,2007)，提出了S-SPEED(Simplified-SPEED)算法，该算法适用于数据本身就非常稀疏的场合，例如暗背景亮信号的MRA应用(Chang Z,Xiang QS,Shen H and Yin FF.Accelerating non-contrast-enhanced MRangiography with inflow inversion recovery imaging by skipped phase encodingand edge deghosting(SPEED).Journal of Magnetic Resonance Imaging.31:757-765,2010)。2006年，常征和向清三将SPEED算法与并行成像技术进一步结合，提出了SPEED-ACE成像法(Chang Z and Xiang QS.Highly accelerated MRI by skipped phase encodingand edge deghosting with array coil enhancement(SPEED-ACE).Med Phys.33:3758-3766,2006)，通过采用多个采集线圈来共同采集k-空间欠采样数据，从而提高成像速度。

以上发表的关于SPEED快速成像方面的文章或已授权的发明专利，重建时是基于采集到的三组欠采样数据来获得重叠鬼影的定位，还未公开过任何基于k空间中心鬼影定位的SPEED快速磁共振成像方法。

发明内容

本发明针对现有SPEED技术的不足，将k空间中心部分的数据用于鬼影定位，从而将SPEED所需采集的三组k空间数据减少到只需采集两组k空间数据，提供了一种新的SPEED数据采集方法，进一步提高了SPEED快速成像的数据采集速度。本发明主要包括九个步骤：k空间数据采集、填零傅立叶重建、差分变换、建立低分辨重叠鬼影图、鬼影阶数的确立、双层鬼影模型求解、鬼影的分离、多个鬼影映射图的配准求和、逆滤波重建。

步骤1：k空间数据采集

在k空间的相位编码方向(即PE方向)每隔N行采集一行数据，共采集两组，分别用S₁和S₂表示。用d₁,d₂表示每组欠采样数据在PE方向上的偏移量，采样方式用N(d₁,d₂)表示。

根据图像大小，在PE方向的k空间中心区域采集16至64行数据，用S_c表示。

步骤2：填零傅立叶重建

对于两组欠采样的数据S₁和S₂，其对应k空间中没有进行数据采集的点用0表示，进行常规的填零傅立叶重建，重建后图像分别用I₁和I₂表示。k空间中每隔N行采集一行数据使得每组数据对应的填零傅立叶重建图像中有N层重叠的鬼影，每个像素点上最多有N层重叠的鬼影。例如，当N＝4时，I₁和I₂上分别有4层重叠的鬼影。

将采集到的k空间中心部分数据S_c也进行填零傅立叶重建，形成一个低分辨率的重建图像I_c。

步骤3：差分变换

对步骤2得到的图像I₁、I₂和I_c分别进行差分变换，得到稀疏的边缘鬼影图像E₁、E₂和E_c。在稀疏的鬼影图像中，每个像素点上通常只有两层重叠的鬼影。

步骤4：建立低分辨重叠鬼影图

在相位编码方向对E_c分别进行长度为N_y×n/N的平移，其中N_y表示沿PE方向的数据矩阵的大小，n表示边缘鬼影的阶数(不同的阶表示鬼影位置不同)，n＝0,1,2,…,N-1。这n个边缘鬼影相加后形成一个重叠的鬼影映射图E_c,n。

步骤5：鬼影阶数的确立

在鬼影映射图E_c,n中，为每个像素点找出两个最强的鬼影，并记录下它们对应的鬼影阶数(n_1s,n_2s)。

步骤6：双层鬼影模型求解

稀疏边缘鬼影图E₁和E₂中，由于每个像素点上通常只有两层鬼影的重叠，因此采用双层稀疏边缘鬼影模型来描述E₁和E₂中的每个像素点。双层稀疏边缘鬼影模型表示为：

公式[1]中为相位因子，G_n1和G_n2分别为每个像素点上需要确定的不同阶的鬼影，n₁和n₂分别表示不同的鬼影阶数。定义为：

公式[2]中d表示每组欠采样数据在PE方向上的偏移量d₁和d₂，n为鬼影阶数。

在公式[1]中，由于E₁、E₂、d和N已知，则基于步骤5得到的鬼影阶数(n_1s,n_2s)，直接解出公式[1]中的两个重叠的鬼影G_n1和G_n2。

步骤7：鬼影的分离

对步骤6得到的G_n1和G_n2中的像素点，按不同的鬼影阶数n进行分类，产生N个分离的鬼影映射图G_n,其中n＝0,1,…,N-1；

步骤8：多个鬼影映射图的配准求和

步骤7得到的N个鬼影映射图G_n，各自对应的鬼影位置不同，通过像素点的移位和对齐来配准。配准后各鬼影图对应的像素点求和后得到没有重叠鬼影的边缘映射图像E₀。

步骤9：逆滤波重建

步骤8得到的边缘映射图E₀经离散傅立叶变换(DFT)到k空间，其对应k空间中实际进行数据采集的点的值用实际采集的数据替代，得到k空间数据R₀。基于逆滤波公式[3]重建出最终的SPEED图像I₀。

公式[3]中IDFT表示离散傅立叶逆变换，k_y表示沿PE方向的k空间位置。

采用本发明方法可以将k空间中心部分数据用于SPEED成像时的鬼影阶数的定位，避免了常规SPEED成像方法的最小平方误差求解鬼影阶数的步骤，从而将所需采集的三组k空间欠采样数据减少到两组，进一步提高了SPEED成像方法的数据采集时间，同时本发明具有以下特点：

(1)本发明只需采集两组k空间欠采样数据，而常规SPEED技术则需采集三组k空间欠采样数据，因此，本发明提高了SPEED成像方法的数据采集速度。

(2)本发明通过利用一部分k空间中心数据来确定鬼影的阶数，避免了常规SPEED方法采用最小平方误差法来求解鬼影阶数的过程，缩短了重建时间。

(3)本发明采取简单有规律方式进行数据的欠采样，无需更改MRI硬件，就能提高数据采集的速度，且易于与常规的临床数据采集方式集成。

附图说明

图1是SPEED数据采集方式的示意图；

图2是SPEED成像过程中部分中间数据的示意图；

图3是采用本发明进行SPEED数据采集和重建实例的结果图。

具体实施方式

以下结合附图对本发明作进一步说明。

本发明主要包括九个步骤：k空间数据采集、填零傅立叶重建、差分变换、建立低分辨重叠鬼影图、鬼影阶数的确立、双层鬼影模型求解、鬼影的分离、多个鬼影映射图的配准求和、逆滤波重建。

步骤1：k空间数据采集

在k空间的相位编码方向(即PE方向)每隔N行采集一行数据，共采集两组，分别用S₁和S₂表示。用d₁,d₂表示每组欠采样数据在PE方向上的偏移量，采样方式用N(d₁,d₂)表示，如图1所示。

根据图像大小，在PE方向的k空间中心区域分别采集16至64行数据，用S_c表示，如图1所示。

步骤2：填零傅立叶重建

如图2所示，对于两组欠采样的数据S₁和S₂，其对应k空间中没有进行数据采集的点用0表示，进行常规的填零傅立叶重建，重建后图像分别用I₁和I₂表示。k空间中每隔N行采集一行数据使得每组数据对应的填零傅立叶重建图像中有N层重叠的鬼影，每个像素点上最多可能有N层重叠的鬼影。例如，当N＝4时，I₁和I₂上分别有4层重叠的鬼影。

如图2所示，将采集到的k空间中心部分数据S_c也进行填零傅立叶重建，形成一个低分辨率的重建图像I_c。

步骤3：差分变换

如图2所示，对步骤2得到的图像I₁、I₂和I_c分别进行差分变换，得到稀疏的边缘鬼影图像E₁、E₂和E_c。在稀疏的鬼影图像中，每个像素点上通常只有两层重叠的鬼影。

步骤4：建立低分辨重叠鬼影图

如图2所示，在相位编码方向对E_c分别进行长度为N_y×n/N的平移，其中N_y表示沿PE方向的数据矩阵的大小，n表示边缘鬼影的阶数(不同的阶表示鬼影位置不同)，n＝0,1,2,…,N-1。这n个边缘鬼影相加后形成一个重叠的鬼影映射图E_c,n。

步骤5：鬼影阶数的确立

如图2所示，在鬼影映射图E_c,n中，为每个像素点找出两个最强的鬼影，并记录下它们对应的鬼影阶数(n_1s,n_2s)。

步骤6：双层鬼影模型求解

如图2所示，稀疏边缘鬼影图E₁和E₂中，由于每个像素点上通常只有两层鬼影的重叠，因此可采用双层稀疏边缘鬼影模型来描述E₁和E₂中的每个像素点。双层鬼影模型可表示为：

公式[1]中为相位因子，G_n1和G_n2分别为每个像素点上需要确定的不同阶的鬼影，n₁和n₂分别表示不同的鬼影阶数。定义为：

公式[2]中d表示每组欠采样数据在PE方向上的偏移量d₁和d₂，n为鬼影阶数。

在公式[1]中，由于E₁、E₂、d和N已知，则基于第5步得到的鬼影阶数(n_1s,n_2s)，可直接解出公式[1]中的两个重叠的鬼影G_n1和G_n2。

步骤7：鬼影的分离

对步骤6得到的G_n1和G_n2中的像素点，按不同的鬼影阶数n进行分类，产生N个分离的鬼影映射图G_n,其中n＝0,1,…,N-1；

步骤8：多个鬼影映射图的配准求和

步骤7得到的N个鬼影映射图G_n，各自对应的鬼影位置不同，可通过像素点的移位和对齐来配准。配准后各鬼影图对应的像素点求和后可得到没有重叠鬼影的边缘映射图像E₀，如图2所示。

步骤9：逆滤波重建

如图2所示，步骤8得到的边缘映射图E₀经离散傅立叶变换(DFT)到k空间，其对应k空间中实际进行数据采集的点的值用实际采集的数据替代，得到k空间数据R₀。基于逆滤波公式[3]重建出最终的SPEED图像I₀。

公式[3]中IDFT表示离散傅立叶逆变换，k_y表示沿PE方向的k空间位置。

以下结合人体膝盖部分的MRI数据，对基于k空间中心鬼影定位的SPEED磁共振成像方法进行实例说明。假设要采集的MRI图像的矩阵大小为k_x×k_y＝256×256。首先以采样方式N(d₁,d₂)＝5(0,1)进行数据采集，在k空间的相位编码PE方向每隔N＝5行采集一行k空间数据，共采集两组，分别得到欠采样的k空间数据S₁和S₂。在信息量集中的k空间中心区域进行全采样，共采集32行相位编码数据。接下来，对两组欠采样数据S₁和S₂和k空间中心数据S_c分别进行常规的填零傅立叶重建，重建后图像分别为I₁、I₂和I_c。然后对图像I₁、I₂和I_c分别进行差分变换，得到稀疏的鬼影图像E₁、E₂和E_c。在相位编码方向对E_c分别进行长度为256×n/5的平移，其中n＝0,1,2,3,4。这5个边缘鬼影相加后生成一个重叠的鬼影映射图E_c,n＝E_c,5。在鬼影映射图E_c,5中，为每个像素点找出两个最强的鬼影，并记录下它们对应的鬼影阶数(n₁,n₂)。基于鬼影阶数(n₁,n₂)和稀疏鬼影图E₁和E₂，可直接解出双层鬼影模型中的两个主要的重叠鬼影G_n1和G_n2。依照鬼影阶数n，对G_n1和G_n2中的像素点分类，产生5个分离的鬼影映射图G_n，n＝0,1,2,3,4。这5个鬼影映射图经过移位、配准和求和后得到无重叠鬼影的边缘映射图E₀。E₀经过离散傅立叶变换(DFT)到k空间，然后用实际采集到的k空间数据替换部分k空间数据，再经过逆滤波公式[3]重建出最终的SPEED图像I₀。如图3所示。图3(a)为人体膝盖的一组欠采样数据的填零重建图I₁，可以看出有5个重叠的鬼影，图3(b)为图3(a)经过差分变换后得到的稀疏鬼影图E₁，图3(c)为去除鬼影后的移位配准图E₀，图3(d)为最终的SPEED重建图I₀，图3(e)为全采样参考图，图3(f)为图3(d)与图3(e)相比的误差图，从图3(f)可以看出，本发明的误差很小。相对均方差误差测量值为(3.65e-4)。可见本发明在减少了SPEED的数据采集量，提高了数据采集速度的同时，获得了较高的成像质量。

Claims (1)

1.基于k空间中心鬼影定位的SPEED磁共振成像方法,其特征在于，该方法包括以下步骤：

步骤1：k空间数据采集

在k空间的相位编码方向每隔N行采集一行数据，共采集两组，分别用S₁和S₂表示；用d₁,d₂表示每组欠采样数据在PE方向上的偏移量，采样方式用N(d₁,d₂)表示，k空间的相位编码方向即PE方向；

根据图像大小，在PE方向的k空间中心区域采集16至64行数据，用S_c表示；

步骤2：填零傅立叶重建

对于两组欠采样的数据S₁和S₂，其对应k空间中没有进行数据采集的点用0表示，进行常规的填零傅立叶重建，重建后图像分别用I₁和I₂表示；k空间中每隔N行采集一行数据使得每组数据对应的填零傅立叶重建图像中有N层重叠的鬼影，每个像素点上最多有N层重叠的鬼影；

将采集到的k空间中心部分数据S_c也进行填零傅立叶重建，形成一个低分辨率的重建图像I_c；

步骤3：差分变换

对步骤2得到的图像I₁、I₂和I_c分别进行差分变换，得到稀疏的边缘鬼影图像E₁、E₂和E_c；

步骤4：建立低分辨重叠鬼影图

在相位编码方向对E_c分别进行长度为N_y×n/N的平移，其中N_y表示沿PE方向的数据矩阵的大小，n表示边缘鬼影的阶数，n＝0,1,2,…,N-1；这n个边缘鬼影相加后形成一个重叠的鬼影映射图E_c,n；

步骤5：鬼影阶数的确立

在鬼影映射图E_c,n中，为每个像素点找出两个最强的鬼影，并记录下它们对应的鬼影阶数(n_1s,n_2s)；

步骤6：双层鬼影模型求解

稀疏边缘鬼影图E₁和E₂中，由于每个像素点上通常只有两层鬼影的重叠，因此采用双层稀疏边缘鬼影模型来描述E₁和E₂中的每个像素点；双层稀疏边缘鬼影模型表示为：

$\begin{array}{c}{E}_1=P_{d1}^{n1}{G}_{n1}+P_{d1}^{n2}{G}_{n2}\\ E_2=P_{d2}^{n1}{G}_{n1}+P_{d2}^{n2}{G}_{n2}\end{array}---\[1\]$

公式[1]中为相位因子，G_n1和G_n2分别为每个像素点上需要确定的不同阶的鬼影，n₁和n₂分别表示不同的鬼影阶数；定义为：

$P_d^n=e^{i(2\mathrm{\π}dn/N)},n=0,1,2,...,N-1---\[2\]$

公式[2]中d表示每组欠采样数据在PE方向上的偏移量d₁和d₂，n为鬼影阶数；

在公式[1]中，由于E₁、E₂、d和N已知，则基于步骤5得到的鬼影阶数(n_1s,n_2s)，直接解出公式[1]中的两个重叠的鬼影G_n1和G_n2；

步骤7：鬼影的分离

对步骤6得到的G_n1和G_n2中的像素点，按不同的鬼影阶数n进行分类，产生N个分离的鬼影映射图G_n,其中n＝0,1,…,N-1；

步骤8：多个鬼影映射图的配准求和

步骤7得到的N个鬼影映射图G_n，各自对应的鬼影位置不同，通过像素点的移位和对齐来配准；配准后各鬼影图对应的像素点求和后得到没有重叠鬼影的边缘映射图像E₀；

步骤9：逆滤波重建

步骤8得到的边缘映射图E₀经离散傅立叶变换到k空间，其对应k空间中实际进行数据采集的点的值用实际采集的数据替代，得到k空间数据R₀；基于逆滤波公式[3]重建出最终的SPEED图像I₀；

$I_0=IDFT\left\{\frac{DFT\[R_0\]}{e^{-j2\mathrm{\π}(k_y/N_y)}-1}\right\}---\[3\]$

公式[3]中IDFT表示离散傅立叶逆变换，k_y表示沿PE方向的k空间位置。