Differences between lsmeans and difflsmeans

Question

I've created an lmer model. The effects of one of my treatments, "Fertilizer" varies depending on whether I use lsmeans or difflsmeans from lmerTest.

Which should I trust, lsmeans or difflsmeans?

I have a split plot experiment. Fertilizer is the whole plot, Harvest is the split plot. Within each plot, there are 10 samples (quad). Data collected over 4 years, which I have made ordinal. The response is total biomass, which I have transformed.

> model <- lmer(logTotal ~ ordYear*Fertilizer*Harvest + (1|(Block:Fertilizer)) + (1|(Plot:Quad)), REML = FALSE, data = df)    


> difflsmeans(model, test.eff="Fertilizer")
Differences of LSMEANS:
                              Estimate Standard Error   DF t-value Lower CI Upper CI p-value   
Fertilizer none-recommended       -0.3         0.0868 16.0   -3.27  -0.4677  -0.0996   0.005 **
Fertilizer none-half              -0.2         0.0868 16.0   -2.19  -0.3742  -0.0060   0.044 * 
Fertilizer none-manure            -0.1         0.0868 16.0   -0.84  -0.2572   0.1110   0.412   
Fertilizer recommended-half        0.1         0.0868 16.0    1.08  -0.0905   0.2777   0.297   
Fertilizer recommended-manure      0.2         0.0868 16.0    2.42   0.0264   0.3946   0.028 * 
Fertilizer half-manure             0.1         0.0868 16.0    1.35  -0.0671   0.3011   0.197   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


> lsmeans(model, list(pairwise~Fertilizer))
NOTE: Results may be misleading due to involvement in interactions
$`lsmeans of Fertilizer`
 Fertilizer    lsmean         SE    df lower.CL upper.CL
 none        4.602075 0.07090356 21.33 4.454764 4.749387
 recommended 4.885723 0.07090356 21.33 4.738411 5.033034
 half        4.792150 0.07090356 21.33 4.644838 4.939462
 manure      4.675189 0.07090356 21.33 4.527877 4.822500

Results are averaged over the levels of: ordYear, Harvest 
Confidence level used: 0.95 

$`pairwise differences of Fertilizer`
 contrast                estimate        SE    df t.ratio p.value
 none - recommended   -0.28364706 0.1002728 21.33  -2.829  0.0455
 none - half          -0.19007464 0.1002728 21.33  -1.896  0.2595
 none - manure        -0.07311311 0.1002728 21.33  -0.729  0.8843
 recommended - half    0.09357242 0.1002728 21.33   0.933  0.7875
 recommended - manure  0.21053396 0.1002728 21.33   2.100  0.1853
 half - manure         0.11696153 0.1002728 21.33   1.166  0.6537

Results are averaged over the levels of: ordYear, Harvest 
P value adjustment: tukey method for a family of 4 means 

Thank you

Answer 1

These are different tests you're performing. LSMEANS calculates least square means for the different groups and tests whether those means differ from 0. DIFFLSMEANS uses these means and can help you to compare group means (e.g. conditions). In your case, it makes more sense to compare conditions with each other, so DIFFLSMEANS would be the test to consider. However, keep in mind that these are post-hoc tests, so first look at the summary/anova table whether you main effect of fertilizer is significant.

Answer 2

This makes me wonder what difflsmeans is doing. By its name, it should be differences of least-squares means, so I don't understand the distinction @Danyou is making (but perhaps overlooked the second part of the lsmeans output for the pairwise differences?). But obviously these are not doing the same thing, the difflsmeans results being very much rounder and with different degrees of freedom.

I do know that in this example, lsmeanscomputes predictions at each combination of ordYear, Fertilizer, and Harvest (assuming all 3 are factors) based on the fixed-effects part of the model. Then it averages together the predictions for each Fertilizer level, giving equal weights to each combination of the other factors. I wonder if the data are unbalanced and if difflsmeans does not give equal weight. You might experiment with the lsmeans call with different weights arguments to see if you can reproduce the difflsmeans results. For example, weights = "prop" or weights = "outer" -- or even weights = "cells" which produces the ordinary marginal means of the data.

If you ask me, I'd trust lsmeans. But I'm biased because I'm its author.

PS -- Another thing I'll note is that the lsmeans results use a Tukey adjustment (this is the defailt for pairwise comparisons, but you can change to some other adjustment), and the difflsmeans intervals and tests appear to be unadjusted for multiplicity.