猩红热是儿童常见流行的呼吸道类传染病，由产毒素的化脓性链球菌（也称为A组链球菌，或GAS）感染引起^[1]，主要临床表现有发热、咽峡炎、弥漫性皮疹，典型的皮疹为红色、针尖状、砂纸质地，从躯干蔓延并固定在弯曲处^[2]。猩红热在全球均有流行，给人群健康和生命安全造成极大威胁。英国1994—2000年，猩红热发病率从250/10万下降到小于5/10万，但2014年发病水平激增，报告了1 500余例猩红热病例，为45年来最高发病水平^[2]。中国在2011年发病率也迅速上升^[3]，同年，香港地区有2例儿童因此而死亡^[4]。目前尚无疫苗能够有效预防猩红热感染，猩红热的流行仍然是全球面临的重要公共卫生问题。

传染病预测有助于及时研判传染病流行趋势，准确预测是制定传染病防控策略的重要依据，对于传染病控制意义重大。近年来，深度学习逐渐被应用在时间序列的分析预测方面，其中循环神经网络 （recurrent neural network，RNN）表现较好；而基于RNN的长期短期记忆网络（long⁃short term memory⁃ neural network，LSTM）解决了RNN在分析预测时序数据中出现梯度消失、梯度爆炸及长期记忆能力不足的问题，能更为有效地处理分析较长的时间序列数据，在环境、医疗等领域的时间序列预测中有着广泛的应用^[5]，且逐渐应用于流感，新型冠状病毒肺炎等传染病的流行趋势预测^[6]。本文以江苏省2005—2019年猩红热逐月发病数据为基础构建LSTM模型，通过与SARIMA模型预测精度的比较，探讨LSTM模型在猩红热发病趋势预测中的应用价值。

1 资料和方法

1.1 资料

2005—2019年江苏省猩红热逐月发病数据来源于传染病报告管理系统。同期气象数据下载自国家气象科学数据中心（http：//data.cma.cn/），包括降雨量、降雨天数、大气压、水汽压、相对湿度、最低气温、最高气温、平均气温、日照时长等。江苏省人口基数大，常住人口数相对稳定，故本文使用猩红热发病数代替发病率以反映流行情况。

1.2 方法

1.2.1 差分自回归滑动平均（seasonal autoregressive integrated moving average，ARIMA）

ARIMA模型即Box⁃Jenkins模型，已经被广泛应用于各类传染病的时间序列分析。多项研究表明，许多传染病时间序列都包含明显的周期性和季节性趋势^[7-8]，猩红热的发病趋势也显示出明显的季节性，在模型中加入明确的季节性趋势会使模型的表现更为良好，因此将ARIMA模型扩展为带季节组分的SARIMA模型。SARIMA模型可以表示为SARI⁃ MA（p，d，q）（P，D，Q）s，其中p、d、q分别是自回归阶数、差分阶数和移动平均阶数；P、D和Q分别是季节自回归阶数、季节差分阶数和季节性移动平均的阶数；S是季节周期的长度^[9]，本研究中定义为12个月。模型拟合过程如下：

首先通过单位根平稳性检验（augmented dickey fuller，ADF），判断猩红热发病序列数据是否平稳。如果不平稳，则采用差分的方法将时间序列转化为平稳序列。然后基于平稳时间序列的自相关函数 （AUTOCORRECTION FUNCtion，ACF）图和偏自相关函数（partial autocorrection function，PACF）图判断模型可能阶数，进而建立备选模型。根据赤池信息准则（Akaike information criterion，AIC）以及贝叶斯信息准则（Bayesian information criterion，BIC），从备选模型中选择AIC和BIC最小的模型作为最优模型^[10]。

1.2.2 LSTM模型

LSTM将传统神经网络的神经元转化为储存器单元的结构，这个结构能够将时间序列中的隐藏信息储存并传递（图1）。储存器单元包含3个门控装置矩阵：输入门、遗忘门和输出门，由一条贯穿时间步的细胞状态参与决定信息的存留^[5]。

第一步，需要决定要从细胞状态中丢弃的信息，由“遗忘门”的sigmoid函数实现（公式1）。它通过 h_t-1（前一个输出）和 x _t（当前输入）合并后，通过sigmoid函数输出值为0和1之间的矩阵，1代表信息被完全保留，而0代表信息被彻底删除。

