# Standard Error of Noise Variance in Least Squares

Suppose I'm doing ordinary least squares with homoskedastic errors, so something like:

$$y = X \beta + \epsilon$$

where $$\epsilon \sim \mathcal{N}(0,\sigma^2)$$.

I know how to estimate the expected value of $$\sigma^2$$, but I'm having trouble identifying how to compute a standard error on that estimate. Something analogous to what's shown on page 47 (pdf page 58) of this book.

• What do you mean by "the expected value of $\sigma^2$"? $\sigma^2$ is a non-random constant. Do you try to say how to compute the standard error of the estimate of $\sigma^2$? Feb 23 at 20:12
• @Zhanxiong it sounds like they are interested in the sampling distribution of the model variance $\sigma^2$. As in stats.stackexchange.com/q/121662/36229 Feb 23 at 20:15
• @shadowtalker I guess so, but a more accurate way to say it is the sampling distribution of $\hat{\sigma}^2$ (which is random), not $\sigma^2$ (which is fixed). Feb 23 at 20:17
• Very closely related questions: stats.stackexchange.com/q/224042/36229 and stats.stackexchange.com/q/20227/36229 and stats.stackexchange.com/q/312337/36229 Feb 23 at 20:25

The other answer here shows you the variance of the estimator of the error variance under the assumption that the errors are normally distributed, which is the specified form in the Gaussian regression model. I would counsel against using this method, since it is often the case that the actual data we use in a regression model does not conform to this assumption. In particular, it is not unusual for the residuals in regression data to exhibit leptokurtosis (or sometimes platykurtosis), contrary to the model assumptions. Many aspects of regression analysis are robust to the normality assumption, but this aspect of the model is not robust --- it is a case where the asserted formula hinges on the Gaussian model assumption rather than being determined by actual analysis of the data used in the model. In general, this to be avoided if we want to use robust methods.

Suppose we are willing to generalise the model slightly by no longer assuming a normal distribution for the error term (we will still assume that the errors are IID with zero mean and fixed variance). In the generalised case we may consider error variables with kurtosis $$\kappa$$ (but still assuming a simple linear regression), in which case the variance of $$\hat{\sigma}^2$$ is actually given by:

$$\mathbb{V}(\hat{\sigma}^2) = \bigg( \kappa - \frac{n-p-4}{n-p-2} \bigg) \frac{\sigma^4}{n-p-1}.$$

In the special case of a mesokurtotic error distribution (e.g., the normal distribution) we have $$\kappa = 3$$ and so the variance formula reduces to the more familiar form:

$$\mathbb{V}(\hat{\sigma}^2) = \frac{2 \sigma^4}{n-p-1}.$$

Now, it is possible to estimate the kurtosis of the error distribution from the residuals in the regression model, so in principle it is possible to estimate the variance of the estimator for the error variance in a way that does not hinge on the assumption of a mesokurtic error distribution. For example, if you have an estimator $$\hat{\kappa}$$ for the kurtosis of the error distribution (e.g., from a kurtosis estimator using the residuals) then you could estimate the standard error of the estimator for the error variance as:

$$\hat{\text{se}}_n = \sqrt{\hat{\mathbb{V}}(\hat{\sigma}^2)} = \sqrt{ \hat{\kappa} - \frac{n-p-4}{n-p-2}} \times \frac{\hat{\sigma}^2}{\sqrt{n-p-1}}.$$

Assuming you use a consistent estimator $$\hat{\kappa}$$, this latter estimator will be a consistent estimator of the true standard deviation of the estimator of the error variance. For more details on the moments of these types of sampling quantities you might find it useful to consult O'Neill (2014) (that paper deals with standard sampling quantities outside of regression, but its results can easily be adapted to the regression setting).

• Thank you, that O'Neill review is an extremely handy reference!
– Ben
Feb 28 at 22:06

$$\DeclareMathOperator{\Var}{Var}$$

Suppose the design matrix is $$N \times (p + 1)$$, I am assuming you want to calculate the variance of $$\hat{\sigma}^2 = \frac{1}{N - p - 1}\epsilon'(I - H)\epsilon$$, which is an unbiased estimator of $$\sigma^2$$, where $$I$$ is the order $$N + 1$$ identity matrix and $$H = X(X'X)^{-1}X'$$. In this thread, it has been shown that \begin{align} (N - p - 1)\hat{\sigma}^2/\sigma^2 \sim \chi_{N - p - 1}^2, \end{align} by which and the variance of a $$\chi_k^2$$ r.v. is $$2k$$ follows that \begin{align} \Var(\hat{\sigma}^2) = \frac{\sigma^4}{(N - p - 1)^2} \times 2(N - p - 1) = \frac{2\sigma^4}{N - p - 1}. \end{align}

Estimate the $$\sigma^4$$ in the numerator by $$\hat{\sigma}^4$$, it follows that the standard error of $$\hat{\sigma}^2$$ is $$\frac{\sqrt{2}\hat{\sigma}^2}{\sqrt{N - p - 1}}$$.