Tutorials › Statistical Genetics › Principal Components and Mixed Models: Complementary Adjustments

Population Stratification Control · Tutorial 191 of 989

Principal Components and Mixed Models: Complementary Adjustments

Learn what each adjustment models, how their roles differ, and when using both in an association model is defensible.

Intermediate 10 min read

What You'll Learn

  • Distinguish fixed-effect principal-component adjustment from a genomic relationship matrix random effect
  • Write a joint association model and state its key assumptions
  • Interpret how a relationship matrix specifies residual covariance
  • Understand why including both adjustments can address different patterns of structure
  • Recognize when joint adjustment may be redundant, unstable, or overly aggressive

Why Consider Two Adjustments?

Population structure can affect a genetic association analysis in more than one way. Ancestry-related differences in phenotype means can create a systematic trend across individuals. In addition, individuals who are genetically similar may have correlated residual outcomes, even after measured covariates and broad ancestry trends have been included. Principal components and genomic relationship matrix mixed models offer different ways to represent these features.

The key distinction is not that principal components always capture “broad” structure while a relationship matrix always captures “fine” structure. A genomic relationship matrix can also reflect broad ancestry, and the two adjustments may overlap. Rather, principal components entered as fixed effects model specified directions of variation in the conditional mean. A relationship-matrix random effect models covariance among residual outcomes. Including both can be useful when both features are plausible and supported, but it is not an automatic improvement.

Definition: A principal-component adjustment includes selected genotype-derived principal-component scores as fixed covariates in the association model. A genomic relationship matrix mixed model includes a random effect whose covariance is proportional to a genomic relationship matrix \(K\). A joint adjustment includes both the fixed principal-component covariates and the relationship-matrix random effect.

The Joint Model and Its Assumptions

For a continuous phenotype \(y_i\) and a genetic dosage \(g_i\) at the variant being tested, a basic joint linear mixed model is

$$ y_i = \beta_0 + \beta_G g_i + \boldsymbol{c}_i^\mathsf{T}\boldsymbol{\gamma} + \boldsymbol{q}_i^\mathsf{T}\boldsymbol{\alpha} + u_i + \varepsilon_i. $$

Here, \(\boldsymbol{c}_i\) contains other fixed covariates, such as age or sex; \(\boldsymbol{q}_i\) contains the selected principal-component scores; and \(\beta_G\) is the association coefficient of interest, measured as the expected phenotype difference per additional copy of the tested allele, conditional on the model. The PC coefficients \(\boldsymbol{\alpha}\) are fixed effects: each included score is assigned its own estimated mean trend. The term \(u_i\) is the random effect, and \(\varepsilon_i\) is the individual-specific residual.

A standard specification assumes

$$ \boldsymbol{u} \sim N(\boldsymbol{0},\tau K), \qquad \boldsymbol{\varepsilon} \sim N(\boldsymbol{0},\sigma^2 I), \qquad \boldsymbol{u} \perp \boldsymbol{\varepsilon}. $$

The matrix \(K\) is an \(n \times n\) positive semidefinite genomic relationship matrix aligned to the analyzed individuals. Its entry \(K_{ij}\) summarizes genomic similarity between individuals \(i\) and \(j\); the exact scale depends on how the matrix was constructed. The variance component \(\tau\) controls the size of the relationship-linked variation, while \(\sigma^2\) is the residual variance not represented by \(K\). Under these assumptions, the outcome covariance is

$$ \operatorname{Var}(\boldsymbol{y}\mid \boldsymbol{g},\boldsymbol{C},\boldsymbol{Q}) = \tau K + \sigma^2 I. $$

The model also assumes that the conditional mean is adequately represented by the included covariates and that the chosen covariance structure is a reasonable description of remaining dependence. The association coefficient is an adjusted association, not necessarily a causal effect. For a binary phenotype, an appropriate generalized mixed model is needed; the equations and worked regression below concern a continuous phenotype.

Worked Example: A Principal Component Changes the Association Estimate

Consider six hypothetical individuals. Let \(g\) be dosage of a tested allele, \(q\) be one centered PC score, and \(y\) be a continuous phenotype in arbitrary phenotype units. Suppose the values are

Individual\(g\)\(q\)\(y\)
10-2-0.4
20-1-2.2
3002.0
4100.8
5115.0
6123.2

First fit an intercept-and-genotype regression without the PC. The mean phenotype at dosage 0 is \((-0.4-2.2+2.0)/3=-0.2\), and at dosage 1 it is \((0.8+5.0+3.2)/3=3.0\). Thus the unadjusted dosage coefficient is \(3.0-(-0.2)=3.2\) phenotype units per allele copy.

