glm.rmap {hdbayes}

Posterior of robust meta-analytic predictive prior (RMAP)

Description

Final step for sampling from the posterior distribution of a GLM using the RMAP by Schmidli et al. (2014) doi:10.1111/biom.12242.

Usage

glm.rmap(
  formula,
  family,
  curr.data,
  probs,
  means,
  covs,
  curr.offset = NULL,
  w = 0.1,
  norm.vague.mean = NULL,
  norm.vague.sd = NULL,
  curr.disp.mean = NULL,
  curr.disp.sd = NULL,
  iter_warmup = 1000,
  iter_sampling = 1000,
  chains = 4,
  ...
)

Arguments

formula

a two-sided formula giving the relationship between the response variable and covariates.

family

an object of class family. See ?stats::family.

curr.data

a data.frame giving the current data.

probs

a vector of mixing proportions in the mixture approximation to the prior induced by the BHM. Obtained from the outputs of the glm.rmap.bhm.approx() function.

means

a matrix with the jth column being the mean vector for the jth component in the mixture approximation to the prior induced by the BHM. Obtained from the outputs of the glm.rmap.bhm.approx() function.

covs

a 3-dimensional array giving the covariance matrices for the mixture approximation to the prior induced by the BHM. Obtained from the outputs of the glm.rmap.bhm.approx() function. the means for the half-normal hyperpriors on the sd hyperparameters of regression coefficients.

curr.offset

a vector whose dimension is equal to the rows of the current data set giving an offset for the current data. Defaults to a vector of 0s.

w

a scalar between 0 and 1 giving how much weight to put on the historical data.

norm.vague.mean

a scalar or a vector whose dimension is equal to the number of regression coefficients giving the means for the vague normal prior on regression coefficients. If a scalar is provided, norm.vague.mean will be a vector of repeated elements of the given scalar. Defaults to a vector of 0s.

norm.vague.sd

a scalar or a vector whose dimension is equal to the number of regression coefficients giving the sds for the vague normal prior on regression coefficients. If a scalar is provided, same as for norm.vague.mean. Defaults to a vector of 10s.

curr.disp.mean

a scalar giving the mean for the half-normal hyperprior on the dispersion parameter for the current data. Defaults to a vector of 0s.

curr.disp.sd

a scalar giving the sd for the half-normal hyperprior on the dispersion parameter for the current data. Defaults to a vector of 10s.

iter_warmup

number of warmup iterations to run per chain. Defaults to 1000. See the argument iter_warmup in sample() method in cmdstanr package.

iter_sampling

number of post-warmup iterations to run per chain. Defaults to 1000. See the argument iter_sampling in sample() method in cmdstanr package.

chains

number of Markov chains to run. Defaults to 4. See the argument chains in sample() method in cmdstanr package.

...

arguments passed to sample() method in cmdstanr package (e.g. seed, refresh, init).

Details

This function samples from the posterior distribution of a GLM using the RMAP. The first component of the RMAP is a prior induced by the Bayesian hierarchical model (BHM). We approximate this component by a mixture of multivariate normal distributions where the parameters are obtained from the outputs of the glm.rmap.bhm.approx() function. The second component is a vague (noninformative) multivariate normal prior. We assume that the covariance matrix of the vague prior is a diagonal matrix.

Value

The function returns an object of class draws_df giving posterior samples.

References

Schmidli, H., Gsteiger, S., Roychoudhury, S., O’Hagan, A., Spiegelhalter, D., and Neuenschwander, B. (2014). Robust meta‐analytic‐predictive priors in clinical trials with historical control information. Biometrics, 70(4), 1023–1032.

See Also

glm.rmap.bhm() for the first step and glm.rmap.bhm.approx() for the second step of implementing RMAP.

Examples

  if (instantiate::stan_cmdstan_exists()) {
    data(actg019) ## current data
    data(actg036) ## historical data
    ## take subset for speed purposes
    actg019 = actg019[1:150, ]
    actg036 = actg036[1:100, ]
    hist_data_list = list(actg036)
    samples_bhm = glm.rmap.bhm(
      formula = outcome ~ scale(age) + race + treatment + scale(cd4),
      family = binomial('logit'),
      hist.data.list = hist_data_list,
      chains = 1, iter_warmup = 1000, iter_sampling = 2000
    )$beta_pred
    res_approx = glm.rmap.bhm.approx(
      samples.bhm = samples_bhm,
      G = 1:5, verbose = FALSE
    )
    glm.rmap(
      formula = outcome ~ scale(age) + race + treatment + scale(cd4),
      family = binomial('logit'),
      curr.data = actg019,
      probs = res_approx$probs,
      means = res_approx$means,
      covs = res_approx$covs,
      chains = 1, iter_warmup = 1000, iter_sampling = 2000
    )
  }