Ecuaciones Lotka-Volterra (Grupo 21)

1 Introducción

Este artículo muestra la resolución del modelo Lotka-Volterra, también conocido como ecuaciones predador-presa, que se usan para el modelado de dos poblaciones que interactuan, una presa y un depredador, usando diferentes métodos en MATLAB (Euler, Runge-Kutta y el método del trapecio). Las ecuaciones fueron propuestas de forma independiente por Alfred Lotka en 1925 y Vito volterra en 1926.

2 Premisas

En nuestro caso la presa es la población de conejos y el depredador la población de zorros, que interactuarán entre ellos en un ecosistema que supondremos aislado, de forma que los conejos sólo serán cazados por los zorros y los zorros sólo se alimentarán de conejos. El modelo Lokta-Volterra parte de unas premisas que aíslan a las dos especies y las interrelacionan entre sí.

La primera premisa expone que si no hay zorros los conejos se reproducen maltusianamente:

[math]R´(t)=αR-βR[/math]

Siendo R´(t): la variación de la población de conejos, R: la población de conejos, α: la constante de natalidad de los conejos y β: la constante de mortalidad de los conejos. Al no haber presencia de zorros los conejos no serán cazados y se multiplicarán exponencialmente siendo la solución:

[math]R(t)=Ke^α[/math]

La segunda premisa dice que si no hubiese conejos que comer la población de zorros disminuiría siguiendo también el modelo de Malthus:

[math]F´=-ϒF+δF[/math]

Igual que en el anterior F: es la población de zorros, F´: es la variación de la población de zorros, ϒ: la constante de mortalidad y δ: la constante de natalidad. Siendo la natalidad nula la población de zorros disminuiría exponencialmente:

[math]F(t)=Ce^{-ϒ}[/math]

Las relaciones entre las dos especies se definen en las premisas 3 y 4 de forma que:

La tercera premisa introduce la relación de los conejos con los zorros, siendo los conejos que son comidos (que mueren) proporcionales a la relación entre conejos y zorros:

[math]βR=bFR[/math]

Siendo c: la constante de proporcionalidad de la relación entre conejos y zorros

Y la cuarta premisa expresa la natalidad de los zorros en función de los conejos que son comidos por los zorros:

[math]δF=dFR[/math]

Definiendo ‘d’ la constante de natalidad de los zorros. De tal forma que nos queda un sistema de ecuaciones:

Igualando:

[math]α=a ; ϒ=c[/math]

Donde a: constante de natalidad de los conejos; b: constante de mortalidad de los conejos; c: constante de mortalidad de los zorros; d: constante de natalidad de los zorros.

3 Método de Euler

Usando el método de Euler y resolviendo (con MATLAB) obtenemos los siguientes gráficas:

En dicho gráfico se puede observar que cuando la población de depredadores (zorros) disminuye, la de presas (conejos) aumentará, y viceversa. En el gráfico se ve como en ningún momento se llegan a igualar las poblaciones, y por lo tanto el ecosistema formado por conejos y zorros se puede considerar estable.

En estas gráficas se observa que asignando los valores propuestos, existirán varios puntos en los cuales las poblaciones de depredador y presa, es decir, zorros y conejos, se igualarán e incluso podemos afirmar que habrá períodos parciales en los cuales la cantidad de zorros será mayor que la de conejos, o que provoca el posterior decrecimiento de la población de zorros.

Si los valores hubieran sido más extremos, podría haber llegado a producirse la extinción de zorros después de un gran decrecimiento de la población de conejos.

3.1 Programa Matlab

a=0.4;
b=0.37;
c=0.3;
d=0.05;
t0=input('Introduzca el valor del tiempo inicial: ');
t1=input('Introduzca el valor deltiempo final: ');        
N=t1*10;           
h=(t1-t0)/N;    
C=zeros(N+1);
Z=zeros(N+1);
C(1)=input('Introduzca el valor de la población inicial de conejos en miles: ');
Z(1)=input('Introduzca el valor de la población inicial de zorros en miles: ');
for n=1:N
C(n+1)=C(n)+ h*(a*C(n)-b*C(n)*Z(n));
Z(n+1)=Z(n)+h*(-c*Z(n)+d*C(n)*Z(n));
end
hold on
x=t0:h:t1;
figure(1)
plot(x,C,'g')
plot(x,Z,'r')
xlabel('tiempo')
ylabel('Número de conejos/zorros')
figure(2)
plot(C,Z,'b')
xlabel('Número de conejos')
ylabel('Número de zorros')
hold off

4 Método Runge-Kutta

Ahora aplicando el método de Runge-Kutta obtenemos:

Comparando los resultados obtenidos en la gráfica mediante ambos métodos numéricos, sea Euler y Runge-Kutta, se observa que las soluciones a través del método Runge-Kutta son más exactos y precisos, pues la diferencia entre máximos es menor, lo que significa que se ajusta más en los resultados a un intervalo mas aproximado.

En el método de Euler la precisión depende de la variable h, cuanto menor sea con respecto al intervalo más convergente es el sistema y por consiguiente más preciso es este método. Sin embargo con el método de Runge-Kutta se pueden obtener buenas aproximaciones a pesar de darle grandes valores a h.

Respecto a la dificultad que se presenta a la hora de usar el método del trapecio para resolver este sistema es que a este método, a pesar de ser más exacto, le surge el inconveniente y la dificultad de despejar Yn+1. También debe tenerse en cuenta que Euler y Runge-Kutta son métodos explícitos, sin embargo el método del trapecio como el método teta, son implícitos y no lineales, y para resolverlos tendríamos que convertirlos en explícitos de la manera ya comentada que tan laboriosa resulta a veces.

4.1 Programa Matlab

t0=input('Introduzca el valor del tiempo inicial: '); 
t1=input('Introduzca el valor deltiempo final: ');
N=t1*10; 
h=(t1-t0)/N;
C=zeros(N+1);
Z=zeros(N+1);
C(1)=input('Introduzca el valor de la población inicial de conejos en miles: ');
Z(1)=input('Introduzca el valor de la población inicial de zorros en miles: ');
a=0.4; 
b=0.37; 
c=0.3; 
d=0.05;
for n=1:N
k1=a*C(n)-b*C(n)*Z(n);
m1=-c*Z(n)+d*C(n)*Z(n);
k2=a*(C(n)+0.5*k1)-b*(C(n)+0.5*k1)*Z(n);
m2=-c*(Z(n)+0.5*m1)+d*C(n)*(Z(n)+0.5*m1);
k3=a*(C(n)+0.5*k2)-b*(C(n)+0.5*k2)*Z(n);
m3=-c*(Z(n)+0.5*m2)+d*C(n)*(Z(n)+0.5*m2);
k4=a*(C(n)+k3)-b*(C(n)+k3)*Z(n);
m4=-c*(Z(n)+m3)+d*C(n)*(Z(n)+m3);
C(n+1)=C(n)+h/6*(k1+2*k2+2*k3+k4);
Z(n+1)=Z(n)+h/6*(m1+2*m2+2*m3+m4);
end
hold on
x=t0:h:t1;
figure(1)
plot(x,C,'g')
plot(x,Z,'r')
xlabel('tiempo')
ylabel('Número de conejos/zorros')
figure(2)
plot(C,Z,'b')
xlabel('Número de conejos')
ylabel('Número de zorros')
hold off

5 Resolución del Problema inverso

También podemos obtener los datos de las poblaciones que había hace un tiempo determinado. Usando matlab obtenemos:

[math] C=3.577,4 conejos // Z=230,4 zorros [/math]

Como se observa en las siguientes gráficas:

Descargas

Descargas

5.1 Programa Matlab

a=0.4;
b=0.37;
c=0.3;
d=0.05;
t0=100;
t1=0;
N=t0*10; 
h=(t0-t1)/N; %valor del tiempo entre cada una de las N iteraciones
t=t0:-h:t1;
C=zeros(N+1);
Z=zeros(N+1);
C=input('Introduzca el valor de la población inicial de conejos en miles: ');
Z=input('Introduzca el valor de la población inicial de zorros en miles: ');
I=[C,Z];
C=I(1);  
Z=I(2);
for n=1:N
k1=[a*I(1)-c*I(1)*I(2), -b*I(2)+d*I(2)*I(1)];
k2=[a*(I(1)+1/2*h*k1(1))-c*(I(1)+1/2*h*k1(1))*(I(2)+1/2*h*k1(2)),-b*(I(2)+1/2*h*k1(2))+d*(I(2)+1/2*h*k1(2))*(I(1)+1/2*h*k1(1))];
k3=[a*(I(1)+1/2*h*k2(1))-c*(I(1)+1/2*h*k2(1))*(I(2)+1/2*h*k2(2)),-b*(I(2)+1/2*h*k2(2))+d*(I(2)+1/2*h*k2(2))*(I(1)+1/2*h*k2(1))];
k4=[a*(I(1)+h*k3(1))-c*(I(1)+h*k3(1))*(I(2)+h*k3(2)),-b*(I(2)+h*k3(2))+d*(I(2)+h*k3(2))*(I(1)+h*k3(1))];
I=I+h*(k1+2*k2+2*k3+k4)/6;   
C(n+1)=I(1);  
Z(n+1)=I(2);
end
hold on
figure(1)
plot(t,C,'g')
plot(t,Z,'b')
figure(2)
plot(C,Z,'r')
hold off