ARIMA季節(jié)模型在預(yù)測新疆地區(qū)丙型肝炎發(fā)病數(shù)中的應(yīng)用

李 澤， 張學(xué)良

(新疆醫(yī)科大學(xué)1公共衛(wèi)生學(xué)院， 2醫(yī)學(xué)工程技術(shù)學(xué)院， 烏魯木齊 830011)

病毒性肝炎是由肝炎病毒引發(fā)的傳染性疾病，主要分為甲、乙、丙、丁、戊共五種類型。其中丙肝(hepatitis C virus，HCV)發(fā)病率逐年上升且死亡率較高，對(duì)人類健康危害較大，主要通過輸血、靜脈毒品注射、血液透析和器官移植等途徑傳播，醫(yī)源性傳播是其主要的擴(kuò)散形式[1]。我國丙肝感染人數(shù)估計(jì)約3 800萬人，50%為病毒攜帶者，是全球感染人數(shù)最多的國家[2]。2005-2014年，新疆維吾爾自治區(qū)累計(jì)報(bào)告丙肝發(fā)病數(shù)83 983例，死亡數(shù)115人。本文使用新疆地區(qū)的丙肝歷史數(shù)據(jù)，結(jié)合時(shí)間序列方法中的ARIMA(Autoregressive Integrated Moving Average)乘積季節(jié)模型建立新疆地區(qū)丙肝月發(fā)病數(shù)模型，在此基礎(chǔ)上進(jìn)行擬合和短期預(yù)測，為丙肝的防控提供一定的依據(jù)。

1 資料與方法

1.1資料來源2005-2014年新疆丙肝月發(fā)病例數(shù)來源于公共衛(wèi)生科學(xué)數(shù)據(jù)中心。

1.2研究方法

1.2.1 ARIMA季節(jié)模型 乘積季節(jié)模型考慮了時(shí)間序列的長期趨勢(shì)、循環(huán)波動(dòng)、季節(jié)變化以及隨機(jī)波動(dòng)之間相互影響[3]，其公式簡記為ARIMA(p, d, q)×(P, D, Q)S。其中p表示自回歸階數(shù)，d表示差分階數(shù)，q表示移動(dòng)平均階數(shù)，對(duì)應(yīng)的參數(shù)P、D和Q分別表示季節(jié)自回歸階數(shù)、季節(jié)差分階數(shù)和季節(jié)移動(dòng)平均階數(shù)。

