朱 丹，王 偉，付東山

（天津醫(yī)科大學(xué)腫瘤醫(yī)院放療科，國家腫瘤臨床醫(yī)學(xué)研究中心，天津市“腫瘤防治”重點(diǎn)實(shí)驗(yàn)室，天津市惡性腫瘤臨床醫(yī)學(xué)研究中心，天津 300060）

高精度、高劑量、高療效、低損失始終是當(dāng)前腫瘤精確放射治療所追求的目標(biāo)。在實(shí)時圖像引導(dǎo)呼吸運(yùn)動跟蹤或者門控放射治療中，從圖像采集、信號處理到射束調(diào)整需要一定的時間，使治療射束照射滯后于腫瘤運(yùn)動，產(chǎn)生系統(tǒng)時間延遲，影響胸腹部腫瘤放療效果，增加正常組織的損傷。呼吸運(yùn)動預(yù)測是實(shí)時跟蹤放療中的一項(xiàng)關(guān)鍵技術(shù)。

國內(nèi)外學(xué)者進(jìn)行了大量研究和探索，已經(jīng)提出了多種呼吸預(yù)測算法。Murphy等[1]探索了神經(jīng)網(wǎng)絡(luò)(Neural Network,NN)預(yù)測呼吸運(yùn)動的能力，研究發(fā)現(xiàn)神經(jīng)網(wǎng)絡(luò)算法的精確性優(yōu)于線性算法。但神經(jīng)網(wǎng)絡(luò)算法存在一些不足[2-4]，其容易發(fā)生局部最小化問題，影響預(yù)測準(zhǔn)確度；也可能因?yàn)檩斎雽雍碗[藏層的節(jié)點(diǎn)過多產(chǎn)生過擬合問題。Sheng等[5]研究出了一種混合預(yù)測算法，應(yīng)用于CyberKnife?放射外科系統(tǒng) (Accuray Incorporated,Sunnyvale,CA,USA)的Synchrony?呼吸同步追蹤系統(tǒng)中，該混合算法結(jié)合了線性算法、模式匹配和模糊邏輯3種方法。線性算法采用最小二乘實(shí)現(xiàn)線性預(yù)測；模式匹配是在歷史呼吸數(shù)據(jù)中尋找和當(dāng)前呼吸周期最相似的一段呼吸信號作為參考，以此來預(yù)測呼吸運(yùn)動；模糊邏輯是將歷史呼吸數(shù)據(jù)按位置信息分成大、中、小3部分，每一部分再按速度信息分為正速度、零速度和負(fù)速度3部分，構(gòu)成3×3模糊控制規(guī)則表進(jìn)行呼吸預(yù)測。Synchrony?呼吸同步追蹤系統(tǒng)已廣泛應(yīng)用于臨床實(shí)踐，具有較好的實(shí)時跟蹤能力。Kakar等[6]采用自適應(yīng)神經(jīng)模糊推理系統(tǒng)（Adaptive Neuro Fuzzy Inference System,ANFIS）預(yù)測呼吸信號，該方法結(jié)合了神經(jīng)網(wǎng)絡(luò)的學(xué)習(xí)能力和模糊邏輯的推理能力，其基于臨床病例數(shù)據(jù)的實(shí)驗(yàn)結(jié)果表明，ANFIS的均方根誤差在亞毫米級別。史少華等[7]提出了一種基于分離有限狀態(tài)模型的呼吸預(yù)測算法，該算法將呼吸基線與呼吸起伏分離，并分別進(jìn)行預(yù)測，預(yù)測誤差不會因?yàn)楹粑\(yùn)動基線的明顯漂移而受到影響。

