段歡歡,李書舟,曹瑛,唐杜,雷明軍,楊振,邱小平
1.南華大學核科學技術學院,湖南衡陽421001;2.中南大學湘雅醫院腫瘤科,湖南長沙410008
相比于常規的適形放療,調強放射治療(Intensity-Modulated Radiation Therapy,IMRT)的最大優勢是通過強度調制實現了靶區內高適形度的劑量分布,同時降低了周圍正常組織的受照量[1]。但是,IMRT也增加了治療過程中計劃劑量傳遞的不確定性,因此,必須制定全面的質量保證(Quality Assurance,QA)和質量控制(Quality Control,QC)程序以評估治療實施的可靠性和安全性[2-4]。目前,傳統的計劃QA 方法主要是在患者治療實施前對計劃進行二維或三維的劑量分布測量,并通過γ 分析方法將其與計劃系統計算的劑量分布進行比較來判斷計劃是否可行[5]。然而,這些測量過程需要耗費大量的時間和人力成本,對于那些沒能通過QA 測量的計劃,也很難確定其失敗的具體原因并加以糾正,這可能會延緩患者的放療進程[6]。因此,臨床上需要一種更精簡、更智能、能預先識別那些不能通過QA 測量計劃的方法。在本研究中,筆者基于深度學習方法建立卷積神經網絡(Convolution Neural Network,CNN)模型來預測腦膠質瘤患者IMRT計劃單個射野的γ 通過率(Gamma Passing Rates,GPR),以研究射野的劑量分布信息用于預測GPR 的可行性并探討改善模型預測精度的方法。這種虛擬QA 的方法將有助于物理師及時優化QA 測量潛在失敗的計劃,有效減少特定患者QA 的工作量并節約臨床資源[7]。
選取2019年1月~2020年5月在本中心接受放療的48例腦膠質瘤患者的IMRT計劃,計劃設計過程都在Eclipse version 13.7(Varian Medical Systems, Palo Alto,CA)治療計劃系統完成,計劃分別采用了5野、6野、7野IMRT照射,共計260個射野,X射線束能量為6 MV,劑量計算采用AAA算法,計算網格為2.5 mm,患者治療實施以及計劃驗證過程都在Varian 23EX直線加速器上進行。
研究中的劑量驗證設備是電子射野影像系統(Electronic Portal Imaging Device,EPID)非晶硅平板探測器,其在40 cm×30 cm 的有效探測平面上分布了1 024×768個探測器,分辨率可達0.392 mm/pixel。在加速器上執行患者的IMRT QA 驗證計劃,并利用EPID 測得射野影像,其灰度值與接收到的劑量直接相關,因此經過校準和刻度后可用于患者治療前的劑量驗證。Eclipse 工作站中的portal dosimetry 劑量測定軟件包會根據患者的計劃及CT 信息推算出射野的劑量分布,并通過γ 分析方法將其與EPID 實測的射野影像進行比較,選取2%(global)/2 mm 標準和10%的最大劑量閾值進行γ 評估并得到每個射野的GPR。為了減小GPR 受加速器性能及其日常輸出波動的影響,所有的計劃驗證工作都在一天之內完成。
每個射野基于portal dosimetry系統計算的劑量分布圖、2%(global)/2 mm 標準以及10%閾值下的GPR將分別作為模型的輸入和輸出。為了減少CNN模型的計算量并提取到更有利的特征,將所有射野的劑量分布圖都減少至227×227個像素,像素值都用對應射野的最大像素值進行標準化,并將小于最大值10%的像素點全部設置為0,像素分辨率為1 mm×1 mm。此外,為了減小EPID 射野影像灰度響應值受射野大小和入射角度的影響,劑量驗證測量時需要保證每個射野的機器跳數(Monitor Unit, MU)都被投射在EPID 平面上,并且機架角度和準直器角度統一設置為0,這也是為了排除多葉準直器(Multileaf Collimators,MLC)重力因素對測量結果的影響。
AlexNet[8]是CNN 領域比較有標志性的一個網絡模型,本文中的CNN 模型則是在AlexNet 的基礎上進行了改動,其結構如圖1所示,我們使用卷積核分別為7×7和5×5的兩個卷積層替換了卷積核為11×11的第一個卷積層,并將步長由4改為2,以及使用兩個卷積核為3×3 的卷積層替換了卷積核為5×5 的第二個卷積層,通過使用更小的濾波器可以實現更深的網絡構造,并得到更好的非線性和降低網絡的權值。此外,使用Relu 作為激活函數可以增加神經網絡模型的非線性,L2 正則化用來防止模型過擬合并實現對模型空間的限制,3 個最大池化層不僅能夠降低特征圖的維度,還可以忽略劑量分布圖傾斜、旋轉等相對位置的變化,有助于提高模型的預測精度。此外,每一個卷積層和全連接層后面都引入了批量標準化(Batch Standardization, BN)以及在全連接層之間使用數據丟失(dropout)技術,可以保證模型有一個更穩定的輸出并在一定程度上防止模型過擬合。最后的全連接層用來對前面所設計的特征做加權和并輸出模型的預測值。

