Parameter estimation

Parameter estimation

• Set up the likelihood equation, solve the maximum likelihood (ML) estimate for a Gaussian mean and variance, and use the invariance property.
• Assess an estimator by its bias, minimum-variance unbiasedness (MVU), and consistency.
• Derive the three Bayes estimates — MMSE, MAVE, MAP — from three different cost functions, and recognise when the three coincide.
• Apply the Cramér–Rao inequality to find the variance lower bound, and determine when an estimator achieves efficiency.
• Extend to multiple parameters using the Fisher information matrix, and set up the best linear unbiased estimator (BLUE) and least-squares estimation (LSE) when the probability density is unavailable.

Use the buttons, the dots or the ← → keys to move between slides.

From "which hypothesis is true" to "what is the parameter"

In chapter 5, the receiver receives a noisy version of the signal and decides which hypothesis is true among M possible hypotheses. In this chapter, assume the receiver has correctly decided the hypothesis, but a parameter associated with the signal may be unknown. The goal is to estimate that parameter optimally, based on a finite number of signal samples. With K independent identically distributed samples Y_1,\ldots,Y_K whose density depends on an unknown parameter \theta, the estimator \hat\theta=g(Y_1,\ldots,Y_K) is a statistic (a function of the data) used to estimate \theta (Barkat, pp. 345 to 346).

A random or nonrandom parameter: the chapter's two branches

The parameter to be estimated may be random or nonrandom. Estimation of a random parameter is called Bayes' estimation (Section 6.5); estimation of a nonrandom parameter is called maximum likelihood estimation (MLE, Section 6.2) (Barkat, p. 345). These two branches run in parallel throughout the chapter: MLE needs no knowledge of a prior density for \theta, while Bayes' estimation exploits that density when available.

The likelihood function and the likelihood equation

With K independent identically distributed samples, the likelihood function is L(\theta)=f_{Y|\Theta}(y|\theta)=\prod_{k=1}^Kf_{Y_k|\Theta}(y_k|\theta). The value \hat\theta maximising L(\theta) is called the maximum likelihood estimate. Because the logarithm is monotonically increasing, maximising L(\theta) is equivalent to maximising \ln L(\theta), and a necessary (not sufficient) condition is to solve the likelihood equation \partial\ln f_{Y|\Theta}(y|\theta)/\partial\theta=0 (Barkat, p. 346). Measured: with K=8 Gaussian samples m=3, \sigma=2, scanning \theta on a grid and locating the likelihood peak numerically gives \hat m = 2.3815, matching the closed form \hat m=\frac1K\sum y_k = 2.3815.

Example 6.1(a): the ML estimate of a Gaussian mean

With Y_k=m+N_k, N_k independent zero-mean Gaussian of variance \sigma^2, and m unknown: \ln f_{Y|M}(y|m)=\text{const}-\sum_{k=1}^K(y_k-m)^2/2\sigma^2. Solving \partial\ln f/\partial m=0 gives \sum_{k=1}^K(y_k-m)/\sigma^2=0, i.e. \hat m_{ml}=\frac1K\sum_{k=1}^Ky_k — the sample mean (Barkat, pp. 347 to 348). Measured: for the specific sample data, \hat m_{ml} = 2.3815, and the numerical derivative of \ln L at that point is 0.00e+00 (approximately zero), confirming it is a stationary point.

Example 6.1(b): the ML estimate of the variance, using invariance

If m is known but \sigma^2 is unknown, maximising \ln L(\sigma^2) is equivalent to minimising g(\sigma^2)=K\ln\sigma+\sum_{k=1}^K(y_k-m)^2/2\sigma^2. Using the invariance property, it is easier to differentiate with respect to \sigma rather than \sigma^2: dg/d\sigma=0 gives \hat\sigma_{ml}=\sqrt{\frac1K\sum(y_k-m)^2}, and hence \hat\sigma^2_{ml}=\frac1K\sum_{k=1}^K(y_k-m)^2 (Barkat, p. 348). Measured: with m known exactly at its true value, \hat\sigma^2_{ml} = 6.3152, against the true \sigma^2 = 4.

Why this is least squares in disguise

