基于記憶學(xué)習(xí)法的放療中呼吸運(yùn)動預(yù)測技術(shù)的研究

萬偉權(quán) 張慧連 徐子海 賀志強(qiáng) 陳超敏*

1(南方醫(yī)科大學(xué)生物醫(yī)學(xué)工程學(xué)院， 廣州 510515)2(廣東新華軟件外包有限公司， 廣州 510515)3(解放軍303醫(yī)院放射治療中心， 南寧 530021)

基于記憶學(xué)習(xí)法的放療中呼吸運(yùn)動預(yù)測技術(shù)的研究

萬偉權(quán)1張慧連2徐子海3賀志強(qiáng)1陳超敏1*

1(南方醫(yī)科大學(xué)生物醫(yī)學(xué)工程學(xué)院， 廣州 510515)2(廣東新華軟件外包有限公司， 廣州 510515)3(解放軍303醫(yī)院放射治療中心， 南寧 530021)

對胸腹部腫瘤進(jìn)行實(shí)時跟蹤放療時，需要通過預(yù)測來補(bǔ)償系統(tǒng)延遲。然而，由于呼吸運(yùn)動的復(fù)雜性，傳統(tǒng)方法難以滿足要求。本文應(yīng)用一種基于記憶學(xué)習(xí)法進(jìn)行呼吸預(yù)測，該方法首先存儲訓(xùn)練數(shù)據(jù)到記憶中，然后查找相關(guān)數(shù)據(jù)應(yīng)答當(dāng)前查詢。在此基礎(chǔ)上，采用“滑窗法”動態(tài)更新訓(xùn)練數(shù)據(jù)集，并針對預(yù)測過程中出現(xiàn)的 “病態(tài)矩陣”采用脊回歸進(jìn)一步改進(jìn)算法，使算法的精確性和魯棒性有了很大提高。實(shí)驗使用POLARIS紅外定位系統(tǒng)采集了10例正常人體表的紅外反射標(biāo)記物的呼吸運(yùn)動數(shù)據(jù)樣本，平均幅度約為20 mm(9.2～37.8 mm)，采用改進(jìn)后的基于記憶學(xué)習(xí)法(預(yù)測步長為1 s)，平均絕對誤差約為0.3 mm(0.08～0.8 mm)，每次估值耗時約1 ms。所提出的方法能夠準(zhǔn)確和實(shí)時捕捉復(fù)雜的呼吸運(yùn)動軌跡。

放射治療；記憶學(xué)習(xí)法；呼吸運(yùn)動；實(shí)時；預(yù)測

引言

精確放射治療技術(shù)基本要求是在腫瘤治療上實(shí)現(xiàn)高精度、高劑量、高療效和低損傷(三高一低)，調(diào)強(qiáng)放療技術(shù)可產(chǎn)生高度適合三維靜態(tài)靶區(qū)形狀的劑量分布。然而由呼吸引起的胸腹部腫瘤運(yùn)動會使目標(biāo)腫瘤移出而正常組織進(jìn)入計劃靶區(qū)，從而極大影響放療的效果，并且會增加發(fā)生正常組織并發(fā)癥的概率。為了減小呼吸運(yùn)動的影響，在傳統(tǒng)的放射治療和IMRT中通常采用腫瘤靶區(qū)射野擴(kuò)邊、呼吸限制和呼吸門控等技術(shù)。但是，這些技術(shù)存在明顯的不足之處，并不能達(dá)到很好的預(yù)期效果。目前，利用圖像引導(dǎo)實(shí)時跟蹤法可以很大程度上減少呼吸運(yùn)動影響，但從腫瘤信息的獲取到射野調(diào)整之間存在一般為幾百毫秒的系統(tǒng)延遲[1]。為補(bǔ)償系統(tǒng)延遲，最好的辦法是通過數(shù)學(xué)建模預(yù)測腫瘤運(yùn)動位置[2]。由于外部呼吸信號和胸腹部腫瘤運(yùn)動具有很好的相關(guān)性[3]，所以利用呼吸信號控制動態(tài)放療過程成為腫瘤實(shí)時跟蹤治療的一個重要研究方向。

