疫情拐点何时出现：基于数学模型的新冠病毒传播预测

　　前言

　　截止到2020年2月2日24:00，我国确诊的新型冠状病毒患者已超过1.4万人，短短一个来月时间确诊病例数已经远超2003年“非典”疫情的全部确诊数，目前每日新增的确诊数仍在攀升，疫情的传播速度超过了大多数人的预期。科学客观地评估新冠肺炎的传染性强弱以及预测患病人数规模和峰值时间，对决策者实施必要的防控措施、评估对经济的影响以及投资者如何应对都具有重要的现实意义。

　　作为数学建模的爱好者，而非流行病学的专业人员，作者通过搜集和学习海内外文献，对经典的流行病传播模型（SEIR）有了较准确的理解，通过适当的简化，较直观地向读者介绍模型的原理和各参数的意义。对于经典且成熟的模型来说，输出结果的可靠性完全由输入参数的准确性决定，因此本文从新冠病毒的具体情况入手，着重从多个角度来评估、检验和校正模型的输入值，提高模型输出结果的可靠性。

　　主要内容&结论

　　主要内容：

　　1）基于SEIR流行病传播模型测算基本传染数R0（也叫基本再生数）来衡量传染性的强弱普遍存在输入参数不易确定、输出结果对参数非常敏感的问题。国内外不同文献对同一疫情（比如SARS）的测算结果可能相差甚远。简单取一组参数得到的输出结果并不可靠。

　　2）本文对SEIR模型做适当简化，减少参数，再针对新冠病毒已明确的具体信息设定参数区间范围，再从多个维度选取校验条件，利用计算机数值模拟的方法求取参数值。通俗的来讲就是让计算机程序在参数区间范围内做大量尝试，试探出较合理的参数值，使得各个校验条件都能较好的相恰。

　　3）通过1月14日到22日武汉迁出人口排名前8的省市随后一周的确诊病例数与武汉迁出人数的关系，来反推武汉潜伏期感染者的比例和数目，获得求解模型的重要约束条件。通过数值求解得到参数值，代入模型对新冠疫情进行评估和预判。

　　主要结论：

　　1）新冠病毒基本传染数明显大于SARS，但单次接触的传染概率可能略低，民众无需恐慌。根据数值模拟确定的参数值，新冠病毒的基本传染数（R0）远高于SARS，即不采取防疫措施的情况下平均每个病人传染的人数更多。但进一步分析表明，新冠病毒R0更高的主要原因可能是潜伏期末期具有传染性和轻症患者比较多，潜伏期末期和轻症患者由于症状不明显，不会及时就医住院或限制自己活动，使得暴露在易感人群中的时间（传染期）更长。而单次接触的传染概率则可能略低于SARS。可见，潜伏期长和轻症患者多增加了疫情防控的难度，对病例的接触者采取更早、更全、更严格的隔离可能是最有效的防控措施。

　　2）当前新增确诊病例数大幅攀升（日新增病例超2500）可能与检测确诊效率提高有关，预计将很快回落（模型中性假设下预测值是1700左右）。1月25日之前，公布新增确诊数大幅低于模型预测值，26日之后则大幅高于模型预测值，可能是因为确诊检测试剂盒供应充分了，将之前已经发病而没有确诊的病例确诊了。

　　3）假定隔离和防控措施能达到或接近SARS后期的效果，则预计未来总患病人数将达到88500人（中性假设下）或58000（乐观假设下）左右，新增确诊人数将在2月17日（中性假设下）或2月11日（乐观假设下）左右出现拐点。

　　4）新冠病毒的重症及致死率低于SARS、禽流感、H1N1等，由于潜伏期长且轻症患者较多，控制疫情传播需要更加严格的隔离防控措施，对短期经济会产生明显影响。

　　5）从2009年美国H1N1大爆发以及2003年“非典”期间美股、港股和A股表现来看，在疫情快速传播期，股市会明显调整，但幅度不会很大，且往往在疫情形势好转之前率先见底。

