Get best of the week

数学家认为COVID-19大流行，或者经典SEIRD模型不起作用的原因

摘要，或关于年轻科学家的闲暇


在过去的几周中，我和我的同事通过使用各种非线性回归方法，对俄罗斯COVID-19流行病发展预测的准确性进行了竞争，从而结束了工作日。如果对明天的预测不可避免地要好，那么对一个多星期的预测只能从总体上反映现实。似乎一切都很清楚：有流行病学模型，有优化方法，有足够详细的数据，足以将它们组合在一起并提前一个月甚至六个月获得准确的预测。在本文中，我将分享我对经典SEIRD模型有什么问题以及如何修复它的想法。而且，当然，我将揭开与您一起笼罩我们未来的秘密面纱。

坐下来，愤怒的男人等待着我们，那些知道微分方程是什么的人（其余的人，都附有漂亮的图画）。


上图以对数标度显示了俄罗斯和三个欧洲国家按感染对数确定的COVID-19确诊病例总数。文本中将进一步解释。


不明飞行物护理分钟


大流行的COVID-19是由SARS-CoV-2冠状病毒（2019-nCoV）引起的潜在严重的急性呼吸道感染，已在全球正式宣布。关于Habré的很多信息都涉及此主题-始终记住，它既可靠又有用，反之亦然。

我们敦促您不要批评任何已发布的信息。


官方资料

, .

洗手，照顾亲人，尽可能呆在家里并远程工作。

阅读有关以下内容的出版物：冠状病毒 | 远程工作

型号SEIRD


SEIRD流行病模型属于所谓的类别 隔间模型，其本质是将总体分为几类（英语隔间），在我们的示例中：$ S $（英语易感人群）-易感人群，$ E $（英，暴露） -那些谁拥有这种疾病的潜伏期，$我$（英语传染性）-生病了，$ R $（已恢复英语）-已恢复，$ D $（英语死了）-死了。然后，将每个组的大小与微分方程系统中的一个变量进行比较，解决该问题，您可以预测该流行病的动态。SEIRD模型有很多修改，例如，SEIR是简化模型，没有单独考虑恢复和死亡的情况。要熟悉其他模型，我可以推荐一篇有关该主题的好文章。

一点理论


首次以三个变量的微分方程组形式的流行病模型 $ S，I，R $出现在1927年W.Kermak和A.McKendrick 的作品中。
这些微分方程的形式为：

$\begin{align} \frac{dS}{dt }&=-\beta \frac{SI}{N},\\ \frac{dI}{dt}&= \beta \frac{SI}{N}-\gamma I,\\ \frac{dR}{dt}&= \gamma I, \end{align}$


除了我们熟悉的变量外，还出现以下常量： $N = S + I + R$ -总人口数， $\beta$ -感染传播率， $\gamma$-恢复速度。

Kermak和McKendrick方程的含义如下：易感人群的数量与他们的数量乘以人口中感染的平均比例成比例地减少$I/N$，受感染人数以相同的速度增长，并根据其中一些事实进行了调整 $\gamma I $感染的数量减少，因此恢复期和疗养次数增加。值得注意的是，SIR模型包含非线性$SI$，因此，方程组的解析解通常变得不可能，但是幸运的是，数值微分方法可以轻松地完成此任务。

在这里添加另一个变量$E$ （潜伏期患病人数），我们得到SEIR模型：

$\begin{align} \frac{dS}{dt }&=-\beta \frac{SI}{N},\\ \frac{dE}{dt}&= \beta \frac{SI}{N}-\kappa E,\\ \frac{dI}{dt}&= \kappa E-\gamma I,\\ \frac{dR}{dt}&= \gamma I, \end{align}$


出现另一个常数的地方 $\kappa$-疾病从潜伏期到开放的转变速度。图取自文章。



显而易见，SEIR模型并不非常适合描述COVID-19，仅因为该模型中存在隐藏的感染携带者$E$不具有传染性。可以通过引入Pengpeng等人的方法纠正这种缺陷。可选参数$\theta$，与患病的相比，表征了潜在的感染携带者的感染程度。我们将尝试将其应用于当前的流行病的修改后的SEIR模型如下所示：

