基于SIQR模型的新冠肺炎期间深圳市防控措施分析

徐韧哲，汪家伟，叶声豪，王 雄

1) 深圳大学高等研究院，广东深圳518060；2)安徽农业大学信息与计算机学院，安徽合肥230036

2019 年12月新型冠状病毒肺炎(coronavirus disease 2019, COVID-19，以下简称新冠肺炎)疫情在中国武汉市出现，并在较短时间内蔓延至其他地区．随着中国各项防控措施的紧急出台，以深圳市为代表的输入型地区疫情很快趋于稳定．世界卫生组织在《中国-世界卫生组织新型冠状病毒肺炎(COVID-19)联合考察报告》中积极评价中国防控措施的影响，并将深圳市作为唯一一个典型城市个案进行防控措施分析．然而，国际上许多其他地区却采取了相比中国更为温和的防疫措施．因此，以量化方法展示不同地区及不同程度的疫情防控措施所带来的效果十分必要，同时也可为平衡经济发展与疫情防控提供重要的理论支持．

经典的传染病模型SEIR(susceptible exposed infectious recovered)早在20世纪便已提出[1]，成为众多传染病扩散研究的基础模型．新冠肺炎疫情爆发以来，已有研究对疫情发展进行分析与预测．JONATHAN等[2]在2020年1月对中国疫情进行初步估计，但研究没有考虑政府防控措施的有效性，致使估计结果远高于病情的实际发展情况．WU等[3]利用考虑了输入输出型病例的SEIR模型对中国实际疫情传染情况进行估计与预测，并根据境外确诊人数反推出2020-01-28中国实际感染人数为6 000人左右．周涛等[4]基于SEIR模型及报告确诊病例数，估算出中国新冠病毒的基本再生数为2.9～3.3，认为及时有效的防控措施能够较快抑制传染病的进一步蔓延．TANG等[5]对SEIR模型进行改进，在模型中引入隔离观察者、潜伏期隔离者及发病期隔离者，并预测不同参数下的可能病例数．通过离散模型也可以动态描述新冠肺炎疫情的发展，并估计所需的流行病学参数，其离散的传播模型和稍作修改的SIR(susceptible infectious recovered)模型在动力学方面具有良好一致性[6]．曹盛力等[7]基于修正SEIR模型对武汉的疫情发展进行回溯分析，为制定未来的疫情干预决策提供一定理论支持．

本研究以疫情输入地为研究对象，考虑疫情发源地与输入地之间的传播关系网络，建立以防疫措施为参数的传播动力学SIQR(susceptible infectious quarantined recovered)模型．以中国人口迁徙数据和深圳市疫情数据为依据，使用基于蒙特卡洛方法的网格搜索算法对SIQR模型的各个参数进行最优化拟合．依据模型对深圳市疫情防控措施进行打分，并模拟不同得分下疫情的发展情况．最后，对比深圳与其他疫情国家，分析其他国家疫情发展的特点．

1 资料与数据

模型考虑的时间范围为2020-01-10至2020-03-10，研究对2020-01-10至2020-03-10期间全国迁徙数据和2020-01-24至2020-03-10期间全国疫情数据进行收集和分析，并假设疫情传播过程中病毒未发生显性突变．

1.1 全国疫情数据

这里只取中国国家卫生健康委员会(国家卫健委)2020-01-24至2020-03-10公布的每日疫情数据[8]，经过分析拟合，最终估算出自1月10日起的全国各省份的真实感染数据，具体步骤如下．

步骤1根据生物界常见的S型生长曲线，使用logistic模型对全国各省确诊人数进行拟合及估计，拟合函数为

(1)

其中，Ci(t)为第t天i地区(省市)的累计确诊病例数；a,b及c为参数，分别影响函数的最值、坡度及位置．

步骤2由于疫情传播过程中存在具有传染性的潜伏病例，还需根据确诊病例数估计各省市的潜伏病例数，从而得到各个时期、各个省市的真实感染人数．这里取平均潜伏期为7.5 d[9]，则潜伏者每日转化为感染者的概率约为1/7.5≈0.13．由此可根据各个省市的每日新增确诊数反推每日存在的潜伏者病例数为

Ei(t)=[Ci(t+1)-Ci(t)]/0.13

(2)

其中，Ei(t)为第t天i地区(省市)的潜伏者病例数．

步骤3将以上两步计算得出的各个省市每日确诊病例数和每日潜伏病例数相加，即可计算出第t天i地区(省市)的感染总人数为

