Dear List,
 
is it possible to run R2OpenBugs and lsoda from desolve together and 'bypass" the differential equation solver in OpenBugs? As suggested, here is a bit of background on my problem: I am trying to use one and two-compartent pharmacokinetic models for individual dose optimization of aminoglycosides and vancomycin using differential equation systems (lsoda) and bayesian feedback (modMCMC). The system works great, but runs for a long time. Currently I use R2OpenBugs and analytical solutions (fast and accurate results), but would like to grow independent of OpenBugs and use diff equations for some dynamic models that can not be solved analytically. 
a simple example of doing so perhaps with a 1 compartment model (f aprorpiate) would be appreciated. this is the code I currently run:
 
require(deSolve)
require(FME)
pars <- list(k = 0.08,v=15)
intimes <- c(0,0.5,12)
input   <- c(800,0,0)
forc <- approxfun(intimes, input, method="constant")
model <- function(pars, times=seq(0, 13, by = 0.01)) {
  derivs <- function(t, state, pars) { 
     with(as.list(c(state, pars)), {
       inp <- forc(t)
       dy1 <- - k* y[1] + inp/v
       return(list(c(dy1)))
 })
}
state <- c(y = 0)
return(lsoda(y = state, times = times, func = derivs, parms = pars))
}
observed <- matrix (nc=2,byrow=2,data=c(0,0,0.5,19.7,11.5,11.3))
colnames(observed) <- c("time", "y")
obj <- function(x, parset = names(x)) {
 pars[parset] <- x
 tout <- seq(0, 13, by = 0.01)
 out <- model(pars, tout)
 return(modCost(obs = observed, model = out))
 }
 
prior <- function(p)
  return( sum(((p-c(0.08,15))/c(0.04,6))^2 ))
final <- modMCMC(p = c(k=0.08,v=15), f = obj, prior=prior,lower=c(0.0001,0.1),jump = NULL, niter = 1000,updatecov=10, wvar0 = 1,var0=NULL,burninlength = 10)
summary(final)
plot(final)
final$bestpar
thanks,
 
Andras
	[[alternative HTML version deleted]]

