1. 引言
2019年底暴发的呼吸道流行传染病迅速传播至全球,对人类健康与社会运行造成了前所未有的冲击[1] [2]。疫情不仅对医疗资源和公共卫生体系提出了巨大挑战,也对社会经济发展和政策决策产生了深远影响。在此背景下,如何准确刻画疫情传播机制、预测疫情发展趋势以及评估防控策略的有效性,成为流行病学与数理建模领域的重要课题。
在诸多方法中,基于常微分方程的传染病动力学模型因其结构清晰、可解释性强,被广泛应用于疾病传播研究和防控决策中。自Kermack与McKendrick [3]提出SIR模型以来,大量扩展模型被提出,如考虑潜伏期的SEIR模型、引入隔离人群或住院人群的多群体模型[4] [5]。这些模型能够揭示疾病在人群中的演化规律,并通过参数估计方法对传播率、恢复率等关键参数进行量化,从而在公共卫生干预中发挥了重要作用[6]。
然而,将上述动力学模型直接用于现实疫情分析仍面临若干挑战。数据稀缺性会导致参数识别困难,传统建模往往假定可以获得潜伏者、感染者、康复者等多类群体的时间序列数据,而在实际防控中,受限于检测能力、报告口径与社会条件,长期稳定可得的数据通常只有每日新增确诊和累计确诊病例[7] [8]。在观测量极为有限的情况下,模型参数之间容易出现可辨识性不清晰,传播率、恢复率等关键参数难以被唯一识别,从而削弱了模型的预测能力和政策指导价值。时变参数估计具有较高复杂性,防控措施调整、病毒变异以及人群行为改变会显著影响传播动力学特征,需要模型能够刻画传播率等参数随时间的动态变化。在观测数据有限且噪声较大的条件下,对时变参数进行稳定、精确的反演尤为困难。模型可解释性与机器学习方法之间存在矛盾。深度学习模型虽然在拟合能力上具有明显优势,但其“黑箱”特性较强,可解释性相对传统机理模型偏低,在公共卫生领域的接受度和推广应用因此受到一定限制。
近年来,物理信息神经网络(Physics-Informed Neural Networks, PINNs)为应对“小样本、高不确定性”的建模难题提供了新途径。PINN的核心思想是将动力学方程残差直接融入损失函数,使神经网络在拟合有限观测数据的同时,也尽可能满足机理模型所描述的动力学约束[9] [10]。该框架由Raissi等[11]首先提出,并已在流体力学、材料科学和生物系统等领域得到广泛应用。随着研究的深入,PINN逐渐被引入传染病动力学建模。Subramanian等[12]利用PINN对SIR模型进行反演,在有限病例数据下恢复了传播率的时变规律;Yang与Karniadakis [13]则在COVID-19场景中展示了PINN在数据稀缺条件下对关键参数的有效估计。尽管如此,多数现有工作仍依赖多类观测量(如确诊、康复与死亡人数等),并常将模型参数视为常数或简单分段常数[14] [15],对单一数据源下的参数识别和复杂时变过程的刻画仍显不足。
针对上述问题,本文在经典SIR模型的基础上提出了一种融合PINN的传染病动力学参数估计新框架(DNN-PINN-SIR)。该方法仅以每日新增和累计确诊病例为观测量,将SIR动力学方程的残差显式嵌入神经网络的训练过程,在保持物理一致性的前提下,实现了对系统隐含状态轨迹以及时变传染率以及常数恢复率的联合估计。通过在损失函数中引入方程约束,该方法一方面有效缓解了数据稀缺导致的参数识别不适定问题,另一方面将传播率和恢复率等关键参数建模为时间函数,更加贴切地反映了防控措施与社会行为变化对传播过程的动态影响。同时,相较于纯“黑箱”的机器学习方法,本框架以内生的SIR机理模型为基础,因而保留了良好的模型可解释性。
本文的主要工作和章节安排如下:第二节介绍所提出模型的构建及PINN网络结构,第三节给出基于南京疫情数据的数值实验过程与结果分析,第四节总结全文并展望后续研究方向。
2. 融合物理信息神经网络的呼吸道流行传染病传播动力学建模
2.1. SIR动力学模型
设研究总人口规模为N,将人群划分为三个舱室:易感者
、感染者
、康复者
,并进一步引入累计感染量
,用于刻画累计感染规模。因本次疫情持续时间较短,模型中忽略出生、自然死亡和其他人口的迁入迁出,总人口在整个研究时段内保持恒定:
(2.1)
动力学模型如下:
(2.2)
其中
为随时间变化的传染率,
表示从疫情开始累计感染人数。在离散观测时刻
模型给出的理论新增病例数可表示为
(2.3)
在模型中,有效再生数定义为
(2.4)
表示在时刻t的实际防控和易感人群水平下,一名感染者在其感染期内平均能够造成的二代病例数。
2.2. DNN-PINN-SIR网络架构
为了在观测数据相对稀缺的条件下同时重建SIR模型的状态轨迹并估计其时间依赖的参数,本文构建了由“状态子网络”和“参数子网络”组成的DNN-PINN-SIR (Daily New Cases-SIR-PINN)框架,并通过自动微分将SIR模型的微分方程显式嵌入损失函数,从而对模型训练过程施加物理约束。设研究时间区间为
,以时间
为唯一输入,构造两个相互独立的前馈神经网络:状态子网络和参数子网络。
状态子网络的输出对应模型中的所有状态变量,用于对观测数据进行拟合,并通过自动微分进入ODE残差项,记为
(2.5)
参数子网络用于识别时变传播率
,其输出进入ODE系统中,记为
(2.6)
其中
是由网络的权重和偏置构成的参数集合。
Figure 1. Architecture of the DNN-PINN-SIR network.
图1. DNN-PINN-SIR网络结构图
在DNN-PINN-SIR框架中(图1),我们将损失函数划分为三部分组成:初值约束损失、数据拟合损失和物理残差损失:
(2.7)
为了保证系统初值的准确性,我们利用归一化后的初始状态向量
对神经网络在
时刻的输出施加约束。初值约束损失定义如下:
(2.8)
其中,
代表各仓室的初始值。
数据拟合部分包括对每日新增病例与累计确诊病例的拟合。记神经网络在时刻
输出的累计感染人数为
,则由神经网络输出得到的新增病例数为
(2.9)
在给定的观测时间节点
上,我们希望网络输出尽可能逼近真实数据。数据拟合损失定义如下:
(2.10)
其中
和
为权重系数,用于平衡两类数据在损失中的贡献。
为将动力学模型信息显式嵌入神经网络,我们要求神经网络在时间节点
上尽量满足给定的ODE系统,物理残差损失定义如下:
(2.11)
其中
为第
个方程在时刻
的残差:
(2.12)
3. 数值实验与结果分析
本章利用南京2021年7月呼吸道传染病区域性暴发为例对所提出的DNN-PINN-SIR模型进行数值验证和性能评估。DNN-PINN-SIR中神经网络结构与BP神经网络基本一致,包括输入层、隐藏层和输出层。具体相关超参数设置见表1,各仓室初始值见表2。
Table 1. Hyperparameter values of the DNN-PINN-SIR network architecture
表1. DNN-PINN-SIR网络结构超参数取值
参数 |
神经
网络深度 |
神经
网络宽度 |
神经
网络深度 |
神经
网络深度 |
迭代次数 |
学习率 |
取值 |
5 |
32 |
3 |
64 |
1 × 106 |
0.001 |
Table 2. Parameter ranges and initial values for each compartment of the DNN-PINN-SIR model
表2. DNN-PINN-SIR参数范围及各仓室初始值
参数 |
|
|
|
|
|
|
含义 |
时变传染率 |
恢复率 |
总人口 |
易感者 初始值 |
感染者 初始值 |
康复者 初始值 |
取值 |
图6 |
0.1 |
9,314,685 |
9,314,676 |
9 |
0 |
3.1. 南京市2021年7月-8月呼吸道传染病区域性暴发及相关数据
本文以南京市2021年7月20日暴发的呼吸道传染病为研究对象,所用数据来自南京市卫生健康委员会发布的官方疫情通报[16]。根据通报信息整理得到连续24天的每日新增确诊病例数及其对应的累计确诊病例数,记为观测序列
(3-1)
分别为每日新增与累计确诊观测序列(图2)。
Figure 2. Observed data of respiratory infectious diseases in Nanjing from July 20 to August 12, 2021
图2. 2021年7月20日至8月12日南京市呼吸道传染病的实际数据
3.2 损失历史与收敛性
通过最小化损失函数,可以优化神经网络中的全部可训练参数。训练过程中采用“两阶段”优化策略:首先使用自适应学习率算法Adam [17]进行预训练,以获得较为理想的参数初值;在此基础上进一步引入高精度二阶拟牛顿优化器L-BFGS-B [18]对网络进行精细化调优。该策略综合了两类优化算法的优点:Adam对初始值相对不敏感、对学习率和初始参数的选取不敏感,但单独使用时收敛到高精度解的速度有限;而L-BFGS-B在已有较好初始点的条件下收敛速度更快、可达到更高的优化精度。网络权重参数采用Xavier随机初始化方案。
在具体训练过程中,当总损失降至预设阈值以下或迭代次数达到设定上限时终止迭代。图3给出了总损失及各组成部分随迭代步数变化的曲线。可以观察到,各项损失整体呈单调下降并最终趋于平稳,说明所采用的训练策略具有良好的收敛性和数值稳定性。
Figure 3. Evolution of the loss terms during training of the DNN-PINN-SIR model as a function of the iteration steps (Total loss denotes the overall loss; Loss U0, Loss U, and Loss F correspond to the losses of the initial-condition constraint, the data-fitting term, and the physics-informed constraint term, respectively)
图3. DNN-PINN-SIR模型训练过程中各项损失随迭代步数的变化情况(Total loss表示总损失,Loss U0、Loss U和Loss F分别对应初始条件约束、观测数据拟合项和物理约束项的损失)
3.3. 最小二乘法与DNN-PINN-SIR的拟合结果对比
鉴于本文仅使用每日新增与累计确诊两类观测数据量且研究窗口较短,若同时将传播率
与恢复率
设为时变函数,容易产生较强的参数相关性并降低可辨识性。为缓解该不适定问题并保持训练稳定性,本文将恢复率
近似为常数。
假设时变传播率为常数,因本轮南京疫情的康复天数为10天,所以两种方法中的恢复率
均设为0.1。通过采用传统最小二乘法对传播率进行分阶段拟合,拟合结果估计得到的每日新增确诊人数和累计确诊人数如图4所示。对比发现传统最小二乘法所估计得到的常数传播率无法精准地刻画疫情的传播状态,对于疫情感染峰值等时间点无法精确捕捉。
在DNN-PINN-SIR模型中,各项损失函数中的权重占比均为1。图5显示了南京疫情数据中每日新增确诊病例曲线及DNN-PINN-SIR模型的拟合结果。从图5(a)可以看出,模型模拟得到的日新增病例轨迹与观测数据整体吻合较好:一方面,模型较准确地捕捉到了疫情暴发初期新增病例曲线的单峰形态及峰值出现时间;另一方面,在疫情后期回落阶段,模型给出的下降趋势平滑、幅度合理。仅在个别病例数较少、数据波动较明显的早期时段,模型存在轻微的高估或低估,但偏差幅度较小,对整体传播趋势的判断影响有限。图5(b)给出了累计确诊病例的观测数据和模型拟合曲线。可以看到,相比新增病例数据,累计确诊曲线与模型结果的重合程度更高,除疫情初期存在极小偏差外,整个研究时段两条曲线几乎完全重合。这表明DNN-PINN-SIR模型在累积感染规模的刻画方面具有较高精度,能够较为真实地再现本轮疫情的整体传播过程。
(a) (b)
Figure 4. Comparison between the observed data and the piecewise least-squares fitting results for daily new and cumulative confirmed cases. (a) daily new confirmed cases; (b) cumulative confirmed cases
图4. 每日新增与累计确诊病例的观测数据及使用最小二乘法分阶段拟合结果对比。(a) 每日新增确诊病例;(b) 累计确诊病例
(a) (b)
Figure 5. Comparison between the observed data and the DNN-PINN-SIR model fitting results for daily new and cumulative confirmed cases. (a) daily new confirmed cases; (b) cumulative confirmed cases
图5. 每日新增与累计确诊病例的观测数据及 DNN-PINN-SIR 模型拟合结果对比。(a) 每日新增确诊病例;(b) 累计确诊病例
图6给出了基于DNN-PINN-SIR模型识别得到的时变传染率
的估计结果。从图中可以看出,疫情暴发初期
处于较高水平,说明人群接触频繁、传播风险较高。随后,
在数日内迅速下降,并在7月下旬至8月上旬呈现出明显的阶梯式衰减趋势,最终趋近于一个较低且稳定的水平。
整体来看,
的变化过程与日新增病例数的演化高度一致:传播率的快速下降对应了新增病例在实施病例隔离、密接追踪和社会面管控等防控措施后出现的明显回落,反映出相关干预在抑制实际有效传播率方面发挥了显著作用。说明基于DNN-PINN-SIR得到的时变传染率能够合理刻画本轮疫情不同时段的传播强度变化,为评估防控措施效果提供了定量依据。
Figure 6. Time-varying transmission rate
identified by the DNN-PINN-SIR model over time
图6. DNN-PINN-SIR模型识别得到的时变传染率
随时间的变化曲线
图7给出了本次疫情的有效再生数
随时间的变化轨迹。可以看到,2021年7月21日左右
达到约7.3的峰值,表明疫情暴发初期传染性极强,单个感染者在其感染期内平均可引起多例继发感染,疫情处于快速扩增阶段。此后,
随时间明显下降,并于7月31日首次降至1以下,标志着本轮疫情由扩张阶段转入相对受控阶段,新发病例的增长动力开始减弱。
Figure 7. Temporal evolution of the effective reproduction number,
, computed from the identified time-varying transmission rate
and the susceptible population size; the dashed line indicates the critical threshold
图7. 基于识别得到的
及易感人群规模计算的有效再生数
随时间的变化曲线,其中虚线表示临界水平
在7月底以后,
始终维持在1以下,说明在持续防控措施作用下,病毒传播链条被有效压缩,疫情整体进入衰减过程。由此可见,基于模型估计得到的有效再生数不仅刻画了疫情从暴发到受控再到消退的阶段性演化特征,也从传播动力学角度印证了各项防控措施的实际干预效果。
3.4. ODE反演验证
为进一步检验DNN-PINN-SIR所估计参数的物理一致性和泛化能力,本节将识别得到的时变传染率
代入SIR动力学模型中,得到相应的状态轨迹和预测新增病例曲线。图8给出了由SIR方程前向模拟得到的每日新增确诊病例曲线与观测数据的对比情况。
Figure 8. Comparison between the observed daily new confirmed cases and the trajectories obtained by solving the SIR equations with the identified time-varying parameters
图8. 利用识别得到的时变参数
代入SIR方程得到的每日新增确诊病例轨迹与观测数据的对比
从图8可以看出,ODE模拟得到的新增确诊病例在峰值大小、峰值出现时间以及后续衰减趋势上均与真实数据高度一致,仅在个别日期存在小幅偏差。总体而言,两条曲线在整个研究时段内拟合良好,表明DNN-PINN-SIR识别出的
不仅在神经网络内部能够有效拟合观测数据,而且在作为参数代入独立的SIR方程进行求解时,同样能够再现实际疫情的传播过程。由此说明,所估计的时变参数与机理模型之间具有良好的一致性,模型具有较强的物理可解释性和外推能力。
4. 结论与展望
本文针对仅掌握每日新增与累计确诊病例、且观测时长较短的呼吸道流行传染病场景,构建了一种融合物理信息神经网络(PINN)的传播动力学分析框架。该方法以每日新增确诊病例为主要观测量,将SIR模型的常微分方程残差显式嵌入深度神经网络的损失函数中,在有限样本条件下同时实现对隐含状态轨迹以及时变传染率
的估计。
基于南京市某轮呼吸道流行传染病的实证分析表明:在仅有约24个观测点的小样本条件下,所提出的PINN模型仍能对每日新增和累计确诊曲线给出高精度拟合。相对比传统最小二乘法估计得到的接触传染率结果,DNN-PINN-SIR模型识别得到的
与时间序列具有清晰且合理的演化规律,能够反映防控措施升级、社会接触模式变化等外部因素对传播动力学的影响,为疫情扩散过程的定量解析提供了新的视角和工具。将PINN反演得到的参数代回SIR方程并进行ODE数值求解,其得到的解与PINN直接输出的状态轨迹高度一致,说明估计参数与机理模型之间具有良好的自洽性。
未来工作可从以下几个方向进一步拓展。模型结构拓展:将当前框架推广至更复杂的多群体或年龄分层模型,引入无症状感染者、住院者等多个舱室,以更精细地刻画现实疫情传播过程。结合人群流动、核酸检测、接触追踪等多源异质数据,构建多模态PINN模型,以提升参数估计精度和中长期预测能力。将该框架应用于其他新发或再发传染病的早期暴发情形,为“数据稀缺但决策紧迫”的公共卫生场景提供具有通用性的建模与分析工具。
NOTES
*通讯作者。