mapbayest {mapbayr}

Estimate parameters (maximum a posteriori)

Description

The main function of the mapbayr package. Performs a maximum a posteriori Bayesian estimation of parameters, from a mrgsolve model object and a dataset containing information about administrations and observed concentrations.

Usage

mapbayest(
  x,
  data = NULL,
  method = c("L-BFGS-B", "newuoa"),
  hessian = stats::optimHess,
  select_eta = NULL,
  lambda = 1,
  lloq = NULL,
  force_initial_eta = NULL,
  quantile_bound = 0.001,
  control = list(),
  check = TRUE,
  verbose = TRUE,
  progress = TRUE,
  reset = 50,
  output = NULL,
  ...
)

Arguments

x

the model object

data

NMTRAN-like data set

method

optimization method; the default is "L-BFGS-B" (from stat::optim()), alternatively "newuoa" for minqa::newuoa()

hessian

function used to compute the Hessian and variance-covariance matrix with (default is stats::optimHess, alternatively use nlmixr::nlmixrHess)

select_eta

a vector of numeric values, the numbers of the ETAs to be estimated (default is NULL, all ETAs non-equal to zero)

lambda

a numeric value, the weight applied to the model prior (default is 1)

lloq

a numeric value, the lower limit of quantification. If not NULL, LLOQ and BLQ (below limit of quantification) variables will be added to the data. The related records will be censored with the M3 method. Ignored if LLOQ already in the data.

force_initial_eta

a vector of numeric values to start the estimation from (default to 0 for "L-BFGS-B")

quantile_bound

a numeric value representing the quantile of the normal distribution admitted to define the bounds for L-BFGS-B (default is 0.001, i.e. 0.1%)

control

a list passed to the optimizer (see stats::optim() or minqa::newuoa() documentation)

check

check model code for mapbayr specification (a logical, default is TRUE)

verbose

print a message whenever optimization is reset (a logical, default is TRUE)

progress

print a progress bar (a logical, default is TRUE)

reset

maximum allowed reset of the optimizer with new initial eta values if numerical difficulties, or with new bounds (L-BFGS-B) if estimate equal to a bound. (a numeric, default is 50)

output

if NULL (the default) a mapbayests object is returned; if df a mapbay_tab dataframe is returned

...

for compatibility (not used)

Value

a mapbayests object. Basically a list containing:

• model: the model object

• arg.ofv.optim, arg.ofv.fix, arg.ofv.id: arguments passed to the optimization function. Useful for debugging but not relevant for a basic usage. Access to the data with get_data(x)

• opt.value: the original output of the optimization function

• final_eta: a list of individual vectors of final estimates. Access it with x$final_eta or get_eta(x).

• covariance: a list of individual variance-covariance matrix of estimation. Access it with x$covariance or get_cov(x).

• mapbay_tab: an output table containing the results of your estimations (data, IPRED, PRED, covariates, captured items, ETA etc...). Access it with x$mapbay_tab, as.data.frame(x) or as_tibble(x).

• information: run times and package versions.

See Also

hist.mapbayests

plot.mapbayests

use_posterior

Examples

# First, code a model
code1 <- "$PARAM ETA1 = 0, ETA2 = 0,
KA = 0.5, TVCL = 1.1, TVV = 23.3
$OMEGA 0.41 0.32
$SIGMA 0.04 0
$CMT DEPOT CENT
$PK
double CL=TVCL*exp(ETA1+ETA(1));
double V=TVV*exp(ETA2+ETA(2)) ;
$ERROR
double DV=CENT/V*(1+EPS(1))+EPS(2);
$PKMODEL ncmt = 1, depot = TRUE
$CAPTURE DV CL
"

my_model <- mrgsolve::mcode("my_model", code1)
# Then, import your data
my_data <- data.frame(ID = 1, TIME = c(0, 1.1, 5.2, 12.3), EVID = c(1,0,0,0), AMT = c(500, 0,0,0),
 CMT = c(1,2,2,2), DV = c(0, 15.1, 29.5, 22.3))
print(my_data)

# And estimate
my_est <- mapbayest(x = my_model, data = my_data)
print(my_est)
# see also plot(my_est) and hist(my_est)

# Use your estimation
get_eta(my_est)
get_param(my_est)
as.data.frame(my_est)
use_posterior(my_est)