Lab 9: Competition between two species

In the logistic model with one species, the dynamical behavior of the system 
in discrete time depends primarily on r, the intrinsic rate of growth (see Box 4.1).

If r<0, the species goes extinct. 
If 0<r<1, the number of individuals increases smoothly to carrying capacity.
If 1<r<2, the number of individuals approaches carrying capacity in an oscillatory manner.
If r>2, then the carrying capacity is no longer a stable equilibrium.
	For 2 < r < 2.569944, population dynamics are oscillatory
	For 2.569944 < r <3, population dynamics are chaotic
If r > 3 the population is driven to extinction.

In this lab, we'll explore how interactions between two species affect these transitions
to answer the following question:

Does the stability of the equilibrium with both species present depend
only on the growth rates of the two species (r1 and r2) in a manner that
is similar to the dependence on r in the logistic model?

First, examine the recursions for the Lotka-Volterra model from class and
come up with your best guess for how a12 and a21 might influence the dynamics.

The Mathematica code for the discrete-time model

Here are the dynamical equations in Mathematica form:

Clear[pop1, pop2]
pop1[r1_, r2_, K1_, K2_, a12_, a21_, n10_, n20_, t_] := pop1[r1, r2, K1, K2, a12, a21, n10, n20, t] = 
   pop1[r1, r2, K1, K2, a12, a21, n10, n20, t - 1]*
    (1 + r1*(1 - (pop1[r1, r2, K1, K2, a12, a21, n10, n20, t - 1] + a12*pop2[r1, r2, K1, K2, a12, a21, n10, n20, 
           t - 1])/K1))
pop2[r1_, r2_, K1_, K2_, a12_, a21_, n10_, n20_, t_] := pop2[r1, r2, K1, K2, a12, a21, n10, n20, t] = 
   pop2[r1, r2, K1, K2, a12, a21, n10, n20, t - 1]*
    (1 + r2*(1 - (pop2[r1, r2, K1, K2, a12, a21, n10, n20, t - 1] + a21*pop1[r1, r2, K1, K2, a12, a21, n10, n20, 
           t - 1])/K2))
pop1[r1_, r2_, K1_, K2_, a12_, a21_, n10_, n20_, 0] := n10
pop2[r1_, r2_, K1_, K2_, a12_, a21_, n10_, n20_, 0] := n20

Make sure when you call either pop1 or pop2 to use the arguments in the order given:
	r1, r2, K1, K2, a12, a21, n10, n20, t
where n10 and n20 are the initial numbers of the first and second species.

To plot the number of individuals in species 1 (pop1), for example

ListPlot[Table[pop1[0.5, 0.2, 100, 100, 0.1, 0.1, 10, 10, t], {t, 0, 100}], PlotJoined -> True, 
  PlotRange -> {Automatic, {0, 150}}]

You can plot both species simultaneously by using "DisplayTogether" after loading the Graphics package
(enter the following line):

<< "Graphics`Graphics`"

Species 1 is shown with a solid curve and species 2 by dashing:

DisplayTogether[ListPlot[Table[pop1[2.1, 0.9, 100, 100, 0.9, 0.1, 10, 10, t], {t, 0, 100}], PlotJoined -> True], 
  ListPlot[Table[pop2[2.1, 0.9, 100, 100, 0.9, 0.1, 10, 10, t], {t, 0, 100}], PlotJoined -> True, 
   PlotStyle -> Dashing[{0.01}]], PlotRange -> {Automatic, {0, 150}}]

Questions

To simplify matters, assume that the carrying capacity of both species is K1 = K2 = 100.
Also, assume that both species start at initial population size n10 = n20 = 10
Note that you can set n20 = 0 to recover the dynamics for species 1 based on the logistic if you want.
You can also increase the time range by altering {t,0,100} if desired.

Question 1:
Set r1 = r2 = 1.95, which is right below the cut-off where oscillations persist with only one species.
Investigate the dynamics under with the following competition coefficients:
(a)	a12 = 0	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(b)	a12 = 0.1 	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(c)	a12 = 0.9 	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(d)	a12 = 0.9 	a21 = 0.8	Do the oscillations persist?	Do both species persist?
(e)	a12 = 1.1 	a21 = 0.8	Do the oscillations persist?	Do both species persist?
(f)	a12 = 1.1 	a21 = 1.2	Do the oscillations persist?	Do both species persist?