Ii(t)=Ei(t)+Ci(t)

(3)

1.2 全国迁徙数据

基于百度地图迁徙大数据[10]，将全国各个地区至深圳的人口迁徙情况分两个阶段：

1) 武汉“封城”前：使用2020-01-10至2020-01-23的各省市迁入迁出深圳人口流量强度日均值，单位流量强度为10 万人，按照平均每日迁入深圳流量强度排名，前6名的省份见表1.

表1 武汉“封城”前各省平均每日迁徙 深圳流量强度(前6)Table 1 The average daily migration flow intensity of each province to Shenzhen before the lockdown of Wuhan (top 6) 10万人

2)武汉“封城”后：使用2020-01-24至2020-03-10各省市迁入迁出深圳人口流量强度日均值，单位流量强度为10 万人，平均每日迁入深圳流量强度前6名见表2.

表2 武汉“封城”后各省平均每日迁徙 深圳流量强度(前6)Table 2 The average daily migration flow intensity of each province to Shenzhen after the lockdown of Wuhan (top 6) 10万人

1.3 深圳疫情数据

根据深圳市政府数据开放平台提供的详细病例数据[11]，对深圳市疫情数据进行分析整合，以便后续模型的拟合与分析．假设深圳常住人口总数为2 500万人，各类病例统计情况如图1．以下模型根据病例情况，使用输入型病例的来深时间和感染型病例的染病时间，作为深圳新增感染人数的时间，可有效避免确诊日期的滞后性．深圳市的新增病例中，82.9%为输入型病例，感染病例中55.2%为家庭传播，每组家庭传播中病源平均传染人数为1.7人[12]．以2020-01-23为分界点，此时间前后的病例自染病至入院的平均时间分别为6.5 d和4.2 d；由百度迁徙大数据可得，分界点时间前后深圳市平均城市内部人口流动规模强度分别为5.0和0.8．

图1 深圳每日病例情况统计Fig.1 Statistics of daily cases in Shenzhen

2 模型与方法

2.1 考虑输入型病例的 SIQR 模型

鉴于新冠肺炎患者在潜伏期和发病期同样具有传染性，而隔离后的患者传染性几乎为0，因此，本研究将全部感染病例划分为感染者I和隔离者Q．感染者为染病后具有传染性的患者(无论是否已经发病)，隔离者为发病就诊后已经被隔离的患者和已经被提前隔离的潜伏期患者．除感染病例外，未被感染但存在被传染可能性的为易感者S．感染后康复且不会再被感染的人群为康复者R．模型中新加入因人口流动所带来的输入型患者，这部分患者被包括在每日新增感染者人数中．模型内部人群转化关系如图2．

图2 模型内部人群转化关系示意图Fig.2 Population transformation diagram in the model

设第t日深圳市内部人口总数为N(t)，

N(t)=S(t)+I(t)+Q(t)+R(t)+In(t)-Out(t)

其中，S(t)、I(t)、Q(t)及R(t)分别为第t日对应的易感者、感染者、隔离者及康复者人数； In(t)和Out(t)为第t日迁入和迁出深圳的总人数．

SIQR模型方程组构建如下

其中，M为每日每个感染者的传染人数；α为感染者每日隔离率；γ为隔离者的每日康复率；δ为输入人口管控力度； InS和InI分别为每日输入易感者和感染者； OutS和OutI分别为每日输出易感者和感染者；Ii(t-1)为第t-1天i地区(省市)的感染者总数；Ni为i地区的人口总数；Fi,in和Fi,out为i地区每日迁入和迁出深圳人数．

对以上模型中的M， 有

M=(k2cβ+l)e-f

(4)

其中，k为城市内部人口流动规模；c为城市内部单位流动人口规模下，每人在单位出行频率下的有效密切接触人数；β为密切接触传染率；f为个人防疫力度；l为因家庭内部接触而造成的平均每日感染家庭成员数，由1.3节中数据可知l≈0.3． 考虑到每人每日的出行频率与城市内部人口流动规模有关，因此，每人每日密切接触人数为k2c． 对于密切接触者，其被感染的概率为βe-f．

2.2 深圳疫情仿真

根据深圳市人口流动数据及相关防控政策出台时间，本研究将深圳市疫情发展分为2个阶段：① 2020-01-10至2020-01-23的武汉“封城”前阶段，深圳市的人口输入指数、人口输出指数、城市内部人口流动规模及防疫力度等都处于无防疫措施时的正常水平；② 2020-01-24开始的武汉“封城”后阶段，深圳市的输入和输出人口指数明显下降，来自武汉地区的人员大幅减少，对入境人群的防疫筛查加大力度，城市内部人口流动较平时有大幅降低，城市内部人群有较高的防疫水平．