$\begin{align} \frac{dS}{dt }&=-\beta \frac{S(I + \theta E)}{N},\\ \frac{dE}{dt}&= \beta \frac{S(I + \theta E)}{N}-\kappa E,\\ \frac{dI}{dt}&= \kappa E-\gamma I,\\ \frac{dR}{dt}&= \gamma I. \end{align}$


乍一看，最终的模型有望完全可信。

SEIR模型的数值实验


对于建模，我们将尝试采用以下参数，重点放在开放数据上。假设该疾病平均持续14天（至少轻度形式持续多长时间，最多占病例的80％），我们发现其价值$\gamma=1/14=0,0714$。会接受$\beta=3/14=0,2143$。大小$\theta=0,6$从Pengpeng等人借来的。假设平均潜伏期为3天，$\kappa=1/3 = 0,33 $。我们接受俄罗斯的人口平等$N = 144,5\cdot10^6$人。

作为初始条件，我们使用4月2日俄罗斯的数据，当时3月底采取的限制感染扩散的措施将产生效果，即：

$\begin{align} &S_0=3548,\\ &I_0= 3283,\\ &E_0=0,5 I_0. \end{align}$


评分 $E_0$我们采取了相对任意的态度，这并不重要，因为如您所知，会出问题。

从4月2日到4月24日（含1天）使用Euler方法进行建模的结果是，我们得到了如下图所示：左侧为线性标度，右侧为对数。



圆形标记标记的是俄罗斯实际病例总数的真实数据，方形标记标记的是患者人数。乍一看，结果看起来不错，但有一个例外：使用模型的参数，我们显然没有猜测。在这里，优化方法对我们有帮助。

优化它


如果读者不熟悉优化方法，那么它们是允许找到某些目标函数最小值的算法。在我们的案例中，我们面临非线性回归问题：如何选择微分方程参数的向量$\mathbf{x} = (\beta, \gamma, \kappa, E_0)^\top$ 这样微分方程的解点集 $F$ 尽可能靠近一组观察点 $\mathcal{F}$。

我们使用标准偏差作为模型误差的度量。目标函数将采用以下形式

$ f(\mathbf{x}) = \frac{1}{M}\sqrt{\sum_{i=1}^{M}(F_{i} - \mathcal{F}_{i})^2 + \sum_{i=1}^{M}(G_{i} - \mathcal{G}_{i})^2}, $


哪里 $M$ -点数 $F$ -模型给出的感染总数， $\mathcal{F}$ -实际案件总数， $G$ -目前的患者人数，即模型 $\mathcal{G}$-当前活动案例的实际总数。

使用MATLAB中的优化工具箱，我们将模型参数调整为观测数据。结果，我们获得了下图所示的解决方案。

image

乍一看，一切都很好。差异竟然相等$f(\mathbf{x}) = 131,98$，并且在“凸海眼”上，解决方案的拟合效果很好。让我们看一下接收到的参数：

$\begin{align} &\beta = 0,374,\\ &\gamma = 0,0117,\\ &E_0 = 7,84\cdot 10^6,\\ &\kappa = 4,81\cdot 10^{-5}. \end{align}$



4月24日 记录的近800万潜在患者的价值（约有60000例病例）令人怀疑。我们还发现，到疾病活跃期的平均过渡时间为$1/\kappa = 2079$天。

为什么会这样呢？如果长时间分析曲线的形状，一切都会变得清晰。为此，请使用带有“合理”参数的SEIR模型并在长时间内对其进行仿真（在本实验中，我采用了新的含义$\beta=0,186$）：



与案例总数相对应的曲线在线性范围内具有特征性S形。优化程序试图使该曲线形成曲线。此外，具有“合理”参数的预测本身是可怕的-据此预测，到9月该国将有近90％的人口生病-如果您查看其他国家的结果（与本文开头的图片相同，只是线性范围内），这显然是不现实的）：



在这里，我比较了案件数量排名前5名的三个欧洲国家和俄罗斯。可以看出，与流行病的流行速度相比，我们落后了大约一个月，而与SEIR模型获得的结果相比，这三个国家的三个月病例总数的增长几乎是线性的（甚至慢于线性）。这就提出了三个问题：

  1. ?
  2. SEIR, ?
  3. , , ?