Krauss比較了目前4種熱門的預(yù)測呼吸運(yùn)動方法：線性回歸法、神經(jīng)網(wǎng)絡(luò)法、核密度估計法和支持向量回歸法，研究發(fā)現(xiàn)使用任何一種方法，其均方根誤差都比無預(yù)測的情況要小[4]。但是，線性回歸法不易選擇最優(yōu)歷史狀態(tài)數(shù)，準(zhǔn)確性較差，且只適用于延時很小的系統(tǒng)；神經(jīng)網(wǎng)絡(luò)法受網(wǎng)絡(luò)結(jié)構(gòu)和樣本復(fù)雜性的影響較大，容易發(fā)生過學(xué)習(xí)或低泛化能力；核密度估計法計算量很大，難以滿足實(shí)時性要求；支持向量回歸法對大規(guī)模訓(xùn)練樣本難以實(shí)施。本研究提出一種基于記憶學(xué)習(xí)法[5]預(yù)測呼吸運(yùn)動。該方法首先存儲訓(xùn)練數(shù)據(jù)到記憶中，然后從中查找相似數(shù)據(jù)回答當(dāng)前查詢。相似性越高的數(shù)據(jù)點(diǎn)其權(quán)值越高。通過相對簡單的模型擬合局部區(qū)域而不是全局模型，該方法能夠正確捕捉極其復(fù)雜的非線性關(guān)系。由于加權(quán)函數(shù)的局部特性，該方法對訓(xùn)練數(shù)據(jù)的異常值有著一定的魯棒性。而且，訓(xùn)練和適應(yīng)新數(shù)據(jù)幾乎都是即時完成的。

1 材料和方法

1.1材料

實(shí)驗使用加拿大NDI生產(chǎn)的POLARIS紅外定位系統(tǒng)，采集了10組正常人體表的紅外反射標(biāo)記物的呼吸運(yùn)動數(shù)據(jù)。受試者采用仰臥位，標(biāo)記物置于受試者胸骨劍突下約5 cm處，采樣頻率設(shè)為30 Hz。由于在采集過程中存在一定的噪聲，實(shí)驗對系統(tǒng)采集的原始數(shù)據(jù)進(jìn)行平滑濾波，從而得到10例樣本。表1是10例樣本在頭腳方向的呼吸運(yùn)動幅度。圖1是其中兩例典型的呼吸運(yùn)動軌跡：(a)樣本1：周期、幅度比較規(guī)則的信號 (b)樣本9：基線漂移、不等幅變化的不規(guī)則信號。算法執(zhí)行平臺是Matlab7.9.0。

表1 各樣本的呼吸運(yùn)動幅度

圖1 兩例典型的呼吸運(yùn)動軌跡。(a)較規(guī)則的信號；(b)不規(guī)則的信號Fig.1 Two typical respiratory movement. (a)Regular breathing；(b)Irregular breathing

1.2基于記憶學(xué)習(xí)法原理

以采樣頻率φ=30 Hz進(jìn)行數(shù)據(jù)采集，得到一組歷史呼吸軌跡對應(yīng)的離散樣本{s1,s2,si…sM}，當(dāng)i≤M-L時，構(gòu)建狀態(tài)向量：xi={si-Dτ,si-(D-1)τ…si}T，其中si表示i時刻測量點(diǎn)的運(yùn)動幅度，D(≥2)是輸入向量維數(shù)，τ是元素間隔長度，L是預(yù)測的步長。yi表示i+L時刻測量點(diǎn)的運(yùn)動幅度。訓(xùn)練集即表示為S={(xi,yi)}。

基于記憶學(xué)習(xí)法(memory-based learning，MBL)直到需要應(yīng)答某查詢時才處理訓(xùn)練數(shù)據(jù)。它包括存儲訓(xùn)練數(shù)據(jù)到記憶中，以及查找相似數(shù)據(jù)回答特定查詢。相似性通常利用距離加權(quán)函數(shù)來判斷，相似性越高的點(diǎn)其權(quán)值越高。最簡單的是最近鄰法，選擇最近點(diǎn)的輸出值，但通常會造成輸入和輸出之間不連續(xù)函數(shù)映射，給出不理想的結(jié)果。本文使用高斯核

