自回歸移動平均模型的建立及在廣州市肺結核發病預測中的應用

雷宇 何立乾 張廣川 賴鏗 謝瑋 杜雨華

近年來，隨著防控力度的進一步加強，廣州市結核病疫情呈現逐年下降趨勢，但由于人口基數大、流動人口多等因素，結核病疫情依然嚴峻[1-2]。自回歸移動平均(autoregressive integrated moving average，ARIMA)模型可對具有季節效應的時間序列進行建模，因其具有模型結構簡單、預測精度高等特點，已被廣泛應用于傳染病預測，而掌握傳染病的發病趨勢及流行特征，對于制定針對性防控策略具有重要意義[3-4]。因此，本研究通過對廣州市2010年1月至2019年12月肺結核月報告發病數構建ARIMA模型，并對2021年廣州市肺結核發病情況進行預測，為衛生行政部門制定結核病防治策略提供參考。

資料和方法

一、資料來源

通過《中國疾病預防控制信息系統》中的《傳染病報告信息管理系統》，按照報告發病時間收集2010年1月至2020年12月廣州市肺結核月報告發病數據。

二、研究方法

ARIMA模型基本結構表示為ARIMA(p，d，q)×(P，D，Q)s，其中p、d、q分別表示非季節性自回歸階數、非季節性差分階數和非季節性移動平均階數，P、D、Q分別表示季節性自回歸階數、季節差分階數和季節移動平均階數，s表示季節周期的長度。

1.序列平穩化：將2010年1月至2019年12月的肺結核報告發病數據定義為按月為單位的時間序列，使該時間序列進行差分處理使其平穩化，并通過增強迪基-福勒(augmented dickey-fuller，ADF)檢驗驗證差分后的時間序列的平穩性，根據ADF檢驗結果的P值判斷數據是否平穩化，P<0.05提示數據為平穩數據。

2.模型識別：根據自相關系數(autocorrelation coefficient，ACF)和偏自相關系數(partial autocorrelation coefficient，PACF)確定p、q、P、Q的可能值，再根據差分階數確定d和D的階數，并建立多個候選模型。

3.模型診斷：通過Ljung-Box殘差檢驗判斷模型的殘差序列是否為白噪聲，若殘差為白噪聲(即Ljung-Box殘差檢驗P>0.05)，通過赤池信息量(Akaike information，AIC)準則選擇適合最優模型，AIC越小的模型擬合的程度越好。

5.模型應用：應用最優ARIMA模型預測2021年1—12月廣州市每月的肺結核發病例數。

三、統計學處理

采用Excel 2007軟件整理2010年1月至2020年12月廣州市肺結核月報告發病數據，采用R3.6.2軟件構建ARIMA模型并進行預測。通過2010年1月至2019年12月的肺結核月報告發病數構建和篩選最優ARIMA模型，并通過2020年實際報告數據驗證最優ARIMA模型的預測效果。為探索納入不同年份的發病數據構建模型的最優參數的影響，使用2010年1月至2018年12月的肺結核月報告發病數重新構建模型，并對比2019年實際報告發病數評價模型參數的穩定性和預測的準確性。最后，因2020年實際報告發病例數受新型冠狀病毒肺炎疫情影響，報告發病數代表性欠佳，故采用2010年1月至2019年12月的月報告發病數據進行擬合，利用該模型預測2021年肺結核月報告發病例數。

結 果

一、廣州市肺結核月報告發病數的發病趨勢

2010—2019年廣州市共報告肺結核患者115 887例。通過時間序列按照總體趨勢、季節趨勢和隨機誤差進行分解后，報告發病例數呈逐年下降趨勢，ADF檢驗顯示原始序列為非平穩序列(t=-1.386，P=0.323)；廣州市肺結核報告發病例數具有明顯的季節性，在每年3—5月份報告發病數出現高峰，2月份出現低谷(圖1)。

二、ARIMA模型構建

1.序列平穩化與模型識別：為使數據平穩化，將原始序列經過1階非季節差分和1階季節性差分后，ADF檢驗顯示差分后的序列為平穩序列(t=-6.210，P=0.010)。ACF和PACF都在1階截尾，由此判斷p=1，Q=1。而P和Q的判斷可分別取0、1、2逐個試驗。根據差分變換的次數，即d=1，D=1，初步確定以12個月為周期的ARIMA(1,1,1)×(P,1,Q)12模型。