　　风险提示：1）计算结果对参数较为敏感；2）政策实施对疫情有较大影响。

　　一、如何定量评估流行病的传染性强弱？——SEIR模型

　　SEIR是经典的流行病传播模型，这里我们只介绍其原理和框架，然后引用一些参考文献的结论计算基本传染数，而不对微分方程的推导过程展开讨论。SEIR模型将整个社群人口分为四类：

　　易感人群S类：对该病毒没有抗体，接触传染期人群可能被传染；

　　潜伏期人群E类：被感染后处于潜伏期的人群，还没有传染性，平均时长记为TE；

　　传染期人群I类：代表潜伏期之后已具有传染性的人群，传染期的平均时长为TI；

　　隔离态人群R类：治愈并获得免疫或被有效隔离或死亡，既不能传染他人也不能被传染。

　　特别说明：

　　1）SEIR模型中潜伏期（E类）是特指已感染但还没有传染性的时期，设这一时期的平均时长为TE。但通常对潜伏期的理解是已感染但还没有症状的这段时间（卫健委专家表示新冠病毒的平均潜伏期在10天左右），本次新冠病毒已发生过在没有症状的时候传染给下一个人的案例。国家医疗专家组成员、北京地坛医院感染性疾病诊疗与研究中心首席专家李兴旺表示，“从一般的呼吸道传染病规律来讲，潜伏期末期，病人将要发病的时候可能具有一定的传染率，早期一般不应该有传染性”。因此，本文假定平均为10天的无症状潜伏期中，没有传染性的潜伏期平均为6.5天（参考了SARS的相关数据），最后的3.5天具有传染性，计入传染期TI。

　　2）传染期的平均时长TI是指病人具备传染性后到不会接触不特定的易感人群为止，比如住院或居家卧床等，而不是从病人发病到痊愈的时间。查询了大量文献，发现在对各类疫情建立SEIR模型时，设置的传染期时长通常都在3到5天，甚至更短，这显然不是病人发病到痊愈的时间，比如SARS很多患者都需要住院治疗很长时间。这样理解下，症状越重的病毒其平均传染期TI就越短，因为患者会很快住院或散失活动能力，而症状越轻的病毒则可能传染期较长。

　　假设一个传染期的人（I类）与易感人群（S类）接触后易感者被感染进入潜伏期（E类）的概率为β，一个处于潜伏期（E类）的个体单位时间内将以的概率转变为传染期（I类）个体，一个传染期（I类）个体单位时间内将以的概率转变为隔离态（R类）个体。

　　图表1：SEIR模型

　　则病毒传播的过程可由这四类人群变化的微分方程来刻画[2]：

　　其中S(t)、E(t)、I(t)和R(t)分别表示t时刻社群中处于易感染态、潜伏期、传染期和隔离态个体数目，N表示社群中个体总数目，且N=S(t)+E(t)+I(t)+R(t)。据参考文献[3],基本传染数可表示为：

　　其中λ=ln(Y(t))/t是指数增长率，对某一t时刻发生的累计病例数Y(t)取对数，再除以t。潜伏期和传染期时长可分别表示为Tg=1/γ1和TI=1/γ2，序列间隔（serial interval）则为Tg=TE+TI。记ρ=TE/Tg为潜伏期占序列间隔的比例，则基本传染数可进一步表示为：

　　计算基本传染数R0需要确定的参数有：

　　λ：其等于ln(Y(t))/t，即需要知道第一例病例发生的时间（确定起点以便计算t），以及之后某一时刻的累计病例数Y(t)。

　　Tg和ρ：只要确定了潜伏期TE和传染期TI，Tg和ρ就都确定了。

　　实际测算中，这三个参数都不好确定。由于基本传染数R0是指没有人为的防疫、隔离等干预的情况下平均每个病人能传染多少个人。在实际情形中，在流行病传播初期没有强力的人为干预，但也就缺乏必要的研究和统计，对病毒的潜伏期、传染期以及病例统计数据都很不精确，但R0的数值却对这些参数却非常敏感。

