Autor: CAChemE.org, original por Charles Hill, No painers blog - Licencia CC-BY
La propagación del virus del ébola en el África occidental sigue, por el momento, imparable. Los modelos epidemiológicos tratan de representar este comportamiento mediante fórmulas matématicas. Este documento intenta ejemplificar con un modelo muy simple, y por tanto aproximado, del impacto de variables como la facilidad de contagio y mortalidad del mismo.
Nuestra formación no es en el ámbito de la medicina ni estadística y los resultados aquí mostrados son un ejercicio aplicado de métodos numéricos de integración con el lenguaje de programación Python. Por ello, cualquier mejora en el modelo es bienvenida :)
Los datos que podemos consultar del Ébola en el África occidental son desesperanzadores:
Fuente: Charles Hill, No painers blog - Datos más actualizados (a octubre de 2014) se pueden consultar en la BBC y Forbes.
Este gráfico representa la cantidad total (acumulada) de los casos y muertes producidos por el virus en el África occidental. El eje vertical es logarítmico y, por tanto, la tendencia lineal que aparece a partir de junio indica que los números están creciendo exponencialmente.
Ateniendo al gráfico, el número de casos y defunciones por Ébola se duplica cada mes.
Siendo más precisos, el nº de casos se se duplica cada 29 días. Si se mantuviera esa velocidad de propagación, en nueve meses tendríamos un millón de casos, y en dos años una fracción significativa de la población mudnial habrá tenido la enfermedad. Obviamente necesitamos un modelo mejor ya que el crecimiento exponencial no puede ser ilimitado y esto es sólo un primer intento.
Existen cuatro variables / categorías de personas en este modelo:
$S$ son personas susceptibles que no han contraído el virus todavía. Si entran en contacto con alguien que es infeccioso, podrían contraer el virus.
$E$ son personas expuestas que tienen el virus , pero que aún no presentan los síntomas y no son infecciosos.
$I$ son personas que han sido contagiados por Ébola, han pasado el tiempo de incubación y son infecciosas. Su cuerpo se encuentra en plena batalla contra el virus y sus síntomas.
$R$ son personas que han ganado la batalla y se han recuperado. Ellos ahora son inmunes al virus.
Si te fijas, el número total de personas aún con vida $N$ corresponde con la suma de todos ellos.
Así el modelo epidemiológico a seguir se llama SEIR y puede ser definido por el siguiente sistema de ecuaciones diferenciales:
\begin{align} \frac{dS} {dt} &= -\beta \frac{S}{N} I \\\\ \frac{dE}{dt} &= +\beta \frac{S}{N} I -\sigma E \\\\ \frac{dI}{dt} &= +\sigma E - \gamma I \\\\ \frac{dR}{dt} &= (1-f) \gamma I \end{align}Para nuestro caso, consideraremos que hay una persona afectada únicamente ($E = 1$)
En esta sección se da una explicación del modelo SEIR y los parámetros que se han elegido. Si quieres, puedes ir directamente a los resultados:
La disminución del número de personas sanas con el tiempo será proporcional a las personas infecciosas (I) y la cantidad (proporcional) de las mismas $S/N$. Además, definimos una factor de probabilidad por día, que llamaremos $\beta$, de infectar a otra persona. La variación del nº personas expuestas ($E$), seguirá la misma relación pero aumentando su número (signo positivo). Esto, matemáticamente queda expresado como:
\begin{align} \frac{dS} {dt} &= -\beta \frac{S}{N} I \\\\ \frac{dE}{dt} &= +\beta \frac{S}{N} I \\\\ \end{align}Una vez que una persona está infectada ($E$), es necesario un tiempo para mostrar los síntomas y llegar a ser infecciosa ($I$). Este tiempo depende de cada caso, pero en los brotes anteriores de Ebola, se estimó que se necesita alrededor de 9 días para ello. Vamos a modelar esto con una velocidad de pasar de $E$ a $I$ de $\sigma = 1/9$ por día:
$$ \begin{align} \frac{dE}{dt} = -\sigma E \\\\ \frac{dI}{dt} = +\sigma E \end{align} $$En este punto, la persona enferma puede recuperarse o fallecer. Una vez más, vamos a modelar esto con una velocidad, que se puede encontrar con los brotes anteriores, de $ \gamma = 1 / 8.5$ por día. La probabilidad de fallecer estimada por la OMS en estos momentos es del 50% (Fuente: OMS). Sin embargo, vamos a ser optimistas y considerar que esta tasa puede reducirse hasta un 40% ($f = 0.4$) si los recursos de los centros hospitalarios no se ven comprometidos y los síntomas del virus se tratan a tiempo. Así, la fracción de personas sobreviven es de $1-f = 0.6$ (personas que sobreviven).
$$ \begin{align} \frac{dI}{dt} &= - \gamma I \\\\ \frac{dR}{dt} &= (1-f) \gamma I \end{align} $$Aún no se ha definido el valor de $\beta$, el cual indica la velocidad de personas que se infectan. Para ello haremos uso de un parámetro importante: el número medio de personas que un paciente con ébola podría llegar a infectar $R_0$.
$$R_0 = \frac{\beta}{\gamma}$$En brotes anteriores, estimaron $R_0$ estar entre 1.3 y 2.7, es decir, un paciente con ébola llega a infectar de 1 a 3 personas. Este parámetro difiere entre países y factores como:
Para ver la importancia que tiene disminuir este factor $R_0$, vamos a interactuar con el modelo y considerar que a mayor nº de recortes en sanidad mayor es el nº de contagios que un paciente con ébola puede llegar a realizar antes de estar controlado, así $ R_0 = recortes \: sanidad $.
El ártículo original en el que se basa este trabajo resolvió el sistema de ecuaciones diferenciale numéricamente con parámetros similares obteniendo una aproximación entre su modelo y los datos reales.
Como se puede ver, el modelo no coincide bien al principio (antes de mayo), pero una vez que la fase de crecimiento exponencial comienza, la simulación obtiene resultados similares. Definitivamente, el modelo captura la velocidad de aumento, la relación entre los casos y los fallecimientos.
Entonces, ¿qué impacto tienen los recortes en la velocidad de propagación del ébola? Veámoslo en este gráfico resuelto y obtenido con Python.
En España seguimos a la espera de comprobar cuantas personas llegaron a contagiarse, por lo que hemos variado el valor de $R_0 = recortes \: sanidad$ entre 1 y 2. A la vista de los resultados, no es de extrañar que expertos y personal sanitario estén pidiendo más medios y formación para evitar contagios durante el tratamiento de pacientes con ébola en hospitales.
Para interactuar con este documento de IPython Notebook/JuPyter, puedes descargar e instalar Anaconda de Continuum (ver instalación) u optar por distribuciones portables de Python como Pyzo.
In [4]:
import numpy as np
%matplotlib inline
import matplotlib.pyplot as plt
from scipy import integrate
from IPython.html.widgets import interact, interactive, fixed
from IPython.html import widgets
# Hasta ahora, hemos cargado nuestras bibliotecas favoritas de Python
# Pasamos a definir nuestro sistema de ec. diferenciales:
def DiffEqs(y, t, sigma, gamma, f, beta):
""" Modelo simple epidemiológico para la
propagacion del virus del Ébola"""
# Recuperamos las variables dependientes de entrada
S, E, I, R, D = y
Nn = S+E+I+R # Total de personas vivas
# Definimos las ecuaciones diferenciales
dSdt = -beta*(S/Nn)*I
dEdt = +beta*(S/Nn)*I - sigma*E
dIdt = sigma*E - gamma*I
dRdt = (1.0-f)*gamma*I
dDdt = f*gamma*I
return [dSdt, dEdt, dIdt, dRdt, dDdt]
def ebola(tiempo, recortes_sanidad):
""" Solve the model and plot a basic graph of the results """
tiempo_simulacion = tiempo # días
t = np.linspace(0., tiempo_simulacion, 1000)
#------------------------------------------------------------------------------
# Constantes del modelo
#------------------------------------------------------------------------------
sigma = 1./9 # Velocidad convertirse en contagioso una vez expuesto
gamma = 1./8.5 # Velocidad de recuperación o fallecimiento
f = 0.4 # Fracción de gente que no supera la enfermedad (40%)
# A mayor nº de recortes, mayor nº de personas contagiadas por infectado
R_0 = recortes_sanidad
beta = R_0*gamma # Velocidad con la que se propaga la infección
# Condiciones iniciales
N = 350 # Población total
y0 = [N, 1, 0, 0, 0] # [S, E, I, R, D]
# Llamamos al ode para resolver las eq. diferenciales
solution = integrate.odeint(DiffEqs, y0, t, args=(sigma, gamma, f, beta))
S = solution[:,0] # S es el nº de personas sin infectar
E = solution[:,1] # E es el nº de portadores del virus no contagiosos (aún)
I = solution[:,2] # I es el nº de personas sintomáticas e infecciosas
R = solution[:,3] # R es el nº de personas recuperadas
Ds = solution[:,4] # D es el nº de personas fallecidas
# Creamos la figura para representar los resultados
fig = plt.figure(figsize=(10, 5))
# Representamos las curvas de evolución con respecto al tiempo
plt.plot(t, (E+I)/N*100 )
plt.plot(t, I/N*100)
plt.plot(t, R/N*100)
plt.plot(t, Ds/N*100)
plt.plot(t, S/N*100)
# Añadimos leyenda a los datos
plt.legend(["% infectados",
"% infectados contagiosos",
"% recuperados",
"% fallecidos",
u"% población sana "], loc=1)
plt.suptitle(u'Variación de la velocidad de propagación del virus del Ébola:', fontsize=14, fontweight='bold')
plt.title(u"Impacto de los recortes en el sistema sanitario.", fontsize = 12)
plt.text(10, 50, u"A mayor nº de recortes, mayor \n nº de contagios (accidentales) \n y mayor velocidad de propagación.", style='italic',
bbox={'facecolor':'red', 'alpha':0.5, 'pad':10})
plt.xlabel(u"tiempo / días")
plt.ylabel(u"porcentaje de población")
plt.show()
In [3]:
interact(ebola, tiempo=(100,600), recortes_sanidad=(1,2,0.01))
In [ ]: