張增鵬，陶寅，劉紅冬，陳志，徐榭，4，裴曦，4

1.中國科學(xué)技術(shù)大學(xué)核科學(xué)技術(shù)學(xué)院，安徽合肥，230025；2.中國科學(xué)技術(shù)大學(xué)核醫(yī)學(xué)物理研究所，安徽合肥230025；3.中山大學(xué)附屬腫瘤醫(yī)院放療科，廣東廣州510060；4.安徽慧軟科技有限公司，安徽合肥230088

前言

相比光子，質(zhì)子放射治療具有獨(dú)特的Bragg 峰，在靶區(qū)給予更高劑量的同時，降低對腫瘤后端器官的輻射損傷［1-2］，因此在全球范圍高速發(fā)展。根據(jù)國際粒子治療協(xié)作組（Particle Therapy Co-Operative Group,PTCOG）的統(tǒng)計，截止2020年3月份，全球正在運(yùn)營的質(zhì)子治療中心已達(dá)91 個，我國正在運(yùn)營的質(zhì)子放療中心僅有2家，仍有7家正在建設(shè)中［3-4］?？梢钥吹?，隨著越來越多質(zhì)子中心的投入使用，質(zhì)子治療在放射治療中的比重會逐漸增大［5］。與此同時，與質(zhì)子治療相配套的治療計劃系統(tǒng)也會越來越受到重視。目前，質(zhì)子治療計劃系統(tǒng)幾乎被國外各大廠商壟斷，如瑞典RaySearch 的RayStation，美國Varian 的Eclipse，瑞典Elekta 的Xio 和荷蘭Philips 的Pinnacle等［6-7］。我國質(zhì)子放療計劃研究起步較晚，國內(nèi)擁有自主知識產(chǎn)權(quán)的放療系統(tǒng)主要有中國科學(xué)技術(shù)大學(xué)放射醫(yī)學(xué)物理中心研發(fā)的DeepPlan 和中國科學(xué)院核能安全技術(shù)研究所研發(fā)的KylinRay［8-9］。兩者質(zhì)子模塊均采用筆形束算法，通過質(zhì)子束在體模中的縱向積分深度劑量和橫線束斑展寬計算體模中各體素的劑量沉積。目前由于缺乏質(zhì)子放療臨床數(shù)據(jù)，國產(chǎn)質(zhì)子放療計劃系統(tǒng)開展的臨床數(shù)據(jù)對比實驗較少。

本研究旨在研究DeepPlan 質(zhì)子模塊中建立佛羅里達(dá)大學(xué)質(zhì)子放療中心（University of Florida Health Proton Therapy Institute，UFHPTI）質(zhì)子加速器筆形束模型，通過對比UFHPTI質(zhì)子加速器在筆形束掃描模式下的臨床實驗數(shù)據(jù)，驗證模型構(gòu)建的準(zhǔn)確性并與UFHPTI配套的商用放療計劃軟件RayStation就前列腺癌臨床數(shù)據(jù)進(jìn)行劑量對比。

1 材料與方法

1.1 DeepPlan質(zhì)子模塊

1.1.1 模塊簡介DeepPlan 是由中國科學(xué)技術(shù)大學(xué)放射醫(yī)學(xué)物理中心自主研發(fā)的放療計劃系統(tǒng)［9-12］。其質(zhì)子模塊采用筆形束算法，公式如下：

其中，D（x,y,z,x0,y0,E）表示能量為E，中心坐標(biāo)為（x0,y0）的質(zhì)子束在模體位置（x,y,z）的劑量沉積，IDD是質(zhì)子束在縱向的積分深度劑量，利用蒙卡軟件TOPAS模擬得到，w表示核散射權(quán)重，Gmcs表示多庫倫散射的高斯分布，Gnuc表示核散射的高斯分布，σmcs是能量為E的質(zhì)子束由于多庫倫散射在x軸或y軸方向上高斯分布的標(biāo)準(zhǔn)差。σnuc是能量為E的質(zhì)子束由于核散射在x軸或y軸方向上高斯分布的標(biāo)準(zhǔn)差，這些參數(shù)根據(jù)Soukup等［13］、Gottschalk等［14-15］的工作計算得到。

