Why Assess Stability After Variant Selection?
A joint model can select different variants when the data or analysis inputs change slightly. This is especially likely when variants are strongly correlated through linkage disequilibrium (LD): several candidates may carry much the same information about a signal, while their conditional estimates are imprecise. A selected variant can therefore change even when the locus continues to show evidence of association. Conversely, a signal that disappears under small, defensible changes deserves more caution than one that persists.
Stability analysis makes these dependencies visible. It repeats a specified selection procedure after controlled perturbations, such as resampling individuals, changing the candidate-variant set, or using plausible alternative LD estimates. The goal is not to find the one perturbation that gives the preferred answer. It is to describe which conclusions persist, which depend on the analysis choices, and what kind of stability is being assessed.
Specify the Analysis Before Perturbing It
Suppose a continuous phenotype \(Y_i\) is analyzed for individuals \(i=1,\ldots,n\). A typical locus-wide linear model is
where \(\mathbf{C}_i\) contains the covariates, \(G_{ij}\) is the genotype dosage for variant \(j\), and \(S\) is the set of variants included in the joint model. Under additive coding, dosage is typically 0, 1, or 2 copies of the counted allele. The coefficient \(\beta_j\) is the conditional association with the phenotype per additional counted allele, holding the other modeled variants and covariates fixed.
For the usual linear-model tests in the examples below, assume the mean model is correctly specified, errors have conditional mean zero and constant variance, and observations are independent with normally distributed errors for exact small-sample \(t\)-tests. If individuals are related or otherwise dependent, the analysis must instead use a model or variance estimate appropriate to that dependence. The same phenotype definition, covariates, genotype coding, selection rule, and significance threshold should be used across runs unless one of those elements is itself the planned perturbation.
A stability result is conditional on the chosen perturbation scheme. For example, resampling individuals assesses sensitivity to the observed sample; it does not assess sensitivity to a different ancestry composition unless that is deliberately examined. Candidate-set perturbations assess dependence on which variants are eligible for selection. LD perturbations assess dependence on the correlation information used by the analysis. These are different questions and should be reported separately.
When resampling people, keep each person’s phenotype, genotypes, and covariates together. If the design includes strata, families, sites, or other dependence, resampling should respect that structure. A candidate-set perturbation might omit variants that fail a prespecified quality criterion or compare plausible candidate sets. An LD perturbation should use defensible alternatives rather than arbitrary correlations; when correlation matrices are used, each should be a valid positive-semidefinite correlation matrix. In every case, rerun the full stated selection procedure, not merely the final test for a variant selected in the original data.
Useful Summaries of Selection Stability
Let \(S_b\) be the selected set in perturbation run \(b\), for \(b=1,\ldots,B\). For variant \(j\), define its availability indicator \(A_{jb}=1\) if it was eligible to be selected in run \(b\), and \(A_{jb}=0\) otherwise. Its selection frequency among runs where it was available is
Here \(I(j\in S_b)\) equals 1 when variant \(j\) is selected and 0 otherwise. Runs in which the variant was unavailable do not count as failures to select it. Report the number of available runs as well as \(\hat{f}_j\); a frequency based on two available runs is less informative than the same frequency based on many. If a variant was never available, its frequency is undefined.
A second summary compares two selected sets, for example the original set and a perturbed set. When their union is nonempty, the Jaccard similarity is the number of shared variants divided by the number of distinct variants in either set. It ranges from 0 for disjoint nonempty sets to 1 for identical sets. If both sets are empty, the formula has a zero denominator. For this tutorial, define their Jaccard similarity to be 1 by convention because the two sets are identical; nevertheless, report that both runs selected no variants. This convention should not be confused with evidence that a signal is stable.
Set comparisons need interpretation alongside inclusion frequencies and signal-level results. When two correlated variants alternate as the selected representative, their set similarity may be low even though at least one variant is selected in nearly every run. Conversely, identical selected sets do not guarantee identical effect estimates or equally strong conditional tests.
Worked Examples
Worked Example: Selection Frequency When Correlated Variants Alternate
In a hypothetical study, an analyst repeats the same joint-model selection rule in 10 sample-resampling runs. In every run, exactly one of two highly correlated variants, \(v_A\) or \(v_B\), is retained. Variant \(v_A\) is selected in runs 1–6, and \(v_B\) is selected in runs 7–10. Both variants are available in all 10 runs.
For \(v_A\), the numerator is 6 and the availability denominator is 10:
For \(v_B\), the corresponding calculation is
Neither variant has a high individual selection frequency, but at least one member of this pair is selected in all 10 runs. If the two variants have \(r^2=0.96\) in the relevant sample, the pattern is consistent with unstable choice of representative among strongly correlated candidates. It does not establish that both variants are independent signals, or that either one is causal. A useful report would give both variant-specific frequencies and the fraction of runs retaining at least one member of the pair, here \(10/10=1.00\).
Worked Example: Comparing Sets, Including Empty Sets
Suppose one run selects \(S_1=\{v_1,v_2\}\), and a second selects \(S_2=\{v_2,v_3\}\). Their intersection contains one variant, while their union contains three:
The sets therefore share one of the three distinct selected variants. Now suppose a third run selects no variants, so \(S_3=\varnothing\). Comparing \(S_3\) with \(S_1\) gives an empty intersection and a union of size 2:
Comparing two runs that both select nothing would give a denominator of zero under the ratio formula. By the stated convention, \(J(\varnothing,\varnothing)=1\), but the report should also state that neither run selected a variant. That distinction matters: identical empty outputs are agreement between procedures, not evidence for a reproducible association signal.
Worked Example: Sensitivity to Alternative LD Estimates
Consider a hypothetical locus with two candidate variants, \(v_1\) and \(v_2\). The sample, candidate set, covariates, and selection rule are held fixed; the analyst reruns the conditional model using three defensible LD estimates. The table gives the estimated conditional effect for \(v_2\), its standard error, and the two-sided \(p\)-value calculated from the reported test statistic. For illustration, use the normal approximation \(z=\hat{\beta}/SE\) and its two-sided tail probability.
| LD estimate | \(\hat{\beta}_{v_2}\) (phenotype units per allele) | \(SE\) | \(z\) | Two-sided \(p\) |
|---|---|---|---|---|
| 1 | 0.10 | 0.04 | 2.50 | 0.0124 |
| 2 | 0.07 | 0.05 | 1.40 | 0.162 |
| 3 | 0.02 | 0.06 | 0.333 | 0.739 |
For the first estimate, the statistic is \(0.10/0.04=2.50\), giving a two-sided normal-approximation \(p\)-value of about 0.0124. For the second it is \(0.07/0.05=1.40\), with \(p\approx0.162\); for the third it is \(0.02/0.06=0.333\), with \(p\approx0.739\). At a prespecified threshold of \(\alpha=0.05\), \(v_2\) meets the threshold only under the first LD estimate. Its conditional evidence is therefore sensitive to the LD input under this analysis.
Suppose the accompanying conditional estimates for \(v_1\) are 0.20, 0.19, and 0.21 phenotype units per allele, each with \(SE=0.05\). Their test statistics are respectively 4.00, 3.80, and 4.20, all beyond the two-sided 0.05 normal critical value of approximately 1.96. This pattern suggests a more stable conditional association for \(v_1\) than for \(v_2\) across the three LD estimates. It does not prove that \(v_1\) is causal, and the results remain conditional on the model and on the plausibility of the alternative LD estimates.
What the Stability Frequencies Mean
A selection frequency is a descriptive summary under a specified perturbation distribution. Its interpretation can be seen directly. Let \(I_{jb}\) equal 1 if variant \(j\) is selected in run \(b\), among runs where it is available. If the perturbation runs are generated from a defined procedure, the average of these indicators is the observed selection frequency:
Each term records one selection outcome; summing counts selections, and dividing by the number of opportunities gives the fraction of eligible runs in which selection occurred. Thus the frequency estimates how often this procedure selects the variant under those perturbations. It is not the probability that the variant is causal, nor is it automatically a confidence level or a \(p\)-value. Changing the perturbation scheme can change the frequency.
Interpret stability at more than one level. A variant-level result asks whether the same named variant is repeatedly selected. A locus-level result asks whether evidence for one or more conditional associations persists. A signal-class result can group strongly correlated alternatives when the scientific question concerns the shared signal rather than a uniquely identified variant. These summaries are complementary: a robust locus-level result can coexist with unstable variant identity.
Also inspect effect direction and magnitude, conditional standard errors, and whether the same model terms remain in the selected set. A variant selected in nearly every run may still have estimates that vary considerably. A variant selected only intermittently may reflect a hard threshold: small changes in estimates can move a test statistic across the selection cutoff. Report the selection rule and threshold so readers can see what “selected” means.
Common Mistakes and Practical Reporting
- Calling a high selection frequency a causal probability. The frequency is calculated across specified perturbations of an analysis procedure. It does not assign a probability that a variant causes the phenotype.
- Equating variant identity with signal stability. Correlated variants can substitute for one another. Report both variant-level selection and whether the broader conditional signal persists.
- Ignoring availability. If a variant is absent from some perturbed candidate sets, count selection frequency among runs where it was available and report that denominator.
- Leaving empty-set comparisons undefined without saying so. When both sets are empty, the Jaccard ratio has a zero denominator. State the adopted convention and report that no variants were selected.
- Changing several analysis choices without tracking them. If sample, variant set, and LD estimate all change at once, it may be impossible to identify which change produced an altered result. Perturb one source at a time when the aim is to diagnose sensitivity.
A clear report names the baseline model and selection rule, describes each perturbation and how many runs were performed, gives variant availability and selection frequencies, and summarizes set overlap or locus-level persistence. It should also identify results that change substantially, including changes in conditional effect estimates or significance. This makes stability useful as a diagnostic of dependence on data and modeling choices, rather than as a substitute for inferential evidence or biological validation.
Check Your Understanding
Try these questions without looking back.
- A variant is eligible in 8 perturbation runs and selected in 5. What is its selection frequency, and why should the availability count be reported?
- Two nonempty selected sets share no variants. What is their Jaccard similarity? What is the similarity when both sets are empty under the convention used here, and what should be reported alongside it?
- Two strongly correlated variants alternate as the selected variant, while at least one is retained in nearly every run. What is stable, and what remains uncertain?
- Why is a selection frequency of 0.80 not an 80% probability that the variant is causal?
- If an association changes when the LD estimate changes but the sample and candidate set stay fixed, what source of sensitivity has been examined?