0
$\begingroup$

I'm having an issue with the emmeans package in R, in which some of the pairwise comparisons on my model have zero degrees of freedom. I'm fitting a negative binomial mixed effects glm in which the abundance of whelks (marine snails) depends on the region and year they were collected in. Each observation is a transect, and there are multiple transects within each site, and multiple sites within both regions (Historic and Expanded). Here is my model code:

library(MASS)
library(emmeans)

whelks_60_glm <- glmmPQL(count ~ Year*Region, 
                     random =  ~ 1 | Region/Site, 
                     family = negative.binomial(theta = 
                     theta(whelks_density_60$count)),
                     data = whelks_density_60)

car::Anova(whelks_60_glm)

emmeans(whelks_60_glm, pairwise ~ Year*Region)$contrasts

Here is the model output from car::Anova:

|            |     Chisq| Df| Pr(>Chisq)|
|:-----------|---------:|--:|----------:|
|Year        | 13.624262|  2|  0.0011003|
|Region      |  2.313283|  1|  0.1282729|
|Year:Region | 22.061868|  2|  0.0000162|

And the pairwise output from emmeans:

|contrast                      |   estimate|        SE|  df|    t.ratio|   p.value|
|:-----------------------------|----------:|---------:|---:|----------:|---------:|
|2015,Expanded - 2016,Expanded | -1.8945578| 0.4044370| 573| -4.6844325| 0.0000514|
|2015,Expanded - 2017,Expanded | -2.0242238| 0.4080372| 573| -4.9608807| 0.0000137|
|2015,Expanded - 2015,Historic | -0.0151271| 0.8971538|   0| -0.0168612|       NaN|
|2015,Expanded - 2016,Historic |  0.1553911| 0.9014864|   0|  0.1723721|       NaN|
|2015,Expanded - 2017,Historic | -0.3585129| 0.8985426|   0| -0.3989938|       NaN|
|2016,Expanded - 2017,Expanded | -0.1296660| 0.3786464| 573| -0.3424461| 0.9993783|
|2016,Expanded - 2015,Historic |  1.8794307| 0.8840423|   0|  2.1259510|       NaN|
|2016,Expanded - 2016,Historic |  2.0499489| 0.8884389|   0|  2.3073606|       NaN|
|2016,Expanded - 2017,Historic |  1.5360449| 0.8854516|   0|  1.7347586|       NaN|
|2017,Expanded - 2015,Historic |  2.0090967| 0.8856627|   0|  2.2684671|       NaN|
|2017,Expanded - 2016,Historic |  2.1796149| 0.8900513|   0|  2.4488645|       NaN|
|2017,Expanded - 2017,Historic |  1.6657109| 0.8870695|   0|  1.8777683|       NaN|
|2015,Historic - 2016,Historic |  0.1705182| 0.2280432| 573|  0.7477450| 0.9757638|
|2015,Historic - 2017,Historic | -0.3433859| 0.2169077| 573| -1.5830964| 0.6100892|
|2016,Historic - 2017,Historic | -0.5139040| 0.2253690| 573| -2.2802777| 0.2038030|

It looks like the degrees of freedom are zero for any comparison between different "Regions" in the model. Is this a product of the nested design? How might I be able to get all pairwise comparisons for the Region*Year interaction in this model?

Thanks!

$\endgroup$
  • $\begingroup$ A variable should not be considered a fixed effect as well as a grouping variable of the random effect. It's one or the other. $\endgroup$ – Roland May 22 '18 at 6:39
  • $\begingroup$ Thanks! Changing ~ 1 | Region / Site to 1 | Site worked. If sites are nested within regions, is there any issue with doing this? No site names occur in both regions, so the sites are explicitly nested. $\endgroup$ – Maurice Goodman May 22 '18 at 18:27

Your Answer

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

Browse other questions tagged or ask your own question.