圖1 CNN模型結構示意圖Fig.1 Structure diagram of convolution neural network(CNN)model
傳統機器學習階段(數據集在萬這個數量級),數據集的劃分比例一般為6:2:2,但由于本文的數據集只有260個樣本,為了避免訓練集過少導致模型欠擬合,具有較差的泛化能力,筆者將數據集按照8:1:1的比例進行訓練集、驗證集和測試集的劃分,即208 個射野作為訓練集供模型訓練,26 個射野作為驗證集用于模型測試和選擇,另外26 個射野則作為測試集來評估模型的泛化性。數據集劃分時采用均勻隨機抽樣的方式來保證它們都處于同一分布中,以減少數據特異性對模型預測效果的影響,此外,筆者還將訓練集中每個射野的劑量分布圖沿水平、豎直方向翻轉,以及旋轉10°來進行數據增強,達到提升模型訓練效果的目的。
訓練過程中,最大迭代次數epoch 設置為1 000,小批量訓練樣本量batch_size 設置為16,即每次學習迭代都從16個隨機選定的樣本數據開始訓練。選用Adam 算法來優化代價函數均方根誤差,是因為Adam 算法在基于梯度的優化問題中占用內存較少、計算效率高,適用于機器學習中的許多優化問題。每個卷積層和全連接層都使用Xavier 初始化器對權重進行賦值,并使用L2 正則化技術,正則化系數為0.05,全連接層之間的dropout參數為0.6,這些參數的合理設置有助于減少模型過擬合的風險。為了實現更穩定的訓練,將初始學習率設置為0.003,并在每兩個epoch 后將學習率按照指數衰減來動態調節學習率,衰減率decay_rate 為0.98。使用平均絕對誤差(Mean Absolute Error,MAE)作為損失函數來評估模型的預測表現并進行參數微調和模型選擇,訓練過程中每次迭代都會改變CNN的參數、權重和偏差值,當驗證集loss 值在50 個epoch 內不再下降達到收斂時,停止網絡訓練以防止模型過擬合情況的發生,最佳的參數設定在驗證集上MAE最小時的epoch,訓練和測試過程都在一臺配置為NVIDIAGeForce GTX 860 Ti的GPU上進行。
如圖2所示,紅色和藍色的點分別表示驗證集和測試集上模型預測值和實際測量值的分布情況,紅線和藍線分別表示測量值比預測值大3%以及小3%的情況,黃線則表示預測值等于測量值的情況。如果模型預測精度是完美的,那么這些點都會出現在黃色直線附近,但在實際應用中,模型的預測值和實際測量值總會存在一定的偏差。圖3則將驗證集和測試集上模型預測和實際測量的GPR 進行了線性擬合,其擬合直線的斜率分別為1.056 和0.867,決定系數R2分別為0.91和0.81,表明模型預測值和實際測量值之間具有較好的線性關系。

圖2 GPR預測值和測量值散點圖Fig.2 Scatter plot of predicted and measured gamma passing rates(GPR)

圖3 GPR預測值和測量值擬合直線Fig.3 Linear fitting of predicted and measured GPR
圖4顯示了驗證集和測試集上每個樣本的GPR預測誤差;圖5則統計了GPR預測誤差在每個區間的樣本分布情況。可以看到,在驗證集和測試集上,96%樣本的GPR預測誤差都在±3%以內,且最大預測誤差分別為3.09%和3.54%,然而,大多數樣本的GPR預測值都大于實際測量值,這在測試集上表現的更加明顯,表明該模型可能具有高估射野GPR的傾向。

圖4 每個射野的GPR預測誤差Fig.4 GPR prediction error of each radiation field

圖5 GPR預測誤差的區間分布Fig.5 Interval distribution of GPR prediction errors
CNN 模型在驗證集和測試集上的預測結果如表1所示,驗證集和測試集上模型預測值和實際測量值的均方根誤差(Root Mean Square Error,RMSE)分別為1.31%和1.44%,MAE分別為0.99%和1.17%,皮爾遜相關性系數r分別為0.96 和0.90,而且其平均值、標準偏差和中位值也都比較接近,表明該模型具有較小的預測誤差和較好的泛化性能。