DeepPlan 質(zhì)子模塊的工作流程為：首先通過數(shù)據(jù)讀取得到模體信息和質(zhì)子束信息，包括模體位置、大小、體素材料、等中心點位置、質(zhì)子能量、入射方向、初始束斑大小和質(zhì)子束權(quán)重等，并導(dǎo)入計算所需要的數(shù)據(jù)庫，如上述計算所需要的IDD 數(shù)據(jù)庫、多庫倫散射數(shù)據(jù)庫、核散射數(shù)據(jù)庫和核散射權(quán)重數(shù)據(jù)庫等。接著光線追蹤通過設(shè)定兩個位置點P1、P2計算質(zhì)子束穿過體素的幾何距離和等效水深度。最后劑量計算根據(jù)等效水深度得到計算所需的數(shù)據(jù)參數(shù)，對每個質(zhì)子束利用公式（1）～（3）得到體模的總劑量。這3個部分實現(xiàn)筆形束掃描下質(zhì)子劑量的快速計算。

1.1.2 模型構(gòu)建

1.1.2.1 數(shù)據(jù)讀取數(shù)據(jù)庫構(gòu)建的準(zhǔn)確性是劑量計算的前提條件，比如IDD 數(shù)據(jù)庫構(gòu)建時不同初始能散和角度發(fā)散質(zhì)子束的IDD 是不一樣的。圖1是用蒙卡軟件TOPAS模擬的200 MeV質(zhì)子在不同能散下的IDD數(shù)據(jù)曲線。由圖可知，相同能量不同能散質(zhì)子束的IDD 曲線在Bragg 峰前端有一定差距，能散越小，IDD 曲線越低，這樣根據(jù)式（1）計算的劑量會比實際情況偏低。本研究結(jié)合UFHPTI 質(zhì)子加速器信息并結(jié)合劉紅冬等［16］的研究，插值計算得到加速器輸出各種能量的參數(shù)信息，構(gòu)建包括IDD在內(nèi)的數(shù)據(jù)庫。

圖1 200 MeV質(zhì)子束在不同能散的積分深度劑量Fig.1 Integrated depth dose of 200 MeV proton beam in different energy spreads

1.1.2.2 光線追蹤光線追蹤是劑量計算的重要一環(huán)，算法程序根據(jù)Siddon［17］、Han 等［18］的工作實現(xiàn)。該算法認(rèn)為構(gòu)成體模的每一個體素不是獨(dú)立的，而是一系列等距平行平面的正交集焦點體積，加上質(zhì)子束穿過體模時的兩個位置點信息，可實現(xiàn)位置點連線上各體素幾何距離和等效水深度的快速計算。需要注意的是在實際筆形束掃描模型中，質(zhì)子加速器通過X 和Y 方向掃描磁鐵實現(xiàn)對質(zhì)子束位置的精確控制，因而有X 和Y 兩個焦點，如圖2所示。以UFHPTI質(zhì)子加速器為例，其X、Y 兩個焦點與等中心平面的距離分別為2 329 和1 934 mm。因此，在設(shè)置兩個位置點進(jìn)行光線追蹤時必須考慮X、Y 方向之間的差異。

1.1.2.3 劑量計算劑量計算過程需要考慮的因素較多，包括空氣處理和散射修正等，公式如下：

圖2 筆形束掃描模型虛擬源示意圖Fig.2 Diagram of virtual source in pencil beam scanning model

其中，式（4）用于在空氣中進(jìn)行劑量計算，主要是因為在空氣中等效水深度為0，且在空氣中只有多庫倫散射沒有核散射。式（5）用于在非空氣體模中進(jìn)行劑量計算，與式（1）相比，該公式增加了多庫倫散射修正項δmcs和核散射修正項δnuc，表達(dá)式如上所示。在（8）～（9）式中多庫倫散射修正因子cmcs和核散射修正因子cnuc分別是與多庫倫散射標(biāo)準(zhǔn)差σmcs和核散射標(biāo)準(zhǔn)差σnuc有關(guān)的常數(shù)項。

