lotka-volterra.Rmd

---
title: "Lotka-Volterra Example"
subtitle: "Numerical approximations"
output: html_document
---

```{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 = .5, beta = 0.5, 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 = 5, 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