ARIMA模型在江西省布魯氏菌病發病數預測中的應用

黃玉萍，傅偉杰，熊長輝，劉曉青，胡國良

人間布魯氏菌病(以下簡稱“布病”)是我國法定報告的乙類傳染病，是由布魯氏菌屬細菌引起的變態反應性人獸共患病。經皮膚、消化道及呼吸道傳播，臨床上以長期發熱、多汗、乏力、關節疼痛、肝脾及淋巴結腫大為特點，若診治不及時，不徹底，易導致慢性感染，可遺留骨關節不可逆的器質性損害，使活動受限或出現中樞神經系統后遺癥[1]，嚴重影響了人群健康及社會生產力，是世界上一項重要的公共衛生問題[2]。2015年江西省被國家確定為南方草地開發試點省之一，將大大促進江西省羊養殖、販運、加工銷售、餐飲等行業的發展，羊產業從業人員也將迅速增加。在目前疫情形勢下，運用時間序列模型對布病疫情進行科學的預測，為布病防控提供科學依據顯得十分重要。

1 材料與方法

1.1材料 發病數資料來自江西省疾病預防控制中心法定傳染病報告管理系統，2014年1月至2017年12月布病月發病數作為擬合集建立模型，對2018年1月的發病數進行預測。

1.2 ARIMA模型的建立

1.2.1 模型識別

1.2.1.1序列平穩化 通過觀察原始序列或差分后序列的時間序列圖(擬合曲線在未來一段時間內的均值和方差不發生明顯變化)、ADF檢驗(P<0.05)獲得一個非白噪聲的平穩序列，提取原序列趨勢效應和季節效應[3-4]。

1.2.1.2定階 (偏)自相關函數分別在q、p階后截尾(幾乎95%(偏)自相關系數的都落在2倍標準差范圍以內)。

1.2.2參數估計 對于平穩序列,利用最小二乘法或最大似然法進行參數估計，參數估計值的顯著性檢驗可以得出模型中的一些項是否需要，一般要求每個參數都有統計學意義。如果擬合的多個模型均通過參數顯著性檢驗，引入赤池信息準則( Akaikes information criterion，AIC)概念，進行模型優化，以該數值小者為優[5-6]。

1.2.3模型檢驗 模型的殘差序列自相關函數和偏相關函數均在置信區間內，且P>0.05，該序列通過殘差白噪聲檢驗，模型可用來預測。

1.3模型預測 模型預測采用均方根誤差(Rood Mean Squared Error， RMSE)、平均絕對誤差(Mean Absolute Error，MAE)、希爾不等系數(Theil Inequality Coefficient，Theil IC)、偏方差(Bias Proportion，BP)、方差率(Variance Proportion，VP)、協變率(Covariance Proportion，CP)等指標衡量，前5個指標越小、最后一個指標越大說明預測效果越好。模型的建立和預測通過Eviews10軟件實現。點預測采用誤差進行預測準確性的比較。

2 結 果

2.1江西省2014-2017年布病月發病數序列圖 2014-2017年間江西省布病發病總數218例，單月最多12例，最少0例，數據整體不算平穩(見圖1)。(偏)自相關圖和ADF檢驗(P<0.05)均顯示序列平穩(見圖2)。用Eviews10對數據進行擬合,得到不做差分時的最優模型AR(1)，R2=0.23，AIC=5.138，RMSE=2.934，MAE=2.429，Theil IC=0.279，BP=0.002，VP=0.339，CP=0.659。

圖1 江西省2014-2017年布病月發病數序列圖Fig.1 Monthly incidence of brucellosis in Jiangxi Province,2014 to 2017

圖2 江西省2014-2017年布病月發病數序列的(偏)自相關圖Fig.2 ACF and PACF diagram of the number of brucellosis cases per month in Jiangxi Province,2014 to 2017

2.2一階差分分析 用Eviews10對一階差分后的數據進行擬合，得到最優模型ARIMA(1,1,3)，R2=0.297，AIC=5.204，RMSE=2.933，MAE=2.376，Theil IC=0.282，BP=0.053，VP=0.145，CP=0.801。

2.3二階差分分析 用Eviews10對二階差分后的數據進行擬合，得到最優模型ARIMA(3,2,0)，R2=0.453，AIC=5.776，RMSE=4.065，MAE=3.089，Theil IC=0.326，BP=0.0003，VP=0.044，CP=0.955。

2.4三階差分分析 用Eviews10對三階差分后的數據進行擬合，得到最優模型ARIMA(3,3,0)，R2=0.690，AIC=6.335，RMSE=5.299，MAE=4.247，Theil IC=0.384，BP=0.0003，VP=0.253，CP=0.747。

