手足口病(hand food and mouth disease，HFMD)是由多种肠道病毒引起的常见传染病，以婴幼儿发病为主。近几年来，中国HFMD疫情呈现发病强度高、高峰持续时间长、疫情分布广泛等特点。2008年5月2日HFMD正式纳入丙类传染病管理，各级医疗机构须按照《中华人民共和国传染病防治法》和《传染病信息报告管理规范》的有关规定进行网络直报。HFMD已不再是一种陌生的疾病，从研究趋势可见，随着社会的发展以及人们健康和保健意识的提高，人们不仅满足于掌握和了解HFMD的传播规律，更期望能准确预测不同情况下HFMD的传播特征和流行趋势。由于传染病在时间、空间上具有特定的分布模式，且其分布特征随自然与社会环境等因素的变化而变化，而在不断出现的各种新因素干扰下，HFMD的传播规律变得更为错综复杂，尽管国内外专家在HFMD治疗和防制取得了一定进展，但不同时期、不同地区发生的肠道病毒感染往往呈现不同的流行特点，包括临床特征都有着或多或少的差异，目前缺少全国范围和更长时间跨度的流行病学研究，在人群中的流行模式也没有形成一套完整的理论，对HFMD进行可靠的预测，对于HFMD的科学防控具有重要的指导意义。HFMD资料是一个与时间有关的序列，利用这一资料，建立时间序列的模型，可以较好地对其发病进行预测和监测。及时发现异常，并采取相应的控制措施，可以达到使HFMD的危害降到最低的目的。本文采用自回归移动平均模型(autoregressive integrated moving average model，ARIMA)对全国HFMD的发病数进行时间序列分析并建立预测模型，为HFMD预防控制提供科学依据。 1材料与方法 1.1资料来源

2009年1月1日至2013年6月30日HFMD监测数据：来源于“中国疾病预防控制信息系统”中的“疾病监测信息报告管理系统”(“传染病疫情信息网络直报系统”)，以发病日期为标准日期，以月为基本时间统计单位计算各年HFMD发病数。 1.2方法

数据存储和数据管理采用Excel 2010、SPSS 19.0、SAS 9.2软件；应用SPSS 19.0软件对资料进行统计分析。 1.3时间序列模型

ARIMA的建模过程按以下阶段进行^{[1, 2]}：①序列特征及平稳化：应用时间序列图、自相关系数函数图(ACF图)、偏自相关系数函数图(PACF图)分析序列特征(趋势性、季节性)及平稳性。若序列平稳，进行平稳序列的时间序列分析； 若序列不平稳，结合本文序列特征，采用了普通差分和季节性差分等合适的数据预处理方法，实现序列平稳化后进行分析。本研究中对发病数序列进行1次差分和1次季节性差分，使最终所分析的序列为平稳化序列。②模型的识别：根据ACF图、PACF图和SPSS模型拟合结果，进行模型的初步识别和定阶。③参数估计和模型诊断：利用非线性最小二乘法估计模型参数。参数估计后，对模型的残差序列进行白噪声检验，用以判断ARIMA的适合性。对于残差而言相互独立就为白噪声，原假设是数据是相互独立分布的，如果P足够大则接受原假设，认为残差相互独立，为随机性误差，说明ARIMA模型比较适合，不需要重新选择模型。依据赤池信息量准则(Akaike information criterion，AIC)和贝叶斯信息准则(Bayesian Information Criterions，BIC)确定模型阶数，如不恰当则返回第二阶段，重新选定模型。通过不断改变模型的阶数并调试比较，确定残差序列为白噪声，AIC与标准化BIC值最小且模型较为简洁的为最佳模型。④预测应用：根据HFMD春季高发，以及ARIMA短期预测效果好的特点，用选定的模型对2013年1 6月HFMD发病数做预测并与实际发病数进行比较，以此验证模型的效果。 2结果 2.1HFMD流行趋势及空间分布分析

对全国资料2009年1月1日至2012年12月31日的HFMD发病数绘制时序分布图(图 1)。由图 1可以看出，HFMD 全年各月均有发病，从2010年开始，发病数呈现上升趋势。每年发病数从4月开始上升，5 6月出现峰值，次年1 2月发病数最低。2011年12月出现小高峰。以春季发病为主，具有明显的季节性。

