Bayesian Model Averaging Across Regression Specifications

Why Average Across Regression Specifications?

A regression analysis can be sensitive to specification: which covariates are included, which functional form is used, or how an outcome is modeled. Selecting one specification and then reporting its predictions or coefficient summaries treats that choice as settled. Bayesian model averaging instead carries uncertainty about the candidate specification into the summaries. Each model contributes in proportion to its posterior probability given the observed data.

The candidate models must refer to the same observed data and outcome for their marginal likelihoods to be compared. The posterior probabilities are conditional on the specified candidate set and on the prior probabilities assigned to its models. They do not establish that one model is objectively true, nor do they remove the need to choose a defensible set of specifications.

Definition: Let $M_1,\ldots,M_K$ be a finite set of candidate regression models, with prior model probabilities $P(M_k)$ summing to one. Given observed data $y$, the posterior model probability is $$ w_k=P(M_k\mid y)=\frac{p(y\mid M_k)P(M_k)}{\sum_{j=1}^{K}p(y\mid M_j)P(M_j)}. $$ Here $p(y\mid M_k)$ is the marginal likelihood, which averages the likelihood over the parameter prior within model $M_k$. The weights satisfy $w_k\geq 0$ and $\sum_{k=1}^{K}w_k=1$.

The model probabilities combine evidence from the data with prior model probabilities. For example, if marginal likelihoods are expressed relative to a reference model, their ratios are Bayes factors; multiplying those ratios by the corresponding prior probabilities and normalizing gives posterior model probabilities. Since marginal likelihoods depend on the parameter priors, model averaging inherits that prior sensitivity.

For any target that is defined in every candidate model, averaging has a direct interpretation. Conditional on $M_k$, the target has a posterior distribution given $y$. Its model-averaged posterior distribution is the mixture of those distributions, with mixture weights $w_k$. A prediction target must be the same across models: for instance, the outcome for the same covariate values, rather than unrelated predictions from different settings.

Model-Averaged Predictions

Let $Y^*$ denote a future outcome at a specified covariate vector. Within model $M_k$, the posterior predictive distribution $p(y^*\mid y,M_k)$ incorporates uncertainty about that model's parameters as well as the future outcome's error. Averaging across models gives

p(y^*\mid y)=\sum_{k=1}^{K}w_k\,p(y^*\mid y,M_k).

This is a mixture distribution, not generally a single distribution from the model with the largest weight. Its mean is the weighted average of the model-specific predictive means. Its variance includes both the uncertainty within each model and differences between model-specific predictive means.

Worked Example: Averaging Forecasts and Predictive Variance

Suppose three hypothetical regression specifications have posterior model probabilities $w_1=5/12$, $w_2=1/2$, and $w_3=1/12$. For the same future outcome and covariate values, their posterior predictive means are 10, 12, and 9 units, and their posterior predictive variances are 4, 9, and 1 square units. The weights sum to $5/12+6/12+1/12=1$.

The model-averaged predictive mean is

E(Y^*\mid y)=\frac{5}{12}(10)+\frac{1}{2}(12)+\frac{1}{12}(9) =\frac{50+72+9}{12} =\frac{131}{12} =10.92\text{ units}.

To find the predictive variance, first average the within-model variances, then add the weighted variation among the predictive means:

\begin{aligned} \operatorname{Var}(Y^*\mid y) &=\sum_{k=1}^{3}w_k\left[\operatorname{Var}(Y^*\mid y,M_k) +\left\{E(Y^*\mid y,M_k)-E(Y^*\mid y)\right\}^2\right]\\ &=\left[\frac{5}{12}(4)+\frac{1}{2}(9)+\frac{1}{12}(1)\right]\\ &\quad+\left[\frac{5}{12}\left(10-\frac{131}{12}\right)^2 +\frac{1}{2}\left(12-\frac{131}{12}\right)^2 +\frac{1}{12}\left(9-\frac{131}{12}\right)^2\right]\\ &=6.25+1.243 =7.493\text{ square units}. \end{aligned}

The model-averaged predictive standard deviation is $\sqrt{7.493}=2.737$ units. The first variance component reflects uncertainty inside the models; the second reflects disagreement among their predictive means. Selecting $M_2$, the most probable model, would instead give a forecast mean of 12 units and predictive variance of 9 square units. Those are summaries conditional on $M_2$, not the model-averaged summaries.

Why the Variance Formula Works

The variance decomposition follows by conditioning on which model generated the posterior distribution. Write $\mu_k=E(Y^*\mid y,M_k)$ and $V_k=\operatorname{Var}(Y^*\mid y,M_k)$, and let $\mu=\sum_k w_k\mu_k$. The conditional variance identity gives

\operatorname{Var}(Y^*\mid y) =E\{\operatorname{Var}(Y^*\mid y,M)\mid y\} +\operatorname{Var}\{E(Y^*\mid y,M)\mid y\} =\sum_k w_k V_k+\sum_k w_k(\mu_k-\mu)^2.

The first term averages uncertainty conditional on each model. The second is uncertainty about the conditional mean associated with the model. Omitting the second term understates uncertainty whenever the candidate models disagree about their predictive means. The same identity applies to a parameter or other scalar target that has a coherent definition across the models.

Averaging Parameter Summaries Carefully

If a scalar parameter $\theta$ has the same definition in every model, let $\mu_k=E(\theta\mid y,M_k)$ and $V_k=\operatorname{Var}(\theta\mid y,M_k)$. Its model-averaged mean and variance are

E(\theta\mid y)=\sum_k w_k\mu_k,\qquad \operatorname{Var}(\theta\mid y)=\sum_k w_k\left[V_k+(\mu_k-E(\theta\mid y))^2\right].

