王延年,雍永強,賈曉燦,李全忠,孫長青+
(1.鄭州大學 信息工程學院,河南 鄭州450001;2.鄭州大學 公共衛生學院,河南 鄭州450001;3.河南省人民醫院 內分泌科,河南 鄭州450003)
CGMS (continous glucose monitoring system)借助 埋在臍周皮下的探頭檢測葡萄糖濃度,能連續獲得患者血糖的變化及波動。國內外學者基于CGMS所提供的數據,提出了多種預測方法,常見的有自回歸 (auto-regressive,AR)模型、自回歸移動平均 (auto-regressive moving average,ARMA)模型[3]、支持向量機[4]、極限學習機[5]、神經網絡[6]等。這些方法都只利用CGMS所提供的數據,沒有考慮外部因素 (飲食、藥物注射、運動等)對人體血糖的影響。
S.Pappada等[7]提出的神經網絡方法和Eleni I.Georga等[8]提出的支持向量回歸方法不僅采用CGMS所提供的數據,而且結合了人體的生理模型。結合生理模型的預測方法考慮比較全面,但需要大量病理學和生理學的知識,所提算法和預測模型復雜龐大,預測有一定的延遲。
本文基于CGMS 所提供的血糖數據,提出一種融合ARIMA 和RBF神經網絡模型預測血糖的綜合預測方法。針對外部事件對血糖的影響,提出一種發現和處理算法,在血糖變化受外部因素影響時能自動調整未來一段時間內的血糖預測值,保障預測精度。
采用ARIMA 模型對未來血糖值進行預測,并計算出誤差項,用RBFNN 算法對輸入值、趨勢值、預測殘差、預測殘差的趨勢值等進行學習、網絡訓練與擬合,最后將ARIMA 計算出的預測值與RBF 算法得到的修正值組合得到準確的預測。同時針對外部因素對血糖變化的影響,提出一種發現和處理算法。
融合ARIMA 和RBFNN 的組合預測模型的具體結構如圖1所示。

圖1 融合ARIMA 和RBFNN 的組合預測模型
基于組合預測模型的血糖預測步驟為:
(1)用ARIMA 模型對從CGMS 采集到的血糖序列y(t)進行建模預測。PH 時間段后的預測結果為(t+PH|t)。
(2)利用 (1)得到結果得到ARIMA 產生的預測殘差:e(t)=y(t)-(t|t-PH)。
(3)將e(T),(1-z-Tm)e(t),y(t),(1-z-Tm)y(t)作為RBFNN 的輸入,得到PH 時間段后的預測結果(t+PH |t)。
從CGMS采集的血糖數據由于存在噪聲干擾需要進行平滑濾波處理。卡爾曼濾波具有計算量小、實時性高、效果好等優點,而且采用Kalman濾波對CGMS采集的數據進行平滑處理能有效地減少信號中的噪聲,同時減少預測值和實際值之間的時延[9]。本文采用Kalman濾波對CGMS采集的數據進行平滑處理。
ARIMA 的一般表達式為

記作ARIMA(p,d,q),其中,xt是血糖值組成的時間序列,i 是自回歸 (AutoregRessive,AR)項的參數,θj是滑動平均 (MovingAverage,MA)項的參數,ω為零均值及方差為δ2的高斯白噪聲。用于對數據進行差分處理,可表示為

ARIMA 處理過程如下:
(1)由于血糖序列是一種非平穩時間序列,首先對它進行差分處理,即


那么血糖序列的ARIMA (p,d,q)模型為

其中,{αt}為白噪聲序列。
(2)模型識別。p,d,q 的選擇是ARIMA 建模的關鍵。首先對血糖數據進行差分處理,每次差分比較其自相關(auto correlation function,ACF)和偏自相關 (partial auto correlation function,PACF)組成的條形圖,選擇最佳的d。p,q的確定采用AIC (Akaiki information criterion)準則和BIC (Bayesian information criterion)準則[10],計算AIC和BIC的方程如下[11]

