Efficiency Gains Under a Correct Working Correlation

Why a Correct Working Covariance Matters

In generalized estimating equations, the working covariance determines how residuals within a cluster are weighted. When it correctly represents the covariance of the response, those weights use the observed dependence rather than treating repeated measurements as unrelated. The resulting regression estimator has a particularly simple first-order covariance: the inverse of its information matrix.

This covariance is not automatically smaller for every coefficient than a variance calculated under independence. The comparison depends on the mean model and the contrast being estimated. A useful way to understand the gain is to derive the covariance first, then examine two-visit models where the relevant contrast is explicit.

Covariance When the Working Covariance Is Correct

Suppose there are $n$ independent clusters. Cluster $i$ has response vector $Y_i$, mean vector $\mu_i(\beta)$, and true covariance matrix $V_i$, where the regression parameter $\beta$ has dimension $p$. Write $D_i=\partial\mu_i(\beta)/\partial\beta^T$, a matrix with one row per response and one column per regression parameter. Assume the working covariance in the estimating equation equals $V_i$.

Definition: The correctly weighted GEE and its information matrix are $$ U_n(\beta)=\sum_{i=1}^n D_i^T V_i^{-1}\{Y_i-\mu_i(\beta)\}, \qquad H_n(\beta)=\sum_{i=1}^n D_i^T V_i^{-1}D_i. $$ The GEE estimator $\hat\beta$ is a solution of $U_n(\hat\beta)=0$.

For the covariance calculation, take $V_i$ as known at the true parameter, or as the true covariance evaluated at that parameter. Assume the clusters are independent; the mean is differentiable in a neighbourhood of the true value $\beta_0$; each $V_i$ is positive definite; and the usual regularity conditions for a central limit theorem and a Taylor expansion hold. More specifically, assume $H_n(\beta_0)/n$ converges to a finite positive-definite matrix $H$, the normalized score satisfies a multivariate central limit theorem, and the Taylor remainder in the estimating equation is negligible at the $n^{-1/2}$ scale. Also assume a consistent solution $\hat\beta$ exists in a neighbourhood of $\beta_0$.

Theorem: Under these conditions, if the working covariance equals the true covariance, then $$ \sqrt{n}(\hat\beta-\beta_0)\xrightarrow{d}N(0,H^{-1}), \qquad H=\lim_{n\to\infty}\frac{H_n(\beta_0)}{n}. $$ Thus the first-order covariance of $\hat\beta$ is $H_n(\beta_0)^{-1}$, or approximately $H^{-1}/n$ when expressed using the limiting information.

To see why, first note that $E[Y_i-\mu_i(\beta_0)]=0$, so $E[U_n(\beta_0)]=0$. Independence across clusters and the correct covariance give $$ \operatorname{Var}\{U_n(\beta_0)\} =\sum_{i=1}^n D_i^T V_i^{-1}\operatorname{Var}(Y_i)V_i^{-1}D_i =\sum_{i=1}^n D_i^T V_i^{-1}D_i =H_n(\beta_0). $$ The derivative of the mean score is $-H_n(\beta_0)$ to first order. Terms involving second derivatives of the mean multiply residuals and are negligible under the stated regularity conditions.

A first-order Taylor expansion around $\beta_0$ therefore yields $$ 0=U_n(\hat\beta) =U_n(\beta_0)-H_n(\beta_0)(\hat\beta-\beta_0)+o_p(\sqrt{n}). $$ Equivalently, after using $H_n(\beta_0)/n\to H$, $$ \sqrt{n}(\hat\beta-\beta_0) = \left\{\frac{H_n(\beta_0)}{n}\right\}^{-1} \frac{U_n(\beta_0)}{\sqrt{n}}+o_p(1). $$ The central limit theorem gives $U_n(\beta_0)/\sqrt{n}\xrightarrow{d}N(0,H)$. Applying the continuous mapping theorem to the preceding expansion gives a limiting covariance $H^{-1}HH^{-1}=H^{-1}$, as claimed. $\square$

In applications, covariance parameters are usually estimated rather than known. Under standard conditions—including consistent estimation of those parameters, smooth dependence of the estimating equation on them, and a nonsingular limiting information matrix—their estimation does not change this first-order covariance when the working covariance is correct. The matrix evaluated at fitted values is then a plug-in estimate of the model-based covariance.

Why Inverse-Covariance Weighting Is Efficient

