Zhao-Staudenmayer-Coull-Wand (2006) General Design Bayesian Generalized Linear Mixed Models

Summary

Demonstrates that a single Bayesian generalized linear mixed model (GLMM) framework with general design matrices for random effects subsumes an extraordinarily wide range of regression models — longitudinal mixed models, penalized spline smoothing, generalized kriging, generalized additive models (GAMs), varying-coefficient models, and crossed/multilevel random effects. Markov chain Monte Carlo (MCMC) via WinBUGS handles the non-conjugate full conditionals via log-concavity (adaptive rejection sampling or slice sampling). Three public health applications illustrate the approach. The paper provides a constructive counterargument to the use of penalized quasi-likelihood (PQL) (which is asymptotically biased) and packages the general-design philosophy for applied statisticians.

Key Claims

• Penalized splines as random effects: spline coefficients $u_k$ are treated as $N(0, \sigma^2)$ random effects with diagonal covariance; the smoothing parameter $\lambda$ is the variance component $\sigma^2$. This maps nonparametric regression directly into the GLMM framework, enabling a unified Gibbs / Metropolis-Hastings (MH) sampler.
• Full model decomposition: the linear predictor $\eta = X_R\beta_R + Z_R u_R + X_G\beta_G + Z_G u_G + Z_C u_C$ separates (R) grouped/longitudinal structure, (G) general design (splines, GAMs), and (C) spatial structure (conditional autoregressive, CAR). Each component uses its own covariance: block-diagonal $\Sigma$ for R, $\sigma^2_\ell I$ for each smooth term in G, intrinsic Gaussian autoregression for C.
• Log-concavity of full conditionals: $[v_k \mid v_{-k}, G, y] \propto \exp\{(c_k^T y)v_k - 1^T b(c_k v_k + C_{-k}v_{-k}) - v_k^2/(2(V)_{kk})\}$ is log-concave for canonical exponential families $\to$ adaptive rejection sampling (Gilks-Wild 1992) and slice sampling (Neal 2003) apply directly.
• Hierarchical centering for grouped random effects: reparameterise $\gamma_i = \beta_R + u_i^R$ so the fixed and random intercept/slope are absorbed into a single block. This eliminates the high posterior correlation between $\beta_R$ and $u_R$ that slows MCMC mixing; not applicable to spline or spatial components.
• WinBUGS recommendation: beats PROC NLMIXED (SAS), glmmML (R), and other off-the-shelf tools on the general design GLMM problem; the user need not write custom MCMC code. Limitation: computation time grows rapidly with model size and is platform-restricted (Windows).
• Variance component prior robustness: IG(0.001,0.001), IG(0.01,0.01), folded-Cauchy(s=12,25) and Uniform(0,100) priors on $\sigma$ all tested; posterior means change $\leq$5% and credible intervals widen $\leq$6.5%. Favours IG(0.01,0.01) as default, with Gelman (2005) folded-Cauchy as a robust alternative.
• PQL bias: Breslow-Clayton (1993) pseudo-likelihood (PQL) is asymptotically biased (Lin-Breslow 1996), especially severely for paired binary data; the Bayesian MCMC approach does not share this deficiency.
• Three applications all converge well (Gelman-Rubin $\hat R < 1.04$): (a) Indonesian children binary logistic (vitamin A deficiency, nonlinear age effect, 9 covariates, n=275); (b) 483-family Poisson stress model (child-age spline, race, income effects; inverse income-stress dose-response); (c) Massachusetts Military Reservation cancer standardized incidence ratios (SIRs) Poisson spatial model (CAR spatial random effects, distance-to-contamination spline; elevated incidence < ~10 km).

Concepts Introduced or Extended

• Generalized Linear Mixed Model
• Bayesian Hierarchical Model
• Conditional Autoregressive Model

Entities Mentioned

• M.P. Wand (UNSW; penalized splines, semiparametric regression)
• B.A. Coull (Harvard Biostatistics)

Quotes

"The key is to allow the random effects design matrix to have general structure; hence our label general design."

"Assuming that computation time is not an issue and that Windows is available, we can report that fitting of general design GLMMs via WinBUGS has a large chance of success."

My Take

This is an important unification paper — the insight that penalized splines are just random effects with a diagonal covariance unlocks the full GLMM machinery (Bayesian, frequentist, or hybrid) for nonparametric regression, making GAMs and kriging "free" once you have a working GLMM sampler. The WinBUGS recommendation was pragmatic and correct for 2006; today Stan/brms handles all these models more efficiently with Hamiltonian Monte Carlo / the No-U-Turn Sampler (HMC/NUTS). The hierarchical centering trick is broadly useful and predates the paper, but this is a clean presentation. The prior sensitivity analysis is exemplary — testing four prior families and reporting that changes are small is exactly the right evidence to cite when defending prior choices. Methodological limitation: for large datasets (n >> 500, many random effects), WinBUGS becomes painfully slow and the authors acknowledge this. The paper also predates the integrated nested Laplace approximation (INLA) revolution (Rue-Martino-Chopin 2009) that now dominates Bayesian spatial/GAM modelling in practice.