Nothing
#' Distribution of chemical steady state concentration with uncertainty and variability
#'
#' @description
#' For a given chemical and fixed dose rate this function determines a
#' distribution of steady-state concentrations reflecting measurement uncertainty
#' an population variability. Uncertainty and variability are simulated via the
#' Monte Carlo method -- many sets of model parameters are drawn according to
#' probability distributions described in
#' Ring et al. (2017) (\doi{10.1016/j.envint.2017.06.004})
#' for human variability and Wambaugh et al. (2019)
#' (\doi{10.1093/toxsci/kfz205})
#' for measurement uncertainty. Monte Carlo samples are generated by the
#' function \code{\link{create_mc_samples}}. To allow rapid application of the
#' Monte Carlo
#' method we make use of analytical solutions for the steady-state concentration
#' for a particular model via a given route (when available) as opposed to
#' solving the model numerically (that is, using differential equations).
#' For each sample of the Monte Carlo method (as specified by argument
#' \code{samples}) the parameters for the analytical solution are varied. An
#' ensemble of steady-state predictions are produced, though by default only
#' the quantiles specified by argument \code{which.quantile} are provided. If
#' the full set of predicted values are desired use set the argument
#' \code{return.samples} to \code{TRUE}.
#'
#' @details
#' The chemical-specific steady-state concentration for a dose rate of
#' 1 mg/kg body weight/day can be used for in \emph{in vitro-in vivo}
#' extrapolation (IVIVE) of
#' a bioactive \emph{in vitro} concentration by dividing the \emph{in vitro}
#' concentration
#' by the steady-state concentration to predict a dose rate (mg/kg/day) that
#' would produce that concentration in plasma. Using quantiles of the
#' distribution (such as the upper 95th percentile) allow incorporation of
#' uncertainty and variability into IVIVE.
#'
#' Reverse Dosimetry Toxicodynamic IVIVE
#'
#' \if{html}{\figure{ivive.png}{options: width="60\%" alt="Reverse Dosimetry Toxicodynamic IVIVE"}}
#' \if{latex}{\figure{ivive.pdf}{options: width=12cm alt="Reverse Dosimetry Toxicodynamic IVIVE"}}
#'
#' Figure from Breen et al. (2021) (\doi{10.1080/17425255.2021.1935867}) shows
#' reverse dosimetry toxicodynamic IVIVE. Equivalent external dose is
#' determined by solving the TK model in reverse by deriving the external dose
#' (that is, TK model input) that produces a specified internal concentration
#' (that is, TK model output). Reverse dosimetry and IVIVE using HTTK relies on
#' the linearity of the models. We calculate a scaling factor to relate
#' \emph{in vitro}
#' concentrations (uM) to administered equivalent doses (AED). The scaling
#' factor is the inverse of the steady state plasma concentration (Css)
#' predicted for a 1 mg/kg/day exposure dose rate. We use Monte Carlo to
#' simulate variability and propagate uncertainty; for example, to calculate an
#' upper 95th percentile Css,95 for individuals who get higher plasma
#' concentrations from the same exposure.
#'
#' The Monte Carlo methods used here were recently updated and described by
#' Breen et al. (submitted).
#'
#'
#' httk-pop is used only for humans. For non-human species biological
#' variability is simulated by drawing parameters from uncorellated log-normal
#' distributions.
#'
#' Chemical-specific httk data are available primarily for human and for a few
#' hundred chemicals in rat. All in silico predictions are for human. Thus,
#' when species is specified as rabbit, dog, or mouse, the user can choose
#' to set the argument \code{default.to.human} to \code{TRUE} so that this
#' function uses the
#' appropriate physiological data (volumes and flows) but substitutes human
#' fraction unbound, partition coefficients, and intrinsic hepatic clearance.
#'
#' If the argument \code{tissue} is used, the steady-state concentration in that
#' tissue, if available, is provided. If that tissue is included in the model
#' used (specified by arguement \code{model}) then the actual tissue
#' concentration is
#' provided. Otherwise, the tissue-specific partition coefficient is used to
#' estimate the concentration from the plasma.
#'
#' The six sets of plausible IVIVE
#' assumptions identified by Honda et al. (2019)
#' (\doi{10.1371/journal.pone.0217564})
#' are: \tabular{lrrrr}{
#' \tab \emph{in vivo} Conc. \tab Metabolic Clearance \tab Bioactive Chemical
#' Conc. \tab TK Statistic Used* \cr Honda1 \tab Veinous (Plasma) \tab
#' Restrictive \tab Free \tab Mean Conc. \cr Honda2 \tab Veinous \tab
#' Restrictive \tab Free \tab Max Conc. \cr Honda3 \tab Veinous \tab
#' Non-restrictive \tab Total \tab Mean Conc. \cr Honda4 \tab Veinous \tab
#' Non-restrictive \tab Total \tab Max Conc. \cr Honda5 \tab Target Tissue \tab
#' Non-restrictive \tab Total \tab Mean Conc. \cr Honda6 \tab Target Tissue
#' \tab Non-restrictive \tab Total \tab Max Conc. \cr } *Assumption is
#' currently ignored because analytical steady-state solutions are currently
#' used by this function.
#'
#' @param chem.cas Chemical Abstract Services Registry Number (CAS-RN) -- if
#' parameters is not specified then the chemical must be identified by either
#' CAS, name, or DTXISD
#'
#' @param chem.name Chemical name (spaces and capitalization ignored) -- if
#' parameters is not specified then the chemical must be identified by either
#' CAS, name, or DTXISD
#'
#' @param dtxsid EPA's DSSTox Structure ID (\url{https://comptox.epa.gov/dashboard})
#' -- if parameters is not specified then the chemical must be identified by
#' either CAS, name, or DTXSIDs
#'
#' @param parameters Parameters from the appropriate parameterization function
#' for the model indicated by argument model
#'
#' @param samples Number of samples generated in calculating quantiles.
#'
#' @param which.quantile Which quantile from Monte Carlo simulation is
#' requested. Can be a vector.
#'
#' @param species Species desired (either "Rat", "Rabbit", "Dog", "Mouse", or
#' default "Human"). Species must be set to "Human" to run httkpop model.
#'
#' @param daily.dose Total daily dose, mg/kg BW.
#'
#' @param suppress.messages Whether or not to suppress output message.
#'
#' @param model Model used in calculation,'gas_pbtk' for the gas pbtk model,
#' 'pbtk' for the multiple compartment model,
#' '3compartment' for the three compartment model, '3compartmentss' for
#' the three compartment steady state model, and '1compartment' for one
#' compartment model. This only applies when httkpop=TRUE and species="Human",
#' otherwise '3compartmentss' is used.
#'
#' @param httkpop Whether or not to use population generator and sampler from
#' httkpop. This is overwrites censored.params and vary.params and is only for
#' human physiology. Species must also be set to 'Human'.
#'
#' @param httkpop.dt A data table generated by \code{\link{httkpop_generate}}.
#' This defaults to NULL, in which case \code{\link{httkpop_generate}} is
#' called to generate this table.
#'
#' @param invitrouv Logical to indicate whether to include in vitro parameters
#' in uncertainty and variability analysis
#'
#' @param calcrb2p Logical determining whether or not to recalculate the
#' chemical ratio of blood to plasma
#'
#' @param censored.params The parameters listed in censored.params are sampled
#' from a normal distribution that is censored for values less than the limit
#' of detection (specified separately for each parameter). This argument should
#' be a list of sublists. Each sublist is named for a parameter in
#' "parameters" and contains two elements: "CV" (coefficient of variation) and
#' "LOD" (limit of detection, below which parameter values are censored. New
#' values are sampled with mean equal to the value in "parameters" and standard
#' deviation equal to the mean times the CV. Censored values are sampled on a
#' uniform distribution between 0 and the limit of detection. Not used with
#' httkpop model.
#'
#' @param vary.params The parameters listed in vary.params are sampled from a
#' normal distribution that is truncated at zero. This argument should be a
#' list of coefficients of variation (CV) for the normal distribution. Each
#' entry in the list is named for a parameter in "parameters". New values are
#' sampled with mean equal to the value in "parameters" and standard deviation
#' equal to the mean times the CV. Not used with httkpop model.
#'
#' @param return.samples Whether or not to return the vector containing the
#' samples from the simulation instead of the selected quantile.
#'
#' @param tissue Desired steady state tissue concentration. Default is of NULL
#' typically gives whole body plasma concentration.
#'
#' @param concentration Desired concentration type: 'blood','tissue', or default
#' 'plasma'. In the case that the concentration is for plasma, selecting "blood"
#' will use the blood:plasma ratio to estimate blood concentration. In the case
#' that the argument 'tissue' specifies a particular tissue of the body,
#' concentration defaults to 'tissue' -- that is, the concentration in the
#' If cocentration is set to 'blood' or 'plasma' and 'tissue' specifies a
#' specific tissue then the value returned is for the plasma or blood in that
#' specific tissue.
#'
#' @param output.units Plasma concentration units, either uM or default mg/L.
#'
#' @param invitro.mc.arg.list List of additional parameters passed to
#' \code{\link{invitro_mc}}
#'
#' @param httkpop.generate.arg.list Additional parameters passed to
#' \code{\link{httkpop_generate}}.
#'
#' @param convert.httkpop.arg.list Additional parameters passed to the
#' convert_httkpop_* function for the model.
#'
#' @param parameterize.arg.list A list of arguments to be passed to the model
#' parameterization function (that is, parameterize_MODEL) corresponding to
#' argument "model". (Defaults to NULL.)
#'
#' @param calc.analytic.css.arg.list Additional parameters passed to
#'
#' @param Caco2.options Arguments describing how to handle Caco2 absorption data
#' that are passed to \code{\link{invitro_mc}} and the parameterize_[MODEL]
#' functions.
#' See \code{\link{get_fbio}} for further details.
#'
#' \code{\link{calc_analytic_css}}.
#'
#' @author Caroline Ring, Robert Pearce, John Wambaugh, Miyuki Breen, and Greg Honda
#'
#' @references
#' Wambaugh, John F., et al. "Toxicokinetic triage for
#' environmental chemicals." Toxicological Sciences 147.1 (2015): 55-67.
#'
#' \insertRef{ring2017identifying}{httk}
#'
#' \insertRef{honda2019using}{httk}
#'
#'\insertRef{rowland1973clearance}{httk}
#'
#' @keywords Monte-Carlo Steady-State
#'
#' @seealso \code{\link{calc_analytic_css}}
#'
#' @seealso \code{\link{create_mc_samples}}
#'
#' @return
#' Quantiles (specified by which.quantile) of the distribution of plasma
#' steady-stae concentration (Css) from the Monte Carlo simulation
#'
#' @examples
#' \donttest{
#' # Set the number of samples (NSAMP) low for rapid testing, increase NSAMP
#' # for more stable results. Default value is 1000:
#' NSAMP = 10
#'
#' # Basic in vitro - in vivo extrapolation with httk, convert 3 uM in vitro
#' # concentration of chemical with CAS 2451-62-9 to mg/kg/day:
#' set.seed(1234)
#' 3/calc_mc_css(chem.cas="2451-62-9", samples=NSAMP, output.units="uM")
#' # The significant digits should give the same answer as:
#' set.seed(1234)
#' calc_mc_oral_equiv(chem.cas="2451-62-9", conc=3, samples=NSAMP)
#'
#' set.seed(1234)
#' calc_mc_css(chem.name='Bisphenol A', output.units='uM',
#' samples=NSAMP, return.samples=TRUE)
#'
#' set.seed(1234)
#' calc_mc_css(chem.name='Bisphenol A', output.units='uM',
#' samples=NSAMP,
#' httkpop.generate.arg.list=list(method='vi'))
#'
#' # The following example should result in an error since we do not
#' # estimate tissue partitioning with '3compartmentss'.
#' set.seed(1234)
#' try(calc_mc_css(chem.name='2,4-d', which.quantile=.9,
#' samples=NSAMP,
#' httkpop=FALSE, tissue='heart'))
#'
#' # But heart will work with PBTK, even though it's lumped since we estimate
#' # a partition coefficient before lumping:
#' set.seed(1234)
#' calc_mc_css(chem.name='2,4-d', model='pbtk',
#' samples=NSAMP,
#' which.quantile=.9, httkpop=FALSE, tissue='heart')
#'
#' set.seed(1234)
#' calc_mc_css(chem.cas = "80-05-7", which.quantile = 0.5,
#' output.units = "uM", samples = NSAMP,
#' httkpop.generate.arg.list=list(method='vi', gendernum=NULL,
#' agelim_years=NULL, agelim_months=NULL,
#' weight_category = c("Underweight","Normal","Overweight","Obese")))
#'
#' params <- parameterize_pbtk(chem.cas="80-05-7")
#' set.seed(1234)
#' calc_mc_css(parameters=params,model="pbtk", samples=NSAMP)
#'
#' set.seed(1234)
#' # Standard HTTK Monte Carlo
#' calc_mc_css(chem.cas="90-43-7", model="pbtk", samples=NSAMP)
#' set.seed(1234)
#' # HTTK Monte Carlo with no measurement uncertainty (pre v1.10.0):
#' calc_mc_css(chem.cas="90-43-7",
#' model="pbtk",
#' samples=NSAMP,
#' invitro.mc.arg.list = list(
#' adjusted.Funbound.plasma = TRUE,
#' poormetab = TRUE,
#' fup.censored.dist = FALSE,
#' fup.lod = 0.01,
#' fup.meas.cv = 0.0,
#' clint.meas.cv = 0.0,
#' fup.pop.cv = 0.3,
#' clint.pop.cv = 0.3))
#'
#' # HTTK Monte Carlo with no HTTK-Pop physiological variability):
#' set.seed(1234)
#' calc_mc_css(chem.cas="90-43-7",model="pbtk",samples=NSAMP,httkpop=FALSE)
#'
#' # HTTK Monte Carlo with no in vitro uncertainty and variability):
#' set.seed(1234)
#' calc_mc_css(chem.cas="90-43-7",model="pbtk",samples=NSAMP,invitrouv=FALSE)
#'
#' # HTTK Monte Carlo with no HTTK-Pop and no in vitro uncertainty and variability):
#' set.seed(1234)
#' calc_mc_css(chem.cas="90-43-7" ,model="pbtk",
#' samples=NSAMP, httkpop=FALSE, invitrouv=FALSE)
#'
#' # Should be the same as the mean result:
#' calc_analytic_css(chem.cas="90-43-7",model="pbtk",output.units="mg/L")
#'
#' # HTTK Monte Carlo using basic Monte Carlo sampler:
#' set.seed(1234)
#' calc_mc_css(chem.cas="90-43-7",
#' model="pbtk",
#' samples=NSAMP,
#' httkpop=FALSE,
#' invitrouv=FALSE,
#' vary.params=list(Pow=0.3))
#'
#' # We can also use the Monte Carlo functions by passing a table
#' # where each row represents a different Monte Carlo draw of parameters:
#' p <- create_mc_samples(chem.cas="80-05-7")
#' # Use data.table for steady-state plasma concentration (Css) Monte Carlo:
#' calc_mc_css(parameters=p)
#' # Using the same table gives the same answer:
#' calc_mc_css(parameters=p)
#' # Use Css for 1 mg/kg/day for simple reverse toxicokinetics
#' # in Vitro-In Vivo Extrapolation to convert 15 uM to mg/kg/day:
#' 15/calc_mc_css(parameters=p, output.units="uM")
#' # Can do the same with calc_mc_oral_equiv:
#' calc_mc_oral_equiv(15, parameters=p)
#' }
#'
#' @import stats
#' @importFrom purrr compact
#' @importFrom data.table is.data.table
#' @export calc_mc_css
calc_mc_css <- function(chem.cas = NULL,
chem.name = NULL,
dtxsid = NULL,
parameters = NULL,
samples = 1000,
which.quantile = 0.95,
species = "Human",
daily.dose = 1,
suppress.messages = FALSE,
model = '3compartmentss',
httkpop = TRUE,
httkpop.dt = NULL,
invitrouv = TRUE,
calcrb2p = TRUE,
censored.params = list(),
vary.params = list(),
return.samples = FALSE,
tissue = NULL,
concentration = "plasma",
output.units = "mg/L",
invitro.mc.arg.list = NULL,
httkpop.generate.arg.list =
list(method = "direct resampling"),
convert.httkpop.arg.list = NULL,
parameterize.arg.list = NULL,
calc.analytic.css.arg.list = NULL,
Caco2.options=NULL
)
{
# We need to describe the chemical to be simulated one way or another:
if (is.null(chem.cas) &
is.null(chem.name) &
is.null(dtxsid) &
is.null(parameters))
stop('Parameters, chem.name, chem.cas, or dtxsid must be specified.')
if (is.null(model)) stop("Model must be specified.")
#Appease R CMD check --as-cran variable binding:
Css <- NULL
# We need to know model-specific information (from modelinfo_[MODEL].R])
# to set up the solver:
model <- tolower(model)
if (!(model %in% names(model.list)))
{
stop(paste("Model",model,"not available. Please select from:",
paste(names(model.list),collapse=", ")))
}
parameterize_function <- model.list[[model]]$parameterize.func
# Error handling for tissue argument:
if (!is.null(tissue))
{
if (is.null(model.list[[model]]$alltissues))
{
stop(paste("Tissues are not available for model", model))
}
if (!(tissue %in% model.list[[model]]$alltissues))
{
stop(paste("Tissue", tissue, "not available for model", model))
}
}
# Error handling for concentration arugment:
if (!(concentration %in% c("blood","tissue","plasma")))
{
stop("Concentration must be one of blood, tissue, or plasma")
}
#
#
# CREATE A TABLE OF PARAMETER VALUES WHERE EACH ROW IS A SEPARATE SET OF
# VALUES FOR WHICH Css SHOULD BE CALCULATED
#
#
if (!(data.table::is.data.table(parameters)))
parameter.dt <- do.call(create_mc_samples,
# we use purrr::compact to drop NULL values from arguments list:
args=purrr::compact(c(list(
chem.cas=chem.cas,
chem.name=chem.name,
dtxsid = dtxsid,
parameters=parameters,
samples=samples,
species=species,
suppress.messages=suppress.messages,
model=model,
httkpop=httkpop,
httkpop.dt=httkpop.dt,
invitrouv=invitrouv,
calcrb2p=calcrb2p,
censored.params=censored.params,
vary.params=vary.params,
return.samples=FALSE,
invitro.mc.arg.list=invitro.mc.arg.list,
httkpop.generate.arg.list=httkpop.generate.arg.list,
convert.httkpop.arg.list=convert.httkpop.arg.list,
parameterize.arg.list=parameterize.arg.list,
Caco2.options=Caco2.options))))
else parameter.dt <- parameters
#
# HERE LIES THE ACTUAL MONTE CARLO STEP:
#
# Calculate CSS for each row in the parameter matrix (each row corresponds to
# a different individual):
#
parameter.dt[,Css:= do.call(calc_analytic_css,
# we use purrr::compact to drop NULL values from arguments list:
args=purrr::compact(c(list(
parameters=.SD,
model=model,
suppress.messages=TRUE,
chem.cas=chem.cas,
chem.name=chem.name,
dtxsid=dtxsid,
tissue=tissue,
concentration=concentration,
output.units=output.units,
daily.dose=daily.dose,
clint.pvalue.threshold=
parameterize.arg.list$clint.pvalue.threshold),
calc.analytic.css.arg.list)))]
css.list <- parameter.dt$Css
if (!return.samples) out <- quantile(css.list,which.quantile,na.rm=TRUE)
else out <- css.list
if (!suppress.messages)
{
if (!return.samples)
{
if (is.null(tissue)) cat(paste(toupper(substr(species,1,1)),
substr(species,2,nchar(species)),sep=''),
"plasma concentration returned in",
output.units,
"units for",
which.quantile,
"quantile.\n")
else cat(paste(toupper(substr(species,1,1)),
substr(species,2,nchar(species)),sep=''),
tissue,
"concentration returned in",
output.units,
"units for",
which.quantile,
"quantile.\n")
} else {
if (is.null(tissue)) cat(paste(toupper(substr(species,1,1)),
substr(species,2,nchar(species)),sep=''),
"plasma concentration returned in",
output.units,
"units.\n")
else cat(paste(toupper(substr(species,1,1)),
substr(species,2,nchar(species)),sep=''),
tissue,
"concentration returned in",
output.units,
"units.\n")
}
}
# Cannot guarantee arbitrary precision:
out <- set_httk_precision(out)
return(out)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.