邵升清，夏桂梅

(太原科技大學(xué) 應(yīng)用科學(xué)學(xué)院，山西 太原 030024)

隨著新冠疫情的大暴發(fā)，傳染病問題再一次出現(xiàn)在公眾的視野里，傳染病的預(yù)防問題也受到人們的高度重視.目前，傳染病發(fā)病率的預(yù)測(cè)方法多種多樣，而且得到了廣泛的應(yīng)用.常用的預(yù)測(cè)模型有時(shí)間序列模型ARIMA[1]、馬爾科夫鏈模型[2]、灰色模型GM(1,1)[3]和趨勢(shì)外推模型等[4].隨著社會(huì)的進(jìn)步和計(jì)算機(jī)的發(fā)展，相關(guān)理論研究也不斷完善，出現(xiàn)了利用機(jī)器學(xué)習(xí)的預(yù)測(cè)模型，有支持向量機(jī)預(yù)測(cè)[5]、人工神經(jīng)網(wǎng)絡(luò)模型預(yù)測(cè)[6]、分割K-最鄰近算法預(yù)測(cè)[7]等方法.

本文根據(jù)已有的病毒性肝炎的發(fā)病率數(shù)據(jù)，從統(tǒng)計(jì)學(xué)角度進(jìn)行研究，比較不同樣本容量下ARIMA(p,d,q)模型和GM(1,1)模型對(duì)傳染病的預(yù)測(cè)效果，并選擇最佳模型預(yù)測(cè)短期內(nèi)病毒性肝炎的發(fā)病率及發(fā)展趨勢(shì)，為制定防治措施提供理論依據(jù).

1 數(shù)據(jù)與方法介紹

1.1 研究數(shù)據(jù)

樣本1 我國(guó)1990-2019年的病毒性肝炎發(fā)病率的數(shù)據(jù)，樣本容量n=30.數(shù)據(jù)均來源于《中國(guó)統(tǒng)計(jì)年鑒》.

1.2 研究方法

(i)ARIMA(p,d,q)模型

時(shí)間序列{Xt}的自回歸滑動(dòng)平均模型[8]定義

Xt=φ0+φ1Xt-1+φ2Xt-2+…+φpXt-p+εt-θ1εt-1-θ2εt-2-…-θqεt-q，

(1)

將差分運(yùn)算與ARMA(p,q)模型結(jié)合后，構(gòu)成ARIMA(p,d,q)模型，其中d為差分階數(shù)，B為延遲算子.則ARIMA(p,1,q)的結(jié)構(gòu)為

φ(B)(1-B)Xt=θ(B)εt.

(2)

ARIMA(p,d,q)的結(jié)構(gòu)為

(3)

(ii)GM(1,1)模型

(4)

稱式(4)為GM(1,1)模型，求解得到

(5)

其中，a和b可通過式(6)用最小二乘估計(jì)得到

(6)

對(duì)累加序列{Yt}累減還原，則原序列的預(yù)測(cè)值為

(7)

1.3 模型評(píng)價(jià)指標(biāo)

2 結(jié)果

2.1 ARIMA模型的建立

2.1.1 原始序列的平穩(wěn)性檢驗(yàn)

通常，通過觀察時(shí)間序列圖的曲線來確定對(duì)每個(gè)變量數(shù)據(jù)進(jìn)行ADF檢驗(yàn)時(shí)使用檢驗(yàn)方程哪一個(gè)[10].觀察時(shí)序圖圖1可知，原始序列有截距，無明顯的時(shí)間趨勢(shì)，故選擇表1中類型2的檢驗(yàn)方程xt=μ+φ1xt-1+φ2xt-2+…+φpxt-p+εt.檢驗(yàn)水平α取0.05，由表1知，原始序列類型2的p值均大于0.05，故序列不平穩(wěn)，需要進(jìn)行差分處理.

圖1 原始序列時(shí)序圖

2.1.2 差分序列平穩(wěn)性檢驗(yàn)與白噪聲檢驗(yàn)

(i)差分序列的平穩(wěn)性檢驗(yàn)

圖2是一階差分序列的時(shí)序圖，由圖2知，差分序列無截距μ，也沒有明顯的時(shí)間趨勢(shì)βt，故應(yīng)選擇表1中類型1的檢驗(yàn)方程xt=φ1xt-1+φ2xt-2+…+φpxt-p+εt.由表1知，差分序列類型1的p值小于顯著水平0.05，故一階差分序列平穩(wěn).

圖2 差分序列時(shí)序圖

表1 原始序列與差分序列在不同情況下ADF檢驗(yàn)的p值

(ii)差分序列的白噪聲檢驗(yàn)

H0序列值之間相互獨(dú)立;

H1序列值之間存在相關(guān)關(guān)系.

差分序列在延遲6階的Q統(tǒng)計(jì)量的值為14.593，p-value值為0.02367；延遲12階的Q統(tǒng)計(jì)量的值為32.93，p-value值為0.0009925；故p-value值均小于0.05，應(yīng)拒絕原假設(shè)H0，接受備擇假設(shè)H1，即差分序列之間存在相關(guān)關(guān)系.

2.1.3 模型識(shí)別與模型診斷

