基于实时再生数的COVID-19预测模型及可视化平台开发

何沅宸 孙传瑞 陈胤孜 陈晓葳 王锡玲△

（1复旦大学公共卫生学院卫生统计学教研室 上海 200032； 2复旦大学公共卫生安全教育部重点实验室 上海 200032）

COVID-19变异株Omicron毒株于2021年11月在南非等地被发现［1］，因传染性和免疫逃避能力显著超过其他变异株［2］，其迅速在全球成为主要流行毒株［3］。实时再生数（time-varying reproduction number，以下简称Rt）即t时刻一个病例可传染所致续发病例的平均数量，在流行病学研究中常被用来快速评估某传染病的实时平均传播水平。当Rt＞1时认为疫情会继续发展，Rt＜1时则认为疫情逐渐受到控制。新冠肺炎疫情Rt的实时估计已在国际上广泛应用，英国、日本和香港等国家和地区均向公众实时展示Rt的变化［4-6］以研判疫情走势和调整干预措施。2022年3月上海市出现新冠本土病例并发生社区传播，流行毒株为Omicron BA.2和BA.2.2型［7］，政府采取了一系列应对措施。目前少有研究评价上海市此轮疫情Rt值的实时变化并基于此预测疫情后续发展［8-9］。

新冠疫情可视化平台的开发具有重要意义，其可将信息图表化以提高复杂数据的可读性，满足大众及科研人员的信息需求［10-11］。目前有多个网站进行新冠数据可视化展示，如丁香园、腾讯新闻等［12-13］。然而，这些平台内容以数据描述为主，相似度高，缺乏对上海市各区的详细描述，没有进行新冠病毒实时传播能力Rt的评估和预测［14］。本研究将基于Rt构建贝叶斯分层模型，评估和预测上海市新冠疫情发生发展并开发新冠疫情可视化平台。

资料和方法

数据来源和病例定义2022年3月1日—5月6日上海市COVID-19疫情每日确诊病例、无症状感染者和新冠疫苗接种数据［15］来源于上海市卫生健康委员会（以下简称“卫健委”）官方网站（https：//wsjkw.sh.gov.cn/）。确诊病例的定义参照国家卫健委修订的《新型冠状病毒肺炎诊疗方案（试行第九版）》［16］。无症状感染者未在第九版诊疗方案中明确定义，一般指除确诊病例以外的病原学检测呈阳性的感染者。

政府所采取的疫情防控措施及实施日期来源于上海市官方微信平台“上海发布”和疫情防控工作新闻发布会。根据所影响的人群和地区，本研究将各类措施划定为低、中、高三级强度。低强度包括3月12日起幼儿园及中小学线下停课和3月15日起的切块式、网格化筛查，中强度包括3月28日起以黄浦江为界分批、分区管控，高强度则包括4月3日起全市范围进行多轮抗原和核酸筛查。常住人口年龄构成来自上海市统计局第七次人口普查数据［17］。

模型构建与主要参数设定本研究使用基于实时再生数构建的贝叶斯分层模型进行预测，模型中假定每日报告的感染者数It服从以It为均数和It+为方差的负二项式分布，考虑感染到报告延迟时间分布π和感染者诊断比例Piar来推导每日实际感染者数it，如公式（1）所示。公式（2）中每日实际感染者数it由该时刻的实时再生数Rt和该时刻前感染者的离散卷积和决定，离散卷积和由此前时刻的感染者数通过代际间隔分布g（即一代病例和二代病例感染的时间间隔分布）加权计算得到。实时再生数Rt则受到主要流行毒株的基本再生数R0、人群免疫水平（包括自然感染和接种新冠疫苗所获得）和非药物干预措施强度的影响，S为经人群免疫水平校正后的易感者比例，pin和pva分别指由自然感染和接种疫苗获得的免疫效果，见公式（3）。

上海市疫情发生社区传播后病例数激增，医疗资源不足，存在病例错分可能，研究中以阳性感染者而非确诊病例作为拟合指标，每日报告的阳性感染者数由每日报告的确诊病例数+无症状感染者数-由无症状感染者转归的确诊病例数计算得到。为减少数据的随机波动性，本研究将每日报告数据进行移动平均处理，取所计算日期前3天至后3天数据（n=7）的算术平均值。研究中所用的主要参数来源于国内外相关研究和合理假定，包括流行病学参数、疫苗接种和自然感染后预防感染的免疫效果及人群年龄别疫苗覆盖率，具体设定见表1。考虑到当病例数激增时报告的延迟性会进一步增加，因此以每日报告感染者数是否超过5 000例为界，使其服从均数不同的Gamma分布。

