1
$\begingroup$

I am currently working on logistic regression and came across some articles stating "generalized estimating equations with robust standard error" or "with robust sandwich estimators for variance". Can someone please tell me what these mean?

Below is how I run my logistic regression. I am not sure if R already took care of those error calculations. If not, what associated R code I should run in order to have this included in my analysis?

model <- glm (disease ~ prescore, data = mydata, family=binomial)
Improve this question
$\endgroup$
4
  • 2
    $\begingroup$ No, by default, glm or summary(model) will not calculate robust sandwich standard errors. You would have to use another package, for example sandwich and lmtest. Then, you could write coeftest(model, vcov = sandwich) or something like this. $\endgroup$
    – COOLSerdash
    Commented Dec 4, 2021 at 19:49
  • $\begingroup$ Thank you for the input! Do you know how I can get an odds ratio & 95% CI? I know how to get those values without sandwich errors exp(cbind(OR = coef(model), confint(model))) but I can't seem to incorporate the coeftest and vcov into this code. $\endgroup$
    – R Beginner
    Commented Dec 4, 2021 at 20:21
  • 1
    $\begingroup$ Use the coefci() function from lmtest() to get confidence intervals. You can supply the functions used to construct the standard errors (e.g., sandwich::vcovHC) to its vcov. argument. $\endgroup$
    – Noah
    Commented Dec 5, 2021 at 19:09
  • $\begingroup$ @Noah Thank you!! How about odds ratio or p-value? Can you pls put those together (odds ratio, CI, and p-value) for me? I am still a bit confused.. $\endgroup$
    – R Beginner
    Commented Dec 5, 2021 at 19:13

1 Answer 1

Reset to default
1
$\begingroup$

Here is an example of how to do this using the lalonde dataset in cobalt, where treat is a binary variable.

data("lalonde", package = "cobalt")

fit <- glm(treat ~ age + educ, data = lalonde, family = binomial)
est <- lmtest::coeftest(fit, vcov. = sandwich::vcovHC)

#Coefficient estimates and CIs
cbind(est, confint(est))
#>                 Estimate  Std. Error    z value   Pr(>|z|)       2.5 %
#> (Intercept) -0.266058710 0.424816044 -0.6262916 0.53112371 -1.09868286
#> age         -0.023963432 0.008351086 -2.8694987 0.00411123 -0.04033126
#> educ         0.006699987 0.031100402  0.2154309 0.82943139 -0.05425568
#>                   97.5 %
#> (Intercept)  0.566565436
#> age         -0.007595603
#> educ         0.067655654

#Odds ratios and CIs
exp(cbind(est[,1], confint(est)))
#>                           2.5 %    97.5 %
#> (Intercept) 0.7663941 0.3333098 1.7622042
#> age         0.9763214 0.9604712 0.9924332
#> educ        1.0067225 0.9471899 1.0699968

Created on 2021-12-06 by the reprex package (v2.0.1)

Here I used HC3 robust SEs as implemented in vcovHC() in the sandwich package. lmtest::coeftest() provides a nice interface to incorporate the robust standard errors, and you can use confint() on its output to extract the Wald confidence intervals. Alternatively, you could have used lmtest::coefci(fit, vcov. = sandwich::vcovHC) to get the confidence intervals directly.

To get the odds ratios and their confidence intervals, I exponentiated the coefficients stored in the first column of est and their confidence intervals. It's not appropriate to exponentiate the standard errors.

Improve this answer
$\endgroup$
2
  • $\begingroup$ This is very useful! I noticed there was no p-value in the output. If I want to report p-values, can you pls show me the R code for it? $\endgroup$
    – R Beginner
    Commented Dec 6, 2021 at 15:28
  • 1
    $\begingroup$ The p-value is in the column Pr(>|z|). $\endgroup$
    – Noah
    Commented Dec 6, 2021 at 15:50

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.