　　以非典为例，查询估算SARS病毒的基本传染数R0的文献，发现估算结果从1.1到4.2之间的都有，在2到3之间的比较多。我们以文献[1]中关于SARS的一组参数展示一下计算R0的过程：

　　文中以香港2003年3月28日累计病例数425为依据，2003年2月15日为香港首例病例，作为起点，即Y(41)=425,取Tg=10，ρ=0.6。

　　λ=ln(Y(41))/41=0.1476,代入到R0的计算公式可以得到R0=3.

　　那么对于本次新冠病毒要怎么确定参数呢？在政府采取强力干预前，只有武汉地区有确诊数据，但有多个证据都表明武汉前期的确诊数据可能失真了。比如有文章统计，武汉新冠病毒患者死亡率是湖北省外其他地方的30倍，这明显不合理。可见根据官方公布的确诊数来定一个Y(t)值是不可取的。

　　二、基于已公开信息构建新冠病毒传播预测模型

　　由于前期的确诊数据可能失真较严重，需要从其他已公开信息来估算早期病例数，从而计算R0并对未来疫情发展做出定量预测。

　　从武汉迁出人群的发病率反推武汉的疫情规模

　　受启发于外省市新冠病毒患者的死亡率明显低于武汉从而推断出武汉确诊数据可能失真，作者想到从武汉迁出人口的发病率倒推武汉的感染人数。

　　根据百度地图发布的迁徙大数据，1月14日到1月23日期间，以武汉为迁出点，取迁入人数最大的8个省市为分析样本，这8个省市期间从武汉迁入人数都超过2万，统计1月21日以来的确诊病例数。可以看到，在1月24日之前，累计确诊的人数跟武汉迁出的人口数没有明显的正比例关系，可以推测在前3天的确诊数据可能受各地确诊检测效率的影响。

　　图表2：从武汉迁出人口最多的8省市累计确诊数据

来源：百度地图，各地卫健委

　　进一步，统计8省市武汉迁入人口的发病率（平均每万武汉迁入人口的确诊数据），由于迁入的时间分布非常接近，如果各地确诊数据都比较可靠的话，发病率应该比较接近。为了反映各地确诊数据的一致性，统计每天8省市发病率的标准差/均值，这一比率越小，代表8省市确诊数据的一致性越高。

　　从图表3可以看到，8省市确诊数据在1月25日之后一致性明显提高，可能是因为确诊检测效率大幅提高了（1月24日首个新冠病毒检测试剂盒通过检验，随后多家企业开始生产）。发病率数据的一致性在1月28日达到最高，随后开始下降。

　　这是因为1月28日之后的确诊数据可能受当地第二代传染病人的影响，各地疫情发展的状况是不同的，之前的确诊病例主要是处于潜伏期的武汉迁入人员。因此取1月28日的发病率（1月22日至28日累计确诊数除以武汉迁入人口数）来反推武汉在1月14日至21日处于潜伏期的新冠病毒感染者比例。之所以取22日至28日一个区间的累计数据来反推，是因为单日的确诊数据可靠性比较低，很容易受确诊检验效率的影响，即有可能前几天发病的患者后几天才确诊。

　　取8省市1月28日发病率的中位数32人/每万武汉迁入人口作为参考（取平均数也很接近）在1月14日至21日期间，迁出的约300万+最终留在武汉的900万人口中，潜伏期的感染者人数期间平均值大约为 1100万*32/7天 = 5257人。

　　图表3：武汉迁入人口的发病率（平均每万武汉迁入人口的累计确诊人数）

来源：百度地图，各地卫健委

　　图表4：8省市武汉迁入人口发病率的标准差/均值

来源：百度地图，各地卫健委

　　简化微分方程组，建立数值模拟模型

　　逐一来理解四个微分方程：

