Variance Misspecification and Robust Covariance

Why Variance Misspecification Matters

A quasi-score uses a specified mean and a variance function to weight observations. The mean model determines the target of estimation; the variance function determines how observations contribute to the estimating equation and to the model-based uncertainty calculation. If the mean is correct but the variance function is wrong, the estimating equation can still have expectation zero at the true parameter. Its usual model-based standard errors, however, need not describe the actual sampling variability.

The sandwich covariance addresses this discrepancy by combining the sensitivity of the estimating equation to the parameter with the actual variability of the estimating equation. It is robust to a wrong variance function under stated conditions, not to every possible model failure.

Working Variance, Score, and Covariance

Definition: Suppose $Y_1,\ldots,Y_n$ are independent, with $E[Y_i]=\mu_i(\beta_0)$, where $\beta_0$ is the true parameter and $\mu_i(\beta)$ is a differentiable mean model. Write $d_i(\beta)=\partial\mu_i(\beta)/\partial\beta$, a column vector. A working variance specification is $\phi V_{w,i}(\beta)$, with $V_{w,i}(\beta)>0$ and dispersion $\phi>0$. The corresponding quasi-score contribution and total estimating function are $$ \psi_i(\beta)=\frac{d_i(\beta)}{\phi V_{w,i}(\beta)}\{Y_i-\mu_i(\beta)\}, \qquad U_n(\beta)=\sum_{i=1}^n\psi_i(\beta). $$ A quasi-estimator $\hat\beta$ solves $U_n(\hat\beta)=0$.

The working variance need not equal the true variance $\operatorname{Var}(Y_i)$. If the mean is correct, then $E[\psi_i(\beta_0)]=0$: the weight is a nonrandom function of the design and parameter, and $E[Y_i-\mu_i(\beta_0)]=0$. This is the key reason that a wrong variance function need not cause inconsistency. It can nevertheless change the weights and the estimator's sampling covariance.

For independent observations, define the sensitivity, often called the bread, and the variability, often called the meat, by their limiting per-observation matrices: $$ A=\lim_{n\to\infty}\frac{1}{n}\sum_{i=1}^n E\left[-\frac{\partial\psi_i(\beta_0)}{\partial\beta^T}\right], \qquad B=\lim_{n\to\infty}\frac{1}{n}\sum_{i=1}^n\operatorname{Var}\{\psi_i(\beta_0)\}. $$ Assume these limits exist and are finite, and that $A$ is nonsingular. The transpose in the covariance formula below is needed in general; for the usual quasi-score with a correctly specified mean, $A$ is symmetric.

Theorem: Suppose the observations are independent, the mean model is correct at an interior parameter $\beta_0$, and $\hat\beta\xrightarrow{p}\beta_0$. Suppose $\psi_i(\beta)$ is differentiable near $\beta_0$, the normalized derivative satisfies a uniform law of large numbers there, $A$ exists and is nonsingular, and the score contributions satisfy a central limit theorem with limiting covariance $B$. Assume also the regularity needed for Taylor expansion of $U_n$ between $\hat\beta$ and $\beta_0$. Then $$ \sqrt{n}(\hat\beta-\beta_0)\xrightarrow{d} N\left(0,A^{-1}BA^{-T}\right). $$

Proof. Since $E[\psi_i(\beta_0)]=0$, the central limit theorem gives $$ \frac{1}{\sqrt n}U_n(\beta_0)\xrightarrow{d}N(0,B). $$ Taylor expansion of the estimating equation around $\beta_0$ gives $$ 0=\frac{1}{\sqrt n}U_n(\beta_0) +\left\{\frac{1}{n}\sum_{i=1}^n \frac{\partial\psi_i(\tilde\beta)}{\partial\beta^T}\right\} \sqrt n(\hat\beta-\beta_0), $$ where $\tilde\beta$ lies between $\hat\beta$ and $\beta_0$ (or the corresponding integral form of the derivative is used). Consistency and the assumed uniform law of large numbers imply that the matrix in braces converges in probability to $-A$. It is therefore invertible with probability approaching one. Solving the expansion yields $$ \sqrt n(\hat\beta-\beta_0) = \left\{-\frac{1}{n}\sum_{i=1}^n \frac{\partial\psi_i(\tilde\beta)}{\partial\beta^T}\right\}^{-1} \frac{1}{\sqrt n}U_n(\beta_0). $$ The first factor converges in probability to $A^{-1}$. Slutsky's theorem then gives the stated normal limit, whose covariance is $A^{-1}BA^{-T}$. $\square$

