Statistics/Point Estimation

Introduction

Usually, a random variable $X$ resulting from a random experiment is assumed to follow a certain distribution with an unknown (but fixed) parameter (vector) $\theta \in {ℝ}^k$ ($k$ is a positive integer, and its value depends on the distribution), taking value in a set $\mathrm{\Theta }$, called the parameter space.

Remark.

• In the context of frequentist statistics (the context here), parameters are regarded as fixed.
• On the other hand, in the context of Bayesian statistics, parameters are regarded as random variables.

For example, suppose the random variable $X$ is assumed to follow a normal distribution $𝒩(\mu ,{\sigma }^2)$. Then, in this case, the parameter vector $\theta =(\mu ,\sigma )\in \mathrm{\Theta }$ is unknown, and the parameter space $\mathrm{\Theta }=\{(\mu ,\sigma ):\mu \in ℝ,\sigma >0\}$. It is often useful to estimate those unknown parameters in some ways to "understand" the random variable $X$ better. We would like to make sure the estimation should be "good" enough, so that the understanding is more accurate.

Intuitively, the (realization of) random sample $X_1,\dots ,X_n$ should be useful. Indeed, the estimators introduced in this chapter are all based on the random sample in some sense, and this is what point estimates mean. To be more precise, let us define point estimation and point estimates.

Definition. (Point estimation) Point estimation is a process of using the value of a statistic to give a single value estimate (which can be interpreted as a point) of an unknown parameter.

Remark.

• Recall that statistics are functions of a random sample.
• We call the unknown parameter as population parameter (since the underlying distribution corresponding to the parameter is called a population).
• The statistics is called a point estimator, and its realization is called a point estimate.

• The notation of point estimator commonly has a $\hat{}$.

• Point estimation will be contrasted with interval estimation, which uses the value of a statistic to estimate an interval of plausible values of the unknown parameter.

Example. Suppose $X_1,\dots ,X_n$ are $n$ random samples from the normal distribution $𝒩(\mu ,{\sigma }^2)$.

• We may use the statistic $\bar{X}=\frac{X_1+\cdots +X_n}{n}$ to estimate $\mu$ intuitively, and $\bar{X}$ is called the point estimator, and its realization $\bar{x}$ is called the point estimate.
• Alternatively, we may simply use the statistic $X_1$ (despite it does not involve $X_2,\dots ,X_n$, it can still be regarded as function of $X_1,\dots ,X_n$) to estimate $\mu$. That is, we use the value of the first random sample from the normal distribution as the point estimate of the mean of the distribution! Intuitively, it may seem that such estimator is quite "bad".

• Such estimator, which just takes one random sample directly, is called a single observation estimator.
• We will later discuss how to evaluate how "good" a point estimator is.

In the following, we will introduce two well-known point estimators, which are actually quite "good", namely maximum likelihood estimator and method of moment estimator.

Maximum likelihood estimator (MLE)

As suggested by the name of this estimator, it is the estimator that maximize some kind of "likelihood". Now, we would like to know what "likelihood" should we maximize to estimate the unknown parameter(s) (in a "good" way). Also, as mentioned in the introduction section, the estimator is based on the random sample in some sense. Hence, this "likelihood" should be also based on the random sample in some sense.

To motivate the definition of maximum likelihood estimator, consider the following example.

Example. In a random experiment, a (fair or unfair) coin is tossed once. Let the random variable $X=1$ if head comes up, and $0$ otherwise. Then, the pmf of $X$ is $f(x;p)=p^x(1-p{)}^{1-x},x\in \{0,1\}$, in which the unknown parameter $p$ represents the probability for head comes up, and $p\in \mathrm{\Theta }=\{p:p\in (0,1)\}$.

Now, suppose you get a random sample $X_1,X_2,\dots ,X_n$ by tossing that coin $n$ independent times (such random sample is called an independent random sample, since the random variables involved are independent), and the corresponding realizations $x_1,x_2,\dots ,x_n$. Then, the probability for $X_1=x_1,X_2=x_2,\dots ,\text{ and }{X}_n=x_n$, i.e., the random sample have these realizations exactly, is $$\begin{array}{ccc}\mathbb{P}(X_1=x_1\cap X_2=x_2\cap \cdots \cap X_n=x_n)& =\mathbb{P}(X_1=x_1)\mathbb{P}(X_2=x_2)\cdots \mathbb{P}(X_n=x_n)& \text{by independence}\\ & =f(x_1;p)f(x_2;p)\cdots f(x_n;p)\\ & =p^{x_1}(1-p{)}^{1-x_1}{p}^{x_2}(1-p{)}^{1-x_2}\cdots p^{x_n}(1-p{)}^{1-x_n}\\ & =p^{x_1+x_2+\cdots +x_n}(1-p{)}^{n-x_1-x_2-\cdots -x_n}.\end{array}$$

Remark.

• Remark on notation: You may observe that there is an additional "$;p$" in the pmf of $X$. Such notation means the pmf is with the parameter value $p$. It is included to emphasize the parameter value we are referring to.
• In general, we write $f(\cdot ;\theta )$ for pmf/pdf with the parameter value $\theta$ ($\theta$ may be a vector).

• There are some alternative notations with the same meaning: $f(\cdot |\theta ),f_{\theta }(\cdot ),\dots$.

• Similarly, we have similar notations, e.g. ${\mathbb{P}}_{\theta }(A),\mathbb{P}(A|\theta ),\mathbb{P}(A;\theta ),\dots$, to mean the probability for event $A$ to happen, with the parameter value $\theta$. (It is more common to use the first notation: ${\mathbb{P}}_{\theta }(A)$.)
• We also have similar notations for mean, variance, covariance, etc., like ${𝔼}_{\theta }[\cdot ],{\mathrm{Var}}_{\theta }(\cdot ),{\mathrm{Cov}}_{\theta }(\cdot ),\dots$

Intuitively, with these particular realizations (fixed), we would like to find a value of $p$ that maximizes this probability, i.e.,, makes the realizations obtained to be the one that is "most probable" or "with maximum likelihood". Now, let us formally define the terms related to MLE.

Definition. (Likelihood function) Let $X_1,\dots ,X_n$ be a random sample with a joint pmf or pdf $f$, and the parameter (vector) $\theta \in \mathrm{\Theta }$ ($\mathrm{\Theta }$ is the parameter space). Suppose $x_1,\dots ,x_n$ are the corresponding realizations of the random sample $X_1,\dots ,X_n$. Then, the likelihood function, denoted by $ℒ(\theta ;x_1,\dots ,x_n)$, is the function $\theta \mapsto f(x_1,\dots ,x_n;\theta )$ ($\theta$ is a variable, and $x_1,\dots ,x_n$ are fixed).

Remark.

• For simplicity, we may use the notation $ℒ(\theta ;𝐱)$ instead of $ℒ(\theta ;x_1,\dots ,x_n)$. Sometimes, we may also just write "$ℒ(\theta ;𝐱)$" for convenience.

• When we replace $x_1,\dots ,x_n$ by $X_1,\dots ,X_n$, then the resulting "likelihood function" becomes a random variable, and we denote it by $ℒ(\theta ;X_1,\dots ,X_n)$ or $ℒ(\theta ;𝐗)$.

• The likelihood function is in contrast with the joint pmf or pdf itself, where $\theta$ is fixed and $x_1,\dots ,x_n$ are variables.
• When the random sample comes from a discrete distribution, then the value of likelihood function is the probability $\mathbb{P}(X_1=x_1\cap \cdots \cap X_n=x_n)$ at the parameter vector $\theta$. That is, the probability for getting this specific realization exactly.
• When the random sample comes from a continuous distribution, then the value of likelihood function is not a probability. Instead, it is only the value of the joint pdf at $(x_1,\dots ,x_n)$ (which can be greater than one). However, the value can still be used to "reflect" the probability for getting "very close to" this specific realization, where the probability can be obtained by integrating the joint pdf over a "very small" region around $(x_1,\dots ,x_n)$.
• The natural logarithm of the likelihood function, $\ln ℒ(\theta ;𝐱)$ (or $\ln ℒ(\theta ;𝐗)$ sometimes), is called the log-likelihood function.
• Notice that the "expression" of the likelihood function is actually the same as that of the joint pdf, and just the inputs are different. So, one may still integrate/sum the likelihood function with respect to $x_1,\dots ,x_n$ (which changes the likelihood function to the joint pdf/pmf in such context in some sense) as if it is the joint pdf/pmf to get probabilities.

Definition. (Maximum likelihood estimate) Given a likelihood function $ℒ(\theta ;𝐱)$, a maximum likelihood estimate of the parameter $\theta$ is a value $\hat{\theta }(𝐱)$ at which $ℒ(\theta ;𝐱)$ is maximized.

Remark.

• The maximum likelihood estimator (MLE) of $\theta$ is $\hat{\theta }(𝐗)$ (obtained by replacing "$x$" in $\hat{\theta }(𝐱)$ by "$X$").

• In some other places, the abbreviation MLE can also mean maximum likelihood estimate depending on the context. However, we will just use the abbreviation MLE when we are talking about maximum likelihood estimator here.

• Since $\frac{d}{dy}\ln y=\frac{1}{y}>0$ (the domain of natural logarithm function is the set of all positive real numbers), the natural logarithm function is strictly increasing, i.e., the output is larger when the input is larger. Thus, when we find a value at which $\ln ℒ(\theta ;𝐱)$ is maximized, $ℒ(\theta ;𝐱)$ is also maximized at the same value.

Now, let us find the MLE of the unknown parameter $p$ in the previous coin flipping example.

Example. (Motivating example revisited) Recall that we use a coin flipping example to motivate maximum likelihood estimation. $X$ follows the Bernoulli distribution with success probability $p$. The pmf of $X$ is $f(x;p)=p^x(1-p{)}^{1-x}$. $X_1,\dots ,X_n$ is a random sample from the distribution.

• The likelihood function $ℒ(p)$ is the joint pmf of $X_1,\dots ,X_n$,

$$\begin{array}{ccc}\mathbb{P}(X_1=x_1\cap \cdots \cap X_n=x_n)& =\prod _{i=1}^{n}f(x_i;p)& \text{by independence}\\ & =\prod _{i=1}^n{p}^{x_i}(1-p{)}^{1-x_i}\\ \end{array}$$

• The log-likelihood function $\ln ℒ(p)$ is thus

$$\begin{array}{cc}\ln ℒ(p)& =\sum _{i=1}^n\ln (p^{x_i}(1-p{)}^{1-x_i})\\ & =\sum _{i=1}^n(\ln (p^{x_i})+\ln ((1-p{)}^{1-x_i}))\\ & =\sum _{i=1}^n(x_i\ln (p)+(1-x_i)\ln (1-p))\\ & =\sum _{i=1}^n(x_i\ln (p))+\sum _{i=1}^n((1-x_i)\ln (1-p))\\ & =\ln (p)\sum _{i=1}^n(x_i)+\ln (1-p)\sum _{i=1}^n(1-x_i)\\ & =\ln (p)\sum _{i=1}^n(x_i)+\ln (1-p)(n-\sum _{i=1}^n(x_i))\\ \end{array}$$

• To find the maximum of the log-likelihood function, we may use derivative test learnt in Calculus. Differentiating $\ln ℒ(p)$ with respect to $p$ gives

$$\begin{array}{cc}\frac{d\ln ℒ(p)}{dp}& =\frac{1}{p}\underset{\text{constant wrt }p}{\underbrace{\sum _{i=1}^n{x}_i}}-\frac{1}{1-p}\underset{\text{constant wrt }p}{\underbrace{(n-\sum _{i=1}^n{x}_i)}}\\ & =\frac{(1-p)\sum _{i=1}^n{x}_i-np+p\sum _{i=1}^n{x}_i}{p(1-p)}\\ & =\frac{(1-p)n\bar{x}-np+pn\bar{x}}{p(1-p)}& (\sum _{i=1}^n{x}_i=n\bar{x}=n\cdot \frac{\sum _{i=1}^n{x}_i}{n})\\ & =\frac{n(\bar{x}-p)}{p(1-p)}\end{array}$$

• To find critical point(s) of $\ln ℒ(p)$, we set $\frac{d\ln ℒ(p)}{dp}=0⟹\frac{n(\bar{x}-p)}{p(1-p)}=0⟹p=\bar{x}$ (we have $p(1-p)\ne 0$)
• To verify that $\ln ℒ(p)$ actually attains maximum (instead of minimum) at $p=\bar{x}$, we need to perform derivative test. In this case, we use first derivative test.
• We can see that $\frac{d\ln ℒ(p)}{dp}>0$ when $p<\bar{x}$, which makes $\bar{x}-p>0$, and thus $\frac{d\ln ℒ(p)}{dp}>0$. On the other hand, when $p>\bar{x}$, this makes $\bar{x}-p<0$, and thus $\frac{d\ln ℒ(p)}{dp}<0$. As a result, we can conclude that $\ln ℒ(p)$ attains its maximum at $p=\bar{x}$. It follows that the MLE of $p$ is $\bar{X}$ (not $\bar{x}$, which is instead maximum likelihood estimate!)

Exercise. Use second derivative test to verify that $\ln ℒ(p)$ attains maximum at $p=\bar{x}$.

Solution

• Since $\frac{d^2\ln ℒ(p)}{d{p}^2}=\frac{-np(1-p)-n(\bar{x}-p)(2p)}{p^2(1-p{)}^2}$, in which the numerator is negative and the denominator is positive. Thus, $\frac{d^2\ln ℒ(p)}{d{p}^2}<0$. By second derivative test, this means $\ln ℒ(p)$ attains maximum at $p=\bar{x}$.

Sometimes, there is constraint imposed on the parameter when we are finding its MLE. The MLE of the parameter in this case is called a restricted MLE. We will illustrate this in the following example.

Example. Continue from the previous coin flipping example. Suppose we have a constraint on $p$ where $0\le p\le \frac{1}{2}$. Find the MLE of $p$ in this case.

Solution: For the steps about deriving likelihood function and log-likelihood function, they are the same in this case. Without the restriction, the MLE of $p$ is $\bar{X}$. Now, with the restriction, the MLE of $p$ is $\bar{X}$ only when $\bar{X}\le \frac{1}{2}$ (we always have $\bar{X}\ge 0$ since $X\ge 0$).

If $\bar{X}>\frac{1}{2}$ (and thus $\bar{x}>1/2$), even though $\ln ℒ(p)$ is maximized at $p=\bar{x}$, we cannot set the MLE to be $\bar{X}$ due to the restriction on $p$: $0\le p\le \frac{1}{2}$. Under this case, this means $\frac{d\ln ℒ(p)}{dp}>0$ when $p\le \frac{1}{2}<\bar{X}$ (we have $\frac{d\ln ℒ(p)}{dp}>0$ when $p<\bar{x}$ from previous example), i.e., $\ln ℒ(p)$ is strictly increasing when $p\le \frac{1}{2}$. Thus, $\ln ℒ(p)$ is maximized when $p=\frac{1}{2}$ with the restriction. As a result, the MLE of $p$ is $\frac{1}{2}$ (the MLE can be a constant, which can still be regarded as a function of $X_1,\dots ,X_n$).

