CN106780372A - 一种基于广义树稀疏的权重核范数磁共振成像重建方法 - Google Patents

Abstract

本发明公开了一种基于广义树稀疏的权重核范数磁共振成像重建方法，首先，获取测试磁共振成像采样数据样本进行傅立叶变换；并根据采样的信号构建出树结构的稀疏信号，利用带权重的核范数逼近约束的目标函数的稀疏表达；再通过增广拉格朗日乘子法优化约束的目标函数和交替方向搜索算法对所述测试数据进行迭代更新，直至得到估计的恢复数据；再通过构建树稀疏反变换得到最后的恢复图像。本发明利用充分挖掘图像信号的内部结构关系,将图像块的广义树稀疏结构特性与权重核范数结合，并利用ADMM算法简化计算过程，降低算法复杂度，提高了部分空间数据重建图像的性能，在更少的扫描测量下更精确地重建图像，减少重建图像的伪影，实现快速磁共振成像。

Description

一种基于广义树稀疏的权重核范数磁共振成像重建方法

技术领域

本发明涉及医学影像处理技术领域，尤其是指一种基于广义树稀疏的权重核范数磁共振成像(MRI)重建方法，主要用于对医学图像的清晰快速恢复，减少重建图像的伪影，恢复更多的图像细节。

背景技术

磁共振成像MRI,因其低损害性及高诊断意义而得到医学界的广泛应用,磁共振MRI建立在磁共振原理的基础上。磁共振涉及的基本物理概念主要包括原子的自旋和磁矩,自旋磁矩在外磁场中的能量状态,产生磁共振的条件,拉莫进动,磁化强度矢量,以及射频场对磁化强度矢量和弛豫过程。

传统的磁共振成像MRI要对原始数据按照奈奎斯特采样定理进行密集采样,然后对所采集密集数据通过逆傅立叶变换重建磁共振MRI图像,这将导致重建磁共振成像所需采样量大的难题。压缩感知理论的提出允许对原始数据空间进行欠采样，使采集的原始数据空间样本可以远小于传统磁共振MRI成像中所需采集的原始数据空间样本数，从而大大减少采样所需数据及节约采样时间。而如何从欠采样的数据样本重构出高清晰的磁共振MRI图像是压缩感知磁共振MRI成像方法成功的一个关键因素,也是近年来研究的热点。

现有磁共振压缩感知CS‐MRI图像重建方法，是利用磁共振MRI图像的稀疏性来重建磁共振MRI图像，现阶段有不少学者通过的选取多利用全变差、小波基、以及自适应字典学习等对磁共振MRI图像进行稀疏表示。相比全变差、以及小波基，基于自适应字典学习的稀疏表示有明显的优越性，能更好地刻画磁共振MRI图像中的边缘结构，但是现有的这种利用自适应稀疏字典表示方法在实际MRI图像重建过程中存在一定的问题，如并没有考虑到信号的内部结构性，它忽略了信号本身内部的结构关系。

此外,从欠采样的原始数据K稀疏恢复磁共振MRI图像的稀疏表示系数是一个病态逆问题。通过对充分挖掘信号内部的结构性可以有效提升图像稀疏分解的精度，这种结构稀疏模型已被证明是一个更加鲁棒、更精确的稀疏模型。但是,现有CS‐MRI图像重建方法由于主要是利用磁共振MRI图像的稀疏性实现对图像重建,并没有考虑到信号的内在结构稀疏性，因而难以精确重建出原始真实的磁共振MRI图像,导致医学诊断困难。

发明内容

本发明的目的在于克服现有技术的缺点和不足，提供一种基于广义树稀疏的权重核范数磁共振成像(MRI)重建方法，利用树稀疏结构的低采样率数据，建立基于权重核范数逼近的MRI重建数学模型，采用交替方向迭代方法对模型进行迭代求解,简化计算过程，降低算法复杂度；本发明对欠采样的数据经过树稀疏结构变换后，可以用更少的稀疏系数和更少的计算量来重构图像，因此进一步地稀疏表示图像，并更好地利用图像的局部信息和非局部相似性，得到更好的图像重建效果；在更少的测量数据下得到更精确的重建图像，减少重建图像的伪影，更加适用于实际场景。