表1 驗證集、測試集上GPR預測值和測量值比較Tab.1 Comparison of predicted and measured GPR in validation set and test set
目前,一些研究已經利用機器學習方法預測了IMRT計劃的GPR。Valdes 等[9-10]基于既往498 例IMRT 計劃,使用Poisson 回歸和Lasso 正則算法預測了3%/3 mm標準下的γ通過率,并在之后的研究中通過多中心的驗證來評估模型的泛化性能。但在他們研究中,模型輸入的是一種低密度的二維二極管陣列探測器矩陣驗證結果,γ 分析標準為3%/3 mm,而該標準通常被認為對臨床相關誤差不敏感[11-13]。Lam 等[14]選取了182 例IMRT 計劃,采用3 種基于樹的機器學習算法(AdaBoost、Random Forest 和XGBoost)預測了EPID 驗證測量的GPR,他們的模型可以準確地預測2%/2 mm 標準下GPR,并且最大誤差小于4%,MAE 小于1%。然而,這種傳統的機器學習算法在數據分析整理階段,需要盡量選擇與預測指標相關性較高的特征參數并排除無關特征參數的干擾,這不僅增加了特征值選擇、提取和計算的難度,而且還可能遺漏掉一些重要的特征。
近年來,深度學習方法預測IMRT計劃的GPR也逐漸受到 人 們的關注。Interian 等[15]利用CNN 對IMRT 計劃的通量圖進行學習并預測了3%/3 mm 標準下的GPR,他們發現深度學習CNN 不需要使用專家設計的特征就可以達到與之前傳統機器學習算法相媲美的預測精度。Tomori等[16]選取了60例前列腺癌IMRT 計劃,并建立了一個結構相對簡單的CNN模型去預測4 種標準下的GPR,結果發現QA 模體上的劑量分布信息是預測GPR 的一種有用數據,而且2%/2 mm 標準下模型預測值與測量值之間具有更好的相關性。但相比其它γ分析標準,2%/2 mm 標準下的GPR 預測效果卻最差,其原因可能是他們預測了整個計劃的GPR,從而削弱了小野以及不規則野對GPR的影響[17-18]。
在本研究中,筆者使用IMRT 計劃每個射野的劑量分布圖作為輸入數據,并建立了一個新的CNN 模型去預測單個射野的GPR。結果表明該CNN模型成功地學習了射野劑量分布圖的特征并取得了較小的預測誤差。然而,模型的預測值和測量值之間的匹配并不是完美的,大部分樣本的GPR 預測值都大于實際測量值,且在測試集上表現的更加明顯,這很可能與數據集較小且樣本的均勻性較差有關,數據集中GPR 小于90%的樣本數量相對較少,導致較低GPR 的樣本獲得了一個較高的預測值。因此,想要進一步提高模型的預測精度,還需要獲取更多的樣本數量并保證樣本數據的均勻性。目前,僅從一家放療機構獲取大量的較低GPR 的計劃用于模型訓練十分困難,這可能需要多機構的協作才能實現[19]。
對于深度學習方法來說,數據集太小很容易造成模型過擬合,從而導致模型具有較差的泛化性能。筆者選擇預測單個射野的γ 通過率而不是整個計劃,這不僅可以擴大訓練樣本的數量,而且單個射野的劑量分布信息也可能含有更多影響劑量傳遞準確性的間接因素[20-21]。通過數據增廣來擴展訓練集,并使用數據丟失、批量標準化、L2正則化等技術來防止模型過擬合,這對實際的臨床工作也十分有利,因為積累大量的臨床數據是非常耗時的,獲取成百上千的計劃來供模型訓練不太切合實際。然而,深度學習CNN 相比傳統的機器學習方法卻更為復雜,并且具有較差的可解釋性[22]。因此,對訓練好的CNN 模型進行數據測試和仔細驗證,并通過各種分析方法來評估模型預測的準確性是至關重要的[23]。此外,這項工作也具有一定的局限性,由于本文的數據都來自于同一機構的同一種病例,并且都在同一臺機器上完成治療和QA 測量,所以研究中所建立的CNN模型可能不適用于其它放療機構。在未來的研究中,可以選取基于不同直線加速器和QA 測量設備進行驗證的計劃來重新建模,進一步提高模型的預測精度和泛化性能[24]。
目前,通過這種虛擬QA 的方法來預測IMRT 計劃的GPR 可能還不足以代替傳統基于測量的QA 方法,但物理師們可以根據虛擬QA 的結果提前識別那些不能通過QA 測量的計劃,以便采取更積極主動的方法來優化計劃[25-26]。而且,在對特定患者的IMRT計劃進行劑量驗證時,可以更專注于那些GPR 預測值較低的計劃,從而減輕物理師的工作負擔。此外,對于那些硬件設施還不是很完善的放療機構來說,虛擬QA 的結果結合簡單的劑量驗證方法也可能成為一種簡便的QA 手段,能進一步節約臨床資源并提高治療計劃實施的質量和效率。