Therefore, the MLE of $p$ can be written as a case defined function: $\hat{\theta }=\{\begin{array}{cc}\bar{X},& \bar{X}\le \frac{1}{2}\\ \frac{1}{2},& \bar{X}>\frac{1}{2}\end{array}$, or it can be written as $\hat{\theta }=min\{\bar{X},\frac{1}{2}\}$

Exercise. Find the MLE of $p$ when $\frac{1}{2}\le p\le 1$.

Solution

• When $\bar{X}<\frac{1}{2}$, we cannot set the MLE to be $\bar{X}$ due to the restriction. In this case, we know that $\frac{d\ln ℒ(p)}{dp}<0$ when $p\ge \frac{1}{2}>\bar{X}$, i.e., $\ln ℒ(p)$ is strictly decreasing when $\frac{1}{2}\le p\le 1$. Thus, $\ln ℒ(p)$ is maximized at $p=\frac{1}{2}$, and so the MLE of $p$ is $\frac{1}{2}$.
• When $\bar{X}\ge \frac{1}{2}$, we can set the MLE to be $\bar{X}$ at which $\ln ℒ(p)$ is maximized, and so $\bar{X}$ is the MLE of $p$ in this case.
• Therefore, the MLE of $p$ is $\hat{\theta }=max\{\bar{X},\frac{1}{2}\}$.

To find the MLE, we sometimes use methods other than derivative test, and we do not need to find the log-likelihood function. Let us illustrate this in the following example.

Example. Let $X_1,\dots ,X_n$ be a random sample from the uniform distribution $𝒰[0,\beta ]$. Find the MLE of $\beta$.

Solution: The pdf of the uniform distribution is $f(x;\beta )=\frac{1}{\beta }\mathrm{𝟏}\{0\le x\le \beta \}$. Thus, the likelihood function is $ℒ(\beta )=\prod _{i=1}^n\frac{1}{\beta }\mathrm{𝟏}\{0\le x_i\le \beta \}=\frac{1}{{\beta }^n}\prod _{i=1}^n\mathrm{𝟏}\{0\le x_i\le \beta \}$.

In order for $ℒ(\beta )$ to attain maximum, first, we need to ensure that $0\le x_i\le \beta$ for each $i\in \{1,\dots ,n\}$, so that the product of the indicator functions in the likelihood function is nonzero (the value is actually one in this case). Apart from that, since $\beta \mapsto \frac{1}{{\beta }^n}$ is a strictly decreasing function of $\beta$ (because $\frac{d}{d\beta }\left(\frac{1}{{\beta }^n}\right)=\frac{-n}{{\beta }^{n+1}}<0$ (we have $n,\beta >0$)), we should pick a $\beta$ that is as small as possible so that $\frac{1}{{\beta }^n}$, and hence $ℒ(\beta )$, is as large as possible.

As a result, we should choose a $\beta$ that is as small as possible, subject to the constraint that $0\le x_i\le \beta$ for each $i\in \{1,\dots ,n\}$, which means that $\beta \ge x_i$ (it is always the case that $x_i\ge 0$, regardless of the choice of $\beta$) for each $i\in \{1,\dots ,n\}$. It follows that $ℒ(\beta )$ attains maximum when $\beta$ is the maximum of $x_1,\dots ,x_n$. Hence, the MLE of $\beta$ is $\hat{\beta }=max\{X_1,\dots ,X_n\}$.

Exercise. Show that the MLE of $\beta$ does not exist if the uniform distribution becomes $𝒰[0,\beta )$.

Solution

Proof. In this case, the constraint from the indicator functions become $0\le x_i<\beta$ for each $i\in \{1,\dots ,n\}$. With similar argument, for the MLE of $\beta$, we should choose a $\beta$ that is as small as possible subject to this constraint, which means $\beta >x_i$ for each $i\in \{1,\dots ,n\}$. However, in this case, we cannot set $\beta$ to be the maximum of $x_1,\dots ,x_n$, or else the constraint will not be satisfied and the likelihood function becomes zero due to the indicator function. Instead, we should set $\beta$ to be slightly greater than the maximum of $x_1,\dots ,x_n$, so that the constraint can still be satisifed, and $\beta$ is quite small. However, for each such $\beta >max\{x_1,\dots ,x_n\}$, we can always chooses a smaller $\beta$ that still satisfies the constraint. For example, for each $\beta$, the smaller beta, ${\beta }^{\prime }$ can be selected as $max\{x_1,\dots ,x_n\}+\frac{\beta -max\{x_1,\dots ,x_n\}}{2}>max\{x_1,\dots ,x_n\}$ . Hence, we cannot find a minimum value of $\beta$ subject to this constraint. Thus, there is no maximum point for $\ln ℒ(p)$, and hence the MLE does not exist.

$◻$

In the following example, we will find the MLE of a parameter vector.

Example. Let $X_1,\dots ,X_n$ be a random sample from the normal distribution with mean ${\theta }_1$ and variance ${\theta }_2$, $𝒩({\theta }_1,{\theta }_2)$. Find the MLE of $({\theta }_1,{\theta }_2)$.

Solution: Let $\theta =({\theta }_1,{\theta }_2)$. The likelihood function is $ℒ(\theta ;𝐱)=\prod _{i=1}^n\frac{1}{\sqrt{2\pi {\theta }_2}}\exp (-\frac{(x_i-{\theta }_1{)}^2}{2{\theta }_2})=(2\pi {\theta }_2{)}^{-n/2}\exp (-\sum _{i=1}^n\frac{(x_i-{\theta }_1{)}^2}{2{\theta }_2})$, and hence the log-likelihood function is $\ln ℒ(\theta ;𝐱)=-\frac{n}{2}\ln (2\pi {\theta }_2)-\sum _{i=1}^n\frac{(x_i-{\theta }_1{)}^2}{2{\theta }_2}$. Since this function is multivariate, we may use the second partial derivative test from multivariable calculus to find maximum point(s). However, in this case, we actually do not need to use such test. Instead, we fix the variables one by one to make the function univariate, so that we can use the derivative test for univariate function to find maximum point (with another variable fixed).

Since $$\frac{\partial \ln ℒ(\theta ;𝐱)}{\partial {\theta }_1}=\frac{1}{{\theta }_2}\sum _{i=1}^n(x_i-{\theta }_1)$$ and $\frac{\partial \ln ℒ(\theta ;𝐱)}{\partial {\theta }_2}=-\frac{2n\pi }{4\pi {\theta }_2}+\frac{1}{2{\theta }_2^2}\sum _{i=1}^n(x_i-{\theta }_1{)}^2=-\frac{n}{2{\theta }_2}+\frac{1}{2{\theta }_2^2}\sum _{i=1}^n(x_i-{\theta }_1{)}^2$.

Also, $$\frac{\partial \ln ℒ(\theta ;𝐱)}{\partial {\theta }_1}=0⟹\sum _{i=1}^n(x_i-{\theta }_1)=0⟹-n{\theta }_1+\sum _{i=1}^n{x}_i=0⟹{\theta }_1=\frac{\sum _{i=1}^n{x}_i}{n}=\bar{x}$$, which is independent from ${\theta }_2$ (this is important for us to use this kind of method) and $\frac{\partial \ln ℒ(\theta ;𝐱)}{\partial {\theta }_2}=0⟹\frac{n}{2{\theta }_2}=\frac{1}{2{\theta }_2^2}(\sum _{i=1}^n(x_i-{\theta }_1{)}^2)⟹n=\frac{1}{{\theta }_2}(\sum _{i=1}^n(x_i-{\theta }_1{)}^2)⟹{\theta }_2=\frac{\sum _{i=1}^n(x_i-{\theta }_1{)}^2}{n}$.

Since $$\frac{{\partial }^2\ln ℒ(\theta ;𝐱)}{\partial {\theta }_1^2}=\frac{\partial }{\partial {\theta }_1}(\frac{1}{{\theta }_2}\sum _{i=1}^n(x_i-{\theta }_1))=\frac{1}{{\theta }_2}\sum _{i=1}^n(-1)=\frac{-n}{{\theta }_2}<0$$, by the second derivative test (for univariate function), $\ln ℒ(\theta ;𝐱)$ is maximized at ${\theta }_1=\bar{x}$, given any fixed ${\theta }_2$.

On the other hand, since $$\frac{{\partial }^2\ln ℒ(\theta ;𝐱)}{\partial {\theta }_2^2}=\frac{n}{2{\theta }_2^2}-\frac{1}{{\theta }_2^3}\sum _{i=1}^n(x_i-{\theta }_1{)}^2$$, and thus $${\frac{{\partial }^2\ln ℒ(\theta ;𝐱)}{\partial {\theta }_2^2}|}_{{\theta }_2=\frac{\sum _{i=1}^n(x_i-{\theta }_1{)}^2}{n}}=\frac{1}{2n{(\sum _{i=1}^n(x_i-{\theta }_1{)}^2)}^2}-\frac{n^3}{{(\sum _{i=1}^n(x_i-{\theta }_1{)}^2)}^2}=\frac{1-2n^4}{2n{(\sum _{i=1}^n(x_i-{\theta }_1{)}^2)}^2}<0$$ (since $2n^4>1$).

Thus, by the second derivative test, $\ln ℒ(\theta ;𝐱)$ is maximized at ${\theta }_2=\frac{\sum _{i=1}^n(x_i-{\theta }_1{)}^2}{n}$, given any fixed ${\theta }_1$.

So, now we fix ${\theta }_1=\bar{x}$, and thus we have $\ln ℒ(\theta ;𝐱)$ is maximized at ${\theta }_2=\frac{\sum _{i=1}^n(x_i-\bar{x}{)}^2}{n}=s^2$, where $s^2$ is the realization of the sample variance $S^2$. Now, fix ${\theta }_2$ to be $s^2$, and we know that $\ln ℒ(\theta ;𝐱)$ attains maximum at ${\theta }_1=\bar{x}$ for each fixed ${\theta }_2$, including this fixed ${\theta }_2=s^2$. As a result, $\ln ℒ(\theta ;𝐱)$ is maximized at $({\theta }_1,{\theta }_2)=(\bar{x},s^2)$. Hence, the MLE of $({\theta }_1,{\theta }_2)$ is $(\bar{X},S^2)$.

Exercise.

(a) Calculate the determinant of the Hessian matrix of $\ln ℒ(\theta ;𝐱)$ at $({\theta }_1,{\theta }_2)=(\bar{x},s^2)$, which can be expressed as $\frac{{\partial }^2\ln ℒ}{\partial {\theta }_1^2}(\bar{x},s^2)\frac{{\partial }^2\ln ℒ}{\partial {\theta }_2^2}(\bar{x},s^2)-{(\frac{{\partial }^2\ln ℒ}{\partial {\theta }_2\partial {\theta }_1}(\bar{x},s^2))}^2$.

(b) Hence, verify that $({\theta }_1,{\theta }_2)=(\bar{x},s^2)$ is the maximum point of $\ln ℒ(\theta ;𝐱)$ using the second partial derivative test.

Solution

(a) First,

• $\frac{{\partial }^2\ln ℒ}{\partial {\theta }_1^2}(\bar{x},s^2)\stackrel{\text{above}}{=}{\frac{-n}{{\theta }_2}|}_{{\theta }_2=s^2}=\frac{-n}{s^2}$
• $\frac{{\partial }^2\ln ℒ}{\partial {\theta }_2^2}(\bar{x},s^2)\stackrel{\text{above}}{=}{\frac{n}{2{\theta }_2^2}-\frac{1}{{\theta }_2^3}\sum _{i=1}^n(x_i-{\theta }_1{)}^2|}_{({\theta }_1,{\theta }_2)=(\bar{x},s^2)}=\frac{n}{2(s^2{)}^2}-\frac{1}{(s^2{)}^3}\cdot n{s}^2=\frac{n}{2(s^2{)}^2}-\frac{n}{(s^2{)}^2}=\frac{-n}{2(s^2{)}^2}$
• $\frac{{\partial }^2\ln ℒ}{\partial {\theta }_2\partial {\theta }_1}(\bar{x},s^2)\stackrel{\text{above}}{=}{\frac{\partial }{\partial {\theta }_2}(\frac{1}{{\theta }_2}\sum _{i=1}^n(x_i-{\theta }_1))|}_{({\theta }_1,{\theta }_2)=(\bar{x},s^2)}={-\frac{\sum _{i=1}^n(x_i-{\theta }_1)}{{\theta }_2^2}|}_{({\theta }_1,{\theta }_2)=(\bar{x},s^2)}=-\frac{\sum _{i=1}^n(x_i-\bar{x})}{(s^2{)}^2}=-\frac{\sum _{i=1}^n(x_i)-n\bar{x}}{(s^2{)}^2}=-\frac{n\bar{x}-n\bar{x}}{(s^2{)}^2}=0$

As a result, the determinant of the Hessian matrix is $\frac{-n}{s^2}\cdot \frac{-n}{2(s^2{)}^2}=\frac{n^2}{2(s^2{)}^3}$.

(b) From (a), the determinant of the Hessian matrix is positive. Also, $\frac{{\partial }^2\ln ℒ}{\partial {\theta }_1^2}(\bar{x},s^2)=-\frac{n}{s^2}<0$. Thus, by the second partial derivative test, $\ln ℒ(\theta ;𝐱)$ attains maximum at $({\theta }_1,{\theta }_2)=(\bar{x},s^2)$.

Exercise. Let $X_1,\dots ,X_n$ be a random sample from the exponential distribution with rate parameter $\lambda$, with pdf $f(x;\lambda )=\lambda e^{-\lambda x},x\ge 0$, where $\lambda >0$. Show that the MLE of $\lambda$ is $\frac{1}{\bar{X}}$.

Solution

Proof. The likelihood function is $ℒ(\lambda )=\prod _{i=1}^n(\lambda e^{-\lambda x_i})={\lambda }^n\exp (-\lambda \sum _{i=1}^n{x}_i)$. Thus, the log-likelihood function is $\ln ℒ(\lambda )=n\ln \lambda -\lambda \sum _{i=1}^n{x}_i$. Differentiating the log-likelihood function with respect to $\lambda$ gives $\frac{d}{d\lambda }\ln ℒ(\lambda )=\frac{n}{\lambda }-\sum _{i=1}^n{x}_i$. Setting the derivative to be zero, we get $\frac{n}{\lambda }-\sum _{i=1}^n{x}_i=0⟹\frac{n}{\lambda }-n\bar{x}=0⟹\frac{1}{\lambda }=\bar{x}⟹\lambda =\frac{1}{\bar{x}}$. It remains to verify that $\ln ℒ(\lambda )$ attains maximum at $\lambda =\frac{1}{\bar{x}}$. Since $\frac{d^2}{d{\lambda }^2}\ln ℒ(\lambda )=-\frac{n}{{\lambda }^2}<0$, this is verified. Hence, the MLE of $\lambda$ is $\frac{1}{\bar{X}}$.

$◻$

Example. (Application of maximum likelihood estimation) Suppose you are given a box which contains four balls, with unknown number of red and black balls. Now, you draw three balls out of the box, and find out that you get two red balls and one black ball. Using maximum likelihood estimation, estimate the number of red and black balls inside the box.

Solution: Given the color of the balls drawn, we know that the box contains at least two red balls and at least one black ball. This means the box contains either two red balls or three red balls. Let $r$ be the number of red balls inside the box. Then, the number of black balls inside the box is $4-r$. The possible values of parameter $r$ are 2 and 3.

Now, we compare the probability of getting such result from drawing three balls when $r=2$ and $r=3$.

• For $r=2$, the probability is $\frac{(\genfrac{}{}{0ex}{}{2}{2})(\genfrac{}{}{0ex}{}{2}{1})}{(\genfrac{}{}{0ex}{}{4}{3})}=\frac{1}{2}$ (consider the pmf of hypergeometric distribution).
• For $r=3$, the probability is $\frac{(\genfrac{}{}{0ex}{}{3}{2})(\genfrac{}{}{0ex}{}{1}{1})}{(\genfrac{}{}{0ex}{}{4}{3})}=\frac{3}{4}$.

Hence, the maximum likelihood estimate of $r$ is 3. Thus, the estimated number of the red balls is 3, and that of the black balls is 1.

Exercise. Suppose the box now contains 100 balls, with unknown number of red and black balls. Now, you draw 99 balls out of the box, and find out that you get 98 red balls and one black ball. Using maximum likelihood estimation, estimate the number of red and black balls inside the box.

Solution

Similarly, the box contains at least 98 red balls and one black ball. We use the same notation as in the above example. Then, the number of the black balls is $100-r$, and the possible values of parameter $r$ are 98 and 99.

• For $r=98$, the probability is $\frac{(\genfrac{}{}{0ex}{}{98}{98})(\genfrac{}{}{0ex}{}{2}{1})}{(\genfrac{}{}{0ex}{}{100}{99})}=0.02$
• For $r=99$, the probability is $\frac{(\genfrac{}{}{0ex}{}{99}{98})(\genfrac{}{}{0ex}{}{1}{1})}{(\genfrac{}{}{0ex}{}{100}{99})}=0.99$

Thus, the maximum likelihood estimate of $r$ is 99. Thus, the estimated number of the red balls is 99 and that of the black balls is 1.

Remark.

• The difference of the probabilities between two possible values of $r$ becomes much larger in this case.
• Intuitively, when you have such draw result, you will think that it is quite unlikely that the box has two black balls inside, i.e., the ball that is not drawn is actually black, and somehow you draw all red balls out, but not the black ball.

Method of moments estimator (MME)

For maximum likelihood estimation, we need to utilize the likelihood function, which is found from the joint pmf of pdf of the random sample from a distribution. However, we may not know exactly the pmf of pdf of the distribution in practice. Instead, we may just know some information about the distribution, e.g. mean, variance, and some moments ($r$th moment of a random variable $X$ is $𝔼[X^r]$, we denote it by ${\mu }_r$ for simplicity). Such moments often contain information about the unknown parameter. For example, for a normal distribution $𝒩(\mu ,{\sigma }^2)$, we know that $\mu ={\mu }_1$ and ${\sigma }^2={\mu }_2-({\mu }_1{)}^2$. Because of this, when we want to estimate the parameters, we can do this through estimating the moments.

Now, we would like to know how to estimate the moments. We let $m_r=\frac{\sum _{i=1}^n{X}_i^r}{n}$ be the $r$th sample moment , where $X_i$'s are independent and identically distributed. By weak law of large number (assuming the conditions are satisified), we have

• $\bar{X}=m_1\stackrel{p}{\to }𝔼[X]={\mu }_1$
• $m_2\stackrel{p}{\to }𝔼[X^2]={\mu }_2$ (this can be seen from replacing the "$X$" by "$X^2$" in the weak law of large number, then the conditions are still satisfied, and so we can still apply the weak law of large number)

In general, we have $m_r\stackrel{p}{\to }{\mu }_r$, since the conditions are still satisfied after replacing the "$X$" by "$X^r$" in the weak law of large number.

Because of these results, we can estimate the $r$-th moment ${\mu }_r$ using the $r$-th sample moment $m_r$, and the estimation is "better" when $n$ is large. For example, in the above normal distribution example, we can estimate $\mu$ by $m_1$ and ${\sigma }^2$ by $m_2-(m_1{)}^2$, and these estimators are actually called the method of moments estimator.

To be more precise, we have the following the definition of the method of moments:

Definition. (Method of moments) Let $X_1,\dots ,X_n$ be a random sample from a distribution with pdf or pmf $f(x;{\theta }_1,\dots ,{\theta }_k)$. Write $k$ moment(s), e.g. ${\mu }_1,\dots ,{\mu }_k$, as function(s) of ${\theta }_1,\dots ,{\theta }_k$: $g_1({\theta }_1,\dots ,{\theta }_k),\dots ,g_k({\theta }_1,\dots ,{\theta }_k)$ respectively. Then, the method of moments estimator (MME) of ${\theta }_1,\dots ,{\theta }_k$, ${\hat{\theta }}_1,\dots ,{\hat{\theta }}_k$ respectively, is given by the solution (in the form of ${\hat{\theta }}_1,\dots ,{\hat{\theta }}_k$ in terms of $m_1,\dots ,m_k$, corresponding to the $k$ moments ${\mu }_1,\dots ,{\mu }_k$) to the following system of equations: $$\{\begin{array}{c}{m}_1=g_1({\hat{\theta }}_1,\dots ,{\hat{\theta }}_k)\\ ⋮\\ m_k=g_k({\hat{\theta }}_1,\dots ,{\hat{\theta }}_k)\end{array}$$

Remark.

• When there are $k$ unknown parameters, we need to solve a system of $k$ equations, involving $k$ sample moments.
• Usually, we select the first $k$ moments for the $k$ moments, as in the definition. But this is not necessary, and we may choose other moments, including fractional moments (e.g. $𝔼[X^{1/2}]$, and we use $m_{1/2}$ in this case).

• Because of this, the method of moment estimator is not unique.

Example. Let $X_1,\dots ,X_n$ be a random sample from the normal distribution $𝒩(\mu ,{\sigma }^2)$. Find the MME of $\mu$ and ${\sigma }^2$.

Solution: First, there are two unknown parameters. Thus, we need to solve aa system of 2 equations, involving 2 sample moments and 2 moments. Since ${\mu }_1=\mu$ and ${\mu }_2={\sigma }^2+(\mu {)}^2$, consider the following system of equations: $$\{\begin{array}{cc}{m}_1=\hat{\mu }& (1)\\ m_2=\widehat{{\sigma }^2}+(\hat{\mu }{)}^2& (2)\end{array}$$ Substituting $(1)$ into $(2)$, we get $m_2=\widehat{{\sigma }^2}+(m_1{)}^2\iff \widehat{{\sigma }^2}=m_2-(m_1{)}^2$. Hence, the MME of $\mu$ is $\hat{\mu }=m_1$ and the MME of ${\sigma }^2$ is $\widehat{{\sigma }^2}=m_2-(m_1{)}^2$.

Remark.

• We can see that the process of finding the MME of $\mu$ and ${\sigma }^2$ is much easier than finding the MLE of $\mu$ and ${\sigma }^2$. This is because the expression of the first and second moment in terms of parameters is simple in this case. However, when the expression is more complicated, finding the MME of the parameters can be quite complicated.

Example. Let $X_1,\dots ,X_n$ be a random sample from the exponential distribution with rate parameter $\lambda$. Find the MME of $\lambda$ and compare it to the MLE of $\lambda$.

Solution: Since ${\mu }_1=\frac{1}{\lambda }$, consider the following equation: $m_1=\frac{1}{\hat{\lambda }}$. We then have $\hat{\lambda }=\frac{1}{m_1}$. Hence, the MME of $\lambda$ is $\hat{\lambda }=\frac{1}{m_1}=\frac{1}{\bar{X}}$, which is somehow the same as the MLE of $\lambda$.

Exercise. Let $X_1,\dots ,X_n$ be a random sample from the uniform distribution $𝒰[a,b]$. Show that the MMEs of $a$ and $b$ are $\hat{a}=m_1-\sqrt{3(m_2-m_1^2)}$ and $\hat{b}=m_1+\sqrt{3(m_2-m_1^2)}$ respectively.

Solution

Proof. Since ${\mu }_1=\frac{a+b}{2}$ and ${\mu }_2=\frac{(b-a{)}^2}{12}+{\left(\frac{a+b}{2}\right)}^2=\frac{b^2-2ab+a^2+3a^2+6ab+3b^2}{12}=\frac{4a^2+4b^2+4ab}{12}=\frac{a^2+b^2+ab}{3}$, consider the following system of equations: $$\{\begin{array}{cc}{m}_1=(\hat{a}+\hat{b})/2& (1)\\ m_2=({\hat{a}}^2+\hat{a}\hat{b}+{\hat{b}}^2)/3& (b)\end{array}$$ From $(1)$, we have $\hat{b}=2m_1-\hat{a}$. Substituting it into $(2)$, we have $m_2=({\hat{a}}^2+\hat{a}(2m_1-\hat{a})+(2m_1-\hat{a}{)}^2)/3-({\hat{a}}^2+2m_1\hat{a}-{\hat{a}}^2+4m_1^2-4m_1\hat{a}+{\hat{a}}^2)/3\iff {\hat{a}}^2-2m_1\hat{a}+4m_1^2=3m_2\iff {\hat{a}}^2-2m_1\hat{a}+4m_1^2-3m_2=0$ Solving this equation by quadratic formula, we get $\hat{a}=\frac{2m_1\pm \sqrt{4m_1^2-4(4m_1^2-3m_2)}}{2}=\frac{2m_1\pm \sqrt{12m_2-12m_1^2}}{2}=\frac{2m_1\pm 2\sqrt{3(m_2-m_1^2)}}{2}=m_1\pm \sqrt{3(m_2-m_1^2)}$.

When $\hat{a}=m_1+\sqrt{3(m_2-m_1^2)}$, $\hat{b}=m_1-\sqrt{3(m_2-m_1^2)}<\hat{a}$. However, from the definition of the uniform distribution, we need to have $\hat{a}<\hat{b}$, and thus this case is rejected.

When $\hat{a}=m_1-\sqrt{3(m_2-m_1^2)}$, $\hat{b}=m_1+\sqrt{3(m_2-m_1^2)}>\hat{a}$, which satisfies the definition of the uniform distribution.

Thus, we have the desired result.

$◻$

Properties of estimator

In this section, we will introduce some criteria for evaluating how "good" a point estimator is, namely unbiasedness, efficienecy and consistency.

Unbiasedness

For $\hat{\theta }$ to be a "good" estimator of a parameter $\theta$, a desirable property of $\hat{\theta }$ is that its expected value equals the value of the parameter $\theta$, or at least close to the value. Because of this, we introduce a value, namely bias, to measure how close is the mean of $\hat{\theta }$ to $\theta$.

Definition. (Bias) The bias of an estimator $\hat{\theta }$ is $\mathrm{Bias}(\hat{\theta })=𝔼[\hat{\theta }]-\theta .$

We will also define some terms related to bias.

Definition. ((Un)biased estimator) An estimator $\hat{\theta }$ is an unbiased estimator of a parameter $\theta$ if $\mathrm{Bias}(\hat{\theta })=0$. Otherwise, the estimator is called a biased estimator.

Definition. (Asymptotically unbiased estimator) An estimator $\hat{\theta }$ is an asymptotically unbiased estimator of a parameter $\theta$ if $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Bias}(\hat{\theta })=0$ where $n$ is the sample size.

Remark.

• An unbiased estimator must be an asymptotically unbiased estimator, but the converse is not true, i.e., an asymptotically unbiased estimator may not be an unbiased estimator. Thus, a biased estimator may be an asymptotically unbiased estimator.
• When we discuss the goodness of estimators in terms of unbiasedness, an unbiased estimator is better than an asymptotically unbiased estimator, which is better than an unbiased estimator.

• However, there are also other criteria for evaluating the goodness of estimators apart from unbiasedness, so when we also account for other criteria, a biased estimator may be somehow "better" than an unbiased estimator overall.

Example. Let $X_1,\dots ,X_n$ be a random sample from the Bernoulli distribution with success probability $p$. Show that the MLE of $p$, $\bar{X}$, is an unbiased estimator of $p$.

Proof. Since $𝔼[\bar{X}]=\frac{1}{n}\cdot 𝔼\left[\sum _{i=1}^n{X}_i\right]=\frac{1}{n}\sum _{i=1}^n𝔼[X_i]=\frac{1}{n}\cdot \sum _{i=1}^{n}p=\frac{np}{n}=p$, the result follows.

$◻$

Exercise. Suppose the Bernoulli distribution is replaced by binomial distribution with $n$ trials and success probability $p$. Show that $\bar{X}$ is a biased estimator of $p$. Modify this estimator such that it is an unbiased estimator of $p$.

Solution

Proof. Since $𝔼[\bar{X}]=\frac{1}{n}\sum _{i=1}^n𝔼[X_i]=\frac{1}{n}\sum _{i=1}^{n}np=np\ne p$, $\bar{X}$ is a biased estimator of $p$.

$◻$

We can modify this estimator to $\frac{\bar{X}}{n}$, and then its mean is $\frac{np}{n}=p$. Alternatively, we may choose the estimator to be $\frac{X_i}{n}$ ($i\in \{1,\dots ,n\}$), whose mean is also $p$ (Other estimators whose mean is $p$ are also fine).

Example. Let $X_1,\dots ,X_n$ be a random sample from the normal distribution $𝒩(\mu ,{\sigma }^2)$. Show that the MLE of $\mu$, $\bar{X}$, is an unbiased estimator of $\mu$, and the MLE of ${\sigma }^2$, $S^2$, is an asymptotically unbiased estimator of ${\sigma }^2$.

Proof. First, since $𝔼[\bar{X}]=\frac{1}{n}\sum _{i=1}^n𝔼[X_i]=\frac{1}{n}\sum _{i=1}^n\mu =\mu$, $\bar{X}$ is an unbiased estimator of $\mu$.

On the other hand, $$\begin{array}{cc}𝔼[S^2]& =\frac{1}{n}\sum _{i=1}^n𝔼[(X_i-\bar{X}{)}^2]\\ & =\frac{1}{n}\sum _{i=1}^n\mathrm{Var}(X_i-\bar{X})& \text{since }𝔼[X_i-\bar{X}]=𝔼[X_i]-𝔼[\bar{X}]=\mu -\mu =0\\ & =\frac{1}{n}\sum _{i=1}^n\mathrm{Var}(X_i-\frac{X_1+\cdots +X_{i-1}+X_i+\cdots +X_n}{n})\\ & =\frac{1}{n}\sum _{i=1}^n\mathrm{Var}(\frac{n{X}_i}{n}-\frac{X_1+\cdots +X_{i-1}+X_i+\cdots +X_n}{n})\\ & =\frac{1}{n}\sum _{i=1}^n\mathrm{Var}\left(\frac{X_1+\cdots +X_{i-1}+(n-1)X_i+\cdots +X_n}{n}\right)\\ & =\frac{1}{n}\sum _{i=1}^n\mathrm{Var}(\frac{(n-1)X_i}{n}+\frac{X_1+\cdots +X_{i-1}+X_{i+1}{X}_n}{n})\\ & =\frac{1}{n}\sum _{i=1}^n[\mathrm{Var}\left(\frac{(n-1)X_i}{n}\right)+\mathrm{Var}\left(\frac{X_1+\cdots +X_{i-1}+X_{i+1}+\cdots +X_n}{n}\right)]& \text{by independence}\\ & =\frac{1}{n}\sum _{i=1}^n[\mathrm{Var}\left(\frac{(n-1)X_i}{n}\right)+\mathrm{Var}\left(\frac{X_1+\cdots +X_{i-1}+X_{i+1}+\cdots +X_n}{n}\right)]\\ & =\frac{1}{n}\sum _{i=1}^n[\frac{(n-1{)}^2}{n^2}{\sigma }^2+\frac{1}{n^2}\mathrm{Var}(X_1+\cdots +X_{i-1}+X_{i+1}+\cdots +X_n)]\\ & =\frac{1}{n}\sum _{i=1}^n[\frac{(n-1{)}^2}{n^2}{\sigma }^2+\frac{n-1}{n^2}{\sigma }^2]& \text{by iid}\\ & =\frac{1}{n}\sum _{i=1}^n[\frac{{\sigma }^2}{n^2}(n^2-2n+1+n-1)]\\ & =\frac{1}{n}\sum _{i=1}^n\left[\frac{(n^2-n){\sigma }^2}{n^2}\right]\\ & =\frac{1}{n}\sum _{i=1}^n\left[\frac{(n-1){\sigma }^2}{n}\right]\\ & =\frac{1}{n}\cdot n\cdot \frac{(n-1){\sigma }^2}{n}\\ & =\frac{n-1}{n}{\sigma }^2\\ \end{array}$$ Thus, $\underset{n\to \mathrm{\infty }}{lim}𝔼[S^2]=\underset{n\to \mathrm{\infty }}{lim}\left(\frac{n-1}{n}{\sigma }^2\right)={\sigma }^2\underset{n\to \mathrm{\infty }}{lim}(1-\frac{1}{n})={\sigma }^2(1-\underset{n\to \mathrm{\infty }}{lim}\frac{1}{n})={\sigma }^2$, as desired.

$◻$

Exercise. Modify the estimator of ${\sigma }^2$ such that it becomes an unbiased estimator.

Solution

The estimator can be modified as $\frac{n}{n-1}{S}^2=\frac{\sum _{i=1}^n(X_i-\bar{X}{)}^2}{n-1}$.

Efficiency

We have discussed how to evaluate the unbiasedness of estimators. Now, if we are given two unbiased estimators, $\hat{\theta }$ and $\tilde{\theta }$, how should we compare their goodness? Their goodness is the same if we are only comparing them in terms of unbiasedness. Therefore, we need another criterion in this case. One possible way is to compare their variances, and the one with smaller variance is better, since on average, the estimator is less deviated from its mean, which is the value of the unknown parameter by the definition of unbiased estimator, and thus the one with smaller variance is more accurate in some deviation sense. Indeed, an unbiased estimator can still have a large variance, and thus deviate a lot from its mean. Such estimator is unbiased since the positive deviations and negative deviations somehow cancel out each other. This is the idea of efficiency.

Definition. (Efficiency) Suppose $\hat{\theta }$ and $\tilde{\theta }$ are two unbiased estimators of an unknown parameter $\theta$. The efficiency of $\hat{\theta }$ relative to $\tilde{\theta }$ is $\mathrm{Eff}(\hat{\theta },\tilde{\theta })=\frac{\mathrm{Var}(\tilde{\theta })}{\mathrm{Var}(\hat{\theta })}$. If $\mathrm{Eff}(\hat{\theta },\tilde{\theta })>1$, then we say that $\hat{\theta }$ is relatively more efficient than $\tilde{\theta }$.

Remark.

• Since $\mathrm{Eff}(\hat{\theta },\tilde{\theta })>1\iff \mathrm{Var}(\hat{\theta })<\mathrm{Var}(\tilde{\theta })$, the estimator with smaller variance is relatively more efficient than the estimator with larger variance.
• Normally, the variance should be nonzero, and thus the efficiency should be defined in normal cases.
• Sometimes, it is also called relative efficiency due to the fact that the efficiency describes $\mathrm{Var}(\tilde{\theta })$ equals "how many" $\mathrm{Var}(\hat{\theta })$.
• One may ask that why we use the ratio of variances is used in the definition to compare variances, instead of using the difference in variances. A possible reason is that the ratio of variances does not have any unit (the unit of the variances (if exists) cancels out each other), but the difference in variances can have an unit. Also, using the ratio of variances allows us to also compare different efficiencies numerically, calculated from different variances.

Actually, for the variance of unbiased estimator, since the mean of the unbiased estimator is the unknown paramter $\theta$, it measures the mean of the squared deviation from $\theta$, and we have a specific term for this deviation, namely mean squared error (MSE).

Definition. (Mean squared error) Suppose $\hat{\theta }$ is an estimator of a parameter $\theta$. The mean squared error (MSE) of $\hat{\theta }$ is $\mathrm{MSE}(\hat{\theta })=𝔼[(\hat{\theta }-\theta {)}^2]$.

Remark.

• From this definition, $\mathrm{MSE}(\hat{\theta })$ is the mean value of the square of error $\hat{\theta }-\theta$, and hence the name mean squared error.

Notice that in the definition of MSE, we do not specify that $\hat{\theta }$ to be an unbiased estimator. Thus, $\hat{\theta }$ in the definition may be biased. We have mentioned that when $\hat{\theta }$ is unbiased, then its variance is actually its MSE. In the following, we will give a more general relationship between $\mathrm{MSE}(\hat{\theta })$ and $\mathrm{Var}(\hat{\theta })$, not just for unbiased estimators.

Proposition. (Relationship between mean squared error and variance) If $\mathrm{Var}(\hat{\theta })$ exists, then $\mathrm{MSE}(\hat{\theta })=\mathrm{Var}(\hat{\theta })+[\mathrm{Bias}(\hat{\theta }){]}^2$.

Proof. By definition, we have $\mathrm{MSE}(\hat{\theta })=𝔼[(\hat{\theta }-\theta {)}^2]$ and $\mathrm{Var}(\hat{\theta })=𝔼[(\hat{\theta }-𝔼[\hat{\theta }]{)}^2]$. From these, we are motivated to write $$\begin{array}{cc}\mathrm{MSE}(\hat{\theta })& =𝔼[(\hat{\theta }-\theta {)}^2]\\ & =𝔼\left[\right((\hat{\theta }-𝔼[\hat{\theta }])+(𝔼[\hat{\theta }]-\theta ){)}^2]\\ & =𝔼[(\hat{\theta }-𝔼[\hat{\theta }]{)}^2+2(\hat{\theta }-𝔼[\hat{\theta }])\underset{\text{constant}}{\underbrace{(𝔼[\hat{\theta }]-\theta )}}+(𝔼[\hat{\theta }]-\theta {)}^2]\\ & =\mathrm{Var}(\hat{\theta })+2(𝔼[\hat{\theta }]-\theta )\underset{=𝔼[\hat{\theta }]-𝔼[\hat{\theta }]=0}{\underbrace{𝔼[\hat{\theta }-𝔼[\hat{\theta }]]}}+[\mathrm{Bias}(\hat{\theta }){]}^2\\ & =\mathrm{Var}(\hat{\theta })+[\mathrm{Bias}(\hat{\theta }){]}^2,\end{array}$$ as desired.

$◻$

Example. Let $X_1,\dots ,X_n$ ($n>1$) be a random sample from $𝒩(\mu ,{\sigma }^2)$.

(a) Show that the single observation estimator $X_1$ is an unbiased estimator for $\mu$.

(b) Calculate the MSE of $X_1$ and $\bar{X}$ respectively.

(c) Which of $X_1$ and $\bar{X}$ is a better estimator of $\mu$ in terms of unbiasedness and efficiency?

Solution:

(a) Since $𝔼[X_1]=\mu$, the result follows.

(b) $\mathrm{MSE}(X_1)=\mathrm{Var}(X_1)+0^2={\sigma }^2$, and $\mathrm{MSE}(\bar{X})=\mathrm{Var}(\bar{X})=\frac{1}{n^2}\sum _{i=1}^n\mathrm{Var}(X_i)=\frac{n{\sigma }^2}{n^2}=\frac{{\sigma }^2}{n}$.

(c) Since $\mathrm{MSE}(\bar{X})<\mathrm{MSE}(X_1)\iff \mathrm{Var}(\bar{X})<\mathrm{Var}(X_1)$, $\bar{X}$ is relatively more efficient than $X_1$. Since both $X_1$ and $\bar{X}$ are unbiased estimators of $\mu$, we conclude that $\bar{X}$ is a better estimator of $\mu$ in terms of unbiasedness and efficiency.

Exercise. In addition to the random sample with sample size $n$ in the example, suppose we take another random sample with sample size $m$. Let ${\bar{X}}^{(n)}$ and ${\bar{X}}^{(m)}$ denote the sample mean for the sample with sample size $n$ and $m$ respectively.

(a) Calculate $\mathrm{Eff}({\bar{X}}^{(n)},{\bar{X}}^{(m)})$.

(b) State the condition on the sample sizes $m$ and $n$ under which ${\bar{X}}^{(m)}$ is relatively more efficient than ${\bar{X}}^{(n)}$.

Solution

(a) Since $\mathrm{Var}\left({\bar{X}}^{(n)}\right)=\frac{{\sigma }^2}{n}$ (from example), and $\mathrm{Var}\left({\bar{X}}^{(m)}\right)=\frac{{\sigma }^2}{m}$ (by similar arguments as in the example), $$\mathrm{Eff}({\bar{X}}^{(n)},{\bar{X}}^{(m)})=\frac{{\sigma }^2/m}{{\sigma }^2/n}=\frac{n}{m}$$.

(b) Since $\frac{n}{m}>1\iff n>m$, the condition is $n>m$.

Remark.

• This shows that the sample mean with a larger sample size is relatively more efficient than the one with smaller sample size.

Proposition. $\underset{n\to \mathrm{\infty }}{lim}\mathrm{MSE}(\hat{\theta })=0$ if and only if $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Var}(\hat{\theta })=0$ and $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Bias}(\hat{\theta })=0$.

Proof.

• "if" part is simple. Assume $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Var}(\hat{\theta })=0$ and $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Bias}(\hat{\theta })=0$. Then, $\underset{n\to \mathrm{\infty }}{lim}(\mathrm{Var}(\hat{\theta })+(\mathrm{Bias}(\hat{\theta }){)}^2)=0\Rightarrow \underset{n\to \mathrm{\infty }}{lim}\mathrm{MSE}(\hat{\theta })=0$.
• "only if" part: we can use proof by contrapositive, i.e., proving that if $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Var}(\hat{\theta })\ne 0$ or $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Bias}(\hat{\theta })=0$, then $\underset{n\to \mathrm{\infty }}{lim}\mathrm{MSE}(\hat{\theta })\ne 0$.

