Maximum likelihood

The estimators note judged estimators once they were given. Maximum likelihood produces them: among all values of the parameter, choose the one under which the observed data are most probable. The recipe works for any model that assigns a probability or a density to the data, from a coin to a language model, and most common loss functions in machine learning are it in disguise. For large samples it is also about as good as an estimator can be. This note defines the likelihood, computes the estimates for the common families, derives squared error and cross-entropy from it, and explains its precision through the Fisher information.

Notation. $x_1, \dots, x_n$ are observed values, assumed i.i.d. from $p_\theta$ (a PMF or a density). $L(\theta) = \prod_i p_\theta(x_i)$ is the likelihood, $\ell(\theta) = \log L(\theta)$ the log-likelihood, $\hat\theta = \arg\max_\theta L(\theta)$ the maximum likelihood estimate (MLE).

The likelihood#

The likelihood is the probability of the observed data, read as a function of the parameter. After 7 heads in 10 independent flips, the probability of that particular sequence under heads-probability $p$ is

$$L(p) = p^7 (1 - p)^3 .$$

The data are fixed; $p$ varies. $L$ is not a probability distribution over $p$: it need not integrate to 1, and only ratios $L(p_1)/L(p_2)$ are meaningful, saying how much better one value explains the data than another. Turning it into a distribution over $p$ needs a prior, the subject of the Bayesian inference note.

Figure

The likelihood of the probability of heads

Maximizing it#

Products of many small numbers are awkward, and the logarithm turns them into sums without moving the maximum. For the coin, $\ell(p) = 7 \log p + 3 \log(1 - p)$, and setting the derivative to zero,

$$\ell'(p) = \frac{7}{p} - \frac{3}{1-p} = 0 \quad \Longrightarrow \quad \hat p = \frac{7}{10} .$$

The derivative of the log-likelihood is called the score. The same calculation for the common families gives estimates that look obvious in hindsight:

Model	MLE	Remark
$\mathrm{Bernoulli}(p)$	$\hat p = k/n$, the share of successes
$N(\mu, \sigma^2)$	$\hat\mu = \bar x$, $\hat\sigma^2 = \frac1n \sum (x_i - \bar x)^2$	$\hat\sigma^2$ is biased by $(n-1)/n$
$\mathrm{Poisson}(\lambda)$	$\hat\lambda = \bar x$
$\mathrm{Exp}(\lambda)$	$\hat\lambda = 1/\bar x$	biased upward: $E[\hat\lambda] = \frac{n}{n-1}\lambda$ for $n > 1$, infinite for $n = 1$
$U(0, \theta)$	$\hat\theta = \max_i x_i$	the maximum is on the boundary, where $\ell' \ne 0$

The last row is the estimator of the uniform end point from the previous note: maximum likelihood does not promise unbiasedness. It does promise invariance: the MLE of any function $g(\theta)$ is $g(\hat\theta)$, so the MLE of the standard deviation is $\sqrt{\hat\sigma^2}$ and the MLE of the odds is $\hat p / (1 - \hat p)$.

Loss functions are likelihoods#

A regression model $f_w(x)$ with Gaussian noise says $y = f_w(x) + \varepsilon$, $\varepsilon \sim N(0, \sigma^2)$. The negative log-likelihood of the training set is

$$-\ell(w) = \sum_{i=1}^n \frac{(y_i - f_w(x_i))^2}{2\sigma^2} + n \log\big(\sigma\sqrt{2\pi}\big) ,$$

and minimizing it over $w$ is minimizing the sum of squared errors. A classifier that outputs probabilities $q_w(k \mid x)$ for the classes assigns the observed label $y_i$ the probability $q_w(y_i \mid x_i)$, so

$$-\ell(w) = -\sum_{i=1}^n \log q_w(y_i \mid x_i) ,$$

the cross-entropy loss. A language model trained on text is the same formula over tokens, with the chain rule of the conditional probability note splitting each sentence into next-token predictions.

Figure

Loss functions are negative log-likelihoods

The view explains choices that otherwise look arbitrary. Absolute error, the negative log-likelihood of Laplace noise, is robust because the Laplace distribution expects occasional large residuals; Poisson regression for counts uses the Poisson likelihood; the noise variance $\sigma^2$, irrelevant for the fitted $w$, matters as soon as the model must report its own uncertainty. And over a large sample, the average negative log-likelihood approaches $E_{p}[-\log p_\theta(X)]$, the cross-entropy between the data distribution and the model, so maximizing likelihood minimizes the KL divergence from the data to the model, the view taken by the generative models course.

