在高強度聚焦超聲（HIFU）腦腫瘤治療過程中需嚴格控制治療溫度。本研究通過調控相控換能器激勵信號實現了焦域溫度均勻分布的目的。首先，本文利用志愿者頭部計算機斷層掃描（CT）數據和 82 陣元相控換能器建立開顱 HIFU 治療腦腫瘤的三維數值仿真模型，通過疊加聚焦于兩個設定目標點的信號并調控兩信號間激勵時間差和幅值，研究其對 HIFU 焦域溫度分布及焦域形狀大小的調控。研究結果表明，兩目標點間距離在一定范圍內通過調節兩激勵信號的激勵時間差和幅值可以實現焦域內溫度均勻分布，同時可調控焦域形狀和體積大小。本文研究的仿真結果或可為 HIFU 安全有效地應用于臨床治療提供理論方法和參考。
引用本文： 常詩卉, 王佩國, 菅喜岐. 高強度聚焦超聲腦腫瘤治療焦域溫度均勻分布調控的數值仿真研究. 生物醫學工程學雜志, 2018, 35(6): 877-886. doi: 10.7507/1001-5515.201802009 復制
引言
高強度聚焦超聲（high-intensity focused ultrasound，HIFU）作為新興腫瘤物理治療手段，以低強度聲能無創透過淺層組織，而以高強度聚焦聲能治療深部腫瘤且可重復施治等優點受到臨床廣泛關注，目前已應用于乳腺癌、子宮肌瘤、前列腺癌等實體腫瘤的臨床治療[1-2]。對于腦腫瘤治療而言，由于腦組織包含眾多功能區，因此需要對治療溫度進行嚴格控制，如焦域中心溫度過高可能引發正常腦組織損傷或嚴重后遺癥[3-5]。
在早期的 HIFU 顱腦聚焦研究中，由于顱骨的密度、聲速等與軟組織相差較大，且對超聲波的衰減和吸收較強，研究者需要去除一部分顱骨，以便超聲波進入腦部[6-7]，如 2006 年 Ram 等[8]對 3 名復發膠質瘤患者行開顱術后進行 HIFU 治療。近年來，隨著相控技術、醫學影像技術的發展及相位校正、幅值補償算法的提出，實現了 HIFU 經顱聚焦并開展了相關腦部疾病的臨床治療試驗。2010 年 McDannold 等[3]利用磁共振成像（magnetic resonance imaging，MRI）引導聚焦超聲治療系統（ExAblate 3000，InSightec，Israel）對 3 位膠質瘤患者進行 HIFU 經顱臨床治療試驗，但由于設備功率的限制，3 位患者顱內靶區未形成凝固型壞死。2013 年 Elias 等[4]利用 MRI 引導的聚焦超聲對 15 名震顫患者進行治療，結果顯示在 4 名患者身上出現了持續性感覺異常等副作用。此外，Chang 等[5]利用 MRI 引導聚焦超聲治療系統（ExAblate 4000，InSightec，Israel）對 11 位震顫患者進行除顫治療，其中 3 名患者由于目標區域未達到治療溫度導致治療失敗。
為了提高 HIFU 治療的安全有效性，國內外研究者利用仿組織體模和動物肌肉組織等進行了一些研究[9-11]。本文針對腦腫瘤治療中焦域溫度分布的調控展開研究，基于一名 46 歲男性健康志愿者的頭部計算機斷層掃描（computed tomography，CT）數據和 82 陣元相控換能器建立 HIFU 開顱治療腦腫瘤三維數值仿真模型。基于時間反轉法在靶區相隔一定距離設定兩目標焦點，并分別獲取聚焦這兩焦點的陣元激勵信號，通過這兩信號的疊加及其激勵時間差和幅值的調控，研究焦域溫度均勻分布及焦域形狀大小的調控方法，為安全有效的 HIFU 腦腫瘤治療方案制定提供理論方法和數據。
1 數值仿真的理論基礎
本研究利用韋斯特維爾特（Westervelt）聲波非線性傳播方程計算聲場分布[12-13]，如式（1）所示：
 |
式中，
為拉普拉斯算子，p 為聲壓，c 和 ρ 分別為介質的聲速和密度，β 為聲波非線性系數，t 為時間，聲波擴散系數
，α 為吸收系數，角頻率
= 2πf，f 為換能器的激勵頻率。
不考慮血流灌注時的佩內斯（Pennes）生物熱傳導方程[14-15]，如式（2）所示：
 |
式中，Cr 為組織比熱，T 為溫度，r 為媒質熱傳導率。單位體積的發熱量 Q = 2αI，其中聲強
，Tp 為激勵函數的周期。HIFU 腦腫瘤治療溫度一般控制在 55～65℃，本研究控制最高溫度不超過 65℃。
等效熱劑量 t43 的計算式如式（3）所示[16]：
 |
