童 夢 丁國榮 余 楠 王文波
(武漢科技大學(xué)理學(xué)院 武漢 430065)
隨著生活條件和水平提高,糖尿病患者的年齡越來越小,尤其在發(fā)展中國家影響越來越廣泛,糖尿病也是引發(fā)死亡的主要原因[1]。從國際糖尿病聯(lián)盟(IDF)發(fā)布的消息來看,2017 年我國糖尿病患者人數(shù)已達1.14 億,居世界首位。到2019 年,國內(nèi)糖尿病患者出現(xiàn)局部年輕化,人數(shù)增至4 億。目前主要通過血糖濃度來判斷是否患糖尿病,根據(jù)預(yù)測的血糖可以提前改變胰島素注射量或者通過改變運動量以及飲食等盡量避免異常血糖值的出現(xiàn),使血糖濃度處于正常。
國內(nèi)外學(xué)者提出各種方法來預(yù)測血糖濃度。Daskalaki等[2]運用機器學(xué)習(xí)的血糖預(yù)測方法,使用遞歸神經(jīng)網(wǎng)絡(luò)(RNN)的血糖和胰島素信息,其模型優(yōu)于自回歸(AR)模型和使用外部胰島素輸入的AR模型(ARX);Orieke等[3]提出了一種神經(jīng)模糊模型,利用胰島素和膳食效應(yīng)聯(lián)合預(yù)測Ⅰ型糖尿病患者的血糖水平。Orieke對模型的誤差網(wǎng)格分析顯示87.5%的預(yù)測落在了A 區(qū)域,而剩下的12.5%的預(yù)測落在了4h 預(yù)測窗口內(nèi)的B 區(qū)域。Yang Jun 等[4]提出了自回歸綜合移動平均模型(ARIMA),該模型采用了模型階數(shù)的自適應(yīng)識別算法。該辨識算法采用Akaike信息準(zhǔn)則和最小二乘估計,自適應(yīng)地確定模型的階數(shù),同時估計出相應(yīng)的參數(shù)。結(jié)果表明,該模型優(yōu)于自適應(yīng)單變量模型和ARIMA 模型。血糖預(yù)測時間的延長,可以為醫(yī)生和患者提供充足時間進行血糖濃度控制,提高糖尿病治療的效果。目前長短時記憶網(wǎng)絡(luò)在時間序列預(yù)測方面的效果較為良好,現(xiàn)已應(yīng)用于預(yù)測股票價格、電力負荷預(yù)測、海洋表面溫度或者學(xué)習(xí)的生理模型血糖行為等。LSTM可以比其他網(wǎng)絡(luò)學(xué)習(xí)的更迅速和解決過去沒有用前反饋網(wǎng)絡(luò)算法解決的困難問題。成驍彬[5]提出了經(jīng)驗?zāi)B(tài)分解與層疊式長短期記憶算法相結(jié)的模型,并運用該模型對短時風(fēng)速進行預(yù)測。
由于血糖數(shù)據(jù)通常是非線性、非平穩(wěn)的數(shù)據(jù),通常需要對血糖數(shù)據(jù)進行預(yù)處理來降低其非平穩(wěn)性。血糖預(yù)測方法有基于數(shù)據(jù)驅(qū)動的、生理模型的和結(jié)合數(shù)據(jù)和生理的預(yù)測方法,運用不同方法各有其優(yōu)點和缺點?,F(xiàn)階段已有許多預(yù)測方法出現(xiàn)在血糖濃度序列預(yù)測比如灰度預(yù)測[6]、徑向基函數(shù)預(yù)測[7]、ARIMA預(yù)測、神經(jīng)網(wǎng)絡(luò)預(yù)測、支持向量機預(yù)測等,同時運用群體智能優(yōu)化算法對預(yù)測模型參數(shù)進行優(yōu)化,主要有粒子群算法、遺傳算法、蟻群算法、人工蜂群算法[8]、花朵授粉算法[9]和煙花算法[10]等。本文選取的是粒子群算法對LSTM 進行結(jié)構(gòu)參數(shù)的優(yōu)化。
本文結(jié)合上述文獻,針對血糖信號的特點,對血糖信號進行多尺度分解降低其非線性以及非平穩(wěn)性的特點,提出一種基于變分模態(tài)分解和經(jīng)粒子群優(yōu)化長短時記憶網(wǎng)絡(luò)的血糖濃度序列預(yù)測模型(VMD-PSO-LSTM)。
變分模態(tài)分解(Variational Mode Decomposition,VMD)是由Dragomiretskiy 等在2014 年提出,VMD 是一種將信號分解轉(zhuǎn)化成為變分分解的模式,其分解的核心思想是構(gòu)建和求解變分問題,對非平穩(wěn)、非線性的信號有較好的效果[11~12]。
假設(shè)原始信號f是由K個固有模態(tài)分量(IMF)uk(t)組成。
式(1)中,相位φk(t)是不減函數(shù);Ak(t)表示包絡(luò)函數(shù)。
設(shè)各IMF 的中心頻率是ωk,以各IMF 的和等于原始信號作為約束條件時,變分模態(tài)分解的分解過程步驟如下。
首先對各個IMF進行希爾伯特變換,并且得出uk(t)的解析信號,通過單邊譜將uk(t)的中心頻帶調(diào)制到相應(yīng)的基頻帶上,可以通過乘以實現(xiàn):
然后估計出各模態(tài)函數(shù)的帶寬,通過計算解調(diào)信號梯度的L2范數(shù)實現(xiàn),最終將原問題變成求解有約束的變分問題,如式(3)。
其中K是將原信號分解成為K個IMF,δ(t)為狄拉克函數(shù),*為卷積運算符。
通過引入Lagrange 乘法算子λ(t)和二次懲罰因子α的方法來求解上述有約束的變分問題的最優(yōu)解,將原問題轉(zhuǎn)化為無約束變分問題。
最后,將上述問題運用乘子交替方向算法進行求解,通過求解得到的無約束模型的解,而原問題與它的無約束模型問題同解,由式(5)得到所有的IMF。
(ω)是通過當(dāng)前剩余量進行Wiener 濾波器處理得到,中心頻率通過算法各分量功率譜的重心進行重新估計,并通過式(6)更新。
粒子群算法(Particle Swarm Optimization,PSO)是一種模擬鳥群覓食的進化算法,是由J.Kennedy和R. C. Eberhart 等開發(fā)的。PSO 是從初始的隨機解開始,通過不斷迭代尋找最優(yōu)解,運用適應(yīng)度函數(shù)來評價解的好壞。不斷更新局部最優(yōu)和全局最優(yōu),最終到達停止條件的出求得的最優(yōu)解。PSO 首先初始化得到一組隨機解,然后迭代不斷搜索當(dāng)前空間最優(yōu)的粒子,最終獲取最優(yōu)解的過程。
在多維搜索空間中,如果一個群體是以m個粒子組成,在第t次迭代中,第i個粒子的位置和速度分別為Xi,t和Vi,t,粒子通過監(jiān)督個體最優(yōu)和全局最優(yōu)這兩個最優(yōu)解來更新自己的位置和速度。在找尋這兩個最優(yōu)解的過程中。
粒子通過式(7)來更新下一時刻速度,通過式(8)來更新t+1時刻位置:
式(7)中:w為速度慣性因子;c1和c2為學(xué)習(xí)因子;rand 為[0,1]之間的隨機數(shù);式(8)中λ為速度系數(shù),本文取0.5。Vi,t被最大速度Vmax和最小速度Vmin限制。
由于循環(huán)神經(jīng)網(wǎng)絡(luò)(RNN)容易出現(xiàn)梯度消失或者梯度爆炸的現(xiàn)象,存在著長期依賴關(guān)系,因此提出一種特殊的RNN,長短期記憶網(wǎng)絡(luò)(Long Short-Term Memory Networks,LSTM)。與RNN 相同,長短期記憶網(wǎng)路也是一種重復(fù)神經(jīng)網(wǎng)絡(luò)模塊的鏈?zhǔn)叫问剑?3]。在RNN 的基礎(chǔ)上,長短期記憶網(wǎng)絡(luò)的核心思想就是增加了一個單元狀態(tài)C(Cell state),并且通過輸入們和遺忘門這兩個門來控制單元狀態(tài)C。LSTM的結(jié)構(gòu)如圖1所示。單個LSTM主要包括以下幾個步驟。
圖1 LSTM的結(jié)構(gòu)圖
首先遺忘門計算公式如式(9)所示。
其次對輸入門來說,確定的是更新的信息,計算式(10)、(11)所示。
更新從t-1 時刻到t時刻的細胞的狀態(tài)的計算公式為式(12)。
最后一步是這一層網(wǎng)絡(luò)的輸出信息,LSTM 最終預(yù)測結(jié)果的輸出ht,是由輸出門和單元狀態(tài)共同確定。計算公式如式(13)以及式(14)。
為了更好的評價和比較各個模型預(yù)測的結(jié)果,本文選取平均絕對誤差(Mean Absolute Error,MAE)、均方根誤差(Root Mean Squared Error,RMSE)和克拉克誤差網(wǎng)格分析法(Clarke Error Grid Analysis,CEGA)作為性能評價指標(biāo),各個評價指標(biāo)的計算方法如下:
克拉克誤差網(wǎng)格分析法利用笛卡爾圖來評價血糖預(yù)測方法的準(zhǔn)確度。CEGA中橫縱坐標(biāo)軸都是[0,400],分為如圖2所示的A、B、C、D、E五個區(qū)域。
圖2 克拉克網(wǎng)格分析
在克拉克網(wǎng)格分析中,預(yù)測效果較好的點位于A區(qū)域,其次落在B區(qū)域,在臨床上落在區(qū)域A和區(qū)域B 是可以接受的預(yù)測值,預(yù)測效果最差的點落在E區(qū)域。
由于血糖信號的特點,一般的預(yù)測方法難以提高血糖濃度的預(yù)測精度以及模型的準(zhǔn)確性[14~15]。根據(jù)前文提出的方法,本節(jié)在此基礎(chǔ)上進行改進,基于變分模態(tài)分解方法對非平穩(wěn)血糖時間序列信號進行處理以及PSO 算法對LSTM 結(jié)構(gòu)參數(shù)的優(yōu)化,本文提出一種基于VMD 和PSO-LSTM 的血糖濃度預(yù)測模型。
基于VMD-PSO-LSTM的血糖濃度預(yù)測流程如下:
1)獲取血糖濃度序列。利用CGM 連續(xù)采集糖尿病患者的血糖序列。
2)基于VMD 方法將血糖序列分解成5 個相對平穩(wěn)的IMF分量{IMF1,…,IMF5} 。
3)分別對各固有模態(tài)分量分別建立PSO-LSTM預(yù)測模型。對患者血糖值序列的第i個固有模態(tài)函數(shù)IMFi,通過粒子群算法對LSTM 參數(shù)尋優(yōu),利用優(yōu)化后的LSTM 對IMFi進行預(yù)測,獲得各分量的預(yù)測結(jié)果。
4)患者血糖最終預(yù)測值等于PSO-LSTM 模型預(yù)測得到各分量預(yù)測值的和。
5)將VMD-PSO-LSTM預(yù)測值分別與實際值和其他預(yù)測模型進行對比,通過本文提出的性能評價指標(biāo)分析驗證所提方法的精確性。
VMD-PSO-LSTM 預(yù)測血糖濃度的流程圖如圖3所示。
圖3 VMD-PSO-LSTM流程圖
在實驗中,利用VMD 對血糖序列進行5 層分解。PSO-LSTM 模型結(jié)構(gòu)由輸入層、隱藏層和輸出層組成,本文選取的損失函數(shù)為預(yù)測的均方誤差。LSTM 的內(nèi)部參數(shù)采用Adam 算法進行訓(xùn)練,訓(xùn)練次數(shù)設(shè)置為200 次,訓(xùn)練目標(biāo)取0.0001,該模型將隱藏層神經(jīng)元數(shù),LSTM 的時間窗口大小設(shè)置為LSTM網(wǎng)絡(luò)模型超參數(shù)運用粒子群算法進行優(yōu)化。
為了減少人為因素對模型的影響,根據(jù)血糖序列的具體情況,對超參數(shù)的取值范圍設(shè)置如下:指定隱藏層單元個數(shù)取值范圍[1,300],時間窗口大小取值范圍[1,30]。同時設(shè)置粒子群粒子個數(shù)為30,最大迭代次數(shù)取值為30,學(xué)習(xí)因子c1=c2=2。速度慣性因子w=0.8,速度系數(shù)λ=1。適應(yīng)度函數(shù)選取的是LSTM網(wǎng)絡(luò)的均方誤差。
為了更好地驗證本文所提方法的有效性,本文將PSO-LSTM 和VMD-LSTM 預(yù)測模型對同樣的血糖濃度時間序列數(shù)據(jù)進行預(yù)測,各個模型的預(yù)測對比圖結(jié)果如圖4,預(yù)測誤差對比圖如圖5所示。
圖4 各個模型的預(yù)測結(jié)果對比圖
圖5 各個模型的預(yù)測誤差對比圖
在克拉克網(wǎng)格誤差分析中,對各個模型運用相對誤差大小進行校驗相對誤差大小,其計算需要用到測量真實值和預(yù)測值,該方法是一種通過直觀的比較每個點來實現(xiàn)的檢驗方法。根據(jù)本文的研究,它通過比較的就是各個時刻患者血糖的預(yù)測值和它的真實值之間的誤差。 血糖預(yù)測VMD-PSO-LSTM 模型的相對誤差為ε(i),計算方法如式(17)。
則平均相對誤差計算公式如式(18):
模型的預(yù)測精度為p,且p=(1-)×100%
各個預(yù)測模型的評價指標(biāo)對比如表1。
表1 各個預(yù)測模型的評價指標(biāo)
實驗結(jié)果顯示,本文所提方法的MAE和RMSE的均小于PSO-LSTM 和VMD-LSTM 模型的誤差,預(yù)測精度高于PSO-LSTM 和VMD-LSTM 模型的預(yù)測精度,證明了本文方法能夠很好地預(yù)測血糖濃度時間序列,是一種比較好的預(yù)測模型。
本文針對血糖濃度序列的特點,提出一種VMD 與PSO 優(yōu)化LSTM 相結(jié)合的短期血糖濃度預(yù)測模型。首先利用VMD 方法將通過CGMS 獲取的患者的血糖濃度時間序列分解成5 個固有模態(tài)函數(shù),從而降低血糖時間序列的非線性和非平穩(wěn)性;然后對各個IMF分量分別用經(jīng)PSO優(yōu)化的LSTM模型來預(yù)測,最后將預(yù)測后的各子序列疊加得到最終預(yù)測結(jié)果。
通過實驗對比PSO-LSTM 和VMD-LSTM 方法,結(jié)果表明:VMD-PSO-LSTM 的平均絕對誤差和均方誤差都最小,并且所有的點都落在克拉克網(wǎng)格的A 區(qū)域,預(yù)測精度最高,證明了本文所提出的方法是一種比較好的預(yù)測模型。