只为梦而战

0%

SARS-CoV-2安徽疫情的微分方程模型预测

发表于

本文pdf：link

本文的结果预测于2月21日。

0.符号约定

  • 某种人群的数量-时间关系函数为$X(t)$，其中 X 表示该人群的名称首字母。
  • 患者传染率为$\beta$，患者治愈率$\gamma$，被传染之后成为患者而不是潜伏者，或者由潜伏者成为患者的概率为$\alpha$（由于潜伏者只有变为患者后才会被医院收治，故忽略潜伏者治愈率）。
  • 患者单位时间内接触人数$r_1$，潜伏者单位时间内接触人数$r_2$（本文的算法中单位为人/天）。
  • 安徽总人数$N=6.3236\times10^7$（单位：人）。

1.假设与约定

  • 附件中的数据真实可靠
  • 安徽区域的人数近似稳定，不会大幅度地向外或向内流动
  • 潜伏者在出现症状前，不会被检查出携带病毒并入院收治
  • 死亡人数对疫情模型几乎没有影响。事实上，安徽地区截至2.19的患者死亡率仅为约$0.6\%$
  • 本模型仅用于计算不超过一个月的情况。由于安徽的疫情发展迅速（事实上根据预测，3月上旬就将结束），这已经足够。

2.传染病模型简述

考虑新冠肺炎疫情复杂，选择$\text{seir}$模型作为基础。
将人群分为如下几类：

  • $\text{S (Susceptible)}$，易感者。指未得病者，但缺乏免疫能力，与感染者接触后容易受到感染。
  • $\text{E (Exposed)}$，潜伏者。指接触过感染者，但暂无能力传染给其他人的人，对潜伏期长的传染病适用。
  • $\text{I (Infectious)}$，确诊患者。指染上传染病的人，可以传播给 S 类成员，将其变为 E 类或 I 类成员。
  • $\text{R (Recovered)}$，治愈者。指被隔离或因病愈而具有免疫力的人。如免疫期有限，R 类成员可以重新变为 S 类。包括病死者。
  • 需要注意的是，在本文中，为了方便调参，I代表累计确诊人数。

其满足恒等关系$S(t)+E(t)+I(t)=N$，与微分方程

由于SARS-CoV-2的潜伏者具有感染能力，我们对微分方程进行调整：

具体计算过程中，我们取天为单位。进行迭代求解：

事实上，$r_1,r_2,\gamma$为定值时，将会导致结果极不准确。这不难理解，因为戒严尺度、医疗经验皆发生了变化。

所以我们令$r_1=r_1(t),r_2=r_2(t),\gamma=\gamma(t)$为关于时间的函数。

模型参数

固有参数

  • 人口：$N=6.3236\times10^7$
  • $S_0=N,I_0=9,R_0=0$（初始状态下E和I极少，故S0大致视作N）。
  • 传染率：$\beta=1.7\%$
  • 潜伏者发病概率：$\alpha=13\%$
  • 固有参数中部分选自文献[1]

待定参数

这是一些较难统计的参数。

  • $E_0$，对模型进行二分拟合后，发现取$E_0$在$660$附近时拟合效果较好。

  • 接触人数

    令$r_1$、$r_2$为$a\exp(-bt)$的形式。对已知数据进行拟合，经过二分尝试，可取值如下：

    $r_1=7\exp(-0.25t),r_2=13.2\exp(-0.48t)$

  • 日治愈率
    令$\gamma(t)=a\exp(bt)$，对已知数据进行拟合，经过二分尝试，取$\gamma(t)=0.004\exp(0.1t)$时拟合结果较好。由于约定本模型并不用于时间过长的情形，而事实上，预测的疫情结束事件前，都有$\gamma(t)\le1$，显然符合“治愈率”的概念。

模型实现

使用Excel计算的结果：

img

总结

预计3月10日左右结束疫情，安徽将连续多日出现零增长，潜伏者数量骤降。

风险提示

  • 本文考虑到的因素有限，包括人口的流动、各地医疗水平的不同等都被忽略。
  • 二分拟合的方法难免存在误差。
  • 由于安徽地区患者较少，偶然因素大，有几率发生突发的转折，也有几率数据中的趋势仅仅是一小段随机的曲线。
  • 众所周知，由于政府对确诊患者文件的变动，可能附件中的数据并不完全准确。

模型检验（3.9日晚补充）

img

可以看出，疫情中后期表现基本与模型完全吻合，并正确地预测了拐点、结束时期。

前期误差原因

  • 一方面即“风险提示”中第四点。
  • 前后期存在一些被本文忽略，但实际上不容忽视的其他因素，它在前期有重要影响，但中后期影响力则淡化。

附件（安徽省疫情历史数据）

日期 确诊 痊愈 Δ TMP 痊愈率
1月23日 9 0 0 1 0
1月24日 15 0 0 9 0
1月25日 39 0 0 15 0
1月26日 60 0 0 39 0
1月27日 70 0 0 60 0
1月28日 106 0 0 70 0
1月29日 152 2 2 106 0.018867925
1月30日 200 2 0 152 0
1月31日 237 3 1 200 0.005050505
2月1日 297 3 0 237 0
2月2日 340 5 2 297 0.006802721
2月3日 408 7 2 340 0.005970149
2月4日 480 14 7 408 0.017456359
2月5日 530 20 6 480 0.012875536
2月6日 591 23 3 530 0.005882353
2月7日 665 34 11 591 0.019366197
2月8日 733 47 13 665 0.020602219
2月9日 779 59 12 733 0.017492711
2月10日 830 73 14 779 0.019444444
2月11日 860 88 15 830 0.019815059
2月12日 889 108 20 860 0.025906736
2月13日 910 128 20 889 0.025608195
2月14日 934 166 38 910 0.04859335
2月15日 950 198 32 934 0.041666667
2月16日 962 232 34 950 0.045212766
2月17日 973 266 34 962 0.046575342
2月18日 982 293 27 973 0.038189533
2月19日 986 361 68 982 0.098693759
2月20日 987 424 63 986 0.1008
2月21日（预测日） 988 500 76 987 0.134991119
2月22日 989 567 67 988 0.137295082
2月23日 989 620 53 989 0.125592417
2月24日 989 648 28 989 0.075880759
2月25日 989 692 44 989 0.129032258
2月26日 989 727 35 989 0.117845118
2月27日 989 764 37 989 0.141221374
2月28日 990 814 50 989 0.222222222
2月29日 990 824 10 990 0.056818182
3月1日 990 868 44 990 0.265060241
3月2日 990 888 20 990 0.163934426
3月3日 990 917 29 990 0.284313725
3月4日 990 941 24 990 0.328767123
3月5日 990 956 15 990 0.306122449
3月6日 990 970 14 990 0.411764706
3月7日 990 979 9 990 0.45
3月8日 990 981 2 990 0.181818182
3月9日（检验日&疫情基本结束） 990 984 3 990 0.333333333

数据来自疫情日报。

参考文献

[1] 花长春、陈礼清、田玉铎.利用SEIR模型推演湖北、非湖北和全国疫情拐点 ——2020年突发风险系列[EB]/[OL].https://finance.sina.cn/2020-02-07/detail-iimxyqvz0930302.d.html .(2020-02-07)[2020-02-21]

[2]央视网“疫情日报”搜索[EB]/[OL].https://search.cctv.com/search.php .[2020-02-21]