CN104825161B - 基于过完备字典与先验知识的高质量肺部磁共振成像方法 - Google Patents

Abstract

本发明公开了基于过完备字典与先验知识的高质量肺部磁共振成像方法。该方法首先建立适用于设定对象，设定图像分辨率的过完备字典；通过小角度激发对肺部进行预成像；随后对预成像中获取的K空间数据进行分析，提取局部大点的位置信息；再根据局部大点的位置信息自适应生成大角度欠采样脉冲序列；使用此脉冲序列再次对肺部进行成像；最后挑选之前建立的设定对象和设定图像分辨率的过完备字典对图像进行重建。和现有技术相比，使用本方法可以获得具有更高分辨率和信噪比的高质量肺部磁共振图像。

Description

基于过完备字典与先验知识的高质量肺部磁共振成像方法

技术领域

本发明涉及肺部的磁共振成像(Magnetic Resonance Imaging，MRI)技术领域，具体涉及基于过完备字典与先验知识的高质量肺部磁共振成像方法。

背景技术

近年来随着经济的快速发展，我国的空气污染状况正变得越来越严峻。肺部疾病对国民健康的影响日益严重。急需一种能够对肺部疾病进行早期诊断的肺部成像方法。

使用超极化气体作为造影剂的磁共振成像技术是一种非侵入性的，无放射性的成像方法。其解决了传统磁共振成像方法无法对肺部进行成像的问题，在对肺部疾病的诊断中有着良好的应用前景。该方法需要采用超极化气体作为造影剂，超极化气体的极化度具有不可恢复的特点。成像过程中提高分辨率或通过大角度脉冲激发获取高信噪比都会消耗超极化气体的极化度。所以，提高成像分辨率与提高信噪比存在矛盾。现有技术在分辨率和信噪比之间进行了折中，主要采用小角度激发的脉冲序列对肺部进行成像。

一种解决这种矛盾的思路是在对肺部进行成像时进行欠采样。如果能够有选择性的对肺部K空间进行采样即欠采样，就可以在提高分辨率的同时采用大角度的脉冲激发，同时获得具有高分辨率和高信噪比的图像。但是，使用欠采样技术也有其自身的问题。欠采样技术能否得到好的重建效果取决于能否找到合适的线性变换对重建图像进行稀疏化分解，另外一方面取决于欠采样轨迹是否能充分获取K空间中的重要信息。

本发明公开的方法，通过建立一种用于肺部磁共振图像的过完备字典库，能够针对性地提供将重建图像充分稀疏化的过完备字典；同时，通过预扫描获取先验知识，指导欠采样轨迹的生成过程，有效获取成像对象K空间中的重要信息。通过这两种技术手段，高效地利用超极化气体的极化度，同时获得具有高分辨率和高信噪比的高质量肺部磁共振图像。

发明内容

本发明的目的是针对现有技术存在的问题，提出的一种基于过完备字典与先验知识的高质量肺部磁共振成像方法。

为了实现上述目的，本发明采用了以下技术措施：

基于过完备字典与先验知识的高质量肺部磁共振成像方法，包括以下步骤：

步骤1、建立适用于设定图像分辨率的设定对象的肺部图像(如小鼠，分辨率256*256像素)的过完备字典的方法为:

步骤1.1，从预定张数的设定对象的肺部的图像中任意选出一幅图像，对图像进行分块处理，将每一块都按照以下方法进行排列：将块的左数第一列提取出来作为一个列向量a1，再提取块左数第二列，放在的a1的下方，这样组成了一个新的列向量a2，以此类推，直到提取到块的最右一列放在当前列向量的下方。这样每个块都会得到一个列向量，组成了学习对象矩阵A的每一列；

步骤1.2，对字典D进行初始化得到字典D₀，初始化包括以下步骤：计算学习对象矩阵A的左奇异矩阵，挑选左奇异矩阵的前设定列组成字典D₀，挑选左奇异矩阵作为初始矩阵可以减少计算量，节省图像重建时间；

