Bayesian Model Averaging in Clinical Research

Introduction

Clinical research frequently requires investigators to make decisions about which statistical model should be used to describe a dataset.

For example, suppose a clinical trial contains several baseline covariates that could potentially predict treatment response. A logistic regression model might include age, sex, disease severity, biomarker status, prior treatment, and baseline laboratory values.

But which variables actually belong in the model?

A conventional analysis often answers this question by selecting a single model. Variables may be chosen using clinical knowledge, a prespecified model, stepwise selection, information criteria, or another model-selection procedure.

The problem is that model selection itself is uncertain. Two models may both be reasonably compatible with the observed data. If the analysis simply chooses one of them and treats it as though it were known to be the correct model, the uncertainty associated with that selection is ignored.

What Is Model Uncertainty?

Suppose we are interested in predicting clinical response using three candidate covariates:

• baseline disease severity;
• biomarker status; and
• age.

There are multiple plausible regression models. For example:

A conventional model-selection procedure might determine that \(M_2\) is the "best" model.

But the data may not strongly distinguish \(M_2\) from \(M_4\).

If we select \(M_2\) and then proceed as though \(M_2\) were known to be the correct model, we have effectively ignored uncertainty about whether age should also be included.

The Bayesian Perspective

Let \(M_k\) denote candidate model \(k\). Bayesian inference assigns a prior probability to each model:

After observing the data \(D\), Bayes' theorem gives the posterior probability of the model:

Here:

• \(P(M_k)\) is the prior probability of model \(M_k\);
• \(P(D\mid M_k)\) is the marginal likelihood of the data under model \(M_k\);
• \(P(M_k\mid D)\) is the posterior probability of model \(M_k\); and
• \(K\) is the total number of candidate models.

The denominator ensures that the posterior model probabilities sum to one.

What Is Bayesian Model Averaging?

Once posterior probabilities have been assigned to the candidate models, Bayesian Model Averaging combines the model-specific inferences using those posterior probabilities as weights.

Suppose the parameter of interest is \(\theta\). Under model \(M_k\), the posterior mean is:

The Bayesian Model Averaging estimate is then:

In words:

Why Not Just Choose the Best Model?

Suppose two models have posterior probabilities of 0.55 and 0.45. A model-selection approach would choose the first model.

But that decision discards almost half of the posterior probability assigned to the second model.

Bayesian Model Averaging instead retains the uncertainty:

This can produce more realistic inference when several models remain plausible.

Posterior Inclusion Probability

One of the most useful quantities in Bayesian variable selection is the posterior inclusion probability.

Suppose \(X_j\) is a candidate predictor. Its posterior inclusion probability is:

This quantity answers:

For example, suppose biomarker status appears in models whose posterior probabilities sum to 0.82. Then:

This means that, under the specified model space and priors, there is substantial posterior support for including the biomarker.

It is not the same statement as saying that there is an 82% probability that the biomarker has a clinically meaningful effect.

Posterior Model Odds and Bayes Factors

For two competing models, posterior odds can be written as:

The first term is the Bayes factor:

Thus:

This provides an intuitive connection between Bayesian model comparison and Bayesian Model Averaging.

A Clinical Research Example

Consider a randomized clinical trial evaluating a new treatment for a chronic disease. The primary endpoint is a binary clinical response at Week 24. Investigators are interested in estimating the probability of response for a future patient and in understanding which baseline characteristics contribute to response.

Three candidate predictors are considered:

Rather than assuming that all three variables must be included, the analyst considers all eight possible subsets of these predictors.

Because there are three candidate predictors:

candidate models are possible.

The Eight Candidate Models

Suppose the prior probability is distributed uniformly across the eight models:

This is one possible prior specification. Other prior structures can be used when there is stronger prior knowledge about which variables are likely to matter.

Posterior Model Probabilities

After fitting the eight models and calculating their marginal likelihoods, suppose the posterior probabilities are:

The probabilities sum to:

The single most probable model is \(M_5\), with posterior probability 0.39.

But \(M_5\) has only a 39% posterior probability. That means selecting \(M_5\) as though it were unquestionably correct would ignore 61% of the posterior probability assigned to other models.

Calculating Posterior Inclusion Probabilities

Consider baseline disease severity. Severity appears in models:

• \(M_2\): 0.18
• \(M_5\): 0.39
• \(M_6\): 0.11
• \(M_8\): 0.18

Therefore:

The posterior inclusion probability for disease severity is therefore:

Now consider biomarker status. It appears in \(M_3\), \(M_5\), \(M_7\), and \(M_8\):

Finally, age appears in \(M_4\), \(M_6\), \(M_7\), and \(M_8\):

The results suggest that disease severity has the strongest posterior support for inclusion, followed by biomarker status. Age has substantially less posterior support.

Model-Averaged Regression Coefficients

