Agricultural Statistics & Data Science · Statistical theory core · also Generalized Linear Models, Analysis of Categorical Data, Log-linear Models
Models for counts, proportions and categories: contingency tables, the exponential-family GLM with its link functions, logistic and Poisson regression, log-linear models, overdispersion, ordinal outcomes and generalised linear mixed models, applied to germination, mortality, pest counts and adoption data.
In the exams: AICE-JRF Ag Stat Unit 2 'contingency tables, log linear models' (heavy); GATE DA logistic regression (present); light in UPSC ISS
A large share of agricultural data are not continuous measurements with normal errors. Seeds germinate or not; insects die or survive at each dose; plants are infested or clean; farmers adopt a practice or not; pests are counted per trap; diseases are graded into ordered classes; survey respondents choose among unordered options. Such responses have variances that depend on their means, bounded ranges and discrete values, and least-squares analysis of raw proportions or counts can give fitted probabilities below 0 or above 1, wrong standard errors and misleading tests. Categorical data analysis and generalized linear models (GLMs) provide the correct framework.
GLMs, introduced by Nelder and Wedderburn in 1972, unify the normal linear model, logistic and probit regression for proportions, Poisson regression for counts, gamma regression for skewed positive data and log-linear models for contingency tables under one theory: an exponential-family response, a linear predictor and a link function, fitted by maximum likelihood through iteratively reweighted least squares. Extensions handle overdispersion, ordinal and nominal outcomes, excess zeros and random effects.
This manual covers contingency tables and odds ratios; the exponential family and the GLM; links, deviance and IRLS; logistic regression; Poisson regression; overdispersion; log-linear models; ordinal and multinomial models; zero-inflated and hurdle models; generalized linear mixed models; and model checking. It builds on Probability Distributions and Sampling Distributions (exponential family, binomial, Poisson, negative binomial), Statistical Inference I and II (maximum likelihood, Fisher information, Newton-Raphson and scoring; likelihood-ratio, Wald and score tests) and Linear Models and Regression (least squares, leverage and residual diagnostics), and refers to them rather than repeating them. Rank-based alternatives for ordinal scores are in Nonparametric and Rank-Based Methods, and random-effect estimation in Mixed Models and Variance Components.
The ICAR AICE-JRF/SRF Agricultural Statistics paper names contingency tables and log-linear models in Unit 2; GATE Data Science asks about logistic regression; and the methods are used in every bioassay, germination study, pest survey and adoption study, so research-institute recruitment interviews and ASRB examinations expect working knowledge of them.
A contingency table cross-classifies units by two or more categorical variables. For an r × c table with counts nᵢⱼ, row totals nᵢ., column totals n.ⱼ and grand total n, independence of the row and column variables gives expected counts Eᵢⱼ = nᵢ.n.ⱼ/n. Pearson's statistic X² = Σ(nᵢⱼ − Eᵢⱼ)²/Eᵢⱼ and the likelihood-ratio statistic G² = 2Σnᵢⱼ log(nᵢⱼ/Eᵢⱼ) are both approximately χ² with (r − 1)(c − 1) degrees of freedom under independence. X² is a score statistic and G² a likelihood-ratio statistic (see the hypothesis-testing manual); they agree closely when counts are large. With small expected counts Fisher's exact test (for 2 × 2) or exact conditional tests are used.
The same arithmetic tests homogeneity of proportions across groups (rows sampled separately) and independence (units classified on both variables after a single sample): the sampling schemes differ but the statistics and df coincide. Strength of association is summarised by Cramér's V = √[X²/(n × min(r − 1, c − 1))], between 0 and 1; significance alone says nothing about strength, since X² grows with n.
For ordinal variables, tests that use the order (a linear-by-linear association test, or a trend test for proportions across ordered doses, the Cochran-Armitage test) have one degree of freedom and are far more powerful than the general (r − 1)(c − 1) df test against monotone alternatives.
For a 2 × 2 table with cells a, b (row 1) and c, d (row 2), the proportions of 'success' in the two rows are p₁ = a/(a + b) and p₂ = c/(c + d). Three measures compare them: the difference p₁ − p₂; the relative risk p₁/p₂; and the odds ratio θ = [p₁/(1 − p₁)]/[p₂/(1 − p₂)] = ad/(bc). The odds ratio is the same whichever variable is treated as the response, is estimable from case-control (retrospective) samples, and is the natural parameter of logistic regression. Its log has approximate standard error √(1/a + 1/b + 1/c + 1/d), so a 95% interval is exp[log θ̂ ± 1.96 × SE]. Example: treated seed 80 germinated, 20 not; untreated 60 and 40. θ̂ = (80 × 40)/(20 × 60) = 2.67, SE(log θ̂) = √(1/80 + 1/20 + 1/60 + 1/40) = 0.323, and the interval is exp(0.981 ± 0.633) = 1.42 to 5.02. The relative risk of germination is 0.8/0.6 = 1.33: odds ratios and relative risks differ substantially when outcomes are common.
With a third variable (location, season), a 2 × 2 × K table can be analysed stratum by stratum. The Cochran-Mantel-Haenszel test checks conditional independence controlling for strata, and the Mantel-Haenszel estimator θ̂_MH = Σ(aₖdₖ/nₖ)/Σ(bₖcₖ/nₖ) estimates a common odds ratio; the Breslow-Day test checks homogeneity of odds ratios across strata. Collapsing over strata can reverse the direction of association (Simpson's paradox) when the strata differ both in the explanatory variable's distribution and in the response rate: a variety can have a higher germination rate in each of two seed lots yet a lower rate overall if most of its seed came from the poorer lot.
A generalized linear model has three components. The random component: independent responses yᵢ from an exponential-dispersion family with density f(y; θ, φ) = exp{[yθ − b(θ)]/a(φ) + c(y, φ)}, where θ is the natural (canonical) parameter and φ a dispersion parameter. The systematic component: a linear predictor ηᵢ = xᵢ′β. The link: a monotone differentiable function g with g(μᵢ) = ηᵢ, connecting the mean μᵢ = E(yᵢ) to the predictor.
The exponential-family form gives the mean and variance directly: E(y) = μ = b′(θ) and Var(y) = b″(θ)a(φ) = V(μ)a(φ), where V(μ) is the variance function. For the normal, V(μ) = 1 and a(φ) = σ²; for the Poisson, V(μ) = μ and φ = 1; for a binomial proportion from m trials, V(μ) = μ(1 − μ) and a(φ) = 1/m; for the gamma, V(μ) = μ² (constant coefficient of variation); for the inverse Gaussian, V(μ) = μ³. The variance function is what distinguishes the families in practice: counts whose variance grows with the mean, proportions whose variance is largest at 0.5, and positive measurements whose standard deviation is proportional to the mean.
The normal linear model is the special case with normal errors and identity link, so everything in the linear-models subject is a GLM. Unlike the transformation approach (analysing log counts or arcsine proportions by least squares), a GLM models the mean on its natural scale through the link while keeping the correct variance structure, and it handles zeros and unequal sample sizes without ad hoc corrections.
The canonical link sets η equal to the natural parameter θ: identity for the normal, logit for the binomial, log for the Poisson, reciprocal (inverse) for the gamma and inverse square for the inverse Gaussian. With the canonical link, X′y is a sufficient statistic for β, the likelihood equations take the simple form X′(y − μ) = 0, and observed and expected information coincide, so Newton-Raphson and Fisher scoring are the same algorithm (see the estimation manual). Non-canonical links (probit and complementary log-log for binomial data, log for gamma) are used when they fit better or have a mechanistic interpretation.
Parameters are estimated by maximum likelihood. Their large-sample covariance is φ(X′WX)⁻¹, with W the diagonal matrix of iterative weights at convergence, giving Wald tests and intervals; likelihood-ratio tests compare nested models through deviances (next chapter). When φ is unknown (normal, gamma, quasi-likelihood), it is estimated by the Pearson statistic divided by its degrees of freedom, and F tests replace χ² tests.
To see the exponential-family form in a familiar case, write the Poisson probability e^(−μ)μʸ/y! as exp[y log μ − μ − log y!]. Matching terms gives θ = log μ, b(θ) = e^θ, a(φ) = 1 and c(y) = −log y!. Then b′(θ) = e^θ = μ and b″(θ) = e^θ = μ, confirming that the mean and variance are equal and that the canonical link is the log. For a binomial proportion from m trials, θ = log[p/(1 − p)], b(θ) = log(1 + e^θ) and a(φ) = 1/m, so the canonical link is the logit and V(μ) = μ(1 − μ). These calculations are the standard examination exercise on GLM theory.
For binary and binomial data the common links are the logit, log[p/(1 − p)], whose coefficients are log odds ratios; the probit, Φ⁻¹(p), which arises from a normally distributed tolerance in dose-response models and gives coefficients about 1/1.6 to 1/1.8 of the logit ones, the exact ratio depending on the data; and the complementary log-log, log[−log(1 − p)], which is asymmetric and arises when the probability of 'at least one event' follows from a Poisson process, as in infection from a Poisson number of propagules or mortality from a Poisson number of hits. For counts the log link makes effects multiplicative (rate ratios); the identity link can be used when effects are additive but may give negative fitted means. For positive skewed data with constant CV, the gamma family with log link is a robust default.
The choice of link is partly substantive (what does the coefficient mean?) and partly empirical (which fits better?); links can be compared by deviance or by goodness-of-link tests. The link transforms the mean, not the data, which is the essential difference from analysing transformed responses.
The deviance compares a model's log-likelihood with that of the saturated model, which fits each observation exactly: D = 2φ[ℓ(saturated) − ℓ(model)], and the scaled deviance D/φ. For the normal it is the residual sum of squares; for the Poisson D = 2Σ[y log(y/μ̂) − (y − μ̂)]; for the binomial with mᵢ trials, D = 2Σ{yᵢ log(yᵢ/μ̂ᵢ) + (mᵢ − yᵢ) log[(mᵢ − yᵢ)/(mᵢ − μ̂ᵢ)]} with yᵢ and μ̂ᵢ as counts.
For nested models, the difference in deviances (scaled by φ) is the likelihood-ratio statistic, approximately χ² with df equal to the difference in numbers of parameters; with estimated φ, the F statistic [(D₀ − D₁)/df]/φ̂ is used. This 'analysis of deviance' replaces the analysis of variance for GLMs. Unlike the residual sum of squares in normal models, the deviance of a single model is only approximately χ²_(n−p), and only when counts are large (grouped binomial data with large mᵢ, Poisson counts that are not small); for binary data it is not a goodness-of-fit statistic at all (Chapter 11). AIC = −2ℓ + 2p compares non-nested models.
Example: a germination experiment with 4 seed treatments in 20 trays gives a residual deviance of 31.2 on 16 df for the treatment model and 52.8 on 19 df for the null model. The drop of 21.6 on 3 df is compared with χ²₃,₀.₀₅ = 7.81, so treatments differ. But the residual deviance 31.2 on 16 df is nearly twice its df, signalling overdispersion; with quasi-likelihood and φ̂ ≈ 1.95, the F statistic is (21.6/3)/1.95 = 3.69 on (3, 16) df, still significant but much less overwhelming than the χ² test suggested.
GLMs are fitted by Fisher scoring, which takes the form of repeated weighted least squares. At each iteration, compute the working response zᵢ = ηᵢ + (yᵢ − μᵢ)g′(μᵢ) and the working weights wᵢ = 1/[V(μᵢ)a(φ)g′(μᵢ)²], regress z on X with weights W, and update β = (X′WX)⁻¹X′Wz. Starting values come from the data (μ = y, with adjustments for zeros). For a Poisson log-linear model, g′(μ) = 1/μ, so wᵢ = μᵢ and zᵢ = ηᵢ + (yᵢ − μᵢ)/μᵢ; for a logit model with mᵢ trials, wᵢ = mᵢpᵢ(1 − pᵢ).
Convergence is fast for well-behaved problems. Failures signal data problems: complete or quasi-complete separation in logistic regression (a predictor that perfectly separates successes from failures), where the MLE is infinite and estimates and standard errors grow without bound; zero cells in log-linear models; and fitted means hitting the boundary. Penalised likelihood (Firth's bias-reduced logistic regression) or exact methods handle separation. Finney's probit analysis, computed by hand with working probits and weights before computers, is this same algorithm.
For a binary response, or yᵢ successes out of mᵢ trials, logistic regression models log[pᵢ/(1 − pᵢ)] = xᵢ′β. A one-unit increase in xⱼ multiplies the odds of success by exp(βⱼ), holding other predictors fixed: exp(βⱼ) is an odds ratio. With a single binary predictor, exp(β₁) equals the sample odds ratio of the 2 × 2 table. The fitted probability is p = exp(η)/(1 + exp(η)), always between 0 and 1. On the probability scale the effect of xⱼ is not constant: the slope of p with respect to xⱼ is βⱼp(1 − p), largest (βⱼ/4) at p = 0.5.
Example: in a germination study, log-odds = −1.2 + 0.8 × (seed treated). The odds ratio for treatment is e^0.8 = 2.23. Untreated seed has p = 1/(1 + e^1.2) = 0.231; treated seed has η = −0.4 and p = 0.401. A 2.23-fold change in odds is here a 17-percentage-point change in probability; at other baseline levels the same odds ratio gives a different difference in probability, which is why results should be reported on both scales.
With log dose as predictor, logit p = β₀ + β₁ log(dose) describes quantal dose-response for insecticides, fungicides and herbicides. The median lethal dose (LD50, or LC50 for concentrations) is the dose at which p = 0.5, i.e. log LD50 = −β₀/β₁. Example: logit p = −4.0 + 2.0 log₁₀(dose) gives log₁₀ LD50 = 2.0 and LD50 = 100 units. Its standard error comes from the delta method or, better, Fieller's theorem, which gives asymmetric limits and fails gracefully when the slope is poorly determined. Other effective doses follow from logit(p) = β₀ + β₁ log dose: the LD90 has log dose (log 9 − β₀)/β₁.
Probit models (Finney) give similar LD50s and slightly different tails; natural (control) mortality is handled by Abbott's correction or by including a background-mortality parameter in the model. Relative potency of two preparations is the horizontal distance between parallel dose-response lines on the log scale, and testing parallelism (equal slopes) comes first.
Model assessment for binary data uses the Hosmer-Lemeshow test (grouping observations by fitted probability), classification tables and the area under the ROC curve (concordance between predicted probabilities and outcomes); for grouped binomial data with adequate mᵢ, the deviance and Pearson statistics are approximately χ². Separation, sparse data and many predictors relative to events are the main practical pitfalls; a common rule of thumb asks for at least about ten events per predictor, though recent work shows the requirement depends on the setting.
For counts yᵢ, Poisson regression models log μᵢ = xᵢ′β. Effects are multiplicative: exp(βⱼ) is the rate ratio for a one-unit increase in xⱼ. Counts collected over different exposures (trap-days, leaf area, plot size, years of observation) are modelled with an offset: log μᵢ = log tᵢ + xᵢ′β, so that the model describes the rate μᵢ/tᵢ. A trap left for 7 days and another for 14 have log 7 and log 14 as offsets, a fixed term with coefficient 1.
Example: aphid counts with log μ = 1.5 + 0.4 × (unsprayed) give a rate ratio of e^0.4 = 1.49: unsprayed plants carry 49% more aphids. With the canonical log link, the MLE satisfies X′(y − μ̂) = 0; with an intercept, Σμ̂ᵢ = Σyᵢ, so fitted totals match observed totals. The equidispersion assumption Var(y) = μ is strong and often violated by field counts, so the first check after fitting is the ratio of the Pearson statistic (or deviance) to its degrees of freedom: values well above 1 indicate overdispersion (next chapter).
Poisson regression and contingency-table analysis are two views of the same model: the cell counts of a table can be treated as independent Poisson variables, and log-linear models for tables (Chapter 7) are Poisson regressions with categorical predictors. Conditioning on the total converts Poisson to multinomial sampling without changing the estimates of association parameters, which is why the same software analyses both.
Typical agricultural applications are counts of insects per trap or plant, weed seedlings per quadrat, lesions per leaf, nematodes per sample, and incidence of disease or deaths in livestock herds per animal-year. When counts are large, a normal model on log counts can be adequate; with small counts and many zeros, the Poisson or its overdispersed extensions are needed, and zero-inflated models (Chapter 9) may be appropriate when zeros arise from a separate process.
A worked comparison shows why the offset matters. Two kinds of sticky trap catch 42 thrips in 6 days and 50 thrips in 10 days. Comparing raw counts suggests the second trap is better; comparing rates (7.0 and 5.0 per day) shows the opposite, a rate ratio of 1.4 in favour of the first. In a Poisson GLM with log(days) as offset, the coefficient for trap type estimates log(1.4) = 0.336, and its standard error √(1/42 + 1/50) = 0.209 gives a Wald z of 1.61. Fitting the same model without the offset would estimate the ratio of counts, confounding trap efficiency with exposure time.
Overdispersion means the variance exceeds what the model implies: Var(y) > μ for Poisson counts, or Var(y) > mp(1 − p) for binomial counts. It arises from unobserved heterogeneity (plants, traps or plots differ in their underlying rates), clustering (seeds in the same tray share conditions), and omitted predictors. Ignoring it leaves the coefficient estimates roughly unbiased but makes standard errors too small and tests too liberal, so spurious 'significant' effects appear.
The simplest remedy is quasi-likelihood: keep the mean model but assume Var(y) = φV(μ) with an unknown dispersion φ. The estimates of β are unchanged, φ is estimated by the Pearson statistic divided by its residual degrees of freedom, φ̂ = X²/(n − p), standard errors are multiplied by √φ̂, and F tests replace χ² tests. With φ̂ = 4, standard errors double and Wald z-values halve. Quasi-Poisson and quasi-binomial models are routine in pest-count and germination analyses. Quasi-likelihood requires only the mean and variance assumptions, not a full distribution, so likelihood-based quantities such as AIC are not defined for it.
The negative binomial distribution, a Poisson whose mean varies according to a gamma distribution (see the distributions manual), gives a full likelihood for overdispersed counts with Var(y) = μ + μ²/k (the NB2 form), where k is the aggregation parameter (some software reports α = 1/k). Its variance grows quadratically in μ, whereas quasi-Poisson's grows linearly (φμ). The two weight observations differently: negative binomial gives relatively less weight to large counts. Plotting squared residuals against fitted means helps choose between them. The negative binomial is the standard model for aggregated insect and weed counts, linking to the k used in sequential sampling plans.
For overdispersed binomial data the beta-binomial model lets p vary across clusters according to a beta distribution, giving Var(y) = mp(1 − p)[1 + (m − 1)ρ], with ρ the intra-cluster correlation; Williams' method gives a quasi-likelihood version. An observation-level random effect in a GLMM (Chapter 10) is another flexible way to absorb overdispersion. The choice should be guided by how the data were generated: when trays, plants or plots are natural clusters, modelling them as random effects is usually the most interpretable route.
Log-linear models describe the structure of association in contingency tables by modelling the logarithm of expected cell counts. For a two-way table, independence is log μᵢⱼ = λ + λᵢ^X + λⱼ^Y, and the saturated model adds the interaction λᵢⱼ^XY, which carries all the association (the log odds ratios). With constraints such as sum-to-zero or a reference category, the independence model has 1 + (r − 1) + (c − 1) parameters and its residual degrees of freedom are (r − 1)(c − 1), the familiar chi-square df.
For three variables X, Y and Z, hierarchical models are named by their highest-order terms: (X, Y, Z) mutual independence; (XY, Z) joint independence of Z from (X, Y); (XZ, YZ) conditional independence of X and Y given Z; (XY, XZ, YZ) homogeneous association, in which each pair is associated with conditional odds ratios that do not depend on the level of the third variable; and (XYZ) the saturated model. Hierarchical models include all lower-order terms of any interaction they contain. Models are compared by differences in G², and the most parsimonious adequate model is chosen.
Degrees of freedom follow from counting parameters. For an I × J × K table, the residual df of (X, Y, Z) is IJK − [1 + (I − 1) + (J − 1) + (K − 1)]; each added two-factor term removes (I − 1)(J − 1)-type parameters. For a 2 × 2 × 2 table the mutual-independence model has 8 − 4 = 4 df, (XZ, YZ) has 2 df and (XY, XZ, YZ) has 1 df, the single df for the three-factor interaction. In a survey of farm type (2) by irrigation (2) by adoption (2), testing (XY, XZ, YZ) against the saturated model asks whether the association between irrigation and adoption is the same for both farm types.
Conditional independence (XZ, YZ) means that X and Y are independent within each level of Z, although they may be associated marginally, which is the log-linear expression of Simpson's paradox. Collapsibility conditions state when association between X and Y can be studied in the marginal table without bias: under (XY, XZ) or (XY, YZ), the XY odds ratios are the same conditionally and marginally.
Log-linear models are fitted by maximum likelihood as Poisson GLMs, or by iterative proportional fitting (IPF), which repeatedly scales fitted counts to match the observed margins corresponding to the model's highest-order terms; the MLEs match those margins exactly. Zero cells can make some parameters inestimable (sampling zeros) or must be fixed at zero (structural zeros, as in cells that cannot occur).
When one variable is a binary response, every log-linear model that includes the full interaction among the explanatory variables corresponds to a logistic regression: for example, (XY, XZ, YZ) with Y the response is equivalent to logit P(Y) = α + βᵢ^X + βₖ^Z. Logistic regression is usually the more direct analysis when there is a clear response; log-linear models are preferred when all variables are on an equal footing, as in describing associations among farm type, crop choice and irrigation source in a survey. AICE-JRF examines the model hierarchy, df and interpretation.
Many agricultural responses are ordered categories: disease severity grades, damage classes, quality grades, Likert-type survey answers. The proportional-odds (cumulative logit) model writes log[P(Y ≤ j)/P(Y > j)] = θⱼ − x′β for j = 1, ..., J − 1, with ordered cut-points θⱼ and a single set of effects β shared across all cumulative splits. exp(βₖ) is a cumulative odds ratio: the factor by which the odds of being in a higher category (with this sign convention) change per unit of xₖ, the same for every cut-point. The model arises from a latent continuous variable with logistic errors observed only in intervals; the cumulative probit model uses normal errors instead.
The proportional-odds assumption should be checked, for example by a score test or by fitting separate binary logits at each cut-point and comparing slopes; when it fails, partial proportional odds models relax it for some predictors. Alternatives are the adjacent-category logit model, comparing each category with the next, and the continuation-ratio model, comparing each category with all higher ones, natural for sequential stages such as germination, emergence and establishment. Rank-based methods (see the nonparametric manual) are a model-free alternative for simple comparisons.
Example: in a fungicide trial, disease is graded none, light, moderate or severe. A proportional-odds model gives a fungicide coefficient of −1.1 in the parameterisation log[P(Y ≤ j)/P(Y > j)] = θⱼ − x′β with x = 1 for treated plots, so treated plots have exp(1.1) = 3.0 times the odds of being in a lower (less diseased) category at every cut-point: 'none' versus the rest, 'none or light' versus the rest, and so on. Reporting this single odds ratio, with fitted category probabilities for treated and untreated plots, is clearer than averaging grade codes. Because sign conventions for β differ between software packages, the direction of an effect should always be confirmed from fitted probabilities.
For an unordered response with J categories, such as the choice among several crops, marketing channels or pest-management options, the baseline-category (multinomial) logit model uses one category as reference and fits J − 1 equations log[P(Y = j)/P(Y = J)] = αⱼ + x′βⱼ, each with its own coefficients. The number of parameters therefore grows with J. Fitted probabilities are P(Y = j) = exp(ηⱼ)/Σₖ exp(ηₖ) with η_J = 0. The model implies independence of irrelevant alternatives (IIA): the odds between any two categories do not depend on which other categories are available, an assumption that can fail when alternatives are close substitutes (two similar improved varieties), in which case nested logit or mixed logit models from econometrics are used.
A worked interpretation: in a model of farmers' choice among maize, sorghum and pigeonpea (reference), a coefficient of 0.7 for irrigation access in the maize equation means irrigated farmers have exp(0.7) = 2.0 times the odds of choosing maize rather than pigeonpea. Applying an ordinal model to unordered categories, or a nominal model to ordered ones, loses validity in the first case and power in the second.
Field counts often contain more zeros than a Poisson or even a negative binomial model predicts: many fields have no pest because it has not arrived, many plants show no lesions because they escaped infection, many farm households report no sales of a crop because they do not grow it. Two model families handle such data.
Zero-inflated models (ZIP, ZINB) are mixtures: with probability π an observation is a 'structural' zero; otherwise it comes from a count distribution, which can itself produce zeros. So P(Y = 0) = π + (1 − π)f(0) and P(Y = y) = (1 − π)f(y) for y > 0. The mixing probability π is modelled with a logistic regression and the count mean with a log-linear model, each with its own predictors. For ZIP, E(Y) = (1 − π)λ and Var(Y) = (1 − π)λ(1 + πλ), so zero inflation itself causes overdispersion. Example: with π = 0.3 and λ = 2, P(0) = 0.3 + 0.7e^(−2) = 0.395, compared with 0.135 for a plain Poisson with λ = 2.
Hurdle (two-part) models separate the process into a binary part, zero versus positive, and a zero-truncated count distribution for the positive values. All zeros come from the binary part. They suit situations in which a threshold must be crossed before counts can occur (a pest must establish before it multiplies; a household must decide to sell before the quantity is chosen). Zero-inflated models suit situations in which zeros come from two sources, structural and sampling.
The models can fit similarly and lead to different interpretations, so the choice should be guided by the data-generating story. Comparison by AIC or likelihood is common; the Vuong test was widely used to compare zero-inflated with ordinary count models, but its use for that purpose has been criticised because the models are nested at a boundary. Before adopting a zero-inflated model it is worth checking whether a negative binomial model already predicts the observed zeros adequately, since apparent excess zeros are often just overdispersion. Plotting observed and predicted frequencies of each count value is the most informative check.
Generalized linear mixed models (GLMMs) add random effects to the linear predictor: g(μᵢ | b) = xᵢ′β + zᵢ′b, with b ~ N(0, G). They handle the clustered and blocked structure of agricultural experiments with non-normal responses: random block, plot or tray effects for germination proportions, random field and trap effects for insect counts, random farm or village effects in adoption surveys, and random location-year effects in multi-environment trials with binary or count traits. They also absorb overdispersion through an observation-level random effect. The variance-component theory is in the Mixed Models and Variance Components subject; this chapter covers what changes when the response is non-normal.
The likelihood requires integrating over the random effects, which has no closed form except in special cases. Methods include the Laplace approximation, adaptive Gauss-Hermite quadrature (accurate, feasible for few random effects per cluster), penalised quasi-likelihood (PQL, fast but biased towards zero for variance components and coefficients with binary data and small clusters), and Bayesian MCMC.
In a logistic GLMM the fixed effects are conditional (cluster-specific): exp(βⱼ) is the odds ratio for units in the same cluster. Averaging over clusters gives population-averaged (marginal) effects, which are smaller in magnitude; for the logit link with random-intercept variance σ², the marginal coefficient is approximately β/√(1 + 0.346σ²), so with σ² = 2 it is about 0.77β. Generalized estimating equations (GEE) estimate marginal effects directly, using a working correlation matrix and robust (sandwich) standard errors, without specifying the random-effect distribution. The choice depends on the question: a breeder interested in the effect of a treatment on a given plot wants conditional effects; a policy analyst interested in the average adoption rate across villages wants marginal effects. For log-link Poisson models with random intercepts, only the intercept differs between the two interpretations; the slopes coincide.
Practical points: random effects need enough levels (variance components from four blocks are poorly estimated); singular fits with variance estimates of zero are common and should be reported, not hidden; and tests of variance components face the boundary problem (see the hypothesis-testing manual).
Residuals for GLMs must account for the non-constant variance. The Pearson residual is rᵢ = (yᵢ − μ̂ᵢ)/√V(μ̂ᵢ) (times 1/√φ when φ ≠ 1); the sum of their squares is Pearson's X². The deviance residual is sign(yᵢ − μ̂ᵢ)√dᵢ, where dᵢ is the observation's contribution to the deviance; their squares sum to the deviance. Standardised versions divide by √(1 − hᵢ), with hᵢ the leverage from the weighted hat matrix. Plots of residuals against fitted values (on the linear-predictor scale) and against each predictor reveal wrong links, missing terms and outliers; half-normal plots of absolute residuals with simulated envelopes check the distributional assumption.
For discrete data with small counts or binary responses, Pearson and deviance residuals fall in bands and are hard to read. Randomised quantile residuals (Dunn and Smyth), which map each observation to a normal quantile with a random component for discreteness, are approximately standard normal when the model is right and are the recommended tool, often computed by simulation. Example: a Poisson observation y = 7 with fitted mean 4 has Pearson residual (7 − 4)/2 = 1.5 and deviance residual √{2[7 log(7/4) − 3]} = √1.83 = 1.35.
For grouped data with large counts, the deviance and X² are approximately χ²_(n−p) and test overall fit; values much larger than the df indicate lack of fit or overdispersion. For binary data and sparse tables these approximations fail, and the Hosmer-Lemeshow test (grouping by fitted probability, typically into ten groups, referred to χ² with g − 2 df), calibration plots and the area under the ROC curve are used. Information criteria (AIC, BIC) compare competing models, including different families and links.
Influence diagnostics extend those of linear models: leverage from the weighted hat matrix W^(1/2)X(X′WX)⁻¹X′W^(1/2), Cook's distance, and the change in coefficients when an observation is deleted (approximated by one-step formulas). In bioassays a single anomalous dose group can dominate the slope and LD50; in log-linear models a few large residual cells can point to a missing interaction. Model checking should end with a substantive check: do fitted probabilities and counts make sense, and do predictions hold up on data not used for fitting?
A typical sequence for pest counts in a blocked trial is: fit a Poisson GLMM with block as random effect and log trap-days as offset; check the ratio of Pearson X² to residual df; if it is well above 1, refit with a negative binomial or an observation-level random effect; compare observed and predicted frequencies of zeros and of large counts; inspect randomised quantile residuals against fitted values and treatments; and only then test treatments and report rate ratios with intervals. Each step guards against a failure that would otherwise make standard errors or conclusions unreliable.
A student needs to compute and interpret odds ratios with confidence intervals, test independence in r × c tables, write a GLM's three components and canonical link for each family, explain deviance and IRLS, interpret logistic and Poisson coefficients and offsets, compute an LD50, recognise overdispersion and choose quasi-likelihood or negative binomial models, read the hierarchy of log-linear models, and fit ordinal and multinomial models. Theses on germination, mortality, pest counts or adoption need all of these.
A statistician five years into work at a research station fits GLMs every week: logistic and probit bioassays with LD50s and Fieller limits, binomial GLMMs for germination in blocked trials, negative binomial and quasi-Poisson models for insect counts with trap-day offsets, zero-inflated models for pest incidence surveys, proportional-odds models for disease scores, and multinomial logits for farmers' choices. Checking overdispersion, separation and model fit with randomised quantile residuals is routine.
A senior statistician with twenty-five years of experience sets analysis standards: which models are used for efficacy trials of pesticides and seed treatments, how dose-response results are reported to regulators, how survey data on adoption and practices are modelled with clustering, and how results are expressed on scales decision makers understand (probabilities and rate ratios, not only log-odds). Such a statistician reviews claims built on proportions and counts and catches analyses that ignore overdispersion or clustering.
A PhD with long experience develops models for categorical and count data in agriculture: GLMMs with spatial and genetic random effects for multi-environment trials of disease scores, zero-inflated and hurdle models for ecological counts, ordinal models for genomic prediction of scored traits, and efficient estimation and testing for these models. Such a researcher proves properties, checks small-sample behaviour by simulation and supervises doctoral students.
Derivations are sketched rather than given. Not covered in depth: exact conditional logistic regression, Bayesian GLMs, generalized additive models, survival and time-to-event models, models for compositional or circular data, and the econometric discrete-choice extensions (nested and mixed logit) beyond a mention. The probit-logit coefficient ratio, the events-per-variable rule and Hosmer-Lemeshow grouping are rules of thumb that vary between sources. Regulatory requirements for bioassay reporting differ by country and are not quoted. Statistical tables are not reproduced; questions state the values they need.
Original practice questions written to the syllabus and pattern of the exams named in each question's examTags. They are not past papers and do not reproduce any board's questions. Every numerical answer was computed with a script; the values needed are stated in the question.
These notes are free to read here, with no account.
What the Laisha app adds:
Free for ten days, then a subscription. We would rather write that here than have you find it later.