Maximum Likelihood Estimation

Maximum Likelihood Estimation (MLE) is a method of estimating the parameters of a statistical model by choosing the values that make the observed data most probable. Key result: the estimate is the argument that maximizes the likelihood function,

\hat{θ} = ar g max_{θ} L (θ ∣ x)

. Applies to: parametric inference across statistics, econometrics, machine learning, physics, and genetics. Who/when: systematized by Ronald A. Fisher in 1912–1922, with roots in Gauss and Laplace.

Maximum Likelihood Estimation is a general-purpose procedure for fitting a probability model with unknown parameters to observed data. Given a model that specifies the probability (or density) of the data as a function of parameters $θ$ , one views this probability — now read as a function of $θ$ with the data held fixed — as the likelihood, and selects the $θ$ that maximizes it. It matters because it yields estimators with deep optimality properties (consistency, asymptotic efficiency, and asymptotic normality) that hold under broad conditions, and it unifies a vast range of statistical techniques under a single principle. A surprising fact: many estimators you already know — the sample mean, the ordinary least squares fit, the sample proportion — are all maximum likelihood estimators for appropriate models, even though they were discovered independently and earlier.

The Core Insight

Probability and likelihood use the same formula but read it in opposite directions. Probability fixes the parameter and asks how likely various data are; likelihood fixes the observed data and asks how plausible various parameter values are. MLE says: believe the parameter value under which what you actually saw was least surprising. Because products of probabilities are awkward, we almost always maximize the logarithm of the likelihood — turning products into sums without moving the location of the maximum.

Formal Statement

Likelihood and Log-Likelihood

Let $X_{1}, \dots, X_{n}$ be data with joint probability mass or density function $f (x_{1}, \dots, x_{n} ∣ θ)$ depending on a parameter $θ \in Θ$ (possibly a vector). The likelihood function is the joint density read as a function of $θ$ :

L (θ ∣ x) = f (x_{1}, \dots, x_{n} ∣ θ) .

⟳simply relabel the joint density, holding the observed

x

fixed and letting

θ

vary over the parameter space

Θ

If the observations are independent and identically distributed (i.i.d.) with common density $f (x ∣ θ)$ , the likelihood factorizes:

L (θ ∣ x) = i = 1 \prod n f (x_{i} ∣ θ) .

⟳independence makes the joint density the product of marginals; "identically distributed" makes each factor use the same

f

The log-likelihood is

ℓ (θ) = lo g L (θ ∣ x) = i = 1 \sum n lo g f (x_{i} ∣ θ) .

⟳take the natural logarithm of [2]; the logarithm of a product is the sum of logarithms.

Symbols:

θ

= unknown parameter(s);

Θ

= parameter space;

x_{i}

= the

i

-th observation;

n

= sample size;

f

= the model's density/mass function;

L

= likelihood;

ℓ

= log-likelihood.

Maximum Likelihood Estimator

The maximum likelihood estimator $\hat{θ}$ is any parameter value maximizing the likelihood (equivalently the log-likelihood):

\hat{θ} = ar g θ \in Θ max L (θ ∣ x) = ar g θ \in Θ max ℓ (θ) .

⟳

lo g

is strictly increasing, so it preserves the location of the maximizer — maximizing

L

and maximizing

ℓ

give the same

\hat{θ}

When $ℓ$ is differentiable and the maximum is interior to $Θ$ , $\hat{θ}$ satisfies the likelihood (score) equation:

\frac{\partial ℓ ( θ )}{\partial θ}_{θ = \hat{θ}} = 0, with \frac{\partial ^{2} ℓ}{\partial θ ^{2}}_{\hat{θ}} < 0.

⟳a smooth interior maximum requires the first derivative to vanish (first-order condition) and the second derivative to be negative (second-order condition for a maximum, not a minimum or inflection).

Asymptotic Properties of the MLE

Under regularity conditions (the true $θ_{0}$ lies in the interior of $Θ$ , the model is identifiable, the support of $f$ does not depend on $θ$ , and $ℓ$ is sufficiently smooth), the MLE $\hat{θ}_{n}$ from $n$ i.i.d. observations satisfies:

Consistency: $\hat{θ}_{n} p θ_{0}$ as $n \to \infty$ .
Asymptotic normality:

n (\hat{θ}_{n} - θ_{0}) d N (0, I_{1} (θ_{0})^{- 1}),

⟳Taylor-expand the score