The same calculation has a useful optimality interpretation. Consider estimating a locally linear change in $\beta$ from the cluster residuals. Among linear unbiased estimating functions, the one that weights by the inverse of the true covariance has the smallest covariance. This gives a precise sense in which the correct working covariance can improve efficiency.

Theorem: Let $D_i$ be the mean derivative at the true parameter, let the independent cluster residuals have mean zero and covariance $V_i$, and suppose $H_n=\sum_i D_i^T V_i^{-1}D_i$ is nonsingular. For any matrices $C_i$ satisfying $\sum_i C_i^T D_i=I_p$, the linear estimator increment $\sum_i C_i^T\{Y_i-\mu_i(\beta_0)\}$ has covariance at least $H_n^{-1}$, in the positive-semidefinite ordering. Equality is attained by $C_i=V_i^{-1}D_iH_n^{-1}$.

To prove the result, set $C_i^*=V_i^{-1}D_iH_n^{-1}$ and $\Delta_i=C_i-C_i^*$. Both sets of matrices satisfy the unbiasedness constraint, so $\sum_i\Delta_i^TD_i=0$. Independence gives the covariance of the estimator using $C_i$ as $\sum_i C_i^TV_iC_i$. Substituting $C_i=C_i^*+\Delta_i$, its difference from the covariance using $C_i^*$ is $$ \sum_i \Delta_i^T V_i\Delta_i +\sum_i \Delta_i^T V_iC_i^* +\sum_i (C_i^*)^T V_i\Delta_i. $$ The first term is positive semidefinite. For the cross terms, $V_iC_i^*=D_iH_n^{-1}$, so $$ \sum_i\Delta_i^T V_iC_i^* =\left(\sum_i\Delta_i^TD_i\right)H_n^{-1}=0. $$ Its transpose is also zero. Finally, $$ \sum_i (C_i^*)^TV_iC_i^* =H_n^{-1}\left(\sum_iD_i^TV_i^{-1}D_i\right)H_n^{-1} =H_n^{-1}. $$ Thus every other linear unbiased choice has covariance $H_n^{-1}+\sum_i\Delta_i^TV_i\Delta_i$, which is no smaller than $H_n^{-1}$. $\square$

This is an optimality result for the stated class of linear unbiased estimating functions. It does not say that every coefficient in every model must have a smaller variance than under independence, nor does it justify comparing a correct covariance estimate with an independence standard error that is calculated as though the observations really were independent.

Worked Two-Visit Comparisons

Worked Example: Using a Correlated Visit to Improve Precision

For each of $n$ independent subjects, suppose the two responses follow $$ Y_{i1}=\beta+\varepsilon_{i1},\qquad Y_{i2}=\varepsilon_{i2}, \qquad \operatorname{Var} \begin{pmatrix}\varepsilon_{i1}\\ \varepsilon_{i2}\end{pmatrix} =\sigma^2 \begin{pmatrix}1&\rho\\ \rho&1\end{pmatrix}, $$ where $\sigma^2$ and $\rho$ are known and $|\rho|<1$. The second visit has known mean zero, but its error is correlated with the error at the first visit. Estimate the single regression coefficient $\beta$.

Here $D_i=(1,0)^T$ and $V_i=\sigma^2R$, where $R$ is the displayed correlation matrix. Direct inversion gives $$ R^{-1}=\frac{1}{1-\rho^2} \begin{pmatrix}1&-\rho\\-\rho&1\end{pmatrix}, \qquad D_i^TV_i^{-1}D_i=\frac{1}{\sigma^2(1-\rho^2)}. $$ The correctly weighted estimating equation is therefore $$ 0=\sum_{i=1}^n D_i^TV_i^{-1} \begin{pmatrix}Y_{i1}-\beta\\Y_{i2}\end{pmatrix} =\frac{1}{\sigma^2(1-\rho^2)} \sum_{i=1}^n\{Y_{i1}-\beta-\rho Y_{i2}\}. $$ Solving gives $\hat\beta=n^{-1}\sum_i(Y_{i1}-\rho Y_{i2})$. Its variance is $$ \operatorname{Var}(\hat\beta) =\frac{1}{n}\operatorname{Var}(\varepsilon_{i1}-\rho\varepsilon_{i2}) =\frac{\sigma^2}{n}(1+\rho^2-2\rho^2) =\frac{\sigma^2(1-\rho^2)}{n}. $$ For $\rho=0.6$, this is $0.64\sigma^2/n$, compared with $\sigma^2/n$ for the estimator that uses only the first visit. The second visit improves precision by helping remove correlated error. Treating the visits as independent would omit that opportunity.

