The bread and the filling

Fit a straight line to eighty rows whose error variance is largest at the ends of the design, and read each candidate standard error against the spread the slope actually has across twenty thousand draws. The model-based variance estimate averages 0.6081 of that spread. The unbiased robust correction averages 0.9961 of it. One of the two is a standard error, and it is not the one the regression prints by default.

The usual account of why stops at the word heteroskedasticity and moves on. It is worth going one step further, because the step contains a result most readers carry backwards: the robust standard error is not systematically the larger of the two. At the setting above it is 1.2806 times the model-based one. Move the error variance from the edges of the design into its middle and the same comparison reads 0.8246 — the model-based error is a fifth larger than the one that is right, and quoting the robust error as a conservative margin would there be quoting the smaller number.

That is a statement about a ratio, and the ratio has a closed form with exactly one design quantity in it. What follows computes it three ways, checks the three against each other, and then sweeps the dial that reverses the sign.

The bread is the design; the filling is what the errors did #

The sandwich estimator is $A^{-1}BA^{-1}$, and the two letters are doing different jobs. The bread $A = E[xx']$ is the design’s own second-moment matrix — the thing least squares inverted in order to produce a fit at all — and it knows nothing about the response. The filling $B = E[\varepsilon^2 xx']$ is the same matrix with each observation weighted by how large its error was.

The model-based variance is what the sandwich collapses to when the error variance is one number. Write it $\sigma^2 A^{-1}$: the errors’ contribution factors out of $B$ entirely, one $A$ cancels, and what is left is the familiar $\sigma^2/\sum(x_i - \bar{x})^2$ for a slope. The model-based estimate is not a different estimator; it is the sandwich under an assumption, which is why the two agree exactly where the assumption holds and why the disagreement can point either way when it does not.

The disagreement’s size is decided by how the error variance lines up with the design. Write the variance as $\sigma^2(x) = \sigma_0^2\left(1 + \gamma\left(u^2/v - 1\right)\right)$, with $u$ the covariate less its mean and $v$ its variance. The average of that is $\sigma_0^2$ for every $\gamma$, so the dial changes the shape of the error variance and never its level: a study with $\gamma = 0.8$ and one with $\gamma = -0.4$ have exactly as much noise in them. What differs is where the noise sits, and $\gamma > 0$ puts it where the leverage is.

Carrying that through the two matrices gives the slope’s variance ratio as

$$\frac{V_{\text{sandwich}}}{V_{\text{model}}} = 1 + \gamma(\kappa - 1),$$

with $\kappa$ the covariate’s own kurtosis. Nothing about the amount of noise survives into that expression. The size of the correction is a property of the covariate’s tail and of one shape parameter, and no part of it is a property of how noisy the data is.

Three routes, and they agree #

A ratio is easy to get right by accident, so the level is computed three ways that share no arithmetic.

The first is the spread itself: the variance of twenty thousand fitted slopes, multiplied by the sample size, which comes to 4.89053. That reads no formula at all. The second is the population sandwich, an expectation of a polynomial in the centred covariate against the design’s own central moments, which comes to 4.92000 and reads no data at all. The third is the model-based variance under the same population, 3.00000. The first two agree to within six parts in a thousand, which is about what twenty thousand draws buys; the third is a different quantity and misses the spread by a third.

Reading each estimator against the spread rather than against a formula is what makes the table interpretable. The model-based estimate averages 0.6081 of the variance the slope has, so its standard error is 0.7798 of the one it should report. HC0 reads 0.9576, HC1 0.9821, HC2 0.9961 and HC3 1.0362 — the four differ by less than a twentieth here, and the field’s account of which correction to use is about a design where they differ by a factor of sixteen.

The ratio of the two standard errors follows: 1.2799 counted and 1.2806 in closed form, which is $\sqrt{1.64} = \sqrt{1 + 0.8 \times 0.8}$ on a uniform covariate’s kurtosis of 1.8000. Two routes, seven digits of input in common — the design and the dial — and no shared arithmetic between them.

The ratio is the covariate’s fourth moment #

The closed form says the correction reads the design’s kurtosis, and four designs are enough to see that it reads nothing else.

At $\gamma = 0.8$ the even spread over the doubled unit interval gives a sandwich variance of 4.92000 against a model-based 3.00000. The same spread shifted gives the same pair. The same spread widened gives 1.23000 against 0.75000 — a quarter of each, because widening the design quarters the variance of a slope — and therefore the identical ratio of 1.6400. The exponential design gives 7.40000 against 1.00000, a ratio of 7.4000, on a kurtosis of 9.0000 against the uniform’s 1.8000. Every one of those four matches $1 + \gamma(\kappa - 1)$ to twelve places.

