A Single Workflow for a Structured Cohort
Population structure can create an association between genotype and outcome even when the association within ancestry groups is weak or absent. It can also distort a genuine within-group association. A useful analysis therefore does not begin and end with one adjusted regression. It checks the structure in the data, compares estimates from approaches with different assumptions, and asks whether any remaining disagreement has a plausible explanation.
This tutorial follows one hypothetical simulated cohort through that process. The outcome is quantitative, and the tested variant is coded additively: \(G=0,1,2\) is the number of effect alleles. The association coefficient is measured in outcome units per additional effect allele. The cohort has two broad ancestry strata, labelled A and B for this example. These labels are simplified analytical groupings; they do not represent the full complexity of genetic ancestry or human identity.
Start With the Data and the Ancestry Plots
Suppose the simulated cohort contains 200 unrelated or partially related participants, 100 in each stratum. The genotype counts are:
| Stratum | Genotype 0 | Genotype 1 | Genotype 2 | Mean \(G\) |
|---|---|---|---|---|
| A | 64 | 32 | 4 | 0.40 |
| B | 16 | 48 | 36 | 1.20 |
A principal-component plot for these same participants shows two visible clusters along PC1, with substantial overlap in PC2. Genotype dosage is overlaid using a colour gradient. The higher average dosage in the cluster associated with stratum B is visible, as is a shift in the outcome distribution between clusters. These plots suggest that genotype and outcome both vary with ancestry structure, so a pooled estimate needs scrutiny.
A plot is diagnostic, not a test of confounding and not a substitute for modelling. A visible cluster does not prove that every relevant source of structure has been captured. Conversely, a cloud without sharp clusters does not prove that ancestry-related confounding is absent. The analysis should also inspect sample and variant quality, missingness, relatedness, and whether the chosen principal components capture the relevant structure.
Worked Example: Reading the Genotype Difference Between Strata
The mean genotype dosage in stratum A is calculated from the genotype counts as
For stratum B,
Thus the mean dosage differs by \(1.20-0.40=0.80\) alleles between the strata. This difference is a reason to examine ancestry confounding: if the outcome also differs between strata, a pooled genotype–outcome association can partly reflect the contrast between groups rather than the within-stratum association.
How Pooling Can Change the Association
For a transparent calculation, suppose the simulated outcome follows \(Y=0.5G+2S+\varepsilon\), where \(S=0\) for stratum A and \(S=1\) for stratum B. The genotype coefficient within each stratum is \(0.5\) outcome units per allele, while the stratum B baseline is 2 units higher. Assume \(E[\varepsilon\mid G,S]=0\). This conditional mean assumption implies the population covariances \(\operatorname{Cov}(G,\varepsilon)=0\) and \(\operatorname{Cov}(S,\varepsilon)=0\). For this finite simulated cohort, we additionally stipulate that the realized sample covariance \(\operatorname{Cov}(G,\varepsilon)\) is zero. Without these assumptions, the covariance calculations below would need additional error-covariance terms.
The stratum-specific genotype variances, computed using divisor 100, are \(0.32\) in A and \(0.48\) in B. Since the two strata have equal size, the overall genotype variance is the average within-stratum variance plus the variance of the stratum means:
Here \(S\) has mean \(0.5\), and the empirical covariance between genotype and stratum indicator is \(\operatorname{Cov}(G,S)=0.20\). Because the conditional mean of the error is zero, the pooled covariance is \(\operatorname{Cov}(G,Y)=0.5\operatorname{Var}(G)+2\operatorname{Cov}(G,S)\). The pooled ordinary least-squares slope is therefore:
Worked Example: Calculating the Unadjusted Pooled Slope
Using the empirical moments above and the stated error assumption,
The pooled regression slope is covariance divided by genotype variance:
The unadjusted estimate is about 1.21 outcome units per allele, despite the stipulated within-stratum slope of 0.50. The extra pooled association arises from the joint between-stratum differences in genotype and outcome. The population covariance is zero under \(E[\varepsilon\mid G,S]=0\); the exact sample calculation also relies on the stipulated zero realized sample covariance. Merely specifying a mean-zero error overall would justify neither.
Principal-Component Adjustment
A common next model includes genotype and selected principal components as covariates. In this simplified example, suppose PC1 separates the two strata and is coded so that it is an affine transformation of \(S\). A regression with an intercept, \(G\), and PC1 then adjusts for the between-stratum outcome shift. The coefficient of \(G\) is identified by genotype variation remaining after accounting for PC1, rather than by the contrast in average genotype dosage between the two clusters.
In actual cohorts, PCs usually vary continuously and may not perfectly encode stratum membership. The number and form of PCs should be chosen using the design and diagnostics, not simply by retaining every available PC. A linear PC adjustment may also fail to control a nonlinear ancestry–outcome relationship; flexible adjustment may be appropriate where diagnostics support it. These choices can change the coefficient and its uncertainty.
Worked Example: Estimating the PC-Adjusted Association
For the hypothetical data, take PC1 as \(S\) after centring and rescaling. The adjusted model is \(Y=\beta_0+\beta_GG+\beta_1\mathrm{PC1}+e\). Since PC1 spans the same two group means as \(S\), the genotype information that estimates \(\beta_G\) is the within-stratum genotype variation.
The within-stratum sum of squared genotype deviations is
The simulation is constructed so that the within-stratum residual cross-product is zero: \(\sum_i(G_i-\bar{G}_{S_i})\varepsilon_i=0\). Thus the adjusted slope is
The adjusted point estimate is 0.500 outcome units per allele, rather than 1.214 in the pooled model. This is a designed illustration of a coefficient, not a claim that every PC-adjusted analysis will recover a true effect. Its interpretation depends on adequate control of the relevant structure, a suitable outcome model, and appropriate treatment of relatedness and other sources of dependence.
Compare With Ancestry-Stratified Estimates
Stratified analysis fits the association within each stratum, then reports stratum-specific estimates or combines them using a suitable method. It makes the within-group comparison explicit and can reveal whether the estimated association differs across strata. Its limitations include reduced precision when strata are small, sensitivity to how strata are defined, and the need to account for multiple stratum-specific estimates when making an overall claim.
Worked Example: Calculating and Combining Within-Stratum Slopes
In stratum A, the genotype sum of squares is \(100(0.32)=32\). The simulated outcome has slope 0.50 and zero residual cross-product with genotype centred within A, so the cross-product is \(0.50(32)=16\). The slope is
In stratum B, the genotype sum of squares is \(100(0.48)=48\), and the corresponding cross-product is \(0.50(48)=24\). Therefore,
For a common-slope fit with a separate intercept for each stratum, the pooled within-stratum slope is the ratio of summed cross-products to summed sums of squares:
The two estimates agree here because the simulation gives both strata the same slope. In a real analysis, a meaningful difference between stratum-specific slopes may signal effect heterogeneity, sparse-data instability, residual confounding, or chance. It should not be hidden automatically by reporting only a common estimate.
Add Relatedness With a Mixed Model
Participants can remain correlated after adjustment for global ancestry, especially when the cohort contains relatives or finer-scale population structure. A linear mixed model represents this dependence through a random effect with a specified covariance structure. A useful form for this workflow is
Here \(\mathbf{C}\) contains selected ancestry PCs and any other prespecified covariates, \(\mathbf{K}\) is a genetic relatedness or structure matrix, and the random effect \(\mathbf{u}\) accounts for covariance represented by that matrix. The model assumes the specified covariance is adequate and that the random effect and residual have the stated distributions and are independent of one another. The fixed-effect genotype coefficient remains an association conditional on the included covariates; adding a random effect does not make the result causal.
Worked Example: Interpreting a Mixed-Model Fit
Suppose a fit to the same simulated cohort includes PC1 as a fixed covariate and a relatedness matrix for the random effect. For illustration, assume the fitted output is \(\hat{\beta}_G=0.490\) outcome units per allele, with \(SE=0.120\). A Wald statistic for \(H_0:\beta_G=0\) is
Using a two-sided standard normal reference distribution gives \(p\approx 4.44\times10^{-5}\). Under the fitted model and its assumptions, this is evidence against a zero conditional association. It is not the probability that the null hypothesis is true, and it does not show that the allele causes the outcome. The mixed-model estimate is close to the PC-adjusted and stratified estimates, while its standard error reflects the fitted covariance structure.
The mixed-model result is an illustrative fitted output, not a quantity determined by the genotype counts alone. Its value and uncertainty depend on the full individual-level outcomes, the relatedness matrix, the covariates, and the estimation method. Reproducible reporting should state these choices and relevant diagnostics.
Putting the Estimates Side by Side
| Approach | Illustrative estimate | What it addresses | Key qualification |
|---|---|---|---|
| Unadjusted pooled regression | 1.214 units per allele | Overall genotype–outcome association | Mixes within- and between-stratum contrasts |
| PC-adjusted regression | 0.500 units per allele | Structure represented by included PCs | Depends on PC choice and model form |
| Common-slope stratified analysis | 0.500 units per allele | Within-stratum association | Depends on stratum definition; check heterogeneity |
| PC-adjusted mixed model | 0.490 units per allele | PC adjustment plus covariance represented by relatedness matrix | Depends on adequacy of the covariance model |
The estimates differ because the methods use different comparisons and assumptions, not because one number is automatically correct. In this constructed example, the large unadjusted estimate is explained by the between-stratum contribution. The close agreement of the other three estimates is reassuring for this particular data-generating setup, but it does not establish that all structure or bias has been removed.
Why the Workflow Matters—and What It Cannot Resolve
The sequence of checks is useful because no single adjustment addresses every threat to a valid association analysis. Principal components summarize major axes of genetic variation, but a small set may miss relevant structure or model ancestry effects inadequately. Stratification makes within-group comparisons clear, but groups can be too broad or too small, and a common-effect summary can conceal heterogeneity. A mixed model can account for covariance represented by its relatedness matrix, but it cannot correct every form of confounding, phenotype error, selection bias, or variant-level technical artifact.
The analysis also depends on the target population and sampling process. Case-control sampling, differential participation, or ancestry-related missingness can alter observed relationships. Local ancestry may matter at a locus in an admixed cohort even when global PCs are included. A marker’s association can differ across groups because of allele-frequency differences, linkage disequilibrium, or genuinely different effects; these explanations are not interchangeable. The result should be described as an association in the analyzed sample, with the model and limitations stated.
Use ancestry plots and quality checks to identify clustering, genotype–ancestry alignment, missingness patterns, and potential relatedness.
Calculate the unadjusted estimate, but treat it as a reference rather than a final result when genotype and outcome vary with ancestry.
Fit PC-adjusted and stratified analyses, checking whether their estimates and stratum-specific patterns are compatible.
Use a mixed model when relatedness or structure-induced covariance warrants it; state its covariance assumptions and describe unresolved sources of bias.
Common pitfall: “The unadjusted association is large, so the variant must have a strong effect; adding PCs only makes the result less significant.” The reasoning confuses a pooled contrast with a within-ancestry association and treats attenuation as a loss of truth. The correction is to examine how genotype and outcome vary across ancestry, calculate or compare within-group estimates, and assess whether an adjusted model is appropriate. A smaller adjusted estimate may be less confounded, but adjustment can also be inadequate or excessive; the interpretation depends on the design, model, and biological question.
Check Your Understanding
Use the simulated cohort and the reasoning in this workflow to answer the following.
- Why can a pooled slope differ from the association within each ancestry stratum even when the within-stratum slopes are equal?
- In the covariance calculation, what assumption justifies setting \(\operatorname{Cov}(G,\varepsilon)=0\)? What changes if that covariance is nonzero?
- Why does adding a principal component that separates the two strata change which genotype variation identifies the slope?
- What information does a stratified analysis provide that a single common-slope estimate does not?
- What kind of dependence is the random effect in a mixed model intended to represent, and what does it not guarantee?
- If adjusted estimates are similar, why is it still necessary to discuss residual structure and other limitations?