馬 霞,王曉佳

（太原工業(yè)學(xué)院理學(xué)系，山西太原 030008）

瘧疾是一種蟲(chóng)媒傳染病，通過(guò)蚊蟲(chóng)叮咬傳播，癥狀主要為貧血和脾腫大后全身發(fā)冷、發(fā)熱、出汗等。瘧疾是我國(guó)法定傳染病之一，預(yù)測(cè)瘧疾的發(fā)病情況對(duì)于我國(guó)的醫(yī)療衛(wèi)生事業(yè)有重要意義。ARI?MA模型是Box和Jenkins[1]在20世紀(jì)70年代初提出的著名的預(yù)測(cè)方法。國(guó)內(nèi)外對(duì)ARIMA的應(yīng)用十分廣泛，樊雯婧[2]應(yīng)用ARIMA模型擬合1991?2011年合肥市每月瘧疾發(fā)病率，并預(yù)測(cè)2012年各月份瘧疾發(fā)病率，取得了較好的預(yù)測(cè)效果。肖丹[3]采用貝葉斯統(tǒng)計(jì)模型研究了索馬里的惡性瘧疾發(fā)病情況。鄭慧敏[4]研究了ARIMA模型在深圳市法定傳染病發(fā)病趨勢(shì)預(yù)測(cè)的應(yīng)用。戴琳琳[5]研究了基于模型的青島市GDP的預(yù)測(cè)分析，為青島市經(jīng)濟(jì)決策提供重要的依據(jù)。近幾年，國(guó)內(nèi)外學(xué)者綜合考慮了序列的趨勢(shì)變化、周期差異等，將ARIMA模型應(yīng)用于GDP的預(yù)測(cè)當(dāng)中，結(jié)果顯示模型預(yù)測(cè)效果良好[6?7]。為減小誤差，提高預(yù)測(cè)精度，很多學(xué)者通常選擇采用多個(gè)模型組合的方法或者將單一模型進(jìn)行修正，從而提高預(yù)測(cè)精度[8?9]。運(yùn)用自回歸求和移動(dòng)平均(模型)對(duì)全國(guó)瘧疾各月發(fā)病趨勢(shì)進(jìn)行預(yù)測(cè)，并探討模型的準(zhǔn)確性，同時(shí)其傳染病預(yù)測(cè)模型的研究提供參考，為瘧疾的防控工作提供依據(jù)。

ARIMA(p，d，q)稱為差分自回歸移動(dòng)平均模型，AR是自回歸，p是自回歸項(xiàng)。MA是移動(dòng)平均，q是移動(dòng)平均項(xiàng)的數(shù)目，d是時(shí)間序列穩(wěn)定時(shí)的差分?jǐn)?shù)目。ARIMA模型的表達(dá)式[1]為：

其中Xt表示在t時(shí)刻的一個(gè)時(shí)間序列，εt表示一個(gè)白噪聲(零均值，方差是一個(gè)常數(shù))，d表示差分階數(shù)，B表示后移算子，即BXt=Xt?1，?=1 ?B，φ(B)表示自回歸算子，自回歸系數(shù)多項(xiàng)式為：

θ(B)表示移動(dòng)平均算子，移動(dòng)平均系數(shù)多項(xiàng)式為：

若模型包含季節(jié)性的成分，還應(yīng)當(dāng)考慮帶季節(jié)性的SARIMA（p，d，q）（P，D，Q）S乘積模型，其中P表示季節(jié)周期的自回歸的階數(shù)，D表示差分次數(shù)，D表示滑動(dòng)平均的階數(shù)，s表示季節(jié)周期。其中φ(B)?dXt表示相同周期內(nèi)各周期點(diǎn)之間的相互關(guān)系，φ(B)?D則表示的是不同周期中對(duì)應(yīng)時(shí)點(diǎn)上的相互關(guān)系。SARIMA（p，d，q）（P，D，Q）S乘積模型的數(shù)學(xué)表達(dá)式為：

1 數(shù)據(jù)的預(yù)處理

選取了全國(guó)2012?2019年各月瘧疾發(fā)病人數(shù)，數(shù)據(jù)如表1（數(shù)據(jù)來(lái)源于中國(guó)疾病控制中心網(wǎng)站[10]），表中數(shù)據(jù)可以看出全國(guó)瘧疾患病總?cè)藬?shù)整體是逐年遞減的，2013年患病總數(shù)偏高，之后呈波浪式下降趨勢(shì)。

表1 全國(guó)2012-2019各月度瘧疾發(fā)病人數(shù)/例

2 模型的建立

2.1 模型的檢驗(yàn)