圖3是劑量計算示意圖，在參數(shù)導(dǎo)入之后首先應(yīng)該判斷各體素是否為空氣，判別方法為計算各體素的質(zhì)量組織本領(lǐng)之和T，若T＜10-6，則認(rèn)為在空氣中，采用式（4）進(jìn)行劑量計算；否則，采用式（5）進(jìn)行劑量計算。

圖3 劑量計算示意圖Fig.3 Diagram of dose calculation

1.1.3 測試環(huán)境本文完成DeepPlan質(zhì)子模塊臨床測試的計算機(jī)配置如下：操作系統(tǒng)為Windows10 專業(yè)版（64 位），CPU 為Intel(R)CoreTMi7-87100K、六核十二線程、16 G內(nèi)存，GPU 為NVIDIA GeForce GTX10606GB。DeepPlan 質(zhì)子模塊有c++和c#兩個版本，本研究基于c#版本在軟件Visual Studio 2017上運(yùn)行。

1.2 UFHPTI質(zhì)子加速器臨床實驗數(shù)據(jù)測量

1.2.1 IDD 數(shù)據(jù)測量UFHPTI 使用IBA 公司的BPC探測器（Bragg Peak Chamber,model 34070,PTW,Freiburg,Germany）在IBA 水箱中進(jìn)行積分深度的測量，共有30 組質(zhì)子能量，分別為：226.35、223.91、219.06、214.12、209.16、204.54、199.63、195.1、190.18、185.74、180.9、176.21、171.89、167.79、163.02、158.09、153.04、149.07、144.81、140.63、136.07、131.51、127.38、122.46、118.69、115.09、110.04、106.57、102.36、99.59 MeV 等。針對測試數(shù)據(jù)，在DeepPlan質(zhì)子模塊中計算上述30組能量的質(zhì)子束在水模中的劑量沉積，體模大小為100 mm×100 mm×400 mm，每個體素的大小為1 mm×1 mm×1 mm，對垂直于質(zhì)子束方向的每個平面進(jìn)行劑量求和，利用最大值進(jìn)行歸一即可得到在縱向的積分深度劑量曲線。

1.2.2 空氣中束斑尺寸測量UFHPTI 使用IBA 公司的Lynx PT 閃爍探測器進(jìn)行束斑尺寸的測量。30 組能量在上面已經(jīng)提過，5 個測量位置分別為距離等中心平面20、10、0、-10、-20 cm。為了將二維數(shù)據(jù)轉(zhuǎn)成一維，需要將測試數(shù)據(jù)導(dǎo)入IBA的MyQA軟件提取X或Y方向上的橫向劑量分布，并做歸一處理。

針對測試數(shù)據(jù)，在DeepPlan 質(zhì)子模塊中計算30組能量的質(zhì)子束在體模為空氣中的劑量沉積，平面大小為40 mm×40 mm，平面精度和Lynx PT 的精度保持一致，為0.25 mm×0.25 mm，將計算結(jié)果在X 或Y方向求和得到一維數(shù)據(jù)，并做歸一處理。

1.2.3 SOBP 和橫向劑量的測量為保證劑量分布在縱向和橫向上與靶區(qū)保持適形，需要使用不同能量、不同權(quán)重、不同位置的質(zhì)子束。在SOBP 測量時UFHPTI 使用IBA 公司的多層電離室Zebra 測量中心軸上的絕對劑量；在橫向劑量測量時UFHPTI 使用Lynx PT 閃爍探測器進(jìn)行測量。針對測量數(shù)據(jù)，在DeepPlan 質(zhì)子劑量模塊根據(jù)提供的質(zhì)子束能量信息、權(quán)重信息和位置信息等，在水模中進(jìn)行劑量計算，水模大小為200 mm×200 mm×400 mm，體素大小為2 mm×2 mm×2 mm。

