# Covariance matrix of the residuals in the linear regression model

I estimate the linear regression model:

$$Y = X\beta + \varepsilon$$

where $$y$$ is an ($$n \times 1$$) dependent variable vector, $$X$$ is an ($$n \times p$$) matrix of independent variables, $$\beta$$ is a ($$p \times 1$$) vector of the regression coefficients, and $$\varepsilon$$ is an ($$n \times 1$$) vector of random errors.

I want to estimate the covariance matrix of the residuals. To do so I use the following formula:

$$Cov(\hat{\varepsilon}) = \sigma^2 (I-H)$$

where $$\hat{\varepsilon}=Y-X\hat{\beta}$$, $$\sigma^2$$ is estimated by $$\hat{\sigma}^2 = \frac{e'e}{n-p}$$, $$I$$ is an identity matrix, and $$H = X(X'X)^{-1}X'$$ is a hat matrix (see Kutner, 2005).

However, in some source I saw that the covariance matrix of the residuals is estimated in other way. The residuals are assumed to follow $$AR(1)$$ process:

$$\varepsilon_t = \rho \varepsilon_{t-1} + \eta_t$$

where $$E(\eta) = 0$$ and $$Var({\eta}) = \sigma^2_{0}I$$.

The covariance matrix is estimated as follows

$$Cov(\varepsilon) = \sigma^2 \begin{bmatrix} 1 & \rho & \rho^2 & ... & \rho^{n-1}\\ \rho & 1 & \rho & ... & \rho^{n-2} \\ ... & ... & ... & ... & ... \\ \rho^{n-1} & \rho^{n-2} & ... & ... & 1 \end{bmatrix}$$

where $$\sigma^2 = \frac{1}{1-\rho^2}\sigma^2_0$$

My question is are there two different specifications of the covariance matrix of residuals or these are somehow connected with each other?

• Connected is a pretty vague term but I would claim that they're dis-connected. The first method assumes that the disturbances are iid and normally distributed. The second method assumes that they are autocorrelated and follow an AR(1) process. How you estimate the regression model depends on which model you assume for the disturbances so they are connected only in that sense. May 24, 2020 at 1:35
• The first one is suitable for independent observations while the second one for serial observations on the same sampling unit. May 24, 2020 at 9:04
• @mlofton: does the first method really assume independance of the disturbances ? If it were the case, an estimator of the covariance matrix of errors should just be $\hat {\sigma} I$, but here, the matrix $I - H$ is not necessarily diagonal.. May 25, 2020 at 11:18
• @CherryGarcia: what assumptions are you willling to make on your errors $\varepsilon$ ? Maybe looking at Feasible Generalized Least Squares methods could help.. May 25, 2020 at 11:20
• @Pohoua: Definitely the assumption is independence but the use of $(I-H)$ has something to do ( waving my hands here because I forget so hopefully someone else can explain more clearly ) wiith the fact that the estimates are not independent in that the sum of the residuals estimates has to equal zero because of the nature of OLS. That still doesn't explain why you get a non diagonal estimate but it's related to that. I'd love to know the answer myself so hopefully someone else can explain. Thanks for great question. May 25, 2020 at 16:05

After some investigation, I think I found a small (but crucial!) imprecision in what your post.

The first formula you wrote : $$var(\varepsilon) = \sigma^2 (I - H)$$ is actually not totally exact. The formula should be $$var(\hat \varepsilon) = \sigma ^2 (I - H)$$ where $$\hat\varepsilon = Y - \hat\beta X$$ considering the OLS estimator $$\hat\beta = (X^TX)^{-1}X^TY$$. Thus $$\hat\sigma(I - H)$$ is an estimator of the variance of the estimated residuals associated with OLS estimator. This formula does not suppose independance of the $$\varepsilon_i$$, just that they all have same variance $$\sigma^2$$. But this is not what you want! You want an estimate of the variance of the true residuals, not the estimated residuals under OLS estimation. OLS estimator corresponds to maximum likelihood estimator under the hypothesis that residuals are i.i.d. and normal. The estimated residuals can thus be very poor estimates of the true residuals if these hypothesis are not met, and there covariance matrix can be very different from the covariance of the true residuals.

The second formula you wrote does correspond to the covariance matrix of the $$\varepsilon_i$$ under the hypothesis that they follow an AR(1) process.

Estimating covariance matrix of the residuals of a linear regression without any asumption cannot easily be done: you would have more unknown than datapoints... So you need to specify some form for the covariance matrix of the residuals. Supposing that they follow an AR(1) process (if this is relevent) is a way of doing so. You can also assume that they have a stationnary parametrized autocorrelation function, whose parameters you can estimate, and use it to deduce the covariance matrix.

• @Pohuoa: Thanks for explaining the $(I-H)$ issue. Now I get it that you don't estimate the covariance matrix of the residuals. The $(I-H)$ matrix is estimating something totally different. Cherry Garcia: The estimates of the residuals in the first case ( where you assume they are iid ) are just the $\hat{\epsilon}$ that come out of the regression. In the second case, you're specifying a model for the residuals, namely an AR(1). So, the estimation method needs to know that so that is why you use feasible generalized least squares. There are other methods also besides FGLS. May 26, 2020 at 3:54
• @Pohuoa: Thank you very much for spotting the mistake; I have adjusted my question accordingly. I think that I understand what you mean, but to be sure, let me summarize: 1) Under Gauss-Markov assumptions the covariance matrix of the residuals is $Cov(\varepsilon) = \sigma^2I$ 2) Assuming residuals follow AR(1) process the covariance matrix of the residuals is the second case I described in the question. 3) The covariance matrix of OLS residuals is: $Cov(\hat{\varepsilon})=\hat{\sigma}^2(I-H)$ May 27, 2020 at 19:01
• @CherryGarcia: Yes, exactly. May 27, 2020 at 21:20
• don't mistake residuals for errors. residuals $r$ (that you call estimates of the residuals) are the difference between the real and the fitted y values, you don't have to estimate them, you have them. the errors $\epsilon$ instead are the random part of the data generating process. those are assumed to be IID, not residuals. a good estimate of errors are the studentized residuals (simple residuals are heteroshedastic and hence biased estimators) May 28, 2020 at 11:30

In basic OLS you don't estimate the covariance matrix of residuals. You assume that errors (not residuals) are spherical, meaning that they're not correlated with each other. Residuals will come out of OLS uncorrelated.

What you described as a second method is a different assumption. When applying basic OLS to time series you run into an issue that its assumptions are not practical. In time series the residuals are often correlated. So, you could assume that they're AR(1) process, and that what that method does: it estimates the model assuming the errors are AR(1). This is called feasible generalized least squares

• well but you do compute the covariance matrix of residuals in OLS, that's $\sigma^2(I-H)$ May 28, 2020 at 11:35
• but nothing of $Y$ goes in there other than the residual variance so I kind of agree. Feb 1 at 6:28