What Bootstrap Stability Measures
A selection procedure returns a model, but a small change in the observed data can sometimes return a different one. This is especially relevant when candidate predictors carry overlapping information: one sample may favor one predictor, while another favors a correlated alternative. Bootstrap stability describes how much a specified selection procedure varies when it is applied repeatedly to bootstrap versions of the observed data.
The object being assessed is the whole selection procedure, not just the final fitted GLM. The procedure includes the candidate models, any screening or search steps, the selection criterion, constraints such as interaction hierarchy, and rules for ties or fitting failures. Changing any of these can change the stability results.
For independent observational units sampled from a common population, the ordinary nonparametric bootstrap draws \(n\) units with replacement from the \(n\) observed units. In a GLM with observed responses and predictors, each resampled unit is the whole response-and-predictor row. The same resampling indices must be used for every variable in a row, so the relationships among that row's measurements are preserved.
A bootstrap replicate is not just a refit of the original selected model. It repeats the selection pipeline: apply any data-dependent screening, search the specified candidate set, enforce the same constraints, and select using the same rule. If a transformation or other preprocessing choice is learned from the data, it must be re-estimated within each replicate. Otherwise, the procedure being assessed is not the procedure originally used.
Frequency Summaries
Suppose \(B\) bootstrap replicates produce selected models \(M_1^*,\ldots,M_B^*\), and let the candidate predictors be \(x_1,\ldots,x_p\). Define \(I_{bj}=1\) if replicate \(b\)'s selected model contains \(x_j\), and \(I_{bj}=0\) otherwise. Then the inclusion frequency for predictor \(x_j\) is
The frequency is between zero and one. A value near one means that this predictor appears in most selected models produced by this particular bootstrap procedure; a value near zero means it rarely appears. Frequencies need not sum to one, because a selected model can contain several predictors, or none.
It is also useful to report how often pairs of predictors occur together. For predictors \(x_j\) and \(x_k\), the co-selection frequency is \(B^{-1}\sum_b I_{bj}I_{bk}\). It distinguishes frequent individual inclusion from frequent inclusion in the same model. For example, two predictors might each be selected often but seldom be selected together, as can happen when the procedure tends to choose between alternatives.
A second summary is the frequency of each distinct selected model. If \(\mathcal{M}\) is the set of possible models, define \(\hat{q}(M)=B^{-1}\sum_b \mathbf{1}\{M_b^*=M\}\). Model frequencies show whether the procedure repeatedly returns one model or spreads its selections across several alternatives. They can be accompanied by model size, the number of predictors in a selected model, to show whether variation is mostly among models of similar complexity or also involves substantial changes in size.
Worked Example: Inclusion and Co-Selection Frequencies
A logistic-regression selection procedure is applied to \(B=100\) bootstrap samples. The candidate predictors are \(x_1,x_2,x_3\). The selected-model counts are: \(\{x_1,x_2\}\), 42; \(\{x_1\}\), 25; \(\{x_2\}\), 18; \(\{x_1,x_3\}\), 10; and \(\{x_3\}\), 5. The counts sum to \(42+25+18+10+5=100\).
Predictor \(x_1\) occurs in the first, second, and fourth models, so its inclusion frequency is \((42+25+10)/100=0.77\). Predictor \(x_2\) occurs in the first and third models, giving \((42+18)/100=0.60\). Predictor \(x_3\) occurs in the fourth and fifth models, giving \((10+5)/100=0.15\).
The only model containing both \(x_1\) and \(x_2\) is \(\{x_1,x_2\}\), so their co-selection frequency is \(42/100=0.42\). Thus \(x_1\) and \(x_2\) are each selected fairly often, but they appear together in fewer than half of the replicates. The most frequent exact model is \(\{x_1,x_2\}\), selected in 42% of replicates; no single model is selected by a majority.
What the Frequencies Do—and Do Not—Mean
Bootstrap frequencies describe the behavior of an algorithm under resampling from the empirical distribution represented by the observed data. They are conditional on that data and on all choices that define the selection procedure. A high inclusion frequency is evidence of stability under this resampling exercise; it is not the probability that the predictor is truly associated with the outcome, nor a confidence level for its coefficient.
These summaries also depend on how predictors are represented and on the candidate set. If two predictors carry similar information, small changes in the resampled data may cause the procedure to choose one rather than the other. Each can then have a moderate inclusion frequency even though the procedure often selects at least one of them. Reporting co-selection or group-level inclusion can help describe this pattern, provided the group was defined in a scientifically meaningful way.
The bootstrap scheme matters too. Resampling rows is appropriate when rows are independent observational units sampled from a population. If the design instead fixes predictor values and models response variation conditional on those values, a response-resampling scheme based on an explicitly stated fitted model may target a different question. The resampling unit must respect the sampling design; resampling individual rows is not automatically appropriate for dependent observations.
Proof. Conditional on \(D_n\), the indicators \(I_{1j},\ldots,I_{Bj}\) are independent and identically distributed Bernoulli random variables with success probability \(\pi_j(D_n)\): each is one precisely when the selected model in that replicate contains \(x_j\). They have finite expectation, \(E[I_{bj}\mid D_n]=\pi_j(D_n)\). The strong law of large numbers, applied under this conditional distribution, gives
The limit concerns Monte Carlo approximation of the conditional bootstrap inclusion probability. It does not claim that the probability itself is a population truth probability or that it converges to a particular quantity as the original sample size changes. \(\square\)
Derivation. Each inclusion indicator has conditional variance \(\pi_j(D_n)\{1-\pi_j(D_n)\}\). Independence makes the variance of their sum \(B\) times this amount. Dividing the sum by \(B\) divides its variance by \(B^2\), so
Substitution of \(\hat{\pi}_j\) for \(\pi_j(D_n)\) yields the stated estimate. This quantifies simulation error from using a finite number of replicates; it does not quantify uncertainty from having observed only \(D_n\).
Worked Example: Choosing the Number of Replicates
In the previous example, \(x_1\) was included in 77 of 100 replicates, so \(\hat{\pi}_1=0.77\). The estimated Monte Carlo standard error is
Thus, even if the conditional bootstrap inclusion probability stays fixed, the estimated frequency from only 100 replicates has appreciable simulation variation. If the analyst increases the number to 400 and the estimated frequency remains near 0.77, the corresponding estimated standard error is approximately \(\sqrt{0.77(0.23)/400}=0.0210\), half as large. Increasing \(B\) improves Monte Carlo precision, but it cannot make the observed data more informative.
A Complete Stability Assessment
A useful assessment reports enough detail for another analyst to understand what was resampled and what was repeated. State the observational unit and resampling scheme, the number of replicates, the candidate set, the selection rule, and any constraints or tie-breaking rules. Also say how nonconvergent fits or other failures were handled. Silently dropping failed replicates changes the reported frequencies to describe only the replicates that succeeded.
Fix the candidate GLMs, selection criterion, search strategy, constraints, data-dependent preprocessing, and rules for ties or failures.
Resample the appropriate independent observational units using the stated scheme, keeping each unit's response and predictors together when resampling rows.
Reapply preprocessing and the complete search and selection procedure; record the selected model, including any prescribed handling of a failure.
Report inclusion, co-selection, and exact-model frequencies, with Monte Carlo error where useful; interpret them as conditional procedure stability.
Worked Example: Stability of a Poisson GLM Search
Suppose a study models a count outcome using Poisson GLMs, with candidate predictors age, treatment, and an age-by-treatment interaction. The procedure uses AIC, requires both main effects whenever the interaction is included, and is applied to 500 row-bootstrap samples. Assume every replicate yields a selected model. The recorded counts are: intercept only, 30; age, 110; treatment, 80; age and treatment, 230; and age, treatment, and their interaction, 50. These counts sum to 500.
Age is included in the age-only, age-and-treatment, and interaction models. Its inclusion frequency is \((110+230+50)/500=390/500=0.78\). Treatment is included in the treatment-only, age-and-treatment, and interaction models, giving \((80+230+50)/500=360/500=0.72\). The interaction frequency is \(50/500=0.10\). Age and treatment are selected together in the two models containing both, so their co-selection frequency is \((230+50)/500=0.56\).
The exact model containing age and treatment without the interaction occurs most often, at \(230/500=0.46\). The results indicate that the two main effects recur frequently, while the interaction is rarely chosen under this particular procedure. They do not establish that the interaction is absent from the population, or that either main effect is causally important. The frequency estimates are also specific to the candidate set, AIC search, hierarchy constraint, resampling scheme, and observed data.
Common Pitfalls
A common mistake is to interpret an inclusion frequency of 0.80 as an 80% probability that a predictor is truly relevant. That reasoning confuses the conditional variability of a selection algorithm with uncertainty about a scientific hypothesis. The correction is to report the frequency as: “this procedure selected the predictor in 80% of these bootstrap replicates,” and to keep substantive evidence and inferential claims separate.
Another mistake is to refit only the model chosen on the original data in each bootstrap replicate. That assesses variation in estimates for a fixed model, not stability of model selection. To assess selection stability, the candidate search and selection rule must be repeated in every replicate. Likewise, a high inclusion frequency alone can hide model variation: report the exact-model frequencies, or other summaries of which predictors tend to appear together.
Check Your Understanding
Use these questions to check your interpretation of bootstrap stability.
- Why must a stability analysis repeat the complete selection pipeline rather than refit only the model selected from the original data?
- If a predictor appears in 72 of 120 selected models, what is its inclusion frequency?
- Two predictors each have high inclusion frequencies but a low co-selection frequency. What kind of selection pattern might this indicate?
- What quantity does the Monte Carlo standard error of an inclusion frequency describe, and what uncertainty does it not describe?
- Why is it incorrect to interpret an inclusion frequency of 0.80 as an 80% probability that the predictor is truly associated with the outcome?