LEVEL 1.1 EXXPONENTIAL GROWTH AND DECAY

library(deSolve)
model <- function(time, state, parameters) {
  with(as.list(c(state, parameters)), {
    dN <- r * state['N'] # exponential growth model, replace 'r' with the rate parameter
    return(list(dN))
  })
}

# Initial conditions
state <- c(N = -1)
state2 <- c(N = 1)
# Parameters
parameters <- c(r = 0.1)

# Time sequence
times <- seq(0, 60, by = 1)

# Solve the ODE
out <- ode(y = state, times = times, func = model, parms = parameters)

# Plot the output
plot(out)

Writeup EXXPONENTIAL GROWTH AND DECAY

The above Graph is an implementation of an Ordinary Differential Equation (ODE) model that demonstrates exponential growth and decay. The workflow starts with loading the deSolve package, which is a robust tool in R for numerically solving differential equations. Next, the model is defined as a function with three arguments: time, state, and parameters. In this case, the model is using the formula for exponential growth: dN = r * N, where ‘N’ is the state (or the current population size in this case), ‘r’ is the growth rate, and ‘dN’ is the change in the population size over time. This equation can also represent decay if the value of ‘r’ is negative. Two different initial states are then defined: one with N = -1, representing exponential decay (since a negative population growth rate implies a shrinking population), and the other with N = 1, representing exponential growth (since a positive population growth rate implies an increasing population). Following that, a single parameter ‘r’ is defined and set to 0.1, which is the rate of growth or decay. Depending on the initial conditions (the initial state ‘N’), this will drive exponential growth (if N is positive) or exponential decay (if N is negative). Next, a sequence of time points is created using the seq() function, ranging from 0 to 60 in steps of 1. This is the time frame over which the ODE will be solved and the population change will be tracked.The deSolve package’s ode() function is then used to solve the defined ODE model given the initial state, the time sequence, the model function, and the parameters. The results are stored in the ‘out’ variable. Finally, the results are plotted with the plot() function, showing the population size (‘N’) over time. This plot allows us to visually inspect the exponential growth or decay of the population over time. This simple workflow allows us to simulate and explore the concepts of exponential growth and decay using ODE models in R. It provides a basis for more complex simulations, where multiple interacting states and parameters might be used.

LEVEL 1.2 Logistic Growth (ODE)

# Define the logistic growth model
M1 <- function(time, state, parameters) {
  with(as.list(c(state, parameters)), {
    dN <- r * N * (1 - N/K)
    list(c(dN))
  })
}

# Set initial state and parameters
initial_state <- c(N = 1) 
parameters <- c(r = 0.1, K = 100) 

# Generate a sequence of time points
T2 <- seq(0, 100, by = 0.1)

# Solve the ODE
O1 <- ode(y = initial_state, times = times, func = M1, parms = parameters)

# Plot the solution
plot(O1, main = "Logistic Growth Model")

# Changing the parameters
parameters1 <- c(r = 0.2, K = 100) 
parameters2 <- c(r = 0.1, K = 200) 

output1 <- ode(y = initial_state, times = times, func = M1, parms = parameters1)
output2 <- ode(y = initial_state, times = times, func = M1, parms = parameters2)

# Plot the solutions with different parameters
plot(O1, main = "Logistic Growth Model with Different Parameters")
lines(output1, col = "red")
lines(output2, col = "blue")
legend("right", legend = c("r = 0.1, K = 100", "r = 0.2, K = 100", "r = 0.1, K = 200"), fill = c("black", "red", "blue"))

writeup Logistic Growth (ODE)

This R code simulates a logistic growth model using Ordinary Differential Equations (ODEs). The first section of the code defines a function ‘M1’ representing the logistic growth model. This function takes in three arguments: time, the current state (or population size ‘N’), and parameters of the model. It returns the change in population size over time ‘dN’ as a result of the logistic growth equation: r * N * (1 - N/K). Here, ‘r’ is the growth rate and ‘K’ is the carrying capacity.Next, initial conditions and parameters for the model are set. The initial population size ‘N’ is 1, the growth rate ‘r’ is 0.1, and the carrying capacity ‘K’ is 100. A sequence of time points from 0 to 100 in increments of 0.1 is then generated. Following that, the ode() function from the deSolve package is used to numerically solve the logistic growth model over the defined time points, and the output is plotted. The code then proceeds to simulate the logistic growth model with two different sets of parameters, namely: a higher growth rate of 0.2 with the same carrying capacity of 100, and the initial growth rate of 0.1 with a higher carrying capacity of 200. The solutions for these parameter sets are plotted on the same graph for comparison.

Some of the challenges I faced when completeing this excerise are as follows. First, understanding the purpose and use of the ‘with’ function might be confusing. This function is used to create a temporary environment where the objects listed can be accessed directly without using the ‘$’ operator. Second, interpreting and applying the logistic growth equation provided to be diffucut due to a lack of experience in in mathematical modeling / ecology. Third, I struggled with understanding the structure and inputs of the ode() function and how it integrates the model over time. Lastly, understanding the plot() and lines() functions, and how they are used to plot the output of the ode() function became difficult after the code stoped worked due to an incorrect syntax error. The parameters and placement of the legend might also need some getting used to.

# Load the necessary library
library(deSolve)