ℓ^{'} (\hat{θ}) = 0

about

θ_{0}

, apply the Central Limit Theorem to the score and the Law of Large Numbers to the (negative) curvature, then divide.

where $I_{1} (θ_{0})$ is the Fisher information per observation,

I_{1} (θ) = E [(\frac{\partial}{\partial θ} lo g f (X ∣ θ))^{2}] = - E [\frac{\partial ^{2}}{\partial θ ^{2}} lo g f (X ∣ θ)] .

⟳differentiate the identity

\int f d x = 1

twice under the integral sign; the score has mean zero, and its variance equals minus the expected curvature.

Asymptotic efficiency: $\hat{θ}_{n}$ attains the Cramér–Rao lower bound asymptotically — no consistent, asymptotically unbiased estimator has smaller asymptotic variance.

Why It Holds / Derivation

Sketch of Consistency and Asymptotic Normality

Step 1 — Why the maximizer lands on the truth. Define the normalized log-likelihood

\frac{1}{n} ℓ (θ) = \frac{1}{n} \sum_{i} lo g f (X_{i} ∣ θ)

. By the Law of Large Numbers this converges to

M (θ) = E_{θ_{0}} [lo g f (X ∣ θ)]

. Step 2 — The truth maximizes the limit. Consider

M (θ_{0}) - M (θ) = E_{θ_{0}} [lo g \frac{f ( X ∣ θ _{0} )}{f ( X ∣ θ )}]

, the Kullback–Leibler divergence, which is

\geq 0

and equals

0

only when

f (\cdot ∣ θ) = f (\cdot ∣ θ_{0})

. Hence

θ_{0}

uniquely maximizes

M

KL (θ_{0} ∥ θ) = E_{θ_{0}} [lo g \frac{f ( X ∣ θ _{0} )}{f ( X ∣ θ )}] \geq 0.

⟳apply Jensen's inequality to the concave

lo g

, using

E_{θ_{0}} [f (X ∣ θ) / f (X ∣ θ_{0})] = 1

Since the empirical objective converges to one maximized at $θ_{0}$ , its maximizer $\hat{θ}_{n}$ converges to $θ_{0}$ .

Step 3 — Expand the score equation. Let

s (θ) = ℓ^{'} (θ)

be the score. By definition

s (\hat{θ}_{n}) = 0

. Taylor-expand about

θ_{0}

0 = s (\hat{θ}_{n}) = s (θ_{0}) + s^{'} (θ_{0}) (\hat{θ}_{n} - θ_{0}) + remainder .

⟳first-order Taylor expansion of the score function around the true parameter; the remainder vanishes faster under regularity.

Step 4 — Rearrange and rescale. Solving and multiplying by

n

n (\hat{θ}_{n} - θ_{0}) = \frac{\frac{1}{n} s ( θ _{0} )}{- \frac{1}{n} s ^{'} ( θ _{0} )} .

⟳divide [9] by

s^{'} (θ_{0})

, isolate

(\hat{θ}_{n} - θ_{0})

, and split the

n

factors so numerator and denominator each have a standard limit.

Step 5 — Apply limit theorems. The numerator

\frac{1}{n} s (θ_{0})

is a sum of i.i.d. mean-zero score terms, so by the CLT it converges to

N (0, I_{1} (θ_{0}))

. The denominator

- \frac{1}{n} s^{'} (θ_{0})

converges by the LLN to

I_{1} (θ_{0})

. By Slutsky's theorem the ratio converges to

N (0, I_{1} (θ_{0}) / I_{1} (θ_{0})^{2}) = N (0, I_{1} (θ_{0})^{- 1})

, establishing [6].

■

Ronald Aylmer Fisher | 1890–1962

The British statistician and geneticist Ronald Fisher introduced the term "likelihood" as distinct from probability and developed maximum likelihood into a coherent theory between 1912 and 1922. His 1922 paper "On the Mathematical Foundations of Theoretical Statistics" defined consistency, efficiency, and sufficiency, and introduced what we now call Fisher information. Fisher argued that likelihood, not posterior probability, was the natural vehicle for inference when no prior is available — a stance that anchored the frequentist tradition. He also founded modern experimental design and population genetics, making him arguably the single most influential figure in 20th-century statistics.

Building Intuition

Imagine you find a coin on the ground and flip it 10 times, getting 7 heads. Ask: "For what bias $p$ would this exact outcome be least surprising?" If $p = 0.1$ , seven heads is wildly unlikely; if $p = 0.9$ , it's also somewhat unlikely; somewhere in between the probability of "7 of 10 heads" peaks. MLE finds that peak — here at $p = 0.7$ , exactly the observed fraction. The likelihood is not the probability that $p = 0.7$ (that question requires a prior); it is a ranking of parameter values by how well each explains the data you hold in your hand. The log turns the unwieldy product of ten coin-flip probabilities into a clean sum, and setting its derivative to zero locates the peak.

Likelihood of a coin's bias given 7 heads in 10 flips, peaking at p = 0.7

Worked Examples

Bernoulli / Binomial Proportion

Given:

n = 10

independent coin flips yield

k = 7

heads. Model each flip as

Bernoulli (p)

. Step 1 — Likelihood.

L (p) = (7 10) p^{7} (1 - p)^{3}

. Step 2 — Log-likelihood (drop the constant binomial coefficient):

ℓ (p) = 7 lo g p + 3 lo g (1 - p) + const .

⟳take

lo g

of the Bernoulli likelihood

p^{k} (1 - p)^{n - k}

with

k = 7

n - k = 3

; the binomial coefficient is constant in

p

Step 3 — Differentiate and set to zero:

\frac{d ℓ}{d p} = \frac{7}{p} - \frac{3}{1 - p} = 0 ⟹ 7 (1 - p) = 3 p ⟹ p = \frac{7}{10} .

⟳differentiate [11], clear denominators by multiplying through by

p (1 - p)

, and solve the linear equation.

Step 4 — Confirm a maximum:

ℓ^{''} (p) = - 7/ p^{2} - 3/ (1 - p)^{2} < 0

everywhere, so

p = 0.7

is the maximizer.

\overset{p}{^} = \frac{k}{n} = 0.70 (dimensionless proportion)

Normal Distribution with Unknown Mean and Variance

Given: Five measurements (in cm) of a rod:

x = {10.1, 9.8, 10.0, 10.3, 9.8}

, modeled as i.i.d.

N (μ, σ^{2})

. Step 1 — Log-likelihood:

ℓ (μ, σ^{2}) = - \frac{n}{2} lo g (2 π) - \frac{n}{2} lo g σ^{2} - \frac{1}{2 σ ^{2}} i = 1 \sum n (x_{i} - μ)^{2} .

⟳take

lo g

of the product of

n

Gaussian densities

\frac{1}{2 π σ ^{2}} e^{- (x_{i} - μ)^{2} /2 σ^{2}}

and collect terms.

Step 2 — Solve for $μ$ :

\partial ℓ / \partial μ = \frac{1}{σ ^{2}} \sum (x_{i} - μ) = 0

gives

\overset{μ}{^} = \overset{x}{ˉ} = \frac{1}{n} i = 1 \sum n x_{i} .

⟳the only

μ

-dependent term is the sum of squares; its derivative vanishes when

μ

equals the sample mean.

Numerically $\overset{μ}{^} = (10.1 + 9.8 + 10.0 + 10.3 + 9.8) /5 = 50.0/5 = 10.00$ cm.

Step 3 — Solve for $σ^{2}$ :

\partial ℓ / \partial σ^{2} = - \frac{n}{2 σ ^{2}} + \frac{1}{2 σ ^{4}} \sum (x_{i} - \overset{μ}{^})^{2} = 0

gives

\overset{σ}{^}^{2} = \frac{1}{n} i = 1 \sum n (x_{i} - \overset{μ}{^})^{2} .

⟳treat

σ^{2}

as a single variable, set the partial derivative to zero, and solve; note the divisor is

n

, not

n - 1

Numerically the squared deviations from 10.00 are $0.01, 0.04, 0.00, 0.09, 0.04$ , summing to $0.18$ , so $\overset{σ}{^}^{2} = 0.18/5 = 0.036 cm^{2}$ .

\overset{μ}{^} = 10.00 cm, \overset{σ}{^}^{2} = 0.036 cm^{2} (\overset{σ}{^} \approx 0.190 cm)

Common Misconceptions

Likelihood is not probability of the parameter.

L (θ ∣ x)

is not the probability that

θ

takes a given value. Within the frequentist framework

θ

is fixed (not random), so it has no probability distribution. The likelihood is the data's density viewed as a function of

θ

; it need not integrate to 1 over

θ

, and a "likelihood of 0.3" carries no probabilistic meaning on its own — only ratios of likelihoods are interpretable. The MLE can be biased. Consistency does not imply unbiasedness. The MLE of normal variance,

\overset{σ}{^}^{2} = \frac{1}{n} \sum (x_{i} - \overset{x}{ˉ})^{2}

, has expectation

\frac{n - 1}{n} σ^{2}

— it systematically underestimates. The familiar unbiased estimator uses

n - 1

in the denominator, but that is not the MLE. MLE optimizes likelihood, not unbiasedness, though the bias vanishes as

n \to \infty

When the MLE Behaves Badly

Consider $X_{1}, \dots, X_{n}$ i.i.d. from the uniform $U (0, θ)$ . The density is $1/ θ$ for $0 \leq x \leq θ$ and $0$ otherwise, so the support depends on $θ$ — violating a key regularity condition. The likelihood $L (θ) = θ^{- n} 1 {θ \geq max_{i} x_{i}}$ is not differentiable at its maximizer; setting a derivative to zero fails entirely. Here the MLE is $\hat{θ} = max_{i} x_{i}$ , found by direct reasoning (smaller $θ$ gives zero likelihood; among feasible $θ$ , $θ^{- n}$ is largest at the smallest feasible value). This estimator is biased and not asymptotically normal — it converges at rate $n$ , not $n$ , and its limiting law is exponential, not Gaussian. The lesson: the smooth asymptotic theory of Theorem [6] assumes a fixed support.

Connections

MLE sits at the center of a web of statistical ideas. Under a flat prior, the MLE coincides with the Maximum A Posteriori Estimation mode, linking it to Bayesian Inference. Minimizing the Kullback-Leibler Divergence between the empirical distribution and the model is equivalent to maximizing likelihood, tying MLE to Information Theory. For the Gaussian linear model, MLE reduces exactly to Ordinary Least Squares. The optimality bound it attains is the Cramér-Rao Lower Bound, expressed through Fisher Information. Likelihood ratios drive the Neyman-Pearson Lemma and the Likelihood Ratio Test, and the broader generalization to non-likelihood objectives is M-Estimation.

Real-World Applications

Logistic regression (epidemiology, credit scoring): parameters fit by maximizing the Bernoulli likelihood of binary outcomes — there is no closed form, so iterative Newton–Raphson (IRLS) is used.
Phylogenetics (RAxML, PhyML software): evolutionary trees and branch lengths inferred by maximizing the likelihood of observed DNA sequences under a substitution model.
GPS positioning: receiver location estimated by maximizing the likelihood of measured satellite signal delays given Gaussian timing noise.
Hidden Markov Models (speech recognition, bioinformatics): the Baum–Welch algorithm is an Expectation–Maximization procedure that ascends the likelihood of observed sequences.
Econometrics (discrete choice, ARCH/GARCH volatility models): Nobel-recognized work by Daniel McFadden and Robert Engle rests on maximum likelihood fitting.
Image reconstruction (PET and SPECT medical scanners): the ML-EM algorithm reconstructs tomographic images by maximizing the Poisson likelihood of detected photon counts.
Machine learning: training a neural network classifier with cross-entropy loss is exactly maximizing the likelihood of the labels under the model's predicted class probabilities.

Historical Development

1755–1760 — Thomas Bayes & Roger Cotes: early reasoning about inverting probability to infer causes from effects foreshadows likelihood thinking.
1809 — Carl Friedrich Gauss: in Theoria Motus, derives the method of least squares by assuming the arithmetic mean is the most probable value, implicitly using a likelihood argument for the normal error law.
1812 — Pierre-Simon Laplace: develops large-sample (asymptotic) approximations and the central role of the normal distribution in estimation.
1821 — Gauss: the Gauss–Markov framework refines least-squares estimation, later seen as Gaussian MLE.
1912 — Ronald Fisher: as an undergraduate, proposes maximizing the "absolute probability" of the data — his first formulation of the method.
1922 — Fisher: "On the Mathematical Foundations of Theoretical Statistics" coins likelihood, consistency, efficiency, sufficiency, and information; establishes the modern theory.
1925 — Fisher: Statistical Methods for Research Workers popularizes the method among practicing scientists.
1940s — Harald Cramér & C. R. Rao: prove the Cramér–Rao bound, rigorously establishing the asymptotic efficiency Fisher had asserted.
1977 — Dempster, Laird & Rubin: formalize the Expectation–Maximization (EM) algorithm, making MLE tractable for incomplete-data and latent-variable models.

Subtleties & Limits

Several deep issues qualify the clean asymptotic theory. Non-existence and non-uniqueness: likelihoods can be unbounded (e.g., a Gaussian mixture where one component's variance shrinks to zero, sending the likelihood to infinity), so a global maximizer may not exist; one must seek a sensible local maximum instead. High-dimensional regimes: when the number of parameters grows with $n$ (the Neyman–Scott problem), the MLE can be inconsistent. Boundary and singular models: when $θ_{0}$ lies on the boundary of $Θ$ or the Fisher information is singular (common in mixtures and neural networks), the standard $n$ -normal limit fails and is replaced by more delicate theory (singular learning theory, Watanabe's work). Computation: finding the global maximum of a non-concave likelihood is generally intractable, so practitioners rely on EM, gradient ascent, or Newton methods that find only local optima — quantifying the gap between local and global solutions in models like deep networks remains an active research frontier.