一种基于特性参数反馈的微波热消融温度场仿真方法
技术领域
本发明属于热消融温度场仿真技术领域,尤其涉及一种基于特性参数反馈的微波热消融温度场仿真方法。
背景技术
肝肿瘤是人类生命安全的最大威胁之一,在中国是肿瘤导致死亡的第二杀手,因此肝肿瘤的有效治疗已成为亟待解决的社会问题。微波热消融技术因微创、治疗效果显著等优点而在临床中得到了广泛应用,已成为治疗人体肝肿瘤行之有效的方法,但热消融的质量仍主要取决于临床医生的经验,缺乏客观依据。在临床手术中,肝肿瘤的热凝固区常采用54℃作为边界阈值。因此,需要热消融温度场的精确表征以提高热消融手术的科学性。
微波消融的基本原理是将微波天线经皮插入肿瘤组织内,通过微波加热使细胞蛋白凝固甚至变性,从而杀灭肿瘤细胞并且达到治疗目的。目前常用的临床微波热消融装置为恒功率微波热消融仪,例如,915MHz和2450MHz微波消融仪,其利用恒定的功率(例如,40W、50W、60W等)获得所需的热凝固效果。
现有技术中,基于电磁方程和生物传热方程的有限元仿真已成为预测微波热消融温度的有效方法,其中特性参数随温度的变化对仿真结果具有显著影响,但特性参数如何影响温度场模型仍是未知的。另外,在温度场仿真模型的常规研究中,通常无法获得模型参数随温度的实际变化,模型参数一般采用固定值或文献参考值,从而导致仿真温度的精度参差不齐。公开号为CN103800075A的专利公开了一种基于肝脏肿瘤消融术的患者特定建模的方法和系统,其在仿真过程中基于细胞坏死程度来更新至少一个组织参数,但在临床手术过程中对细胞坏死的评价以及如何调节特定参数仍存在问题;公开号为CN102008351A的专利公开了一种热消融温度场分布的获取方法,其利用温度场分布函数表征各点的温度,具有一定的便捷性,但未考虑到消融过程中的参数变化,从而不能精确地获得温度场分布。因此,温度场模型中特性参数的特异性分析和反馈研究仍是温度场仿真领域亟待解决的技术。
发明内容
针对现有技术中因组织特性参数的不确定性而导致无法精确获得温度场分布的问题,本发明提出了一种基于特性参数反馈的微波热消融温度场仿真方法,通过参数敏感性分析和单针实测反馈获得组织特性参数的精确表征形式并进而通过仿真技术获得恒功率微波热消融的精确温度场模型。现有技术中多采用线性形式表征特性参数随温度的变化,因此本发明基于线性函数来反馈组织特性参数。
应用本发明的前提是:通过COMSOL Multiphysics软件(COMSOL Inc.,Palo Alto,CA,USA)建立温度场仿真模型;利用Minitab软件(State College,PA,USA)中的因子分析(DOE)研究各特性参数(比热容(Cp)、对流交换系数(h)、电导率(Sigma)、导热率(K)和相对介电常数(Epsilon))的敏感性,从而为参数的特异性反馈提供依据;利用水冷式MWA装置(KY-2000;Kangyou Microwave Energy Sources Institute,Nanjing,China)和仿肝组织体模获得热消融实验数据以验证温度场仿真模型的精确性。
本发明采取的技术方案是:首先在COMSOL软件中采用经典的麦克斯韦电磁方程和Pennes传热方程建立热消融仿真模型;基于现有技术来确定组织特性参数的正常取值范围,利用Minitab软件中的DOE对各特性参数进行敏感性分析,获得参数对温度的影响贡献率SS%;在微波热消融实验中利用测温针获得七个点的实测温度,其中选取一根测温针作为实测反馈针;基于参数在不同时刻的SS%对参数进行特异性反馈调节,由此获得各参数的最优表达形式;最后通过对比仿真数据和实验数据,验证此技术的可行性。
一种基于特性参数反馈的微波热消融温度场仿真方法,包括以下步骤:
步骤1,在COMSOL Multiphysics软件中设置预设初始条件和边界条件,建立微波热消融温度场仿真模型;
步骤2,基于仿真模型对体模组织特性参数进行敏感性分析;
步骤2.1,确定仿肝组织特性参数的取值范围,并在Minitab软件中通过DOE设计有关特性参数的16组实验,在Comsol软件中进行仿真;
步骤2.2,基于方差分析理论对仿真结果进行特性参数敏感性分析,获得各参数的方差贡献率SS%;
步骤3,执行微波热消融实验,利用测温针对消融区的七个点进行实测,获得实测温度,其中选取一根针作为反馈针,另外六根针作为验证针;
步骤4,基于参数敏感性分析结果,利用反馈针实测数据获得特性参数的精确表征形式;
步骤5,利用特性参数的精确表征形式进行温度场仿真,并与实测数据进行对比以验证温度场分布模型的精确性。
作为优选,步骤1中的仿真模型为:
其中,ρ为组织密度,SAR为比吸收率,T为温度,t为时间,c为比热容,k为导热率。
作为优选,步骤2.2中,选取微波发射点P0(1mm,0mm)、水冷侧近场点P1(10mm,-10mm)、无水冷侧近场点P2(10mm,10mm)和远场点P3(25mm,0mm)进行敏感性分析,获得各参数的方差贡献率SS%以研究模型参数对不同位置的温度分布的影响;基于公式(2)计算SS%以量化各参数对消融温度场的影响显著性:
其中,i为参数序号;Adj_SSi表示源自第i参数的误差方差和。
作为优选,步骤4中,
所述比热容的反馈表征函数为
其中,Cpat 25℃,Cpat 100℃分别表示25℃和100℃下的比热容。
作为优选,步骤4中,
所述介电常数的反馈表征函数为
其中,εat 25℃,εat 100℃分别表示25℃和100℃下的介电常数。
作为优选,步骤4中,
所述导热率的反馈表征函数为
其中,kat 25℃,kat 100℃分别表示25℃和100℃下的导热率。
作为优选,步骤4中,
所述电导率的反馈表征函数为
其中,σat 25℃,σat 100℃分别表示25℃和100℃下的电导率。
与现有技术相比,本发明具有如下优点和有益效果:
1.本发明能够基于敏感性分析获得各组织特性参数对温度场模型的影响程度,从而能够解决现有技术中无法对参数进行特异性分析的难题。
2.本发明提出了基于单针实测的分段参数反馈方法,能够获得不同消融组织参数的精确表征形式,因此在温度场仿真中获得的温度场分布更符合实际情况,具有更强的特异性和临床适用性。
3.本发明提供了比热容Cp的快速反馈方法,易于操作,能够精确地反馈出比热容的最优表达形式,对温度场仿真的研究具有重要意义。
附图说明
图1为根据本发明的实施方案的微波热消融实验装置的详细视图;
图2为根据本发明的实施方案的基于特性参数反馈的微波热消融温度场仿真方法的概述流程图;
图3为根据本发明的实施方案的特性参数SS%变化趋势的示意图;
图4为根据本发明的实施方案的54℃仿真等温面的示意图;
图5为根据本发明的实施方案的各验证针处的仿真值与实验值的对照示意图。
具体实施方式
下面结合附图和具体实施方案对本发明作进一步说明。
在本实施方案中,微波热消融实验装置包括微波消融仪KY-2000、数据采集装置(34970A;Agilent Technologies Inc.,Santa Clara,CA,USA)、测温针、蠕动泵和热消融体模。图1示出了根据本发明的实施方案的微波热消融实验装置的详细视图,其中微波消融仪的频率为2450MHz,功率设置为60W;水冷式微波天线具有1.9mm的直径和180mm的长度;微波能量采用缝隙发射,缝隙宽度为1.0mm,距尖端11mm。
在本实施方案中,微波热消融体模是基于江汉保等人研制的仿肝组织体模配方制备的并且由3.5重量份羧甲基纤维素钠、27.15重量份聚氯乙烯、0.35重量份氯化钠和69重量份蒸馏水构成。所制备的体模的尺寸为100×100×100mm3,其在25℃下的特性参数如表1所示。
表1 MW体模的基本特性参数表
本发明的基于特性参数反馈的微波热消融温度场仿真方法的概述流程图如图2所示,包括以下步骤:
步骤1,在COMSOL Multiphysics软件中,建立几何模型,设置物理参数,定义边界条件和初始条件,利用麦克斯韦电磁方程和Pennes传热方程获得温度场仿真模型。体模组织中因不存在血液灌注,其传热方程表示如下:
其中ρ为组织密度(kg/m3),SAR为比吸收率(W/m3),T为温度(℃),t为时间(s),c为比热容(J/(Kg·℃)),k为导热率(J/(s·m·℃))。
步骤2,基于仿真模型对体模组织特性参数进行敏感性分析。
步骤2.1,基于现有技术确定仿肝组织体模特性参数的取值范围,如表2所示,其中参数高水平取1,低水平取-1;在Minitab软件中通过DOE设计有关特性参数的16组实验,在Comsol软件中进行仿真;
表2模型参数和取值范围
步骤2.2,选取微波发射点P0(1mm,0mm)、水冷侧近场点P1(10mm,-10mm)、无水冷侧近场点P2(10mm,10mm)和远场点P3(25mm,0mm)进行敏感性分析,获得各参数的方差贡献率SS%以研究模型参数对不同位置的温度分布的影响;基于公式(2)计算SS%以量化各参数对消融温度场的影响显著性:
式中,i为参数序号;Adj_SSi表示源自第i参数的误差方差和。
图3示出了仿真过程中各特性参数SS%变化趋势的示意图,从该图可以看出,水冷侧近场点P1点模型参数对于不同时间的仿真温度具有特异性影响,因此选取P1作为反馈点。基于P1点敏感性分析结果,可知在前100秒内,比热容影响非常显著,可对该参数进行反馈,另外介电常数具有一定影响,因此在反馈出比热容之后对介电常数进行反馈;在400秒和500秒时刻,导热率和电导率对温度分布具有显著影响,因此在此时间对导热率和电导率进行反馈;水冷效果在整个消融过程中对温度的影响几乎为零,故不考虑。
步骤3,执行微波热消融实验,消融时间取10分钟,利用测温针对消融区的七个点进行实测,获得实测温度,其中选取P1位置的针作为反馈针,V1-V6位置的六根针作为验证针,测温点的坐标如表3所示;
表3测温点的坐标
步骤4,基于参数敏感性分析结果,利用反馈针实测数据获得特性参数的精确表征形式;
步骤4.1,在消融前100秒,对比热容进行反馈;基于传热方程,在该时间段内,由于电特性参数和导热率k对温度几乎无影响,因此温升斜率的误差主要源自比热容随温度的变化。设定比热容随温度变化的公式为Cp=Cpat25℃+kc(T-25),根据实验结果,Cpat25℃为3700(J/(kg·k)),100秒时的温度TS为48.61896℃,温升斜率ks为0.242,而仿真过程中温升斜率kM为0.375,因此基于公式Cpat25℃/(Cpat25℃+kc(T-25))=kS/kM,可求得kc的值为55.56。Cp的反馈表达形式如公式(3)所示:
步骤4.2,基于比热容的反馈表征函数,对介电常数进行反馈;由于前100秒时间内,比热容对温度结果具有一定影响,因此在此时间段对介电常数ε进行反馈分析。设定介电常数随温度变化的公式为ε=εat25℃+kε(T-25),调节kε的值,计算反馈点的仿真温度(TS)与实测温度(TM)之间的温度差值其中,取TM-TS在100秒时刻的绝对和,通过1stOpt(7D-Soft High Technology Inc.,Beijing,China)拟合获得与kε的函数关系最小化获得最优的介电常数变化率kε为0.99876。因此ε的反馈表达形式如公式(4)所示:
步骤4.3,基于比热容和介电常数的反馈表征函数,对导热率和电导率进行反馈调节;根据敏感性分析结果,在400秒-500秒时间内,导热率k和电导率σ对温度结果具有显著影响,因此在此时间段对k和σ进行反馈分析。设定导热率和电导率随温度变化的公式分别为k=kat25℃+kk(T-25)和σ=σat25℃+kσ(T-25),调节kk和kσ的值,计算反馈点的仿真温度(TS)与实测温度(TM)之间的温度差值其中,取TM-TS在400秒和500秒时刻的绝对和,通过1stOpt拟合获得与kk和kσ的函数关系最小化获得最优的导热率变化率kk为0.0066,最优的电导率变化率kσ为-0.01026。因此k和σ的反馈表达形式如公式(5)和公式(6)所示:
步骤5,利用特性参数的精确表征形式进行温度场仿真,并与实测数据进行对比以验证温度场分布模型的精确性。图4示出了基于反馈参数的54℃等温面示意图,其中黑色轮廓线为54℃消融边界。
利用验证针的实测温度数据验证温度场分布模型的精确性。图5示出了本发明方法中获得的仿真值与实验值的对照示意图。为了显示仿真值与实测值之间的误差,分别基于公式(7)、(8)、(9)获得最大误差(α)、平均误差(β)和标准偏差(δ)。表5列出了V1-V6点的相应误差。
其中,i表示时间采样点的序号,n表示时间采样点的总数,Ti S为第i采样点的仿真温度,Ti M为第i采样点的实测温度。
表4验证点的仿真值与实测值之间的误差
从表4中可以看出,反馈仿真结果与实测结果具有很好的一致性:各点最大误差的平均值为3.47℃;平均误差的平均值为1.55℃;各点标准偏差的平均值为0.95℃,满足临床需要的误差范围(小于5℃)。本发明提出的基于特性参数反馈的仿真方法能够针对不同组织进行特异性参数反馈,具有极好的实用性,可有效地解决微波消融仿真中特性参数无法精确获得的难题。