模拟新型冠状病毒肺炎无症状感染者传播风险*

孙廷哲，崔隽

1. 安庆师范大学生命科学学院，安徽安庆 246133

2. 中山大学生命科学学院，广东广州 510275

新型冠状病毒肺炎（COVⅠD-19，coronavirus disease 2019） 是对人类健康的巨大威胁［1］。研究新型冠状病毒肺炎的传播特征将有助于揭示有效的防控措施，遏制病毒扩散［2］。中国政府在疫情暴发初期即迅速采取有效防疫措施，成功遏制了疫情在国内的扩散。国内多省市在2020 年3 月已经实现新增确诊病例为零，表明疫情已被遏制。而在国外，由于没有采取积极的措施，疫情已近乎失控。截至欧洲中部时间2021年4月27日下午4时，全球已报道确诊病例147 539 302 例，死亡病例3 116 444 （https://covid19.who.int/）。而全球疫情最为严重的美国累计确诊病例达到31 742 914例，死亡566 842例。

新型冠状病毒无症状感染者是一类特殊群体，其无明显症状，但病毒载量与确诊病例相当，同时具备病毒传播能力［3］。另有研究表明年轻病人表现为重症的机会较低［4-5］。随着冬季来临和人员流动的增加，病毒传播风险加大，对精准防控提出了不小的挑战。数学模型对新型冠状病毒传播特征的研究起到了有效的辅助作用［6］。许多学者建立了流行病学模型并且对COVID-19传播特征进行了分析，对疫情发展进行了预测并提出了具体调控策略［7-9］。一些学者也通过模拟手段探索了病毒通过中间宿主的传播特征［10-11］。我们之前的研究发现，较之境外输入病例，无症状感染者诱发的疫情暴发速度更快、规模更大［12］，无症状感染者也具有传播风险。虽然确诊病例的清零动态得到了较为广泛的关注，但无症状感染者在人群中的清零时间不易被准确测量，因此需要用模型模拟手段予以估计。

为此，基于之前构建的包含无症状感染者的流行病学微分方程模型［12］，对国内流行病数据进行定量拟合，结果显示模型与实际流行病学数据拟合良好。模拟结果表明群体中无症状感染者的清零时间显著迟于确诊病例首次 “清零” 时间。若精准防控措施不到位、不能有效抑制人员流动和相互接触，那么人群中潜在的无症状感染者将很有可能引发第二波疫情。因此，群体中无症状感染者的延迟消除特征为疫情防控提供了有益的参考。

1 SCIRA模型构建

1.1 模型基础

在之前的工作中，我们已构建了用于描述新型冠状病毒传播的 “SCIRA” 模型［12］，此SCIRA模型基于之前的SCIR 模型［13］引入了无症状感染者。我们发现经典四变量模型无法很好地拟合国内流行病学数据，存在较大误差［12］。而引入无症状感染者后，模型可以较好地拟合国内流行病学数据［12］。模型构建简要阐述如下：基于经典SEIR模型基本假设［14］和SCIR 模型［13］对人群的划分，将人群分为易感人群（未被感染的健康人群，susception population，S），暴露人群（需要被密切健康监控人群，包括潜在的密切接触者，close contacts，C），确诊病例（diagnosed，I），康复人群（recovered，R） 和无症状感染者（asymptomatic patients，A）。由于国内多地进入重大突发公共卫生事件Ⅰ级响应（如江西省宣布紧急状态时间为2020 年1 月24 日），各省市间人员流动得到了显著抑制，因此不考虑各省之间输入病例的影响。无症状感染者表示个体虽被感染，但无明显症状，且具备传播能力［15］。COVID-19传播的重要途径是感染人群和易感人群之间的相互接触，易感人群（S） 与密切接触者（C）、确诊病例（I）、无症状感染者（A） 之间的接触都有可能造成病毒的传播［13］。模型用如下微分方程进行描述

