Independence GEE as a Marginal GLM Score

Why Compare Independence GEE With a GLM?

A generalized estimating equation (GEE) uses a marginal mean model for outcomes that may be grouped into clusters, such as repeated measurements on the same patient. The working covariance describes how the estimating equation weights residuals; it need not assert the true dependence structure. With an independence working correlation, each residual is weighted as if it were separate from the others. This makes the coefficient equation coincide with the ordinary marginal GLM score, provided the mean, link, variance function, and observation weights agree.

That coefficient equivalence does not force the uncertainty estimates to agree. The ordinary GLM model-based covariance uses the assumed observation-level covariance, whereas a cluster-robust estimate uses variation among whole-cluster score contributions. If observations within a cluster are dependent, those two calculations can differ even when they use the same fitted coefficients.

The Independence-Working Equation

Let $Y_{ij}$ be response $j$ in cluster $i$, for $i=1,\ldots,m$ and $j=1,\ldots,n_i$. Suppose the marginal mean is specified by $\mu_{ij}=E[Y_{ij}]$ and $g(\mu_{ij})=x_{ij}^{T}\beta$, where $x_{ij}$ is a fixed covariate vector, $g$ is the link, and $\beta$ is the coefficient vector. An offset can be included in the linear predictor without changing the argument below. Define $D_i=\partial\mu_i/\partial\beta^T$, the matrix of derivatives of the cluster mean vector, and let $A_i$ be diagonal with entries $\phi v(\mu_{ij})$. Here $v$ is the GLM variance function and $\phi$ is a common dispersion parameter. For simplicity, there are no prior weights; with such weights, the corresponding weighted variances must be used in both procedures.

Definition: The independence-working GEE is $$ U_{\mathrm{ind}}(\beta)=\sum_{i=1}^{m}D_i^T A_i^{-1}(Y_i-\mu_i)=0, $$ where $Y_i=(Y_{i1},\ldots,Y_{in_i})^T$ and $\mu_i=(\mu_{i1},\ldots,\mu_{in_i})^T$. Equivalently, its working covariance for cluster $i$ is $V_i=A_i$, which corresponds to setting the working correlation matrix to the identity.

For an individual observation, the contribution to this equation is \[ \frac{\partial\mu_{ij}}{\partial\beta}\, \frac{Y_{ij}-\mu_{ij}}{\phi v(\mu_{ij})}. \] Thus the equation adds the usual GLM residual contribution over all observations. Grouping observations into clusters does not change this sum when the working covariance is diagonal.

Why the Coefficient Equations Match

The equality can be seen directly from the GLM score. Suppose each response has an exponential-dispersion density or mass function with log likelihood, up to terms not depending on its mean, $\{Y_{ij}\theta_{ij}-b(\theta_{ij})\}/\phi$. Then $\mu_{ij}=b'(\theta_{ij})$, and $v(\mu_{ij})=b''(\theta_{ij})$. Differentiating with respect to $\mu_{ij}$ gives \[ \frac{\partial\ell_{ij}}{\partial\mu_{ij}} =\frac{Y_{ij}-\mu_{ij}}{\phi v(\mu_{ij})}. \] The chain rule then gives the coefficient score.

Theorem: For the same marginal mean model, link, variance function, dispersion convention, and observation weights, the independence-working GEE coefficient equation equals the ordinary GLM coefficient score equation. In particular, their coefficient solutions coincide when the equations have the same unique solution. This is an equality of estimating equations; it does not by itself assert that the joint distribution of responses in a cluster is independent.

