ARIMA及Holt-Winters指數(shù)平滑模型在河南省肺結(jié)核流行趨勢(shì)預(yù)測(cè)中的應(yīng)用

任嘉豪，徐 潔，楊海燕

鄭州大學(xué)公共衛(wèi)生學(xué)院流行病學(xué)教研室 鄭州 450001

我國(guó)是肺結(jié)核病的主要流行區(qū)域。河南省人口基數(shù)大，農(nóng)村居民占比高，肺結(jié)核疾病負(fù)擔(dān)仍然較重[1-2]。整合移動(dòng)平均自回歸(auto regressive integrated moving average，ARIMA)和Holt-Winters指數(shù)平滑(Holt-Winters exponential smoothing)模型較多應(yīng)用于疾病相關(guān)時(shí)間序列的分析，特別是對(duì)于具有季節(jié)性、趨勢(shì)性等分布特征的疾病(如流感、手足口病等)的流行趨勢(shì)具有良好的預(yù)測(cè)功能[3-5]，在我國(guó)傳染病疫情防控及醫(yī)療衛(wèi)生資源分配方面應(yīng)用廣泛且進(jìn)展顯著[6-7]。本研究利用河南省2013年1月至2022年2月肺結(jié)核發(fā)病監(jiān)測(cè)數(shù)據(jù)建立ARIMA及Holt-Winters指數(shù)平滑預(yù)測(cè)模型，評(píng)估模型的預(yù)測(cè)效能，并利用所建立的模型對(duì)2022年3月至12月的發(fā)病數(shù)進(jìn)行預(yù)測(cè)，為指導(dǎo)河南省肺結(jié)核防控措施的制定提供決策依據(jù)。

1 資料與方法

1.1 資料來(lái)源在線訪問(wèn)河南省衛(wèi)生健康委員會(huì)官方網(wǎng)站，收集河南省2013年1月至2022年2月法定傳染病報(bào)告數(shù)據(jù)構(gòu)建數(shù)據(jù)集，其中2013年1月至2020年12月的肺結(jié)核發(fā)病數(shù)據(jù)用于模型構(gòu)建，2021年1月至2022年2月的數(shù)據(jù)用于模型預(yù)測(cè)效能的評(píng)估。

1.2 模型構(gòu)建

1.2.1序列特征分析及模型選擇 對(duì)于2013年1月至2020年12月的肺結(jié)核發(fā)病數(shù)據(jù)，以時(shí)間(月)為橫坐標(biāo)、以發(fā)病數(shù)為縱坐標(biāo)繪制時(shí)序分布圖，分析序列的趨勢(shì)性及季節(jié)性分布特征。Holt-Winters指數(shù)平滑模型可用于同時(shí)包含水平項(xiàng)、趨勢(shì)項(xiàng)及季節(jié)效應(yīng)的時(shí)間序列。ARIMA模型分為非季節(jié)性和乘積季節(jié)性模型[6]。由于肺結(jié)核發(fā)病具有典型的季節(jié)性及趨勢(shì)性分布特征，故構(gòu)建Holt-Winters指數(shù)平滑模型及ARIMA乘積季節(jié)模型。

1.2.2ARIMA模型[7]首先，構(gòu)建肺結(jié)核發(fā)病數(shù)時(shí)序分布圖，并通過(guò)ADF法分析序列的平穩(wěn)。由于序列非平穩(wěn)，nsdiff()函數(shù)輸出結(jié)果為1，故對(duì)原始數(shù)據(jù)進(jìn)行一次季節(jié)性差分，將其轉(zhuǎn)換為平穩(wěn)序列。繪制平穩(wěn)序列的自相關(guān)函數(shù)(auto correlation function，ACF)圖及偏自相關(guān)函數(shù)(partial auto correlation function，PACF)圖。ARIMA模型中，參數(shù)p、d、q分別表示為普通性的自回歸項(xiàng)數(shù)、差分階數(shù)及移動(dòng)平均數(shù)；P、D、Q分別表示季節(jié)性的自回歸項(xiàng)數(shù)、差分階數(shù)及移動(dòng)平均數(shù)，s為周期[6]。根據(jù)序列特征，參數(shù)值設(shè)定如下：d取值0，D取值1，s取值12，p取值0、1，q取值0、1、2，P取值0、1、2，Q取值0、1；據(jù)此共設(shè)置10個(gè)備選模型，以AICc值最小為標(biāo)準(zhǔn)，選擇最佳模型。