SIQR模型对参数的敏感性较高，除β为差别较小的参数外，其余参数在不同时间点和不同地区取值不同．因此，本研究使用基于蒙特卡洛方法的网格搜索算法，对深圳市不同防疫阶段进行最适参数的数值求解，从而得到防疫措施对模型参数的影响．武汉“封城”前后深圳市各阶段参数对比见表3．其中，参数δ、f及c由拟合得出；k由百度迁徙大数据获取；α取深圳病例自发病至入院平均时间的倒数．

表3 深圳各阶段参数对比Table 3 Comparison of parameters in different stages of Shenzhen

根据深圳市自2020-01-10起的累计实际感染者人数，进行模型拟合，结果如图3．由于深圳市人口基数大，感染人数相对较少，结果中不再展示S、 Q及R人数，仅展示衡量城市管控措施最为关键的指标——每日感染人数．在防疫前深圳市的每日新增感染者人数持近似指数型上升，且新增感染者中80%左右都是输入型病例．在武汉“封城”后，由于深圳的输入病例，尤其是来自湖北的输入病例大幅减少，因此，每日新增感染者人数开始下降．且在深圳内部由于人口流量管控及戴口罩等个人防疫措施的作用，传染率始终保持较低水平．同时，政府的高效行动使感染者的隔离率进一步提升，从而有效避免疫情在城市内部爆发．由图3可知，拟合优度R2=98.2%， 表明模型很好地模拟了深圳疫情的发展规律．

图3 深圳市感染数拟合结果Fig.3 Fitting result of the number of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) infected people in Shenzhen

3 分析与评估

基于深圳市疫情传播模型，使用控制变量法对各项防控措施进行量化评估，以解释其在疫情防控中所起到的真实作用．评估参数包括武汉“封城”时间T、δ、f及k． 为了便于比较，根据各参数有效范围对各项防疫措施进行打分，无防疫措施时为0分，最佳理想防疫状态下为100分．

3.1 武汉“封城”时间

武汉“封城”时间与输入人口管控都直接影响深圳的外来输入型病例．武汉“封城”时间是2020-01-23，即为自2020-01-10起的第13 d，模拟假设T分别为18、23及28 d，并在深圳其他防疫措施不变的情况下估计疫情发展，结果见图4．

图4 自2020-01-10起武汉“封城”时间对深圳疫情影响Fig.4 Effect of Wuhan lockdown time on Shenzhen epidemic since January 10th, 2020

由图4可见，若武汉“封城”时间为自2020-01-10起的第18 d，深圳的染病人数将从419人增加至1 895 人；而若为28 d后，深圳的染病人数会大幅提高至28 698人，是当前防疫措施下染病人数的69.5倍．一旦处置不当便可能面临疫情失控的风险．中国政府所采取的“封城”时间点，从经济和人民生产生活的角度分析较为理想，及时避免了疫情在其他地区的爆发．

3.2 输入人口管控

武汉进行“封城”之后，来自疫情发源地的威胁大幅减少，但由于疫情已出现扩散态势，作为疫情输入地的其他地区仍不能懈怠，需对地区的输入人口进行大力度的监察和管控．假设深圳自2020-01-23起全面对外管控，规定输入人口管控力度得分为100δ/2.5， 比较得分分别为0、20、40及80分时的感染人数，模型估计结果如图5．

图5 输入人口管控力度对深圳疫情影响Fig.5 Effect of input population control strength on Shenzhen epidemic

由图5可见，虽然自2020-01-23起武汉已经进行“封城”，但若深圳的管控就此松懈，当δ=0时，深圳的染病人数依然会大幅提升至2 147人；当输入人口管控力度减少至20分时，染病人数为950人；当输入人口管控力度减少至40分时，染病人数为571人．因此，武汉“封城”后深圳的管控措施在很大程度上减弱了外来病例输入，输入人口的管控力度效果已近乎达到理想情况，深圳市政府很好地完成了作为疫情输入地区的输入性病例防控．

3.3 个人防疫力度与城市内部人口流动管控

