R help - Jun 2004 - olesolve: stepsize

Hi, 

I am doing a project on the simulation of glucose metabolism based on a
pharmacokinetic modeling in which we have 4 differential equations. I did
this in R by using the odesolve package. It works very well, but I have two
questions:

 

Here is the odemodel function

_________________________________________________

Ogtt.Odemodel <- function(t, y, p) { 

    absx <- c(-60, -45, -30, -15, 0, 5, 10, 15, 20, 25, 30, 45, 60, 75, 90,
105, 120, 140, 160, 180, 210, 240, 270, 301) 

    absx <- absx + 60

    absy <- c(0, 0, 0, 0, 0, 156.1225165, 266.8509934, 432.5794702,
609.8807946, 759.0364238, 824.7649004, 636.0562915, 482.7450332,
395.1870861, 329.1953641, 591.2408942, 214.6274834, 180.316, 123.3885,
111.0497, 88.83972, 53.30383, 21.32153, 0) 

    sx   <- c(-60, -37.5, -22.5, -7.5, 0, 2.5, 7.5, 12.5, 17.5, 22.5, 27.5,
37.5, 52.5, 67.5, 82.5, 97.5, 112.5, 130, 150, 170, 195, 225, 255, 285, 315)


    sx   <- sx + 60

    sy   <- c(81.32, 81.32, 81.32, 81.32, 81.32, 218.24, 329.41, 448.82,
527.06, 588.82, 568.24, 613.53, 564.12, 415.88, 354.12, 325.29, 267.65,
222.35, 115.29, 98.31, 89.81, 85.57, 81.32, 81.32, 81.32) 

    abs  <-  approxfun(absx, absy, method = "linear") 

    s    <- approxfun(sx, sy, method = "linear") 

    Rtb  <- 81.32 

    ib   <-
((1-p["f"])*Rtb*p["kg"]*(p["hgo0"]/p["vg"])/(p["k21"]
+ p["k01"]
+ p["kl"] - p["k21"]*p["k12"]/(p["k12"]
+ p["k02"]))/(p["ki"]*p["vi"]))

    # G1(t)-y(1); G2(t)-y(2);I(t)-y(3); X(t)-y(4) 

    yd1 <- (abs(t) + p["hgo0"] - (p["kl"] +
p["fx"]*y[4])*y[1]*p["vg"] +
p["k12"]*y[2])/p["vg"] - (p["k21"] +
p["k01"])*y[1]

    yd2 <- p["k21"]*y[1]*p["vg"] - (p["k12"] +
p["k02"] +
(1-p["fx"])*y[4])*y[2] 

    yd3 <- (1-p["f"])*s(t)*p["kg"]*y[1]/p["vi"]
- p["ki"]*y[3]

    yd4 <- -p["p2"]*y[4] + p["p3"]*(y[3]-ib) 

    list(c(yd1, yd2, yd3, yd4)) 

}

 

__________________________________________________

 

First is that The computing time for each individual simulation is about 6
sec if I set up time from 0 - 360 min and calculate on each minute.  In
order to reduce computing time, I want to calculate for every two minutes,
but I don't know how to do it.  I try to make t to be t <- seq(0, 360, by
2), however, the computing time is still about 6 sec instead of 3.  I also
found that there is a hmin parameter in odesolve. But again, it doesn't
work.

 

The second question is that I found that some time, the last two or three
points of calculation is NA.   

 

Thanks.

 

Guan

 


	[[alternative HTML version deleted]]

"Zhixin Guan" <zguan at usc.edu> writes:
> Hi, 
> 
> I am doing a project on the simulation of glucose metabolism based on a
> pharmacokinetic modeling in which we have 4 differential equations. I did
> this in R by using the odesolve package. It works very well, but I have two
> questions:
> 
>  
> 
> Here is the odemodel function
> 
> _________________________________________________
> 
> Ogtt.Odemodel <- function(t, y, p) { 
.....
> }
> First is that The computing time for each individual simulation is about 6
> sec if I set up time from 0 - 360 min and calculate on each minute.  In
> order to reduce computing time, I want to calculate for every two minutes,
> but I don't know how to do it.  I try to make t to be t <- seq(0,
360, by > 2), however, the computing time is still about 6 sec instead of 3. 
I also
> found that there is a hmin parameter in odesolve. But again, it doesn't
> work.
The thing you're not telling us is what solver you are using! lsoda(),
I suppose. In that case, the times you specify are just read-off
times, whereas the actual integration proceeds in potentially much
smaller time step. If you want the solver to take larger steps, you
should probably loosen the tolerances and maybe also increase hmin. 
  
> 
> The second question is that I found that some time, the last two or three
> points of calculation is NA.   
This might be due to using approxfun with values beyond the end of
your sx data.

-- 
   O__  ---- Peter Dalgaard             Blegdamsvej 3  
  c/ /'_ --- Dept. of Biostatistics     2200 Cph. N   
 (*) \(*) -- University of Copenhagen   Denmark      Ph: (+45) 35327918
~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk)             FAX: (+45) 35327907