CN108577840A

CN108577840A - 一种动态磁共振稳健pca成像方法

Info

Publication number: CN108577840A
Application number: CN201810141292.9A
Authority: CN
Inventors: 杨敏; 周宝来; 荆晓远; 晏士友
Original assignee: Nanjing Post and Telecommunication University
Current assignee: Nanjing Post and Telecommunication University; Nanjing University of Posts and Telecommunications
Priority date: 2018-02-11
Filing date: 2018-02-11
Publication date: 2018-09-28

Abstract

本发明公开了一种动态磁共振稳健PCA成像方法，属于图像重建领域。本方法首先根据压缩感知理论，进行部分k‑t空间测量的动态磁共振图像采样，然后引入低秩加稀疏模型，将其应用到欠采样动态磁共振成像，最后运用交替方向法（ADMM）解决稳健PCA优化问题。本方法不但保持了图像的结构信息，同时低秩加稀疏的特性，在同样欠采样的情况下提高了图像重建质量。

Description

一种动态磁共振稳健PCA成像方法

技术领域

本发明涉及图像重建的方法，特别涉及一种动态磁共振稳健PCA成像方法，属于图像重建技术领域。

背景技术

动态磁共振成像，是有着空间和时间信息的磁共振信号，被用于多项临床应用，例如心血管影像，动态对比增强磁共振成像等。然而，磁共振成像是一个内在的慢过程，由于核松弛，外周神经刺激，功率消耗和噪声信号的约束，它的空间和时间的核磁共振分辨率受到限制，另外，长时间扫描会影响病人的舒适度，因而也增加运动伪影的产生。

因此，随着磁共振成像的发展，许多方法已经被提出，来降低采样时间。经典的方法包括回波平面成像，快速低角度镜头成像和使用多个接收线圈的并行磁共振成像。

压缩感知(CS)在MRI应用中能够提高成像速度和效率。CS理论要求采集空间和空间之间的图像稀疏表示没有关联。幸运的是，MR图像序列常常在空间和时间域中提供冗余信息，这为CS的应用提供了有利条件。此外，由于广泛的时空相关性，结果图像的稀疏表示这个想法是很容易扩展到动态MRI重建(DMRI)，K-T FOCUSS是一个成功的方法，它用了FOCUSS算法在时空变换域产生了稀疏约束，延长FOCUS技术与运动估计和运动补偿来压缩感知心脏MRI框架。但当运动是非周期性的时候，往往效果并不理想，近年来，研究者们致力于利用矩阵的低秩性质而不是单纯的向量稀疏性。Lingala et al.提出了k-t SLR算法，利用在KLT域的低秩和全局稀疏特性进行MRI重建。然而，该算法没有考虑到MRI图像的结构稀疏性。这种局限性阻碍了进一步的改进。有些学者提出了基于块的字典学习的DMRI重建技术。然而，基于块的学习不能有效地用于DMRI重建，主要是因为DMRI序列的大小往往很大，学习这样的大型数据集的字典一般是低效的。即使我们不考虑计算的局限性，它获取学习稀疏词典这样巨大的DMRI训练序列也是不现实的。目前，稳健主成分分析(RPCA)已被用于恢复的动态图像数据的低等级结构的探讨。

发明内容

本发明所要解决的技术问题是：为了克服以上现有算法的缺陷，提高成像质量及加快成像速度,本发明提出一种动态磁共振稳健PCA成像方法。

本发明为解决上述技术问题采用以下技术方案：

一种动态磁共振稳健PCA成像方法，其特征在于，步骤包括：

步骤1、根据压缩感知理论，进行部分k-t空间测量的动态磁共振图像采样，得到k-t空间测量矩阵y；

步骤2、将低秩加稀疏模型引入步骤1所述动态磁共振图像，具体如下：

步骤201、基于低秩矩阵填充的方法，用矩阵X表示要恢复的动态磁共振图像数据，有其中矩阵的每一列代表一个N_x×N_y的二维图像，共采集N_t张在时间轴上排列；

步骤202、由有限维度的时空磁共振成像模型，矩阵y与X间有：

y＝E(X)+n

其中y∈C^P表示步骤1所得k-t空间测量矩阵，是MRI的编码操作，P＜＜N_xN_y×N_t，X是要恢复出的动态磁共振图像矩阵，n∈C^P是噪声向量；

步骤203、利用低秩和稀疏惩罚组合优化关于矩阵X的问题如下：

其中，ψ表示秩先验信息，φ表示稀疏先验信息，α,β分别是对应的惩罚系数；