1.2.3Holt-Winters指數(shù)平滑模型[8]該利用時(shí)間序列中的歷史觀測(cè)值，對(duì)序列加權(quán)，權(quán)重由近及遠(yuǎn)按等比級(jí)數(shù)減少，故預(yù)測(cè)值實(shí)際是全部歷史平均值的加權(quán)平均。其模型可近似表示為：

Y=level+slope×t+s_t+irregular_t

其中，Y代表某一特定時(shí)間點(diǎn)肺結(jié)核預(yù)測(cè)發(fā)病數(shù)；level為常數(shù)水平項(xiàng)，用平滑參數(shù)α控制；slope代表趨勢(shì)項(xiàng)，用斜率參數(shù)β控制；s_t代表t時(shí)刻的季節(jié)效應(yīng)，用光滑參數(shù)γ控制，irregulart_t代表t時(shí)刻的隨機(jī)項(xiàng)。

1.3 模型評(píng)價(jià)用所構(gòu)建的模型對(duì)2021年1月至2022年2月河南省肺結(jié)核發(fā)病數(shù)進(jìn)行預(yù)測(cè)。對(duì)模型殘差的自相關(guān)性進(jìn)行Box-Ljung檢驗(yàn)。與實(shí)際值比較，計(jì)算預(yù)測(cè)值的相對(duì)誤差及平均相對(duì)誤差，計(jì)算模型的均方根誤差(root mean squared error，RMSE)和平均絕對(duì)百分誤差(mean absolute percentage error，MAPE)[9]。用模型對(duì)2022年3月至12月河南省肺結(jié)核發(fā)病數(shù)進(jìn)行預(yù)測(cè)。

1.4 統(tǒng)計(jì)學(xué)處理使用R軟件(4.1.2版本)構(gòu)建時(shí)間序列、擬合模型和預(yù)測(cè)，檢驗(yàn)水準(zhǔn)α=0.05。

2 結(jié)果

2.1 建模序列的分布特征2013年1月至2020年12月河南省共報(bào)告肺結(jié)核557 980例，月平均發(fā)病5 166例，其中2013年3月發(fā)病數(shù)最多(7 478例)，2020年2月最少(2 720例)。時(shí)序分布圖見(jiàn)圖1，由圖1可知，肺結(jié)核發(fā)病呈現(xiàn)多峰分布，且具有典型的季節(jié)性分布特征，每年3、4月為發(fā)病高峰。

圖1 建模序列的時(shí)序分布圖

2.2 模型構(gòu)建結(jié)果

2.2.1ARIMA模型 序列的ACF和PACF圖見(jiàn)圖2。備選模型AICc的比較見(jiàn)表1，最終選擇的最佳模型為ARIMA(1,0,1)(2,1,0)[12]。

圖2 ACF(A)與PACF(B)圖

表1 備選模型AICc的比較

2.2.2Holt-Winters指數(shù)平滑模型 用Holt-Winters指數(shù)平滑模型對(duì)原始數(shù)據(jù)集進(jìn)行時(shí)間序列分析，參數(shù)α、β、γ分別為0.218 1，0.000 3和0.000 4，基線level為647 3，slope為-25.82。

2.3 模型評(píng)價(jià)經(jīng)Box-Ljung檢驗(yàn)，ARIMA模型和Holt-Winters指數(shù)平滑模型預(yù)測(cè)殘差均服從正態(tài)分布(χ2=0.219、0.007，P=0.640、0.931)，提示兩個(gè)模型均擬合良好。兩個(gè)模型預(yù)測(cè)效能評(píng)價(jià)見(jiàn)表2及圖3。ARIMA模型預(yù)測(cè)值的RMSE為374.62，MAPE為5.32%，平均相對(duì)誤差為10.07%；Holt-Winters指數(shù)平滑模型分別為322.65、4.87%和8.12%，說(shuō)明兩個(gè)模型預(yù)測(cè)準(zhǔn)確性均較好，適用于河南省肺結(jié)核流行趨勢(shì)的預(yù)測(cè)。兩個(gè)模型對(duì)2022年3月至12月河南省肺結(jié)核發(fā)病數(shù)的預(yù)測(cè)結(jié)果見(jiàn)表3。

表2 兩種模型預(yù)測(cè)效能的評(píng)價(jià)

續(xù)表2

圖3 Holt-Winters指數(shù)平滑模型(A)和ARIMA模型(B)預(yù)測(cè)結(jié)果

表3 兩個(gè)模型對(duì)2022年3月至12月河南省肺結(jié)核發(fā)病數(shù)的預(yù)測(cè)

3 討論

