Deriving the HC0 Robust Covariance Estimate

Why Derive an HC0 Estimate?

Ordinary least squares (OLS) estimates regression coefficients by minimizing squared residuals. When the errors have unequal variances, OLS can still estimate the coefficients appropriately under a correct conditional mean model and mean-zero errors, but the usual covariance calculation based on one common error variance is no longer appropriate. The sandwich covariance estimator provides a way to account for unequal variances without specifying their pattern.

The sandwich formula depends on the unknown variances of the individual errors. HC0 makes the formula usable by replacing each unknown variance with the corresponding squared OLS residual. The construction is straightforward, but its meaning depends on the model assumptions, and its matrix dimensions are worth checking carefully.

Model, Assumptions, and the Exact Covariance

Write the multiple linear regression model in matrix form as \(y=X\beta+\varepsilon\). Here \(y\) is an \(n\times 1\) vector of responses, \(X\) is an \(n\times p\) design matrix, \(\beta\) is a \(p\times 1\) vector of coefficients, and \(\varepsilon\) is an \(n\times 1\) vector of errors. Assume \(X\) has full column rank, so \(X'X\) is invertible. The \(i\)th row of \(X\) is \(x_i'\), where \(x_i\) is a \(p\times 1\) vector.

Model assumptions for the HC0 derivation: Conditional on the observed design matrix \(X\), the conditional mean is \(E(\varepsilon\mid X)=0\), and the errors are uncorrelated with variances that may differ by observation: \[ \operatorname{Var}(\varepsilon\mid X)=\Omega =\operatorname{diag}(\sigma_1^2,\ldots,\sigma_n^2). \] The conditional mean of the response is therefore \(E(y\mid X)=X\beta\). No common error variance is assumed.

