Predator-prey models by Prof. Matt Miller, Department of Mathematics, University of South Carolina email: miller@math.sc.edu
We use Maple's DEtools to study solutions of the Lotka-Volterra system and its refinements as described in section 4.1 of Calculus in Context and chapter 6 of E-K. As you play with the models, keep these questions in mind:
1. What is the long term behavior of the system? 2. In the case of oscillations, what is the period (time interval from peak to peak or trough to trough), and what is the amplitude? 3. How does changing the initial conditions affect your answers to qustions 1 and 2? 4. Does the system have any steady states (equilibria)? Do these appear to be stable or unstable? 5. If there are steady states, are they in any way related to the long term behavior? We first consider problem 6 on page 190 of CIC; see also E-K pp 218-220, equations (7ab). Here h represents the hare population and u represents 60 times the lynx population (since the lynx population is numerically much smaller than the hare population, we scale it up to fit on the same graph). We are going to work with three different initial conditions.
> restart: with(plots): with(DEtools):
Warning, the name changecoords has been redefined
> rate_eqn1:= diff(h(t),t)=(0.1)*h-(0.005)*h*(1/60)*u; rate_eqn2:=diff(u(t),t)=(0.00004)*h*u-(0.04)*u; vars:= [h(t), u(t)];
> init1:=[h(0)=2000,u(0)=600]; init2:=[h(0)=2000,u(0)=1200]; init3:=[h(0)=2000, u(0)=3000]; domain := 0 .. 320;
We plot the hare and lynx populations jointly against time using the first of the given initial conditions. You should repeat this with the other initial conditions. Get a feeling for the accuracy of the computations by changing the step size, and for the long term behavior by changing the time interval. Keep a record of the results with questions 1-5 in mind!
> L:= DEplot({rate_eqn1, rate_eqn2}, vars, domain,{init1 }, stepsize=0.5, scene=[t, u], arrows=NONE): H:= DEplot({rate_eqn1, rate_eqn2}, vars, domain,{init1 }, stepsize=0.5, scene=[t, h], arrows=NONE):
> display( {L,H} , title = `Hares and 60 * Lynxes vs. time` );
>
Which graph is which? You may want to inset options such as linecolor= or thickness= to distinguish them.
Next we plot the hare and lynx populations against one another in what is called a PHASE PORTRAIT. We do this for three different initial conditions. ***Can you identify which curve goes with which initial condition? How is the independent variable t showing up in these pictures? (Hint: try it again with time interval t = 0 .. 20.) ***
> DEplot({rate_eqn1, rate_eqn2}, vars, t= 0 .. 160, {init1, init2, init3}, stepsize=0.5, scene=[h,u], title=`Hares vs. 60 * Lynxes for t = 0 .. 160`, arrows=NONE);
What is the significance of the next calculation? (Hint: try using these values of h and u as initial conditions.)
> equil:= solve( {rhs(rate_eqn1), rhs(rate_eqn2)}, {h , u });
***Discuss the answers to questions 1-5 above in light of your examination of the model.***
Next we study solutions of the Lotka-Volterra system of CIC pp 184-186; see also E-K pp 220-222 equations (8ab). In this model the prey is assumed to grow logistically in the absence of any predators. Can you see how the rate equations have been changed from the original L-V model to incorporate this assumption? This time h represents the hare (rabbit) population and u represents 100 times the lynx (fox) population. We are going to work with three different initial conditions.
> rate_eq1:= diff(h(t),t) = (0.1)*h-.00001*h^2-(0.005)*h*(1/100)*u ; rate_eq2:= diff(u(t),t) = (0.00004)*h*u-(0.04)*u ; vars:= [h(t), u(t)];
> init1:=[h(0)=2000,u(0)=500]; init2:=[h(0)=2000,u(0)=1000]; init3:=[h(0)=2000,u(0)=5000]; domain:= 0 .. 380;
First we plot the hare and lynx populations jointly against time using the first of the given initial conditions. As above you should repeat this with the other initial conditions, different time intervals, different step sizes, etc.
> L:= DEplot([rate_eq1, rate_eq2], vars, domain,{init1 },stepsize=0.5, scene=[t, u], arrows=NONE): H:= DEplot([ rate_eq1, rate_eq2], vars, domain, {init1 }, stepsize=0.5, scene=[t, h], arrows=NONE):
> display( {L,H}, title = `Hares and 100 * Lynxes vs. time` );
Next we plot the hare and lynx populations against one another in what is called a PHASE PORTRAIT. We do this for two different initial conditions. Can you identify which curve goes with which initial condition? How is the independent variable t showing up in these pictures? (Hint: try it again with time interval t = 0 .. 20, or use the option arrows = SLIM.)
> DEplot([rate_eq1, rate_eq2], vars, domain,{init2, init3 }, stepsize=0.5, scene=[h,u], title = `Hares vs. 100 * Lynxes for t = 0 .. 320`, arrows= SLIM ) ;
> equil:= solve( {rhs(rate_eq1), rhs(rate_eq2)}, {h , u});
What is the significance of this last calculation? ***Answer questions 1-5 for this model.***
Finally we study solutions of the May system (problem 5 on page 189 of CIC; see also E-K page 265 problem 32 where the predators are whales and the prey are krill). Here x is measured in units of "hectohares" (i.e., the number of hares in units of 100) and y is the number of lynxes. Choose a variety of initial conditions, time intervals, stepsizes, and so forth. ***Are there any QUALITATIVE similarities and/or differences that you notice between the May model and the two L-V models?***
> eq1:=diff( x(t), t) = 0.6 * x *(1 -(x / 10)) - 0.5 * x * y /(x + 1);
> eq2:=diff( y(t), t) = 0.1 * y * ( 1 - y / (2 * x) );
> vars:= [x(t), y(t)];
> init1:=[0, 10, 10]; init2:=[0, 10, 15]; init3:=[0, 10, 5]; domain:= 0 .. 120;
> X:= DEplot([eq1, eq2], vars, domain, {init1 }, stepsize=0.5, scene=[t,x], linecolor=blue):
> Y:= DEplot([eq1, eq2], vars, domain, {init1 }, stepsize=0.5, scene=[t,y], linecolor=green):
> display( {X, Y}, title = `May model: Rabbits/100 and Foxes vs. time` );
> DEplot([eq1, eq2], vars, domain, {init1, init2, init3}, stepsize=0.25, scene=[x,y], title = `May model phase portrait`, arrows = SLIM );
> equil:= solve( {rhs(eq1), rhs(eq2)}, {x, y});
What is the significance of the last calculation? *** In Biology 301 you may have heard the term "limit cycle". Can you explain where this term comes from? Predict what will happen if you use initial conditions close to, but not at, the equilibrium values. Test your prediction!***
***At last it's time to judge the models. What sort of field measurements would you want to have in order to choose one over another? Are there any purely mathematical features of the predictions of the models that might help? (One can critique the models on the basis of unreasonable assumptions that might go into their construction--this is a separate matter.) For example, you have surely observed that each of the three models (L-V with unbounded prey growth, L-V with bounded prey growth, and May) exhibits a different possibility for long term behavior. Why is this considered to be such an important aspect of the model? What about sensitivity to initial conditions--how does long term behaviour depend on initial conditions, and what does this mean in terms of actual observation of populations?***
> restart:
