Fjedsen - 4 months ago 24

R Question

I am working through the "Practical Guide to Ecological Modeling" by Soetaert and Herman. I arrived at chapter 6 (Download) with the example "6.7.2 Numerical Solution of a Nutrient-Algae Chemostat Model - Euler Integration". So I tried to replicate the modelling with R using the deSolve-package and the documentation. I worked through the first part and derived the following code:

`require(deSolve)`

# Setting parameter values

parameters <- c(Nin=10, pmax=1, ks=1, dil=0.24)

# Setting state variable values

state <- c(DIN=0.1, Phyto=1)

# Defining function

Lorenz<-function(t, state, parameters) {

with(as.list(c(state, parameters)),{

# Rate of change

dDIN <- -pmax * (DIN/DIN+ks) * Phyto - dil * (DIN-Nin)

dPhyto <- pmax * (DIN/DIN+ks) * Phyto - dil * Phyto

# Return the rate of change

return(list(c(dDIN, dPhyto)))

})

}

# Define time frame

times<-seq(0,20,by=0.1)

# Solving equations

out<-ode(y=state,times=times,func=Lorenz,parms=parameters,method="euler")

head(out)

# Plotting results

par(mfrow=c(1,3), oma = c(0, 0, 3, 0))

plot(out, xlab = "time", ylab = "-")

plot(out[, "Phyto"], out[, "DIN"], pch = ".")

mtext(outer = TRUE, side = 3, "Lorenz model", cex = 1.5)

So I get no error, but I know as a fact that this solution is wrong. When I use the example code of the documentation everything works. My idea was that the error is within the defining function part. It seems strange not to include a time step operation as ...[i+1]... but the documentation works without them.

Does anybody have any experience with

`deSolve`

Answer

Your error are wrong set brackets, you have to use this:

```
dDIN <- -pmax * DIN/(DIN+ks) * Phyto - dil * (DIN-Nin)
dPhyto <- pmax * DIN/(DIN+ks) * Phyto - dil * Phyto
```

I think it´s clear why, isn´t it?