根據(jù)2012?2019年全國(guó)瘧疾的發(fā)病數(shù)制作時(shí)序圖(圖1)，可看出2013年的數(shù)據(jù)明顯高于其他年份，瘧疾的發(fā)病情況具有明顯的季節(jié)性，序列是非平穩(wěn)的，因此，需要進(jìn)行一階季節(jié)性差分，去除季節(jié)的影響，最終獲得平穩(wěn)序列（圖2）。選用圖形檢測(cè)法和單位根檢驗(yàn)方法檢驗(yàn)平穩(wěn)性。由表2可知，ADF單位根檢驗(yàn)的統(tǒng)計(jì)量為?6.219 312，該結(jié)果明顯小于各顯著水平臨界值，符合平穩(wěn)性條件，即該序列不存在單位根，發(fā)病數(shù)Yt是平穩(wěn)的，因此，ARIMA（p，d，q）中的d=0，故可以構(gòu)建ARIMA（p，q）模型。

圖1 2012-2019年全國(guó)瘧疾發(fā)病數(shù)序列圖

圖2 經(jīng)處理后的瘧疾月度發(fā)病數(shù)序列圖

表2 處理后瘧疾月度發(fā)病數(shù)序列的ADF單位根檢驗(yàn)結(jié)果

2.2 模型的識(shí)別與參數(shù)估計(jì)

根據(jù)序列的平穩(wěn)性和季節(jié)特征，通過(guò)觀察序列自相關(guān)系數(shù)圖和偏相關(guān)系數(shù)圖（圖3和圖4），ACF在1后隨參數(shù)P的增大而突然減小，PACF在2后隨參數(shù)P的增大而突然減小，從而初步確定連續(xù)模型為ARI?MA(1，0，2)。通過(guò)觀測(cè)季節(jié)自相關(guān)系數(shù)圖和偏相關(guān)系數(shù)圖（圖5和圖6），季節(jié)模型的參數(shù)P和Q，一般不超過(guò)二階，將0，1，2分別代入從低階到高階進(jìn)行調(diào)試，根據(jù)殘差分析得出模型ARIMA(1，0，2)(0，1，1)12、模型ARIMA(1，0，2)(1，1，0)12和模型ARIMA(2，0，2)(0，1，1)12為備選模型。BIC值最小的是最優(yōu)預(yù)測(cè)模型，顯示擬合優(yōu)度最好的是ARIMA(1，0，2)(0，1，1)12(如表3)。

圖3 0階差分自相關(guān)圖

圖4 0階差分偏自相關(guān)圖

圖5 1階差分自相關(guān)圖

圖6 1階差分偏自相關(guān)圖

表3 備選模型參數(shù)比較

2.3 模型的診斷

通過(guò)驗(yàn)證殘差序列，判斷模型ARIMA(1，0，2)(0，1，1)12的適用性。由AFC圖和PACF圖所給信息，根據(jù)殘?jiān)肼暀z驗(yàn)圖(圖7)可知相關(guān)系數(shù)均分布在95%的置信區(qū)間內(nèi)，判斷ARIMA(1，0，2)(0，1，1)12模型是適用的。

圖7 殘?jiān)肼暀z驗(yàn)圖

2.4 模型的預(yù)測(cè)

用ARIMA(1，0，2)(0，1，1)12對(duì)2012?2019年全國(guó)瘧疾各月發(fā)病人數(shù)進(jìn)行擬合，得出其擬合值與瘧疾實(shí)際值進(jìn)行對(duì)比如圖8所示?？芍懠驳陌l(fā)病數(shù)均在擬合值的置信區(qū)間內(nèi)，通過(guò)該模型預(yù)測(cè)出2020年1?12月全國(guó)瘧疾發(fā)病的人數(shù)如表4。把每個(gè)月的數(shù)據(jù)進(jìn)行累加求和便可得2020年瘧疾的發(fā)病數(shù)為2 148。雖然預(yù)測(cè)的模型預(yù)測(cè)的數(shù)據(jù)和實(shí)際值不完全吻合，但是從表5可知，預(yù)測(cè)值與真實(shí)值的相對(duì)誤差很小。進(jìn)一步說(shuō)明模型的擬合效果較好，因此，該可以用來(lái)預(yù)測(cè)我國(guó)瘧疾的流行發(fā)病情況。

圖8 2020年全國(guó)每月的瘧疾發(fā)病數(shù)預(yù)測(cè)圖

表4 預(yù)測(cè)2020年我國(guó)瘧疾每月發(fā)病數(shù)與95%置信區(qū)間

表5 2014-2020年我國(guó)瘧疾發(fā)病數(shù)的預(yù)測(cè)值與實(shí)際值的對(duì)比

3 結(jié)語(yǔ)

根據(jù)2012?2019年的各月瘧疾發(fā)病數(shù)，利用ARI?MA復(fù)合季節(jié)模型，可看出全國(guó)瘧疾發(fā)病數(shù)呈明顯的以年為周期的震動(dòng)現(xiàn)象。根據(jù)殘差分析，自相關(guān)和偏自相關(guān)系數(shù)，正態(tài)化的BIC值確定模型的擬合優(yōu)度，確定了模型ARIMA（1，0，2）（0，1，1）12的適用性，預(yù)測(cè)了2020年各月度瘧疾發(fā)病數(shù)，均在預(yù)測(cè)值的95%置信區(qū)間內(nèi)，因此，該模型可以用來(lái)預(yù)測(cè)瘧疾的流行發(fā)病情況。