模型拟合与预测性能评估本研究以2022年3月1日—4月29日为模型拟合段，以4月30日—5月6日为预测性能评估段。模型拟合段以最大似然法［25］估计不同强度的防控措施对Rt的影响，并选择最佳参数预测疫情的后续发展，拟合偏差比例的计算如公式（4）所示，其中kt指实际观察值，Pt，median是t时刻预测概率分布的中位数。鉴于感染到报告的延迟性，未纳入4月23—29日的数据用于预测，假定此时间段内的Rt值和4月22日的数值相同。

通过对比预测结果与官方报告数据来评估模型预测性能。预测所得的每日感染者数是由非负整数集合的概率分布，其评估首先需要验证预测概率分布和实际观察值在统计学上具有一致性，再衡量集中度［26-28］。本研究中选择随机化概率分布变换（probability integral transform，PIT）和分级概率评分（ranked probability score，RPS）来分别对预测结果统计学的一致性和集中度进行评估。

PIT常被用来探究实际值和连续性的预测概率分布的统计学一致性，而对于每日新增感染数这样的离散性资料，随机化PIT是一种新的替代方法［27-28］，计算方法如公式（5）所示，其中kt指实际观察值，Pt(x)是t时刻的预测概率分布，根据定义Pt(-1)=0，v是标准均匀分布且与k无关。若Pt确为实际观察值kt的预测概率分布，那么ut将是标准均匀分布，再用Anderson-Darling检验来验证该假设在统计学上是否具有意义。当检验结果P≥0.05时认为没有充分证据表明实际值和预测概率分布在统计学上不具有一致性；当0.01＜P＜0.05时认为有证据表明两者可能不具有统计学的一致性，而当P≤0.01时认为有足够的证据表明两者不具有一致性，即预测概率分布相对实际观察值是失准的。RPS计算如公式（6）所示，其中k，Pt(x)含义同公式（5），当计算结果越接近于0，则认为实际观察值越有可能来自于预测概率分布且预测概率分布的集中度越好。研究分析均用R 4.1.3软件完成，其中模型构建使用RStan 2.29［29］。

可视化平台开发可视化平台的开发主要包括数据提取与整理、数据可视化、预测模型构建和平台搭建等4个步骤，基本架构如图1所示。使用Python 3.10的Request、Re和Datetime第三方库提取数据。使用R 4.1.3软件完成数据整理、分析和可视化。使用网页开发语言和R shiny 1.7.0框架［30］交互实现平台搭建。R shiny框架已将网页开发中的Ajax技术、Bootstrap框架和jQuery功能框架进行封装以实现数据通信，而为满足页面布局的美化需求，后续再分别使用HTML进行网页前端设计，CSS进行网页布局优化以及Java Script提高网页互动效果。平台使用Echarts实现条形图、线图、地图、热图等数据可视化，并支持鼠标悬浮或点击等交互。该平台的系统架构搭建并运行于Ubuntu服务器。

图1 新冠肺炎疫情预测预警平台基本架构Fig 1 Architecture of the COVID-19 prediction and early warning platform

结 果

疫情概况2022年3月1日，上海市卫健委报告首2例本土新冠感染者。截至4月29日，累计报告确诊病例52 576例（其中18 686例由无症状感染者转归而来），无症状感染者530 642例（图2A）。上海市本土疫情发展呈现先缓慢增加后迅速上升的态势，在4月13日达到疫情高峰，当日新增27 605例感染者。此后疫情放缓，每日报告的感染者数波动下降，4月27日起每日新增感染者数已小于10 000例。

图2 上海市2022年COVID-19疫情发展趋势、实时再生数及预测结果曲线Fig 2 COVID-19 epidemic curve， Rt values and forecast result of Shanghai in 2022

根据疫情发展，上海市的非药物干预措施逐渐从“精准防控”策略过渡到3月15日开始的切块式、网格化防控筛查，又调整为3月28日起以黄浦江为界分批、分区实行封控，再到4月初对全市范围进行核酸筛查和防控管理。随着防控措施升级，Rt曲线总体呈现波动下降，在4月5日其均值首次降至1以下（0.94，95%CI：0.87～1.03），结合延迟时间分布，与4月13日疫情达到高峰相符（图2B）。在4月13—14日Rt的95%CI上限再次超过1，达到1.01，与4月22日感染者新增数出现新的小高峰的趋势贴合。考虑感染到报告的延迟时间，4月23—29日（即图中灰色区域）的Rt估计值仍存在较大不确定性。

模型拟合结果拟合时间段中，经移动平均处理后的每日新增感染者数与预测概率分布中位数拟合良好，拟合偏差比例为14.41%±15.20%。在4月20—26日模型的预测中位数明显高于对应的移动平均值（如图2C蓝线和黄线所示），其原因可能是4月22日新增感染者数又出现一个小高峰而模型拟合未能体现。各级防控措施的调整对Rt值的相对效应如表2所示，上海市所有防控措施可将Rt值降低至1以下（0.30，95%CI：0.28～0.32），而在放松部分区域的防控后，Rt值增加至0.43（95%CI：0.41～0.46），仍未超过关键阈值水平，但使得本轮疫情中每日新增感染者数的下降速度在一定程度上放缓。

