Cluster Size and Information in Random-Effects Integration

Why Cluster Size Is Not the Whole Story

In a generalized linear mixed model, integrating over a cluster’s random effects combines its conditional response likelihood with the random-effects distribution. The shape of that product determines how difficult the integration is. A larger cluster often makes the conditional integrand more concentrated, which can make a local approximation effective. But cluster size alone does not guarantee this: response probabilities, predictor values, and the directions in which the data inform the random effects all matter.

This tutorial focuses on that link between information and integration. It does not treat a large cluster as automatically easy or a small cluster as automatically difficult. Instead, it identifies the curvature and shape that determine whether a local Gaussian approximation is plausible, and explains how those features affect Laplace approximation and quadrature.

Definition: For cluster $i$, let $b_i$ be its $q$-dimensional random-effects vector, with density $f(b_i)$, and let $y_i$ denote its responses. Conditional on $b_i$, suppose the responses are independent with the specified GLMM conditional distributions. At fixed model parameters, define the log conditional integrand as $$ \ell_i(b_i)=\log f(b_i)+\sum_{j=1}^{n_i}\log f(y_{ij}\mid b_i), $$ where $n_i$ is the number of observations in the cluster. The conditional integrand is $\exp\{\ell_i(b_i)\}$; its integral over $b_i$ contributes to the marginal likelihood. If $\tilde b_i$ maximizes $\ell_i$, the negative Hessian $J_i=-\nabla^2\ell_i(\tilde b_i)$ describes local curvature at that mode.

A large eigenvalue of $J_i$ means the log integrand falls quickly in the corresponding direction, so the integrand is locally narrow in that direction. A small eigenvalue means it is relatively broad. When $J_i$ is positive definite, the local Gaussian approximation has covariance matrix $J_i^{-1}$. This is a description of local integrand scale, not a claim that the exact conditional distribution of the random effects is Gaussian.

Where the Information Comes From

The log integrand has two sources of curvature. The random-effects density contributes curvature even if the cluster contains little response information. The conditional likelihood adds curvature that depends on the observed design and response model. For Gaussian random effects $b_i\sim N(0,D)$, the prior contribution to $J_i$ is $D^{-1}$. Thus the total curvature includes a baseline from the random-effects distribution as well as information from the cluster.

For example, in a Bernoulli logit random-intercept model, write $\eta_{ij}=x_{ij}^{T}\beta$, so that $P(Y_{ij}=1\mid b_i)=p_{ij}(b_i)=\operatorname{logit}^{-1}(\eta_{ij}+b_i)$. For a scalar random intercept $b_i\sim N(0,\sigma^2)$, the log integrand, up to terms constant in $b_i$, is $$ \ell_i(b_i)=\sum_{j=1}^{n_i}\{y_{ij}(\eta_{ij}+b_i)-\log(1+\exp(\eta_{ij}+b_i))\}-\frac{b_i^2}{2\sigma^2}. $$ Differentiation gives $$ -\ell_i''(b_i)=\sum_{j=1}^{n_i}p_{ij}(b_i)\{1-p_{ij}(b_i)\}+\frac{1}{\sigma^2}. $$

Each response contributes at most $1/4$ units of likelihood curvature, attained when its conditional probability is $1/2$. When a probability is near zero or one, its contribution is smaller. Thus $n_i$ responses do not necessarily provide $n_i/4$ units of information: their probabilities may differ, and the sum of their contributions can be much lower. In a random-slope model, information also depends on predictor values. Responses can inform an intercept while providing little information about a slope if the relevant predictor barely varies.

Worked Example: More Balanced Responses, Greater Intercept Curvature

Consider a hypothetical Bernoulli logit random-intercept cluster with $\beta=0$, $\sigma^2=1$, and all fixed predictors equal to zero. Compare a cluster of $n_i=8$ responses, of which 4 are successes, with a cluster of $n_i=32$ responses, of which 16 are successes.

For either cluster, at $b_i=0$ the conditional probability is $p=0.5$. The derivative of the log integrand there is $\sum_j(y_{ij}-p)-b_i/\sigma^2=0$, so $b_i=0$ is the mode. The curvature is therefore $$ J_i=\sum_{j=1}^{n_i}p(1-p)+\frac{1}{\sigma^2} =n_i(0.5)(0.5)+1. $$ For $n_i=8$, $J_i=8(0.25)+1=3$, and the local Gaussian scale is $1/\sqrt{3}=0.5774$ on the random-intercept scale. For $n_i=32$, $J_i=32(0.25)+1=9$, giving $1/\sqrt{9}=0.3333$. Under these balanced conditions, the larger cluster produces a narrower integrand around the mode. These scales are local curvature summaries, not exact posterior standard deviations.

Worked Example: The Same Cluster Size Can Give Different Information

Consider two hypothetical clusters, each with $n_i=20$, in a Bernoulli logit random-intercept model with $\sigma^2=1$. In the first, set the fixed linear predictor to $\eta=\log(0.1/0.9)$ and observe 2 successes. In the second, set $\eta=0$ and observe 10 successes.

In the first cluster, $p(0)=0.1$, and the likelihood score at $b_i=0$ is $2-20(0.1)=0$. The prior score is also zero, so $b_i=0$ is the mode. Its curvature is $J_i=20(0.1)(0.9)+1=2.8$, giving local scale $1/\sqrt{2.8}=0.5976$.

In the second cluster, $p(0)=0.5$, and the likelihood score is $10-20(0.5)=0$. Again the prior score is zero and the mode is $b_i=0$. Now $J_i=20(0.5)(0.5)+1=6$, giving local scale $1/\sqrt{6}=0.4082$. Both clusters have 20 observations, but the second has greater likelihood curvature at its mode because its responses are less saturated. This comparison holds the mode at zero in both cases; it does not imply that response probabilities alone determine integrand shape in every cluster.

From Curvature to Laplace Approximation

Laplace approximation replaces the log integrand near its mode by a quadratic function. Its justification is that a sufficiently concentrated integrand can place most of its mass close to the mode, where the quadratic expansion is useful. The approximation for a $q$-dimensional random-effects integral is $$ \int \exp\{\ell_i(b_i)\}\,db_i \ \approx\ \exp\{\ell_i(\tilde b_i)\}(2\pi)^{q/2}|J_i|^{-1/2}. $$

The formula follows directly from the local Taylor expansion. Let $u=b_i-\tilde b_i$. At an interior mode, $\nabla\ell_i(\tilde b_i)=0$, so $$ \ell_i(\tilde b_i+u)=\ell_i(\tilde b_i)-\frac{1}{2}u^TJ_i u+R_i(u), $$ where $R_i(u)$ contains the terms of order three and higher. Dropping this remainder gives an integrand proportional to $\exp(-u^TJ_i u/2)$. Since $J_i$ is positive definite, its Gaussian integral is $(2\pi)^{q/2}|J_i|^{-1/2}$, yielding the stated expression.

The role of concentration can be made more precise under regularity conditions. In fixed dimension, suppose the cluster’s log-likelihood contribution grows with $n_i$, its curvature grows with $n_i$ in every relevant direction, the mode remains interior, and the derivatives are smooth and well behaved near the mode. Then local displacements are typically of order $n_i^{-1/2}$. The quadratic term remains order one at that scale, while higher-order terms shrink. In regular settings, the leading relative error from the local approximation is typically of order $n_i^{-1}$. This is a qualitative asymptotic result, not a universal error bound: it can fail when some directions receive little information, the mode is poorly behaved, or the integrand is strongly nonquadratic over the region carrying appreciable mass.

Worked Example: Cluster Size Does Not Guarantee Slope Information

Consider a hypothetical Bernoulli logit random-intercept and random-slope model with $b_i=(b_{0i},b_{1i})^T\sim N(0,I)$, fixed effects equal to zero, and linear predictor $b_{0i}+x_{ij}b_{1i}$. For this illustration, compare two clusters of eight responses, balanced at 4 successes and 4 failures. Use the curvature at $b_i=0$, where each response has probability $0.5$.

First let four predictor values be $-1$ and four be $1$, with two successes at each value. The score for both random effects is zero at $b_i=0$, so that point is the mode. The likelihood curvature is $0.25X^TX$, where $X$ has rows $(1,x_{ij})$. Here $X^TX=\begin{pmatrix}8&0\\0&8\end{pmatrix}$, so adding the prior curvature $I$ gives $J_i=\begin{pmatrix}3&0\\0&3\end{pmatrix}$. Both local scales are $1/\sqrt{3}=0.5774$.