(1)

(2)

(3)

(4)

1.3算法改進(jìn)

1.3.1動態(tài)更新訓(xùn)練集

上述的基于記憶學(xué)習(xí)法的訓(xùn)練集是靜態(tài)不變的，預(yù)測準(zhǔn)確性隨著訓(xùn)練集的增大而提高。然而訓(xùn)練集太大會導(dǎo)致搜索速度過慢。另外，呼吸運(yùn)動隨著時間變化而不規(guī)則變化，查詢點(diǎn)與訓(xùn)練集相隔時間越長，其相似程度越低。因此，先采集一段數(shù)據(jù)構(gòu)成初始訓(xùn)練集，隨著測量的進(jìn)行，采用“滑窗法”持續(xù)進(jìn)行更新，使得訓(xùn)練集始終包含最新數(shù)據(jù)[6]。

1.3.2脊回歸[7]

采用動態(tài)更新訓(xùn)練集方法后，對應(yīng)不同訓(xùn)練集長度，預(yù)測準(zhǔn)確性都有了很大的提高。當(dāng)訓(xùn)練集大小為5 s(一個典型的呼吸周期結(jié)束的最小持續(xù)時間)時，絕大部分?jǐn)?shù)據(jù)擬合程度最好，但某些時刻會出現(xiàn)偏差非常大的情況。分析原因，是因為過擬合而使信息矩陣XTWTWX奇異或接近奇異成為“病態(tài)”[8]，造成其逆矩陣的極度膨脹，進(jìn)而大大增加回歸系數(shù)的誤差均方，影響回歸擬合的準(zhǔn)確性和穩(wěn)健性。本研究嘗試使用脊回歸進(jìn)行結(jié)果優(yōu)化。

為簡單描述，令X=WX,Y=WY，b=(XTX)-1XTY=(b0,b1…bD+1)T，引入一個非負(fù)因子σ，則b(σ)=(XTX+σθI)-1XTY，I為單位矩陣，θ為奇異臨界值。這樣det(b(σ))的值就能迅速變大，從而改進(jìn)估計值質(zhì)量。脊回歸分析通常要先對X變數(shù)做中心化和標(biāo)量化處理，以使不同變數(shù)處于同樣數(shù)量級上而便于比較。這就是引入新變數(shù)Z，令

(5)

于是局部模型函數(shù)變?yōu)?/p>

(6)

上述βZ表示回歸系數(shù)β是由Z變數(shù)估計，它們在統(tǒng)計學(xué)又稱為標(biāo)準(zhǔn)化回歸系數(shù)。βZ的最小平方估計為

(7)

所采用的是Hoerl等建議的計算方法

(8)

式中，D+1為回歸模型的參數(shù)數(shù)目(不包括β0)；s2為式(7)的離回歸均方；biZ為對于βiZ的最小平方估計數(shù)。式(8)實(shí)際上是離回歸均方對回歸系數(shù)平方平均值的一個比率。

1.4算法評估

1.4.1算法比較

1.4.2評價標(biāo)準(zhǔn)

采用下述指標(biāo)進(jìn)行模型預(yù)測效果的檢驗：

(1)均方誤差(root mean squared error, RMSE)

(9)

(2)平均絕對誤差(mean absolute error, MAE)

(10)

實(shí)驗中還將對RDMBL與其他比較算法的平均絕對值進(jìn)行配對t檢驗，以驗證算法改進(jìn)效果的統(tǒng)計顯著性。

2 結(jié)果

設(shè)定M=150，即先采取5 s的數(shù)據(jù)構(gòu)建初始訓(xùn)練集，在預(yù)測階段隨著新數(shù)據(jù)的采集按照先進(jìn)先出原則動態(tài)更新。選定τ=12(0.4 s)以實(shí)現(xiàn)呼吸運(yùn)動特征與測量噪聲之間的平衡；L=30(1 s)已經(jīng)足夠應(yīng)對延時長的情況；D=2以最小化預(yù)測所花費(fèi)的時間。各種方法的平均絕對誤差和均方誤差如表2所示。