在放療過程中，患者呼吸運(yùn)動的基線、幅度、頻率經(jīng)常發(fā)生變化，增加了呼吸運(yùn)動的復(fù)雜度及其預(yù)測的難度。對非規(guī)則呼吸的預(yù)測仍是跟蹤或門控放療研究的重點(diǎn)和難點(diǎn)。呼吸運(yùn)動的基線變化會導(dǎo)致患者解剖結(jié)構(gòu)發(fā)生顯著的時空變化，引起胸腹部腫瘤放療中幾何和劑量學(xué)的誤差[8]。本文提出一種基于小波分解和自適應(yīng)神經(jīng)模糊推理系統(tǒng)的呼吸運(yùn)動預(yù)測算法(WANFIS)。該WANFIS算法利用小波分解將呼吸信號分成基線、低頻、高頻3部分，并分別采用3種不同方法對3部分信號進(jìn)行預(yù)測?；€部分較平滑，采用線性擬合方法進(jìn)行預(yù)測；低頻部分是呼吸信號去除基線后的呼吸幅度部分，采用ANFIS進(jìn)行預(yù)測，ANFIS模型采用位置信息作為輸入?yún)?shù)，通過位置的模糊劃分構(gòu)造模糊集合，建立模糊規(guī)則控制表；高頻部分主要是噪聲，采用簡單移動平均法進(jìn)行預(yù)測。最后綜合3部分預(yù)測值作為呼吸運(yùn)動的預(yù)測結(jié)果。通過30例非規(guī)則呼吸運(yùn)動臨床數(shù)據(jù)的回顧性分析，比較 NN、Synchrony、ANFIS、WANFIS 4種算法的歸一化均方根誤差（Normalized Root Mean Square Error，nRMSE）、平均值（Mean）和標(biāo)準(zhǔn)偏差（STDEV）。

1 資料與方法

1.1 呼吸樣本數(shù)據(jù)來源及預(yù)測過程 本文所采用的病例數(shù)據(jù)來源于天津醫(yī)科大學(xué)腫瘤醫(yī)院射波刀中心CyberKnife治療胸腹部腫瘤患者光學(xué)傳感器記錄的體表呼吸運(yùn)動數(shù)據(jù)，采樣頻率為26 Hz。選擇30例有基線漂移或者振幅突變的非規(guī)則呼吸運(yùn)動數(shù)據(jù)對照比較NN、Synchrony、ANFIS以及WANFIS 4種呼吸預(yù)測算法的精確性。

根據(jù)已報(bào)道的Synchrony呼吸同步追蹤系統(tǒng)的延遲時間[9]，將預(yù)測步長τs定為3個采樣點(diǎn)。實(shí)時預(yù)測時，先對N個歷史呼吸數(shù)據(jù){St-N+1,…,St-1,St}進(jìn)行歸一化預(yù)處理，再應(yīng)用預(yù)測算法進(jìn)行預(yù)測，最后對預(yù)測值進(jìn)行反歸一化得到τs個采樣點(diǎn)后的呼吸運(yùn)動預(yù)測位置St+τs。每新采集一個呼吸數(shù)據(jù)，按照先進(jìn)先出原則更新歷史呼吸數(shù)據(jù)。

1.2 小波分解結(jié)合自適應(yīng)神經(jīng)模糊推理系統(tǒng)的預(yù)測方法（WANFIS） WANFIS算法流程圖如圖1所示。該算法首先通過小波分解將歷史呼吸信號分解成基線低頻高頻部分，再采用不同方法分別進(jìn)行預(yù)測?；€部分反映呼吸信號的基線漂移情況，比較平滑，采用線性擬合算法進(jìn)行預(yù)測；低頻部分是去除基線后的呼吸幅度部分，采用ANFIS算法進(jìn)行預(yù)測；高頻部分主要是噪聲干擾，采用簡單移動平均法進(jìn)行預(yù)測。最后綜合3部分預(yù)測值作為呼吸運(yùn)動的預(yù)測結(jié)果。

圖1 WANFIS算法流程圖Fig 1 Flow chart of WANFIS algorithm

1.2.1 小波分解 小波變換是一種時間頻率局部化分析方法，通過伸縮平移運(yùn)算對信號逐步進(jìn)行多尺度細(xì)化。目前利用小波函數(shù)進(jìn)行分解和重構(gòu)的主流方法是Mallat算法[10]，分解示意圖如圖2所示：