　　1）第一个微分方程是刻画易感人群（S类）数量变化的。理论上讲，如果某一病毒基本传染数R0大于1，在不加防控的情况下它会传染几乎所有人群，但随着被传染的人和已治愈有抗体的人越来越多，易感人群数量越来越小，传播的速度就会下降，因为传染期人群与易感人群接触的机会会减少，因此理论模型中需要对易感人群的数量进行刻画。但是在现实的疫情中，除了极少数失控的疫情，绝大部分情况下，已感染的人数（含已治愈隔离的）占总人口的比例非常小，易感人群数量可以近似认为是不变的,等于总人口（S(t)=N）。因此这个微分方程可以省去。

　　2）第二个微分方程的意思是单位时间内潜伏期人群（E）的变化等于新增加的感染者减去由潜伏期过渡到传染期的数量。表达成离散的形式如下：

　　其中β是平均每个病人单位时间内的传染人数，则，整个传染期能传染多少个人，就是基本传染数R0的定义。

　　3）第三个微分方程是单位时间内传染期人群（I）的变化等于新由潜伏期过渡到传染期的人数，减去从传染期进入到隔离态的人数。表达成离散形式：

　　4)第四个微分方程是单位时间隔离态的人数等于新增的从传染期进入隔离态的人数。离散形式：

　　由此，我们得到3个离散的方程式，反映的是三个状态量的迭代关系，即如果知道所有参数值，就可以从t-1时刻推导出t时刻的值，确定一个起始点就可以得到整个序列。

　　数值模拟求解

　　三个迭代方程涉及到3个参数：β、γ1、γ2，根据它们的定义可以转化为另外三个参数：R0、TE、TI，其中:

　　1）R0、TI作为优化的目标参数，参考各类文献两者的数值模拟区间分别为[2,7]、[3,8].

　　2）对于潜伏期时间，各参考文献中大多认为SARS的潜伏期为5到6.5天，但明确表明SARS潜伏期内没有传染性，国家卫健委专家表示新冠病毒潜伏期为10天左右，但潜伏期可能有传染性，考虑到两者都是冠状病毒，专家也表示两者有同源性和很多相似的地方，将新冠病毒无传染性的潜伏期时间TE假定为6.5天（SEIR模型中潜伏期不具有传染性，假定新冠病毒10天左右的无症状潜伏期里后面3.5天有传染性，计入到传染期）。

　　3）约束条件：a.

　　b.在1月14日至21日期间潜伏期人群的平均人数在5257左右。两个待定参数，两个约束条件，理论上是可以求解的。

　　4）模拟优化过程：参数R0、TI各自模拟区间内随机取值，通过3个迭代方程可以计算出E(t)、I(t)、R(t)的时间序列（以2019年12月8日为起点出现首个确诊患者），通过时间序列可以计算λ以及1月14日至21日期间潜伏期人数，判断约束条件是否接近。大量模拟后取约束条件最接近的参数值，作为近似解。

　　参数优化结果：

　　基本传染数R0比之前大多估计都高一些，传染期在不加干预的情况下长达7天，可能与有大量的轻症患者以及发病前夕无症状时就具有传染性有关。因为没有症状或者症状很轻会使得病人不会停止活动、更不会住院，平均传染期就变长了。

　　三、高传染性下防控措施能达效果吗？

　　自然传染数是指在前期没有额外的防疫干预措施的情况下，平均每个病人能传染的人数。根据传播动力学[1]：

　　其中k是一个有传染能力的患者平均每天与易感人群的接触次数，b是每次接触传染成功的概率，D是可以传播的时间，可以认为等同于TI。

　　新冠病毒R0明显高于SARS主要是因为平均每个病人传播的时间更长

　　比较新型冠状病毒和SARS，考虑到新冠病毒传播初期刚好遇到春运，人员流动性高，假设单位时间接触次数比SARS时期高10%。

　　图表5：新型冠状病毒与SARS对比

　　可见，新冠病毒明显高于SARS主要是因为平均每个病人的传染期时间更长，而单次接触的传染概率可能略低于SARS。新冠病毒传染期格外长可能主要是因为其有大量的轻症患者，不会及时就医或采取措施，而会较长时间处于活跃的传染状态。

