# Understanding variance estimators

I am having trouble understanding the following. Let $\mu$ and $\sigma^2$ be the true mean and variance, $\bar{x}$ and $s^2$ the measured mean and variance for a random variable $x$, where $$\displaystyle s^2 = \frac{1}{N+k}\sum_i (x_i-\bar{x})^2.$$

• If $s^2$ is an unbiased estimate of variance then $k=-1$.
• If $s^2$ has the smallest mean square spread from true variance then $k=1$.
• Are you asking for a derivation of these results or are you asking for a more intuitive explanation? – assumednormal Jul 5 '11 at 5:10
• @Max Preferably both :) Intuitively, I see them as equivalent because for most data $k\ll N$ so $\displaystyle \frac{1}{N+k}\approx\frac{1}{N}$ – yayu Jul 5 '11 at 6:41
• I have edited the problem. Someone would check it. I hope it is in the spirit of what you want. – Theta30 Jul 5 '11 at 8:11
• I have gone even further -- as long as I hopefully haven't changed the original idea it seems clear now. – user88 Jul 5 '11 at 8:37

The idea of the unbaised variance estimate, is to have $E(s^{2})=\sigma^{2}$ where the expectation is with respect to the sampling distribution of $s^{2}$ or equivalently, with respect to the sampling distribution of $x_1,\dots,x_N$. So if we knew the true mean and the true variance, but not the value of $s^{2}$, then $s^2$ would have expected value of $\sigma^{2}$.

Now we have:

$$E(s^{2})=E\left[ \frac{1}{N+k}\sum_i (x_i-\bar{x})^2\right] =\frac{1}{N+k}E\left[\sum_i x_i^2-N\bar{x}^{2}\right]$$

$$=\frac{1}{N+k}\left[\sum_i E(x_i^2)-NE(\bar{x}^{2})\right]$$

We can use the identity $E(Y^{2})=V(Y)+[E(Y)]^{2}$, and the fact that we know $E(x_i)=\mu$ and $V(x_i)=\sigma^{2}$, and the summation becomes:

$$\sum_i E(x_i^2)=\sum_i V(x_i)+[E(x_i)]^2=\sum_i (\sigma^2+\mu^2)=N(\sigma^2+\mu^2)$$

Now to calculate the second expectation, we re-write $NE(\bar{x}^{2})$ as follows:

$$NE(\bar{x}^{2})=NE\left(\left[\frac{1}{N}\sum_i x_i\right]^{2}\right)=\frac{1}{N}E\left(\sum_i \sum_j x_jx_i\right)=\frac{1}{N}\sum_i \sum_j E(x_jx_i)$$ $$=\frac{1}{N}\left(\sum_i E(x_i^2)+\sum_{i\neq j} E(x_jx_i)\right)=(\sigma^2+\mu^2)+\frac{1}{N}\sum_{i\neq j} E(x_jx_i)$$

Now we have another identity we can use $E(YZ)=Cov(Y,Z)+E(Y)E(Z)$. Your problem statement does not specify whether or not the sample is independent, but it does say that they have the same distribution. So we can take $Cov(x_ix_j)=\rho\sigma^2$ for some correlation $-\frac{1}{N-1}\leq\rho\leq 1$ (lower bound required for positive variance). We then get:

$$NE(\bar{x}^{2})=(\sigma^2+\mu^2)+\frac{1}{N}\sum_{i\neq j} (\rho\sigma^2+\mu^2)=(\sigma^2+\mu^2)+\frac{N(N-1)}{N}(\rho\sigma^2+\mu^2)$$

$$=\sigma^2(1+(N-1)\rho)+N\mu^{2}$$ Putting this together we have:

$$E(s^2)=\frac{1}{N+k}\left[N(\sigma^2+\mu^2)-\left(\sigma^2(1+(N-1)\rho)+N\mu^{2}\right)\right]=\frac{N-1}{N+k}\sigma^2(1-\rho)$$

So if we choose $k=-1$ and we assume that $\rho=0$ (i.e. independence), then we have $E(s^2)=\sigma^2$, and we get an unbiased estimate for $s^2$. However, if we assume that $\rho=-\frac{1}{N-1}$ (i.e. the sum is fixed), then we get $E(s^2)=\frac{N-1}{N+k}\sigma^2(1+\frac{1}{N-1})=\frac{N}{N+k}\sigma^2$ and we should set $k=0$ for an unbiased estimate. So one can interpret the intuition behind $N-1$ as accounting for the fact that the true mean has been estimated by the sample mean (and hence is not "fixed").

For minimum mean square error, we require $MSE(s^2)=E\left[(s^2-\sigma^2)^2\right]$ to be a minimum for some choice of $k$. Expanding the mean square error we get:

$$MSE(s^2)=E[s^4]-2\sigma^2E[s^2]+\sigma^4$$

$E(s^2)$ has been calculated already, now to calculate $E(s^4)$. squaring $s^2$ gives us:

$$s^4=\frac{1}{(N+k)^2}\left[\sum_i x_i^2-N\bar{x}^{2}\right]^2$$ $$=\frac{1}{(N+k)^2}\left(\left[\sum_i x_i^2\right]^2-2\left[\sum_i x_i^2\right]\left[N\bar{x}^{2}\right]+\left[N\bar{x}^{2}\right]^2\right)$$ $$=\frac{\left[\sum_i x_i^4+\sum_{i\neq j} x_j^2x_i^2\right]-2\frac{1}{N}\left[\sum_i x_i^2\right]\left[\sum_i x_i^2+\sum_{i\neq j} x_jx_i\right]+\frac{1}{N^2}\left[\sum_i x_i^2+\sum_{i\neq j} x_jx_i\right]^2}{(N+k)^2}$$ $$=\frac{(\frac{1}{N}-1)^2\left[\sum_i x_i^4+\sum_{i\neq j} x_j^2x_i^2\right]+2\frac{1}{N}(\frac{1}{N}-1)\left[\sum_i x_i^2\right]\left[\sum_{i\neq j} x_jx_i\right]+\frac{1}{N^2}\left[\sum_{i\neq j} x_jx_i\right]^2}{(N+k)^2}$$ $$=\frac{f(x_1,\dots,x_N)}{(N+k)^2}$$

And you can see that without some assumptions that the expectation will in general be a function of the fourth order moments $E(x_ix_jx_kx_l),E(x_i^2x_jx_k),E(x_i^3x_j),E(x_i^2x_j^2),E(x_i^4)$ (which are not given in the question). However, its dependence on $k$ is quite simple, so we can still solve the variation problem algebraically with $F=E[f(x_1,\dots,x_N)]$. So we have:

$$MSE(s^2)=\frac{F}{(N+k)^2}-2\sigma^2\left[\frac{N-1}{N+k}\sigma^2(1-\rho)\right]+\sigma^4$$

Taking derivative with respect to $k$, set to zero, and solve for $k$:

$$-2\frac{F}{(N+k)^3}+2\frac{N-1}{(N+k)^2}\sigma^4(1-\rho)=0$$ $$\implies k=\frac{F}{(N-1)\sigma^4(1-\rho)}-N$$

This shows that unless $F=c\sigma^4(1-\rho)$, where $c$ depends only on the sample size, the optimum value of $k$ will be a function of the parameters, and hence you have no "solution" per se because it depends on things you don't know. You can show that if you assume an independent normal distribution for $x_i$ (so $\rho=0$), then $F=(N^2-1)\sigma^4$ and you get $k=+1$ as the optimum value.

• (+1) Nicely done. The generality of this approach reveals the essence of the problem. Note that one common situation is when the $x_i$ are independent, $\sigma$ is a scale parameter, and the fourth moment exists. Then, $F$ necessarily is a constant times $\sigma^4$ (and we can easily recover that constant from the kurtosis, the mgf, the cf, or the cumulant generating function). – whuber Jul 5 '11 at 16:27
• @whuber - that was my intention with my answer, to show that the rule is not as general as many people suppose. I didn't know that about $F$ though. – probabilityislogic Jul 6 '11 at 0:53

What's the problem with just dividing by $N$? You don't take into account that you aren't subtracting the true population mean off of each $x_i$, but rather than estimate of it.

One way that I like to think about it is, suppose that I gave you the sample mean $\bar{x}$. How many data points $N-k$ would I have to give you so that you could tell me the exact values of the remaining $k$? Well, a mean gives one equation involving every observation and, if there was only one unknown observation, we could solve the equation. In summary, knowing the mean and $N - 1$ data points is the same as knowing every data point.

In calculating the sample variance, I know the mean, so it's like I only have $N - 1$ effective data points giving me information; the last one I could guess using the others plus the mean. We divide by the number of effective data points that we have, $N-1$. This is known as a degrees-of-freedom correction (we only have $N-1$ degrees of freedom, parameters that we don't know, given that we know the mean).

Now, here's the math: Let $$\begin{equation*} s^2 = \frac{1}{N}\sum_{i=1}^N{(y_i - \bar{y})^2}. \end{equation*}$$ Then, \begin{align*} E[s^2] &= E\left[ \frac{1}{N}\sum_{i=1}^N{(y_i - \bar{y})^2} \right] \\ &= E\left[\frac{1}{N}\sum_{i=1}^N{\left((y_i - \mu) - (\bar{y} - \mu)\right)^2} \right] \\ &= \frac{1}{N}\sum_{i=1}^NE\left[(y_i - \mu)^2\right] - 2\frac{1}{N}E\left[(\bar{y} - \mu)\sum_{i=1}^N{(y_i - \mu)}\right] \\ &\qquad + E\left[(\bar{y} - \mu)^2\right] \\ &= \frac{1}{N}\sum_{i=1}^NE\left[(y_i - \mu)^2\right] - E\left[(\bar{y} - \mu)^2\right] \\ &= \text{Var}(y_i) - \frac{\text{Var}(y_i)}{N} = \frac{N-1}{N}\text{Var}(y_i). \end{align*}

This uses the fact that the variance of the sample mean is the variance of $y_i$ divided by $N$.

Hence, an unbiased estimator requires multiplying $s^2$ by $N/(N-1)$, giving the equation that you sought.

As you mention, dividing by $N-1$ is close to dividing by $N$; the two get close as $N$ gets big. Hence, $s^2$ is a consistent estimator---its bias goes to 0 as $N$ gets big.

• The penultimate paragraph of this reply answers the first question. Can you also address the second question? – whuber Jul 5 '11 at 13:40
• Looks like @probabilityislogic beat me to it. – Charlie Jul 5 '11 at 15:52