I have an experimental design in which I measured the abundance of a set of invertebrates in 10 sites, grouped in two categories depending on the presence of an invasive species (5 invaded sites vs. 5 uninvaded sites). Each of these 10 sites was visited three times across the year to account for potential seasonal effects. Now, after performing a GLMM in which the interaction invasion status x season was significant, I explored posthoc comparisons using emmeans.
Since I am solely interested in exploring how invaded and uninvaded sites differ in each season --I have no interest in testing how invertebrate abundance varies across the year in uninvaded (or invaded sites)-- I explored posthoc comparisons usins
pairs(emmeans(model,~invasion status|season))
However, this means invaded and uninvaded sites are compared for each season separately, with a single test in each case. Consequently, I cannot apply p-value adjustments for multiple comparisons, even though three tests are actually being considered (one for each season). I assume this must be correct, or the package would otherwise run differently. However, I would like to know if this is safe to assume and whether there is a way to control for these multiple comparisons (e.g., by grouping tests, or some other way I am unaware). Thank you in advance.
The rbind() method for emmGrid objects allows you to specify a p-value adjustment to combine comparisons that are made across multiple factor levels. See documentation and vignette. You can do it like in this example which uses the Bonferroni correction but other options are possible. To see all the adjustment options look at the "p-value adjustments" section in the documentation for ?summary.emmGrid.
comparisons <- pairs(emmeans(model,~invasionstatus|season))
rbind(comparisons, adjust = 'bonferroni')