表2 不同方法結(jié)果比較(單位：mm)

從表2可以清楚比較各種方法的優(yōu)劣。從總體性能表現(xiàn)來看，最近采樣法

從表2和圖2結(jié)果來看，采用動態(tài)更新訓(xùn)練集方法后，DMBL的精確性和抗噪聲魯棒性有了很大的提高，但結(jié)合圖2(b)可以看出某些時刻仍存在因病態(tài)矩陣而出現(xiàn)大的誤差。從圖2(c)中RDMBL的性能表現(xiàn)可以看出，進(jìn)一步使用脊回歸改進(jìn)后消減了病態(tài)矩陣的影響，使得預(yù)測誤差降到了最低。實(shí)驗中，在預(yù)測樣本1,2,3和4呼吸信號時出現(xiàn)了嚴(yán)重病態(tài)矩陣問題，因此使用脊回歸后預(yù)測性能改進(jìn)明顯；而其他信號沒有或出現(xiàn)趨良態(tài)的病態(tài)矩陣，脊回歸的作用就沒有體現(xiàn)出來。圖3中(a)和(b)表明當(dāng)DMBL出現(xiàn)嚴(yán)重病態(tài)矩陣的情況下，RDMBL能消除大的誤差，從而提高性能。圖3(c)則表明當(dāng)病態(tài)矩陣問題不明顯時，相比DMBL，RDMBL無法進(jìn)一步改進(jìn)性能。

為評估實(shí)測性能及差異的統(tǒng)計顯著性，表3為RDMBL與其他算法的t檢驗結(jié)果，其中RDMBL與最近采樣法、線性回歸法和MBL比較的P值均小于0.05，這表明相比傳統(tǒng)算法，所提出的改進(jìn)方法對性能提升具有統(tǒng)計顯著性。而RDML與DMBL的t檢驗結(jié)果P=0.082，在a=0.1的顯著性水平下，說明脊回歸在一定程度上提升了算法性能。

圖2 不同算法對樣本3的預(yù)測誤差曲線圖。(a) MBL；(b) DMBL；(c)RDMBL。Fig.2 The prediction errors on sample 3 using different algorithms. (a) MBL; (b) DMBL; (c) RDMBL.

表3 RDMBL與其他算法的t檢驗結(jié)果

圖3 DMBL和RDMBL的預(yù)測結(jié)果比較。(a)樣本1；(b)樣本4；(c)樣本8。Fig.3 The comparison between DMBL and RDMBL. (a)Sample 1; (b)Sample 4; (c)Sample 8

圖4 絕對誤差直方圖。(a)樣本1；(b)樣本4；(c)樣本8；(d)樣本9Fig.4 The histograms of the absolute error.(a)Sample 1; (b)Sample 4; (c)Sample 8; (d)Sample 9

最后，圖4是RDMBL預(yù)測值與真實(shí)值的絕對誤差直方圖分布，結(jié)合圖1、圖3和表2可以看出當(dāng)呼吸運(yùn)動曲線較為規(guī)則時，擬合效果好，誤差小，當(dāng)曲線不規(guī)則時誤差比較大，這是因為當(dāng)前查詢點(diǎn)與訓(xùn)練數(shù)據(jù)集相似性低，但動態(tài)更新訓(xùn)練集法使得算法能迅速適應(yīng)新數(shù)據(jù)，在圖4中的直方圖表現(xiàn)為不管呼吸運(yùn)動軌跡是否規(guī)則，絕大部分時間點(diǎn)的預(yù)測值與實(shí)際值的誤差集中在很小值范圍，充分說明了算法具有很好的穩(wěn)健性和精確性。另外，執(zhí)行一次RDMBL預(yù)測算法僅需要約1 ms，完全適用于實(shí)時預(yù)測。

3 討論和結(jié)論

