核磁共振驰豫信号反演问题
- 格式:pdf
- 大小:427.42 KB
- 文档页数:21
核磁共振弛豫信号的多指数反演核磁共振弛豫信号的多指数反演王武蕾王梦佳刘⽂英2010年8⽉24⽇摘要本⽂针对核磁共振信号处理问题,系统地研究了多指数反演的⽅法。
⾸先,本⽂构造了⼀个具有双峰特征的模拟弛豫时间谱,并在不同信噪⽐下,进⾏多指数正演计算,得到反演所需的原始回波数据,正演模拟结果如图2所⽰。
为求解2T时间谱,本⽂从第⼀类积分⽅程的反演求解⼊⼿,在预先知道2T弛豫时T分布的情况下,推导出适合于核磁共振弛豫时间多指数反演的两种算法——奇异间2j值分解⽅法和阻尼最⼩⼆乘算法。
对算法的具体实现过程进⾏了详细的论述,从理论和实例处理两⽅⾯分析讨论了两种算法的优缺点。
针对这两种算法,分析了噪声对其解谱的影响,并确定了反演的条件。
在进⾏核磁共振信号的多指数反演处理时,应优先≥时,可选⽤奇异值分解选⽤阻尼最⼩⼆乘算法,只有当原始数据的信噪⽐SNR40反演算法。
若预先给定弛豫时间分布,由于布点间断、不连续导致所得结果分辨率较低,结果不是很理想。
本⽂使⽤基于差分进化算法对核磁共振弛豫信号进⾏多指数反演,在将反演问题转化为带⾮负约束的⾮线性拟合优化问题的基础上,根据⾮负最⼩⼆乘⽅法(NNLS)所确定的弛豫时间的⼤致组分数,直接利⽤差分进化算法进⾏反演进化,从⽽避免了对弛豫时间2T分布的确定。
仍利⽤所构造的模拟时间谱,在不同的信噪⽐下使⽤该⽅法进⾏指数反演,从计算精度、抗噪能⼒和计算速度上说,所得结果都⽐上述算法优越,尤其是在信噪⽐较低的情况下,仍保证了弛豫时间谱的真实性。
在对反演算法进⾏研究的基础上,针对提供的实际核磁共振数据,分别使⽤上述三种算法对原始信号和去噪信号进⾏多指数反演,反演结果见图11。
根据结果分析,阻尼最⼩⼆乘算法和差分进化算法具有较好的抗噪能⼒和较⾼的计算精度,奇异值分解算法只有在信噪⽐较⾼时才能得到精确的2T时间谱,这与问题⼀、⼆中的分析结果保持⼀致。
同时,本⽂使⽤对数均匀分布和2的幂指数分布对弛豫时间布点,⼆者所得结果基本相同,对反演计算没有产⽣重要影响,具体数值见表4。
承诺书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.我们完全明白,在竞赛开始后参赛队员不能以任何方式(包括电话、电子邮件、网上咨询等)与队外的任何人(包括指导教师)研究、讨论与赛题有关的问题。
我们知道,抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料(包括网上查到的资料),必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。
我们郑重承诺,严格遵守竞赛规则,以保证竞赛的公正、公平性。
如有违反竞赛规则的行为,我们将受到严肃处理。
我们参赛选择的题号是(从A/B/C/D中选择一项填写): C我们的参赛报名号为(如果赛区设置报名号的话):所属学校(请填写完整的全名):中南大学参赛队员(打印并签名) :1. 杨金风2. 胡静3. 周盛泽指导教师或指导教师组负责人(打印并签名):日期: 2010年 08月 23 日赛区评阅编号(由赛区组委会评阅前进行编号):编号专用页赛区评阅编号(由赛区组委会评阅前进行编号):赛区评阅记录(可供赛区评阅时使用):评阅人评分备注全国统一编号(由赛区组委会送交全国前编号):全国评阅编号(由全国组委会评阅前进行编号):核磁共振驰豫信号反演问题摘要核磁共振在测井技术和岩心分析中已经得到广泛应用。
本文主要讨论核磁共振驰豫信号反演问题。
针对横向驰豫时间2T 的不同布点方式,建立数学模型,通过NNLS 算法和差分进化算法进行2T 谱的多指数反演,并对结果进行相关讨论。
考虑到原始数据含有噪声,在前两问中均采用去噪信号进行求解,从而使所得2T 谱更符合要求;第三问中,给定一个双峰2T 谱,并利用matlab7.0产生随机噪声生成模拟信号,通过对比反演结果,总结测量误差对前两问算法结果的影响。
在预先给定弛豫信号分布的情况下,在不同布点方式中,根据2T 的时间分布区间,分别合理确定m 值进行取点,此时j T 2均确定。
求解j f 的过程转化为曲线拟合问题,且需保证0>j f ,建立以曲线拟合残差为目标函数的非负优化模型,采用非负最小二乘(NNLS )算法,利用matlab7.0求解出2T 谱和其反演曲线与原曲线对比图(见图4.1至4.11)。
磁力共振成像中的反演和重建技术研究磁共振成像(MRI)是一项主要用于医学诊断的技术,它通过利用医学梯度磁场、射频波及计算机技术来制造成像。
MRI成像具有非侵入式、无放射性等特点,在目前医学中得到越来越广泛的使用。
其中,磁共振成像中的反演和重建技术,是MRI技术的核心。
一、MRI成像原理MRI成像的原理是基于磁共振现象。
磁共振是指物体在磁场中的微小振荡,其振荡频率与物体本身的化学成分、形态、内部运动状态等相关。
利用磁共振现象,在MRI成像中,通过在人体内部产生一个高强度的医学磁场,再加上一个高频的医学射频波,可以让人体内的氢原子(含于水分子中的质子核)发生共振,进而产生信号。
这些信号被检测出来,并进行处理,最后呈现出一个人体内部的影像,来反映不同组织在磁场下的磁共振现象。
二、MRI成像中的反演和重建技术在MRI成像中,由于检测系统的限制,信号的采集是有限制的。
因此,在信号采集之后,进一步对数据进行反演和重建技术处理是很重要的。
通常,MRI成像中的反演和重建技术可以分为以下几个步骤:1、信号反演在信号反演中,需要利用前期采集的数据,通过数学模型来推算出人体内部的磁共振信号强度分布,这称为反演问题。
反演问题的数学模型通常是一组偏微分方程,例如亥姆霍兹方程、泊松方程和拉普拉斯方程等。
不同的场景和物理过程需要使用不同的方程组,以便更准确地描述磁共振信号的产生和传播过程。
2、图像重建反演问题求解出信号分布后,就可以进行图像重建。
图像重建是指将反演求解的磁共振信号强度分布,通过某种算法,产生一个人体内部的影像。
常用的图像重建算法包括快速傅里叶变换算法、Hough变换算法、哈尔变换算法等等。
这些算法都有其各自的优缺点和适用范围,需要在具体场景下进行选择。
3、图像后处理经过图像重建后,我们得到了一幅初步的影像。
但是这幅影像常常会受到许多因素的影响,例如机器噪音、磁场不均、各种电磁干扰等。
因此,需要进行一些去噪、增强和校正等后处理操作,以便得到更加准确和清晰的影像。
利用改进截断奇异值分解法反演核磁共振弛豫时间高阳;肖立志【摘要】截断奇异值分解(TSVD)法已应用于核磁共振弛豫时间T2谱反演,而奇异值截止值的选取和非负约束的实现是其关键.本文提出一种新的确定奇异值截止值方法,同时利用联合迭代反演法(SIRT)实现非负约束,使本文的改进TSVD法适用于低信噪比资料;该方法具有只需进行一次奇异值分解、计算速度快、保持T2谱分布连续等优点.数值模拟结果表明,改进TSVD法对信噪比仅为5的测量数据依然可得到优于其他方法的较好反演结果;实际岩心资料反演结果进一步证明该方法适用于页岩核磁共振岩心分析和核磁共振测井.【期刊名称】《石油地球物理勘探》【年(卷),期】2015(050)002【总页数】6页(P376-381)【关键词】核磁共振;截断奇异值分解;弛豫时间;多指数反演【作者】高阳;肖立志【作者单位】中国石油大学(北京)理学院,北京102249;中国石油大学(北京)地球物理与信息工程学院,北京102249【正文语种】中文【中图分类】P6311 引言储层岩石核磁共振(NMR)横向弛豫时间(T2谱)包含大量信息,在岩心分析和测井解释中具有重要作用。
因此T2谱反演结果的好坏直接影响到岩石储集物性及流体特性等参数的计算,如岩石孔隙度和孔隙尺寸分布、渗透率、饱和度、流体扩散系数和黏度等[1,2]。
奇异值分解(SVD)法是一种经典的核磁共振多指数反演方法[3]。
王才志等[4]通过对保留的奇异值加上正则化方式拓展了SVD 法;肖立志等[5]和张恒荣[6]则利用 Backus-Gilbert(BG)理论[7,8]确定奇异值截止值。
这些方法均从不同角度改进了SVD法,使其更适用于核磁共振横向弛豫时间的多指数反演,然而这些反演方法对低信噪比数据反演效果并不理想。
本文综合BG折中准则和截止值最小原则,提出一种确定截止值的新方法。
当噪声水平较低时,此方法与BG折中准则一致,保持了在高信噪比条件下TSVD法的优越性;当噪声水平较高时,据此方法确定的截止值比据BG折中准则确定的值小,但尚能较好地保持信号信息,保证在低信噪比条件下TSVD法也能正确地反演T2谱,特别是短弛豫部分与传统TSVD法相比有明显的优越性。
文章编号:1008-0058(2001)04-0395-03核磁测井数据弛豫谱反演软件研究翁爱华,李舟波,王雪秋(吉林大学勘查技术学院,吉林长春 130026)摘要:简单介绍了核磁测井数据弛豫谱反演方法原理,用混合编程方法和动态连接库技术,将用FOR TRAN 语言编写的核心处理模块转换为动态连接库,并在C ++语言中动态调用。
应用快速应用程序开发工具完成了图形用户界面的核磁测井数据处理软件设计。
关键词:核磁测井;弛豫时间谱;混合编程;软件设计中图分类号:P631.8 文献标识码:A收稿日期:2001-08-02基金项目:国家自然科学基金项目(49874028)资助作者简介:翁爱华(1969-),男,安徽省天长市人,博士,主要从事电(磁)法勘探与核磁测井教学与科研工作1 目前国内使用的核磁测井数据处理软件都是随国外仪器绑定销售的,且只作数据反演处理和通用的较为简单的储层评价,很难在其中嵌入地区适应新的储层解释评价方法。
因此,开发独立于仪器的核磁测井数据反演与处理分析软件很有必要。
目前软件开发的语言有多种。
但工程计算上,FO R TRAN 语言具有不可替代的优势,而在进行图形用户界面程序开发时,采用C 语言,尤其是最新发展起来的C ++语言非常方便。
结合两者的优点是工程处理软件最有效的设计与开发途径之一[1]。
FO R TRAN 程序有多种被C 语言调用的方式,如OBJ 文件嵌入,动态连接库(DLL )动态连接等。
由于DLL 具有如下优点:①DLL 独立于主调用程序,可随时更新维护,②占用系统资源小,③程序模块可在进程间共享,④语言独立,以及⑤可随时装载与卸载连接库,故在程序设计中得到广泛应用。
本文的FO R TRAN 语言核心处理模块将被编译成DLL 形式,并由C ++图形用户程序动态调用。
整个软件设计工作关键有三点:一是核心数据处理模块,即T 2-谱反演计算程序的实现;二是动态连接库的生成与调用;三是用快速应用程序开发工具制作图形用户界面,重点是核磁测井自由感应衰减信号(F ID )与弛豫时间谱的图形显示技术。
一种提高核磁共振快弛豫信号测量精度的方法
马莹颖;张宫;廖中枢
【期刊名称】《波谱学杂志》
【年(卷),期】2024(41)2
【摘要】核磁共振技术能够直接探测样品中水、烃等含氢流体信号并对其含量及所处环境进行定量表征,广泛应用在油气勘探领域.在对页岩油气、致密油气等含有大量超快弛豫组分样品进行测量时,由于CPMG脉冲序列观测到的前端信号数据量过少导致测量信号不稳定,快弛豫组分结果可靠性较低.提出了多次微变回波间隔测量(MCTEM)的方法,在测量阶段加密了对核磁共振前端信号数据的采集,反演阶段采用了多次测量回波数据的联合反演,可以得到更加精确的快弛豫T2谱.通过数值模拟和岩心物理实验对MCTEM方法进行了验证,结果显示MCTEM方法能够改善由回波间隔TE偏高引起的T2谱快弛豫组分峰消失的问题,可以有效提高快弛豫信号的测量精度.
【总页数】11页(P162-172)
【作者】马莹颖;张宫;廖中枢
【作者单位】长江大学油气资源与勘探技术教育部重点实验室
【正文语种】中文
【中图分类】TE112.6
【相关文献】
1.利用Phillips-Twomey方法对核磁共振弛豫测量结果的反演及其效果分析
2.核磁共振弛豫信号多指数反演新方法
3.稳态核磁共振实验中弛豫时间T2的测量方法
4.Books
5.一种面向核磁共振短弛豫弱信号的反演方法
因版权原因,仅展示原文概要,查看原文内容请购买。
分层多指数磁共振弛豫信号反演方法研究林婷婷;慧芳;蒋川东;林君【摘要】磁共振测深技术传统反演方法包括平滑反演和分块反演,通过分别获取初始振幅和平均弛豫时间构建地层含水量及有效孔隙度.然而,这些方法局限于单指数拟合方式,损失了大部分有效采集信息,受限于多孔地质环境解释,并在某些情况下无法刻画含水层清晰分界面.针对上述问题,本文建立了基于MRS全数据的多指数反演方法,依据全部采集时间下的有效信息,通过弛豫时间e指数分解,推导出新的磁共振正演核函数,结合泛函极小值方程,直接反演建立含水量,弛豫时间及地层深度三个重要参数关系,适用于复杂地质环境解释.为得到快速稳定的反演结果和更清晰的含水层分界面,本文借鉴分块反演思想,进一步构建了新的反演目标函数,利用基于不等式约束的空间信赖域算法进行优化,最终实现了一种基于分层反演与多指数结合的磁共振弛豫信号反演方法.模型数据以及实测算例表明该方法的效果和优势,并具备较高的计算效率,本研究为地面磁共振反演提供了一种新的思路与方法.【期刊名称】《地球物理学报》【年(卷),期】2013(056)008【总页数】13页(P2849-2861)【关键词】磁共振测深;弛豫数据;多指数拟合;分层反演【作者】林婷婷;慧芳;蒋川东;林君【作者单位】吉林大学地球信息探测仪器教育部重点实验室/仪器科学与电气工程学院,长春 130021;吉林大学地球信息探测仪器教育部重点实验室/仪器科学与电气工程学院,长春 130021;吉林大学地球信息探测仪器教育部重点实验室/仪器科学与电气工程学院,长春 130021;吉林大学地球信息探测仪器教育部重点实验室/仪器科学与电气工程学院,长春 130021【正文语种】中文【中图分类】P6311 引言磁共振测深(Magnetic Resonance Sounding,简称MRS)是一种应用核磁共振现象探测地下水的地球物理方法[1-2],与其他物探方法相比具有非侵入性和能够直接定量探测地下水的特点.在磁共振测深方法[3-4]中,地下水中氢质子在地磁场的环境下处于平衡状态,在地面铺设线圈中通入当地拉莫尔频率的交变电流产生激发磁场,氢质子受磁场激发产生能级跃迁,激发停止后在线圈中接收到氢质子由激发态恢复到平衡态过程中产生的MRS信号.磁共振测深的反演是通过MRS信号初始振幅和平均弛豫时间来反映地下含水量的大小及含水层平均有效孔隙大小[5].Legchenko[6]于1998年提出利用MRS 信号不同脉冲矩下的初始振幅进行反演,仅能得到随深度变化的含水量结果.随着技术的发展,Mohnke[7]在2002年提出了平滑反演和分块反演方法,能快速得到地下含水量和平均弛豫时间.基于此,人们开展了大量的线性反演和非线性反演算法研究.翁爱华等[8]于2005年采用非线性反演技术得到层状大地的MRS 含水量垂向分布.戴苗[9]于2009年采用奇异值分解和模拟退火方法实现了MRS平滑反演.随着研究的不断深入和应用的日益广泛,这两种传统反演方法固有的缺点也逐渐显现,并成为制约MRS技术进一步发展的瓶颈.首先,平滑和分块均需要通过获取初始振幅信息进行反演,由于仪器死区时间的客观存在,初始振幅需通过单一指数拟合后延拓得到,含水量反演结果受噪声干扰及拟合程度影响严重;其次,反演得到的平均弛豫时间仅反映含水层平均孔隙度情况,不利于多孔隙地质环境的解释.地面MRS多弛豫时间反演起步于2005年.Mohnke[10]提出了时间步长反演方法,利用信号随时间变化的测深曲线,得到不同时间下含水量与深度的关系,实现了多弛豫时间反演.但该方法过度依赖于含水量时间曲线,结果极不稳定.Mike等[11]考虑反演模型垂向平滑问题,研究了一种稳定的QT算法,不过由于平滑约束的固有特点,出现含水层边界面刻画模糊,反演结果不清晰、迭代耗时长等缺陷.本文针对目前国内外磁共振反演方法的不足,结合MRS信号整体性与多指数特性,借鉴分层反演的思想,首次提出并实现了磁共振全弛豫数据分层多指数反演方法.采用多指数拟合方式将MRS信号所有采集数据一次性带入反演,直接建立起地下探测深度、弛豫时间、含水量三个重要参数的关系.该方法通过构建泛函极小值方程,利用不等式约束的空间信赖域算法进行优化,可快速、精确的反演获得各含水层位置、弛豫时间及其对应含水量等信息.为全面验证本文提出的反演方法,文中首先建立了多含水层单一孔隙和多含水层多孔隙两种不同的地下水模型.反演分别采用平滑、分块和分层多指数三种方法依次进行,并从模型逼近程度,目标值拟合误差等方面展开对比分析.最后,笔者对白城四方坨子测试点的磁共振野外实测数据进行反演解释,估算了涌水量[12-13]的大小,与钻井资料进行一致性比对.本研究通过模型数据与实测数据,充分证实了该反演方法的优越性.2 磁共振数据反演方法2.1 磁共振信号初始振幅反演磁共振是具有非零磁矩的原子核在交变磁场的作用下发生共振吸收的现象.计算地下水中氢原子核的磁共振角频率ω0公式为其中,B0是测试地点的静地磁场;γp=0.267518 Hz/nT是氢质子的旋磁比;f0是当地氢质子的拉莫尔频率[14].为测量MRS信号,在测试地点铺设线圈,线圈的直径大小取决于所要探测的最大深度.线圈中通入当地拉莫尔频率的交变电流,在地下形成激发磁场.公式为其中,I0为电流强度,激发持续时间0<t<τ,由脉冲矩q=I0τ体现激发强度.激发停止后,地面接收线圈感应到磁共振响应信号,方程为其中,ωr和φr分别代表接收信号的角频率和相位;平均弛豫时间T*2反映地下平均有效孔隙大小;信号的初始振幅E0表示在激发停止时刻磁共振信号的幅值大小,它与地下含水量直接相关.磁共振信号的包络称为自由感应衰减(Free Induction Decay,简称FID)信号,呈指数规律衰减.Weichman[15]等人提出通过方程(4)描述FID信号,建立信号初始振幅与含水量的关系:其中,L代表了地下水探测的深度范围;K1 D(q,z)是一维核函数[16],描述不同脉冲矩下不同深度处核磁共振信号的响应,反映了激发磁场的穿透能力及分辨率,体现地下磁共振响应的灵敏度,如图1所示,横轴代表激发脉冲矩,纵轴代表地下深度,颜色代表了响应强度;M0是氢原子核的磁化强度,代表平衡条件下单位体积内的磁矩;B⊥是激发磁场垂直于地磁场的分量;0≤m(z)≤1代表地下含水量情况.按深度z将方程(4)离散化,得到了初始振幅与含水量之间的矩阵关系为平滑反演方法是将地下从地表至最大探测深度通过等间距或递增间距的方式使地下含水层呈层状分布,公式为采用吉洪诺夫正则化[17-18]方法计算每个含水层对应的含水量m(zi)及平均弛豫时间(zi),i=1,…,b.其吉洪诺夫泛函极小值方程为其中,是经过仪器死区时间td后的信号幅值;α>0为正则化参数,Ω(z)是正则化因子.平滑反演通过单一指数拟合得到不同激发脉冲矩下的信号幅值E1(td)及,根据公式(9)延拓计算初始振幅E0,求解方程(7)和方程(8),获得地下总含水量及平均弛豫时间的垂向分布.公式为分块反演的整个过程与平滑反演类似,不同之处在平滑反演预设地下含水层按深度z呈层状分布,而分块反演是已知地下含水层个数,各含水层深度未知,需要反演获得,其结果使得各含水层之间分界明显,易于判断含水层位置.图1 传统一维核函数图Fig.1 Traditional one dimensional kernel matrix2.2 磁共振弛豫信号的多指数描述目前,磁共振平滑反演与分块反演方法通过平均弛豫时间仅将地下环境假想为单一孔隙.然而,核磁共振测井及实验室核磁共振[19]技术指出,每一地层中探测到的MRS信号可视为不同孔隙产生的弛豫信号叠加.根据磁共振表面弛豫扩散理论,孔隙内流体的弛豫时间与孔隙大小以及形状密切相关,较小孔隙对应信号的弛豫时间短,反之,较大孔隙对应信号的弛豫时间长[20].由上述分析可知,地面磁共振信号应具有多指数特性,而多指数FID信号可通过方程(10)表述:其中,T2i(i=1,2,…n)表示地下环境中第i类孔隙物质对应的弛豫时间,Ei 表示第i类弛豫单元对总MRS信号的贡献.分析多指数方式描述磁共振信号对获取信号初始振幅的意义.图2直观给出单指数方式和多指数方式针对含有噪声的FID仿真信号的拟合情况,黑色点即以信噪比以10dB为例的FID信号;空心圆圈为初始振幅值;黑色及灰色实线分别为采用单指数拟合及多指数拟合FID信号的结果.由图可见,单指数拟合后得到的信号初始振幅低于模型设定值;而多指数方式获得的初始振幅与模型设定值吻合.从趋势而言,采用多指数拟合能够更好的逼近FID信号整个衰减过程.图2 磁共振FID仿真信号的单指数和多指数拟合对比Fig.2 The comparison of mono-exponential and multi-exponential fitting of FID model data表1分析了不同信噪比下分别采用单一指数拟合和多指数拟合的初始振幅拟合误差.通过对比可以发现,在不同信噪比下,多指数拟合信号初始振幅的准确度均高于单一指数拟合,尤其在低信噪比(1~2dB)的情况下仍能比较理想的进行初始振幅拟合.初始振幅的大小直接决定了反演的准确性,因此,采用多指数方式描述磁共振信号对反演解释将具有重要意义.2.3 分层多指数反演方法由前所述,平滑反演和分块反演方法均通过初始振幅-脉冲矩曲线实现地下水垂向分布反演,忽略了大部分弛豫信号,导致各脉冲矩下信息独立,造成较大反演偏差.考虑到核函数的积分特性,即信号在相邻时间相邻脉冲矩下的信息并非完全独立,本文由此提出了分层模式与多指数结合的磁共振全弛豫信号反演方法,能够直接建立地下含水层位置、各层含水量大小及其映射弛豫时间的关系.表1 不同信噪比下初始振幅拟合误差Table 1 Initial amplitude fitting error in different SNR信噪比/dB 30 20 10 5 2 1初始振幅拟合误差/%单指数方式33.83 34.46 39.72 41.64 41.78 45.77多指数方式1.69 4.15 3.39 14.24 12.9 24.182.3.1 弛豫信号多指数的正演推导为实现弛豫信号的反演方法,需要重新推导磁共振弛豫信号的正演方程,建立弛豫信号与地下含水量之间的关系.公式为其中,mi(T2i,z)>0,i=1,2,…,n代表地层z处不同弛豫时间 T2i对应的含水量大小.K1 D(q,z,t,T2i)是本文建立的弛豫信号的一维核函数矩阵,宏观上说,新核函数通过不同弛豫时间T2的e指数因子刻画了在多指数前提下脉冲矩随着时间变化的磁共振响应.根据弛豫时间与含水地层对应关系[21],确定30ms<T2<1000ms,在此范围内选取合适的T2谱分布作为先验信息.为使该正演过程更加清晰,笔者在图3绘制了新的核函数,横轴代表弛豫时间T2,纵轴代表信号时间t.m行n列的子核函数矩阵能够精细描述出不同信号采集时间下、不同弛豫时间T2对应的随深度及脉冲矩变化的MRS响应.图3中可见,随着时间增加,同一T2下对MRS的响应逐渐衰减,小T2对应的响应衰减快,大T2对应的响应衰减慢,这与本文2.2中的磁共振表面弛豫扩散理论吻合.将方程(11)离散化得到多指数下磁共振FID信号与含水量的关系,用矩阵形式表示为2.3.2 分层控制下的目标函数磁共振分层多指数反演之目的是要获得从地表至最大探测深度的地下水分布情况及赋存状态,即求解方程(12)中的mT2,z,对该矩阵的求解过程也是决定反演收敛速度的关键.由于采用全弛豫信号进行反演,涉及的数据量庞大,同时为避免光滑连续模型中出现的迭代缓慢及水层分界模糊现象,文中采取了分层控制方式突出含水边界,经有限步的计算就可以得到精确解.将待反演的含水层个数设为控制量,分层数为k.若地下实际存在的含水层个数为θ,则地层分层数至少为2θ+1.仅当含水层个数小于激发脉冲矩个数时能获取较好的反演效果,固而确定k的取值范围为2θ+1≤k<g,g代表激发脉冲矩个数.对含水层层顶深度zt和层底深度zb求解,就能够确定含水层的位置及厚度,它们是构成本文分层反演的核心,方程(13)给出了分层多指数反演中含水量与深度及弛豫时间的关系:在方程(12)中,反演参数表现为,这是n×k维列向量,磁共振全弛豫信号数据是p×g维列向量,核函数是(n×k)×(p×g)矩阵.由于方程(12)是不适定的,需通过其近似解使函数f)达到极小值,从而得到含水量结果.公式为2.3.3 空间信赖域反演人们求解含水量反演这类非线性规划问题的方法主要有拟牛顿法和共轭梯度法等.但这些方法过于依赖初值选取,有时还会限于局部极值,无法找到全局极值.鉴于传统方法的缺点,我们引入空间信赖域算法[21-22],这是现代应用较广泛的地球物理反演方法,基本思想在于:通过一系列信赖域子问题的最优值逼近最优化问题的解.该方法稳定的数值性能,适合于求解病态最优化问题,并具备全局收敛性.对于无约束或等式约束的优化问题,该方法已被证实是可靠而强适的.笔者通过构造信赖域子问题,得到不等式约束优化问题的信赖域算法.建立目标函数的最小二乘形式,方程为其中,,首先需给定初始信赖域半径Δ0,它刻画了多数情况下相信该信赖域模型的广义球,解被限定于球体内,一定程度上增加了计算的稳定性.每次迭代的信赖域步ξ可由求解子问题(16)得到:其中,grad(J)d 和 Hess(J)d 分别代表泛函J 在迭代点~md的梯度和Hessian阵,分别记为其中,是的一阶导数;是F′)的转置;是F)二阶导数的转置.由于Hessian阵的第二项很难精确求得,且计算量庞大,固而将其进行近似处理.其次,确定信赖域步ξd,再根据比值rd(方程(18))及算法设置条件确定是否接受并调节信赖域半径,得到下一步的信赖域半径Δd+1.接受后,获得,带入方程(15)计算函数值,若满足精度要求则输出结果,否则在新的信赖域半径下计算ξd+1,进行下一次迭代直到满足条件为止.(1)不等式约束的施加为了减小解的范围,根据反演变量的物理意义对其进行不等式约束.zt与zb代表地下深度,定义地表深度为零,则有zt≤0,zb≤0.含水量指单位体积内含水体积占总体积的百分比,确定一个宽泛的取值范围,则有第j层含水层弛豫时间为T2i对应的含水量m(T2i,zj)>0,(i=1,…,n;j=1,…,k),第j层含水层总含水量此外,也可以通过钻井、地质资料等其他方式获得更准确的含水量变化范围(如m(T2i,zj)<0.3).上述变量的取值范围共同构成了解空间M,则建立起不等式约束下的目标函数为将不等式约束施加于反演方程能够有效的约束信赖域空间,提高反演效率. (2)正则化参数的确定图3 分层多指数一维核函数图Fig.3 Layered multi-exponential one dimensional kernel matrix其中,I为单位矩阵.为保证算法顺利进行,调节正则化参数α使Hess(J)d+αI是正定的,则ξ唯一的由方程(21)决定,有然后判断ξd是否落在信赖域半径Δd内,满足‖ξd‖≤Δd则接受该点,并转入下一步迭代,否则按αd+1=αd·γ(γ>1)调节参数,使每次迭代的目标函数为一单调下降序列,直至(15)式迭代收敛到一个稳定数值解.图4为分层多指数反演流程图.图5是本文归纳的平滑反演、分块反演和分层多指数反演的总体流程框图.3 模型对比3.1 模型实验一构造单一弛豫时间T2多含水层模型,分别设置含水层弛豫时间及含水量,比对平滑、分块、分层多指数三种方法反演效果.模型数据是在如下条件下通过正演计算获得:(1)边长100m的单匝方形线圈,采用收发共线方式;(2)假设地下为均匀半空间,电阻率为500μm.地磁场强度B0=46970nT,地磁倾角60°;图6 模型一Fig.6 Model one(3)拉莫尔频率为2000Hz,激发脉冲矩从0.4As到10As分为20次激发地下0m到100m;(4)假设地下存在三个含水层(图6),第一层设置于地下浅层5m至15m,含水量为10%,弛豫时间T2为100ms;第二层设置于地下中层25m到35m,含水量为10%,弛豫时间T2为300ms;第三层设置于地下深层55m到75m,含水量为10%,弛豫时间T2为200ms;(5)由于磁共振方法采集信号水平在纳伏级,易受到噪声的干扰,导致反演精度下降.因此,为使仿真结果接近实际并验证该方法的适用性,在正演模型中加入5%Guass白噪声.设置反演地层个数k=10;初始正则化参数α=100.图7给出了模型一的反演结果.(a)为含噪磁共振数据,蓝色实线为不同脉冲矩下的信号衰减曲线,红色圈为不同脉冲矩下接收时刻信号的幅值E1(td).图(b)、(c)中左图分别为平滑反演和分块反演得到的含水量垂向分布,右图为平均弛豫时间垂向分布,红色虚线为模型值,蓝色实线为反演值.(d)中左图为分层多指数反演得到的地下含水量分布,横轴代表弛豫时间T2,纵轴代表深度,颜色代表含水量;右图为总含水量结果,红色虚线为模型值,蓝色实线为反演值.表2列出了采用平滑、分块以及分层多指数三种反演方法处理模型一数据的反演结果,对比了反演含水层的层顶深度、层底深度、含水量及弛豫时间,并对总含水量的拟合情况进行了误差分析.图7 模型一数据反演结果(a)模型一数据;(b)平滑反演结果;(c)分块反演结果;(d)分层多指数反演结果.Fig.7 Inversion results of model one data(a)Data set of model one;(b)Result of smooth inversion;(c)Result of block inversion;(d)Result of layered multi-exponential inversion.通过图7及表2可以得出结论:平滑反演和分块反演能够得到地下含水量及平均弛豫时间的垂向分布,但由于方法的固有缺陷,无法定位深层含水量的准确位置,总含水量拟合数值偏小,误差较大.分层多指数反演能够得到地下随深度及弛豫时间T2的含水量分布情况.同时它精确区分了浅层、中层、深层的含水构造,与模型设置基本一致,总含水量反演结果准确,拟合误差仅为2.59%.3.2 模型实验二构造多弛豫时间T2含水层模型.模型数据是在如下条件下通过正演计算获得:表2 模型一数据反演结果Table 2 Inversion results of model one层顶深度/m/%模型一-5 -15 100 10 10-25 -35 300 10 10-55 -75 200 10 10平滑反演-3.5-15.5 185(T*2)8.26-21.5-37.5 389(T*2)11.29-43.5-100 188.5(T*2)6.852.81-15.5 211.2(T*2)5.84-23.5-33.5 410.7(T*2)10.64-49.5-100 200(T*2)4.780分块反演4.16分层多指数反演-4.5 -14.5 100 9.91 10.36-23.5 -34.5 300 7.91 9.8-53.5 -78.5 200 6.797.812.59(1)边长为100m的单匝方形线圈,采用收发共线方式;(2)假设地下为均匀半空间,电阻率为500μm.地磁场强度B0=46970nT,地磁倾角为60°;(3)拉莫尔频率为2000Hz,激发脉冲矩从0.4As到10As分为20次激发地下0m到100m;(4)假设地下存在两层含水层(图8),第一层设置于地下浅层20m至40m,含有两种弛豫时间成分,T21为100ms,对应含水量为5%,T22为400ms,对应含水量为5%;第二层设置于地下深层60m至80m,含有两种弛豫时间成分,T21为200ms,对应含水量为5%,T22为600ms,对应含水量为5%;(5)模型中加入5%Guass白噪声.设置反演地层个数k=8,初始正则化参数α=100.图8 模型二Fig.8 Model two图9 模型二数据反演结果(a)模型二数据;(b)平滑反演结果;(c)分块反演结果;(d)分层多指数反演结果.Fig.9 Inversion results of model two data(a)Data set of model two;(b)Result of smooth inversion;(c)Result of block inversion;(d)Result of layered multi-exponential inversion.表3 模型二数据反演结果Table 3 Inversion results of model two层顶深度/m/%-20 -40100 5400 510模型二-60 -80200 5600 510平滑反演分块反演-18.5-50.5 362(T*2)8.62-58.5-100 455(T*2)6.45-19.5-38.5 364(T*2)8.73-60.5-100 491(T*2)5.324.793.59-20.5 -39.51004.98400 4.4110.5分层多指数反演2.86-63.5 -84.5200 2.13650 1.798.24图9和表3给出了三种方法处理模型二数据的反演结果.结论如下:平滑反演和分块反演无法准确获得深部含水层位置信息,并且仅能获得平均弛豫时间,无法区分同一含水层的不同弛豫时间分布,含水量拟合误差分别达到4.79%及3.59%.而分层多指数反演能精确实现浅层、深层含水定位,与模型基本相符.同时,在浅层和深层区分出给出不同弛豫时间T2及其映射的含水量分布情况,含水量拟合误差仅为2.86%.表3中给出最大局部含水量所对应的弛豫时间T2.此外,图9(d)中右图所示的总含水量则代表每一含水层多个弛豫时间所映射的含水量之和,与模型设置吻合.需要指出的是,随着探测深度的增加,磁共振响应的分辨率逐渐降低,T2反演精度下降,向相邻的T2值均匀扩散,这是导致深层两个弛豫时间映射含水量偏小的本质原因.但从总含水量的数值而言,分层多指数反演更接近于模型,显著优于传统反演方法.4 野外数据实例经模型数据验证分层多指数方法后,本文对野外实测MRS数据进行反演处理,证明其适用性.研究组利用吉林大学自主研发的JL-MRS磁共振地下水探测仪于白城市镇赉县四方坨子嫩江岸边进行了150m地下水探测试验.测试点地磁场55936nT,采用边长150m方形线圈,激发频率为2383Hz.钻探结果表明含水层见于井深3m到21m、63m到106m和126m到148m处,含水层厚度共83m,为孔隙承压水,静止水位4m,当水位下降至7m时的涌水量为320t/h.图10a为实测数据,平均信噪比达34.68dB(相当于信号中仅含有1.85%的噪声,在核磁共振反演中,信噪比高于6.02dB数据即为真实可信[23]).图10(b—d)是分别采用平滑、分块和分层多指数方法对数据进行反演的结果,图10e为钻井成图.表4给出了白城测点钻井地层结构与三种反演方法的结果对比情况.平滑反演得到地下存在两层含水层,浅层位于地下5.25m至38.25m,深层存在含水层但无法确定其具体位置,其平均弛豫时间T*2分布平缓过渡,无法与含水层准确对应.分块反演显示地下主要存在两层含水层,浅层位于地下6.75m至29.25m,同样无法准确体现深层含水层的位置,与平滑反演相较,T*2有明确划分,然而深层含水层仅有较大弛豫时间247.2ms的反演解释结果与野外实际地层情况有一定偏差.上述两种方法反演的含水层均为连续变化,未见相对隔水层出现.分层多指数方法的初始正则化参数α=100,反演设置地层个数k=10(参数优化后,E0-q曲线拟合误差最小).结果显示地下主要存在两个含水层,浅层位于地下3.75m至33.75m,从弛豫时间上看地表土对应的T2为53ms;隔水构造出现在33.75m与41.25m之间,该隔水构造对应右图所示的较小的含水量分布;中间层亚粘土地层位于地下41.25m至65.25m,对应的弛豫时间T2为85ms比理论值略大,因实际探测中,亚粘土属于束缚水构造(一般认为磁共振探测为自由水),反演拟合易造成偏差,而成为上下两个含水层弛豫时间值的过渡.深层含水层层顶位置反演准确,从弛豫时间上看,地下63.75m至150m主要存在两个弛豫时间成分,分别是35ms和196ms,能够对应于钻井结果63m到148m的岩砾石结构与亚粘土.随着深度的增加,激发场逐渐减弱,导致了磁共振响应分辨率随着深度的增加而逐渐降低,使得分层多指数方法亦不能明确的给出第三层出水层,但其弛豫时间的分布差异提示在主要含水层中或有隔水介质出现.图10 白城测点数据反演结果(a)白城测点数据;(b)平滑反演结果;(c)分。
核磁共振弛豫信号多指数反演新方法
陈华;潘克家;谭永基
【期刊名称】《测井技术》
【年(卷),期】2009(033)001
【摘要】孔隙介质核磁共振(NMR)弛豫信号的多指数反演在NMR岩心分析与测井解释中起着关键作用.通过将核磁共振弛豫信号多指数反演问题转化为一个带非负约束的非线性优化问题,先用非负最小二乘方法(NNLS)确定弛豫时间组分数,然后由差分进化(DE)算法求解优化问题.新方法计算非常简单,不依赖于初值,而且容易实现非负约束.计算机模拟和实际数据反演结果表明,该方法对核磁共振弛豫信号分析非常有效,尤其是在信噪比较低的情况下,仍能保证弛豫谱的真实性.
【总页数】5页(P37-41)
【作者】陈华;潘克家;谭永基
【作者单位】中国石油大学,山东,东营,257061;复旦大学数学科学学院,上
海,200433;复旦大学数学科学学院,上海,200433
【正文语种】中文
【中图分类】P631.84
【相关文献】
1.核磁共振弛豫信号的多指数反演及应用 [J], 苏俊磊;孙建孟;张守伟
2.利用改进截断奇异值分解法反演核磁共振弛豫时间 [J], 高阳;肖立志
3.低场核磁共振扩散-横向弛豫二维反演算法 [J], 李新军;聂生东;王远军;杨培强
4.分层多指数磁共振弛豫信号反演方法研究 [J], 林婷婷;慧芳;蒋川东;林君
5.孔隙介质核磁共振弛豫测量多指数反演影响因素研究 [J], 廖广志;肖立志;谢然红;付娟娟
因版权原因,仅展示原文概要,查看原文内容请购买。
基于差分进化算法的核磁共振!!谱多指数反演!潘克家!)"陈华#)谭永基!)!)(复旦大学数学科学学院,上海#$$%&&)#)(中国石油大学(华东),东营#’($)!)(#$$(年!!月!’日收到;#$$*年!月!+日收到修改稿)提出了一种基于差分进化(,-)算法的核磁共振弛豫信号多指数反演新方法.将核磁共振!#谱反演问题转化为带非负约束的非线性优化问题,不需要预先给定横向弛豫时间!#分布,直接利用差分进化算法进行反演计算.在测量信号低信噪比的情况下,计算机模拟和实验数据反演都表明了该方法在分析处理/01弛豫信号中的有效性.关键词:核磁共振,多指数反演,差分进化,岩心分析"#$$:())$!国家自然科学基金(批准号:!$%&!$&$)资助的课题."-23456:7894:;<4=>>.?:!@引言在核磁共振测井中,一般采用A0BC 方法测量自旋回波串[!,#],纵向弛豫时间!!和横向弛豫时间!#都是用来描述自旋回波信号的弛豫特征.由于和横向弛豫时间!#相对应的!#谱能够提供许多岩石物性和流体特性的信息,越来越受到人们的关注[&—(].而!#谱的求解本质上是一个多指数反演的问题[*,+],对于多指数反演,近年来先后出现了罚函数法(D1,)、奇异值分解法(EF,)、非负最小二乘(//GE )以及联合迭代重建算法(EH1I ),并且在罚函数法和奇异值分解法的基础上又发展了一些改进算法,这些算法从不同角度来解决多指数反演问题[!$].但这些算法除了受到横向弛豫时间布点数、原始回波采集个数等因素影响外,还易受信噪比的影响,在不做任何参数校正的情况下,这些算法都要求信噪比E/1J %$@对于低信噪比(E/1K &$)情况,上述算法会出现基线偏离的情况,导致孔隙度估计不准[!!],另外,在处理!#谱的非负约束条件时,它们均采用迭代法修正,大大增加了计算时间.差分进化(L5MMNON:P546NQ>6RP5>:,,-)算法[!#,!&]是由EP>ON 和BO5?N 于!++)年共同提出的一种采用浮点矢量编码在连续空间中进行随机搜索的优化算法,由于其原理简单,受控参数少,无需设置初值和进行导数计算,且易于理解和实现,是一种比较有发展前途的算法.核磁共振弛豫信号多指数反演问题,可以转化为一个带非负约束的非线性优化问题,然后利用差分进化算法求解.#@核磁共振!#弛豫信号数学描述根据核磁共振理论,核磁共振信号与横向弛豫时间谱满足第一类SONL=>63积分方程:"(#)T "!#34U !#35:$(!#)NV #!#L !#W !(#),(!)其中横向弛豫时间谱$(!#)为横向弛豫时间!#对应的幅度值,从物理意义上讲,$(!#)为非负值.!#35:与!#34U 是测量的衰减信号所能分辨的最短和最长弛豫时间,一般!#35:T $@!3X ,!#34U T !$$$$3X ,!(#)为随机噪声.把(!)式离散化后,所测的核磁共振信号"(#%)可以认为是由一系列满足单指数衰减规律的回波信号叠加,即"(#%)T #&’T !$’·N V #%!#’W !(#%),!$%$(,(#)其中&,(分别为测量到的回波个数和横向弛豫时间!#分布个数,#%为采集时刻(通常为回波间隔I-的整数倍),$(#%)为#%时刻测得的回波幅度.第’(卷第+期#$$*年+月!$$$2&#+$Y#$$*Y’(($+)Y’+’)2$)物理学报ZAIZ B[\EHAZ EH/HAZF>6.’(,/>.+,EN9PN3]NO ,#$$*!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%#$$*A=5:.B=<X.E>?.!"多指数反演优化模型横向弛豫时间!#分布跨度大,通常范围从$"%&’到%$$$$&’(因仪器不同而有所差异),导致方程(#)的求解是一个典型的不适定问题(测量数据很小的误差也会导致最后的计算结果偏差很大,这也增大了反演过程的难度(现有的各种!#谱反演方法[%)—%*],!#分布都是事先采取某种对数均匀布点得到,然后由方程(#)得到一个线性方程组,再利用+,-,+./0,112+等方法求解(而对"(!#)非负约束的处理,都是采用某种迭代修正的思想(由于!#分布事先给定,导致这些方法都存在一个普遍的问题:就是在测量!#组分离散且分布较宽时,反演得到的!#谱分辨率较低(并且,这些基于迭代处理非负约束的方法,有时收敛得比较慢,从而大大增加了问题的计算量(文献[%3]提出了一种反演与拟合相结合处理核磁共振弛豫数据的方法,求解使用了245467489:;<8=><8?@方法,对初值的依赖性比较强,可能陷入局部极值(为了避免这些问题,方程(#)的求解转化为一个带非负约束的优化问题,同时反演出!#和"(!#),得到横向弛豫时间!#和"(!#)的关系图,即!#谱图(考虑到原问题的病态性质和测量的误差,通过如下两种正则化方法[%A ,%B ],将原问题转化为如下带非负约束的优化问题(范数平滑&C6!#$D %%$E!&’D %"’·4E ($!#()’#F !!&’D %"#’,(!)’(@("’"$,%#’#&(())曲率平滑&C6!#$D %%$E!&’D %"’·4E ($!#()’#F !!&E%’D #"’F%E #"’F "’E ()%#,(G )’(@("’"$,%#’#&((*)其中,!为正则化参数,它的选取与原始数据的误差水平相关(对上述非线性优化问题,采取差分进化算法求解(计算结果表明,基于曲率光滑的正则化优化算法能更好地抑制"(!#)的振荡,对正则化参数的选取不如范数平滑时敏感,对此类多指数反演问题更加有效(和文献[%3]采取范数平滑不同,最后的反演实例都采用曲率平滑,求解优化问题(G ),(*)反演!#谱()"差分进化算法差分进化算法是一种基于实数编码的具有保优思想的全局优化进化算法(在%BB*年举行的第一届国际.HHH 进化优化竞赛上,对提出的各种方法进行了现场验证,-H 被证明是最快的进化算法(与传统的进化算法如遗传算法等相较,除了具有收敛速度快,控制变量少,更易于理解和编程实现等特点之外,其主要区别在于:%)传统的IJ 算法采取二进制编码,而-H 算法采取实数编码;#)在IJ 中,通过两个父代个体交叉变异产生两个子个体,而-H 中,通过三个父代个体交叉变异产生一个子个体;!)传统IJ 产生的子个体通过某种概率选择直接取代父个体,而-H 中新产生的子个体,只有当它比相应的父代个体优良时才替换父个体(!"#"基本差分进化算法为求非线性函数"()%,)#,…,)*)的最小值,进化过程中的第+代利用,-个*维参数向量构成种群{)$,+,$D %,#,…,,-},其中种群大小,-在整个进化过程中保持不变()"%"%"种群初始化在*维空间,随机产生满足约束条件的,-个向量)’$,$D 8<6?’·(.’,/E .’,0)F .’,0构成初始种群(其中,8<6?’为区间($,%)上的随机数,.’,/和.’,0分别为第’个变量)’的上界和下界()"%"#"变异对每一个目标向量)$,+,在{%,#,…,,-}范围内随机选择!个互异整数1%,1#,1!,且使得1’$$(%#’#!),得到变异向量2$,+F%D )1%,+F 3·()1#,+E )1!,+),其中,变异因子3%($,#)为一常数,控制偏差向量的放大程度()"%"!"交叉交叉是为了增加种群的多样性,通过如下方式得到试验向量:4’$,+F%D2’$,+F%,8<6?’#56或’D 1#.1$,)’$,+F%,8<6?’K 56和’$1#.1${,3G B G B 期潘克家等:基于差分进化算法的核磁共振!#谱多指数反演!!",#,…,"$式中%&’(!为(),")之间的随机数,%’*%!!{",#,…,"}为随机选择的序列,用它来确保#$,%+"从&$,%+"至少获取一个分量$交叉因子’(![),"]$,-"-,-选择差分进化的选择方案基于局部竞争机理,如果#$,%+"的目标函数值小于)$,%+"的目标函数值,则令)$,%+"等于#$,%+";反之,令)$,%+"等于)$,%$,-"-.-边界条件处理在有边界约束的问题中,必须确保产生新个体的参数值落在问题的可行域中,一个简单方法是将不符合边界约束的新个体用在可行域中随机产生的参数向量代替$!"#"差分进化算法其他形式上面介绍的是最基本的/0操作程序,实际应用中还发展了/0的多个变形形式,并用符号/01)1 *加以区分,其中)限定当前被变异的向量是“随机的”或“最佳的”;*是所利用的差向量的个数$利用这个表示方法,前面叙述的基本/0策略可描述为/01%&’(1"-其他多种形式,主要体现在变异操作的不同,如:")/01*2341"&$,%+"!)*234,%++·(),",%5),#,%)$#)/01*2341#&$,%+"!)*234,%++·(),",%5),#,%+),6,%5),,,%)$6)/01%&’(1#&$,%+"!),",%++·(),#,%5),6,%+),,,%5),.,%)$,)/01%&’(7487*2341"&$,%+"!)$,%+!·()*234,%5)$,%)++·(),",%5),#,%)$实际反演计算表明,/01%&’(7487*2341"算法既保持了种群的多样性,又具有较快的收敛速度,对此类9:;多指数反演问题最为有效$!"$"%&算法控制参数选择控制参数对于一个全局进化优化算法的影响是比较大的,/0算法的控制参数选择有如下经验规则["#,"6]:")种群数量$根据问题规模,-.的合理选择在."—")"之间$过大的-.导致算法收敛很慢,过小的-.可能导致/0算法过早收敛,得不到目标函数的全局极小值$#)变异因子$在理论上也没有一个最佳的选择准则,比较好的容许区间为[)-,,"],不落在此区间的变异因子只在极少的情况下有效$通常取初始值+!)-.,若种群过早收敛,则增大+或-.$ 6)交叉因子$通常/0算法如果收敛,则交叉因子’(越大,收敛速度越快$故为得到较快的收敛速度,一般先选取’(!"-)或者’(!)-<$ ,)加权因子$加权因子!!(),"),通常取作)-=$.-反演数值结果这一节,采用差分进化算法反演核磁共振/#谱分布$第一小节中用理论模型验证算法的有效性,并研究了信噪比>9;和正则化参数"对反演结果的影响$通过计算实例和对比说明了/0反演算法抗噪能力强,收敛速度比较快$第二小节给出了运用差分进化算法反演/#谱的一个实例$’"("理论模型试算结果与分析在这一节的计算实例中,先给定一个具有双峰特性的/#谱分布,由(#)式计算出不同时刻0$的信号强度值*$,再对*加上一定大小>9;的高斯白噪音得到"*,然后由数据(0$,"*$)用/0算法反演/#分布$信噪比>9;按如下方式定义:>9;!#)·?8@")#*###*5"*##$当/#组分1!")时,构造了."#个测量数据进行反演;当/#组分1!#)时,我们构造了")#,个测量数据进行反演$从表"可以看出,对于反演#)个未知参数的情形,不加噪声时反演结果非常精确,相对误差2#范数为)-)6A;即使加入信噪比为=)的噪声后,反演结果也还很精确,相对误差仅为"-=.A$这说明了/0反演算法对于比较小的1,高信噪比数据反演非常有效,精确$下面研究1!#)(反演,)个参数)时,正则化参数对反演结果的影响$反演结果如图"所示:当"取")5=时计算结果相当好;但不加正则化项(1取作))时,计算出来的3值振荡得比较厉害;而正则化参数"过大取")5B时又会把解磨得太光,甚至失去函数本该有的双峰特性$所以,选择合适的正则化参数对此类问题的求解起着至关重要的作用$目前还没有非常好的正则化参数选取方法,下一步=.<.物理学报.C卷希望能研究出一种自适应的正则化参数选取算法,来更好地解决此类实际问题!表"!#"$,%&’分别为()和*$时反演结果"+""++"+,"+-"+."+/"+0"+*"+1"+"$23456"7$$$8"$9,"711.8"$9,,71*"8"$9,071-,8"$9,"7.*.8"$9+,7"/+8"$9+/7,"$8"$9+"7+.18"$9"+7."+8"$9".7$"+8"$9"%&’#)"7$$$8"$9,"711.8"$9,,71*"8"$9,071-,8"$9,"7.*.8"$9+,7"/,8"$9+/7,""8"$9+"7+.18"$9"+7.",8"$9".7$"08"$9"%&’#"$$171*08"$9-+7$",8"$9,-7$*,8"$9,*7+,+8"$9,"7/.18"$9+,7,./8"$9+/7.,.8"$9+"7+018"$9"+7.,.8"$9"-711*8"$9"#"#+#,#-#.#/#0#*#1#"$23456.7*"18"$9+07//$8"$9""7.0/+7,0$+7-$0+7$,*,7,-0.7$+*+7,0-+7,*,8"$9"%&’#).7*+"8"$9+07//$8"$9""7.0/+7,0"+7-$0+7$,1,7,-*.7$+1+7,0++7,0,8"$9"%&’#"$$.7***8"$9+071++8"$9""7/-++7-+.+7-+0+7$.+,7,"1-71-/+7,$-+7,/+8"$9"图"正则化参数对反演结果的影响图+信噪比对反演结果的影响图+为在不同信噪比下选择合适的正则化参数进行反演,得到的一组反演结果!从图中可以看出,:2算法对测量数据的信噪比不敏感,即使在噪声比较大,信噪比仅为"$时,仍能反演出比较好的结果!而传统的反演算法,如传统的%;:,&&<%等,由于没有引进正则化项,对测量数据的误差比较敏感!在文献["$]中给出了几种传统算法的信噪比稳定区间,其中%;:算法对信噪比的稳定区间大约为-$=%&’="$$$$,>’:算法的稳定区间大约为,$=%&’,%?’@算法的稳定区间约为-$=%&’,并且该文中也专门针对传统方法对信噪比比较敏感,发展了相应的校正方法提高反演质量!文献["0]中提到&&<%方法得到的"+谱可信度较高,为了与:2反演算法进行比较,下面采用&&<%方法和改进的%;:方法[".]反演!#+$时的理论模型!并且,还假设"+弛豫时间布点已知,在程序中布点数为,$,且其中的+$个点正好与理论模型中给定的数据相同,然后利用&&<%方法和改进%;:方法反演"+谱,计算结果分别如图,和图-所示!从图中可以看出,&&<%算法和%;:算法在没有误差,或者信噪比较高的情况下计算非常准确,但在有测量误差比较大的情况下,反演得到的结果振荡比较厉害!当信噪比为/$时,&&<%方法已经得不到满意的反演结果,当信噪比为,$时,改进的%;:算法得到的"+谱振荡也比较厉害;而基于:2的优化算法由于加入了曲率平滑项,能够很好地抑制"+谱的振荡,即使在信噪比为"$的情况下,也能得到如图+所示的比较可信的"+谱!前面所有反演计算,都是采用:2AB4CDE6FEGHI6A"求解优化问题!:2控制参数选取如下:变异因子$取为$7*,交叉因子%&取为"7$!加权因子!取为$7*,’(取为未知参数个数)的"$倍!图+中在各种信噪比情况下的计算时间和目标函数值如表+所示!表中计算时间为用"$+-个数据反演-$个参数值,JFB6B4C 程序在KHC6LMN (’):,7+O 电脑上的运行时间!从表+中可以看出,尤其在信噪比比较小的情况下,:2算法的收敛得非常快!当%&’#"$1.1.1期潘克家等:基于差分进化算法的核磁共振"+谱多指数反演图!信噪比对""#$方法反演结果的影响图%信噪比对改进$&’方法反演结果的影响时,仅用了()*,就得到了如图+中比较好的反演结果,由此也说明了’-反演算法非常适合低信噪比核磁共振弛豫信号反演,表+’-算法收敛速度$".()+)%)/)())计算时间0*()+(12()%+23目标函数值(4+5(43)6()7((4%)6()7!!45)6()73%4%)6()75!"#"应用实例使用英国.8*9:;:<8=:*>?@A8:>*公司的BC.C"D图3实际数据反演结果比较+#!型核磁共振岩心分析仪,对长庆某井岩心进行岩石核磁共振横向弛豫时间!+谱测量,试验参数设定为:回波间隔时间E-F )45A*,回波采样个数()+%个,经数据压缩后,实际回波采样个数为(+/个,’-反演算法正则化参数!取为54)6()7(),!+布点数"取为+);而G.’和$&’算法!+均采取对数均匀布点,布点数为!)4计算结果如图3中所示,’-算法在不预先给定!+的情况下,得到的结果与已有算法取得的结果非常符合,14结论通过计算机模拟和实际算例,得到如下结论:(4和传统算法如$&’相比,基于’-的优化算法不依赖于初值选取,计算稳定;更容易实现非负约束,计算速度快,+4优化算法能同时反演出!+分布,得到分辨率比较高的!+谱,!4优化算法引入正则化项,抗噪能力比较强,可用于低信噪比!+谱反演,%4优化算法非常适合!+组分离散且分布较宽("H %)),信噪比较低的".B 弛豫信号反演;当!+组分较多甚至连续分布时,可采用$&’,$=.E ,""#$等方法求解,[(]I98JA;:I ,K?;AA8?B I (223#,$%&’,()*,$$%(([+]LM;:N I O ,PQ?;*;RQ I S +))!#,$%&’,()*,$&%335)123物理学报35卷[!]"#$$%&’(,)#*+,"-..-!"#$/%&’/()"#/!0![1]234,53%6789-..-*+,"-./’/01$/"#1:;(<6=><6+?+)[许峰、黄永仁-..-物理学报"#1:;][;]"+67@(,53=":AAB *+,"-./’/01$/$%:AC:(<6=><6+?+)[孟庆安、胡传民:AAB 物理学报$%:AC:][C ]D>+67E F ,=>+6D ,=>+6D G ,D>#67,5-..:2.1$/-./’/#&;;0[B ]H<I 8,232=:AA;*+,"-./’/01$/$$:01B (<6=><6+?+)[李鲠颖、徐学诚:AA;物理学报$$:01B ][0]’#J7<%I =,’J#K69,E ,4%6&%LL<6<):AA03/!"#$/%&’/#’!C;[A ]’#J7<%I =,’J#K69,E ,4%6&%LL<6<)-...3/!"#$/%&’/#$(-B![:.]H<%#I D ,2<%#H D ,2<+95,43,,-..B 2.1$&’&3/4&56./’/"&A!-(<6=><6+?+)[廖广志、肖立志、谢然红、付娟娟-..B 地球物理学报"&A!-][::]G+67(5,H<D ’,"M 2G ,G%672@-..!3/3171$8$19/(N%J&>EO</NP/)’’-!-(<6=><6+?+)[翁爱华、李舟波、莫修文、王学秋-..!吉林大学学报(地球科学版)’’-!-][:-]E&#J69,)J<O+F :AAB 3/475:/;6,1)/##!1:[:!]F%+Q#),(Q<""-..C <=>/3/;6&>/%&’/#%)::BC[:1]G%67D R ,2<%#H D ,H<3S 8-..!0+1&$+&1$2.1$"(E+J/I )’’!-!(<6=><6+?+)[王忠东、肖立志、刘堂宴-..!中国科学I ’’!-!][:;],<%679D ,8%#8),"<%#E ,D>%67=E -..;*+,"-&,>5/01$/!%;B (<6=><6+?+)[姜瑞忠、姚彦平、苗盛、张春生-..;石油学报!%;B ][:C ]G%67G ",H<),8+D 5-..:0+1&$+&1$2.1$"(E+J/()’#B!.(<6=><6+?+)[王为民、李培、叶朝辉-..:中国科学(’#B!.][:B ]G%675,H<I 8-..;*+,"-./’/01$/"$:1!:(<6=><6+?+)[王鹤、李鲠颖-..;物理学报"$:1!:][:0]’3&Q+J ,),9++P?,(,R%K?#6E T :A0:0(*!3/?=)&>/*$"7/#*!0:[:A ]R366F ,,H%S#JJ%O%I (,G%J6+J ,H ,’+J7U%6R ,:AA10-<-0!CB 1;!"#$%+&’()*&*$%+#%*,&-.%)*)/!!.(&0$-"1%*2!34+.&5)*5%//&-&*$%+#&,)#"$%)*+#6)-%$71!)%6F+V,<%:)W =>+653%-)S%68#67V,<:):)(0+.5575@!",.&)",1+"70+1&$+&’,A=B"$8$19&>’1,/,0."$#."1-..1!!,2.1$")-)(2.1$"8$19&>’1,/5@-&,>57&=),C5$#/1$#-;B.C:,2.1$")(9+O+<X+P :;Y#X+UZ+J -..B ;J+X<?+P U%63?OJ<*&J+O+<X+P :A ,%63%J[-..0)(Z?&J%O&(6+K U3Q&<V+\*#6+6&<%Q <6X+J?<#6U+&>#P $#J Y"9J+Q%\%&<#6?<76%Q?<?*J+?+6&+P %6P &+?&+P ,K><O><?Z%?+P #6P<$$+J+6&<%Q +X#Q3&<#6(RN )%Q7#J<&>U/S>+<6X+J?<#6#$D -?*+O&J3U <6Y"9<?J+P3O+P &#%6#6V6+7%&<X+O#6?&J%<6&#*&<U<L%&<#6*J#ZQ+U ,%6P ?#QX+P P<J+O&Q[Z[RN %Q7#J<&>U K<&>#3&*J+V%??<76<67&>+P<?&J<Z3&<#6#$J+Q%\%&<#6&<U+D -/S>+X%Q<P<&[%6P +$$+O&<X+6+??#$#3J U+&>#P <?P+U#6?&J%&+P Z[O#U*3&+J ?<U3Q%&<#6?%6P &>+<6X+J?<#6#$*J%O&<O%Q Y"9P%&%36P+J Q#K EY9/,-./0123:63OQ+%J U%76+&<O J+?#6%6O+,U3Q&<V+\*#6+6&<%Q <6X+J?<#6,P<$$+J+6&<%Q +X#Q3&<#6,O#J+%6%Q[?<?4566:BCC.!)J#]+O&?3**#J&+P Z[&>+Y%&<#6%Q Y%&3J%Q EO<+6O+4#36P%&<#6#$=><6%(IJ%6&Y#/:.1!:.!.)/W NVU%<Q :^]*%6_[%>##/O6:C A ;A 期潘克家等:基于差分进化算法的核磁共振D -谱多指数反演。