So I am trying to get an understanding of causal inference and how it differs from the usual contrasts. I regularly use the emmeans package in R, and I am wondering when the function emmeans() mentions it has averaged over the covariates is this essentially performing G-computation? At least for regular OLS or identity link GLM models it seems like that

One possible difference I see is that G-computation takes place on the response scale, so I wonder if you use the transform argument in emmeans when working with GLMs with a non-identity link, then would it be performing the G-computation? It seems like it does the delta method to convert the SEs to the response scale from the link scale.

And then if I modeled the treatment probability, inverted it, and used it as a covariate in a model, and then used emmeans -- would this be doubly robust estimation?

As my reference I am using the book here https://www.hsph.harvard.edu/miguel-hernan/causal-inference-book/ pages 163-167.

Edit 09/07/22: Anyone wondering about the connection --it is written about here https://vincentarelbundock.github.io/marginaleffects/articles/gformula.html. The marginaleffects package can do this

  • $\begingroup$ I don't really understand what emmeans does, but I was able to replicate g-computation results for linear models (you need to specify type = "prop" in the call to contrasts(), though). For logistic regression models, I could not get the results to agree. I would recommend not using emmeans for g-computation. The stdReg package does a good job for most purposes and provides correct standard errors. $\endgroup$
    – Noah
    Apr 20, 2021 at 5:12
  • 1
    $\begingroup$ As emmeans developer, I do understand what it does but am not an expert on causal inference so not so sure of that. It does seem that you might want to look at the cov.reduce argument (see help for ref_grid), which does provide some flexibility in setting covariates at their predicted values from another model. See also the "framing" example in the vignettes. I have tried to document things pretty clearly and the vignettes n general should help. $\endgroup$
    – Russ Lenth
    Apr 20, 2021 at 13:36
  • 2
    $\begingroup$ There is nothing inherently causal about the output of a linear, generalized linear, mixed or any other kind of model. The thing that makes an estimate causal comes from assumptions about the setting that are external to the statistical model. This means that there is no way that R (or anything else) could do causal inference without you knowing about it. $\endgroup$
    – dimitriy
    Apr 22, 2021 at 22:02

1 Answer 1


Re-reading your question, my understanding is that you are asking if emmeans() does G-computation as part of what it ordinarily does. And based on my very limited understanding of causal models and G-computation, I would say the answer is NO. That is simply because we don't treat covariates in any special way. For a numerical covariate, the default action is to compute its mean and use that as a reference value for all subsequent estimates, regardless of whether it is regarded as a mediator or not. We just treat it as a direct effect.

There may be some options in emmeans() that do allow the user to treat covariates in a different way. For example, we can fit a model y ~ treat + M where treat is a treatment and M is a mediator. Then suppose we subsequently do

emmeans(model, "treat", cov.reduce = M ~ treat)

This instructs emmeans to not use the average value of M, but rather to use lm() to fit the model M ~ treat (with the same dataset) and use its predictions for the value of M. In that way, the reference value of M is different for each treatment level. This is equivalent to creating a covariate C that is equal to the residuals of the M ~ treat model, fitting the model y ~ treat + C, and using emmeans() in the ordinary way by using C's mean (which is zero) as the reference value. Perhaps this is similar to what G-computation does -- I am not sure, but perhaps someone else can shed some light on this. But at least it does something special with covariates thought to be mediators, and that seems more akin to what is needed in causal inference.


A comment to this answer suggests doing something like emmeans(model, "treat", cov.reduce = FALSE, weights = "prop") but that this is very inefficient as it creates a huge reference grid. I believe that the following may do the same thing:

emmeans(model, "treat", submodel = ~ treat)

The above puts a linear constraint on the estimates whereby all the effects other than those of treat are replaced by predictions of those effects from the given submodel. See vignette("xplanations", "emmeans") for the gory details. But in words, what happens is that we are trying to obtain the predictions we would have obtained from the submodel, while still accounting for the reduction in error variance achieved by including the covariate in the model. I think this in fact does relate to some causal-inference methods, but I lack the depth of knowledge in that are to be sure.

In the case of mixed models and generalized linear models, the submodel constraint will not be quite the same as would be obtained by fitting the submodel with the same method. To accomplish this (or at least get closer), one can use a new feature in version 1.6.0 of emmeans to bring in covariate predictions from an external model. Suppose model was fitted using something like model <- lmer(y ~ M + treat + (1|SUBJ), ...)

Mmod <- lmer(M ~ treat + (1|SUBJ, ...)
Mpred <- function(grid)
             list(M = predict(Mmod, newdata = grid, re.form = ~ 0))

emmeans(model, "treat", cov.reduce = extern ~ Mpred)

This is like the original cov.reduce = M ~ treat, except it uses Mmod instead of lm(M ~ treat) to do the predictions of M.

  • $\begingroup$ Thanks, based on the causal inf reference above, page 164, it seems like the G computation (also called standardization) method in the binary treatment case essentially makes 2 additional copies of the dataset but deletes the Ys and changes the treatment to all 0 and all 1 respectively. Then it fits the model on the original dataset and gets predictions on the 2 copies, and the predictions are then averaged and then a difference is computed. So it seems like G-comp is giving E(Y|T=1,X) - E(Y|T=0,X) which looks very similar to a contrast to me, with covariates X averaged out. $\endgroup$
    – Vattaka
    Apr 23, 2021 at 0:52
  • $\begingroup$ So I found out that using the cov.keep argument in ref_grid and keeping the continuous variables at their value, then using weights=proportional in the call to emmeans(grid,pairwise~T,weights=proportional) seems to give the effect estimate as G-comp. But it can be computationally intensive because looks like the ref grid keeps unique values for the continuous variables. To do G-comp, what I need is to average the expected value after, and not set the covariate to its average before (so use the empirical distribution of the covariates basically). Wonder if there is a faster way? $\endgroup$
    – Vattaka
    May 15, 2021 at 16:52
  • $\begingroup$ Did you try your proportional weighting withcov.reduce = M ~ treat? That does not reduce M to its mean, but instead to different predicted values for each treat. $\endgroup$
    – Russ Lenth
    May 15, 2021 at 21:18
  • 1
    $\begingroup$ @Vattaka Please see my addendum. I had to make a couple of corrections to the external cov.reduce part as there were errors $\endgroup$
    – Russ Lenth
    May 18, 2021 at 19:13
  • 1
    $\begingroup$ @Vattaka FYI I have been looking at this recently and have added a counterfactuals argument to emmeans (actually ref_grid) whereby counterfactual predictions are computed for the factors named there. I think that is a step forward in the kind of analysis you are asking about. It's only on GitHub now but will be included in the next CRAN update. $\endgroup$
    – Russ Lenth
    Oct 25, 2022 at 4:38

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.