Why Vary the Random-Effects Variance?
A GLMM’s fixed effects describe conditional relationships: they relate covariates to the response while holding the random effects fixed. A population-averaged quantity instead integrates the response over a distribution of random effects. Consequently, the fitted random-effects variance can affect marginal means and contrasts even when the fixed-effect coefficients do not change. The effect is especially important with nonlinear links, because averaging a nonlinear response function is generally not the same as applying it to an average linear predictor.
A variance-sensitivity analysis asks a focused question: if the random-effects variance were set to other plausible values, how much would the integrated response mean or marginal contrast change? This is a model-based assessment, not a new estimate of the variance and not, by itself, a confidence interval. It helps reveal whether a population-averaged conclusion depends strongly on a variance component that may be uncertain or hard to estimate.
Definition and Setup
For a random-intercept model, take \(z=1\) and \(b_j\sim N(0,\tau^2)\). In a variance-sensitivity calculation, one evaluates \(m(x;\tau^2,\beta)\) across selected values of \(\tau^2\). A common first analysis holds the fixed effects and all other model components fixed, changing only \(\tau^2\). This isolates the mathematical influence of the variance on the integrated quantity. A separate analysis can refit or update the other parameters under alternative assumptions, but that answers a broader question.
This exercise relies on the specified conditional response model, link, and random-effects distribution being appropriate for the calculation. In particular, a normal random intercept is an assumption, not a consequence of integrating. The values in the sensitivity grid should be scientifically or statistically plausible—for example, values around the fitted variance or values supported by a profile-likelihood analysis. Varying the variance over arbitrary values can produce a numerical comparison without a useful interpretation.
Worked Example: Binary Probabilities Across a Variance Grid
Worked Example: Binary Probabilities Across a Variance Grid
Suppose a hypothetical binary-outcome GLMM has a logit link, a random intercept \(b_j\sim N(0,\tau^2)\), and conditional linear predictor \(-1+\beta_T T+b_j\), where \(T=0\) denotes the reference condition and \(T=1\) the comparison condition. Let the fitted treatment coefficient be \(\beta_T=1\). Hold these fixed effects constant and compare \(\tau=0\), \(0.5\), and \(1\).
When \(\tau=0\), there is no random-intercept variation in the calculation. The two response probabilities are $$ m(0;0)=\operatorname{logit}^{-1}(-1)=0.2689,\qquad m(1;0)=\operatorname{logit}^{-1}(0)=0.5000. $$ The marginal risk difference is \(0.5000-0.2689=0.2311\), or 23.11 percentage points.
For positive \(\tau\), approximate integration using the three-point normal quadrature rule $$ E\{f(\tau Z)\}\approx \frac{2}{3}f(0)+\frac{1}{6}f(\sqrt{3}\tau)+\frac{1}{6}f(-\sqrt{3}\tau),\qquad Z\sim N(0,1), $$ where \(f(u)=\operatorname{logit}^{-1}(\eta+u)\). At \(\tau=0.5\), the nonzero nodes are approximately \(0.866\) and \(-0.866\). For \(T=0\), $$ m(0;0.25)\approx \frac{2}{3}(0.2689)+\frac{1}{6}\operatorname{logit}^{-1}(-0.134)+\frac{1}{6}\operatorname{logit}^{-1}(-1.866) \approx 0.2794. $$ For \(T=1\), symmetry gives \(m(1;0.25)=0.5000\): the two nonzero-node probabilities sum to 1, and the center probability is 0.5. Thus the risk difference is approximately \(0.2206\).
At \(\tau=1\), the nodes are approximately \(1.732\) and \(-1.732\). The reference probability is $$ m(0;1)\approx \frac{2}{3}(0.2689)+\frac{1}{6}\operatorname{logit}^{-1}(0.732)+\frac{1}{6}\operatorname{logit}^{-1}(-2.732) \approx 0.3020. $$ The comparison probability remains \(0.5000\), so the risk difference is approximately \(0.1980\). Under these illustrative assumptions, increasing the random-intercept standard deviation from 0 to 1 raises the integrated reference probability and reduces the risk difference from 23.11 to 19.80 percentage points. These are approximate quadrature results; higher-order quadrature or adaptive integration should be used when precision matters.
The example does not establish that greater variance must reduce a binary marginal contrast. It shows the result for a particular intercept and treatment coefficient. The direction and size of change depend on where the linear predictors lie on the response curve and on how the contrast is defined.
Why the Response Curve Matters
A short derivation makes the role of curvature precise. Let \(f\) be a twice-differentiable inverse-link response function and let \(Z\sim N(0,1)\). Define \(M(v)=E\{f(\eta+\sqrt{v}Z)\}\), where \(v=\tau^2\). Differentiating the normal convolution with respect to its variance gives $$ \frac{dM(v)}{dv}=\frac{1}{2}E\{f''(\eta+\sqrt{v}Z)\}. $$ At \(v=0\), this reduces to \(M'(0)=f''(\eta)/2\). Thus, for a small increase in random-effects variance, the integrated mean initially moves according to the curvature of the response function near the conditional predictor. Where the response curve is locally convex, added variation tends to raise the mean; where it is concave, added variation tends to lower it. If curvature changes over the range of random effects, the net change is an average and can vary with \(v\).
For the logistic inverse link, \(f''(\eta)=p(1-p)(1-2p)\), where \(p=\operatorname{logit}^{-1}(\eta)\). At the reference predictor \(-1\), this second derivative is positive, so adding a small amount of variance initially increases the integrated probability. At predictor 0 it is zero, and symmetry keeps the integrated probability at 0.5 for a symmetric random intercept. The treatment contrast therefore changes because the two covariate conditions are affected differently.
Worked Example: Count Rates Under a Log Link
Worked Example: Count Rates Under a Log Link
Consider a hypothetical Poisson random-intercept GLMM with a log link, \(b_j\sim N(0,\tau^2)\), and conditional mean \(\exp(\eta+b_j)\). Suppose the conditional means at two covariate settings, before adding the random intercept, are 2 and 3 events per unit of exposure. Equivalently, their linear predictors are \(\log(2)\) and \(\log(3)\).
For a normal random intercept, the moment-generating function gives \(E\{\exp(b_j)\}=\exp(\tau^2/2)\). Therefore $$ m(x;\tau^2)=\exp(\eta)\exp(\tau^2/2). $$ At \(\tau=0\), the marginal rates are 2 and 3 events per exposure unit, and their rate difference is 1 event per exposure unit. At \(\tau=0.5\), the variance is \(0.25\), so the multiplier is \(\exp(0.25/2)=\exp(0.125)\approx1.133\). The marginal rates are \(2(1.133)=2.266\) and \(3(1.133)=3.399\) events per exposure unit; the difference is 1.133 events per exposure unit.
At \(\tau=1\), the multiplier is \(\exp(0.5)\approx1.649\). The marginal rates are approximately 3.297 and 4.946 events per exposure unit, and the difference is 1.649 events per exposure unit. In contrast, the marginal rate ratio is $$ \frac{m(x_1;\tau^2)}{m(x_0;\tau^2)} =\frac{\exp(\eta_1)\exp(\tau^2/2)}{\exp(\eta_0)\exp(\tau^2/2)} =\exp(\eta_1-\eta_0), $$ so it does not change with \(\tau^2\) in this random-intercept model when the same random-intercept distribution applies to both settings. The absolute rate difference changes, even though this rate ratio does not.
This invariance has specific conditions: the log link, normal random intercept, and a common random-intercept variance and distribution across the compared settings. It does not imply that all marginal contrasts from count GLMMs are insensitive to variance. For example, random slopes can make the variance of the linear predictor depend on the covariate setting.
How to Conduct and Interpret the Analysis
Name the response mean, probability, rate, or contrast, and state the covariate settings and population over which random effects are integrated.
Use the fitted variance and a justified range of plausible alternatives. Keep clear whether the grid contains variances \(\tau^2\) or standard deviations \(\tau\).
For each grid value, integrate the response mean and then form the contrast. Use numerical integration appropriate to the random-effects distribution and model.
Report the grid, the integrated means, and the contrast with units. Describe whether the substantive conclusion is stable over the chosen range, without presenting the range as a confidence interval.
For a multivariate random effect \(b_j\sim N(0,\Sigma)\), a sensitivity analysis may vary one variance component at a time or vary a coherent set of entries in \(\Sigma\). A covariance matrix must remain positive semidefinite: not every independent combination of variances and covariances is valid. If more than one component is changed, state the scenario clearly; otherwise it is difficult to tell which change explains the result.
Worked Example: A Marginal Odds Ratio Also Changes
Worked Example: A Marginal Odds Ratio Also Changes
Use the binary example’s integrated probabilities. At \(\tau=0.5\), the reference probability is approximately 0.2794 and the comparison probability is 0.5000. The corresponding marginal odds ratio is $$ \frac{0.5000/(1-0.5000)}{0.2794/(1-0.2794)} =\frac{1}{0.2794/0.7206} \approx 2.579. $$ At \(\tau=1\), using probabilities 0.3020 and 0.5000, it is $$ \frac{0.5000/(1-0.5000)}{0.3020/(1-0.3020)} =\frac{1}{0.3020/0.6980} \approx 2.311. $$ The conditional odds ratio is \(\exp(1)\approx2.718\) at both variance settings, because the conditional treatment coefficient was held fixed. The marginal odds ratios differ because integrating over random intercepts changes the marginal probabilities. These calculations illustrate why a conditional coefficient should not be reported as though it were a marginal effect.
Diagnostics, Comparisons, and Limitations
Check the numerical integration as well as the model. Repeat calculations with more quadrature points, or a suitable simulation-based integration, and confirm that the reported means and contrasts are stable to numerical accuracy. For simulation, use enough draws and, when comparing grid values, common random draws can reduce simulation noise in differences. Also verify that the random-effects distribution and variance range match the intended target. A variance estimate near a boundary, a wide profile likelihood, or sensitivity to model specification can make the chosen range particularly consequential.
A fixed-effects-only calculation at \(\tau^2=0\) is a useful reference point, but it is not generally an adequate substitute for integrating over group heterogeneity. Conversely, comparing conditional coefficients alone does not answer how population-averaged probabilities or rates change. A variance grid is also different from a confidence interval: it displays how a target changes under selected parameter values, while an interval quantifies uncertainty under a specified inferential procedure. To include parameter uncertainty, one must also account for uncertainty in the fixed effects and variance components, for example through an appropriate simulation or resampling procedure.
A one-at-a-time grid has an important limitation: it holds other parameters fixed, even though estimates of fixed effects and variance components may be statistically dependent. It isolates one mechanism but does not describe every plausible joint parameter change. Results also depend on the random-effects distribution, link function, covariate settings, and target population. State these conditions and avoid claiming that a sensitivity pattern is universal.
Common Mistake
Check Your Understanding
Answer each question, then compare your reasoning with the solutions.
- In a random-intercept GLMM, what does an integrated response mean average over?
- For a Poisson log-link model with a normal random intercept, what multiplicative factor relates the integrated mean to the conditional mean \(\exp(\eta)\)?
- Why can a binary marginal risk difference change when the variance changes even if the conditional treatment coefficient is held fixed?
- Does a grid of marginal contrasts over plausible variances constitute a confidence interval? Explain.
Solutions. (1) It averages the conditional response mean over the specified random-effects distribution. (2) The factor is \(\exp(\tau^2/2)\), because \(E\{\exp(b_j)\}=\exp(\tau^2/2)\) for \(b_j\sim N(0,\tau^2)\). (3) The inverse link is nonlinear, and the integrated probabilities at the two covariate settings can respond differently to variance because the response curve has different curvature at their linear predictors. (4) No. It shows how the estimates vary over selected parameter values; a confidence interval requires an inferential procedure that accounts for sampling uncertainty.