1.2.2 ARIMA季節(jié)模型建模步驟 (1)平穩(wěn)性檢驗(yàn)：常見的平穩(wěn)性檢驗(yàn)有圖檢法和單位根檢驗(yàn)，如ADF檢驗(yàn)、DFGLS檢驗(yàn)、KPSS檢驗(yàn)和NP檢驗(yàn)，其中ADF檢驗(yàn)和KPSS檢驗(yàn)運(yùn)用較多[4]；(2)數(shù)據(jù)變換：非平穩(wěn)序列在經(jīng)過Box-Cox變換和差分處理后可轉(zhuǎn)換為平穩(wěn)序列，它是一種將倒數(shù)變換、指數(shù)變換、對(duì)數(shù)變換結(jié)合起來的變換方法[5]，同時(shí)能實(shí)現(xiàn)方差齊性并消除異方差[6]，數(shù)據(jù)變換后需重新做平穩(wěn)性檢驗(yàn)；(3)純隨機(jī)性檢驗(yàn)：純隨機(jī)性檢驗(yàn)選用QBP或QLB統(tǒng)計(jì)量，當(dāng)P<0.05時(shí)認(rèn)為此時(shí)間序列為非白噪聲序列，說明此平穩(wěn)序列中包含值得提取的信息；(4)確定模型結(jié)構(gòu)：繪制自相關(guān)圖ACF和偏自相關(guān)圖PACF，根據(jù)表1中的規(guī)則，估算模型ARIMA(p, d, q)×(P, D, Q)S中參數(shù)p、d、P和Q的范圍，從而確定候選模型；(5)估計(jì)模型參數(shù)：使用矩估計(jì)作為最大似然估計(jì)和最小二乘法迭代的初始值，并估計(jì)各個(gè)候選模型的參數(shù)；(6)模型和參數(shù)顯著性檢驗(yàn)：若模型殘差通過白噪聲檢驗(yàn)且滿足方差齊性，說明此最優(yōu)模型的殘差為白噪聲，否則選擇其他次優(yōu)候選模型，其次還需對(duì)模型中的參數(shù)做顯著性檢驗(yàn)，如果有任何一個(gè)參數(shù)不顯著，則不再選擇此模型，而重新選擇其他候選模型再次檢驗(yàn)；(7)尋找最小信息準(zhǔn)則模型：為了選擇其中最合理的模型，還需要計(jì)算其信息準(zhǔn)則函數(shù)值，常見有AIC、AICc、BIC、DIC、HQC，因AIC/AICc在理論上比BIC更有優(yōu)勢(shì)[7]，且當(dāng)樣本量足夠大時(shí)AICc會(huì)收斂于AIC[8]，同時(shí)AICc更適用于時(shí)間序列模型，因此本文選用AICc作為最優(yōu)模型的評(píng)價(jià)指標(biāo)；(8)模型的交叉驗(yàn)證和預(yù)測：考慮到時(shí)間序列的特點(diǎn)，不宜采用K-fold交叉驗(yàn)證，選用Hold-Out較為合適，把時(shí)序數(shù)據(jù)劃分為訓(xùn)練集和驗(yàn)證集，在訓(xùn)練集上建立模型并估計(jì)參數(shù)，再將候選模型的預(yù)測值和驗(yàn)證集進(jìn)行比較從而判斷誤差，常見的擬合效果評(píng)價(jià)指標(biāo)有MSE、MAPE和SMAPE。

1.3數(shù)據(jù)處理軟件使用R語言3.3.3，預(yù)測包forecast 8.0，時(shí)間序列包tseries 0.10-38，單元根檢驗(yàn)包fUnitRoots 3010.78。

2 結(jié)果

2.1平穩(wěn)性檢驗(yàn)使用R語言繪制2005-2014年新疆地區(qū)丙肝月發(fā)病數(shù)時(shí)序圖，見圖1。對(duì)數(shù)據(jù)做ADF和KPSS平穩(wěn)性檢驗(yàn)，前者P=0.231，后者P<0.01，說明該時(shí)序是非平穩(wěn)的，需要進(jìn)行數(shù)據(jù)變換。

2.2數(shù)據(jù)變換為減少結(jié)果出現(xiàn)異方差的可能性，直接對(duì)原始數(shù)據(jù)做λ=0的Box-Cox變換，即自然對(duì)數(shù)變換。為得到季節(jié)差分和非季節(jié)差分項(xiàng)，對(duì)原始序列做非平穩(wěn)序列的確定性分析，圖2中可以看出明顯的季節(jié)性變化，因此需要做1階12步季節(jié)差分，即D=1，S=12。隨后對(duì)差分?jǐn)?shù)據(jù)再做平穩(wěn)性檢驗(yàn)，發(fā)現(xiàn)依然是非平穩(wěn)的，所以嘗試1階非季節(jié)差分，即d=1,檢驗(yàn)后發(fā)現(xiàn)此時(shí)序已平穩(wěn)。

圖1 2005-2014年新疆地區(qū)丙肝月發(fā)病數(shù)時(shí)序圖

圖2 丙肝月發(fā)病數(shù)的確定性分析

2.3純隨機(jī)性檢驗(yàn)采用QLB統(tǒng)計(jì)量進(jìn)行白噪聲檢驗(yàn)，差異有統(tǒng)計(jì)學(xué)意義(P<0.01)，說明此變換后的平穩(wěn)序列不是白噪聲，序列中包含值得提取的信息。

