Dear All,
please help with some thoughts on overcoming the following issues, if possible:
#R Code
require(deSolve)
require(FME)
pars <- list(k = 0.06,v=18)
intimes <- c(0,0.5,12,12.5,50)
input <- c(800,0,800,0,0)
forc <- approxfun(intimes, input, method="constant")
model <- function(pars, times=seq(0, 50, by = 1)) {
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,3,17.4,5.5,15.3,8,13.5,11.5,11.3,13,29.8,15.5,26.3,16.5,25,18,23.2,24,17.2,30,12.7,40,7.7,49.5,4.8))
colnames(observed) <- c("time", "y")
obj <- function(x, parset = names(x)) {
pars[parset] <- x
tout <- seq(0, 50, by = 1)
out <- model(pars, tout)
return(modCost(obs = observed, model = out))
}
prior <- function(p)
return( sum(((p-c(0.06,18))/c(0.04,6))^2 ))
final <- modMCMC(p = c(k=0.06,v=18), f = obj,
prior=prior,lower=c(0.0001,0.1),jump = NULL, niter = 1000,updatecov=100, wvar0 =
1,var0=NULL,burninlength = 50)
summary(final)
plot(final)
the values for "observed" were generated with k=0.05, v=20, which are
the retunrs I expect from the modMCMC run. Please comment on:
1. is there any way to get modMCMC to run multiple times faster with this data?
2. my accepted number of runs is very low. Is there a way to increase that?
3. any thoughts on improving convergence?
4. at times I only have one or two "observed" data points to work
with, which often results in "suboptimal" outputs from the modMCMC
runs. any thoughts on how to improve identifiability of parameter values with
less "observed" data points being available?
thank you for the help,
Andras
[[alternative HTML version deleted]]