当已有染病患者进入城市，患者会在城市内部传播疫情．因此，在深圳市具有较好对外防控基础上，也需要分析比较不同内部管控对疫情的影响．城市内部个人防疫力度f代表市民因戴口罩及保持社交距离等个人措施，使因密切接触而导致传染概率的降低程度；城市内部人口流动管控k代表市民因减少外出及保持居家隔离等措施所导致有效密切接触次数的减少程度．规定个人防护力度得分计算方式为100f/2.5． 根据百度迁徙大数据，深圳市正常情况下的k=5， 防疫状态下k减小至0.8．规定城市内部人口流动管控得分为100(5-k)/5． 由于这两项防控措施相辅相成，同时分析当f=0、 20、 60及100和k=0、 40及84时的疫情传播情况，结果如图6．

由图6可见，若深圳人口流动管控和个人防疫力度同时为0时，将会出现疫情大爆发，感染人数将暴增至1 346万人．当放松人口流动管控至40分时，若个人防疫力度为60分，感染人数会增加到955人；若个人防疫力度为20分，感染人数会增加到468万人．因此，当个人防护资源不充足、难以保证个人防疫力度时，城市人口流动管控可以大幅度有效降低疫情风险．当个人防护力度足够大时，足以有效控制疫情，此时可以适当放松人口流动管控，以满足人民生产生活需要．

图6 不同个人防疫力度下，不同城市内部人口流动规模对疫情的影响Fig.6 The effect of population flow scale in city on epidemic situation under different individual epidemic prevention efforts

3.4 深圳防疫措施与其他国家对比

选择疫情发展特点明显的国家与深圳的疫情发展情况比较，以当地的疫情初始报告时间为起点(意大利：2020-01-31；美国：2020-02-01；新加坡：2020-01-27；深圳：2020-01-10)，将感染人数改为以感染人数占该地区总人数百分比的形式表示，结果如图7．同时，根据各地区对应的模型参数进行各项防控措施打分，结果见表4．

意大利疫情发展的特点是在疫情初期几乎完全没有对外防控措施，但随着疫情爆发，国家内部及时“封城”并大力防疫管控，此后疫情增长有所放缓．新加坡疫情恰好相反，在疫情初期一直有着较为严格的输入人口管控，但城市内部的防疫力度不强，虽然初期保持相对较低的感染率，但随着时间推进，逐渐出现疫情爆发．美国的输入人口管控和城市内部防疫力度在3地中相对比较均衡，但管控和防疫力度仍在较低水平，前期的输入病例导致感染人数快速上升，在没有完善的内部防疫措施下，美国很快成为新的疫情震中．对比之下，深圳在综合打分上高于以上3个国家，具有全方位、力度大的疫情防控措施，这些措施有效避免了疫情在深圳的爆发．

图7 深圳疫情发展与其他国家对比Fig.7 Comparison of epidemic development in Shenzhen with that in other countries

表4 深圳与其他国家防疫措施打分
Table 4 Comparison of the scores of epidemic prevention measures between Shenzhen and other countries

输入人口管控力度得分个人防疫力度得分城市内部人口流动管控得分深 圳809284意大利105080美 国404525新加坡853530

结 语

本研究基于新型冠状病毒肺炎输入型病例地区特点，建立考虑输入型病例的 SIQR 传播动力学模型对深圳市的疫情发展进行模拟，使用基于蒙特卡洛方法的网格搜索算法对疫情不同阶段的参数进行拟合，并对各项防控措施进行量化评估，以解释其在疫情防控中所起到的实际作用．分析结果如下：

1)在保持其他防疫措施不变的情况下，若武汉“封城”时间再推迟15 d，深圳的染病人数会大幅提高至28 698人，是当前防疫措施下染病人数的69.5倍．考虑到当时处于春运期间，综合经济发展和人民生产生活，武汉“封城”时间点较为理想．

2)在武汉“封城”后，深圳市政府对来自其他地区的输入人口管控使深圳可能的染病人数下降80%，很好地完成了作为疫情输入地区的输入性病例防控工作．

3)在城市内部防控方面，若深圳人口流动管控和个人防疫力度同时为0分，将会出现疫情大爆发，感染人数将暴增至1 346万人．当放松人口流动管控至40分时，若个人防疫力度为60分，感染人数会增加至955人；但若个人防疫力度为20分，感染人数会增加到468万人．因此，在疫情期间适当放松人口流动管控时，应时刻保持高水平的个人防疫力度，否则可能会出现新一轮的疫情爆发．

与意大利、美国及新加坡防疫情况的对比分析发现，在全球出现多个疫情震中的情况下，只有外部防控和内部防控全面同时进行，才能有效抑制疫情的发展．