Gaussian Likelihoods for Quantitative Traits by Genotype

Why Model a Quantitative Trait by Genotype?

A quantitative phenotype, such as a blood concentration or measured height, may differ across genotype groups. A Gaussian model describes the phenotype distribution within each group by a mean and a residual variance. The likelihood then expresses how well proposed values of those parameters account for the observed measurements, conditional on the genotypes recorded for the individuals.

There are two common variance specifications. A shared-variance model allows the mean to differ by genotype but assumes the residual variance is the same in every group. A separate-variance model allows both the mean and residual variance to differ by genotype. The choice is a modeling assumption, not a conclusion guaranteed by observing different group means.

Definition: Let \(Y_i\) be the quantitative phenotype for individual \(i\), and let \(G_i\) denote that individual's observed genotype group. For each genotype group \(g\), write \(n_g\) for the number of observed individuals and \(y_{ig}\) for the phenotype of individual \(i\) in that group. A Gaussian likelihood conditional on the observed genotype groups models each \(Y_i\mid G_i=g\) as normal with mean \(\mu_g\) and a positive residual variance. The likelihood is the product of the corresponding phenotype densities over the observed individuals.

Writing the Conditional Likelihood

Assume that, conditional on the observed genotypes, phenotype measurements from different individuals are independent. In the shared-variance model, assume that each phenotype in group \(g\) has distribution \(N(\mu_g,\sigma^2)\), with a genotype-specific mean \(\mu_g\) and a common variance \(\sigma^2>0\). The normal density for one observation is

f(y_{ig}\mid G_i=g)=\frac{1}{\sqrt{2\pi\sigma^2}}\exp\left\{-\frac{(y_{ig}-\mu_g)^2}{2\sigma^2}\right\}.

Independence makes the likelihood the product of these densities. Equivalently, its log-likelihood is the sum of their logarithms:

L(\{\mu_g\},\sigma^2)=\prod_g\prod_{i=1}^{n_g}\frac{1}{\sqrt{2\pi\sigma^2}}\exp\left\{-\frac{(y_{ig}-\mu_g)^2}{2\sigma^2}\right\},

\ell(\{\mu_g\},\sigma^2)=-\frac{n}{2}\log(2\pi\sigma^2)-\frac{1}{2\sigma^2}\sum_g\sum_{i=1}^{n_g}(y_{ig}-\mu_g)^2,\qquad n=\sum_g n_g.

The observed genotype values identify which mean applies to each measurement. In this conditional likelihood, genotype probabilities are not included: the likelihood concerns the phenotype densities given the observed groups. A model with a different variance for each genotype group replaces \(\sigma^2\) in group \(g\)'s densities with \(\sigma_g^2>0\).

Estimating Means and Variances

For a nonempty group, define its sample mean and within-group residual sum of squares by

\bar{y}_g=\frac{1}{n_g}\sum_{i=1}^{n_g}y_{ig},\qquad S_g=\sum_{i=1}^{n_g}(y_{ig}-\bar{y}_g)^2.

For any proposed mean \(\mu_g\), the identity

\sum_{i=1}^{n_g}(y_{ig}-\mu_g)^2=S_g+n_g(\bar{y}_g-\mu_g)^2

shows that the residual sum of squares is minimized at \(\mu_g=\bar{y}_g\). Since the Gaussian log-likelihood decreases as the residual sum of squares increases when variance is fixed, this is also the likelihood-maximizing mean. In the shared-variance model, all group means are estimated this way, and the total residual sum of squares is \(S=\sum_g S_g\).

Theorem: If every modeled group is nonempty and the total within-group residual sum of squares \(S=\sum_g S_g\) is positive, the shared-variance likelihood has maximum-likelihood estimates \(\hat{\mu}_g=\bar{y}_g\) and \(\hat{\sigma}^2=S/n\). In the separate-variance model, if every group is nonempty and \(S_g>0\) in every group, the estimates are \(\hat{\mu}_g=\bar{y}_g\) and \(\hat{\sigma}_g^2=S_g/n_g\) for each group.

