How Random-Effect Covariance Shapes Marginal Predictions

Why Covariance Matters at a Chosen Covariate Value

In a random-slope GLMM, clusters can differ both in their baseline linear predictor and in how that predictor changes with a covariate. These two sources of variation need not be independent. Their covariance determines whether clusters with higher-than-average baselines tend to have steeper or flatter slopes, and that relationship changes the amount of between-cluster variation at different covariate values.

The focus here is the distribution of the cluster-specific linear predictor at a chosen value of the covariate. This is a latent-scale quantity: its variation describes how the linear predictor differs across clusters before the inverse link is applied. It is not the variance of the observed outcome, nor does it by itself give a response-scale prediction.

Definition: For cluster \(j\) at covariate value \(x\), write the linear predictor as \(\eta_j(x)=\beta_0+\beta_1x+b_{0j}+b_{1j}x\). Here, \(\beta_0+\beta_1x\) is the fixed part, \(b_{0j}\) is the cluster-specific random intercept, and \(b_{1j}\) is the cluster-specific random slope. Assume the random-effects vector \((b_{0j},b_{1j})^\mathsf{T}\) has a bivariate normal distribution with mean zero and covariance matrix \(G\), independently across clusters.

In a GLMM, the model also specifies the conditional distribution of the response given the random effects, with the conditional mean linked to \(\eta_j(x)\). For the results below, the key assumptions are the stated linear predictor, mean-zero random effects, and their covariance matrix. The normality assumption gives a normal distribution for \(\eta_j(x)\); its mean and variance follow directly from the random-effects moments.

Deriving the Distribution at \(x\)

Let the random-intercept variance be \(\tau_{00}\), the random-slope variance be \(\tau_{11}\), and their covariance be \(\tau_{01}\). The covariance matrix is

G= \begin{pmatrix} \tau_{00} & \tau_{01}\\ \tau_{01} & \tau_{11} \end{pmatrix}.

The fixed component is constant across clusters at a given \(x\). Since both random effects have mean zero, the expected linear predictor across clusters equals this fixed component. For the variance, apply the variance formula for a sum: the variance of \(b_{0j}+xb_{1j}\) includes twice their covariance multiplied by \(x\).

\begin{aligned} E\{\eta_j(x)\}&=\beta_0+\beta_1x,\\ \operatorname{Var}\{\eta_j(x)\} &=\operatorname{Var}(b_{0j}+xb_{1j})\\ &=\tau_{00}+x^2\tau_{11}+2x\tau_{01}. \end{aligned}

Because a linear combination of jointly normal random variables is normal, the full distribution is

\eta_j(x)\sim N\!\left(\beta_0+\beta_1x,\; \tau_{00}+2x\tau_{01}+x^2\tau_{11}\right).

The covariance contribution is \(2x\tau_{01}\), not merely \(x\tau_{01}\): the cross-product appears twice when the variance of \(b_{0j}+xb_{1j}\) is expanded. A positive covariance increases the variance when \(x>0\) and decreases it when \(x<0\); a negative covariance has the reverse pattern. This describes a change in spread at a particular covariate value, not a change in the fixed mean.

Assumption check: The covariance matrix must be positive semidefinite: \(\tau_{00}\geq 0\), \(\tau_{11}\geq 0\), and \(\tau_{01}^2\leq\tau_{00}\tau_{11}\). These conditions ensure the random-effects distribution is valid and the predictor variance is nonnegative for every \(x\).

Worked Examples

Worked Example: Predictor Distribution at a Specific Value

Suppose a hypothetical model uses \(x\) measured in years from a meaningful reference point. Its fixed coefficients are \(\beta_0=0.4\) and \(\beta_1=0.3\) per year. The fitted random-effects covariance estimates are \(\tau_{00}=0.36\), \(\tau_{11}=0.04\) per year squared, and \(\tau_{01}=-0.06\) per year. Find the distribution of the cluster-specific predictor at \(x=2\) years.

First check the covariance estimate is admissible. Its determinant is \(0.36(0.04)-(-0.06)^2=0.0144-0.0036=0.0108>0\), and both variances are positive, so the matrix is positive definite.