　　国家卫健委的专家：“说它是新型冠状病毒肺炎，不如说是新型冠状病毒感染，因为检测阳性的患者不见得一定出现肺炎。就像流感、重感冒，有的轻症患者不需要住院便能治愈。病例数虽多，但大多数都是轻症。”而恰恰是轻症患者较多增加了不小的防控难度。

　　对防控措施应该有信心

　　采取防控措施就是要降低k、b、D这几个参数值，从而降低传染数。虽然新增确诊数还在不断攀升，但从1月23日之后采取的一系列措施来看，我们应该对控制疫情有信心。在流行病传播模型里，只要将基本传染数R0降到1以下，疫情就不会失控。

　　假定R0就如前文测算的5.38，远高于SARS。首先，交通管制、取消聚会等措施使得人员流动大幅下降，单位时间的接触次数k下降70%甚至更多都有可能。

　　其次，通过宣传和提醒使得民众对疫情具有极高的警惕性，有轻微症状的患者大多也会尽快就医，从而减少传染期的长度，假设传染期的长度由7天降为4天，则D下降43%。

　　第三，通过佩戴口罩、消毒等防护措施，单次接触传染的概率假如能下降30%。则传染数就下降为5.38*（1-70%）（1-43%）（1-30%）=0.64，明显小于1，疫情就能较快控制住了。事实上，WHO专家曾表示，中国在2003年非典期间后期，严格的防疫措施使得SARS的基本传染数下降到0.4。

　　防控措施的效果体尚未体现到新增确诊数

　　为什么采取了这么严格的防控措施，新增确诊数还在不断攀升呢？首先，可能由于确诊检测试剂盒上市较晚，一些地方的基层医院缺乏确诊的条件或者确诊效率很低，造成确诊数据失真。可能存在前些天发病的病例过了几天才有条件确诊，这就可能造成近日新增确诊数大幅增加。

　　其次，1月24日采取了进出武汉的交通管制，之前大量的人口从武汉外流到各地，按照10天左右的潜伏期，现在是发病的高峰期，随着他们的确诊，跟他们密切接触者大多会进入隔离观察，隔离观察者中还是会有比较多的潜伏期感染者，因此接下来新增确诊数可能还会保持高位，直到这批隔离观察者的潜伏期过去，被确诊，之后的新增确诊数就会大幅下降。从时间上看，可能要到2月10日之后。

　　以上海为例，其近日公布的新增确诊病例中，大部分还是曾湖北居住或旅游，说明目前新增的确诊病例主要还是原来潜伏期的感染者，严格的管控措施产生的效果还没有反映到新增确诊数据上来。

　　图表6：上海每日新增确诊病例中有湖北居住或旅行史的占比

　　来源：上海卫健委

　　四、基于数学模型的全国新冠病毒患病人数预测：我们能打赢这场防疫战

　　利用前文建立的迭代方程组和估算的参数，可以对未来疫情发展进行预测。但在采取严格防控措施后，基本传染数下降到什么水平目前还没有数据可以估算，因为新冠病毒的潜伏期比较长，当前的新增确诊数据还更多地受之前潜伏期感染人数的影响。

　　疫情发展预测：中性与乐观假设

　　中性假设：

　　1）假设自1月23日后，民众减少出行和聚会，单位时间接触次数下降50%，由于提高了对疫情的警惕，有轻微症状的患者会及时就医，使得传染期的时长从7天降为4.5天，由于佩戴口罩、消毒等，单次接触传染的概率下降15%。则基本传染数将从5.38降为1.47。

　　2）假设自2月15日起，由于隔离观察制度的严格实施，大部分新增病例都在隔离观察名单里，基本失去了再次传染的可能，基本传染数下降到0.5，稍高于非典后期（WHO研究认为非典后期在采取严控措施后R0降到了0.4）。

