Modelo de Lotka-Volterra

Anterior: Jogo da Vida | Índice: Ecologia | Próximo: Modelo de Lotka-Volterra amortecido

Versão tradicional

No modelo de Lotka-Volterra temos as seguintes considerações:

• Na ausência de predadores, a população de presas aumenta a uma taxa proporcional à população atual;
• Na ausência de presas, os predadores irão à extinção;
• O número de encontro entre presas e predadores é proporcional a produto das duas populações.
  • Estes encontros beneficiam os predadores em detrimento das presas.

Dessa forma, as equações são:

• $\frac{dx}{dt}=x(a-\alpha y)$
• $\frac{dy}{dt}=y(-c+\gamma x)$

Onde:

• $a:$ taxa de crescimento de presas sem predadores;
• $\alpha :$taxa de decréscimo da população de presas devido a predação;
• $c:$ taxa de mortalidade da população de predadores sem presas;
• $\gamma$: taxa de crescimento de predadores devido a predação.

Separação de variáveis

Utilizando a separação de variáveis, temos:

$$\frac{dy}{dt}\frac{dt}{dx}=\frac{y(-c+\gamma x)}{x(a-\alpha y)}$$

Logo:

$$\frac{dy}{dx}=\frac{y(-c+\gamma x)}{x(a-\alpha y)}$$$$\frac{(a-\alpha y)}{y}dy=\frac{(-c+\gamma x)}{x}dx$$$$(\frac{a}{y}-\alpha )dy=(-\frac{c}{x}+\gamma )dx$$Integrando ambos os lados:

$$a\ln y-\alpha y=-c\ln x+\gamma x+C$$$$a\ln y-\alpha y+c\ln x-\gamma x=C$$

Onde $C$ é uma constante de integração. Para plotarmos um gráfico, considerando apenas $a=\alpha =\gamma =c=1$ Temos então:

$$\ln y+\ln x-(x+y)=C$$

Um ponto de equilíbrio fora da origem é obtido quando:

$$\frac{dx}{dt}=x(a-\alpha y)=0\to y=\frac{a}{\alpha }=1$$$$\frac{dy}{dt}=y(-c+\gamma x)=0\to x=\frac{c}{\gamma }=1$$

Então neste caso, o sistema oscila em torno de $(1,1)$ e a constante $C$ é definida pelas condições iniciais $(x_0,y_0)$. Para a condição em que $x_0=y_0=1$, então:

$$\ln 1+\ln 1-(1+1)=C$$$$-2=C$$

Então para este conjunto de parâmetros e condições iniciais: $$\ln y+\ln x-(x+y)+2=0$$

Neste caso o sistema vai permanecer no ponto de equilíbrio. Para outras condições iniciais, o sistema vai oscilar em torno do ponto de equilíbrio. Obviamente além do ponto $(1,1)$, temos um ponto de equilíbrio em $(0,0)$. Vamos analisar a dinâmica na vizinhança dos pontos através de um processo simples de linearização.

$f(x,y)=\ln y+\ln x-(x+y)-C$ com as condições $a=\alpha =c=\gamma =1$ e condição inicial arbitrária, plotado no GeoGebra.

Linearização em torno do ponto de equilíbrio

Primeiro podemos perceber que o sistema é quase-linear em torno de $(0,0)$, verificando que satisfaz:

$$\underset{(x,y)\to (0,0)}{\lim }\left[\frac{\text{parte não linear}}{\text{parte linear}}\right]=0$$

Então lembrando as equações:

• $\frac{dx}{dt}=\left[xa\right]-\left(\alpha xy\right)=\left[\text{linear}\right]-\left(\text{não linear}\right)$
• $\frac{dy}{dt}=-\left[yc\right]+\left(\gamma yx\right)=-\left[\text{linear}\right]+\left(\text{não linear}\right)$

Logo:

$$\underset{(x,y)\to (0,0)}{\lim }-\frac{\alpha xy}{xa}=\underset{(x,y)\to (0,0)}{\lim }-\frac{\alpha }{a}y=0$$$$\underset{(x,y)\to (0,0)}{\lim }-\frac{\gamma xy}{cy}=\underset{(x,y)\to (0,0)}{\lim }-\frac{\gamma }{c}x=0$$

Desprezando então os termos não lineares podemos escrever o seguinte sistema linearizado em torno da origem:

