Tutorials › High-Dimensional Genomic Data Analysis › Subsampling Frequencies for Selection Stability

Penalized Regression for p Much Larger Than n · Tutorial 542 of 795

Subsampling Frequencies for Selection Stability

You will learn to define a reproducible subsampling-and-fitting rule, compute selection frequencies, and interpret their Monte Carlo variation without treating them as probabilities that genes are truly associated.

Advanced 10 min read

What You'll Learn

  • Define a selection indicator for each feature and subsample.
  • Construct selection frequencies from repeated penalized model fits.
  • Distinguish a fixed fitting rule from one that retunes the penalty in every subsample.
  • Interpret frequencies relative to the chosen subsampling and fitting procedure.
  • Quantify Monte Carlo variation from a finite number of repeated subsamples.
  • Diagnose unstable or incomparable frequency calculations.

Why Repeat the Fit?

A penalized model fitted once to genomic data gives one selected-feature set. That set can depend on which individuals happened to be observed, as well as on modeling choices such as the penalty strength. Repeating the fit on different subsets of the observed individuals makes one aspect of that sensitivity visible: how often does each feature enter the fitted model under a specified procedure?

This tutorial develops that calculation. It builds on the earlier material on selection stability as a feature-level quantity, but focuses on how to construct the frequencies reproducibly and how to interpret the numbers. The frequency is a property of a defined subsampling-and-fitting procedure. It is not, by itself, a probability that a gene is biologically relevant or a test of association.

Definition: Let the observed data be \(D=\{(y_i,x_i):i=1,\ldots,n\}\), with \(p\) candidate features. Specify a subsample size \(m\), a rule for drawing a subset \(S_b\) of \(m\) individuals, and a penalized fitting rule. For feature \(j\), let \(I_{bj}=1\) if feature \(j\) is selected in the fit on subsample \(b\), and \(I_{bj}=0\) otherwise. With \(B\) repeated subsamples, its selection frequency is $$ \hat{\pi}_j=\frac{1}{B}\sum_{b=1}^{B} I_{bj}. $$ It is a number between 0 and 1; for example, \(\hat{\pi}_j=0.70\) means that feature \(j\) was selected in 70% of these \(B\) fits.

Specify the Procedure Before Counting

A selection frequency is only interpretable when the rule that produced each indicator is clear. A complete specification should state the response and candidate predictors, the model family, the penalty type, the subsample size and sampling scheme, the preprocessing, the penalty-selection rule, and the criterion for calling a coefficient selected. These details define what is being repeated.

For example, in a penalized linear model, a common criterion is that feature \(j\)'s fitted coefficient is nonzero. A numerical implementation may return tiny values near zero, so the analysis should use a stated numerical tolerance if exact zeros are not reliable. In penalized logistic regression, the same inclusion rule can be applied to the fitted predictor coefficients. A nonzero intercept is not counted as selection of a genomic feature.

Preprocessing must also be part of the rule. If predictors are standardized, decide whether their means and scales are calculated separately within each subsample or once using the full data. Calculating them within each subsample makes each fit depend only on that subsample. Using full-data summaries instead is a different procedure; it can be reasonable in a defined analysis pipeline, but must be stated and applied consistently. Any outcome-informed filtering or feature screening should likewise be repeated within each subsample if the goal is to include that screening step in the stability assessment.

The penalty can be handled in at least two ways. One can fix its value in advance for every subsample, or reselect it in each subsample using a specified rule, such as cross-validation within that subsample. These choices answer different questions. With a fixed value, frequencies describe selection under one common penalty setting. With retuning, they describe the combined behavior of subsampling and the retuning procedure. A penalty chosen once on the full data and then reused is also a defined option, but its frequencies do not account for variation in choosing the penalty.

1
Define the target procedure.
Choose the model, penalty, subsample size and sampling scheme, preprocessing, tuning rule, and selection criterion before computing frequencies.
2
Draw and fit.
For each repetition, draw a subset of individuals and apply the specified fitting rule to those data.
3
Record feature-level indicators.
For every candidate feature, record 1 if selected and 0 otherwise. Keep features that are not selected in any fit in the record as well.
4
Average and report.
Sum each feature's indicators and divide by \(B\); report the procedure and number of repetitions alongside the frequencies.

Worked Frequency Calculations

Worked Example: Counting Selections Across Eight Fits

Suppose a hypothetical genomic data set has \(n=40\) individuals and three candidate features, \(G_1,G_2,G_3\). The analyst draws \(B=8\) subsets of \(m=20\) individuals without replacement within each draw. The same prespecified penalized fitting rule is applied in every subset. The selected features are:

Fit 1: \(G_1,G_2\); fit 2: \(G_1\); fit 3: \(G_1,G_2,G_3\); fit 4: \(G_1,G_3\); fit 5: \(G_1,G_2\); fit 6: \(G_2\); fit 7: \(G_1,G_2,G_3\); fit 8: \(G_1,G_3\).

For \(G_1\), the indicator is 1 in fits 1, 2, 3, 4, 5, 7, and 8, so its count is 7: $$ \hat{\pi}_1=\frac{7}{8}=0.875. $$ For \(G_2\), the count is 5, giving \(\hat{\pi}_2=5/8=0.625\). For \(G_3\), the count is 4, giving \(\hat{\pi}_3=4/8=0.500\). Thus \(G_1\) was selected in 87.5% of these fits, \(G_2\) in 62.5%, and \(G_3\) in 50.0%.

These results describe the eight fits under the stated procedure. The observed ordering does not establish that \(G_1\) has a stronger biological effect than the other features, nor does it show that any feature is truly associated with the outcome.

Worked Example: Interpreting Frequencies Under Two Tuning Rules

Consider a second hypothetical analysis with 30 individuals and four candidate genes. Six subsets of 18 individuals are drawn. Under a rule that fixes the penalty before subsampling, the recorded selection counts are 5, 3, 1, and 0 for genes \(A,B,C,D\), respectively. The frequencies are: $$ \hat{\pi}_A=\frac{5}{6}=0.833,\quad \hat{\pi}_B=\frac{3}{6}=0.500,\quad \hat{\pi}_C=\frac{1}{6}=0.167,\quad \hat{\pi}_D=\frac{0}{6}=0. $$

Now suppose the analyst instead retunes the penalty by cross-validation separately in each subset. The counts become 4, 4, 2, and 1. The corresponding frequencies are \(4/6=0.667\), \(4/6=0.667\), \(2/6=0.333\), and \(1/6=0.167\). The arithmetic is the same, but the procedures differ: the second set incorporates variation from both the sampled individuals and the within-subset tuning.

Neither set of frequencies is universally preferable. The appropriate one depends on the question. If the aim is to describe selection at a fixed penalty, use the first rule. If the intended modeling workflow retunes on the available data, the second rule may better represent that workflow. A report should not present one as though it were the result of the other.

What Variation Does the Frequency Summarize?

For a fixed observed data set and a fully specified random subsampling-and-fitting procedure, define \(\pi_j\) as the probability that feature \(j\) is selected when a subset is drawn according to that procedure. This probability is conditional on the observed data and the chosen algorithm. The calculated frequency \(\hat{\pi}_j\) approximates it through repeated draws. It does not automatically describe how results would vary across new studies or a new population.

A short derivation clarifies the role of the number of repetitions. If the \(B\) subsamples are drawn independently from the specified subsampling mechanism, then for a fixed feature \(j\), the indicators \(I_{1j},\ldots,I_{Bj}\) are independent Bernoulli random variables with selection probability \(\pi_j\). Therefore \(E(I_{bj})=\pi_j\) and \(\operatorname{Var}(I_{bj})=\pi_j(1-\pi_j)\). By linearity of expectation and independence, $$ E(\hat{\pi}_j)=\frac{1}{B}\sum_{b=1}^{B}E(I_{bj})=\pi_j,\qquad \operatorname{Var}(\hat{\pi}_j)=\frac{\pi_j(1-\pi_j)}{B}. $$

Replacing the unknown \(\pi_j\) with \(\hat{\pi}_j\) gives the estimated Monte Carlo standard error \(SE\approx\sqrt{\hat{\pi}_j(1-\hat{\pi}_j)/B}\). This quantifies simulation noise from using a finite number of repetitions under the specified procedure. It is not a confidence interval for a population-level selection probability and does not include uncertainty from observing a particular genomic study.

Worked Example: Checking Monte Carlo Precision

Suppose a hypothetical analysis uses \(B=100\) independently drawn subsamples, and gene \(H\) is selected in 72 fits. Its frequency is \(\hat{\pi}_H=72/100=0.72\). The estimated Monte Carlo standard error is: $$ SE\approx\sqrt{\frac{0.72(1-0.72)}{100}} =\sqrt{0.002016} =0.0449. $$

This indicates that repeating the random subsampling experiment with 100 draws can produce noticeable Monte Carlo variation around the procedure's conditional selection probability. If the number of independent draws were increased to 400 while the underlying probability remained similar, the estimated standard error would be approximately \(\sqrt{0.72(0.28)/400}=0.0224\), half as large. More repetitions improve the numerical precision of the frequency, but cannot repair a poorly chosen model or make the result generalize beyond the specified data and procedure.

Diagnostics, Alternatives, and Limitations

Several checks make frequency calculations more informative and help identify implementation errors:

  • Verify the denominator. Every feature frequency must use the same intended number \(B\) of valid fits. If a fit fails, document how it is handled; silently dropping failed repetitions can change the procedure.
  • Check the selection rule. Confirm that exact zeros or the stated numerical tolerance are applied consistently. In a model with several coefficient representations, define precisely which coefficient counts as inclusion.
  • Audit the pipeline. Confirm that subsamples contain the intended individuals, preprocessing and tuning follow the prespecified rule, and no outcome-dependent operation is inadvertently carried out on the full data when it was meant to be repeated within subsets.
  • Look for Monte Carlo instability. Compare estimates after additional independent repetitions or across batches. A frequency that changes materially as repetitions are added is not yet numerically settled.
  • Report design choices. State \(m\), \(B\), whether sampling is with or without replacement, the penalty and tuning rule, and the selection criterion. Without these, the number is difficult to reproduce or compare.

A single fit to all observations is the simplest alternative: it gives one selected set, not a frequency, and conceals sensitivity to the subsampling perturbations being studied. Repeated subsampling makes that sensitivity visible but requires additional computation and depends on the chosen subset size and fitting rule. Repeated cross-validation primarily evaluates predictive performance across folds; its fold-wise prediction errors are not interchangeable with feature selection frequencies. Bootstrap resampling is another possible perturbation scheme, but it changes the sampling mechanism and therefore the quantity being summarized.

Frequencies can be affected by correlated genomic predictors: different correlated features may be selected in different subsets even when the model's predictive behavior is similar. A low frequency does not prove that a feature has no signal, and a high frequency does not prove that its effect is causal or independently reproducible. The frequencies also depend on the candidate feature set, preprocessing, penalty, and sample size. They are best reported as a descriptive summary of a specified fitting procedure, with those dependencies made explicit.

Common Mistakes

  • Wrong: “A frequency of 0.8 means there is an 80% probability this gene is truly associated.” Correction: It means the gene was selected in 80% of the fits in the stated subsampling experiment. No probability of biological truth follows from this calculation alone.
  • Wrong: “All frequencies are comparable, even if we changed the tuning rule between analyses.” Correction: Frequencies summarize the procedure used; changing the rule changes what is being summarized.
  • Wrong: “A frequency of zero proves the feature is irrelevant.” Correction: It was not selected in the observed repetitions. A different data set, feature set, or fitting rule could produce a different result.
  • Wrong: “Increasing \(B\) addresses every source of uncertainty.” Correction: Increasing \(B\) reduces Monte Carlo error in approximating the procedure's selection probability. It does not remove uncertainty from the observed study or correct model misspecification.
Key takeaway: A selection frequency is the fraction of repeated fits that include a feature under a defined subsampling and modeling procedure. Its meaning depends on that procedure, and its Monte Carlo precision depends on the number of repetitions.

Check Your Understanding

Try each problem, then check the solutions below.

  1. A feature is selected in 18 of 30 fits. Calculate its selection frequency and interpret it in context.
  2. Why can’t a frequency of 0.60 be interpreted as a 60% probability that the feature is truly associated with the outcome?
  3. In 100 independent subsampling repetitions, a feature has frequency 0.50. Estimate its Monte Carlo standard error.
  4. What changes in interpretation if the penalty is retuned within every subsample rather than fixed across subsamples?

Solutions. (1) \(\hat{\pi}=18/30=0.60\); the feature was selected in 60% of the 30 fits under the stated procedure. (2) The frequency describes selection by an algorithm on subsets of the observed data; it is not a probability model for whether the feature is truly associated. (3) \(SE\approx\sqrt{0.50(1-0.50)/100}=\sqrt{0.0025}=0.050\). This measures Monte Carlo variation, not population-level uncertainty. (4) Retuning adds variation from the penalty-selection step, so the frequency summarizes the combined subsampling-and-retuning procedure rather than selection at one fixed penalty.

Clinical Biostats