2.4確定模型結(jié)構(gòu)繪制該平穩(wěn)序列的ACF和PACF，見圖3和圖4。根據(jù)表1的判斷方法，非季節(jié)參數(shù)q可能取值0、1、2，季節(jié)參數(shù)Q可能取值為0、1，非季節(jié)參數(shù)p可能取值為0、1、2、3，季節(jié)參數(shù)P可能取值為0、1、2。因此共有3×2×4×3=72個(gè)候選模型。

圖3 平穩(wěn)序列的自相關(guān)圖ACF

2.5估計(jì)模型參數(shù)使用R語言構(gòu)建了72個(gè)候選模型，每個(gè)模型的參數(shù)均會(huì)被自動(dòng)估計(jì)。

2.6模型和參數(shù)顯著性檢驗(yàn)對(duì)72個(gè)候選模型的殘差做統(tǒng)計(jì)量的白噪聲檢驗(yàn)。隨后做顯著性檢驗(yàn)，自由度為2005年1月-2014年6月訓(xùn)練集的月數(shù)總數(shù)114減去當(dāng)前候選模型的參數(shù)數(shù)量。通過計(jì)算得到72個(gè)候選模型中，有12個(gè)模型呈現(xiàn)顯著性。

圖4 平穩(wěn)序列的偏自相關(guān)圖PACF

2.7尋找最小信息準(zhǔn)則模型計(jì)算上述12個(gè)模型的AICc，見表2。ARIMA(2,1,0)×(1,1,0)12即為最優(yōu)模型。圖5是該模型的殘差平方圖，可以看出沒有明顯的趨勢(shì)，并未呈現(xiàn)出異方差性。表3為模型的參數(shù)顯著性檢驗(yàn)，顯示所有參數(shù)均顯著非零。

表2 通過顯著性檢驗(yàn)的候選模型的AICc值

圖5 最優(yōu)模型的殘差平方圖

參數(shù)模型回歸系數(shù)標(biāo)準(zhǔn)誤t值Par1-0.7720.097-7.985<0.001ar2-0.2740.096-2.8630.003sar1-0.4430.096-4.610<0.001

2.8模型的交叉驗(yàn)證和預(yù)測為了驗(yàn)證模型ARIMA(2,1,0)×(1,1,0)12的外推能力，將2005-2014年的月時(shí)序數(shù)據(jù)劃分為兩部分，2005年1月-2014年6月的月數(shù)據(jù)為訓(xùn)練集，2014年7月-2014年12月的月數(shù)據(jù)為驗(yàn)證集。做Hold-Out交叉驗(yàn)證，訓(xùn)練集MAPE=1.44%，驗(yàn)證集MAPE= 4.80%，驗(yàn)證集SMAPE=2.37%，擬合與預(yù)測效果均較好，擬合情況見圖6?？梢钥闯?，模型ARIMA(2,1,0)×(1,1,0)12在驗(yàn)證集上的外推能力較好。表4給出了驗(yàn)證集上的誤差，平均誤差為4.67%。圖6預(yù)測部分顯示出2015年的丙肝發(fā)病總數(shù)為11 788例，略高于2014年的11 715例，預(yù)測數(shù)據(jù)見表5，發(fā)病數(shù)峰值1 154例,出現(xiàn)在3月。

圖6 最優(yōu)模型的擬合、驗(yàn)證和預(yù)測圖

時(shí)間實(shí)際值預(yù)測值絕對(duì)誤差相對(duì)誤差2014.79271013860.082014.8967923-44-0.052014.9842857150.022014.10753764110.012014.11913973600.062014.121010956-54-0.06

表5 2015年丙肝預(yù)測月發(fā)病數(shù)

3 討論