For this quasi-score, the expected sensitivity has a particularly useful form. Differentiating $\psi_i$ produces a term proportional to the residual $Y_i-\mu_i(\beta_0)$; its expectation is zero. The remaining term gives $$ E\left[-\frac{\partial\psi_i(\beta_0)}{\partial\beta^T}\right] = \frac{d_i(\beta_0)d_i(\beta_0)^T}{\phi V_{w,i}(\beta_0)}. $$ Thus $A$ is determined by the working weights and the mean-model derivatives. In contrast, $B$ depends on the actual variance: $$ \operatorname{Var}\{\psi_i(\beta_0)\} = \frac{d_i(\beta_0)d_i(\beta_0)^T} {\{\phi V_{w,i}(\beta_0)\}^2}\operatorname{Var}(Y_i). $$ The model-based covariance uses the working variance to stand in for the actual variance. The sandwich keeps the sensitivity but estimates the score variability separately.

What the Model-Based Standard Error Assumes

When the working variance is correct, $\operatorname{Var}(Y_i)=\phi V_{w,i}(\beta_0)$. The two displayed matrices then agree: each score contribution has variance equal to its expected sensitivity, so $B=A$. The asymptotic covariance becomes $A^{-1}/n$. If the variance is misspecified, that equality need not hold, and the model-based covariance $A^{-1}/n$ can be too large or too small.

A useful special case is a common variance multiplier. Suppose the true variance is $c\phi V_{w,i}(\beta_0)$ for every $i$, with $c>0$. The expression for the score variance shows that $B=cA$. Therefore the actual asymptotic covariance is $cA^{-1}/n$, while the model-based covariance remains $A^{-1}/n$. This identity depends on the same multiplier applying throughout; general variance misspecification does not reduce to a single scale adjustment.

Worked Example: A Poisson Mean With Twice the Working Variance

Let $Y_1,\ldots,Y_n$ be independent with $E[Y_i]=\mu$ and $\operatorname{Var}(Y_i)=2\mu$. Fit an intercept-only quasi-Poisson mean using the working variance $V_w(\mu)=\mu$ and $\phi=1$. Find the estimate and compare its model-based and actual standard errors when $\mu=4$ and $n=100$.

Estimate. Here $d_i(\mu)=1$, so the score is $$ U_n(\mu)=\sum_{i=1}^n\frac{Y_i-\mu}{\mu} =\frac{1}{\mu}\left(\sum_{i=1}^nY_i-n\mu\right). $$ The root is $\hat\mu=\bar Y$, which has expectation $\mu$.

Model-based standard error. The working variance says $\operatorname{Var}(Y_i)=\mu$, so the model-based variance of the sample mean is $\mu/n=4/100=0.04$. Its standard error is $\sqrt{0.04}=0.20$.

Actual and robust standard error. Under the stated distribution, $\operatorname{Var}(\bar Y)=2\mu/n=8/100=0.08$, giving standard error $\sqrt{0.08}\approx0.283$. The sandwich obtains the same result: for $\psi_i=(Y_i-\mu)/\mu$, the sensitivity per observation is $A=1/\mu$, while the meat per observation is $B=2/\mu$. Hence $A^{-1}BA^{-1}=\mu^2(2/\mu)=2\mu=8$, and the covariance of $\hat\mu$ is $8/n=0.08$. The model-based standard error is too small because the actual variance is twice the working variance.

Worked Regression Comparisons

Worked Example: Heteroskedastic Linear Regression

Consider independent errors in the model $Y_i=\beta_0+\beta_1x_i+\varepsilon_i$, with $E[\varepsilon_i]=0$. Use a working variance of 1 for every observation. For the four fixed design values $x=(-1,0,1,2)$, suppose the actual error variances are $(1,1,4,4)$. Calculate the model-based covariance and the actual sandwich covariance of the least-squares estimator.

Let $X$ have rows $(1,x_i)$, and let $\Omega$ be the diagonal matrix of the four error variances. The estimating function is $X^T(Y-X\beta)$. Its sensitivity and variability are $$ A=X^TX= \begin{pmatrix}4&2\\2&6\end{pmatrix}, \qquad B=X^T\Omega X= \begin{pmatrix}10&11\\11&21\end{pmatrix}. $$ For example, the off-diagonal entry of $B$ is $\sum_i x_i\operatorname{Var}(\varepsilon_i)=(-1)(1)+0(1)+1(4)+2(4)=11$, and its lower-right entry is $\sum_i x_i^2\operatorname{Var}(\varepsilon_i)=1+0+4+16=21$.

