Standard error in estimated marginal means are all the same

Question

I use the package emmeans to calculate estimated marginal means and I don't know why the standard errors are equal within the factors:

> warp.lm <- lm(breaks ~ wool * tension, data = warpbreaks)
> emmeans (warp.lm,  ~ wool | tension)
tension = L:
 wool emmean   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 emmean   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 emmean   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 
> # or equivalently emmeans(warp.lm, "wool", by = "tension")
> 
> emmeans (warp.lm, poly ~ tension | wool)
$`emmeans`
wool = A:
 tension emmean   SE df lower.CL upper.CL
 L         44.6 3.65 48     37.2     51.9
 M         24.0 3.65 48     16.7     31.3
 H         24.6 3.65 48     17.2     31.9

wool = B:
 tension emmean   SE df lower.CL upper.CL
 L         28.2 3.65 48     20.9     35.6
 M         28.8 3.65 48     21.4     36.1
 H         18.8 3.65 48     11.4     26.1

Confidence level used: 0.95 

$contrasts
wool = A:
 contrast  estimate   SE df t.ratio p.value
 linear      -20.00 5.16 48 -3.878  0.0003 
 quadratic    21.11 8.93 48  2.363  0.0222 

wool = B:
 contrast  estimate   SE df t.ratio p.value
 linear       -9.44 5.16 48 -1.831  0.0733 
 quadratic   -10.56 8.93 48 -1.182  0.2432 

Why is every SE equal to 3.65? For the contrasts I get the same SE for each wool but not for linear or quadratic. How can I calculate these values "by hand"?

Answer 1

Because it is a balanced experiment, and you are using a model that presumes the error variance is homogeneous. Accordingly, the SE of each cell mean is $s/\sqrt n$ where $s$ is the estimated error SD and $n$ is the number of observations in each mean.

In this particular example, there are 6 cell means with equal counts, and you can get $s$ from the model:

> nrow(warpbreaks) / 6
[1] 9

> sigma(warp.lm) / sqrt(9)
[1] 3.646761
```