Proof. For each observation, the exponential-dispersion form gives $\partial\ell_{ij}/\partial\theta_{ij}=(Y_{ij}-b'(\theta_{ij}))/\phi$ and $d\mu_{ij}/d\theta_{ij}=b''(\theta_{ij})=v(\mu_{ij})$. Therefore, by the chain rule, \[ \frac{\partial\ell_{ij}}{\partial\mu_{ij}} =\frac{\partial\ell_{ij}}{\partial\theta_{ij}} \frac{d\theta_{ij}}{d\mu_{ij}} =\frac{Y_{ij}-\mu_{ij}}{\phi v(\mu_{ij})}. \] Also, $\partial\mu_{ij}/\partial\beta$ is the observation’s row in $D_i$. Applying the chain rule from $\mu_{ij}$ to $\beta$ and summing gives \[ \frac{\partial\ell}{\partial\beta} =\sum_{i=1}^{m}\sum_{j=1}^{n_i} \left(\frac{\partial\mu_{ij}}{\partial\beta}\right) \frac{Y_{ij}-\mu_{ij}}{\phi v(\mu_{ij})} =\sum_{i=1}^{m}D_i^T A_i^{-1}(Y_i-\mu_i). \] The last expression is $U_{\mathrm{ind}}(\beta)$. Hence the equations are identical. If the GLM has a common unknown dispersion, its factor $1/\phi$ does not change the coefficient equation’s roots. $\square$

The same algebra defines the usual GLM quasi-score when a full likelihood is not being used. In that case, “same score” means the same estimating equation, not necessarily the score of a correctly specified joint probability model.

Worked Example: Poisson Counts With Exposure

Suppose four observations have Poisson means $\mu_j=t_j\exp(\beta)$, where the exposures are $t=(3,2,1,2)$, and the observed counts are $y=(1,3,2,4)$. The log link with $\log t_j$ as an offset gives $\partial\mu_j/\partial\beta=\mu_j$, while the Poisson variance function is $v(\mu_j)=\mu_j$ and $\phi=1$.

The independence GEE contribution for observation $j$ is \[ \frac{\partial\mu_j}{\partial\beta}\frac{y_j-\mu_j}{\mu_j} =y_j-\mu_j. \] Consequently its equation is \[ 0=\sum_{j=1}^{4}(y_j-\mu_j) =10-\exp(\beta)(3+2+1+2) =10-8\exp(\beta). \] Solving gives $\exp(\hat\beta)=10/8=5/4$ and $\hat\beta=\log(5/4)\approx0.2231$. The ordinary Poisson GLM score is the same equation, so it gives the same estimate. The exposure sum is $8$; using any other total would change the fitted rate.

Same Coefficients, Different Covariance

The coefficient equation alone does not determine how its solution varies across repeated samples. To see the distinction, let $\Sigma_i=\operatorname{Cov}(Y_i)$ denote the actual covariance within cluster $i$. With fixed covariates, the variance of the independence-working score at the true coefficient is \[ \operatorname{Var}\{U_{\mathrm{ind}}(\beta)\} =\sum_{i=1}^{m}D_i^T A_i^{-1}\Sigma_i A_i^{-1}D_i, \] assuming clusters are independent. This follows by writing the score as a sum of cluster contributions $U_i=D_i^T A_i^{-1}(Y_i-\mu_i)$, taking each contribution’s variance, and using independence to make cross-cluster covariances zero.

If the working covariance $A_i$ is also the actual covariance $\Sigma_i$, this score variance reduces to $\sum_i D_i^T A_i^{-1}D_i$. Its inverse is the familiar model-based covariance, under the usual regularity and correctly specified variance assumptions. If within-cluster covariances differ from the working ones, the score variance generally does not reduce to that expression. The cluster-robust covariance estimates the score variance using the observed cluster contributions, with the appropriate sensitivity matrix on either side. This is why the two covariance estimates can differ even though the coefficient equation is unchanged.

Key distinction: Independence working correlation specifies the weights used to solve the coefficient equation. It does not establish that the observations are truly independent. Model-based uncertainty uses the working covariance as the response covariance; cluster-robust uncertainty allows the observed cluster scores to reflect within-cluster dependence.

Worked Example: Logistic Intercept With Three Clusters

Consider three clusters, each with two binary responses: $(0,0)$, $(0,1)$, and $(1,1)$. Fit a marginal logistic model with intercept $\beta$, so $p=\mu=\exp(\beta)/(1+\exp(\beta))$. The independence GEE is \[ \sum_{i=1}^{3}\sum_{j=1}^{2}(y_{ij}-p)=0. \] There are three successes among six responses, so this becomes $3-6p=0$, giving $\hat p=1/2$ and $\hat\beta=\log\{(1/2)/(1/2)\}=0$. This is also the ordinary logistic GLM coefficient equation.

For the model-based calculation, the derivative of the mean with respect to the intercept is $p(1-p)=1/4$, and the Bernoulli variance is also $1/4$. The sensitivity is therefore $H=6(1/4)=3/2$, so the model-based variance estimate is $H^{-1}=2/3$.

At the fitted value, each cluster’s score is its number of successes minus $2p$. The three cluster scores are $-1,0,1$, so their squared values sum to $2$. The uncorrected cluster-robust estimate is $H^{-1}(2)H^{-1}=(2/3)^2(2)=8/9$. It differs from $2/3$, although the fitted coefficient is identical. This small-sample calculation is illustrative; it applies no finite-sample correction and is not a claim that either estimate is precise with only three clusters.

Worked Example: Gaussian Identity-Link Regression

Take four responses with covariates $x=(0,0,1,1)$ and values $y=(1,3,2,4)$. Use the marginal model $E[Y_j]=\beta_0+\beta_1x_j$, with working variance $\sigma^2$ for every observation. The independence GEE is the ordinary linear-model score: \[ \frac{1}{\sigma^2}\sum_{j=1}^4 \begin{pmatrix}1\\x_j\end{pmatrix} (y_j-\beta_0-\beta_1x_j)=0. \] After multiplying by $\sigma^2$, the equations are $\sum_j r_j=0$ and $\sum_jx_jr_j=0$, where $r_j=y_j-\beta_0-\beta_1x_j$.

Expanding gives $4\beta_0+2\beta_1=10$ and $2\beta_0+2\beta_1=6$. Subtracting the second equation from the first yields $2\beta_0=4$, so $\hat\beta_0=2$; substituting gives $\hat\beta_1=1$. The fitted values in observation order are $2,2,3,3$, and the residuals are $-1,1,-1,1$. Their total is zero, and the two residuals with $x=1$ also sum to zero, verifying both equations. These are the ordinary least-squares coefficients as well as the independence-GEE coefficients.

Common Mistakes and Practical Interpretation

Assuming the independence working correlation proves independence. It is a choice of working weights. Repeated outcomes can remain dependent, which is precisely why the cluster-robust covariance may be useful.
Expecting robust and model-based standard errors to match. They agree only under conditions that make the working covariance suitable for the actual response covariance (or by coincidence in a particular sample). Identical coefficient estimates do not imply identical covariance estimates.
Comparing unlike GLMs. The score equations match only when mean model, link, variance function, dispersion convention, and weights are aligned. Changing an offset or variance assumption can change the equation, not merely the standard error.
Calling every independence-GEE fit an ordinary likelihood fit. The coefficient equation may match a GLM score or quasi-score even when a joint likelihood for the correlated cluster outcomes has not been specified.

In practice, the independence working structure is a useful reference point: it reveals which coefficient equation is being solved and makes clear that clustering enters the robust uncertainty calculation through cluster-level contributions. The distinction matters whenever an analysis reports both fitted marginal effects and their standard errors.

Check Your Understanding

Use the definitions and calculations above to check the distinction between the coefficient equation and its covariance estimate.

What assumptions must be aligned for the independence-working GEE equation to equal the ordinary GLM coefficient score?
In the Poisson example, why does the exposure total enter the equation, and what is its value?
For the logistic example, state the three cluster scores and explain why their squared sum is not the model-based sensitivity.
Does choosing an independence working correlation imply that repeated responses are independent? Explain.
Give one reason the cluster-robust and model-based covariance estimates can differ even when the coefficient estimates match.