diff --git a/R/poped.R b/R/poped.R index 7697d17c..a8c08c15 100644 --- a/R/poped.R +++ b/R/poped.R @@ -3060,3 +3060,68 @@ babel.poped.database <- function(popedInput, ..., optTime=NA) { call.=FALSE) } } + +#' Get the bpop_idx by variable name for a poped database created by `babelmixr2` +#' +#' This may work for other poped databases if the population parameters are named. +#' +#' @param popedInput The babaelmixr2 created database +#' @param var variable to query +#' @return index of the variable +#' @export +#' @author Matthew L. Fidler +#' +#' @examples +#' +#' if (requireNamespace("PopED", quietly=TRUE)) { +#' +#' f <- function() { +#' ini({ +#' tV <- 72.8 +#' tKa <- 0.25 +#' tCl <- 3.75 +#' tF <- fix(0.9) +#' pedCL <- 0.8 +#' +#' eta.v ~ 0.09 +#' eta.ka ~ 0.09 +#' eta.cl ~0.25^2 +#' +#' prop.sd <- fix(sqrt(0.04)) +#' add.sd <- fix(sqrt(5e-6)) +#' +#' }) +#' model({ +#' V<-tV*exp(eta.v) +#' KA<-tKa*exp(eta.ka) * (pedCL**isPediatric) # add covariate for pediatrics +#' CL<-tCl*exp(eta.cl) +#' Favail <- tF +#' +#' N <- floor(t/TAU)+1 +#' y <- (DOSE*Favail/V)*(KA/(KA - CL/V)) * +#' (exp(-CL/V * (t - (N - 1) * TAU)) * +#' (1 - exp(-N * CL/V * TAU))/(1 - exp(-CL/V * TAU)) - +#' exp(-KA * (t - (N - 1) * TAU)) * (1 - exp(-N * KA * TAU))/(1 - exp(-KA * TAU))) +#' +#' y ~ prop(prop.sd) + add(add.sd) +#' }) +#' } +#' +#' e <- et(c( 1,8,10,240,245)) +#' +#' babel.db <- nlmixr2(f, e, "poped", +#' popedControl(m = 2, +#' groupsize=20, +#' bUseGrouped_xt=TRUE, +#' a=list(c(DOSE=20,TAU=24,isPediatric = 0), +#' c(DOSE=40, TAU=24,isPediatric = 0)))) +#' +#' babelBpopIdx(babel.db, "pedCL") +#' +#' } +babelBpopIdx <- function(popedInput, var) { + .w <- which(rownames(popedInput$parameters$bpop) == var) + if (length(.w) == 1L) return(.w) + stop("cannot find '", var, "' in the baelmixr2 poped model", + call.=FALSE) +} diff --git a/src/poped.cpp b/src/poped.cpp index 9ca203fe..13aaf5d6 100644 --- a/src/poped.cpp +++ b/src/poped.cpp @@ -467,7 +467,9 @@ void popedSolveFidMat(arma::mat &matMT, NumericVector &theta, int id, int nrow, setIndIdx(ind, j); kk = getIndIx(ind, j); curT = getTime(kk, ind); - isMT = getIndEvid(ind, kk) >= 10 && getIndEvid(ind, kk) <= 99; + int evid = getIndEvid(ind, kk); + isMT = evid >= 10 && evid <= 99; + // REprintf("curT: %f; evid: %d; isMT: %d; kk: %d\n", curT, evid, isMT, kk); if (isDose(getIndEvid(ind, kk))) { rxInner.calc_lhs(id, curT, getOpIndSolve(op, ind, j), lhs); continue;