式中:k——參數模型的個數,n——樣本數。
(3)估計參數并檢驗模型。模型中的所有參數的估計都采用極大似然估計法,然后進行檢驗。方法是檢驗誤差序列是否具有隨機性;如果檢驗通不過就繼續進行識別。
(4)利用選中的最佳參數組成預測模型。ARIMA 模型采用已經采集到的在xt,xt-1,…xt-n時刻的血糖值,預測未來xt+m時刻的血糖值。方程為xnn+m=E[xn+m|xn,xn-1,…,x1]。詳細的過程可參考文獻 [12]。
對于線性時間序列,ARIMA 的預測簡單精準。從CGMS收集到的血糖數據是線性的,但影響血糖變化的因素復雜多樣。只使用線性預測模型會帶來時延,影響預測的準確性。
RBF神經網絡在處理非線性問題時能力很強。不僅學習速度快,且避免局部最小問題,只有少量權值需要調整。
RBF神經網絡學習過程:
(1)利用高斯函數激勵函數,輸入和輸出的非線性映射如下

其中,xp=(xp1,xp2,…xpn)T為第p個輸入樣本,ci為隱含層高斯函數的中心,σ為高斯函數的方差,為歐式范數。
(2)隱含層和輸出層之間的線性映射,函數為

式中:h——隱含層的節點總數,ωij——隱含層到輸出層的線性加權值,yi——網絡中第j 個輸出節點的實際輸出,n——輸出節點數。d為期望輸出值,隱含層的σ可表示為

(3)參照樣本數據,校正輸出層和隱含層之間的參數。本文中RBF的輸入包括4個部分:
(1)ARIMA 產生的預測殘差

(2)ARIMA 預測殘差的趨勢(1-z-Tm)e(t)。
(3)目前的血糖值y(t),即從CGMS讀取的數據。
(4)在時間步Tm目前的血糖趨勢值(1-z-Tm)y(t)。
本文所使用的血糖數據來自河南省人民醫院,是對病人連續監測得到的,采用動態血糖監測系統 (CGMS)采集。CGMS每10秒中接收一次電流信號,每5分鐘將獲得的平均值轉換成葡萄糖值儲存起來,24小時共儲存288個葡萄糖值,共持續監測72個小時。任意選取其中的24 個小時,其中前258個數據作為訓練集,對血糖數據進行建模,采用最后30個樣本作為測試集。24小時所采取的某病人的血糖值如圖2所示。
將288個數據樣本輸入到MATLAB R2012A 中,采用ARIMA 的數據分析功能,得到血糖值自相關圖和偏自相關圖,如圖3所示。從圖3可知,血糖值序列有較高的自我相關性。需要進行差分處理,使血糖值序列變成平穩的時間序列。取差分階數d=1時,血糖序列基本平穩,所以最佳差分階數為d=1。其1階偏相關圖和自相關圖如圖4所示。

圖2 24小時采取的病人的血糖值

圖3 血糖值的ACF圖和PACF圖

圖4 1階偏相關和自相關圖
采用從低階到高階逐步試探法來識別模型的參數,嘗試采用不同的p和q值,比較AIC和BIC取得最小值。
從表1的模型指標對比可以看出,模型ARIMA (0,1,1)比較適合預測血糖值。采用一步預測法進行預測,預測結果如圖5所示:ARIMA 模型用于對血糖的預測,能夠好地把握血糖的變化趨勢,但精度不高。

表1 血糖序列的AIC和BIC檢測值

圖5 ARIMA 模型的預測結果
根據ARIMA (0,1,1)預測結果和血糖實際值得到血糖值的殘差序列,計算出ARIMA 預測誤差的趨勢和當前血糖值的趨勢,結合當前的血糖值對RBF進行建模和訓練,得到預測值。組合模型的預測結果如圖6所示。
從圖6可知,ARIMA-RBFNN 組合預測模型不僅能夠準確把握血糖的變化趨勢,而且能比較準確地預測出短時間能的血糖值。
為了能評價和比較預測結果,本文使用兩個性能指標
(1)均方誤差 (mean squared error,MSE)


圖6 ARIMA-RBFNN 組合預測結果
(2)平均絕對相對誤差 (mean absolute error,MAE)

從表2的對比結果可知,ARIMA-RBFNN 的血糖預測精度要高于單一的ARIMA 模型的預測精度,預測誤差有較大的降低,是一種有效的血糖預測方法。用ARIMARBFNN 模型進行血糖值的預測,只需要利用CGMS提供的歷史數據,相對于結合各種生理模型的預測方法[5,6],ARIMA-RBF神經網絡組合預測算法簡單可靠。