为实现上述目的，本发明所提供的技术方案为：一种基于广义树稀疏的权重核范数磁共振成像重建方法，首先，获取测试磁共振成像采样数据样本进行傅立叶变换；然后根据采样的信号构建出树结构的稀疏信号，利用带权重的核范数逼近约束的目标函数的稀疏表达；接着再通过增广拉格朗日乘子法优化约束的目标函数和交替方向搜索算法对所述测试数据进行迭代更新，直至得到估计的恢复数据；最后再通过构建树稀疏反变换得到最后的恢复图像；其包括以下步骤：

1)采样数据是根据临床志愿者提供的，输入磁共振MRI的笛卡尔欠采样数据或非笛卡尔欠采样数据，这个由具体型号的设备参数和被采集对象决定；

2)将欠采样数据进行相应的傅立叶变换，并通过变换后的数据进行相应的树稀疏结构模型构建；

3)采用自适应奇异值分解求解得到每一分类对应的权重核范数逼近；

4)采用交替方向搜索(ADMM)算法，迭代更新直到迭代次数达到预设的阈值，得到最终每一分类对应的稀疏向量；

5)对树稀疏向量进行反变换最终求解得到对应的稀疏向量重建MRI图像。

所述步骤2)包括以下步骤：

2.1)根据MRI设备获取MR欠采样数据y；

2.2)进行相应的傅立叶变换，即：θ＝Fy；

2.3)构建树稀疏结构模型G，即：Gθ。

在步骤3)中，所述权重核范数逼近是一种新的范数模式，既具有l₀范数的全局性，又具有l₀范数的准确性，它的目标是寻找混合范数最优解，同时满足采样约束条件，其相应的模型：

其中，σ_i(x)表示x的第i个奇异值，w_i为对应的权重；

原子范数优化问题：

这里σ_i(x)为对应的奇异值，λ为辅助参量；

基于树结构的权重核范数最小化的求解能够转换成一个无约束的最小化问题：

其中，y为欠采样的数据，w表示权重，*表示原子范数以及A为欠定观测矩阵；

该优化问题等价于增广拉格朗日乘子法无约束问题：

其中，z，ρ是辅助参量，d为拉格朗日乘子。

在步骤4)中，结合交替方向搜索算法，优化分解为：

d^k+1＝d^k-(Gθ-z+d^k)

轮换地更新一个变量，同时固定其它变量：固定辅助变量z^k，求解x^k+1；固定图像x^k，通过迭代软阈值方法更新z^k+1；

对于估算出的数据z的子问题，通过迭代软阈值方法可解，如下：

S_W(z^k+1)＝max(σ_i(x_i)-w_i,0)

其中，S_W(z^k+1)表示软阈值，σ_i(x)表示x的第i个奇异值，w_i为对应的权重；

而在步骤5)，通过迭代优化得到的数据，进行稀疏信号树结构反变换就求解重构图像。

本发明与现有技术相比，具有如下优点与有益效果：

1、本发明由于在图像重建中,利用构建广义树稀疏的信号结构，能更加精确恢复磁共振MRI图像。

2、本发明由于结合ADMM算法，该算法对正则参数自适应的选取,使得重构方法更加鲁棒性，恢复速度更快。

3、本发明提出了利用权重核范数逼近的方法来求解以更优于1-范数问题，使得重构图像更加精确，呈现更少的混叠伪影,恢复的图像有更好的保真度。

附图说明

图1为本发明方法的流程示意图。

图2为本发明方法与其他五种方法的恢复精度对比效果图。

图3为不同的迭代次数变化情况下，方法恢复出来的信噪比值示意图。

具体实施方式

下面结合具体实施例对本发明做进一步的说明。

如图1所示，本实施例所述的基于广义树稀疏的权重核范数磁共振成像(MRI)重建方法，具体是：首先，获取测试磁共振成像(MRI)采样数据样本进行傅立叶变换；然后根据采样的信号构建出树结构的稀疏信号，利用带权重的核范数逼近约束的目标函数的稀疏表达；再通过增广拉格朗日乘子法优化约束的目标函数和交替方向搜索(ADMM)算法对所述测试数据进行迭代更新，直至得到估计的恢复数据；最后，通过构建树稀疏反变换得到最后的恢复图像。其包括以下步骤：