The variance denominator in these maximum-likelihood estimates is the number of observations in the relevant likelihood, not a degrees-of-freedom adjustment. The estimate \(S_g/(n_g-1)\), when defined, is commonly used as an unbiased sample variance under the model, but it is not the Gaussian maximum-likelihood estimate.

Why the Variance Estimate Uses the Likelihood Denominator

After substituting the group means, the shared-variance log-likelihood is

\ell(\sigma^2)=-\frac{n}{2}\log(2\pi\sigma^2)-\frac{S}{2\sigma^2}.

Differentiating with respect to the positive parameter \(\sigma^2\) gives

\frac{d\ell}{d\sigma^2}=-\frac{n}{2\sigma^2}+\frac{S}{2(\sigma^2)^2}.

Setting this derivative to zero yields \(n\sigma^2=S\), so the stationary value is \(S/n\). When \(S>0\), the derivative is positive for \(0<\sigma^2<S/n\) and negative for \(\sigma^2>S/n\), proving that this value maximizes the log-likelihood. In the separate-variance model, the log-likelihood is a sum of group-specific terms; applying the same argument to group \(g\) replaces \(n\) and \(S\) with \(n_g\) and \(S_g\).

The condition that the relevant residual sum of squares is positive matters. If it is zero, then every observation in the relevant likelihood lies exactly at its fitted group mean. As the positive variance approaches zero, the density factors for those exact fits increase without bound. There is then no attained maximum at a positive variance.

Worked Example: Estimating a Shared Residual Variance

Suppose a hypothetical dataset contains phenotypes measured in three genotype groups. The values are 8, 10, and 12 units for group 0; 11, 13, and 15 units for group 1; and 14, 16, and 18 units for group 2. Assume independent Gaussian phenotypes conditional on genotype, a separate mean in each group, and one shared positive variance.

The group means are \(\bar{y}_0=10\), \(\bar{y}_1=13\), and \(\bar{y}_2=16\), each in phenotype units. In each group, the squared deviations from its mean sum to \(4+0+4=8\) units squared. Thus \(S=8+8+8=24\) units squared, and \(n=9\).

Because \(S>0\), the shared-variance MLE exists and is

\hat{\sigma}^2=\frac{S}{n}=\frac{24}{9}=2.667\ \text{units}^2.

The fitted residual standard deviation is \(\sqrt{2.667}=1.633\) units. The likelihood at these fitted values is

L=(2\pi\cdot 2.667)^{-9/2}\exp\left\{-\frac{24}{2\cdot 2.667}\right\} \approx (2\pi\cdot 2.667)^{-9/2}\exp(-4.5).

This example estimates a common within-genotype variance while permitting different group means. It does not assert that the three observed sample means are known population values; they are the fitted means under this model.

Allowing a Different Variance in Each Group

In the separate-variance model, conditional independence gives a product of group likelihoods, each with its own \(\mu_g\) and \(\sigma_g^2\). Maximizing the group-specific contribution gives \(\hat{\mu}_g=\bar{y}_g\) and, when \(S_g>0\), \(\hat{\sigma}_g^2=S_g/n_g\). This specification is more flexible, but it estimates a separate variance from each group's observations, so the amount of data within each group is especially relevant.

Worked Example: Fitting Separate Variances

Consider hypothetical measurements in two genotype groups. Group A has values 2, 4, and 6 units. Group B has values 5, 7, 9, and 11 units. Assume independent normal phenotypes conditional on genotype, with a separate mean and positive variance for each group.

For group A, \(\bar{y}_A=(2+4+6)/3=4\) units and \(S_A=(2-4)^2+(4-4)^2+(6-4)^2=8\) units squared. Since \(S_A>0\), \(\hat{\sigma}_A^2=8/3=2.667\) units squared.

For group B, \(\bar{y}_B=(5+7+9+11)/4=8\) units and \(S_B=9+1+1+9=20\) units squared. Since \(S_B>0\), \(\hat{\sigma}_B^2=20/4=5.000\) units squared. The fitted standard deviations are approximately 1.633 units for A and 2.236 units for B.