表2 兩種預測方法誤差分析
通過預測結果的觀察,某些時刻預測值與真實值差距較大,經過分析研究,差距較大的時間段大多有飲食、胰島素的注入、運動等事件的發生。如圖7 所示。預測時間段受飲食的影響,飲食后的預測精度大大降低。
定義滿足如下條件的點為奇異點

式中:threshold(k)——判定某點是否為奇異點的門檻,大小可調,當k時刻真實血糖值y(k)與預測的血糖值(k|k-PH)滿足上述條件時為奇異點。
通過實驗和觀察發現,奇異點一旦出現,將會持續3到5個采集點。可以在預測過程中不斷加入實際血糖值來判斷奇異點,為了提高預測的精度,提出了如下的算法:
(1)根據公式判斷采集點是否是奇異點;
(2)如果是奇異點,將當前時刻未來4個時間段血糖預測更新為


圖7 受飲食影響的ARIMA-RBFNN 預測結果
(3)在發現奇異點時,其5個時間段后繼續執行判斷奇異點,如果再次發現,則回到 (1);如果沒有發現,則不做任何處理。
采用ARIMA-RBF預測模型對奇異點進行預測后的效果對比如圖8所示,進行預測后的MSE 和MAE 為0.0069和0.0616,而不是改進前的0.0609和0.1416。

圖8 改進后與改進前的對比
本文基于CGMS,提出一種融合ARIMA 和RBFNN 的組合預測模型,并與ARIMA進行對比。在組合預測的基礎上,針對外部因素對血糖的影響,提出一種發現和處理算法,能夠在血糖受外部因素影響時調整未來一段時間的血糖預測值。實例結果表明,該組合預測模型結合奇異點發現和處理算法能夠顯著提高預測能力和預測精度,很好地反應血糖上升或下降的趨勢,具有很好的推廣和應用價值。
[1]Ceriello A,Novials A,Ortega E,et al.Evidence that hyperglycemia after recovery from hypoglycemia worsens endothelial function and increases oxidative stress and inflammation in healthy control subjects and subjects with type 1diabetes[J].Diabetes,2012,61 (11):2993-2997.
[2]Georga E I,Protopappas V C,Polyzos D,et al.A predictive model of subcutaneous glucose concentration in type 1diabetes based on random forests [C]//Engineering in Medicine and Biology Society.Annual International Conference of the IEEE,2012:2889-2892.
[3]Eren-Oruklu M,Cinar A,Quinn L,et al.Estimation of future glucose concentrations with subject-specific recursive linear models[J].Diabetes Technology & Therapeutics,2009,11(4):243-253.
[4]Wang Z,Lai L,Xiong D,et al.Study on predicting method for acute hypotensive episodes based on wavelet transform and support vector machine[C]//3rd International Conference on Biomedical Engineering and Informatics.IEEE, 2010:1041-1045.
[5]Mo X,Wang Y,Wu X.Hypoglycemia prediction using extreme learning machine(ELM)and regularized ELM [C]//Control and Decision Conference.IEEE,2013:4405-4409.
[6]Allam F,Nossair Z,Gomma H,et al.Prediction of subcutaneous glucose concentration for type-1diabetic patients using a feed forward neural network [C]//International Conference on Computer Engineering &Systems.IEEE,2011:129-133.
[7]Pappada S M,Cameron B D,Rosman P M,et al.Neural network-based real-time prediction of glucose in patients with insulin-dependent diabetes[J].Diabetes Technology & Therapeutics,2011,13 (2):135-141.
[8]Georga E I,Protopappas V C,Ardigo D,et al.Multivariate prediction of subcutaneous glucose concentration in type 1diabetes patients based on support vector regression [J].IEEE Journal of Biomedical and Health Informatics,2013,17 (1):71-81.
[9]Facchinetti A,Sparacino G,Cobelli C.An Online SelfTunable Method to Denoise CGM Sensor Data [J].IEEE Trans Bio-Med Eng,2010,57 (3):634-641.
[10]Ye R,Suganthan P N,Srikanth N,et al.A hybrid ARIMADENFIS method for wind speed forecasting [C]//IEEE International Conference on Fuzzy Systems.IEEE,2013:1-6.
[11]Shumway R H,Stoffer D S.Time series analysis and its applications:With R examples[M].Springer,2010.
[12]Box G E P,Jenkins G M,Reinsel G C.Time series analysis:Forecasting and control[M].John Wiley &Sons,2013.