Count Models
Regression for count data, threaded by one problem — extra-Poisson variation, the near-universal fact that real counts vary more than a plain Poisson allows, and the different ways to model it. It starts with the two core GLMs, Poisson and Negative Binomial regression, where the NB's dispersion parameter restores calibrated inference under overdispersion. The hierarchical Poisson GLMM adds group random effects on the log-rate, and the hierarchical binomial GLMM is its twin on grouped proportions — together completing the Gaussian / Poisson / binomial hierarchical-GLM grid. Zero-inflated Poisson and Negative Binomial handle counts with far more zeros than any Poisson can produce, splitting structural from sampling zeros via a data-augmentation Gibbs sampler. And the Besag–York–Mollié spatial model closes the set — disease mapping where an intrinsic-CAR prior pools each area toward its map neighbours rather than a global mean. Each is built from scratch and cross-checked across PyMC and R.
Bayesian Count Regression — Poisson & Negative Binomial
A from-scratch implementation of the two core count GLMs — Poisson
(equidispersion) and Negative Binomial (NB2, with a dispersion parameter
r) — each fit both by maximum likelihood (Newton/IRLS) and Bayesian
random-walk Metropolis. The throughline is overdispersion: on equidispersed
data MLE, Bayes, PyMC, and R glm coincide exactly (Bernstein–von Mises), but on
overdispersed data the Poisson keeps β roughly right while its standard errors run ~1.7× too
small — and NB recovers the dispersion and restores calibrated inference. Applied to the
Epil epilepsy trial (Thall & Vail 1990; Pearson dispersion ≈4.71), where
Poisson reports no treatment effect but NB finds a 21% seizure reduction (RR 0.79) — a gap
traced to a single influential high-count patient that NB down-weights. Validated across four
engines (from-scratch MLE & RW-Metropolis, PyMC, R glm/MASS::glm.nb),
with the correlated repeated measures motivating a hierarchical-count follow-up
(Breslow–Clayton 1993).
Bayesian Hierarchical Poisson
The hierarchical extension of count regression — a Poisson GLMM with group
random effects on the log-rate (log μi = xi'β + bgi),
a log-normal mixing that absorbs within-group correlation and overdispersion. Implemented
from scratch via Metropolis-within-Gibbs (no Pólya–Gamma for Poisson), across three
applications spanning the conjugacy spectrum: Epil (normal random intercept,
non-conjugate) — where adding a patient random effect for the 4 correlated visits
overturns the single-level treatment conclusion, pushing the progabide effect onto
the significance boundary (the classic Breslow–Clayton pseudo-replication lesson);
Pumps (gamma rate, conjugate) — an almost-pure-Gibbs sampler that is the
explicit-random-effects form of the Negative Binomial;
and Rugby (two crossed attack/defence effects) — recovering the 2014 Six
Nations standings from 15 matches. Each application is cross-checked three ways — from-scratch
sampler, PyMC NUTS, and R (lme4::glmer or the NB marginal).
Bayesian Hierarchical Binomial
A binomial GLMM — logistic regression on grouped counts with a group random
effect that absorbs over-dispersion — the binomial twin of the
hierarchical Poisson (Epil), completing the
hierarchical-GLM grid (Gaussian / Poisson / binomial). Grouped binomial counts routinely vary
more than a single logistic curve allows; a per-group random effect bi
on the logit soaks up that extra-binomial spread — the logistic-scale analogue of the Negative
Binomial's gamma mixing. Implemented from scratch via the same Metropolis-within-Gibbs as the
Poisson GLMM, binomial likelihood swapped in. Three applications: Seeds (BUGS;
21 plates, 2×2 factorial) — over-dispersion in a designed experiment (the Breslow–Clayton
interval-widening point); Surgical (BUGS; 12 hospitals) — the intercept-only
"league table" problem, where partial pooling stabilises uneven-volume mortality rates (a raw
0/47 = 0% pulled to ~6%) and shows the ranking is genuinely uncertain; and Baseball
(Efron–Morris) — James–Stein shrinkage two ways, contrasting the logit-normal GLMM with a
conjugate beta-binomial (the binomial mirror of the gamma-Poisson/Pumps vs
Poisson-GLMM/Epil conjugacy split). Cross-checked against PyMC and R
(lme4::glmer, VGAM).
Bayesian Zero-Inflated Count Regression
A mixture count model for data with far more zeros than any Poisson can produce —
zero-inflated Poisson (ZIP) and Negative Binomial (ZINB) —
separating structural zeros (never at risk) from sampling zeros (at risk but
caught none). Two linked regressions (a log-link count model and a logit zero-inflation model)
are fit from scratch by a data-augmentation Gibbs sampler: drawing a latent
regime indicator si decouples the awkward mixture into two ordinary
GLMs — a logistic regression and a count regression, reusing the
logit and
Poisson/NB machinery from earlier examples. Applied
to the UCLA/PyMC fish example (250 park groups, 57% zeros, variance ≈40× the mean): the data
demand zero-inflation, and ZINB absorbs both the excess zeros and the overdispersion
(r ≈ 0.35). Children roughly halve the catch, a camper raises it ~2.3×. Cross-checked against
PyMC (built-in zero-inflated likelihoods) and R pscl::zeroinfl.
Bayesian Hierarchical Spatial Poisson Model
The spatial extension of the hierarchical Poisson — the classic
Besag–York–Mollié (BYM) disease-mapping model (Clayton & Kaldor 1987;
the WinBUGS "Lip" example), fit from scratch by Metropolis-within-Gibbs. Each area carries two
random effects: an unstructured term θi and a
structured spatial term φ with an intrinsic CAR prior
that pulls each district toward the average of its map neighbours. Where the earlier
hierarchical-Poisson models pooled toward a global
mean, the ICAR prior is "pooling with a map" — borrowing strength from adjacent areas.
Applied to Scotland lip cancer (56 districts): outdoor work raises risk (RR ≈ 1.4 per +10%
agriculture/fishing/forestry), and the variation is almost entirely spatial
(σφ ≫ σθ) — turning a noisy raw SMR map into a coherent north-vs-central-belt
risk surface. Cross-checked against PyMC (pm.ICAR) and R
(CARBayes::S.CARbym).