Fisher information#

How precise is $\hat\theta$? The figure suggests the answer: a sharp peak of the log-likelihood pins the parameter down, a flat one does not. The expected curvature per observation is the Fisher information

$$I(\theta) = -E_\theta\big[\tfrac{\partial^2}{\partial\theta^2} \log p_\theta(X)\big] = \operatorname{Var}_\theta\big(\tfrac{\partial}{\partial\theta} \log p_\theta(X)\big) ,$$

the variance of the score. For a Bernoulli observation $I(p) = 1/(p(1-p))$. The curvature that the figure shows is a different, data-dependent quantity, the observed information $J_n(\hat\theta) = -\ell_n''(\hat\theta)$ of the whole sample; at the true $\theta$ its expectation is $n I(\theta)$, and under the conditions below the two agree for large $n$. For the coin they coincide at the estimate: $J_n(\hat p) = n/(\hat p(1 - \hat p)) = n I(\hat p)$. Under regularity conditions (among them: the model is correctly specified and identifiable, so that different $\theta$ give different distributions, the support does not depend on $\theta$, the log-likelihood is smooth, the true $\theta$ is not on the boundary), two results make the MLE the default estimator:

• Asymptotic normality. $\hat\theta \approx N\big(\theta, \frac{1}{n I(\theta)}\big)$ for large $n$. For the coin, the standard error is $\sqrt{p(1-p)/n}$: $0.145$ after 10 flips and $0.046$ after 100, and a simulation of $\hat p$ from 100 flips gives a spread of $0.0458$.
• Efficiency. The Cramér–Rao bound says that no unbiased estimator has variance below $1/(n I(\theta))$. The MLE reaches it asymptotically, so for large samples no reasonable estimator is more precise.

With $d$ parameters, $I(\theta)$ is a $d \times d$ matrix, and its inverse, divided by the sample size, is the asymptotic covariance of the estimate. Its role in machine learning goes beyond error bars: natural gradient methods precondition by it, and continual-learning methods such as elastic weight consolidation use it to measure which weights matter.

When maximum likelihood goes wrong#

• Small samples. After 3 heads in 3 flips, $\hat p = 1$: tails are declared impossible. A word never seen in training gets probability zero in an n-gram model, which is why such models are smoothed, as in the n-gram language models note. Priors fix this systematically.
• Unbounded likelihoods. In a mixture of two Gaussians, a component centred on one data point with variance shrinking to zero makes the likelihood grow without limit: for five points the log-likelihood is $-7.64$, $-3.06$ and $3.84$ for component widths $10^{-1}$, $10^{-3}$, $10^{-6}$. The maximum does not exist, and fitting procedures must bound the variances from below.
• Too many parameters. When the number of parameters grows with the data, the MLE can be inconsistent. With many pairs of measurements, each pair with its own mean and a common variance, the MLE of the variance converges to half its true value (the Neyman–Scott problem); overparametrized models without regularization fit noise in the same way.
• The wrong model. If no $p_\theta$ is the true distribution, the MLE converges to the member closest in KL divergence, which can be a poor description of what matters, for instance of the tails.

What changed#

	Before	With maximum likelihood
Finding an estimator	a new idea per problem	one recipe: maximize $\ell(\theta) = \sum \log p_\theta(x_i)$
Loss functions	chosen by habit	negative log-likelihoods of a noise model: MSE, MAE, cross-entropy
Precision	estimated case by case	$\hat\theta \approx N(\theta, 1/(nI(\theta)))$; Cramér–Rao bound
Failure modes		small samples, unbounded likelihoods, too many parameters, wrong models

Limitations#

• The model must be specified. Maximum likelihood finds the best parameter within a family; choosing the family, and checking it, is outside the method.
• Asymptotic guarantees need regularity and large $n$. With few observations, boundary estimates, or parameters that change the support (the uniform end point), the normal approximation and the Cramér–Rao bound do not apply.
• A point estimate hides uncertainty. The MLE alone does not say how confident to be; the Fisher information, the likelihood curve, a confidence interval or a posterior distribution are needed.

What you should remember#

