Los modelos pueden ser realmente útiles a la hora hacer predicciones, pero hay que ser consciente de sus limitaciones a fin de saber en qué situaciones se pueden utilizar.
La primera pregunta que nos deberíamos hacer en nuestro caso es qué modelo entre el SIR y el SEIR es el más adecuado. Se elegirá uno u otro en función de si el periodo de incubación 1/γ es lo suficientemente corto como para poder despreciarlo (SIR) o si, por el contrario, su larga duración afecta a la dinámica de la epidemia (SEIR).
El siguiente punto a considerar es si se resuelve el caso mediante el modelo determi-nista o uno de los estocásticos. Para ello, es necesario plantearse si la dinámica de cada individuo puede tener una importante influencia en la evolución del brote (casi siempre es así) o si se puede despreciar, lo que depende del tamaño de la población N. Para poblaciones pequeñas, la aleatoriedad gana importancia y, por tanto, el modelo determi-nista (poblaciones infinitas) deja de ser válido. En la Fig. 11 se muestra cómo aumenta la incertidumbre en el caso de una población con pocos individuos.
11.1:N = 103 11.2:N = 106
Figura 11: Comparación entre dos simulaciones del SIR hechas con el modelo de las Binomiales con el fin de ver el aumento de las fluctuaciones en una población con una N pequeña. Parámetros utilizados: µ= 0,16; β = 2µ.
Por último, si se decide utilizar el modelo estocástico a fin de tener en cuenta las posibles fluctuaciones en la solución, se ha de optar entre el Gillespie o el de las Binomiales.
El segundo es una aproximación del primero y se ha mostrado en la Sección 4 las posibles diferencias con respecto a la solución teórica a las que puede dar lugar. Sin embargo, la enorme carga computacional que supone utilizar el Gillespie con una población grande (los diagramas de fase requerían del orden de 8 horas de ejecución, mientras que los del modelo de las Binomiales menos de 5 min) vuelve realmente complicado llevar a cabo determinados análisis. Por tanto, puede valer la pena sacrificar cierta precisión a fin de ganar eficacia computacional. No obstante, el modelo de las Binomiales deja de ser válido si las probabilidades Pβ, Pµ o Pγ son demasiado grandes, ya que podría darse el caso de que un individuo se infectara y se recuperara dentro del intervalo de tiempo discreto (el modelo sólo considera un evento por individuo en cada intervalo de tiempo).
6. Conclusión
En este trabajo, motivado por la pandemia del COVID-19, se ha buscado explicar de manera detallada los modelos epidemiológicos SIR y SEIR sin demografía, así como también su resolución tanto determinista como estocástica. Se ha planteado la situación de un confinamiento parcial modificando los modelos originales y se han visto los efectos que tiene sobre la evolución de la epidemia. El lector podrá reproducir cualquier apartado si así lo desea.
Se han comparado los tres casos (determinista, Gillespie y Binomiales) destacando las posibles diferencias que podían aparecer entre cada uno de ellos. Hemos visto que la mediana de las soluciones del Gillespie da lugar a la solución determinista, mientras que su aproximación (el de las Binomiales) retrasa el pico de infectados y aumenta la duración del brote, aunque predice correctamente el valor del máximo de incidencia y la cantidad de recuperados final.
Se ha analizado la aparición del régimen crítico en los tres modelos dependiendo de los parámetros del confinamiento Xs y Pth que eligiéramos en un caso supercrítico con R0 = 2 para una población de un millón de habitantes, con el resultado de que tanto en el SIR como en el SEIR la criticidad aparece en Xs = 0,5 para los tres casos, comprobando así la validez de la Ec. (2.26). Por otro lado, la gran similitud entre los diagramas de fase de los dos modelos estocásticos nos indica que el de las Binomiales puede ser una buena alternativa al Gillespie (teniendo en cuenta sus limitaciones), debido a que el trabajo computacional que requiere éste último puede dificultar estudios en los que se requieran un gran número de simulaciones con un tamaño de población grande.
Destacar que es importante el estudio de este régimen, ya que al caracterizarse por una evolución constante de la curva de infectados durante un largo periodo, significa que el relajamiento de las medidas de contención antes de tiempo darían lugar a otro brote.
En la situación ideal planteada, se podría pensar que un confinamiento parcial del 50 % de la población sería una buena opción para acabar con la epidemia: no se confina a todos los individuos consiguiendo no paralizar del todo la economía y, al reducir drásticamente el número de susceptibles, se logra erradicar al patógeno. Sin embargo, como hemos visto en la Fig. 10, eso sería un error, pues daría lugar al régimen crítico alargando la duración de la epidemia y provocando un rebrote en cuanto las autoridades decidieran desconfinar a la población.
El presente trabajo ha tratado una situación idealizada a la que se le puede añadir complejidad con el fin de ser más realista. Por tanto, puede servir como punto de par-tida para futuros estudios de la criticidad en los que, por ejemplo, se tenga en cuenta la demografía (nacimientos y muertes) o se considere un sistema de varias poblaciones
interconectadas. Por otro lado, también se puede introducir el efecto de la vacunación o utilizar modelos con otros compartimentos (SIS, SEIS, MSEIR...).
7. Apéndices
7.1. Cálculo de R
0mediante la Next Generation Matrix
El número reproductivo básico de un modelo epidemiológico se obtiene utilizando la llamada Next Generation Matrix [9]. En este apartado explicaremos su significado aplicándola en el modelo SEIR (fuentes que se han consultado [10],[11] y [12]).
Sean las ecuaciones de expuestos (2.28) e infectados (2.29) linealizadas en torno al punto de absorción (S ≈1,I << 1,E <<1):
dE
dt ≈βI−γE (7.1)
dI
dt =γE−µI (7.2)
Si definimos el vector X~ = (E, I)0, el sistema de ecuaciones anterior se puede escribir de la siguiente manera: DondeT es la parte de transmisión y S la de transición.
Si resolvemos la Ec. (7.3), nos queda:
X~ =X~0exphT S−1−1tSi (7.4) T S−1es laNext Generation Matrix yR0se define como su autovalor dominante, ya que éste será el que determinará la dinámica de las curvas de los distintos compartimentos de la clase ”portadores del patógeno” (expuestos e infectados en nuestro caso). Así, el autovalor más grande de la matriz T S−1 será:
R0 = β
µ (7.5)
Quedando demostrado que el SIR y el SEIR tienen el mismo número reproductivo básico.
Para el SIR el cálculo deR0 es trivial, véase [11].
En el caso del modelo de las Binomiales, estamos discretizando el sistema de ecuaciones en intervalos en los que las transiciones ya no vienen determinadas por las tasas, sino por las probabilidades:
Por ello, elR0 en este modelo viene dado por:
R0 = Pβ
Pµ (7.8)
7.2. Método de ”binning” aplicado a los datos del Gillespie
El método de ”binning” permite agrupar datos continuos en ”bins” donde cada ”bin”
corresponde a un tiempo discreto.
A modo de ejemplo se mostrará cómo se ha aplicado en los datos del SIR:
1. Primeramente, modificaremos el algoritmo de la sección 3.1.1 añadiendo el cálculo de sbin, ibin y rbin ( inicialmente sbin =s0, ibin = i0 y rbin = r0). En cada iteración se han de redefinir e ir guardando los nuevos valores:
sbin →sbin+s ibin →ibin+i rbin →rbin+r
2. Al acabar la simulación, guardamos la duración máxima del brote tmáx, la lista de tiempo continuo y las de los valores de sbin, ibin y rbin .
3. Definimos los ”bins” como una lista de tiempos discretos: [0,1,2, ..., tmáx].
4. El dato del ”bin” n para una variable concreta ( sbin, ibin o rbin) será la suma de todos los valores de la variable en cuestión, cuyos tiempos continuos estén entre el
”bin”ny eln+1, dividida entre el número de estos valores. Por ejemplo, imaginemos que tenemos unos valores de sbin = [s1, s2..., sf] cuyos tiempos continuos son t = [0,5; 0,6;...; 0,9]. Como estos están entre 0 y 1, corresponden al ”bin” (0). Entonces los susceptibles en t= 0 serán:
s(0) =
Pf j=1sj
f (7.9)
De esta manera se consiguen datos de s, i y r discretizados con los que es posible calcular la mediana entre distintas simulaciones del Gillespie.
Referencias
[1] M. J. Keeling y P. Rohani, Modeling Infectious Diseases in Humans and Animals (Princeton University Press, 2008).
[2] Número reproductivo básico, https://www.ncbi.nlm.nih.gov/pmc/articles/
PMC6291769/.
[3] COVID-19 en España, https://cnecovid.isciii.es/covid19/#ccaa.
[4] F. Radicchi y G. Bianconi, «Epidemic plateau in critical susceptible-infected-removed dynamics with nontrivial initial conditions», Physical Review E102 (2020).
[5] W. O. Kermack y A. G. McKendrick, «A Contribution to the Mathematical Theory of Epidemics», Proceedings of the Royal Society of London. Series A, Containing Papers of a Mathematical and Physical Character115, 700-721 (1927).
[6] P. L. DeVries, A First Course in Computational Physics, 1st (John Wiley & Sons, Inc., USA, 1993).
[7] D. T. Gillespie, «Exact Stochastic Simulation of Coupled Chemical Reactions»,The Journal of Physical Chemistry81, 2340-2361 (1977).
[8] D. Balcan, V. Colizza, B. Gonçalves, H. Hu, J. J. Ramasco y A. Vespignani, «Multis-cale mobility networks and the spatial spreading of infectious diseases»,Proceedings of the National Academy of Sciences106, 21484-21489 (2009).
[9] O. Diekmann y J. Heesterbeek, Mathematical epidemiology of infectious diseases:
model building, analysis and interpretation, English, Wiley series in mathematical and computational biology (John Wiley y Sons, United States, 2000).
[10] Notes on R0, http : / / web . stanford . edu / ~jhj1 / teachingdocs / Jones on -R0.pdf.
[11] Compute R0 using Next Generation Operator, https://qrlssp.asu.edu/sites/
default/files/brn_mtbi_2016.pdf.
[12] The construction of next-generation matrices for compartmental epidemic models, https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2871801/#RSIF20090386C2.