一种加权冗余稀疏约束磁共振图像的迭代重建方法与流程

本发明涉及磁共振图像重建方法，尤其是涉及利用压缩感知技术加快磁共振成像时的一种加权冗余稀疏约束磁共振图像的迭代重建方法。

背景技术：

磁共振成像被广泛应用于医学临床诊断中，但在传统的磁共振成像中，获取完整的k空间数据耗时较长，对受采集者不友好，其次成像速度也制约着图像的质量。非均匀采样和压缩感知满足快速成像的迫切需求，压缩感知通过约束图像的稀疏性，只需要采集少量的k空间数据，就可以恢复出完整的k空间数据。设计基于压缩感知的磁共振图像的重建方法关键在于磁共振图像的稀疏表示。磁共振图像的稀疏表示主要有两种类型：正交稀疏表示(m.lustig,d.donoho,andj.m.pauly,"sparsemri:theapplicationofcompressedsensingforrapidmrimaging,"magneticresonanceinmedicine,vol.58,no.6,1182-1195,2007)和冗余稀疏表示(x.qu,d.guo,b.ning,y.hou,y.lin,s.cai,andz.chen,"undersampledmrireconstructionwithpatch-baseddirectionalwavelets,"magneticresonanceimaging,vol.30,no.7,964-977,2012)、(x.qu,y.hou,f.lam,d.guo,j.zhong,andz.chen,"magneticresonanceimagereconstructionfromundersampledmeasurementsusingapatch-basednonlocaloperator,"medicalimageanalysis,vol.18,no.6,843-856,2014)、zhan等(z.zhan,j.cai,d.guo,y.liu,z.chen,andx.qu,"fastmulticlassdictionarieslearningwithgeometricaldirectionsinmrireconstruction,"ieeetransactionsonbiomedicalengineering,vol.63,no.9,1850-1861,2016)；并且，y.liu等人在文献(y.liu,z.zhan,j.-f.cai,d.guo,z.chen,andx.qu,"projectediterativesoft-thresholdingalgorithmfortightframesincompressedsensingmagneticresonanceimaging,"ieeetransactionsonmedicalimaging,vol.35,no.9,2130-2140,2016)中证明了正交稀疏表示有利于理论分析、快速重建算法设计但是稀疏表示能力有限，冗余稀疏表示能更好的捕捉图像的细节，从而更好的移除图像噪声和伪影。针对上述两种类型的稀疏表示，目前有两种不同的磁共振图像稀疏重建模型：综合型重建模型和分解型重建模型，两种重建模型均是用l1范数作为l0范数的一个凸近似，约束图像的稀疏性(y.liu,j.-f.cai,z.zhan,d.guo,j.ye,z.chen,andx.qu,“balancedsparsemodelfortightframesincompressedsensingmagneticresonanceimaging,”plosone,vol.10,no.4,p.e0119584,2015)。m.salmanasif等人(m.salmanasifandjustinromberg,"fastandaccuratealgorithmsforre-weightedl1-normminimization,"ieeetransactionsonsignalprocessing,vol.61,no.23,5905-5916,2013)提出在正交稀疏表示下用加权l1范数来约束图像的稀疏性，得到了比l1范数稀疏约束更好的重建结果，但是在冗余稀疏表示下的加权l1范数稀疏约束重建还未见有报道。

技术实现要素：

本发明的目的在于提供利用压缩感知技术加快磁共振成像时的一种加权冗余稀疏约束磁共振图像的迭代重建方法。

本发明包括以下步骤：

1)提出加权冗余稀疏约束重建模型；

