Point Estimation and Maximum Likelihood

Point estimation is the silent engine under every quoted number in experimental physics. When the LIGO and Virgo collaborations report the masses and spins of a merging black-hole binary, they are performing maximum-likelihood (and Bayesian) parameter estimation against a Gaussian-noise model, and the Fisher information matrix forecasts how tightly each parameter can be measured before the full analysis is run. At the Large Hadron Collider the discovery significance of a new particle is a likelihood-ratio statistic, and its mass and cross-section are extracted by maximizing a likelihood built from Poisson-distributed bin counts — the very $\sqrt N$ rule of Worked Example 2, scaled to billions of events. Cosmologists fit the six parameters of the standard cosmological model to the microwave background by maximizing a likelihood over millions of pixels, and the Cramér–Rao bound, in its Fisher-matrix guise, is the standard tool for designing the next survey before it collects a single photon. The same asymptotic-normality result that makes the MLE Gaussian underlies the error ellipses in every fit, and the bias–variance trade-off is exactly the principle behind the regularized estimators that let machine-learning models fit millions of parameters from finite data. Wherever a physicist writes a central value with an error bar, Fisher's log-likelihood is doing the work.

Flip the likelihood on its head. Once the data are in hand they are fixed constants; what varies is the parameter, and $L(\theta)=\prod_i p(x_i\mid\theta)$ becomes a function over the space of possible truths. Each data point casts a vote for the $\theta$ values that could plausibly have produced it, and the votes multiply. The place where the product peaks is the value the data endorse most strongly — the maximum-likelihood estimate. What matters just as much as where the peak sits is how sharp it is: a needle-thin peak pins $\theta$ down tightly, a broad hump leaves many values nearly as plausible. That sharpness — the curvature of $\ln L$ at its summit — is precisely the Fisher information, and its reciprocal is the variance you should expect in $\hat\theta$.

Now separate the two ways an estimator can go wrong, because they are genuinely independent. A biased estimator is a rifle whose sights are misaligned: every shot lands consistently off the bullseye no matter how steady your hand. A high-variance estimator is a shaky hand: even with perfect sights the shots scatter widely. These failures do not trade against each other automatically — you can have either, both, or neither — and the mean-squared error simply adds their squares, $\mathrm{MSE}=\mathrm{Var}+b^2$. The useful fact is that buying a little bias can sometimes purchase a large reduction in variance — the controlled bargain behind shrinkage, regularization, and estimation with too little data.

The Cramér–Rao bound is the no-free-lunch theorem that closes the picture. It says the sharpness of the likelihood peak — set purely by how strongly the distribution's shape depends on $\theta$ — caps the precision any unbiased estimator can achieve, and no cleverer statistic or smarter algorithm can beat it. It is the statistical cousin of an uncertainty relation, and maximum likelihood is remarkable precisely because it saturates the bound as the data pile up.

The likelihood is the joint probability of the observed data read as a function of the parameter; independence turns the joint density into a product, and the logarithm turns the product into a sum that is far easier to differentiate and to reason about asymptotically. The log-likelihood $\ell$ is the object every later formula acts on.

where: $x_i$ — the $i$-th observed data point (independent draws); $p(x\mid\theta)$ — the model probability density or mass at $x$ for parameter $\theta$; $\theta$ — the unknown parameter(s) to be estimated; $N$ — sample size. Reading: multiplying probabilities for independent data becomes adding log-probabilities, so information from separate observations accumulates term by term. The maximum-likelihood estimate is the location of the peak. In the smooth (regular) case it is found by setting the first derivative — the score — to zero and checking that the curvature is negative.

where: $\hat\theta_{\mathrm{ML}}$ — the maximum-likelihood estimate, the stationary point that is a maximum. Reading: the recipe is universal but presumes a differentiable interior maximum — a caveat that bites when the support of $p$ depends on $\theta$. To judge an estimator we split its total error into a systematic part and a random part. The bias is the offset of its average from the truth; the mean-squared error decomposes cleanly, with no cross term, into variance plus squared bias.

where: $b(\theta)$ — bias, the systematic deviation of the estimator's mean from the truth; $\mathrm{Var}(\hat\theta)$ — sampling variance of the estimator; $\mathbb{E}_\theta$ — expectation over data drawn from $p(\cdot\mid\theta)$. Reading: because the two contributions add in quadrature, driving one to zero can raise the other, and minimizing their sum — not either alone — is the real goal. The Fisher information quantifies how much a single observation tells you about $\theta$. It is the variance of the score, and equivalently — under regularity — the expected curvature of the log-likelihood; the two forms are the workhorse identity of the whole theory.