圖2 Mallat小波分解示意圖Fig 2 Schematic diagram of Mallat wavelet decomposition

當(dāng)j=0 時，cA0即為原信號 s；當(dāng)j≥1 時，cAj為第j層分解得到的低頻系數(shù)。cAj與低通濾波器卷積并下采樣可得到第j+1層分解后的近似系數(shù)cAj+1：

cAj與高通濾波器卷積并下采樣可得到第j+1層分解后的細(xì)節(jié)系數(shù)cDj+1：

式(1)和式(2)中：h（m）和g（m）分別是低通濾波器L_D和高通濾波器H_D的沖擊響應(yīng)序列。每分解一次，k的取值范圍減小一半。

下采樣過程的作用是使信號長度減半，減少小波變換的復(fù)雜度。原信號s的長度選擇2的整數(shù)次冪，本文取N=1 024。對小波分解得到的高低頻系數(shù)進(jìn)行單支重構(gòu)可以得到原信號在該系數(shù)對應(yīng)尺度下的高低頻分量，其長度與原信號一致。原信號經(jīng)過一次分解得到第1層的低頻分量和高頻分量，分解出的高頻分量不再分解，對低頻分量繼續(xù)分解得到第2層的低頻分量和高頻分量，依次類推，從而得到信號多個尺度層次上的分量。

頻帶為0～ω的信號經(jīng)過一次小波分解成高低頻兩部分，高頻部分頻帶為ω/2～ω，低頻部分頻帶為0～ω/2[11-12]。根據(jù)Nyquist采樣定理，本文原信號s的頻帶為0～13 Hz，經(jīng)過一次小波分解后低頻分量對應(yīng)的頻帶為0～7.5 Hz，高頻分量對應(yīng)的頻帶為7.5～13 Hz。多尺度分解依次類推，Mallat塔式分解示意圖如圖3所示。

圖3 Mallat塔式分解示意圖Fig 3 Schematic diagram of Mallat pyramidal decomposition

在小波分解中，存在兩個影響因素：一是小波基的選取，本文選取5階Daubechies小波；二是分解層數(shù)的確定，本文首先將原信號分解到第8層，將這層的低頻分量作為基線部分，然后將原信號去除基線的其余部分進(jìn)行一次小波分解，分解得到的低頻分量作為低頻部分，高頻分量作為高頻部分。圖4和圖5是兩個呼吸運(yùn)動信號經(jīng)過小波分解成了基線、低頻和高頻3部分的示例。

圖4 小波分解呼吸信號示例1Fig 4 First example of wavelet decomposition of respiratory signal

圖5 小波分解呼吸信號示例2Fig 5 Second example of wavelet decomposition of respiratory signal

在小波分解過程中，為了消除邊界效應(yīng)，需要選擇合適的邊界延拓方式[13]。在利用小波分解得到基線部分時，采用0階平滑延拓；在利用小波分解得到低頻部分和高頻部分時，采用1階平滑延拓。

1.2.2 線性擬合和簡單移動平均 基線部分采用線性擬合方法進(jìn)行預(yù)測。線性擬合預(yù)測是根據(jù)信號過去n個已知抽樣值序列，通過線性方程，預(yù)測τs個采樣點(diǎn)以后的信號位置的方法：

式（3）中：系數(shù)α通過最小二乘法優(yōu)化，使擬合誤差平方和最小。

簡單移動平均法是將信號過去n個已知抽樣值序列的平均值作為τs個采樣點(diǎn)以后的信號位置預(yù)測值的方法。利用該方法預(yù)測信號分解得到的高頻部分，以減小噪聲對呼吸預(yù)測的干擾。本文取n=6。