步骤1.3，以公式(1)所描述的目标进行计算，获得稀疏系数矩阵Γ；

其中，λ为稀疏系数矩阵Γ的L1范数的加权系数，m为字典被更新的次数(m≥0)，m＝0时D_m为初始化得到的字典D₀，采用的算法为FISTA(Fast Iterative Shrinkage-Thresholding Algorithm)迭代算法，计算方法为：

(1)：设y₁＝Γ₀∈Rⁿ,t₁＝1

……

(k)：(k>1)，计算

其中L＝2maxλ(D^TD),λ(D^TD)为D^TD的特征值，通过设定k来决定算法计算的迭代次数，取k>20可以获得较高的图像质量。这种迭代算法可以通过解决凸优化问题找到全局最优解，同时以1/k²的速度快速收敛。

判断稀疏系数矩阵Γ是否满足公式(2)，若满足公式(2)，则对单幅图像建立字典的过程完成，进入步骤1.5，若不满足公式(2)，则进入步骤1.4

其中λ为矩阵Γ的L1范数的加权系数，ε₁为设定的优化值；

步骤1.4，按照K-SVD(Kernel-Singular Value Decomposition)方法对D_m进行逐列更新，获得新的字典D_m+1，将D_m+1代入步骤1.3进行迭代计算。其中K-SVD方法是一种经典的字典学习算法，在矩阵维数较高时具有良好的效果，可以通过迭代得到优化的字典；

步骤1.5，重复步骤1.1—1.3，直到对步骤1.1中所有预定张数的设定对象的肺部的图像的建立字典过程完成，则每一幅图像都有一个对应的字典。将所有获得的字典按列组合获得新的字典，对新字典中的任意两列相互计算点积，若点积大于设定的阈值(如0.6)，则删除其中任意一列，这样可以保证字典中每一列的相对独立性，去掉多余的列来减少计算量，节省重建时间，最终获得适用于设定对象的肺部的过完备字典；

步骤2、通过变角FLASH(Fast Low Angle Shot)序列，设定小角度激发脉冲序列，对设定对象肺部进行预成像，获得相应的K空间数据，其中小角度激发可以获得高分辨率的肺部图像信息；

步骤3、对步骤2中获得的K空间数据进行分析，获取其中局部大点的位置信息，挑选局部大点的标准为数据幅值大于其四邻域或八邻域的点。这样挑选得到的点包含了图像中较多地轮廓和细节信息，与挑选其他点相比可以减少重建图像的伪影，提高信噪比；

步骤4、预先在计算机中随机生成多个的固定欠采样率的笛卡尔采样矩阵，通过计算这些固定欠采样率的笛卡尔采样矩阵的点扩散函数，按点扩散函数的平均值从小到大进行排序，挑选出排名靠前(如前10)的一些笛卡尔采样矩阵，储存在计算机中。之后根据步骤3中获得的局部大点的位置信息，计算各个笛卡尔采样矩阵中与局部大点位置的相干性，得到其中相干性最大的笛卡尔采样矩阵，将该相干性最大的采样矩阵作为采样轨迹。通过变角FLASH(Fast Low Angle Shot)序列，生成大角度激发的欠采样脉冲序列。其中预先获得伪影较小的采样矩阵虽然计算量巨大，但是可以在成像之前完成并反复使用，提高了重建图像的质量又能降低重建时间，而大角度激发的脉冲序列可以提高信噪比；

步骤5、根据步骤4中生成的采样轨迹，对设定对象的肺部进行超极化磁共振成像，获得设定对象的肺部K空间数据y；

步骤6、挑选步骤1中建立的设定图像分辨率的设定对象的肺部的过完备字典，对步骤5中获得的K空间数据进行图像迭代重建的步骤为：

步骤6.1将得到的K空间数据y填零后进行傅里叶变换得到初始图像x_n＝0，即迭代次数n＝0；

步骤6.2字典D选用步骤1已经得到的用于肺部磁共振图像的过完备字典，；

