经颅聚焦超声治疗脑血栓的数值仿真研究

孙天宇，张艳秋，潘 婷，张默涵，凌子超，菅喜岐

(天津医科大学生物医学工程与技术学院，天津 300070)

0 引言

脑卒中是我国成年人致残和致死的首要原因，且血栓性缺血脑卒中居多[1]。经颅聚焦超声溶栓治疗是近年兴起的新型疗法，是将超声波聚焦在目标区域，利用超声波的机械效应、空化效应等进行溶栓的疗法，具有无创、颅内出血风险低、可重复施治的优势。

1999年Westermark等[2]使用频率为1.1 MHz的碎石机产生的冲击波、压电换能器连续超声和脉冲超声对体外模型的血凝块辐照30 s，观察到脉冲超声辐照后血凝块重量显著降低(99%±2%)。2009年Maxwell等[3]使用1 MHz的换能器进行体外实验，当脉冲重复频率1 kHz峰值负压达到-6 MPa时出现血栓溶解，峰值负压在-8 MPa及以上时溶栓效率显著。2012年Burgess等[4]使用高强度聚焦超声(High Intensity Focused Ultrasound,HIFU)对兔栓塞性卒中模型进行体内溶栓实验，结果表明溶栓率和输入声功率正相关。2013年Ahadi等[5]使用ExAblate4000半球形换能器在占空比50%和脉冲长度200 ms的条件下进行体外经离体人体颅骨实验，结果表明随声输入功率的增加，血凝块的裂解也将增加。同一年Hölscher等[6]使用ExAblate4000半球形换能器进行体外经离体人体颅骨实验，发现颅骨厚度越大血凝块的溶解率越低。2015年Xu等[7]研究了占空比分别为2.3%、9%和18%时，对脉冲聚焦超声诱导的溶栓效率的影响，结果发现使用9%的占空比时溶栓效率最高并且血凝块碎片的平均直径最小，且在这一占空比下的空化活性最强。

在临床治疗时由于人体颅骨的非均质结构，超声波穿过颅骨时可能出现散焦、焦点偏移和颅骨处热量沉积等问题，须通过对多阵元相控换能器的相位控制才能实现在颅内目标靶区精准聚焦。2003年Aubry等[8]在0.8～2.0 MHz范围内利用CT图像的亨氏值对组织的声学参数进行计算。2015年张千等[9]在0.5MHz的条件下利用CT图像的亨氏值对人体颅骨及脑组织的声学参数进行数值仿真研究。本文基于46岁健康男性志愿者头颅CT图像和82阵元随机分布相控阵换能器建立声场数值仿真模型，基于时间反转法和时域有限差分法(Finite Difference Time Domain,FDTD)数值解析Westervelt声波非线性传播方程获得阵元的激励信号，并基于该信号进行声压场的数值仿真。探究不同超声激励频率下所需的输入声功率，筛选效果最佳的参数，为临床试验中聚焦超声经颅溶栓治疗提供参考数据。

1 基本方程式

Westervelt声波非线性传播方程式[10]为

2 数值仿真模型及仿真参数

图1为切除大小为L×L(mm×mm)的颅骨作为声窗的开颅辐照的数值仿真模型。图2为与图1相对应的经颅辐照的数值仿真模型，在颅内治疗深度均为46 mm。82阵元凹球面换能器的曲率半径R=80 mm、开口直径为100 mm。FDTD法数值仿真区域为150 mm×150 mm×180 mm的立方体，数值仿真的时间步长为10 ns，数值仿真的空间步长为0.3 mm，声轴为z轴。

图1 开颅辐照数值仿真模型(单位：mm)Fig.1 Numerical simulation model of craniotomy irradiation(unit:mm)

图2 经颅辐照数值仿真模型(单位：mm)Fig.2 Numerical simulation model of transcranial irradiation(unit:mm)

基于志愿者CT图像的亨氏值H计算头颅组织的密度、声速、衰减系数等声学参数，其计算公式为[8]