$$\left(\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right)=\left(\begin{array}{cc}a& 0\\ 0& -c\end{array}\right)\left(\begin{array}{c}x\\ y\end{array}\right)$$Calculando os autovalores da matriz, obtemos então:

$$-(a-\lambda )(-c-\lambda )=0$$$$(a-\lambda )(c+\lambda )=0$$

os seguintes autovalores $\lambda =\{a,-c\}$. Como um dos valores tem parte real positiva, então é um ponto instável, especificamente devido aos sinais opostos é um ponto de sela. Como é instável significa que se a condição inicial for próxima de $(0,0)$, a evolução do sistema vai se afastar do ponto de equilíbrio. Essa aproximação também indica que próximo do ponto de equilíbrio, a dinâmica pode ser descrita tanto pelo conjunto de equações não lineares, como pelo sistema linear.

Agora o segundo ponto de equilíbrio, de maneira geral é $(\frac{c}{\gamma },\frac{a}{\alpha })$. Primeiro reescrevemos o sistema em torno do ponto de equilíbrio, isto é, fazemos um deslocamento $u=x-\frac{c}{\gamma }$ e $v=y-\frac{a}{\alpha }$. Então temos $dx=du$ e $dv=dy$ e substituindo, para $\dot{x}$:$$\frac{du}{dt}=(u+\frac{c}{\gamma })a-\alpha (u+\frac{c}{\gamma })(v+\frac{a}{\alpha })$$$$\frac{du}{dt}=ua+\frac{c}{\gamma }a-\alpha uv-\frac{\alpha c}{\gamma }v-ua-\frac{ca}{\gamma }$$$$\frac{du}{dt}=-\alpha uv-\frac{\alpha c}{\gamma }v$$E para $\dot{y}$:

$$\frac{dv}{dt}=-(v+\frac{a}{\alpha })c+\gamma (v+\frac{a}{\alpha })(u+\frac{c}{\gamma })$$$$\frac{dv}{dt}=-cv-\frac{ca}{\alpha }+\gamma vu+\frac{\gamma a}{\alpha }u+cv+\frac{ca}{\alpha }$$$$\frac{dv}{dt}=\gamma vu+\frac{\gamma a}{\alpha }u$$Podemos analisar o comportamento em torno do ponto de equilíbrio:

$$\underset{(u,v)\to (0,0)}{\lim }-\frac{\alpha uv}{\frac{\alpha vc}{\gamma }}=\underset{(u,v)\to (0,0)}{\lim }-\frac{\gamma }{c}u=0$$$$\underset{(u,v)\to (0,0)}{\lim }\frac{\gamma vu}{\frac{\gamma au}{\alpha }}=\underset{(u,v)\to (0,0)}{\lim }\frac{\alpha }{a}v=0$$Desprezando os termos não lineares então:

$$\left(\begin{array}{c}\dot{u}\\ \dot{v}\end{array}\right)=\left(\begin{array}{cc}0& -\frac{\alpha c}{\gamma }\\ \frac{\gamma a}{\alpha }& 0\end{array}\right)\left(\begin{array}{c}u\\ v\end{array}\right)$$Então os autovalores correspondentes:

$$-{\lambda }^2-\frac{\gamma a}{\alpha }\frac{\alpha c}{\gamma }=0$$$$\lambda =\pm \sqrt{-ac}=\pm \sqrt{ac}i$$

Como temos raízes puramente imaginárias e ${\lambda }_1={\lambda }_2^{*}$, temos um centro, ponto de estabilidade. Isto é, se a condição inicial for próxima de $(\frac{c}{\gamma },\frac{a}{\alpha })$ o sistema evoluirá de forma que o estado do sistema permanecerá próximo do ponto de equilíbrio.

Classificação dos pontos de estabilidade de acordo com os autovalores.

Segundo método de Lyapunov

Para avaliar o ponto $(x_1,y_1)=(0,0)$, podemos usar de maneira análoga ao exemplo do segundo critério de Lyapunov:

$$V\left(\mathit{x}\right)=\frac{x^2}{\alpha }-\frac{y^2}{\gamma }$$

Como já discutimos $V\left({\mathit{x}}_0\right)=0$ e a região $W^{+}\{(x,y)|\left|x\right|>\left|y\right|\}$ onde $V\left(\mathit{x}\right)>0$ para $\mathit{x}\ne {\mathit{x}}_0$ , sendo ${\mathit{x}}_0$ um ponto de acumulação em $W^{+}$. Então:

$$\begin{array}{rl}\dot{V}\left(\mathit{x}\right)& =\left[\mathrm{\nabla }V\right]\cdot \left[\mathit{f}\left(\mathit{x}\right)\right]\\ & =(\frac{2x}{\alpha },-\frac{2y}{\gamma })(\dot{x},\dot{y})\\ & =2x^2\frac{a}{\alpha }-2x^{2}y+2y^2\frac{c}{\gamma }-2y^{2}x\\ & =2x^2(\frac{a}{\alpha }-y)+2y^2(\frac{c}{\gamma }-x)\end{array}$$Lembrando do nosso segundo ponto de equilíbrio $(x_2,y_2)=(\frac{c}{\gamma },\frac{a}{\alpha })$:

$$\dot{V}\left(\mathit{x}\right)=2x^2(y_2-y)+2y^2(x_2-x)$$

Então se estamos próximos suficiente do ponto de equilíbrio em análise $(x_1,y_1)=(0,0)$, temos então uma instabilidade local pois $\dot{V}\left(\mathit{x}\right)>0$ é positivo definido em $W^{+}$, uma vez que $\left|y\right|<\left|y_2\right|$ , $\left|x\right|<\left|x_2\right|$. Olhando o segundo ponto de equilíbrio, $(x_2,y_2)=(\frac{c}{\gamma },\frac{a}{\alpha })$, podemos manipular as equações da seguinte forma:

• $\frac{dx}{dt}=x(a-\alpha y)=x\alpha (\frac{a}{\alpha }-y)=x\alpha (y_2-y)$
• $\frac{dy}{dt}=y(-c+\gamma x)=y\gamma (-\frac{c}{\gamma }+x)=y\gamma (-x_2+x)$

Definindo então a seguinte função de Lyapunov:

$$V(x,y)=x-x_2[1+\ln \left(\frac{x}{x_2}\right)]+\frac{\alpha }{\gamma }(y-y_2[1+\ln \left(\frac{y}{y_2}\right)])$$

Analisando no ponto de e equilíbrio, temos:

$$\begin{array}{rl}V(x_2,y_2)& =x_2-x_2[1+\ln \left(\frac{x_2}{x_2}\right)]+\frac{\alpha }{\gamma }(y_2-y_2[1+\ln \left(\frac{y_2}{y_2}\right)])\\ & =x_2-x_2+\frac{\alpha }{\gamma }(y_2-y_2)\\ & =0\end{array}$$

Agora precisamos que para $(x,y)\ne 0$ tenhamos $V>0$, na região próxima ao ponto de equilíbrio. Separando em dois termos:

$$\begin{array}{rl}V(x,y)& =[x-x_2(1+\ln \left(\frac{x}{x_2}\right))]+\frac{\alpha }{\gamma }[y-y_2(1+\ln \left(\frac{y}{y_2}\right))]\\ & =V\left(x\right)+\frac{\alpha }{\gamma }V\left(y\right)\end{array}$$

De forma geral temos $V\left(z\right)=z-z_2(1+\ln \left(\frac{z}{z_2}\right))$, e precisamos que $V\left(z\right)>0$ quando $z\ne z_2$. Além de ser facilmente visto via gráfico para todos os valores aceitáveis $z,z_2\in [0,1]$, também podemos analisar a seguinte desigualdade:

$$\begin{array}{rl}z-z_2(1+\ln \left(\frac{z}{z_2}\right))& >0\\ z& >z_2(1+\ln \left(\frac{z}{z_2}\right))\\ \frac{z}{z_2}& >1+\ln \left(\frac{z}{z_2}\right)\\ e^{\frac{z}{z_2}}& >e\frac{z}{z_2}\\ e^u& >eu\end{array}$$Podemos ver quer a desigualdade desigualdade é válida exceto se $x=1$. Mas como fizemos a seguinte substituição $u=\frac{z}{z_2}$ então $u=1\to z=z_2$, e de fato que queremos que seja positiva definida fora do ponto de equilíbrio. Uma vez que sabemos que $V(x,y)$ é positivo definido, calculamos então:

$$\begin{array}{rl}\dot{V}\left(\mathit{x}\right)& =\left[\mathrm{\nabla }V\right]\cdot \left[\mathit{f}\left(\mathit{x}\right)\right]\\ & =[\frac{\partial V\left(x\right)}{\partial x},\frac{\alpha }{\gamma }\frac{\partial V\left(y\right)}{\partial y}]\cdot [x\alpha (y_2-y),y\gamma (-x_2+x)]\\ & =[1-\frac{x_2}{x},\frac{\alpha }{\gamma }(1-\frac{y_2}{y})]\cdot [x\alpha (y_2-y),y\gamma (-x_2+x)]\\ & =\left(\frac{x-x_2}{x}\right)\left(x\alpha (y_2-y)\right)+\frac{\alpha }{\gamma }\left(\frac{y-y_2}{y}\right)\left(y\gamma (-x_2+x)\right)\end{array}$$

Então:$$\dot{V}\left(\mathit{x}\right)=\alpha (x-x_2)(y_2-y)-\alpha (y_2-y)(x-x_2)=0$$Temos então a condição de estabilidade $\dot{V}\le 0$ concordando como que já havíamos obtidos anteriormente.

Solução numérica

Um exemplo resolvido numericamente pode ser visto em Modelo de Lotka-Volterra amortecido, onde foi aproveitado os códigos desenvolvidos para este mesmo.

Versão adimensional

Originalmente temos 4 parâmetros, mas podemos realizar uma série de manipulações visando uma redução da quantidade de parâmetros. Escrevendo o sistema como:

$$\begin{array}{rl}\frac{dN}{dt}& =aN-bNP\\ \frac{dP}{dt}& =cPN-dP\end{array}$$

Podemos definir então $\hat{t}=at$. Multiplicando ambas equações por $1/a$:

$$\begin{array}{rl}\frac{1}{a}\frac{dN}{dt}& =N-\frac{b}{a}NP\\ \frac{1}{a}\frac{dP}{dt}& =\frac{c}{a}PN-\frac{d}{a}P\end{array}$$

Se definimos $p=(b/a)P$ e multiplicamos a segunda equação por $b/a$:

$$\begin{array}{rl}\frac{dN}{d\hat{t}}& =N-Np\\ \frac{d}{d\hat{t}}\left(\frac{bP}{a}\right)& =\frac{c}{a}\left(\frac{b}{a}P\right)N-\frac{d}{a}\left(\frac{b}{a}P\right)\end{array}$$

Agora se multiplicamos a primeira linha e definimos $n=(c/d)N$:

$$\begin{array}{rl}\frac{d}{d\hat{t}}\left(\frac{c}{d}N\right)& =\frac{c}{d}N-\left(\frac{c}{d}N\right)p\\ \frac{dp}{d\hat{t}}& =\frac{c}{a}pN-\frac{d}{a}p\end{array}$$

Definindo então $\alpha =\frac{d}{a}$:

$$\begin{array}{rl}\frac{dn}{d\hat{t}}& =n-np\\ \frac{dp}{d\hat{t}}& =\frac{c}{a}\frac{d}{c}pn-\alpha p\end{array}$$

Ou então, ficamos apenas com um único parâmetro.

$$\begin{array}{rl}\frac{dn}{d\hat{t}}& =n(1-p)\\ \frac{dp}{d\hat{t}}& =\alpha p(n-1)\end{array}$$

Separação de variáveis

Aplicando a separação de variáveis, temos então: $$\begin{array}{rl}\frac{dp}{dn}& =\frac{\alpha p(n-1)}{n(1-p)}\\ \frac{(1-p)}{p}dp& =\alpha \frac{(n-1)}{n}dn\\ (\frac{1}{p}-1)dp& =\alpha (1-\frac{1}{n})dn\\ \ln p-p+K& =\alpha (n-\ln \left(n\right))\\ K& =\alpha n+p-\alpha \ln \left(n\right)-\ln p\end{array}$$

Ou ainda, apenas:

$$K=\alpha n+p+\ln \left(n^{\alpha }p\right)$$

Referências

Principais materiais utilizados

1. A survey of constructing Lyapunov functions for mathematical models in population biology (Sze-Bi, Revista Taiwanesa de Matemática )
2. Estabilidade de pontos de equilíbrio e existência de soluções periódicas em alguns modelos bidimensionais (Salvador Tavares de Oliveira, UNESP)
3. Modelagem Matemática e estabilidade de sistemas predador-presa (Paulo Laerte Natti e outros, UEL)
4. Modelo de Lotka-Volterra: a dinâmica predador-presa (Rafael Biasi Pata e Elisa Regina Cara, UNIPAMPA)

Citações

Anterior: Jogo da Vida | Índice: Ecologia | Próximo: Modelo de Lotka-Volterra amortecido