模型参数中，λ为（潜在） 密切接触者（C）转化为易感者（S） 的平均速率。由于C包含潜在未被监测到的密切接触者人群，因此λ取值需要被重新估计，而不能简单设置为1/14 （14 d为平均隔离时间）。β1、β2和β3分别表示S和C、A或I群体之间的平均接触速率。S和C、A或I群体的接触将增加人群中的潜在密切接触者数量。ε表示无症状感染者自我康复速率［16-17］。ν1和ν2分别描述密切接触者（C） 转化为无症状感染者（A） 和确诊病例（I） 的平均速率。基于前模型假设，无症状感染者（A） 通过速率ν3转化为确诊病例（I）［16-17］。μ为确诊病例平均康复速率。考虑隔离留观及确诊病例的增多对有限医疗资源的限制作用，在康复过程中引入一种阈值行为（threshold behavior）。确诊病例的暴发式增长将导致医疗资源的严重不足或饱和运行，从而降低确诊病例的治疗和康复速率。由于各省份医疗承载能力不同，因此确诊病例数量对康复速率μ的影响亦不同。参数K描述此过程的阈值特征，n为协同性。

初始确诊感染者和康复者数量参照江西省卫健委（https://hc.jiangxi.gov.cn/） 发布数据设置（I=18，R=0）。在2020 年1 月21 日、22 日和23日，江西省分别报告新增确诊2 例、1 例和4 例。此7例确诊病例和1月24日报告确诊病例合并为初始值。2020 年1 月24 日之前，江西省未报告病例治愈或死亡，因此初始值R=0。初始时间设定为2020 年1 月24 日，模拟时间步长为1 d。初始密切接触者通过模型拟合取整设置为453 246。初始易感人群“S” 根据江西省常住人口数量设定（46 207 736，依2019 年数据）。依据之前的模型，死亡病例和康复病例进行合并处理［13，18］。动力学参数通过拟合流行病学数据确定（表1）。由于群体数目较大，故采用τ-leap 随机模拟算法［19］，随机模拟过程中保持总人口数恒定。常微分方程组通过MATLAB 软件（版本号R2018b） 的ode23s 算子进行数值求解。

表1 模型参数Table 1 Model parameters

1.2 参数拟合

参数拟合通过MATLAB 工具箱PottersWheel实现［20］。江西省的流行病学数据取自江西省卫健委官方网站（https://hc.jiangxi.gov.cn/），时间跨度 为2020 年1 月24 日～3 月31 日。数据记录 于Excel 表格中作为PottersWheel 工具箱的标准输入用于拟合模型参数。通过模型拟合确诊病例数量（数据个数N=68） 和康复者数量（N=68）。模型χ2值用于评估模型，χ2为模型和实际流行病数据之间的误差［20］。 依据Maiwald 等［20］的定义，若χ2/N< 1 表明模型可靠。 随机参数组拟合使用PottersWheel 工具箱中的Fit sequence analysis 功能实现（共400 组）。通过随机参数组拟合，选取最优拟合参数组作为当前模型参数，参数取值如表1 所示。

2 结果与分析

2.1 模型拟合流行病学数据

首先用模型拟合江西省新型冠状病毒肺炎流行病学数据（2020 年1 月24 日～2020 年3 月31 日）。模型参数拟合使用PottersWheel 工具箱［20］。自3 月12 日起，江西省未再次发现本土确诊病例，因此，确诊病例（I）全部治愈后其数量一致设置为0。由于严格的管控措施，境外输入性新型冠状病毒肺炎确诊病例（如3 月21 日和3 月28 日各一例） 不与易感人群形成有效接触、不与本土人群形成流行病学关联，因此不计入此模型。为避免参数拟合落入局部极小区域，通过400组随机初始参数作为拟合起点，获得最优拟合参数组（参数值见表1）。由于无症状感染者数量的初始拟合值接近0（0. 062 1），因此取整设置为0。模型拟合结果参见图1A（χ2/N=0. 086）和图1B（χ2/N=0. 041）。对确诊病例和累计康复患者的拟合χ2/N< 1，验证了模型的可靠性。 而使用经典SEIR 模型在拟合COVID-19确诊病例峰值时具有较大偏差［12］，因此引入无症状感染者可以显著提高模型对流行病学数据的拟合。注意到最优拟合参数βi（i=1，2，3）的取值的数量级维持在10-10～10-9，这表明在严密防控措施下（如社交隔离、居家隔离、正确佩戴口罩等），不同群体之间的有效接触显著降低（表1）。另外，n的最优拟合值接近1，表明医疗承载能力对康复速率μ的影响协同性较低。为研究COVID-19传播特性，根据Driessche等［21］的结论，计算模型的基本生成数R0，得

图1 2020年江西省流行病学数据拟合结果Fig. 1 The ODE model was trained by epidemic data from Jiangxi province in 2020

