Skip to content

deSolve ODE Output Question

2 messages · walke554, David Winsemius

Original
walke554
# 
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 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.
David Winsemius
#
On May 20, 2015, at 6:45 AM, walke554 wrote:

            

      
      
      
    

        
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.)

            

      
      
      
    

        
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.