在步骤1)中，所述加权冗余稀疏约束重建模型的具体方法为：测量得到的欠采样的k空间数据y表示为y＝ufx，其中，x表示待重建的磁共振图像经元素重排后的列向量，f表示傅里叶变换，u表示欠采样矩阵；用ψ表示对重排后的磁共振图像的冗余稀疏变换，ψ^*表示ψ的复共轭，φ＝ψ^*表示ψ的反变换；待重建的磁共振图像x在冗余稀疏变换ψ的表示下的稀疏变换域系数表示为α＝ψx，此外通过引入一个权值矩阵w＝diag(w1,w2,...,wn)衡量稀疏变换域系数α中每个元素的重要性，其中diag(w1,w2,...,wn)表示主对角线元素为wi(i＝1,2,...,n)，的对角阵，wi(i＝1,2,...,n)是权值，用来衡量稀疏变换域系数α中第i(i＝1,2,...,n)个元素的重要性；权值可以通过求解未加权的冗余稀疏约束重建模型先得到预重建结果在冗余稀疏变换ψ下的稀疏表示系数其中表示的第i个元素，于是wi可由计算得到，其中表示的第i个元素的模值，ε表示一个阈值；提出的加权冗余稀疏约束重建模型为：

其中，||·||1表示向量的1范数，表示向量的2范数的平方，正则化参数λ用于权衡||wψx||1和两项的重要性。

2)通过迭代运算求解重建模型获得重建磁共振图像；

在步骤2)中，所述通过迭代运算求解重建模型获得重建磁共振图像的具体方法可为：待重建的磁共振图像x用稀疏反变换φ与变换域系数α的乘积来表示，即x＝φα；x按照公式(2)的迭代运算来获得：

其中，xn和xn+1是在第n和n+1次迭代时得到的图像，u^t表示u的转置，f^*表示反傅里叶变换，γ是步长，其取值范围为0＜γ≤1；步长γ与正则化参数λ、权值wi的乘积为γλwi，表示对给定向量z＝(ψ(xn+γf^*u^t(y-ufxn)))的所有元素进行软阈值操作，软阈值操作的定义为：给定向量z，若向量z的第i个元素zi的绝对值|zi|≤γλwi，则zi＝0，若|zi|＞γλwi，则zi＝sgn(zi)(|zi|-γλwi)，其中sgn(zi)是符号函数；

3)在每次迭代中利用相邻两次迭代的重建结果对磁共振图像进行修正，用于加快磁共振图像的重建速度。

在步骤3)中，所述加快磁共振图像的重建速度的具体方法可为：引入一组因子{tn}，其中n是大于等于0的整数，初始化的因子t0＝1，已知这组因子中的第n项因子tn，通过表达式可以计算第n+1项因子tn+1；在第n+1次迭代重建图像后，利用第n和n+1次迭代时得到的图像xn和xn+1对迭代的图像进行修正，从而加快图像重建的收敛速度。修正的方法是：

公式(3)中的符号表示将公式(3)中右边的结果赋值给公式(3)左边的xn+1中；

迭代过程中，当相邻两次的重建图像xn和xn+1的2范数误差||xn+1-xn||2/||xn||2缩小到指定的范围时，停止迭代运算；将最后一步迭代的重建结果xn+1输出为最终的重建结果图像。

本发明的有益效果是：图像重建的参数少，仅引入步长γ、正则化参数λ和阈值ε；重建的速度快，重建图像质量高。

附图说明

图1为迭代过程中的人脑图像和全采样人脑图像的相对l2范数误差变化曲线。

图2为全采样的人脑图像、采样率为30％的高斯采样模板、前沿方法的重建图像、本发明重建的人脑图像以及相应的与全采样的人脑图像之间的误差图。在图2中，(a)是全采样的人脑图像，(b)是前沿的投影快速迭代阈值算法重建出来的人脑图像(y.liu,z.zhan,j.-f.cai,d.guo,z.chen,andx.qu,"projectediterativesoft-thresholdingalgorithmfortightframesincompressedsensingmagneticresonanceimaging,"ieeetransactionsonmedicalimaging,vol.35,no.9,2130-2140,2016)，(c)是本发明重建出来的人脑图像，(d)是采样率为30％的高斯采样模板，(e)是前沿的投影快速迭代阈值算法重建出来的人脑图像与全采样的人脑图像之间的误差图，(f)是本发明重建出来的人脑图像与全采样的人脑图像之间的误差图。

具体实施方式

以下实施例将结合附图对本发明作进一步的说明。

本发明实施例使用磁场强度为3特斯拉的磁共振成像仪对志愿者的大脑进行成像。本实施例使用的采样序列为t2加权快速自旋回波序列，序列的回波时间te＝99ms，重复时间tr＝6100ms，选层厚度3mm，视野fov＝220mm×220mm。

