Al inicio de una epidemia el aumento de casos es exponencial, una función que sorprende a la mayoría por su explosivo crecimiento, sin embargo esto solo es la fase inicial de un modelo de crecimiento denominado sigmoide.
El modelo de crecimiento sigmoide más conocido fue propuesto por Pierre François Verhulst entre 1838 y 1847 y por ello se le suele llamar función de Verhulst o función logística. Sin embargo para multitud de procesos de crecimiento biológico se suele utilizar el modelo propuesto en 1825 por Benjamin Gompertz para estudiar la relación entre la edad y el incremento de la mortalidad en humanos. Lo que él denominaba "el agotamiento promedio del poder de un hombre para evitar la muerte", o la "porción de su poder restante para oponerse a la destrucción". La función de Gompertz, rápidamente adoptada por la industria de los seguros para proyectar el riesgo de muerte, es ampliamente utilizada en biología como modelo de regresión para estudiar el crecimiento de las poblaciones de animales, bacterias, crecimiento de tumores y supervivencia de pacientes de cáncer, y por supuesto modelos de infección.
Fig 1. Función de Gompertz (azul) y función de Verhulst (logística, rojo) utilizando la capacidad de calculadora gráfica del buscador de Google introduciendo las fórmulas tal y como están escritas.
1. La función de Gompertz
La función de Gompertz es muy similar a la función logística, pero se diferencia de ella en que es asimétrica. El crecimiento es más rápido al principio de la curva que al final. Esta pequeña diferencia la hace más apropiada para describir el crecimiento biológico.
La función de Gompertz es la solución a la ecuación diferencial que describe los cambios en la población (P) con el paso del tiempo (t) en función de su capacidad de crecimiento intrínseca (c, constante) y la máxima población que el ecosistema puede soportar, lo que se define como capacidad de carga (K, carrying capacity). Podéis consultar mi artículo sobre "El problema de la población humana", donde este concepto se explica con gran detalle. La ecuación diferencial que describe el comportamiento de las poblaciones es:
dP/dt = c ln(K/P)P
Esta ecuación cuando se integra y resuelve nos da la función de Gompertz, que describe los cambios en la población en función del tiempo:
P(t) = K e^-ln(K/P0) e^-ct
Donde P0 es el tamaño de la población inicial. Esta ecuación la podéis ver escrita de varias formas equivalentes. En Excel se puede introducir como:
=K*EXP((-LN(K/P0))*EXP((-C)*XX))
Donde K es el número que refleja el máximo de casos que se van a alcanzar, la asíntota superior.
P0 es el número de casos iniciales.
C es el número que corresponde a la constante de crecimiento.
XX es el primer número de celda de donde debe tomar el dato de tiempo (número de días desde el inicio).
Una de las propiedades de la función de Gompertz es que su segunda derivada se vuelve cero cuando la población alcanza el valor de la capacidad de carga dividido por e, P = K/e, (e = 2,718). Esto quiere decir que antes de alcanzar la mitad del total de casos que va a haber, el número de nuevos casos alcanza su máximo y comienza a descender. Como he indicado el descenso en el número de casos es más lento que el ascenso (figura 2). Esta es una de las principales diferencias con la función de Verhulst (logística) que es simétrica.
Fig. 2. Características de la función de Gompertz para una capacidad de carga K de 200.000 y una constante de crecimiento c de 0,06, con una P0 de 2 en el día 1. El pico de nuevos casos se alcanza en el día 42 con 4.400 casos al día.
2. La función de Gompertz para la pandemia de coronavirus
El crecimiento de una enfermedad infecciosa se puede ajustar a la función de Gompertz. Para ello hacen falta dos datos, uno es la capacidad de carga, que corresponde al máximo número de casos de infección que se van a alcanzar, el numero de población que habrá sido infectada al final del episodio epidémico. Este dato puede estimarse o calcularse porque cuando el número de nuevos casos alcanza su máximo basta multiplicar el número total de infectados en ese punto por el número e: K = Pe. El segundo dato necesario es la constante de crecimiento, que es la que determina la velocidad de crecimiento en el eje de las ordenadas (Y), o lo que es lo mismo la inclinación de la curva. Este valor se puede hallar mediante ajuste a los datos.
Veamos un ejemplo. Tomemos los datos de muertes por coronavirus COVID19 en China, que son más homogéneos que los de casos por los cambios que hubo en los criterios de diagnóstico, de la página de Worldometer de China. El análisis de las nuevas muertes cada día muestra un perfil un tanto irregular, pero se puede decir que el máximo de nuevas muertes se alcanzó el 11 de febrero con 146 muertos. Aunque hay dos días que tuvieron más muertos (6 y 22 de febrero), el perfil muestra claramente que se trataron de desviaciones puntuales, no de la tendencia principal, con la curva comenzando el 22 de enero con 17 muertos. El valor de K se puede estimar a partir del total de muertos el 11 de febrero, 1.259 multiplicándolo por el número e (2,718) lo que da una K de 3.422. El valor de c se ajusta para que la curva de la función coincida con la de los datos.
Tomando un valor de K de 3450 y un valor de c de 0,088 la función de Gompertz en Excel toma la forma:
=3450*EXP((-LN(3450/17))*EXP((-0.088)*B3))
Donde la columna B tiene el número de días desde el comienzo con B3=1 correspondiendo al 23 de enero partiendo de un valor inicial P0=17 el día anterior.
La representación de la función comparada con los datos reales es:
Fig. 3. Muertes en China ajustadas a la función de Gompertz.
Es importante destacar que las medidas que se toman para frenar al coronavirus afectan a los valores de K y de c, y por lo tanto a la curva final de Gompertz a la que se ajusta la epidemia. En el caso de las muertes en China, las medidas tomadas redujeron el valor de K para las muertes al entorno de las 3.500, y redujeron el valor de c al disminuir la velocidad de propagación.
3. Las matemáticas de la epidemia en España
La gente se pregunta cuanto va a durar la necesidad de mantener las actuales medidas (o incluso medidas más estrictas), hasta que remita la epidemia. Es algo a lo que todavía no se puede responder. Estimar la K, el número de casos que va a haber, sin haber llegado al máximo de los nuevos casos está lleno de incertidumbre. No ayuda que desde el 11 de marzo haya cambiado el criterio para hacer los tests y ahora ya no se hagan a quienes presentan síntomas leves. Según lo observado en China, el pico de nuevas muertes tiene lugar unos 12-13 días después del pico de nuevos casos, y sería deseable no tener que esperar tanto.
Se puede hacer un ajuste de la función de Gompertz al número de casos de coronavirus en España, pero está sujeto a un error muy grande, porque hay un gran número de combinaciones de K y c que en estos momentos pueden dar un ajuste razonablemente bueno. Sin embargo sirve para hacernos una idea de las cifras que podemos manejar.
Por ejemplo, utilizando los datos de casos detectados en España hasta ayer 15 de marzo, podemos estimar K en 200.000 casos y c en 0,06 lo que produce un ajuste razonablemente bueno. Bajo este supuesto los casos seguirían aumentando hasta llegar a los 200.000, con una cifra de muertos de en torno a 4.000-5.000. El número de casos nuevos alcanzaría su máximo el 4 de abril, y el número de nuevos muertos en torno al 16-17 de abril. Para finales de junio habría menos de 100 nuevos casos diarios. Podemos estimar que en torno a 1-4 millones de personas habrían pasado la enfermedad, por lo que más del 90% seguiría siendo susceptible. Llegados a este punto supongo que el gobierno relajaría las medidas y pasaría a hacer un control mucho más agresivo sobre los focos y nuevos casos, porque es lo que yo haría.
Fig.4. Posible ajuste de los casos de coronavirus en España a una entre muchas funciones de Gompertz, que en este caso presupone una incidencia máxima K de 200.000 casos y una constante de crecimiento c de 0,06. Es altamente improbable que valores estimados tan pronto en la epidemia sean correctos, pero sirven para hacerse una idea de la progresión de la epidemia. Estos supuestos son los mismos que los de la figura 2.
Este escenario no constituye una predicción, ni tan siquiera una proyección. Los datos son insuficientes para tener la más mínima fiabilidad, pero cuando se aproxime el pico de nuevos casos, y sobre todo cuando lo confirme el pico de nuevas muertes tendremos una idea razonable de la evolución del coronavirus en España. Este artículo no trata de predecir el futuro, sino de mostrar como funciona una epidemia, algo que los responsables de nuestra seguridad no parecen haber entendido con la necesaria prontitud, y trata de dar una idea de los tiempos que se manejan. Por ello cuando el gobierno ha dicho hoy que será necesario prolongar el estado de alarma más allá de los 15 días establecidos por la ley es porque ya sabe que incluso aunque todo vaya bien y el pico de nuevos casos tuviera lugar a principios de abril, esto tiene cuerda al menos hasta junio.