I want to prove that $V$ is an unbiased estimator of the covariance matrix $$(X'X)^{-1}(X'DX)(X'X)^{-1},$$ where $D=diag(\sigma^2,...,\sigma^2)=E(ee'|X)$ in a linear model.

$$V = \frac{n}{n-k}(X'X)^{-1}\left(\sum_{i=1}^n{X_iX'_i\hat{e}^2_i}\right)(X'X)^{-1}$$

To do so, I first find the conditional expectation of $V$.

$$E[V|X]= \frac{n}{n-k}(X'X)^{-1}\left(\sum_{i=1}^n{X_iX'_iE(\hat{e}^2_i}|X)\right)(X'X)^{-1}$$

However, I am not sure how to proceed from here.


1 Answer 1


I do not think your claim is correct. You are analyzing what MacKinnon and White (1985) refer to as the heteroskedasticity-robust variance estimator $HC_1$. I will argue that, instead, the estimator they call $HC_2$ is unbiased when in fact no heteroskedasticity is present, while $HC_1$ only is in a special case. Here is more detail:

Notice that the residuals satisfy $\hat e=My$ for $M=I-H$ and hat matrix $$H=\{h_{ij}\}_{i,j=1,\ldots,n}=\{x_i'(X'X)^{-1}x_j\}_{i,j=1,\ldots,n}=X(X'X)^{-1}X'.$$ Hence, $$ \hat e_i=(1-h_{ii})y_i-\sum_{j\neq i}h_{ij}y_j, $$ so that $$\hat e_i^2=\left[(1-h_{ii})y_i-\sum_{j\neq i}h_{ij}y_j\right]^2$$ In conditional expectation and under random sampling (i.i.d.) with conditional homoskedasticity - assumptions you do not state but that I need - $$ \begin{align*} E(\hat e_i^2|X)&=Var(\hat e_i|X) \end{align*}$$ in view of $E(\hat e|X)=ME(y|X)=MX\beta=0$. Now, again by the i.i.d. assumption,

\begin{align*} Var(\hat e_i|X)&=(1-h_{ii})^2Var(y_i|X)+\sum_{j\neq i}h_{ij}^2Var(y_j|X)\\ &=(1-h_{ii})^2\sigma^2+\sum_{j\neq i}h_{ij}^2\sigma^2, \end{align*}

as covariances are zero and all conditional variances are assumed to be the same. Next, notice that $$ \sum_{j=1}^nh_{ij}^2=h_{ii},$$ as, by symmetry of $H$, \begin{align*} \sum_{j=1}^nh_{ij}^2&=\sum_{j=1}^nh_{ij}h_{ji}\\ &=\sum_{j=1}^nx_i'(X'X)^{-1}x_jx_j'(X'X)^{-1}x_i\\ &=x_i'(X'X)^{-1}\sum_{j=1}^nx_jx_j'(X'X)^{-1}x_i\\ &=x_i'(X'X)^{-1}(X'X)(X'X)^{-1}x_i\\ &=x_i'(X'X)^{-1}x_i\\ &=h_{ii}\\ \end{align*}

Hence, by multiplying out and rearranging, \begin{align*} E(\hat e_i^2|X)&=(1-2h_{ii}+h_{ii}^2)\sigma^2+\sum_{j\neq i}h_{ij}^2\sigma^2\\ &=(1-h_{ii})\sigma^2-h_{ii}\sigma^2+\sum_{j=1}^nh_{ij}^2\sigma^2\\ &=(1-h_{ii})\sigma^2, \end{align*} Hence, under homoskedasticity, unbiasedness would call for adjusting squared residuals by $1-h_{ii}$, not multiplying by $n/(n-k)$.

Only in a "balanced design" (note that $\sum_ih_{ii}=k$, so that $k/n$ is the average value of the $h_{ii}$) in which $h_{ii}=k/n$ for all $i$ would the two coincide, as we would then have $$ \frac{1}{1-h_{ii}}=\frac{1}{1-k/n}=\frac{n}{n-k} $$ So all in all, for $$ HC_2=(X'X)^{-1}\left(\sum_{i=1}^n{x_ix'_i\frac{\hat{e}^2_i}{1-h_{ii}}}\right)(X'X)^{-1} $$ we obtain \begin{align*} E(HC_2|X)&=(X'X)^{-1}\left(\sum_{i=1}^nx_ix'_i\sigma^2\right)(X'X)^{-1}\\&=\sigma^2(X'X)^{-1} \end{align*} Hence, $HC_2$'s conditional expected value equals the conditional variance of $\hat\beta$ under homoskedasticity, which is different from the unconditional variance, see here.

Here is a little simulation to illustrate. I take the regressor to be fixed in repeated samples here to avoid the distinction between conditional and unconditional variance discussed above.

I find that the bias of $HC_2$ tends to be an order of magnitude smaller, although both biases are very small for the designs considered here.


n <- 50
x <- rnorm(n, sd=3)
x <- rt(n, df=2)
x <- runif(n, -10, 10)

sigma <- .2

mc.function <- function(n){
  u <- rnorm(n, sd=sigma)
  y <- 2*x + u
  limo <- lm(y~x-1)
  return(c(vcovHC(limo, "HC1"), vcovHC(limo, "HC2")))

true.cond.var <- sigma^2/sum(x^2)

vcovs <- replicate(1000, mc.function(n))
(bias <- rowMeans(vcovs-true.cond.var))
abs(bias[1]) > abs(bias[2])
  • $\begingroup$ I believe that when you sum the squared residuals, though, you obtain a quantity closely related to the trace of $H,$ which depends only on $n$ and $p,$ regardless of the design, thereby leading to the claim in the question. $\endgroup$
    – whuber
    Commented Nov 16, 2021 at 17:31
  • $\begingroup$ Yeah, $tr(H)=tr(X(X'X)^{-1}X')tr((X'X)^{-1}X'X)=k$ only depends on $k$, which gives rise to the unbiasedness of the standard estimator $\sum_i\hat e_i^2/(n-k)$ under homoskedasticity, eg stats.stackexchange.com/questions/76738/…. I am not sure how that applies here, though, since the heteroskedasticity robust estimator operates on the product $x_i\hat e_i$, so that we need to work on individual squared residuals, not their sum. See also econstor.eu/bitstream/10419/189084/1/qed_wp_0537.pdf $\endgroup$ Commented Nov 16, 2021 at 18:30
  • 1
    $\begingroup$ The part I did not catch--and now I see you asked for clarification--is whether the model is heteroscedastic or not. If, at the outset, you were to stipulate that you are discussing the heteroscedastic model, there would be less chance of misunderstanding. (+1) $\endgroup$
    – whuber
    Commented Nov 16, 2021 at 19:39

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.