The mean is \(\beta_0+\beta_1x=0.4+0.3(2)=1.0\). The variance is

$$ 0.36+2(2)(-0.06)+(2^2)(0.04) =0.36-0.24+0.16=0.28. $$

Thus \(\eta_j(2)\sim N(1.0,0.28)\), where the second argument is the variance. The across-cluster standard deviation is \(\sqrt{0.28}=0.5292\) linear-predictor units. The negative covariance reduces the variance at this positive \(x\), relative to the \(0.36+0.16=0.52\) that would result if the random intercept and slope were uncorrelated.

Worked Example: Where Is Predictor Variation Smallest?

Consider a second hypothetical model with \(\tau_{00}=0.25\), \(\tau_{11}=0.09\), and \(\tau_{01}=0.06\). The covariate \(x\) is measured in months from its chosen reference point. Find the value of \(x\) that minimizes the across-cluster variance, and interpret it over the observed range \(-3\leq x\leq 3\).

The covariance matrix is valid because \(0.25(0.09)-(0.06)^2=0.0225-0.0036=0.0189>0\), and both variances are positive. The variance function is

$$ V(x)=0.25+0.12x+0.09x^2. $$

For a quadratic \(a+bx+cx^2\) with \(c>0\), the minimum occurs at \(-b/(2c)\). Here, \(x=-0.12/(2(0.09))=-0.6667\) months. Substitution gives

V(-0.6667)\approx0.25+0.12(-0.6667)+0.09(-0.6667)^2 =0.2100.

This value lies within the stated covariate range, so the model’s between-cluster predictor variance is smallest near \(x=-0.667\) months in that range. The positive covariance makes the variance larger for positive \(x\) and smaller for negative \(x\), relative to the variance function without the cross-term. This is a statement about the fitted latent-scale spread, not about which covariate value has the largest mean response.

Worked Example: Centering Changes Components but Not the Predictor

Use the covariance values from the first example and define a centered covariate \(u=x-1\), so \(x=u+1\). Find the random-effect covariance entries under this parameterization and verify the predictor variance at \(x=2\), which corresponds to \(u=1\).

The random contribution is \(b_{0j}+b_{1j}x=(b_{0j}+b_{1j})+b_{1j}u\). Thus the new random intercept is \(b'_{0j}=b_{0j}+b_{1j}\), while the new random slope remains \(b'_{1j}=b_{1j}\). Their variances and covariance are

\begin{aligned} \operatorname{Var}(b'_{0j})&=0.36+0.04+2(-0.06)=0.28,\\ \operatorname{Var}(b'_{1j})&=0.04,\\ \operatorname{Cov}(b'_{0j},b'_{1j})&=-0.06+0.04=-0.02. \end{aligned}

At \(u=1\), the variance is \(0.28+2(1)(-0.02)+1^2(0.04)=0.28\), matching the variance at \(x=2\) in the first example. Centering changes the reported intercept variance and covariance because it changes the reference point for the intercept. It does not change the random contribution, or its variance, at the same physical covariate value.

Variance Shape, Reference Points, and Model Choices

The variance function is quadratic in \(x\). If \(\tau_{11}>0\), it can be written by completing the square:

\operatorname{Var}\{\eta_j(x)\} =\tau_{11}\left(x+\frac{\tau_{01}}{\tau_{11}}\right)^2 +\left(\tau_{00}-\frac{\tau_{01}^2}{\tau_{11}}\right).

The minimum occurs at \(x^*=-\tau_{01}/\tau_{11}\), with minimum variance \(\tau_{00}-\tau_{01}^2/\tau_{11}\). The covariance-matrix condition makes this minimum nonnegative. If \(\tau_{11}=0\), positive semidefiniteness requires \(\tau_{01}=0\), and the variance is constant at \(\tau_{00}\); the formula for \(x^*\) does not apply. Even when \(\tau_{11}>0\), the mathematical minimum may fall outside the data’s covariate range, so it should not be described as an observed-range feature unless it lies within that range.