表2 上海市疫情防控措施对Rt值的效应估计Tab 2 The estimated relative effectiveness of nonpharmaceutical interventions on Rt in Shanghai

疫情走势预测模型预测结果如图2C所示，本轮疫情预计累积感染者数将达到625 040例（95%CI：600 465～670 211），预计将于5月17日（5月12—22日）单日新增感染者数＜1 000例，将于6月5日（5月31日—6月11日）单日新增感染者数＜100例。对4月30日—5月6日的预测结果进行评估，随机化PIT的结果是P=0.12，尚无足够证据表明结果概率分布和实际观察值无统计学一致性，分级概率评分指标为382.53±105.29。

可视化平台可视化平台实现了上海市疫情变化趋势描述和疫情预测两大功能。疫情变化趋势描述包括上海市及上海市各区新增、累积阳性感染者和死亡人数变化趋势、各区累积感染率和感染者来源（隔离管控或风险人群排查）。疫情预测包括实时再生数、新增阳性感染者预测曲线和预测感染人数关键节点。为满足用户端的不同需求，平台内置图表适应电脑和手机端多种主流浏览器。该平台于2022年4月19日起开放访问，访问地址：http：∥124. 222.178.35：3838/COVID19-tracker/。

讨 论

本研究使用基于实时再生数构建的模型，评估和预测了2022年3月以来上海市本土疫情的发展态势，并完成了新冠肺炎预测预警平台的搭建。研究中考虑上海市现有免疫水平并设置随时间变化的模型参数，包括感染者从感染到报告的延迟时间分布和各级非药物干预措施的相对效应，预测结果更贴合疫情的实际发展。而在预测段中虽假定Rt值不变，模型评估指标和后续实际报告数据依然证实了预测结果的准确性。2022年底，我国新冠疫情防控策略已有所调整，但调整相关参数后，此模型也可用于预测其他新发呼吸道传染病以评估公共卫生风险和干预措施有效性。

模型拟合结果显示，上海市所采取的非药物干预措施将Rt值降低并保持在1以下，且本轮疫情预测将于6月中下旬实现“动态清零”目标。另有研究通过上海市每日报告的新增感染者数据估计Rt值，同样证实了上海市防控措施的效果［31］。刘可伋等［32］使用时滞动力学模型预测疫情发展并估计实时再生数，但与本研究不同的是未充分考虑感染到报告的延迟时间，预测结果具有滞后性。疫情后期，连续14天无新增感染者区域的管控放松，而存在传播风险的地区仍在严格管控中，防疫政策的变化对结果影响不大。但当跨区及跨省市人员流动增加时，新冠传播风险增加，感染规模可能进一步扩大。

研究中假定上海市自2022年3月1日报告首2例感染者后已暂停全市的新冠疫苗接种，而部分区域在4月下旬已开始恢复对老年人等重点人群的疫苗接种［33］，这部分人群的免疫水平未在模型中体现，考虑到日均接种千余剂次和疫苗产生有效抗体水平需要2～3周的时间，可认为其对预测结果影响有限。但若老年人群的接种率进一步增高时，其所提供的保护效果将使得疫情控制速度快于预测结果。鉴于上海市内各区域人员流动频繁，且非药物干预措施的调整多是基于全市范围，研究使用全市数据进行预测，区级疫情预测需要更多数据来支撑分析。

预测预警平台使用Python进行公开数据实时爬取，基于实时再生数进行疫情走势预测，使用Echarts和R shiny实现数据可视化和交互，具有实时性、科学性、页面友好等优点。通过数据可视化，高效展示了上海市及各区新冠疫情变化的时空差异，降低公众捕获并理解公开疫情数据的时间成本。目前已有预测预警平台无法提供Rt值以估算疫情实时传播速度［34-35］。本研究基于实时再生数的预测模型可视化，填补了已有平台对疫情数据深层次挖掘的空白，提供更加科学的防控措施效果评价和疫情预警预判。该可视化平台可快速扩展到其他城市的疫情数据解读及预测预警，乃至其他可能发生的具有类似流行病学特征的感染性疾病，具有较好的应用前景。

作者贡献声明何沅宸 研究设计，数据分析，论文撰写和修改。孙传瑞 数据分析，论文撰写和修改。陈胤孜 数据整理和分析。陈晓葳 数据采集。王锡玲 研究设计和指导，论文修改，经费支持。

利益冲突声明所有作者均声明不存在利益冲突。