Diferencia entre revisiones de «Modelos epidemiológicos (Grupo 19)»

Revisión del 23:10 4 mar 2013

1 Interpretación de los diferentes parámetros

Tenemos un modelo epidemiológico “SI” definido mediante las ecuaciones:

[math]\frac{\mathrm{d} S}{\mathrm{d} T}=-aSI [/math]

[math]\frac{\mathrm{d} I}{\mathrm{d} T}=aSI-bI-cI [/math]

Siendo:

• S: población susceptible de infectarse por la enfermedad

• I: población infectada por la enfermedad

• a: Parámetro de la interacción entre el número de individuos en ambas clases.
• b: Parámetro que representa el número de fallecimientos por la enfermedad.
• c: parámetro que representa el numero de curas de la enfermedad.

La primera ecuación representa la variación de la población susceptible de contraer la enfermedad en función del tiempo. Vemos que es una función decreciente debido al signo negativo que acompaña al parámetro a. Esto implica que la población susceptible de este modelo siempre va a disminuir siendo su máximo el valor inicial [math]S_0[/math] . La función es dependiente tanto de la población susceptible como de la infectada. Al estar multiplicadas ambas variables se convierte en una ecuación no lineal, factor determinante a la hora de resolver el sistema.

La segunda ecuación representa la variación de la población infectada con respecto del tiempo. Podemos estudiarla en dos partes dividiéndola en una no lineal (análoga a la primera ecuación) y la parte lineal que tan solo depende de los infectados (I). La primera parte es dependiente de nuevo del parámetro a con signo positivo lo que hace crecer el número de infectados. Por el contrario, la parte lineal varía con los parámetros b y c. En este caso el signo es negativo ajustándose a la realidad pues representan el número de muertes y de curas.

El número de infectados vemos que aumentará o disminuirá según la relación que exista entre los parámetros a, b y c y el número variable de infectados y susceptibles.

2 Resolución del sistema mediante Euler

Para resolver el sistema lo aproximamos mediante el método de Euler con los valores iniciales [math]S_0,I_O[/math] y los valores de los parámetros [math]a=0.003,b=0.3,c=0.2[/math]

1).[math]S_0=700,I_0=1[/math]

clear all
t0=0; tN=30; %Intervalo de tiempo tomado
S0=700; I0=1; %Valores iniciales de infectados(I) y susceptibles(S) 
N=300; h=(tN-t0)/N; %Determinación del paso 
a=0.003; %valores de los parametros
b=0.3;
c=0.2; 
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
t=t0:h:tN; %Desarrollo del tiempo desde t0 hasta tN tomando intervalos de paso h

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

Obtenemos diferentes aproximaciones según los pasos empleados:

centro

Comprobamos que cuanto menor es el paso a utilizar mayor es la precisión de la aproximación, de forma que el valor de la población susceptible a los 30 días toma un valor que tiende a 11.2

• [math]h=10^-1 → S=9.75[/math]
• [math]h=10^-2 → S=11.02[/math]
• [math]h=10^-3 → S=11.14[/math]
• [math]h=10^-4 → S=11.16[/math]

2.1 Interpretación

Vemos que la población susceptible va disminuyendo conforme aumenta la población infectada. Durante los dos primeros días la enfermedad presenta un crecimiento suave mientras que los días posteriores tiene un crecimiento muy rápido hasta el quinto dia, donde el número de infectados alcanza su máximo siendo de 303 personas. En este momento los susceptibles son sólo 166 personas.

Llegado este punto el número de infectados disminuye por las muertes y curas hasta la desaparición de la enfermedad ([math]I=0[/math]), quedando un valor de población susceptible de 10 personas,valor que se mantiene constante :

[math]I=0[/math]

[math]\frac{\mathrm{d} S}{\mathrm{d} T}=-aSI=0 [/math]

[math]S=cte[/math]

2).[math]S_0=5000,I_0=5[/math]

Al emplear el mismo código de matlab que en el caso anterior pero tan solo modificando las condiciones iniciales también se altera el proceso a lo largo del tiempo.

