Point estimation: bias, variance and MSE · Probability & Statistics

Everything so far has gone one direction: given a distribution, work out how the data behaves. Statistical inference runs the other way — given data, work out the distribution.

The simplest version is point estimation: produce a single number as your best guess at an unknown parameter, and be honest about how good it is.

Estimators are random variables

An estimator \hat\theta is a function of the data used to estimate a parameter \theta. Since the data is random, the estimator is a random variable — it has its own distribution, called the sampling distribution.

This is the central idea of the section. \theta is a fixed unknown constant; \hat\theta is random and varies from sample to sample. Every question about "how good is this estimate" is a question about the distribution of \hat\theta.

Notation: \hat\theta is the estimator (a rule); the number it produces on a particular dataset is an estimate.

Bias

$\operatorname{Bias}(\hat\theta) = E[\hat\theta] - \theta$

An estimator is unbiased if E[\hat\theta] = \theta — on average, over many samples, it lands on the truth.

The sample mean is unbiased for \mu:

$E[\bar X] = \frac{1}{n}\sum E[X_i] = \frac{n\mu}{n} = \mu$

Where n-1 comes from

Here is the classic result, and it finally answers §0's question.

Consider the "obvious" variance estimator with divisor n:

$\hat\sigma^2 = \frac{1}{n}\sum(X_i - \bar X)^2$

It is biased downward:

$E[\hat\sigma^2] = \frac{n-1}{n}\sigma^2$

The reason: the deviations are measured from \bar X, not the true \mu. And \bar X is, by construction, the point that minimises the sum of squared deviations for this sample. So the sum comes out too small — the estimator has been flattered by using the data twice.

Correcting the factor gives the sample variance:

$s^2 = \frac{1}{n-1}\sum(X_i - \bar X)^2, \qquad E[s^2] = \sigma^2$

The n-1 is the number of degrees of freedom: once \bar X is fixed, only n-1 of the deviations can vary freely — the last is determined, since they must sum to zero.

