-
¿Quieres apoyar a nuestro foro haciendo una donación?, entra aquí.
Exterminio de Plagas, Coronavirus: La pandemia del COVID-19
- Autor de tema Catrillanca
- Fecha de inicio
-
- Etiquetas
- coronavirus virus virus chino virus corona
- Estado
ABusheadorsh
Culiad@
- Registrado
- 2018/04/27
- Mensajes
- 3.338
- Sexo
No soy epidemiólogo, pero imparto clases de métodos matemáticos para ingenieros. En autocuarentena en casa, he decidido implementar en Matlab un modelo sencillo de la epidemia de coronavirus. En lugar de los conocidos modelos SIR y SEIR, me he decantado por el modelo SEIRV con tasas de transmisión variables, que es un poquito más realista. Por supuesto, estimar el valor de sus parámetros es difícil, así que he tomado los valores que se acaban de publicar para la epidemia en Wuhan, China. Quizás sea de interés para algunos de los lectores de este blog. Pero, cuidado, siendo un modelo matemático tan sencillo, que tiene en cuenta tan pocos factores, sus resultados no deben ser extrapolados sin ton ni son. La epidemiología rigurosa requiere usar modelos más complejos que incluyan toda la información disponible.
El modelo de juguete que voy a presentar es solo eso, un juguete que me servirá para ilustrar el uso de Matlab a la hora de resolver problemas de valores iniciales en ecuaciones diferenciales ordinarias. Por ello me centraré en lo fácil que es programar en este lenguaje este tipo de modelos matemáticos. He seleccionado el lenguaje Matlab porque es el que uso en mis clases para ingenieros (uso R en mis clases para bioquímicos); quienes prefieran una alternativa open source pueden usar Octave de GNU, que es un lenguaje casi idéntico. El código se puede adaptar de forma sencilla a otros lenguajes matemáticos, por ejemplo, a Mathematica usando NDSolve; pero en otros lenguajes hay que usar librerías o paquetes específicos, por ejemplo, deSolve en R, o diffeqpy en Python.
El artículo (open access) que me ha servido en bandeja el modelo y sus parámetros es Chayu Yang, Jin Wang, “A mathematical model for the novel coronavirus epidemic in Wuhan, China,” Mathematical Biosciences and Engineering 17: 2708-2724 (11 Mar 2020), doi: https://doi.org/10.3934/mbe.2020148.
% componentes del vector y = SEIRV
S=1;E=2;I=3;R=4;V=5;
dSEIRVdt = @(t,y) [
% dS/dt =
Lambda-betaE(y(E))*y(S)*y(E)-betaI(y(I))*y(S)*y(I)-betaV(y(V))*y(S)*y(V)-mu*y(S)
% dE/dt =
betaE(y(E))*y(S)*y(E)+betaI(y(I))*y(S)*y(I)+betaV(y(V))*y(S)*y(V)-(alfa+mu)*y(E)
% dI/dt =
alfa*y(E)-(w+gamma+mu)*y(I)
% dR/dt =
gamma*y(I)-mu*y(R)
% dV/dt =
xi1*y(E)+xi2*y(I)-sigma*y(V)
];
Introducir las ecuaciones diferenciales dY/dt = F(t,Y) en Matlab es muy fácil, basta definir la función F(t,Y). En este caso tenemos cinco ecuaciones acopladas, luego Y(t) será un vector de cinco componentes. He llamado a la función dSEIRVdt(t,Y) donde Y(1) es S, Y(2) es E, etc. En este código he usado una función anónima @(t,y) en Matlab. A mis estudiantes les recomiendo que no abusen de las funciones anónimas, porque el resultado es un código mucho más difícil de entender que cuando se usan funciones. Te muestro la diferencia.
function dYdt = dSEIRVdt (t,SEIRV)
% function dYdt = dSEIRVdt (t,SEIRV) % modelo SEIRV que usa constantes globales
%%
global Lambda betaE betaI betaV mu alfa w gamma xi1 xi2 sigma;
S=SEIRV(1); E=SEIRV(2); I=SEIRV(3); R=SEIRV(4); V=SEIRV(5);
%%
dSdt = Lambda-betaE(E)*S*E-betaI(I)*S*I-betaV(V)*S*V-mu*S;
dEdt = betaE(E)*S*E+betaI(I)*S*I+betaV(V)*S*V-(alfa+mu)*E;
dIdt = alfa*E-(w+gamma+mu)*I;
dRdt = gamma*I-mu*R;
dVdt = xi1*E+xi2*I-sigma*V;
%%
dYdt = [dSdt ; dEdt; dIdt; dRdt; dVdt ]
end
En el modelo SEIRV se usan cuatro poblaciones de humanos: susceptibles (S) son las personas que pueden ser infectadas; expuestas (E) al virus son las personas que lo contagian, pero que no presentan síntomas durante el periodo de incubación, por lo que no requieren hospitalización y aislamiento riguroso; infectadas (I) son las personas que presentan síntomas y están hospitalizadas, sin contacto con las expuestas; y recuperadas (R) son las personas que se han inmunizado contra el virus tras superar la infección; además, se incluye la concentración de virus (V) en el entorno debida a las personas expuestas (pues las infectadas influyen poco al estar hospitalizadas, sin contacto directo con las susceptibles).
Un modelo tan sencillo tiene 12 parámetros (los 11 mostradosmás un parámetro adicional común a todas las funciones beta) cuyos valores deben ser ajustados con datos epidemiológicos. Un ajuste completo es inviable (explorar un espacio con 12 dimensiones con 5 muestras por dimensión requiere unos 250 mil puntos; pero aún no hay tantos infectados). En la práctica, se estiman algunos parámetros de forma razonable y razonada, y se ajustan unos poquitos parámetros.
El artículo en que me baso estima 9 parámetros a partir de la literatura previa y ajusta 3 parámetros usando los datos epidemiológicos de Wuhan entre el 23 de enero de 2020 y el 10 de febrero de 2020, es decir, en la fase de crecimiento exponencial. Por tanto, dicho ajuste debe ser tomado como un simple ejemplo; en ningún caso será realista si se extrapola a otros países o a toda la pandemia.
% parametros del modelo SEIRV http://doi.org/10.3934/mbe.2020148
Lambda = 271.23 ;% flujo de poblacion hacia Wuhan (por dia)
betaE0 = 3.11e-8;% tasa de transmision entre S y E (/personas/dia)
betaI0 = 0.62e-8;% tasa de transmision entre S e I (/personas/dia)
betaV0 = 1.03e-8;% tasa de transmision entre S y V (/personas/dia)
c = 1.01e-4;% coeficiente de ajuste de la transmision
mu = 3.01e-5;% tasa natural de fallecimientos (por día)
alfa = 1/7 ;% 1/alfa = periodo de incubacion (/dias)
w = 0.01 ;% tasa de fallceimientos por la infeccion (por dia)
gamma = 1/15 ;% tasa de recuperacion de la infeccion (por dia)
sigma = 1 ;% tasa de eliminacion del virus del ambiente (por dia)
xi1 = 2.30 ;% dispersion del virus por infectados asintomaticos
xi2 = 0 ;% y por infectados aislados (por persona por dia por ml)
El parámetro Λ representa el flujo de entrada de personas a Wuhan; μ es la mortalidad de la población sin relación con la infeccion; 1/α es el periodo de incubación en el que la infección es asintomática (7 días); w es la mortalidad debida a la infección; 1/γ es el número de días para la recuperación de un infectado (15 días); ξ1 y ξ2 es la contribución a la concentración ambiental del virus debida a las personas expuestas e infectadas, resp.; y σ es el ritmo de degradación del virus en el entorno (tanto por causas ambientales, como por desinfecciones programadas). Las funciones βE(E) y βI(I) representan las tasas de contagio de persona a persona entre las expuestas y las susceptibles, y entre las infectadas y las susceptibles, resp.; la función βV(V) representa la tasa de contagio de una persona debida a la concentración de virus en el entorno. Estas tres funciones decrecen de forma monótona, pues a mayores valores de E, I y V se tomarán medidas de control más drásticas para reducir el contagio.
% funciones de contagio (positivas y decrecientes con su argumento)
betaE = @(E) betaE0/(1+c*E); % para expuestos
betaI = @(I) betaI0/(1+c*I); % para infectados
betaV = @(V) betaV0/(1+c*V); % para el ambiente
En el artículo se usa la misma expresión matemática para las tres funciones que describen el contagio; para no introducir tres nuevos parámetros, solo se usa un único parámetro c común a las tres, que depende de las medidas de control que tratan de evitar los contagios. Para c=0 se obtiene el modelo SEIRV con parámetros constantes.
Esta figura del artículo muestra que para c=0 los resultados no describen la evolución de la epidemia en Wuhan desde el 23 de enero de 2020 (día 0 en la figura). El modelo SEIRV con parámetros constantes predice casi tres millones de contagiados y todos sabemos que en Wuhan ni siquiera se han alcanzado los cien mil. Esto nos recuerda que ajustar el modelo con parámetros estimados en fecha muy temprana (antes del 10 de febrero) conlleva una predicción muy exagerada del pico máximo de infectados. En redes sociales mucha gente usa una simple función exponencial con el número de reproducción y obtienen valores anuméricos sin ningún sentido. Ten mucho cuidado cuando alguien te muestre ese tipo de cálculos; son simplemente basura anumérica.
Esta figura muestra la predicción del modelo para c estimado usando los datos entre 23 de enero y 10 de febrero. Predice un pico de personas infectadas de unas 50 mil, lo que es mucho más realista. Sin embargo, predice una caída del número de personas infectadas mucho más lenta de la observada en la epidemia en Wuhan (300 días es casi un año). Además, predice un número de personas expuestas (infectadas asintomáticas) inferior a las que acaban siendo infectadas. Habrás leído en redes sociales a mucha gente estimar que tiene que haber entre 10 y 30 veces más personas expuestas que infectadas; este modelo no da cuenta de ese hecho (incluso aunque inicialmente haya más personas expuestas que infectadas).
El modelo SEIRV del artículo es muy sencillo. Sus conclusiones deben ser tomadas con mucha cautela. El artículo me ha gustado porque incluye un cálculo de los puntos de equilibrio y demuestra que para una número de reproducción R0>1 siempre existe un punto de equilibrio no trivial. Esta figura te muestra el plano de fases para las personas expuestas (E) e infectadas (I). El punto fijo calculado para R0 = 4.25 es E* = 1353, e I* = 2735, con S* = 2 583 683, R* = 6 528 015, y V* = 3111. Este punto fijo implica que la infección se vuelve endémica; siempre habrá un reservorio de virus entre la población. Aunque el objetivo de las medidas de la OMS y de muchos Gobiernos ha sido impedir este hecho, cada día parece más difícil evitar que así ocurra con la pandemia de COVID-19. Pero, repito de nuevo, el modelo SEIRV es muy sencillo y sus conclusiones debe ser puestas en cuarentena. Una vez pase la pandemia, gracias a modelos mucho más complicados, se podrá saber si el nuevo coronavirus se ha convertido en un nuevo coronavirus humano estacional.
% condicion inicial
S0 = 8998505 ;% susceptibles de infectar
E0 = 1000 ;% expuestos (infectados asintomaticos)
I0 = 475 ;% infectados (con sintomas y aislados)
R0 = 10 ;% recuperados
V0 = 10000 ;% concentracion de virus en el ambiente
Volvamos a mi código en Matlab. La condición inicial que se usa en el artículo está basada en información oficial de Wuhan (fuente). Se considera una población de susceptibles de casi 9 millones de habitantes en Wuhan, 475 personas infectadas, 10 recuperadas y 1000 personas expuestas, todo ello a fecha de 23 de enero de 2020. La carga viral inicial en el entorno es muy difícil de estimar, pero los autores proponen un valor de 10 mil partículas virales por mililitro; no hay ningún razonamiento en el artículo que justifique este número, que me parece que a sido elegido más a ojo que otra cosa.
% resolucion numerica
SEIRV0 = [S0; E0; I0; R0; V0];
[t,SEIRV] = ode23t(dSEIRVdt,[0,3000],SEIRV0);
La ventaja de usar Matlab (o lenguajes matemáticos como Mathematica, R, Python, etc.) es que la resolución del sistema de ecuaciones diferenciales ordinarias solo requiere una única línea. En este código uso la función ode23t, que se suele usar en problemas moderadamente rígidos (moderate stiffness). ¿Por qué la uso? Porque he probado con todas…
% probar toda la odesuite
odesuite = {“ode23″,”ode45″,”ode15s”,”ode23s”,”ode23t”,”ode23tb”};
for ii=1:length(odesuite)
time = cputime;
eval(strcat(‘[t,SEIRV]=’,odesuite{ii},'(dSEIRVdt,[0,30000],SEIRV0);’));
cost(ii) = cputime-time,
end
He probado con todos los métodos de la odesuite de Matlab. El coste en tiempo de CPU ha sido de 0.58, 1.0, 0.047, 0.016, 0.016, y 0.016 segundos. Así que los métodos ode23s, ode23t, y ode23tb tardan unas 66 veces menos que ode45; por ello he elegido ode23t.
clf;
SEIRV0 = [S0; E0; I0; R0; V0];
[t,SEIRV] = ode23t(dSEIRVdt,[0,300],SEIRV0);
subplot(2,2,1); plot(t,SEIRV,E),t,SEIRV
[t,SEIRV] = ode23t(dSEIRVdt,[0,3000],SEIRV0);
subplot(2,2,2); plot(t,SEIRV
Este código y esta figura muestra los resultados de Matlab para el modelo SEIRV. El pico de la curva es razonable pero la caída del número de personas infectadas ha sido mucho más rápida (300 días es casi un año). Sin lugar a dudas, el método requiere una nueva estimación de sus parámetros.
% dibuja el plano de fase
I0v = [ 5 10 15 15 15 15 15 10 6 3 0 0 0 0 ]*1e4;
E0v = [ 0 0 0 1 2 3 4 4 4 4 3.5 2.5 1.5 0.5 ]*1e4;
clf; box on;
for iv=1:length(I0v)
SEIRV0 = [S0; E0v(iv); I0v(iv); R0; V0];
[t,SEIRV] = ode23t(dSEIRVdt,[0,50000],SEIRV0);
subplot(2,2,1); plot(SEIRV
subplot(2,2,2); semilogx(t,SEIRV
end;
subplot(2,2,1); xlabel(‘E’); ylabel(‘I’); title(‘phase portrait’);
subplot(2,2,2); xlabel(‘t’); title(‘E, I’); axis([0.1 5e4 0 1.5e5]);
Este código y esta figura muestran el plano de fases (E, I) a la derecha, así como las soluciones E(t) e I(t) en función del tiempo (en días). Se observa que las soluciones tienden a un punto intermedio alrededor de t = 100 días con E(100) ≈ 21 000, e I(100) = 39 000; tras este punto, todas las soluciones en el pano de fase evolucionan hacia el mismo punto de equilibrio E* = 1347, e I* = 2509. Este punto no coincide con el que aparece en el artículo, pero se puede calcular directamente en Matlab.
% calcula punto de equilibrio
SEIRVstar = fsolve(@(y) dSEIRVdt(0,y), SEIRV0)
SEIRVstar([E, I])
Recuerda que un punto de equilibrio del sistema de ecuaciones diferenciales ordinarias dY/dt = F(Y) es un punto Y* tal que F(Y*)=0; por tanto dicho punto es estacionario (usado como condición inicial el sistema se mantiene en ese punto como solución constante). El cálculo de Matlab ofrece el mismo resultado que se ha determinado con ode23t, que no coincide con el artículo. Supongo que la razón es el error introducido por el método numérico usado por los autores (que no aclaran cuál es en su artículo).
En resumen, sirva esta pieza como ejemplo del uso de Matlab para resolver ecuaciones diferenciales ordinarias. Los resultados epidemiológicos mostrados son poco rigurosos, pues el modelo SEIRV es muy sencillo (y yo creo que habría que reestimar sus parámetros a partir de los resultados actuales de la epidemia en Wuhan para obtener un mejor ajuste). Sin embargo, es un modelo mucho mejor que usar una ley exponencial, aunque la hayamos visto tantas veces en redes sociales que duele a los ojos tanto anumerismo. Si conoces Matlab (u Octave) , te recomiendo ejecutar el código con diferentes parámetros; usa los comentarios si logras resultados que merezca la pena remarcar.
hola señor ha subido bastante el nivel de discusión
Se que los modelos matemáticos que son usados para simular el crecimiento de microorganismos son muy fiables ,pero mas allá del programa Matlab u otro software de simulación,un complemento muy importante es contar con un computador capaz que de modelar del crecimiento de esos seres o incluso como se comportara la infeccion en el globo.
Acá por ejemplo la mayoría de los centros de investigación avanzada en el área de biotecnología no usan o no tiene simuladores, principalmente por el costo de adquirir y operar un supercomputador es muy elevado, a excepcion de la U chile que tiene varios Superquipos hay uno que tiene la capacidad de calculo de 25.mil cpu con ellos han logrado me parece algo relacionado con la biolixiviacion.
Es un hecho que las potencias los esten usando para predecir varios aspectos biologicos hasta los economicos que dejara esta pandemia, asi tambien variaciones de infectividad en diversos climas, analizar su estructura y encontrar una ventana para neutralizarla.
saludos
Clonosaurio
Babos@
- Registrado
- 2016/06/01
- Mensajes
- 1.320
- Sexo
todas las universidades suspendiendo las clases, los academicos rogando por medidas mas extremas y el saco wea que esta en la moneda como si nada
El brazos cortos debe estar como loco revisando y comprando acciones a la baja, de compañías con potencial de recuperarse cuando pase la emergencia.
Al pueblo de chile simplemente se lo pasa por el pico , y no va a tomar ninguna medida efectiva contra el coronavirus hasta que la presión sea insostenible.
Es su venganza por haber quedado en ridículo, cuando su megalomanía lo hacia querer mostrarse como un líder mundial, en la COP25 y en la APEC del año pasado
Ojalá las crisis del último tiempo nos hagan replantearnos, y a futuro dejemos de tener un régimen de gobierno presidencialista, porque como pais elegimos puros weones. Primero chanchelet y ahora este pelmazo. Dos weones que odian a Chile a mas no poder, y nos transformarán en un paisaje de ruinas y cementerios.
CptPrice
Vac@
- Registrado
- 2019/11/11
- Mensajes
- 1.981
- Sexo
Nunca te había leído que putearas a alguien hermanito
Egregore
Animal
- Registrado
- 2019/11/30
- Mensajes
- 2.187
- Sexo
Depende para qué la quieres.
Si es para filtrar partículas debe tener sello completo en la cara, o sea respiras a través del material de la mascarilla (por ejemplo para cortar fibra de vidrio, nubes de aserrín, cemento, etc, y en este caso aerosol de saliva)
Si es para otros usos menos importantes te podría servir, como para eso que dice hacer, o para ocultar qué tipo de mascarilla llevas debajo.
Sixpool
Plasta
- Registrado
- 2019/10/12
- Mensajes
- 767
¿Te refieres a esa tecnología 5G que es la verdadera guerra de las comunicaciones hoy en día? ¿De las antenas 5G que está plagada la ciudad china de Wuhan y que otros países que han instalado muchas de esas antenas son Italia e Irán? ¿La tecnología 5G cuestionada por bajar las defensas de las personas cerca de sus antenas?
Egregore
Animal
- Registrado
- 2019/11/30
- Mensajes
- 2.187
- Sexo
Yo la simplifiqué para que fuera factible fabricarla con el mínimo viable.Como que es ligeramente diferente a la receta del @Egregore con unos cuantos ingredientes más
Prueba comprar los insumos que aparecen listados en esa.
Si tienes acceso a insumos químicos, de más te sirve esa otra.
Pero el "principio activo" del alcohol gel es: alcohol final a 70% y gel para que el daño "celular" repetitivo en tus manos no sea tan fuerte.
Chupallita
Maricon
- Registrado
- 2008/04/04
- Mensajes
- 4.450
- Sexo
Wena la propaganda!!!
Bob Virutilla
Exiliado
- Registrado
- 2020/02/28
- Mensajes
- 467
Mas pasta de estadista que el esperpento que tenemos por presidente
UsuarioFeliz
Come Mierda
- Registrado
- 2019/04/16
- Mensajes
- 19.384
- Sexo
Te acordai cuando el piñera en una reunion de conchesumadres encaro a los otros lideres, dandoselaxde eclectico el culiao
Oye wn. Recapitulando la gorda les metio el pico supremo en el ojo a todos los chilenos wn. No hay wea en que no haya dejado la caga......
lmfao
Bestia
- Registrado
- 2018/12/23
- Mensajes
- 2.691
- Sexo
comercios, pymes, independientes o cualquier personas que vivan de trabajo relacionado esta en riesgo.
En estos casos es el estado el mejor candidato para ayudar, con dinero y de forma directa (cash is king). Es eso o el colapso del país.
Focusom
Exiliado
- Registrado
- 2019/10/20
- Mensajes
- 1.427
- Sexo
Lo unico que me alegra es que después de la gran cagasa que ya se mandó este gobierno. Solo que aun no lo saben, el sin bolas quedará como el peor presidente de la historia incluso peor que el borracho chicho.
Ojalá que a futuro tenga que pagar con carcel su abandono de labores
Ojalá que a futuro tenga que pagar con carcel su abandono de labores
Minsal detalla preparación para eventual fase 4 por Covid-19: FF.AA. dispondrán instalaciones para atención de pacientes
Fuente: Emol.com - https://www.emol.com/noticias/Nacional/2020/03/14/979851/Minsal-coronavirus-apoyo-FFAA.html
https://www.emol.com/noticias/Nacional/2020/03/14/979851/Minsal-coronavirus-apoyo-FFAA.html
CAGAMOS......
Fuente: Emol.com - https://www.emol.com/noticias/Nacional/2020/03/14/979851/Minsal-coronavirus-apoyo-FFAA.html
https://www.emol.com/noticias/Nacional/2020/03/14/979851/Minsal-coronavirus-apoyo-FFAA.html
CAGAMOS......
LordIzziNosnam00
¡Colaborador!
- Registrado
- 2006/02/26
- Mensajes
- 3.489
- Sexo
Eso amerita la apertura de un nuevo thread en una sección como el area 51
¿No has escuchado hablar de los "filósofos", gente que se dedica a darle vueltas a asuntos como las formas de pensar ideas y su relación con la realidad, o de "místicos" o "ascetas" que se aislan del mundanal ruido en un acto introspectivo de meditación profunda? Quien te haya dicho eso, y convencido de ello, no tiene idea del poder de la propia mente.
No soy epidemiólogo, pero imparto clases de métodos matemáticos para ingenieros. En autocuarentena en casa, he decidido implementar en Matlab un modelo sencillo de la epidemia de coronavirus. En lugar de los conocidos modelos SIR y SEIR, me he decantado por el modelo SEIRV con tasas de transmisión variables, que es un poquito más realista. Por supuesto, estimar el valor de sus parámetros es difícil, así que he tomado los valores que se acaban de publicar para la epidemia en Wuhan, China. Quizás sea de interés para algunos de los lectores de este blog. Pero, cuidado, siendo un modelo matemático tan sencillo, que tiene en cuenta tan pocos factores, sus resultados no deben ser extrapolados sin ton ni son. La epidemiología rigurosa requiere usar modelos más complejos que incluyan toda la información disponible.
El modelo de juguete que voy a presentar es solo eso, un juguete que me servirá para ilustrar el uso de Matlab a la hora de resolver problemas de valores iniciales en ecuaciones diferenciales ordinarias. Por ello me centraré en lo fácil que es programar en este lenguaje este tipo de modelos matemáticos. He seleccionado el lenguaje Matlab porque es el que uso en mis clases para ingenieros (uso R en mis clases para bioquímicos); quienes prefieran una alternativa open source pueden usar Octave de GNU, que es un lenguaje casi idéntico. El código se puede adaptar de forma sencilla a otros lenguajes matemáticos, por ejemplo, a Mathematica usando NDSolve; pero en otros lenguajes hay que usar librerías o paquetes específicos, por ejemplo, deSolve en R, o diffeqpy en Python.
El artículo (open access) que me ha servido en bandeja el modelo y sus parámetros es Chayu Yang, Jin Wang, “A mathematical model for the novel coronavirus epidemic in Wuhan, China,” Mathematical Biosciences and Engineering 17: 2708-2724 (11 Mar 2020), doi: https://doi.org/10.3934/mbe.2020148.
% componentes del vector y = SEIRV
S=1;E=2;I=3;R=4;V=5;
dSEIRVdt = @(t,y) [
% dS/dt =
Lambda-betaE(y(E))*y(S)*y(E)-betaI(y(I))*y(S)*y(I)-betaV(y(V))*y(S)*y(V)-mu*y(S)
% dE/dt =
betaE(y(E))*y(S)*y(E)+betaI(y(I))*y(S)*y(I)+betaV(y(V))*y(S)*y(V)-(alfa+mu)*y(E)
% dI/dt =
alfa*y(E)-(w+gamma+mu)*y(I)
% dR/dt =
gamma*y(I)-mu*y(R)
% dV/dt =
xi1*y(E)+xi2*y(I)-sigma*y(V)
];
Introducir las ecuaciones diferenciales dY/dt = F(t,Y) en Matlab es muy fácil, basta definir la función F(t,Y). En este caso tenemos cinco ecuaciones acopladas, luego Y(t) será un vector de cinco componentes. He llamado a la función dSEIRVdt(t,Y) donde Y(1) es S, Y(2) es E, etc. En este código he usado una función anónima @(t,y) en Matlab. A mis estudiantes les recomiendo que no abusen de las funciones anónimas, porque el resultado es un código mucho más difícil de entender que cuando se usan funciones. Te muestro la diferencia.
function dYdt = dSEIRVdt (t,SEIRV)
% function dYdt = dSEIRVdt (t,SEIRV) % modelo SEIRV que usa constantes globales
%%
global Lambda betaE betaI betaV mu alfa w gamma xi1 xi2 sigma;
S=SEIRV(1); E=SEIRV(2); I=SEIRV(3); R=SEIRV(4); V=SEIRV(5);
%%
dSdt = Lambda-betaE(E)*S*E-betaI(I)*S*I-betaV(V)*S*V-mu*S;
dEdt = betaE(E)*S*E+betaI(I)*S*I+betaV(V)*S*V-(alfa+mu)*E;
dIdt = alfa*E-(w+gamma+mu)*I;
dRdt = gamma*I-mu*R;
dVdt = xi1*E+xi2*I-sigma*V;
%%
dYdt = [dSdt ; dEdt; dIdt; dRdt; dVdt ]
end
En el modelo SEIRV se usan cuatro poblaciones de humanos: susceptibles (S) son las personas que pueden ser infectadas; expuestas (E) al virus son las personas que lo contagian, pero que no presentan síntomas durante el periodo de incubación, por lo que no requieren hospitalización y aislamiento riguroso; infectadas (I) son las personas que presentan síntomas y están hospitalizadas, sin contacto con las expuestas; y recuperadas (R) son las personas que se han inmunizado contra el virus tras superar la infección; además, se incluye la concentración de virus (V) en el entorno debida a las personas expuestas (pues las infectadas influyen poco al estar hospitalizadas, sin contacto directo con las susceptibles).
Un modelo tan sencillo tiene 12 parámetros (los 11 mostradosmás un parámetro adicional común a todas las funciones beta) cuyos valores deben ser ajustados con datos epidemiológicos. Un ajuste completo es inviable (explorar un espacio con 12 dimensiones con 5 muestras por dimensión requiere unos 250 mil puntos; pero aún no hay tantos infectados). En la práctica, se estiman algunos parámetros de forma razonable y razonada, y se ajustan unos poquitos parámetros.
El artículo en que me baso estima 9 parámetros a partir de la literatura previa y ajusta 3 parámetros usando los datos epidemiológicos de Wuhan entre el 23 de enero de 2020 y el 10 de febrero de 2020, es decir, en la fase de crecimiento exponencial. Por tanto, dicho ajuste debe ser tomado como un simple ejemplo; en ningún caso será realista si se extrapola a otros países o a toda la pandemia.
% parametros del modelo SEIRV http://doi.org/10.3934/mbe.2020148
Lambda = 271.23 ;% flujo de poblacion hacia Wuhan (por dia)
betaE0 = 3.11e-8;% tasa de transmision entre S y E (/personas/dia)
betaI0 = 0.62e-8;% tasa de transmision entre S e I (/personas/dia)
betaV0 = 1.03e-8;% tasa de transmision entre S y V (/personas/dia)
c = 1.01e-4;% coeficiente de ajuste de la transmision
mu = 3.01e-5;% tasa natural de fallecimientos (por día)
alfa = 1/7 ;% 1/alfa = periodo de incubacion (/dias)
w = 0.01 ;% tasa de fallceimientos por la infeccion (por dia)
gamma = 1/15 ;% tasa de recuperacion de la infeccion (por dia)
sigma = 1 ;% tasa de eliminacion del virus del ambiente (por dia)
xi1 = 2.30 ;% dispersion del virus por infectados asintomaticos
xi2 = 0 ;% y por infectados aislados (por persona por dia por ml)
El parámetro Λ representa el flujo de entrada de personas a Wuhan; μ es la mortalidad de la población sin relación con la infeccion; 1/α es el periodo de incubación en el que la infección es asintomática (7 días); w es la mortalidad debida a la infección; 1/γ es el número de días para la recuperación de un infectado (15 días); ξ1 y ξ2 es la contribución a la concentración ambiental del virus debida a las personas expuestas e infectadas, resp.; y σ es el ritmo de degradación del virus en el entorno (tanto por causas ambientales, como por desinfecciones programadas). Las funciones βE(E) y βI(I) representan las tasas de contagio de persona a persona entre las expuestas y las susceptibles, y entre las infectadas y las susceptibles, resp.; la función βV(V) representa la tasa de contagio de una persona debida a la concentración de virus en el entorno. Estas tres funciones decrecen de forma monótona, pues a mayores valores de E, I y V se tomarán medidas de control más drásticas para reducir el contagio.
% funciones de contagio (positivas y decrecientes con su argumento)
betaE = @(E) betaE0/(1+c*E); % para expuestos
betaI = @(I) betaI0/(1+c*I); % para infectados
betaV = @(V) betaV0/(1+c*V); % para el ambiente
En el artículo se usa la misma expresión matemática para las tres funciones que describen el contagio; para no introducir tres nuevos parámetros, solo se usa un único parámetro c común a las tres, que depende de las medidas de control que tratan de evitar los contagios. Para c=0 se obtiene el modelo SEIRV con parámetros constantes.
Esta figura del artículo muestra que para c=0 los resultados no describen la evolución de la epidemia en Wuhan desde el 23 de enero de 2020 (día 0 en la figura). El modelo SEIRV con parámetros constantes predice casi tres millones de contagiados y todos sabemos que en Wuhan ni siquiera se han alcanzado los cien mil. Esto nos recuerda que ajustar el modelo con parámetros estimados en fecha muy temprana (antes del 10 de febrero) conlleva una predicción muy exagerada del pico máximo de infectados. En redes sociales mucha gente usa una simple función exponencial con el número de reproducción y obtienen valores anuméricos sin ningún sentido. Ten mucho cuidado cuando alguien te muestre ese tipo de cálculos; son simplemente basura anumérica.
Esta figura muestra la predicción del modelo para c estimado usando los datos entre 23 de enero y 10 de febrero. Predice un pico de personas infectadas de unas 50 mil, lo que es mucho más realista. Sin embargo, predice una caída del número de personas infectadas mucho más lenta de la observada en la epidemia en Wuhan (300 días es casi un año). Además, predice un número de personas expuestas (infectadas asintomáticas) inferior a las que acaban siendo infectadas. Habrás leído en redes sociales a mucha gente estimar que tiene que haber entre 10 y 30 veces más personas expuestas que infectadas; este modelo no da cuenta de ese hecho (incluso aunque inicialmente haya más personas expuestas que infectadas).
El modelo SEIRV del artículo es muy sencillo. Sus conclusiones deben ser tomadas con mucha cautela. El artículo me ha gustado porque incluye un cálculo de los puntos de equilibrio y demuestra que para una número de reproducción R0>1 siempre existe un punto de equilibrio no trivial. Esta figura te muestra el plano de fases para las personas expuestas (E) e infectadas (I). El punto fijo calculado para R0 = 4.25 es E* = 1353, e I* = 2735, con S* = 2 583 683, R* = 6 528 015, y V* = 3111. Este punto fijo implica que la infección se vuelve endémica; siempre habrá un reservorio de virus entre la población. Aunque el objetivo de las medidas de la OMS y de muchos Gobiernos ha sido impedir este hecho, cada día parece más difícil evitar que así ocurra con la pandemia de COVID-19. Pero, repito de nuevo, el modelo SEIRV es muy sencillo y sus conclusiones debe ser puestas en cuarentena. Una vez pase la pandemia, gracias a modelos mucho más complicados, se podrá saber si el nuevo coronavirus se ha convertido en un nuevo coronavirus humano estacional.
% condicion inicial
S0 = 8998505 ;% susceptibles de infectar
E0 = 1000 ;% expuestos (infectados asintomaticos)
I0 = 475 ;% infectados (con sintomas y aislados)
R0 = 10 ;% recuperados
V0 = 10000 ;% concentracion de virus en el ambiente
Volvamos a mi código en Matlab. La condición inicial que se usa en el artículo está basada en información oficial de Wuhan (fuente). Se considera una población de susceptibles de casi 9 millones de habitantes en Wuhan, 475 personas infectadas, 10 recuperadas y 1000 personas expuestas, todo ello a fecha de 23 de enero de 2020. La carga viral inicial en el entorno es muy difícil de estimar, pero los autores proponen un valor de 10 mil partículas virales por mililitro; no hay ningún razonamiento en el artículo que justifique este número, que me parece que a sido elegido más a ojo que otra cosa.
% resolucion numerica
SEIRV0 = [S0; E0; I0; R0; V0];
[t,SEIRV] = ode23t(dSEIRVdt,[0,3000],SEIRV0);
La ventaja de usar Matlab (o lenguajes matemáticos como Mathematica, R, Python, etc.) es que la resolución del sistema de ecuaciones diferenciales ordinarias solo requiere una única línea. En este código uso la función ode23t, que se suele usar en problemas moderadamente rígidos (moderate stiffness). ¿Por qué la uso? Porque he probado con todas…
% probar toda la odesuite
odesuite = {“ode23″,”ode45″,”ode15s”,”ode23s”,”ode23t”,”ode23tb”};
for ii=1:length(odesuite)
time = cputime;
eval(strcat(‘[t,SEIRV]=’,odesuite{ii},'(dSEIRVdt,[0,30000],SEIRV0);’));
cost(ii) = cputime-time,
end
He probado con todos los métodos de la odesuite de Matlab. El coste en tiempo de CPU ha sido de 0.58, 1.0, 0.047, 0.016, 0.016, y 0.016 segundos. Así que los métodos ode23s, ode23t, y ode23tb tardan unas 66 veces menos que ode45; por ello he elegido ode23t.
clf;
SEIRV0 = [S0; E0; I0; R0; V0];
[t,SEIRV] = ode23t(dSEIRVdt,[0,300],SEIRV0);
subplot(2,2,1); plot(t,SEIRV
[t,SEIRV] = ode23t(dSEIRVdt,[0,3000],SEIRV0);
subplot(2,2,2); plot(t,SEIRV
Este código y esta figura muestra los resultados de Matlab para el modelo SEIRV. El pico de la curva es razonable pero la caída del número de personas infectadas ha sido mucho más rápida (300 días es casi un año). Sin lugar a dudas, el método requiere una nueva estimación de sus parámetros.
% dibuja el plano de fase
I0v = [ 5 10 15 15 15 15 15 10 6 3 0 0 0 0 ]*1e4;
E0v = [ 0 0 0 1 2 3 4 4 4 4 3.5 2.5 1.5 0.5 ]*1e4;
clf; box on;
for iv=1:length(I0v)
SEIRV0 = [S0; E0v(iv); I0v(iv); R0; V0];
[t,SEIRV] = ode23t(dSEIRVdt,[0,50000],SEIRV0);
subplot(2,2,1); plot(SEIRV
subplot(2,2,2); semilogx(t,SEIRV
end;
subplot(2,2,1); xlabel(‘E’); ylabel(‘I’); title(‘phase portrait’);
subplot(2,2,2); xlabel(‘t’); title(‘E, I’); axis([0.1 5e4 0 1.5e5]);
Este código y esta figura muestran el plano de fases (E, I) a la derecha, así como las soluciones E(t) e I(t) en función del tiempo (en días). Se observa que las soluciones tienden a un punto intermedio alrededor de t = 100 días con E(100) ≈ 21 000, e I(100) = 39 000; tras este punto, todas las soluciones en el pano de fase evolucionan hacia el mismo punto de equilibrio E* = 1347, e I* = 2509. Este punto no coincide con el que aparece en el artículo, pero se puede calcular directamente en Matlab.
% calcula punto de equilibrio
SEIRVstar = fsolve(@(y) dSEIRVdt(0,y), SEIRV0)
SEIRVstar([E, I])
Recuerda que un punto de equilibrio del sistema de ecuaciones diferenciales ordinarias dY/dt = F(Y) es un punto Y* tal que F(Y*)=0; por tanto dicho punto es estacionario (usado como condición inicial el sistema se mantiene en ese punto como solución constante). El cálculo de Matlab ofrece el mismo resultado que se ha determinado con ode23t, que no coincide con el artículo. Supongo que la razón es el error introducido por el método numérico usado por los autores (que no aclaran cuál es en su artículo).
En resumen, sirva esta pieza como ejemplo del uso de Matlab para resolver ecuaciones diferenciales ordinarias. Los resultados epidemiológicos mostrados son poco rigurosos, pues el modelo SEIRV es muy sencillo (y yo creo que habría que reestimar sus parámetros a partir de los resultados actuales de la epidemia en Wuhan para obtener un mejor ajuste). Sin embargo, es un modelo mucho mejor que usar una ley exponencial, aunque la hayamos visto tantas veces en redes sociales que duele a los ojos tanto anumerismo. Si conoces Matlab (u Octave) , te recomiendo ejecutar el código con diferentes parámetros; usa los comentarios si logras resultados que merezca la pena remarcar.
Pudo haber mostrado los resultados de su modelo matemático po cipadrito
Claro claro, lo que pasa es que según me cuenta un amigo químico farmaceuta, a tu receta y a las 2 de la OMS que se compartieron aquí, les falta el espesante, que puede ser colapez o gelatina sin sabor (aparte de otros químicos difíciles de conseguir, que se me olvidaron en el instante que me los dijo
), porque en caso contrario, estas recetas van a durar poco y se evaporará fácilmente entre apertura-y-apertura del recipiente en uso.
Chupallita
Maricon
- Registrado
- 2008/04/04
- Mensajes
- 4.450
- Sexo
Mmmmmmmm... un régimen parlamentario dices tú? Ándate a la....
- Estado