Bayesian Covariate Selection With Shrinkage Priors

Why Use Bayesian Shrinkage for Covariates?

A covariate search can involve many plausible effects, while the available pharmacokinetic data may provide limited information about each one. Estimates from a large search can be unstable: some coefficients may appear substantial because of sampling variation, and correlated covariates can make it difficult to distinguish their separate contributions. Bayesian shrinkage addresses this by combining the likelihood with prior distributions that favor small effects, while allowing the data to support larger effects when warranted.

For example, let \(x_i\) be a covariate for individual \(i\), and let \(\beta\) describe its effect on typical clearance:

Here \(CL_{\mathrm{pop}}\) is typical clearance at the reference value of the covariate, and \(\eta_i\) is the individual random effect. The coefficient \(\beta\) is a change in log clearance per one-unit increase in \(x_i\); \(\exp(\beta)\) is the corresponding multiplicative change in typical clearance. Bayesian shrinkage specifies a prior for \(\beta\), or for a collection of such coefficients, and uses it together with the model likelihood to obtain a posterior distribution.

Spike-and-Slab Priors and Inclusion Probability

For covariate \(j\), define an indicator \(\gamma_j\), where \(\gamma_j=1\) means the coefficient is included and \(\gamma_j=0\) means it is excluded. A basic spike-and-slab specification is

The inclusion probability \(\pi_j\) is the prior probability of inclusion. The slab variance \(\tau_j^2\) controls how broadly nonzero effects are distributed before observing the data. A smaller \(\tau_j\) expresses a stronger prior preference for modest effects among included covariates. The posterior inclusion probability, abbreviated PIP, is \(P(\gamma_j=1\mid D)\), where \(D\) denotes the data.

A useful way to calculate the PIP is through a Bayes factor. The Bayes factor \(BF_{10}\) compares the marginal likelihood of the data under inclusion, with the slab prior integrated over \(\beta_j\), to the likelihood under exclusion, where \(\beta_j=0\). Posterior odds equal prior odds multiplied by this Bayes factor:

This identity follows directly from Bayes' rule: the posterior probability of each model is proportional to its prior probability times its marginal likelihood. Dividing the inclusion-model expression by the exclusion-model expression gives the posterior-odds formula. The marginal likelihood under inclusion averages the likelihood over the slab; it is not simply the likelihood evaluated at the best-fitting coefficient.

Coefficient Estimates Within a Spike-and-Slab Model

The PIP answers a model-inclusion question; it does not by itself describe the size of the effect conditional on inclusion. In the normal summary-likelihood example, the posterior for \(\beta\) conditional on inclusion is normal. Completing the square in the product of the normal likelihood and normal slab gives posterior variance \(V\) and mean \(m\):

The factor multiplying \(\hat{\beta}\) in the posterior mean is \(\tau^2/(\tau^2+s^2)\), which lies between zero and one. This shows explicitly how the slab shrinks the estimate toward zero. The full posterior is a mixture: it has a point mass at zero with probability \(1-\mathrm{PIP}\), and the conditional slab posterior with probability \(\mathrm{PIP}\).

Continuous Shrinkage Priors

A continuous shrinkage prior replaces the exact spike with a distribution that can concentrate strongly near zero. A simple example is \(\beta\sim N(0,\tau^2)\); its posterior mean under the normal summary likelihood is \(\tau^2/(\tau^2+s^2)\hat{\beta}\). More elaborate global-local priors, such as the horseshoe, give each coefficient a local scale and share a global scale across coefficients. In one common form, conditional on the scales, \(\beta_j\sim N(0,\tau^2\lambda_j^2)\); the local scales \(\lambda_j\) have heavy-tailed priors. A small global scale encourages broad shrinkage, while large local scales can let some effects escape much of that shrinkage.

Because a continuous prior has no point mass at zero, it does not by itself select an exact subset of covariates. A posterior probability such as \(P(\beta_j>0\mid D)\) or the probability that \(\beta_j\) lies in a prespecified practical-effect interval can be reported, but neither is a PIP for exact inclusion. If a decision requires a selected subset, a rule must be stated separately, such as a prespecified threshold on a coefficient's posterior magnitude. Such a threshold is a decision rule, not a property that makes the continuous prior a spike-and-slab prior.

Interpretation, Diagnostics, and Limitations

A PIP is the posterior probability that a covariate is included, conditional on the observed data, the candidate model set, and the specified prior. For instance, a PIP of 0.80 says that the posterior probability of inclusion is 80% under those assumptions. It is not a p-value, a frequentist type I error rate, or proof that the covariate causes a change in clearance. A PIP also does not measure effect magnitude; report coefficient summaries alongside it when the scientific question concerns direction or size.

The result depends on the prior inclusion probability and the slab distribution. A low prior inclusion probability can reduce PIPs when many candidate covariates are considered; a diffuse slab can also change the marginal likelihood because it spreads prior mass over a wider range of effects. Check whether conclusions persist across scientifically defensible prior settings. For continuous shrinkage, inspect posterior coefficient distributions and the sensitivity of conclusions to global and local scale choices.

Correlated covariates require particular care. If two predictors carry similar information, the data may support their joint contribution without clearly identifying which one should receive inclusion probability. Individual PIPs can then be sensitive to prior structure, while the joint posterior pattern may be more informative. Bayesian shrinkage does not remove confounding or make poorly identified effects identifiable. Model adequacy, covariate measurement, and the pharmacokinetic structure still matter.

For complex nonlinear mixed-effects models, posterior computation can be challenging. Indicators may mix poorly when models are similar or when coefficients and random effects are strongly dependent; continuous shrinkage can also yield difficult posterior geometry. Examine convergence and effective sample sizes for coefficients, inclusion indicators, and important derived quantities, and compare chains or runs. Prior predictive checks help reveal whether proposed priors imply implausibly large changes in typical pharmacokinetic parameters. A well-mixing computation does not establish that the model or prior is scientifically appropriate.

The main advantage of Bayesian shrinkage is that it represents uncertainty about both effect sizes and, for spike-and-slab models, model inclusion. Its limitations are equally important: conclusions are conditional on candidate models and prior choices, computation can be demanding, and posterior inclusion is not a universal measure of biological importance. Report the prior specification, the model used, coefficient units, and the posterior summaries that support the conclusion.