The mean alone can hide important features of the mixture. In particular, if models give quite different parameter values, the posterior can be multimodal or otherwise non-normal. A normal approximation based only on its mean and variance may then misrepresent posterior probabilities or credible intervals.

A parameter may not be present in every specification. In that case, do not average coefficient labels without first defining a common target. If “the effect of $x$” is interpreted as zero in specifications that exclude $x$, the model-averaged distribution includes that point value for those models. If instead the question concerns the coefficient only when $x$ is included, report the posterior distribution conditional on inclusion, renormalizing the weights over models that include $x$. These answer different questions and should be labeled accordingly.

Worked Example: A Coefficient Absent From One Model

Consider two hypothetical specifications for an outcome measured in log points. The first includes $x$, and has posterior model probability 0.65; its coefficient for $x$, measured per one-unit increase in $x$, has posterior mean $-0.80$ and variance $0.09$. The second excludes $x$, has probability 0.35, and assigns the coefficient for this common-target comparison the value zero. This zero is part of the specification being compared; it is not an estimated coefficient from the excluded model.

Under this definition, the model-averaged coefficient mean is

E(\beta_x\mid y)=0.65(-0.80)+0.35(0)=-0.520\text{ log points per unit of }x.

Its variance combines the included-model variance and the between-model variation:

\begin{aligned} \operatorname{Var}(\beta_x\mid y) &=0.65\{0.09+(-0.80+0.52)^2\} +0.35\{0+(0+0.52)^2\}\\ &=0.65(0.1684)+0.35(0.2704)\\ &=0.2041\text{ squared log points per squared unit of }x. \end{aligned}

The resulting posterior standard deviation is $\sqrt{0.2041}=0.4518$ log points per unit of $x$. By contrast, conditional on the model including $x$, the posterior mean remains $-0.80$ and the posterior variance remains $0.09$. The overall mean of $-0.520$ answers a model-averaged question with zero assigned when $x$ is excluded; it is not the conditional-on-inclusion coefficient mean.

Mixtures, Probabilities, and Model Selection

Posterior probabilities of events are also averaged across models. If $A$ is an event about a prediction or a common parameter, then

P(A\mid y)=\sum_k w_k P(A\mid y,M_k).

This event-based calculation uses the full model-specific posterior distributions; it does not generally follow from averaging point estimates. For example, the probability that a prediction exceeds a threshold can differ substantially from the probability calculated under the single most probable model.

Worked Example: A Threshold Probability Under Model Averaging

Suppose a hypothetical future outcome has posterior predictive distribution $N(0,1)$ under $M_1$, and $N(2,1)$ under $M_2$, where the second argument is the variance. Let the posterior model probabilities be 0.70 and 0.30. Find the model-averaged probability that the outcome exceeds 1 unit.

Under $M_1$, standardizing gives $P(Y^*>1\mid y,M_1)=P(Z>1)=0.1587$, where $Z\sim N(0,1)$. Under $M_2$, it gives $P(Y^*>1\mid y,M_2)=P(Z>-1)=0.8413$. Therefore

P(Y^*>1\mid y)=0.70(0.1587)+0.30(0.8413) =0.1111+0.2524 =0.3635.

Thus the model-averaged probability of exceeding 1 unit is about 0.364. Choosing $M_1$ alone would report 0.159, while choosing $M_2$ alone would report 0.841. Neither selected-model probability is the model-averaged result.

Model selection and model averaging therefore answer different questions. Selecting the model with the greatest posterior probability conditions subsequent summaries on that one specification. Averaging combines summaries or distributions across the candidate set. Selection may be appropriate for a decision that explicitly requires one specification, but a high posterior probability does not make the other models' uncertainty vanish.

Practical Computation and Common Pitfalls

For a finite set of models, the computation follows a clear sequence. Establish the candidate set and prior model probabilities; obtain marginal likelihoods on a comparable basis; normalize to posterior model probabilities; then combine the model-specific posterior distributions for the target of interest. Predictions must use the same future outcome and covariate values in every specification. Parameter summaries must refer to the same estimand or clearly state how models omitting a parameter are handled.

If posterior simulations are available within each model, the same mixture can be approximated by weighting model-specific draws according to $w_k$, or by drawing a model according to those probabilities and then drawing a posterior quantity conditional on that model. For numerical summaries, preserve the mixture structure where possible: averaging credible-interval endpoints does not in general produce a credible interval for the mixture. Mixture quantiles should be calculated from the combined distribution.

Pitfall: reporting the best-supported model as if it were the average. The largest $w_k$ identifies the most probable candidate in the set, but does not equal one and does not erase the remaining weights. Use the weighted mixture for model-averaged conclusions.
Pitfall: averaging only conditional variances. This ignores disagreement among model-specific means. Include the between-model variance term.
Pitfall: averaging coefficients with different meanings. A coefficient that refers to different transformations or conditioning sets may not be a common target. Define a comparable estimand first, or report model-specific summaries.
Pitfall: treating posterior model probabilities as prior-free truth probabilities. They depend on the candidate set, prior model probabilities, and parameter priors used in the marginal likelihoods. State the model-averaging scope accordingly.

Key takeaway: Bayesian model averaging combines each model's posterior distribution using its posterior model probability. The resulting uncertainty includes both uncertainty within specifications and disagreement across them; it is not the same as conditioning on the single most probable specification.

Check Your Understanding

Use the definitions and calculations in this tutorial to answer the following.

Why must posterior model probabilities sum to one, and what quantities determine each weight?
For a prediction, what two components make up the model-averaged predictive variance?
If a coefficient is absent from some candidate models, what additional decision is needed before reporting a model-averaged coefficient mean?
Why is averaging the endpoints of model-specific credible intervals generally not a valid way to obtain a model-averaged credible interval?
How does a prediction probability under the most probable model differ from the model-averaged probability of the same event?