其中：Φ为介质孔隙率，ρwater、cwater、αwater分别为水的密度、声速和衰减系数，ρbone、cbone、αbone分别为皮质骨的密度、声速和衰减系数。其他各项数值仿真参数如表1所示[11]。激励阵元的脉冲信号如图3所示，其中T为脉冲重复周期，T1为正弦信号所占时间长度(即脉冲宽度)，T1/T为占空比。基于时间反转法将该激励信号以虚拟点声源的形式设置在拟聚焦的目标靶点，数值仿真虚拟点声源的声压信号传播到每个阵元的声压信号，通过自相关和互相关法获取对应阵元的相位和输入声强的幅值补偿系数，将获得对应阵元相位差的幅值补偿系数与点声源信号相结合作为对应阵元的激励信号[12]。

表1 数值仿真参数表Table 1 Simulation parameters

图3 脉冲信号波形Fig.3 Pulse signal waveform

3 结果

3.1 声窗尺寸的影响

以图1所示的数值仿真模型为例，在脉冲重复周期T=0.2ms、脉冲宽度T1=0.02ms、占空比T1/T=10%、辐照时间t=0.2 s、f=0.8 MHz和输入声功率为60W的条件下，不同声窗尺寸L×L时形成焦点处负压如表2所示。由表2可知，当声窗尺寸大于等于85 mm×85 mm时，声窗基本对超声波的传播没有影响。

表2 不同开口大小及焦点处负压Table 2 Negative pressures at focus for different opened windows

3.2 颅骨的影响

以图1中声窗尺寸为85 mm×85 mm的开颅模型和图2中的经颅模型为例，在相同脉冲信号和不同超声激励频率的条件下形成焦点处-6 MPa时焦域声压分布、声轴(z轴)上的声压分布、声压最大的焦点处垂直z轴的x轴上的声压分布如图4所示，图5为焦点处-8 MPa时的声压分布图。图4(c)、4(d)及图5(c)、5(d)中，黑色实线为经颅轴向声压，灰色虚线为开颅轴向声压。由图4、5可知，随着超声激励频率的增大，负压焦域面积变小，焦域处的旁瓣增多，相同频率条件下经颅和开颅模型中分别形成-6 MPa和-8 MPa的焦域长轴和短轴长度几乎相等，经颅聚焦形成的旁瓣较强。形成焦点处最大负压分别为-6 MPa和-8 MPa所需输入声功率曲线和形成最大负压位置曲线如图6所示。由图6可知，相同频率情况下焦点处形成的负压越大所需输入声功率越大，经颅所需输入声功率远远大于开颅的输入声功率，约为开颅的1.5倍；随着频率的增加，焦点处负压达到-6 MPa或-8 MPa时所需输入声功率先减少后增加，激励频率为0.8 MHz时所需输入声功率均最小，开颅时达到-6 MPa和-8 MPa所需输入声功率分别为60 W和107 W，经颅时达到-6 MPa和-8 MPa所需输入声功率分别为95 W和170 W；0.5～0.6 MHz时焦点向换能器方向前移，其他频率时均在偏离设定焦点80 mm处向换能器方向前移1.1 mm的78.9 mm处形成焦点。

图4 焦点声压为-6 MPa时经颅和开颅的焦域声压分布图Fig.4 Sound pressure distribution diagrams in the focal areas of transcranial and craniotomy models when the sound pressure at focus is-6 MPa

图5 焦点声压为-8 MPa时经颅和开颅的焦域声压分布图Fig.5 Sound pressure distribution diagrams in the focal areas of transcranial and craniotomy models when the sound pressure at focus is-8 MPa

图6 焦点负声压达到要求所需的输入声功率和焦点位置随频率变化曲线Fig.6 Variation curves of the required input acoustic power for the special negative pressure level at focus and the focus position with frequency

3.3 辐照时间的影响

在开颅和经颅模型中超声激励频率为0.8 MHz的条件下，分别形成-6 MPa和-8 MPa对应输入声功率分别为60 W和107 W和95 W和170 W时，不同辐照时间时形成的焦点位置在-6～-8 MPa的焦域面积如表3所示。由表3可知，焦点位置和焦域面积与辐照时间无关，开颅模型中形成的焦域面积略大于经颅模型中形成的焦域面积。

