How bias arises in using maximum likelihood to determine the variance of a Gaussian?

Question

While reading Christopher M. Bishop Pattern Recognition and Machine Learning, I ran into the following explanation for why there is an error in estimating the variance of Gaussian distribution using Maximum Likelihood. On page 48,

Illustration of how bias arises in using maximum likelihood to determine the variance of a Gaussian. The green curve shows the true Gaussian distribution from which data is generated, and the three red curves show the Gaussian distributions obtained by fitting to three data sets, each consisting of two data points shown in blue, using the maximum likelihood results (1.55) and (1.56). Averaged across the three data sets, the mean is correct, but the variance is systematically under-estimated because it is measured relative to the sample mean and not relative to the true mean.

I'm pasting figure 1.15 and the required equations (1.55, 1.56) here for convenience.

$\mu_{ML} = \frac{1}{N}\sum_{n=1}^{N}x_n$

$\sigma_{ML}^2 = \frac{1}{N}\sum_{n=1}^{N}(x_n - \mu_{ML})^2$

I can imagine how maximum likelihood estimation would look pictorially. However, I'm not able to understand the figure or the explanation. I would appreciate if someone could label the figure or explain what points are considered in each part of the figure and why the curve could look like that in each case.

I looked up exercise 1.12, and I'm able to understand why there is a factor of $\frac{N-1}{N}$ between the true $\sigma^2$ and $\sigma_{MLE}^2$.

${\rm I\!E}[\sigma_{MLE}^2]= (\frac{N-1}{N})\sigma^2 $

Answer 1

It appears that a very basic statistical result gets obfuscated here.

The MLE for the variance in the Gaussian case is, as the OP writes,

$$\sigma_{ML}^2 = \frac{1}{N}\sum_{n=1}^{N}(x_n - \mu_{ML})^2 = \frac{1}n\sum_{i=1}^nx_i^2 - (\bar x)^2$$

where $\bar x$ stands for the sample mean. Then

$$E[\sigma_{ML}^2] = E(X^2) - E[(\bar x)^2]$$

By Jensen's inequality,

$$E[(\bar x)^2] > [E(\bar x)]^2 = \mu ^2 \implies -E[(\bar x)^2] < -\mu ^2$$

and therefore

$$E[\sigma_{ML}^2] < E(X^2) - \mu^2 = \text{Var}(X)$$

Answer 2

If one wants to realise the impact of using the MLE estimate on the estimation of the model density, one needs more than three replicas: here is an experiment with 10⁴ simulations of a N$(0,1)$ sample of size 10:

where the red curve is the N$(0,1)$ density and the gray curves are the estimated ones. The R code associated with the picture is

samplz=matrix(rnorm(1e5),1e4,1e1)
curve(dnorm,-4,4,col="sienna",lwd=2,ylab="f(x)",ylim=c(0,1))
for (i in 1:1e4)
curve(dnorm(x,mean=mean(samplz[i,]),sd=sqrt(9/10)*sd(samplz[i,])),
      from=-4,to=4,col="gray",add=TRUE)
curve(dnorm,-4,4,col="sienna",lwd=2,ylab="f(x)",ylim=c(0,1),add=TRUE)

If instead one uses the unbiased estimator in

dnorm(x,mean=mean(samplz[i,]),sd=sd(samplz[i,]))

the curves are slightly less peaked but the difference is hard to spot:

since

 > dnorm(0)
 [1] 0.3989423
 > dnorm(0,sd=sqrt(9/10))
 [1] 0.4205221

means that at its apex the difference in the densities is about 0.02.

$\qquad\qquad$

From a more global perspective, the fact that the MLE is biased is unsurprising and not a cause for concern. Since the MLE of a transform of the parameter $\theta$ is the (same) transform of the MLE of $\theta$, the MLE is necessarily biased for most transforms of the parameter $\theta$ or, equivalently, for most parameterisations. This does not clash with the efficiency of the MLE (in the Cramèr-Rao sense) in regular settings.