A subtle issue arises when averaging regression coefficients. Suppose the biomarker coefficient is estimated under each model containing the biomarker. A model-averaged coefficient can be constructed by assigning a value of zero to models that exclude the variable and weighting the coefficient by the posterior model probability.

For a predictor \(X_j\), this gives:

where:

for models that do not contain \(X_j\).

This is sometimes called model averaging with zero coefficients. It naturally reflects uncertainty about whether the variable belongs in the model.

A Worked Coefficient Example

Suppose the biomarker coefficient is estimated as follows:

The model-averaged coefficient is:

which gives approximately:

Notice that the model-averaged coefficient is smaller than the coefficient in any individual model containing the biomarker.

That occurs because some posterior probability is assigned to models in which the biomarker is absent and its coefficient is therefore zero.

Conditional vs. Unconditional Effects

This distinction is particularly important when reporting BMA results. A coefficient conditional on inclusion answers:

An unconditional model-averaged coefficient answers a different question:

Both can be useful, but they should not be confused.

Bayesian Model Averaging for Prediction

Prediction is one of the most natural applications of BMA. Suppose we want to predict the probability of response for a future patient with covariate vector \(x\). Under model \(M_k\), the posterior predictive probability is:

BMA averages the prediction over the posterior distribution of models:

This is an important advantage of BMA. Rather than making a prediction conditional on one selected model, the prediction accounts for uncertainty about which model generated the data.

Clinical Trial Applications

Bayesian Model Averaging can arise in several areas of clinical research.

1. Prognostic modeling

When numerous baseline characteristics could predict an outcome, BMA can incorporate uncertainty about which predictors are genuinely useful.

2. Treatment-effect heterogeneity

Suppose treatment effect may depend on disease severity, biomarker status, baseline risk, and previous treatment. Rather than committing prematurely to a single interaction structure, BMA can average across plausible models.

3. Survival analysis

Different covariate structures can be considered in proportional hazards, accelerated failure-time, or other survival models. BMA can account for uncertainty about the prognostic model while producing model-averaged inference or prediction.

4. Dose-response modeling

Dose-response analyses frequently involve uncertainty about the appropriate functional form. Potential models may include:

• linear dose-response;
• quadratic dose-response;
• Emax models;
• sigmoid Emax models; and
• other nonlinear structures.

BMA provides a framework for averaging predictions or derived quantities across these candidate models.

5. Historical borrowing and prior-data models

When historical control data are available, analysts may consider multiple models for the degree of similarity between historical and current patients. BMA can be used when several borrowing structures are considered plausible.

6. Biomarker development

In exploratory biomarker analyses, many candidate predictors or interactions may be plausible. BMA can provide posterior inclusion probabilities rather than relying entirely on a single selected model.

BMA vs. Stepwise Selection

The important conceptual difference is that stepwise selection attempts to identify one model, whereas BMA treats the model itself as uncertain.

BMA vs. Bayesian Variable Selection

The terms Bayesian variable selection and Bayesian Model Averaging are closely related but are not identical.

Bayesian variable selection is the broader problem of determining which predictors should be included. BMA is a way of performing inference while averaging over the resulting model uncertainty.

For a small number of predictors, it may be possible to enumerate every possible model. With \(p\) candidate variables, the number of possible subsets is:

For example:

This exponential growth quickly makes exhaustive enumeration impractical.

Bayesian Model Averaging and MCMC

When the number of candidate models becomes very large, MCMC can be used to explore the model space. One approach is Markov Chain Monte Carlo Model Composition (MC³).

The algorithm moves between candidate models, proposing additions, removals, or changes to predictors. Models that have higher posterior probability are visited more frequently.

R Implementation

For a small model space, the conceptual workflow can be implemented by fitting the candidate models and combining their posterior probabilities. A simplified illustration is:

# Candidate models
models <- list(
  M1 = response ~ 1,
  M2 = response ~ severity,
  M3 = response ~ biomarker,
  M4 = response ~ age,
  M5 = response ~ severity + biomarker,
  M6 = response ~ severity + age,
  M7 = response ~ biomarker + age,
  M8 = response ~ severity + biomarker + age
)

# Each model would be fit under a Bayesian framework.
# Posterior model probabilities would then be obtained.

posterior_model_prob <- c(
  0.03, 0.18, 0.04, 0.02,
  0.39, 0.11, 0.05, 0.18
)

sum(posterior_model_prob)
# 1

# Posterior inclusion probability for severity
pip_severity <- sum(
  posterior_model_prob[c(2, 5, 6, 8)]
)

pip_severity
# 0.86

In practice, specialized Bayesian model-selection packages can perform the model fitting, posterior model exploration, and averaging. The exact implementation depends on the likelihood, prior distributions, model space, and whether the analysis concerns regression coefficients, prediction, survival outcomes, dose-response models, or another clinical endpoint.

