🔄 🔆 👨🏿‍✈️ La pandemia de COVID-19 a través de los ojos de un matemático, o por qué el clásico modelo SEIRD no funciona ✌🏾 📯 🤷🏿

Resumen, o sobre el ocio de los jóvenes científicos.

En las últimas semanas, mis colegas y yo hemos terminado el día de trabajo compitiendo en la precisión del pronóstico para el desarrollo de la epidemia de COVID-19 en Rusia utilizando varios métodos de regresión no lineal. Y si el pronóstico para mañana inevitablemente resulta ser bueno, entonces el pronóstico para un período de más de una semana refleja la realidad solo en términos generales. Parece que todo está claro: hay modelos epidemiológicos , hay métodos de optimización , hay datos suficientemente detallados : es suficiente combinar esto y obtener un pronóstico preciso para un mes, o incluso seis meses, por adelantado. En este artículo compartiré mis pensamientos sobre lo que está mal con el modelo clásico de SEIRD y cómo solucionarlo. Y, por supuesto, abriré el velo del secreto que envuelve nuestro futuro contigo.

Siéntese, un matan furioso nos espera para aquellos que saben cuáles son las ecuaciones diferenciales (por lo demás, se adjuntan bellas imágenes).

La figura anterior muestra el número total de casos confirmados de COVID-19 en una escala logarítmica para Rusia y los tres países europeos en el top 5 en términos de la cantidad de infectados. La explicación está más adelante en el texto.

Minuto de atención ovni

La pandemia COVID-19, una infección respiratoria aguda potencialmente grave causada por el coronavirus SARS-CoV-2 (2019-nCoV), se ha anunciado oficialmente en el mundo. Hay mucha información sobre Habré sobre este tema; recuerde siempre que puede ser confiable / útil, y viceversa.

Le instamos a que sea crítico con cualquier información publicada.

Fuentes oficiales
C
C
()

, .

Lávese las manos, cuide a sus seres queridos, quédese en casa siempre que sea posible y trabaje de forma remota.

Leer publicaciones sobre: coronavirus | trabajo remoto

Modelo SEIRD

El modelo epidémico SEIRD pertenece a la clase de los llamados modelos compartimentales , cuya esencia es dividir a la población en varios grupos ( compartimentos ingleses ), en nuestro caso:

$inline$ ( Inglés susceptible) - susceptible,

$inline$ ( Eng. Expuestos) - aquellos que tienen la enfermedad en el período de incubación,

$inline$ ( Inglés infeccioso) - enfermo,

$inline$ ( Inglés recuperado) - recuperado,

$inline$ ( Inglés muerto) - los muertos. Luego, el tamaño de cada grupo se compara con una variable en el sistema de ecuaciones diferenciales, resolviendo qué, usted puede predecir la dinámica de la epidemia. Hay muchas modificaciones del modelo SEIRD, por ejemplo, SEIR es un modelo simplificado que no considera por separado los casos de recuperación y muerte. Para conocer otros modelos, puedo recomendar un buen artículo sobre este tema.

Poco de teoría

Por primera vez, un modelo epidémico en forma de un sistema de tres ecuaciones diferenciales para variables

$inline$ apareció en el trabajo de W. Kermak y A. McKendrick en 1927.
Estas ecuaciones diferenciales tienen la forma:

$\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}$

donde, además de las variables que nos son familiares, aparecen las siguientes constantes:

$inline$ - tamaño total de la población,

$\beta$ - tasa de transmisión de infecciones,

$\gamma$ - velocidad de recuperación.

El significado de la ecuación de Kermak y McKendrick es el siguiente: el número de susceptibles disminuye en proporción a su número multiplicado por la proporción promedio de infectados en la población

$inline$ , el número de infectados crece al mismo ritmo, ajustado por el hecho de que algunos de ellos

$\gamma I$ recuperación, y el número de convalencias aumenta debido a una disminución en el número de infectados. Vale la pena señalar que el modelo SIR contiene no linealidad