步骤6.3对图像x_n矩阵进行分块处理，将每一块都按照以下方法进行排列：将块的左数第一列提取出来作为一个列向量a1，再提取块左数第二列，放在的a1的下方，这样组成了一个新的列向量a2，以此类推，直到提取到块的最右一列放在当前列向量的下方。这样每个块都会得到一个列向量，组成了学习对象矩阵A_n的每一列；

步骤6.4以获得公式(3)的最优解为目标进行计算，获得稀疏系数矩阵Γ；

其中，λ为稀疏系数矩阵Γ的L1范数的加权系数，采用的算法为FISTA(FastIterative Shrinkage-Thresholding Algorithm)迭代算法，计算方法为：

(1)：设y₁＝Γ₀∈Rⁿ,t₁＝1

……

(k)：(k>1)，计算

其中L＝2maxλ(D^TD),λ(D^TD)为D^TD的特征值，通过设定k来决定算法计算的迭代次数，取k>20可以获得较高的图像质量。这种迭代算法可以通过解决凸优化问题找到全局最优解，同时以1/k²的速度快速收敛。

步骤6.5将步骤6.4中所得的DΓ还原成块后求平均值得到图像矩阵x’，再对x’求傅里叶变换得到K空间数据S，按照公式(4)、公式(5)计算得到图像矩阵x_n+1，

x_n+1＝F^-1(Fx_n+1) (5)

其中Fx_n+1(k_i,k_j)为对图像矩阵x_n+1进行傅立叶变换得到的矩阵，F^-1为傅立叶逆变换，S(k_i,k_i)为K空间数据S中某个点的像素值，S₀为对肺部K空间数据y填零得到的矩阵，γ为S₀的加权系数，Ω为采样轨迹所覆盖的区域。对图像矩阵x_n+1用步骤6.3中的方法得到新的学习对象矩阵A_n+1，判断x_n+1和A_n+1是否满足公式(6)或公式(7)，满足其中一个则计算结束，得到最终重建的图像矩阵x_n+1，x_n+1和A_n+1同时不满足公式(6)和公式(7)，则将A_n+1迭代为步骤6.4中的A_n并返回步骤6.4进行处理。

n≥Ν (7)

其中F_u为部分傅立叶变换，y为肺部K空间数据，λ为稀疏系数矩阵Γ的L1范数的加权系数，N为设定的迭代次数的上限，ε₂为设定的优化值。

本发明具有以下优点和积极效果：

本发明可以提高对超极化气体极化度的利用率，获得具有高分辨率和高信噪比的肺部磁共振图像，为肺部疾病的早期诊断和治疗提供良好的技术支持。

附图说明

图1发明方法流程图；

图2过完备字典建立流程图

图3肺部K空间大点分布情况；

图4 K空间欠采样扫描轨迹；

图5利用过完备字典与先验知识进行图像重建流程图；

图6(a)利用过完备字典与先验知识方法获得肺部磁共振图像效果。

图6(b)利用小角度激发的脉冲序列获得肺部磁共振图像效果。

具体实施方式

下面结合附图和实施示例对本发明进一步说明(方法流程图如图1所示)：

步骤1、挑选出100幅大小为256*256像素(以下的矩阵及向量单位均为像素)的小鼠肺部磁共振图像数据，按照以下步骤建立过完备字典(流程图如图2所示)：

步骤1.1，从100幅图像中任意选出一幅图像，对图像进行分块处理，块的大小为8*8，每个块都相互重叠，相邻的块之间相差一行或一列，即重叠7行或者7列，这样可以获得62001个块，将每一块都按照以下方法进行排列：将块的左数第一列提取出来作为一个列向量a1，大小为8*1，再提取块左数第二列，放在的a1的下方，这样组成了一个新的列向量a2，大小为16*1，以此类推，直到提取到块的最右一列放在当前列向量的下方。这样每个块都会得到一个大小为64*1的列向量，组成了学习对象矩阵A的每一列，学习对象矩阵A的大小为64*62001；