The practical reading is that the amount a robust standard error can differ from a model-based one is set before the response is looked at. A covariate with a heavy tail admits a large correction; a covariate spread evenly admits at most a modest one, since $\kappa - 1 = 0.8$ bounds the whole effect at $1 + 0.8\gamma$. It is the same quantity in the same position as the leverage in the essay on what one observation does to a slope, seen as a moment rather than as a point.

The direction is the sign of one dial #

Sweep $\gamma$ from $-0.4$ to $0.8$ at eighty rows, five thousand draws a setting, and the ratio of the two standard errors traces out $\sqrt{1 + \gamma(\kappa - 1)}$: 0.8246, 0.9165, 1.0000, 1.0770, 1.1489, 1.2166, 1.2806. The counted ratios are 0.8271, 0.9173, 0.9996, 1.0760, 1.1477, 1.2154 and 1.2798 — within a thousandth at every setting. Two of the seven sit below one and four above.

The two below one are the point of the sweep. At $\gamma = -0.4$ the error variance is largest in the middle of the design, where the leverage is least, so the observations that move the slope most are the quietest ones and the slope is more stable than a constant-variance calculation says. The model-based standard error is then 1.2222 times the one the slope’s own spread calls for — it over-states the uncertainty by a fifth — while the robust one is 1.0109 times it.

Across all seven settings the robust standard error stays between 1.0054 and 1.0109 of the right one while the model-based one runs from 1.2222 down to 0.7856. That is the honest summary of what a robust standard error buys: not a larger number, but the same number under a change of shape that moves the other one by fifty per cent of its own value.

Two readings follow, and both cut against habit. A robust standard error is not conservative. It is an estimate of a different quantity that happens to be the right one, and on a study whose noisy observations sit in the middle of its design it will report a narrower interval than the default. And an analyst who switches to robust errors and finds the numbers barely moved has learned something specific: not that the model is fine, but that $\gamma(\kappa - 1)$ is small on this design, which can happen because the error variance is nearly constant or because the covariate has no tail for it to lean on.

There is a diagnostic buried in that, and it is a poor one. The ratio of the two standard errors on a single dataset is an estimate of $\sqrt{1 + \gamma(\kappa - 1)}$, so an analyst who computes both is holding a noisy reading of exactly the quantity that decides whether the default was safe. The trouble is that both numerator and denominator are noisy at the sample sizes where the answer matters, and the ratio is near one under a wide band of departures — at $\gamma = 0.2$ the closed form is 1.0770, which a reader would round to “no difference” while the model-based standard error is already 0.9354 of the one it should report. A test built on that comparison is the check before the standard error in a different subject, and it has the same defect there: the departures large enough to detect reliably are the ones large enough to see, and the ones that quietly cost coverage are the ones the check misses.

Where they agree exactly, and what that costs #

Set the dial to nought and the two estimators are estimating the same thing. The sweep confirms it at every level: the scaled counted variance is 2.97826 against a population sandwich of 3.00000, the counted ratio of the standard errors is 0.9999 against a closed form of 1.0000, and the model-based estimate averages 1.0069 of the spread against HC2’s 1.0066.

The last pair is the one worth staring at, because it is where the case for a default comes from. Under a constant error variance the two estimates are equally right on average and they are not equally precise. The model-based estimate’s own coefficient of variation across the twenty thousand draws is 0.1608; HC0’s is 0.2156 and HC3’s is 0.2175. The robust estimate is right and noisy; the model-based one, in the one case where it is unbiased, is right and steady. That is the whole of the trade, and the essay that prices it where the risk is absent counts what the noise costs.

The same column at $\gamma = 0.8$ says the trade does not reverse. The model-based estimate’s coefficient of variation is 0.1976 against HC2’s 0.2615, so it is still the steadier estimate — it is simply steady around the wrong number. A precise estimate of the wrong quantity is the failure mode this whole field is about, and it is not visible in any diagnostic that reads only the estimate’s own variability.

One squared residual for every variance #

There is an objection to the whole construction that deserves an answer rather than a shrug, and it is usually raised as a rhetorical question: the filling asks for the error variance at every row, and there is exactly one observation at every row, so how can anything be estimated at all?