Since $\det(A)=4(6)-2(2)=20$, $$ A^{-1}=\frac{1}{20} \begin{pmatrix}6&-2\\-2&4\end{pmatrix}. $$ The model-based covariance is $A^{-1}$. For the sandwich, first multiply $$ A^{-1}B=\frac{1}{20} \begin{pmatrix}38&24\\24&62\end{pmatrix}. $$ Multiplying on the right by $A^{-1}$ gives $$ A^{-1}BA^{-1}= \frac{1}{400} \begin{pmatrix}180&20\\20&200\end{pmatrix} = \begin{pmatrix}0.45&0.05\\0.05&0.50\end{pmatrix}. $$ The model-based variances of the intercept and slope are $0.30$ and $0.20$; the actual variances are $0.45$ and $0.50$. Thus the corresponding standard errors are approximately $0.548,0.447$ under the working model and $0.671,0.707$ under the actual covariance. Unequal error variances have changed both marginal variances and their covariance.

Worked Example: Checking the Sandwich Matrix Multiplication

Suppose an estimating equation has limiting sensitivity and meat $$ A=\begin{pmatrix}3&3\\3&9\end{pmatrix}, \qquad B=4qq^T,\qquad q=\begin{pmatrix}1\\1\end{pmatrix}. $$ Compute the sandwich covariance and compare it with the model-based covariance $A^{-1}$.

The determinant of $A$ is $27-9=18$, so $$ A^{-1}=\frac{1}{18}\begin{pmatrix}9&-3\\-3&3\end{pmatrix} =\begin{pmatrix}1/2&-1/6\\-1/6&1/6\end{pmatrix}. $$ Direct multiplication gives $A^{-1}q=(1/3,0)^T$. Since $B=4qq^T$, the sandwich can be calculated without skipping a matrix product: $$ A^{-1}BA^{-1} =4(A^{-1}q)(A^{-1}q)^T =4\begin{pmatrix}1/3\\0\end{pmatrix} \begin{pmatrix}1/3&0\end{pmatrix} = \begin{pmatrix}4/9&0\\0&0\end{pmatrix}. $$ By contrast, the model-based covariance is $$ A^{-1}=\begin{pmatrix}1/2&-1/6\\-1/6&1/6\end{pmatrix}. $$ The meat here has rank one, so the sandwich has zero variance in one parameter direction. This is a degenerate illustration of why the two matrices cannot be assumed equal; in regular applications, one checks the conditions that support the desired inference as well as the covariance calculation.

Estimating Robust Covariance and Knowing Its Limits

The population matrices $A$ and $B$ are unknown. A common empirical sandwich estimate replaces the sensitivity by the observed derivative and the meat by the sum of outer products of fitted score contributions: $$ \widehat A=-\frac{1}{n}\sum_{i=1}^n \frac{\partial\psi_i(\hat\beta)}{\partial\beta^T}, \qquad \widehat B=\frac{1}{n}\sum_{i=1}^n \psi_i(\hat\beta)\psi_i(\hat\beta)^T, \qquad \widehat{\operatorname{Cov}}(\hat\beta) =\frac{1}{n}\widehat A^{-1}\widehat B\widehat A^{-T}. $$ Under the theorem's regularity conditions, the mean model's correctness, and laws of large numbers for these empirical matrices, $\widehat A$ and $\widehat B$ converge to $A$ and $B$. Continuity of matrix inversion at nonsingular $A$ then gives consistency for the asymptotic covariance. The formula uses independence: when observations have dependence, the meat must account for covariance between contributions rather than treating them as independent.

A sandwich estimate does not make an incorrect mean model correct. If $E[Y_i]\ne\mu_i(\beta_0)$, the score may not have expectation zero at the intended parameter, so the consistency and centering argument above can fail. Nor does “robust” mean exact in a small sample: the sandwich estimates a large-sample covariance and can be unstable with limited data or influential observations.

Common Mistake

Common mistake: “The estimating equation is unbiased, so the model-based standard error must be right.” Unbiasedness of $U_n(\beta_0)$ concerns its expected value, not its variance. The model-based calculation also assumes the working variance correctly describes the variability entering the score. If that assumption fails, use the actual score variability through the sandwich, provided the mean and regularity conditions hold.

When reviewing a quasi-likelihood analysis, keep three questions separate: Is the mean model plausible? Does the working variance describe the data well enough for the intended weighting and model-based uncertainty? If not, do the independence and regularity assumptions justify a sandwich covariance? These questions distinguish a robust standard error from a general cure for model misspecification.

Check Your Understanding

Use the definitions and calculations above to test your understanding.

Why can a quasi-score have expectation zero at $\beta_0$ even when its working variance is wrong?
In the sandwich formula, what does the sensitivity matrix measure, and what does the meat matrix measure?
If the true variance is $c$ times the working variance for every observation, how does the asymptotic covariance compare with the model-based covariance?
For the heteroskedastic regression example, verify the intercept-slope entry of $B=X^T\Omega X$ from the four design points and variances.
Give one reason a sandwich covariance might fail to provide valid inference even if it has been calculated correctly.