Worked Example: Comparing the Correct Covariance With Independence

Use the same two-visit model, but now compare the correct estimator with the independence-working estimator. Under independence weights, the estimating equation uses $D_i^T(\sigma^2I)^{-1}$, so it gives $\hat\beta_{\mathrm{ind}}=n^{-1}\sum_iY_{i1}$. Its actual variance under the stated correlated model is $\sigma^2/n$. The correct estimator from the previous example has variance $\sigma^2(1-\rho^2)/n$.

At $\rho=0.6$, the actual variance ratio is $$ \frac{\operatorname{Var}(\hat\beta)} {\operatorname{Var}(\hat\beta_{\mathrm{ind}})} =\frac{0.64\sigma^2/n}{\sigma^2/n}=0.64. $$ Thus the correct estimator has 36% lower variance in this example. Its standard error is smaller by a factor of $\sqrt{0.64}=0.8$. This comparison concerns the estimators’ actual variances under the specified model. The independence model-based variance for the independence estimator happens to be $\sigma^2/n$ here, but in general a model-based variance based on an incorrect covariance need not equal an estimator’s actual variance.

Worked Example: A Common Mean Can Have Less Information Under Positive Correlation

Now suppose each subject has the same mean $\mu$ at both visits, with the same covariance matrix $\sigma^2R$ and known $\rho$. Then $D_i=(1,1)^T$. The correct information per subject is $$ D_i^TV_i^{-1}D_i =\frac{1}{\sigma^2(1-\rho^2)} \begin{pmatrix}1&1\end{pmatrix} \begin{pmatrix}1&-\rho\\-\rho&1\end{pmatrix} \begin{pmatrix}1\\1\end{pmatrix} =\frac{2}{\sigma^2(1+\rho)}. $$ Consequently, the correct covariance of the estimated common mean is $\sigma^2(1+\rho)/(2n)$. For positive correlation, this is larger than $\sigma^2/(2n)$, the variance obtained by incorrectly treating all $2n$ responses as independent. Indeed, the estimator is the grand mean, whose actual variance is $$ \operatorname{Var}\left\{\frac{1}{2n}\sum_{i=1}^n(Y_{i1}+Y_{i2})\right\} =\frac{1}{4n}\,2\sigma^2(1+\rho) =\frac{\sigma^2(1+\rho)}{2n}. $$ The increase is not a failure of correct weighting: positive correlation means the two visits contain less independent information about a shared mean. The independence calculation understates that uncertainty.

Interpretation and a Common Pitfall

The examples distinguish two claims that are easy to confuse. First, when the working covariance is correct, inverse-covariance weighting gives the model-based covariance $H_n^{-1}$ and is optimal among the linear unbiased estimating functions described above. Second, whether this is numerically smaller than an independence-based quantity depends on what the mean model asks the data to estimate. Correlation can be useful when one visit helps predict and remove noise relevant to a contrast, as in the auxiliary-visit example. For a shared mean with positive correlation, repeated measurements provide less information than independent measurements would.

Common mistake: “The measurements are positively correlated, so using the correct working correlation must always reduce the variance.” This reasoning treats correlation as if it affected every mean contrast in the same way. The common-mean calculation shows the problem: positive correlation increases the actual variance relative to the independence model-based value. The correction is to identify the mean derivative $D_i$, calculate $D_i^TV_i^{-1}D_i$, and compare covariances for the specific estimator and model rather than infer a gain from the sign of the correlation alone.

Key takeaway: With a correct working covariance, the GEE score variance equals its information matrix, so the first-order regression covariance is the inverse information. Efficiency is a property of the weighting and the target mean contrast together—not of the correlation sign in isolation.

Check Your Understanding

Use the definitions and two-visit calculations above to answer the following.

Why does the variance of the correctly weighted score equal $H_n(\beta_0)$ when the working covariance is correct?
In the auxiliary-visit model, what is the coefficient multiplying $Y_{i2}$ in the correctly weighted estimator, and why does it depend on $\rho$?
For $\rho=0.6$, calculate the variance ratio of the correct estimator to the estimator using only the first visit.
For a common mean and positive $\rho$, why is the independence model-based variance too small?
What conditions on the matrices $C_i$ express that a linear estimating function is unbiased for a local change in $\beta$?