丙肝逐漸成為突出的公共衛(wèi)生問題，給社會(huì)造成了一定的經(jīng)濟(jì)負(fù)擔(dān)。本研究結(jié)果顯示,2005年新疆地區(qū)丙肝病例數(shù)較少，然而在2006年之后丙肝病例數(shù)持續(xù)增加，2008年之后相對(duì)穩(wěn)定，但發(fā)病數(shù)仍緩慢上升。2004-2010年新疆地區(qū)的法定傳染病發(fā)病率中，丙肝的發(fā)病率平均為26.45/10萬，死亡率平均為0.04/10萬[9]，已成為影響新疆傳染病發(fā)病率的主要原因之一。此外，新疆地區(qū)地域遼闊，各地區(qū)間的經(jīng)濟(jì)發(fā)展水平、衛(wèi)生意識(shí)和習(xí)慣差距大，易導(dǎo)致貧窮和疾病的惡性循環(huán)，因此通過現(xiàn)有數(shù)據(jù)尋找適用于新疆地區(qū)的丙肝預(yù)測模型，將會(huì)為丙肝的防控提供一定幫助。

時(shí)間序列主要研究事物發(fā)展和變化的規(guī)律并預(yù)測未來趨勢(shì)，ARIMA是較為常用的平穩(wěn)時(shí)間序列擬合模型。本文中針對(duì)新疆地區(qū)丙肝發(fā)病數(shù)建立了AICc最小的ARIMA乘積季節(jié)模型ARIMA(2,1,0)×(1,1,0)12用于預(yù)測新疆地區(qū)丙肝發(fā)病數(shù)。但由于ARIMA模型更加適合短期預(yù)測，在做長期預(yù)測時(shí)，最好可以更多地考慮歷史數(shù)據(jù)，從而獲得精確的預(yù)測結(jié)果。本研究結(jié)果顯示個(gè)別候選模型在驗(yàn)證集上的MAPE小于最優(yōu)模型，如模型ARIMA(0,1,2)×(1,1,0)12在驗(yàn)證集上的MAPE=4.71%，優(yōu)于本研究選定模型的外推能力，但是最小信息函數(shù)受到模型的極大似然函數(shù)值和模型中未知參數(shù)個(gè)數(shù)的影響，說明它充分提取了數(shù)據(jù)的信息，對(duì)數(shù)據(jù)的建模更加充分。

[1] 陳兆云, 劉繼文, 孟存仁，等. 新疆地區(qū)漢族、維吾爾族、哈薩克族丙肝患者基因型研究[J]. 新疆醫(yī)科大學(xué)學(xué)報(bào), 2015, 38(7): 855-857.

[2] 王曉軍, 張榮珍, 胡苑笙，等. 我國病毒性肝炎流行現(xiàn)狀研究[J]. 疾病監(jiān)測, 2004, 19(8): 290-292.

[3] 王燕. 時(shí)間序列分析：基于R[M]. 北京: 中國人民大學(xué)出版社, 2015：158.

[4] 陳雙金. 時(shí)間序列單位根檢驗(yàn)方法比較[D]. 成都: 電子科技大學(xué), 2013.

[5] 吳劉倉, 黃麗, 戴琳. Box-Cox變換下聯(lián)合均值與方差模型的極大似然估計(jì)[J]. 統(tǒng)計(jì)與信息論壇, 2012, 27(5): 3-8.

[6] 崔玫意, 張玉虎, 陳秋華. Box-Cox正態(tài)分布及其在降雨極值分析中的應(yīng)用[J]. 數(shù)理統(tǒng)計(jì)與管理, 2017, 36(1): 8-15.

[7] ANDERSON B. Model Selection and Multimodel Inference: A Practical Information-Theoretic Approach[M].2nd ed.Berlin: Springer-Verlag, 2002.

[8] BURNHAM K P, ANDERSON D. Multimodel inference: understanding AIC and BIC in model selection[J].Soc Methods Res,2004,33(2):261-304.

[9] 鄭強(qiáng),王新旗,曹巖,等.2004～2010年新疆法定傳染病流行趨勢(shì)分析[J]. 疾病預(yù)防控制通報(bào),2011,26(6):6-10.