Now fit \(y_i=\beta_0+\beta_Gg_i+\alpha q_i+\varepsilon_i\). The data were constructed so that \(y_i=1+0.8g_i+1.2q_i+e_i\), where \(e=(1,-2,1,-1,2,-1)\). The residuals sum to \(1-2+1-1+2-1=0\), their sum within each dosage group is zero, and their PC cross-product is \((-2)+2+0+0+2-2=0\). Consequently, they are orthogonal to the intercept, dosage, and PC columns, and the fitted coefficients are exactly \(\hat{\beta}_0=1\), \(\hat{\beta}_G=0.8\), and \(\hat{\alpha}=1.2\). The adjusted association is therefore \(0.8\) phenotype units per allele copy in this hypothetical sample. The change from 3.2 to 0.8 illustrates how a PC-associated phenotype trend can contribute to an unadjusted association when dosage and ancestry vary together.

What the Relationship Matrix Adds

A fixed PC coefficient and a relationship-matrix random effect play different mathematical roles. The fixed coefficient changes the expected phenotype along a selected PC direction. The random effect instead says that residual deviations for genetically similar individuals tend to covary. This distinction can be seen directly from the covariance model.

Worked Example: Translating a Relationship Matrix Into Covariance

Suppose a hypothetical two-person subset has

$$ K= \begin{pmatrix} 1 & 0.4\\ 0.4 & 1 \end{pmatrix}, \qquad \tau=0.5,\qquad \sigma^2=0.75. $$

Under the mixed-model assumptions, the outcome covariance is \(\tau K+\sigma^2I\). Substitution gives

$$ 0.5 \begin{pmatrix} 1 & 0.4\\ 0.4 & 1 \end{pmatrix} +0.75 \begin{pmatrix} 1 & 0\\ 0 & 1 \end{pmatrix} = \begin{pmatrix} 1.25 & 0.20\\ 0.20 & 1.25 \end{pmatrix}. $$

Each person has modeled residual variance \(1.25\) phenotype-units-squared, and the two residual outcomes have covariance \(0.20\) phenotype-units-squared. This covariance is not a claim that either person has a mean shift of \(0.20\). It represents dependence around the model’s conditional means. A fixed PC effect would instead alter those conditional means according to the individuals’ PC scores.

The two approaches can overlap because PCs and \(K\) are often derived from related genotype data. To make the distinction more explicit, suppose \(K\) has eigendecomposition \(K=U\Lambda U^\mathsf{T}\). In a direction given by an eigenvector with eigenvalue \(\lambda_j\), the random effect has variance \(\tau\lambda_j\). If that direction is also included as a PC fixed effect, the model estimates a mean coefficient in that direction as well. One parameter describes a mean trend; the other describes variation and covariance around the mean. The presence of both is not, by itself, a mathematical contradiction, but overlapping directions can make the model less informative or less stable.

Why Fixed Mean Adjustment and Random Covariance Can Complement Each Other

The basic reason for combining the adjustments is that mean structure and residual dependence are distinct properties. A PC fixed effect can account for a systematic outcome trend along the included ancestry directions, without assuming that every pair of individuals with similar scores has correlated errors. A relationship-matrix random effect can represent correlated outcomes according to genomic similarity, without requiring the mean to follow a particular slope along every PC.

This complementarity is most defensible when there is a substantive reason to expect both patterns. For example, a measured ancestry axis may correspond to a mean phenotype difference that should be adjusted explicitly, while additional genomic similarity among individuals may induce residual correlation. Using both then allows the fixed effects to represent the specified mean trend and \(K\) to represent residual covariance. The random effect may also capture structure not fully represented by a small set of fixed PCs.

A short derivation shows why the two components are not interchangeable. In matrix notation, write the fixed-effect part as \(X\boldsymbol{\beta}\), where \(X\) includes the intercept, tested genotype, other covariates, and selected PCs. Under the model, \(E(\boldsymbol{y}\mid X)=X\boldsymbol{\beta}\), while \(\operatorname{Var}(\boldsymbol{y}\mid X)=V=\tau K+\sigma^2 I\). Thus

$$ E(\boldsymbol{y}\mid X)=X\boldsymbol{\beta} \quad\text{and}\quad \operatorname{Var}(\boldsymbol{y}\mid X)=\tau K+\sigma^2 I. $$

Changing a PC coefficient changes the first expression, the conditional mean. Changing \(\tau\) changes the second expression, the covariance. In a generalized least squares fit, the coefficient estimate has the form