2.5模型的擬合度和預測誤差分析 從擬合圖(見圖3)來看，擬合曲線最吻合的為ARIMA(3,3,0)，其次為ARIMA(3,2,0)；AR(1)和ARIMA(1,1,3)擬合效果較差。就各模型擬合的指標(見表1)來看，AR(1)的AIC值(5.138)和希爾不等指數(0.279)最小；ARIMA(1,1,3)的均方根誤差(2.933)和平均絕對誤差(2.376)最小；ARIMA(3,2,0)的偏差率(0.0003)和方差率(0.044)最小，協變率(0.955)最大；ARIMA(3,3,0)的R2值(0.69)最大。

2018年1月江西省布病報告數為5例，4個最優模型中，AR(1)的預測誤差最小為0.21，其次是ARIMA(3,3,0)。拋開預測誤差較大的ARIMA(1,1,3)和ARIMA(3,2,0)，比較AR(1)和ARIMA(3,3,0)，可以發現ARIMA(3,3,0)有著更吻合的趨勢圖和最大的R2值，而AR(1)有著更小的AIC值和希爾不等指數。

圖3 四種模型擬合圖的比較Fig.3 Comparison of fitting graphs of four models

表1 四種模型擬合參數的比較
Tab.1 Comparison of fitting indicators of four models

模型模型擬合指標預測值誤差R2AICRMSEMAETheil ICBPVPCPAR(1)0.2305.1382.9342.4290.2790.0020.3390.6595.210.21ARIMA(1,1,3)0.2975.2042.9332.3760.2820.0530.1450.8016.301.3ARIMA(3,2,0)0.4535.7764.0653.0890.3260.00030.0440.9556.371.37ARIMA(3,3,0)0.6906.3355.2994.2470.3840.00030.2530.7474.10-0.9

3 討 論

江西省2014-2017年布病報告總數為218例，病例多集中于5-9月，與北方地區不一致[7]。發病高峰不穩定。原因可能是布病近年來才在江西省內出現，相關醫務工作者及患者對其較陌生，對于發現疾病的靈敏度不夠，從而出現了診斷和報告的滯后現象。往年1月份江西省布病報告數一般為0～2例，2018年1月報告數為5例，明顯超過以往水平。報告單位都是省內市級及以上醫療單位，2017年底江西省內布病病例分布覆蓋所有地級市[8]，這說明各市級醫療單位對于布病有了一定的警惕性，但同時要加強縣(區)級醫療機構的人員培訓，以實現布病的早發現、早治療，早處理、早報告。值得關注的是，流調報告顯示該5例病例中，3例有明顯的羊接觸史；1例為1歲兒童，食用新鮮蒸煮過的羊奶；1例為家務待業者，無明確傳染源。結合病例中非職業人群的占比由2015年的30%增長至2017年的48%，表明布病的傳播形式多樣，疾病的溯源和防控依然存在著較大難度。

ARIMA模型是根據數據序列偏、自相關函數建立起線性的數據間相互依賴的定量模型，整合了趨勢因素、周期因素和隨機誤差等因素，既吸取了回歸分析的優點又發揮了移動平均的長處，且只需發病資料就可進行預測，是短期預測最為成熟的時間序列預測方法之一[9-10]。

本研究探索在不同差分次數條件下創建預測模型，得出了幾個最優模型。其中，AR(1)和ARIMA(3,3,0)的預測值最接近于實際值，而ARIMA(3,2,0)和ARIMA(1,1,3)的預測誤差較大，說明赤池信息準則和R2值是確定最優模型的兩個較好指標。另一方面，ARIMA(3,2,0)的R2、協方率和擬合圖的表現都好于ARIMA(1,1,3)，預測值的誤差卻較大，猜測可能是由于預測誤差本身存在的影響。同時，一般認為ARIMA模型法至少需要50個以上數據建模[11]，本研究只有48個數據且含有3個零值的小樣本。

模型的選擇很大程度上是研究者根據ACF圖和PACF圖再結合自身經驗綜合確定的相對最優模型，所以研究者給出的并不是最優模型而是基于自身經驗做出的最優選擇[12]。專家指出p、q的值一般不超過5，然而有研究跳出這個限制嘗試建立了幾個可能的模型，發現ARIMA(7,1,0)(1,0,1)12有著更小的絕對誤差和相對誤差(平均值<5%)，并且殘差檢驗大于0.8提示擬合效果更好[13]。盡管ARIMA模型在很大程度上提取了時間序列信息，但原理上依然是一個線性模型。而布病的發生受不同宿主動物與傳播途徑等因素的影響，因此建立的預測模型并不唯一，在今后的研究中考慮在條件允許的情況下盡可能的全面收集布病發病的影響因素, 將相關因素納入模型中提高模型精度，以達到最佳的預測效果。

利益沖突：無