The fitted likelihood factors by group:

L=(2\pi\cdot 2.667)^{-3/2}\exp\left(-\frac{8}{2\cdot2.667}\right) (2\pi\cdot5)^{-4/2}\exp\left(-\frac{20}{2\cdot5}\right).

The exponent terms are \(-1.5\) and \(-2\), respectively. These estimates describe different fitted residual spreads; the data and model alone do not establish why the spreads might differ.

Boundary Cases: When a Positive-Variance MLE Does Not Exist

The variance formulas are not unconditional recipes. An empty genotype group has no phenotype observations from which to estimate its mean or its group-specific variance. More subtly, a nonempty group can have zero residual sum of squares, for example when it contains one observation or when all its observations are identical. In a separate-variance model, this makes the likelihood for that group's variance unbounded as the variance approaches zero. In the shared-variance model, the corresponding issue arises if the total within-group residual sum of squares is zero.

Worked Example: A Zero Group Residual Sum of Squares

Suppose hypothetical group A measurements are 5 and 5 units, while group B measurements are 7 and 9 units. First fit the separate-variance model. Group A has mean 5 units and \(S_A=0\) units squared. Its likelihood contribution, as a function of its positive variance \(v=\sigma_A^2\), is

L_A(v)=(2\pi v)^{-2/2}\exp(0)=\frac{1}{2\pi v}.

As \(v\) decreases toward zero, \(1/(2\pi v)\) increases without bound. Therefore there is no positive-variance MLE for group A; the formal expression \(S_A/n_A=0/2=0\) is a boundary limit, not an allowed or attained estimate. Group B has mean 8 units and \(S_B=(7-8)^2+(9-8)^2=2\) units squared, so its separate-variance MLE is \(2/2=1\) unit squared. That finite estimate does not prevent the full separate-variance likelihood from being unbounded because of group A.

Under the shared-variance model, however, the total residual sum of squares is \(S=0+2=2\) units squared and \(n=4\). The shared-variance estimate is \(2/4=0.500\) units squared, which is positive. The distinction arises because group B's residual variation contributes to the one common variance, whereas the separate model lets group A's variance shrink on its own.

Common Mistakes and Interpretation

Using a phenotype likelihood as if it included genotype probabilities. The likelihood here is conditional on the observed genotype groups. It contains phenotype densities, not probabilities for observing those genotypes.
Using \(n_g-1\) as the MLE denominator. The Gaussian likelihood estimate is \(S_g/n_g\) for a separate variance, or \(\sum_g S_g/n\) for a shared variance. The unbiased sample-variance denominator answers a different estimation criterion.
Assuming every variance formula gives an attained estimate. A zero residual sum of squares leads to an unbounded likelihood as the positive variance approaches zero. Do not report zero as a positive-variance MLE.
Calling the shared-variance model a common-mean model. Sharing the variance does not require sharing the means; each genotype group can have its own \(\mu_g\).
Interpreting fitted means causally. A difference between genotype-specific fitted means is a modeled association in these data. It is not, by itself, proof that genotype caused the phenotype difference.

Key takeaway: A Gaussian likelihood for quantitative phenotypes by genotype multiplies normal densities conditional on the observed groups. Its mean estimates are group sample means; its variance estimates use within-group residual sums of squares divided by the corresponding observation counts, provided the relevant residual sum of squares is positive.

Check Your Understanding

Try these questions by identifying the model assumption before applying a formula.

In a shared-variance model, what does the likelihood condition on, and what quantities are allowed to differ by genotype?
A genotype group has four observations and within-group residual sum of squares 12 units squared. What is its separate-variance MLE, if the group has a separate positive variance?
Why does the Gaussian likelihood MLE for a shared variance divide the total residual sum of squares by \(n\), rather than by \(n\) minus the number of genotype groups?
In a separate-variance model, one nonempty group has identical measurements. What happens to its variance likelihood as the variance approaches zero?
If one genotype group is empty, which group-specific parameters cannot be estimated from the conditional phenotype data?