Search code examples
rplotdifferential-equations

how to solve and plot a differential equation in R?

I want to solve and graph a differential equation for exponential growth but I can't quite understand how to use the deSolve library. My equation is N = N_0 * e^(rt) and the code that I tried is

library(deSolve)

## Time
t <- seq(0, 5, 1)

## Initial population
N0 <- 2

## Parameter values
r = 1

fn <- function(t, N0, r) with(r, list(N0 * exp(r*t)))

## Solving and ploting 
out <- ode(N0, t, fn, params)
plot(out, lwd=2, main="exp")

but the output that I hope is not what I want. The graphs that I want to obtain are the following:

enter image description here

I hope you can help me. Thank you

Solution

  • The model function fn should contain the derivative, the integration is then done by the solver. First order growth can of course be solved analytically, but this is not always possible for more complex models.

    library(deSolve)

## == derivative ==
fn <- function(t, N, r) {
  # dN/dt = r * N
  list(r * N)
}

r <- 1       # Parameter value
N <- 0:100   # sequence of N
t <- 0       # dummy as the derivative is not time dependent

plot(N, fn(t, N, r)[[1]], type="l")

## == integration ==
t <- seq(0, 5, .1)  # time
N0 <- 2             # initial state

## numerical solver
out <- ode(N0, t, fn, r)
plot(out, lwd=2, main="exp")

## for comparison: analytical integration
lines(t, N0*exp(r*t), lwd=2, lty="dotted", col="red")