2.模型診斷：根據P、Q的取值，由于P、Q的取值一般在0、1和2之間[6]，因此可初步確定9個ARIMA(1,1,1)×(P,1,Q)12。根據模型總體的顯著性、擬合優度等指標進行比較，選出AIC最小的3個模型。擬合的3個備選模型參數估計值及檢驗結果見表1。用Ljung-Box統計量對3個備選模型的殘差值進行檢驗，殘差均為白噪聲，由AIC統計量可知，3個模型中ARIMA(1,1,1)×(0,1,1)12擬合效果最好。

3.模型評價：采用ARIMA(1,1,1)×(0,1,1)12模型對廣州市2020年1—12月肺結核發病例數進行預測，除3月份和4月份外，實際值均在預測值的95%CI內；除2—6月份，預測值與實際值較為接近，MAPE為7.29%，提示該ARIMA模型預測精度較高。2020年實際報告發病總例數為8111例，預測發病總例數為8642例。見圖2，表2。

為進一步驗證該模型參數的穩定性及預測的準確性，通過2010年1月至2018年12月的肺結核月報告數據重新擬合AMRMA模型。結果顯示，AMRMA(1,1,1)×(0,1,1)12仍為最優預測模型，提示該模型適合廣州市肺結核報告發病數據的預測，見表3。對比預測值和實際值發現，2019年報告肺結核患者的實際例數均在ARIMA模型預測數的95%CI值內，且MAPE為4.91%，AMRIMA模型的預測效果較好，見表4。

注 原始序列圖展示2010—2019年肺結核報告發病例數的變化情況；趨勢圖展示2010—2019年肺結核報告發病數據的中心化移動平均值求得的趨勢項；季節變化圖展示2010—2019年肺結核報告發病數據的季節趨勢，縱坐標為季節指數，季節指數>1表示當月報告發病例數>平均水平，季節指數=1表示沒有季節性，季節指數<1表示當月報告發病例數<平均水平；隨機分布圖展示原始序列的平穩程度，縱坐標數值越接近1表示原始序列越平穩；橫坐標各個小刻度代表月份圖1 2010—2019年廣州市肺結核發病數時間序列

表1 2010—2019年廣州市肺結核報告發病備選預測模型的參數值

表2 2020年1—12月廣州市肺結核報告發病例數預測值與實際報告值的比較

圖2 2020年1—12月廣州市肺結核報告發病例數預測值與實際值的比較

4.模型預測：應用ARIMA(1,1,1)×(0,1,1)12模型對廣州市2021年1—12月肺結核報告發病例數進行預測。結果顯示，2021年廣州市肺結核預測發病數為8270例，月預測發病平均值為689例，較2020年略有上升，見表5。

討 論

準確掌握發病趨勢對于肺結核的早期暴發流行具有十分重要的預警作用，可及早提示結核病防治機構(簡稱“結防機構”)開展相關的防控工作，如重點人群篩查、健康教育和藥品、防護物品采購等[7-8]。ARIMA模型能較好地將時間序列的依存性與隨機干擾因素相結合，已陸續應用在肺結核的發病預測上，但目前尚未見該模型應用在廣州市肺結核發病預測上。因此，本研究通過構建適合廣州市預測的最優模型，預測2021年肺結核發病情況，為制定肺結核應急及防控措施提供參考。

表3 2010—2018年廣州肺結核報告發病備選預測模型的參數值

表4 2019年1—12月廣州市肺結核報告發病例數預測值與實際報告值的比較

表5 2021年1—12月廣州市肺結核預測發病例數

本研究發現，2010—2019年廣州市肺結核報告發病例數整體呈現逐年下降趨勢，這可能與廣州市近年來切實落實結核病防治規劃工作有關，疫情穩步下降，與既往研究相符[1, 9]。此外，廣州市肺結核發病具有明顯的季節性，發病高峰為3—5月份，發病低谷為2月份。冬春季是呼吸道傳染病的高發季節，而我國的農歷春節基本上在2月份，由于風俗習慣、大量流動人口的遷出、2月份天數少等原因，導致2月份就醫的患者例數減少，確診患者也隨之減少，形成“春節效應”[10]。春節過后的3—5月份，大量的流動人口返回廣州務工，加上2月份積累的就診延誤患者，導致3—5月份報告發病數急劇升高，與全國報告發病情況相似[11]。