我将首先回答第三个问题。当我们做出预测时，我们面临着一个相当不愉快的任务：构建模型所依据的数据并不理想-它们包含误差，噪声，并且基于它们建立的模型也包含一些误差。当我们使用模型点继续时间序列时，如果我们预测随时间增加的函数，则误差会累积，而且会很快累积。这就是我们的情况。而且，该模型是反映真实情况的模型，非常有限。在一个新的大城市中，这种流行病的突然发展，使用更有效的治疗方法，改变收集信息的方法-所有这些都可能在真实数据中造成很多错误，以至于长期的预测将完全脱离现实。

拐杖和自行车：我们修改了SEIR模型


让我们尝试回答这个问题，为什么流行病的增长会放缓到线性的趋势。鉴于我们现在拥有的受感染人数众多，与人与人之间沟通速度有限相关的规模经济在其中发挥着重要作用。

更准确地说，我们记得：SEIR模型中的病例数与总体中的平均病例数成正比$I/N$。此规则在少数人群中可以很好地起作用，每个人都可以与每个人进行交流，并且患者分布均匀。在现实中，尤其是在成千上万的人的规模上，如果您带两个随机生病的人，事实证明他们不仅彼此之间从未交流过，也没有看到彼此，他们甚至没有乘坐同一辆地铁。实际上，他们生活在不同的城市。使他们团结在一起的是一系列社会纽带，这导致该病毒已传播给他们。

例如，我以细胞自动机的形式构建了一个流行病模型，其中每个细胞仅与4个相邻细胞相互作用。这等同于人口中的每个人都有4个社会联系这一事实-对于人类来说这是一个很小的数目，但是限制社会联系的效果越快地显现出来。每次迭代的概率为0.1，可以感染被感染细胞的4个邻居中的每个邻居。该病平均持续14天。下图显示了一个200x200单元池的仿真结果，其中$k$是迭代次数。



蓝色表示易感，黄色-有病，绿色-已恢复。最有趣的是案件数量图表的外观。而且它们看起来像是计划的：在短暂的次指数级增长阶段之后，就像在SEIR模型中一样，线性增长处于延长阶段，就像在现实中一样。



我的目标不是要获得与现实情况定量相似的照片。如果您想要更高的信誉，我可以推荐最近在Habré上发表的Sergey Potekhin项目。对于细心的读者，下面给出线性增长的更严格证明。

大量人口中线性增长流行病定理的证明
: $d$- . . , 20 ( ) $d \approx 4$. $d$- . $n^{\frac{1}{d}}$, , $P$, $2P$. , ,

$n_{k+1}=\left(n_k^{\frac{1}{d}} + 2P \right)^d,$



$n$ : $n_k = n(t_k); n_{k+1} = n(t_k + 1)$. , :

$n'_{k+1}=n'_{k}\left(1 + \frac{2P}{n_k^{\frac{1}{d}}} \right)^{d-1}.$



, $n_k$ $n'_{k+1}=n'_{k}$, , .

4% 16- .



关于全球流行病的统计数据，我们将看到同一件事：尽管经典模型向我们保证了病例数呈指数级增长的事实，但在一个月的时间里，增长一直是线性的。



好奇的任务
  1. COVID-19, offline .
  2. .1 ?
  3. .1 2, .


现在介绍SEIR模型的修改。我们可以做的最简单的事情是将非线性分量乘以某个函数，具体取决于患者的数量。在少数情况下，此函数应接近于1；在大量情况下，应渐近趋于零。最简单的合适人选是

$\mathcal{\varphi}(I,E) = e^{-\alpha(I +\theta E)^{K_0}}.$


参数选择 $\alpha$ 和 $K_0$可以补偿原始模型中的指数增长。

添加到模型（以获得更多信息）和组件$D$-死亡人数。获取SEIRD模型的修改版本：

$\begin{align} \frac{dS}{dt }&=-\beta \frac{S(I + \theta E) \mathcal{\varphi}(I,E)}{N},\\ \frac{dE}{dt}&= \beta \frac{S(I + \theta E) \mathcal{\varphi}(I,E)}{N}-\kappa E,\\ \frac{dI}{dt}&= \kappa E-\gamma I-\mu D,\\ \frac{dR}{dt}&= \gamma I,\\ \frac{dD}{dt}&= \mu I. \end{align} $



仿真结果如下图所示。