• Case 1: when $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Var}(\hat{\theta })\ne 0$, it means $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Var}(\hat{\theta })>0$ since the variance is nonnegative. Also, $\underset{n\to \mathrm{\infty }}{lim}(\mathrm{Bias}(\hat{\theta }){)}^2\ge 0$. It follows that $\underset{n\to \mathrm{\infty }}{lim}\mathrm{MSE}(\hat{\theta })>0$, i.e., the MSE does not equal zero.
• Case 2: when $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Bias}(\hat{\theta })\ne 0$, it means $\underset{n\to \mathrm{\infty }}{lim}(\mathrm{Bias}(\hat{\theta }){)}^2>0$. Also, $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Var}(\hat{\theta })\ge 0$. It follows that $\underset{n\to \mathrm{\infty }}{lim}\mathrm{MSE}(\hat{\theta })>0$, i.e., the MSE does not equal zero.

$◻$

Remark.

• As a result, if we know that $\underset{n\to \mathrm{\infty }}{lim}\mathrm{MSE}(\hat{\theta })=0$, then we know that $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Bias}(\hat{\theta })=0$, i.e., $\hat{\theta }$ is an asymptotically unbiased estimator (in addition to $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Var}(\hat{\theta })=0$) ($\hat{\theta }$ may be an unbiased estimator).