Maximising the Gaussian likelihood function is equivalent to minimising the sum of squared deviations \sum_k(y_k-m)^2 — this is exactly the least-squares principle (met again explicitly in part 5 of this module, Section 6.9). For a Gaussian model, ML and least squares coincide; for a different model (e.g. Example 6.5 in part 3), the two principles can give different results. Measured: scanning \sum_k(y_k-m)^2 against m directly and locating the numerical minimum gives 2.3815, matching \hat m_{ml} = 2.3815.

Self-check, part 1

• What quantity does the ML estimate \hat\theta maximise?
• Is the likelihood equation \partial\ln f_{Y|\Theta}/\partial\theta=0 a necessary or a sufficient condition?
• What does the invariance property of ML say about g(\hat\theta)?

Hints: the likelihood function L(\theta) (or \ln L(\theta)); only necessary, not sufficient (there may be several stationary points); g(\hat\theta) is the ML estimate of g(\theta) when g is one-to-one.

When a composite hypothesis meets an unknown parameter

In Example 5.9 (chapter 5), a hypothesis-testing problem had one hypothesis containing an unknown parameter, making that hypothesis composite. If that parameter is nonrandom, a practical way to proceed is to replace the parameter by its ML estimate directly inside the likelihood ratio, instead of averaging over a prior density (as module 21 did when the parameter was random). The resulting test, using \hat\theta_{ml} in place of the unknown \theta in \Lambda(y), is called the generalized likelihood ratio test (GLRT) (Barkat, pp. 349 to 351).

GLRT: practical, but no longer strictly optimal

GLRT does not minimise the Bayes risk in a strict sense (there is no prior density to average over), and is not necessarily the uniformly most powerful test (UMP, met in module 21) — but it typically performs well in practice and is the most common approach when the composite parameter has no reasonable prior probability model. The core idea: estimate first, test second — use the best available estimate of the unknown parameter, then proceed with the test as though that parameter were known (Barkat, pp. 349 to 351).

Unbiased: the expected value equals the true parameter

Since \hat\theta is a random variable (a function of random data), it can take different values across trials — criteria are needed to assess its \"goodness.\" The first criterion: \hat\theta is unbiased if E[\hat\theta]=\theta for every \theta. More generally, write E[\hat\theta]=\theta+b(\theta); if b(\theta) does not depend on \theta (a known, constant bias), one can correct \hat\theta-b to get an unbiased estimate; if b(\theta) depends on \theta (an unknown bias), no correction is possible since \theta is itself the unknown (Barkat, pp. 353 to 354).

Unbiased is not enough: small variance is also needed