A Simple BMA Calculation in R

Suppose model-specific posterior means for a parameter are available.

model_prob <- c(
  0.03, 0.18, 0.04, 0.02,
  0.39, 0.11, 0.05, 0.18
)

beta_biomarker <- c(
  0,
  0,
  0.80,
  0,
  0.72,
  0,
  0.61,
  0.69
)

bma_beta <- sum(
  model_prob * beta_biomarker
)

bma_beta
# approximately 0.54

The calculation is simply the weighted posterior average:

The important point is not the arithmetic. The important point is that the inference explicitly carries forward uncertainty about the model.

Model Averaging for Treatment Effects

BMA can also be used when the quantity of interest is a treatment effect. Suppose the treatment effect may depend on baseline disease severity. We could consider models such as:

• treatment only;
• treatment + severity;
• treatment + biomarker;
• treatment + severity + biomarker;
• treatment + severity × treatment; and
• other prespecified clinically plausible structures.

Each model produces a posterior distribution for the treatment effect. BMA then produces a model-averaged posterior distribution:

This equation is particularly important. BMA does not merely average posterior means. It can average the entire posterior distribution.

As a result, the model-averaged posterior can incorporate both:

• within-model parameter uncertainty; and
• between-model uncertainty.

Posterior Predictive Distribution

For clinical prediction, the model-averaged posterior predictive distribution is often more relevant than a model-averaged regression coefficient. For a future outcome \(Y_{\text{new}}\):

This equation says:

Advantages of Bayesian Model Averaging

• Accounts for model uncertainty. Inference does not depend entirely on one selected model.
• Provides posterior inclusion probabilities. These quantify the posterior support for individual predictors.
• Improves predictive uncertainty assessment. Model uncertainty can be incorporated into predictions.
• Provides a coherent Bayesian framework. Model uncertainty and parameter uncertainty can be propagated through the same probabilistic framework.
• Can reduce dependence on arbitrary model-selection thresholds. There is no requirement to declare one model definitively "correct."

Limitations and Practical Considerations

BMA is powerful, but it is not automatically superior in every clinical analysis.

1. The model space matters

BMA averages over the models that you specify. If important models are excluded from the candidate set, BMA cannot average over them.

2. Prior model probabilities matter

Posterior model probabilities depend on both the data and the prior model probabilities. When data are weakly informative, prior assumptions can have a meaningful influence.

3. Computational complexity can be substantial

With \(p\) predictors, the number of possible subset models is \(2^p\). Large model spaces may require MCMC or other stochastic search strategies.

4. Interpretation can be more complicated

A model-averaged coefficient does not have exactly the same interpretation as a coefficient from a prespecified regression model. Clinical interpretation should distinguish conditional effects, unconditional effects, posterior inclusion probabilities, and predictive quantities.

5. Causal inference requires additional care

BMA is fundamentally a model-averaging framework. It does not automatically convert an observational regression analysis into a valid causal analysis. Variables should be selected and modeled according to the causal question, study design, and relevant assumptions.

When Should a Clinical Statistician Consider BMA?

BMA is particularly attractive when all of the following are true:

What BMA Does Not Mean

Several common interpretations of BMA are incorrect.

A Practical Clinical Workflow

For a real clinical analysis, the workflow can be summarized as follows.

Key Takeaways

The central equation is:

The most important practical concepts are:

• Posterior model probability measures how much posterior support a candidate model receives.
• Posterior inclusion probability measures how much posterior probability is assigned to models containing a particular predictor.
• Model-averaged inference incorporates uncertainty about which model is appropriate.
• Model-averaged prediction incorporates model uncertainty into predictions for future patients.
• BMA does not remove modeling assumptions. The candidate model space and prior distributions remain important.

For clinical statisticians, the most important conceptual shift is simple:

Bayesian Model Averaging turns that second question into a formal statistical analysis.

References

Hoeting, J.A., Madigan, D., Raftery, A.E., & Volinsky, C.T. (1999). Bayesian model averaging: A tutorial. Statistical Science, 14(4), 382–417.

Madigan, D. & Raftery, A.E. (1994). Model selection and accounting for model uncertainty in graphical models using Occam's window. Journal of the American Statistical Association, 89, 1535–1546.

Raftery, A.E., Madigan, D., & Hoeting, J.A. (1997). Bayesian model averaging for linear regression models. Journal of the American Statistical Association, 92, 179–191.

Hinne, M., Gronau, Q.F., van den Bergh, D., & Wagenmakers, E.J. (2020). A conceptual introduction to Bayesian model averaging. Advances in Methods and Practices in Psychological Science.

George, E.I. & McCulloch, R.E. (1993). Variable selection via Gibbs sampling. Journal of the American Statistical Association, 88, 881–889.