图 1 2009年1月31日至2012年12月31日 全国HFMD发病数时序分布 Figuer 1 Time distribution of HFMD cases in China， 1 January 2009-31 December 2012

图选项

2009年1月1日至2012年12月31日全国HFMD发病数序列存在明显的季节性趋势，季节性周期为12个月。为不平稳序列，需进行差分处理使其平稳化。对发病数序列进行1次差分和1次季节性差分，序列做自相关分析图(即ACF图，图 2)和偏相关分析图(即PACF图，图 3)，ACF图显示相关系数在延迟数目>2的时候均落入置信区间并逐渐趋于0，说明该时间序列具有平稳性。PACF图显示偏自相关系数序列均接近0，只有两个相对偏离，可初步判断该序列用于二阶移动平均模型。

图 2 原始时间序列经一阶差分和一阶季节性差分后的ACF图 Figuer 2 Time series after the first order difference and first-order seasonal difference of ACF

图选项

图 3 原始时间序列经一阶差分和一阶季节性差分后的PACF图 Figuer 3 Time series after the first order difference and first-order seasonal difference of PACF

图选项

根据数据季节性的特征以及对PCAF图的分析初步判断模型为ARIMA(p，1，2)(P，1，Q)₁₂，模型的参数采取从低阶到高阶逐个进行尝试的办法。备选模型的拟合优度检验经过多次试验，根据模型的检验结果和参数之间的相关性进行调试、检验，选取AIC值和标准化BIC值最小和R²最大的模型为最佳模型。模型ARIMA(0，1，2)(0，1，0)₁₂的R²=0.784表明方程拟合度很好。AIC值=10.288标准化BIC值=21.955，在拟合的所有模型中最小。采用Ljung-Box方法检验残差白噪声，残差序列的自相关系数及偏相关系数均在95%CI以内(图 4)，残差为白噪声；Ljung-Box统计量为11.074，P=0.805，表明残差为随机性误差；且模型中的参数检验具有意义(表 1)。因此，ARIMA(0，1，2)(0，1，0)₁₂模型为最优模型。

图 4 模型拟合优度检验的残差自相关、偏自相关系数 Figuer 4 Coefficient of residual ACF and residual PACF by model test of goodness of fit

图选项

按照ARIMA建模方法，对全国2009年1月1日至2012年12月31日的HFMD发病数建模，并用该模型对全国2013年1月1日至2013年6月30日发病数进行预测，得到预测数据(表 2)，预测的平均 相对误差为0.52。模型拟合图(图 5)可 见预测值与观测值相比虽然趋势相同但是具体预测值与观测值相差较大。可见ARIMA(0，1，2)(0，1，0)₁₂模型的预测效果一般。

图 5 ARIMA(0，1，2)(0，1，0)12模型拟合图 Figuer 5 Fitting of ARIMA(0，1，2)(0，1，0)12 model

图选项

对全国31个省直辖市的2009年1月1日至2012年12月31日的HFMD资料分别绘制发病数时序分布图，发现2种明显的发病模式。北京、天津、河北、山东、辽宁等15个北方省市呈现明显的单峰模式(即不明显双峰模式，见图 6)，而上海、江苏、浙江、安徽、福建等16个南方省市则呈现明显的双 峰模式(图 7)。由图 6可见除了2009年3 7月间出现2个明显高峰外，2010 2012年3年均在5 6 月呈现一个明显的高峰，这与全国情况相同。图 7则 显示2009 2012年每年均在4 5月和11月出现2个明显发病高峰，特别是2011年和2012年越发明显。

图 6 全国HFMD发病数时序分布的单峰模式(2009.1.1-2012.12.31) Figuer 6 Single peak pattern of HFMD incidence

图选项

图 7 全国HFMD发病数时序分布的双峰模式 (2009.1.1-2012.12.31) Figuer 7 Double peak pattern of HFMD incidence

图选项

以这两种发病模式分层，按照上文叙述的分析过程分别建立ARIMA并对全国2013年1月1日至2013年6月30日的HFMD发病数进行预测。结果得出2个模型：北方城市的单峰模式，即不明显双峰模式最优模型为ARIMA(0，1，0)，(1，1，0)₁₂(模型1)，南方城市双峰模式最优模型为ARIMA(0，0，0)，(1，1，0)₁₂(模型2)，见表 3。结果可见，将全国资料分为两层分别得出2个模型，综合对全国2013年1月1日至2013年6月30日的HFMD发病进行预测。平均相对误差为0.12，比全国只用一个模式拟合时的相对误差明显减小。证明因地区发病的变异和模式不同在中国应用分层的ARIMA预测HFMD发病数效果更好。