基于記憶學(xué)習(xí)法主要依賴歷史數(shù)據(jù)描述自變量與因變量的關(guān)系，它尋找歷史狀態(tài)中相似的“近鄰”并利用這些“近鄰”預(yù)測下一個時刻的狀態(tài)。呼吸運(yùn)動有復(fù)雜的非線性特征和遲滯現(xiàn)象，構(gòu)造適當(dāng)?shù)淖宰兞靠梢圆蹲狡溥\(yùn)動特性。相比直接由當(dāng)前觀測值構(gòu)成狀態(tài)向量，本文采用一種擴(kuò)展?fàn)顟B(tài)向量法，即自變量由當(dāng)前觀測值與最近兩個反映遲滯時間的歷史樣本生成一個擴(kuò)展的狀態(tài)，能更準(zhǔn)確地把握感興趣時間點(diǎn)的局部動態(tài)特征。采用高斯核距離加權(quán)函數(shù)判斷當(dāng)前數(shù)據(jù)與歷史數(shù)據(jù)的相似性，在選擇帶寬時，采用西爾弗曼法則而不是常用的交叉驗證，可以減少運(yùn)算量。基于記憶學(xué)習(xí)法對訓(xùn)練數(shù)據(jù)集要求較高，由表2和圖2(a)可以看出，MBL在訓(xùn)練數(shù)據(jù)過少時難以捕捉呼吸運(yùn)動的動態(tài)特征,預(yù)測誤差非常大。由圖2(b)可知，采用“滑窗法”動態(tài)更新訓(xùn)練集后，由于始終包含最新數(shù)據(jù)，絕大部分時間能夠跟蹤系統(tǒng)的動態(tài)變化特征。從表2看出，DMBL相比MBL，平均絕對誤差和均方誤差均大幅度下降，說明選擇合適的動態(tài)更新訓(xùn)練集方法可以有效減小訓(xùn)練集規(guī)模并提高估值的準(zhǔn)確性。

但從圖2(a)和圖3可以看出，DMBL在某些時刻仍存在較大的誤差，原因是回歸分析中信息矩陣通常會因為“過擬合”導(dǎo)致病態(tài)，而且病態(tài)矩陣問題還未有很好的解決方法。回歸分析在放療中呼吸預(yù)測和構(gòu)建內(nèi)-外關(guān)系都有運(yùn)用，當(dāng)樣本過小尤其內(nèi)部腫瘤運(yùn)動數(shù)據(jù)樣本數(shù)量稀疏時，很可能會遇到病態(tài)矩陣問題，一般的處理方法是擴(kuò)大樣本，這顯然不適用于實(shí)際要求。本研究首次嘗試在放療中運(yùn)用脊回歸給病態(tài)矩陣引入一個適當(dāng)小的偏差來提高估計數(shù)的穩(wěn)定性和準(zhǔn)確性。由圖3和表2中數(shù)據(jù)的對比結(jié)果可以看出，脊回歸可以消除因嚴(yán)重病態(tài)矩陣引起的大誤差，但對趨良態(tài)的病態(tài)矩陣的作用并不明顯。采用有效的方法減小甚至消除病態(tài)矩陣帶來的消極影響也是回歸分析的重要課題，本文的嘗試有借鑒意義但不完善，可以跟蹤該方向的最新文獻(xiàn)尋找脊回歸的優(yōu)化方法或更為有效的新方法。

本研究結(jié)合動態(tài)更新訓(xùn)練集方法和脊回歸算法，提出一種改進(jìn)的基于記憶學(xué)習(xí)法，對實(shí)驗中所使用的十組數(shù)據(jù)，實(shí)現(xiàn)了用高效的方式(每次估值耗時約1 ms)動態(tài)跟蹤并且預(yù)測復(fù)雜的呼吸運(yùn)動，平均絕對誤差僅約0.3 mm(0.08～0.8 mm)。并且，用t檢驗證明了實(shí)驗結(jié)果的統(tǒng)計顯著性。實(shí)驗結(jié)果表明，改進(jìn)后的算法精確性高，實(shí)時性好，由于加權(quán)函數(shù)的局部特性和脊回歸的作用，該算法有著很好的穩(wěn)健性，特別適合處理誤差大或被異常值污染的情形。而且，該模型適用于延時較長的系統(tǒng)。

