forked from allisonhorst/lotka-volterra-example
-
Notifications
You must be signed in to change notification settings - Fork 0
/
lotka-volterra.Rmd
86 lines (60 loc) · 3.52 KB
/
lotka-volterra.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
---
title: "Lotka-Volterra Example"
subtitle: "Numerical approximations"
output: html_document
---
Now making a change in R.
# Lokta-Volterra Example
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
library(tidyverse)
library(deSolve)
library(kableExtra)
```
**Credit:** This code is closely based on the article [Numerically solving differential equations with R](https://rstudio-pubs-static.s3.amazonaws.com/32888_197d1a1896534397b67fb04e0d4899ae.html)
### The Lotka-Volterra equations
As described in the lecture, the Lotke-Volterra models have been used to describe predator-prey populations.
#### Prey equation:
$$\frac{dx}{dt}=\alpha x-\beta xy$$
From Wikipedia: "The prey are assumed to have an unlimited food supply and to reproduce exponentially, unless subject to predation; this exponential growth is represented in the equation above by the term $\alpha x$. The rate of predation upon the prey is assumed to be proportional to the rate at which the predators and the prey meet, this is represented above by $\beta xy$. If either x or y is zero, then there can be no predation."
#### Predator equation:
$$\frac{dy}{dt}=\delta xy - \gamma y$$
From Wikipedia: "In this equation, $\delta xy$ represents the growth of the predator population. (Note the similarity to the predation rate; however, a different constant is used, as the rate at which the predator population grows is not necessarily equal to the rate at which it consumes the prey). The term $\gamma y$ represents the loss rate of the predators due to either natural death or emigration, it leads to an exponential decay in the absence of prey.
Where:
- $x$ is prey population (e.g. rabbits)
- $y$ is predator population (e.g. wolves)
- $\alpha, \beta, \gamma, \delta$ are positive parameters
To find an approximate solution in R, we will need four things: - Parameter values - A sequence of times over which we'll approximate predator & prey populations - An initial condition (initial populations of predator & prey at t = 0) - The differential equations that need to be solved
Solving the Lotke-Volterra equation:
```{r}
# Create a sequence of times (days):
time <- seq(0, 25, by = 0.05)
# Set some parameter values (these can change - keep it in mind):
parameters <- c(alpha = .75, beta = 0.8, delta = 0.5, gamma = 1)
# Set the initial condition (prey and predator populations at t = 0).
# Recall: x = prey, y = predator
init_cond <- c(x = 7, y = 4)
# Prepare the series of differential equations as a function:
lk_equations <- function(time, init_cond, parameters) {
with(as.list(c(init_cond, parameters)), {
dxdt = alpha * x - beta * x * y
dydt = delta * x * y - gamma * y
return(list(c(dxdt, dydt)))
})
}
# Find the approximate the solution using `deSolve::ode()`:
approx_lk <- ode(y = init_cond, times = time, func = lk_equations, parms = parameters)
# Check the class:
class(approx_lk)
# We really want this to be a data frame, and we want both prey (x) and predator (y) to be in the same column -- we'll learn why in EDS 221 (tidy data)
approx_lk_df <- data.frame(approx_lk) %>%
pivot_longer(cols = c(x,y), names_to = "species", values_to = "population")
# Plot it!
ggplot(data = approx_lk_df, aes(x = time, y = population)) +
geom_line(aes(color = species))
```
## Updating parameter values
What happens as you change the different parameters (and re-run the entire code chunk)? How does that align with what you see in the graph? Some things to keep in mind:
- $\alpha$ is a growth rate for prey
- $\gamma$ is a mortality rate for predator
## End