$$ \hat{\boldsymbol{\beta}}=(X^\mathsf{T}V^{-1}X)^{-1}X^\mathsf{T}V^{-1}\boldsymbol{y}. $$

The covariance model therefore affects how observations are weighted when estimating the fixed effects, while the PC columns specify which mean trends are estimated. This is why including both is not simply “adjusting twice” in every case: they act on different parts of the model. But if a fixed PC trend and the covariance structure encode nearly the same information, estimates can become imprecise, and the model may gain little over a simpler alternative.

Worked Example: Interpreting a Joint-Model Association

Suppose a hypothetical cohort analysis tests a dosage effect on a continuous phenotype measured in points. The joint model includes two PC scores as fixed effects and a genomic relationship matrix random effect. Assume the fitted output is \(\hat{\beta}_G=0.18\) points per allele copy with \(SE=0.10\) points per allele copy. For illustration, use a large-sample normal test of \(H_0:\beta_G=0\).

The test statistic is

$$ z=\frac{\hat{\beta}_G}{SE}=\frac{0.18}{0.10}=1.80. $$

The two-sided normal p-value is \(2P(Z\geq1.80)\approx 0.0719\), where \(Z\sim N(0,1)\). The corresponding approximate 95% confidence interval is

$$ 0.18\pm1.96(0.10) =0.18\pm0.196 =(-0.016,\ 0.376) $$

points per allele copy. At a 0.05 significance level, this illustrative test does not reject \(H_0\). The interval is compatible with a small negative association, no association, or a positive association up to about 0.376 points per copy under the model. It is not a 95% probability statement about the fixed unknown coefficient. The estimate is conditional on both the included PC trends and the covariance assumptions; it does not establish that the adjustment is correct or that the association is causal.

Choosing Whether to Use Both

The decision should follow the structure of the analysis, rather than a rule that more adjustments are always safer. First specify which fixed covariates and PC directions are justified for the conditional mean. Then ask whether residual correlation linked to genomic similarity remains plausible and whether the relationship matrix is appropriate for the sample. Check that the genotype data used to construct \(K\) and the outcome data refer to the same individuals in the same order, and that the matrix has been prepared for the intended model.

A useful comparison is between a PC-adjusted fixed-effect model, a relationship-matrix mixed model, and a joint model. Compare the estimated dosage effect and its uncertainty, not only whether a p-value crosses a threshold. Large changes can indicate that the approaches are handling overlapping structure, that important dependence is being modeled differently, or that the association is sensitive to modeling choices. They do not by themselves identify which model is correct. The variance-component estimate and model diagnostics can help assess whether the random effect contributes, but a near-zero estimate does not prove there is no population structure.

Key takeaway: PCs entered as fixed effects specify selected ancestry-related mean trends; a genomic relationship matrix random effect specifies residual covariance linked to genomic similarity. A joint model is useful when both are needed, but overlap, weak information, or inappropriate assumptions can make the combination redundant or unstable.

Common Mistakes

  • Wrong reasoning: “The PCs and relationship matrix come from genotypes, so including both necessarily adjusts for the same thing twice.” Correction: They can overlap, but fixed PC effects model conditional means, whereas the random effect models covariance. Assess whether both roles are justified.
  • Wrong reasoning: “The relationship matrix random effect replaces all ancestry adjustment.” Correction: It can account for covariance associated with genomic similarity, but it does not automatically represent every mean trend of interest. Fixed PCs may still be appropriate.
  • Wrong reasoning: “Adding both must improve calibration.” Correction: The joint model depends on correctly specified mean and covariance structures. Redundancy, weak data support, or model misspecification can undermine precision or inference.
  • Wrong reasoning: “A significant adjusted coefficient proves the variant causes the phenotype.” Correction: The coefficient is a conditional association under the specified model and assumptions; adjustment alone does not establish causation.

Check Your Understanding

Answer these questions using the distinction between conditional means and residual covariance.

  1. In a joint model, what does a PC fixed-effect coefficient describe, and what does the relationship-matrix random effect describe?
  2. If \(K_{ij}\) is large for two individuals, what feature of the mixed model does that affect?
  3. Why can a relationship matrix capture some broad ancestry structure without making fixed PC adjustment conceptually identical?
  4. Name one circumstance in which including both PCs and a relationship-matrix random effect may be reasonable, and one concern that should be checked.
  5. In the worked association example, why does the confidence interval not mean there is a 95% probability that \(\beta_G\) lies inside it?
Clinical Biostats