1.2.4 前列腺癌的劑量計算RayStation 是由RaySearch 公司于2008年開始研發(fā)，2011年市場發(fā)布的治療計劃系統(tǒng)［19］。RayStation 的技術(shù)主要有輪廓的勾畫、形變配準(zhǔn)及定量評估、自適應(yīng)治療、劑量追蹤、多標(biāo)準(zhǔn)優(yōu)化和后備計劃等。作為與UFHPTI質(zhì)子加速器相配套的治療計劃系統(tǒng)，RayStation 在筆形束模型下的質(zhì)子劑量計算有兩種算法，一種是筆形束解析算法，另一種是蒙卡算法。筆形束算法比蒙卡算法計算速度快，但準(zhǔn)確性低于蒙卡算法：對于大量非均勻介質(zhì)體模的計算和所有質(zhì)子束使用射程移位器的情況等，主要以蒙卡算法為主。本文以RayStation 蒙卡算法的計算結(jié)果進(jìn)行對比。在UFHPTI 的兩個前列腺癌臨床放療計劃中，相關(guān)信息如下：病例1的體模大小為512 mm×512 mm×361 mm，體素大小為1.171 9 mm×1.171 9 mm×1.000 0 mm，270°和90°方向的質(zhì)子束總數(shù)分別為1 039 和977；病例2的體模大小為512 mm×512 mm×293 mm，體素大小為1.171 9 mm×1.171 9 mm×1.000 0 mm，270°和90°方向的質(zhì)子束總數(shù)分別為1 307 和1 303。由于RayStation 中給出的三維劑量信息體素大小為2 mm×2 mm×2 mm，為了進(jìn)行g(shù)amma 通過率對比，需要對病例1 和2 的體模按相同體素大小進(jìn)行重采樣。最后使用PTW 公司的VeriSoft 軟件對二者的計算結(jié)果進(jìn)行g(shù)amma分析。

2 結(jié)果與討論

2.1 IDD結(jié)果比較

為展示方便，圖4所示將30 組IDD 數(shù)據(jù)從高到低分為6 組（a～f 組），每組5 個能量進(jìn)行比較，每組能量的IDD 曲線都利用該能量下最大值進(jìn)行歸一，深度方向進(jìn)行適當(dāng)截取，DeepPlan 質(zhì)子模塊計算結(jié)果用DeepPlan 表示，UFHPTI質(zhì)子加速器臨床實驗數(shù)據(jù)用Measure 表示。從圖中可以看到，DeepPlan 質(zhì)子模塊計算結(jié)果與UFHPTI 質(zhì)子加速器臨床實驗數(shù)據(jù)在Bragg 峰前后曲線相吻合。兩者平均相對誤差為0.01%，少數(shù)能量如圖4f 組中的99.59、102.36 MeV等，DeepPlan 計算結(jié)果與實驗數(shù)據(jù)在最大值位置附近出現(xiàn)一個0.23%的相對誤差。出現(xiàn)這個情況的原因可能是測試數(shù)據(jù)最小深度間隔為3 mm，精度較低，導(dǎo)致在測量峰值最大位置時出現(xiàn)一定誤差。

圖4 30組不同能量下IDD結(jié)果對比Fig.4 Comparison of integrated depth dose under different energies

2.2 空氣中束斑尺寸結(jié)果對比

如圖5所示，以能量為199.63 MeV 的質(zhì)子束為例，藍(lán)色曲線表示UFHPTI質(zhì)子加速器的臨床測試結(jié)果，紅色曲線表示DeepPlan 質(zhì)子模塊的計算結(jié)果，圖中包含5個不同位置的對比結(jié)果，束斑寬度由小到大分別為距離等中心平面20、10、0、-10、-20 cm 處的結(jié)果。對于某一位置的計算結(jié)果和測量數(shù)據(jù)，導(dǎo)入MATLAB R2018b 版本中，使用軟件自帶的高斯擬合函數(shù)求出束斑大小。圖6是不同能量質(zhì)子束在等中心平面束斑大小的對比結(jié)果，大部分相對誤差在0.5%以內(nèi)，只有163.02 MeV 比實際測試略低一些，相對誤差為1.13%。除此之外，表1展示距離等中心平面20、10、-10 和-20 cm 處的計算結(jié)果。從表中可知，DeepPlan 質(zhì)子模塊的計算結(jié)果與臨床實驗數(shù)據(jù)大部分相對誤差在0.3%以內(nèi)，平均相對誤差為0.15%?？梢钥吹?，由于質(zhì)子束在空氣和體模中的輸運(yùn)不同，在劑量計算時將二者分開可使結(jié)果更加準(zhǔn)確。