The OLS estimator is \(\hat{\beta}=(X'X)^{-1}X'y\). Substituting \(y=X\beta+\varepsilon\) and using \((X'X)^{-1}X'X=I\) gives \[ \hat{\beta}-\beta=(X'X)^{-1}X'\varepsilon. \] Conditional on \(X\), its covariance is therefore \[ \operatorname{Var}(\hat{\beta}\mid X) =(X'X)^{-1}X'\Omega X(X'X)^{-1}. \] This is the unequal-variance sandwich covariance formula. It is exact under the stated conditional mean and covariance assumptions, but it cannot be directly calculated unless \(\Omega\) is known.

Worked Example: From Known Variances to a Residual-Based Estimate

Consider a hypothetical regression with an intercept and one predictor, whose values are \(0,1,2,3,4\). Suppose the response values are \(11,9,14,11,15\). For illustration, also suppose the conditional error variances are known to be \(1,1,4,4,9\). First calculate the exact sandwich covariance that would apply if those variances were known, then calculate the HC0 estimate from the fitted residuals.

The design matrix has five rows and two columns: \[ X= \begin{pmatrix} 1&0\\ 1&1\\ 1&2\\ 1&3\\ 1&4 \end{pmatrix}, \qquad X'X= \begin{pmatrix} 5&10\\ 10&30 \end{pmatrix}. \] The determinant of \(X'X\) is \(5(30)-10(10)=50\), so \[ (X'X)^{-1} =\frac{1}{50} \begin{pmatrix} 30&-10\\ -10&5 \end{pmatrix} = \begin{pmatrix} 0.6&-0.2\\ -0.2&0.1 \end{pmatrix}. \]

The variance matrix is \(5\times5\), while the middle matrix \(X'\Omega X\) is \(2\times2\). Its entries are the sums of the variances, predictor values times variances, and squared predictor values times variances: \[ X'\Omega X= \begin{pmatrix} 1+1+4+4+9&0(1)+1(1)+2(4)+3(4)+4(9)\\ 0(1)+1(1)+2(4)+3(4)+4(9)&0^2(1)+1^2(1)+2^2(4)+3^2(4)+4^2(9) \end{pmatrix} = \begin{pmatrix} 19&57\\ 57&197 \end{pmatrix}. \] Thus the exact covariance is \[ \begin{pmatrix} 0.6&-0.2\\ -0.2&0.1 \end{pmatrix} \begin{pmatrix} 19&57\\ 57&197 \end{pmatrix} \begin{pmatrix} 0.6&-0.2\\ -0.2&0.1 \end{pmatrix} = \begin{pmatrix} 1.04&-0.52\\ -0.52&0.45 \end{pmatrix}. \]

In practice, the variances are not known. For these hypothetical observations, OLS gives \(\hat{\beta}=(10,1)'\), and the residual vector is \(\hat{\varepsilon}=(1,-2,2,-2,1)'\). Squaring its entries gives \(1,4,4,4,1\), so the HC0 middle matrix is \[ X'\operatorname{diag}(1,4,4,4,1)X = \begin{pmatrix} 14&28\\ 28&72 \end{pmatrix}. \] The resulting estimate is \[ \widehat{\operatorname{Var}}_{\mathrm{HC0}}(\hat{\beta}) = \begin{pmatrix} 0.6&-0.2\\ -0.2&0.1 \end{pmatrix} \begin{pmatrix} 14&28\\ 28&72 \end{pmatrix} \begin{pmatrix} 0.6&-0.2\\ -0.2&0.1 \end{pmatrix} = \begin{pmatrix} 1.20&-0.32\\ -0.32&0.16 \end{pmatrix}. \] The exact covariance and its HC0 estimate differ here because the squared residuals are estimates, not the actual error variances.

Constructing HC0

Let \(\hat{\varepsilon}=y-X\hat{\beta}\) be the vector of OLS residuals, with \(i\)th entry \(\hat{\varepsilon}_i\). HC0 uses the diagonal matrix of their squares as an estimate of the unknown diagonal matrix \(\Omega\):

Definition: The HC0 estimate of the OLS coefficient covariance matrix is \[ \widehat{\operatorname{Var}}_{\mathrm{HC0}}(\hat{\beta}) =(X'X)^{-1}X'\operatorname{diag}(\hat{\varepsilon}_1^2,\ldots,\hat{\varepsilon}_n^2)X(X'X)^{-1}. \] Equivalently, writing \(x_i\) as a column vector, \[ \widehat{\operatorname{Var}}_{\mathrm{HC0}}(\hat{\beta}) =(X'X)^{-1} \left(\sum_{i=1}^{n}\hat{\varepsilon}_i^2x_ix_i'\right) (X'X)^{-1}. \]

The two expressions are identical: the \(i\)th diagonal entry in the middle matrix weights the \(i\)th row of \(X\), giving the outer product \(\hat{\varepsilon}_i^2x_ix_i'\). The middle matrix is often called the “meat” of the sandwich; the two copies of \((X'X)^{-1}\) form its sides.

The dimensions confirm that the product is a covariance matrix for the \(p\) coefficients. The matrix \(X'\) is \(p\times n\), the diagonal residual-square matrix is \(n\times n\), and \(X\) is \(n\times p\); hence the middle product is \(p\times p\). Each \((X'X)^{-1}\) is also \(p\times p\), so the complete product is \(p\times p\). No \(n\times n\) matrix is being reported as the coefficient covariance.

Worked Example: Checking the HC0 Matrix Dimensions

Suppose a hypothetical model has \(n=8\) observations and \(p=3\) coefficients, including an intercept. For observation \(i\), let the residual be \(\hat{\varepsilon}_i=2\) and the predictor vector be \(x_i=(1,3,-1)'\). Calculate that observation’s contribution to the HC0 middle matrix.

The contribution is the residual square times the predictor outer product: \[ \hat{\varepsilon}_i^2x_ix_i' =4 \begin{pmatrix} 1\\3\\-1 \end{pmatrix} \begin{pmatrix} 1&3&-1 \end{pmatrix} = 4 \begin{pmatrix} 1&3&-1\\ 3&9&-3\\ -1&-3&1 \end{pmatrix} = \begin{pmatrix} 4&12&-4\\ 12&36&-12\\ -4&-12&4 \end{pmatrix}. \] This is \(3\times3\), as required. Adding the corresponding \(3\times3\) contribution from each of the other seven observations gives the full middle matrix; multiplying on both sides by \((X'X)^{-1}\) gives the \(3\times3\) HC0 estimate. The full design matrix is \(8\times3\), and \(X'X\) and its inverse are \(3\times3\).

Why Squared Residuals Are a Plug-In Estimate

For the true errors, the \(i\)th diagonal entry of \(\Omega\) is \(\sigma_i^2\). The \(i\)th squared OLS residual is not generally an unbiased estimate of that individual variance. To see why, define the hat matrix \(H=X(X'X)^{-1}X'\) and write \(M=I-H\). Since \(\hat{\varepsilon}=M\varepsilon\), the residuals are influenced by the fitted regression as well as by the errors. Conditional on \(X\), \[ E(\hat{\varepsilon}_i^2\mid X) =\sum_{j=1}^{n}M_{ij}^2\sigma_j^2, \] which is generally not equal to \(\sigma_i^2\). OLS projection reduces the total squared residual norm relative to the error vector, but individual residuals need not be smaller than their corresponding errors; as the preceding formula shows, their variances depend on the design and on all the error variances.

So HC0 is a plug-in estimate of the sandwich covariance, not an exact finite-sample identity and not a claim that each residual square recovers its own error variance. Under standard large-sample regularity conditions—such as a fixed number of coefficients, suitable finite error moments, a well-behaved design, and no observation dominating the fit—the residual-based middle matrix can estimate the needed variance contribution. The exact strength of such an approximation depends on the design and error conditions.

Worked Example: HC0 in a One-Coefficient Model

For a compact arithmetic check, consider a hypothetical regression through the origin, \(y_i=\beta x_i+\varepsilon_i\), with predictor values \(x=(1,2,3)'\) and responses \(y=(1,-2,1)'\). This model has no intercept; that restriction would need a substantive justification in an application. Calculate the OLS coefficient, residuals, and HC0 covariance estimate.

Here \(X\) is \(3\times1\), \(X'X=1^2+2^2+3^2=14\), and \(X'y=1(1)+2(-2)+3(1)=0\). Therefore \[ \hat{\beta}=(X'X)^{-1}X'y=\frac{0}{14}=0, \] and the fitted values are all zero, leaving residuals \((1,-2,1)'\). The HC0 middle term is a scalar: \[ X'\operatorname{diag}(1,4,1)X =1^2(1)+2^2(4)+3^2(1) =1+16+9=26. \] Thus \[ \widehat{\operatorname{Var}}_{\mathrm{HC0}}(\hat{\beta}) =\frac{1}{14}(26)\frac{1}{14} =\frac{26}{196} =\frac{13}{98} \approx 0.1327. \] A one-coefficient model has a \(1\times1\) covariance matrix, so the scalar result is dimensionally consistent.

What HC0 Does—and Does Not—Address

HC0 changes the estimated coefficient uncertainty to allow the error variance to differ across observations. It does not change the OLS coefficient estimates: the same \(\hat{\beta}\) is used. Nor does it correct a misspecified conditional mean. If the expected response is curved but the fitted mean is forced to be linear, substituting squared residuals into the sandwich does not make that mean model correct.

The diagonal form of \(\Omega\) also matters. The derivation assumes errors are uncorrelated conditional on \(X\). If errors are correlated across observations, the conditional covariance matrix is not diagonal, and this HC0 expression does not estimate the full sandwich middle term. The assumptions must match the data structure; the word “robust” does not mean robust to every modeling failure.

Common Mistakes

Using residuals instead of squared residuals. The middle matrix estimates variance contributions, so its diagonal weights are \(\hat{\varepsilon}_i^2\), not \(\hat{\varepsilon}_i\).
Putting the residual-square matrix in the wrong dimensions. It is \(n\times n\), while the resulting coefficient covariance is \(p\times p\).
Calling HC0 the exact covariance. The exact sandwich contains the unknown \(\Omega\); HC0 substitutes residual squares and is an estimate.
Assuming “robust” means the mean model can be wrong. HC0 addresses unequal error variances under the stated regression mean and error assumptions; it does not repair curvature or omitted structure in the conditional mean.

Key takeaway: HC0 is obtained by replacing \(\Omega=\operatorname{diag}(\sigma_1^2,\ldots,\sigma_n^2)\) in the unequal-variance OLS sandwich with \(\operatorname{diag}(\hat{\varepsilon}_1^2,\ldots,\hat{\varepsilon}_n^2)\). The result is a \(p\times p\) estimate of coefficient covariance.

Check Your Understanding

Check that you can connect the assumptions, substitution, and matrix structure.

State the conditional mean and conditional covariance assumptions used in deriving the HC0 estimate.
If \(X\) is \(n\times p\), what are the dimensions of \(\operatorname{diag}(\hat{\varepsilon}_1^2,\ldots,\hat{\varepsilon}_n^2)\) and of the final covariance estimate?
Why is \(\hat{\varepsilon}_i^2\) generally not an unbiased estimate of \(\sigma_i^2\) in a fitted regression?
In the outer-product expression, what is the \(i\)th observation’s contribution to the middle matrix?
Name one important failure HC0 does not fix, even if its covariance estimate is calculated correctly.