時(shí)間序列分析是基于歷史觀測(cè)數(shù)據(jù)對(duì)未來(lái)數(shù)據(jù)演化進(jìn)行預(yù)測(cè)的一種方法，廣泛應(yīng)用于與傳染病密切相關(guān)的公共衛(wèi)生領(lǐng)域[10-11]。本研究對(duì)河南省2013年1月至2022年2月肺結(jié)核發(fā)病數(shù)據(jù)的分析結(jié)果顯示，肺結(jié)核發(fā)病具有季節(jié)性，每年3、4月發(fā)病數(shù)最多，提示疾病的流行與氣候因素密切相關(guān)。因肺結(jié)核發(fā)病具有季節(jié)性，本研究基于2013年1月至2020年12月河南省肺結(jié)核發(fā)病數(shù)，構(gòu)建了ARIMA乘積季節(jié)模型和Holt-Winters指數(shù)平滑模型，用以預(yù)測(cè)河南省肺結(jié)核流行趨勢(shì)，并利用2021年1月至2022年2月的數(shù)據(jù)對(duì)模型預(yù)測(cè)效能進(jìn)行了評(píng)價(jià)。經(jīng)統(tǒng)計(jì)學(xué)檢驗(yàn)，構(gòu)建的兩個(gè)模型均擬合良好，預(yù)測(cè)值的平均相對(duì)誤差較小，預(yù)測(cè)效能較為理想。

預(yù)測(cè)模型的主要局限性在于遠(yuǎn)期預(yù)測(cè)能力不足。有研究[12]應(yīng)用ARIMA模型，以周為單位對(duì)流感樣疾病發(fā)病趨勢(shì)進(jìn)行了預(yù)測(cè)，49周后預(yù)測(cè)值與實(shí)際值的絕對(duì)誤差較大，認(rèn)為該模型不適用于長(zhǎng)遠(yuǎn)期預(yù)測(cè)。也有研究[13]指出，ARIMA模型對(duì)流感樣疾病44周后的預(yù)測(cè)效果與10周無(wú)統(tǒng)計(jì)學(xué)差異，在一定程度上肯定了ARIMA模型稍長(zhǎng)時(shí)期(不超過(guò)建模時(shí)間跨度的1/5)的預(yù)測(cè)效能。張魯玉等[14]基于2013至2018年我國(guó)丙型病毒性肝炎的發(fā)病數(shù)據(jù)，使用這兩類(lèi)模型預(yù)測(cè)了2019年的流行趨勢(shì)，結(jié)果表明，Holt-Winters指數(shù)平滑模型的預(yù)測(cè)精度相對(duì)較高。

本研究結(jié)果顯示，Holt-Winters指數(shù)平滑模型與ARIMA乘積季節(jié)模型預(yù)測(cè)的2021年前4個(gè)月河南省肺結(jié)核發(fā)病數(shù)的平均相對(duì)誤差較小，分別為3.72%和3.60%；2021年5月、6月及8月預(yù)測(cè)的相對(duì)誤差較大，且兩類(lèi)模型預(yù)測(cè)值均高于實(shí)際值；后4個(gè)月(2021年9月至12月)預(yù)測(cè)的平均相對(duì)誤差均較小，分別為2.78%和2.51%。這可能與2021年5月、6月及8月為河南省部分地區(qū)新型冠狀肺炎流行期有關(guān)，因此時(shí)期居民多采取居家封閉措施，人群流動(dòng)減少，從而一定程度上限制了肺結(jié)核的流行。 2021年各月肺結(jié)核實(shí)際發(fā)病數(shù)的降低提示該年度河南省內(nèi)對(duì)肺結(jié)核流行的防控措施效果良好。2022年1月和2月，兩種模型預(yù)測(cè)的相對(duì)誤差又再度升高，提示兩類(lèi)模型長(zhǎng)期預(yù)測(cè)效能可能并不理想。兩類(lèi)模型對(duì)2022年3至12月的預(yù)測(cè)結(jié)果提示，河南省肺結(jié)核發(fā)病例數(shù)將進(jìn)一步降低，當(dāng)下的措施對(duì)于肺結(jié)核防控仍有重要作用。

綜上，ARIMA乘積季節(jié)模型與Holt-Winters指數(shù)平滑模型對(duì)河南省肺結(jié)核流行趨勢(shì)的預(yù)測(cè)較為準(zhǔn)確，在短期預(yù)測(cè)上有較好的推廣價(jià)值，兩種模型長(zhǎng)期預(yù)測(cè)能力有待進(jìn)一步評(píng)價(jià)。