Being unbiased only guarantees the average of the estimate is correct; if the variance is large, the specific error in any one trial can still be sizable (the book's Fig. 6.4). The second criterion is small variance. An estimate \hat\theta is called minimum variance unbiased (MVU) if, among all unbiased estimates \theta' of \theta, \text{var}[\hat\theta]\le\text{var}[\theta'] for every \theta' — i.e. it has the smallest variance in the class of unbiased estimators (Barkat, p. 354).

Consistent: the estimate improves as more data arrives

An estimator \hat\theta (based on K samples) is called consistent if \lim_{K\to\infty}P(|\hat\theta-\theta|>\varepsilon)=0 for every \varepsilon>0 (Barkat, p. 354). Checking this definition directly is often difficult; a more convenient theorem: if \hat\theta is unbiased and \lim_{K\to\infty}\text{var}[\hat\theta]=0, then \hat\theta is consistent (Barkat, pp. 354 to 355). Measured: for \hat m_{ml}=\frac1K\sum y_k, the variance \sigma^2/K decreases as K grows: at K=4,16,64, the measured variances are 1.0151, 0.2433, 0.0641, decreasing exactly as 1/K.

Example 6.3: checking the bias of \hat m_{ml} and \hat\sigma^2_{ml}

Part (a): E[\hat m_{ml}]=E[\frac1K\sum Y_k]=\frac1K\cdot Km=m, so \hat m_{ml} is unbiased. Part (b), with m known: E[\hat\sigma^2_{ml}]=E[\frac1K\sum(Y_k-m)^2]=\sigma^2, so \hat\sigma^2_{ml} (when m is known) is also unbiased (Barkat, pp. 355 to 356). Measured: a Monte Carlo simulation of 20,000 trials with m=3, \sigma=2, K=8: the average of \hat m_{ml} across trials is 3.0086 (matching m=3), and the average of \hat\sigma^2_{ml} (using the true m) is 4.0016 (matching \sigma^2=4).

The trap: when m is also estimated, \hat\sigma^2_{ml} becomes biased

The unbiased result of Example 6.3(b) holds only when m is known exactly. In practice one usually must estimate both m and \sigma^2 simultaneously (met in part 4, Example 6.6): using \hat m_{ml} in place of the true m in \hat\sigma^2_{ml}=\frac1K\sum(y_k-\hat m_{ml})^2 costs one \"degree of freedom,\" making E[\hat\sigma^2_{ml}]=\frac{K-1}{K}\sigma^2 — biased below the true value (this is exactly why the familiar sample-variance estimator usually divides by K-1 instead of K to become unbiased). Measured: the same Monte Carlo simulation but using \hat m_{ml} estimated from the sample itself instead of the true m: the average of \hat\sigma^2_{ml} is 3.4997, against the biased prediction \frac{K-1}{K}\sigma^2 = 3.5.

Self-check, part 2

• What does the generalized likelihood ratio test (GLRT) replace the unknown parameter with?
• Is an unbiased estimator with a large variance always a \"good\" estimator?
• Why does \hat\sigma^2_{ml} become biased when m must also be estimated?

Hints: the maximum likelihood estimate \hat\theta_{ml} of that parameter; no, small variance is also needed; because using \hat m_{ml} in place of the true m costs one degree of freedom.

The Bayes risk for a random parameter

In Bayes estimation, a nonnegative cost C(\theta,\hat\theta) is assigned to each pair (\theta,\hat\theta), and the risk \Re=E[C(\theta,\hat\theta)] is minimised. Only costs that are functions of the error \tilde\theta=\theta-\hat\theta are considered: (1) squared error C=(\theta-\hat\theta)^2; (2) absolute value of error C=|\theta-\hat\theta|; (3) a uniform cost function, zero if |\theta-\hat\theta|<\Delta/2 and one otherwise (Barkat, pp. 355 to 357).

MMSE: minimising squared error gives the conditional mean

With squared cost, minimising the risk \Re_{ms} leads to \hat\theta_{ms}=\int\theta f_{\Theta|Y}(\theta|y)\,d\theta=E[\theta|y] — the minimum mean-square (MMSE) estimate is exactly the conditional mean of \theta given Y (Barkat, pp. 357 to 358). The corresponding minimum risk is exactly the conditional variance of \theta given Y, averaged over all values of Y: \Re_{ms}=E_Y[\text{var}(\theta|Y)].

MAVE: minimising absolute error gives the median

With absolute-value cost, minimising the risk leads to \int_{-\infty}^{\hat\theta_{abs}}f_{\Theta|Y}(\theta|y)\,d\theta=\int_{\hat\theta_{abs}}^{\infty}f_{\Theta|Y}(\theta|y)\,d\theta — meaning \hat\theta_{abs}=\hat\theta_{mave} is exactly the median of the posterior density f_{\Theta|Y}(\theta|y) (Barkat, pp. 358 to 359).

MAP: minimising the uniform cost gives the mode (density peak)

With the uniform cost, when \Delta is small enough, minimising the risk is equivalent to maximising the probability P(\hat\theta-\Delta/2\le\Theta\le\hat\theta+\Delta/2|y), i.e. maximising the posterior density f_{\Theta|Y}(\theta|y) itself at \theta=\hat\theta. This is the maximum a posteriori (MAP) estimate, solved from the MAP equation \partial\ln f_{\Theta|Y}(\theta|y)/\partial\theta=0, or equivalently (using Bayes' rule) \partial\ln f_{Y|\Theta}(y|\theta)/\partial\theta+\partial\ln f_\Theta(\theta)/\partial\theta=0 (Barkat, pp. 359 to 361).

Example 6.4: a random Gaussian parameter, the three estimates coincide

With Y_k=M+N_k, M and N_k independent, both zero-mean Gaussian and sharing the same variance \sigma^2: the posterior density f_{M|Y}(m|y) turns out still Gaussian, of mean \frac1{K+1}\sum y_k, standard deviation \sigma'=\sigma/\sqrt{K+1} (Barkat, pp. 360 to 362). Because Gaussian is symmetric, mean = median = mode, so \hat m_{ms}=\hat m_{map}=\hat m_{mave}=\frac1{K+1}\sum_{k=1}^Ky_k. Measured: with K=5 and specific sample data, the numerical posterior density (direct integration) gives conditional mean 0.6059, matching the closed form 0.6059.

Example 6.5: when the three estimates genuinely differ

Let Y=X+N, with X taking only two discrete values \{0,1\} equally likely, f_X(x)=\tfrac12\delta(x)+\tfrac12\delta(x-1), and N having Laplace density f_N(n)=\tfrac12e^{-|n|}. Since the posterior density is nonzero only at x=0,1, the MAP estimate simply picks whichever value is closer to n=y: \hat x_{map}=1 if n\ge1/2, else \hat x_{map}=0 (Barkat, pp. 362 to 364). Meanwhile \hat x_{ms}=e^{-n-1}/(e^{-n}+e^{-n-1}) turns out to be a constant independent of n — a surprising result due to the special structure of the Laplace density combined with two Dirac impulses.

Measured: \hat x_{map} jumps, \hat x_{ms} stays constant

Measured: at n=-1,0,0.5,1,2, \hat x_{map} is respectively 0,0,1,1,1 — jumping exactly at n=0.5. Meanwhile \hat x_{ms}=e^{-n-1}/(e^{-n}+e^{-n-1})=1/(1+e) at every n, a constant equal to 0.2689 — entirely independent of n, strikingly different from MAP's jump behaviour; measured at five different values of n, it gives the same value 0.2689,0.2689,0.2689,0.2689,0.2689 every time, confirming it is indeed constant. This is clear evidence that the choice of cost function (and hence the estimation criterion) can substantially affect the result.

Self-check, part 3

• What is the MMSE estimate \hat\theta_{ms} exactly, in terms of the posterior density?
• When do the three MMSE, MAVE, MAP estimates coincide?
• In Example 6.5, what is surprising about \hat x_{ms}?

Hints: the conditional mean E[\theta|y]; when the posterior density is symmetric about its peak (like Gaussian); it does not depend on the observation n, always the constant 1/(1+e).

When the exact variance is hard to find: a lower bound is needed

If the parameter \theta is nonrandom (the prior density is zero for every \theta), the MAP equation (6.31) reduces to the likelihood equation (6.3) — so ML can be viewed as a special case of MAP. To assess the \"goodness\" of an unbiased estimator in this case, the error variance is needed, but its exact expression is often very hard to find. Instead, a lower bound on the error variance is studied, called the Cramér–Rao bound (Barkat, pp. 364 to 365).

The Cramér–Rao theorem

Given the observation vector Y of K samples, and \hat\theta an unbiased estimator of \theta. Then \text{var}[(\hat\theta-\theta)|\theta]\ge1/E[(\partial\ln f_{Y|\Theta}(y|\theta)/\partial\theta)^2], with E[(\partial\ln f_{Y|\Theta}/\partial\theta)^2]=-E[\partial^2\ln f_{Y|\Theta}/\partial\theta^2] (the Fisher information) (Barkat, pp. 365 to 366).

A proof sketch: the Schwarz inequality

The proof starts from E[(\hat\theta-\theta)|\theta]=0 (by unbiasedness), differentiates with respect to \theta, and uses \partial f_{Y|\Theta}/\partial\theta=f_{Y|\Theta}\cdot\partial\ln f_{Y|\Theta}/\partial\theta to get \int(\hat\theta-\theta)f_{Y|\Theta}\cdot(\partial\ln f_{Y|\Theta}/\partial\theta)\,dy=1. Applying the Schwarz inequality [\int xy\,dt]^2\le\int x^2dt\int y^2dt with x=(\hat\theta-\theta)\sqrt{f_{Y|\Theta}} and y=(\partial\ln f_{Y|\Theta}/\partial\theta)\sqrt{f_{Y|\Theta}} gives exactly the Cramér–Rao inequality (Barkat, pp. 366 to 367). Equality holds if and only if \partial\ln f_{Y|\Theta}/\partial\theta=c(\theta)[\hat\theta-\theta] — the efficiency condition.

An efficient estimator: attaining the bound exactly, and equalling ML

Any unbiased estimator satisfying equality in the Cramér–Rao inequality is called efficient. If an efficient estimator exists, it can be shown to coincide with the ML estimator — since both solve the same stationarity equation \partial\ln f_{Y|\Theta}/\partial\theta=0 at \hat\theta=\hat\theta_{ml} (Barkat, pp. 367 to 368). Note: not every ML estimator is efficient — the Schwarz-equality condition is not always satisfied; but when it is, ML is exactly the efficient choice.

Example 6.6: ML jointly for m and \sigma^2

With Y_k=m+N_k, both m and \sigma^2 unknown: solving \partial\ln f/\partial m=0 and \partial\ln f/\partial\sigma^2=0 simultaneously gives \hat m_{ml}=\frac1K\sum y_k and \hat\sigma^2_{ml}=\frac1K\sum(y_k-\hat m_{ml})^2 (Barkat, pp. 368 to 369) — exactly matching the part-2 observation that \hat\sigma^2_{ml} becomes biased once m is also estimated. \hat m_{ml} remains unbiased, since E[\hat m_{ml}]=m regardless of whether \sigma^2 is known.

Measured: \hat m_{ml} attains the Cramér–Rao bound exactly

Using the equality condition (6.50), \partial\ln f/\partial m=K(\bar y-m)/\sigma^2 has exactly the form c(\theta)[\hat\theta-\theta] with c=K/\sigma^2 — so \hat m_{ml} is an efficient estimator, and the Cramér–Rao bound is attained exactly, equal to \sigma^2/K (Barkat, pp. 369 to 370). Measured: a Monte Carlo simulation of 50,000 trials with \sigma=1.5, K=10: the empirical variance of \hat m_{ml} is 0.2249, against the theoretical Cramér–Rao bound \sigma^2/K = 0.225 — matching almost exactly, confirming efficiency.

Self-check, part 4

• Why use the Cramér–Rao bound instead of computing the error variance directly?
• What is the equality condition in the Cramér–Rao inequality?
• If an efficient estimator exists, which estimator does it coincide with?

Hints: because the exact variance is often very hard to compute; \partial\ln f_{Y|\Theta}/\partial\theta=c(\theta)[\hat\theta-\theta]; the maximum likelihood (ML) estimator.

When several parameters must be estimated at once

In many radar and communication applications, several parameters must be estimated simultaneously — for example a radar target's range and velocity, or a communication signal's amplitude, arrival time, and carrier frequency. With a nonrandom parameter vector \theta=[\theta_1\;\theta_2\;\cdots\;\theta_K]^T, the likelihood equation (6.3) becomes a system of K simultaneous equations \partial\ln f_{Y|\theta}/\partial\theta_i=0, i=1,\ldots,K, written compactly with the gradient operator as \nabla_\theta[\ln f_{Y|\Theta}(y|\theta)]=0 (Barkat, pp. 371 to 372).

The Fisher information matrix and the matrix Cramér–Rao bound

The Cramér–Rao bound extends to vectors: if \hat\theta is absolutely unbiased, the error covariance is bounded below by the inverse Fisher information matrix J: E[(\theta-\hat\theta)(\theta-\hat\theta)^T|\theta]\ge J^{-1}, with elements J_{ij}=E[(\partial\ln f_{Y|\Theta}/\partial\theta_i)(\partial\ln f_{Y|\Theta}/\partial\theta_j)|\theta] (Barkat, pp. 373 to 375). The variance of \theta_i alone is bounded below by the ith diagonal element of J^{-1} — a direct generalisation of the single-parameter case in part 4.

When the density is unknown: the best linear unbiased estimator (BLUE)

In many practical problems, the probability density of the data is unknown, only the first and second moments are known — then MMSE or the Cramér–Rao bound cannot be applied. A practical solution: restrict the estimator to a linear form \hat\theta_{lms}=aY+b, and find a,b minimising the mean-square risk. Solving the two zero-derivative equations gives a=(E[\theta Y]-E[\theta]E[Y])/(E[Y^2]-E^2[Y]) and b=E[\theta]-aE[Y] — expressions needing only the first and second moments (Barkat, pp. 378 to 380).

Measured: BLUE matches a numerical grid minimum

Measured: simulating \theta and Y linearly correlated with noise (200,000 samples), computing a,b from the closed form above, then comparing the mean-square risk at that (a,b) with neighbouring points on a \pm0.01 grid. Measured: (a,b) = (0.7508, 0.9472), the risk there is 1.0735, smaller than or equal to the risk at every neighbouring grid point tried (the smallest margin is 2.25e-04), confirming the closed-form (a,b) is indeed the minimum.

Least squares: no statistics needed at all

Unlike the methods above (which rest on a minimum-variance unbiased estimator), the least-squares estimator (LSE) only minimises the squared deviation between the observed data and the assumed signal model, requiring no knowledge of the noise's statistics whatsoever. With the linear model Y=H\theta+N (H a known K\times M matrix, K\ge M), minimising J(\theta)=(Y-H\theta)^T(Y-H\theta) gives the normal equations H^TH\theta=H^TY, with solution \hat\theta_{ls}=(H^TH)^{-1}H^TY if H has full rank (Barkat, pp. 387 to 389).

Measured: the closed form matches a direct numerical solve

Measured: with a random 6\times2 matrix H and a known true \theta, generating Y=H\theta+N with small Gaussian noise, then comparing the closed form (H^TH)^{-1}H^TY with a direct numerical least-squares solve (no closed form, pure numerical optimisation). The two components of \hat\theta_{ls} by the closed form are 1.3883,-1.7178, and by the direct numerical solve are 1.3883,-1.7178 — agreeing to within 8.88e-16, confirming the closed form is correct.

Looking back at the chapter: a decision tree of estimation methods

The whole chapter suggests a practical decision tree: if the parameter is nonrandom and the density is known, use ML (checking efficiency via Cramér–Rao); if the parameter is random and a prior density is known, use Bayes (MMSE for smallest variance, MAP for the posterior peak, MAVE for a median robust to outliers); if the density is unknown but the first and second moments are known, use BLUE; if there is no statistical information at all, only a linear model, use least squares. This is the overall picture linking chapter 5 (decision) and chapter 6 (estimation): the same Bayesian framework, differing only in what is sought — discrete (a hypothesis) or continuous (a parameter).

Self-check, part 5

• What does the matrix Cramér–Rao bound bound from below?
• What information about the data does the BLUE estimator need?
• How does least squares (LSE) differ from the earlier methods?

Hints: the error covariance of the vector estimator; only the first and second moments, no probability density needed; it needs no noise statistics at all, only minimising the squared deviation from a linear model.

Takeaways

• The ML estimate \hat\theta maximises the likelihood function L(\theta), solved from the likelihood equation \partial\ln f_{Y|\Theta}/\partial\theta=0; the invariance property allows a convenient change of variable before differentiating.
• A good estimator needs to be unbiased, have small variance (ideally MVU), and be consistent; \hat\sigma^2_{ml} becomes biased when the mean must also be estimated simultaneously.
• Three cost functions — squared, absolute, uniform — give three Bayes estimates MMSE (E[\theta|y]), MAVE (median), MAP (mode); they coincide when the posterior density is symmetric.
• The Cramér–Rao inequality gives the lower bound 1/E[(\partial\ln f_{Y|\Theta}/\partial\theta)^2] on the error variance; an estimator attaining the bound exactly is called efficient, and coincides with ML.
• Density unknown but moments known: use BLUE \hat\theta=aY+b; no statistics at all: use least squares \hat\theta_{ls}=(H^TH)^{-1}H^TY from the linear model Y=H\theta+N.