centro

Al tomar una población susceptible mucho mayor observamos que en menos de un día el número de infectados alcanza un valor máximo superior a 4000 personas, descendiendo drásticamente el número de susceptibles a 150.

Vemos que en menos de dos días toda la población susceptible queda infectada desapareciendo la enfermedad por completo antes de las dos semanas.

En ambos casos, tomando valores iniciales diferentes, las ecuaciones del sistema tienden al equilibrio y en general desaparece la enfermedad al llegar a un valor de [math]I=0[/math].

3 Interpretación con nuevos datos iniciales

Decidimos tomar unos nuevos datos iniciales [math]S_0=250,I_0=5[/math]

Para contrastar con los casos anteriores hemos tomado una población inicial menor. Podemos apreciar un cambio significativo en la evolución de la infección en la población debido a que la enfermedad se desarrolla más lentamente y el número máximo de infectados es de 21 a los 10 días.

centro

Se deduce que al haber menos población susceptible de ser infectada la enfermedad tiene más dificultad en propagarse. En este caso la relación entre la población susceptible al final y al principio es mucho mayor que en los supuestos anteriores, quedando 100 personas sin contraer la enfermedad.

4 Resolución del sistema mediante Runge-Kutta

Aproximamos el mismo sistema por medio del método Runge-Kutta, empleando los mismo datos iniciales que en los casos anteriores resueltos por el método de Euler (apartado 2). .

clear all
t0=0; tN=30; %Intervalo de tiempo tomado
y0=[700;1]; %Valores de población iniciales (Susceptibles e Infectados)
N=300; h=(tN-t0)/N; %Determinación del paso
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.2; %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

%Interpretación gráfica
t=t0:h:tN;
figure(1)
hold on
plot(t,S,'b')
plot(t,I,'r')
hold off
figure(2)
plot(S,I,'b')

1). [math]S_0=700,I_0=1[/math]

centro

2). [math]S_0=5000,I_0=5[/math]

centro

3). [math]S_0=250,I_0=5[/math]

centro

Tanto el método de Euler como el Runge-Kutta son métodos basados en iteraciones. Como Euler se basa en buscar la pendiente de la recta tangente necesita un paso muy pequeño para obtener una mejor aproximación a la solución. Por el contrario podemos obtener un buen resultado empleando un paso mayor con el método de Runge-Kutta. Si se tienen en cuenta estos factores se pueden lograr aproximaciones con gran exactitud.

Otra alternativa a la hora de resolver el sistema sería emplear un método implícito como lo es el Método Trapezoidal. Sin embargo, el inconveniente de emplear dicho método aparece a la hora de despejar las variables dependientes. Este proceso puede resultar muy trabajoso dependiendo de las ecuaciones que se planteen, como es nuestro caso.

5 Cálculo de los valores iniciales

En este caso se nos plantea calcular los valores iniciales de nuestro sistema sabiendo los valores para un tiempo determinado.

Conocemos [math]S_0=20000,I_0=3000[/math] para un [math]t=15[/math]

clear all
t0=0; tN=15; %Intervalo de tiempo tomado
S15=2; I15=0.3; %Valores iniciales de infectados(I) y susceptibles(S) 
N=1500; h=(tN-t0)/N; %Determinación del paso 
a=0.003; %valores de los parametros
b=0.3;
c=0.2; 
S(15)=S15; %Asignacion del valor incial para la primera componente
I(15)=I15;

%Resolucion del sistema de forma matricial
for n=14:-1:1
   A=[S(n+1);I(n+1)]-h*[-a*S(n+1)*I(n+1);a*S(n+1)*I(n+1)-(b+c)*I(n+1)];
   S(n)=A(1); %Asignacion de los distintos valores de S como la primera componente de la matriz A
   I(n)=A(2); %Asignacion de los distintos valores de I como la primera componente de la matriz A
end
t=t0:h:tN; %Desarrollo del tiempo desde t0 hasta tN tomando intervalos de paso h

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