Simulation and cases
By: Fabio Pinto
(Susceptible, exposed, infected, recovered)
One of the most widely used epidemiological models is the so-called SIR model, which was proposed by W. O. Kermack and A. G. McKendrick in 1927.
In a population of fixed size N in which an epidemic has been unleashed that spreads by contagion, in a time t individuals can be in three different states:
By: Fabio Pinto
(Susceptible, exposed, infected, recovered)
One of the most widely used epidemiological models is the so-called SIR model, which was proposed by W. O. Kermack and A. G. McKendrick in 1927.
In a population of fixed size N in which an epidemic has been unleashed that spreads by contagion, in a time t individuals can be in three different states:
- susceptible S(t),
- infected I(t),
- recovered R(t)
Wang, and others (2020) have opted for an evolved CRS model, SEIR , which gave acceptable results in the case of Ebola in 2004 (Chowell, Hengartnet, Castillo-Chavez, Fenimore, & Hyman, 2004) in which the population infected but not infecting others (Exposed) is included, where now N=S+E+I+R.
Referencias:
- Qianying Lin, Shi Zhao, Daozhou Gao, Yijun Lou, Shu Yang, Salihu S. Musa, Maggie H. Wang, Yongli Cai, Weiming Wang, Lin Yang, Daihai He, A conceptual model for the coronavirus disease 2019 (COVID-19) outbreak in Wuhan, China with individual reaction and governmental action, International Journal of Infectious Diseases 93 (2020), 211–216. Publicado: March 04, 2020. https://www.ijidonline.com/article/S1201-9712(20)30117-X/fulltext
- Gutierrez and Varona, 2020. Universidad De La Rioja. [online] Análisis de la posible evolución de la epidemia de coronavirus COVID-19 por medio de un modelo SEIR. Publicado: March 15, 2020. https://www.unirioja.es/apnoticias/servlet/Archivo?C_BINARIO=12051 [Accessed 13 April 2020].
In [23]:
from IPython.display import Image Image(filename='data/sir.png')
Out[23]:
In [24]:
Image(filename='data/seir.png')
Out[24]:
In [25]:
from IPython.display import Image Image(filename='data/edosir.png')
Out[25]:
In [26]:
Image(filename='data/edoseir.png')
Out[26]:
- Beta: Tasa de infeccion
- Gamma: Tasa de recuperación
- Ommega: Tasa de incubación
In [27]:
import scipy.integrate as spi import numpy as np import matplotlib.pyplot as plt import pandas as pd
Resolución de las ecuaciones diferenciales
In [28]:
## Definimos las ecuaciones diferenciales del modelo
def odeSIR(y, t, N, beta, gamma):
S, I, R = y
dSdt = -beta * S * I / N
dIdt = beta * S * I / N - gamma * I
dRdt = gamma * I
return dSdt, dIdt, dRdt
In [29]:
def odeSEIR(y, t, N, beta, gamma, omega):
S, E, I, R = y
dSdt = -beta * S * I / N
dEdt = beta * S * I / N - omega*E
dIdt = omega*E - gamma * I
dRdt = gamma * I
return dSdt, dEdt, dIdt, dRdt
Modelado y simulación del modelo SIR
In [30]:
def modeloSIR(N, beta, gamma, days):
I0, R0 = 1, 0
## S susceptibles
S0 = N - I0 - R0
## beta = ratio de contacto entre personas
## gamma = inverso del tiempo de recuperación
# Vector de tiempo
t = range(0,days)
# Vector de condiciones iniciales
y0 = S0, I0, R0
# Evaluación del sistema de ecuaciones diferenciales
ret = spi.odeint(odeSIR, y0, t, args=(N, beta, gamma))
S, I, R = ret.T
return S, I, RModelado y simulación del modelo SEIR
In [31]:
def modeloSEIR(N, beta, gamma, omega, days):
"""N = población total """
I0, R0, E0 = 1, 0, 0
## S susceptibles
S0 = N - I0 - R0- E0
## beta = ratio de contacto entre personas
## gamma = inverso del tiempo de recuperación
## omega = tasa de incubación
# Vector de tiempo
t = range(0,days)
# Vector de condiciones iniciales
y0 = S0, E0,I0, R0
# Evaluación del sistema de ecuaciones diferenciales
ret = spi.odeint(odeSEIR, y0, t, args=(N, beta, gamma, omega))
S, E, I, R = ret.T
return S, E, I, R
Función para hacer un breve resumen de los datos
In [32]:
def informeSIR(d):
mean = int(np.mean(d))
maxx = int(np.max(d))
max_day = d.argmax()
print('Promedio Diario = ' + str(mean))
print('Máxima tasa = ' + str(maxx))
print('Dia del pico = ' + str(max_day))
'Aplanar la curva': Modelo SIR con diferentes medidas de contención
La presencia o no de medidas de contención (aislamiento, cuarentena, etc) determinan la variación de la tasa de contagio
Para el ejemplo estamos tomando un beta de 2.5 inicial para ausencia de medidas de contención y 0.5 para medidas de contención serias.
In [54]:
S_1, I_1, R_1 = modeloSIR(1e5,2.5,1/5,100) # sin ninguna medida de contención S_2, I_2, R_2 = modeloSIR(1e5,1.5,1/5,100) # con medidas de contención leves S_3, I_3, R_3 = modeloSIR(1e5,0.5,1/5,100) # con medidas de contención fuertes
In [77]:
t = range(0,len(S_1))
plt.subplots(figsize=(10, 10))
plt.subplot(311)
plt.plot(t, S_1, 'g', label = 'Susceptibles')
plt.plot(t, I_1, 'r', label = 'Infectados')
plt.plot(t, R_1, 'b', label = 'Recuperados')
plt.ylabel('Conteo')
plt.title('Modelo SIR con beta = 2.5')
plt.legend()
plt.subplot(312)
plt.plot(t, S_2, 'g', label = 'Susceptibles')
plt.plot(t, I_2, 'r', label = 'Infectados')
plt.plot(t, R_2, 'b', label = 'Recuperados')
plt.ylabel('Conteo')
plt.title('Modelo SIR con beta = 1.5')
plt.legend()
plt.subplot(313)
plt.plot(t, S_3, 'g', label = 'Susceptibles')
plt.plot(t, I_3, 'r', label = 'Infectados')
plt.plot(t, R_3, 'b', label = 'Recuperados')
plt.xlabel('Dias')
plt.ylabel('Conteo')
plt.title('Modelo SIR con beta = 0.5')
plt.legend()
Out[77]:
<matplotlib.legend.Legend at 0x108d9f71ac8>
Comparación de curvas de infectados con diferentes medidas de contención
In [39]:
t = range(0,len(I_1)) plt.plot(t, I_1, 'g', label = 'Infectados sin medidas') plt.plot(t, I_2, 'r', label = 'Infectados con medidas leves') plt.plot(t, I_3, 'b', label = 'Infectados con medidas serias') plt.legend() plt.show()
In [40]:
print("---------Infectados sin medidas-----------")
informeSIR(I_1)
print("---------Infectados con medidas leves-----------")
informeSIR(I_2)
print("---------Infectados con medidas serias-----------")
informeSIR(I_3)
---------Infectados sin medidas-----------
Promedio Diario = 5000
Máxima tasa = 71489
Dia del pico = 6
---------Infectados con medidas leves-----------
Promedio Diario = 4997
Máxima tasa = 58608
Dia del pico = 10
---------Infectados con medidas serias-----------
Promedio Diario = 4462
Máxima tasa = 23330
Dia del pico = 40
Comparación de curvas de recuperados con diferentes medidas de contención
In [41]:
t = range(0,len(I_1)) plt.plot(t, R_1, 'g', label = 'Recuperados sin medidas') plt.plot(t, R_2, 'r', label = 'Recuperados con medidas leves') plt.plot(t, R_3, 'b', label = 'Recuperados con medidas serias') plt.legend() plt.show()
Comparación de curvas de susceptibles con diferentes medidas de contención
In [42]:
t = range(0,len(S_1)) plt.plot(t, S_1, 'g', label = 'Susceptibles sin medidas') plt.plot(t, S_2, 'r', label = 'Susceptibles con medidas leves') plt.plot(t, S_3, 'b', label = 'Susceptibles con medidas serias') plt.legend() plt.show()
'Aplanar la curva': Modelo SEIR con diferentes medidas de contención
In [73]:
S1, E1, I1, R1 = modeloSEIR(1e5,2.5,1/5,1/7,200) # sin ninguna medida de contención S2, E2, I2, R2 = modeloSEIR(1e5,1.5,1/5,1/7,200) # con medidas de contención leves S3, E3, I3, R3 = modeloSEIR(1e5,0.5,1/5,1/7,200) # con medidas de contención fuertes
In [75]:
t = range(0,len(S1))
plt.subplots(figsize=(10, 10))
plt.subplot(311)
plt.plot(t, S1, 'g', label = 'Susceptibles')
plt.plot(t, E1, 'c', label = 'Expuestos')
plt.plot(t, I1, 'r', label = 'Infectados')
plt.plot(t, R1, 'b', label = 'Recuperados')
plt.ylabel('Conteo')
plt.title('Modelo SIR con beta = 2.5')
plt.legend()
plt.subplot(312)
plt.plot(t, S2, 'g', label = 'Susceptibles')
plt.plot(t, E2, 'c', label = 'Expuestos')
plt.plot(t, I2, 'r', label = 'Infectados')
plt.plot(t, R2, 'b', label = 'Recuperados')
plt.ylabel('Conteo')
plt.title('Modelo SIR con beta = 1.5')
plt.legend()
plt.subplot(313)
plt.plot(t, S3, 'g', label = 'Susceptibles')
plt.plot(t, E3, 'c', label = 'Expuestos')
plt.plot(t, I3, 'r', label = 'Infectados')
plt.plot(t, R3, 'b', label = 'Recuperados')
plt.xlabel('Dias')
plt.ylabel('Conteo')
plt.title('Modelo SIR con beta = 0.5')
plt.legend()
Out[75]:
<matplotlib.legend.Legend at 0x108d842ec48>
In [47]:
print("---------Infectados sin medidas-----------")
informeSIR(I1)
print("---------Infectados con medidas leves-----------")
informeSIR(I2)
print("---------Infectados con medidas serias-----------")
informeSIR(I3)
---------Infectados sin medidas----------- Promedio Diario = 2499 Máxima tasa = 26478 Dia del pico = 31 ---------Infectados con medidas leves----------- Promedio Diario = 2498 Máxima tasa = 23033 Dia del pico = 44 ---------Infectados con medidas serias----------- Promedio Diario = 2227 Máxima tasa = 9578 Dia del pico = 119
In [48]:
print("---------Expuestos sin medidas-----------")
informeSIR(E1)
print("---------Expuestos con medidas leves-----------")
informeSIR(E2)
print("---------Expuestos con medidas serias-----------")
informeSIR(E3)
---------Expuestos sin medidas----------- Promedio Diario = 3499 Máxima tasa = 49266 Dia del pico = 28 ---------Expuestos con medidas leves----------- Promedio Diario = 3498 Máxima tasa = 39267 Dia del pico = 40 ---------Expuestos con medidas serias----------- Promedio Diario = 3121 Máxima tasa = 14002 Dia del pico = 114
In [49]:
t = range(0,len(I1)) plt.plot(t, I1, 'g', label = 'Infectados sin medidas') plt.plot(t, I2, 'r', label = 'Infectados con medidas leves') plt.plot(t, I3, 'b', label = 'Infectados con medidas serias') plt.legend() plt.show()
In [50]:
t = range(0,len(S1)) plt.plot(t, S1, 'g', label = 'Susceptibles sin medidas') plt.plot(t, S2, 'r', label = 'Susceptibles con medidas leves') plt.plot(t, S3, 'b', label = 'Susceptibles con medidas serias') plt.legend() plt.show()