A caution: s^2 is unbiased for \sigma^2, but s is not unbiased for \sigma. Taking a square root does not preserve unbiasedness, because E[\sqrt{Y}] \ne \sqrt{E[Y]} (Jensen's inequality). This is a good example of unbiasedness being a fragile property.

Variance, and the trade-off

Unbiasedness alone is not enough. An estimator can be centred correctly and still be useless if it's wildly variable.

Consider estimating \mu by "just use X_1 and ignore the rest". It's perfectly unbiased — E[X_1] = \mu — but has variance \sigma^2 rather than \sigma^2/n. Nobody would use it.

So we need a criterion combining both.

Mean squared error

$\operatorname{MSE}(\hat\theta) = E\big[(\hat\theta - \theta)^2\big]$

And the decomposition that makes it useful:

$\boxed{\;\operatorname{MSE}(\hat\theta) = \operatorname{Var}(\hat\theta) + \big[\operatorname{Bias}(\hat\theta)\big]^2\;}$

Proof. Add and subtract E[\hat\theta]:

$E\big[(\hat\theta - \theta)^2\big] = E\Big[\big((\hat\theta - E[\hat\theta]) + (E[\hat\theta] - \theta)\big)^2\Big]$

Expanding gives three terms. The first is \operatorname{Var}(\hat\theta), the third is \text{Bias}^2 (a constant), and the cross term vanishes because E[\hat\theta - E[\hat\theta]] = 0. ∎

Total error splits into being wrong on average, and being unstable. That's the bias–variance decomposition, and it's the same one that governs overfitting in machine learning.

Biased can beat unbiased

Since MSE weighs both, an estimator with a little bias and much less variance can win outright.

The clearest example is estimating \sigma^2 for a Normal sample. Comparing divisors:

Divisor	Bias	MSE
n - 1	0	\dfrac{2\sigma^4}{n-1}
n	-\sigma^2/n	\dfrac{(2n-1)\sigma^4}{n^2}
n + 1	-2\sigma^2/(n+1)	\dfrac{2\sigma^4}{n+1}

The n+1 divisor has the lowest MSE of the three, despite being the most biased. Unbiasedness is a convention, not a law — and shrinkage methods like ridge regression (§6) deliberately trade bias for variance.

Consistency

$\hat\theta_n \xrightarrow{p} \theta \quad \text{as } n \to \infty$

A consistent estimator converges to the truth with enough data. A sufficient condition, straight from the MSE decomposition: if both bias and variance $\to 0$, then MSE \to 0, which implies consistency.

Note that consistency and unbiasedness are independent properties. \hat\sigma^2 with divisor n is biased at every n, but its bias \to 0, so it is consistent. Conversely X_1 is unbiased for every n but never converges — so it is not consistent.

Worked example

X_1, \dots, X_n \sim \text{Uniform}(0, \theta). Compare two estimators of \theta: \hat\theta_1 = 2\bar X and \hat\theta_2 = X_{(n)} (the maximum).

Estimator 1. Since E[X_i] = \theta/2:

$E[\hat\theta_1] = 2 \cdot \frac{\theta}{2} = \theta \quad \text{— unbiased}$

With \operatorname{Var}(X_i) = \theta^2/12:

$\operatorname{Var}(\hat\theta_1) = 4 \cdot \frac{\theta^2}{12n} = \frac{\theta^2}{3n} = \operatorname{MSE}$

Estimator 2. From §2's order statistics, f_{X_{(n)}}(x) = nx^{n-1}/\theta^n, giving

$E[X_{(n)}] = \frac{n}{n+1}\theta \quad \text{— biased downward}$

which makes sense: the maximum of a sample can never exceed \theta.

$\operatorname{Var}(X_{(n)}) = \frac{n\theta^2}{(n+1)^2(n+2)}$

$\operatorname{MSE}(\hat\theta_2) = \frac{n\theta^2}{(n+1)^2(n+2)} + \frac{\theta^2}{(n+1)^2} = \frac{2\theta^2}{(n+1)(n+2)}$

The comparison is dramatic. MSE of the first is O(1/n); MSE of the second is O(1/n^2). At n = 10:

$\operatorname{MSE}(\hat\theta_1) = \frac{\theta^2}{30} \approx 0.0333\theta^2, \qquad \operatorname{MSE}(\hat\theta_2) = \frac{2\theta^2}{132} \approx 0.0152\theta^2$

The biased estimator is more than twice as good, and the gap widens with n. And the bias is trivially removable: \frac{n+1}{n}X_{(n)} is unbiased and keeps the O(1/n^2) rate, beating both.

The reason the maximum wins is that it uses the structure of the problem — for a uniform, the largest observation is enormously informative about the upper endpoint, while the mean throws that information away. That intuition becomes the formal notion of sufficiency, three lessons from now.

Doing it in Python

Sampling distributions are what simulation is for — generate many samples, compute the estimator on each, and look at the resulting spread:

import numpy as np

rng = np.random.default_rng(0)
theta, n, trials = 10.0, 10, 200_000

samples = rng.uniform(0, theta, size=(trials, n))
est1 = 2 * samples.mean(axis=1)          # 2 * sample mean
est2 = samples.max(axis=1)               # the maximum
est3 = (n + 1) / n * est2                # bias-corrected maximum

for name, est in [("2 * mean", est1), ("max", est2), ("(n+1)/n * max", est3)]:
    bias = est.mean() - theta
    var = est.var()
    print(f"{name:>16}  bias {bias:>8.4f}  var {var:>8.4f}  MSE {var + bias**2:>8.4f}")

Verifying the n-1 correction directly:

import numpy as np

rng = np.random.default_rng(1)
mu, sigma, n, trials = 0.0, 3.0, 5, 400_000

x = rng.normal(mu, sigma, size=(trials, n))

biased = x.var(axis=1, ddof=0)      # divide by n
unbiased = x.var(axis=1, ddof=1)    # divide by n-1

print("true variance      :", sigma**2)
print("E[divide by n]     :", round(biased.mean(), 4),
      f"  (theory {(n-1)/n * sigma**2:.4f})")
print("E[divide by n-1]   :", round(unbiased.mean(), 4), "  <- unbiased")
print()
print("but the square root is NOT unbiased:")
print("  E[s]   =", round(np.sqrt(unbiased).mean(), 4))
print("  sigma  =", sigma)

And the MSE comparison across divisors, showing n+1 winning:

import numpy as np

rng = np.random.default_rng(2)
sigma, n, trials = 2.0, 8, 400_000
x = rng.normal(0, sigma, size=(trials, n))
ss = ((x - x.mean(axis=1, keepdims=True)) ** 2).sum(axis=1)

print(f"{'divisor':>10} {'bias':>10} {'variance':>12} {'MSE':>12}")
for d, label in [(n - 1, "n-1"), (n, "n"), (n + 1, "n+1")]:
    est = ss / d
    bias = est.mean() - sigma**2
    print(f"{label:>10} {bias:>10.4f} {est.var():>12.4f} {est.var() + bias**2:>12.4f}")
print("\nLowest MSE goes to n+1 -- the most biased of the three.")

Your turn

1. \hat\theta has bias 2 and variance 9. Find its MSE.

2. Is \bar X consistent for \mu? Prove it.

3. \hat\theta_1 is unbiased with variance 10; \hat\theta_2 has bias 1 and variance 4. Which would you prefer?

Solutions

$\operatorname{MSE} = \operatorname{Var} + \text{Bias}^2 = 9 + 4 = 13$

2. Yes. Two routes.

Directly: \operatorname{Bias} = 0 and $\operatorname{Var}(\bar X) = \sigma^2/n \to 0$, so \operatorname{MSE} \to 0, which implies convergence in probability.

By Chebyshev:

$P(|\bar X - \mu| > \varepsilon) \le \frac{\sigma^2}{n\varepsilon^2} \to 0$

Either way \bar X \xrightarrow{p} \mu — which is precisely the Weak Law of Large Numbers. Consistency of the sample mean is the LLN.

3. Prefer \hat\theta_2, on MSE:

$\operatorname{MSE}_1 = 10 + 0 = 10, \qquad \operatorname{MSE}_2 = 4 + 1 = 5$

Half the expected squared error, despite the bias.

That said, MSE isn't automatically the right criterion. If systematic error matters more than random error — a scale that always reads 1 kg heavy is worse than one that's noisy but centred — unbiasedness may be worth paying for. And if the bias direction is known, it can often just be corrected, as with the uniform maximum above. The decomposition tells you what you're trading; it doesn't tell you what to value.

Check yourself in code

Compare three estimators of \theta for a Uniform(0,\theta) sample and confirm the biased maximum beats the unbiased 2\bar X on MSE.

Print exactly this:

2*mean   bias 0.0026 MSE 3.3263
max      bias -0.9075 MSE 1.5128
corrected bias 0.0018 MSE 0.8341
max beats 2*mean: True

Use default_rng(0), \theta = 10, n = 10, 200000 trials. Round biases to 4 decimal places and MSEs to 4.

import numpy as np

rng = np.random.default_rng(0)
theta, n, trials = 10.0, 10, 200_000
samples = rng.uniform(0, theta, size=(trials, n))

est1 = 2 * samples.mean(axis=1)
est2 = samples.max(axis=1)
est3 = (n + 1) / n * est2

def report(name, est):
    bias = est.mean() - theta
    mse = est.var() + bias**2
    print(f"{name} bias {round(bias, 4)} MSE {round(mse, 4)}")
    return mse

m1 = report("2*mean  ", est1)
# Report "max     " and "corrected", then whether max beats 2*mean on MSE.

import numpy as np

rng = np.random.default_rng(0)
theta, n, trials = 10.0, 10, 200_000
samples = rng.uniform(0, theta, size=(trials, n))

est1 = 2 * samples.mean(axis=1)
est2 = samples.max(axis=1)
est3 = (n + 1) / n * est2

def report(name, est):
    bias = est.mean() - theta
    mse = est.var() + bias**2
    print(f"{name} bias {round(bias, 4)} MSE {round(mse, 4)}")
    return mse

m1 = report("2*mean  ", est1)
m2 = report("max     ", est2)
m3 = report("corrected", est3)
print("max beats 2*mean:", bool(m2 < m1))

An estimator is a random variable, so judging it means describing its distribution. Bias measures whether it's centred correctly; variance measures how stable it is; MSE combines them as \operatorname{Var} + \text{Bias}^2. Unbiasedness is a convention rather than a requirement — a little bias often buys a lot of stability.

Next: two general recipes for constructing estimators, rather than guessing them.

26. Point estimation: bias, variance and MSE