1.2.3 自適應(yīng)神經(jīng)模糊推理系統(tǒng)（ANFIS） ANFIS是一種基于Takagi-Sugeno模型的模糊推理系統(tǒng)，它將模糊控制的模糊化、模糊推理和反模糊化三個基本過程全部用神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)實(shí)現(xiàn)，利用神經(jīng)網(wǎng)絡(luò)的學(xué)習(xí)機(jī)制自動地從輸入輸出樣本數(shù)據(jù)中生成ifthen模糊規(guī)則，構(gòu)成自適應(yīng)神經(jīng)模糊控制器，通過在線學(xué)習(xí)不斷調(diào)整模糊推理控制規(guī)則，實(shí)現(xiàn)系統(tǒng)的自適應(yīng)能力。

圖6是ANFIS的一種模型結(jié)構(gòu)[14]，有3個輸入?yún)?shù)，一個輸出參數(shù)f。輸入?yún)?shù)來自呼吸信號低頻部分，輸出f是對應(yīng)時間τs個采樣點(diǎn)以后的位置預(yù)測值A(chǔ)NFIS模型結(jié)構(gòu)中每一層的結(jié)點(diǎn)具有相同的函數(shù)。

圖6 ANFIS模型結(jié)構(gòu)圖Fig 6 Model structural diagram of ANFIS

式（4）中：{a，b，c}是隸屬函數(shù)參數(shù)，該層每個結(jié)點(diǎn)輸出相應(yīng)模糊集的隸屬度。

第2層：這一層每個結(jié)點(diǎn)的輸出是所有輸入?yún)?shù)對模糊集隸屬度的代數(shù)積，代表一條規(guī)則的激勵強(qiáng)度wr。模型構(gòu)成2×2×2模糊集合，共生成8條模糊規(guī)則（r=1，2，...，8）。

第3層：將各條規(guī)則的激勵強(qiáng)度歸一化：

第4層：這一層每個結(jié)點(diǎn)的輸出為每條規(guī)則的輸出，是每條規(guī)則的歸一化激勵強(qiáng)度與對應(yīng)函數(shù)的乘積（r=1，2，...，8），其中：

……

fr是規(guī)則r對應(yīng)的函數(shù)，{br，1，br，2，br，3，br，4} 是規(guī)則函數(shù)參數(shù)（r=1，2，...，8）。

第5層：計(jì)算所有規(guī)則的總輸出，得到模型輸出：

依據(jù)模型輸入輸出參數(shù)構(gòu)造訓(xùn)練集：

矩陣共k+1行（k=N-3τs-1），前三列為模型的輸入數(shù)據(jù)，第4列為模型的輸出數(shù)據(jù)。應(yīng)用訓(xùn)練集訓(xùn)練預(yù)測模型，采用梯度下降法調(diào)節(jié)隸屬函數(shù)參數(shù)，采用最小二乘法調(diào)節(jié)規(guī)則函數(shù)參數(shù)。進(jìn)行預(yù)測時，將最新數(shù)據(jù)｝輸入已訓(xùn)練好的模型，得到τs個采樣點(diǎn)以后的預(yù)測值：

2 結(jié)果

本文基于MATLAB平臺對NN、Synchrony、ANFIS、WANFIS 4種算法進(jìn)行對照比較。圖7表示了呼吸信號的基線出現(xiàn)明顯漂移狀態(tài)下四種算法的預(yù)測結(jié)果，WANFIS 與 NN、ANFIS、Synchrony相比，預(yù)測結(jié)果更加平滑，更接近實(shí)際呼吸運(yùn)動。WANFIS減小了基線漂移的影響。圖8表示了呼吸信號發(fā)生振幅突變情況下4種算法的預(yù)測結(jié)果。WANFIS的預(yù)測結(jié)果優(yōu)于其他3種算法。下面采用nRMSE、Mean、STDEV 3種不同的指標(biāo)定量比較分析四種算法預(yù)測呼吸運(yùn)動的精確性和魯棒性。每例呼吸信號的預(yù)測誤差統(tǒng)計(jì)長度都為2 000個采樣點(diǎn)。

nRMSE代表了呼吸運(yùn)動預(yù)測曲線與實(shí)際呼吸信號間的距離，并且排除了呼吸信號本身幅度大小對結(jié)果的影響：