其中，f_t表示 t 时刻遗忘门矩阵；σ表示sigmoid函数（公式2）；h_t表示 t 时刻LSTM神经网络神经元的输出；x_t表示t时刻LSTM神经网络神经元的输入； W表示权重矩阵，b是门的偏差，其下标表示由这个矩阵连接的变量。

第二步，需要决定在细胞状态中存储和更新的信息。首先，“输入门”的sigmoid函数决定将要更新的值（公式3）。然后用tanh函数将创建候选细胞状态向量${\tilde{C}}_t$，该向量用于细胞状态的更新（公式4）。最后需要将上一个细胞状态值C_t-1更新为C_t，将上一个状态值乘以f_t，代表需要忘记的部分，之后将得到的值加上i_t×${\tilde{C}}_t$，得到新的细胞状态C _t（公式5）。

其中，i_t表示t时刻输入门矩阵；${\tilde{C}}_t$表示t时刻候选神经元细胞状态矩阵；C_t表示t时刻的神经元细胞状态矩阵；tanh表示tanh函数（公式6）。

第三步，需要决定输出的值。首先，运行sig⁃ moid函数，建立“输出门”，它决定细胞状态输出的部分（公式7）。然后，将细胞状态通过tanh函数将值规范化到-1和1之间，并乘以“输出门”矩阵（公式8），至此只输出了模型决定输出的部分^[11]。

其中，o_t表示t时刻输出门矩阵。

训练LSTM神经网络模型：

首先，根据实际问题确定输入输出张量的维数，分别构造训练集、测试集；对训练集、测试集的样本输入进行归一化以消除不同指标量纲的影响，具体计算方法如下（公式9）：

其中，x表示时间序列X中的一个数据；x为x进行归一化后的数据；X_max为时间序列X中的最大值； X_min为时间序列X中的最小值。

其次，确定网络结构：设定LSTM网络层数、全连接层数及各层神经元个数，设定损失函数以及全连接层随机舍弃率，设定模型迭代的批量以及周期数。

最后，用训练集训练LSTM神经网络模型。

1.2.3 模型评价

为了评价模型预测精度以及其可扩展性，使用平均绝对误差百分比（mean absolute percentage er⁃ ror，MAPE）和均方根误差（Root Mean Squared Error， RMSE）作为评价标准。本研究以2005—2017年数据为训练集，2018—2019年数据为测试集，比较SARIMA模型及LSTM模型在拟合优度以及在外推预测效果上的优劣。MAPE和RMSE计算公式如下 （公式10、公式11）：

式中，y_i为第i个实际值，x_i为第i个预测值。

1.3 统计学方法

本文采用R4.0完成描述性统计分析，采用py⁃ thon3.8完成SARIMA和LSTM神经网络模型拟合和预测分析。

2 结果

2.1 江苏省2005—2019年猩红热流行特征概述

2005年1月—2019年12月，江苏省共计报告猩红热30 482例。猩红热月发病数存在明显季节性，且呈双峰，第一高峰出现在3—7月，第二高峰出现在11月至次年1月。2005—2013年间，江苏省猩红热发病趋势平稳，2014年出现异常高发状态，之后猩红热发病水平呈现逐渐上升的长期趋势（图2）。

2.2 SARIMA模型

经过Ljung⁃Box检验（P< 0.05），2005—2019年江苏省猩红热逐月发病数据是非白噪声序列，能够使用SARIMA模型提取其中的线性关系。

对2005—2017年江苏省猩红热逐月发病数序列进行ADF检验，Dickey⁃Fuller=-2.5436，P=0.349，存在单位根，该序列为非平稳序列。对原始序列进行一阶一般差分以及一阶季节性差分，差分后的序列经过ADF检验，Dickey⁃Fuller=-4.9266，P=0.01，提示差分处理后的序列为平稳序列。

通过猩红热发病时间序列数据平稳化的过程， d=1，D=1。绘制平稳化之后序列的自相关函数图 （ACF）和偏自相关函数图（PACF）（图3）。通过观察偏自相关图，可取p=2或3，P=1；通过观察自相关图，可取q=2或3，Q=0或1。由此可建立8个备选模型（表1）。

通过比较各个模型的AIC、BIC，并结合模型显著性检验结果，确定最优模型为SARIMA（3，1，2）（1，1，1）₁₂，AIC=1 511.45，BIC=1 535.15，白噪声检验结果表明模型的残差为白噪声序列。

2.3 LSTM神经网络

