自回归移动平均模型的建立及在广州市肺结核发病预测中的应用

雷宇 何立乾 张广川 赖铿 谢玮 杜雨华

近年来，随着防控力度的进一步加强，广州市结核病疫情呈现逐年下降趋势，但由于人口基数大、流动人口多等因素，结核病疫情依然严峻[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模型对广州市肺结核发病例数的拟合效果较好，可用于广州市肺结核的短期预测和动态分析，具有良好的应用价值。

本研究存在一定的局限性，因拟合模型的数据来源于《中国疾病预防控制信息系统》中的《传染病报告信息管理系统》，在实际工作中不同地区存在不同程度的漏报现象，从而导致《传染病报告信息管理系统》反映的疫情情况被低估。但从目前公开的研究结果中，尚未获得广州市相关漏报数据，故无法对数据进行相应的校正，因此本研究预测结果较真实情况有所低估，但对于制定结核病防控策略仍有积极的指导意义。