Competition - Solutions
LV in discrete time
NO QUESTIONS WERE ASKED IN THE FIRST TAB. THIS TAB THEREFORE REPEATS THE FIRST TAB FROM THE EXERCISE.
The Lotka-Volterra competition model in discrete time can be written as the following system of two difference equations: \[ \begin{align} N_1(t+1)&=N_1(t)+N_1(t)r_1\bigg(1-\frac{N_1(t)+\alpha N_2(t)}{K_1}\bigg)\\ N_2(t+1)&=N_2(t)+N_2(t)r_2\bigg(1-\frac{N_2(t)+\beta N_1(t)}{K_2}\bigg) \end{align} \]
In this exercise, you will investigate the behavior of this Lotka-Volterra competition model by iterating this system of equations. To do this, we first define a function that describes how the population size of two competing species one time step in the future depends on their population sizes at the current time step.
LVcomp <- function(N, r1, r2, K1, K2, alpha, beta) {
N1.t1 <- N[1] + r1 * N[1] * (1 - (N[1] + alpha * N[2])/K1)
N2.t1 <- N[2] + r2 * N[2] * (1 - (N[2] + beta * N[1])/K2)
c(N1.t1, N2.t1)
}This user-defined function has seven arguments. The first is the
two-dimensional population vector N, then comes the
intrinsic growth rates, the carrying capacities, and the competition
coefficients for the each species. In the next code-block, all parameter
values and the initial population size vector N.0
specifying the initial densities for the two species are determined.
The following code produces the time series for the two competing
species. We use code that you have seen in a similar way in previous
exercises. First, we create an empty data frame, called N,
in which we store the time series. It has two columns, one for each
species, and t+1 rows, one row for the initial population vector, and t
rows for the future population vectors.
t <- 50 # we set the length of the time series, which determines the number of rows in the data frame
N <- data.frame(Species1=rep(NA, t+1), Species2=rep(NA, t+1))Let us define an initial population vector at time \(t=0\) and store it in the first row of the data frame.
N0 <- c(5, 20)
N[1, ] <- N0 # we store the initial population vector in the first rowNext, we use a for-loop to iterate the time series. In each iteration
of the for-loop we use the function LVcomp to produce the
population vector at the next time step and store the outcome in a new
row in the data frame N. Note, that in this example we use
the following parameters:
r1=.9, r2=.9, K1=100, K2=100, alpha=0.8, beta=0.9.
for (i in 1:t){
N[i + 1, ] <- LVcomp(N[i, ], r1=.9, r2=.9, K1=100, K2=100, alpha=0.8, beta=0.9)
} # this for-loop iterates the population dynamics t-times
head(N) # let's have a sneak peak at the data frame, just to make sure that everything went correctly## Species1 Species2
## 1 5.00000 20.00000
## 2 8.55500 33.59000
## 3 13.52680 51.33877
## 4 19.05411 68.19760
## 5 23.57928 77.19172
## 6 26.69189 78.29420Finally, let us plot the time series. To plot the columns in a data
frame as separate lines, we can use the function
matplot.
matplot(0:t, N, type = "p", pch = c(1,16), lty = 1:2, col = 1, ylim = c(0, 90), ylab="Number of individuals", xlab="Time")
legend("top", c(expression("species 1 "), expression("species 2 ")), pch = c(1,16), bty = "n")
Logistic growth in discrete time
Before we actually look at competition between two species, let us go back one step. We can use the above code to investigate the population dynamics of a single species that growths according to logistic growth in discrete time. Recall the lecture about the discrete time logistic model. All that we have to do for that is to set the initial population size of one species, say species 2, equal to zero, that is \(N_2(t=0)=0\). This works because in a population without immigration and emigration, a species that is not present at the beginning will never be present, and the model behaves exactly as a single species model.
Here, you can see the population dynamics as resulting from the above parameter values. The plot also contains a horizontal line indicating the value of the carrying capacity of species 1.
t <- 30 # we set the length of the time series, which determines the number of rows in the data frame
N <- data.frame(Species1=rep(NA, t+1), Species2=rep(NA, t+1)) # empty data frame
N0 <- c(10, 0) # starting population vector; note that we set N2=0!
N[1, ] <- N0 # we store the initial population vector in the first row
for (i in 1:t){
N[i + 1, ] <- LVcomp(N[i, ], r1=0.9, r2=0.9, K1=100, K2=100, alpha=0.8, beta=0.9)
} # this for-loop iterates the population dynamics t-times
# let's plot the time series
matplot(0:t, N, type = "p", pch =c(1,16), col = 1, ylim = c(0, 120), ylab="Number of individuals", xlab="Time")
abline(h = 100, lty = 3) # plots a horizontal line indicating the carrying capacity of species 1
legend("top", c(expression("species1 "), expression("species2")), pch =c(1,16), bty = "n")
Vary the intrinsic growth rate r1 to investigate the
effect this has on the population dynamics. For example, you could try
the values r1=0.5, 0.9, 1.5, 1.9, 2.1, 2.4, 2.5, 2.6.
What do you observe? Feel free to choose additional values to complete your understanding.
The below figure shows the time series for the eight different \(r\)-values. From these you can observe the following: (i) the carrying capacity is approached monotonically for 0<r<1, (ii) the carrying capacity is approached in an oscillating manner for 1<r<2 (alternating between over- and undershooting), (iii) for r>2 the population does not reach the carrying capacity but undergoes sustained cycles, the amplitude of these cycles increases with r, (iv) for 2<r<2.449 this limit cycle is a 2-cycle, which means that the same population density is reached every second time step, (v) for 2,45<r<2.57 the limit cycle is a 4-cycle (by looking at more \(r\)-values one could observe that the length of the cycle doubles with increasing r to become cycle of length 8, 16, 32 and so forth), (vi) this regular cyclic pattern stops once r>2.57, from then on the population dynamics show deterministic chaos which means that the dynamics never visit the same value more than once.
Does the behavior you observe change when you alter the initial population size at time zero? Investigate this in particular for \(r=2.6\). How do you interpret your observation?
The long-term outcome of the population dynamics is independent of the initial population size. In other words, the exact same cyclic pattern is reached for different initial population sizes. This is true as long as the population dynamics is not chaotic. Once r is so large that the population dynamics shows deterministic chaos, small changes in the initial population size result in completely different time series, as shown in the figure below.
length <- 30
N1 <- data.frame(Species1=rep(NA, t+1), Species2=rep(NA, t+1))
N0 <- c(10, 0)
N1[1, ] <- N0
for (i in 1:t){
N1[i + 1, ] <- LVcomp(N1[i, ], r1=2.6, r2=0.9, K1=100, K2=100, alpha=0.8, beta=0.9)
}
length <- 30
N2 <- data.frame(Species1=rep(NA, t+1), Species2=rep(NA, t+1))
N0 <- c(50, 0)
N2[1, ] <- N0
for (i in 1:t){
N2[i + 1, ] <- LVcomp(N2[i, ], r1=2.6, r2=0.9, K1=100, K2=100, alpha=0.8, beta=0.9)
}
length <- 30
N3 <- data.frame(Species1=rep(NA, t+1), Species2=rep(NA, t+1))
N0 <- c(90, 0)
N3[1, ] <- N0
for (i in 1:t){
N3[i + 1, ] <- LVcomp(N3[i, ], r1=2.6, r2=0.9, K1=100, K2=100, alpha=0.8, beta=0.9)
}
plot(0:t, N1$Species1, type = "l", pch =c(1,16), col = 1, ylim = c(0, 130), ylab="N", xlab="Time", main="r=2.6")
abline(h = 100, lty = 3)
lines(0:t, N2$Species1, type = "l", pch =c(1,16), col = 2)
lines(0:t, N3$Species1, type = "l", pch =c(1,16), col = 3)
legend("bottom", c(expression("N0=10"), expression("N0=50"), expression("N0=90")), col =c(1, 2, 3), lty=c(1,1,1), bty = "n")
Competition between two species
In order to properly investigate the competition between two species, both need to have a non-zero density at t=0.
r1 <- 1 # intrinsic growth rate of species 1
r2 <- 1.5 # intrinsic growth rate of species 2
K1 <- 100 # carrying capacity of species 1
K2 <- 100 # carrying capacity of species 2
alpha <- 0.8 # competitive effect of species 2 on species 1
beta <- 0.9 # competitive effect of species 1 on species 2
N0 <- c(10, 5) # vector of initial population sizes, N1(0) and N2(0)
t <- 30 # number of time steps the simulation should runFor the remainder of this exercise, investigate the effect of the various parameters on the outcome of competition between two species. Specifically, investigate the following questions.
- First focus on the case that both carrying capacities are equal,
K1=K2, and that the intrinsic growth rates are relatively small, say,r1, r2 <1. Under this scenario, look at the following cases: \[ \begin{align} \mathrm{i)\,\,} &\alpha<1<\beta\\ \mathrm{ii)\,\,} &\alpha>1>\beta\\ \mathrm{iii)\,\,} &1<\alpha, \beta\\ \mathrm{iv)\,\,} &\alpha, \beta<1 \end{align} \] In summary, all cases behave exactly as predicted based on the theoretical results you have seen in the lecture. But let us nevertheless go through these case one by one. Let us first look at the case \(\alpha<1<\beta\). The simulation below shows that in this case species 1 wins and species 2 goes extint. This is because species 2 suffers heavily from competition with species 1 while species 1 suffers to a lesser extend from competition with species 2.
r1 <- 0.9 # intrinsic growth rate of species 1
r2 <- 0.9 # intrinsic growth rate of species 2
K1 <- 100 # carrying capacity of species 1
K2 <- 100 # carrying capacity of species 2
alpha <- 0.9 # competitive effect of species 2 on species 1
beta <- 1.1 # competitive effect of species 1 on species 2
N0 <- c(10, 5) # vector of initial population sizes, N1(0) and N2(0)
t <- 30 # number of time steps the simulation should run
N <- matrix(NA, nrow = t + 1, ncol = 2) # defines a matrix in which the population vectors of the time series will be stored
N[1, ] <- N0 # enters the initial population vector as a column into this matrix
for (i in 1:t) N[i + 1, ] <- LVcomp(N[i, ], r1, r2, K1, K2, alpha, beta) # this for-loop iterates the population dynamics t-times
matplot(0:t, N, type = "p", pch = c(1,16), col = 1, ylim = c(0, 120), ylab="Number of individuals", xlab="Time")
abline(h = K1, lty = 3) # plots a line indicating the carrying capacity of species 1
text(0, K1+1, "K1", adj = c(0, 0))
legend("top", c(expression("species1 wins"), expression("species2 goes extinct")), lty = 1:2, bty = "n")
In the reverse case, \(\alpha>1>\beta\), the simulations show that species 2 wins and species 1 goes extinct.
alpha <- 1.2 # competitive effect of species 2 on species 1
beta <- 0.9 # competitive effect of species 1 on species 2
N0 <- c(10, 5) # vector of initial population sizes, N1(0) and N2(0)
t <- 30 # number of time steps the simulation should run
N <- matrix(NA, nrow = t + 1, ncol = 2) # defines a matrix in which the population vectors of the time series will be stored
N[1, ] <- N0 # enters the initial population vector as a column into this matrix
for (i in 1:t) N[i + 1, ] <- LVcomp(N[i, ], r1, r2, K1, K2, alpha, beta) # this for-loop iterates the population dynamics t-times
matplot(0:t, N, type = "p", pch = c(1,16), col = 1, ylim = c(0, 120), ylab="Number of individuals", xlab="Time")
abline(h = K1, lty = 3) # plots a line indicating the carrying capacity of species 1
text(0, K1+1, "K1", adj = c(0, 0))
legend("top", c(expression("species1 goes extinct"), expression("species2 wins")), pch = c(1,16), bty = "n")
Next, let us look at the case that both competition coefficients are larger than 1, \(\alpha, \beta>1\). In this case, the simulations show that one species goes extinct and the other one reaches its carrying capacity. Interestingly, which one goes extinct and which one wins depends on the initial population densities. Below, two simulations are shown. In the first one, species 1 is more common and it turns out to be the winner of the competition. In the second one, species 2 is more common and it turns out to be the winner of the competition.
N0 <- c(2, 1) # vector of initial population sizes, N1(0) and N2(0)
alpha <- 1.1 # competitive effect of species 2 on species 1
beta <- 1.1 # competitive effect of species 1 on species 2
N <- matrix(NA, nrow = t + 1, ncol = 2) # defines a matrix in which the population vectors of the time series will be stored
N[1, ] <- N0 # enters the initial population vector as a column into this matrix
for (i in 1:t) N[i + 1, ] <- LVcomp(N[i, ], r1, r2, K1, K2, alpha, beta) # this for-loop interates the population dynamics t-times
matplot(0:t, N, type = "p", pch = c(1,16), col = 1, ylim = c(0, 120), ylab="Number of individuals", xlab="Time")
abline(h = K1, lty = 3) # plots a line indicating the carrying capacity of species 1
text(0, K1+1, "K1", adj = c(0, 0))
legend("top", c(expression("species1 wins"), expression("species2 goes extinct")), pch = c(1,16), bty = "n")
N0 <- c(1, 2) # vector of initial population sizes, N1(0) and N2(0)
alpha <- 1.1 # competitive effect of species 2 on species 1
beta <- 1.1 # competitive effect of species 1 on species 2
N <- matrix(NA, nrow = t + 1, ncol = 2) # defines a matrix in which the population vectors of the time series will be stored
N[1, ] <- N0 # enters the initial population vector as a column into this matrix
for (i in 1:t) N[i + 1, ] <- LVcomp(N[i, ], r1, r2, K1, K2, alpha, beta) # this for-loop iterates the population dynamics t-times
matplot(0:t, N, type = "p", pch = c(1,16), col = 1, ylim = c(0, 120), ylab="Number of individuals", xlab="Time")
abline(h = K1, lty = 3) # plots a line indicating the carrying capacity of species 1
text(0, K1+1, "K1", adj = c(0, 0))
legend("top", c(expression("species1 goes extinct"), expression("species2 wins")), pch = c(1,16), bty = "n")
Finally, let us look at the case that both competition coefficients are less than 1, \(\alpha, \beta < 1\). In this case, the simulations show that the two species can coexist. This is true independent of the initial population vector. Note that the sum of the density of the two species together is larger than the carrying capacity of a single species.
alpha <- 0.8 # competitive effect of species 2 on species 1
beta <- 0.9 # competitive effect of species 1 on species 2
N <- matrix(NA, nrow = t + 1, ncol = 2) # defines a matrix in which the population vectors of the time series will be stored
N[1, ] <- N0 # enters the initial population vector as a column into this matrix
for (i in 1:t) N[i + 1, ] <- LVcomp(N[i, ], r1, r2, K1, K2, alpha, beta) # this for-loop interates the population dynamics t-times
matplot(0:t, N, type = "p", pch = c(1,16), col = 1, ylim = c(0, 120), ylab="Number of individuals", xlab="Time")
abline(h = K1, lty = 3) # plots a line indicating the carrying capacity of species 1
text(0, K1+1, "K1", adj = c(0, 0))
legend("top", c(expression("species1 "), expression("species2")), pch = c(1,16), bty = "n")
- The fact that the single-species logistic equation allows for
overshooting and complicated population dynamical cycles when
r>2 carries over to the LV-competition model in discrete
time (note, that this is not possible in the continuous time versions of
both these models). Investigate this possibility for the special case
that the intrinsic growth rates and carrying capacities of both species
are equal to each other,
r1 = r2andK1 = K2. What happens if both growth rates are close to 2 but slightly less and what happens if you increase both growth rates slightly above 2? You could combine this investigation with the different configurations that the competition competition coefficients can take as listed under (1).
As an example, let us study the case that
both competition coefficients are less than 1, that is we expect the two
species to coexist. First, we run a simulation with
r1 = 1.9 =r2 and then a simulation with
r1 = 2.1 = r2. We see that in the first case the
coexistence equilibrium is reached in an oscillating manner. In the
second case, interestingly, the two species coexist but both species
fluctuate in a 2-cycle.
r1 <- 1.9 # intrinsic growth rate of species 1
r2 <- 1.9 # intrinsic growth rate of species 2
K1 <- 100 # carrying capacity of species 1
K2 <- 100 # carrying capacity of species 2
alpha <- 0.9 # competitive effect of species 2 on species 1
beta <- 0.8 # competitive effect of species 1 on species 2
N0 <- c(10, 5) # vector of initial population sizes, N1(0) and N2(0)
t <- 30 # number of time steps the simulation should run
N <- matrix(NA, nrow = t + 1, ncol = 2)
N[1, ] <- N0
for (i in 1:t) N[i + 1, ] <- LVcomp(N[i, ], r1, r2, K1, K2, alpha, beta)
matplot(0:t, N, type = "l", col = 1, ylim = c(0, 120))
abline(h = K1, lty = 3)
text(0, K1, "K1", adj = c(0, 0))
legend("top", c(expression("species1 "), expression("species2")), lty = 1:2, bty = "n")
r1 <- 2.1 # intrinsic growth rate of species 1
r2 <- 2.1 # intrinsic growth rate of species 2
K1 <- 100 # carrying capacity of species 1
K2 <- 100 # carrying capacity of species 2
alpha <- 0.9 # competitive effect of species 2 on species 1
beta <- 0.8 # competitive effect of species 1 on species 2
N0 <- c(10, 5) # vector of initial population sizes, N1(0) and N2(0)
t <- 30 # number of time steps the simulation should run
N <- matrix(NA, nrow = t + 1, ncol = 2)
N[1, ] <- N0
for (i in 1:t) N[i + 1, ] <- LVcomp(N[i, ], r1, r2, K1, K2, alpha, beta)
matplot(0:t, N, type = "l", col = 1, ylim = c(0, 120))
abline(h = K1, lty = 3)
text(0, K1, "K1", adj = c(0, 0))
legend("top", c(expression("species1 "), expression("species2")), lty = 1:2, bty = "n")