步骤3、运用交替方向法解决稳健PCA优化问题，具体为：

步骤301、用稳健PCA描述凸极小化问题为

其中X是一大小为N_xN_y×N_t的MR数据矩阵，λ_p是分解参数，矩阵X被分解成一个低秩矩阵L和稀疏矩阵S的叠加，L对应低秩矩阵成分即缓慢变化的部分，S对应稀疏矩阵成分即快速变化的部分；

步骤302、通过交替方向乘子法求解稳健PCA：交替方向乘子法的迭代方法是分别极小化L和S，然后再更新拉格朗日乘子，其增广拉格朗日函数为：

其中<·,·>表示内积操作，Z是线性约束的拉格朗日乘子，δ是惩罚参数；

步骤303、求解极小化问题min_L,SJ(L,S)，其中凸目标函数定义如下：

其中E表示磁共振成像的编码操作，是下奈奎斯特采样和傅里叶变换，F_t代表沿着时间方向的傅里叶变换，μ是正则化参数；

步骤304、通过变量分裂法结合增广拉格朗日来极小化问题，算法基于交替方向法；变量分裂如下：

相应的增广拉格朗日函数是:

其中Z_i是拉格朗日乘子，表示真实部分；

忽略与优化不相关的常量，则表示成：

步骤305、采用交替方向乘子法ADMM对模型进行求解。

进一步地，步骤305所述采用交替方向乘子法ADMM对模型进行求解，具体过程如下：

(a)固定Q,L,S，更新P：

(b)固定P,L,S，更新Q：

(c)固定P,Q,S，更新L：

(d)固定P,Q,L，更新S：

将问题转化为求解P,Q,L,S各自的子问题，分别求闭解。

进一步地，步骤305所述惩罚参数δ₁,δ₂设置成固定值1，当算法迭代次数达到100或者停止准则||X^k+1-X^k||_F/||X^k||_F≤10^-6时，算法结束。

将问题转化为求解P,Q,L,S各自的子问题，分别求闭解。用交替方向法极小化L_A是将问题转化为求解P,Q,L,S各自的子问题，这些子问题都有闭解。核范数极小化问题通过奇异值阈值法求解；L1范数问题通过软阈值求解，其他子问题二次项成一个线性等式。

本发明采用以上技术方案与现有技术相比，具有以下技术效果：

本方法不但保持了图像的结构信息，同时低秩加稀疏的特性，在同样欠采样的情况下提高了图像重建质量。

附图说明

图1是本发明的流程图；

图2是幻影模拟数据图像对比。

具体实施方式

下面结合附图对本发明的技术方案做进一步的详细说明：

本技术领域技术人员可以理解的是，除非另外定义，这里使用的所有术语(包括技术术语和科学术语)具有与本发明所属领域中的普通技术人员的一般理解相同的意义。还应该理解的是，诸如通用字典中定义的那些术语应该被理解为具有与现有技术的上下文中的意义一致的意义，并且除非像这里一样定义，不会用理想化或过于正式的含义来解释。

本发明提供一种动态磁共振稳健PCA成像方法，流程如图1所示，包括：

步骤1、根据压缩感知理论，进行部分k-t空间测量的动态磁共振图像采样：

已知动态磁共振成像的过程中存在关系：

d(k,t)＝∫γ(x,t)e^-j2π(k·x)dx+n(k,t)

其中d(k,t)表示测量信号，γ(k,t)表示成像函数，n(k,t)表示添加的噪声。在MRI中，噪声通常是均匀分布的复数数据。动态磁共振成像目的是在高时空分辨率下，从部分测量值d(k,t)中重建出γ(k,t)。

步骤2、引入低秩加稀疏模型，将其应用到欠采样动态磁共振成像：

假设N_t个N_x×N_y维度的图像，基于低秩矩阵填充的方法，可表示成一个矩阵X，矩阵的每一列代表一个向量化的MR图像序列，由于图像之间的高度相关性，故此矩阵非常近似低秩的。有限维度的时空磁共振成像模型可写成：

y＝E(X)+n

其中y∈C^P表示k-t空间测量矩阵，是MRI的编码操作，P＜＜N_xN_y×N_t，X是要恢复出的矩阵，n∈C^P是噪声向量。

利用低秩和稀疏惩罚组合优化问题如下：

其中，ψ表示秩先验信息，φ表示稀疏先验信息。

步骤3、运用交替方向法(ADMM)解决稳健PCA优化问题：

稳健主成分分析(RPCA)是将一个已知矩阵分解成一个低秩成分和一个稀疏成分的一种数学方法。已知一个MR数据矩阵大小为N_xN_y×N_t，RPCA描述了如下凸极小化问题：

通过交替方向乘子法(ADMM)能有效求解RPCA，其增广拉格朗日函数是：

其中<·,·>表示内积操作，Z是线性约束的拉格朗日乘子，δ是惩罚参数。ADMM的迭代方法是分别极小化L和S，然后再更新拉格朗日乘子。

图像重建就是求解极小化问题min_L,SJ(L,S)，其中凸目标函数定义如下：

其中E表示磁共振成像的编码操作，是下奈奎斯特采样和傅里叶变换，F_t代表沿着时间方向的傅里叶变换，μ是正则化参数，λ_p是分解参数。

重建问题中所包含的的那些先验信息假定了动态图像数据可以被分解成近似低秩和近似稀疏成分。附加操作算子F_t说明了，本文方法处理的是欠采样动态图像数据的重建，它是用沿着时间方向的傅里叶变换作为稀疏变换以提高稀疏性。

本文提出的方法称为k-t L+S，如下所示：

由于目标函数的第一项代表数据拟合项，保证了重建数据与k-t空间采样的部分数据之间的一致性，第二项围绕的是低秩加稀疏分解。

通过变量分裂法结合增广拉格朗日来极小化问题，算法基于交替方向法。变量分裂如下：

相应的增广拉格朗日函数是:

其中Z_i是拉格朗日乘子，表示真实部分，忽略与优化不相关的常量，还可以被写成：

采用交替方向乘子法ADMM对模型进行求解，具体过程如下：

(e)固定Q,L,S，更新P：

(f)固定P,L,S，更新Q：

(g)固定P,Q,S，更新L：

(h)固定P,Q,L，更新S：

用交替方向法极小化L_A是将问题转化为求解P,Q,L,S各自的子问题，这些子问题都有闭解。核范数极小化问题通过奇异值阈值法求解；L1范数问题通过软阈值求解，其他子问题二次项成一个线性等式。惩罚参数δ₁,δ₂设置成固定值1，当算法迭代次数达到100或者停止准则||X^k+1-X^k||_F/||X^k||_F≤10^-6时，k-t L+S算法结束，以确保算法收敛。

如图2所示的模拟数据图像对比，本方法不但保持了图像的结构信息，同时低秩加稀疏的特性，在同样欠采样的情况下提高了图像重建质量。

以上所述仅是本发明的部分实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims

1.一种动态磁共振稳健PCA成像方法，其特征在于，步骤包括：

步骤201、基于低秩矩阵填充的方法，用矩阵X表示要恢复的动态磁共振图像数据，有其中矩阵的每一列代表一个N_x×N_y的二维图像，共采集N_t张在时间轴上排列的图像；

y＝E(X)+n

其中y∈C^P表示步骤1所得k-t空间测量矩阵，E:是MRI的编码操作，P＜＜N_xN_y×N_t，n∈C^P是噪声向量；

步骤3、运用交替方向法解决稳健PCA优化问题，具体为：

步骤301、用稳健PCA描述凸极小化问题为

其中λ_p是分解参数，矩阵X被分解成一个低秩矩阵L和稀疏矩阵S的叠加，L对应低秩矩阵成分即缓慢变化的部分，S对应稀疏矩阵成分即快速变化的部分；

步骤303、求解极小化问题min_L,S J(L,S)，其中凸目标函数J(L,S):定义如下：

步骤304、基于交替方向法，通过变量分裂法结合增广拉格朗日来极小化问题；变量分裂如下：

相应的增广拉格朗日函数是:

其中Z_i是拉格朗日乘子，表示真实部分；

忽略与优化不相关的常量，则表示成：

步骤305、采用交替方向乘子法ADMM对模型进行求解，获得重建图像。

2.如权利要求1所述的一种动态磁共振稳健PCA成像方法，其特征在于，步骤305所述采用交替方向乘子法ADMM对模型进行求解，具体过程如下：

(a)固定Q,L,S，更新P：

(b)固定P,L,S，更新Q：

(c)固定P,Q,S，更新L：

(d)固定P,Q,L，更新S：

将问题转化为求解P,Q,L,S各自的子问题，分别求闭解。

3.如权利要求2所述的一种动态磁共振稳健PCA成像方法，其特征在于，步骤305所述惩罚参数δ₁,δ₂设置成固定值1，当算法迭代次数达到100或者停止准则||X^k+1-X^k||_F/||X^k||_F≤10^-6时，算法结束。