Linear Models and Regularisation

Generalised linear models

GLMs are one recipe with three swappable parts — a linear score, a link that bends it to the target's range, and a noise family — uniting linear, logistic and Poisson regression.

On this page 5
  1. Why it exists
  2. How it works
  3. A real example you have seen
  4. Remember this
  5. What to learn next

One lesson, three depths. Pick the one that fits you today — you can switch any time.

Beginner — No maths. Plain English.

A generalised linear model is one recipe with three swappable parts.

A straight-line score. A bend that maps that score into the right range. And a choice of how the data wobbles around it.

Think of a kitchen mixer with attachments. One machine; swap the whisk for the dough hook and it makes bread instead of cream. You do not buy a new machine per dish.

Linear regression and logistic regression look like different machines. They are the same mixer with different attachments — and seeing that unlocks a whole shelf of further attachments.

Why it exists

Plain linear regression predicts any number from minus infinity to plus infinity. Real targets often refuse that range. A probability lives between 0 and 1. A count of shop visits cannot be negative. Predicting "-2.3 customers" is not a small error; it is a category error.

The generalised linear model (GLM) fixes this with three choices:

  1. The linear score. Multiply each feature by a coefficient and add up — the familiar part, kept in every GLM.
  2. The link. A bending function that maps the score onto the target's legal range. An S-curve squeezes scores into 0-to-1 for probabilities. An exponential bend keeps counts positive.
  3. The family. A description of how observations wobble around the prediction. Yes/no outcomes wobble like coin flips; counts wobble like arrivals at a counter; measurements wobble like a bell curve.

Pick bell curve + no bend: linear regression. Coin flip + S-curve: logistic regression. Arrivals + exponential: Poisson regression. One machine, many attachments.

How it works

features ──→ linear score ──→ bend (link) ──→ prediction in legal range
                                                ▲
              family: how real data wobbles ────╯

measurement (any number) : no bend        + bell-curve wobble
probability (0 to 1)     : S-curve bend   + coin-flip wobble
count (0, 1, 2, ...)     : exponential    + arrivals wobble

A real example you have seen

Insurance runs on GLMs to this day. How many claims will this driver file (a count)? Will this policy lapse (a yes/no)? How large will a claim be (a positive amount)? Three questions, three attachments, one modelling machine that regulators can audit line by line.

Remember this

  • A GLM = linear score + link + family — three parts, each swappable.
  • Linear, logistic and Poisson regression are the same machine with different attachments.
  • Choose the family by asking: what range does my target live in, and how does it wobble?

What to learn next

Developer — Code and libraries.

Setup

bash
pip install statsmodels

Outputs verified with statsmodels 0.14.6 on CPU.

One machine, the logistic attachment

glm_driving_test.py
import numpy as np
import statsmodels.api as sm

# hours of driving practice -> passed the test (1) or not (0)
hours = np.array([2., 4., 5., 6., 8., 10., 12., 14., 16., 18., 20., 24.])
passed = np.array([0, 0, 0, 0, 0, 1, 0, 1, 1, 1, 1, 1])

X = sm.add_constant(hours)
model = sm.GLM(passed, X, family=sm.families.Binomial()).fit()
print("intercept, slope:", model.params.round(3))

new = sm.add_constant(np.array([6., 11., 18.]))
for h, p in zip((6, 11, 18), model.predict(new)):
    print(f"{h:2d} hours of practice -> P(pass) = {p:.2f}")
Output
intercept, slope: [-7.396  0.67 ]
 6 hours of practice -> P(pass) = 0.03
11 hours of practice -> P(pass) = 0.49
18 hours of practice -> P(pass) = 0.99

The walkthrough

family=sm.families.Binomial() is the attachment. It selects coin-flip wobble, and brings its default link — the logit, the S-curve's inverse — along with it. Swap that one argument for sm.families.Poisson() and the identical code fits counts. That swap-one-word property is the entire practical payoff of the GLM view.

sm.add_constant supplies the intercept. statsmodels does not add one silently. Forgetting it forces the S-curve through a fixed point and quietly wrecks the fit — the most common statsmodels bug in existence.

Reading the coefficients. The linear score here is -7.396 + 0.67 × hours, in log-odds — the link's units. Exponentiate to speak human: exp(0.67) ≈ 1.95, so each practice hour nearly doubles the odds of passing. GLM coefficients are read through their link: multiplicative odds for logit, multiplicative rates for Poisson's log link.

The predictions respect the range. 0.03, 0.49, 0.99 — all legal probabilities, the S-curve saturating at both ends. A straight line through this data would cheerfully predict 1.2 for a 30-hour student.