　　乐观假设：

　　1）假设自1月23日后单位时间平均接触次数下降55%，其余与中性假设1）不变。

　　2）假设自2月10日起，隔离观察制度开始起作用，基本传染数下降到0.4。

　　在这两种假设下，可以对未来累计病例数、新增病例数等作出预测，结果如下：

　　1）当前（2月1日）估算的累计病例人数为22732人，与公布的累计确诊病例相差较大，可能主要是因为前期缺乏高效的确诊技术手段，叠加上流感季，武汉医疗资源严重不够，加上新冠病毒的轻症患者较多，症状跟感冒类似，可能未确诊前就已经痊愈了。

　　2）预计未来新冠病毒总患病人数最高将达到88500人（中性假设下）或58000（乐观假设下）人左右，新增确诊人数将在2月17日（中性假设下）或2月11日（乐观假设下）左右出现拐点。

　　3）当前新增确诊病例数大幅攀升（新增病例超2500）可能与检测确诊效率提高有关，预计将很快回落（模型中性假设下预测值是1700左右）。1月25日之前，公布新增确诊数大幅低于模型预测值，26日之后则大幅高于模型预测值，可能是因为确诊检测试剂盒供应充分了，将之前已经发病而没有确诊的病例确诊了。如果看1月24到2月1日期间总的新增确诊量，则模型预估和官方公布的数据很接近，分别为13469人和13582人。

　　图表7：中性假设下模型预测的新冠病毒确诊数和实际公布的新增确诊数

数据来源：模型预测，国家卫健委

　　图表8：乐观假设下模型预测的新冠病毒确诊数和实际公布的新增确诊数

数据来源：模型预测，国家卫健委

　　疫情影响分析：需要足够重视，但无需恐慌

　　比较过去20年发生的比较大的疫情，这次新型冠状病毒的传染性比较强，可能只弱于2009年发源于美国的H1N1流感。但从死亡率来看（以湖北省外的数据为参考），或许只高于普通流感，比其他几类都明显更低。但即使采用湖北省外的数据，新冠病毒的致死率也是普通流感的10倍左右。

　　因此，政府和民众采取足够重视的防控和防护措施是必要的，但也无需恐慌，2009年席卷全球的H1N1疫情，感染人数和死亡人数都非常之高，但由于致死率较低，也没有引起很大的恐慌。

　　图表9：过去20年全球较大的几次疫情情况

数据来源：WHO，中国疾控中心，美国疾病控制和预防中心（CDC）

　　说明：2009年H1N1美国的数据是CDC事后建模预估的，当时由于疫情发展失控，没有统计和确认病例数了，据事后的研究，当时H1N1在墨西哥境内死亡率为2%，在其他地区的死亡率为0.1%左右。

　　股市短期可能反应过度，若超跌是买入机会

　　2009年美国爆发H1N1流感疫情，最后完全失控。期间美国股市仅在2009年5月、6月有小幅回调，幅度-6%左右。美国股市对疫情冲击反应较小。

　　2003年SARS疫情期间，港股和A股做出了明显反应，下跌幅度-14%左右，但随后都反弹明显，且反弹时间领先于疫情拐点。可见，股市可能在情绪影响下短期反应过度，随后会修复。

　　在1月23日发布的《关注新型肺炎“概念股”的估值偏差》文中提到，“事件性的、一次性的利好，即使对当期业绩影响很大（盈利翻倍），对股价的理论影响也很小”，其实反过来也一样，事件性、一次性的利空影响，对公司价值的影响也是较小的，股价不应该反应过度。当然利空与利好不同的是：短期利空（如收入剧减、资金回笼困难等）可能引起公司短期周转困难从而引发更大风险。因此也建议政府相关部门出台一些措施，对疫情中遭遇短期困难的企业给予支持和帮助。

　　图表10：2009年H1N1流感爆发期间美股表现

数据来源：Wind，WHO

　　图表11：SARS期间美股、港股、A股表现

　　数据来源：Wind，中国网