Centering \(x\) can make the intercept easier to interpret, for example by placing \(x=0\) at a substantively meaningful value. It also changes the intercept variance and intercept-slope covariance, as the third example showed. Therefore, the sign or size of the intercept-slope covariance must be interpreted with its covariate origin in view. The variance at any fixed original covariate value remains invariant when the model is re-expressed correctly.

The main alternatives are different assumptions about the random-effects structure. A random-intercept-only model sets the slope variation to zero and has constant predictor variance across \(x\). An independent random-intercept-and-slope model sets \(\tau_{01}=0\), leaving variance \(\tau_{00}+x^2\tau_{11}\), symmetric around \(x=0\). A correlated random-slope model allows the full quadratic, including asymmetry around the reference point. These structures answer different modeling questions; a simpler structure should not be treated as equivalent when the data support meaningful intercept-slope association.

Diagnostics, Limitations, and a Common Mistake

Inspect the estimated covariance matrix and the implied variance function over the covariate range relevant to the analysis. Check that the matrix is positive semidefinite, that the variance is not implausibly small or large in important regions, and that the estimated covariance is not being interpreted without regard to centering. Convergence warnings, estimates at a boundary, or substantial sensitivity to small changes in model specification can indicate that the available data provide weak information about one or more covariance components. Sparse observations at the extremes of \(x\) are a particular concern for estimating slope variation.

A fitted covariance matrix is an estimate, not known truth. The variance curve calculated from it does not include uncertainty in the estimated fixed effects or covariance parameters. It also relies on the chosen random-effects distribution and linear random-slope form; unusual cluster distributions or nonlinear cluster-specific trajectories may not be well represented. Finally, normality of \(\eta_j(x)\) follows from normal random effects, but applying an inverse link does not generally produce a normal response-scale distribution. Response-scale marginal predictions require additional integration and are a separate calculation.

Common mistake: “The random intercept and slope have negative covariance, so the predictor variance must be smaller at every \(x\).”
Correction: The covariance contribution is \(2x\tau_{01}\). With negative covariance, it reduces variance for positive \(x\) but increases it for negative \(x\). The total variance also includes both variance components and must be evaluated at the chosen covariate value.

Another frequent error is to call \(\tau_{00}\) the predictor variance everywhere. It is the predictor variance at \(x=0\); away from that reference point, slope variation and the covariance contribution also matter.

Check Your Understanding

Work through the questions, then compare your answers with the solutions below.

Suppose \(\tau_{00}=0.49\), \(\tau_{11}=0.16\), and \(\tau_{01}=-0.14\). Find the predictor variance at \(x=1\). Is the covariance matrix admissible?
For \(\tau_{11}=0.25\) and \(\tau_{01}=0.10\), at what value of \(x\) is the predictor variance minimized? What additional condition on \(\tau_{00}\) ensures a valid covariance matrix?
If \(\tau_{01}=0\), compare the variance at \(x=2\) with the variance at \(x=-2\). What does this tell you about the variance profile around zero?
Explain why the covariance between random intercept and slope can change after centering \(x\), even though the predictor variance at a particular original covariate value does not change.

Solutions. (1) The variance is \(0.49+2(1)(-0.14)+1^2(0.16)=0.37\). The determinant is \(0.49(0.16)-(-0.14)^2=0.0784-0.0196=0.0588>0\), and both variances are positive, so the matrix is admissible. (2) The minimum is at \(x^*=-\tau_{01}/\tau_{11}=-0.10/0.25=-0.4\). Validity requires \(\tau_{00}\geq0\) and \(\tau_{01}^2\leq\tau_{00}\tau_{11}\), so here \(\tau_{00}\geq0.01/0.25=0.04\). (3) With zero covariance, the variance is \(\tau_{00}+4\tau_{11}\) at either \(x=2\) or \(x=-2\). The quadratic variance profile is symmetric around zero. (4) Centering replaces the random intercept by a linear combination of the old intercept and slope. That changes its variance and its covariance with the slope, but the random linear combination at a fixed original \(x\) remains the same, so its variance there is unchanged.