Diferencia entre revisiones de «Modelos epidemiológicos Grupo 6C»

Revisión actual del 14:17 6 mar 2015

Trabajo realizado por estudiantes
Título	Modelos epidemiológicos. Grupo 6-C
Asignatura	Ecuaciones Diferenciales
Curso	Curso 2014-15
Autores	Manuel Morales López 1175 David Toledo Menéndez 1228 Sergio Rodríguez Torcal 994 Jose María Rodríguez Vicente 1213 Lourdes Sánchez-Ocaña Merino 1248 Jorge Villa Lobo 1237
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

Contenido

1 Enunciado
2 Interpretación de parámetros
3 Método de Euler y Trapecio con S=0
4 Método de Euler y Trapecio con S=100
5 Modelo completo
6 Comparación Runge-Kutta con Euler
7 Tasa de infectados por contagio dependiente del tiempo. Método de Heun
8 Calibración del parámetro a

1 Enunciado

En el desarrollo de una epidemia se distinguen dos tipos de individuos: los que ya han contraido la enfermedad o infectados I, y los que son susceptibles de contraerla por encontrarse en zona de riesgo S. Supongamos que se dan las siguientes hip´otesis: 1. La poblaci´on de personas infectadas se altera por el fallecimiento o la cura de las mismas. En ambos casos, la tasa de cambio depende del n´umero de personas infectadas; 2. La tasa de individuos que pasan de ser susceptibles a contraer la enfermedad a estar infectados es proporcional a la interacci´on entre el n´umero de individuos en ambas clases. Consideramos las variables: t tiempo, S(t) poblaci´on de individuos susceptibles a contraer la enfermedad, I(t) poblaci´on de individuos infectados; y el sistema: dS dt = −aSI dI dt = aSI − bI − cI donde a, b, c son parametros.

2 Interpretación de parámetros

En el problema: el coeficiente "a" es la tasa de infectados por contagio, "b" la de muertos y "c" la de curados.

3 Método de Euler y Trapecio con S=0

