Maximum Likelihood for Linear Mixed Models

Why the Marginal Likelihood Matters

A linear mixed model represents genetic similarity through covariance: individuals with related genotypes can have correlated phenotypes even after accounting for measured covariates. Maximum likelihood (ML) fits the model by asking which mean and covariance parameters make the observed phenotype vector most plausible under its marginal distribution. The covariance is therefore part of the likelihood, not merely an adjustment made after fitting the mean.

The key calculation uses the distribution of the observed data after the random genetic effects have been integrated out. This tutorial constructs that likelihood, derives its score equations, and uses them in a small hypothetical genetic example. The focus is likelihood fitting itself; the likelihood of error contrasts and its use for variance-component estimation are separate topics.

Marginal Gaussian Likelihood

For example, a genetic random effect \(u\) and residual error \(\varepsilon\) can be modeled as \(u\sim N(0,\sigma_g^2K)\) and \(\varepsilon\sim N(0,\sigma_e^2I)\), independently. Here \(K\) is a known genetic relationship matrix, \(\sigma_g^2\) is the genetic variance component, \(\sigma_e^2\) is the residual variance, and \(I\) is the identity matrix. Since \(y=X\beta+u+\varepsilon\), integrating out \(u\) and \(\varepsilon\) gives \(V=\sigma_g^2K+\sigma_e^2I\). If \(K\) is positive semidefinite and \(\sigma_e^2>0\), this \(V\) is positive definite.

Write \(r=y-X\beta\) for the residual vector. The log likelihood, including its normalizing constant, is

The quadratic term measures the residual relative to its covariance: directions with greater modeled variance contribute less per unit residual. The log-determinant term, \(\log|V|\), penalizes covariance matrices that spread probability across a larger volume. Both terms depend on covariance parameters, so maximizing only the quadratic term would not be ML.

Score Equations

For a fixed covariance matrix, maximizing the log likelihood over \(\beta\) is generalized least squares. Differentiating with respect to \(\beta\) gives the fixed-effect score:

Setting this score to zero yields the generalized least-squares estimate at a given \(\theta\), provided \(X^\mathsf{T}V^{-1}X\) is invertible:

This is a conditional fit: it estimates the fixed effects for the current covariance parameters. In ML fitting, the covariance parameters must still be chosen by maximizing the likelihood after accounting for this dependence of \(\hat\beta\) on \(\theta\).

For a scalar covariance parameter \(\theta_j\), define \(V_j=\partial V/\partial\theta_j\). When \(V\) is differentiable and positive definite, the covariance-parameter score is

The first term reflects the change in the log determinant; the second reflects the change in the residual quadratic form. At an interior ML solution, the relevant scores are zero. If parameters are constrained, a maximum can instead occur at the edge of the allowed parameter space, where a score need not be zero.

Why These Scores Have This Form

The score equations follow by differentiating the log likelihood while holding the other parameters fixed. The fixed-effect derivative follows from differentiating the quadratic form in \(y-X\beta\). For a covariance change \(dV\), the standard matrix differentials are \(d\log|V|=\operatorname{tr}(V^{-1}dV)\) and \(dV^{-1}=-V^{-1}(dV)V^{-1}\). Applying them to the log likelihood gives

Taking \(dV=V_j\,d\theta_j\) and dividing by \(d\theta_j\) gives the covariance score above. The derivation also clarifies why the determinant cannot be omitted: it contributes the entire trace term. In practice, the score equations may not have a closed-form joint solution, so numerical optimization of the likelihood is common.

Fitting a Genetic Variance Component

For a genetic mixed model, ML fitting evaluates the full marginal likelihood at candidate covariance parameters. A common strategy is to compute \(\hat\beta(\theta)\) for each candidate, substitute it into \(\ell(\beta,\theta)\), and maximize the resulting profile log likelihood over the allowed covariance parameters. “Profile” here means that the fixed effects are optimized for each covariance choice; it does not remove the determinant term or replace the marginal likelihood.

Interpretation and Common Pitfalls

ML is useful because it gives a single objective for estimating fixed effects and covariance parameters under the assumed marginal distribution. In a genetic model, the relationship matrix determines how genetic variance induces covariance among individuals; the likelihood uses that pattern along with the observed phenotypes. The estimates remain conditional on the model being appropriate, including its mean structure, covariance structure, and Gaussian assumptions. An estimated covariance parameter is not, by itself, evidence of a causal genetic effect.

• Dropping the log determinant. The reasoning “larger variance always makes residuals less costly, so choose the largest variance” ignores \(\log|V|\). The correction is to maximize the full log likelihood, including both the determinant and residual terms.
• Treating the fixed-effect estimate as independent of covariance parameters. The GLS estimate generally changes with \(V(\theta)\). Recompute \(\hat\beta(\theta)\) at each candidate covariance value when profiling the likelihood.
• Assuming every score-zero point is a maximum. A zero score can identify a minimum or another stationary point. Check the likelihood’s behavior and the parameter constraints, as in the worked fit.
• Claiming a positive-variance maximum when the observation is zero. In the one-observation example, if \(y=0\), then \(\ell(v)=-\tfrac12[\log(2\pi)+\log v]\), which increases without bound as \(v\to0^+\). There is no maximizer over \(v>0\); the likelihood has a supremum approached at the boundary. The formula \(\hat v=y^2\) gives an interior maximum only when \(y\ne0\).