Why Use Inclusion Indicators?
Covariate selection in a nonlinear mixed-effects (NLME) model is not only a question of estimating the size of an effect. It is also a question of whether a candidate effect should be present in the model at all. A Bayesian inclusion indicator represents that choice explicitly: the indicator is zero when the effect is excluded and one when it is included. The resulting posterior distribution can describe uncertainty across candidate models rather than returning only one selected model.
This approach is useful when a set of scientifically motivated candidate effects is available and the analyst wants to quantify how support for each effect depends on both the data and prior assumptions. It does not remove the need to specify a plausible NLME model, examine its fit, or make an informed choice of candidate covariates.
Model and Prior Definition
Suppose subject \(i\) has observations \(y_i\), covariates \(x_{ij}\), and subject-specific random effects \(\eta_i\). A candidate covariate effect \(j\) may enter a population parameter through a specified relationship. For example, a covariate could affect the log of a typical clearance value. The particular functional form must be chosen before interpreting an inclusion indicator: the indicator concerns that proposed effect in that proposed form.
A common prior specification is \(\gamma_j \sim \mathrm{Bernoulli}(\pi_j)\). Thus, before observing the data, the prior probability of including effect \(j\) is \(\pi_j\). Conditional on inclusion, one might use a normal slab:
Here \(\beta_j\) is the covariate coefficient, and \(\tau_j\) controls the prior spread of included effects. A prior with a point mass at zero is different from a single continuous prior concentrated near zero: the former assigns positive prior probability to an exactly absent effect. In practice, scaling the covariate and choosing a slab scale that allows clinically plausible effects can make the prior easier to interpret.
The NLME likelihood accounts for the nonlinear observation model, residual variation, and random effects. If \(f(y_i\mid\eta_i,\theta,\gamma)\) is the conditional likelihood contribution for subject \(i\), the marginal likelihood for a model averages over that subject’s random effects:
Here \(\theta\) denotes the model’s other parameters, and \(y\) denotes the complete data. The model also requires priors for its continuous parameters, including relevant structural, random-effect, and residual parameters. Posterior inclusion probabilities are conditional on the specified candidate set, NLME likelihood, and prior distributions. They do not establish that the model is correctly specified.
Worked Example: Defining an Indicator for a Clearance Effect
Consider a hypothetical NLME model in which a subject’s typical clearance depends on body weight. Let \(WT_i\) be weight in kilograms, centered at 70 kg and scaled in 10-kg units. Specify
For this example, \(\gamma\) is the inclusion indicator, \(\beta\) is the weight effect, and \(\eta_{CL,i}\) is the subject-specific random effect on log clearance. Suppose the included-effect slab is \(N(0,0.25^2)\) and the prior inclusion probability is \(\pi=0.30\). If \(\gamma=0\), the covariate contribution is zero at every weight. If \(\gamma=1\), a one-unit increase in \(x_i\), corresponding to 10 kg, changes typical clearance by a factor of \(\exp(\beta)\), holding the other model terms fixed.
For illustration, if the coefficient were \(\beta=\log(1.25)=0.2231\), the ratio in typical clearance between 80 kg and 70 kg would be \(\exp(0.2231\times 1)=1.25\). That is a 25% higher typical clearance at 80 kg under this specified model. This calculation interprets a hypothetical coefficient; it does not imply that the posterior probability of inclusion is 25%, nor does it show that the data support the effect.
Posterior Inclusion Probabilities
The posterior inclusion probability (PIP) for effect \(j\) is the posterior probability that its indicator equals one, given the observed data and the specified model and priors:
If \(M_\gamma\) denotes a model defined by a complete vector of indicators, its posterior model probability is proportional to its marginal likelihood times its prior model probability. The PIP sums the posterior probabilities of all models that include effect \(j\):
A PIP of 0.80 means that, under the stated candidate set, likelihood, and priors, the posterior probability assigned to models containing the effect is 0.80. It is not the probability that the null hypothesis is true, a frequentist p-value, or the probability that the coefficient is positive. The last question concerns the posterior distribution of \(\beta_j\) conditional on inclusion.
Worked Example: Summing Posterior Model Probabilities
Suppose a hypothetical analysis considers two candidate effects, A and B, and reports the following posterior model probabilities:
| Model | Included effects | Posterior model probability |
|---|---|---|
| \(M_{00}\) | Neither A nor B | 0.10 |
| \(M_{10}\) | A only | 0.55 |
| \(M_{01}\) | B only | 0.15 |
| \(M_{11}\) | A and B | 0.20 |
The probabilities sum to \(0.10+0.55+0.15+0.20=1.00\). Effect A appears in \(M_{10}\) and \(M_{11}\), so its PIP is \(0.55+0.20=0.75\). Effect B appears in \(M_{01}\) and \(M_{11}\), so its PIP is \(0.15+0.20=0.35\).
The posterior gives more model probability to models containing A than to models containing B. This is not, by itself, a conclusion that A has a clinically important effect or that B has no effect. Those judgments also require the effect-size distribution, the model’s adequacy, and the scientific context.
Why Prior Inclusion Odds Matter
For a single candidate effect, compare an included model \(M_1\) with an excluded model \(M_0\). The Bayes factor \(BF_{10}\) is the marginal likelihood under \(M_1\) divided by the marginal likelihood under \(M_0\). Bayes’ rule gives posterior odds as prior odds multiplied by the Bayes factor:
To see the result directly, write \(P(M_k\mid y)=P(y\mid M_k)P(M_k)/P(y)\) for \(k=0,1\). Taking the ratio cancels the common denominator \(P(y)\), leaving the displayed expression. In a spike-and-slab model, the included-model marginal likelihood averages over the slab and the other uncertain parameters. Consequently, the Bayes factor depends on the slab prior as well as on the data and likelihood.
Worked Example: Changing the Prior Inclusion Probability
Suppose a hypothetical analysis produces \(BF_{10}=4\) for one covariate effect under a fixed slab prior. First use prior inclusion probability \(\pi=0.50\). The prior odds are \(0.50/0.50=1\), so the posterior odds are \(4\times1=4\). Converting odds to probability gives
Now set \(\pi=0.10\), keeping the same Bayes factor. The prior odds are \(0.10/0.90=0.1111\), and the posterior odds are \(4\times0.1111=0.4444\). Therefore,
With \(\pi=0.80\), the prior odds are \(0.80/0.20=4\); the posterior odds are \(4\times4=16\), giving \(P(M_1\mid y)=16/(1+16)=0.9412\). The data evidence, represented here by the fixed Bayes factor, is the same in all three calculations. The posterior probabilities differ because the prior inclusion odds differ. These values should therefore be reported with the relevant prior choices, not treated as prior-free summaries.
Prior Sensitivity, Diagnostics, and Alternatives
Sensitivity analysis should examine more than one prior inclusion probability when the choice is not strongly determined by the scientific context. It should also examine plausible slab scales. A very narrow slab can make the included model imply effects close to zero; a very diffuse slab can spread prior probability over implausibly large effects. Because the marginal likelihood averages fit over the slab, changing its scale can change the Bayes factor and therefore the PIP. The covariate’s units and scaling must be stated so that the slab’s implications are interpretable.
Useful checks include verifying that the candidate effects and functional forms were specified clearly, that the fitted NLME model describes the data adequately, and that the posterior calculations are numerically reliable. When posterior computation uses sampling, check convergence and mixing for both continuous parameters and inclusion indicators, and assess Monte Carlo error in reported PIPs. Poorly explored indicators can make an apparently precise PIP unreliable. Separately, residual patterns or implausible parameter estimates can reveal model inadequacy that a high PIP does not repair.
In a model space with several candidates, individual PIPs summarize inclusion across models; they do not necessarily identify one uniquely supported model. In particular, when candidate effects convey overlapping information, posterior support can be distributed among them. The PIP is also conditional on the candidate set: an omitted alternative cannot receive posterior support. These are reasons to present sensitivity results and model assumptions alongside the probabilities.
Stepwise likelihood-ratio or Wald procedures make inclusion decisions through repeated tests and specified thresholds. L1-penalized screening instead shrinks coefficients and can set some estimates to zero through a penalty. A spike-and-slab approach differs by assigning explicit posterior probability to inclusion and exclusion, and by expressing uncertainty across models. None of these approaches automatically resolves confounding, ensures predictive performance, or establishes a causal effect. A PIP threshold such as 0.5 is not a universal decision rule; the threshold and consequences of retaining or omitting an effect depend on the purpose of the model.
Common Pitfall
Correction: A PIP of 0.90 is a posterior probability of inclusion conditional on the candidate model space, likelihood, and prior choices. It neither establishes a causal effect nor removes dependence on the prior and model specification. Report the specified covariate relationship, prior inclusion probability, slab scale, and relevant sensitivity results.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions.
- In a spike-and-slab prior, what does \(\gamma_j=0\) mean, and what does \(\gamma_j=1\) mean?
- A candidate effect has posterior inclusion probability 0.65. What probability does this describe, and what does it not establish?
- For \(BF_{10}=3\) and prior inclusion probability \(\pi=0.25\), calculate the posterior inclusion probability for a single effect.
- Why can changing the slab scale change the posterior inclusion probability even when the observed data are unchanged?
Solutions. (1) \(\gamma_j=0\) fixes the coefficient at zero; \(\gamma_j=1\) assigns it the continuous slab prior. (2) It is the posterior probability that the indicator is one, conditional on the data, candidate set, likelihood, and priors. It is not a p-value, a probability that the null is true, or evidence of causation. (3) Prior odds are \(0.25/0.75=1/3\). Posterior odds are \(3\times(1/3)=1\), so the posterior probability is \(1/(1+1)=0.50\). (4) The slab scale changes the distribution of included effect sizes and thus the included model’s marginal likelihood, which changes its Bayes factor and posterior inclusion probability.