bamlss

Fit Bayesian Additive Models for Location Scale and Shape (and Beyond)

This is the main model fitting function of the package. Function bamlss() is a wrapper function that parses the data and the model formula, or extended bamlss.formula, as well as the bamlss.family into a bamlss.frame. The bamlss.frame then holds all model matrices and information that is needed for setting up estimation engines. The model matrices are based on mgcv infrastructures, i.e., smooth terms are constructed using smooth.construct and smoothCon. Therefore, all mgcv model term constructors like s, te, t2 and ti can be used. Identifiability conditions are imposed using function gam.side.

After the bamlss.frame is set up function bamlss() applies optimizer and/or sampling functions. These functions can also be provided by the user. See the details below on how to create new engines to be used with function bamlss().

Finally, the estimated parameters and/or samples are used to create model output results like summary statistics or effect plots. The computation of results may also be controlled by the user.

Usage

bamlss(formula, family = "gaussian", data = NULL,
  start = NULL, knots = NULL, weights = NULL,
  subset = NULL, offset = NULL, na.action = na.omit,
  contrasts = NULL, reference = NULL, transform = NULL,
  optimizer = NULL, sampler = NULL, samplestats = NULL,
  results = NULL, cores = NULL, sleep = NULL,
  combine = TRUE, model = TRUE, x = TRUE,
  light = FALSE, ...)

Details

The main idea of this function is to provide infrastructures that make it relatively easy to create estimation engines for new problems, or write interfaces to existing software packages.

The steps that are performed within the function are:

• First, the function parses the data, the formula or the extended bamlss.formula as well as the bamlss.family into a model frame like object, the bamlss.frame. This object holds all necessary model matrices and information that is needed for subsequent model fitting engines. Per default, all package mgcv smooth term constructor functions like s, te, t2 and ti can be used (see also function smooth.construct), however, even special user defined constructors can be included, see the manual of bamlss.frame.

• In a second step, the bamlss.frame can be transformed, e.g., if a mixed model representation of smooth terms is needed, see function randomize.

• Then an optimizer function is started, e.g., a function that finds posterior mode estimates of the parameters. A convention for model fitting engines is that such functions should have the following arguments:

  optimizer(x, y, family, start, weights, offset, ...)

  Internally, function bamlss() will send the x object that holds all model matrices, the response y object, the family object, starting values for the parameters, possible weights and offsets of the created bamlss.frame to the optimizer function (see the manual of bamlss.frame for more details on the x, y and other objects). The job of the optimizer is to return a named numeric vector of optimum parameters. The names of the parameters should be such that they can be uniquely mapped to the corresponding model matrices in x. See function parameters for more details on parameter names. The default optimizer function is bfit. The optimizer can return more information than only the optimum parameters. It is possible to return a list, the convention here is that an element named "parameters" then holds the named vector of estimated parameters. Possible other return values could be fitted values, the Hessian matrix, information criteria or information about convergence of the algorithm, etc. Note that the parameters are simply added to the bamlss.frame in an (list) entry named parameters.

• After the optimization step, a sampler function is started. The arguments of such sampler functions are the same as for the optimizer functions

  sampler(x, y, family, start, weights, offset, ...)

  Sampler functions must return a matrix of samples, each row represents one iteration and the matrix can be coerced to mcmc objects. The function may return a list of samples, e.g., if multiple chains are returned each list entry then holds one sample matrix of one chain. The column names of the sample matrix should be the same as the names of estimated parameters. For a possible naming convention see function parameters, which ensures unique mapping of samples with the model matrices in the x object of the bamlss.frame. The samples are added to the bamlss.frame in an (list) entry named samples.

• Next, the samplestats function is applied. This function can compute any quantity from the samples and the x object, the arguments of such functions are

  samplestats(samples, x, y, family, ...)

  where argument samples are the samples returned from the sampler function, and x, y and family are the same objects as passed to the optimizer and or sampler functions. For example, the default function in bamlss() for this task is also called samplestats and returns the mean of the log-likelihood and the log-posterior computed of all samples, as well as the DIC.

• The last step is to compute more complex information about the model using the results function. The arguments of such results functions are

  results(bamlss.frame, ...)

  here, the full bamlss.frame including possible parameters and samples is passed to the function within bamlss(). The default function for this task is results.bamlss.default which returns an object of class "bamlss.results" for which generic plotting functions are and a summary function is provided. Hence, the user can control the output of the model, the plotting and summary statistics, too.

Note that function transform(), optimizer(), sampler(), samplestats() and results() can be provided from the bamlss.family object, e.g., if a bamlss.family object has an element named "optimizer", which represents a valid optimizer function such as bfit, exactly this optimizer function will be used as a default when using the family.

References

Umlauf N, Klein N, Zeileis A (2019). BAMLSS: Bayesian Additive Models for Location, Scale and Shape (and Beyond). Journal of Computational and Graphical Statistics, 27(3), 612--627. 10.1080/10618600.2017.1407325

Umlauf N, Klein N, Simon T, Zeileis A (2019). “bamlss: A Lego Toolbox for Flexible Bayesian Regression (and Beyond).” arXiv:1909.11784, arXiv.org E-Print Archive. https://arxiv.org/abs/1909.11784

See Also

bamlss.frame, family.bamlss, bamlss.formula, randomize, bamlss.engine.setup, bfit, GMCMC, continue, coef.bamlss, parameters, predict.bamlss, plot.bamlss

Aliases

• bamlss

Examples

# NOT RUN {
## Simulated data example.
d <- GAMart()
f <- num ~ s(x1) + s(x2) + s(x3) + te(lon, lat)
b <- bamlss(f, data = d)
summary(b)
plot(b)
plot(b, which = 3:4)
plot(b, which = "samples")

## Use of optimizer and sampler functions:
## * first run optimizer,
b1 <- bamlss(f, data = d, optimizer = bfit, sampler = FALSE)
print(b1)
summary(b1)

## * afterwards, start sampler with staring values,
b2 <- bamlss(f, data = d, start = coef(b1), optimizer = FALSE, sampler = GMCMC)
print(b2)
summary(b2)

## Continue sampling.
b3 <- continue(b2, n.iter = 12000, burnin = 0, thin = 10)
plot(b3, which = "samples")
plot(b3, which = "max-acf")
plot(b3, which = "max-acf", burnin = 500, thin = 4)
# }