步骤1.2，设定字典D的大小为64*64，对字典D进行初始化得到字典D₀，初始化包括以下步骤：计算学习对象矩阵A的左奇异矩阵，挑选左奇异矩阵的前64列组成字典D₀；

步骤1.3，以获得公式(1)的最优解为目标进行计算，获得稀疏系数矩阵Γ；

其中，λ＝140，m为字典被更新的次数(m≥0)，m＝0时D_m为初始化得到的字典D₀，采用的算法为FISTA(Fast Iterative Shrinkage-Thresholding Algorithm)迭代算法，计算方法为：

(1)：设y₁＝Γ₀∈Rⁿ,t₁＝1

……

(k)：(k>1)，计算

其中L＝2maxλ(D^TD),λ(D^TD)为D^TD的特征值，通过设定k来决定算法计算的迭代次数，取k>20可以获得较高的图像质量。这种迭代算法可以通过解决凸优化问题找到全局最优解，同时以1/k²的速度快速收敛。

判断稀疏系数矩阵Γ是否满足公式(2)，若满足条件，则对单幅图像建立字典的过程完成，进入步骤1.5，不满足公式(2)则进入步骤1.4

其中ε₁＝20000，λ＝140；

步骤1.4，按照K-SVD(Kernel-Singular Value Decomposition)方法对字典D_m进行逐列更新，获得新的字典D_m+1，将D_m+1迭代入步骤1.3的公式(1)中的D_m进行迭代计算；

步骤1.5，重复步骤1.1到步骤1.3，完成对这100幅小鼠肺部图像的建立字典过程，共获得100个大小为64*64字典，将所有获得的字典按列组合获得新的字典，大小为64*6400，对新字典中的任意两列相互计算点积，若点积大于设定的阈值0.6，则删除其中任意一列，最终获得了大小为64*3556的，用于小鼠肺部磁共振图像的过完备字典；

步骤2、通过选择变角的FLASH(Fast Low Angle Shot)序列，设定翻转角为8度，对一只小鼠的肺部进行预成像，设定TR＝2s，TE＝8.667ms，采样轨迹为笛卡尔采样方式，进行满采样扫描，得到K空间数据，获得相应的大小为256*256的K空间数据，通过傅里叶变换重建图像；

步骤3、对步骤2中获得的K空间数据进行分析，提取局部大点的位置信息。挑选局部大点的标准为数据幅值大于其八邻域的点，经过统计，K空间中八邻域局部大信号点数量占总点数的12.8％，局部大点的位置分布信息，如图3所示；

步骤4、预先在计算机中随机生成1000万个的欠采样率为25％的笛卡尔采样矩阵，通过计算1000万个的欠采样率为25％的笛卡尔采样矩阵的点扩散函数，按各个点扩散函数的平均值从小到大进行排序，挑选出排名靠前的20个笛卡尔采样矩阵，储存在计算机中。根据步骤3中获得的局部大点的位置分布信息，计算笛卡尔采样矩阵中与局部大点位置的相干性，得到其中相干性最大为71％的采样矩阵，将该相干性最大的采样矩阵作为采样轨迹，如图4所示。通过变角FLASH(Fast Low Angle Shot)序列，设定TR＝2s，TE＝8.667ms，生成翻转角为20度的欠采样脉冲序列；

步骤5、根据步骤4中生成的欠采样脉冲序列，对步骤2中的小鼠肺部进行超极化磁共振成像，获得肺部K空间数据y，大小为64*256；

步骤6、挑选步骤1中建立的小鼠肺部的过完备字典，对步骤5中获得的K空间数据y进行图像迭代重建，流程图如图5所示，图像迭代重建的步骤为：

步骤6.1，将得到的肺部K空间数据y填零后进行傅里叶变换得到初始图像x_n＝0，初始图像x_n中的n＝0，即迭代次数n＝0；

步骤6.2，字典D选用步骤1已经得到的用于肺部磁共振图像的过完备字典，大小为64*3556；