上述模型為實(shí)時跟蹤放療中系統(tǒng)延遲問題提供了一種有吸引力的解決方案，下一步的工作應(yīng)用到大量臨床數(shù)據(jù)進(jìn)行測試[9]，考慮采用多點(diǎn)跟蹤，并且建立起外部呼吸標(biāo)識物運(yùn)動和腫瘤運(yùn)動之間的精確關(guān)系，從而與腫瘤靶區(qū)跟蹤技術(shù)相結(jié)合，那么在胸腹部腫瘤放射治療中將會有更加實(shí)際的臨床應(yīng)用。

[1] Jin JY, Yin FF. Time delay measurement for linac based treatment delivery in synchronized respiratory gating radiotherapy [J]. Medical Physics, 2005,32(5): 1293-1296.

[2] Murphy MJ. Tracking moving organs in real time [J]. Semin Radiat Oncol, 2004,14(1): 91-100.

[3] Gierga DP, Brewer J, Sharp G,etal. The correlation between internal and external markers for abdominal tumors: Implications for respiratory gating [J]. Radiation Oncology Biol Phys, 2005,61(5): 1551-1558.

[4] Krauss A, Nill S, Oelfke U. The comparative performance of four respiratory motion preditors for real-time tumour tracking [J]. Phys Med Biol, 2011,56(16): 5303-5317.

[5] Li Ruijiang, Lewis JH, Berbeco RI,etal. Real-time tumor motion estimation using respiratory surrogate via memory-based learning [J]. Phys Med Biol, 2012,57(15): 4771-4786.

[6] Ruan D, Fessler JA, Balter JM. Real-time prediction of respiratory motion based on local regression methods [J]. Phys Med Biol, 2007,52(23): 7137-7152.

[7] 莫惠棟. 脊回歸技術(shù)及其應(yīng)用 [J]. 作物學(xué)報, 2002,28(4): 433- 438.

[8] 莫惠棟. 回歸分析中的病態(tài)矩陣及其改進(jìn) [J]. 作物學(xué)報, 2006,32(1): 1- 6.

[9] Ernst F, Dürichen R, Schlaefer A,etal. Evaluating and Comparing Algorithms for Respiratory Motion Prediction [J]. Phys Med Biol, 2013,58(11): 3911-3929.

TheStudyonPredictingRespiratoryMotionviaMemory-BasedLearninginRadiotherapy

WAN Wei-Quan1ZHANG Hui-Lian2XU Zi-Hai3HE Zhi-Qiang1CHEN Chao-Min1*

1(InstituteofBiomedicalEngineering，SouthernMedicalUniversity，Guangzhou510515,China)2(GuangdongSunwahTechConsultingGroup，Guangzhou510515,China)2(RadiotherapyCenterofPLA303Hospital,Nanning530021,China)

Prediction is necessary to compensate the system latency in the real-time tracking radiation therapy for thoracic and abdominal cancers. However, because of the complexity of the breathing motion, conventional methods are far from clinical requirements. This paper proposed a memory-based learning method to predict respiratory motion. The method stores the training data in memory, then finds relevant data to answer a particular query. Furthermore, the paper adopts dynamic update the training data method and ridge regression aimed at “ill-condition matrix” to greatly improve the accuracy and robustness of the algorithm. Our experiment collected ten respiratory motion data with average amplitude of 20 mm (9.2～37.8 mm) from humans’ body surface using POLARIS infrared positioning system. Using our methods (prediction horizon is 1s), mean absolute error (MAE) was reduced to 0.3 mm (0.08～0.8 mm), per estimate takes 1 ms. The results confirm that the proposed method is able to capture highly complex breathing movement accurately in real time.

radiotherapy；memory-based learning；respiratory motion；real-time；prediction

10.3969/j.issn.0258-8021. 2014. 02.003

2013-06-20， 錄用日期:2014-02-10

廣東省重大科技專項 (2012A080104010)

R814

A

0258-8021(2014) 02-0148-07

*通信作者。E-mail: gzccm@fimmu.com