Where is scikit-learn? LogisticRegression and PoissonRegressor fit the same models, geared for prediction pipelines. statsmodels gives the statistician's console — standard errors, confidence intervals, deviance — via model.summary(). Choose by what you need to read.

Common mistakes

Forgetting add_constant on new data too. Predicting on raw [6, 11, 18] without the constant column mismatches shapes or, worse, misaligns columns. Every design matrix, training or prediction, gets the same treatment.

Reading logit coefficients as probability changes. "0.67 means +67% chance per hour" is wrong twice over. It is log-odds, and its probability effect varies along the curve — steep near 0.5, flat near the ends.

Regularisation by default. scikit-learn's LogisticRegression penalises unless told penalty=None; statsmodels' GLM does not penalise at all. Coefficients from the two will differ until you align the settings — a perennial "why don't these match" mystery.

Choosing the family by habit rather than by target. Gaussian on counts, Gaussian on strictly-positive skewed amounts — both fit and both lie. The target's range and wobble pattern choose the family; convenience does not.

Try it yourself

Call model.summary() and find the slope's confidence interval — does it exclude zero? Then refit the same data with family=sm.families.Gaussian() (the wrong attachment) and compare predictions at 2 and 24 hours with the logistic ones. Where does the straight line leave the legal range?

What to learn next

Researcher — Mathematics and papers.

The formal recipe

A GLM (Nelder and Wedderburn, 1972, Generalized linear models) has three components:

  1. Random: $y_i$ drawn from an exponential-family distribution with mean $\mu_i$.
  2. Systematic: linear predictor $\eta_i = x_i^\top \beta$.
  3. Link: $g(\mu_i) = \eta_i$, for a smooth monotone $g$.

The exponential family in canonical form:

$$ f(y; \theta, \phi) = \exp!\left( \frac{y\,\theta - b(\theta)}{a(\phi)} + c(y, \phi) \right) $$

Where:

  • $\theta$ — the natural parameter; $b(\theta)$ — the log-partition (cumulant) function.
  • $\phi$ — the dispersion parameter; $a, c$ — family-specific functions.
  • Mean and variance follow from $b$: $\mu = b'(\theta)$, $\operatorname{Var}(y) = a(\phi)\, b''(\theta) = a(\phi) V(\mu)$.

The variance function $V(\mu)$ is each family's signature: constant for Gaussian, $\mu(1-\mu)$ for Bernoulli, $\mu$ for Poisson, $\mu^2$ for Gamma. The canonical link is $g = (b')^{-1}$ — identity, logit, log respectively — which makes $X^\top y$ sufficient and the log-likelihood concave in $\beta$.

Fitting: IRLS

Maximum likelihood proceeds by iteratively reweighted least squares: at each step, form working responses $z_i = \eta_i + (y_i - \mu_i)\, g'(\mu_i)$ and weights $w_i = [V(\mu_i)\, g'(\mu_i)^2]^{-1}$, then solve the weighted least-squares problem:

$$ \beta \leftarrow (X^\top W X)^{-1} X^\top W z $$

This is exactly Newton–Raphson with the Fisher information (Fisher scoring); cost $O(nd^2 + d^3)$ per iteration, typically < 10 iterations. Every logistic regression you have ever fitted converged through some variant of this loop.

Inference and fit quality

Asymptotics: $\hat\beta \sim \mathcal{N}(\beta, (X^\top W X)^{-1} a(\phi))$, giving Wald tests; nested models compare via the deviance $D = 2[\ell_{sat} - \ell_{model}]$, with $\Delta D \sim \chi^2$ under the null. Overdispersion — variance exceeding $V(\mu)$'s promise, endemic in real counts — is handled by quasi-likelihood (Wedderburn, 1974), which needs only a mean-variance relationship, or by richer families (negative binomial); see the Poisson lesson.

Descendants

  • GAMs (Hastie and Tibshirani, 1990): replace $x^\top\beta$ with smooth functions $\sum_j f_j(x_j)$ — nonlinearity, link and family kept.
  • GLMMs: add random effects for grouped data (Breslow and Clayton, 1993).
  • Penalised GLMs: the ridge and lasso penalties graft directly onto the likelihood (glmnet).
  • A neural network with a sigmoid output and cross-entropy loss is a GLM whose linear predictor grew hidden layers — the family/link discipline survives intact in every modern classifier's final layer.

Reference text: McCullagh and Nelder (1989), Generalized Linear Models, 2nd ed.

What to learn next