Aim: I am trying to numerically solve a Lotka-Volterra ODE in R, using de sde.sim() function in the sde package. I would like to use the sde.sim() function in order to eventually transform this system into an SDE. So initially, I started with a simple ODE system (Lotka-Volterra model) without a noise term.

The Lotka-Volterra ODE system: $$ \left\{\begin{matrix}\frac{dx}{dt}= \alpha x -\beta xy. \\ \frac{dy}{dt}= \delta xy -\gamma y \end{matrix}\right. $$

with initial values for x = 10 and y = 10.

The parameter values for alpha, beta, delta and gamma are 1.1, 0.4, 0.1 and 0.4 respectively (mimicking this example).

Attempt to solve the problem:

d <- expression((1.1 * x[0] - 0.4 * x[0] * x[1]), (0.1 * x[0] * x[1] - 0.4 * x[1]))
s <- expression(0, 0)
X <- sde.sim(X0=c(10,10), T = 10, drift=d, sigma=s) 

However, this does not seem to generate a nice cyclic behavior of the predator and prey population.


What is going wrong in the attempt to solve this system of ODEs in sde.sim()?


I did not yet succeed to use the sde.sim() function of the SDE package, however, I succeeded to solve the system (with and without noise) using the suggestions of Chris and the diffeqr package in R.


# Lotka-Volterra Model (SDE without noise)
f <- function(u,p,t) {
  du1 = p[1]*u[1]-p[2]*u[1]*u[2]
  du2 = p[3]*u[1]*u[2]-p[4]*u[2]

g <- function(u,p,t) {

u0 = c(10.0, 10.0)
tspan <- list(0.0,100.0)
p = c(1.1, 0.4, 0.1, 0.4)
sol = diffeqr::sde.solve(f,g,u0,tspan,p=p,saveat=0.005)
udf = as.data.frame(sol$u)

plot_ly(udf, x = sol$t, y = ~V1, name = 'Prey', type = 'scatter', mode = 'lines') %>%
  add_trace(y = ~V2, name = 'Predator', mode = 'lines')

So, not the solution to the problem, but a step in the right direction.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.