%Apartado 2 trabajo1 2015
clear all
t0=0;
tN=20;
y0=2000;
h=0.1;
t=t0:h:tN;
N=(tN-t0)/h;
y=zeros(1,N+1); %euler 
y(1)=y0;
z=zeros(1,N+1); %trapecio 
z(1)=y0;
w=zeros(1,N+1);
w(1)=y0;
x=zeros(1,N+1);
x(1)=y0;
for i=1:N
y(i+1)=y(i)-h*(0.3+0.001)*y(i);
z(i+1)=(1/(1-h/2*(-0.31)))*(z(i))+h/2*(-0.31*z(i));
end
[t',y',z'];
i=1;
while abs(y0-w(i))<3/4*y0
  w(i+1)=(1/(1-h/2*(-0.31)))*(w(i))+h/2*(-0.31*w(i)); %trapecio, evitamos la variable auxiliar yy
  t(i+1)=t(i)+h;
  i=i+1;
end
disp('Tiempo final con trapecio:')
disp(t(i))
%gráfico
i=1;
while abs(y0-x(i))<3/4*y0
  x(i+1)=x(i)-h*(0.3+0.001)*x(i); %euler, evitamos la variable auxiliar yy
  t(i+1)=t(i)+h;
  i=i+1;
end
disp('Tiempo final con Euler:')
disp(t(i))
hold on
plot(t,y,'b')
plot(t,z,'r')
plot(t,w,'g')
plot(t,x,'k')
legend('Euler','Trapecio','Tiempo que tarda por Trapecio','Tiempo que tarda por Euler','Location','best'); 
hold off

Esta gráfica nos muestra que en t=4.60 se alcanza la condición final, que el número de infectados se reduzca a 500 mientras que con el trapecio el tiempo final es t=4.50

Tiempo en llegar a 500 infectados

4 Método de Euler y Trapecio con S=100

Se puede interpretar como que a partir de un valor limite de "S" entre 100 y 200, el numero de infectados se mantiene constante en el tiempo. De igual forma, si "S" es menor que este valor limite, el numero de infectados desciende en el tiempo, y si es mayor que el valor limite "S", asciende.

%Apartado 3 trabajo 6-C 2015
clear all
t0=0; tN=40; %Intervalo de tiempo tomado
S0=input('introduce valor de población susceptible: '); 
I0=2000; %Valores iniciales de infectados(I) y susceptibles(S) 
h=0.1; %Determinación del paso 
t=t0:h:tN; %Desarrollo del tiempo desde t0 hasta tN tomando intervalos de paso h
a=0.003; %valores de los parametros
b=0.3;
c=0.01;
N=(tN-t0)/h;
S(1)=S0; %Asignacion del valor incial para la primera componente
I(1)=I0;
 
%Resolucion del sistema de forma matricial
for n=1:N
   A=[S0;I(n)]+h*[-a*S0*I(n);a*S0*I(n)-(b+c)*I(n)];
   S(n+1)=S0;%Asignacion de los distintos valores de S como la primera compnente de la matriz A
   I(n+1)=A(2); %Asignacion de los distintos valores de I como la primera componente de la matriz A
end

%Interpretación gráfica 
hold on
figure(1) 
plot(t,S,'r')
plot(t,I,'b')
hold off

Gráfica para 100 supceptibles

Gráfica para 200 supceptibles

Podemos llegar a la conclusión de que el número de personas infectadas varía notablemente cuando el número de susceptibles se mantiene constante a lo largo del tiempo. Esto se debe a que cuanto mayor sea el número de susceptibles, mayor será el número de infectados, es decir, con 200 susceptibles, mantiene un ritmo constante hasta que llega un momento en el cual el número de infectados se dispara, debido a que la tasa de personas susceptibles es mayor que las personas curadas o fallecidas. En conclusión, debe de haber un valor de personas susceptibles para que la función de infectados sea lineal, es decir, sea constante en el tiempo.

5 Modelo completo

%Apartado 4 trabajo 6-C 2015
clear all
t0=0; tN=40; %Intervalo de tiempo tomado
S0=input('introduce valor de población susceptible: '); 
I0=input('introduce valor de población infectada inicial: '); %Valores iniciales de infectados(I) y susceptibles(S) 
h=input('introduce valor del paso de tiempo: '); %Determinación del paso 
t=t0:h:tN; %Desarrollo del tiempo desde t0 hasta tN tomando intervalos de paso h
a=0.003; %valores de los parametros
b=0.3;
c=0.01;
N=(tN-t0)/h;
S(1)=S0; %Asignacion del valor incial para la primera componente
I(1)=I0;
 
%Resolucion del sistema de forma matricial
for n=1:N
   A=[S(n);I(n)]+h*[-a*S(n)*I(n);a*S(n)*I(n)-(b+c)*I(n)];
   S(n+1)=A(1);%Asignacion de los distintos valores de S como la primera compnente de la matriz A
   I(n+1)=A(2); %Asignacion de los distintos valores de I como la primera componente de la matriz A
end

%Interpretación gráfica 
hold on
figure(1) 
plot(t,S,'r')
plot(t,I,'b')
hold off
maximo=max(I)
posicion=find(I==maximo)
tiempo=t(posicion)

Como podemos apreciar en las gráficas el número máximo de enfermos esperados con 800 susceptibles y 20 infectados inicialmente, es aproximadamente 500, mientras que para 10000 susceptibles y 40 infectados inicialmente, observamos un aumento de personas enfermas, aproximadamente de 9500.

left

right

left

right

left

right

left

right

6 Comparación Runge-Kutta con Euler

%Apartado 5 trabajo 6-C 2015
clear all
t0=0; tN=30; %Intervalo de tiempo tomado
y1=input('introduce valor de población susceptible inicial: ');
y2=input('introduce valor de población infectada inicial: ');
y0=[y1;y2]; %Valores de población iniciales (Susceptibles e Infectados)
h=input('introduce valor del paso del tiempo: ');
N=(tN-t0)/h;
y=y0;
S(1)=y(1); %Asignación del valor incial para la primera componente de S
I(1)=y(2); %Asignación del valor incial para la primera componente de I
a=0.003; % Valores de los parámetros
b=0.3; c=0.01; %Resolución empleando el método Runge Kutta
 
for n=1:N  
    k1=[-a*y(1)*y(2);a*y(1)*y(2)-(b+c)*y(2)];    
    k2=[-a*(y(1)+1/2*h*k1(1))*(y(2)+1/2*h*k1(2));a*(y(1)+1/2*h*k1(1))*(y(2)+1/2*h*k1(2)-(b+c)*(y(2)+1/2*h*k1(2)))];         
    k3=[-a*(y(1)+1/2*h*k2(1))*(y(2)+1/2*h*k2(2));a*(y(1)+1/2*h*k2(1))*(y(2)+1/2*h*k2(2)-(b+c)*(y(2)+1/2*h*k2(2)))];
    k4=[-a*(y(1)+h*k3(1))*(y(2)+h*k3(2));a*(y(1)+h*k3(1))*(y(2)+h*k3(2)-(b+c)*(y(2)+h*k3(2)))];    
    y=y+h/6*(k1+2*k2+2*k3+k4);    
    S(n+1)=y(1);    
    I(n+1)=y(2);
end
t=t0:h:tN;
figure(1)
hold on
plot(t,S,'b')
plot(t,I,'r')
hold off

for n=1:N
    A=[S(n);I(n)]+h*[-a*S(n)*I(n);a*S(n)*I(n)-(b+c)*I(n)];
    S(n+1)=A(1);%Asignacion de los distintos valores de S como la primera compnente de la matriz A
    I(n+1)=A(2); %Asignacion de los distintos valores de I como la primera componente de la matriz A
end
%Interpretación gráfica
figure(2)
hold on
plot(t,S,'b')
plot(t,I,'r')
hold off

Texto de la leyenda

La principal dificultad a la hora de emplear un método implícito radica en que estos métodos, al contrario que los explícitos no emplean la información del punto anterior para calcular el siguiente y requieren por tanto algoritmos mas complejos. Sin embargo presentan mejores aproximaciones que los métodos explícitos.

7 Tasa de infectados por contagio dependiente del tiempo. Método de Heun

clear all
%Datos
t0 = 0;
tN = 40;
h = 0.1;
N = round((tN-t0)/h); 
t = t0:h:tN;

S = zeros(1,N+1);
I = zeros(1,N+1);
S(1) = 1600;
I(1) = 40;
%Valores de los parámetros

a = zeros(1,N+1);
a(1) = 0.003/(1+t0);
b = 0.3;
c = 0.01;
 
%Calcular los vectores infectados y susceptibles.

for i = 1:N
    a(i+1) = 0.003/(1+t(i));
 
    K1 = -a(i)*S(i)*I(i);  
    K2 = -a(i)*S(i)*I(i) + K1*h;
    S(i+1) = S(i) + h/2*(K1+K2); 
    Z1 = (a(i)*S(i)*I(i)) -I(i)*(b+c);
    Z2 = (a(i)*S(i)*I(i)) -I(i)*(b+c) + Z1*h;
    I(i+1) = I(i) + h/2*(Z1+Z2);
 
end
 
%Gráficas
hold on
plot (t,S,'-')
plot (t,I,'-r')
hold off

Con el paso del tiempo el valor de la variable "a" disminuye considerablemente, puesto que en este apartado la "a" varia segun el tiempo, por esta razon, se ve amortiguada la caida de gente susceptible en la funcion.

left

8 Calibración del parámetro a

clear all
%Datos

t0 = 0;
tN = 8;
h = 0.0001;

N = (tN-t0)/h;
t = t0:h:tN;

S = zeros(1,N+1);
I = zeros(1,N+1);

S(1) = 1600;
I(1) = 40;


%Valores de los parámetros

a = 0.0005:h:0.0020;
b = 0.3;
c = 0.01;

%Definimos los vectores M, T y D, donde guardaremos respectivamente
%los máximos de infectados, el tiempo en el que se producen y la 
%diferencia entre ese tiempo y los 5 días  

M=zeros(1,length(a));
T=zeros(1,length(a));
D=zeros(1,length(a));

%Calcular los vectores infectados y susceptibles

for j=1:length(a)
A=a(j);
  for i = 1:N
    K1 = -A*S(i)*I(i);  
    K2 = -A*S(i)*I(i) + K1*h;
    S(i+1) = S(i) + h/2*(K1+K2); 
    Z1 = (A*S(i)*I(i)) -I(i)*(b+c);
    Z2 = (A*S(i)*I(i)) -I(i)*(b+c) + Z1*h; 
    I(i+1) = I(i) + h/2*(Z1+Z2);
  end
M(j)=max(I)
posicion=find(I==M(j));
tiempo=t(posicion);
T(j)=tiempo
D(j)=abs(T(j)-5);
end
[M']
[D']
[T']
T1=min(D);
posicion=find(D==T1);

%soluciones obtenidas

A_buscado=a(posicion)
Maximo_en_ese_caso=M(posicion)

Cuando calibramos el coeficiente de a con una experiencia anterior en la que el máximo de personas infectadas se alcanza a los 5 días obtenemos el resultado de que a=0.0008 es el valor que mas se ajusta a esa experiencia, tomando como paso h=10^-4.

Diferencia entre revisiones de «Modelos epidemiológicos Grupo 6C»

Revisión actual del 14:17 6 mar 2015

Contenido

1 Enunciado

2 Interpretación de parámetros

3 Método de Euler y Trapecio con S=0

4 Método de Euler y Trapecio con S=100

5 Modelo completo

6 Comparación Runge-Kutta con Euler

7 Tasa de infectados por contagio dependiente del tiempo. Método de Heun

8 Calibración del parámetro a

Menú de navegación

Herramientas personales

Espacios de nombres

Variantes

Vistas

Más

Buscar

Navegación

Herramientas

@@ Línea 1: / Línea 1: @@
 {{ TrabajoED | Modelos epidemiológicos. Grupo 6-C | [[:Categoría:Ecuaciones Diferenciales|Ecuaciones Diferenciales]]|[[:Categoría:ED14/15|Curso 2014-15]] |
 Manuel Morales López 1175
-David Toledo Menéndez 1228
-Sergio Rodríguez Torcal 994
+David Toledo Menéndez 1228
-Jose María Rodríguez Vicente 1213
-Lourdes Sánchez-Ocaña Merino 1248
+Sergio Rodríguez Torcal 994
+Jose María Rodríguez Vicente 1213
+Lourdes Sánchez-Ocaña Merino 1248
 Jorge Villa Lobo 1237 }}
 ==Enunciado==
@@ Línea 365: / Línea 370: @@
+La principal dificultad a la hora de emplear un método implícito radica en que estos métodos, al contrario que los explícitos no emplean la información del punto anterior para calcular el siguiente y requieren por tanto algoritmos mas complejos. Sin embargo presentan mejores aproximaciones que los métodos explícitos.
 ==Tasa de infectados por contagio dependiente del tiempo. Método de Heun==
-Aquí se suele poner algo...
 {{matlab|codigo=
@@ Línea 379: / Línea 383: @@
 N = round((tN-t0)/h);
 t = t0:h:tN;
-S(1) = 1600;
-I(1) = 40;
 S = zeros(1,N+1);
 I = zeros(1,N+1);
+S(1) = 1600;
+I(1) = 40;
 %Valores de los parámetros
@@ Línea 399: / Línea 404: @@
      S(i+1) = S(i) + h/2*(K1+K2);
      Z1 = (a(i)*S(i)*I(i)) -I(i)*(b+c);
-     Z2 = (a(i)*S(i)*I(i)) -I(i)*(b+c) + L1*h;
+     Z2 = (a(i)*S(i)*I(i)) -I(i)*(b+c) + Z1*h;
      I(i+1) = I(i) + h/2*(Z1+Z2);
@@ Línea 408: / Línea 413: @@
 plot (t,S,'-')
 plot (t,I,'-r')
-hold off
+hold off}}
-}}
+Con el paso del tiempo el valor de la variable "a" disminuye considerablemente, puesto que en este apartado la "a" varia segun el tiempo, por esta razon, se ve amortiguada la caida de gente susceptible en la funcion.
+[[Archivo:archivoapartado62015.jpg|360px|thumb|left]]
@@ Línea 415: / Línea 449: @@
+{{matlab|codigo=
+clear all
+%Datos
+t0 = 0;
+tN = 8;
+h = 0.0001;
+N = (tN-t0)/h;
+t = t0:h:tN;
+S = zeros(1,N+1);
+I = zeros(1,N+1);
+S(1) = 1600;
+I(1) = 40;
+%Valores de los parámetros
+a = 0.0005:h:0.0020;
+b = 0.3;
+c = 0.01;
+%Definimos los vectores M, T y D, donde guardaremos respectivamente
+%los máximos de infectados, el tiempo en el que se producen y la
+%diferencia entre ese tiempo y los 5 días
+M=zeros(1,length(a));
+T=zeros(1,length(a));
+D=zeros(1,length(a));
+%Calcular los vectores infectados y susceptibles
+for j=1:length(a)
+A=a(j);
+  for i = 1:N
+    K1 = -A*S(i)*I(i);
+    K2 = -A*S(i)*I(i) + K1*h;
+    S(i+1) = S(i) + h/2*(K1+K2);
+    Z1 = (A*S(i)*I(i)) -I(i)*(b+c);
+    Z2 = (A*S(i)*I(i)) -I(i)*(b+c) + Z1*h;
+    I(i+1) = I(i) + h/2*(Z1+Z2);
+  end
+M(j)=max(I)
+posicion=find(I==M(j));
+tiempo=t(posicion);
+T(j)=tiempo
+D(j)=abs(T(j)-5);
+end
+[M']
+[D']
+[T']
+T1=min(D);
+posicion=find(D==T1);
+%soluciones obtenidas
+A_buscado=a(posicion)
+Maximo_en_ese_caso=M(posicion)
+}}
+Cuando calibramos el coeficiente de a con una experiencia anterior en la que el máximo de personas infectadas se alcanza a los 5 días obtenemos el resultado de que a=0.0008 es el valor que mas se ajusta a esa experiencia, tomando como paso h=10^-4.