与原始模型相比，标准误差几乎没有变化。参数值已经很实际了。为了方便起见，我将活动阶段的初始感染人数指定为$I_0$。

$ \ begin {align}＆\ beta = 0.219，\\＆\ gamma = 0.0102，\\＆E_0 = 0.13 \ cdot I_0，\\＆\ kappa = 1/3，\\＆\ mu = 1 ，13 \ cdot 10 ^ {-3}。 \结束{align} $



该模型也很好地对导数进行插值-病例数和死亡人数每天增加的值。



让我们尝试做出预测。我们将预测期设为2个月，并继续使用优化程序找到的参数对解决方案进行建模。



乍一看，还不错，但是您不能希望对您的祖国有这样的预测：新病例的数量将继续减少，但总数将继续增加。在这种情况下，只有借助疫苗或等到几乎整个人口都病了才能停止流行。每天新增死亡人数约为200。这清楚地表明，如果我们不加强抗击流行病的措施，将会发生什么。那是等待我们的吗？为了这个不光明的未来，我们许多人都辛苦地坐在家里，买了荞麦和厕纸？

下面，我将考虑两种情况，再看看2020年4月28日未来几个月的雾蒙蒙的距离，我不能肯定地说哪种会进一步发展。现在，在打破新案例曲线的时刻，我们正处于预测某些事物的双重问题。

美国脚本


世界霸主处于极其令人羡慕的地位。他在做出关键决定以至于一开始就减慢了疫情的发展速度时被耽搁了，但他仍然无法应对新病例的自然增长。

修改后的SEIRD模型在3月2日开始的前33点接受了培训，从正负实际预测了4月份的流行过程。



如您所见，四月份的增长几乎是线性的。该模型稍微夸大了4月份的死亡率，但总体情况是正确的。



此图显示了美国新病例和死亡人数的每日增加。它与模型对俄罗斯的预测非常相似。

德国情景


有纪律的德国人设法扭转了曲线，他们的增长比线性慢。此外，为使模型适用，我必须在4月6日手动增加恢复率$ \伽玛$是SEIRD模型的1.7倍，否则无法解释这种情况下的急剧下降。



从3月10日开始，在前27个点对模型进行了训练。我还更改了非线性函数。对于德国，与时间有关的指数更好：

$ \ mathcal {\ varphi}（t）= K_0 e ^ {-\ alpha t}。


这种类型的功能表明社会联系中断的数量不断增加，从而导致感染扩散。在这里，您可以清楚地看到自我隔离的好处。



上面显示了新病例和死亡人数的每日增加。与美国一样，真实数据包含明显的波动，为期7天。这意味着在周末接触的数量增加，因此，被感染的数量增加。[UPD：相反，它正在下降，因为Yandex自隔离指数可能间接证明。]

结论


做出短期和长期的预测不仅是对好奇心的致敬。如果发生流行病，您需要知道应该准备多少张床，要生产多少台呼吸机以及要在几个月内为医生存放个人防护设备。官员们必须了解所采取的措施是否足够或是否应实行新的禁止和限制措施。理想情况下，该模型应该很好地反映现实情况，以便可以看到每个新采用的措施的作用力，然后可以加强有用的措施并取消原来无用的决定。

尽管有些保留，但仍可以用字母描述尚未获得COVID-19免疫力的每个人$ S $，各种各样的病人-信 $我$等而且，SEIRD模型甚至可以解释某些内容。但是她可以很粗略地预测遥远的未来。

在本文中，我故意引用了当我们根据流行病的动态重复美国的命运时的负面情景。如果这种情况真的成立，到6月底，我们将有30万多例登记的疾病病例和1万多例死亡。尽管这种情况没有实现是有先决条件的，但我还是建议您根据以下原则与之联系：“希望最好，为最坏的情况做准备”。俗话说，如果美国宇航局的最好的头脑已经采取了对付美国流行病的斗争，那确实是一件坏事。

到目前为止，剩下的就是减少对公共场所的访问，使用使用合适的呼吸器，洗手，在街上外出后别忘了用酒精擦拭智能手机，并遵循其他简单建议。

但是仍然可以做出更准确的预测吗？当然是。但是大约在其他时间。

如果您想使用源代码并自己提出开发方案，那么这里是github的链接。

More articles:


All Articles