式中，Tt 為輻照 t 秒時組織溫度，t0 和 tfinal 分別指輻照開始和結束時間，R 的取值隨溫度而變化。溫度大于等于 43℃ 時，取值 0.5，溫度低于 43℃ 時，取值 0.25。本研究將等效熱劑量 90 min 以上的區域定義為可治療焦域[17-18]。
2 仿真模型及調控方法
2.1 志愿者信息采集
本研究利用一名 46 歲男性健康志愿者的頭部 CT 圖像進行三維重建，建立了數值仿真模型，CT 圖像由天津醫科大學附屬腫瘤醫院提供（掃描參數為 120 kV，100 mA；掃描層厚 3.0 mm）。志愿者簽署知情同意書之后進行了 CT 掃描。
2.2 仿真模型及陣元激勵信號
如圖 1 所示為開顱狀態下 HIFU 腦腫瘤治療數值仿真模型以及擬聚焦目標點設定位置示意圖，其中仿真模型由換能器、水及腦組織構成，換能器為 82 陣元半隨機分布相控換能器，內表面曲率半徑為 80 mm、開口直徑 100 mm、陣元直徑 8.0 mm，激勵頻率為 0.7 MHz；腦組織為 CT 圖像三維重建后去除聲窗處顱骨后的腦組織；換能器內表面與腦組織表面在聲軸（z 軸）上的距離為 55 mm；數值仿真區域為邊長 100 mm 的正方體，仿真過程中空間步長為 0.25 mm，時間步長為 10 ns。
圖1
數值仿真模型和設定焦點位置圖
Figure1.
The simulation model and the preset targets position schematic
在聲軸方向或垂直聲軸方向設定兩個目標焦點 F1 和 F2，其位置圖如圖 1 右圖所示。首先在 F1 處設置點聲源
，時域有限差分（finite difference in time domain，FDTD）數值仿真獲得換能器各陣元上的聲壓信號，并經時間反轉法處理后得到聚焦于 F1 的陣元激勵信號，如式（4）所示：
 |
其中，I1i 為換能器各陣元激勵信號聲強，
為各陣元激勵信號初始相位，下標 i 為陣元標號。同理在 F2 處也設置點聲源
，獲得聚焦于 F2 處的陣元激勵信號，如式（5）所示：
 |
同時，實現 F1、F2 處聚焦的陣元激勵信號如式（6）所示：
 |
Δt 為
相對于
觸發的時間延遲，M 為激勵信號幅值調控系數，S1i 與 S2i 在 82 個陣元上的激勵總功率相等時，M = 0.50。
2.3 仿真參數
腦組織的密度（ρ）、聲速（c）、衰減系數（α）由志愿者頭部 CT 圖像的 CT 值（以亨氏單位 H 表示）計算得到，計算公式如式（7）～式（10）所示[19]：
 |
 |
 |
 |
式中，Φ 為孔隙率，ρbone、cbone、αbone 分別為皮質骨的密度、聲速、衰減系數，在仿真計算中分別取 1 600 kg/m3、3 200 m/s、8 dB/mm；ρwater、cwater、αwater 分別為水的密度、聲速、衰減系數，仿真計算中分別取 998 kg/m3、1 500 m/s、0.2 dB/mm。其它參數如表 1 所示[17, 19-20]，T0 為開始輻照前的初始溫度。
表1
仿真中用到的常數參數
Table1.
Constant parameters used in simulation
3 仿真結果
3.1 聲軸方向的調控
3.1.1 時間延遲 Δt 的影響
如圖 1 所示，以仿真模型和聲軸方向設定焦點為例，當設定焦點坐標為 F1（0，0，75.00）、F2（0，0，85.00）即 L = 10.0 mm 時，在 82 陣元換能器的輸入總功率 Pw = 41.22 W，激勵信號幅值調控參數 M = 0.50，焦域內最高溫度為 65℃ 的條件下，聲波激勵時間差 Δt 在激勵信號的一個周期內變化時形成的溫度場分布如圖 2 所示，當 Δt 在一個聲波周期內變化時，均在聲軸上形成了 55℃ 以上焦域；其中 Δt = 0～600 ns 時，焦域為圍繞兩個設定聚焦點的一個連續區域；Δt = 800～1 200 ns 時，焦域為分別圍繞兩個設定聚焦點的兩個區域。
圖2
最高溫度為 65℃ 時不同 Δt（0～1 400 ns）的溫度場分布（L = 10.0 mm，M = 0.50，Pw = 41.22 W）
Figure2.
The temperature distribution with different Δt (0～1 400 ns) when the maximum focal temperature was 65℃(L = 10.0 mm,M = 0.50,Pw = 41.22 W)
3.1.2 平臺式溫度的調控
在 L = 10.0 mm、Δt = 700 ns 條件下，最高溫度達到 65℃ 時的溫度分布圖如圖 3 所示。當 M = 0.50 時的溫度分布如圖 3 左上圖所示，對應的聲軸溫度為曲線圖上灰色曲線，該曲線上有距離為 Δz 的兩個溫度峰值 T1、T2，且最大溫度峰值與溫度低谷之間的最大溫差 ΔTm = 17.7℃。當調控 M = 0.38 時，形成如圖 3 左下圖所示的溫度分布，對應聲軸溫度曲線為黑色曲線所示，此時兩溫度峰值 T1、T2 相等，溫差 ΔTm 減小至 5.2℃；調控 M 不會改變溫度峰值的位置。
圖3
焦域溫度分布及聲軸溫度曲線圖（L = 10.0 mm，Δt = 700 ns，Pw = 41.22 W）
Figure3.
Temperature distribution and the temperature along the acoustic axis (L = 10.0 mm, Δt = 700 ns, Pw = 41.22 W)
由圖 3 可知，通過調控 M 使 T1 = T2 時，ΔTm 會相應減小。在忽略熱傳導的情況下通過公式（2）可以推導使得 T1 = T2 的 M 值計算式如式（11）所示：
 |
