| bnec {bayesnec} | R Documentation |
bnec
Description
Fits a variety of NEC models using Bayesian analysis and provides a model averaged predictions based on WAIC model weights
Usage
bnec(
formula,
data,
x_range = NA,
resolution = 1000,
sig_val = 0.01,
loo_controls,
x_var = NULL,
y_var = NULL,
trials_var = NULL,
model = NULL,
random = NULL,
random_vars = NULL,
...
)
Arguments
formula |
Either a |
data |
A |
x_range |
A range of predictor values over which to consider extracting ECx. |
resolution |
The length of the predictor vector used for posterior predictions, and over which to extract ECx values. Large values will be slower but more precise. |
sig_val |
Probability value to use as the lower quantile to test significance of the predicted posterior values against the lowest observed concentration (assumed to be the control), to estimate NEC as an interpolated NOEC value from smooth ECx curves. |
loo_controls |
A named |
x_var |
Removed in version 2.0. Use formula instead. Used to be a
|
y_var |
Removed in version 2.0. Use formula instead. Used to be a
|
trials_var |
Removed in version 2.0. Use formula instead. Used to be a
|
model |
Removed in version 2.0. Use formula instead. Used to be a
|
random |
Removed in version 2.0. Use formula instead. Used to be a
named |
random_vars |
Removed in version 2.0. Use formula instead. Used to be a
|
... |
Further arguments to |
Details
Overview
bnec serves as a wrapper for (currently) 23 (mostly) non-linear
equations that are classically applied to concentration(dose)-response
problems. The primary goal of these equations is to provide the user with
estimates of No-Effect-Concentration (NEC),
No-Significant-Effect-Concentration (NSEC), and Effect-Concentration
(of specified percentage 'x', ECx) thresholds.
These in turn are fitted through the brm function from
package brms and therefore further arguments to brm
are allowed. In the absence of those arguments, bnec makes
its best attempt to calculate distribution family, priors and initialisation
values for the user based on the characteristics of the data. Moreover, in
the absence of user-specified values, bnec sets the number of
iterations to 1e4 and warmup period to floor(iterations / 5) * 4.
The chosen models can be extended by the addition of brms special
"aterms" as well as group-level effects. See ?bayesnecformula
for details.
The available models/equations/formulas
The available equations (or models) can be found via the models
function. Since version 2.0, bnec requires a specific formula
structure
which is fully explained in the help file of bayesnecformula.
This formula incorporates the information regarding the chosen model(s). If
one single model is specified, bnec will return an object of
class bayesnecfit; otherwise if model is either a concatenation
of multiple models and/or a string indicating a family of models,
bnec will return an object of class
bayesmanecfit, providing they were successfully fitted. The
major difference being that the output of the latter includes Bayesian model
averaging statistics. These classes come with multiple associated
methods such as plot, autoplot,
summary, print,
model.frame and formula.
model may also be one of "all", meaning all of the available models
will be fit; "ecx" meaning only models excluding a specific NEC step
parameter will be fit; "nec" meaning only models with a specific NEC step
parameter will be fit; "bot_free" meaning only models without a "bot"
parameter (without a bottom plateau) will be fit; "zero_bounded" are models
that are bounded to be zero; or "decline" excludes all hormesis models, i.e.,
only allows a strict decline in response across the whole predictor range.
Notice that if one of these group strings is provided together with a
user-specified named list for the brm's argument
prior, the list names need to contain
the actual model names, and not the group string , e.g. if
model = "ecx" and prior = my_priors then
names(my_priors) must contain models("ecx"). To check
available models and associated parameters for each group,
use the function models or to check the parameters of a
specific model use the function show_params.
No-effect toxicity estimates
Regardless of the model(s) fitted, the resulting object will contain a
no-effect toxicity estimate. Where the fitted model(s) are NEC models (threshold
models, containing a step function - all models with "nec" as a
prefix) the no-effect estimate is a true
no-effect-concentration (NEC, see Fox 2010). Where the fitted model(s) are
smooth ECx models with no step function (all models with "ecx" as a
prefix), the no-effect estimate is a no-significant-effect-concentration
(NSEC, see Fisher and Fox 2023).
In the case of a bayesmanecfit that contains a mixture of both
NEC and ECx models, the no-effect estimate is a model averaged combination of
the NEC and NSEC estimates, and is reported as the N(S)EC
(see Fisher et al. 2023).
Further argument to brm
If not supplied via the brm argument family, the
appropriate distribution will be guessed based on the characteristics of the
input data. Guesses include: "bernoulli" / bernoulli / bernoulli(), "Beta" /
Beta / Beta(), "binomial" / binomial / binomial(), "beta_binomial" /
"beta_binomial", "Gamma" / Gamma / Gamma(), "gaussian" / gaussian /
gaussian(), "negbinomial" / negbinomial / negbinomial(), or "poisson" /
poisson / poisson(). Note, however, that "negbinomial" and "betabinomimal2"
require knowledge on whether the data is over-dispersed. As
explained below in the Return section, the user can extract the dispersion
parameter from a bnec call, and if they so wish, can refit the
model using the "negbinomial" family.
Other families can be considered as required, please raise an issue on the GitHub development site if your required family is not currently available.
As a default, bnec sets the brm argument
sample_prior to "yes" in order to sample draws from the priors in
addition to the posterior distributions. Among others, these samples can be
used to calculate Bayes factors for point hypotheses via
hypothesis.
Model averaging is achieved through a weighted sample of each fitted models
posterior predictions, with weights derived using functions
loo and loo_model_weights from
brms. Argument to both these functions can be passed via the
loo_controls argument. Individual model fits can be pulled out
for examination using function pull_out.
Additional technical notes
As some concentration-response data will use zero concentration
which can cause numerical estimation issues, a small offset is added (1 /
10th of the next lowest value) to zero values of concentration where
x_var are distributed on a continuous scale from 0 to infinity, or
are bounded to 0, or 1.
NAs are thrown away
Stan's default behaviour is to fail when the input data contains NAs. For
that reason brms excludes any NAs from input data prior to fitting,
and does not allow them back in as is the case with e.g. stats::lm and
na.action = exclude. So we advise that you exclude any NAs in your
data prior to fitting because if you so wish that should facilitate merging
predictions back onto your original dataset.
Value
If argument model is a single string, then an object of class
bayesnecfit; if many strings or a set,
an object of class bayesmanecfit.
References
Fisher R, Fox DR (2023). Introducing the no significant effect concentration (NSEC).Environmental Toxicology and Chemistry, 42(9), 2019–2028. doi: 10.1002/etc.5610.
Fisher R, Fox DR, Negri AP, van Dam J, Flores F, Koppel D (2023). Methods for estimating no-effect toxicity concentrations in ecotoxicology. Integrated Environmental Assessment and Management. doi:10.1002/ieam.4809.
Fox DR (2010). A Bayesian Approach for Determining the No Effect Concentration and Hazardous Concentration in Ecotoxicology. Ecotoxicology and Environmental Safety, 73(2), 123–131. doi: 10.1016/j.ecoenv.2009.09.012.
See Also
bayesnecformula,
check_formula,
models,
show_params
Examples
## Not run:
library(bayesnec)
data(nec_data)
# A single model
exmp_a <- bnec(y ~ crf(x, model = "nec4param"), data = nec_data, chains = 2)
# Two models model
exmp_b <- bnec(y ~ crf(x, model = c("nec4param", "ecx4param")),
data = nec_data, chains = 2)
## End(Not run)