其中N=S + C + I + R + A。代入各参数取值得R0=0. 462 5 < 1。在严密防控措施状态下，模型的R0< 1，仅存在无病平衡点（DFE，disease free equilibrium），且DFE是局部渐近稳定的［21］。表明防控措施有效，感染者群体将最终消失。综上所述，模型可较好地拟合江西省新型冠状病毒肺炎流行病学数据且表明江西省防控措施有效遏制了病毒传播。

2.2 群体中无症状感染者完全康复时间有明显时滞

进一步采用τ-leap 随机模拟算法对模型进行模拟。结果显示确诊病例和无症状感染者的清零时间具有较大的差异（图2A）。注意到无症状感染者的完全治愈（即“清零”）时间较之确诊病例首次完全治愈时间具有显著的时滞（图2B）。1 000 组随机模拟结果显示，此时滞显著大于0 （Wilcoxon 符号秩检验，p=1. 106 2×10-24）。注意，核密度估计只在所有时滞项的最大和最小值之间进行，因此可发现最小时滞为0 d。部分时滞可以超过50 d（图2B）。为排除江西省流行病学数据的特异性影响，安徽和江苏省COVID-19流行病学数据也用模型进行了拟合（参数值详见文献［12］）。结果显示，无症状感染者完全治愈时间亦显著迟于确诊病例完全康复时间（图3A 和3C）。在安徽省和江苏省的模型中［12］，无症状感染者清零的时滞也显著大于0（ 安徽p=4. 421 7×10-24， 江苏p=6. 919 1×10-54） 中位数时滞分别达到15 d 和9 d（图3B和3D，虚线）。因此，无症状感染者的完全治愈和确诊病例首次清零时间相比具有一定时滞，也就是说群体需要更多时间 “清除” 无症状感染者。后续使用江西模型做进一步讨论。

图2 2020年江西省模型中确诊病例和无症状感染者动态变化Fig. 2 Temporal differences in complete recovery time for asymptotic and diagnosed patients in Jiangxi model in 2020

图3 2020年无症状感染者和确诊病例完全治愈时间对比图，模拟基于江苏省和安徽省模型Fig. 3 Temporal differences in complete recovery time for asymptotic and diagnosed patients in Jiangsu and Anhui model in 2020

2.3 基于模型估计无症状感染者的流行病学特征

基于江西模型的时序模拟显示，无症状感染者达到峰值的时间早于确诊病例（图4A）。基于1 000 组随机模拟结果显示，此时滞的中位数约为12 d （图4B）。进一步对群体中无症状感染者所占比例进行估计（即无症状感染者占总感染者比例，总感染者数量= 无症状感染者数量+ 确诊患者数量）。结果显示群体中无症状感染者比例随时间发生变化（图4C）。注意到确诊病例和无症状感染者在持续严格防控措施下最终降为零，因此晚期无症状感染者比例将不再计算。依照模型计算，无症状感染者或转归为确诊病例（ν3A），或经历自愈过程（εA）。因此，在模型中比较此两参数大小，即可对确诊转归和自愈过程进行比较。基于400组随机初始参数拟合流行病学数据，取前200组拟合较好的参数组进行比较。结果显示，ν3的取值显著大于ε（图4D）。因此，多数无症状感染者将最终转归为确诊病例，相对而言仅有少数无症状感染者可以通过自愈途径康复。以上分析揭示了无症状感染者较之确诊病例的传播特性差异。

图4 2020年无症状感染者新型冠状病毒肺炎传播特征分析Fig. 4 Epidemic features of asymptomatic patients in 2020

2.4 无症状感染者可诱导第二波疫情暴发

