```
# NOT RUN {
## Two way anova:
data(warpbreaks)
## An additive model
m0 <- lm(breaks ~ wool + tension, data=warpbreaks)
## Estimate mean for each wool type, for tension="M":
K <- LE_matrix(m0, at=list(wool=c("A", "B"), tension="M"))
K
## Vanilla computation:
K %*% coef(m0)
## Alternative; also providing standard errors etc:
linest(m0, K)
esticon(m0, K)
## Estimate mean for each wool type when averaging over tension;
# two ways of doing this
K <- LE_matrix(m0, at=list(wool=c("A", "B")))
K
K <- LE_matrix(m0, effect="wool")
K
linest(m0, K)
## The linear estimate is sometimes called to "least squares mean"
# (LSmeans) or popupulation means.
# Same as
LSmeans(m0, effect="wool")
## Without mentioning 'effect' or 'at' an average across all
#predictors are calculated:
K <- LE_matrix(m0)
K
linest(m0, K)
## Because the design is balanced (9 observations per combination
#of wool and tension) this is the same as computing the average. If
#the design is not balanced, the two quantities are in general not
#the same.
mean(warpbreaks$breaks)
## Same as
LSmeans(m0)
## An interaction model
m1 <- lm(breaks ~ wool * tension, data=warpbreaks)
K <- LE_matrix(m1, at=list(wool=c("A", "B"), tension="M"))
K
linest(m1, K)
K <- LE_matrix(m1, at=list(wool=c("A", "B")))
K
linest(m1, K)
K <- LE_matrix(m1, effect="wool")
K
linest(m1, K)
LSmeans(m1, effect="wool")
K <- LE_matrix(m1)
K
linest(m1, K)
LSmeans(m1)
# }
```

Run the code above in your browser using DataLab