← All examples

Variable Selection

Which predictors belong in the model? The Bayesian answer is not a single winner but a posterior over the whole model space — inclusion probabilities and a ranked list of models, with honest uncertainty. This arc builds that machinery from scratch: stochastic search variable selection (SSVS) with George–McCulloch spike-and-slab priors; exact model enumeration with conjugate g-priors and Bayesian model averaging; VAR restriction search that selects a whole system's lag structure and error covariance at once; log-linear model selection that chooses which associations belong in a contingency table; and, when even the model's size is unknown, a trilogy of trans-dimensional samplers — a reversible-jump change-point model, a Dirichlet-process mixture, and free-knot splines — inferring the number of regimes, of components, and of knots. Further nodes extend the reach across outcome types — a logistic selector (Pólya–Gamma augmentation carries the whole Gaussian machinery over to a binary outcome) and a survival selector that handles censored time-to-event data — and across scale: a p ≫ n genomics horseshoe, where continuous global-local shrinkage replaces discrete indicators and a specialised sampler makes selection over thousands of predictors routine. Each is cross-checked across Python and R.

What an inclusion probability is, and where its zero point sits. Every method below reports the same headline number: the posterior inclusion probability, the share of the posterior weight sitting on models that contain a given predictor. If half the probability lands on models including weight, weight’s inclusion probability is 0.5. Three things about it are easy to get wrong. Its neutral point is not zero — it is whatever the prior said, and these pages give each predictor an even chance of being in or out, so 0.5 is the number a predictor starts with: above it the data have raised the odds, below it the data have lowered them, and 0.49 is not weak evidence for a predictor but evidence that went the other way. It is not an effect size and not a p-value — a coefficient can be large and excluded, or small and kept. And it is conditional: it asks whether a predictor earns its place given the others, not whether it moves with the outcome on its own. In the heart-disease example below, adiposity tracks the disease about as closely as LDL cholesterol does (r = +0.25 against +0.26) and is dropped anyway, at inclusion 0.08, because whatever it knows is already in the model through the predictors it shares its information with. That gap between marginal association and conditional necessity is the whole reason to run a selection at all.

Marginal likelihood, and the Occam penalty nobody writes down. To weigh one model against another, each needs a score: the marginal likelihood, the probability that model assigned to the data actually observed — averaged over every parameter value its prior allowed, before seeing anything. The averaging is what turns it into a selection criterion. A model with more free parameters can accommodate more datasets it might have seen, so it must spread the same total probability more thinly across them, and any single dataset receives less of it. Complexity is therefore penalised as a consequence of the integral, not added by hand the way AIC and BIC bolt on a parameter-count term — the automatic Occam balance the pages below refer to. The ratio of two marginal likelihoods is a Bayes factor. The price is that the integral runs over the entire prior, which makes it far more sensitive to that prior than any posterior mean would be; that is why the priors here are picked for tractability as much as realism, conjugate ones giving the integral in closed form and a Laplace approximation standing in where the likelihood is not conjugate.

Why the answer is a distribution over models, not a winner. Once every model carries a probability there is no need to choose one. Bayesian model averaging forms the quantity of interest — a coefficient, a prediction — as a weighted average across models, each weighted by its posterior probability, so a predictor present in half the credible models contributes half as much. Why bother is visible in Part 2: the single best of 512 models holds about 5% of the posterior, so a procedure reporting it alone discards 95% of what the data support. Stepwise selection compounds that by quoting the winner’s standard errors as though the search that produced the winner had never happened. Model averaging keeps the uncertainty about which model inside the answer instead of dropping it at the last step.

Collapsing, and why a change of dimension needs no Jacobian. One phrase recurs below — the coefficients are integrated out, or collapsed. It is the same operation each time: when the prior is conjugate the continuous parameters can be removed analytically, leaving a probability that depends only on the discrete structure — which predictors are in, where the change-points fall, how the points are grouped. Two things follow. The sampler then explores a far smaller space, which is why these purpose-built engines mix well where general-purpose ones labour. And when the model’s size changes there is no continuous parameter vector to be matched across the two dimensions, so the change-of-variables term — the Jacobian — that a reversible jump normally has to supply is not needed at all. That is what no Jacobian means in Parts 5 to 7, and it is why those samplers are as short as they are.