where: $I(\theta)$ — Fisher information for one observation; for $N$ iid points it is additive, $I_N=N\,I_1$; $\partial_\theta\ln p$ — the score function, whose expectation is zero. Reading: more information means a sharper peak and, by the next equation, a lower attainable variance. The Cramér–Rao inequality converts curvature into a precision floor. For a possibly biased estimator the bound carries the factor $(1+b')^2$; for an unbiased one it collapses to the reciprocal of the total Fisher information.

where: $b'(\theta)$ — derivative of the bias with respect to $\theta$; $I_1(\theta)$ — single-observation Fisher information. Reading: the floor falls as $1/N$, so precision improves only as $\sqrt{N}$ — quadrupling the data halves the error, a law of diminishing returns baked into statistics itself. Applied to $N$ Gaussian samples, maximum likelihood returns exactly the intuitive estimators — but with a twist in the variance. The sample mean is unbiased and efficient; the ML variance divides by $N$ and is biased low, which is why the unbiased version divides by $N-1$.

where: $\bar x$ — sample mean, unbiased and efficient for $\mu$; $\hat\sigma^2_{\mathrm{ML}}$ — ML variance estimate, biased by $(N-1)/N$; the unbiased sample variance uses the divisor $N-1$. Reading: estimating $\mu$ from the same data consumes one degree of freedom, so the residuals scatter a touch less than the true errors. Maximum likelihood has a property no competitor shares so cleanly: it commutes with reparametrization. The MLE of any function of $\theta$ is that function of the MLE — you may estimate in whichever coordinate is convenient and transform afterward.

where: $g$ — any measurable reparametrization. Reading: estimate a rate and report a lifetime without redoing the optimization — though this invariance concerns the estimate, not its bias. Finally, maximum likelihood is the flat-prior limit of Bayesian estimation. Adding the log-prior to the log-likelihood and maximizing gives the maximum-a-posteriori estimate; when the prior is uniform the extra term is constant and MAP coincides with MLE.

where: $\pi(\theta)$ — prior density encoding pre-data belief. Reading: a prior acts as a regularizer that pulls the estimate toward favored values, and its influence fades as the likelihood sharpens with more data.

1 · Measuring g: a Gaussian mean and its honest uncertainty (fully worked)

Eight repeated timings of a pendulum yield the following values for the local gravitational acceleration (in m/s²): 9.81, 9.79, 9.85, 9.78, 9.83, 9.80, 9.84, 9.76. Treating them as iid Gaussian samples, find the ML estimate of $g$, the ML and unbiased variance estimates, and the standard error of the mean you would quote.

Step 1. The MLE of the mean is the sample average (equation 6). Summing gives $\sum_i x_i=78.46$, so $\hat g=\bar x=78.46/8=\mathbf{9.8075\ m/s^2}$. The sample mean is unbiased and efficient for a Gaussian mean, so no bias correction is needed here.

Step 2. Form the residual sum of squares: $\sum_i(x_i-\bar x)^2=6.75\times10^{-3}\ \mathrm{(m/s^2)^2}$. The ML variance divides by $N$: $\hat\sigma^2_{\mathrm{ML}}=6.75\times10^{-3}/8=8.44\times10^{-4}$, i.e. $\hat\sigma_{\mathrm{ML}}=0.0290\ \mathrm{m/s^2}$ — what the likelihood equation literally returns.

Step 3. Apply Bessel's correction for an unbiased spread: $s^2=6.75\times10^{-3}/7=9.64\times10^{-4}$, so $s=\mathbf{0.0311\ m/s^2}$. The ratio $\hat\sigma^2_{\mathrm{ML}}/s^2=(N-1)/N=7/8=0.875$ is exactly the bias factor from Derivation 1 — the ML variance is 12.5% low.

Step 4. The uncertainty on the mean is the standard error, $\mathrm{SEM}=s/\sqrt N=0.0311/\sqrt 8=\mathbf{0.0110\ m/s^2}$, shrinking as $1/\sqrt N$ just as Cramér–Rao demands (for a Gaussian mean $I_1=1/\sigma^2$, so the floor is $\sigma^2/N$ and $\bar x$ attains it). Report g = 9.808 ± 0.011 m/s²; quoting the ML variance for $s^2$ would understate the spread by 0.875.

2 · A counting experiment: Poisson rate and the √N rule (fully worked)

A particle detector is run four times under identical conditions and records 18, 25, 21, and 24 counts. Assuming Poisson statistics, estimate the mean rate per run, its uncertainty, and verify that the estimator saturates the Cramér–Rao bound.

Step 1. The Poisson MLE is the mean count (Derivation logic: $\partial_\lambda\ell=-N+\lambda^{-1}\sum_i n_i=0$). With $\sum_i n_i=88$ over $N=4$ runs, $\hat\lambda=\bar n=88/4=\mathbf{22.0\ counts/run}$. This is unbiased because $\mathbb{E}[n_i]=\lambda$.

Step 2. The single-count Fisher information is $I_1(\lambda)=1/\lambda$ (from $\partial_\lambda^2\ln p=-n/\lambda^2$ and $\mathbb{E}[n]=\lambda$). The Cramér–Rao floor is therefore $1/(N I_1)=\lambda/N$. Evaluated at the estimate, $\mathrm{Var}(\hat\lambda)=22.0/4=5.50$, so the standard error is $\sqrt{5.50}=\mathbf{2.35\ counts/run}$.

Step 3. Check saturation directly. For a Poisson variable $\mathrm{Var}(n)=\lambda$, so the sampling variance of the mean is $\mathrm{Var}(\bar n)=\mathrm{Var}(n)/N=\lambda/N$, identical to the Cramér–Rao floor — the Poisson rate estimator is efficient for every $N$, not merely asymptotically, because the score is linear in $\hat\lambda$.

Step 4. Connect to the everyday "$\sqrt N$ rule". The total count is $88$, whose statistical uncertainty is $\sqrt{88}=9.4$, a relative error of $9.4/88=\mathbf{10.7\%}$; this is why a histogram bin holding $N$ events carries an error bar of $\sqrt N$. Report 22.0 ± 2.3 counts/run, with the fractional precision improving only as $1/\sqrt{\text{total counts}}$.

3 · Reparametrization invariance: from decay rate to lifetime

Decay times $t_1,\dots,t_N$ are drawn from $p(t\mid\tau)=\tau^{-1}e^{-t/\tau}$. Show the MLE of the mean lifetime is the sample mean, and use invariance to write down the MLE of the decay rate $\gamma=1/\tau$ without any further calculus.

Step 1. The log-likelihood is $\ell(\tau)=\sum_i(-\ln\tau-t_i/\tau)=-N\ln\tau-\tau^{-1}\sum_i t_i$. Setting $\partial_\tau\ell=-N/\tau+\tau^{-2}\sum_i t_i=0$ gives $\hat\tau=\tfrac1N\sum_i t_i=\bar t$, the average observed decay time.

Step 2. The Fisher information for the exponential is $I_1(\tau)=1/\tau^2$, so the asymptotic variance is $\tau^2/N$ and the fractional error on a lifetime measured from $N$ decays is $1/\sqrt N$ — a thousand decays give a 3% lifetime.

Step 3. By functional invariance (equation 7), the MLE of the rate is simply the reciprocal of the MLE of the lifetime: $\hat\gamma=1/\hat\tau=N/\sum_i t_i$. No re-optimization is needed — estimate in whichever variable is convenient and transform, remembering only that unbiasedness does not survive the nonlinear map.

Maximum likelihood, Fisher information, efficiency, sufficiency — nearly the entire vocabulary of point estimation was minted by a single person in the span of a few years, and then sharpened into a hard inequality by two more a generation later.

Texts. Stone & Goldbart, Mathematics for Physics (Cambridge, 2009), Ch. 11 (Probability and Statistics) · Arfken, Weber & Harris, Mathematical Methods for Physicists, 7th ed., chapter on Probability and Statistics · Riley, Hobson & Bence, Mathematical Methods for Physics and Engineering, 3rd ed., Ch. 31 (Statistics) · Cowan, Statistical Data Analysis (Oxford, 1998), Ch. 5–6 · Kendall & Stuart, The Advanced Theory of Statistics, Vol. 2. Key results to keep in mind: the Gaussian MLE and the $(N-1)/N$ bias follow from the identity $\sum(x_i-\bar x)^2=\sum(x_i-\mu)^2-N(\bar x-\mu)^2$; the Fisher-information equivalence $\mathbb{E}[(\partial_\theta\ln p)^2]=-\mathbb{E}[\partial_\theta^2\ln p]$ follows from differentiating the normalization twice; the Cramér–Rao inequality follows from Cauchy–Schwarz; and the asymptotic-normality result $\hat\theta_{\mathrm{ML}}\to\mathcal{N}(\theta_0,1/(N I_1))$ follows from Taylor-expanding the score. The numerical anchor values are: Gaussian $\bar x=9.8075$, $s=0.0311$, $\mathrm{SEM}=0.0110$, bias ratio $0.875$ at $N=8$; Poisson $\hat\lambda=22.0$, $\mathrm{SE}=2.35$, saturating $\lambda/N=1/(N I_1)$.

Prerequisites. Random Variables and Probability Distributions · Expectation, Variance and Moments · Conditional Probability and Bayes' Theorem · Central Limit Theorem · Gaussian (Normal) Distribution. Leads to. Multivariate Gaussian Distribution · Moment-Generating and Characteristic Functions · Least-Squares Fitting and Regression · Bayesian Inference and Posterior Estimation.

Ask the tutor: Ask me to walk through why the Gaussian variance picks up its $(N-1)/N$ bias, or to Fisher-forecast the error bar on a parameter before you take any data.