The hardware and bandwidth for this mirror is donated by METANET, the Webhosting and Full Service-Cloud Provider.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]metanet.ch.

A comparison of geex and sandwich for robust covariance estimation

Bradley Saul

2022-07-24

This examples uses the vaccinesim dataset from the inferference package to compare the estimated covariance matrix obtained from geex and sandwich. An example \(\psi\) function written in R.

This function computes the score functions for a GLM.

eefun <- function(data, model){
  X <- model.matrix(model, data = data)
  Y <- model.response(model.frame(model, data = data))
  function(theta){
    lp  <- X %*% theta
    rho <- plogis(lp)

    score_eqns <- apply(X, 2, function(x) sum((Y - rho) * x))
    score_eqns
  }
}

Compare sandwich variance estimators to sandwich treating individuals as units:

library(geex)
library(inferference)
mglm    <- glm(A ~ X1, data = vaccinesim, family = binomial)
estimates <- m_estimate(
  estFUN = eefun,
  data = vaccinesim,
  root_control = setup_root_control(start = c(-.35, 0)),
  outer_args = list(model = mglm))

# Compare point estimates
coef(estimates) # from GEEX

## [1] -0.36869683 -0.02037916

coef(mglm) # from the GLM function

## (Intercept)          X1 
## -0.36869683 -0.02037916

# Compare variance estimates
vcov(estimates)

##               [,1]          [,2]
## [1,]  0.0028345579 -0.0007476536
## [2,] -0.0007476536  0.0003870030

sandwich::sandwich(mglm)

##               (Intercept)            X1
## (Intercept)  0.0028345579 -0.0007476536
## X1          -0.0007476536  0.0003870030

Pretty darn good! Note that the geex method is much slower than sandwich (especially using method = 'Richardson' for numDeriv), but this is because sandwich uses the closed form of the score equations, while geex compute them numerically. However, geex’s real utility comes when you have more complicated estimating equations. Also, the analyst has the ability to code faster \(\psi\) functions by optimizing their code or using Rccp, for example.

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.