式（9）中：xi是信號x第i個實(shí)際值是第i個預(yù)測值是x的平均值。

表1列出了30例呼吸信號4種算法的nRMSE。其中，28個呼吸病例的WANFIS的nRMSE小于Synchrony，24個小于ANFIS。

圖7 基線漂移呼吸信號4種算法的預(yù)測結(jié)果Fig 7 Prediction results of four algorithms for baseline drift respiratory signal

圖8 振幅突變呼吸信號4種算法的預(yù)測結(jié)果Fig 8 Prediction results of four algorithms for amplitude jump respiratory signal

表1 4種算法的nRMSE比較Tab 1 nRMSE comparison of four algorithms

表2統(tǒng)計(jì)了4種算法預(yù)測誤差的平均值和標(biāo)準(zhǔn)偏差。WANFIS預(yù)測誤差的平均值和標(biāo)準(zhǔn)偏差明顯小于NN和Synchrony；對于30個病例中的25個，WANFIS小于ANFIS。

表2 四種算法預(yù)測誤差的平均值和標(biāo)準(zhǔn)偏差比較Tab 2 Comparison of Mean and STDEV of prediction error among four algorithms

3 討論

本研究提出的WANFIS算法不僅具有神經(jīng)網(wǎng)絡(luò)的學(xué)習(xí)能力和模糊邏輯的推理能力，還結(jié)合了小波變換的分解能力，將呼吸信號分解成基線、低頻、高頻3部分，根據(jù)每部分的特征分別進(jìn)行預(yù)測?；?0例臨床數(shù)據(jù)回顧性分析，將WANFIS算法與NN、Synchrony、ANFIS這3種典型預(yù)測算法進(jìn)行對照比較。30例臨床數(shù)據(jù)是有基線漂移或者振幅突變的非規(guī)則呼吸運(yùn)動。本文提出的WANFIS算法的nRMSE平均值為0.09，小于NN的0.17、Synchrony的0.11以及ANFIS的0.11。

Kakar等[6]比較了ANFIS對經(jīng)過呼吸訓(xùn)練和不經(jīng)過呼吸訓(xùn)練的呼吸信號的預(yù)測結(jié)果，發(fā)現(xiàn)ANFIS更適合預(yù)測規(guī)則穩(wěn)定的呼吸信號。對于不規(guī)則呼吸信號，ANFIS預(yù)測誤差過大。因此本文將小波分解和ANFIS結(jié)合起來提出WANFIS算法，將呼吸信號的基線分離出來，從不規(guī)則呼吸中提取出相對規(guī)則的信號成分采用ANFIS進(jìn)行預(yù)測，減小了基線漂移對預(yù)測精度的影響，使預(yù)測結(jié)果更加準(zhǔn)確。史少華等[7]利用求若干個周期呼吸數(shù)據(jù)平均值的方法提取基線，來解決基線漂移問題。但該方法不能立即檢測到基線的變化，無法應(yīng)用于臨床。WANFIS通過小波分解可以即時探測到基線的變化，能夠應(yīng)用于臨床實(shí)時跟蹤中。另外，如果預(yù)測結(jié)果含有太多高頻成分，可能會在放療呼吸運(yùn)動補(bǔ)償中導(dǎo)致機(jī)器顫動[5]。WANFIS將呼吸信號高頻部分分離出來通過簡單移動平均進(jìn)行預(yù)測，相當(dāng)于對呼吸信號濾波去噪，減小了噪聲對呼吸預(yù)測的影響，使預(yù)測的呼吸曲線更加平滑。這是小波分解的另一個作用。

WANFIS算法的精確性和魯棒性都優(yōu)于NN、Synchrony以及ANFIS 3種典型算法，其預(yù)測結(jié)果更接近實(shí)際呼吸運(yùn)動，能更有效地補(bǔ)償系統(tǒng)時間延遲。本研究的下一步工作是將WANFIS算法結(jié)合臨床實(shí)踐，應(yīng)用于呼吸實(shí)時跟蹤放療和呼吸門控放療。