具体实施过程如下：

第一步：提出加权冗余稀疏约束重建模型

对于实施例中待重建的人脑磁共振图像经元素重排后的列向量x，使用如图2.(d)所示的采样模板对待重建的人脑磁共振图像x对应的k空间数据fx进行欠采样得到欠采样的k空间数据y。欠采样的k空间数据y与待重建的人脑磁共振图像x之间的关系可以表示为：y＝ufx。其中f表示傅里叶变换，u表示欠采样矩阵。在重建过程中采用sidwt(移不变小波变换)ψ(在本实施例中选用多贝西小波基，分解级数选为4)，对待重建的人脑磁共振图像x进行冗余稀疏表示，得到的稀疏变换域系数α可以表示为α＝ψx。此外通过引入一个权值矩阵w＝diag(w1,w2,…,wn)来衡量稀疏变换域系数α中每个元素的重要性,其中diag(w1,w2,…wn)表示主对角线元素为wi(i＝1,2,…,n)，的对角阵，wi(i＝1,2,…,n)是权值，用来衡量稀疏变换域系数α中第i(i＝1,2,…,n)个元素的重要性。权值可以通过求解未加权的冗余稀疏约束重建模型先得到预重建结果在冗余稀疏变换ψ下的稀疏表示系数其中表示的第i个元素，于是wi可由计算得到，其中表示的第i个元素的模值，在本实例中阈值ε取值为0.001。提出的加权冗余稀疏约束重建模型为：

其中，||·||1表示向量的1范数，表示向量的2范数的平方，正则化参数λ用于权衡||wψx||1和两项的重要性。实例人脑图像的重建中，设定正则化参数λ＝0.0001。

第二步：重建图像

人脑图像x通过公式(5)的迭代运算来获得

其中xn和xn+1是在第n和n+1次迭代时得到的图像，u^t表示u的转置，f^*表示反傅里叶变换，γ是迭代步长，实施例中设置迭代步长γ＝1。步长γ与正则化参数λ、权值wi的乘积为γλwi。表示对给定向量z＝(ψ(xn+γf^*u^t(y-ufxn)))的所有元素进行软阈值操作，软阈值操作的定义为：给定向量z，若向量z的第i个元素zi的绝对值|zi|≤γλwi，则zi＝0，若|zi|＞γλwi，则zi＝sgn(zi)(|zi|-γλwi)，其中sgn(zi)是符号函数。

第三步：加快图像重建的速度

在第n+1次迭代时，利用第n和n+1次迭代时得到的图像xn和xn+1对迭代的图像进行修正，加快迭代重建的收敛速度。修正的方法是

第n项因子tn和第n+1项因子tn+1满足关系式其中初始的第n项的因子t0＝1。公式(6)中的符号表示将公式(6)中右边的结果赋值给公式(6)左边的xn+1中。

迭代过程中，当相邻两次的重建图像xn和xn+1的2范数误差||xn+1-xn||2/||xn||2小于10^-4时重建图像稳定，停止迭代运算。将最后一步迭代的重建结果xn+1输出为最终的重建结果图像。

迭代重建过程中的重建人脑图像和全采样人脑图像的相对l2范数误差变化如图1所示，从图1可以看出4s后，相对l2范数误差的值趋于稳定，不再下降。前沿的投影快速迭代阈值算法重建出来的人脑图像如图2(b)所示，与全采样的人脑图像(如图2(a)所示)之间的误差图如图2(e)所示，相对l2范数误差为0.1044；本发明重建出的人脑图像如图2(c)所示，与全采样的人脑图像之间的误差图如图2(f)所示，相对l2范数误差为0.0807。

可以看出，利用本发明可以快速重建出高质量的磁共振图像，且重建结果优于前沿的投影快速迭代阈值算法(y.liu,z.zhan,j.-f.cai,d.guo,z.chen,andx.qu,"projectediterativesoft-thresholdingalgorithmfortightframesincompressedsensingmagneticresonanceimaging,"ieeetransactionsonmedicalimaging,vol.35,no.9,2130-2140,2016)。