步骤6.3，对图像x_n矩阵进行分块处理，块的大小为8*8，每个块都相互重叠，相邻的块之间相差一行或一列，即重叠7行或者7列，这样可以获得62001个块，将每一块都按照以下方法进行排列：将块的左数第一列提取出来作为一个列向量b1，大小为8*1，再提取块左数第二列，放在的b1的下方，这样组成了一个新的列向量b2，大小为16*1，以此类推，直到提取到块的最右一列放在当前列向量的下方。这样每个块都会得到一个大小为64*1的列向量，组成了学习对象矩阵A_n的每一列，学习对象矩阵A_n的大小为64*62001；

步骤6.4，以获得公式(1)的最优解为目标进行计算，获得稀疏系数矩阵Γ；

其中，λ＝140，采用的算法为FISTA(Fast Iterative Shrinkage-ThresholdingAlgorithm)迭代算法，计算方法为：

(1)：设y₁＝Γ₀∈Rⁿ,t₁＝1

……

(k)：(k>1)，计算

其中k＝25，L＝2maxλ(D^TD),λ(D^TD)为D^TD的特征值；

步骤6.5，将步骤6.4中所得的DΓ还原成块后求平均值得到图像矩阵x’，再对x’求傅里叶变换得到K空间数据S，按照公式(4)、公式(5)计算得到图像矩阵x_n+1，

x_n+1＝F^-1(Fx_n+1) (5)

其中Fx_n+1(k_i,k_j)为对图像矩阵x_n+1进行傅立叶变换得到的矩阵，F^-1为傅立叶逆变换，S(k_i,k_i)为K空间数据S中某个点的像素值，S₀为对肺部K空间数据y填零得到的矩阵，γ＝14，Ω为采样轨迹所覆盖的区域。对图像矩阵x_n+1用步骤6.3中的方法得到新的学习对象矩阵A_n+1，判断是否满足公式(6)或公式(7)，满足则计算结束，得到最终重建的图像矩阵x_n+1，不满足则将A_n+1迭代为步骤6.4中的A_n并返回步骤6.4进行处理。

n≥Ν (7)

其中F_u为部分傅立叶变换，y为肺部K空间数据，λ＝140，N＝15，ε₂＝25000，γ＝14。

最终在迭代了15次后得到了小鼠的肺部图像，如图6(a)。6(b)为采用小角度激发的脉冲序列方法得到的图像。

本文中所描述的具体实施方式仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施方式做各种各样的修改或补充或采用类似的方式替代，但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。

Claims (3)

1.基于过完备字典与先验知识的高质量肺部磁共振成像方法，其特征在于，包括以下步骤：

步骤1、建立适用于设定图像分辨率的设定对象的肺部图像的过完备字典；

步骤2、通过变角FLASH序列，设定小角度激发脉冲序列，对设定对象的肺部进行预成像，获得相应的K空间数据；

步骤3、对步骤2中获得的K空间数据进行分析，获取其中局部大点的位置信息，挑选局部大点的标准为数据幅值大于其四邻域或八邻域的点；

步骤4、预先在计算机中随机生成多个的固定欠采样率的笛卡尔采样矩阵，通过计算这些固定欠采样率的笛卡尔采样矩阵的点扩散函数，按点扩散函数的平均值从小到大进行排序，挑选出排名靠前的若干个笛卡尔采样矩阵，储存在计算机中，根据步骤3中获得的局部大点的位置信息，计算各个笛卡尔采样矩阵中与局部大点位置的相干性，得到其中相干性最大的笛卡尔采样矩阵，将该相干性最大的笛卡尔采样矩阵作为采样轨迹，通过变角FLASH序列，生成大角度激发的欠采样脉冲序列；

步骤5、根据步骤4中生成的采样轨迹，对设定对象的肺部进行超极化磁共振成像，获得设定对象的肺部的K空间数据y；