传染病受到众多不确定因素的影响，这使得在传染病预测预警上存在着诸多的困难，导致了传染病早期预防控制工作一直较为滞后。因此如何及时有效地预测预警传染病的暴发流行一直是传染病预测控制工作 的重点^{[3, 4, 5]}。目前，已有诸多数学 模型运用到传染病的预测预警当中。但当在预测处于变化中的传染病疫情时，疫情数据呈现出较大的波动，降低了其反映真实暴发流行情况的能力，进一步使模型的预测效能降低^[6]。时间序列模型克服了影响预测对象的因素错综复杂或有关影响因素数据资料无法获得的难题^[7]，将复杂因素的综合效应统一蕴含到时间变量中，建立时间序列模型，以达到预测的目的。ARIMA是时间序列分析中重要而基本的模型之一，能够较好地描述时序的变动情况，特别适用于时序规律不明显，或有明显季节性和周期性的情况^[8]，广泛应用于传染病发病率的预测，特别是针对有季节性变动的时间序列提出的建模方法^{[9, 10]}/sup>。

本研究通过对2009 2012年HFMD发病情况的分析来看，全国每年HFMD发病数从4月开始上升，5 6月出现峰值，次年1 2月发病数最低。以春季发病为主，具有明显的季节性。北方省市和南方省市出现两种不同的发病模式。北京、天津、河北、山东、辽宁等15个北方省市呈现明显的单峰模式，即不明显双峰模式。北方省市 2009年3 7月见出现2个明显高峰，2010 2012年均在5 6月呈现一个明显的高峰。而上海、江苏、浙江、安徽、福建等16个南方省市则呈现明显的双峰模式，即2009 2012年每年均在4 5月和11月出现2个明显发病高峰，特别是2011年和2012年越发明显。可能与南北方气候有关^[11]。本研究根据实际发现的全国HFMD发病情况进行了分层ARIMA分析方法的应用，针对北方省市和南方省市分别建立模型并预测，平均相对误差仅为0.12，比全国只用一个模式拟合时的相对误差(0.52)明显减小。不同模型预测值与实际观测值的比较也可以明显地看出，全国分层建立ARIMA后的预测值比全国ARIMA预测值更接近于实际的观测值。分层建立ARIMA的分析方法不仅对HFMD的研究有重要的意义，而且对其他传染病的时间序列分析甚至对整个传染病防控工作有重要的提示作用。

ARIMA适用于具有典型趋势特征变化数据的预测^[12]。ARIMA要求时间序列长度在30个数据以上^[13]，本文以2009 2012年共48个月的发病数为基础建立模型，并对2013年上半年6个月的发病数进行预测并与实际发病数进行对比。ARIMA对2013年1 6月进行预测发现2013年5 6月仍为发病高峰期；与2013年实际发病情况一致。说明ARIMA对未来6月内的预测值在HFMD发病情况未发生显著变化时，能较好地预测HFMD未来的走势及发病数；在发病情况发生显著变化时(如发生暴发或流行时)，即可直观地判定出疫情走势的变更，并通过相关的信息系统发出预警信息，为尽早介入并采取针对性的防控措施提供有效的科学依据^[14]。以往文献多在地区范围进行ARIMA建模^{[15, 16]}，或者HFMD以外的其他传染病在全国范围进行建模^[17]，而本研究在全国范围，以4年共48个月的HFMD发病数为基础建立ARIMA预测预警模型，在国内众多基于ARIMA的文献中实属罕见。本研究另一个重要的创新点在于首次应用分层建立ARIMA的方法对全国水平进行预测，这对其他传染病在全国水平上的建模及预测预警起到了重要的提示和示范作用。综上所述，对于全国HFMD甚至其他传染病发病的预测可以根据不同地区的发病模式分层后分别估计出最优ARIMA，再综合几种模型共同对全国进行预测，并根据预测发出预警，具有加强传染病防控工作的实际应用价值。