在最高溫度達到 65℃、T1 = T2 的條件下，M 隨 Δt 的變化曲線如圖 4 左上圖所示，其中黑色點線為|T1-T2| ≤ 0.05℃ 的數值仿真 M 值，灰色點線為利用式（11）的計算值，由該圖可知仿真值與計算值幾乎一致。調節 M 值后，ΔTm 隨 Δz 的變化散點圖及擬合曲線如圖 4 右上圖所示，其中 ΔTm 隨 Δz 的擬合變化關系式，如式（12）所示（擬合度
= 0.94）：
 |
若將 ΔTm ≤ 1℃ 作為平臺式溫度分布的條件，則根據式（12）計算得 0 < Δ z ≤ 8.0 mm，再對照圖 4 下圖可知，要調控出平臺式均勻溫度分布的 Δt 取值范圍為 50 ns ≤ Δt ≤ 150 ns。如圖 5 所示為在該條件范圍內形成的平臺式溫度分布圖，可治療焦域體積為 48.3～33.6 mm3。
圖5
聲軸方向均勻溫度分布圖（L = 10.0 mm，Pw = 41.22 W）
Figure5.
The uniform temperature distribution along acoustic axis (L = 10.0 mm, Pw = 41.22 W)
圖4
65℃ 時 M、ΔTm 及 Δz 的變化曲線（L = 10.0 mm，Pw = 41.22 W）
Figure4.
The variation of M, ΔTm and Δzwhen Δt was changed with the maximum temperature reached 65℃(L = 10.0 mm, Pw = 41.22 W)
3.1.3 L 對焦域溫度分布的影響
將 F1 固定，F2 向 F1 相對于靠近換能器一側調節時，F1 與 F2 之間的距離 L 設定為 Lnear，F2 位于 F1 遠離換能器一側時設定為 Lfar。在 Pw = 41.22 W 和 M = 0.50 的條件下，不同 L 時最高溫度達到 65℃ 所需輻照時間 t 隨 Δt 的變化曲線如圖 6 所示，其中'o'表示形成一個連續焦域，'×'表示形成兩個焦域。
圖6
不同 L 時輻照時間隨 Δt 的變化曲線（Pw = 41.22 W，M = 0.50；‘o’表示形成一個 55℃ 以上焦域，‘×’形成兩個焦域）
Figure6.
The variation of irradiation time with Δt under different L (Pw = 41.22 W, M = 0.50; ‘o’ stand for one continue focal region upon 55℃ was formed and ‘×’stand for two focal regions were formed)
由圖 6 可知，當 L ≤ 12.5 mm 時，通過 Δt 的調控可實現一個焦域的聚焦，且 L 越小可實現一個焦域聚焦的 Δt 的范圍越大；當 L ≥ 15.0 mm 時，通過調節 Δt 無法實現一個焦域聚焦。
在不同 L 的條件下，Pw = 41.22 W、焦域最高溫度達到 65℃ 時，溫度分布中兩溫度峰值間的距離 Δz 隨 Δt 的變化曲線如圖 7 所示。
圖7
不同 L 時 Δz 隨 Δt 的變化曲線圖（Pw = 41.22 W）
Figure7.
The variation of Δz with Δt under different L(Pw = 41.22 W)
由圖 7 可知，當 7.5 mm ≤ L ≤ 12.5 mm 時，通過調節 Δt 可以使得 0 < Δ z ≤ 8 mm，在此條件下調節 M 可以形成 T1 = T2 且 ΔTm ≤ 1℃ 的平臺式溫度分布。如表 2 所示為調控形成平臺式溫度分布的 Δt、M 及可治療焦域體積。由表 2 可知，相同的 L 條件下，通過 Δt 和 M 的調節，可在設定擬聚焦的目標位置實現平臺溫度分布的焦域；不同的 L 條件下，形成的平臺式溫度分布可治療體積也不同，L 越大形成的可知治療焦域體積越大。L 在 7.5～12.5 mm 的范圍內平臺式溫度分布可治療焦域的體積調控范圍為 30.09～71.56 mm。
表2
不同 L 時得到沿聲軸方向平臺溫度分布的參數及可治療焦域體積
Table2.
Parameters for generating the uniform temperature field along acoustic axis with different L and the volume of treatable focal region
3.2 垂直于聲軸方向的溫度分布調控
3.2.1 Δt 對溫度分布的影響
如圖 1 所示，在聲軸兩側設定 F1（0，1.00，75.00）和 F2（0，–1.00，75.00）（L = 2.0 mm），當最高溫度達到 65℃ 時，形成的溫度分布（Pw = 41.22 W，M = 0.50）如圖 8 所示。
圖8
Δt = 0～1 400 ns、最高溫度為 65℃ 時溫度分布圖（L = 2.0 mm，M = 0.50，Pw = 41.22 W）
Figure8.
Δt = 0～1 400 ns, the temperature distribution when the maximum focal temperature reached 65℃ (L = 2.0 mm,M = 0.50, Pw = 41.22 W)
由圖 8 可知，隨著 Δt 由 0 增大到 1 400 ns，55℃ 以上焦域由一個連續區域逐漸分離為兩個區域，隨后又逐漸融合為一個連續區域。
3.2.2 垂直于聲軸方向平臺溫度分布的調控
以圖 8 中 Δt = 400 ns 的溫度分布為例，當調控 M 時沒有形成沿 y 軸方向的均勻溫度分布。在這里當設定 F1（0，1.00，75.00）、F2（0，–1.00，77.50），調控 M = 0.46，t = 10.9 s 時形成最高溫度 65℃ 的溫度分布及 55℃ 以上區域如圖 9 的上半部分所示，其中左圖為 y-z 平面溫度分布圖，右圖的白色區域分別為 55℃ 以上焦域的 y-z 及 x-z 截面圖。由圖 9 可知，y 方向上形成了 ΔTm = 0.99℃ 的平臺溫度分布；55℃ 以上焦域在 y-z 平面上為 8 × 4 mm （底 × 高）的平行四邊形，x-z 平面上為 8 × 2 mm （長軸 × 短軸）的橢圓形。而在 F1 固定不變，F2 設定為（0，–1.00，80.00）的條件下，M = 0.47 和輻照時間 t = 23.0 s 時形成溫度場分布及 55℃ 以上焦域如圖 9 的下半部分所示，由圖可知，y 方向上溫度分布均勻（ΔTm = 0.61℃），55℃ 以上焦域在 y-z 平面上近似 10 × 5 mm （長 × 寬）的四邊形，x-z 平面上為 10 × 3 mm （長軸 × 短軸）的橢圓形。
圖9
y 軸方向均勻溫度分布圖（L = 2.0 mm，Pw = 41.22 W）
Figure9.
The uniform temperature distribution on y direction (L = 2.0 mm, Pw = 41.22 W)
3.2.3 L 對平臺溫度分布焦域的調控
在垂直聲軸方向設定不同距離 L 的目標焦點，即當 L = 2.0～3.0 mm 時，通過調節目標焦點相對位置及 Δt 和 M 的值，實現了 y 軸方向均勻的溫度分布，調控參數如表 3 所示。
表3
不同 L 時得到 y 軸方向均勻溫度分布的參數及可治療焦域大小
Table3.
Parameters and size of treatable focal region along with uniform temperature field on y direction under different L
由表 3 可知，均勻溫度分布區域的溫度差 ΔTm 均小于 1℃，可治療焦域在 y 方向的短軸長度可實現 4.25～6.25 mm 內的調控，可治療焦域的體積可在 58.7～201.8 mm3內調控。
4 討論與分析
本研究在靶區設定兩個間隔一定的擬聚焦目標點，采用時間反轉法獲取其激勵信號并進行疊加激勵換能器陣元，然后通過改變觸發時間延遲和幅值來調節焦域的溫度分布。時間反轉方法是目前 HIFU 顱腦聚焦較普遍采用的方法，其能夠通過相位矯正在目標區域形成焦域[17-19]。本研究調控過程中，調控兩激勵信號的觸發時間延遲會改變焦域溫度分布，但最終形成焦域是圍繞兩個設定目標焦點位置的區域；調控激勵信號的幅值比例系數不會改變焦域位置。為實現焦域內均勻的溫度分布，本研究采用雙激勵信號疊加激勵換能器的方法，不對激勵信號進行調控，得到的焦域內會有一個或兩個能量集中的溫度峰值；通過調控激勵信號，得到了溫度起伏 ΔTm 不超過 1℃ 的溫度分布均勻的焦域。本研究以 ΔTm 不超過 1℃ 作為焦域溫度均勻分布的標準；在實際治療過程中，需要根據治療需求，對焦域內溫度起伏 ΔTm 進行調節。在傳統的 HIFU 單焦點輻照過程中，往往通過改變輻照功率或輻照時間改變焦域的溫度分布和體積，這樣易引發淺部正常組織的傷害[21-22]。本研究在輻照功率、輻照時間一定的前提下實現了焦域溫度分布的調控，并形成了沿聲軸和垂直聲軸的均勻溫度分布，同時可調控治療焦域長、短軸長度，從而改變治療焦域體積大小和形狀。
上述仿真結果表明，調控得到垂直聲軸方向均勻溫度分布的 L 和 Δt 的范圍均較小，這是由于 HIFU 聚焦形成的焦域垂直聲軸方向長度（短軸）較短。本研究通過兩個激勵信號疊加并調控激勵信號間時間延遲和幅值調控系數，實現了單次輻照焦域均勻溫度分布，與 Fan 等[23]提出的理想 HIFU 治療焦域相一致；且在實現焦域溫度均勻分布的同時，可根據臨床實際的需要調控焦域體積、形狀和焦域內溫度分布。當治療較大體積腫瘤時，需要將多次輻照形成損傷進行堆疊，對輻照路徑及間隔時間的優化有待于進一步深入研究。
本研究采用一名 46 歲男性健康志愿者頭部 CT 數據建立三維數值仿真模型，研究 HIFU 腦腫瘤治療焦域溫度分布的調控，該模型與實際腦腫瘤治療有一定差異，但腦腫瘤組織與正常腦組織的差異性非常小，幾乎不會對 HIFU 焦域溫度場產生影響[24]。而在制定 HIFU 腦腫瘤治療方案時，需要準確定位腦腫瘤的位置、大小以及形狀，并根據這些參數制定安全、有效的治療計劃，這部分研究正在進行之中。此外，本研究通過調控雙激勵信號的觸發延遲時間和幅值比例系數實現對焦域溫度分布的調控，其中信號間的觸發延遲時間為納秒（ns）級，這期間組織溫升很小，不會影響組織的物理及生物特性；在整個治療過程中，靶區組織會產生變性，考慮組織變性后聲學及熱力學影響的研究正在進行中。HIFU 腦部疾病治療的理想方式是無需開顱術的 HIFU 經顱治療。目前，在本研究結果及方法的基礎上，考慮顱骨的 HIFU 經顱聚焦焦域溫度分布的調控及相應的實驗驗證研究正在進行之中。
引言
高強度聚焦超聲（high-intensity focused ultrasound，HIFU）作為新興腫瘤物理治療手段，以低強度聲能無創透過淺層組織，而以高強度聚焦聲能治療深部腫瘤且可重復施治等優點受到臨床廣泛關注，目前已應用于乳腺癌、子宮肌瘤、前列腺癌等實體腫瘤的臨床治療[1-2]。對于腦腫瘤治療而言，由于腦組織包含眾多功能區，因此需要對治療溫度進行嚴格控制，如焦域中心溫度過高可能引發正常腦組織損傷或嚴重后遺癥[3-5]。
在早期的 HIFU 顱腦聚焦研究中，由于顱骨的密度、聲速等與軟組織相差較大，且對超聲波的衰減和吸收較強，研究者需要去除一部分顱骨，以便超聲波進入腦部[6-7]，如 2006 年 Ram 等[8]對 3 名復發膠質瘤患者行開顱術后進行 HIFU 治療。近年來，隨著相控技術、醫學影像技術的發展及相位校正、幅值補償算法的提出，實現了 HIFU 經顱聚焦并開展了相關腦部疾病的臨床治療試驗。2010 年 McDannold 等[3]利用磁共振成像（magnetic resonance imaging，MRI）引導聚焦超聲治療系統（ExAblate 3000，InSightec，Israel）對 3 位膠質瘤患者進行 HIFU 經顱臨床治療試驗，但由于設備功率的限制，3 位患者顱內靶區未形成凝固型壞死。2013 年 Elias 等[4]利用 MRI 引導的聚焦超聲對 15 名震顫患者進行治療，結果顯示在 4 名患者身上出現了持續性感覺異常等副作用。此外，Chang 等[5]利用 MRI 引導聚焦超聲治療系統（ExAblate 4000，InSightec，Israel）對 11 位震顫患者進行除顫治療，其中 3 名患者由于目標區域未達到治療溫度導致治療失敗。
為了提高 HIFU 治療的安全有效性，國內外研究者利用仿組織體模和動物肌肉組織等進行了一些研究[9-11]。本文針對腦腫瘤治療中焦域溫度分布的調控展開研究，基于一名 46 歲男性健康志愿者的頭部計算機斷層掃描（computed tomography，CT）數據和 82 陣元相控換能器建立 HIFU 開顱治療腦腫瘤三維數值仿真模型。基于時間反轉法在靶區相隔一定距離設定兩目標焦點，并分別獲取聚焦這兩焦點的陣元激勵信號，通過這兩信號的疊加及其激勵時間差和幅值的調控，研究焦域溫度均勻分布及焦域形狀大小的調控方法，為安全有效的 HIFU 腦腫瘤治療方案制定提供理論方法和數據。
1 數值仿真的理論基礎
本研究利用韋斯特維爾特（Westervelt）聲波非線性傳播方程計算聲場分布[12-13]，如式（1）所示：
|
式中，
為拉普拉斯算子，p 為聲壓，c 和 ρ 分別為介質的聲速和密度，β 為聲波非線性系數，t 為時間，聲波擴散系數
，α 為吸收系數，角頻率
= 2πf，f 為換能器的激勵頻率。
不考慮血流灌注時的佩內斯（Pennes）生物熱傳導方程[14-15]，如式（2）所示：
|
式中，Cr 為組織比熱，T 為溫度，r 為媒質熱傳導率。單位體積的發熱量 Q = 2αI，其中聲強
，Tp 為激勵函數的周期。HIFU 腦腫瘤治療溫度一般控制在 55～65℃，本研究控制最高溫度不超過 65℃。
等效熱劑量 t43 的計算式如式（3）所示[16]：
|
式中，Tt 為輻照 t 秒時組織溫度，t0 和 tfinal 分別指輻照開始和結束時間，R 的取值隨溫度而變化。溫度大于等于 43℃ 時，取值 0.5，溫度低于 43℃ 時，取值 0.25。本研究將等效熱劑量 90 min 以上的區域定義為可治療焦域[17-18]。
2 仿真模型及調控方法
2.1 志愿者信息采集
本研究利用一名 46 歲男性健康志愿者的頭部 CT 圖像進行三維重建，建立了數值仿真模型，CT 圖像由天津醫科大學附屬腫瘤醫院提供（掃描參數為 120 kV，100 mA；掃描層厚 3.0 mm）。志愿者簽署知情同意書之后進行了 CT 掃描。
2.2 仿真模型及陣元激勵信號
如圖 1 所示為開顱狀態下 HIFU 腦腫瘤治療數值仿真模型以及擬聚焦目標點設定位置示意圖，其中仿真模型由換能器、水及腦組織構成，換能器為 82 陣元半隨機分布相控換能器，內表面曲率半徑為 80 mm、開口直徑 100 mm、陣元直徑 8.0 mm，激勵頻率為 0.7 MHz；腦組織為 CT 圖像三維重建后去除聲窗處顱骨后的腦組織；換能器內表面與腦組織表面在聲軸（z 軸）上的距離為 55 mm；數值仿真區域為邊長 100 mm 的正方體，仿真過程中空間步長為 0.25 mm，時間步長為 10 ns。
圖1
數值仿真模型和設定焦點位置圖
Figure1.
The simulation model and the preset targets position schematic
在聲軸方向或垂直聲軸方向設定兩個目標焦點 F1 和 F2，其位置圖如圖 1 右圖所示。首先在 F1 處設置點聲源
，時域有限差分（finite difference in time domain，FDTD）數值仿真獲得換能器各陣元上的聲壓信號，并經時間反轉法處理后得到聚焦于 F1 的陣元激勵信號，如式（4）所示：
|
其中，I1i 為換能器各陣元激勵信號聲強，
為各陣元激勵信號初始相位，下標 i 為陣元標號。同理在 F2 處也設置點聲源
，獲得聚焦于 F2 處的陣元激勵信號，如式（5）所示：
|
同時，實現 F1、F2 處聚焦的陣元激勵信號如式（6）所示：
|
Δt 為
相對于
觸發的時間延遲，M 為激勵信號幅值調控系數，S1i 與 S2i 在 82 個陣元上的激勵總功率相等時，M = 0.50。
2.3 仿真參數
腦組織的密度（ρ）、聲速（c）、衰減系數（α）由志愿者頭部 CT 圖像的 CT 值（以亨氏單位 H 表示）計算得到，計算公式如式（7）～式（10）所示[19]：
|
|
|
|
式中，Φ 為孔隙率，ρbone、cbone、αbone 分別為皮質骨的密度、聲速、衰減系數，在仿真計算中分別取 1 600 kg/m3、3 200 m/s、8 dB/mm；ρwater、cwater、αwater 分別為水的密度、聲速、衰減系數，仿真計算中分別取 998 kg/m3、1 500 m/s、0.2 dB/mm。其它參數如表 1 所示[17, 19-20]，T0 為開始輻照前的初始溫度。
表1
仿真中用到的常數參數
Table1.
Constant parameters used in simulation
3 仿真結果
3.1 聲軸方向的調控
3.1.1 時間延遲 Δt 的影響
如圖 1 所示，以仿真模型和聲軸方向設定焦點為例，當設定焦點坐標為 F1（0，0，75.00）、F2（0，0，85.00）即 L = 10.0 mm 時，在 82 陣元換能器的輸入總功率 Pw = 41.22 W，激勵信號幅值調控參數 M = 0.50，焦域內最高溫度為 65℃ 的條件下，聲波激勵時間差 Δt 在激勵信號的一個周期內變化時形成的溫度場分布如圖 2 所示，當 Δt 在一個聲波周期內變化時，均在聲軸上形成了 55℃ 以上焦域；其中 Δt = 0～600 ns 時，焦域為圍繞兩個設定聚焦點的一個連續區域；Δt = 800～1 200 ns 時，焦域為分別圍繞兩個設定聚焦點的兩個區域。
圖2
最高溫度為 65℃ 時不同 Δt（0～1 400 ns）的溫度場分布（L = 10.0 mm，M = 0.50，Pw = 41.22 W）
Figure2.
The temperature distribution with different Δt (0～1 400 ns) when the maximum focal temperature was 65℃(L = 10.0 mm,M = 0.50,Pw = 41.22 W)
3.1.2 平臺式溫度的調控
在 L = 10.0 mm、Δt = 700 ns 條件下，最高溫度達到 65℃ 時的溫度分布圖如圖 3 所示。當 M = 0.50 時的溫度分布如圖 3 左上圖所示，對應的聲軸溫度為曲線圖上灰色曲線，該曲線上有距離為 Δz 的兩個溫度峰值 T1、T2，且最大溫度峰值與溫度低谷之間的最大溫差 ΔTm = 17.7℃。當調控 M = 0.38 時，形成如圖 3 左下圖所示的溫度分布，對應聲軸溫度曲線為黑色曲線所示，此時兩溫度峰值 T1、T2 相等，溫差 ΔTm 減小至 5.2℃；調控 M 不會改變溫度峰值的位置。
圖3
焦域溫度分布及聲軸溫度曲線圖（L = 10.0 mm，Δt = 700 ns，Pw = 41.22 W）
Figure3.
Temperature distribution and the temperature along the acoustic axis (L = 10.0 mm, Δt = 700 ns, Pw = 41.22 W)
由圖 3 可知，通過調控 M 使 T1 = T2 時，ΔTm 會相應減小。在忽略熱傳導的情況下通過公式（2）可以推導使得 T1 = T2 的 M 值計算式如式（11）所示：
|
在最高溫度達到 65℃、T1 = T2 的條件下，M 隨 Δt 的變化曲線如圖 4 左上圖所示，其中黑色點線為|T1-T2| ≤ 0.05℃ 的數值仿真 M 值，灰色點線為利用式（11）的計算值，由該圖可知仿真值與計算值幾乎一致。調節 M 值后，ΔTm 隨 Δz 的變化散點圖及擬合曲線如圖 4 右上圖所示，其中 ΔTm 隨 Δz 的擬合變化關系式，如式（12）所示（擬合度
= 0.94）：
|
若將 ΔTm ≤ 1℃ 作為平臺式溫度分布的條件，則根據式（12）計算得 0 < Δ z ≤ 8.0 mm，再對照圖 4 下圖可知，要調控出平臺式均勻溫度分布的 Δt 取值范圍為 50 ns ≤ Δt ≤ 150 ns。如圖 5 所示為在該條件范圍內形成的平臺式溫度分布圖，可治療焦域體積為 48.3～33.6 mm3。
圖5
聲軸方向均勻溫度分布圖（L = 10.0 mm，Pw = 41.22 W）
Figure5.
The uniform temperature distribution along acoustic axis (L = 10.0 mm, Pw = 41.22 W)
圖4
65℃ 時 M、ΔTm 及 Δz 的變化曲線（L = 10.0 mm，Pw = 41.22 W）
Figure4.
The variation of M, ΔTm and Δzwhen Δt was changed with the maximum temperature reached 65℃(L = 10.0 mm, Pw = 41.22 W)
3.1.3 L 對焦域溫度分布的影響
將 F1 固定，F2 向 F1 相對于靠近換能器一側調節時，F1 與 F2 之間的距離 L 設定為 Lnear，F2 位于 F1 遠離換能器一側時設定為 Lfar。在 Pw = 41.22 W 和 M = 0.50 的條件下，不同 L 時最高溫度達到 65℃ 所需輻照時間 t 隨 Δt 的變化曲線如圖 6 所示，其中'o'表示形成一個連續焦域，'×'表示形成兩個焦域。
圖6
不同 L 時輻照時間隨 Δt 的變化曲線（Pw = 41.22 W，M = 0.50；‘o’表示形成一個 55℃ 以上焦域，‘×’形成兩個焦域）
Figure6.
The variation of irradiation time with Δt under different L (Pw = 41.22 W, M = 0.50; ‘o’ stand for one continue focal region upon 55℃ was formed and ‘×’stand for two focal regions were formed)
由圖 6 可知，當 L ≤ 12.5 mm 時，通過 Δt 的調控可實現一個焦域的聚焦，且 L 越小可實現一個焦域聚焦的 Δt 的范圍越大；當 L ≥ 15.0 mm 時，通過調節 Δt 無法實現一個焦域聚焦。
在不同 L 的條件下，Pw = 41.22 W、焦域最高溫度達到 65℃ 時，溫度分布中兩溫度峰值間的距離 Δz 隨 Δt 的變化曲線如圖 7 所示。
圖7
不同 L 時 Δz 隨 Δt 的變化曲線圖（Pw = 41.22 W）
Figure7.
The variation of Δz with Δt under different L(Pw = 41.22 W)
由圖 7 可知，當 7.5 mm ≤ L ≤ 12.5 mm 時，通過調節 Δt 可以使得 0 < Δ z ≤ 8 mm，在此條件下調節 M 可以形成 T1 = T2 且 ΔTm ≤ 1℃ 的平臺式溫度分布。如表 2 所示為調控形成平臺式溫度分布的 Δt、M 及可治療焦域體積。由表 2 可知，相同的 L 條件下，通過 Δt 和 M 的調節，可在設定擬聚焦的目標位置實現平臺溫度分布的焦域；不同的 L 條件下，形成的平臺式溫度分布可治療體積也不同，L 越大形成的可知治療焦域體積越大。L 在 7.5～12.5 mm 的范圍內平臺式溫度分布可治療焦域的體積調控范圍為 30.09～71.56 mm。
表2
不同 L 時得到沿聲軸方向平臺溫度分布的參數及可治療焦域體積
Table2.
Parameters for generating the uniform temperature field along acoustic axis with different L and the volume of treatable focal region
3.2 垂直于聲軸方向的溫度分布調控
3.2.1 Δt 對溫度分布的影響
如圖 1 所示，在聲軸兩側設定 F1（0，1.00，75.00）和 F2（0，–1.00，75.00）（L = 2.0 mm），當最高溫度達到 65℃ 時，形成的溫度分布（Pw = 41.22 W，M = 0.50）如圖 8 所示。
圖8
Δt = 0～1 400 ns、最高溫度為 65℃ 時溫度分布圖（L = 2.0 mm，M = 0.50，Pw = 41.22 W）
Figure8.
Δt = 0～1 400 ns, the temperature distribution when the maximum focal temperature reached 65℃ (L = 2.0 mm,M = 0.50, Pw = 41.22 W)
由圖 8 可知，隨著 Δt 由 0 增大到 1 400 ns，55℃ 以上焦域由一個連續區域逐漸分離為兩個區域，隨后又逐漸融合為一個連續區域。
3.2.2 垂直于聲軸方向平臺溫度分布的調控
以圖 8 中 Δt = 400 ns 的溫度分布為例，當調控 M 時沒有形成沿 y 軸方向的均勻溫度分布。在這里當設定 F1（0，1.00，75.00）、F2（0，–1.00，77.50），調控 M = 0.46，t = 10.9 s 時形成最高溫度 65℃ 的溫度分布及 55℃ 以上區域如圖 9 的上半部分所示，其中左圖為 y-z 平面溫度分布圖，右圖的白色區域分別為 55℃ 以上焦域的 y-z 及 x-z 截面圖。由圖 9 可知，y 方向上形成了 ΔTm = 0.99℃ 的平臺溫度分布；55℃ 以上焦域在 y-z 平面上為 8 × 4 mm （底 × 高）的平行四邊形，x-z 平面上為 8 × 2 mm （長軸 × 短軸）的橢圓形。而在 F1 固定不變，F2 設定為（0，–1.00，80.00）的條件下，M = 0.47 和輻照時間 t = 23.0 s 時形成溫度場分布及 55℃ 以上焦域如圖 9 的下半部分所示，由圖可知，y 方向上溫度分布均勻（ΔTm = 0.61℃），55℃ 以上焦域在 y-z 平面上近似 10 × 5 mm （長 × 寬）的四邊形，x-z 平面上為 10 × 3 mm （長軸 × 短軸）的橢圓形。
圖9
y 軸方向均勻溫度分布圖（L = 2.0 mm，Pw = 41.22 W）
Figure9.
The uniform temperature distribution on y direction (L = 2.0 mm, Pw = 41.22 W)
3.2.3 L 對平臺溫度分布焦域的調控
在垂直聲軸方向設定不同距離 L 的目標焦點，即當 L = 2.0～3.0 mm 時，通過調節目標焦點相對位置及 Δt 和 M 的值，實現了 y 軸方向均勻的溫度分布，調控參數如表 3 所示。
表3
不同 L 時得到 y 軸方向均勻溫度分布的參數及可治療焦域大小
Table3.
Parameters and size of treatable focal region along with uniform temperature field on y direction under different L
由表 3 可知，均勻溫度分布區域的溫度差 ΔTm 均小于 1℃，可治療焦域在 y 方向的短軸長度可實現 4.25～6.25 mm 內的調控，可治療焦域的體積可在 58.7～201.8 mm3內調控。
4 討論與分析
本研究在靶區設定兩個間隔一定的擬聚焦目標點，采用時間反轉法獲取其激勵信號并進行疊加激勵換能器陣元，然后通過改變觸發時間延遲和幅值來調節焦域的溫度分布。時間反轉方法是目前 HIFU 顱腦聚焦較普遍采用的方法，其能夠通過相位矯正在目標區域形成焦域[17-19]。本研究調控過程中，調控兩激勵信號的觸發時間延遲會改變焦域溫度分布，但最終形成焦域是圍繞兩個設定目標焦點位置的區域；調控激勵信號的幅值比例系數不會改變焦域位置。為實現焦域內均勻的溫度分布，本研究采用雙激勵信號疊加激勵換能器的方法，不對激勵信號進行調控，得到的焦域內會有一個或兩個能量集中的溫度峰值；通過調控激勵信號，得到了溫度起伏 ΔTm 不超過 1℃ 的溫度分布均勻的焦域。本研究以 ΔTm 不超過 1℃ 作為焦域溫度均勻分布的標準；在實際治療過程中，需要根據治療需求，對焦域內溫度起伏 ΔTm 進行調節。在傳統的 HIFU 單焦點輻照過程中，往往通過改變輻照功率或輻照時間改變焦域的溫度分布和體積，這樣易引發淺部正常組織的傷害[21-22]。本研究在輻照功率、輻照時間一定的前提下實現了焦域溫度分布的調控，并形成了沿聲軸和垂直聲軸的均勻溫度分布，同時可調控治療焦域長、短軸長度，從而改變治療焦域體積大小和形狀。
上述仿真結果表明，調控得到垂直聲軸方向均勻溫度分布的 L 和 Δt 的范圍均較小，這是由于 HIFU 聚焦形成的焦域垂直聲軸方向長度（短軸）較短。本研究通過兩個激勵信號疊加并調控激勵信號間時間延遲和幅值調控系數，實現了單次輻照焦域均勻溫度分布，與 Fan 等[23]提出的理想 HIFU 治療焦域相一致；且在實現焦域溫度均勻分布的同時，可根據臨床實際的需要調控焦域體積、形狀和焦域內溫度分布。當治療較大體積腫瘤時，需要將多次輻照形成損傷進行堆疊，對輻照路徑及間隔時間的優化有待于進一步深入研究。
本研究采用一名 46 歲男性健康志愿者頭部 CT 數據建立三維數值仿真模型，研究 HIFU 腦腫瘤治療焦域溫度分布的調控，該模型與實際腦腫瘤治療有一定差異，但腦腫瘤組織與正常腦組織的差異性非常小，幾乎不會對 HIFU 焦域溫度場產生影響[24]。而在制定 HIFU 腦腫瘤治療方案時，需要準確定位腦腫瘤的位置、大小以及形狀，并根據這些參數制定安全、有效的治療計劃，這部分研究正在進行之中。此外，本研究通過調控雙激勵信號的觸發延遲時間和幅值比例系數實現對焦域溫度分布的調控，其中信號間的觸發延遲時間為納秒（ns）級，這期間組織溫升很小，不會影響組織的物理及生物特性；在整個治療過程中，靶區組織會產生變性，考慮組織變性后聲學及熱力學影響的研究正在進行中。HIFU 腦部疾病治療的理想方式是無需開顱術的 HIFU 經顱治療。目前，在本研究結果及方法的基礎上，考慮顱骨的 HIFU 經顱聚焦焦域溫度分布的調控及相應的實驗驗證研究正在進行之中。