表3 不同辐照时间的条件下经颅和开颅模型的焦点位置和焦域面积Table 3 Focus position and focal area in the transcranial and craniotomy models under different irradiation times

3.4 占空比的影响

在开颅和经颅模型中超声激励频率为0.8 MHz的条件下，分别形成-6 MPa和-8 MPa对应输入声功率分别为60 W和107 W和95 W和170 W时，不同占空比时形成的焦点位置和焦域声压在-6～-8 MPa的焦域面积如表4所示。由表4可知，焦点位置和焦域面积与占空比无关，开颅模型中形成的焦域面积略大于经颅模型中形成的焦域面积。

表4 不同占空比的条件下经颅和开颅模型的焦点位置和焦域面积Table 4 Focus position and focal area in the transcranial and craniotomy models under different pulse duty factors

4 结论

本文基于82阵元凹球面随机分布相控换能器，结合志愿者头颅CT图像扫描数据和水体、脑组织的相关参数，建立聚焦超声经颅溶栓治疗的数值仿真模型，得到以下结果：

(1)在颅内治疗深度为46 mm、颅骨处的开口直径为73.68 mm且声窗尺寸≥85 mm×85 mm时，声窗对超声波的传播几乎无影响。

(2)相同频率情况下焦点处形成的负压越大所需输入声功率越大，经颅所需输入声功率大于开颅的输入声功率，约为开颅的1.5倍；随着频率的增加，焦点处负压达到-6 MPa或-8 MPa时所需输入声功率先减少后增加，激励频率为0.8 MHz时所需输入声功率均最小。此时颅骨的平均厚度接近其一个波长(4 mm)，透声系数最大，与声学原理相一致。开颅时达到-6 MPa和-8 MPa所需输入声功率分别为60 W和107 W，经颅时达到-6 MPa和-8 MPa所需输入声功率分别为95 W和170 W。

(3)0.5～0.6 MHz时焦点向换能器方向前移，其他频率时均在偏离设定焦点80 mm处向换能器方向前移1.1 mm的78.9 mm处形成焦点。

(4)随着超声激励频率的增大，焦域面积变小，焦域处的旁瓣增多，相同频率条件下经颅和开颅模型中形成的焦域形状、大小相近，经颅聚焦形成的旁瓣较强。

(5)占空比为10%～50%且辐照时间小于0.6 s时，焦点位置和焦域面积与辐照时间和占空比无关，且开颅模型中形成的焦域面积略大于经颅模型中形成的焦域面积。

基于上述结果得到以下结论：

(1)频率相同时焦点处形成的负压越大所需输入声功率越大，经颅所需输入声功率约为开颅时的1.5倍。

(2)频率越大焦域面积越小但焦域处的旁瓣增多，频率相同时经颅和开颅的焦域形状、大小相近但经颅时的旁瓣较强。

(3)占空比为10%～50%且辐照时间小于0.6 s时，占空比和辐照时间对焦点位置和焦域面积没有影响。

(4)经颅聚焦治疗所需输入声功率约为开颅的1.5倍。

5 讨论

在实际临床应用中由于颅骨的缺失、病变等无法进行经颅治疗的情况发生，为此，研究并建立开颅和经颅两个模型，对比数值仿真分析了两种情况下形成的焦域特性。临床治疗时热损伤也是需要考虑的问题，在这里以经颅骨模型中超声频率最大1.0 MHz、辐照时间最长0.6 s、输入功率最大270 W及占空比50%的情况为例，基于Pennes生物热传导方程[12-15]数值仿真结果发现焦点处的最高温升约2℃，不会出现热损伤。

本文仅以一名志愿者的头颅CT数据为例建立相控阵超声换能器开颅和经颅数值仿真模型，发现激励频率在0.8 MHz颅骨厚度接近其波长时所需输入声功率最低，针对不同志愿者头颅的差异性及其颅骨厚度与曲率等的相关性和不同治疗深度的研究正在进行之中。