The answer is that the filling is never estimated. What is estimated is a single weighted average of it — $\hat{B} = \frac{1}{n}\sum e_i^2 x_i x_i'$ is a sum of $n$ terms whose expectation is the population filling, and a sum of $n$ noisy terms converges even though not one of the terms individually estimates anything. Each $e_i^2$ is a hopeless estimate of $\sigma^2(x_i)$ — it has one degree of freedom and a coefficient of variation of about $\sqrt{2}$ — and the average of eighty of them is a perfectly good estimate of the one number the slope’s variance needs.

That is why the robust estimate’s own coefficient of variation is 0.2604 at eighty rows rather than the $\sqrt{2}$ of a single squared residual: the averaging has done its work. It is also why the noise does not vanish the way the model-based estimate’s does. The model-based estimate pools every residual into one $s^2$ with $n - 2$ degrees of freedom and then multiplies by a known matrix; the robust one lets each residual keep its own position in the design, so the rows with the most leverage — the rows whose squared residuals matter most to the sum — are exactly the rows there are fewest of. The estimate is right on average and is carried by its extremes, which is the same structure the essay on which resampling keeps a non-constant error variance finds when it asks what a bootstrap must preserve.

And the residual is not the error. The fit removes part of every error, and it removes more of it where the leverage is high, so $E[e_i^2] = \sigma^2(1 - h_{ii})$ rather than $\sigma^2$ — the exact statement of the essay that separates residuals from errors. The uncorrected filling is therefore short by construction before any noise is considered, which is the whole content of the four corrections and the reason the choice between them is a question about the design.

What could have produced this without the claim being true #

The estimand could have moved. The sweeps in this essay hold the fitted model’s mean correct: the truth is a straight line and only the error variance departs from what least squares assumes. That is deliberate. The projection a wrong mean converges to moves the target as well as the spread, and mixing the two would leave a reader unable to say which defect produced the number. Here the target is the truth’s own slope, the coverage question is clean, and every ratio above is a statement about variance alone.

The residuals could be carrying the misfit rather than the errors. Under the mean model used here they are not, because there is no misfit. Under a curved truth they would be: the filling of a robust estimate is the squared residual, and a squared residual contains the systematic arch as well as the error. That arch is not small — on the even design it runs to 0.3234 at both ends — and it enters the filling as though it were noise, which means a robust standard error computed on a misspecified mean is estimating the spread of an estimator around a projection, not the spread of anything around the truth.

The count could be an artefact of eighty rows. It is not, in the sense that the closed forms are population quantities and agree with the count, and it is in the sense that eighty rows is a size at which the estimator has a small-sample defect of its own: HC0 reads 0.9576 rather than one here, and at twenty rows the field’s coverage sweep shows that shortfall growing to a fifth. The ratio in this essay is asymptotic and the estimate of it is not, and keeping the two apart is the whole subject of the coverage measurement.

The dial could be doing two things at once. It is not: $\sigma^2(x) = \sigma_0^2(1 + \gamma(u^2/v - 1))$ integrates to $\sigma_0^2$ against the design for every admissible $\gamma$, which is checked rather than asserted, and the machinery refuses a $\gamma$ steep enough to take the variance negative somewhere inside the design. A dial that changed the level as well as the shape would produce the same picture for a reason nobody could disentangle.

Where the arithmetic runs out #

The collapsed form $1 + \gamma(\kappa - 1)$ is exact for a symmetric design under a straight-line mean, and it is the form the sweeps here are checked against. The general sandwich is $A^{-1}BA^{-1}$ with no such collapse, and on a skewed design under a curved truth the slope element picks up cross-terms between the misfit and the variance shape that no single moment summarises. The full computation still runs — it is the same polynomial expectation — but there is no two-symbol answer to quote, and reporting one would be a tidiness that the object does not have.

It is also stated for a variance that is a stated function of the covariate. A variance that depends on something unmeasured — a group, a time, a latent stratum — is not covered by any of this, and the correction it needs is not a different weight on a row but a different unit of independence altogether. That is the count that is not the rows, where a sandwich whose filling is summed over clusters behaves like a sample of the size of the cluster count rather than of the row count.

And the whole of this essay is about expectations. Every number above is an average over twenty thousand draws or an integral against a population, and an analyst has one dataset. What the average of a variance estimate is worth to somebody holding a single draw is a question about its spread and about the reference distribution it is read against, and both of those get worse exactly where the sample gets small. That is what the coverage measurement counts, and its answer at twenty rows is that an estimate averaging 0.9576 of the truth delivers an interval covering 88.73% rather than 95% — which is a shortfall no amount of being right on average repairs.