步骤6、挑选步骤1中建立的设定图像分辨率的设定对象的肺部的过完备字典，对步骤5中获得的K空间数据y进行图像迭代重建。

2.根据权利要求1所述的基于过完备字典与先验知识的高质量肺部磁共振成像方法，其特征在于，所述的步骤1包括以下步骤：

步骤1.1、从预定张数的设定对象肺部图像中任意选出一幅图像，对图像进行分块处理，每一块的各列按列序号纵向拼接为单列的列向量，每个块对应的列向量构成了学习对象矩阵A的每一列；

步骤1.2、对字典D进行初始化得到字典D₀，初始化包括以下步骤：计算学习对象矩阵A的左奇异矩阵，挑选左奇异矩阵的前设定列组成字典D₀；

步骤1.3、以获得以下公式1的最优解为目标进行计算，获得稀疏系数矩阵Γ；

其中，λ为稀疏系数矩阵Γ的L1范数的加权系数，m为字典被更新的次数，m≥0，m＝0时D_m为初始化得到的字典D₀，采用的算法为FISTA迭代算法，

判断稀疏系数矩阵Γ是否满足以下公式2，若满足以下公式2，则对单幅图像建立字典的过程完成，进入步骤1.5；若不满足以下公式2，进入步骤1.4

其中λ为矩阵Γ的L1范数的加权系数，ε₁为设定的优化值；

步骤1.4、按照K-SVD方法对字典D_m进行逐列更新，获得新的字典D_m+1，将D_m+1代入步骤1.3进行迭代计算；

步骤1.5、重复步骤1.1～1.3，直到对步骤1.1中所有预定张数的设定对象肺部图像的建立字典过程完成，则每一幅图像都有一个对应的字典，将每一幅图像对应的字典按列组合获得新的字典，对新的字典中的任意两列相互计算点积，若点积大于设定的阈值，则删除其中任意一列，最终获得适用于设定对象的肺部的图像的过完备字典。

3.根据权利要求1所述的基于过完备字典与先验知识的高质量肺部磁共振成像方法，其特征在于，所述的步骤6包括以下步骤：

步骤6.1、将得到的K空间数据y填零后进行傅里叶变换得到初始图像x₀，即初始图像为x_n＝0，n为迭代次数；

步骤6.2、字典D选用步骤1获得的过完备字典；

步骤6.3、对图像x_n矩阵进行分块处理，每一块的各列按列序号纵向拼接为单列的列向量，每个块对应的列向量构成了学习对象矩阵A_n的每一列；

步骤6.4、以获得以下公式3的最优解为目标进行计算，获得稀疏系数矩阵Γ；

其中，λ为稀疏系数矩阵Γ的L1范数的加权系数，采用的算法为FISTA迭代算法，

步骤6.5、将步骤6.4中所得的DΓ还原成块后求平均值得到图像矩阵x’，再对x’求傅里叶变换得到K空间数据S，按照下面公式4和公式5计算得到图像矩阵x_n+1，

x_n+1＝F^-1(Fx_n+1) 公式5

其中Fx_n+1(k_i,k_j)为对图像矩阵x_n+1进行傅立叶变换得到的矩阵，F^-1为傅立叶逆变换，S(k_i,k_j)为K空间数据S中某个点的像素值，S₀为对肺部K空间数据y填零得到的矩阵，γ为S₀的加权系数，Ω为采样轨迹所覆盖的区域，对图像矩阵x_n+1用步骤6.3中的方法得到新的学习对象矩阵A_n+1，若x_n+1和A_n+1满足下列公式6和公式7中任意一个公式则计算结束，得到最终重建的图像矩阵x_n+1；若x_n+1和A_n+1不满足下列公式6和公式7中任意一个公式，则将A_n+1迭代为步骤6.4中的A_n并返回步骤6.4进行处理，

n≥Ν 公式7

其中F_u为部分傅立叶变换，y为肺部K空间数据，λ为稀疏系数矩阵Γ的L1范数的加权系数，N为设定的迭代次数的上限，ε₂为设定的优化值。