本研究通過篩選模型，最終確定廣州市肺結核預測最優模型為ARIMA(1,1,1)×(0,1,1)12，模型及其參數均具有統計學意義，但與其他地區研究發現的最優模型參數存在差異，這可能與不同地區的經濟、醫療、人口特征等差異有關[12]。因此，建立預測模型需要根據每個地區特定的報告發病數據來擬合不同的最優模型。此外，為進一步驗證模型參數的穩定性和預測的準確性，本研究納入2010—2018年的數據構建模型進而驗證，結果發現最優模型仍為ARIMA(1,1,1)×(0,1,1)12，其MAPE為4.91%，低于10%，提示該模型適合廣州市肺結核發病數據的預測，且參數穩定，預測精度良好。

ARIMA模型因其結構簡單、易于實現且預測精度高，可將傳染病發展的復雜因素(如人口流動、氣溫、濕度等因素)綜合蘊含在時間變量中，借助其趨勢變化、周期變化等特點，實現對疾病未來走勢進行模擬預測，因而在結核病預測研究中使用較多[13-15]。同時，由于使用ARIMA進行預測時，一般需要數個周期性變化的數據，方可較好地擬合模型，故在急性傳染病的早期階段或未出現周期性變化趨勢前，該方法擬合效果較差。此外，ARIMA模型進行長期預測，會出現一定的偏差，因此只適合用于短期預測。組合利用不同預測模型的優劣來彌補單一模型的缺陷，可以達到提高預測精度和增加穩定性的效果，目前已有學者開始探索不同預測模型組合應用在結核病預測上的效果[16]。

通過比較2020年肺結核報告發病實際值與預測值，發現2020年2—6月實際值與預測值差異較大，其原因可能與2019年底我國暴發了新型冠狀病毒肺炎密切相關。2020年1月廣州發現本地新型冠狀病毒肺炎患者，廣州政府隨即制定了一系列嚴格有效的防控措施，如居家自我隔離、保持社交距離、加強核酸檢測、出門需佩戴口罩等，新型冠狀病毒肺炎疫情迅速得到有效控制，4月份起本地患者基本已經消除。同時，因部分結防機構人員被抽調支援抗擊新型冠狀病毒肺炎疫情、可疑或確診患者擔心在就診期間感染新型冠狀病毒而出現的就診延誤等原因[17]，導致肺結核患者發現力度明顯下降，出現ARIMA模型在短周期內肺結核患者數與預測數出現較大差異的情況，但構建的ARIMA模型的MAPE為7.29%，低于10%，提示該模型具有良好的預測精度，從整體上仍能反映2020年肺結核報告發病數據，可用于廣州市肺結核發病的短期預測。

本研究結果顯示，2021年廣州市肺結核預測發病例數為8270例，較2020年模型發病預測值(8642例)呈持續下降趨勢，而較2020年實際報告發病值(8111例)則略有上升，這可能與模型構建采用了2010—2019年數據有關。有研究發現，與2019年同期相比，2020年1月23日至6月30日期間全國結核病患者發現數下降了20%，湖北省下降了29%，武漢市下降了44%[18]。由此可見，受新型冠狀病毒肺炎疫情影響，全國不同地區的結核病發現力度均受到不同程度的影響，從而導致2020年結核病的報告發病水平出現明顯下降，故未納入2020年報告發病數據擬合ARIMA模型。此后，除個別地區存在散發疫情外，全國新型冠狀病毒肺炎疫情得到有效控制，隨著復工復產的推進，2020年7月起，每月實際報告發病數與預測發病數已十分接近，提示結核病患者發現水平已逐漸恢復至正常水平。但由于結核病屬于慢性傳染性疾病，人群發病可能需要一定的周期才能從報告發病例數上體現出來。同時，2021年廣州市將繼續推進分級診療和綜合防治服務，加強重點人群主動篩查及在疑似患者中應用分子生物學檢查等措施，進一步加強患者的發現力度。因此，預計2021年患者例數可能會較2020年出現小幅反彈。此外，本研究對2021年肺結核預測發病例數的置信區間可用于廣州市肺結核疫情的預警，當超過上限值時，需警惕部分地區可能出現肺結核疫情的暴發流行，需要重點關注。

綜上，ARIMA(1,1,1)×(0,1,1)12模型對廣州市肺結核發病例數的擬合效果較好，可用于廣州市肺結核的短期預測和動態分析，具有良好的應用價值。

本研究存在一定的局限性，因擬合模型的數據來源于《中國疾病預防控制信息系統》中的《傳染病報告信息管理系統》，在實際工作中不同地區存在不同程度的漏報現象，從而導致《傳染病報告信息管理系統》反映的疫情情況被低估。但從目前公開的研究結果中，尚未獲得廣州市相關漏報數據，故無法對數據進行相應的校正，因此本研究預測結果較真實情況有所低估，但對于制定結核病防控策略仍有積極的指導意義。