Dear All, I like to simulate a physiologically based pharmacokinetics model using R but am having a problem with the daspk routine.
The same problem has been implemented in Berkeley madonna and Winbugs so that I know that it is working. However, with daspk it is not, and the numbers are everywhere! Please see the following and let me know if I am missing something... Thanks a lot in advance, In-Sun #--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------- library("deSolve") y <- c(Agi = 0,Alu = 0, Abr = 0, Ah = 0, Ali = 0, Ak = 0, Am = 0, Ask = 0, Aad = 0, Apa = 0, Asp = 0, Aar = 0, Ave = 0) times = seq(0, 100, length=100) pars <- c( dose = 80 * 0.26, doseduration = 10, Vmax = 1.44, Km = 0.96, s = 1.33, fp = 0.236, Kpfgi=0.324, Kpflu = 1.092, Kpfbr= 0.155 , Kpfh=0.767, Kpfli = 0.551, Kpfk=0.537, Kpfm=0.339, Kpfsk=0.784, Kpfad=0.465, Kpfpa=0.595, Kpfsp=0.410, Qar = 51.9, Qve = 51.9, Qgi = 12.3, Qlu = 51.9, Qbr = 3.2, Qh = 6.4, Qli = 16.5, Qk = 12.8, Qm = 7.6, Qsk = 5.0, Qad = 0.4, Qpa = 1.0, Qsp = 1.0, Var = 7.0, Vve = 14.1, Vgi = 12.4, Vlu = 1.3, Vbr = 1.3, Vh = 1.2, Vli = 12.4, Vk = 2.2, Vm = 140.0, Vsk = 49.0, Vad = 11.2, Vpa = 1.0, Vsp = 1.0 ) Fun_ODE <- function(t,y, pars){ with (as.list(c(y, pars)), { It <- dose/doseduration Car <- Aar/Var Cve <- Ave/Vve Clu <- Alu/Vlu Cli <- Ali/Vli Cbr <- Abr/Vbr Ch <- Ah/Vh Cpa <- Apa/Vpa Csp <- Asp/Vsp Cgi <- Agi/Vgi Ck <- Ak/Vk Cm <- Am/Vm Cad <- Aad/Vad Csk <- Ask/Vsk Kpbbr <- s*fp*Kpfbr Kpbli <- s*fp*Kpfli Kpbh <- s*fp*Kpfh Kpbpa <- s*fp*Kpfpa Kpbsp <- s*fp*Kpfsp Kpbgi <- s*fp*Kpfgi Kpbk <- s*fp*Kpfk Kpbm <- s*fp*Kpfm Kpbad <- s*fp*Kpfad Kpbsk <- s*fp*Kpfsk Kpblu <- s*fp*Kpflu dAar <- (Clu/Kpblu - Car)*Qar dAve <- if (t < 10) It + Cbr*Qbr/Kpbbr + Ch *Qh/Kpbh + Cli*Qli/Kpbli + Ck*Qk/Kpbk + Cm*Qm/Kpbm + Csk * Qsk /Kpbsk + Cad*Qad/Kpbad - Cve*Qve else Cbr*Qbr/Kpbbr + Ch *Qh/Kpbh + Cli*Qli/Kpbli + Ck*Qk/Kpbk + Cm*Qm/Kpbm + Csk * Qsk /Kpbsk + Cad*Qad/Kpbad - Cve*Qve dAlu <- (Cve-Clu/Kpblu)*Qlu dAli <- ((Qli - Qgi- Qpa-Qsp)*Car + Cgi*Qgi/Kpbgi + Csp*Qsp/Kpbsp + Cpa*Qpa/Kpbpa - Cli*Qli/Kpbli) - Vmax*Cli/Kpfli/(Km + Cli/Kpfli) dAbr <- (Car - Cbr/Kpbbr)*Qbr dAh <- (Car - Ch/Kpbh)*Qh dApa <- (Car - Cpa/Kpbpa)*Qpa dAsp <- (Car - Csp/Kpbsp)*Qsp dAgi <- (Car - Cgi/Kpbgi)*Qgi dAk <- (Car - Ck/Kpbk)*Qk dAm <- (Car - Cm/Kpbm)*Qm dAad <- (Car - Cad/Kpbad)*Qad dAsk <- (Car - Csk/Kpbsk)*Qsk return(list(dy = c(dAar, dAve, dAlu, dAli, dAbr, dAh, dApa, dAsp, dAgi, dAk, dAm, dAad, dAsk), Car = Car, Cve=Cve, Clu=Clu, Cli=Cli, Cbr=Cbr, Ch=Ch, Cpa=Cpa, Csp=Csp, Cgi=Cgi, Ck=Ck, Cm=Cm, Cad=Cad, Csk=Csk)) }) } ODE <- as.data.frame(daspk(y = y, times = times, func = Fun_ODE, parms = pars, atol = 1e-10, rtol = 1e-10)) -- Dr In-Sun Nam Knutsson Research Associate The Centre for Applied Pharmacokinetic Research (CAPKR) School of Pharmacy and Pharmaceutical Sciences University of Manchester Stopford Building Oxford Road Manchester U.K. Phone: +44 161 275 2355 Email: in-sunnam.knuts...@manchester.ac.uk [[alternative HTML version deleted]] ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.