江門市疾病預(yù)防控制中心(518000)
黃 國 朱宇平 黃煥鶯
手足口病(hand-foot-mouth disease,HFMD)是人腸道病毒引起的一種常見傳染病,嬰幼兒及學(xué)齡前兒童多發(fā),患兒以發(fā)熱,手、足、口腔等部位皮疹或皰疹為主,傳染性強(qiáng),在人群密集的地方短時間內(nèi)可引起暴發(fā)或流行[1-2]。近年來江門市手足口病呈逐年高發(fā)趨勢,為有效預(yù)防控制手足口病暴發(fā)流行,本研究構(gòu)建自回歸積分滑動平均模型(autoregressive integrated moving average model,ARIMA)預(yù)測江門市手足口病發(fā)病趨勢,探討該模型在預(yù)測手足口病發(fā)病率中的應(yīng)用。
1.資料來源 江門市手足口病月報告病例數(shù)來自中國疾病預(yù)防控制信息系統(tǒng),人口數(shù)來源于江門市統(tǒng)計局。
2.研究方法 利用SPSS統(tǒng)計軟件中ARIMA模型分析方法,首先根據(jù)江門市2009年1月-2017年6月手足口病月發(fā)病率建立時間序列。手足口病月發(fā)病率時間序列為季節(jié)性時間序列,故采用乘積季節(jié)模型,即ARIMA(p,d,q)×(P,D,Q)s。其中d為平穩(wěn)化過程中差分的階數(shù),p、q為自回歸和移動平均階數(shù)。P、Q為季節(jié)性自回歸和移動平均階數(shù),D為季節(jié)差分階數(shù),S為季節(jié)周期。通過數(shù)據(jù)平穩(wěn)化處理、模型識別、參數(shù)估計與檢驗(yàn)等步驟,探索建立模型,將2017年7-12月實(shí)際發(fā)病率與預(yù)測發(fā)病率相對誤差(相對誤差=|實(shí)際值-預(yù)測值|/實(shí)際值)進(jìn)行比較作為外推驗(yàn)證,評價模型預(yù)測效果。
(1)序列平穩(wěn)化 根據(jù)2009年1月-2017年6月江門市手足口病發(fā)病率序列圖、自相關(guān)系數(shù)(ACF)圖和偏相關(guān)系數(shù)(PACF)圖判斷序列平穩(wěn)性。若序列為非平穩(wěn)序列,對原序列進(jìn)行非季節(jié)差分或季節(jié)差分,消除序列長期趨勢和周期性變化的影響,使序列平穩(wěn)[3]。
(2)模型識別 根據(jù)差分后序列自相關(guān)函數(shù)(ACF)圖和偏自相關(guān)函數(shù)(PACF)圖,為模型進(jìn)行初步識別和定階。
(3)模型參數(shù)估計和檢驗(yàn) 利用非線性最小二乘法估計模型參數(shù),在參數(shù)有統(tǒng)計學(xué)意義基礎(chǔ)條件上,用擬合優(yōu)度比較模型優(yōu)劣。模型的擬合優(yōu)度采用標(biāo)準(zhǔn)化的貝葉斯準(zhǔn)則比較,標(biāo)準(zhǔn)化BIC值最小,Ljung-BoxQ統(tǒng)計量P值>0.05的模型為最優(yōu)。
(4)評價模型預(yù)測效果 比較2017年7-12月的實(shí)際發(fā)病率與預(yù)測發(fā)病率相對誤差,驗(yàn)證模型預(yù)測效果。
1.繪制序列圖 繪制2009年1月-2017年6月江門市手足口病發(fā)病率時間序列圖(圖1)。由圖可見江門市手足口病發(fā)病有明顯的季節(jié)性,以12個月為流行周期。每個流行周期出現(xiàn)2個流行高峰,大高峰出現(xiàn)在5-7月,小高峰出現(xiàn)在9-10月,低谷出現(xiàn)在11月份到次年2月份。該序列既有季節(jié)周期性波動特點(diǎn),又有逐年上升趨勢,故采用乘積ARIMA模型。
圖1 2009年1月-2017年6月江門市HFMD月發(fā)病率時間序列圖
2.序列平穩(wěn)化 2009年1月-2017年6月江門市手足口病月發(fā)病率時序圖呈周期性波動趨勢,不能滿足平穩(wěn)化的要求,根據(jù)時序圖季節(jié)性波動特征,對原序列進(jìn)行自然對數(shù)轉(zhuǎn)換和一階季節(jié)差分,差分后的時間序列自相關(guān)(ACF)和偏自相關(guān)(PACF)函數(shù)無明顯截尾和拖尾現(xiàn)象(圖2、3),也不呈線性衰減趨勢,差分后的時間序列圖(圖4)接近平穩(wěn),提示差分后序列適合時間序列模型。
圖2 2009年1月~2017年6月江門市HFMD月發(fā)病率季節(jié)差分ACF函數(shù)圖
3.模型的參數(shù)估計和診斷 確定模型類型后,需要確定p、d、q和P、D、Q的值,對模型定階。根據(jù)序列季節(jié)化特征和平穩(wěn)化處理過程,d=0,D=1。根據(jù)自相關(guān)函數(shù)圖和偏自相關(guān)函數(shù)圖,p=1,q=1。季節(jié)模型P、Q值較難判斷,根據(jù)文獻(xiàn)[4-5],參數(shù)P、Q很少超過2階,可分別取0~2由低階到高階摸索試驗(yàn),結(jié)合模型的擬合優(yōu)度、殘差以及系數(shù)間的相關(guān)性進(jìn)行估計。采用Ljung-Box方法檢驗(yàn)殘差白噪聲,非白噪聲模型排除。經(jīng)試驗(yàn),模型ARIMA(1,0,1)(0,1,1)12標(biāo)準(zhǔn)化BIC值(9.87)最小,平穩(wěn)R2=0.73,殘差序
列的自相關(guān)系數(shù)及偏相關(guān)系數(shù)均在95%置信區(qū)間內(nèi)(圖5),Ljung-Box=21.76,P=0.11。由此,ARIMA(1,0,1)(0,1,1)12模型被選為最優(yōu)模型。
圖3 2009年1月-2017年6月江門市HFMD月發(fā)病率季節(jié)差分PACF函數(shù)圖
圖4 2009年1月-2017年6月江門市HFMD月發(fā)病率季節(jié)差分序列圖
圖5 模型ARIMA(1,0,1)(0,1,1)12殘差序列ACF、PACF函數(shù)圖
4.模型預(yù)測 按照ARIMA建模方法,對2009年1月-2017年6月江門市手足口病月發(fā)病率時間序列建模,再以2017年7-12月全市手足口病月發(fā)病率為驗(yàn)證數(shù)據(jù)進(jìn)行驗(yàn)證,并繪制實(shí)際值和預(yù)測值序列圖,見圖6。根據(jù)預(yù)測值與實(shí)際值相對誤差來判斷模型的預(yù)測效果(表1)。
圖6 ARIMA(1,0,1)(0,1,1)12模型擬合圖
表1 2017年7-12月江門市HFMD月發(fā)病率預(yù)測值與實(shí)際值比較
手足口病2008年5月起被納入丙類傳染病[6],其傳染性強(qiáng),病原學(xué)復(fù)雜,傳播途徑多,可多次重復(fù)感染,手足口病預(yù)防控制工作難度大,早防早控工作一直難以落實(shí)到位。及時有效地預(yù)測預(yù)警發(fā)病趨勢,是該病預(yù)防控制工作的重點(diǎn)和難點(diǎn)[7-8]。時間序列模型將復(fù)雜因素的綜合效應(yīng)統(tǒng)一蘊(yùn)含到時間變量中,克服了疾病發(fā)病影響因素錯綜復(fù)雜,或有關(guān)數(shù)據(jù)資料無法獲得的難題[9-10],在具有典型趨勢特征變化的數(shù)據(jù)預(yù)測上適用性好[11-12]。
本研究利用2009年1月-2017年6月江門市手足口病月發(fā)病率資料,通過序列平穩(wěn)化、模型識別、參數(shù)估計及診斷、模型預(yù)測效果評價等步驟,建立了ARIMA(1,0,1)(0,1,1)12模型。該模型對江門市手足口病發(fā)病率進(jìn)行了較好地擬合,說明在短時間、實(shí)際發(fā)病趨勢無較大波動時,ARIMA模型可以對發(fā)病情況進(jìn)行較好的預(yù)測,特別是季節(jié)性模型可以對手足口病季節(jié)性特征做很好的擬合,提前判斷疫情走勢,為防控策略的制定提供科學(xué)依據(jù)[13]。
從區(qū)間估計看,本研究預(yù)測數(shù)據(jù)與實(shí)際發(fā)病情況區(qū)間估計一致,實(shí)際發(fā)病率全部落入預(yù)測值95%CI內(nèi)。在預(yù)測的精度上,2017年7-12月手足口病月發(fā)病率預(yù)測最小相對誤差為7.53%,最大相對誤差為22.38%,平均相對誤差為18.14%,但預(yù)測值的95%CI寬度偏大。在手足口病實(shí)際防控工作中,受社會因素、氣候因素、人群免疫水平等影響[14],手足口病發(fā)病情況復(fù)雜多變,特別是江門市作為珠三角地區(qū),人口密集且流動性大,需更進(jìn)一步探索手足口病發(fā)病預(yù)測模型研究工作,使其更具有實(shí)際的指導(dǎo)意義與價值。