Uniformly minimum-variance unbiased estimator

Now, we know that the smaller the variance of an unbiased estimator, the more efficient (and "better") it is. Thus, it is natural that we want to know what is the most efficient (i.e., the "best") unbiased estimator, i.e., the unbiased estimator with the smallest variance. We have a specific name for such unbiased estimator, namely uniformly minimum-variance unbiased estimator (UMVUE) . To be more precise, we have the following definition for UMVUE:

Definition. (Uniformly minimum-variance unbiased estimator) The uniformly minimum-variance unbiased estimator (UMVUE) is an unbiased estimator with the smallest variance among all unbiased estimators.

Indeed, UMVUE is unique, i.e., there is exactly one unbiased estimator with the smallest variance among all unbiased estimators, and we will prove it in the following.

Proposition. (Uniqueness of UMVUE) If $W$ is an UMVUE of a function of a parameter $\tau (\theta )$, then $W$ is unique.

Proof. Assume that $W$ is an UMVUE of $\tau (\theta )$, and $W^{\prime }$ is another UMVUE of $\tau (\theta )$. Define the estimator $W^{\ast }=\frac{1}{2}(W+W^{\prime })$. Since $𝔼[W^{\ast }]=\frac{1}{2}(𝔼[W]+𝔼[W^{\prime }])=\frac{1}{2}(\tau (\theta +\theta )=\tau (\theta )$, $W^{\ast }$ is an unbiased estimator of $\tau (\theta )$.

Now, we consider the variance of $W^{\ast }$. $$\begin{array}{cc}\mathrm{Var}(W^{\ast })& =\frac{1}{4}\mathrm{Var}(W+W^{\prime })\\ & =\frac{1}{4}[\mathrm{Var}(W)+\mathrm{Var}(W^{\prime })+2\mathrm{Cov}(W,W^{\prime })]\\ & \le \frac{1}{4}\mathrm{Var}(W)+\frac{1}{4}\mathrm{Var}(W^{\prime })+\frac{1}{2}\sqrt{\mathrm{Var}(W)\mathrm{Var}(W^{\prime })}& (\text{covariance inequality})\\ & =\frac{1}{4}\mathrm{Var}(W)+\frac{1}{4}\mathrm{Var}(W)+\frac{1}{2}\sqrt{(\mathrm{Var}(W){)}^2}& (\mathrm{Var}(W)=\mathrm{Var}(W^{\prime })\text{ since }W\text{ and }{W}^{\prime }\text{ are both UMVUE})\\ & =\frac{1}{2}\mathrm{Var}(W)+\frac{1}{2}\mathrm{Var}(W)& (\mathrm{Var}(W)>0)\\ & =\mathrm{Var}(W).\end{array}$$ Thus, we now have either $\mathrm{Var}(W^{\ast })<\mathrm{Var}(W)$ or $\mathrm{Var}(W^{\ast })=\mathrm{Var}(W)$. If the former is true, then $W$ is not an UMVUE of $\tau (\theta )$ by definition, since we can find another unbiased estimator, namely $W^{\ast }$, with smaller variance than it. Hence, we must have the latter, i.e., $$\mathrm{Var}(W^{\ast })=\mathrm{Var}(W).$$ This implies when we apply the covariance inequality, the equality holds, i.e., $$\mathrm{Cov}(W,W^{\prime })=\sqrt{\mathrm{Var}(W)\mathrm{Var}(W^{\prime })}⟺\rho (W^{\prime },W)=1,$$ which means $W^{\prime }$ is increasing linearly with $W$, i.e., we can write $$W^{\prime }=aW+b$$ for some constants $a>0$ and $b$.

Now, we consider the covariance $\mathrm{Cov}(W,W^{\prime })$. $$\mathrm{Cov}(W,W^{\prime })\stackrel{\text{ above }}{=}\mathrm{Cov}(W,aW+b)\stackrel{\text{ properties }}{=}a\mathrm{Cov}(W,W)\stackrel{\text{ property }}{=}a\mathrm{Var}(W).$$ On the other hand, since the equality holds in the covariance inequality, and $\mathrm{Var}(W)=\mathrm{Var}(W^{\prime })$ (since they are both UMVUE), $$\mathrm{Cov}(W,W^{\prime })=\sqrt{\mathrm{Var}(W)\mathrm{Var}(W^{\prime })}=\sqrt{(\mathrm{Var}(W){)}^2}=\mathrm{Var}(W).$$ Thus, we have $a=1$.

It remains to show that $b=0$ to prove that $W=W^{\prime }$, and therefore conclude that $W$ is unique.

From above, we currently have $W^{\prime }=W+b⟹𝔼[W^{\prime }]=𝔼[W]+b⟹\tau (\theta )=\tau (\theta )+b⟹b=0$, as desired.

$◻$

Remark.

• Thus, when we are able to find an UMVUE, then it is the unique one, and the variance every other possible unbiased estimator is strictly greater than the variance of the UMVUE.

Cramer-Rao lower bound

Without using some results, it is quite difficult to determine the UMVUE, since there are many (perhaps even infinitely many) possible unbiased estimator, so it is quite hard to ensure that one particular unbiased estimator is relative more efficient than every other possible unbiased estimators.

Therefore, we will introduce some approaches that help us to find the UMVUE. For the first approach, we find a lower bound on the variances of all possible unbiased estimators. After getting such lower bound, if we can find an unbiased estimator with variance to be exactly equal to the lower bound, then the lower bound is the minimum value of the variances, and hence such unbiased estimator is an UMVUE by definition.

Remark.

• There are many possible lower bounds, but when the lower bound is greater, it is closer to the actual minimum value of the variances, and hence "better".
• An unbiased estimator can still be an UMVUE even if its variance does not achieve the lower bound.

A common way to find such lower bound is to use the Cramer-Rao lower bound (CRLB), and we get the CRLB through Cramer-Rao inequality. Before stating the inequality, let us define some related terms.

Definition. (Fisher information) The Fisher information about a parameter $\theta$ with sample size $n$ is $${ℐ}_n(\theta )=𝔼\left[{\left(\frac{\partial \ln ℒ(𝜽;𝐗)}{\partial \theta }\right)}^2\right]$$ where $\ln ℒ(𝜽;𝐗)=\ln ℒ(𝜽;X_1,\dots ,X_n)$ is the log-likelihood function (as a random variable).

Remark.

• $\frac{\partial \ln ℒ(𝜽;𝐗)}{\partial \theta }$ is called the score function, and is denoted by $S(\theta ;𝐗)$.
• The "$𝜽$" may or may not be a parameter vector. If it is just a single parameter (usually the case here), then it is the same as "$\theta$". We use "$𝜽$" instead of "$\theta$" to emphasize that the "$\theta$" in ${ℐ}_n(\theta )$ and $S(\theta ;𝐗)$ is referring to the "$\theta$" in "$\frac{\partial }{\partial \theta }$"
• It is possible to define "Fisher information about a parameter vector", but in this case the Fisher information takes the form of a matrix instead of a single number, and it is called Fisher information matrix. However, since it is more complicated, we will not discuss it here.
• Since the expected value of the score function

$$𝔼[S(\theta ;𝐗)]𝔼\left[\frac{\partial \ln ℒ(𝜽;𝐗)}{\partial \theta }\right]={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\cdots {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{\partial \ln ℒ(𝜽;𝐱)}{\partial \theta }\cdot ℒ(𝜽;𝐱)d{x}_n\cdots d{x}_1={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\cdots {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{\frac{\partial ℒ(𝜽;𝐱)}{\partial \theta }}{ℒ(𝜽;𝐱)}\cdot ℒ(𝜽;𝐱)d{x}_n\cdots d{x}_1={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\cdots {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{\partial ℒ(𝜽;𝐱)}{\partial \theta }d{x}_n\cdots d{x}_1,$$

and under some regularity conditions which allow interchange of derivative and integral, this equals $\frac{\partial }{\partial \theta }{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\cdots {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}ℒ(𝜽;𝐱)d{x}_n\cdots d{x}_1=\frac{\partial }{\partial \theta }(1)=0$, the Fisher information about $\theta$ is also the variance of the score function, i.e., $\mathrm{Var}(S(\theta ;𝐗))=\mathrm{Var}\left(\frac{\partial \ln ℒ(𝜽;𝐗)}{\partial \theta }\right)$.

For the regularity conditions which allow interchange of derivative and integral, they include

1. the partial derivatives involved should exist, i.e., the (natural log) of the functions involved is differentiable
2. the integrals involved should be differentiable
3. the support does not depend on the parameter(s) involved

We have some results that assist us to compute the Fisher information.

Proposition. Let $X_1,\dots ,X_n$ be a random sample from a distribution with pdf or pmf $f$. Also, let $ℐ(\theta )=𝔼\left[{\left(\frac{\partial \ln f(X;𝜽)}{\partial \theta }\right)}^2\right]$, the Fisher information about $\theta$ with sample size one. Then, under some regularity conditions which allow interchange of derivative and integral, ${ℐ}_n(\theta )=nℐ(\theta )$.

Proof. $$\begin{array}{cc}{ℐ}_n(\theta )& =𝔼\left[{\left(\frac{\partial \ln ℒ(𝜽;𝐱)}{\partial \theta }\right)}^2\right]\\ & =\mathrm{Var}\left(\frac{\partial \ln ℒ(𝜽;𝐱)}{\partial \theta }\right)& \text{by above remark}\\ & =\mathrm{Var}\left(\frac{\partial }{\partial \theta }(\ln \prod _{i=1}^{n}f(X_i;𝜽))\right)& (ℒ(\theta ;𝐱)=\prod _{i=1}^{n}f(x_i;\theta ))\\ & =\mathrm{Var}\left(\frac{\partial }{\partial \theta }(\sum _{i=1}^n\ln f(X_i;𝜽))\right)\\ & =\mathrm{Var}(\sum _{i=1}^n\frac{\partial }{\partial \theta }\ln f(X_i;𝜽))& \text{by linearity of differentiation}\\ & =\sum _{i=1}^n\mathrm{Var}(\frac{\partial }{\partial \theta }\ln f(X_i;𝜽))& \text{by independence}\\ & =n\mathrm{Var}(\frac{\partial }{\partial \theta }\ln f(X_i;𝜽))& \text{by identically distributed property}\\ & =n𝔼\left[{\left(\frac{\partial \ln f(X;𝜽)}{\partial \theta }\right)}^2\right]& \text{by above remark}\\ & =nℐ(\theta ).\end{array}$$

$◻$

Proposition. Under some regularity conditions which allow interchange of derivative and integral, $ℐ(\theta )=-𝔼\left[\frac{{\partial }^2\ln f(X;𝜽)}{\partial {\theta }^2}\right]$.

Proof. $$\begin{array}{cc}𝔼\left[\frac{{\partial }^2\ln f(X;𝜽)}{\partial {\theta }^2}\right]& =𝔼\left[\frac{\partial }{\partial \theta }\left(\frac{\partial \ln f(X;𝜽)}{\partial \theta }\right)\right]\\ & =𝔼\left[\frac{\partial }{\partial \theta }(\frac{1}{f(X;𝜽)}\cdot \frac{\partial f(X;𝜽)}{\partial \theta })\right]\\ & =𝔼[\frac{1}{f(X;𝜽)}\cdot \frac{{\partial }^{2}f(X;𝜽)}{\partial {\theta }^2}-\frac{\partial f(X;𝜽)}{\partial \theta }\cdot \frac{1}{(f(X;𝜽){)}^2}\cdot \frac{\partial f(X;𝜽)}{\partial \theta }]\\ & =𝔼[\frac{1}{f(X;𝜽)}\cdot \frac{{\partial }^{2}f(X;𝜽)}{\partial {\theta }^2}-{\left(\frac{\partial f(X;𝜽)}{\partial \theta }\right)}^2\cdot \frac{1}{(f(X;𝜽){)}^2}]\\ & =𝔼[\frac{1}{f(X;𝜽)}\cdot \frac{{\partial }^{2}f(X;𝜽)}{\partial {\theta }^2}]-𝔼\left[{\left(\frac{\partial \ln f(X;𝜽)}{\partial \theta }\right)}^2\right]\\ & =𝔼[\frac{1}{f(X;𝜽)}\cdot \frac{{\partial }^{2}f(X;𝜽)}{\partial {\theta }^2}]-ℐ(\theta )\\ \end{array}$$ Now, it suffices to prove that $𝔼[\frac{1}{f(X;𝜽)}\cdot \frac{{\partial }^{2}f(X;𝜽)}{\partial {\theta }^2}]=0$, which is true since $$\begin{array}{cc}𝔼[\frac{1}{f(X;𝜽)}\cdot \frac{{\partial }^{2}f(X;𝜽)}{\partial {\theta }^2}]& ={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{1}{f(x;𝜽)}\cdot \frac{{\partial }^{2}f(x;𝜽)}{\partial {\theta }^2}\cdot f(x;𝜽)dx\\ & ={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{{\partial }^{2}f(x;𝜽)}{\partial {\theta }^2}dx\\ & =\frac{{\partial }^2}{\partial {\theta }^2}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}f(x;𝜽)dx\\ & =\frac{{\partial }^2}{\partial {\theta }^2}(1)\\ & =0.\\ \end{array}$$

$◻$

Remark.

• This proposition can be quite useful, since after partially differentiating $\ln f(X;𝜽)$, it is likely that many $X$'s will vanish, and thus the computation of the expectation will be easier.

Theorem. (Cramer-Rao inequality) Let $X_1,\dots ,X_n$ be a random sample from a distribution, and let $W$ be an unbiased estimator of $\tau (\theta )$ (a function of $\theta$). Then, under some regularity conditions which allow interchange of derivative and integral, $\mathrm{Var}(W)\ge \frac{({\tau }^{\prime }(\theta ){)}^2}{{ℐ}_n(\theta )}$.

Proof. Since $W$ is an unbiased estimator of $\tau (\theta )$, we have by definition $𝔼[W]=\tau (\theta )$. By definition of expectation, we have $𝔼[W]={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\cdots {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}wℒ(\theta ;𝐱)d{x}_n\cdots d{x}_1$ where $ℒ(\theta ;𝐱)$ is the likelihood function. Thus, $$\begin{array}{cccc}& & {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\cdots {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}wℒ(\theta ;𝐱)d{x}_n\cdots d{x}_1& =\tau (\theta )\\ & \Rightarrow & \frac{\partial }{\partial \theta }{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\cdots {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}wℒ(\theta ;𝐱)d{x}_n\cdots d{x}_1& =\frac{\partial }{\partial \theta }\tau (\theta )\\ & \Rightarrow & {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\cdots {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{\partial }{\partial \theta }(wℒ(\theta ;𝐱))d{x}_n\cdots d{x}_1& ={\tau }^{\prime }(\theta )\\ & \Rightarrow & {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\cdots {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}w\frac{\partial }{\partial \theta }(ℒ(\theta ;𝐱))\cdot \frac{1}{ℒ(\theta ;𝐱)}\cdot ℒ(\theta ;𝐱)d{x}_n\cdots d{x}_1& ={\tau }^{\prime }(\theta )\\ & \Rightarrow & {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\cdots {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}w\frac{\partial \ln ℒ(\theta ;𝐱)}{\partial \theta }ℒ(\theta ;𝐱)d{x}_n\cdots d{x}_1& ={\tau }^{\prime }(\theta )\\ & \Rightarrow & 𝔼[W\cdot \frac{\partial \ln ℒ(\theta ;𝐱)}{\partial \theta }]& ={\tau }^{\prime }(\theta )\\ & \Rightarrow & 𝔼[WS(\theta ;𝐗)]& ={\tau }^{\prime }(\theta )& (S(\theta ;𝐗)=\frac{\partial \ln ℒ(\theta ;𝐱)}{\partial \theta })\\ & \Rightarrow & 𝔼[WS(\theta ;𝐗)]-𝔼[W]\underset{=0}{\underbrace{𝔼[S(\theta ;𝐗)]}}& ={\tau }^{\prime }(\theta )& (𝔼[S(\theta ;𝐗)]=0\text{ by remark about Fisher information})\\ & \Rightarrow & \mathrm{Cov}(W,S(\theta ;𝐗))& ={\tau }^{\prime }(\theta )\\ \end{array}$$ Consider the covariance inequality: $(\mathrm{Cov}(X,Y){)}^2\le \mathrm{Var}(X)\mathrm{Var}(Y)$. We have $$(\mathrm{Cov}(W,S(\theta ;𝐗)){)}^2\le \mathrm{Var}(W)\mathrm{Var}(S(\theta ;𝐗))⟹({\tau }^{\prime }(\theta ){)}^2\le \mathrm{Var}(W)\mathrm{Var}(S(\theta ;𝐗))⟹\mathrm{Var}(W)\ge \frac{({\tau }^{\prime }(\theta ){)}^2}{\mathrm{Var}(S(\theta ;𝐗))}=\frac{({\tau }^{\prime }(\theta ){)}^2}{{ℐ}_n(\theta )}.$$ (${ℐ}_n(\theta )=\mathrm{Var}(S(\theta ;𝐗))$ by remark about Fisher information)

$◻$

Remark.

• $\frac{({\tau }^{\prime }(\theta ){)}^2}{{ℐ}_n(\theta )}$ is called the Cramer-Rao lower bound (CRLB).
• When $\tau (\theta )=\theta$, meaning that $W$ is an unbiased estimator of $\theta$, since $({\tau }^{\prime }(\theta ){)}^2=1^2=1$, the CRLB becomes $\frac{1}{{ℐ}_n(\theta )}$.

Example. Let $X_1,\dots ,X_n$ be a random sample from the normal distribution $𝒩(\mu ,{\sigma }^2)$. Show that the MLE of $\mu$, $\bar{X}$, is an UMVUE of $\mu$.

Proof. First, we can see that the regularity conditions are satisfied in this case. So, we may consider the CRLB of $\mu$ as follows. Since $$\ln f(X;\mu ,{\sigma }^2)=\ln \frac{1}{\sqrt{2\pi {\sigma }^2}}\exp (-\frac{(X-\mu {)}^2}{2{\sigma }^2})=-\frac{1}{2}(\ln 2\pi {\sigma }^2)-\frac{(X-\mu {)}^2}{2{\sigma }^2}=-\frac{1}{2}(\ln 2\pi {\sigma }^2)-\frac{X^2-2\mu X+{\mu }^2}{2{\sigma }^2},$$ we have $$ℐ(\mu )=-𝔼\left[\frac{{\partial }^2}{\partial {\mu }^2}(-\frac{1}{2}(\ln 2\pi {\sigma }^2)-\frac{X^2-2\mu X+{\mu }^2}{2{\sigma }^2})\right]=-𝔼\left[\frac{\partial }{\partial \mu }(-\frac{-2X+2\mu }{2{\sigma }^2})\right]=-𝔼[-\left(\frac{2}{2{\sigma }^2}\right)]=-𝔼[\underset{\text{constant wrt }X}{\underbrace{-\frac{1}{{\sigma }^2}}}]=\frac{1}{{\sigma }^2}.$$ Thus, the CRLB of $\mu$ is $\frac{1}{nℐ(\mu )}=\frac{1}{n(1/{\sigma }^2)}=\frac{{\sigma }^2}{n}$.

On the other hand, the variance of $\bar{X}$ is $\mathrm{Var}(\bar{X})=\frac{{\sigma }^2}{n}$ (which is shown in a previous example), which equals the CRLB of $\mu$. It follows that $\bar{X}$ is an UMVUE of $\mu$.

$◻$

Exercise. A student claims that $\frac{X_1}{\sqrt{n}}$ is another UMVUE of $\mu$, since $\mathrm{Var}\left(\frac{X_1}{\sqrt{n}}\right)=\frac{{\sigma }^2}{(\sqrt{n}{)}^2}=\frac{{\sigma }^2}{n}$, which equals the CRLB of $\mu$ as well. Is the claim correct? Why?

Solution

Recall that UMVUE is an unbiased estimator.

The claim is wrong, since $\frac{X_1}{\sqrt{n}}$ is not an unbiased estimator in general. This is because $𝔼\left[\frac{X_1}{\sqrt{n}}\right]=\frac{\mu }{\sqrt{n}}\ne \mu$ unless $n=1$. But if $n=1$, then this estimator is simply $X_1$, which is exactly the same as $\bar{X}=\frac{X_1}{1}=X_1$. So, the estimator is not another UMVUE in this case.

Sometimes, we cannot use the CRLB method for finding UMVUE, because

• the regularity conditions may not be satisfied, and thus we cannot use the Cramer-Rao inequality, and
• the variance of the unbiased estimator may not be equal to the CRLB, but we cannot conclude that it is not an UMVUE, because it may be the case that the CRLB is not attainable at all, and the smallest variance among all unbiased estimators is actually the variance of that estimator, which is larger than the CRLB.

We will illustrate some examples for these two cases in the following.

Example. Let $X_1,\dots ,X_n$ be a random sample from the uniform distribution $𝒰[0,\beta ]$. If we want to find the UMVUE of $\beta$, we cannot use the Cramer-Rao inequality to find it, since the support $[0,\beta ]$ depends on the parameter $\beta$.

Example. Let $X_1,\dots ,X_n$ be a random sample from the normal distribution $𝒩(\mu ,{\sigma }^2)$. It is given that in this case, $\frac{n{S}^2}{{\sigma }^2}\sim {\chi }_{n-1}^2$, where ${\chi }_k^2$ is the chi-squared distribution with $k$ degrees of freedom, and its variance is $2k$. Calculate $\mathrm{Var}(\frac{n}{n-1}\cdot S^2)$, and the CRLB of ${\sigma }^2$.

Solution: By the given information, we have $$\mathrm{Var}\left(\frac{n{S}^2}{{\sigma }^2}\right)=2(n-1)⟹\frac{n^2}{{\sigma }^4}\mathrm{Var}(S^2)=2(n-1)⟹\mathrm{Var}(S^2)=\frac{2(n-1){\sigma }^4}{n^2}.$$ Hence, $\mathrm{Var}(\frac{n}{n-1}\cdot S^2)=\frac{n^2}{(n-1{)}^2}\cdot \frac{2(n-1){\sigma }^4}{n^2}=\frac{2{\sigma }^4}{n-1}$.

On the other hand, since $$\begin{array}{cc}ℐ(\theta )& =-𝔼\left[\frac{{\partial }^2}{\partial (\sigma {)}^2}(-\frac{1}{2}(\ln 2\pi {\sigma }^2)-\frac{(X-\mu {)}^2}{2{\sigma }^2})\right]\\ & =-𝔼\left[\frac{\partial }{\partial \sigma }(-\frac{4\pi \sigma }{2(2\pi {\sigma }^2}+\frac{(X-\mu {)}^2}{{\sigma }^3})\right]\\ & =-𝔼\left[\frac{\partial }{\partial \sigma }(-\frac{1}{\sigma }+\frac{(X-\mu {)}^2}{{\sigma }^3})\right]\\ & =-𝔼[\frac{1}{{\sigma }^2}-\frac{3(X-\mu {)}^2}{{\sigma }^4}]\\ & =-\frac{1}{{\sigma }^2}+\frac{3}{{\sigma }^4}𝔼[(X-\mu {)}^2]\\ & =-\frac{1}{{\sigma }^2}+\frac{3}{{\sigma }^4}\cdot {\sigma }^2\\ & =\frac{2}{{\sigma }^2},\\ \end{array}$$ and ${\left(\frac{d}{d\sigma }{\sigma }^2\right)}^2=(2\sigma {)}^2=4{\sigma }^2$, the CRLB of ${\sigma }^2$ is $$\frac{4{\sigma }^2}{2n/{\sigma }^2}=\frac{2{\sigma }^4}{n}.$$

Remark.

• $\frac{n}{n-1}\cdot S^2$ is an unbiased estimator of ${\sigma }^2$, since $𝔼[\frac{n}{n-1}\cdot S^2]=\frac{n}{n-1}\left(\frac{n-1}{n}{\sigma }^2\right)={\sigma }^2$.
• We can observe that $\mathrm{Var}(\frac{n}{n-1}\cdot S^2)$ is greater than CRLB. But does this mean $\frac{n}{n-1}\cdot S^2$ is not the UMVUE of ${\sigma }^2$? We do not know, since we are not sure that whether there is another unbiased estimators with variance less than $\mathrm{Var}(\frac{n}{n-1}\cdot S^2)$, and it is possible that the CRLB is not attainable.

Since the CRLB is sometimes attainable and sometimes not, it is natural to question that when can the CRLB be attained. In other words, we would like to know the attainment conditions for the CRLB, which are stated in the following corollary.

Corollary. (Attainable condition for the CRLB) Let $X_1,\dots ,X_n$ be a random sample from a distribution, and let $W$ be an unbiased estimator of $\tau (\theta )$. Suppose the regularity conditions in the Cramer-Rao inequality are satisfied. Then, the CRLB can be attained, i.e., there exists some $W$ such that $\mathrm{Var}(W)=\frac{({\tau }^{\prime }(\theta ){)}^2}{{ℐ}_n(\theta )}$, if and only if $$k(W-\tau (\theta ))=S(\theta ;𝐗)$$ where $S(\theta ;𝐗)=\frac{\partial \ln ℒ(\theta ;𝐗)}{\partial \theta }$ is the score function, and $k$ is a constant.

Proof. Considering the proof for Cramer-Rao inequality, we have $$\mathrm{Var}(W)=\frac{({\tau }^{\prime }(\theta ){)}^2}{{ℐ}_n(\theta )}⟺(\mathrm{Cov}(W,S(\theta ;𝐗)){)}^2=\mathrm{Var}(W)\mathrm{Var}(S(\theta ;𝐗))$$ We can write $\mathrm{Cov}(W,S(\theta ;𝐗))$ as $\mathrm{Cov}(W-\underset{\text{constant}}{\underbrace{\tau (\theta )}},S(\theta ;𝐗))$ (by result about covariance). Also, $\mathrm{Var}(W)=\mathrm{Var}(W-\underset{\text{constant}}{\underbrace{\tau (\theta )}})$ (by result about variance). Thus, we have $$\begin{array}{cccc}& & (\mathrm{Cov}(W-\tau (\theta ),S(\theta ;𝐗)){)}^2& =\mathrm{Var}(W-\tau (\theta ))\mathrm{Var}(S(\theta ;𝐗))\\ & \iff & \frac{(\mathrm{Cov}(W-\tau (\theta ),S(\theta ;𝐗)){)}^2}{\mathrm{Var}(W-\tau (\theta ))\mathrm{Var}(S(\theta ;𝐗))}& =1\\ & \iff & \frac{(\mathrm{Cov}(S(\theta ;𝐗),W-\tau (\theta )){)}^2}{\mathrm{Var}(W-\tau (\theta ))\mathrm{Var}(S(\theta ;𝐗))}& =1\\ & \iff & (\rho (S(\theta ;𝐗),W-\tau (\theta )){)}^2& =1\\ & \iff & \rho (S(\theta ;𝐗),W-\tau (\theta ))& =\pm 1\end{array}$$ where $\rho (\cdot ,\cdot )$ is the correlation coefficient between two random variables. This means $S(\theta ;𝐗)$ increases or decreases linearly with $W-\tau (\theta )$, i.e., $$S(\theta ;𝐗)=k(W-\tau (\theta ))+c$$ for some constants $c,k$. Now, it suffices to show that the constant $c$ is actually zero.

We know that $𝔼[W]=\tau (\theta )$ (since $W$ is an unbiased estimator of $\tau (\theta )$), and $𝔼[S(\theta ;𝐗)]=0$ (from remark about Fisher information). Thus, applying expectations on both side gives $$𝔼[S(\theta ;𝐗)]=k𝔼[W-\tau (\theta )]+c⟺𝔼[S(\theta ;𝐗)]=k(\underset{=0}{\underbrace{𝔼[W]-\tau (\theta )}})+c⟺0=0+c⟺c=0.$$ Then, the result follows.

$◻$

Remark.

• Considering the proof, we know that if we have such attainable condition satisfied, the variance of the unbiased estimator $W$ equals the CRLB of $\tau (\theta )$, i.e., that estimator is the UMVUE of $\tau (\theta )$.

Example. We have shown that the log-likelihood function of a random sample $X_1,\dots ,X_n$ from the normal distribution $𝒩(\mu ,{\sigma }^2)$ is $\ln ℒ(\mu ,{\sigma }^2)=-\frac{n}{2}\ln (2\pi {\sigma }^2)-\sum _{i=1}^n\frac{(x_i-\mu {)}^2}{2{\sigma }^2}$. Prove that the CRLB of $\mu$ is attainable using the attainable conditions for the CRLB.

Proof. The score function is $$S(\mu )=\frac{\partial \ln ℒ(\mu ,{\sigma }^2)}{\partial \mu }=\frac{\partial }{\partial \mu }(-\frac{n}{2}\ln (2\pi {\sigma }^2)-\sum _{i=1}^n\frac{(X_i-\mu {)}^2}{2{\sigma }^2})=-\sum _{i=1}^n\frac{-2X_i+2\mu }{2{\sigma }^2}=\sum _{i=1}^n\frac{X_i-\mu }{{\sigma }^2}=\frac{1}{{\sigma }^2}(\sum _{i=1}^n{X}_i-\sum _{i=1}^n\mu )=\frac{n}{{\sigma }^2}(\bar{X}-\mu ).$$ Since we have $\tau (\mu )=\mu$ and $\hat{\mu }=\bar{X}$ (which is an unbiased estimator of $\mu$), the attainable conditions for the CRLB are satisfied (the constant "$k$" is $\frac{n}{{\sigma }^2}$ in this case), and hence the CRLB of $\mu$ is attainable.

$◻$

Remark.

• Indeed, we know that the CRLB of $\mu$ is attainable before this proof since we have already found an unbiased estimator of $\mu$, namely $\bar{X}$, whose variance is exactly equal to the CRLB previously.

Example. Continue from the previous example. Prove that the CRLB of ${\sigma }^2$ is not attainable using the attainable conditions for the CRLB.

Proof. The score function in this case is $$S(\sigma )=\frac{\partial \ln ℒ(\mu ,{\sigma }^2)}{\partial \sigma }=\frac{\partial }{\partial \sigma }(-\frac{n}{2}\ln (2\pi {\sigma }^2)-\sum _{i=1}^n\frac{(X_i-\mu {)}^2}{2{\sigma }^2})=-\frac{n}{2}\frac{4\pi \sigma }{2\pi {\sigma }^2}+\frac{2}{2{\sigma }^3}\sum _{i=1}^n(X_i-\mu {)}^2=-\frac{n}{\sigma }+\frac{1}{{\sigma }^3}\sum _{i=1}^n(X_i-\mu {)}^2=\underset{\text{constant}}{\underbrace{\frac{n}{{\sigma }^3}}}(\sum _{i=1}^n\frac{(X_i-\mu {)}^2}{n}-{\sigma }^2).$$ Taking the constant $k=\frac{n}{{\sigma }^3}$, a potential candidate for the unbiased estimator $\hat{\sigma }$ that achieves the CRLB is $\sum _{i=1}^n\frac{(X_i-\mu {)}^2}{n}$. However, we notice that $\sum _{i=1}^n\frac{(X_i-\mu {)}^2}{n}$ is not calculable since $\mu$ is unknown. It follows that there does not exist some $W$ such that $S(\sigma )=k(W-\tau (\sigma ))$, where $k$ is some constant and $\tau (\sigma )={\sigma }^2$.

$◻$

Remark.

• Even if we know that the CRLB of ${\sigma }^2$ is not attainable, we still do not know whether $\frac{n}{n-1}\cdot S^2$ is the UMVUE, since it is possible that some unbiased estimator with smaller variance (but not achieving the CRLB).

We have discussed MLE previously, and MLE is actually a "best choice" asymptotically (i.e., as the sample size $n\to \mathrm{\infty }$) according to the following theorem.

Theorem. Suppose $\hat{\theta }$ is the MLE of an unknown parameter $\theta$ from a distribution. Then, under some regularity conditions, as $n\to \mathrm{\infty }$, $$\frac{\hat{\theta }-\theta }{\sqrt{1/{ℐ}_n(\theta )}}\stackrel{d}{\to }𝒩(0,1).$$

Proof. Partial proof: we consider the Taylor series of order 2 for $\frac{d}{d\theta }\ln ℒ(\theta )$, and we will get $$\frac{d}{d\theta }\ln ℒ(\hat{\theta })=\frac{d}{d\theta }\ln ℒ(\theta )+(\hat{\theta }-\theta )\frac{d^2}{d{\theta }^2}\ln ℒ(\theta )+\frac{1}{2}(\hat{\theta }-\theta {)}^2\frac{d^3}{d{\theta }^3}\ln ℒ(\theta ){|}_{\theta ={\theta }^{\ast }}$$ where ${\theta }^{\ast }$ is between $\theta$ and $\hat{\theta }$. Since $\hat{\theta }$ is the MLE of $\theta$, from the derivative test, we know that $\frac{d}{d\theta }\ln ℒ(\hat{\theta })=0$ (we apply regularity condition to ensure the existence of this derivative). Hence, we have $$\begin{array}{cccc}& & \frac{d}{d\theta }\ln ℒ(\theta )+(\hat{\theta }-\theta )\frac{d^2}{d{\theta }^2}\ln ℒ(\theta )+\frac{1}{2}(\hat{\theta }-\theta {)}^2\frac{d^3}{d{\theta }^3}\ln ℒ(\theta ){|}_{\theta ={\theta }^{\ast }}& =0\\ & \Rightarrow & -\sqrt{n}(\hat{\theta }-\theta )\frac{d^2}{d{\theta }^2}\ln ℒ(\theta )-\frac{\sqrt{n}}{2}(\hat{\theta }-\theta {)}^2\frac{d^3}{d{\theta }^3}\ln ℒ(\theta ){|}_{\theta ={\theta }^{\ast }}=\sqrt{n}\frac{d}{d\theta }\ln ℒ(\theta )\\ & \Rightarrow & \sqrt{n}(\hat{\theta }-\theta )=\frac{\frac{d}{d\theta }\ln ℒ(\theta )/\sqrt{n}}{-n^{-1}\frac{d^2}{d{\theta }^2}\ln ℒ(\theta )-(2n{)}^{-1}(\hat{\theta }-\theta )\frac{d^3}{d{\theta }^3}\ln ℒ(\theta ){|}_{\theta ={\theta }^{\ast }}}.\end{array}$$ Since $$\mathrm{Var}\left(\sum _{i=1}^n\frac{\partial \ln f(X_i;\theta )}{\partial \theta }\right)=\sum _{i=1}^n\mathrm{Var}\left(\frac{\partial \ln f(X_i;\theta )}{\partial \theta }\right)=\sum _{i=1}^n𝔼\left[{\left(\frac{\partial \ln f(X_i;\theta )}{\partial \theta }\right)}^2\right]=nℐ(\theta )(1),$$ by central limit theorem, $$\frac{\frac{d}{d\theta }\ln ℒ(\theta )}{\sqrt{n}}=\frac{1}{\sqrt{n}}\sum _{i=1}^n\frac{\partial \ln f(X_i;\theta )}{\partial \theta }\stackrel{d}{\to }𝒩(0,(1/n)nI(\theta ))\equiv 𝒩(0,ℐ(\theta )).$$ Furthermore, we apply the weak law of large number to show that $$-n^{-1}\frac{d^2}{d{\theta }^2}\ln ℒ(\theta )=-\frac{1}{n}\sum _{i=1}^n\frac{{\partial }^2\ln f(X_i;\theta )}{\partial {\theta }^2}\stackrel{p}{\to }-𝔼\left[\frac{{\partial }^2\ln f(X_i;\theta )}{\partial {\theta }^2}\right]=ℐ(\theta )(2).$$ It can be shown in a quite complicated way (and using regularity conditions) that $$-(2n{)}^{-1}(\hat{\theta }-\theta )\frac{d^3}{d{\theta }^3}\ln ℒ(\theta ){|}_{\theta ={\theta }^{\ast }}\stackrel{p}{\to }0.(3).$$ Considering $(2)$ and $(3)$, using property of convergence in probability, we have $$-n^{-1}\frac{d^2}{d{\theta }^2}\ln ℒ(\theta )-(2n{)}^{-1}(\hat{\theta }-\theta )\frac{d^3}{d{\theta }^3}\ln ℒ(\theta ){|}_{\theta ={\theta }^{\ast }}\stackrel{p}{\to }ℐ(\theta )+0=ℐ(\theta )(4).$$ Considering $(1)$ and $(4)$, and using Slutsky's theorem, we have $$\sqrt{n}(\hat{\theta }-\theta )=\frac{\frac{d}{d\theta }\ln ℒ(\theta )/\sqrt{n}}{-n^{-1}\frac{d^2}{d{\theta }^2}\ln ℒ(\theta )-(2n{)}^{-1}(\hat{\theta }-\theta )\frac{d^3}{d{\theta }^3}\ln ℒ(\theta ){|}_{\theta ={\theta }^{\ast }}}\stackrel{d}{\to }\frac{Y}{ℐ(\theta )}$$ where $Y\sim 𝒩(0,ℐ(\theta ))$, and hence $\frac{Y}{ℐ(\theta )}\sim 𝒩(0,\frac{ℐ(\theta )}{[ℐ(\theta ){]}^2})\equiv 𝒩(0,1/ℐ(\theta ))$. It follows that $$\sqrt{n}(\hat{\theta }-\theta )\stackrel{d}{\to }𝒩(0,1/ℐ(\theta )).$$ This means $$\hat{\theta }-\theta \stackrel{d}{\to }𝒩(0,1/(nℐ(\theta )))\equiv 𝒩(0,1/{ℐ}_n(\theta )),$$ and thus $$\frac{\hat{\theta }-\theta }{\sqrt{1/{ℐ}_n(\theta )}}\stackrel{d}{\to }𝒩(0,\frac{1/(nℐ(\theta ))}{1/\underset{=nℐ(\theta )}{\underbrace{{ℐ}_n(\theta )}}})\equiv 𝒩(0,1)$$ as desired.

$◻$

Remark.

• Equivalently, we can write $\hat{\theta }\stackrel{d}{\to }𝒩(\theta ,1/{ℐ}_n(\theta ))$. Thus, the variance of MLE of $\theta$ achieves the CRLB of $\theta$ asymptotically. This means the MLE of $\theta$ is the UMVUE of $\theta$ asymptotically.
• The regularity conditions are basically similar to the regularity conditions mentioned in Cramer-Rao inequality.

Since we are not able to use the CRLB to find UMVUE in some situations, we will introduce another method to find UMVUE in the following, which uses the concepts of sufficiency and completeness.

Sufficiency

Intuitively, a sufficient statistic $T(X_1,\dots ,X_n)$, which is a function of a given random sample $X_1,\dots ,X_n$, contains all information needed for estimating the unknown parameter (vector) $\theta$. Thus, the statistic $T(X_1,\dots ,X_n)$ itself is "sufficient" for estimating the unknown parameter (vector) $\theta$.

Formally, we can define and describe sufficient statistic as follows:

Definition. (Sufficient statistic) A statistic $T=T(X_1,\dots ,X_n)$ is a sufficient statistic for the unknown parameter (vector) $\theta$ if the conditional distribution of the random sample $X_1,\dots ,X_n$ given $T$ does not depend on $\theta$.

Remark.

• The definition may be expressed as

$$f(x_1,\dots ,x_n|T;\theta )=f(x_1,\dots ,x_n|T)$$

where $f$ is the joint pdf or pmf of $X_1,\dots ,X_n$.

• The equation means the joint conditional pmf or pdf of $X_1,\dots ,X_n$ given (the value of) $T$ is the same as the joint conditional pmf or pdf of $X_1,\dots ,X_n$ given (the values of $T$), and with the parameter value $\theta$.

• This means the pmf of pdf is not changed even if the parameter value $\theta$ is provided, which in turn means the joint conditional pmf or pdf of $X_1,\dots ,X_n$, given the value of $T$, actually does not depend on $\theta$.

• $f(x_1,\dots ,x_n|T)$ refers to $f_{X_1,\dots ,X_n|T}(x_1,\dots ,x_n|t)$ before the realization $T=t$, and it is a random variable (the randomness comes from $T$).
• After the realization $T=t$, the equation still holds ($T$ is modified to $T=t$).

Example. Consider a random sample $X_1,\dots ,X_n$ from $𝒩(\mu ,{\sigma }^2)$. It can be shown that $\bar{X}$ is a sufficient statistic for $\mu$, but not a sufficient statistic for ${\sigma }^2$.

This can be shown by applying the definition. However, we will later give an alternative and often more convenient method to check the sufficiency of a statistic, and find sufficient statistics. We will explain informally why it is true here.

• $\bar{X}$ contains the information of the central tendency of the distribution, which should be the information needed to estimate the mean $\mu$. Hence, it is a sufficient statistic for $\mu$.
• However, $\bar{X}$ does not contain the information of the dispersion of the distribution (it only tells that the "central location", but for a particular central location, the dispersion can be very different), which should be the information needed to estimate the variance ${\sigma }^2$. Hence, it is not a sufficient statistic for ${\sigma }^2$.

Remark.

• From here, we can also expect that the sufficient statistic is not unique, since, for example, $2\bar{X}$ should also contain the information of the central tendency (since we can divide it by 2 to get the value of $\bar{X}$, and thus get the information).
• Indeed, in general, given $T$ is a sufficient statistic for $\theta$, then $v(T)$ is also a sufficient statistic for $\theta$, provided that $v$ is a bijective function (also known as invertible function, one-to-one correspondence, or bijection), so that its inverse exists.

Let us state the above remark about transformation of sufficient statistic formally below.

Proposition. Let $T$ be a sufficient statistic of the unknown parameter (vector) $\theta$. Then, $v(T)$ is also a sufficient statistic of $\theta$ for each bijective function $v$.

Now, we discuss a theorem that helps us to check the sufficiency of a statistic, namely (Fisher-Neyman) factorization theorem.

Theorem. (Factorization theorem) Let $f(x_1,\dots ,x_n;\theta )$ be the joint pdf of pmf of a random sample $X_1,\dots ,X_n$. A statistic $T=T(X_1,\dots ,X_n)$ is a sufficient statistic of $\theta$ if and only if there exist functions $g$ and $h$ such that $$f(x_1,\dots ,x_n;\theta )=g(T(x_1,\dots ,x_n);\theta )h(x_1,\dots ,x_n)$$ where $g$ depends on $x_1,\dots ,x_n$ only through $T(x_1,\dots ,x_n)$, and $h$ does not depend on $\theta$.

Proof. Since the proof for continuous case is quite complicated, we will only give a proof for the discrete case. For simplicity of presentation, let $𝐗=(X_1,\dots ,X_n)$, $T=T(X_1,\dots ,X_n)$, $𝐱=(x_1,\dots ,x_n)$, and $t=T(x_1,\dots ,x_n)$, and hence there are notations for different types of pmfs from these. By definition, $f_{𝐗|T}(𝐱|t;\theta )=f_{𝐗|T}(𝐱,t)$. Also, we have $𝐗=𝐱⟺𝐗=𝐱\cap T(𝐗)=T(𝐱)⟺𝐗=𝐱\cap T=t$. Thus, we can write $$f_{𝐗,T}(𝐱,t;\theta )=f_{𝐗}(𝐱;\theta )(\ast )$$.

"only if" ($\Rightarrow$) direction: Assume $T$ is a sufficient statistic. Then, we choose $g(t;\theta )=f_T(t;\theta )$ and $h(𝐱)=f_{𝐗|T}(𝐱|t)$, which does not depend on $\theta$ by the definition of sufficient statistic. It remains to verify that the equation actually holds for this choice.

Hence, $$f_{𝐗}(𝐱;\theta )=f_{𝐗,T}(𝐱,t;\theta )\stackrel{\text{ def }}{=}{f}_{𝐗|T}(𝐱|t;\theta )f_T(t;\theta )\stackrel{\text{ sufficiency }}{=}{f}_{𝐗|T}(𝐱|t)f_T(t;\theta )=h(𝐱)g(t;\theta ).$$

"if" ($\Leftarrow$) direction: Assume we can write $f_{𝐗}(𝐱;\theta )=g(t;\theta )h(𝐱)$. Then, $$f_T(t;\theta )\stackrel{\text{ marginal pmf }}{=}\sum _{𝐱}^{}{f}_{𝐗,T}(𝐱,t;\theta )\stackrel{\text{ (*) }}{=}\sum _{𝐱}^{}{f}_{𝐗}(𝐱;\theta )\stackrel{\text{ assumption }}{=}\sum _{𝐱}^{}g(t;\theta )h(𝐱)=\underset{\text{independent from }𝐱}{\underbrace{g(t;\theta )}}\sum _{𝐱}^{}h(𝐱).$$ Now, we aim to show that $f_{𝐗|T}(𝐱|t)$ does not depend on $\theta$, which means $T$ is a sufficient statistic for $\theta$. We have $$f_{𝐗|T}(𝐱|t)\stackrel{\text{ def }}{=}\frac{f_{𝐗,T}(𝐱,t;\theta )}{f_T(t;\theta )}\stackrel{\text{ (*) }}{=}\frac{f_{𝐗}(𝐱;\theta )}{f_T(t;\theta )}=\frac{\stackrel{\text{assumption}}{\overbrace{g(t;\theta )h(𝐱)}}}{\underset{\text{above}}{\underbrace{g(t;\theta )\sum _{𝐱}^{}h(𝐱)}}}=\frac{h(𝐱)}{\sum _{𝐱}^{}h(𝐱)},$$ which does not depend on $\theta$, as desired.

$◻$

Remark.

• $h(x_1,\dots ,x_n)$ can also be a constant, which does not depend on $\theta$ clearly.

Example. Consider a random sample $X_1,\dots ,X_n$ from $𝒩(\mu ,{\sigma }^2)$. Find a sufficient statistic of $\theta =(\mu ,{\sigma }^2)$.

Solution: The joint pdf of $X_1,\dots ,X_n$ is $$\begin{array}{cc}f(x_1,\dots ,x_n;\theta )& =\prod _{i=1}^n\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp (-\frac{(x_i-\mu {)}^2}{2{\sigma }^2})\\ & =(2\pi {\sigma }^2{)}^{-n/2}\exp \left(\sum _{i=1}^n\frac{(x_i-\mu {)}^2}{2{\sigma }^2}\right)\\ & =(2\pi {\sigma }^2{)}^{-n/2}\exp \left(\sum _{i=1}^n\frac{(x_i-\bar{x}+\bar{x}-\mu {)}^2}{2{\sigma }^2}\right)\\ & =(2\pi {\sigma }^2{)}^{-n/2}\exp \left(\sum _{i=1}^n\frac{(x_i-\bar{x}{)}^2+2(x_i-\bar{x})(\bar{x}-\mu )+(\bar{x}-\mu {)}^2}{2{\sigma }^2}\right)\\ & =(2\pi {\sigma }^2{)}^{-n/2}\exp \left(\sum _{i=1}^n\frac{(x_i-\bar{x}{)}^2+(\bar{x}-\mu {)}^2}{2{\sigma }^2}\right)& (\sum _{i=1}^n(x_i-\bar{x})(\bar{x}-\mu )=(\bar{x}-\mu )\sum _{i=1}^n(x_i-\bar{x})=(\bar{x}-\mu )(\sum _{i=1}^n{x}_i-\sum _{i=1}^n\bar{x})=(\bar{x}-\mu )(n\bar{x}-n\bar{x})=0)\\ & =(2\pi {\sigma }^2{)}^{-n/2}\exp \left(\frac{1}{2{\sigma }^2}(\sum _{i=1}^n(x_i-\bar{x}{)}^2+\sum _{i=1}^n(\bar{x}-\mu {)}^2)\right)\\ & =\underset{h(x_1,\dots ,x_n)}{\underbrace{(2\pi {)}^{-n/2}}}\underset{g(T(x_1,\dots ,x_n);\theta )}{\underbrace{{\sigma }^{-n}\exp \left(\frac{1}{2{\sigma }^2}(n{s}^2+n(\bar{x}-\mu {)}^2)\right)}}& ((\bar{x}-\mu {)}^2\text{ is independent from }i).\\ \end{array}$$ Notice that the function $g$ depends on $x_1,\dots ,x_n$ only through $T(x_1,\dots ,x_n)=(\bar{x},s^2)$, so we can conclude that $T(X_1,\dots ,X_n)=(\bar{X},S^2)$.

Remark.

• We can also write $(\bar{X},S^2)$ as $(S^2,\bar{X})$, which is also a sufficient statistic for $\theta$.

• Intuitively, this is because the latter one also contains the same statistics, and thus contains the same information.
• Alternatively, we can define the function $v$ as $(z_1,z_2)\mapsto (z_2,z_1)$, which is a bijective function, so $v(\bar{X},S^2)=(S^2,\bar{X})$ is also a sufficient statistic for $\theta$.

• We need to separate "${\sigma }^{-n}$" out from "$(2\pi {\sigma }^2{)}^{-n/2}$", since for the function $h(x_1,\dots ,x_n)$, it cannot depend on $\theta =(\mu ,{\sigma }^2)$. Hence, we cannot include "${\sigma }^{-n}$" in the definition of the $h(x_1,\dots ,x_n)$ function.
• There are many ways to define the $g$ and $h$ functions in this case.

For some "nice" distributions, which belong to exponential family, sufficient statistics can be found using another alternative method easily and more conveniently. This method works because of the "nice" form of the pdf or pmf of those distributions, which can be characterized as follows:

Definition. (Exponential family) The distribution of a random variable $X$ belongs to the exponential distribution if the pdf or pmf of $X$ has the form of $$f(x;\theta )=h(x)g(\theta )\exp (\sum _{i=1}^s{\eta }_i(\theta )T_i(x))$$ where $\theta =({\theta }_1,\dots ,{\theta }_s)\in \mathrm{\Theta }\subseteq {ℝ}^s$, for some functions $h,g,{\eta }_i,T_i$ ($i=1,2\dots ,s$).

Remark.

• The value of $s$ depends on the number of unknown parameters.

• Notice that $s$ can be 1, and in this case the "$\theta$" is just a single parameter.

• Exponential family includes many common distributions, e.g. normal, exponential, gamma, chi squared, beta, Bernoulli, Poisson, geometric, etc.

• However, some common distributions do not belong to exponential family, e.g. Student's $t$-distribution, $F$-distribution, Cauchy distribution and hypergeometric distribution.

Example. The normal distribution belongs to the exponential family, where $\theta =(\mu ,{\sigma }^2)\in {ℝ}^2$ (so the "$s$" is 2 in this case), since its pdf can be expressed as $$f(x;\theta )=\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp (-\frac{(x-\mu {)}^2}{2{\sigma }^2})=\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp (-\frac{x^2-2\mu x+{\mu }^2}{2{\sigma }^2})=[\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp (-\frac{{\mu }^2}{2{\sigma }^2})]\exp (-\frac{x^2-2\mu x}{2{\sigma }^2})=\underset{h(x)}{\underbrace{(1)}}\underset{g(\theta )}{\underbrace{[\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp (-\frac{{\mu }^2}{2{\sigma }^2})]}}\exp [\underset{{\eta }_1(\theta )}{\underbrace{\frac{\mu }{{\sigma }^2}}}\cdot \underset{T_1(x)}{\underbrace{x}}+\underset{{\eta }_2(\theta )}{\underbrace{(-\frac{1}{2{\sigma }^2})}}\cdot \underset{T_2(x)}{\underbrace{x^2}}]$$

Theorem. (Sufficient statistic for exponential family) Let $X_1,\dots ,X_n$ be a random sample from a distribution belonging to the exponential family, with pdf or pmf $f(x;\theta )$ where $\theta \in {ℝ}^s$. Then, a sufficient statistic for $\theta$ is $$T(X_1,\dots ,X_n)=(\sum _{j=1}^n{T}_1(X_j),\dots ,\sum _{j=1}^n{T}_s(X_j)).$$

Proof. Since the distribution belongs to the exponential family, the joint pdf or pmf of $X_1,\dots ,X_n$ can be expressed as $$\begin{array}{cc}f(x_1,\dots ,x_n;\theta )& =\prod _{j=1}^n[h(x_j)g(\theta )\exp (\sum _{i=1}^s{\eta }_i(\theta )T_i(x_j))]\\ & =[\prod _{j=1}^{n}h(x_j)](g(\theta ){)}^n\exp (\sum _{j=1}^n\sum _{i=1}^s{\eta }_i(\theta )T_i(x_j))\\ & =[\prod _{j=1}^{n}h(x_j)](g(\theta ){)}^n\exp (\sum _{i=1}^s\sum _{j=1}^n{\eta }_i(\theta )T_i(x_j))& (\text{changing summation order, where the upper bounds are constants})\\ & =[\prod _{j=1}^{n}h(x_j)](g(\theta ){)}^n\exp (\sum _{i=1}^s\underset{\text{independent from }j}{\underbrace{{\eta }_i(\theta )}}\sum _{j=1}^n{T}_i(x_j))\\ & =[\prod _{j=1}^{n}h(x_j)](g(\theta ){)}^n\exp ({\eta }_1(\theta )\sum _{j=1}^n{T}_1(x_j)+\cdots +{\eta }_s(\theta )\sum _{j=1}^n{T}_s(x_j)).\\ \end{array}$$ From here, for applying the factorization theorem, we can identify the purple part of the function as "$h(x_1,\dots ,x_n)$", and the red part of the function as "$g(T(x_1,\dots ,x_n);\theta )$". We can notice that the red part of the function depends on $x_1,\dots ,x_n$ only through $(\sum _{j=1}^n{T}_1(x_j),\dots ,\sum _{j=1}^n{T}_s(x_j))$. The result follows.

$◻$

Example. Consider a random sample $X_1,\dots ,X_n$ from $𝒩(\mu ,{\sigma }^2)$. Show that a sufficient statistic of $\theta =(\mu ,{\sigma }^2)$ is $(\bar{X},S^2)$, using the result for finding sufficient statistic for exponential family.

Proof. From previous example, we have shown that normal distribution belongs to the exponential family, and from the expression there, we can see that a sufficient statistic of $\theta$ is $T=(\sum _{j=1}^{n}X,\sum _{j=1}^n{X}^2)=(n\bar{X},n\overline{X^2})$.

Since $S^2=\frac{1}{n}\sum _{j=1}^n(X_j-\bar{X}{)}^2=\frac{1}{n}\sum _{j=1}^n(X_j^2-2X_j\bar{X}+(\bar{X}{)}^2)=\frac{\sum _{j=1}^n{X}_j^2}{n}-\frac{2\bar{X}}{n}\sum _{j=1}^n{X}_j+(\bar{X}{)}^2=\overline{X^2}-2(\bar{X}{)}^2+(\bar{X}{)}^2=\overline{X^2}-(\bar{X}{)}^2$, we can define function $v$ as $$(z_1,z_2)\mapsto (z_1/n,z_2/n-(z_1/n{)}^2),$$ which can be shown to be an bijective function.

Thus, $v(T)=(\bar{X},S^2)$ is also an sufficient statistic of $\theta$.

$◻$

Now, we will start discussing how is sufficient statistic related to UMVUE. We begin our discussion by Rao-Blackwell theorem.

Theorem. (Rao-Blackwell theorem) Let $W$ be an arbitrary unbiased estimator of $\tau (\theta )$, and $T$ be a sufficient statistic of $\theta$. Define $\phi (T)=𝔼[W|T]$. Then, $\phi (T)$ is an unbiased estimator of $\tau (\theta )$ and $\mathrm{Var}(\phi (T))\le \mathrm{Var}(W)$.

Proof. Assume $W$ is an arbitrary unbiased estimator of $\tau (\theta )$, and $T$ is a sufficient statistic of $\theta$.

First, we prove that $\phi (T)$ is an unbiased estimator of $\tau (\theta )$. Before proving the unbiasedness, we should ensure that $\phi (T)$ is actually an estimator, i.e., it is a statistic, which is a function of random sample, and needs to be independent from $\theta$ (so that it is calculable): since $W$ is a function of random sample, and $T$ is a sufficient statistic, which make the conditional distribution of $W$, given $T$, independent of $\theta$. Also, $\phi (T)=𝔼[W|T]$ is a function of $W$, and thus is also a function of random sample.

Now, we prove that $\phi (T)$ is an unbiased estimator of $\tau (\theta )$: since $𝔼[\phi (T)]=𝔼[𝔼[W|T]]\stackrel{\text{ law of total expectation }}{=}𝔼[W]\stackrel{\text{ unbiasedness }}{=}\tau (\theta )$, $\phi (T)$ is an unbiased estimator of $\tau (\theta )$.

Next, we prove that $\mathrm{Var}(\phi (T))\le \mathrm{Var}(W)$: by law of total variance, we have $$\mathrm{Var}(W)=\mathrm{Var}(𝔼[W|T])+𝔼[\mathrm{Var}(W|T)]\stackrel{\text{ def }}{=}\mathrm{Var}(\phi (T))+\stackrel{\ge 0}{\overbrace{𝔼[\underset{\ge 0}{\underbrace{\mathrm{Var}(W|T)}}]}}\ge \mathrm{Var}(\phi (T)),$$ as desired.

$◻$

Remark.

• The random variable $\phi (T)=𝔼[W|T]$ is determined by first finding $\phi (t)=𝔼[W|T=t]$, and then replacing $t$ by $T$. Here, $\phi (t)$ is the realization of $\phi (T)$.
• From Rao-Blackwell theorem, we know that $\phi (T)=𝔼[W|T]$ is a better (or at least "same quality") estimator than $W$ in efficiency sense. Notice that the theorem does not state that $\phi (T)$ is the best estimator in efficiency sense (i.e., UMVUE). Instead, it only states that $\phi (T)$ is better than $W$ in efficiency sense.
• After applying this theorem once, the $\phi (T)$ can act as the "arbitrary unbiased estimator of $\tau (\theta )$", and we can apply this theorem again, and so on. This means that after applying this theorem many times, the "$\phi (T)$" we get will be an UMVUE.

• We can interpret this process as keep "improving" the unbiased estimator $W$, until it becomes the best one (in efficiency sense), i.e., it is an UMVUE.
• Since UMVUE is unique, UMVUE must be the conditional expectation of a random variable given a sufficient statistic $T$, which is a function of $T$.
• Thus, we now can narrow the candidates for UMVUE to functions of sufficient statistic $T$.

To actually determine the UMVUE, we need another theorem, called Lehmann-Scheffé theorem, which is based on Rao-Blackwell theorem, and requires the concept of completeness.

Completeness

Definition. (Complete statistic) Let $X_1,\dots ,X_n$ be a random sample from a distribution with a parameter (vector) $\theta$ lying in the parameter space $\mathrm{\Theta }$. A statistic $T$ is a complete statistic if $𝔼[g(T)]=0$ for each $\theta \in \mathrm{\Theta }$ implies $\mathbb{P}(g(T)=0)=1$ for each $\theta \in \mathrm{\Theta }$.

When a random sample $X_1,\dots ,X_n$ is from a distribution in exponential family, then a complete statistic can also be founded easily, similar to the case for sufficient statistic.

Theorem. (Complete statistic for exponential family) If $X_1,\dots ,X_n$ is a random sample from a distribution in exponential family where the unknown parameter (vector) $\theta \in \mathrm{\Theta }\subseteq {ℝ}^s$, Then, $$T(X_1,\dots ,X_n)=(\sum _{j=1}^n{T}_1(X_j),\sum _{j=1}^n{T}_2(X_j),\dots ,\sum _{j=1}^n{T}_s(X_j))$$ is a complete statistic, given that the parameter space $\mathrm{\Theta }$ contains an open set in ${ℝ}^s$.

Proof. Omitted.

$◻$

Remark.

• Open sets are a generalization of open intervals. Indeed, open sets in $ℝ$ is simply open intervals.
• Intuitively, open sets refers to the sets that, for each point in the set, the set contain all points that are sufficiently near to that point.
• For example, a line in ${ℝ}^2$ (which can be interpreted as a set) is not an open set since for each point in the line, the line does not contain all points that are sufficiently near to that point (there are some points "above" and "below" the line that are not contained in the set).
• Also, a disk (a region in a plane bounded by a circle) in ${ℝ}^3$ is not an open set since for each point in the disk, the disk does not contain all points that are sufficiently near to that point (there are some points "above" and "below" the disk that are not contained in the disk).
• From the previous theorem about sufficient statistic for exponential family, we know that $T(X_1,\dots ,X_n)$ is also a sufficient statistic of $\theta$ under such conditions.

• When a statistic is sufficient for a parameter (vector) $\theta$ and is also a complete statistic, we call such statistic as a complete and sufficient statistic for $\theta$.

Theorem. (Lehmann-Scheffé theorem) If $T$ is a complete and sufficient statistic for $\theta$ and $𝔼[\phi (T)]=\tau (\theta )$, then $\phi (T)$ is the unique UMVUE of $\tau (\theta )$ (with probability 1).

Proof. Assume $T$ is a complete and sufficient statistic for $\theta$ and $𝔼[\phi (T)]=\tau (\theta )$.

Since $T$ is a sufficient statistic for $\theta$, we can apply the Rao-Blackwell theorem. From Rao-Blackwell theorem, if $W$ is an arbitrary unbiased estimator of $\tau (\theta )$, then $\phi (T)$ is another unbiased estimator where $\mathrm{Var}(\phi (T))\le \mathrm{Var}(W)$.

To prove that $\phi (T)$ is the unique UMVUE of $\tau (\theta )$, we proceed to show that regardless of the choice of the unbiased estimator $W$ of $\tau (\theta )$, we get the same $\phi (T)$ from the Rao-Blackwell theorem (with probability 1). Then, we will have for every possible unbiased estimator $W$ of $\tau (\theta )$, $\mathrm{Var}(\phi (T))\le \mathrm{Var}(W)$ (with probability 1) , which means $\phi (T)$ is the UMVUE, and is also the unique UMVUE since we always get the same $\phi (T)$ .

Assume that $W^{\prime }$ is another unbiased estimator of $\tau (\theta )$ ($W^{\prime }\ne W$). By Rao-Blackwell theorem again, there is an unbiased estimator $\psi (T)=𝔼[W^{\prime }|T]$ ($\psi (T)\ne \phi (T)$) where $\mathrm{Var}(\psi (T))\le \mathrm{Var}(W^{\prime })$. Since both $\phi (T)$ and $\psi (T)$ are unbiased estimators of $\tau (\theta )$, we have for each $\theta \in \mathrm{\Theta }$, $$𝔼[\phi (T)]=𝔼[\psi (T)]⟹𝔼[\phi (T)-\psi (T)]=0.$$ Since $T$ is a complete statistic, we have $$\mathbb{P}(\phi (T)-\psi (T)=0)=1⟹\mathbb{P}(\phi (T)=\psi (T))=1,$$ which means $\phi (T)=\psi (T)$ (with probability 1), i.e., we get the same $\phi (T)$ from the Rao-Blackwell theorem in this case (with probability 1).

$◻$

Remark.

• The "$\phi (T)$" in this theorem is a function of $T$, and we know from the proof and Rao-Blackwell theorem that it is actually $𝔼[W|T]$ where $W$ is an arbitrary unbiased estimator of $\tau (\theta )$.

• Thus, when we apply this theorem, as long as we can find a function of $T$, $\varphi (T)$, (perhaps by some inspections) such that $𝔼[\varphi (T)]=\tau (\theta )$, we know that $\varphi (T)$ is the unique UMVUE of $\tau (\theta )$. Also, due to the uniqueness of UMVUE, the $\varphi (T)$ is actually $\phi (T)=𝔼[W|T]$ where $W$ is an arbitrary unbiased estimator of $\tau (\theta )$.
• We can find $\phi (T)$ by some inspections, as in above, in simple cases. However, in more complicated cases, it may not be immediately transparent that what should be the explicit form of $\phi (T)$ such that $𝔼[\phi (T)]=\tau (\theta )$. In such case, we need to find an unbiased estimator of $\tau (\theta )$ and evaluate $𝔼[W|T]$ to get the explicit form of $\phi (T)$.

Example. Consider a random sample $X_1,\dots ,X_n$ from $𝒩(\mu ,{\sigma }^2)$. Let the unknown parameter vector $\theta =(\mu ,{\sigma }^2)$.

(a) Show that the sufficient statistic for $\theta$, namely $(\bar{X},S^2)$, is also a complete statistic.

(b) Hence, show that $\bar{X}$ and $\frac{n}{n-1}\cdot S^2$ is the UMVUE of $\mu$ and ${\sigma }^2$ respectively.

Solution:

(a)

Proof. It suffices to show that the parameter space $\mathrm{\Theta }=\{(\mu ,{\sigma }^2):\mu \in ℝ,{\sigma }^2>0\}$ contains an open set in ${ℝ}^2$. This is true since the parameter space $\mathrm{\Theta }$ is the whole region above $x$-axis if we represent it using the Cartesian coordinate system, and thus contains an open set.

$◻$

(b)

Proof. Since $𝔼[\bar{X}]=\mu$ and $𝔼[\frac{n}{n-1}\cdot S^2]={\sigma }^2$ (we have shown these before), and $\bar{X}$ and $\frac{n}{n-1}\cdot S^2$ is function of complete and sufficient statistic $\bar{X}$ (of $\mu$) and $S^2$ (of ${\sigma }^2$) respectively, by Lehmann-Scheffé theorem, we have the desired result.

$◻$

Remark.

• We have shown that $\frac{n}{n-1}\cdot S^2$ does not attain the CRLB of ${\sigma }^2$, and the CRLB of ${\sigma }^2$ is actually unattainable. Thus, we were not able to determine whether $\frac{n}{n-1}\cdot S^2$ is the UMVUE of ${\sigma }^2$ before. Now, we know that $\frac{n}{n-1}\cdot S^2$ is actually the UMVUE of ${\sigma }^2$ with the help of Lehmann-Scheffé theorem.

Example. Consider a random sample $X_1,\dots ,X_n$ from the Bernoulli distribution with success probability $p$, i.e., $\mathrm{Ber}(p)$, with pmf $f(x;p)=p^x(1-p{)}^{1-x},x=0,1$.

(a) Find a complete and sufficient statistic $T$ for $p$.

(b) Hence, find the UMVUE of $p$.

(c) Show that $\mathrm{𝟏}\{X_1=1\}$ is an unbiased estimator of $p$, and $𝔼[\mathrm{𝟏}\{X_1=1\}|T]$ is the UMVUE of $p$.

Solution:

(a) The pmf $f(x;p)=p^x(1-p{)}^{1-x}=(1-p){\left(\frac{p}{1-p}\right)}^x=\underset{h(x)}{\underbrace{(1)}}\underset{g(\theta )}{\underbrace{(1-p)}}\exp \left(\underset{T(x)}{\underbrace{x}}\underset{\eta (p)}{\underbrace{\ln \left(\frac{p}{1-p}\right)}}\right)$. This means Bernoulli distribution belongs to the exponential family. Also, the parameter space $\mathrm{\Theta }=\{p:0\le p\le 1\}$ contains an open set in $ℝ$. Hence, $T=\sum _{j=1}^n{X}_j$ is a complete and sufficient statistic for $p$.

(b) Notice that $𝔼[T/n]=𝔼[\bar{X}]=\frac{np}{n}=p$. Hence, $\bar{X}$ (which is a function of $T$) is the UMVUE of $p$.

(c)

Proof. Since $𝔼[\mathrm{𝟏}\{X_1=1\}]=(1)\mathbb{P}(X_1=1)=p$, $\mathrm{𝟏}\{X_1=1\}$ is an unbiased estimator of $p$.

Now, we consider $𝔼[\mathrm{𝟏}\{X_1=1\}|T]=𝔼[\mathrm{𝟏}\{X_1=1\}|\sum _{j=1}^n{X}_j]$. We denote $\sum _{j=1}^n{X}_j$ by $S_n$. Then, this expectation becomes $𝔼[\mathrm{𝟏}\{X_1=1\}|S_n]$. In the following, we evaluate $𝔼[\mathrm{𝟏}\{X_1=1\}|S_n=s_n]$. $$\begin{array}{ccc}𝔼[\mathrm{𝟏}\{X_1=1\}|\sum _{j=1}^n{X}_j=s_n]& =(1)\mathbb{P}(\mathrm{𝟏}\{X_1=1\}=1|\sum _{j=1}^n{X}_j=s_n)& (\text{definition})\\ & =\mathbb{P}(X_1=1|\sum _{j=1}^n{X}_j=s_n)\\ & =\frac{\mathbb{P}(\sum _{j=1}^n{X}_j=s_n|X_1=1)\mathbb{P}(X_1=1)}{\mathbb{P}(\sum _{j=1}^n{X}_j=s_n)}& (\text{Bayes' theorem})\\ & =\frac{\mathbb{P}(\sum _{j=2}^n{X}_j=s_n-1)\cdot p}{\mathbb{P}(\sum _{j=1}^n{X}_j=s_n)}\\ \end{array}$$ Notice that $\sum _{j=1}^n{X}_j$ follows the binomial distribution with $n$ trials with success probability $p$, i.e., $\mathrm{Binom}(n,p)$, and $\sum _{j=2}^n{X}_j\sim \mathrm{Binom}(n-1,p)$. Hence, $$\begin{array}{ccc}\frac{\mathbb{P}(\sum _{j=2}^n{X}_j=s_n-1)\cdot p}{\mathbb{P}(\sum _{j=1}^n{X}_j=s_n)}& =\frac{(\genfrac{}{}{0ex}{}{n-1}{s_n-1})p^{s_n-1}(1-p{)}^{n-1-s_n+1}\cdot p}{(\genfrac{}{}{0ex}{}{n}{s_n})p^{s_n}(1-p{)}^{n-s_n}}& (\text{binomial distribution pmf's})\\ & =\frac{\frac{(n-1)!}{(s_n-1)!(n-s_n)!}}{\frac{n!}{s_n!(n-s_n)!}}\\ & =\frac{(n-1)!s_n(s_n-1)!}{n(n-1)!(s_n-1)!}& (s_n!=s_n(s_n-1)!\text{ and }n!=n(n-1)!)\\ & =\frac{s_n}{n}.\end{array}$$ Now, replacing $s_n$ by $S_n=\sum _{j=1}^n{X}_j$ gives $$𝔼[\mathrm{𝟏}\{X_1=1\}|\sum _{j=1}^n{X}_j]=\frac{\sum _{j=1}^n{X}_j}{n}=\bar{X},$$ which is the UMVUE of $p$, as desired.

$◻$

Exercise. Can we find the UMVUE of $p$ using the CRLB of $p$? If yes, find it using this way. If no, explain why.

Solution

No. This is because the log-likelihood function is not differentiable (it has nonzero value only when $x=0,1$), and thus the Fisher information is undefined. Hence, the CRLB does not exist.

Exercise. Consider a random sample $X_1,\dots ,X_n$ from the Poisson distribution with rate parameter $\lambda$, with pmf $f(x;\lambda )=\frac{e^{-\lambda }{\lambda }^x}{x!}$ .

(a) Find a complete and sufficient statistic for $\lambda$.

(b) Find the UMVUE of $\lambda /n$.

Solution

(a) The pmf is $$f(x;\lambda )=\frac{e^{-\lambda }{\lambda }^x}{x!}=\frac{e^{-\lambda }}{x!}\underset{={\lambda }^x}{\underbrace{\exp (x\ln \lambda )}}=\underset{h(x)}{\underbrace{\frac{1}{x!}}}\cdot \underset{g(\lambda )}{\underbrace{e^{-\lambda }}}\exp (\underset{T(x)}{\underbrace{x}}\underset{\eta (\lambda )}{\underbrace{\ln \lambda }}).$$ Hence, Poisson distribution belongs to the exponential family, and a complete and sufficient statistic for $\lambda$ is $T=\sum _{j=1}^n{X}_j$.

(b) Take $\tau (\lambda )=\lambda /n$. Since $𝔼[T]=𝔼\left[\sum _{j=1}^n{X}_j\right]=n\lambda$, we have $$𝔼[T/n^2]=\lambda /n=\tau (\lambda ).$$ Thus, the UMVUE of $\tau (\lambda )=\lambda /n$ is $\frac{\sum _{j=1}^n{X}_j}{n^2}=\frac{\bar{X}}{n}$ (which is a function of $T$).

Consistency

In the previous sections, we have discussed unbiasedness and efficiency. In this section, we will discuss another property called consistency.

Definition. (Consistent estimator) $\hat{\theta }$ is a consistent estimator of the unknown parameter $\theta$ if $\hat{\theta }\stackrel{p}{\to }\theta$.

Remark.

• By the definition of convergence in probability, $\hat{\theta }\stackrel{p}{\to }\theta$ means $\mathbb{P}(|\hat{\theta }-\theta |>\epsilon )\to 0$ as $n\to \mathrm{\infty }$, for each $\epsilon >0$.

Proposition. If $\hat{\theta }$ is an (asymptotically) unbiased estimator of an unknown parameter $\theta$ and $\mathrm{Var}(\hat{\theta })\to 0$ as $n\to \mathrm{\infty }$, then $\hat{\theta }$ is a consistent estimator of $\theta$.

Proof. Assume $\hat{\theta }$ is an (asymptotically) unbiased estimator of an unknown parameter $\theta$ and $\mathrm{Var}(\hat{\theta })\to 0$ as $n\to \mathrm{\infty }$. Since $\hat{\theta }$ is an (asymptotically) unbiased estimator of $\theta$, we have $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Bias}(\hat{\theta })=0$ (this is true for both asymptotically unbiased estimator and unbiased estimator of $\theta$). In addition to this, we have by assumption that $\underset{n\to \mathrm{\infty }}{lim}\mathrm{Var}(\hat{\theta })=0$. By definition of mean squared error, these imply that $\underset{n\to \mathrm{\infty }}{lim}\mathrm{MSE}(\hat{\theta })=0\Rightarrow \underset{n\to \mathrm{\infty }}{lim}𝔼[(\hat{\theta }-\theta {)}^2]=0$. Thus, as $n\to \mathrm{\infty }$, we have by Chebyshov's inequality (notice that $\mathrm{MSE}(\hat{\theta })=𝔼[(\hat{\theta }-\theta {)}^2]$ exist from above), for each $\epsilon >0$, $$\mathbb{P}(|\hat{\theta }-\theta |>\epsilon )\le \frac{𝔼[(\hat{\theta }-\theta {)}^2]}{{\epsilon }^2}\to \frac{0}{{\epsilon }^2}=0.$$ Since probability is nonnegative ($\ge 0$), and this probability is less than or equal to an expression that tends to be 0 as $n\to \mathrm{\infty }$, we conclude that this probability tends to be zero as $n\to \mathrm{\infty }$. That is, $\hat{\theta }$ is a consistent estimator of $\theta$.

$◻$

Remark.

• Unbiasedness alone does not imply consistency.

Example. Let $X_1,\dots ,X_n$ be a random sample from $𝒩(\mu ,{\sigma }^2)$. Then, $X_1$ is an unbiased estimator of $\mu$ since $𝔼[X_1]=\mu$. However, there exist some $\epsilon >0$ such that $$\mathbb{P}(|X_1-\mu |>\epsilon )\nrightarrow 0$$ as $n\to \mathrm{\infty }$, i.e., $\underset{n\to \mathrm{\infty }}{lim}\mathbb{P}(|X_1-\mu |>\epsilon )\ne 0$ for some $\epsilon >0$. Since $|X_1-\mu |$ is independent from $n$, this means $$\mathbb{P}(|X_1-\mu |>\epsilon )\ne 0$$ for some $\epsilon >0$, which is true. Hence, $X_1$ is not a consistent estimator of $\mu$.

Exercise.

(a) Propose a consistent estimator of $\mu$, and show that it is actually a consistent estimator of $\mu$ (Hint: consider weak law of large number).

(b) Propose a consistent estimator of the coefficient of variation (or relative standard deviation) $\frac{\sigma }{\mu }$ (assuming $\mu \ne 0$ so that it is defined), and show that it is actually a consistent estimator of $\frac{\sigma }{\mu }$ (Hint: consider weak law of large number, and properties for convergence in probability. You can use the fact that normal distribution has a finite fourth moment).

Solution

(a) $\bar{X}$ is a consistent estimator of $\mu$.

Proof. By weak law of large number (notice that mean $\mu$ and variance ${\sigma }^2$ are finite for normal distribution), $\bar{X}\stackrel{p}{\to }\mu$ as desired.

$◻$

(b) $\frac{\sqrt{S^2}}{\bar{X}}$ is a consistent estimator of $\frac{\sigma }{\mu }$.

Proof. By weak law of large number (variance is finite, and fourth moment is finite), $\overline{X^2}\stackrel{p}{\to }𝔼[X^2]$. Also, by continuous mapping theorem, since $(\bar{X}{)}^2\stackrel{p}{\to }{\mu }^2$. Thus, by properties about convergence in probability and result about sample variance, $$S^2=\overline{X^2}-(\bar{X}{)}^2\stackrel{p}{\to }𝔼[X^2]-{\mu }^2={\sigma }^2.$$ By continuous mapping theorem again, $$\sqrt{S^2}\stackrel{p}{\to }\sqrt{{\sigma }^2}=\sigma$$ (since $\sigma >0$). Hence, by properties about convergence in probability (we assume that $\mu \ne 0$) again, $$\frac{\sqrt{S^2}}{\bar{X}}\stackrel{p}{\to }\frac{\sigma }{\mu }$$ as desired.

$◻$