圖5 199.63 MeV質(zhì)子束在不同位置結(jié)果對比Fig.5 Comparison of 199.63 MeV proton beams at different positions

圖6 不同能量質(zhì)子束在等中心平面的束斑尺寸Fig.6 Spot size of proton beams with different energies in the isocenter plane

表1 30組不同能量質(zhì)子束在空氣不同位置的束斑比較Tab.1 Comparison of the spot size of 30 groups of proton beams with different energies at different positions in the air

2.3 SOBP和橫向劑量結(jié)果對比

圖7是縱向和橫向的結(jié)果對比，黑色曲線表示UFHPTI 質(zhì)子加速器的臨床測量結(jié)果，紅色曲線表示DeepPlan 質(zhì)子模塊的計算結(jié)果。由圖7a 在縱向擴(kuò)展布拉格峰的對比中，可以看到實現(xiàn)射程為200 mm，展寬100 mm 目的。DeepPlan 劑量計算結(jié)果與臨床實驗數(shù)據(jù)平均相對誤差1.07%，最大相對誤差3.91%。圖7b 展示的是等中心平面橫向劑量分布的對比，等中心平面在距離水模表面15 cm 的位置，可以看到平坦區(qū)域達(dá)到100 mm 目的。DeepPlan 劑量計算結(jié)果與臨床實驗數(shù)據(jù)平均相對誤差1.92%，最大相對誤差4.09%。

圖7 縱向和橫向結(jié)果對比Fig.7 Comparison of longitudinal and transverse results

2.4 前列腺癌結(jié)果對比

DeepPlan與RayStation的結(jié)果在該平面二維gamma通過率達(dá)到了97.1%，但在末端右下角未通過，出現(xiàn)一些冷點和熱點，見圖8?？赡茉蚴窍啾菵eepPlan解析算法，RayStation采用蒙卡算法，對非均勻介質(zhì)部分計算結(jié)果更好，這是后續(xù)工作需要改進(jìn)的地方。表2具體展示了兩個前列腺癌病例在兩個射野方向的三維gamma通過率，測試標(biāo)準(zhǔn)為3 mm/3%。由表2可知每個射野方向在大于10%處方劑量區(qū)域里的gamma通過率都達(dá)到95%以上，其中病例1在兩個方向的gamma通過率分別為96.4%和97.5%，病例2在兩個方向的gamma通過率分別為99.3%和98.9%。

圖8 DeepPlan與RayStation關(guān)于前列腺癌的二維劑量結(jié)果比較Fig.8 Comparison of two-dimensional dose calculation results of DeepPlan and RayStation on prostate cancer

表2 DeepPlan與RayStation關(guān)于前列腺癌的gamma分析結(jié)果Tab.2 Gamma analysis results of DeepPlan and RayStation on prostate cancer

3 結(jié)語

本文對DeepPlan 質(zhì)子模塊在數(shù)據(jù)讀取、光線追蹤和劑量計算等方面建模，通過與UFHPTI質(zhì)子加速器在30組IDD對比、30組空氣中質(zhì)子束斑發(fā)散對比、SOBP 和橫向劑量分布對比，驗證DeepPlan 質(zhì)子模塊模型構(gòu)建的準(zhǔn)確性。同時以兩個前列腺癌的放療計劃為實例，與商用放療計劃系統(tǒng)RayStation 進(jìn)行g(shù)amma 分析，表明DeepPlan 質(zhì)子模塊在筆形束掃描模式下可以初步應(yīng)用于臨床前列腺癌的劑量計算。后續(xù)還將繼續(xù)對治療計劃系統(tǒng)DeepPlan 質(zhì)子模塊進(jìn)行改進(jìn)，包括對非均勻介質(zhì)的修正和加入射程移位器的劑量計算等，對更多類型的病例進(jìn)行測試分析，進(jìn)一步增強(qiáng)DeepPlan 質(zhì)子模塊在臨床應(yīng)用的可行性。

感謝UFHPTI提供的臨床實驗數(shù)據(jù)。