Generate matrix specifying linear estimate.
LE_matrix(object, effect = NULL, at = NULL)
# Default S3 method
LE_matrix(object, effect = NULL, at = NULL)
aggregate_linest_list(linest_list)
get_linest_list(object, effect = NULL, at = NULL)
Model object
A vector of variables. For each configuration of these the estimate will be calculated.
Either NULL, a list or a dataframe. 1) If a list, then the list must consist of covariates (including levels of some factors) to be used in the calculations. 2) If a dataframe, the dataframe is split rowwise and the function is invoked on each row.
Linear estimate list (as generated by get_linest_list
).
Check this
## 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"))
#> List of 2
#> $ new.fact.lev:List of 2
#> ..$ wool : chr [1:2] "A" "B"
#> ..$ tension: chr "M"
#> $ grid.data :'data.frame': 2 obs. of 2 variables:
#> ..$ wool : chr [1:2] "A" "B"
#> ..$ tension: chr [1:2] "M" "M"
K
#> (Intercept) woolB tensionM tensionH
#> [1,] 1 0 1 0
#> [2,] 1 1 1 0
## Vanilla computation:
K %*% coef(m0)
#> [,1]
#> [1,] 29.27778
#> [2,] 23.50000
## Alternative; also providing standard errors etc:
linest(m0, K)
#> wool tension estimate std.error statistic df p.value lwr upr
#> 1 A M 29.3 3.16 9.26 50 2.00e-12 22.9 35.6
#> 2 B M 23.5 3.16 7.43 50 1.27e-09 17.1 29.9
esticon(m0, K)
#> estimate std.error statistic p.value beta0 df
#> [1,] 2.9278e+01 3.1618e+00 9.2599e+00 2.0017e-12 0.0000e+00 50
#> [2,] 2.3500e+01 3.1618e+00 7.4325e+00 1.2661e-09 0.0000e+00 50
## 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")))
#> List of 2
#> $ new.fact.lev:List of 1
#> ..$ wool: chr [1:2] "A" "B"
#> $ grid.data :'data.frame': 2 obs. of 1 variable:
#> ..$ wool: chr [1:2] "A" "B"
K
#> (Intercept) woolB tensionM tensionH
#> [1,] 1 0 0.3333333 0.3333333
#> [2,] 1 1 0.3333333 0.3333333
K <- LE_matrix(m0, effect="wool")
#> List of 2
#> $ new.fact.lev:List of 1
#> ..$ wool: chr [1:2] "A" "B"
#> $ grid.data :'data.frame': 2 obs. of 1 variable:
#> ..$ wool: chr [1:2] "A" "B"
K
#> (Intercept) woolB tensionM tensionH
#> [1,] 1 0 0.3333333 0.3333333
#> [2,] 1 1 0.3333333 0.3333333
linest(m0, K)
#> wool estimate std.error statistic df p.value lwr upr
#> 1 A 31.0 2.24 13.9 50 8.81e-19 26.5 35.5
#> 2 B 25.3 2.24 11.3 50 2.26e-15 20.8 29.7
## The linear estimate is sometimes called to "least squares mean"
# (LSmeans) or popupulation means.
# Same as
LSmeans(m0, effect="wool")
#> List of 2
#> $ new.fact.lev:List of 1
#> ..$ wool: chr [1:2] "A" "B"
#> $ grid.data :'data.frame': 2 obs. of 1 variable:
#> ..$ wool: chr [1:2] "A" "B"
#> wool estimate std.error statistic df p.value lwr upr
#> 1 A 31.0 2.24 13.9 50 8.81e-19 26.5 35.5
#> 2 B 25.3 2.24 11.3 50 2.26e-15 20.8 29.7
## Without mentioning 'effect' or 'at' an average across all
#predictors are calculated:
K <- LE_matrix(m0)
K
#> (Intercept) woolB tensionM tensionH
#> [1,] 1 0.5 0.3333333 0.3333333
linest(m0, K)
#> estimate std.error statistic df p.value lwr upr
#> 1 28.1 1.58 17.8 50 2.74e-23 25 31.3
## 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)
#> [1] 28.14815
## Same as
LSmeans(m0)
#> estimate std.error statistic df p.value lwr upr
#> 1 28.1 1.58 17.8 50 2.74e-23 25 31.3
## An interaction model
m1 <- lm(breaks ~ wool * tension, data=warpbreaks)
K <- LE_matrix(m1, at=list(wool=c("A", "B"), tension="M"))
#> List of 2
#> $ new.fact.lev:List of 2
#> ..$ wool : chr [1:2] "A" "B"
#> ..$ tension: chr "M"
#> $ grid.data :'data.frame': 2 obs. of 2 variables:
#> ..$ wool : chr [1:2] "A" "B"
#> ..$ tension: chr [1:2] "M" "M"
K
#> (Intercept) woolB tensionM tensionH woolB:tensionM woolB:tensionH
#> [1,] 1 0 1 0 0 0
#> [2,] 1 1 1 0 1 0
linest(m1, K)
#> wool tension estimate std.error statistic df p.value lwr upr
#> 1 A M 24.0 3.65 6.58 48 3.23e-08 16.7 31.3
#> 2 B M 28.8 3.65 7.89 48 3.22e-10 21.4 36.1
K <- LE_matrix(m1, at=list(wool=c("A", "B")))
#> List of 2
#> $ new.fact.lev:List of 1
#> ..$ wool: chr [1:2] "A" "B"
#> $ grid.data :'data.frame': 2 obs. of 1 variable:
#> ..$ wool: chr [1:2] "A" "B"
K
#> (Intercept) woolB tensionM tensionH woolB:tensionM woolB:tensionH
#> [1,] 1 0 0.3333333 0.3333333 0.0000000 0.0000000
#> [2,] 1 1 0.3333333 0.3333333 0.3333333 0.3333333
linest(m1, K)
#> wool estimate std.error statistic df p.value lwr upr
#> 1 A 31.0 2.11 14.7 48 1.91e-19 26.8 35.3
#> 2 B 25.3 2.11 12.0 48 4.71e-16 21.0 29.5
K <- LE_matrix(m1, effect="wool")
#> List of 2
#> $ new.fact.lev:List of 1
#> ..$ wool: chr [1:2] "A" "B"
#> $ grid.data :'data.frame': 2 obs. of 1 variable:
#> ..$ wool: chr [1:2] "A" "B"
K
#> (Intercept) woolB tensionM tensionH woolB:tensionM woolB:tensionH
#> [1,] 1 0 0.3333333 0.3333333 0.0000000 0.0000000
#> [2,] 1 1 0.3333333 0.3333333 0.3333333 0.3333333
linest(m1, K)
#> wool estimate std.error statistic df p.value lwr upr
#> 1 A 31.0 2.11 14.7 48 1.91e-19 26.8 35.3
#> 2 B 25.3 2.11 12.0 48 4.71e-16 21.0 29.5
LSmeans(m1, effect="wool")
#> List of 2
#> $ new.fact.lev:List of 1
#> ..$ wool: chr [1:2] "A" "B"
#> $ grid.data :'data.frame': 2 obs. of 1 variable:
#> ..$ wool: chr [1:2] "A" "B"
#> wool estimate std.error statistic df p.value lwr upr
#> 1 A 31.0 2.11 14.7 48 1.91e-19 26.8 35.3
#> 2 B 25.3 2.11 12.0 48 4.71e-16 21.0 29.5
K <- LE_matrix(m1)
K
#> (Intercept) woolB tensionM tensionH woolB:tensionM woolB:tensionH
#> [1,] 1 0.5 0.3333333 0.3333333 0.1666667 0.1666667
linest(m1, K)
#> estimate std.error statistic df p.value lwr upr
#> 1 28.1 1.49 18.9 48 6.98e-24 25.2 31.1
LSmeans(m1)
#> estimate std.error statistic df p.value lwr upr
#> 1 28.1 1.49 18.9 48 6.98e-24 25.2 31.1