I would appreciate any help regarding emmeans package. I am fitting dummy-variable regression model (ANCOVA) with follow-up post hoc test in emmeans.

My data includes the following variables produced in experimental setting. Scale is dependent (outcome) variable and Condition, BMI, Sex, Age are independent (predictor) variables. I am truly interested in the effect of Condition on Scale (dependent variable / outcome). Condition is factor with eight levels (e.g. Condition 1, 2, ... 8.) and the effect of interest that was experimentally controlled, BMI, Sex, Age are considered confounders in the model and I am including them but not interpreting.

Given the experimental nature, there are specific levels in variable Condition that I am interested in. Therefore, the post-hoc analysis is a planned contrast comparison with Holm adjustment.

The specified model looks like this:

model1 <- lm(Scale ~ Condition + BMI + Sex + Age, data)

Which results in:

             Estimate Std. Error t value Pr(>|t|)    
(Intercept)  0.199016   0.402865   0.494 0.621790    
Condition2  -0.107303   0.166441  -0.645 0.519787    
Condition3  -0.259154   0.165977  -1.561 0.119844    
Condition4   0.001307   0.166773   0.008 0.993755    
Condition5  -0.572274   0.167330  -3.420 0.000744 ***
Condition6  -0.180062   0.170184  -1.058 0.291176    
Condition7  -0.411715   0.164065  -2.509 0.012799 *  
Condition8  -0.422063   0.177381  -2.379 0.018179 *  
BMI          0.003585   0.006122   0.586 0.558692    
SexMale      0.015283   0.084327   0.181 0.856345    
Age          0.010749   0.004628   2.323 0.021096 *   
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

I am not really interested in the output above that much because Conditions need to be compared between themselves. Therefore, I am fitting the post-hoc; however, I am confused how to make sense of cov.reduce argument in the function.

Help page around this is sparse and it states: Using cov.reduce = FALSE specifies that, instead of using the mean, the reference grid should use all the unique values of each covariate.

Leaving cov.reduce as TRUE (default) results in this:

emm_model1 <- emmeans(model1, "Condition", cov.reduce = TRUE)
 Condition  emmean    SE  df lower.CL upper.CL
 1          0.2386 0.125 224  -0.0068   0.4841
 2          0.1313 0.112 224  -0.0886   0.3512
 3         -0.0205 0.112 224  -0.2418   0.2007
 4          0.2399 0.105 224   0.0336   0.4462
 5         -0.3336 0.112 224  -0.5548  -0.1125
 6          0.0586 0.116 224  -0.1692   0.2863
 7         -0.1731 0.109 224  -0.3881   0.0419
 8         -0.1834 0.126 224  -0.4318   0.0650

Results are averaged over the levels of: Sex 
Confidence level used: 0.95 

I don't understand why conditions are averaged only on the Sex predictor though. I've also tried cov.reduce = FALSE which averages across all unique levels but as there are numerical variables the covariance matrix is large and it may crash.

Intuitively I feel that it should average across all variables, given that, the best result was achieved with cov.reduce = range resulting in:

 Condition  emmean    SE  df lower.CL upper.CL
 1          0.3033 0.158 224 -0.00712   0.6137
 2          0.1960 0.148 224 -0.09487   0.4869
 3          0.0442 0.145 224 -0.24071   0.3290
 4          0.3046 0.145 224  0.01920   0.5900
 5         -0.2690 0.151 224 -0.56730   0.0294
 6          0.1232 0.137 224 -0.14725   0.3937
 7         -0.1084 0.153 224 -0.40923   0.1924
 8         -0.1188 0.158 224 -0.42940   0.1919

Results are averaged over the levels of: BMI, Sex, Age

However, I still don't fully understand what is the impact of this argument. The help page doesn't do it justice and it seems people are not really concerned with it. I've noticed it has only impact on emm_model1; however contrasts are not affected at all (see below brief example).

My next step to produce the contrast with (contrast_testing is list of pre defined contrasts and emm_model1 was specified above):

contrast(object = emm_model1, contrast_testing, adjust = "holm")

    contrast           estimate    SE  df t.ratio p.value
     Condition 1 vs 7     0.4117 0.164 224  2.509  0.1536
     etc. etc. etc. etc. etc. etc. (I've cut out the rest)

    Results are averaged over the levels of: BMI, Sex, Age

No matter how I fiddle with cov.reduce, results are always the same.

My questions are "What is the purpose of cov.reduce?", "How should I understand it and interpret it?". "If someone asks why does it average on Sex or other variables only, what would be the answer?"

Thank you for any help.


cov.reduce only affects how numeric predictors (in your case, BMI and Age) are handled. If you don't specify cov.reduce at all (which I recommend, because you don't want to examine these further), predictions are made at the averages of those predictors. So cov.reduce is not an issue in your analysis at all, in my opinion.

If you will look at the "basics" vignette -- vignette("basics") or on the CRAN page, you will read that estimates are based on a regular grid (called the reference grid) consisting of all possible combinations of the levels of the factors. If more than one value of any covariate is specified, then those come into play too; but we're not doing that; so in your example, we get a grid with 16 nodes, 8 conditions and 2 sexes. If you ask for EMMs for Condition, then Sex is averaged over.

Your model has only linear terms forf BMI and Age, so holding values of those predictors fixed, the predictions for Condition and Age will be exactly the same relative to one another -- they'll just be shifted up or down by a fixed amount, the same at all 16 nodes of the grid. So you don't need to account for any more than the means of those variables.

One concern I'll mention, though, is if the model is adequate. Have you tried fitting a model with the interaction of Condition and Age? And have you looked at residual plots? It's important to have the model right before going into emmeans or other post hoc analyses.

I invite you to read some of the vignettes associated with the package. There is also an index of topics.

  • $\begingroup$ Thank you for your explanation. I will read more vignettes associated with the package but you've provided the answer I was looking for. $\endgroup$ – gofraidh Jun 6 '19 at 12:17

Please remember this is a statistics site, not a code one. However, a quick look at the emmeans documentation package reveals the following:

Using cov.reduce

cov.reduce may be a function, logical value, formula, or a named list of these. If a single function, it is applied to each covariate.

If logical and TRUE, mean is used. If logical and FALSE, it is equivalent to specifying ‘function(x) sort(unique(x))’, and these values are considered part of the reference grid; thus, it is a handy alternative to specifying these same values in at.

If a formula (which must be two-sided), then a model is fitted to that formula using lm; then in the reference grid, its response variable is set to the results of predict for that model, with the reference grid as newdata. (This is done after the reference grid is determined.) A formula is appropriate here when you think experimental conditions affect the covariate as well as the response.

If cov.reduce is a named list, then the above criteria are used to determine what to do with covariates named in the list. (However, formula elements do not need to be named, as those names are determined from the formulas’ left-hand sides.) Any unresolved covariates are reduced using "mean".

Any cov.reduce specification for a covariate also named in at is ignored.

  • $\begingroup$ Hello thank you for answering. I apologise if my question came across as coding question, it is statistical in the sense that I am interested in what cov.reduce does and how it impacts the analysis. The output above states: "Results are averaged over the levels of: Sex" if cov.reduce = TRUE; however, I do not understand the meaning of this. If mean of Sex is used (presumably mean of male and females), what happens with the other covariates? $\endgroup$ – gofraidh May 31 '19 at 12:07

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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