# Define the quadratic growth model
model <- function(time, state, parameters) {
  with(as.list(c(state, parameters)), {
    dv <- g # dv/dt = g (acceleration due to gravity)
    dx <- v # dx/dt = v (velocity)
    return(list(c(dv, dx)))
  })
}

# Set initial state and parameters
initial_state <- c(v = 0, x = 0) # At t=0, v=0 (object at rest), x=0 (start from origin)
parameters <- c(g = 9.8) # Gravity on earth's surface

# Generate a sequence of time points
times <- seq(0, 10, by = 0.1)

# Solve the ODE
output <- ode(y = initial_state, times = times, func = model, parms = parameters)

# Plot the position
plot(output[, "time"], output[, "x"], type = "l", col = "blue", xlab = "Time", ylab = "Position", main = "Quadratic behavior")

# Change the parameter and resolve
parameters <- c(g = 9.8 / 2) # Half the gravity
output2 <- ode(y = initial_state, times = times, func = model, parms = parameters)

# Add the new solution to the plot
lines(output2[, "time"], output2[, "x"], col = "red")

# Add legend
legend("topleft", legend = c("g = 9.8", "g = 9.8/2"), col = c("blue", "red"), lty = 1)

Quadratic Models Writeup

Source Code Origin: Reddit

This R script demonstrates how to model quadratic behavior using Ordinary Differential Equations (ODE) in the context of an object falling under the force of gravity. The deSolve library is used to solve the differential equations. The model is defined through a function that takes the current state of the system and its parameters, here gravity’s acceleration ‘g’, and returns the derivative of state variables velocity ‘v’ and position ‘x’. The initial state is set with the object at rest and at the origin.The time sequence generated serves as the framework on which the ODE is solved using the ode function from deSolve. The result is then plotted showing the position of the object over time. Next, the gravity parameter is halved, simulating a scenario as if the object was in a place with half the Earth’s gravity. The ODE is solved again under these new conditions. The resultant solution is added to the existing plot using the lines function, illustrating the differences in the object’s position over time under different gravitational accelerations.

This code provides a practical demonstration of how changing parameters can affect the solutions of an ODE and how these changes can be visualized on the same plot. For someone new to R, understanding the workflow of this code can offer valuable insights into numerical solution of differential equations and their visualization. The possibility of errors might occur in understanding the ‘with’ function, differentiating between position and velocity in the model function, and understanding the structure required by the ode function for the model, parameters, and initial state.

Level 1.4: Predator Prey

# Load required package
library(deSolve)

# Define the Lotka-Volterra model
lotka_volterra_model <- function(time, state, parameters) {
  with(as.list(c(state, parameters)), {
    dPrey <- growth_rate_prey * Prey - predation_rate * Prey * Predator
    dPredator <- predation_efficiency * predation_rate * Prey * Predator - death_rate_predator * Predator
    return(list(c(dPrey, dPredator)))
  })
}

# Set initial state and parameters
initial_state <- c(Prey = 20, Predator = 10)
parameters <- c(growth_rate_prey = 0.1, predation_rate = 0.02, predation_efficiency = 0.6, death_rate_predator = 0.3)

# Generate a sequence of time points
times <- seq(0, 200, by = 0.1)

# Solve the ODE
output <- ode(y = initial_state, times = times, func = lotka_volterra_model, parms = parameters)

# Convert output to data frame
output_df <- as.data.frame(output)
names(output_df) <- c("Time", "Prey", "Predator")

# Plot the populations over time
plot(output_df$Time, output_df$Prey, type = 'l', ylim = c(0, max(output_df$Prey, output_df$Predator)), ylab = "Population", xlab = "Time", col = 'blue')
lines(output_df$Time, output_df$Predator, col = 'red')
legend("topright", legend = c("Prey", "Predator"), fill = c("blue", "red"), bty = "n")

# Plot the phase plot
plot(output_df$Prey, output_df$Predator, type = 'l', xlab = "Prey population", ylab = "Predator population", col = 'darkgreen')

Predator Prey Write up

The script in question utilizes R and the ‘deSolve’ library to simulate the dynamics of predator-prey interactions as represented by the Lotka-Volterra model. This model encapsulates the cyclical relationship between two species, typically a predator and its prey. The ‘lotka_volterra_model’ function defines the system of differential equations representing the Lotka-Volterra model. It comprises of the growth rate of the prey population (‘dPrey’), which depends on the natural growth of prey and their rate of predation, and the growth rate of the predator population (‘dPredator’), which hinges on the rate of predation and the natural death rate of predators. Initial populations of both prey and predator are defined, alongside other model parameters such as growth rates, predation rates, predation efficiency, and predator death rates. A sequence of time points is then generated, spanning a suitable range for observing population dynamics. The ode function from ‘deSolve’ is used to numerically solve the system of equations across the specified time points, yielding the population sizes of prey and predator at each point in time. The output is plotted to visualize the population dynamics over time. Adjustments to the initial population values of predators and prey, or alterations to the parameters of the model, will affect these dynamics.For example, increasing the initial predator population or the predation rate might cause a quicker decline in the prey population. On the other hand, increasing the initial prey population or their growth rate can support a larger predator population.The script also includes a phase plot, which is a plot of the predator population against the prey population. This provides another perspective on the predator-prey dynamics, illustrating the cyclical nature of their relationship independent of time. In conclusion, this script provides a hands-on exploration of the Lotka-Volterra model and the dynamics of predator-prey interactions, and reveals how various parameters influence these dynamics.