PQL In Sparse Binary Data

Why Sparsity Matters for Binary PQL

Binary data can be sparse in more than one sense. A study may have few observations within each cluster, few events overall, or clusters in which every observed outcome is the same. These features can occur together, but they are not interchangeable. For PQL, they matter because the method replaces the binary response model with a locally linear, Gaussian working model. When the data provide little within-cluster variation, or fitted event probabilities are close to zero or one, that approximation has less support from the observed outcomes.

The key question is not simply whether a data set is small. It is whether the observed outcome pattern gives a plausible basis for the local Gaussian approximation and for separating cluster-level variation from the binary sampling variation. Sparse patterns are a reason to scrutinize the approximation, not a guarantee that a particular estimate is biased.

Model and Working Approximation

Consider a binary random-intercept logistic model. For cluster \(i\) and observation \(j\), let \(Y_{ij}\) be zero or one, let \(x_{ij}\) be a covariate vector, and let \(b_i\) be a cluster-specific random intercept. The model is

The standard assumptions here are conditional independence of observations given the random effects and covariates, a Bernoulli conditional distribution, a correctly specified logit mean model, and independent normally distributed random intercepts with common variance \(\tau^2\). PQL constructs an approximate Gaussian working model by linearizing the link around current fitted values. Thus, the pseudo-response used in the working fit is not the original zero-or-one outcome.

For a current linear predictor \(\eta_{ij}\), write \(p_{ij}=\operatorname{logit}^{-1}(\eta_{ij})\). The derivative of the logit link \(g(p)=\log\{p/(1-p)\}\) is \(g'(p)=1/[p(1-p)]\). A first-order Taylor expansion around \(p_{ij}\) gives the working response \(z_{ij}\) and its conditional working variance:

The corresponding working weight is \(w_{ij}=p_{ij}(1-p_{ij})\). To see the variance expression, the linearized response is a constant plus \((Y_{ij}-p_{ij})/[p_{ij}(1-p_{ij})]\). Since a Bernoulli variable has variance \(p_{ij}(1-p_{ij})\), scaling that random part gives variance \[ \frac{p_{ij}(1-p_{ij})}{\{p_{ij}(1-p_{ij})\}^2} =\frac{1}{p_{ij}(1-p_{ij})}. \] This is a local approximation: it does not turn the original binary response into a genuinely Gaussian variable.

What Sparse Patterns Do to the Approximation

A cluster with one or two observations contains little information about its own response pattern. If the outcomes are all zero, for example, the data do not reveal whether that pattern reflects a low cluster probability, chance variation from a moderate probability, or some combination. With more observations, within-cluster variation can be more informative, although it does not automatically resolve every modeling issue.

Near-uniform outcomes also matter through the fitted probabilities. If a fitted probability is near zero or one, \(p(1-p)\) is small. The working weight is then small, while the working variance and the adjustment from \(Y-p\) to \(z\) can be large. Consequently, a single outcome that is surprising under the current fitted probability can produce a very large pseudo-response, even though the original observation is only zero or one.

This relationship explains why sparse data and extreme fitted probabilities can reinforce each other. Few observations or few events may lead to uncertain cluster-level patterns; a fitted model may then assign very small probabilities to some outcomes. The resulting working data can be dominated by a few large adjustments and weakly weighted observations. That is precisely where a Gaussian approximation based on local linearization deserves caution.

A Practical Assessment of Sparsity

Before interpreting a binary PQL fit, describe the data at the cluster level as well as overall. Useful summaries include the number of clusters, the distribution of cluster sizes, the total number and proportion of events, and the number of clusters with no events, all events, or mixed outcomes. When cluster sizes vary, report these counts by size range as well: a uniform outcome in a cluster of size two carries a different amount of information than one in a much larger cluster.

A simple audit should distinguish data sparsity from fitted-probability extremity. The first is visible in the outcome table; the second is visible only after fitting. For the fitted model, examine the range of conditional fitted probabilities and the working weights \(p_{ij}(1-p_{ij})\). A concentration of very small weights, large working responses, or both indicates that the working fit may be relying heavily on local approximations in regions where the Bernoulli distribution is strongly non-Gaussian.

These checks diagnose vulnerability, not failure by themselves. A data set can contain many uniform clusters because the event rate is low, while PQL estimates may still be relatively stable in a particular design. Conversely, a moderate overall event proportion does not rule out sparsity if events are concentrated in only a few clusters. The cluster-level distribution is essential context.

Comparison, Limitations, and a Common Mistake

PQL is computationally convenient because its working step resembles a weighted Gaussian mixed-model fit. Its limitation here is that the Gaussian response is constructed by local linearization, while the original observations remain binary. Likelihood-based integration is a main alternative for fitting GLMMs: it targets the likelihood after accounting for random effects rather than relying on this same working-response approximation. It can provide a useful comparison, but it is not automatically immune to sparse-data problems, model misspecification, or limited information.

The assessment in this tutorial is also limited. Uniform clusters and extreme fitted probabilities identify conditions under which the PQL approximation is least comfortable; they do not quantify the bias or establish that a particular variance estimate is wrong. The same observed pattern can occur by chance, and the impact depends on the full design, covariates, number of clusters, and underlying probabilities. Sparse-data summaries should therefore inform interpretation and further assessment, not replace it.

Common mistake: “Most clusters are all zero, so the random-intercept variance must be large.” This reasoning confuses an outcome pattern with its cause. If events are rare, many all-zero clusters can occur even with a common event probability; for instance, the worked benchmark gave a 72.9% chance of an all-zero cluster of size three at \(p=0.10\). Uniform outcomes provide limited within-cluster contrast, but they do not by themselves identify a large random-effect variance. Check event frequency, cluster sizes, and event distribution before drawing conclusions.