[R] deSolve ODE Output Question

David Winsemius dwinsemius at comcast.net
Wed May 20 18:34:19 CEST 2015


On May 20, 2015, at 6:45 AM, walke554 wrote:

> Hello,
> 
> I am working on a simple ODE problem with the deSolve package, and I was
> hoping that someone could answer a question about how the deSolve package
> does integration.  
> 
> Here is my program:
> 
> #The function
> STDMod<-function(t,y,p){
> 	IH = y[1];
> 	IL = y[2];
> 	with(as.list(p), {
> 		dIH.dt = Beta[1,1]*(nH-IH)*IH + Beta[1,2]*(nH-IH)*IL - gamma*IH;
> 		dIL.dt = Beta[2,1]*(nL-IL)*IH + Beta[2,2]*(nL-IL)*IL - gamma*IL;
> 		return(list(c(dIH.dt,dIL.dt)));
> 	})
> }
> 
> #giving the parameters
> 
> Beta = matrix(data=c(10, 0.1, 0.1, 1.0), ncol=2, nrow=2)
> nH = 0.2
> nL = 0.8
> IH0 = 1e-5
> IL0 = 0
> 
> gamma = 1
> 
> p = list(Beta=Beta, gamma=gamma, nH=nH, nL=nL)
> 
> y0 = c(IH0, IL0)
> 
> #Running the ode integrator
> 
> steps= 10;
> t = seq(from=0, to=30, by=.01);
> out = ode(y=y0, times=t, func=STDMod, parms=p);
> 
> My understanding is that the 'out' matrix would be the values of the STDmod
> function at each time step, given the initial values of the state parameters
> IH and IL.  However, for the very first time step, I am getting a value
> different than what the derivative and the initial values add up to.

My reading of the returned value from a call to ode is that it is not the values returned by the 'func'-function which only contain the derivative vectors, but rather the integrated values. That is what you are seeing in any case. You need to remember that integrators calculate the derivative and then further multiply by the "dt" increment which in your case will make the incremental values one-hundredth of the value you calculated by hand below.  (Basic calculus.)

> 
> My output:
> <http://r.789695.n4.nabble.com/file/n4707448/ROutput.png> 
> 
> but 'IH0' + 'the value of dIH.dt = Beta[1,1]*(nH-IH)*IH +
> Beta[1,2]*(nH-IH)*IL - gamma*IH for timestep one' = 1.999e-5.  Is there
> something that I am missing?  I am hoping to use this for more complicated
> derivatives, but I want to make sure I am using the package correctly first.
> 
> Thanks!
> 
> 
> 
> View this message in context: http://r.789695.n4.nabble.com/deSolve-ODE-Output-Question-tp4707448.html
> Sent from the R help mailing list archive at Nabble.com.

I normally also respond to the list and the questioner but have gotten an annoying number of bounces from Nabble users so am not responding directly.

-- 
David Winsemius
Alameda, CA, USA



More information about the R-help mailing list