Now let all eight predictor values equal zero, with four successes. The intercept score is zero at $b_i=0$, and the slope score is zero because every $x_{ij}=0$; this is again the mode. Now $X^TX=\begin{pmatrix}8&0\\0&0\end{pmatrix}$, so $J_i=\begin{pmatrix}3&0\\0&1\end{pmatrix}$. The intercept scale remains $0.5774$, but the slope scale is $1/\sqrt{1}=1$: the data add no slope curvature, and that direction is informed only by the random-effects distribution. Both clusters have the same size, but only the first has within-cluster predictor variation that informs the random slope.

Implications for Quadrature and Practical Assessment

Quadrature approximates the integral using weighted evaluations of the integrand. When a cluster’s integrand is narrow, a fixed, nonadaptive set of nodes may place too few evaluations in the region carrying most of its mass. Adaptive quadrature addresses this scale problem by locating and scaling nodes using the integrand’s mode and local curvature. It can therefore use evaluations more effectively when clusters have different modes or widths. Neither approach makes a poor local shape Gaussian: a sharply concentrated but skewed integrand can still challenge Laplace approximation, and quadrature accuracy also depends on node placement and the dimension of integration.

Cluster information is useful when deciding what to watch for in an analysis. A few considerations are especially important:

Check which directions are informed. A large cluster can still have weak information for a random slope if its predictor varies little.
Account for response probabilities. In binary models, probabilities near zero or one contribute little local likelihood curvature, even when the cluster is not small.
Distinguish sharpness from Gaussian shape. Curvature describes the local width, not the integrand’s skewness or its behavior away from the mode.
Consider dimension. Information can grow in each random-effect direction while integration remains demanding when the dimension $q$ is large.

These considerations are not a substitute for evaluating the integration used in a particular model. They explain why approximation quality can vary across clusters and why the number of observations alone is an incomplete guide. Laplace approximation is attractive when the important mass is concentrated and locally close to quadratic. Quadrature can better represent departures from that local approximation when its node scheme adequately covers the important region, but it can require substantial computation, especially in multiple dimensions.

Common Mistake and Limitations

A common mistake is to reason, “This cluster has many observations, so its random-effects integral must be easy.” The conclusion does not follow. The observations may contribute little information because their conditional probabilities are extreme, a random-slope predictor may have little variation, or the integrand may remain asymmetric around its mode. The correction is to consider the information and shape in each random-effect direction, not merely $n_i$.

The local-curvature argument also has limits. It assumes a well-defined, interior mode and a positive-definite curvature matrix there. Curvature is local; it cannot by itself rule out substantial mass away from the mode or describe a long tail. In addition, a direction informed mainly by the random-effects density does not gain the same concentration from a larger cluster as a direction informed by the likelihood. These are reasons to interpret the Laplace formula as an approximation whose suitability depends on the model and data, not as an automatic consequence of having many within-cluster observations.

Key takeaway: Within-cluster sample size matters through the information it supplies about the random effects. Response probabilities and predictor design determine how that information accumulates; the resulting curvature and integrand shape, rather than cluster size alone, govern the prospects for Laplace and quadrature accuracy.

Check Your Understanding

Answer each question, then compare with the solutions below.

In a Bernoulli logit random-intercept model with $\sigma^2=1$, a balanced cluster has $n_i=12$ and $p=0.5$ at its mode. What is the curvature $J_i$ and the corresponding local scale?
With $n_i=20$, why does $p=0.1$ contribute less binary likelihood curvature than $p=0.5$, even though both are valid conditional probabilities?
A cluster has many responses but every value of a slope predictor is zero. What can its responses inform about the random intercept and random slope?
Does large curvature at the mode prove that Laplace approximation is accurate? Explain.

Solutions. (1) $J_i=12(0.5)(0.5)+1=4$; the local scale is $1/\sqrt{4}=0.5$. (2) The response contribution is $p(1-p)$: it is $0.09$ at $p=0.1$, versus $0.25$ at $p=0.5$, so the total likelihood curvature for 20 responses is $1.8$ versus $5$, before adding prior curvature. (3) The responses can inform the intercept, but with all slope predictors zero the likelihood has no slope information; the slope direction is informed by the random-effects distribution. (4) No. Curvature describes local width, not skewness or the integrand away from its mode, so a sharply concentrated integrand can still depart from a quadratic shape.