通过对模型的进一步模拟，评估确诊病例“清零” 后回归正常生活（即无严密防控措施） 的安全性。江西省确诊病例在3 月12 日 “清零”，实现全部治愈。3月12日以后，在模型中变化实施严格防控措施的持续时间。由于严格防控措施可有效降低人群间相互接触，因此，在无防控或弱防控措施条件下，人群间相互接触将显著增加。模型中βi（i=1，2，3）表示不同群体间相互作用，因此，基于之前模型假设［12］，将βi（i=1，2，3）乘以系数F。因此，βi′=F·βi（i=1，2，3）。以βi′替换βi代入模型来模拟调控措施放松条件下新型冠状病毒肺炎传播特征（此时F> 1）。以F=3为例模拟正常生活即未实施严密防控措施的状态。确诊病例首次归零时间记为D1（图5A，“确诊病例清零”）。之后，严格防控措施持续至D2日（D2日之前，F=1，图5A）。D2日后，防控措施削弱，各群体间相互接触显著增加（F=3，图5A）。因此，确诊病例清零后严格防控措施持续时间为D（D=D2-D1）。变化严密防控持续天数D（图5A），对D的每次取值，皆进行500 组随机模拟（每次模拟记录300 d数据） 并监测新型冠状病毒肺炎的二次暴发频率。模拟结果显示当确诊病例清零后即刻放松防控措施将不可避免地触发第二波疫情（D =0，图5B）。即使严密防控在确诊病例清零后持续30 d，仍有相当比例的随机模拟会触发第二波疫情（74/500，图5B）。当D＞40 d 以后，COVID-19 再次暴发的频率急剧下降；直至D≥55 d后，第二波疫情的触发频率降至4‰（2/500） 或以下（图5B）。这些模拟结果显示，确诊病例清零后，实施持续的严密防控措施对于控制新型冠状病毒肺炎的二次暴发是至关重要的。

图5 2020年确诊病例首次清零后持续严密防控对新型冠状病毒肺炎疫情二次暴发的影响Fig.5 Effects of duration for strict interventions on potentially secondary COVID-19 outbreak after diagnosed patients are all cured for the first time in 2020

3 讨 论

建立包含无症状感染者的微分方程模型，较之传统的SEIR 模型可以更好地拟合流行病学数据［12］。结果显示群体中的无症状感染者对新型冠状病毒肺炎的传播至关重要。最近的研究表明，新型冠状病毒肺炎无症状感染者的病毒载量甚至超过有症状感染者［22］，且无症状感染者也具有一定的病毒传播能力［23］。通过拟合新型冠状病毒肺炎流行病学数据，本文建立了数据驱动模型（datadriven model）。在之前的模型中，我们比较了无症状感染者和输入病例对COVID-19传播的影响，指出无症状感染者触发的疫情再次暴发速度更快、规模更大［12］。本研究中我们发现人群中无症状感染者的完全清除时间要显著迟于确诊病例清零时间，表明即使确诊病例清零，人群中无症状感染者仍大概率存在。另外，模拟结果显示无症状感染者的达峰时间（time to peak） 更短，相当一部分无症状感染者将转归为确诊病例。因此，无症状感染者的严密管控对遏制疫情暴发显得尤为重要。中华人民共和国国家卫生健康委员会自2020 年4月1日起将无症状感染者情况列入每日通报，这充分表明了监测无症状感染者的重要意义。另外，模拟显示无症状感染者占总感染者的比例呈动态变化，此结果可部分解释不同研究中所估算的无症状感染者比例差异较大的现象［24-25］。

在确诊病例清零之后，严密防控措施是否必要也在模型中予以探讨。随机模拟结果显示，在确诊病例清零后，即使严密管控措施持续30 d，群体中仍然有疫情二次暴发风险。在有效疫苗大规模接种前，只有进行长期的精准管控，才能有效地遏制新型冠状病毒肺炎疫情的蔓延。

另外，模型中使用参数F调控不同人群的相互接触速率。由于无症状感染者自愈和无症状感染者转归为确诊病例是一个动态过程，各群体的具体数量无法得到静态的、准确的估计。因此，对任意指定两群体间相互接触（即某一指定βi，i=1，2，3）的调控是不易实现的。较为可能的措施是同时对所有群体的相互接触进行严格调控，包括大规模COVID-19检测、集中隔离、社区管控以及科学佩戴口罩。这些措施将有效降低不同人群间相互接触即βi值。而管控放松的模拟采用增大F值的方法［12］，旨在强调严密管控的重要性。另外，当前模型仍具有一定的不足：模型并未区分潜伏期病例、隐性感染无症状感染者和潜伏期无症状感染者［26］；模型未描述公共卫生事件Ⅰ级响应之前COVID-19 传播特征；群体的自然出生和死亡亦未予考虑；进一步构建精细的模型将有助于解决上述缺陷。综上所述，确诊病例首次清零并非表征一种无病毒（coronavirus-free） 状态，即 “零病例不等于零风险”，持续的科学管控对于遏制疫情再次暴发仍然十分必要。