Reference

Least-squares means with interaction

Without interaction

at a special value

# Compute least-squares means for specified factors or factor combinations in a linear model

require(lsmeans)
## Loading required package: lsmeans
## Loading required package: emmeans
## The 'lsmeans' package is now basically a front end for 'emmeans'.
## Users are encouraged to switch the rest of the way.
## See help('transition') for more information, including how to
## convert old 'lsmeans' objects and scripts to work with 'emmeans'.
# Use dataset 'fiber'

### model 1 without interaction, Covariance example 
fiber.lm <- lm(strength ~ diameter + machine, data = fiber)

# adjusted means and comparison
lsmeans(fiber.lm,   list (pairwise ~ machine ) ) 
## $`lsmeans of machine`
##  machine lsmean    SE df lower.CL upper.CL
##  A         40.4 0.724 11     38.8     42.0
##  B         41.4 0.744 11     39.8     43.1
##  C         38.8 0.788 11     37.1     40.5
## 
## Confidence level used: 0.95 
## 
## $`pairwise differences of machine`
##  1     estimate   SE df t.ratio p.value
##  A - B    -1.04 1.01 11  -1.024  0.5781
##  A - C     1.58 1.11 11   1.431  0.3596
##  B - C     2.62 1.15 11   2.283  0.1005
## 
## P value adjustment: tukey method for comparing a family of 3 estimates
# adjusted means and contrast, treating machine C as control
lsmeans (fiber.lm, "machine", contr = "trt.vs.ctrlk")
## $lsmeans
##  machine lsmean    SE df lower.CL upper.CL
##  A         40.4 0.724 11     38.8     42.0
##  B         41.4 0.744 11     39.8     43.1
##  C         38.8 0.788 11     37.1     40.5
## 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast estimate   SE df t.ratio p.value
##  A - C        1.58 1.11 11   1.431  0.3062
##  B - C        2.62 1.15 11   2.283  0.0796
## 
## P value adjustment: dunnettx method for 2 tests
# adjusted means and comparison, after setting "diameter" as a specific value
lsmeans(fiber.lm,   list (pairwise ~ machine ), at = list(diameter =20) ) 
## $`lsmeans of machine`
##  machine lsmean    SE df lower.CL upper.CL
##  A         36.4 0.928 11     34.4     38.5
##  B         37.5 0.988 11     35.3     39.7
##  C         34.9 0.726 11     33.3     36.5
## 
## Confidence level used: 0.95 
## 
## $`pairwise differences of machine`
##  1     estimate   SE df t.ratio p.value
##  A - B    -1.04 1.01 11  -1.024  0.5781
##  A - C     1.58 1.11 11   1.431  0.3596
##  B - C     2.62 1.15 11   2.283  0.1005
## 
## P value adjustment: tukey method for comparing a family of 3 estimates

With interaction

### model 2 with interaction
fiber.lm2 <- lm(strength ~ diameter*machine, data=fiber)

# means when diameter is its mean
(lsmeans(fiber.lm2, ~ machine | diameter))
## diameter = 24.1:
##  machine lsmean    SE df lower.CL upper.CL
##  A         40.2 0.777  9     38.5     42.0
##  B         41.6 0.858  9     39.7     43.5
##  C         38.5 0.966  9     36.3     40.7
## 
## Confidence level used: 0.95
# means when setting diameter as its mean
lsmeans(fiber.lm2,   list (pairwise ~ machine ), at = list(diameter =24.1) ) 
## NOTE: Results may be misleading due to involvement in interactions
## $`lsmeans of machine`
##  machine lsmean    SE df lower.CL upper.CL
##  A         40.2 0.779  9     38.4     41.9
##  B         41.6 0.861  9     39.6     43.5
##  C         38.5 0.962  9     36.3     40.7
## 
## Confidence level used: 0.95 
## 
## $`pairwise differences of machine`
##  1     estimate   SE df t.ratio p.value
##  A - B    -1.39 1.16  9  -1.194  0.4858
##  A - C     1.68 1.24  9   1.357  0.4017
##  B - C     3.07 1.29  9   2.374  0.0956
## 
## P value adjustment: tukey method for comparing a family of 3 estimates
### model 3 for category covariates
warp.lm <- lm(breaks ~ wool * tension, data = warpbreaks)