$inline$ , por lo que la solución analítica del sistema de ecuaciones se vuelve generalmente imposible, pero, afortunadamente, los métodos de diferenciación numérica pueden hacer frente fácilmente a esta tarea.

Agregar otra variable aquí

$inline$ (el número de personas con la enfermedad en el período de incubación), obtenemos el modelo 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}$

donde aparece otra constante

$\kappa$ - la tasa de transición de la enfermedad desde la etapa de incubación a la abierta. Figura tomada del artículo .

Está claro de inmediato que el modelo SEIR no es muy adecuado para describir COVID-19, aunque solo sea porque en este modelo hay portadores ocultos de infección.

$inline$ No contagioso. Esta deficiencia se puede corregir introduciendo, siguiendo a Pengpeng et al. parámetro opcional

$\theta$ , caracterizando el grado de infectividad de los portadores latentes de infección en comparación con los enfermos. El modelo SEIR modificado, que intentaremos aplicar a la epidemia actual, se verá así:

$\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}$

A primera vista, el modelo resultante promete ser completamente creíble.

Experimento numérico con el modelo SEIR

Para el modelado, intentaremos tomar los siguientes parámetros, centrándonos en los datos abiertos . Suponiendo que la enfermedad dura 14 días en promedio (al menos cuánto dura la forma leve, que representa hasta el 80% de los casos), encontramos el valor

$\gamma=1/14=0,0714$ . Va a aceptar

$\beta=3/14=0,2143$ . Magnitud

$\theta=0,6$ tomado de Pengpeng et al . Dado el período de incubación promedio de 3 días, tome

$\kappa=1/3 = 0,33$ . Aceptamos la población de Rusia igual

$N = 144,5\cdot10^6$ persona.

Como condiciones iniciales, utilizamos los datos para Rusia el 2 de abril, cuando las medidas introducidas a fines de marzo para limitar la propagación de la infección tendrían su efecto, a saber:

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

Clasificación

$inline$ lo tomamos de manera relativamente arbitraria, y no importa, porque, como saben, algo saldrá mal.

Como resultado del modelado por el método de Euler con un paso de 1 día del 2 al 24 de abril inclusive, obtenemos los gráficos que se muestran a continuación: a la izquierda en la escala lineal, a la derecha en el logarítmico.

Los marcadores redondos marcan datos reales sobre el número total de casos en Rusia, marcadores cuadrados sobre el número de pacientes. A primera vista, los resultados se ven bien, excepto uno: con los parámetros del modelo, claramente no lo adivinamos. Y aquí los métodos de optimización nos ayudan.

Optimizarlo

Los métodos de optimización, si el lector no está familiarizado con ellos, son algoritmos que permiten encontrar el mínimo de alguna función objetivo. En nuestro caso, tenemos ante nosotros el problema de la regresión no lineal: cómo elegir el vector de parámetros de la ecuación diferencial

$\mathbf{x} = (\beta, \gamma, \kappa, E_0)^\top$ para que el conjunto de puntos de solución de la ecuación diferencial

$inline$ estaba lo más cerca posible del conjunto de puntos de observación

$\mathcal{F}$ .

Usamos la desviación estándar como una medida del error del modelo. La función objetivo tomará la forma

$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},$

Dónde

$inline$ - el número de puntos

$inline$ - el número total de casos de infección que da el modelo,

$\mathcal{F}$ - el número total real de casos,

$inline$ - el número de pacientes en este momento, que da el modelo,

$\mathcal{G}$ - el total real de casos activos actuales.

Usando la Caja de herramientas de optimización en MATLAB, ajustaremos los parámetros del modelo a los datos de observación. Como resultado, obtenemos la solución que se muestra en la figura a continuación.

A primera vista, todo está bien. La discrepancia resultó igual

$f(\mathbf{x}) = 131,98$ , y "en el ojo de mar convexo", el ajuste de la solución salió bien. Veamos los parámetros recibidos:

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

El valor de casi 8 millones de pacientes latentes, con alrededor de 60 mil casos registrados
el 24 de abril, es algo dudoso. También encontramos que el tiempo de transición promedio a la fase activa de la enfermedad es

$1/\kappa = 2079$ dias.

¿Por qué pasó esto? Todo se aclarará si analizamos la forma de la curva en una escala de tiempo larga. Para hacer esto, tome nuestro modelo SEIR con parámetros "plausibles" y simúlelo durante un largo período de tiempo (en este experimento adquirí un nuevo significado

$\beta=0,186$ ):

La curva correspondiente al número total de casos tiene una forma de S característica en una escala lineal. El programa de optimización trató de darle a esta forma una curva. Además, el pronóstico en sí mismo con parámetros "plausibles" es horrible - según él, casi el 90% de la población del país estará enfermo en septiembre - obviamente no es realista si se miran los resultados de otros países (la misma imagen que al principio del artículo, solo en una escala lineal ):

Aquí comparo los tres países europeos en el top 5 en el número de casos, y Rusia. Se puede ver que estamos aproximadamente un mes detrás del ritmo de la epidemia, y que durante un mes, el crecimiento en el número total de casos en los tres países es casi lineal (e incluso más lento que lineal), en contraste con los resultados obtenidos en el modelo SEIR. Esto plantea tres preguntas:

?
SEIR, ?
, , ?

Comenzaré respondiendo la tercera pregunta. Cuando predecimos algo, nos enfrentamos a una tarea bastante desagradable: los datos sobre los que construimos el modelo no son ideales: contienen errores, ruido, y el modelo construido sobre su base también contiene algunos errores. Cuando continuamos la serie temporal con nuestros puntos de modelo, el error se acumula, y bastante rápido, si predecimos una función que aumenta con el tiempo. Y este es el caso en nuestro caso. Además, el modelo es un modelo que refleja la situación real es muy limitada. El repentino desarrollo de la epidemia en una nueva gran ciudad, el uso de un método de tratamiento más efectivo, un cambio en el método de recopilación de información: todo esto puede cometer tantos errores en los datos reales que un pronóstico a largo plazo estará completamente lejos de la realidad.

Muletas y bicicletas: modificamos el modelo SEIR

Intentemos responder a la pregunta de por qué el crecimiento de la epidemia se está desacelerando a uno lineal. Con la cantidad de personas infectadas que tenemos ahora, las economías de escala desempeñan un papel importante asociado con la velocidad limitada de comunicación entre las personas.

Más precisamente, recordamos: el número de casos en el modelo SEIR es directamente proporcional al número promedio de casos en la población

$inline$ . Esta regla funciona bien en poblaciones pequeñas, donde todos pueden comunicarse con todos, y los pacientes están distribuidos de manera uniforme. En realidad, especialmente en la escala de decenas y cientos de miles de personas, si tomas a dos personas enfermas al azar, resulta que no solo nunca se comunicaron entre sí y no se vieron, ni siquiera viajaron en el mismo vagón del metro. Y de hecho viven en diferentes ciudades. Todo lo que los une es una cadena de lazos sociales, lo que llevó al hecho de que el virus se les transmitió.

Como ejemplo, construí un modelo epidémico en forma de autómata celular, donde cada célula interactúa con solo 4 vecinas. Esto es equivalente al hecho de que cada individuo en la población tiene 4 contactos sociales: este es un número muy pequeño para la población humana, pero cuanto más rápido se manifieste el efecto de limitar las conexiones sociales. En cada iteración con una probabilidad de 0.1, cada uno de los 4 vecinos de la célula infectada puede infectarse. La enfermedad dura un promedio de 14 días. Los resultados de la simulación para un grupo de celdas de 200x200 se presentan en la figura a continuación, donde

$inline$ Es el número de iteración.

