Resampling, Likelihood and Bayes

Maximum likelihood estimation

MLE picks the parameter values under which your data would have been least surprising — the single fitting principle behind regression, classification and deep learning losses.

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.

Maximum likelihood estimation picks the explanation under which your data would have been least surprising.

You come home to a torn cushion, scattered stuffing, and a dog looking guilty. Two suspects: the dog, or a burglar. Under "dog did it", this scene is completely ordinary. Under "burglar", it is bizarre — burglars take things; they do not chew cushions.

You pick the dog. Not with certainty — because the evidence fits that story better. That reasoning, made into arithmetic, is maximum likelihood.

Why it exists

Data never announces its own settings. You see the wait times between support tickets; nobody tells you the true average rate. You see coin flips; nobody certifies the coin's bias.

Ronald Fisher, around 1922, proposed one universal recipe: try each candidate setting, score how well it explains the observed data, and keep the top scorer. That score is the likelihood — the probability of your actual data, viewed as a function of the unknown setting.

One recipe now fits nearly everything: coins, wait times, regression lines, neural networks.

How it works

candidate setting  →  "how probable was MY data under this setting?"

rate 0.05/min  →  data probability: tiny      (tickets too spread out)
rate 0.20/min  →  data probability: highest   ★ keep this one
rate 0.90/min  →  data probability: tiny      (tickets too bunched)

slide through the candidates, keep the peak

Two practical notes hide in that picture. Scoring a whole dataset means combining the probability of every single point. Computers need a careful way to do that, or the running total collapses to nothing. The developer tab shows the standard trick.

And the answer often matches common sense. For a coin, the maximum-likelihood bias is the fraction of heads you saw. The recipe earns its keep on problems where common sense has no formula ready.

A real example you have seen

Every time a spam filter is trained, maximum likelihood (or a close relative) is choosing its numbers. It picks the settings under which the pile of real emails and spam looks least surprising. Insurance companies fit accident-rate settings to claims data the same way. So did the model behind your phone's autocorrect.

Remember this

  • The likelihood scores a setting by how well it explains your actual data.
  • MLE keeps the setting at the peak of that score.
  • Common sense often agrees with it; its value is the problems where common sense has no answer ready.

What to learn next

Developer — Code and libraries.

Setup

bash
pip install numpy scipy

Outputs verified with numpy 1.26 and scipy 1.14, CPU.

Fitting a rate to waiting times

Minutes between successive support tickets. Model: waiting times follow an exponential distribution — the standard model for "events arriving at random". One unknown: the rate.

mle_rate.py
import numpy as np
from scipy import optimize

gaps = np.array([3.2, 11.5, 0.8, 6.1, 2.4, 9.7, 1.5, 4.8, 7.3, 2.9])

def neg_log_likelihood(rate):
    # log of each gap's probability density under this rate, summed.
    # Minimising the negative = maximising the likelihood.
    return -np.sum(np.log(rate) - rate * gaps)

result = optimize.minimize_scalar(neg_log_likelihood,
                                  bounds=(0.01, 5), method="bounded")
print(f"MLE rate: {result.x:.4f} tickets per minute")
print(f"closed-form answer, 1/mean: {1 / gaps.mean():.4f}")
Output
MLE rate: 0.1992 tickets per minute
closed-form answer, 1/mean: 0.1992

The walkthrough

np.log(rate) - rate * gaps is the log of the exponential density, written directly. For each gap, "how plausible is this wait under the candidate rate", logged. The sum scores the whole dataset. Logging first is the standard trick the beginner tab pointed at: multiplying ten densities gives a small number, and multiplying ten thousand underflows to exactly 0.0 in floating point. Logs turn that product into a sum, which stays representable, and the log is increasing, so the peak does not move.

Why the negative? Optimisers minimise by convention. Minimising negative log-likelihood is maximising likelihood. You have met this object before under an alias: cross-entropy loss is the negative log-likelihood of a classification model. Training a classifier is MLE at industrial scale.

The two printed lines agree, and that is the lesson. For the exponential, calculus solves the peak exactly: one over the sample mean. The optimiser found the same answer numerically, without the calculus. The numeric route is what generalises — most real models have no closed form, and then minimize with a gradient (or autograd) takes over.

Multiple parameters change the tool, not the idea. For a two-parameter fit, write neg_log_likelihood(params) unpacking both, and call optimize.minimize with a starting guess.

Common mistakes

Maximising raw likelihood instead of log-likelihood. Ten probabilities near 0.01 multiply to 1e-20; a thousand underflow to exactly 0.0, and the optimiser sees a flat landscape. Sum logs, never multiply probabilities.

Letting the optimiser wander into forbidden territory. A negative rate makes np.log(rate) blow up with RuntimeWarning: invalid value encountered in log. Fix: bounded methods as above, or optimise the log of the parameter so any real value is legal.

Trusting the MLE from tiny samples. MLE's guarantees are large-sample guarantees. Famous small-sample bug: the MLE of variance divides by n, understating spread — the reason the n−1 version exists. With 10 points, treat the fit as a rough sketch; bootstrap it for an honest wobble.

Comparing likelihoods across different models carelessly. More parameters never hurt the training likelihood — the same trap as overfitting. Penalise complexity (AIC/BIC) or validate on held-out data.

Try it yourself

Fit a normal distribution to gaps instead: two parameters, optimize.minimize, negative log density 0.5*np.log(2*np.pi*sig**2) + (x-mu)**2/(2*sig**2) summed. Check the fitted mu against gaps.mean(). Then ask yourself which model an ML engineer should prefer for waiting times, and why.

What to learn next

Researcher — Mathematics and papers.

Definition and mechanics

Given i.i.d. data $x_1, \dots, x_n$ from density $p(x \mid \theta)$:

$$ \hat\theta_{MLE} = \arg\max_\theta \; \ell(\theta), \qquad \ell(\theta) = \sum_{i=1}^{n} \log p(x_i \mid \theta) $$

Where:

  • $\theta$ — the parameter vector; $\ell$ — the log-likelihood.
  • The score is $s(\theta) = \nabla_\theta \ell$; the MLE solves $s(\hat\theta) = 0$ in interior, smooth cases.

The Fisher information measures the likelihood's curvature:

$$ I(\theta) = \mathbb{E}\left[ s(\theta) s(\theta)^\top \right] = -\mathbb{E}\left[ \nabla^2_\theta \log p(x \mid \theta) \right] $$

(the equality holding under regularity conditions).

Asymptotic theory

Under regularity (identifiability, smoothness, true $\theta_0$ interior):

  • Consistency: $\hat\theta \xrightarrow{p} \theta_0$.
  • Asymptotic normality: $\sqrt{n}(\hat\theta - \theta_0) \xrightarrow{d} \mathcal{N}!\left(0, I_1(\theta_0)^{-1}\right)$ with $I_1$ the per-observation information.
  • Efficiency: the asymptotic variance attains the Cramér–Rao bound — no consistent, asymptotically normal estimator does better.
  • Invariance: $\widehat{g(\theta)} = g(\hat\theta)$ for any function $g$ — a property Bayesian posteriors and unbiased estimators lack.

Standard errors follow from the inverse observed information $\left[-\nabla^2 \ell(\hat\theta)\right]^{-1}$; the three classical tests (Wald, score, likelihood-ratio) are asymptotically equivalent expansions of the same object, with Wilks' theorem giving $2[\ell(\hat\theta) - \ell(\theta_0)] \xrightarrow{d} \chi^2_k$.

MLE as KL projection

Maximising $\ell$ is minimising the empirical Kullback–Leibler divergence:

$$ \hat\theta = \arg\min_\theta \; \mathrm{KL}!\left(\hat{F}n \,|\, p\theta\right) $$

with $\hat F_n$ the empirical distribution. Under misspecification (no $\theta$ makes $p_\theta$ true), the MLE converges to the KL-closest model — the quasi-MLE — and honest standard errors need the sandwich covariance $A^{-1} B A^{-1}$, with $A$ the expected Hessian and $B$ the score's variance (White, 1982, Maximum likelihood estimation of misspecified models).

The bridge to ML losses

  • Cross-entropy loss = negative log-likelihood of a categorical model; logistic regression is Bernoulli MLE.
  • Squared error = Gaussian negative log-likelihood with fixed variance; absolute error = Laplace.
  • Language-model pretraining objectives are conditional MLE over token sequences.
  • Regularised losses correspond to maximum a posteriori (MAP) estimation — the next lessons' territory: L2 penalty = Gaussian prior, L1 = Laplace prior.
  • Latent-variable likelihoods (mixtures, HMMs) are typically non-convex and optimised by expectation–maximisation (Dempster, Laird and Rubin, 1977) or variational bounds.

Pathologies worth knowing

Unbounded likelihoods (Gaussian mixture with one component's variance $\to 0$); non-identifiable parameters (label switching); boundary parameters breaking the $\chi^2$ asymptotics (mixture order testing); and the Neyman–Scott problem — infinitely many nuisance parameters make MLE inconsistent, motivating profile and marginal likelihoods.

Foundational papers: Fisher (1922), On the mathematical foundations of theoretical statistics; Le Cam (1953) for the asymptotic optimality theory.

What to learn next