# means after adjusting the interaction effect
(lsmeans(warp.lm, ~ wool | tension))
## tension = L:
##  wool lsmean   SE df lower.CL upper.CL
##  A      44.6 3.65 48     37.2     51.9
##  B      28.2 3.65 48     20.9     35.6
## 
## tension = M:
##  wool lsmean   SE df lower.CL upper.CL
##  A      24.0 3.65 48     16.7     31.3
##  B      28.8 3.65 48     21.4     36.1
## 
## tension = H:
##  wool lsmean   SE df lower.CL upper.CL
##  A      24.6 3.65 48     17.2     31.9
##  B      18.8 3.65 48     11.4     26.1
## 
## Confidence level used: 0.95
# set tension ="M"
lsmeans(warp.lm,   list (pairwise ~ wool ), at = list(tension ="M")) 
## NOTE: Results may be misleading due to involvement in interactions
## $`lsmeans of wool`
##  wool lsmean   SE df lower.CL upper.CL
##  A      24.0 3.65 48     16.7     31.3
##  B      28.8 3.65 48     21.4     36.1
## 
## Confidence level used: 0.95 
## 
## $`pairwise differences of wool`
##  1     estimate   SE df t.ratio p.value
##  A - B    -4.78 5.16 48  -0.926  0.3589
lsmeans(warp.lm,   list (pairwise ~ wool ), at = list(tension ="H")) 
## NOTE: Results may be misleading due to involvement in interactions
## $`lsmeans of wool`
##  wool lsmean   SE df lower.CL upper.CL
##  A      24.6 3.65 48     17.2     31.9
##  B      18.8 3.65 48     11.4     26.1
## 
## Confidence level used: 0.95 
## 
## $`pairwise differences of wool`
##  1     estimate   SE df t.ratio p.value
##  A - B     5.78 5.16 48   1.120  0.2682

Compare slopes and test

### Trends, obtain slopes  
fiber.lst=( lstrends(fiber.lm, "machine", var="diameter") )
fiber.lst
##  machine diameter.trend    SE df lower.CL upper.CL
##  A                0.954 0.114 11    0.703     1.21
##  B                0.954 0.114 11    0.703     1.21
##  C                0.954 0.114 11    0.703     1.21
## 
## Confidence level used: 0.95
# pairwise comparisons thereof slopes 
pairs(fiber.lst)
##  contrast estimate SE df t.ratio p.value
##  A - B           0  0 11     NaN     NaN
##  A - C           0  0 11     NaN     NaN
##  B - C           0  0 11     NaN     NaN
## 
## P value adjustment: tukey method for comparing a family of 3 estimates
fiber.lst2=( lstrends(fiber.lm2, "machine", var="diameter") )
fiber.lst2
##  machine diameter.trend    SE df lower.CL upper.CL
##  A                1.104 0.194  9    0.666     1.54
##  B                0.857 0.224  9    0.351     1.36
##  C                0.864 0.208  9    0.394     1.33
## 
## Confidence level used: 0.95

Least-squares means with interaction for LMM

library(lme4)
## Loading required package: Matrix
##############
m1 <- lmer(Reaction ~ Days + (1 | Subject), sleepstudy)
lsmeans(m1,    ~ Days  , at = list(Days =20)  ) 
##  Days lsmean   SE  df lower.CL upper.CL
##    20    461 15.4 103      430      491
## 
## Degrees-of-freedom method: kenward-roger 
## Confidence level used: 0.95
lsmeans(m1,    ~ Days    ) 
##  Days lsmean   SE df lower.CL upper.CL
##   4.5    299 9.05 17      279      318
## 
## Degrees-of-freedom method: kenward-roger 
## Confidence level used: 0.95