El color azul indica susceptible, amarillo - enfermo, verde - recuperado. Lo más interesante es cómo se ven los gráficos de la cantidad de casos. Y se ven como planeados: después de una breve fase de crecimiento sub-exponencial, al igual que en el modelo SEIR, hay una fase prolongada de crecimiento lineal, como en la realidad.

Mi objetivo no era obtener una imagen que se pareciera cuantitativamente a la realidad. Si desea más credibilidad, puedo recomendar el proyecto de Sergey Potekhin, que se publicó recientemente en Habré . Para lectores meticulosos, la prueba más estricta de crecimiento lineal se da a continuación.

Prueba del teorema de la epidemia de crecimiento lineal en una gran población.

$inline$ - . . , 20 ( )

$d \approx 4$ .

$inline$ - .

$n^{\frac{1}{d}}$ , ,

$inline$ ,

$inline$ . , ,

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

$inline$ :

$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}.$

$inline$

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

4% 16- .

En cuanto a las estadísticas mundiales sobre la epidemia, veremos lo mismo: durante un mes, el crecimiento ha sido lineal, a pesar del hecho de que los modelos clásicos nos prometieron una continuación del aumento exponencial en el número de casos.

Tareas para curiosos

COVID-19, offline .
.1 ?
.1 2, .

Ahora sobre la modificación del modelo SEIR. Lo más simple que podemos hacer es multiplicar el componente no lineal por alguna función, dependiendo del número de pacientes. Con un pequeño número de casos, esta función debería estar cerca de 1, con un gran número, debería asintóticamente tender a cero. El candidato adecuado más simple es

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

Selección de parámetros.

$\alpha$ y

$inline$ puede compensar el crecimiento exponencial en el modelo original.

Agregar al modelo, para mayor información, y el componente

$inline$ - El número de muertes. Obtenga una versión modificada del modelo 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}$

Los resultados de la simulación se muestran en la figura a continuación.

El error estándar en comparación con el modelo original casi no ha cambiado. Los valores de los parámetros ya son realistas. Por conveniencia, designé el número inicial de infectados en la fase activa como

$inline$ .

$\ begin {align} & \ beta = 0.219, \\ & \ gamma = 0.0102, \\ & E_0 = 0.13 \ cdot I_0, \\ & \ kappa = 1/3, \\ & \ mu = 1 , 13 \ cdot 10 ^ {- 3}. \ end {align}$

El modelo también interpola muy bien los derivados: los valores del aumento diario en el número de casos y el número de muertes.

Intentemos hacer un pronóstico. Tomamos un horizonte de pronóstico de 2 meses y continuamos modelando soluciones con los parámetros encontrados por el programa de optimización.

A primera vista, no está mal, pero no puede desear tal pronóstico para su amada patria: el número de casos nuevos seguirá disminuyendo, pero el número total seguirá creciendo. En este caso, una epidemia se puede detener solo con la ayuda de una vacuna, o esperando hasta que casi toda la población haya estado enferma. El número de nuevas muertes se establece en aproximadamente 200 por día. Esta es una ilustración clara de lo que sucederá si no fortalecemos las medidas para combatir la epidemia. ¿Es eso lo que nos espera? Y en aras de este futuro no muy brillante, muchos de nosotros nos sentamos diligentemente en casa, ¿hemos comprado trigo sarraceno y papel higiénico?

A continuación consideraré dos escenarios, y mirando la distancia nublada de los próximos meses del 28 de abril de 2020, no puedo decir con certeza cuál desarrollará más los eventos. Ahora, en el momento de romper la curva de nuevos casos, estamos en un punto donde es doblemente problemático predecir algo.

Guión estadounidense

El hegemón mundial estaba en una posición extremadamente poco envidiable. Retrasado en la toma de decisiones clave que retrasarían el crecimiento de la epidemia al principio, todavía no puede hacer frente al aumento natural en los nuevos casos.

El modelo SEIRD modificado, entrenado en los primeros 33 puntos, a partir del 2 de marzo, más o menos pronostica de manera realista el curso de la epidemia en abril.