1)采样数据是根据临床志愿者提供的，输入磁共振MRI的笛卡尔欠采样数据，也可以是非笛卡尔采样，这个由具体型号的设备参数和被采集对象决定。

2)将欠采样数据进行相应的傅立叶变换，并通过变换后的数据进行相应的树稀疏结构模型构建(由于是在欠采样数据上进行，它的规模远小于全采样数据,其数据量取决于采样率和采样模式)，具体如下：

2.1)根据MRI设备获取MR欠采样数据y；

2.2)进行相应的傅立叶变换，即：θ＝Fy；

2.3)构建树稀疏结构模型G，即：Gθ。

3)采用自适应奇异值分解求解得到每一分类对应的权重核范数逼近，如下：

在CS‐MRI重建中，l₁范数最小化是求解非线性凸优化问题,它能确保重构解的稳定性,但得到的解不是最优解。l₀范数解才是我们要求的理想最优解，然而l₀范数问题是NP‐hard问题，很难直接求解。所以,如何确保解的准确性和收敛性是重构面临的重要问题。而本发明提出的权重核范数逼近是一种新的范数模式，既具有l₀范数的全局性，又具有l₀范数的准确性，它的目标是寻找混合范数最优解，同时满足采样约束条件，其相应的模型：

其中，σ_i(x)表示x的第i个奇异值，w_i为对应的权重；

原子范数优化问题：

这里σ_i(x)为对应的奇异值；

基于树结构的权重核范数最小化的求解能够转换成一个无约束的最小化问题：

其中，y为欠采样的数据，w表示权重，*表示原子范数以及A为欠定观测矩阵；

该优化问题等价于增广拉格朗日乘子法无约束问题：

其中，z，ρ是辅助参量，d为拉格朗日乘子。

4)采用交替方向搜索算法，迭代更新直到迭代次数达到预设的阈值，得到最终每一分类对应的稀疏向量，如下：

结合ADMM算法，优化分解为：

d^k+1＝d^k-(Gθ-z+d^k)

轮换地更新一个变量，同时固定其它变量：固定辅助变量z^k，求解x^k+1；固定图像x^k，通过迭代软阈值方法更新z^k+1。

对于估算出的数据z的子问题，通过迭代软阈值方法可解，如下：

S_W(z^k+1)＝max(σ_i(x_i)-w_i,0)

其中，S_W(z^k+1)表示软阈值，σ_i(x)表示x的第i个奇异值，w_i为对应的权重；

5)对树稀疏向量进行反变换最终求解得到对应的稀疏向量重建MRI图像。其中通过迭代优化，并将得到的稀疏信号树结构反变换就求解重构图像。

本发明效果通过以下实验进一步说明：

实验条件：本发明采用标准核共振MRI上采样数据进行算法比较；所用编程平台为MatlabR2010b；其欠采样率为30％，经过恢复后的核共振成像的的质量通过使用由下式来衡量:

其中,为x原图像,为恢复出来的图像。

图2中(a)‐(e)分明RecPF方法、FCSA方法,WaTMRI方法、,NESTA Tree方法和本发明对人脑图像重构仿真,(f)‐(i)分别为对应RecPF方法、FCSA方法,WaTMRI方法、,NESTA Tree方法和本发明对人脑图像重构后丢失的像素。

图2所显示的人脑图像的恢复结果可以看出,本发明的基于广义树稀疏的权重核范数磁共振成像(MRI)重建方法恢复出来的图像比其他方法恢复出来的图像更干净,清晰,图像边缘更锐利,恢复更多的图像细节。

如图3所示，从图中可以看出不同的迭代次数变化情况下，几种同类方法恢复出来的信噪比值随迭代次数的变化情况。特别地，随着迭代次数的增加，本发明所采用的方法具有良好的恢复稳定性。进一步说明本发明方法重建结果可以获得更清晰的图像对比度以及更精确的解剖结构描述,恢复的图像有更好的保真度。

以上所述实施例只为本发明之较佳实施例，并非以此限制本发明的实施范围，故凡依本发明之形状、原理所作的变化，均应涵盖在本发明的保护范围内。

Claims (3)

