磁共振波谱分析——基于 范数的阈值迭代算法在MRS去噪中的应用---毕业论文

磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用摘要磁共振波谱（magnetic resonance spectroscopy，MRS）是一种利用核磁共振现象和化学位移作用，进行系列特定原子核及其他合物定量分析的方法。成为目前唯一无创性研究人体内部器官、组织代谢、生理生化改变的定量分析方法。而分析这些波谱的特征、变化和相互关系等，据此推测人体的病理生理信息，构成了医学波谱分析的主要内容。但是，由于各种原因，实际得到的波谱经常含有噪声和畸变。在考虑谱图实际意义的条件下，对采集得到的信号如何进行优化，去除噪声，提高谱图质量，成为该领域备受关注的问题。本文首先对核磁共振(nuclear magnetic resonance，NMR)、MRS、信号与信号处理等重要概念及其应用领域做出诠释。接着，介绍了时域采集、频域观测的概念，从信号的稀疏性角度引入了范数及基追踪去噪(basis pursuit denoising ,BPDN)。基于上述理论，本文将稀疏性与一致性作为优化的目标，建立去噪问题的数学模型，即转化为基追踪（basis pursuit,BP）问题的求解。在此基础上，实现了基于范数的阈值迭代算法。利用磁共振谱仪实际采集到的一维、二维数据进行仿真实验，并详细比较了仿真结果。另外，对重要参数的设置及影响进行实验与讨论。使用小波与傅立叶作为变换域进行实验，后者的去噪效果更好。通过指标衡量以及对二维、间接维谱图的比较表明，基于范数的MRS去噪方法能有效抑制噪声，恢复出较好的一维和二维谱峰。关键词：MRS；L1范数；阈值迭代；去噪AbstractMagnetic resonance spectroscopy (MRS) is a method by using the phenomenon of nuclear magnetic resonance (NMR) and chemical shifts to do quantitative analysis of a series of specific nuclei and other complex materials. It becomes the only non-invasive study of human internal organs, tissues metabolic, physiological and biochemical changes for quantitative analysis. The analysis of the spectrum characteristics, relationship, changes and so on, which can speculate the pathophysiology information of the human body, constitutes a medical analysis of the main contents of spectroscopy.However, due to magnetic resonance equipment, radio frequency excitation and other reasons, the spectrum acquired from the equipment frequently contains noise and distortion. Considering practical application of spectroscopy, how to remove the noise and improve the quality of spectroscopy is in great concern.In this thesis, concept of NMR, MRS, signals and signal processing with their application background is briefly presented first. Then, from the point view of sparse signals, time and frequency signal processing is introduced. This provides firm foundation for L1 norm and basis pursuit De-Noising.In this paper, we present a MRS denoising method based on fast L1 minimization by applying iterative thresholding. Setting sparsity and consistency as the goal of optimization, denoising mathematical model is introduced into solving basis pursuit problems. By comparing the denoised result of two-dimensional spectroscopy and indirect spectroscopy, simulation results demonstrate that proposed method can suppress the noise well.Key words: MRS; L1-norm; Iterative thresholding; Denoising第一章 引言1.1 选题背景快速多维NMR在异构位移和蛋白质构成的研究中起着重要的作用。对于磁共振波谱的分析和研究，对一些疾病的病理生理变化、早期诊断、预后和疗效的判断都有非常重要的意义。但是，由于磁共振设备和射频激发等原因，实验中得到的波谱经常含有噪声和畸变。通过对设备的改良、采集方式的优化可以降低噪声。更多得是通过对得到的频谱进行去噪处理，从而降低噪声。常用的去噪方法有：傅立叶去噪、小波去噪等。基于稀疏性理论，提出的基追踪去噪方法，常常与阈值迭代、剃度下降等数值算法共同应用于图像去噪领域。但在MRS领域，仍欠缺通用且效果良好的去噪理论与方法。1.2 研究内容及意义本论文研究的基础是一维、二维频谱在采集的过程中会产生噪声，影响谱的质量。如何在现有的理论基础上，结合实际应用情况，设计出合理、有效、易行的去噪方案，是该领域急需解决的问题，也是本文研究的核心内容。结合稀疏性原理、适当的数值算法，在MRS领域进行去噪、重建是当前热门的课题。因而，本论文的研究内容在理论上紧跟国际前沿，研究成果为课题组的长远工作奠定了基础。此外，研究内容在核磁共振的应用领域也有其实际意义。成功地对于设备采集后的信号进行去噪、恢复原谱质量，能更好地进行医学波谱分析，据此推测人体的病例生理信息。从而避免由于噪声导致频谱无法辨识、有效峰被淹没，以至无法准确推断病人病情的情况。3第一章 引言1.3 研究思路采用理论探讨与仿真实验相结合的方式对去噪问题进行研究。对设备中采回的实际信号添加合理的噪声，结合稀疏性理论对其特性进行分析、观察。建立MRS去噪平台，通过仿真实验的结果对算法及其参数、变换域进行调试，不断优化去噪结果，使其能适应各种类型的噪声模型。图1-1为研究思路框图： 图11：研究思路框图2磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用1.4 论文组织结构本文的写作思路从相关知识的介绍，到算法的设计，最后是实验结果的分析。据此，论文的组织结构如下：第二章，介绍核磁共振（NMR，Nuclear Magnetic Resonance）、磁共振波谱（MRS, Magnetic Resonance Spectroscopy）、信号与信号处理等基础概念。第三章，引入在MRS仿真实验中的有关问题，包括时域采集、频域观测、稀疏性等，并简要介绍了两种常用的去噪方法：傅立叶变换与小波变换。第四章，依据该领域的理论知识与实际应用情况，建立基于范数的阈值迭代算法去噪的数学模型： subject to ，定义阈值迭代公式：。在此基础上设计算法，并从理论上做出评价。第五章，对一维、二维信号分别进行实验，调试参数，比较实验结果。分析二维谱图、间接维谱图，分辨有效信息与噪声，使用NMSE指标进行评价。分别使用傅立叶与小波作为变换域进行实验，比较实验结果，判定稀疏性。第六章，总结与展望。通过以上的分析与研究，证明基于范数的阈值迭代算法无论在理论上还是实际的仿真实验中，对于MRS的去噪问题，都十分有效。本论文的意义在于，在特定的应用领域及其背景下，充分应用现有的知识，设计算法，进行实验。基于实验结果的分析和讨论，对今后的工作有一定的启发作用。 3第二章 磁共振波谱分析的基本概念第二章 磁共振波谱分析的基本概念本章介绍了磁共振波谱分析的相关概念、原理及其应用领域，充分阐述了该领域的基础知识，为下文去噪问题的研究提供了背景和依据。2.1 核磁共振（NMR）核磁共振（NMR，Nuclear Magnetic Resonance）是基于原子尺度的量子磁物理性质。具有奇数质子或中子的核子，具有内在的性质：核自旋，自旋角动量。核自旋产生磁矩。NMR观测原子的方法，是将样品置于外加强大的磁场下，现代的仪器通常采用低温超导磁铁。核自旋本身的磁场，在外加磁场下重新排列，大多数核自旋会处于低能态。我们额外施加电磁场来干涉低能态的核自旋转向高能态，再回到平衡态便会释放出射频，这就是NMR讯号。利用这样的过程，可以进行分子科学的研究，如分子结构，动态等。NMR技术即核磁共振谱技术，是将核磁共振现象应用于分子结构测定的一项技术。对于有机分子结构测定来说，核磁共振谱扮演了非常重要的角色，核磁共振谱与紫外光谱、红外光谱和质谱一起被有机化学家们称为“四大名谱”。目前对核磁共振谱的研究主要集中在1H和13C两类原子核的图谱。早期的核磁共振谱主要集中于氢谱，这是由于能够产生核磁共振信号1H原子在自然界丰度极高，由其产生的核磁共振信号很强，容易检测。随着傅立叶变换技术的发展，核磁共振仪可以在很短的时间内同时发出不同频率的射频场，这样就可以对样品重复扫描，从而将微弱的核磁共振信号从背景噪音中区分出来，这使得人们可以收集13C核磁共振信号。核磁共振适合于液体、固体。如今的高分辨技术，还将核磁用于了半固体及微量样品的研究。核磁谱图已经从过去的一维（1D）谱图发展到如今的二维（2D）、三维（3D）甚至四维（4D）谱图，陈旧的实验方法被放弃，新的实验方法迅速发展，它们将分子结构和分子间的关系表现得更加清晰。4磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用在世界的许多大学、研究机构和企业集团，都可以听到核磁共振这个名词，包括我们在日常生活中熟悉的大集团。而且它在化工、石油、橡胶、建材、食品、冶金、地质、国防、环保、纺织及其它工业部门用途日益广泛。在中国，其应用主要在基础研究方面，企业和商业应用普及率不高，主要原因是产品开发不够、使用成本较高。但在石油化工、医疗诊断方法应用较多。12.2 磁共振波谱（MRS）磁共振波谱分析（MRS）是测定活体内某一特定组织区域化学成分的唯一的无损伤技术，是磁共振成像和磁共振波谱技术完美结合的产物，是在磁共振成像的基础上又一新型的功能分析诊断方法。70年代中期MRS就已应用于人和动物组织器官的活体检测，近年来随着MRS技术的迅速发展和美国FDA对MRS技术的认可，MRS已从实验室研究转入临床应用阶段，磁共振波谱对一些疾病的病理生理变化、早期诊断、预后和疗效的判断都有非常重要的意义。现在用于MRS检测的核素有1H、13C、19F、23Na、31P。磁共振波谱分析目前主要应用于脑部、心脏、骨骼肌和肝脏等方面的研究，以脑部最为广。脑部磁共振波谱研究较多的有脑梗死、脑肿瘤、脑白质和脑灰质疾病、癫痫和代谢性疾病等，尤其是颅脑肿瘤研究较多，对脑肿瘤与非肿瘤性病变鉴别、脑肿瘤良恶性鉴别、恶性肿瘤分级、肿瘤术后复发与坏死的鉴别、原发与转移瘤的鉴别等均有很大的临床应用价值，此外，还能鉴别颅咽管瘤与垂体瘤，脑内肿瘤与脑外肿瘤，确定脑室内的中枢神经细胞瘤等。在心脏方面的应用主要是在心肌缺血、心肌病等心肌代谢方面的研究。肝脏31P-MRS主要研究包括肝代谢性疾病、肝炎肝硬化及肝肿瘤等。MRS能提供前列腺组织的代谢信息有助于鉴别前列腺癌和前列腺增生。MRS还能无创性地检测骨骼肌磷脂代谢和能量代谢的代谢产物及细胞内pH值，研究骨及软组织肿瘤的磷脂和能量代谢的异常变化。因此，在医学上对磁共振波谱进行分析和研究有重大意义。分析这些波谱的特征、变化和相互关系等，据此推测人体的病理生理信息，构成了医学波谱分析的主要内容。25第二章 磁共振波谱分析的基本概念2.3 信号与信号处理 信号是信息的载体，任何携带信息的物理量都可以称为信号。高速度、高精度的信号处理已经成为了一套完备的理论体系，其中包括各种快速的和优良的算法。随着新理论和新算法的提出，计算机技术、超大规模集成电路技术等的飞速发展，信号处理的理论和技术还在不断地丰富、完善和提高。从信号中提取尽可能多的有用信息，增强信号的有用分量，估计信号的特征参数，识别信号的特征，同时抑制或消除不需要的，甚至有害的信号分量等，是信号处理的根本目的。为此，需要对信号进行分析和变换、扩展和压缩、滤波、参数估计、特性识别等加工，通称为信号处理。信号处理的实现方法分为软件实现、硬件实现和硬软件结合实现。一，软件实现方法。按照信号处理的原理和算法，编写程序（如用Matlab语言），在通用计算机上实现信号处理。这种方法的优点是灵活方便，程序中的参数容易改变，程序也便于修改，缺点是计算速度较慢，达不到高速实时处理的目的。二，硬件实现方法。按照信号处理的原理和算法，设计硬件实现方案，从硬件角度完成信号处理。其优点是运算速度快，实时性强，但不够灵活。三，硬软件集合实现方法。3一般来说，信号可以表示为独立自变量的函数。磁共振波谱就是一种信号，对于磁共振波谱的分析和研究过程，可以看作信号处理过程。磁共振波谱作为信号，按照信号中自变量地数目可分为一维信号和多维信号。在应用领域中用到的主要是多维信号，一维信号处理可以作为研究多维信号处理的基础。研究的谱是离散时间信号，其自变量时间在信号的定义域内是离散的。可以通过对连续时间信号的采样获得。信号强度是核磁共振谱关注的重要信息。处于相同化学环境的原子核在核磁共振谱中会显示为同一个信号峰，通过解析信号峰的强度可以获知这些原子核的数量，从而为分子结构的解析提供重要信息。在此，研究的主要问题是磁共振波谱的去噪问题，也就是一维、多维信号的去噪问题，采用的方式是软件实现。6磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用第三章 频谱的去噪问题本章以磁共振波谱为背景，通过阐述时域频域、稀疏性等问题，说明频谱去噪的必要性及实施条件。并介绍了两种基本的去噪方式：傅立叶变换与小波变换。3.1 时域采集与频域观测对信号的刻画往往采取两种最基本的形式:时域形式和频域形式。从信号的时域描述，可分辨出信号的某些特征，如变化的快慢、有值范围、连续与离散等。然而信号频域描述更是人们所关注的，有些在时域内看不清的问题，在频域则一目了然。 3.1.1 时域采集时域（时间域）自变量是时间,即横轴是时间,纵轴是信号的变化。其动态信号是描述信号在不同时刻取值的函数。对信号进行时域分析时，有时一些信号的时域参数相同，但并不能说明信号就完全相同。因为信号不仅随时间变化，还与频率、相位等信息有关，这就需要进一步分析信号的频率结构，并在频率域中对信号进行描述。动态信号从时间域变换到频率域主要通过傅立叶级数和傅立叶变换实现。周期信号采用傅立叶级数，非周期信号采用傅立叶变换。对于MRS波谱，信号采集在时域进行的。但是，由于实际物理设备的限制，采集过程会造成信号含有噪声和畸变。这些都将直接影响对谱峰信息的化学分析，这也就是去噪工作的意义所在。在仿真实验中，对于MRS的时域信号主要考察两种噪声来源，一维的信号截断和二维笛卡儿模板采集。 7第三章 频谱的去噪问题图31：一维信号截断图 图32：二维信号笛卡儿采集模板 图31为一维信号的截断图。对于12016的一维的信号，进行了30的截断。即保留原信号长度，将信号尾部数据全部置零。此种方式，常用于衰减的MRS信号。图32为二维信号的采集模板。模板大小256256，适用于同样大小的二维数据。采集方式是笛卡儿采样。 3.1.2 频域观测频域（频率域）自变量是频率,即横轴是频率,纵轴是该频率信号的幅度,也就是通常说的频谱图。频谱图描述了信号的频率结构及频率与该频率信号幅度的关系。 将时域信号进行傅立叶变换，转换到频域。在MRS领域，频谱是感兴趣的研究对象。图33和图34分别为一维信号时域与频域在截断前后的对比图、二维信号采集前后的频谱图： 8磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用Fourier 图33（b）图33（a）Fouriersampling图33（d）图33（c）图33：一维信号截断前后的时域和频谱。（a）截断前的时域图，（b）截断后的时域图，（c）截断前的频域图，（d）截断后的频域图。图33（a）为实际设备中采回的完整的12016的一维信号，信号呈衰减趋势。经过傅立叶变换到频域，得到谱图33（c）。4个明显的谱峰是应关注的重要信息，包括峰的高度与平滑度。原信号经过截断，得到图33（b）。截断是仿真中对于一维信号常用的采集方式。由于信号是衰减的，因此去掉信号尾部的信息，将其置零，保留长度。截断后的时域信号经过傅立叶变换，得到谱图33（d）。可以看到，由于截断的影响，谱峰明显降低，且在底部出现多余的小峰，曲线变得不光滑，这些都是必须去除的噪声。 9第三章 频谱的去噪问题图34（a） 图34（b）图34（c） 图34（d） 图34：二维信号采集前后，二维谱、间接维谱。（a）截断前的频域二维谱图,（b）截断后的频域二维谱图，（c）截断前的频域间接维投影图，（d）截断后的频域间接维投影图。 在二维MRS中，关注的是二维谱图及间接维投影。图34（a）为原二维信号的二维谱图，显示了重要信息的分布与强弱。图34（b）为经过笛卡儿采集后的二维谱图，增加了很多多余信息且中心能量减弱。图34（c）为原二维信号在间接维上的谱图，有一个很明显的主峰。图34（d）为经过笛卡儿采样后的间接维谱图，主峰的高度大大降低。 综上所述，无论是一维还是二维信号，都形成了噪声,淹没了有效信号。下面，我们对这些带噪谱进行稀疏去噪。 10磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用3.2 稀疏性3.2.1 稀疏分解与表示基于稀疏分解的信号去噪，将被噪声污染的信号分解成信号的稀疏成分和其他成分。稀疏成分对应于信号中的有用信息，其他成分对应于信号中的噪声。由的稀疏成分重建信号，可达到去除信号噪声的目的。,X是信号，是稀疏变换的库。X的效果取决于库中的原子。假设研究的信号为f，信号大小为。若将信号分解在一组完备正交的基上，则这组基的数目应为。由于基的正交性，在由信号所组成的空间中的分布是稀疏的。信号的能量在分解以后将分散分布在不同的基上，这种能量分布的分散将导致用基的组合来表示信号时，表达不简洁。即信号表示不是稀疏的，这不利于信号的处理，如识别和压缩等。为了得到信号的稀疏表示，基的构造必须使基在信号组成的空间中足够密，基的正交性不再得到保证，所以此时的基也不再是真正意义上的基了，而改称为原子。由这些原子组成的集合是过完备的，称为过完备库。信号在过完备库上的分解结果一定是稀疏的。基于信号稀疏分解的去噪方法近似认为，信号中的稀疏成分对应于有用信息，而其他成分则对应于信号中的噪声。传统的去噪是直接从信号中去除噪声。与此不同，基于稀疏分解的去噪方法是先从信号中提取出稀疏成分，然后利用稀疏成分重建信号，则重建的信号即为去除噪声后的信号。式中为信号的稀疏成分；为其他成分，即信号中去除稀疏成份后的残差部分。4 11第三章 频谱的去噪问题3.2.2 基追踪（Basis Pursuit问题）通过变换来表示信号由于在许多应用中具有简单性和有效性而一直受到重视。传统的信号表示理论大多基于非冗余的正交基函数的变换，如傅立叶变换、小波变换等。但傅立叶变换无法表达信号的时频局域性质，而小波变换虽在处理具有点状奇异性的信号时具有良好的特性,但由一维小波张成的可分离小波只具有有限的方向,不能“最优”表示含线或者面奇异的高维函数。基于这些变换，普遍采用的是超完备(overcomplete)冗余表示，其基本思想是基函数用称之为字典的超完备的冗余函数系统取代,字典中的元素被称为原子，信号由原子的线性组合来表示。其中原子的数目比信号的维数大得多，由此产生了冗余。由于这种超完备性，就有很多表示信号的方法，其中具有最少系数(最稀疏)的表示是最简单的，也是我们认为最优的一种表示方法。找到最好的、也就是最稀疏的信号的信号表示，等同于解决下述问题：minimize subject to 其中，是序列中非零项的个数。从一个随机冗余原子库中寻找信号的稀疏扩展是一个NP难问题。由此，提出解决下述稍有差别的问题：minimize subject to 这一微妙的修改，转化了问题的性质。最小化范数，被称为Basis Pursuit(BP),是较为简单的问题。即把最小化范数转化为最小化 范数,使得问题简单化。对于最小化范数,我们可以通过线性规划方法来解决。基追踪去噪(Basis Pursuit Denoising ,BPDN)是BP算法的推广,演变成无约束优化问题：其中参数控制着表示信号s时的失真度。该式具有二次规划结构，可以通过很多方法进行求解。 12磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用基追踪去噪方法是信号稀疏表示领域的一种新方法。在分解信号的同时，也自动地提取出噪声成分,完成了对信号的去噪处理。它寻求从完备的（过完备）函数（基）集合中得到信号的最稀疏的表示，即用尽可能少的基精确地表示原信号，从而获得信号的内在本质特性。基追踪方法采用表示系统的范数作为信号稀疏性的度量，通过最小化范数将信号稀疏表示问题定义为一类有约束的极值问题，进而转化为线性规划问题进行求解。5基追踪方法采用了一种新的线性规划算法内点算法，使得在所有的字典向量中极小化一个全局目标函数。BPDN为本文所研究的去噪算法提供了理论基础，是其数学模型的前身。3.3 传统去噪方法3.3.1 快速傅立叶变换长期以来，傅立叶变换是信号去噪的主要手段，利用傅立叶变换给信号去噪等价于信号通过一个低通或高通滤波器。它利用傅立叶变换把信号映射到频域内加以分析,但是傅立叶变换在给信号去噪的同时，也模糊了信号的位置信息，并且在实际应用中，大多数信号都是非平稳的，非平稳信号谱沿时间轴无限扩展，此时傅立叶变换的基函数很难与其匹配，这给傅立叶变换对非平稳信号去噪带来了困难。设原始信号为，其傅立叶变换为：其傅立叶逆变换为：确定了在整个时间域上的频谱特性。傅立叶变换将信号的时域特征和频域特征联系起来，能分别从时域和频域上观察信号，但不能把二者有机结 13第三章 频谱的去噪问题合起来。无论在时域或是在频域，傅立叶变换都是定义在R上的全局积分。在用傅立叶变换对原始信号进行处理时，识别出的频率在什么时候产生并不知道，因此傅立叶变换是一种全局的变换，不能反映某一局部时间内信号的频谱特性，也即傅立叶变换在时间域上没有任何分辨率。这样在利用傅立叶变换做信号分析时就面临着一对矛盾：时域和频域的局部化矛盾。离散傅里叶变换（Discrete Fourier Transform,简称DFT），它使数字信号处理可以在频域采用数字运算的方法。DFT有多种快速算法，统称为快速傅里叶变换（Fast Fourier Transform,简称FFT）,能极大地提高数字信号处理的实时性和便捷性，为实时频谱分析提供了可能，成为分析离散信号和系统的有力工具。傅里叶频谱分析的基础是傅里叶变换，用傅里叶变换及反变换可以求的傅里叶谱。单相检波的核磁共振，经采样和累加后得到的正y恢复到它们的基线，可以由下式计算基线纵坐标的平均值：值减去值得到值： 值接着被转换成频域中相应的傅立叶正弦变换式和余弦变换式。可用下列方程描述傅立叶变换: 和 它们的积分实际上被定义为对离散的、截断了的系列逐项求和。利用FFT法在频域处理噪声的优点在于节省时间及大量的内存地址。63.3.2 小波变换小波分析是把一个信号分解成由原始小波经过移位和缩放后的一系列小波,因此小波是小波变换的基函数,即小波可以表示一些函数的基函数。小波(wavelet)是定义在有限间隔且平均值为零的函数。 14磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用设（表示平方可积的实数空间，即能量有限的信号空间），其傅立叶变换为。当满足容许性条件： 时，我们称为一个基小波或母小波（Mother Wavelet）。将母小波经伸缩、平移后，就可以得到一个小波序列：其中称为尺度因子，称为平移因子。若尺度因子a和平移因子b的取值连续变化，是依赖于参数a和b的连续小波基函数，这样就得到了连续小波变换的定义。设是基小波，是基小波经伸缩和平移后得到的连续小波序列，则信号的连续小波变换（Continuous Wavelet Transform,简记为CWT），可以定义为：其中是的复共轭。进一步有，在基小波的Fourier变换满足容许性条件时，若信号在处连续，可以从连续小波变换中恢复出原始信号，即连续小波变换的逆变换存在：在中变量a反映函数的尺度，变量b检测沿时间轴或位置轴平移的位置。一般情况下，母小波函数的能量集中在原点，连续小波函数的能量集中在b 点。 15第三章 频谱的去噪问题小波变换就是将任意信号表示为小波的叠加,信号的小波叠加表示将信号分解为不同的尺度级。在每一尺度级,该信号又在这一尺度级对应的分辨率下被分解。尺度级对应着频率,频率越高,对应的分辨率越高,换句话说,我们将信号分解成了不同的频率通道成份,并将每一频率成份又按相位进行了分解频率越高者,相位划分越细;反之则越疏。小波等多分辨率分析的方法，能够克服传统傅立叶变换只能得到频率信息，但不能定位频率在时域/空域上定位频率发生的位置的缺点。在时域、频率域都具有表征信号局部特征的能力,一些在原时间域上因混叠而不易观察的信号特征,能够在频率域的某个尺度上进行分离,并且得到显著地体现,达到获取有效信号的目的，从而在剔除噪声的同时，保留原波谱的结构和细节信息。小波去噪通常通过以下3个步骤予以实现：1小波分解：选定小波和分解的层数，对信号进行小波分解，在各层上得到逼近部分（低频带）和细节部分（高频带）。噪声部分通常包含在高频中。2设定各层细节的阈值，取出各层中的细节部分，根据选定的阈值，对得到的小波系数进行阈值处理。3小波逆变换重构信号，得到去噪后的频谱。12 16磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用第四章 基于l1范数的阈值迭代算法本章在前两章的基础上，提出了应用于MRS去噪的基于l1范数的阈值迭代算法。写作思路依据整个流程，从数学模型的建立到算法设计，遵循从抽象到具体的过程。4.1 建立数学模型4.1.1 优化目标该算法将稀疏性与时域上的一致性作为优化目标。稀疏性，即找到稀疏域，满足范数最小值约束；一致性，即恢复后的信号还原到时域上要与原信号保持一定范围内的一致。在上一章稀疏性一节中，已经提到关于最稀疏问题的欠定方程，在这里，重新定义为： subject to (1)y是时域上得到的观测值，带噪。A是算子，下文将有具体解释。x是在K空间上未知的向量。代表x中非零元素的个数，即范数。这是一个非凸组合最优问题，即为一个找最稀疏方案的NP难问题。通常，有三种方法解决NP难问题：如果输入值非常小，一个时间复杂度非常高的算法可以符合要求，尽管这种方案并不完善；分离在多项式中有较高灵活性的一些特殊案例；找到一个近似方程，用以代替原问题。在这里，采用第三种方案，将范数修正为范数。问题转化为： subject to （2） 17第四章 基于范数的阈值迭代算法其中，是x中非零元素的绝对值之和，即范数。根据（2），提出由作为约束条件的求最小值的目标函数： subject to （3）这个方程也可以写成拉格朗日形式： （4）公式和本质上是相同的，都是将稀疏性和一致性同时作为优化的标准。已经在信号处理领域被提出，称为基追踪去噪（basis pursuit denoising,BPDN）。基于实际算法操作过程中的易用性，将公式变形为： subject to （5）式（5）与等同，优化问题和可以被视为线性方程。中，优化目标的意义为：在稀疏域上经过去噪处理的信号，返回到时间域上与带噪的观测值的残差必须控制在以内。约束条件：利用范数判定稀疏性，稀疏域上的非零值的绝对值必须尽可能小。 4.1.2 阈值迭代阈值迭代在信号处理的应用已经有很长的一段历史。它的优点在于简单易操作，且有启发作用。阈值迭代又分为软阈值迭代和硬阈值迭代。将转换为一个二次方程，其形式类似于线性方程： subject to （6）是两个向量的内积。是一个矩阵向量的积。，代表向量x中的项必须为非零值。 18磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用如果A是正交的，就等同于： （7）由此，得到阈值非线性特性，用以解决最小值问题： （8）软阈值操作定义为： （9）硬阈值操作定义为： （10）软阈值与硬阈值在每次迭代中，都将低于阈值的值置零。二者的区别在于，硬阈值操作完全保留大于阈值的项，而软阈值迭代对原值进行一定的变换。 图4-1：软阈值迭代示意图 19第四章 基于范数的阈值迭代算法4.2 算法设计基于以上数学模型和分析，为了去除NMR频谱在采集过程中产生的噪声，将问题转换为解决优化问题： subject to （11）y是给定的的观察值向量。,F作为部分傅立叶变换，将原不带噪的时域信号转换为带噪的频域信号。作为其反变换，是F的逆过程。S是从频域到变换域的转变，稀疏域的选择是任意的。如果是小波，则进行正向的小波变换；甚至可以是频域本身，此时S操作即为None操作。是S的逆过程。相当于一个给定的的矩阵，pm。使用如下迭代方式： （12）是第k次迭代得到的谱，是阈值为t的阈值迭代函数，阈值随着迭代次数k的增加而降低。4.2.1 算法的操作流程 对于受到噪声污染的信号，转换到稀疏域上进行阈值迭代，每次迭代改变阈值，并将迭代后的结果返回到时域上做一致性检验。迭代的起点是零，终点至残差到达一定比例。 20磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用算法的伪代码： Input: matrix A, observation vector y,Output: solution of Algorithm:Initialization:Residual Correlation Threshold Solution of length nWhile : update solution compute residual compute correlation update threshold end while 21第四章 基于范数的阈值迭代算法4.2.2 算法的评价 高维度的NMR波谱分析非常耗时，需要数日甚至数周的时间。本文提出的去噪方案从两个角度诠释去噪问题： （1）稀疏性。使用频域和小波域这两个稀疏空间来对信号进行稀疏表示。（2）非线性去噪。作为优化标准的应用，产生一组稀疏的系数表示。对于算法的限制，主要来自于计算的复杂度。该算法的核心即以稀疏性和信号保真度为优化目标，以阈值迭代为操作方式，将信号变换到一个可以只用少数有意义系数就能表示从而去噪的空间。因此，积的选择，越稀疏越好。基于范数的阈值迭代算法可以作为多种用途的解决方案，尤其是原本需要耗费大量时间的应用领域。对于类似的问题，潜在的解决方案是利用少数非零值的稀疏性。在大型的应用系统中，一个近似的最稀疏解决方案是通过阈值迭代实现的，这种方案的运行速度远远快于其他方式。因此，结合上述内容，提出了基于范数的阈值迭代算法，解决快速多维NMR波谱问题，用于去除实际应用中可能产生的噪声。这是一种快速而有效的去噪算法。该算法的关键，在于算子A的定义，其中包含了部分傅立叶及变换域。在时域上模拟加噪，在频域上观测谱图，在稀疏域上进行迭代操作。A算子很好地连接了以上三个操作，将作为优化目标的范数及信号的一致性，与阈值迭代算法无缝地结合起来。每次的迭代过程，既实现了一致性的检验，又对信号的稀疏性做出了衡量。 22磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用第五章 实验结果分析运用基于l1范数的阈值迭代算法，进行MRS去噪得到了很有效且快速的效果。图5-1是一维信号的去噪效果。 图5-1(a) 图5-1(b) 图5-1(c) 图5-1：一维去噪效果对比图。（a）是不带噪的原谱，（b）是带噪的谱，（c）是去噪后的谱图5-1(a)(c)，依次是不带噪的谱、截断后带噪的谱、去噪后的谱。去噪的效果非常明显,在谱底部由于截断产生的大量噪声被去除，主要的峰被恢复。根据算法设计思路，影响去噪效果的参数，有以下四个：噪声的形式,变换域的选择,阈值的衰减速率，迭代的终止条件即残差终值。噪声的模式不是本文讨论的重点，在以下实验中，固定这一参数。因此，主要研究阈值衰减速率、残差终值以及变化域的选择对于去噪效果、时间复杂度、收敛程度等因素的影响。下文主要对二维信号进行实验讨论。 23第五章 实验结果分析5.1 阈值衰减率与残差终值本实验的基本设置为：源数据大小：128128噪声模式：固定变换域：频域 图52（a） 图52（b） 图52（c） 图52（d） 图52：参数对比实验的原数据及噪声数据。（a）是原二维谱，（b）是带噪声的二维谱，（c）是原间接维投影，（d）是带噪的间接维投影由图5-2，可以观察到，带噪声的二维谱、间接维谱与原谱有很大差别：二维上出现了许多噪点，有效的点能量降低、位置发生偏移；间接维上有效峰的高度降低，出现了大量的旁瓣。下面的实验就是对以上噪声进行去噪处理的参数探讨，分别从迭代收敛趋势、二维谱、间接维谱、NMSE指标的角度对去噪程度进行评价。 24磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用5.1.1 阈值衰减速率残差终值设为10-6，将阈值衰减速率依次设置维0.1，0.6，0.9 图53（a）仅迭代6次，残差衰减速度非常快。图53（a）图53（b）迭代进行22次，曲线较上图平滑一些。图53（b） 图53（c）迭代98次，衰减速率比较小。图53（c）图53：不同阈值衰减速率的迭代趋势对比。(a)阈值衰减速率0.1,(b)阈值衰减速率0.6,(c)阈值衰减速率0.9 25第五章 实验结果分析 图54（a）二维谱去噪效果很差，仍存在很多噪点。图54（a）图54（b）二维谱去噪效果较好。噪点较上图减少许多，但与原谱相比，仍残余许多噪点。 图54（b）图54（c）二维谱较上图，更接近原值，噪点更少。 图54（c）图54：不同阈值衰减速率的二维谱对比。(a)阈值衰减速率0.1,(b)阈值衰减速率0.6,(c)阈值衰减速率0.9 26磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用图55（a）间接维投影有效谱谱峰低，去噪效果差。图5-5(a)图55（b）间接维投影去噪效果比上图略好。主要的峰有被突出，但幅度很低，峰的底部不平滑。图5-5(b)图55（c）间接维投影去噪效果较上图更好一些，有效峰的高度大幅度增加。图5-5(c)图55：不同阈值衰减速率的间接维投影对比(a)阈值衰减速率0.1,(b)阈值衰减速率0.6,(c)阈值衰减速率0.9 27第五章 实验结果分析采用规一化均方差NMSE (Normalized Mean Square Error)衡量含噪谱与原谱之间的差距。一般在01之间，值越小，说明差异越少，质量越好。NMSE的定义：其中，为原始信号，为带噪的信号或去噪后的信号。13表5-1：不同阈值衰减速率的NMSE值对比 单位：100%阈值衰减速率噪声谱去噪后的谱0.10.60.9NMSE0.58720.53690.44270.3743综上所述，阈值衰减速率的取值在01之间。值越小，迭代次数越少，程序的运行速度越快。但衰减曲线不平滑，去噪效果差。当其值非常接近1时，时间复杂度非常大，运行时间长。因此，将值设置为0.9。5.1.2 残差终值 迭代的终止条件，即最终残差所占的比例，是一个非常重要的参数，控制去噪的程度。阈值衰减速率设为0.9，将残差终值分别设为10-1,10-3,10-6,10-10进行实验。 图56显示了不同残差终值取值下的，迭代次数与当次残差比例的关系，由于阈值衰减速率值固定，因此曲线的衰减趋势是固定的。 28磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用 图56（a） 图56（b）图56（c） 图56（d）图56：不同残差终值的迭代衰减对比。(a)残差终值10-1,(b)残差终值10-3,(c)残差终值10-6,(d)残差终值10-10图56（a）(d)的残差终值分别为10-1,10-3,10-6, 10-10。曲线的下降趋势相同，但终止条件不同，因此迭代次数逐渐增加。迭代次数的增加，意味着耗费更多的时间。因此，随着残差终值的减小，去除的信息减少，保留的信息赠多。去噪不足，显然会造成大量噪声残余；去噪过度，可能会意外去除非噪信息。 图57比较了不同残差终值取值下的二维谱。二维谱的观察，重点在于各点的分布及强弱，注意与原谱进行对比，分辨有效信息与噪点。 29第五章 实验结果分析图57（a） 图57（b）图57（c） 图57（d）图57：不同残差终值的二维谱对比。(a)残差终值10-1,(b)残差终值10-3,(c)残差终值10-6,(d)残差终值10-10图57（a）(d)的残差终值分别为10-1,10-3,10-6, 10-10。（b）(d)差别不明显，分布情况基本一致。比较（a）与（c）,在图中标出的是有明显差别的位置。（c）比（a）多出一些细节，而这些细节是原谱上有的，可以认为是有效的而非噪点。因此，可以认为（a）去噪过度。 30磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用 图58（a） 图58（b）图58（c） 图58（d）图58：不同残差终值的间接维投影对比。(a)残差终值10-1,(b)残差终值10-3,(c)残差终值10-6,(d)残差终值10-10 图58比较了不同残差终值取值下的间接维投影。间接维投影是实际应用中对于二维信号考察的重点。图58（a）(d)的残差终值分别为10-1,10-3,10-6, 10-10。（b）(d)差别不明显，分布情况、幅度、平滑度基本一致。比较（a）与（c）,在图中标出的是有明显差别的位置。显然，（c）的整体谱峰都要高于（a），也与原间接维谱更接近。说明（a）的能量损失较大。 31第五章 实验结果分析最后，通过NMSE指标进行比较，结果如下：表5-2：不同残差终值的NMSE值对比 单位：100%残差终值噪声谱去噪后的谱10-110-310-610-10NMSE0.58720.39190.37440.37430.3744从迭代次数、二维谱、间接维谱、NMSE指标，对不同残差终值去噪能力进行考察，可以得出以下结论：残差终值是非常重要的参数，控制去噪的程度。值过大，会遗失重要信息；值过小，会残留噪声，且耗费较高的时间复杂度。当低于某个临界值时，与原谱的一致性基本处于稳定状态。因此，一般取10-6左右。5.2 小波稀疏与傅立叶稀疏对比稀疏域的选择，直接影响去噪效果。小波和傅立叶是两种常用的变换域，究竟哪种去噪效果更好，哪种更稀疏。在本节，对比它们的收敛趋势、去噪情况及稀疏性，指出哪种方案更优。参数选择：阈值衰减速率设为0.9 残差终值设为10-6 固定噪声模式源数据： 512512 32磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用5.2.1 收敛趋势对比 图59（a） 图59（b）图59：小波与傅立叶对比收敛趋势。(a) 采用傅立叶变换去噪，(b)采用小波变换去噪 由于迭代的终点是相同的，由图59可以看出，使用小波作为稀疏域即变换域比傅立叶迭代次数略少一些，但小波变换本身比傅立叶变换复杂，耗时多。因此可以认为，二者所耗费的时间复杂度是相近的。5.2.2 去噪效果对比与上节实验相同，通过对二维谱、间接维谱、NMSE值的考察，来判别去噪效果的优劣。 33第五章 实验结果分析5.2.2.1 二维谱图510（a） 图510（b） 图510（c） 图510（d） 图510：小波与傅立叶对比二维谱。(a)原谱，(d)带噪谱，(c) 采用傅立叶变换去噪，(d)采用小波变换去噪显然，含噪的(b)比(a)多出了许多噪点，且有效点的能量被削弱，位置发生偏移。在频域进行去噪的(c)去除了几乎所有噪点，基本恢复到了(a)的分布情况；在小波域去噪的(d)在去除噪点的同时，也去除了大量的有效点，二维谱上只剩下少量信息。 34磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用5.2.2.2 间接维投影 图511（a） 图511（b）图511（c） 图511（d） 图511：小波与傅立叶对比间接维投影。(a)原谱，(d)带噪谱，(c) 采用傅立叶变换去噪，(d)采用小波变换去噪 图511（a）（d）分别为不带噪的原间接维谱、带噪的间接维谱、在频域去噪后的间接维谱、在小波域去噪后的间接维谱。得到的结果与二维谱类似。加噪后的谱，噪声与有效峰混合在一起，难以分辨。在频域去噪的效果较好，基本恢复原谱分布。在小波域去噪有一定效果，但没有完全分离噪声和有效峰：噪声的幅度仍然较高，部分有效峰被淹没。 35第五章 实验结果分析5.2.2.3 NMSE值 表5-3：傅立叶与小波的NMSE值对比 单位：100%噪声谱去噪后的谱傅立叶小波NMSE0.82730.19250.5670NMSE值的比较同样说明了，使用傅立叶变换即在频域去噪，得到的谱与原谱的一致性更高。 5.2.3 稀疏性对比两种变换的比较，实质上是稀疏性的比较。即信号在小波域上更稀疏，还是在频域上更稀疏，决定了去噪效果的优劣。这是由于，本文提出的去噪算法是基于范数的阈值迭代算法，也就是基于稀疏性的算法。因此，变换域的选择必须以稀疏性为前提。选择的域越稀疏，去噪效果就越好。下图是归一化后，傅立叶系数与小波系数的稀疏程度的对比。显然，傅立叶的曲线在几乎所有位置，幅度都低于小波的曲线。表明了，傅立叶能用更少的系数来表示原信号，频域更稀疏。这也就解释了为什么前者的去噪效果更好。图5-12：小波与傅立叶稀疏性对比 36磁共振波谱分析基于范数的阈值迭代算法在MRS去噪中的应用第六章 总结与展望 磁共振波谱的去噪一直是核磁共振领域备受关注的一项课题。如何运用现有的去噪原理，并结合实际应用领域，将其转化为切实可行且有效率的算法，是十分有意义的工作。综合全文，可以得出以下结论：第一，基于范数