Como puede ver, el crecimiento en abril es casi perfectamente lineal. El modelo exagera ligeramente la mortalidad en abril, pero el panorama general es correcto.

Esta imagen muestra el aumento diario de nuevos casos y muertes en los Estados Unidos. Es muy similar a lo que el modelo predijo para Rusia.

Escenario de Alemania

Los alemanes disciplinados lograron girar la curva a su favor, y su crecimiento es más lento que el lineal. Además, para que el modelo sea relevante, tuve que agregar manualmente un aumento en la tasa de recuperación el 6 de abril

$\ gamma$ 1,7 veces, de lo contrario, no se puede explicar una caída tan pronunciada en el número de casos en términos del modelo SEIRD.

El modelo fue entrenado en los primeros 27 puntos, a partir del 10 de marzo. También cambié la función no lineal. Para Alemania, el exponente dependiente del tiempo es mejor:

$\ mathcal {\ varphi} (t) = K_0 e ^ {- \ alpha t}.$

Este tipo de función indica un aumento acumulativo en el número de lazos sociales interrumpidos y, en consecuencia, la propagación de la infección. Aquí tiene una ilustración clara de los beneficios del autoaislamiento.

Lo anterior muestra el aumento diario de nuevos casos y muertes. Como en el caso de los Estados Unidos, los datos reales contienen fluctuaciones pronunciadas con un período de 7 días. Esto significa que los fines de semana aumenta la cantidad de contactos y, en consecuencia, aumenta la cantidad de infectados. [UPD: por el contrario, está disminuyendo, ya que el índice de autoaislamiento Yandex puede testificar indirectamente . ]

Conclusión

Hacer predicciones, tanto a corto como a largo plazo, no es solo un tributo a la curiosidad. En caso de una epidemia, debe saber cuántas camas se deben preparar, cuántos ventiladores se deben producir, cuántos meses se deben abastecer de equipos de protección personal para los médicos. Los funcionarios deben entender si las medidas tomadas son suficientes o si se deben introducir nuevas prohibiciones y restricciones. Idealmente, el modelo debería reflejar la realidad tan bien que sería posible ver la fuerza de acción de cada medida recientemente adoptada, y luego sería posible fortalecer medidas útiles y cancelar decisiones que resultaron inútiles.

Aunque con algunas reservas, todos los que aún no han adquirido la inmunidad de COVID-19 realmente se pueden describir con la carta

$inline$ , toda la variedad de personas enfermas - la carta

$inline$ etc. Además, el modelo SEIRD puede incluso explicar algo. Pero ella puede predecir algo en el futuro distante de manera muy aproximada.

Deliberadamente cité en el artículo solo un escenario negativo cuando repetimos el destino de los Estados Unidos en términos de la dinámica de la epidemia. Si este escenario resulta ser cierto, para fines de junio tendremos más de 300 mil casos registrados de la enfermedad y más de 10 mil muertes. Aunque existen requisitos previos para el hecho de que este escenario no se traduzca en realidad, le aconsejaría que se relacione con él de acuerdo con el principio: "espere lo mejor, prepárese para lo peor". Como dice el dicho, si las mejores mentes de la NASA han emprendido la lucha contra la epidemia en los Estados Unidos , entonces esto es realmente algo malo.

Hasta ahora, todo lo que nos queda es visitar menos lugares públicos, usarcon respiradores adecuados , lávese las manos, no olvide limpiar el teléfono inteligente con alcohol después de sacarlo a la calle y siga otras recomendaciones simples.

Pero aún así, ¿es posible hacer predicciones más precisas? Sí, por supuesto. Pero sobre esto en otro momento.

Si desea jugar con el código fuente y sugerir un escenario de desarrollo usted mismo, aquí está el enlace al github .

La pandemia de COVID-19 a través de los ojos de un matemático, o por qué el clásico modelo SEIRD no funciona