如图4所示，各气象因素之间存在很强的相关性。因子分析的充分性检验（Bartlett检验，χ²=4 461.37， P ＜0.05）表明气象因素适合做因子分析。为降低数据冗余对模型效率的影响，对气象因素进行因子分析，提取出众多气象因素的公因子。根据因子特征值及各个因子累积贡献率，选择前3个因子，特征值分别为6.14、1.92、0.33，累积贡献率达到93.19%。因子载荷热图见图5。

本研究以前3个月的逐月发病数及气象公因子作为模型输入，以当月发病数为模型输出。最终建立的LSTM模型由4层LSTM层和1层全连接层构成，每层含32个神经元。

2.4 模型比较

SARIMA（3，1，2）（1，1，1）₁₂在训练集和测试集上的MAPE分别为22.47%、35.97%，RMSE分别为43.64、227.85，LSTM模型在训练集和测试集上的MAPE分别为16.94%、22.46%，RMSE分别为19.30、 152.46。LSTM对原始数据的拟合优度以及外推预测精度均明显优于SARIMA（3，1，2）（1，1，1）₁₂。模型预测结果见图6。

3 讨论

猩红热是江苏省常见的呼吸道类传染病，近年来发病水平不断上升，且在全国范围内均有类似的流行趋势^[12]，猩红热防控仍面临重大挑战。建立准确的预测模型对于猩红热疫情防控以及医疗资源配置具有重要的指导意义。目前猩红热流行趋势预测的模型主要是ARIMA模型^[13]。ARIMA模型是针对一系列随时间变化，而又互相关联的动态数据建立的预测模型。为了能够发现数据之间的关联关系，ARIMA模型对数据的数量和数据的自相关性要求较高。SARIMA模型能够在ARIMA模型的基础上，针对发病具有季节性特点的疾病，提取季节滞后的自相关信息。但是，猩红热的发病受到多种因素，如平均温度、相对湿度等气象因素的影响^[14]，这会造成单因素模型预测效果的不理想，同时，在较长的时间跨度上，气象因素的变化将会对模型带来更大的影响，导致模型的长期预测效能下降。

神经网络模型以其优秀的非线性拟合能力、能够纳入多种预测因素等优势，在不同领域的时间序列预测中表现出明显的优越性。随着计算机软硬件技术的发展，深度学习技术也得到了飞速发展， RNN作为深度学习的网络模型之一，能够很好地处理序列数据，但是在较长时间序列的预测上，容易造成长期相关信息的遗漏。在RNN的基础上， LSTM神经网络模型引入门控系统，选择性地保留或者剔除在长时间步上的相关信息，能够捕获长期依赖关系，在时间序列的预测上表现出显著性能，因此在许多实际应用，如天气预报以及新冠肺炎的预测中发挥着重要作用^[15-16]。

本研究使用SARIMA模型以及LSTM神经网络模型对2005—2019年江苏省猩红热月发病数进行拟合预测。SARIMA模型表现较差，前瞻性预测误差超过了35%。LSTM模型前瞻性预测精度接近80%，预测效果明显优于SARIMA模型，预测结果能够很好的反映猩红热实际发病水平的变化趋势。建立多因素预测模型时，如何选择合适的预测指标是模型构建的难点之一。相关性分析结果表明各气象指标与猩红热流行之间均存在显著相关性，可以作为构建LSTM的预测指标。但是，气象因素之间存在高度相关性，如果全部纳入预测模型中，将降低模型的性能，影响模型预测效果。本研究通过因子分析的方法提取气象因素的公因子，以公因子作为预测指标建立LSTM模型，去除气象因素之间的多重共线性影响，对模型预测效果起到了很好的优化作用。相较于SARIMA模型，LSTM神经网络模型的拟合效果及可扩展性均优于SARIMA模型，主要体现在LSTM神经网络模型对数据条件要求较低，能够结合多种影响因素数据，有效拟合历史发病数据，并且能够根据即时输入数据的波动，准确地预测2018年和2019年发病序列中的异常高发值，这也要求LSTM模型在实际应用中，需要不断获取最新的数据资料，加入到已经建立的模型中，进行不断地学习修正，才能保证模型在各种复杂因素的影响下保持着较高的预测精度。

本研究结果提示LSTM神经网络模型在江苏省猩红热发病时间序列预测中表现良好，可以用于流行趋势研判和风险评估，为优化和调整猩红热监测、防控策略和措施提供依据。

参考文献