Conclusions:

DisplayTogether[ListPlot[Table[pop1[1.95, 1.95, 100, 100, 0, 0, 10, 10, t], {t, 0, 100}], PlotJoined -> True], 
  ListPlot[Table[pop2[1.95, 1.95, 100, 100, 0, 0, 10, 10, t], {t, 0, 100}], PlotJoined -> True, 
   PlotStyle -> Dashing[{0.01}]], PlotRange -> {Automatic, {0, 150}}]

Question 2:
Set r1 = r2 = 2.05, which is right above the cut-off where oscillations persist with only one species.
Investigate the dynamics under with the following competition coefficients:
(a)	a12 = 0	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(b)	a12 = 0.1 	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(c)	a12 = 0.9 	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(d)	a12 = 0.9 	a21 = 0.8	Do the oscillations persist?	Do both species persist?
(e)	a12 = 1.1 	a21 = 0.8	Do the oscillations persist?	Do both species persist?
(f)	a12 = 1.1 	a21 = 1.2	Do the oscillations persist?	Do both species persist?

Conclusions:

Question 3:
Set r1 = 2.05 and r2 = 0.5 (i.e. repeat Question 2 when the second species has a lower intrinsic growth rate).
Investigate the dynamics under with the following competition coefficients:
(a)	a12 = 0	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(b)	a12 = 0.1 	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(c)	a12 = 0.9 	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(d)	a12 = 0.9 	a21 = 0.8	Do the oscillations persist?	Do both species persist?
(e)	a12 = 1.1 	a21 = 0.8	Do the oscillations persist?	Do both species persist?
(f)	a12 = 1.1 	a21 = 1.2	Do the oscillations persist?	Do both species persist?

Conclusions:

Question 4:
Set r1 = 3.05 and r2 = 0.5, so that you would predict that species 1 would go extinct in the absence of species 2.
Investigate the dynamics under with the following competition coefficients:
(a)	a12 = 0	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(b)	a12 = 0.1 	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(c)	a12 = 0.9 	a21 = 0.2	Do the oscillations persist?	Do both species persist?
(d)	a12 = 0.9 	a21 = 0.8	Do the oscillations persist?	Do both species persist?
(e)	a12 = 1.1 	a21 = 0.8	Do the oscillations persist?	Do both species persist?
(f)	a12 = 1.1 	a21 = 1.2	Do the oscillations persist?	Do both species persist?

Conclusions:

Question 5:
Make up your own predictions based on the above simulations and test it with additional simulations.

Conclusions:

The continuous model

(Do this if you have time.)

From your experience with the logistic model, would you predict that the oscillations
would get more pronounced or would dampen if the model were set in continuous
time rather than in discrete time?

To answer this, choose a parameter set from above that gives oscillations and
use NDSolve to solve the system of differential equations numerically:

Clear[system, n1, n2]
system[r1_, r2_, K1_, K2_, a12_, a21_, n10_, n20_, tmin_, tmax_] := 
  system[r1, r2, K1, K2, a12, a21, ninit1, ninit2, tmin, tmax] = 
   NDSolve[{Derivative[1][n1][t] == n1[t]*r1*(1 - (n1[t] + a12*n2[t])/K1), 
     Derivative[1][n2][t] == n2[t]*r2*(1 - (n2[t] + a21*n1[t])/K2), n1[0] == n10, n2[0] == n20}, {n1[t], n2[t]}, 
    {t, tmin, tmax}]

("system" is the numerical solution to the two sets of equations.  
It contains information about both the numbers of species 1 and 2).

To plot the number of individuals in species 1 (n1), for example

Plot[Evaluate[n1[t] /. system[3.1, 3.5, 100, 100, 0.2, 0.2, 10, 10, 0, 100]], {t, 0, 100}, 
  PlotRange -> {{0, 100}, {0, 100}}]