1.一种基于广义树稀疏的权重核范数磁共振成像重建方法，其特征在于：首先，获取测试磁共振成像采样数据样本进行傅立叶变换；然后根据采样的信号构建出树结构的稀疏信号，利用带权重的核范数逼近约束的目标函数的稀疏表达；接着再通过增广拉格朗日乘子法优化约束的目标函数和交替方向搜索算法对所述测试数据进行迭代更新，直至得到估计的恢复数据；最后再通过构建树稀疏反变换得到最后的恢复图像；其包括以下步骤：

1)采样数据是根据临床志愿者提供的，输入磁共振MRI的笛卡尔欠采样数据或非笛卡尔欠采样数据，这个由具体型号的设备参数和被采集对象决定；

2)将欠采样数据进行相应的傅立叶变换，并通过变换后的数据进行相应的树稀疏结构模型构建；

3)采用自适应奇异值分解求解得到每一分类对应的权重核范数逼近；

4)采用交替方向搜索算法，迭代更新直到迭代次数达到预设的阈值，得到最终每一分类对应的稀疏向量；

5)对树稀疏向量进行反变换最终求解得到对应的稀疏向量重建MRI图像。

2.根据权利要求1所述的一种基于广义树稀疏的权重核范数磁共振成像重建方法，其特征在于，所述步骤2)包括以下步骤：

2.1)根据MRI设备获取MR欠采样数据y；

2.2)进行相应的傅立叶变换，即：θ＝Fy；

2.3)构建树稀疏结构模型G，即：Gθ。

3.根据权利要求1所述的一种基于广义树稀疏的权重核范数磁共振成像重建方法，其特征在于：在步骤3)中，所述权重核范数逼近是一种新的范数模式，既具有范数的全局性，又具有范数的准确性，它的目标是寻找混合范数最优解，同时满足采样约束条件，其相应的模型：

$\left|\right|x|{|}_{w,*}=\underset{i}{\Σ}{w}_i{\mathrm{\σ}}_i\left(x\right)$

其中，σ_i(x)表示x的第i个奇异值,，w_i表示对应的权重；

原子范数优化问题：

这里σ_i(x)为对应的奇异值；

基于树结构的权重核范数最小化的求解能够转换成一个无约束的最小化问题：

$\underset{x\∈R^N}{min}\frac{1}{2}\left|\right|y-Ax|{|}_2^2+|\left|G\mathrm{\θ}\right|{|}_{w,*}$

其中，y为欠采样的数据，w表示权重，*表示原子范数以及A为欠定观测矩阵；

该优化问题等价于增广拉格朗日乘子法无约束问题：

$L(x,z,d)=\frac{1}{2}\left|\right|y-Ax|{|}_2^2+|\left|G\mathrm{\θ}\right|{|}_{w,*}+d(G\mathrm{\θ}-z)+\frac{\mathrm{\ρ}}{2}\left|\right|G\mathrm{\θ}-z|{|}_2^2$

其中，z，ρ是一个辅助参量，d为拉格朗日乘子。

在步骤4)中，结合交替方向搜索算法，优化分解为：

$\underset{x\∈R^N}{x^{k+1}=\arg \min }=\frac{\mathrm{\ρ}}{2}\left|\right|G\mathrm{\θ}-z^k+d^k|{|}_2^2+\frac{1}{2}||{\mathrm{Ax}}^k-y|{{|}^2}_{w,*}$

$\underset{z\∈R^N}{z^{k+1}=\mathrm{argmin}}\frac{\mathrm{\ρ}}{2}\left|\right|G\mathrm{\θ}-z+d^k|{|}_2^2+\mathrm{\λ}|\left|G\mathrm{\θ}\right|{|}_{w,*}$

d^k+1＝d^k-(Gθ-z+d^k)

轮换地更新一个变量，同时固定其它变量：固定辅助变量z^k，求解x^k+1；固定图像x^k，通过迭代软阈值方法更新z^k+1；

对于估算出的数据z的子问题，通过迭代软阈值方法可解，如下：

S_W(z^k+1)＝max(σ_i(x_i)-w_i,0)

其中，S_W(z^k+1)表示软阈值，σ_i(x)表示x的第i个奇异值，w_i为对应的权重；

而在步骤5)，通过迭代优化得到的数据，进行稀疏信号树结构反变换就求解重构图像。