• The likelihood $L(\theta) = \prod_i p_\theta(x_i)$ is the probability of the data as a function of the parameter; only ratios matter.
• The MLE maximizes $\ell(\theta) = \sum_i \log p_\theta(x_i)$; for the common families it is the share, the mean, the $1/n$ variance, $1/\bar x$, the maximum.
• Squared error, absolute error and cross-entropy are negative log-likelihoods of Gaussian, Laplace and categorical models; maximum likelihood minimizes KL from the data to the model.
• Fisher information $I(\theta) = \operatorname{Var}(\text{score})$; $\hat\theta \approx N(\theta, 1/(nI(\theta)))$; no unbiased estimator beats $1/(nI(\theta))$.
• Small samples, unbounded likelihoods, growing parameter counts and misspecified models break it.

Code#

The coin’s likelihood, least squares as maximum likelihood, small-sample bias and an unbounded likelihood:

import numpy as np, math
from scipy import stats, optimize

rng = np.random.default_rng(0)

# 7 heads in 10 and 70 in 100: the MLE, the likelihood ratio of p = 0.5, the standard error from Fisher information
for k, n in ((7, 10), (70, 100)):
    loglik = lambda p: k * math.log(p) + (n - k) * math.log(1 - p)
    p_hat = optimize.minimize_scalar(lambda p: -loglik(p), bounds=(1e-6, 1 - 1e-6), method='bounded').x
    print(n, round(p_hat, 3), f'{math.exp(loglik(0.5) - loglik(0.7)):.2g}', round(math.sqrt(0.7 * 0.3 / n), 3))
# 10 0.7 0.44 0.145 | 100 0.7 0.00027 0.046
print(round(float((rng.binomial(100, 0.7, 200_000) / 100).std()), 4))   # 0.0458: the simulated spread of p-hat

# least squares is maximum likelihood under Gaussian noise: both fits of a line give the same coefficients
x = rng.uniform(0, 1, 50)
y = 2 + 3 * x + 0.5 * rng.normal(size=50)
ls = np.polyfit(x, y, 1)
nll = lambda w: -stats.norm.logpdf(y, w[1] + w[0] * x, 0.5).sum()
ml = optimize.minimize(nll, [0.0, 0.0]).x
print(ls.round(4), ml.round(4))                                   # [2.7233 2.1203] [2.7233 2.1203]

# small-sample bias of MLEs, n = 5: the normal variance (true 1) and the exponential rate (true 1)
print(round(float(rng.normal(size=(200_000, 5)).var(axis=1).mean()), 3), round(float((1 / rng.exponential(1, (200_000, 5)).mean(axis=1)).mean()), 3))   # 0.798 1.249

# an unbounded likelihood: a two-component Gaussian mixture with one component shrinking onto a data point
data = np.array([-1.2, 0.3, 0.5, 1.1, 2.4])
for s in (1e-1, 1e-3, 1e-6):
    print(s, round(float(np.log(0.5 * stats.norm.pdf(data, data[0], s) + 0.5 * stats.norm.pdf(data, data.mean(), data.std())).sum()), 2))
# 0.1 -7.64 | 0.001 -3.06 | 1e-06 3.84

Next#

Maximum likelihood returns one value and treats the parameter as fixed but unknown. After 3 heads in 3 flips it declares tails impossible, because it has no way to say that a fair-looking coin is a priori likely to be roughly fair. The next note treats the parameter itself as random, with a prior distribution updated by Bayes’ formula into a posterior, and shows that common regularizers are priors in disguise.

References#

• R. A. Fisher. On the Mathematical Foundations of Theoretical Statistics. Philosophical Transactions of the Royal Society A, 222, 1922.
• L. Wasserman. All of Statistics. Springer, 2004. Chapter 9 (maximum likelihood, Fisher information, asymptotic normality, efficiency).
• G. Casella, R. L. Berger. Statistical Inference, 2nd ed. Duxbury, 2002. Section 7.2 (maximum likelihood), Section 7.3 (the Cramér–Rao inequality), Chapter 10 (asymptotics).
• C. M. Bishop. Pattern Recognition and Machine Learning. Springer, 2006. Sections 1.2.4–1.2.5 (maximum likelihood for the Gaussian, curve fitting as maximum likelihood), Section 9.2.1 (singularities of the mixture likelihood).
• J. Neyman, E. L. Scott. Consistent Estimates Based on Partially Consistent Observations. Econometrica, 16(1), 1948.