What Makes a Population-Averaged Result Reproducible?
A population-averaged GLMM result is not fully described by a coefficient and a standard error. A reader needs to know which population was averaged over, how the random effects were integrated out, which response-scale quantity was estimated, and how uncertainty was calculated. Without these details, two analysts can report different numbers from the same fitted model while both use the phrase “marginal effect.”
The key is to report the estimand before reporting its estimate. An estimand specifies the population quantity of interest; an estimate is its value calculated from fitted model parameters and the chosen target covariate distribution. The introductory material on marginal effects establishes why integrated predictions can differ from conditional predictions. Here the focus is how to make a report explicit enough to interpret and reproduce.
A report should identify both distributions in this definition. The covariate distribution \(Q\) might be the empirical distribution of individuals in a target sample, an explicitly weighted distribution of covariate strata, or a distribution representing a named target population. The random-effects distribution might be a normal distribution with a stated variance or covariance matrix. It also matters whether the target concerns new clusters drawn from that distribution, or some other random-effects target.
A standard GLMM specification also states the conditional response model and link. For example, a binary random-intercept model may have \(Y_{ij}\mid b_j,x_{ij}\sim\operatorname{Bernoulli}(p_{ij})\), with \(\operatorname{logit}(p_{ij})=\beta_0+\beta_1a_{ij}+\beta_2x_{ij}+b_j\) and \(b_j\sim N(0,\sigma_b^2)\). Conditional independence of responses given covariates and random effects, the chosen response distribution, and the random-effects distribution are model assumptions. The marginal quantity is obtained by integrating the conditional mean under those assumptions, not by simply reporting a fixed-effect coefficient.
A Reporting Recipe
A concise report can still be reproducible if it answers the following questions in a consistent order. State the response and its units; identify the target population and the covariate distribution used to represent it; give the fitted conditional model and random-effects distribution; specify the quantity integrated over; define the response-scale estimand and contrast; and state the uncertainty method and interval.
Define who or what the average represents and how covariate values are weighted. If the target is an empirical distribution, state whose records form it and whether the weighting is individual- or cluster-based.
Report the response distribution, link, relevant fixed effects, random-effects structure, and fitted random-effects distribution. Clarify whether the integration represents new clusters from that distribution.
Give the marginal mean or rate and the contrast (difference, ratio, or another specified comparison), including units and the values being compared.
Name the interval procedure, what sources of uncertainty it includes, and any numerical integration or model checks relevant to the result.
In general, if a contrast is \(C(\theta)\), uncertainty concerns the sampling variation of its estimate \(\widehat{C}=C(\hat{\theta})\), and may also concern the estimated target covariate distribution. A confidence interval is a procedure that, under its assumptions and repeated sampling, covers the fixed estimand at its stated rate. It is not a probability statement about the parameter after seeing the particular interval.
Worked Examples
Worked Example: A Standardized Binary-Outcome Difference
Hypothetical model and target. A binary response is modeled with a logistic random-intercept GLMM, \(\operatorname{logit}\{P(Y_{ij}=1\mid A_{ij},b_j)\}=\beta_0+\beta_1A_{ij}+b_j\), where \(b_j\sim N(0,0.64)\). Suppose the fitted intercept is \(-1\) and the fitted treatment coefficient is \(0.6\). The target is a new individual in a new cluster, with \(A\) set to 0 or 1, and the random intercept integrated over the stated normal distribution. The response-scale estimand is the marginal probability of response; the contrast is the treatment-minus-control risk difference.
Integration calculation. For illustration, approximate the normal integral with three-point Gauss–Hermite quadrature. For a normal random intercept with standard deviation \(0.8\), the nodes are approximately \(-1.3856,0,1.3856\), with weights \(1/6,2/3,1/6\). At each node, calculate the inverse-logit of the linear predictor. The resulting approximate control probabilities are \(0.0843,0.2689,0.5952\). Thus
For treatment, add \(0.6\) to each node-specific linear predictor. The corresponding probabilities are approximately \(0.1436,0.4013,0.7282\), so
The estimated risk difference is \(0.4128-0.2925=0.1203\), or about 12.03 percentage points. This is an approximate marginal difference for the stated new-cluster target, not the conditional coefficient \(0.6\). In an analysis, the quadrature approximation should be checked by increasing the number of integration nodes until the reported estimates are stable. A report should also give the uncertainty procedure and interval; the calculation above alone does not quantify uncertainty.
Worked Example: A Standardized Count-Rate Difference
Hypothetical model and target. Suppose a Poisson GLMM has a log link and a random intercept \(b_j\sim N(0,0.25)\). The exposure offset is included so that the conditional mean count is the exposure multiplied by the rate. Assume the fitted conditional rates at exposure \(1\) are 2 events per person-year for control individuals with covariate value \(x=0\), and 3 for those with \(x=1\). Treatment multiplies each conditional rate by \(1.4\). The target covariate distribution assigns weight \(0.25\) to \(x=0\) and \(0.75\) to \(x=1\), with individual weights. The estimand is the new-cluster marginal rate standardized to this distribution.
For a normal random intercept, \(E\{\exp(b_j)\}=\exp\{0.5\operatorname{Var}(b_j)\}\). This follows from the moment-generating function of a normal variable: if \(b_j\sim N(0,\sigma_b^2)\), then \(E\{\exp(b_j)\}=\exp(0+\sigma_b^2/2)\). Here the integration multiplier is \(\exp(0.125)\approx1.1331\).
The covariate-standardized conditional control rate before random-effects integration is \(0.25(2)+0.75(3)=2.75\) events per person-year. After integration it is \(2.75(1.1331)\approx3.116\) events per person-year. For treatment, the two conditional rates are \(2(1.4)=2.8\) and \(3(1.4)=4.2\). Therefore the standardized marginal treatment rate is
The marginal rate difference, calculated from the unrounded rates, is approximately \(4.3626-3.1162=1.2465\) events per person-year. The marginal rate ratio is \(4.3626/3.1162\approx1.400\). Both contrasts are meaningful, but they answer different questions and should not be interchanged. The report must retain the exposure definition, target covariate weights, and random-effects distribution so that the rate has a clear interpretation.
Worked Example: Reporting an Interval for a Marginal Contrast
Hypothetical analysis. Consider a binary-outcome GLMM whose target is all individuals in a specified eligible-population sample, weighted equally at the individual level. For each target record, the fitted probability is integrated over a normal new-cluster random intercept; those integrated probabilities are averaged over the target records. Suppose this procedure gives marginal probabilities of \(0.31\) under control and \(0.42\) under treatment. The estimated risk difference is
or 11 percentage points. Suppose a cluster bootstrap, with 2,000 resamples, gives a percentile 95% confidence interval of \(0.03\) to \(0.19\). In this procedure, clusters are resampled, the model is refitted in each resample, and the integration and standardization are repeated. This allows the interval to reflect model-estimation uncertainty and cluster sampling variation, subject to the bootstrap design and model assumptions. If the target covariate distribution is estimated from the resampled clusters, its sampling variation is also represented; if it is held fixed by design, the report should say so.
A suitable interpretation is: “In the specified eligible-population target, the estimated treatment-minus-control difference in response probability was 11 percentage points (cluster-bootstrap 95% confidence interval, 3 to 19 percentage points).” The interval describes uncertainty in the procedure for estimating the target contrast; it does not mean that there is a 95% probability that this fixed contrast lies between 3 and 19 percentage points. Nor does the model-based contrast alone establish a causal effect: a causal interpretation requires a design and assumptions that justify it.
Uncertainty, Diagnostics, and Alternatives
There are two common approaches to uncertainty for a marginal contrast. A delta-method interval uses an estimated covariance matrix for the model parameters and a derivative of the contrast with respect to those parameters. It can be efficient when the estimator is sufficiently regular and the approximation is adequate. A simulation or bootstrap approach repeatedly evaluates the full marginal estimand under draws or resamples, which can better accommodate nonlinear contrasts and complex standardization. The chosen method should match the fitted model and sampling design; merely attaching the standard error of a fixed-effect coefficient to a marginal contrast is not justified.
For a reproducible report, state whether uncertainty includes fixed-effect estimates, variance and covariance parameters for the random effects, and estimation of the target covariate distribution. For a cluster-sampled study, resampling individuals independently can understate uncertainty because it breaks the clustering structure. A cluster-level resampling scheme is often more appropriate when clusters are the independent sampling units, but the precise scheme must reflect the actual design. Numerical integration should also be checked: repeat it with more nodes or tighter tolerances and confirm that estimates and interval endpoints do not materially change.
Compare a marginal result with a conditional coefficient only after naming the scale and target of each. A conditional coefficient describes a model parameter conditional on the random effects and covariates, on the link scale. A marginal contrast averages response-scale predictions over the stated distributions. These are different summaries, not competing estimates of a single automatically defined number. A prediction conditional on a particular estimated cluster effect is also distinct from a prediction integrated over the random-effects distribution for a new cluster.
Several reporting failures are diagnostic warning signs: no named target population; no indication of individual or cluster weighting; an unspecified random-effects distribution; an interval method with no resampling or approximation details; a contrast reported without units; or an estimate that changes noticeably as integration accuracy increases. Check that treatment or setting values are clearly stated, the response scale is named, and the same target distribution is used for each side of the contrast unless a different comparison is explicitly intended. Also distinguish a confidence interval from a prediction interval: the latter concerns variation in future outcomes or clusters, rather than uncertainty only about a population mean or contrast.
Common Reporting Mistake
A frequent mistake is to report “the treatment effect was 0.6” because the fitted GLMM treatment coefficient is \(0.6\), while the stated scientific question concerns a population-averaged difference in response probability. The coefficient alone does not give that risk difference: it is on the model’s link scale and is conditional on the model structure. The correction is to define the response-scale marginal means for the chosen target, calculate the requested contrast, and report its uncertainty. If the conditional coefficient is also scientifically relevant, report it separately and label its scale and conditional interpretation.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions below.
- A report gives a marginal response probability but does not say which covariate distribution was used. What important information is missing?
- For a target population, the estimated marginal probabilities are \(0.28\) under control and \(0.37\) under treatment. Calculate and interpret the treatment-minus-control risk difference.
- A Poisson random-intercept model has \(b_j\sim N(0,0.36)\). What random-effects integration multiplier applies to the conditional rate?
- Why is a confidence interval for a fixed treatment coefficient not automatically an interval for a marginal risk difference?
Solutions. (1) The report must identify the target population and how its covariates were weighted or distributed; otherwise the averaging target is not reproducible. (2) The difference is \(0.37-0.28=0.09\), or 9 percentage points. It is the estimated marginal treatment-minus-control response-probability difference for the stated target, not necessarily a causal effect. (3) The multiplier is \(E\{\exp(b_j)\}=\exp(0.36/2)=\exp(0.18)\approx1.197\). Thus the marginal rate integrates the conditional rate over the specified normal random-intercept distribution. (4) The coefficient is on the link scale and is a conditional model parameter, whereas the marginal risk difference is a nonlinear response-scale contrast integrated over random effects and covariates. Its uncertainty must be calculated for that contrast, accounting for the relevant estimated quantities.