由于樣本數(shù)據(jù)的自相關(guān)函數(shù)圖和偏自相關(guān)函數(shù)圖沒有很好的截尾性質(zhì)，故對(duì)p和q分別取遍0、1、2 、3、4的不同的階數(shù)的模型進(jìn)行檢驗(yàn)，同時(shí)，也利用截尾性質(zhì)和AIC/BIC準(zhǔn)則確定一些備選模型.對(duì)備選模型進(jìn)行檢驗(yàn)，結(jié)果見表2.

表2 不同模型下殘差的正態(tài)性、獨(dú)立性、參數(shù)顯著性檢驗(yàn)結(jié)果匯總表

說明：數(shù)字“1”表示通過檢驗(yàn)，即殘差滿足正態(tài)性、獨(dú)立性或參數(shù)估計(jì)值顯著，數(shù)字“0”表示未通過檢驗(yàn).“*”號(hào)表示顯著性比較好.

從表2知，模型ARIMA(0,1,3)的參數(shù)估計(jì)值顯著，殘差服從正態(tài)分布，且相互獨(dú)立.又由標(biāo)準(zhǔn)化殘差圖知，模型的殘差序列基本落入(-2,2)內(nèi)，滿足零均值等方差的特點(diǎn).

綜上，模型ARIMA(0,1,3)通過了模型檢驗(yàn).

2.1.4 參數(shù)估計(jì)

由(3)式計(jì)算知，ARIMA(0,1,3)的表達(dá)式為

xt=xt-1+εt-θ1εt-1-θ2εt-2-θ3εt-3，

(8)

經(jīng)計(jì)算，模型參數(shù)的估計(jì)值為MA(1)=0.7487482，MA(2)=0.4133138，MA(3)=0.6645649,且參數(shù)估計(jì)值的p值都小于0.05，說明參數(shù)估計(jì)值顯著不為零.

將參數(shù)估計(jì)值代入(8)式中，得到擬合模型ARIMA(0,1,3)的表達(dá)式為

xt=xt-1+εt-0.7487482εt-1-0.4133138εt-2-0.6645649εt-3.

2.1.5 模型預(yù)測(cè)

用樣本1(n=30)建立的模型ARIMA(0,1,3)預(yù)測(cè)2020—2022年病毒性肝炎發(fā)病率，分別為90.7401 /10萬、88.4158 /10萬、88.9763/10萬，則我國(guó)病毒性肝炎未來的發(fā)病率趨勢(shì)呈下降狀態(tài)，預(yù)測(cè)值誤差的標(biāo)準(zhǔn)差為5.2880，10.5109，15.4407.從圖3看出，原始序列值在均擬合模型的置信區(qū)間內(nèi)，說明模型擬合效果比較好.

圖3 ARIMA(0,1,3)模型下原始序列的擬合與預(yù)測(cè)圖

說明：豎虛線右側(cè)為預(yù)測(cè)值，左側(cè)為擬合曲線.“*”為原始序列值，黑色實(shí)線為原始序列的擬合曲線，黑色虛線為95%的置信線.

2.3 GM(1,1)模型的建立

2.3.1 模型建立

樣本1(n=30)的參數(shù)估計(jì)值a= -0.0021，b= 306.53，則擬合方程為

Yt+1=146312.25e0.0021t-145966.67.

2.3.2 模型擬合度檢驗(yàn)

使用后驗(yàn)差比值C檢驗(yàn)法和小誤差概率P檢驗(yàn)法對(duì)GM(1,1)模型進(jìn)行檢驗(yàn)[11-12].對(duì)于樣本1(n=30)，GM(1,1)模型的C值為0.899，P值為0.5，模型預(yù)測(cè)精度等級(jí)不合格.

2.3.3 模型預(yù)測(cè)

用樣本1(n=30)建立的模型GM(1,1)預(yù)測(cè)2020年—2022年病毒性肝炎發(fā)病率，分別為97.6385 /10萬、98.4129 /10萬、99.1934/10萬.該模型預(yù)測(cè)精度等級(jí)不合格，預(yù)測(cè)效果較差.

3 結(jié)論

對(duì)于樣本容量n=30的序列，ARIMA(0,1,3)擬合效果較好，MRE值為4.68%(見表3)；GM(1,1)模型的擬合度檢驗(yàn)不合格，且MER值為15.00%，不適合進(jìn)行預(yù)測(cè).

表3 兩種模型的比較

灰色系統(tǒng)GM(1,1)模型和ARIMA(p,d,q)模型都可以進(jìn)行預(yù)測(cè)，但每個(gè)模型各有利弊，故實(shí)際應(yīng)用中須使用最優(yōu)模型進(jìn)行預(yù)測(cè).

ARIMA(p,d,q)模型常用于有時(shí)間特性(如季節(jié)性，周期性)的樣本序列，對(duì)大樣本數(shù)據(jù)擬合效果較好.ARIMA(p,d,q)模型是提取時(shí)間序列中的相關(guān)信息，并以此建模進(jìn)行預(yù)測(cè)，故其局限性是只適用于平穩(wěn)非白噪聲時(shí)間序列.

灰色系統(tǒng)GM(1,1)模型常用于已知信息少，且數(shù)據(jù)規(guī)律性差的樣本.該模型將無規(guī)律的數(shù)據(jù)轉(zhuǎn)化為有規(guī)律的數(shù)據(jù)序列，根據(jù)此規(guī)律進(jìn)行后期預(yù)測(cè).