基于压缩感知的磁共振图像重构算法研究

《基于压缩感知的磁共振图像重构算法研究》由会员上传分享，免费在线阅读，更多相关内容在学术论文-天天文库。

分类号：密级：UDC：编号：河北工业大学硕士学位论文基于压缩感知的磁共振图像重构算法研究论文作者：王晓云学生类别：全日制学科门类：工学学科专业：通信与信息系统指导教师：马杰职称：教授资助基金项目：国家自然科学基金(61203245)、河北省自然科学基金(F2012202027) DissertationSubmittedtoHebeiUniversityofTechnologyforTheMasterDegreeofElectronicsandCommunicationEngineeringTheMagneticResonanceImageReconstructionMethodsbasedonCompressedSensingbyWangXiaoyunSupervisor:Prof.MaJieMarch2016ThisworkwassupportedbyNationalNaturalScienceFoundation(61203245),HebeiProvinceNaturalScienceFoundation(No.F2012202027) 摘要压缩感知是近年来提出的一种新的信号处理理论，其特点是在采样率很低的情况下仍然能够精确的重建信号，自提出以来就得到了广泛的关注并成为研究的热点。目前压缩感知理论应用于许多领域，其中最为经典的应用之一就是将压缩感知理论应用于磁共振成像，尤其是在动态磁共振成像的应用中。磁共振成像是对通过K空间数据进行下采样然后经过一系列数据处理最后呈现的加权图像。在医疗实际中扫描获取K空间数据时，为了得到准确的图像被检测者往往需要保持身体的静止，而长时间的保持一个状态会产生疲劳甚至抖动。这些不自主的抖动以及器官的运动都会导致运动伪影。因此如何缩短扫描时间快速得到数据是一个重点。将压缩感知理论与磁共振成像相结合，只利用极少量的K空间数据便可以重建图像，这样就大大减少了扫描时间。对于动态磁共振成像，帧与帧之间高度相关，背景部分变化较少，运动的部分变化较为剧烈，这就类似于视频图像，因此将用于视频背景建模的低秩矩阵恢复理论应用于动态磁共振成像。首先将下采样得到的不同时间帧的K空间数据展开成一个列向量，再由这些列向量组合成为一个矩阵。这就构造出了一个满足低秩特性的矩阵。对这个矩阵进行低秩稀疏分解，便可以分别重建出原图像的稀疏部分和低秩部分。利用低秩稀疏分解模型进行动态磁共振图像重建，本文首先介绍了应用迭代软阈值算法(IST)进行重建，并针对其存在重建精度一般，重建速度慢的问题，一方面引入了分别利用加速近端梯度法（APG）和非精确增广拉格朗日算法（IALM）求解，推导出求解模型和步骤，并通过实验对比验证三种算法的优劣。另一方面，为了进一步提高重建精度，本文在矩阵低秩稀疏分解模型的基础上引入全变分(TV)正则项，达到进一步去噪声和去伪影的目的。通过对心脏灌注动态MRI成像和心电影MRI成像的仿真实验表明：对于低秩稀疏矩阵分解模型的重建，APG算法和IALM算法比IST算法速度更快，精度更高，IALM算方更优于APG算法；模型引入TV正则项后再利用IALM算法重建，重建速度虽然比之前的IALM算法有所降低，但依然远远优于IST算法，并且重建精度高于之前的IALM算法。关键词：压缩感知；低秩矩阵恢复；稀疏表示；动态MRIⅠ ABSTRACTCompressedsensingisanewsignalprocessingtheorywhichisputforwardrecentyears.Itischaracterizedbyaccuratelyreconstructingsignalundertheconditionoflowsamplingrate.Sinceproposedithasbeenreceivedwidespreadattentionandbecomehotspots.Recently,compressedsensinghasappliedtomanyfields;oneofthemostclassicapplicationsisthatitwasappliedtomagneticresonanceimaging,especiallyinthedynamicmagneticresonanceimaging.MagneticresonanceimagingwasconductedbyKspacedatasamplingandthroughaseriesofdataprocessingfinallypresentingtheweightedimages.Inmedicalpractice,inordertogetaccurateimage,patientsoftenneedtokeepbodystaticwhilescanning.Butmaintainastateforalongtimewillproducefatigueorevenshaking.Thisinvoluntaryjitterwillleadstoartifacts.Sohowtoshortenthescantimeandacceleratethedataacquisitionismostimportant.ThecombinationofcompressedsensingandmagneticresonanceimagingwhichusingonlyatinyamountofKspacedatatoreconstructedcangreatlyreducethescanningtime.Fordynamicmagneticresonanceimaging,therearehighlycorrelatedamongframes.Thebackgroundcomponenthaslittlechange,whilethedynamiccomponentisrapidlychangingovertime.Thisissimilartothevideosequence,sothelowrankmatrixtheorywhichisusedforvideomodelingissuitablefordynamicmagneticresonanceimagingreconstruction.Consideringeachtemporalframeasacolumnofaspacetimematrix,throughwhichweproducealow-rankmatrix.Bytakinglow-rankandSparsematrixdecompositionwecanreconstructtheoriginalimageaswellasthelow-rankandsparsecomponents.Usinglow-rankandSparsematrixdecompositionmodelfordynamicmagneticresonanceimagereconstruction,thispaperfirstintroducetheapplicationofIterativeSoftThreshold(IST)algorithmfordynamicMagneticResonanceImagingreconstruction.Consideringthismethodexiststwoproblems:thepoorreconstructionaccuracyandconvergespeed,twomeasureshavebeentaken.Ontheonehand,Acceleratetheproximalgradientmethod(APG)andInexactAugmentedLagrangianMethod(IALM)havebeenintroducedforfastreconstruction.Ontheotherhand,inordertoachievehighprecision,TotalVariation(TV)regularizationbasedonlow-rankandsparsedecompositionmodelhasbeenintroducedtofurtherremovenoiseandartifacts.Inordertoverifytheeffectivenessoftheproposedmethods,thesimulationaboutreconstructionofcardiacperfusionMRIandcardiaccineMRIhasbeendone.Theresultofthesimulationdemonstratesthat:usingAPGandIALMforLow-rankandSparsematrixdecompositionreconstructionresultinhigherreconstructionⅢ performancethanISTmethod,besides,IALMismoresuperiortotheAPGalgorithm.AlthoughintroducingTVregularizationreducesthereconstructionspeedthanIALM,itstillfastthanISTandthereconstructionaccuracyishigherthanIALM.KEYWORDS:compressedsensing,low-rankmatrixcompletion,sparsity,dynamicMRIⅣ 目录第一章绪论.....................................................................................................................................-1-1.1CS-MRI研究的目的及意义.....................................................................................................-1-1.2国内外研究现状及发展趋势..................................................................................................-2-1.3本文主要研究内容及结构安排.............................................................................................-4-1.3.1主要研究内容.....................................................................................................-4-1.3.2全文结构安排.....................................................................................................-5-第二章CS-MRI图像重建理论........................................................................................-7-2.1压缩感知理论...............................................................................................................-7-2.2CS-MRI图像重建原理................................................................................................-8-2.2.1稀疏表示.............................................................................................................-9-2.2.2采样模式............................................................................................................-11-2.2.3重建算法...........................................................................................................-14-2.3本章小结....................................................................................................................-16-第三章低秩矩阵恢复理论.............................................................................................-17-3.1低秩矩阵恢复理论....................................................................................................-17-3.2低秩矩阵恢复算法....................................................................................................-17-3.2.1迭代阈值算法..................................................................................................-18-3.2.2加速近端梯度法..............................................................................................-19-3.2.3增广拉格朗日法..............................................................................................-21-3.3低秩矩阵恢复的应用................................................................................................-24-3.4本章小结....................................................................................................................-25-第四章基于低秩稀疏分解模型的动态CS-MR图像重建.............................................-27-4.1动态CS-MRI图像重建低秩稀疏分解模型............................................................-27-4.2基于L+S分解模型的CS-MRI重建算法...............................................................-28-4.2.1IST算法............................................................................................................-28-4.2.2APG算法..........................................................................................................-30-Ⅴ 4.2.3IALM算法.......................................................................................................-32-4.3实验结果及分析........................................................................................................-36-4.4本章小结....................................................................................................................-41-第五章改进的L+S分解模型及重建算法.....................................................................-43-5.1全变分去噪................................................................................................................-43-5.2联合TV正则化MRI图像重建...............................................................................-45-5.2.1引入TV正则项的L+S重建模型及重建算法..............................................-45-5.3实验结果及分析........................................................................................................-46-5.4本章小结..................................................................................................................-50-第六章总结与展望.........................................................................................................-51-参考文献...........................................................................................................................-53-攻读学位期间所取得的相关科研成果...........................................................................-57-致谢...........................................................................................................................-59-Ⅵ 第一章绪论1.1CS-MRI研究的目的及意义磁共振成像(MagneticResonanceImaging,MRI)是上世纪80年代发展起来的一种利用核磁共振现象进行医疗诊断的影像检查技术。利用磁场共振原理对物体的内部结构进行成像，使得人们可以不经过手术而直接取得生物器官和组织的内部详细图像。磁共振成像对人体软组织和脑组织有较高的分辨力，且多方位成像，它能够直接扫描出器官的管状面，横截面等各种体层的图像，这就减少了病人不必要的手术和其对病人的损伤。与传统的CT和X射线相比，它以射频脉冲作为成像的能量源，不存在电离辐射，对人体的伤害小。由于MRI的无创性以及图像分辨率高，并且对人体任意层面都能成像，因此目前它已经发展成为生物医学工程领域重要的技术之一。虽然MRI具有很多优点，但其仍有有待提高的地方即扫描速度慢。磁共振成像过程中时间消耗主要包括扫描获取K空间数据的时间和基于采集数据进行重建的时间。针对此问题，一方面可以通过提高硬件设备的性能来加快扫描速度。通常是通过提高梯度场的切换速度以及主磁场的场强来加快K空间数据的采集。然而高速切换的梯度线圈使得硬件系统的成本提高很多，并且如果梯度场的切换速度太快会刺激人体的肌肉和神经；此外当射频脉冲序列的密集度很高时，大量的射频能量会聚集在人体内，从而伤害到人体。另一种办法是并行成像技术[1]，通过利用多通道接收阵列线圈的不同的空间敏感度，在不改变磁场梯度的条件下进行空间编码，目的是减少相位编码的次数，此外通过设计相应的解码滤波器来恢复图像。理论上，并行成像的加速倍数等于接收线圈的个数，但由于有限的空间约束使得线圈的数量不能无限增加，并且随着线圈数量的增加，线圈之间的干扰也会增加从而导致线圈之间变得高度相关，进而造成对噪声极不稳定的结果。因此在不改变硬件设备的基础上，寻求一种快速精确重建MRI的方法很有必要。根据压缩感知理论可知，如果信号是稀疏的，则可以利用某种算法对这些极少量的信息进行计算处理进而恢复原信息。而医学图像通常满足稀疏条件，为CS-MRI提供了理论基础。本文基于压缩感知理论，对K空间的数据快速扫描，只获取少量的K空间数据，就能利用欠采样重建的方法对MRI进行重建。并且选取合适的稀疏表示、采样模式、重建模型及重构算法，就能保证精确快速的MRI重建，这克服了提高硬件性能带来的成本消耗和潜在伤害，具有广阔的研究价值和应用前景。-1- 1.2国内外研究现状及发展趋势压缩感知（CompressedSensing，CS）是由Donoho和Candes等在2004年提出的一种新的信息获取理论[2][3]。该理论突破了Nyquist采样定理，充分考虑了信号的稀疏性，实现以较低的采样率进行采样；一步到位的完成信号的采样和压缩，然后对信息进行存储和传输，而在重构端利用合适的重构算法完成对信号的近乎精确的重构。CS理论提出以来，得到了广泛的应用。其中，磁共振成像是最成功的应用之一。MRI可以看作是CS理论的一个特殊情况，因为K空间实质上是经傅里叶变换后的频率域空间，每一个K空间数据都包含了所有图像域的信息。在这种条件下，利用CS理论就可以从一小部分K空间数据中精确恢复出原图像。2007年M.Lusting,D.L.Donoho等人[4][5]结合CS理论和MRI成像，将CS理论成功应用于脑成像、心脏成像、快速三维血管造影，并取得良好的重建效果。自此M.Lusting等人的开创性研究打开了利用CS理论进行MRI图像重建的新局面。各种改进的基于CS理论的利用部分K空间数据进行重建的方法被提出并进行研究。在基于CS理论的MRI重建中，主要围绕三个问题进行研究，即稀疏性、采样模式和重构算法。为了减小重构误差即更加精确地重构信号，通常希望原始信号越稀疏越好。通常需要对原始信号进行稀疏变换，传统的稀疏变换利用预先设定好的稀疏变换模型进行稀疏变换，主要包括全差分[6]和小波变换[7][8]等。全差分变换可以控制重建图像的局部差异，小波变换强调的是奇异点的各向同性，它们都可以很好抑制由于K空间欠采样引入的伪影。但全差分和小波变换都存在一个缺点，即对于含有丰富曲线和边缘的图像，不能对其进行有效的稀疏表示，由此重建的图像会出现边缘模糊。此外，许多MRI图像在一定的变换域内是稀疏的，例如血管造影在图像域就已经是稀疏的，可以通过限差分变换进行稀疏表示；脑部图像可以在小波域稀疏表示等。因此针对含有不同特征的图像，应该用不同的稀疏变换对其进行有效表示。为了充分有效的对图像稀疏表示，大量学者进行了研究。Ravishankar等人提出一种基于K-SVD自适应字典学习的MRI重建算法DLMRI[9]。当下采样率相同时，利用DLMRI算法重建的结果要比Lustig等人的SparseMRI方法重建的结果好很多。文献[10]将Contourlet变换运用到稀疏表示中，取得较好的曲线和边缘效果。Contourlet[11]变换提供了一种局域的、多方向的、多分辨率的图像表示法，通常看作是小波变换的一种新扩展。虽然利用Contourlet变换可以保留图像的曲线和边缘信息，但是它的一个缺点是计算复杂度要比小波变换的高。同时在对图像块稀疏表示方面也进行了研究。文献[12]中提出了基于块的方向性小波变换(Patch-BasedDirectionWavelets,PBDW)来稀疏重建MRI图像。文献[13]中提出了基于非局部图像块算子的欠采样MRI图像重建算法(Patch-BasedNonlocal-2- Operator,PANO)，利用图像块的自相似性或者不同对比度图像之间的互相似性对图像进行有效的稀疏表示。总之，由于很多MRI图像在一个变换域内不是充分稀疏的。为了提高重建图像的质量，如何运用多种变换对图像的不同的结构信息进行充分的稀疏表示是一个热点研究方向。对K空间数据的采样模式设计，与CS理论中确定测量矩阵的过程相对应。K空间数据的分布情况与MRI图像质量紧密相关：MRI图像的对比度由K空间中心数据分布情况决定，而图像的空间分辨力则对应于边缘数据的分布。K空间中每个数据点到K空间中心的距离代表了空间频率的大小，距离越大则对应的空间频率就越高，相应的信号空间分辨力也越高。当采样模式是全采样时，直接对得到的数据进行标准傅立叶逆变换就可以获得待重建的图像。而当是下采样时，通过补零傅立叶逆变换后得到的是含噪声的图像。采样模式不同则影响噪声分布不同。常见的采样模式有变密度笛卡尔采样[14]、二维K空间欠采样、射线状采样轨迹[15]、螺旋采样轨迹等[16][17]。相较于笛卡尔采样，射线状采样对运动伪影的稳定性更高，当加速因子很大时，仍可进行采样，尤其对对比度较高的物体效果更好。螺旋采样可以有效的利用系统的硬件，通常应用于快速或实时成像的场合，不足之处是这种非笛卡尔采样的重建更为复杂。笛卡尔采样由于可以利用FFT快速算法进行图像恢复并且具有鲁棒性，所以目前实际应用中用的最广泛的非相干采样是变密度笛卡尔采样。在重建算法方面，通常的做法是把目标函数做近似处理，凸松弛得到一个凸函数，然后进行优化求解。目前对CS-MRI重建的算法大致有贪婪算法[18]、非凸优化、凸优化算法[19-21]等。贪婪算法的重建速度很好，缺点的是其重建精度比较低，并且精确重建的理论保证不足。并且大多数贪婪算法需要明确待重构信号稀疏度，而在实际应用中这个条件难以达到。非凸优化算法，虽然能得到最稀疏的解，但其代价是最高的算法复杂度。因此，目前研究的热点主要集中在凸优化算法上。解决这一问题的常见算法有迭代阈值法[22][23]、梯度投影法、全变分法[24][25]、同伦算法、交替方向法[26]、增广拉格朗日法[27]等。在动态磁共振成像（dynamicMRI，dMRI）中，成像过程是对同一片层不同时间的数据进行采样，以此来追踪反映感兴趣区域的图像，例如心脏灌注MRI、心电影MRI等。扫描数据时由于物理或生理上的不可抑制抖动，会造成运动伪影，因此dMRI成像往往对扫描速度要求更高。在早期的研究中，文献[28]将快速dMRI成像分成三类：根据K空间数据之间的相关性，从部分下采样的数据中恢复出缺失的数据；根据时间帧之间的相关性，在已得到的高分辨率数据的基础上进行插值估算其他时间帧缺失的数据；同时结合时间帧的相关性以及K空间数据的相关性，例如k-tBLAST/SENSE[28]便是应用了时间和空间的相关性来加速dMRI采样的。运用CS理论-3- 可以使采集数据的过程进一步加快。文献[29]介绍的“k-tSPARSE”技术，根据动态磁共振图像在图像空间满足小波变换的稀疏性，并且在K空间满足傅立叶变换的稀疏性，基于这个性质Lustig将CS理论运用到dMRI重构上，取得良好的重建效果。此后大量针对dMRI重建的研究都是基于CS理论展开的。分析dMRI图像的特点可知，时间帧之间相关性很大，一组动态MRI图像与一组视频监控图像非常类似。一般监控视频都是由固定相机拍摄的图像序列，这些图像序列的背景成分基本不变，变化较大的是活动的物体或人等。在视频背景建模中，将一组视频图像分别排成列向量并组成一个矩阵，再对此进行矩阵低秩稀疏分解，达到分离静止的背景和活动前景的目的。相应的，一组dMRI图像也具有这些特征，可以通过将同一片层不同时间帧的图像拉成一个列向量，作为空间-时间矩阵的一个列。依次抽取多个时间帧，构建一个低秩矩阵，然后将低秩矩阵恢复问题应用于动态MRI重建中[30][31]。由此，将压缩感知和低秩矩阵恢复相结合的思想为进一步提高dMRI成像速度提供了新方向。起初，学者们研究的重点是找到一种途径，使矩阵同时满足低秩性和稀疏性。不同以往，文献[32][33]采用新思路，提出将数据矩阵分解成为低秩成分和稀疏成分相叠加（Low-rankplusSparse，L+S）的形式，文献[33]中分解模型用于鲁棒主成分分析（RPCA），即从一个缺失或破损的矩阵中恢复出原始信息，可以用于视频监控、背景分离、图像对齐等，目前已成功的应用于计算机视觉领域。GaoH等人在文献[34]中初步研究了将低秩稀疏分解（L+S）运用到动态MRI重建，并在文献[35]中提出基于加速弥散加权序列的图像重建。在此基础上RicardoOtazo[36]等人将L+S模型运用到加速动态对比增强MRI重建上，并在文献[37]中基于L+S分解模型，利用IST算法对临床医学图像进行了重建，取得了很好的效果，但在重建速度和精度上仍有可观的提升空间。压缩感知理论运用到磁共振成像中，克服了传统技术在空间分辨率、时间分辨率和收敛速度以及器官运动的敏感度方面的局限性，具备高精度快速成像的特点。虽然将CS理论应用于MRI图像重建取得了很大的成功，但仍存在一些问题，目前的动态CS-MRI图像重建都是在基于某种模型、某种稀疏表示、某种采样模式的基础上重建的。如何有效的将重建模型和重建算法结合以达到快速精确重建的目的是研究的热点。1.3本文主要研究内容及结构安排1.3.1主要研究内容本文将压缩感知和低秩矩阵恢复相结合，运用矩阵低秩稀疏分解模型进行动态MRI图像重建。首先通过几个简单的例子直观的说明压缩感知理论在信号重建过程中-4- 的应用，对其基本原理进行介绍。进而从三方面展开讨论CS-MRI重建原理，总结运用压缩感知进行MRI重建的原理流程图。其次，在基于压缩感知的基础上，结合低秩矩阵恢复理论，对原始数据矩阵处理，并运用低秩稀疏分解模型进行重建。对比于原有的迭代软阈值算法（IST），本文引入了加速近端梯度法（APG）、非精确增广拉格朗日法（IALM），并通过Matlab仿真实验和评价指标比较各个算法的优劣，寻找出最优算法。最后，为了得到更好的重建精度，在矩阵低秩稀疏分解模型的基础上，引入全变分正则项对图像进一步的去噪声和去伪影处理。通过Matlab软件仿真，比较改进后的模型的重建效果，得出结论。1.3.2全文结构安排本文一共分为六章。第一章：绪论。介绍了课题研究的目的和意义，探讨了基于压缩感知的磁共振图像重建的研究现状和发展趋势。在此基础上引出本文的研究内容，并安排文章结构。第二章：CS-MRI图像重建理论。围绕压缩感知理论，根据其在磁共振图像重建中的重要作用，分别从稀疏表示、采样模式以及重建算法三方面进行介绍。并用一些简单的例子直观的描述了CS-MRI的基本原理，为后续研究内容的展开打下理论基础。第三章：低秩矩阵恢复理论。首先介绍了什么是低秩矩阵以及矩阵低秩分解的模型，接着讨论了低秩矩阵恢复的重建算法：迭代阈值算法（IT）、近端梯度算法（APG）以及增广拉格朗日算法（ALM）等。在此基础上简单介绍了低秩矩阵恢复有哪些应用，为第四章的内容做铺垫。第四章：基于低秩稀疏分解模型的动态CS-MRI图像重建。在前三章的理论基础上，分析动态磁共振图像的特点，给出了适用于动态磁共振图像重建的矩阵L+S分解模型。接着讨论了针对这个模型的重建算法。先介绍已经成功应用的迭代软阈值算法（IST），在此基础上结合低秩矩阵恢复理论，本文引入了近端梯度算法（APG）和非精确的增广拉格朗日算法（IALM），经过理论推导分别得到了利用各个算法求解该模型的过程和步骤。最后利用Matlab软件，通过对心脏灌注MRI和心电影MRI图像的重建实验来检测算法的性能。第五章：改进的L+S分解模型及重建算法。本章在第四章的基础上，为了进一步的提高dMRI图像的重建精度，提出在模型上进行改进。首先介绍了全变分（TV）去噪的特征、原理和基本数学模型，并给出了常用的算法。进而基于L+S分解模型，本文引入了TV正则项，并推导新模型下的求解流程。为了验证新模型的性能，用IALM算法分别对新模型和原模型进行求解，通过对心脏灌注MRI和心电影MRI图像的重-5- 建结果比较两者的优劣，得出结论。第六章：总结和展望。归纳总结本文的方法和结论，并对今后研究工作的重点和方向进行展望。-6- 第二章CS-MRI图像重建理论2.1压缩感知理论压缩感知理论是由Candes等人提出的信号处理方法，它可以利用某种算法从非常少量的采样数据中重建得到原来的信号。压缩感知的主要内容是：如果一个高维信号是可压缩的，或者在某个变换域上是稀疏的，则通过一个与稀疏变换基不相关的测量矩阵可以将信号投影到一个低维空间上，然后基于少量的投影通过最优化求解来重建原始信号。CS理论指出从非常少量的采样数据重构信号需要满足三个条件：①图像是稀疏的或变换稀疏的；②下采样引起的失真是不相关的，它们类似于噪声；③非线性重建，利用非线性算法，加强图像稀疏性和一致性。如图2.1所示，高维信号x本身是可压缩的，或者在某个变换域ψ上具有稀疏性，其稀疏度用向量a表示，则其稀疏变换的表达式如下：x=ψa（2-1）图2.1稀疏变换映射如图2.2所示，将信号x投影到一个低维空间y上是通过一个与变换基ψ不相关的矩阵ϕ实现的。其表达是如下：y=ϕx=ϕψa=Φa（2-2）图2.2高维信号映射为低维信号-7- 式（2-2）中ψ表示稀疏变换基，ϕ表示测量矩阵，a表示稀疏系数，Φ=ϕψ表示感知矩阵。要想从低维信号y中重建原信号x，只需要计算出稀疏系数a，然后带入到（2-2）式并求解便可以得到原来的信号。其求解模型如下所示：minas.ty=Φa（2-3）0指矩阵中有几个非零元素，即零范数。求解上式，为了得到最稀疏的形式，0需要找到所有可能存在的稀疏情况，这是一个NP难问题。文献[3]指出当感知矩阵满足约束等距性[38]（RIP）时，上式可以凸松弛为下面的表达式：minas.ty=Φa（2-4）1式中表示L1范数，它是指将向量内所有的元素取绝对值后再求和。当感知1矩阵Φ满足RIP条件时，式（2-4）的解等效于式（2-3）的解。下面考虑有噪声影响时，则式（2-4）变形为（2-5），表达式如下：minas.ty−Φa≤ε（2-5）12ε是与噪声相关的常数。对于上式的求解其本质就是找到信号在稀疏变换域的表示中非零元素的位置以及数值大小，即解决定位定值的问题。求解L1范数最小化的方法有很多，例如Bregman迭代法、内点法、梯度投影法等。2.2CS-MRI图像重建原理将CS理论运用到MRI中，为进一步减少扫描时间提供了新途径。文献[4]指出，要想应用CS理论进行精确的图像重建，必须要满足下面两个条件：图像是稀疏的；测量矩阵和稀疏矩阵是不相关的。MRI恰好满足这些条件。首先，医学图像在一定的稀疏变换基下是可压缩的，例如小波变换、有限差分、学习字典等。其次，MRI数据是来自K空间（频域），而不是时域，这有利于通过随机采样生成不相关的伪随机噪声，便于信号的分离提取。另外，类似于视频图像通常比静态图像的可压缩性更强的特征，动态MRI具有广泛的时空相关性，因此CS理论更适用于动态MRI成像。类似于压缩感知中一般信号的重建模型，将CS运用于MRI的重建模型如下所示：minTms.t.Fm−y<ε（2-6）s12m：待重建图像；T：线性算子，将图像压缩到稀疏域；F：欠采样的F变换，对应于K空间的欠采样机制；sy：扫描得到的K空间数据；-8- ε：控制重建的准确性；在基于CS理论的MRI重建中，前面提到的CS理论精确重建的三个原则依然很重要，下面将从采样方式，稀疏表示和重建算法三方面进行介绍。2.2.1稀疏表示稀疏性是指信号在某种变换基下，大部分的变换系数是零，而只有很少一部分系数是不为零的。并且这些少量的不为零的大系数包含了大部分的原始信息。根据CS理论可知，重建信号的过程实质是求解稀疏系数矩阵的过程，即对稀疏系数矩阵中的非零元素进行定位定量的计算。因此，想要更快更好的重建原始信号，则希望矩阵越稀疏越好，这就对稀疏变换提出了要求。通常医学图像不是稀疏的，利用一些稀疏变换对其进行处理已达到稀疏的目的，然后在此基础上进行稀疏图像重建。常见的稀疏变换有小波变换、傅立叶变换、全变分变换、离散余弦变换等。下面通过一幅T2加权的大脑图像经小波正交变换重建的例子说明利用CS-MRI理论对稀疏图像进行重建的应用。如图2.3所示，图（a）是一幅大脑MRI图像，它本身不是稀疏的，经过小波变换后得到了图（b），小波系数的每个频带表示一个尺度的图像，带内小波系数的位置表示其在空间中的位置。图（b）表示的是图像在不同分辨率和不同方向上的边缘信息。（a）原始大脑图像（b）小波变换后的图2.3大脑MRI图像经过小波稀疏变换将变换后得到的小波系数按照降序排列，然后分别选取前1%、10%、20%和40%的小波系数对图像重建，重建结果如图2.4所示。不同的比例系数的重建结果与原图像的差异如图2.5所示。-9- （a）1%（b）10%（c）20%（d）40%图2.4利用不同比例的小波系数进行重建的结果（a）1%（b）10%（c）20%（d）40%图2.5利用不同比例的小波系数进行重建的图像与原图像之间的差异结合图2.4和图2.5可知，采用1%的小波稀疏重建得到图像比较模糊，并且跟原图像相比差异较大。进而逐步增加小波系数的比例进行重建，从图2.5的（b）-（d）可以看出，利用10%的小波系数便可以精确的进行图像重建。由上述内容可知，将图像稀疏编码，就可以利用压缩感知理论从少量的观测信号中较为精确的重构出原图像。另一方面，对于不同的医学图像，其最佳的稀疏表示是不同的，因此运用不同的稀疏表示重建结果也会不同。如图2.6和图2.7所示，分别利用离散余弦变换、小波变换和有限差分变换对T1加权的大脑图像和轴向三维腿部外围的对比增强血管造影进行稀疏重建，选取重构系数比例为1%、5%、10%、25%、35%、50%时的重建效果。1%,5%,15%,25%,35%,50%FiniteDiff.WaveletDCT图2.6基于不同稀疏表示方法进行重建的T1加权大脑图像-10- 1%,5%,15%,25%,35%,50%FiniteDiff.WaveletDCT图2.7不同稀疏表示方法下的轴向三维腿部外围的对比增强血管造影通过分析可知对于血管造影图像，由于含有丰富的边缘信息，所以用限差分做稀疏变换较好。而对于脑部的稀疏表示，限差分变换的效果则不如离散余弦变换或者小波变换的效果好。综上可知，根据不同的医学图像，要具体问题具体分析，选取有效的稀疏表示方式。2.2.2采样模式MRI图像重建是指通过下采样获得的K空间数据来重建原始图像的技术。K空间是指经Fourier变换后得到的频域空间，将K空间数据进行逆傅立叶变换就可以得到图像空间。K空间数据的分布情况与MRI图像紧密相关：MRI图像的对比度是由分布在K空间中心部分的数据决定的，而图像的空间分辨力是由K空间的边缘数据决定的。因此用不同的方式对K空间进行采样，则重建出来的图像的质量与特性也是大不相同的。由CS理论可知，希望由下采样引起的失真是不相关的，类似噪声。这对采样模式提出了要求。下面通过一个简单的一维稀疏信号重建的例子来讨论不同采样模式的重建效果。如图2.8所示，一个稀疏的一维信号（a），经Fourier变换到频域，模拟K空间得到图像（b），在频域进行采样率为四分之一的均匀采样后见（c），对采样样本进行标准的反傅立叶变换得到如（d）所示的图像。根据图（a）和图（d）可以看出重建后的信号出现混叠现象，从图（d）中无法确定原信号。由此可知针对均匀下采样方式是无法有效的重构出原信号的。-11- SparseSignalk-spaceofSparseSignal0.181.20.1610.140.80.120.60.1xX0.40.080.060.20.0400.02-0.20020406080100120020406080100120140nf（a）原一维稀疏信号（b）Fourier变换到频域空间（模拟K空间）equispacedundersamplingequispaced4-foldundersampling0.181.20.1610.140.80.120.60.1Xx0.080.40.060.20.0400.02-0.20020406080100120140020406080100120fn（c）频域（K空间）均匀采样（d）标准反傅立叶变换重建图2.8对K空间均匀采样后进行重建的结果下面讨论对K空间进行随机采样时情况会是怎样。如图2.9所示，（a）表示在频域进行随机采样，采样率仍为四分之一。对采样信号同样进行傅立叶反变换后得到图像（b），从（b）中看出除了两个大信号可以分辨出来外，另一个小信号被大信号处泄漏的能量所湮没，整体呈现类噪声的现象。在这个基础上进行非线性重建，利用阈值函数筛选出较大的信号得到图（c），计算在大信号处的能量泄漏得到图（d），再从重建的图像（b）中消除泄漏的能量，则小信号便可以显现出来，最后达到重建图像的目的。randomundersamplingrandom4-foldundersampling0.181.20.1610.140.80.120.10.6Xx0.080.40.060.20.0400.02-0.20020406080100120140020406080100120fn（a）频域（K空间）随机采样（b）利用随机采样反Fourier变换重建-12- 1.21.2110.80.80.60.6xx0.40.40.20.200-0.2-0.2020406080100120020406080100120nf（c）利用阈值函数筛选出大信号成分（d）大信号成分造成的能量泄漏图2.9对K空间随机采样并重建由此可知CS-MRI理论要求采样模式是不相关的。通常采样模式的不相关性是由点分布函数[4](PointSpreadFunction,PSF)来衡量的。令F表示欠采样的傅立叶变换，s*F表示补零的反傅里叶变换，则点分布函数定义为：s*PSF(i,j)=(FF)(i,j)（2-7）ss点分布函数是衡量零填充线性重建时，第j个像素点处受到的来自于第i个像素点的能量泄漏的影响的大小的函数。这种能量泄漏反映到图像中就是噪声干扰。在奈奎斯特采样定理下，像素之间不存在干扰因此PSF(i,j)=0，而下采样时像素之间i≠j存在干扰，点扩散函数PSF(i,j)是一个非零值。随机下采样的目的就是为了展开像i≠j素间的干扰噪声，使其均布在图像中，从而使得最大的能量泄漏趋向较小。因此用点分布函数的最大非对角元素衡量相关性的大小。一个衡量非相关性的简单方法是求旁瓣峰值比（Sidelobe–to-PeakRatio，SPR）的最大值：PSF(i,j)max（2-8）i≠jPSF(i,i)此外，有一些图像是在变换域稀疏而不是图像域稀疏的，这种情况下利用变换点分布函数（TransformPointSpreadFunction,TPSF）衡量其相关性。令T表示稀疏变换，则变换点分布函数的定义如下所示：**TPSF(i,j)=(ψFFψ)(i,j)（2-9）ssmaxTPSF(i,j)i≠j其相关性通过计算TPSF的最大非对角元素进行衡量。目前用的最广泛的采样模式是笛卡尔网格直线采样轨迹，临床中使用较多的脉冲序列是笛卡尔。由此重建很简单，只要进行IFFT变换，更重要的是笛卡尔采样模型重建精度和对系统的鲁棒性很好。其他常见的采样模式有径向采样、螺旋采样等，见-13- 图2.10。径向采样相对于笛卡尔采样，运动伪影较小，并且可以充分的欠采样，尤其对于高对比度的图像。螺旋采样充分运用了硬件梯度系统，比较适用于实时以及快速成像系统中。（a）笛卡尔采样（b）径向采样（c）螺旋采样（d）二维随机采样图2.10K空间采样模式2.2.3重建算法利用CS-MRI理论进行图像重建是求解稀疏矩阵的过程，其重建模型重写如下：minTms.t.Fm−y<ε（2-10）s12求解上式是一个解最优化的问题。求解方法主要包括贪婪算法和凸优化算法等，贪婪算法的重建速度很好，缺点的是其重建精度低，精确重建的理论保证不足。并且大多数贪婪算法需要明确知道待重构信号的稀疏度，而在实际应用中这个条件难以达到。因此常用的是凸优化算法。目前较为常见的凸优化算法主要有迭代阈值法、梯度投影法、全变分法、同伦算法、交替方向法、增广拉格朗日法等。下面简单介绍比较常用的收敛速度较快、稳定性较高的非线性共轭梯度下降法[4]。为了将带有约束条件的凸优化问题转变成没有约束的优化问题来进行求解，这里我们引入了拉格朗日乘子。其模型如下：2minFm−y+λTm（2-11）sm21其中，λ表示正则化参数，被用来权衡数据的一致性与稀疏性之间的关系。通过选取不同的λ计算式（2-11），以使Fm−y≈ε。将式（2-11）对m求导得到：s2*∇f(m)=2F(Fm−y)+λ∇Tm（2-12）ss1因为绝对值函数不是光滑的，所以直接求解式（2-12）并不能很好的找到m的解。*因此式（2-12）中的绝对值函数可以通过用光滑函数来代替：x=xx+μ，其中μ是正数，表示光滑因子。在这种近似下有：dxx≈（2-13）dxx*x+μ-14- *令W表示对角矩阵，它的对角元素表达式是w=(Tm)(Tm)+μ，结合（2-12）iii式得到如下表达式：**−1∇f(m)≈2F(Fm−y)+λTWTm（2-14）ss非线性共轭梯度的算法流程如表2.1所示。表2.1非线性共轭梯度算法流程图非线性共轭梯度算法1.输入：y：采样得到的k空间数据；F：欠采样的Fourier变换；sT：稀疏变换算子；λ：权重因子；2.初始化：ε：停止精度；maxIter:最大迭待次数；α、β：线性搜索因子；k=0；m=0；)g=∇f(m；Δm=−g；00003.While(k>maxIterorg<ε)k2t=1;*4.while(f(m+tΔm)>f(m)+αt⋅real(gΔm)){t=βt;}kkkkk更新线性搜索的步长；5.m=m+tΔm；更新图像k+1kk6.g=∇f(m)；根据式（2-15）求解；k+1k+1227.γ=gk+1gk；228.Δm=−g+γΔm；k+1k+1k9.k=k+1;10.endwhile;11.输出m。重建算法是CS-MRI理论的重要内容，它与前两部分讨论的内容紧密相关。许多重建算法是基于某种稀疏表示或采样模式下，确立某种数学模型，根据一定的先验知识进行重建，从而达到提高精度和速度的目的。因此根据实际情况选择合适的稀疏变换方式和采样模式，分析系统特点确立正确的数学模型，相应的选择合适的重建算法是快速精确重建的关键。-15- 2.3本章小结这一章首先讨论了压缩感知的基本原理，给出了利用压缩感知精确重建信号的三个基本条件。此基础上从稀疏表示、采样矩阵和重建算法三方面介绍了将压缩感知理论应用于磁共振图像重建的基本原理，并分别用简单直观的例子进行了说明。为后续内容打下理论基础。综合本章内容，下面给出利用压缩感知对MRI图像进行重建的流程图如图2.11所示：开始医学图像gfourier变换全采样K空间数据D采样模式Fs欠采样K空间数据d一逆fourier变换K空间数据d’fourier变换稀疏变换T含噪声图像M取大系数，重建图像M’逆fourier变换一更新M否满足精度要是输出M图2.11CS-MRI图像重建流程图-16- 第三章低秩矩阵恢复理论3.1低秩矩阵恢复理论低秩矩阵恢复是由JohnWright[39]等人提出的，它是指识别出矩阵中被遮挡或被破坏的元素并对其进行修复以得到原矩阵。若想实现这个目的，则待恢复的矩阵要满足一个条件：矩阵是低秩的，即矩阵的信息分布在线性的低维度子空间中，因此又称为低秩稀疏分解（Sparseandlow-rankmatrixdecomposition）或者鲁棒主成分分析（RobustPrincipalComponentsAnalysis,RPCA）。此外，矩阵被损坏的元子应只是一小部分，即噪声的大小不限但分布应该是稀疏的。对矩阵进行低秩稀疏分解的模型见式（3-1）：minrank(A)+λEs.tA+E=D（3-1）0式中A表示低秩矩阵，E表示低秩矩阵，D表示被破坏的矩阵，λ表示低秩矩阵的权重。因为rank(A)和E都是非凸的非线性优化函数，所以对于式（3-1）的求解0是一个NP难问题。根据数学知识可知，矩阵的秩可以用该矩阵中不为零的奇异值的数目来等价表示，因此对于秩的求解可以用该矩阵的核范数（即矩阵的奇异值之和）来代替，并且根据压缩感知理论可知，在一定的条件下，L0范数是可以用L1范数来近似的。松弛后的数学模型如下：minA+λEs.tA+E=D（3-2）*1其中的意义是指对矩阵中所有的不为零的奇异值进行求和，也就是核范数。*表示的是L1范数，它指对所有的矩阵元素先取绝对值然后再求和。经过近似替代，1就将原问题转化成了一个多变量的凸优化求解问题。对式（3-2）的求解有多种方法，将在下一节进行介绍。3.2低秩矩阵恢复算法在Recht、Candes[3]等人对低秩矩阵恢复的基本理论进行证明的基础上，学者们围绕求解式（3-2）展开了大量研究，目前常用的优化算法有很多种，下面将简要的介绍几种主要的算法。-17- 3.2.1迭代阈值算法对式（3-2）进行正则化就能得到下列优化函数：μ22minA+λE+(A+E)s.tA+E=D（3-3）A,E*12FF式中μ是一个较小的正数。与上式对应的拉格朗日函数为：μ22L(A,E,Y)=A+λE+(A+E)+Y,D−A−E（3-4）*12FF于是μ2μ2argminL(A,E,Y)=argminA+A−Y+λE+E−Y（3-5）A,EA,E*2F12F使用迭代阈值算法[40-42]（IterativeThresholding，IT）交替更新矩阵A、E和Y。更新其中一个时，固定另外两个的值。表达式如下：A=argminL(A,E,Y)k+1kkA2（3-6）11YYkk=argminA+A−=D()1μAμ*2μμFE=argminL(A,E,Y)k+1k+1kE2（3-7）λ1YYkk=argminE+E−=S()λμEμ12μμFY=Y+δ(D−A−E)（3-8）k+1kkk+1k+1δ式中k表示步长，满足0<δ<1。D(W)表示奇异值收缩算子：k1μTD(W)=US（W）V（3-9）1μ1μ其中S（W）表示阈值收缩算子：λμS（W）=max(W−1μ,0)sgn(W)（3-10）1μsgn(W)是符号函数，本文后续章节出现的D(W)、S（W）的意义和表达式与1μλμ此相同。综上得出IT算法的计算步骤见表3.1。-18- 表3.1低秩矩阵恢复的迭代阈值算法实现步骤低秩矩阵恢复迭代阈值法（IT）1.输入：D、λ、μ；2.whilenotconvergeddo3.(U,S,V)=svd(Y)kT4.A=US[S]Vk+11μ5.E=S[]Yk+1λμk6.Y=Y+δ(D−A−E)；k+1kkk+1k+17.endwhile8.输出：A、E。k+1k+1迭代阈值算法迭代过程比较简单，不足之处是收敛的速度慢，并且步长的大小不容易确定，因此它的应用范围比较有限。3.2.2加速近端梯度法经凸松弛处理将公式（3-2）中的等式约束加入到目标函数中，得到下面的拉格朗日函数：12L(A,E,μ)=μ(A+λE)+D−A−E（3-11）*12F针对式（3-11），文献[43]中介绍了用加速近端梯度算法（AcceleratedProximalGradient,APG）求解优化问题。令g(A,E,μ)=μ(A+λE)（3-12）*112f(A,E)=D−A−E（3-13）2F代入式（3-11）则得到的拉格朗日函数变形为L(A,E,μ)=g(A,E,μ)+f(A,E)。上式中g(A,E,μ)是不可微的，f(A,E)是光滑的并且具有李普希兹连续梯度，所以存在0L>使下式成立：f∇f(A,E)−∇f(A,E)≤L(A−A,E−E)（3-14）fFF其中∇f(A,E)表示矩阵函数f(A,E)关于矩阵变量A和E的Frechet梯度。在此L=2。下面通过与D同型的两个矩阵Y、Y对拉格朗日函数L(A,E,μ)进行部分二fAB次逼近。-19- Q(A,E,μ,Y,Y)=g(A，E,μ)+f(Y,Y)+∇f(Y,Y),(A−Y,E−Y)AEAEAEAEL2f+(A−Y,E−Y)AE2F=g(A,E,μ)+f(Y,Y)+Y−(D−Y),(A−Y)（3-15）AEAEA+Y−(D−Y),(E−Y)EAEL2L2ff+A−Y+E−YAE2F2Fkk当E=E,Y=Y,Y=Y,μ=μ时，kAAEEkkkA=argminQ(A,E,μ,Y,Y)k+1kkAEAL2kkfk=argminμA+Y−(D−Y),A+A−YkAEAA*2F（3-16）L2fkkk=argminμA+A−(Y+(D−Y−Y)L)kAAEfA*2Fkkk=D(Y+(D−Y−Y)L)μkLfAAEfkk类似的，当A=A,Y=Y,Y=Y,μ=μ时，k+1AAEEkkkE=argminQ(A,E,μ,Y,Y)k+1k+1kAEE（3-17）kkk=S(Y+(D−Y−Y)L)λμkLfEAEfY和Y的更新表达式如下：AEkY=A+(t−1)(A−A)t（3-18）Akk-1kk-1kkY=E+(t−1)(E−E)t（3-19）Ekk-1kk-1k2其中t=(1+1+4t)2。参数的更新如下：k+1kμ=max(ημ,μ)k+1k（3-20）其中μ是事先定义的正数，η满足0<η<1。k+1相比于IT算法，APG算法降低了迭代次数，其重建速度在一定程度上得到了提高。利用近端梯度算法求解优化问题的步骤见表3.2。-20- 表3.2APG算法求解低秩矩阵恢复的步骤APG算法求解低秩矩阵恢复步骤1.输入：目标矩阵D；初始化：A=A=0,E=E=0,t=t=0，μ>0；0−10−10−12.whilenotconvergeddok3.Y=A+(t−1)(A−A)tAkk-1kk-1kkY=E+(t−1)(E−E)t；Ekk-1kk-1kkkkk4.G=Y+(D−Y−Y)2；AAAE5.更新A：k+1kT（U,S,V）=svd(G)，A=US(S)V；Ak+1μk26.更新E：k+1kkkkkG=Y+(D−Y−Y)2；E=S[G]；EEAEk+1λμk2E27.t=(1+1+4t)2；k+1kμ=max(ημ,μ)8.k+1k；9.k=k+1；10.endwhile;11.输出：A、E。3.2.3增广拉格朗日法在介绍运用增广拉格朗日算法求解式（3-2）前，这里首先介绍一般的增广拉格朗日乘子法（AugmentedLagrangeMultiplier，ALM）。对于一个约束优化问题：minf(X)s.th(X)=0（3-21）nnm式中f:R→R,h:R→R，写成增广拉格朗日函数的形式为：μ2L(X,Y,μ)=f(X)+Y,h(X)+h(X)（3-22）2F式中μ是一个正数。与普通的拉格朗日函数不同的是，增广拉格朗日函数将约束条件作为一个惩罚项加入到了等式中。ALM算法的求解过程是通过不断迭代以达到增广拉格朗日函数的最小化的目的。先更新X，再根据这个新的X和μ更新Y，然kkk后再更新下一个X，如此往复经过若干次迭代后X将收敛到某一值，根据设定的误kk差范围得到原问题最优解。利用ALM算法求解低秩矩阵恢复问题时，式（3-2）的增广拉格朗日函数的表达-21- 式为如下形式：μ2L(A,E,Y,μ)=A+λE+Y,D−A−E+D−A−E（3-23）*12F上式中的Y表示拉格朗日乘子。通过配方合并后两项，得到下式2μYL(A,E,Y,μ)=A+λE+D−A−E+（3-24）*12μF这里对（3-24）式的求解可以通过两种方法实现，它们分别为精确的增广拉格朗日函数法（ExactALM，EALM）和非精确的增广拉格朗日函数法（InexactALM，IALM）。使用EALM算法进行重建的步骤如下，交替的迭代求解矩阵A和E直到满j足停止条件，首先固定E，当E=E时有下式：k+1j+1jA=argminL(A,E,Y,μ)k+1k+1kkA2μYkjk=argminA+D−A−E+)（3-25）k+1A*2μkFYjk=D(D−E+)1μkk+1μkj+1根据新的A更新矩阵E：k+1j+1j+1E=argminL(A,E,Y,μ)k+1k+1kkE2j+1Yk=argminλE+μD−A−E+(3-26)kk+1E1μkFj+1Yk=S(D−A+)λμkk+1μkj+1j+1**设A和E分别收敛于A、E，据此更新Y如下式：k+1k+1k+1k+1**Y=Y+μ(D−A−E)(3-27)k+1kkk+1k+1最后更新参数μ：ρμ,μE*−E*D<εkkk+1kμk+1=F(3-28)μ，其他k其中ρ是大于1的常数，ε是较小的正数。利用EALM算法对矩阵恢复进行求解的步骤流程见表3.3。-22- 表3.3利用EALM算法解低秩矩阵恢复步骤低秩矩阵恢复的EALM算法1.输入：D，λ，初始化Y，E=0，μ，k=0；0002.whilenotconvergeddo03.E=E,j=0；k+104.whilenotconvergeddoj−15.(U,S,V)=svd(D−E+μY)；k+1kkj+1T6.Ak+1=USμ−1[S]V；kj+1j+1−17.Ek+1=Sλμ−1[D−Ak+1+μkYk]；k8.j=j+1；9.endwhile**10.Y=Y+μ(D−A−E)；μ=ρμ；k+1kkk+1k+1k+1k11.k=k+1；12.endwhile13.输出：A、E。EALM算法在更新迭代过程中，每一步都需要求出子问题的精确解，这比较耗时。事实上只需要对A和E各自更新一次便可以得到该子问题的近似解，并且此解足以使算法收敛到最优解，据此EALM算法便可以简化为收敛速度更快的算法，称之为非精确的增广拉个朗日乘子法（IALM）。IALM算法的目标函数仍为式（3-24），只是更新A和E的表达式略有不同，其更新表达式如下：−1Ak+1=argminL(A,Ek+1,YK,μk)=Dμ−1(D−Ek+1+μkYk)(3-29)Ak−1Ek+1=argminL(Ak+1,E,YK,μk)=Sλμ−1(D−Ak+1+μkYk)(3-30)Ek利用IALM算法进行矩阵恢复的求解步骤见表3.4。-23- 表3.4矩阵恢复的IALM算法矩阵恢复的IALM算法1.输入：D，λ，初始化Y，E=0，μ，k=0；0002.whilenotconvergeddo−13.(U,S,V)=svd(D−E+μY)；kkkT4.Ak+1=USμ−1[S]V；k−15.Ek+1=Sλμ−1[D−Ak+1+μkYk]；k6.Y=Y+μ(D−A−E)；k+1kkk+1k+17.μ=ρμ；k+1k8.k=k+1；9.endwhile10.输出：A、E。相较于其他算法，利用非精确增广拉格朗日算法进行矩阵恢复的速度较快并且能达到更高的精度。3.3低秩矩阵恢复的应用首先在图像处理领域低秩矩阵恢复有广泛的应用。例如Candes等人就在监控视频的背景建模中成功的运用了RPCA。一般监控视频都是由固定相机拍摄的图像序列，这些图像序列之间时间相关性很强，所以将每一个图像展开成一个列向量。因此多个图像序列便能组成一个矩阵，并且这个矩阵是满足低秩特征的。矩阵中的低秩部分描述了视频中的变化缓慢的、近乎静止的背景成分，而稀疏部分则对应于变换较为剧烈的、活动的前景部分。利用低秩矩阵恢复理论便可以将背景和前景分离[33]，如图3.1所示。此外，在文献[44]中介绍的鲁棒联合图像对齐的应用，主要是将一组本质上线性相关的图像通过利用低秩矩阵恢复进行对齐，并去掉图像中存在的损毁、遮挡等。在文献[45][46]中提出了用于提取二维图像中低秩纹理的具有变换不变性的低秩纹理模型。另外，低秩矩阵恢复还可以用于旧电影的修复以及去除视频中的雨线等。低秩矩阵恢复模型与算法的应用研究方兴未艾，它在图像处理、计算机图形学、计算机视觉等领域的应用有巨大的前景。-24- （a）监控视频（b）低秩部分（即背景）（c）稀疏部分（即前景）图3.1低秩矩阵恢复在视频背景建模中的应用3.4本章小结这一章主要介绍了低秩矩阵恢复理论的基本内容。第一小节首先介绍了什么是低秩矩阵恢复以及该理论对应的数学模型。对模型进行了分析，经过凸松弛变换将原来直接求解时存在的NP难问题转换成为一个凸优化求解问题。紧接着在第二小节对求解优化问题的几种方法进行了介绍。围绕凸优化求解模型，分别介绍了IT算法、APG算法以及ALM算法（包括EALM和IALM算法），并且简要讨论了各个算法的优点和不足。最后在第三小节介绍了低秩矩阵恢复的应用，阐述了在实际生活中其巨大的研究价值和应用前景。综上，本章主要目的是引入低秩稀疏分解模型的原理，结合其重构算法以及其在实际应用中的特点，为后续章节中其在CS-MRI重建中的应用做铺垫。-25- -26- 第四章基于低秩稀疏分解模型的动态CS-MRI图像重建压缩感知理论的不断发展为MRI图像的快速重建开辟了一个新途径。尤其是近来学者们将压缩感知和低秩矩阵恢复相结合，为进一步提高动态MRI成像速度提供了新方向。在第三章中，介绍了低秩恢复理论在视频建模中的应用，而动态MRI成像与视频图像非常类似，因此考虑将低秩稀疏矩阵分解模型（LowrankandSparse，L+S）用于动态MRI重建。其中，低秩部分L描述背景成分的时间相关性，稀疏部分S描述前景的动态成分。GaoH等人在文献[34]中初步研究了将L+S分解运用到动态MRI重建，并在文献[35]中提出基于加速弥散加权序列的图像重建。在此基础上RicardoOtazo[36]等人将L+S模型运用到加速动态对比增强MRI重建上，并在文献[37]中基于L+S分解模型，利用IST算法对临床医学上的MRI图像进行重建，并且取得了可观的结果。大量研究表明运用L+S分解模型进行重建，将背景成分L分离后，使S部分的稀疏性更强，整体的重建精度较好，但仍有改进的空间，并且在重建速度上有待提高。下面就对基于L+S分解模型的动态CS-MRI图像重建的重建模型和算法进行介绍。4.1动态CS-MRI图像重建低秩稀疏分解模型类似于视频图像中背景去除的处理方法，Otazo等人基于GaoH的工作将L+S分解模型应用于动态MRI重建。因为动态MRI图像可表示为背景和动态成分的叠加，所以对于一组动态MRI图像可以分解成低秩和稀疏两部分：M=L+S，L和S分别对应背景和动态部分。背景部分对应图像帧之间的高度相关信息，且随时间变化缓慢，是低秩的，可通过核范数求解。动态部分包含了每一帧的不同信息，随时间变化较快，因为连续帧之间实质性的差异通常限定在相对较少的像素点上，所以认为是稀疏的，通过范数增强其稀疏性。需要注意的是，为了更加稳定的区分背景成分和动态成分，1将L+S模型应用于动态MRI重建，需要满足以下三个“不相关”原则[37]：①K空间采样机制E与低秩矩阵L互不相关；②K空间采样机制E和稀疏变换基T互不相关；③低秩矩阵L跟稀疏矩阵S互不相关（即低秩部分不稀疏，且稀疏部分不低秩）。前两个原则的意义在于去噪；第三个的意义在于分离L和S的稳定性。动态MRI基本满足这三条原则。首先，动态MRI的K空间采样机制是随机的对-27- 每个时间点进行不同的变密度下采样。在这中采样模式下，低空间频率通常是充分采样的，并且随着远离K空间中心欠采样因子会增加。低秩部分的列向量无法由随机选择的K空间的高空间频率子集近似表示，其行向量也无法由随机选择时间函数的子集来近似表示。并且当使用傅立叶变换时，由于K空间与时域空间的傅立叶关系，使得他们的不相关性最大。因此，动态MRI满足前两条不相关原则。其次，第三个不相关原则只与重建时用到的稀疏变换有关。在实际应用中，医学图像通常不满足低秩分量和稀疏分量的严格不相关：低秩部分可能在变换域稀疏，稀疏分量可能具有低秩的特征。但是低秩分量的秩远远小于稀疏分量的秩，并且低秩部分的奇异值远大于稀疏部分的，因此大部分背景信息仍保留在低秩部分。此外，在实际应用中，运用低秩稀疏分解模型的目的是重建原图像，因此不需要将低秩分量和稀疏分量完全分离。文献[37]将L+S分解模型应用于动态MRI图像重建，模型为：minL+λTSs.t.E(L+S)=d（4-1）*1其中，λ是调优参数，T是代表对S进行稀疏变换的算子；E表示编码算子或采样算子，作用是对图像进行傅立叶变换，再对傅立叶变换系数欠采样，也就是对图像进行空间编码[47][48]；d表示经过欠采样后得到的K空间的数据，L+S=M表示原图像。类似于低秩矩阵恢复的求解过程，求解（4-1）式通常是将非凸优化问题转化成凸优化问题，接着用无约束优化问题来代替原来的约束优化问题进行求解，然后基于不同的数学模型选择不同的重建方法。下面将一一介绍。4.2基于L+S分解模型的CS-MRI重建算法这一节将介绍求解式（4-1）的凸优化算法，主要包括已有的迭代软阈值算法（IterativeSoftThreshold，IST），以及受低秩矩阵恢复重建算法的启发，本文引入近端梯度算法（APG）和非精确增广拉格朗日算法（IALM）进行重建。并通过对心脏灌注动态MRI成像和心电影MRI成像的仿真实验验证三种算法的重建速度和精度。4.2.1IST算法对式（4-1）进行正则化处理得到下式：12minE(L+S)−d+λL+λTS（4-2）LSL,S22*1对于上式可以通过交替方向法[49]，分裂Bregman方法[50][51]或其他凸优化算法进行求解。这里介绍迭代简单的IST算法。IST算法通过交替对低秩矩阵L和稀疏矩阵TS进行软阈值迭代，根据每次迭代的-28- 结果更新图像M，直到满足停止准则：L+S−(L+S)≤εL+S（4-3）k+1k+1kkkk22其中ε是很小的正数，控制重建的精度。更新低秩矩阵是对矩阵L的奇异值设置阈值并进行收缩处理，更新TS是对其元素设置阈值并进行处理。两者的更新表达式如下：Lk+1=DλL(Mk−Sk)(4-4)−1S=T(S(T(M−L)))(4-5)k+1λSkk其中D、S分别是奇异值收缩算子和软阈值收缩算子，其表达式在第三章式λLλS（3-9）、（3-10）中给出。由L、S得到更新图像M的表达式：kkHM=L+S−E(E(L+S)−d)(4-6)k+1k+1k+1k+1k+1式中E(L+S)−d表示K空间去噪。k+1k+1综上，给出动态CS-MRI图像重建的IST算法计算流程见表4.1。表4.1IST算法计算流程IST算法计算流程1、输入：d：多线圈欠采样的K空间数据；E：编码算子或采样算子；T：稀疏变换；λL：奇异值阈值；λS：稀疏变换部分的阈值；H2、初始化：M=Ed，S=0，ε；003、whilenotconvergeddo4、更新L：Lk+1=DλL(Mk−Sk)−15、更新S：S=T(S(T(M−L)))k+1λSkkH6、更新M：M=L+S−E(E(L+S)−d)k+1k+1k+1k+1k+17、endwhile;8、输出：L，S。-29- 4.2.2APG算法IST算法迭代过程比较简单易懂，能达到一定的精度。但是由于迭代的次数较多，导致重建时间变长。根据动态MRI图像与视频建模之间的相似性，以及第三章中讨论的APG和IALM低秩矩阵恢复方法，设想将APG算法和IALM算法应用于模型为L+S的动态MRI的重建。基于这个设想，本文将分别推导利用APG算法和IALM算法进行动态MRI重建的模型表达式和计算步骤，并在4.3节中进行实验验证。下面首先介绍APG算法。将L+S重建模型重写如下：minL+λTSs.t.E（L+S）=d（4-7）*1其中L+S=M，通过低秩稀疏分解最终是要求出原图像M。由于算子E含有H(n)(n)(n)H下采样过程，所以根据式（4-7）的约束条件可知Ed=L+S=M，E是对（n）(n)(n)偶算子，M是含有噪声的图像，L、S是相应的低秩部分和稀疏部分。因为整个过程中实际可掌握的只有K空间下采样的数据d，因此用下列表达式来代替式（4-7）。(n)(n)(n)(n)H(n)minL+λTSs.t.L+S=Ed=M(4-8)*1通过上式求解原图像M是一个图像去噪的过程，将后面的等式约束凸松弛到目标函数中得到下列拉格朗日函数：12(n)(n)(n)(n)(n)L(L,S,μ)=μ(L+λTS)+M−L−S(4-9)*12F为了书写简洁，将式（4-9）写为：12L(L,S,μ)=μ(L+λTS)+M−L−S(4-10)*12F令g(L,S,μ)=μ(L+λTS)(4-11)*112f(L,S)=M−L−S(4-12)2F则L(L,S,μ)=g(L,S,μ)+f(L,S)(4-13)(n)参考3.2.2节中的内容，下面通过与M同型的两个矩阵Y、Y对拉格朗日函数LS-30- L(L,S,μ)进行部分二次逼近。Q(L,S,μ,Y,Y)=g(L,S,μ)+f(Y,Y)+∇f(Y,Y),(L−Y,S−Y)LSLSLSLSL2f+(L−Y,S−Y)LS2F=g(L,S,μ)+f(Y,Y)+Y−(M−Y),(L−Y)LSLSL+Y−(M−Y),(S−Y)SLSL2L2ff+L−Y+S−YLS2F2F(4-14)kk当S=S,Y=Y,Y=Y,μ=μ时，kLLSSkkkL=argminQ(L,S,μ,Y,Y)k+1kkLSLL2kkfk=argminμL+Y−(M−Y),L+L−YkLSLL*2F(4-15)L2fkkk=argminμL+L−(Y+(M−Y−Y)L)kLLSfL*2Fkkk=D(Y+(M−Y−Y)L)μkLfALSfkk类似的，当L=L,Y=Y,Y=Y,μ=μ时，k+1LLSSkkkS=argminQ(L,S,μ,Y,Y)k+1k+1kLSS(4-16)−1kkk=T(S(T(Y+(M−Y−Y)L)))λμkLfSLSfY和Y的更新表达式如下：LSkY=L+(t−1)(L−L)t(4-17)Lkk-1kk-1kkY=S+(t−1)(S−S)t(4-18)Skk-1kk-1k2t=(1+1+4t)2其中k+1k。参数的更新如下：μ=max(ημ,μ)k+1k(4-19)其中μ是事先定义的正数，η满足0<η<1。利用APG算法求解基于L+Sk+1分解模型的dMRI图像重建的计算流程见表4.2。-31- 表4.2利用APG算法进行dMRI图像重建的计算流程APG算法用于dMRI重建的流程1、输入：d：多线圈欠采样的K空间数据E：编码算子或采样算子；T：稀疏变换；λ：稀疏部分的权重系数；H2、初始化：M=Ed、k=0；0L=L=0,S=S=0,t=t=0，μ>0；0−10−10−13.whilenotconvergeddok4.Y=L+(t−1)(L−L)tLkk-1kk-1kkY=S+(t−1)(S−S)t；Skk-1kk-1kkkkk5.G=Y+(M−Y−Y)2；LLLS6.更新L：k+1kT（U,S,V）=svd(G)，L=US(S)V；Lk+1μk27.更新S：k+1kkkk−1kG=Y+(M−Y−Y)2；S=T(S[T(G)])；SSLSk+1λμk2SH8.更新M：M=L+S−E(E(L+S)−d)；k+1k+1k+1k+1k+129.t=(1+1+4t)2；k+1kμ=max(ημ,μ)10.k+1k；11.k=k+1；12.endwhile;13.输出：L、S。4.2.3IALM算法根据式（4-8）运用IALM算法进行变形处理得到下列无约束优化问题：μ2(n)(n)(n)(n)(n)(n)(n)(n)minL+λTS+Y,M−L−S+M−L−S（4-20）L,S*12F其中Y表示拉格朗日乘子，通常是一个未知的标量，用来表示约束方程的梯度的线性组合里每个向量的系数。通过Y的引入，便将原来的带有约束条件的最优化问题变换成没有约束条件的最优化问题。利用配方将后两项合并，进一步变形得到拉格朗日函数如下所示：-32- 2(n)(n)(n)(n)μ(n)(n)(n)YL(L,S,Y,μ)=L+λTS+M−L−S+（4-21）*12μF式中μ是一个很小的正数，Y是拉格朗日乘子。拉格朗日函数式（4-21）的求解(n)(n)(n)(n)利用交替更新的方式进行，首先固定Sk和Yk，更新求解Lk+1使)L(Lk+1,Sk,Yk,μk最(n)(n)(n)(n)小；然后固定Lk+1和Yk，更新求解Sk+1使L(Lk+1,Sk+1,Yk,μk)最小，在此基础上更新（n）(n)(n)噪声图像Mk+1。更新L时运用奇异值阈值方法，更新S时运用软阈值迭代方法，如此迭代数次就可以收敛到目标函数的最优解。其收敛性的证明在后面给出。下面给(n)(n)（n）出更新L、S、Y、μ和M的表达式。2(n)(n)μ(n)(n)(n)YkL=argminL+M−L−S+k+1(n)k+1*kk+1kLk+12μkF（4-22）(n)(n)Yk=Dμ−1(Mk−Sk+)μk2(n)(n)μ(n)(n)(n)YkS=argminλTS+M−L−S+k+1(n)k+11kk+1k+1S2μ（4-23）kF−1−1=T(Sλμ−1(T(Mk−Lk+1+μYk)))(n)(n)(n)Y=Y+μ(M−L−S)（4-24）k+1kkkk+1k+1μ=ρμ（4-25）k+1k(n)(n)(n)H(n)(n)M=L+S−E(E(L+S)−d)（4-26）k+1k+1k+1k+1k+1基于L+S分解模型的动态MRI重建算法计算流程见表4.3。-33- 表4.3L+S模型动态MRI重建的IALM算法L+S重建的IALM算法流程1、输入：d：多线圈欠采样的K空间数据E：编码算子或采样算子；T：稀疏变换；λ：稀疏部分的权重系数；Hμ2、初始化：M=Ed、Y、S=0、、ε、iter=0，k=0；03、whilenotconvergeddoYk4、更新L：Lk+1=Dμ−1(Mk−Sk+)；μk5、更新S：−1(((−1)))；Sk+1=TSλμ−1TMk−Lk+1+μYk6、更新Y：Yk+1=Yk+μk(Mk−Lk+1−Sk+1)；7、更新μ：μ=ρμ；k+1k+1kH8、更新M：M=L+S−E(E(L+S)−d)；k+1k+1k+1k+1k+19、判断停止条件：L+S-（L+S）≤εL+S；k+1k+1kkkk2210、更新迭代次数：iter=iter+1，k=k+1；11、endwhile12、输出：L，S。下面给出利用IALM算法进行L+S分解的动态MRI重建的收敛性证明。对于基于L+S模型重建的IALM算法，如果每一次更新{}μ满足是非减的，并且k+∞−1(n)(n)μk=+∞，则Lk、Sk可以收敛到问题的最优解。其收敛性证明如下：k=1文献[52]中的引理4：若μ是非减的，则式（4-27）中各项是非负的，并且其和k是有限的。+∞∧μ−1Y−Y,S−S+L−L*,Y−Y*+S−S*,Y−Y*<+∞（4-27）kk+1Kk+1Kk+1k+1k+1k+1k=1∧***其中L、S、Y表示目标函数的最优解，Yk+1∈∂Lk+1*，Yk+1∈∂(λSk+1)，∂1表示求梯度。参考[52]引理4的证明有下式：+∞2−2μkYk+1−Yk<+∞（4-28）Fk=1根据Y=Y+μ(M−L−S)，得下式：kk−1k−1kk-34- −1M−L−S=μY−Y→0（4-29）kkkkk−1FF(n)(n)由此可知(L,S)的任何一个聚点都可能是式（4-21）的最优解。另一方面，假kk∧*设优化问题的最优解表示为f，而Yk∈∂Lk，Yk∈∂(λSk)，因此有下式：*1L+λTS≤L+λSkkkk*1*1∧****≤L+λTS−Yk,L−Lk−Yk,S−Sk*1∧********=f+Y-Yk,L−Lk+Y−Yk,S−Sk−Y,L−Lk+S−Sk∧******=f+Y-Yk,L−Lk+Y−Yk,S−Sk−Y,M−Lk−Sk（4-30）由引理4可知+∞∧−1****μkLk−L,Yk−Y+Sk−S,Yk−Y<+∞（4-31）k=1+∞−1因为μk=+∞，所以一定存在(Lkj,Skj)满足：k=1∧****L−L,Yk−Y+S−S,Y−Y→0（4-32）kjjkjkj带入式（30）得到：*limL+λTS≤f（4-33）j→+∞kj*kj1**由式（33）可知(L,S)比较接近于最优解(L，S)。因为μ→+∞，并且{Y}kjkjkk22*−2*是有界的（证明见[52]中的引理1），所以有Skj−SF+μkjYkj−YF→0。22*−2*又因为S−S+μY−Y是非增的（原因可参考文献[27]中引理4的证明），kkkFF22*−2*所以可知S−S+μY−Y→0，由此得到kkkFF*limS=S（4-34）kk→+∞**由于limM−L-S=0，并且M=L+S，所以得到kkk→+∞*limL=L（4-35）kk→+∞-35- 至此证毕得出结论，即利用IALM算法求解L+S模型可以收敛到目标函数的最优解。4.3实验结果及分析对于前一部分介绍的L+S模型的重建算法，这一节将通过仿真实验对其进行验证。分别利用IST算法、APG算法和IALM算法分别对心脏灌注MRI成像和心电影MRI成像进行重建。然后比较根据不同算法重建效果的均方误差（MSE），差错率（err），峰值信噪比（PSNR）和信噪比（SNR）。若用g表示原始图像，f表示重构图像，则各个准则的数学表达式分别为：2均方误差：MSE=()1M*N(g−f)（4-36）i,j2峰值信噪比：PSNR=10lg(255MSE)（4-37）差错率：err=(g−f)(g)（4-38）FF22信噪比：SNR=10lgg(g−f)（4-39）i,ji,j仿真是在物理内存为4GB的Window7系统的电脑上进行的，仿真软件是Matlab2013。实验采用的采样模式为笛卡尔采样，稀疏变换是三维傅立叶变换。实验数据来自于通过用12元线圈新型梯度回波序列对健康志愿者进行3T扫描。心脏灌注图像数据取自于中心室短轴位舒张中期的采样，采样得到40帧大小为128x128的图像矩阵。采样加速因子为8。相关的图像参数包括：FOV=320x320mm2，切片厚度为8mm，翻转角10度，TE/TR=1.2/2.4ms，空间分辨率3.2x3.2mm2，时间分辨率307ms。心电影MRI图像数据通过24帧大小为256x256的矩阵获取，FOV=320x320mm2，加速因子为8。心脏灌注MRI结合心电影MRI是临床用来辅助诊断心肌壁肥厚和心肌壁梗死等疾病的重要工具。心脏灌注成像如图4.1所示，抽取第5、8、11、14帧数据重建结果。最左边是利用IST算法进行的重建，中间是利用APG算法，最右边是利用IALM算法。图4.1中白色高亮的部分反映的是心脏中的血液，高亮部分周围的黑色圆环反应的是心肌壁，根据图像作为辅助工具来判断心肌壁的健康情况。-36- L+SL+SL+SLLLSSS（a）IST算法(b)APG算法（c）IALM算法图4.1IST、APG、IALM对心脏灌注MRI的重建从图4.1中看出，三种算法的重建效果在L部分和S部分存在一些差异，这跟初始参数的设置有关，并且我们的目标是通过L+S分解模型来重建原图像M，所以这里L和S部分的差异可以忽略，着重考虑L+S重建的部分。图4.1中从最上面一行可以看出三种方法都能清晰的看出心肌壁的形状，重建效果基本相同，从主观上看不相上下。接下来将分别从重构时间t、均方误差MSE、峰值信噪比PSNR、差错率err以及信噪比SNR等四方面针对不同算法进行比较。实验抽取时间帧为2、5、8、11、14、17、23、26的重建结果进行比较，实验数据见表4.4，为了更直观的反映不同算法的重建效果，根据实验数据画出折线图，如图4.2所示。表4.4心脏灌注MRI不同时间帧重建结果Time2581114172326IST8.91279.99399.673010.6827.32838.82358.55147.6804MSEAPG4.92506.35845.52106.03693.29564.79604.72514.3957/10-7IALM4.86886.19895.38265.83683.23554.71314.55884.2610108.630108.133108.275107.844109.480108.674108.810109.277IST74248440PSNR111.206110.097110.710110.322112.951111.322111.386111.700APG（dB）84764176111.256110.207110.820110.469113.031111.397111.542111.835IALM67904737IST0.10080.08520.08480.07970.05850.04640.04850.0517errAPG0.07350.06690.06320.05920.03880.03370.03550.0386IALM0.07300.0660.06230.05820.03840.03340.03490.0380SNRIST19.934321.392521.432621.97224.650326.674326.292325.7282(dB)APG22.678423.487923.986124.553428.222729.442228.984628.2700IALM22.73923.604524.103824.707428.305729.52529.145928.4116ReconstructionIST：t=69.855967APG：t=59.410806IALM:t=54.544524time（S）-37- 113.511.010.5IST113.0IST10.0APGAPG112.59.5IALMIALM9.0112.08.5111.58.07.5111.0)7.0-7110.56.56.0110.0MSE(x105.5PSNR(dB)109.55.04.5109.04.0108.53.53.0108.02.52.0107.502468101214161820222426280246810121416182022242628timetime（a）IST、APG、IALM重建图像的MSE（b）IST、APG、IALM重建图像的PSNR320.110IST0.10531ISTAPG0.10030APGIALM0.095IALM290.090280.0850.080270.075260.07025err0.065SNR(dB)240.0600.055230.050220.045210.040200.035190.03002468101214161820222426280246810121416182022242628timetime（c）IST、APG、IALM重建图像的err（d）IST、APG、IALM重建图像的SNR图4.2IST、APG、IALM对心脏灌注MRI重建结果折线图根据表4.4和图4.2可知，相对于IST算法重建结果的均方误差，APG算法和IALM算法重建结果的均方误差更低，近似降低了50%，并且IALM算法性能上更加优于APG算法。差错率上本文引入的两种算法也是远远的小于IST算法：APG算法比IST算法平均降低了26.6%；IALM算法比IST平均降低了27.5%。相应的，APG算法与IALM算法的PSNR以及SNR都比IST算法的大很多：APG算法的峰值信噪比平均增长了2.3%，信噪比平均增长了11.5%；IALM算法的峰值信噪比平均增长了2.5%，信噪比平均增长了12%。同时，在重建速度上，APG算法比IST算法提高了近15%，IALM算法比IST算法提高了近22%。总的来说，APG算法和IALM算法在重建精度和速度上较IST有很大的提升，并且IALM算比APG算法更好。图4.3所示为心电影MRI成像，依次抽取24帧数据中的第2、8、14、20帧数据重建的图像。从左到右对应的重建算法分别是IST算法、APG算法、IALM算法。-38- L+SL+SL+SLLLSSS（a）IST算法（b）APG算法（c）IALM算法图4.3IST、APG、IALM算法对心电影MRI的重建从图4.3中可以看出，根据三种方法重建的图像除了低秩和稀疏两部分的图像的分配比例不同外，在L+S部分重建效果是很近似的，根据APG算法和IALM算法重建的图像比IST算法更清晰一些。为了进一步说明本文引入的两种算法的优势，下面将从客观标准上进行比较。分别抽取时间帧为2、5、8、11、14、17、20、23的重建图像，重构结果数据如表4.5所示，为了更直观的比较，这里根据实验所得的数据画出折线图，如图4.4所示。表4.5心电影MRI不同时间帧的重建结果Time2581114172023IST1.94702.12091.95952.14921.84511.96811.98801.7545MSEAPG1.91642.07371.96462.13501.83361.92551.96461.7409/10-7IALM1.89142.05251.93592.11561.81571.90551.95491.7106115.237114.865115.209114.808115.470115.190115.146115.689IST26406473PSN115.306114.963115.198114.836115.497115.285115.198115.723RAPG03097521(dB)115.363115.007115.262114.876115.540115.330115.219115.799IALM09044752IST0.02740.02990.02780.02830.02880.03020.02820.0299errAPG0.02710.02950.02780.02820.02860.02980.02800.0297IALM0.02690.02930.02760.02800.02850.02970.02790.0294IST31.248530.495531.115230.955530.816930.390230.987030.4815SNRAPG31.339030.616231.126131.004630.863830.506331.055730.5425(dB)IALM31.398030.663131.192031.045930.908730.553331.078730.6212Reconstructiontime(S)IST：t=235.602168APG:t=152.936865IALM:t=139.994591-39- 2.20115.92.15IST115.8ISTAPGAPG2.10IALM115.7IALM115.62.05115.52.00)115.4-71.95115.31.90115.2MSE(x101.85PSNR(dB)115.11.80115.01.75114.91.70114.81.65114.70246810121416182022242602468101214161820222426timetime（a）IST、APG、IALM重建图像的MSE（b）IST、APG、IALM重建图像的PSNR0.031031.5ISTIST31.40.0305APGAPGIALM31.3IALM0.030031.20.029531.131.00.029030.9err0.028530.8SNR(dB)0.028030.730.60.027530.50.027030.40.026530.30246810121416182022242602468101214161820222426timetime（b）IST、APG、IALM重建图像的err（d）IST、APG、IALM重建图像的SNR图4.4心电影MRI不同方法重建结果折线图根据表4.5和图4.4可知，本文引入的APG算法和IALM算法在在重建精度上有很大的提高。与IST算法相比较，运用APG算法得到的峰值信噪比平均增长了0.04%，信噪比平均增长了2.3%，均方误差平均降低了一个百分点，差错率降低较小平均为近一个百分点，重建速度提高了35%；而运用IALM算法得到的峰值信噪比平均提高了0.08%，信噪比平均提高了4%，均方误差平均降低了两个百分点，差错率有一个百分点的下降，在重建速度上加快了40%。此外，根据四副折线图可以看出，在t=8时，APG算法的均方误差稍微大于IST算法，差错率和IST算法基本相同，峰值信噪比稍微小于IST算法，信噪比大于IST算法。而IALM算法在整体效果上都是好于IST算法和APG算法的。综合心脏灌注MRI和心电影MRI图像的重建可知，针对L+S分解模型的dMRI重建，本文引入的APG算法和IALM算法在重建精度和速度上都比原来的IST算法好，同时IALM算法的性能要比APG算法的更好，且稳定性更强。-40- 4.4本章小结本章首先给出了利用低秩稀疏分解模型进行动态MRI图像重建的基本理论，介绍了求解该模型已有的IST算法，根据模型的特点，将低秩矩阵恢复算法中的APG算法和IALM算法推广应用到L+S分解模型的MRI重建中，并推导了求解的公式，以及求解流程。第三部分基于IST算法、APG算法和IALM算法进行matlab软件仿真，分别用心脏灌注MRI和心电影MRI进行实验。从主观观察和客观指标两方面进行比较，经实验验证，本文引入的APG算法和IALM算法无论在重建精度还是重建速度上，都优于IST算法，并且相较于APG算法，IALM算法的重建效果更好，且性能更稳定。-41- -42- 第五章改进的L+S分解模型及重建算法虽然基于L+S分解模型的MRI图像重建虽然具有一定的重建精度，但仍可以进一步的提高。特别是在对动态MRI图像进行重建的过程中，由于K空间数据的下采样造成重建图像中存在伪影，为了去除伪影更好的重建图像，提高重建精度，本文基于全变分去噪声去模糊的特点，提出在原有的L+S模型的基础上引入全变分（TotalVariation）正则项。下面首先介绍一般全变分去噪声去伪影的基本原理，并给出求解方法，在此基础上将TV正则项引入到L+S模型中，利用第四章的IALM方法对新重模型的重建步骤进行推导。最后利用Matlab仿真实验验证在模型引入TV正则项后的重建效果。5.1全变分去噪作为图像去噪与图像复原领域最为成功的方法之一，全变分图像去噪模型最早是由Rudin[53]等人提出的。它的优点是充分利用了图像本身的正则性，便于将真实图像的几何正则性从被污染的噪声图像的解中反映出来[54]。用全变分方法处理图像可以在平滑噪声的同时，保留图像的边缘以及纹理信息，达到强化边缘的目的。TV去噪模型[55]如下：2x=argminx−y+αTV(x)（5-1）x上式中x用来表示原始图像，y用于表示含有噪声的图像，α则表示TV项的权重。其中TV的表达式如下：h2v2TV（x）=(Δix)+(Δjx)（5-2）i,jhvhv其中Δix=xi−xl(i)，Δix=xi−xa(i)。一阶微分算子Δi、Δi表示分别对行和列求一阶微分，下标l(i)指第i行的左边的最近的元素，下标a(i)指第i列的上边的最近的元素。TV最小化的方法有多种，其中优化最小化算法（Majorization-Minimization，MM）的收敛速度较快，所以在此利用MM算法进行求解。观察可知虽然TV(x)是x的凸函数，但并不是处处可微分的，因此直接求解式（5-1）通常难度较大，所以寻找一个简单的优化问题替代原问题。将式（5-1）表示为函数F(x)，设计一个函数G(x)，使其满足G(x,x)≥F(x),∀x,x，并且G(x,x)=F(x)。Figueiredo[55]等人成功构造出了kkkkk函数G(x)，其表达式如下：-43- TT(k)TG(x,x)=x(DΛD+I)x−2xy（5-3）k(k)(k)(k)其中Λ=diag(w,w)（5-4）−1(k)(h2v2)w=α2（Δx）+（Δx）（5-5）ikikT[]hTvTD=(D)(D)（5-6）hv式（5-6）中D、D分别表示对行向量和列向量的一阶微分。式（5-3）是一个二次函数，因此求解x相当于求解下式：T(k)x=solution{(DΛD+I)x=y}（5-7）k+1x考虑到矩阵的维数太大，直接求解式（5-7）非常困难。并且，当h2v2(k)T(k)（Δx）+（Δx）=0时，相应的w和矩阵(DΛD+I)中的一些元素将趋近于无ikik(k)穷大。而当Λ不断趋于无穷大时，找不到求解矩阵的数据稳定性方法。下面通过矩阵求逆引理来处理此问题，令：−1T(k)()TT(k)−1−1DΛD+I=I−D(DD+(Λ))D（5-8）T(k)则x=y−Dz（5-9）k+1(k)T(k)−1其中z=solution{(DD+(Λ))z=Dy}（5-10）经过变换处理后，式（5-9）可以利用数值稳定性方法求解。求解式（5-10）这里采用阈值收缩算法，表达式如下：(k)T(k)−1z=Sα(Dy(DD+(Λ)))（5-11）综上，利用MM算法进行TV去噪的流程如表5.1所示，算法的收敛性证明参见文献[56]。表5.1MM算法求解TV去噪TV去噪MM算法1、初始化：x=y，iter=0，停止准则ε；2、whilenotconvergeddo(k)(k)(k)(k)3、求解：Λ=diag(w,w)，（w由式（5-5）给出）；T(k)(k)4、更新：x=y−Dz，（z由式（5-10）给出）；k+15、判断停止准则：x−x≤εx；k+1kk226、更新迭代次数：iter=iter+1；7、endwhile；8、输出：x。-44- 5.2联合TV正则化MRI图像重建由欠采样重建的图像通常包含噪声和伪影，文献[4]中指出可以通过加入TV正则项进行去噪声和去伪影。受此启发，本文在L+S分解模型的基础上，联合TV正则项，并利用第四章介绍的较优的IALM算法进行重建，最后通过心脏灌注MRI和心电影MRI重建验证新模型的性能。5.2.1引入TV正则项的L+S重建模型及重建算法将动态MRI重建的L+S分解模型重写如下：(n)(n)(n)(n)H(n)minL+λTSs.t.L+S=Ed=M（5-12）*1在原来的L+S分解模型中引入TV正则项，得到新的重建模型表达式如下：(n)(n)(n)(n)(n)HminL+λTS+αTV(M)s.t.L+S=Ed（5-13）*1其中λ、α分别表示稀疏项和TV正则项的权重，其值的选择是通过比较重建效果，在一个范围中选定的值。对上式运用IALM算法处理，得到没有约束条件的求最优问题：μ2(n)(n)(n)H(n)(n)H(n)(n)minL+λTS+αTV(M)+Y,Ed−L−S+Ed−L−SL,S*12F（5-14）利用配方将后两项合并，进一步变形得到拉格朗日函数如（5-15）所示：2(n)(n)(n)(n)(n)μH(n)(n)YL(L,S,Y,μ)=L+λTS+αTV(M)+Ed−L−S+*12μF（5-15）其中μ是一个较小的正数，Y是拉格朗日乘子。λ、α分别代表稀疏部分L和TV正则项的权重的大小。对式（5-15）求解，通过4.2.3节中公式（4-22）-（4-26）更新(n)(n)L、S、Y、μ以及M。TV项的求解利用MM算法，算法涉及的表达式以及参数-45- 等都已经在5.1节中做了详细说明，其更新表达式如下：(n)(n)T(n)T(k)−1M=M−D(Sα(DM(DD+(Λ))))（5-16）综合以上内容，下面给出利用IALM算法对L+S+TV模型的重建流程，如表5.2所示。表5.2L+S+TV重建的IALM算法流程L+S+TV重建的IALM算法流程1、输入：d：多线圈的欠采样K空间数据E：编码算子或采样算子；T：稀疏变换；λ：稀疏部分的权重系数；α：TV项的权重系数；Hμ2、初始化：M=Ed、Y、S=0、、ε、iter=0；3、whilenotconvergeddoTT(k)−14、更新TV项：M=M−D(S(DM(DD+(Λ))))kα5、更新L：(-1)；L=SVTM−S+μYk+1-1kkkμ−1−16、更新S：Sk+1=T(Sλμ−1(T(Mk−Lk+1+μYk)))；7、更新Y：Y=Y+μ(M−L−S)；kkkk+1k+1μμ=ρμ8、更新k+1：k+1k；H9、更新M：M=L+S−E(E(L+S)−d)；k+1k+1k+1k+110、判断停止条件：L+S-（L+S）≤εL+S；k+1k+1kkkk2211、更新迭代次数：iter=iter+1；12、endwhile4、输出：L，S。5.3实验结果及分析下面通过Matlab仿真实验验证利用新模型进行MRI重建的效果。实验数据仍是来-46- 自对健康志愿者的3T扫描得到的采样得到40帧大小为128x128的心脏灌注图像矩阵以及24帧大小为256x256的心电影图像矩阵。实验中采用的采样模式是笛卡尔采样，稀疏变换是三维傅立叶变换。实验通过用IALM算法对L+S模型重建和对引入TV项的新模型重建进行对比，记IALM算法对原模型的重建为IALM；记IALM对新模型的重建为TV+IALM。分别抽取心脏灌注MRI扫描数据的第5、8、11、14帧数据进行重建，结果如图5.1所示。L+SL+SLLSS（a）IALM（b）TV+IALM图5.1IALM、TV+IALM对心脏灌注MRI的重建图5.1中高亮部分表示血液，周围的黑色圆环表示心肌壁，从图中可以看出，利用新模型重建的图像，高亮部分的血液和心肌壁的对比度更明显，比用原模型的重建结果更清晰明显。以上从主观判断上说明了新模型的优点，下面从客观标准来对比新模型的改进方面。分别抽取第2、5、8、11、14、17、23、26帧数据的重建结果，计算它们的MSE、PSNR、err和SNR，并记录在表5.3中，然后根据表5.3中的数据画出折线图如图5.2所示，以便更直观的反映出两种模型下的重建效果。表5.3心脏灌注MRI不同时间帧重建结果时间帧2581114172326MSEIALM4.86886.19895.38265.83683.23554.71314.55884.2610/10-7TV+IALM4.70126.07735.37165.64653.09294.63164.45974.1629111.256110.207110.820110.469113.031111.397111.542111.835IALMPSNR67904737(dB)111.408110.293110.829110.613113.227111.473111.637111.936TV+IALM77802578IALM0.07300.0660.06230.05820.03840.03340.03490.0380errTV+IALM0.07180.06540.06240.05720.03760.03310.03450.0375IALM22.73923.604524.103824.707428.305729.52529.145928.4116SNR(dB)TV+IALM22.879123.690624.100924.849228.50129.600329.243928.5097重构时间（S）IALM:t=54.54452TV+IALM:t=65.498277-47- 113.56.46.2IALMIALM6.0TV+IALM113.0TV+IALM5.85.6112.55.45.25.0112.0)4.8-74.6111.54.44.2MSE(x104.0PSNR(dB)111.03.83.6110.53.43.23.0110.02.82.6109.5024681012141618202224262830323436024681012141618202224262830323436timetime（a）IALM、TV+IALM重建图像的MSE（b）IALM、TV+IALM重建图像的PSNR0.08031.030.50.075IALM30.0IALMTV+IALM29.5TV+IALM0.07029.028.50.06528.027.50.06027.026.50.05526.0err25.50.050SNR(dB)25.024.50.04524.023.50.04023.022.50.03522.021.50.03021.0024681012141618202224262830323436024681012141618202224262830323436timetime（c）IALM、TV+IALM重建图像的err（d）IALM、TV+IALM重建图像的SNR图5.2IALM、TV+IALM对心脏灌注MRI重建结果折线图结合表5.3和图5.2可以看出，新模型对心脏灌注MRI重建效果比原模型有所改进。但是从上面的折线图中可以明显看出，两条折线是非常接近的，点与点之间差距很小，这是由于心脏灌注MRI图像边缘信息明显，纹理细节较少，因此引入TV正则项后改进效果不太显著。此外，由于在原模型的基础上加入了TV项的更新，所以重建速度略有降低，但仍比IST算法高很多。下面将新模型应用于心电影MRI重建，比对重建效果。将新模型应用与心电影MRI成像，依次抽取24帧数据中的2、8、14、20帧数据重建的图像如图5.3所示。-48- L+SL+SLLSS（a）IALM（b）TV+IALM图5.3IALM和TV+IALM算分对心电影MRI的重建从图5.3中可以看出，新模型重建的心电影MRI在低秩部分和原图部分都更加清晰，直观上效果更好，下面从实验数据进一步分析两种模型的重建效果。分别抽取心电影MRI图像的第2,、5、8、11、14、17、20、23帧重建数据进行比对，仍从MSE、PSNR、err和SNR四个标准进行衡量，实验数据如表5.4所示，根据实验得到的数据画出折线图形式如图5.4所示。表5.4心电影MRI不同时间帧的重建结果时间帧2581114172023MSEIALM1.89142.05251.93592.11561.81571.90551.95491.7106/10-7TV+IALM1.87142.04701.90742.10811.79091.87911.93481.6824PSNR115.363115.007115.262114.876115.540115.330115.219115.799IALM（dB）09044752115.409115.019115.326114.891115.600115.391115.264115.871TV+IALM16581456errIALM0.02690.02930.02760.02800.02850.02970.02790.0294TV+IALM0.02680.02930.02740.02800.02830.02950.02780.0292SNRIALM31.398030.663131.192031.045930.908730.553331.078730.6212(dB)TV+IALM31.444730.675231.257431.061630.969330.614531.124230.6946重构时间(S)IALM:t=139.994591TV+IALM:t=187.163462-49- 2.20116.02.15IALM115.9IALMTV+IALMTV+IALM2.10115.82.05115.72.00115.6)1.95115.5-71.90115.41.85115.3MSE(x10PSNR(dB)1.80115.21.75115.11.70115.01.65114.91.60114.8024681012141618202224024681012141618202224timetime（a）IALM、TV+IALM重建图像的MSE（b）IALM、TV+IALM重建图像的PSNR0.029831.5IALM0.0296TV+IALM31.4IALM0.0294TV+IALM0.029231.30.029031.20.02880.028631.10.028431.00.0282err0.028030.90.0278SNR(sB)30.80.027630.70.02740.027230.60.027030.50.02680.026630.4024681012141618202224024681012141618202224timetime（c）IALM、TV+IALM重建图像的err（d）IALM、TV+IALM重建图像的SNR图5.4心电影MRI不同方法重建结果折线图结合表5.4和图5.4可知，新模型对心电影MRI的重建效果要比心脏灌注MRI中好。新模型在均方误差上比原模型平均降低了一个百分点，差错率上降低了0.5%。在峰值信噪比上，新模型比原模型提高了0.04%，在信噪比上提高了0.15%，在虽然重建速度有所降低，但从重建精度和速度上都是优于原始模型和原始的IST算法的。5.4本章小结为了进一步的去除图像伪影和噪声，本章在L+S分解模型的基础上引入了TV正则项。首先介绍了TV去噪的模型以及求解方法，进而给出了联合TV正则化的L+S重建模型，推导出求解方法，利用第四章中引入的IALM算法进行求解。为了验证新模型的性能，将新模型的重建效果与第四章中的IALM算法进行对比，实验数据同样是心脏灌注MRI和心电影MRI。实验结果表明，相对于L+S模型的IALM算法重建，新模型的重建精度得到了进一步的提升，虽然在重建速度上有所减缓，但仍远快于IST算法，在重建精度和速度上达到了较好的平衡。-50- 第六章总结与展望近年来磁共振成像作为一种重要的医学辅助工具得到了广泛的研究和应用。研究的热点主要是如何有效的加快K空间扫描速度，并保证成像质量。将压缩感知理论应用于磁共振称成像可以达到这个目的，即可以从极少量的数据中根据一定的重构算法重建原信号。此外，对于动态磁共振成像，其与视频图像有很大的相似性，都是背景成分变化极小，前景成分变化较为剧烈，这就为动态磁共振图像重建提供了一个新方法，即在压缩感知理论的基础上，再与低秩矩阵恢复理论结合进行动态MRI重建。类似于视频图像，可以用主成分分析的方法通过低秩稀疏矩阵分解模型对动态MRI进行重建。本文首先从三方面介绍了CS-MRI重建的基本原理的，接着介绍了低秩矩阵恢复理论的模型及其求解方法，在此基础上，将L+S分解模型应用于dMRI图像重建。针对该模型，先介绍了已有的IST算法，针对其存在的不足提出将APG算法和IALM算法应用到dMRI重建，并通过实验验证了IST算法、APG算法、IALM算法的性能。另一方面，为了得到更好的重建精度，本文提出在L+S分解模型的基础上，引入TV正则项，并利用IALM算法进行重建。最后通过实验与原模型的IALM算法重建效果进行对比，验证新模型的性能。通过matlab仿真实验表明，针对于L+S分解模型，本文引入的APG算法和IALM方法在重建速度和精度上都远远高于IST算法，尤其是IALM算法，它的速度和稳定性更胜于APG算法。另一方面，在L+S分解模型的基础上引入TV正则项，利用IALM算法进行重建，相较于原模型在速度上有所降低，在重建精度上得到了进一步的提升。综合考虑，本文提出的改进模型以及引入的APG算法和IALM算法在整体的重建精度和速度上都远远优于原模型的IST算法。虽然本文提出的模型和算法在重建精度以及速度上都有所改进，但是由于在重建过程中，迭代更新M项时涉及到多次傅立叶变换，需要耗费大量时间，因此减少这一部分的时间消耗是进一步加快重建速度的关键，也是下一步要研究的内容。-51- -52- 参考文献[1]UeckerM,LaiP,MurphyMJ,etal.ESPIRiT—aneigenvalueapproachtoautocalibratingparallelMRI:whereSENSEmeetsGRAPPA[J].MagneticResonanceinMedicine,2014,71(3):990-1001.[2]DonohoDL.Compressedsensing[J].InformationTheory,IEEETransactionson,2006,52(4):1289-1306.[3]CandesEJ,RombergJ,TaoT.Robustuncertaintyprinciples:Exactsignalreconstructionfromhighlyincompletefrequencyinformation[J].InformationTheory,IEEETransactionson,2006,52(2):489-509.[4]LustigM,DonohoD,PaulyJM.SparseMRI:TheapplicationofcompressedsensingforrapidMRimaging[J].Magneticresonanceinmedicine,2007,58(6):1182-1195.[5]LustigM,DonohoDL,SantosJM,etal.CompressedsensingMRI[J].SignalProcessingMagazine,IEEE,2008,25(2):72-82.[6]CaiJF,DongB,OsherS,etal.Imagerestoration:Totalvariation,waveletframes,andbeyond[J].JournaloftheAmericanMathematicalSociety,2012,25(4):1033-1089.[7]DemirelH,AnbarjafariG.Discretewavelettransform-basedsatelliteimageresolutionenhancement[J].GeoscienceandRemoteSensing,IEEETransactionson,2011,49(6):1997-2004.[8]ZhangY,DragottiPL.Onthereconstructionofwavelet-sparsesignalsfrompartialFourierinformation[J].SignalProcessingLetters,IEEE,2015,22(9):1234-1238.[9]RavishankarS,BreslerY.MRImageReconstructionfromHighlyUndersampledk-spaceDatabyDictionaryLearning[J].IEEETransonMedicalImaging,2011,30:1028-1041.[10]QuXB,ZhangWR,GuoD,CaiCB,etal.Iterativethres-holdingcompressedsensingMRIbasedoncontourlettransform[J].InverseProblemsinScienceandEngineering,2010,18(6):737-758.[11]MinhN.Do,MartinVetterli.Thecontourlettransform:Anefficientdirectionalmultiresolutionimagerepresentation.IEEEtransactionsonimageprocessing.2005,14(12):2091-2106.[12]QuX,GuoD,NingB,et.al.UndersampledMRIreconstructionwithpatch-baseddirectionalwaveltets[J].Magneticresonanceimaging,2012,30:964-977.[13]QuXB,HouYK,LamF,etal.Megneticresonanceimagereconstructionfromundersampledmeasurementsusingapatch-basednonlocaloperator[J],MedicalImageAnalysis,2014,18:843-856.[14]SongT,StainsbyJA,HoVB,etal.FlexiblecardiacT1mappingusingamodifiedlook–lockeracquisitionwithsaturationrecovery[J].MagneticResonanceinMedicine,2012,67(3):622-627.[15]FengL,GrimmR,BlockKT,etal.Golden‐angleradialsparseparallelMRI:Combinationofcompressedsensing,parallelimaging,andgolden‐angleradialsamplingforfastandflexibledynamicvolumetricMRI[J].Magneticresonanceinmedicine,2014,72(3):707-717.[16]ShinT,NayakKS,SantosJM,etal.Three‐dimensionalfirst‐passmyocardialperfusionMRIusingastack‐of‐spiralsacquisition[J].MagneticResonanceinMedicine,2013,69(3):839-844.[17]FieldenSW,MuglerJP,HagspielKD,etal.NoncontrastperipheralMRAwithspiralechotrainimaging[J].MagneticResonanceinMedicine,2015,73(3):1026-1033.[18]LiuE,TemlyakovVN.Theorthogonalsupergreedyalgorithmandapplicationsincompressedsensing[J].InformationTheory,IEEETransactionson,2012,58(4):2040-2047.[19]CandèsEJ,RechtB.Exactmatrixcompletionviaconvexoptimization[J].FoundationsofComputationalmathematics,2009,9(6):717-772.-53- [20]YuanJ,QiuW,UkwattaE,etal.Anefficientconvexoptimizationapproachto3DprostateMRIsegmentationwithgenericstarshapeprior[J].ProstateMRImageSegmentationChallenge,MICCAI,2012,7512.[21]RajchlM,YuanJ,WhiteJ,etal.Afastconvexoptimizationapproachtosegmenting3Dscartissuefromdelayed-enhancementcardiacMRimages[J].MedicalImageComputingandComputer-AssistedIntervention–MICCAI2012,2012:659-666.[22]MalekiA,DonohoDL.Optimallytunediterativereconstructionalgorithmsforcompressedsensing[J].SelectedTopicsinSignalProcessing,IEEEJournalof,2010,4(2):330-341.[23]VoroninS,ChartrandR.Anewgeneralizedthresholdingalgorithmforinverseproblemswithsparsityconstraints[C].Acoustics,SpeechandSignalProcessing(ICASSP),2013IEEEInternationalConferenceon.IEEE,2013:1636-1640.[24]ChenC,LiY,HuangJ.CalibrationlessparallelMRIwithjointtotalvariationregularization[M].MedicalImageComputingandComputer-AssistedIntervention–MICCAI2013.SpringerBerlinHeidelberg,2013:106-114.[25]LiaoCS,ChoiJH,ZhangD,etal.DenoisingstimulatedRamanspectroscopicimagesbytotalvariationminimization[J].TheJournalofPhysicalChemistryC,2015,119(33):19397-19403.[26]BoydS,ParikhN,ChuE,etal.Distributedoptimizationandstatisticallearningviathealternatingdirectionmethodofmultipliers[J].FoundationsandTrends®inMachineLearning,2011,3(1):1-122.[27]WuC,TaiXC.AugmentedLagrangianMethod,DualMethods,andSplitBregmanIterationforROF,VectorialTV,andHighOrderModels[J].SIAMJ.ImagingSciences,2010,3(3):300-339.[28]TsaoJ,BoesigerP,PruessmannKP.k‐tBLASTandk‐tSENSE:DynamicMRIwithhighframerateexploitingspatiotemporalcorrelations[J].MagneticResonanceinMedicine,2003,50(5):1031-1042.[29]LustigM,SantosJM,DonohoDL,etal.ktSPARSE:HighframeratedynamicMRIexploitingspatio-temporalsparsity[C]//Proceedingsofthe13thAnnualMeetingofISMRM,Seattle.2006,2420.[30]MajumdarA,WardRK,AboulnasrT.Non-convexalgorithmforsparseandlow-rankrecovery:ApplicationtodynamicMRIreconstruction[J].Magneticresonanceimaging,2013,31(3):448-455.[31]MajumdarA,WardR.FastSVDfreelow-rankmatrixrecovery:ApplicationtodynamicMRIreconstruction[C].MedicalImaging,m-HealthandEmergingCommunicationSystems(MedCom),2014InternationalConferenceon.IEEE,2014:24-29.[32]ChandrasekaranV,SanghaviS,ParriloPA,etal.Rank-sparsityincoherenceformatrixdecomposition[J].SocietyforIndustrialandAppliedMathematicsJournalonOptimization,2011,21(2):572-596.[33]CandèsEJ,LiX,MaY,etal.Robustprincipalcomponentanalysis?[J].JournaloftheACM(JACM),2011,58(3):1-37.[34]GaoH,RapacchiS,WangD,etal.Compressedsensingusingpriorrank,intensityandsparsitymodel(PRISM):applicationincardiaccineMRI[C].InProceedingsofthe20thAnnualMeetingofISMRM,Melbourne,Australia,2012.P.2242.[35]GaoH,LiL,HuX.CompressivediffusionMRI-part1:whylow-rank?[C].InProceedingsofthe21thAnnualMeetingofISMRM,SaltLakeCity,Utah,USA,2013.P.610.[36]OtazoR,CandesE,SodicksonDK.Low-rankandsparsematrixdecompositionforacceleratedDCE-MRIwithbackgroundandcontrastseparation[C].ProceedingsoftheISMRMWorkshoponDataSamplingandImageReconstruction,Sedona,Arizona,USA.2013.[37]OtazoR,CandèsE,SodicksonDK.Low‐rankplussparsematrixdecompositionforaccelerated-54- dynamicMRIwithseparationofbackgroundanddynamiccomponents[J].MagneticResonanceinMedicine,2015,73(3):1125-1136.[38]CandesEJ,TaoT.Decodingbylinearprogramming[J].InformationTheory,IEEETransactionson,2005,51(12):4203-4215.[39]WrightJ,GaneshA,RaoS,etal.Robustprincipalcomponentanalysis:Exactrecoveryofcorruptedlow-rankmatricesviaconvexoptimization[C].Advancesinneuralinformationprocessingsystems.2009:2080-2088.[40]BeckA,TeboulleM.Afastiterativeshrinkage-thresholdingalgorithmforlinearinverseproblems[J].SIAMjournalonimagingsciences,2009,2(1):183-202.[41]CaiJF,CandèsEJ,ShenZ.Asingularvaluethresholdingalgorithmformatrixcompletion[J].SIAMJournalonOptimization,2010,20(4):1956-1982.[42]CandèsEJ,LiX,MaY,etal.Robustprincipalcomponentanalysis?[J].JournaloftheACM(JACM),2011,58(3):11.[43]LinZ,GaneshA,WrightJ,etal.Fastconvexoptimizationalgorithmsforexactrecoveryofacorruptedlow-rankmatrix[J].ComputationalAdvancesinMulti-SensorAdaptiveProcessing(CAMSAP),2009,61:1-18.[44]PengY,GaneshA,WrightJ,etal.RASL:Robustalignmentbysparseandlow-rankdecompositionforlinearlycorrelatedimages[J].PatternAnalysisandMachineIntelligence,IEEETransactionson,2012,34(11):2233-2246.[45]ZhangZ,LiangX,GaneshA,etal.TILT:transforminvariantlow-ranktextures[M]//ComputerVision–ACCV2010.SpringerBerlinHeidelberg,2011:314-328.[46]ZHANGZ,GANESHA,LIANGX.TILT:TransformInvariantLow-RankTextures[J].Internationaljournalofcomputervision,2012,99(1):1-24.[47]WangH,WangX,ZhouY,etal.Smoothedrandom-liketrajectoryforcompressedsensingMRI[C].EngineeringinMedicineandBiologySociety(EMBC),2012AnnualInternationalConferenceoftheIEEE.IEEE,2012:404-407.[48]DekaB,DattaS.APracticalUnder-SamplingPatternforCompressedSensingMRI[M].AdvancesinCommunicationandComputing.SpringerIndia,2015:115-125.[49]GoldsteinT,O'DonoghueB,SetzerS,etal.Fastalternatingdirectionoptimizationmethods[J].SIAMJournalonImagingSciences,2014,7(3):1588-1623.[50]ZouJ,FuY,ZhangQ,etal.SplitBregmanalgorithmsformultiplemeasurementvectorproblem[J].MultidimensionalSystemsandSignalProcessing,2015,26(1):207-224.[51]NienH,FesslerJA.AconvergenceproofofthesplitBregmanmethodforregularizedleast-squaresproblems[J].arXivpreprintarXiv:1402.4371,2014.[52]LinZ,ChenM,MaY.Theaugmentedlagrangemultipliermethodforexactrecoveryofcorruptedlow-rankmatrices[J].arXivpreprintarXiv:1009.5055,2010.[53]RudinLI,OsherS,FatemiE.Nonlineartotalvariationbasednoiseremovalalgorithms[J].PhysicaD:NonlinearPhenomena,1992,60(1):259-268.[54]QinZ,GoldfarbD,MaS.Analternatingdirectionmethodfortotalvariationdenoising[J].OptimizationMethodsandSoftware,2014(ahead-of-print):1-22.[55]FigueiredoM,DiasJB,OliveiraJP,etal.Ontotalvariationdenoising:Anewmajorization-minimizationalgorithmandanexperimentalcomparisonwithwavaletdenoising[C].ImageProcessing,2006IEEEInternationalConferenceon.IEEE,2006:2633-2636.-55- [56]Bioucas-DiasJM,FigueiredoMAT,OliveiraJP.Totalvariation-basedimagedeconvolution:amajorization-minimizationapproach[C].Acoustics,SpeechandSignalProcessing,2006.ICASSP2006Proceedings.2006IEEEInternationalConferenceon.IEEE,2006,2:II861-II864.-56- 攻读学位期间所取得的相关科研成果一、发表论文[1]马杰,王晓云,张志伟,刘雅莉.一种基于全变分正则化低秩稀疏分解的动态MRI重建方法[J].光电子激光,2016,(1):87-96.[2]刘雅莉,马杰,王晓云,张志伟.一种改进的K-SVD字典学习算法[J].河北工业大学学报,已录用,拟于2016年45卷2期发表。二、参与基金[1]国家自然科学基金(61203245)[2]河北省自然科学基金(F2012202027)-57- -58- 致谢随着论文的撰写接近尾声，研究生生活也马上就要结束了。在这短短的两年半时间中我成长学习了很多。研究生的学习需要有耐心和恒心，需要发挥自身的主动性，从门外汉到慢慢的了解一些东西是一个循序渐进厚积薄发的过程。从刚开始接触课题时的手足无措，不知所云，到后来似乎明白了一些，再到论文的撰写，这些都需要有效时间的付出。每一个阶段都需要时间与精力的投入，而这一切都离不开老师和同学对我的指导、关怀和陪伴。在此论文即将完成之际，我要诚挚的感谢他们。感谢我的导师马杰老师。感谢他督促我专心学习,教导我端正学习态度，认真对待研究生生活；感谢他给了我很多学习锻炼的机会,并在我一筹莫展时给出的宝贵意见，鼓励我继续前行。他治学严谨，对我们要求严格，不仅教会我们做学问，还教会我们做人。他是我们的良师益友，在学习和生活上不遗余力的帮助我们，言传身教鞭策我们前进。感谢我的室友兼闺蜜刘蓓蕾同学。本科加研究生，七年，我们是同学、是室友、是老乡、是闺蜜。七年来我们在“相看两生厌”中培养了深厚的感情与吓人的默契，感谢每一次的调侃斗嘴让单调的生活不那么枯燥；感谢七年来彼此的陪伴和鼓励让我们在最迷茫的时候从没想过要放弃。感谢她对我的支持，对我的论文写作提供了很大的帮助。感谢我的同门刘雅莉、郭海红和陈风华同学，以及506的同学们。是她们的鼓励和交流给了我很多的启发和灵感。我们交流心得、经验，互相鼓励、帮助，郁闷的时候互相吐吐槽，然后满血复活又是活蹦乱跳继续搞科研。是她们的陪伴与鼓励让实验室充满了阳光与正能量。感谢我的父母，是他们的无条件的支持、理解和鼓励，我才有机会跻身于研究生，才有机会结识可敬、可爱的老师和同学们。学校给了我广阔的蓝天，是他们为我插上了梦想的翅膀。感谢他们，对他们的爱是我风雨中搏击翱翔的力量。-59-