testFactors: Evaluate and Test Combinations of Factor Levels
In phia: Post-Hoc Interaction Analysis
Description
Usage
Arguments
Details
Value
Note
Author(s)
See Also
Examples
Description
Calculates and tests the adjusted mean value of the response and other terms of a fitted model, for specific linear combinations of factor levels and specific values of the covariates. This function is specially indicated for post-hoc analyses of models with factors, to test pairwise comparisons of factor levels, simple main effects, and interaction contrasts. In multivariate models, it is possible to define and test contrasts of intra-subjects factors, as in a repeated-measures analysis.
Usage
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
## Default S3 method:
testFactors(model, levels, covariates, offset, terms.formula=~1,
inherit.contrasts=FALSE, default.contrasts=c("contr.sum","contr.poly"),
lht=TRUE, ...)
## S3 method for class 'lm'
testFactors(model, ...)
## S3 method for class 'glm'
testFactors(model, ..., link=FALSE)
## S3 method for class 'mlm'
testFactors(model, levels, covariates, offset, terms.formula=~1,
inherit.contrasts=FALSE, default.contrasts=c("contr.sum","contr.poly"),
idata, icontrasts=default.contrasts, lht=TRUE, ...)
## S3 method for class 'lme'
testFactors(model, ...)
## S3 method for class 'mer'
testFactors(model, ..., link=FALSE)
## S3 method for class 'merMod'
testFactors(model, ..., link=FALSE)
## S3 method for class 'testFactors'
summary(object, predictors=TRUE, matrices=TRUE, covmat=FALSE, ...)
Arguments
model
fitted model. Currently supported classes include "lm", "glm", "mlm", "lme", and "mer" or "merMod" (excluding models fitted by nlmer).
levels
list that defines the factor levels that will be contrasted; see the Details for more information.
covariates
optional numeric vector with specific values of the covariates (see Details).
offset
optional numeric value with a specific value of the offset (if any).
terms.formula
formula that defines the terms of the model that will be calculated and tested. The default value is ~1, and stands for the adjusted mean value. See the Details for more information.
inherit.contrasts
logical value: should the default contrasts of model factors be inherited from the model data frame?
default.contrasts
names of contrast-generating functions to be applied by default to factors and ordered factors, respectively, if inherit.contrasts is FALSE (the default); the contrasts must produce an intra-subjects model matrix in which different terms are orthogonal. The default is c("contr.sum", "contr.poly").
lht
logical indicating if the adjusted values are tested (via linearHypothesis).
link
for models fitted with glm or glmer, logical indicating if the adjusted mean values should represent the link function (FALSE by default, i.e. represent the adjusted means of the response variable).
idata
an optional data frame giving a factor or factors defining the intra-subjects model for multivariate repeated-measures data, as defined in Anova or linearHypothesis.
icontrasts
names of contrast-generating functions to be applied in the within-subject “data”. The default is the same as default.contrasts.
object
object returned by testFactors.
predictors
logical value: should summary return the values of the predictor values used in the calculations?
matrices
logical value: should summary return the matrices used for testing by linearHypothesis?
covmat
logical value: should summary return the covariance matrix of the adjusted values?
...
other arguments passed down to linearHypothesis.
Details
The only mandatory argument is model, which may include any number of factor or numeric predictors, and one offset. The simplest usage of this method, where no other argument is defined, calculates the adjusted mean of the model response variable, pooling over all the levels of factor predictors, and setting the numeric predictors (covariates and offset, ifany) to their average values in the model data frame.
The calculations will be done for the linear combinations of factor levels defined by levels. This argument must be a list, with one element for each factor of the model that has to be manipulated (including factors of the intra-subjects design, if suitable). The factors that are not represented in this list will be pooled over, and elements that do not correspond to any factor of the model will be ignored with a warning. levels may be a named list, where the name of each element identifies the represented factor, and its contents may be one of the following:
A character string of length 1 or 2, with the name of one or two factor levels. In the former case, the calculations will be restricted to this level of the factor, as in a simple main effects analysis; in the latter, a pairwise contrast will be calculated between both factor levels.
A numeric vector without names, as long as the number of levels in the factor. This will create a linear combination of factor levels, with the elements of the vector as coefficients. For instance, if the factor f has three levels, an element f=c(0.5, 0.5, 0) will average the two first levels, and f=c(0.5, 0.5, -1) will contrast the average of the two first levels against the third one.
A numeric vector with names equal to some or all the levels of the factor. This is a simplification of the previous option, where some levels can be omitted, and the coefficient of each level is determined by the names of the vector, which do not have to follow a specific order. Omitted levels will automatically be set to zero.
A numeric matrix, as an extension of the two previous options for calculating several combinations at a time. Combinations are defined in columns, so if the matrix does not have row names, the number of rows must be equal to the number of levels in the factor, or if the matrix does have row names, they must coincide with the levels of the factor.
Alternatively, levels may be a single formula or an unnamed list of formulas, of the type factorname ~ K1*level1 + K2*level2 ... (see contrastCoefficients for further details). Both types of lists (named list of string or numeric vectors and matrices, and unnamed lists of formulas) may be mixed.
The argument covariates may be used for setting specific values of the model numeric predictors (covariates). It must be a vector of numeric values. If the elements of this vector have names, their values will be assigned to the covariates that match them; covariates of the model with names not represented in this vector will be set to their default value (the average in the model data frame), and elements with names that do not match with covariates will be ignored. On the other hand, if covariates has no names, and its length is equal to the number of covariates of the model, the values will be assigned to those covariates in the same order as they occur in the model. If it has a different length, the vector will be trimmed or reclycled as needed to fit the number of covariates in the model; this feature may be used, for instance, to set all covariates to the same value, e.g. covariates = 0. The argument offset can likewise be used to define a specific value for the offset of the model.
To analyse terms other than the adjusted mean value of the response, use the argument terms.formula. For instance, if the model has the covariates var1, var2, ..., the slopes of the response with respect to them may be added to the analysis by defining terms.formula as ~var1 + var2 .... This formula may be used more generally, for analysing interactions, omitting the mean response, adding the main effects of factors, etc. A different analysis is done for each term of this formula, that must also be contained in the formula of model. For instance, if terms.formula is equal to ~ var1*var2, the function will analyse the adjusted intercept, plus the terms var1, var2, and var1:var2. The intercept stands for the mean value of the response, and terms formed by one or more covariates stand por the slope of the response with respect to the product of those covariates.
If any of the variables in the term is a factor, the function analyses a full set of contrasts for that factor of the remaining part of the term; for instance if var1 were a factor, the term var1 would stand for the contrasts of the intercept, and var1:var2 would stand for the contrasts of the slope var2, across the levels of var1. The set of contrasts used in the analysis is normally defined by the argument default.contrasts: by default, if the factor is ordered it will be a set of “polynomial contrasts”, and otherwise “sum contrasts”; however, if inherit.contrasts is TRUE the contrasts will directly be copied from the ones used to define the model. Factors that have explicit contrasts defined in the model data frame will use those contrasts, regardless of the values defined for default.contrasts and inherit.contrasts. The analysis assumes that the contrasts are orthogonal to the intercept, which is the usual case if the default arguments are provided, and a warning will be issued if non-orthogonal contrasts are used; take special care of not using “treatment contrats” if inherit.contrasts is set to TRUE or default.contrasts is changed.
In generalized linear models, the adjusted means represent the expected values of the response by default, but the expected value of the link function may be shown by setting the argument link=FALSE. On the other hand, slope values and standard errors always refer to the link function.
For multivariate models, the arguments idata, and icontrasts may be used to define an intra-subjects model for multivariate repeated-measures data, as described for Anova or linearHypothesis in package car. Note, however, that the combinations of intra-subjects factor levels are defined in levels, and other arguments defined in those functions like idesign, imatrix or iterms will have no effect in testFactors.
The significance of adjusted values is tested by a call to linearHypothesis for each term, unless lht is set to FALSE. Extra arguments may be passed down to that function, for instance to specify the test statistic that will be evaluated.
Value
An object of class "testFactors", that contains the adjusted values and their standard errors for each term, and the otuput of the test, plus other variables used in the calculations. The summary method for this object will display those variables, unless they be omitted by setting the optional arguments predictors, matrices or covmat to FALSE. The argument predictors refers to the coefficients of specified combinations of factor levels, the values of covariates, and the contrast matrices used for terms that include factors; matrices refers to the “linear hypothesis matrix” used by linearHypothesis, and in multivariate linear models, to the “response transformation matrix” as well — if it exists; covmat refers to the variance-covariance matrix of the adjusted values.
Moreover, summary groups the results of the tests for all terms in one table. By default this table shows the test statistics, their degrees of freedom, and the p-values. If the model is of class "lm", it also shows the sums of squares; and if it is of class "mlm", only the first type of test statistic returned by linearHypothesis (by default “Pillai”) is shown. This variable shape of the ANOVA table is controlled by additional classes assigned to the object (either "testFactors.lm" or "testFactors.mlm", as suitable).
Note
The tests of mixed models are done under the assumption that the estimation of the random part of the model is exact.
Author(s)
Helios De Rosario-Martinez, helios.derosario@gmail.com
See Also
linearHypothesis in package car.
interactionMeans, and testInteractions as useful wrappers of testFactors.
Examples
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
# Example with factors and covariates
# Analyse prestige of Canadian occupations depending on
# education, income and type of occupation
# See ?Prestige for a description of the data set
prestige.mod <- lm(prestige ~ (education+log2(income))*type, data=Prestige)
# Pairwise comparisons for factor "type", to see how it influences
# the mean value of prestige and interacts with log.income
# 1: "white collar" vs "blue collar"
wc.vs.bc <- list(type=c("wc", "bc"))
testFactors(prestige.mod, wc.vs.bc, terms.formula=~log2(income))
# 2: "professional" vs. "blue collar"
prof.vs.bc <- list(type=c("prof", "bc"))
testFactors(prestige.mod, prof.vs.bc, terms.formula=~log2(income))
# 3: "professional" vs. "white collar"
prof.vs.wc <- list(type=c("prof", "wc"))
testFactors(prestige.mod, prof.vs.wc, terms.formula=~log2(income))
# Interaction contrasts in a repeated-measures experiment
# See ?OBrienKaiser for a description of the data set
mod.ok <- lm(cbind(pre.1, pre.2, pre.3, pre.4, pre.5,
post.1, post.2, post.3, post.4, post.5,
fup.1, fup.2, fup.3, fup.4, fup.5) ~ treatment*gender,
data=OBrienKaiser)
# intra-subjects data:
phase <- factor(rep(c("pretest", "posttest", "followup"), each=5))
hour <- ordered(rep(1:5, 3))
idata <- data.frame(phase, hour)
Anova(mod.ok, idata=idata, idesign=~phase*hour)
# Contrasts across "phase", for different contrasts of "treatment"
# Using different definitions of the argument "levels"
# 1: Treatment "A" vs. treatment "B".
A.vs.B <- list(treatment=c("A", "B"))
# The following are equivalent:
# A.vs.B <- list(treatment=c(A=1, B=-1, control=0))
# A.vs.B <- list(treatment=c(A=1, B=-1))
# A.vs.B <- list(treatment ~ A - B)
# A.vs.B <- treatment ~ A - B
testFactors(mod.ok, A.vs.B, idata=idata, terms.formula=~0+phase)
# 2: Controls vs. treatments
control.vs.AB <- list(treatment=c(A=0.5, B=0.5, control=-1))
# The following is equivalent:
# control.vs.AB <- treatment ~ (A+B)/2 - control
testFactors(mod.ok, control.vs.AB, idata=idata, terms.formula=~0+phase)
# Shortcut to get only the adjusted values and simplified ANOVA tables
contr <- list(A.vs.B=A.vs.B, control.vs.AB=control.vs.AB)
anovaTables <- function(contrast) summary(testFactors(mod.ok, contrast,
idata=idata, terms.formula=~0+phase),
predictors=FALSE, matrices=FALSE)
lapply(contr,anovaTables)
Example output
Loading required package: car
Loading required package: carData
Call: testFactors(model = prestige.mod, levels = wc.vs.bc, terms.formula = ~log2(income))
Term (Intercept)
Adjusted mean:
-1.44478
Std. Error:
2.842545
Linear hypothesis test
Hypothesis:
typewc + 10.7951020408163 education:typewc + 12.5561779064803 log2(income):typewc = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 3666.007
2 89 3655.397 1 10.61044 0.25834 0.61252
------
Term log2(income)
Adjusted slope for log2(income):
-5.653038
Std. Error:
3.051886
Linear hypothesis test
Hypothesis:
log2(income):typewc = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 3796.316
2 89 3655.397 1 140.9197 3.43105 0.067296 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------
Call: testFactors(model = prestige.mod, levels = prof.vs.bc, terms.formula = ~log2(income))
Term (Intercept)
Adjusted mean:
10.62698
Std. Error:
3.784129
Linear hypothesis test
Hypothesis:
typeprof + 10.7951020408163 education:typeprof + 12.5561779064803 log2(income):typeprof = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 3979.313
2 89 3655.397 1 323.916 7.88657 0.0061205 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------
Term log2(income)
Adjusted slope for log2(income):
-6.535559
Std. Error:
2.616708
Linear hypothesis test
Hypothesis:
log2(income):typeprof = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 3911.609
2 89 3655.397 1 256.212 6.23814 0.014341 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------
Call: testFactors(model = prestige.mod, levels = prof.vs.wc, terms.formula = ~log2(income))
Term (Intercept)
Adjusted mean:
12.07176
Std. Error:
3.42693
Linear hypothesis test
Hypothesis:
typeprof - typewc + 10.7951020408163 education:typeprof - 10.7951020408163 education:typewc + 12.5561779064803 log2(income):typeprof - 12.5561779064803 log2(income):typewc = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 4165.051
2 89 3655.397 1 509.6539 12.40883 0.00067667 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------
Term log2(income)
Adjusted slope for log2(income):
-0.8825204
Std. Error:
3.104149
Linear hypothesis test
Hypothesis:
log2(income):typeprof - log2(income):typewc = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 3658.716
2 89 3655.397 1 3.319776 0.08083 0.77684
------
Type II Repeated Measures MANOVA Tests: Pillai test statistic
Df test stat approx F num Df den Df Pr(>F)
(Intercept) 1 0.96954 318.34 1 10 6.532e-09 ***
treatment 2 0.48092 4.63 2 10 0.0376868 *
gender 1 0.20356 2.56 1 10 0.1409735
treatment:gender 2 0.36350 2.86 2 10 0.1044692
phase 1 0.85052 25.61 2 9 0.0001930 ***
treatment:phase 2 0.68518 2.61 4 20 0.0667354 .
gender:phase 1 0.04314 0.20 2 9 0.8199968
treatment:gender:phase 2 0.31060 0.92 4 20 0.4721498
hour 1 0.93468 25.04 4 7 0.0003043 ***
treatment:hour 2 0.30144 0.35 8 16 0.9295212
gender:hour 1 0.29274 0.72 4 7 0.6023742
treatment:gender:hour 2 0.57022 0.80 8 16 0.6131884
phase:hour 1 0.54958 0.46 8 3 0.8324517
treatment:phase:hour 2 0.66367 0.25 16 8 0.9914415
gender:phase:hour 1 0.69505 0.85 8 3 0.6202076
treatment:gender:phase:hour 2 0.79277 0.33 16 8 0.9723693
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Call: testFactors(model = mod.ok, levels = A.vs.B, terms.formula = ~0 + phase, idata = idata)
Term phase
Adjusted mean at contrasts of phase:
phase1 phase2
-0.8750000 -0.9583333
Std. Error:
phase1 phase2
0.7010657 0.9763701
Response transformation matrix:
phase1 phase2
pre.1 -0.2 -0.2
pre.2 -0.2 -0.2
pre.3 -0.2 -0.2
pre.4 -0.2 -0.2
pre.5 -0.2 -0.2
post.1 0.0 0.2
post.2 0.0 0.2
post.3 0.0 0.2
post.4 0.0 0.2
post.5 0.0 0.2
fup.1 0.2 0.0
fup.2 0.2 0.0
fup.3 0.2 0.0
fup.4 0.2 0.0
fup.5 0.2 0.0
Sum of squares and products for the hypothesis:
phase1 phase2
phase1 1.934211 2.118421
phase2 2.118421 2.320175
Sum of squares and products for error:
phase1 phase2
phase1 12.41667 12.41667
phase2 12.41667 24.08333
Multivariate Tests:
Df test stat approx F num Df den Df Pr(>F)
Pillai 1 0.1359042 0.707756 2 9 0.51824
Wilks 1 0.8640958 0.707756 2 9 0.51824
Hotelling-Lawley 1 0.1572791 0.707756 2 9 0.51824
Roy 1 0.1572791 0.707756 2 9 0.51824
------
Call: testFactors(model = mod.ok, levels = control.vs.AB, terms.formula = ~0 + phase, idata = idata)
Term phase
Adjusted mean at contrasts of phase:
phase1 phase2
2.604167 2.145833
Std. Error:
phase1 phase2
0.6177004 0.8602678
Response transformation matrix:
phase1 phase2
pre.1 -0.2 -0.2
pre.2 -0.2 -0.2
pre.3 -0.2 -0.2
pre.4 -0.2 -0.2
pre.5 -0.2 -0.2
post.1 0.0 0.2
post.2 0.0 0.2
post.3 0.0 0.2
post.4 0.0 0.2
post.5 0.0 0.2
fup.1 0.2 0.0
fup.2 0.2 0.0
fup.3 0.2 0.0
fup.4 0.2 0.0
fup.5 0.2 0.0
Sum of squares and products for the hypothesis:
phase1 phase2
phase1 22.06921 18.18503
phase2 18.18503 14.98446
Sum of squares and products for error:
phase1 phase2
phase1 12.41667 12.41667
phase2 12.41667 24.08333
Multivariate Tests:
Df test stat approx F num Df den Df Pr(>F)
Pillai 1 0.6473884 8.261917 2 9 0.0091798 **
Wilks 1 0.3526116 8.261917 2 9 0.0091798 **
Hotelling-Lawley 1 1.8359816 8.261917 2 9 0.0091798 **
Roy 1 1.8359816 8.261917 2 9 0.0091798 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------
$A.vs.B
Adjusted values for factor combinations in model:
lm(formula = cbind(pre.1, pre.2, pre.3, pre.4, pre.5, post.1, post.2, post.3, post.4, post.5, fup.1, fup.2, fup.3, fup.4, fup.5) ~ treatment * gender, data = OBrienKaiser)
------
Adjusted values
Term phase
Adjusted mean at contrasts of phase:
phase1 phase2
-0.8750000 -0.9583333
Standard error:
phase1 phase2
0.7010657 0.9763701
---
------
Multivariate Test: Pillai test statistic
Df test stat approx F num Df den Df Pr(>F)
phase 1 0.1359 0.70776 2 9 0.5182
------------
$control.vs.AB
Adjusted values for factor combinations in model:
lm(formula = cbind(pre.1, pre.2, pre.3, pre.4, pre.5, post.1, post.2, post.3, post.4, post.5, fup.1, fup.2, fup.3, fup.4, fup.5) ~ treatment * gender, data = OBrienKaiser)
------
Adjusted values
Term phase
Adjusted mean at contrasts of phase:
phase1 phase2
2.604167 2.145833
Standard error:
phase1 phase2
0.6177004 0.8602678
---
------
Multivariate Test: Pillai test statistic
Df test stat approx F num Df den Df Pr(>F)
phase 1 0.64739 8.2619 2 9 0.00918 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------------
phia documentation built on May 2, 2019, 8:23 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Description Usage Arguments Details Value Note Author(s) See Also Examples
Description
Calculates and tests the adjusted mean value of the response and other terms of a fitted model, for specific linear combinations of factor levels and specific values of the covariates. This function is specially indicated for post-hoc analyses of models with factors, to test pairwise comparisons of factor levels, simple main effects, and interaction contrasts. In multivariate models, it is possible to define and test contrasts of intra-subjects factors, as in a repeated-measures analysis.
Usage
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 | ## Default S3 method:
testFactors(model, levels, covariates, offset, terms.formula=~1,
inherit.contrasts=FALSE, default.contrasts=c("contr.sum","contr.poly"),
lht=TRUE, ...)
## S3 method for class 'lm'
testFactors(model, ...)
## S3 method for class 'glm'
testFactors(model, ..., link=FALSE)
## S3 method for class 'mlm'
testFactors(model, levels, covariates, offset, terms.formula=~1,
inherit.contrasts=FALSE, default.contrasts=c("contr.sum","contr.poly"),
idata, icontrasts=default.contrasts, lht=TRUE, ...)
## S3 method for class 'lme'
testFactors(model, ...)
## S3 method for class 'mer'
testFactors(model, ..., link=FALSE)
## S3 method for class 'merMod'
testFactors(model, ..., link=FALSE)
## S3 method for class 'testFactors'
summary(object, predictors=TRUE, matrices=TRUE, covmat=FALSE, ...)
|
Arguments
model |
fitted model. Currently supported classes include |
levels |
list that defines the factor levels that will be contrasted; see the Details for more information. |
covariates |
optional numeric vector with specific values of the covariates (see Details). |
offset |
optional numeric value with a specific value of the offset (if any). |
terms.formula |
formula that defines the terms of the model that will be calculated and tested. The default value is ~1, and stands for the adjusted mean value. See the Details for more information. |
inherit.contrasts |
logical value: should the default contrasts of model factors be inherited from the model data frame? |
default.contrasts |
names of contrast-generating functions to be applied by default to factors and ordered factors, respectively, if |
lht |
logical indicating if the adjusted values are tested (via |
link |
for models fitted with |
idata |
an optional data frame giving a factor or factors defining the intra-subjects model for multivariate repeated-measures data, as defined in |
icontrasts |
names of contrast-generating functions to be applied in the within-subject “data”. The default is the same as |
object |
object returned by |
predictors |
logical value: should |
matrices |
logical value: should |
covmat |
logical value: should |
... |
other arguments passed down to |
Details
The only mandatory argument is model, which may include any number of factor or numeric predictors, and one offset. The simplest usage of this method, where no other argument is defined, calculates the adjusted mean of the model response variable, pooling over all the levels of factor predictors, and setting the numeric predictors (covariates and offset, ifany) to their average values in the model data frame.
The calculations will be done for the linear combinations of factor levels defined by levels. This argument must be a list, with one element for each factor of the model that has to be manipulated (including factors of the intra-subjects design, if suitable). The factors that are not represented in this list will be pooled over, and elements that do not correspond to any factor of the model will be ignored with a warning. levels may be a named list, where the name of each element identifies the represented factor, and its contents may be one of the following:
A character string of length 1 or 2, with the name of one or two factor levels. In the former case, the calculations will be restricted to this level of the factor, as in a simple main effects analysis; in the latter, a pairwise contrast will be calculated between both factor levels.
A numeric vector without names, as long as the number of levels in the factor. This will create a linear combination of factor levels, with the elements of the vector as coefficients. For instance, if the factor
fhas three levels, an elementf=c(0.5, 0.5, 0)will average the two first levels, andf=c(0.5, 0.5, -1)will contrast the average of the two first levels against the third one.A numeric vector with names equal to some or all the levels of the factor. This is a simplification of the previous option, where some levels can be omitted, and the coefficient of each level is determined by the names of the vector, which do not have to follow a specific order. Omitted levels will automatically be set to zero.
A numeric matrix, as an extension of the two previous options for calculating several combinations at a time. Combinations are defined in columns, so if the matrix does not have row names, the number of rows must be equal to the number of levels in the factor, or if the matrix does have row names, they must coincide with the levels of the factor.
Alternatively, levels may be a single formula or an unnamed list of formulas, of the type factorname ~ K1*level1 + K2*level2 ... (see contrastCoefficients for further details). Both types of lists (named list of string or numeric vectors and matrices, and unnamed lists of formulas) may be mixed.
The argument covariates may be used for setting specific values of the model numeric predictors (covariates). It must be a vector of numeric values. If the elements of this vector have names, their values will be assigned to the covariates that match them; covariates of the model with names not represented in this vector will be set to their default value (the average in the model data frame), and elements with names that do not match with covariates will be ignored. On the other hand, if covariates has no names, and its length is equal to the number of covariates of the model, the values will be assigned to those covariates in the same order as they occur in the model. If it has a different length, the vector will be trimmed or reclycled as needed to fit the number of covariates in the model; this feature may be used, for instance, to set all covariates to the same value, e.g. covariates = 0. The argument offset can likewise be used to define a specific value for the offset of the model.
To analyse terms other than the adjusted mean value of the response, use the argument terms.formula. For instance, if the model has the covariates var1, var2, ..., the slopes of the response with respect to them may be added to the analysis by defining terms.formula as ~var1 + var2 .... This formula may be used more generally, for analysing interactions, omitting the mean response, adding the main effects of factors, etc. A different analysis is done for each term of this formula, that must also be contained in the formula of model. For instance, if terms.formula is equal to ~ var1*var2, the function will analyse the adjusted intercept, plus the terms var1, var2, and var1:var2. The intercept stands for the mean value of the response, and terms formed by one or more covariates stand por the slope of the response with respect to the product of those covariates.
If any of the variables in the term is a factor, the function analyses a full set of contrasts for that factor of the remaining part of the term; for instance if var1 were a factor, the term var1 would stand for the contrasts of the intercept, and var1:var2 would stand for the contrasts of the slope var2, across the levels of var1. The set of contrasts used in the analysis is normally defined by the argument default.contrasts: by default, if the factor is ordered it will be a set of “polynomial contrasts”, and otherwise “sum contrasts”; however, if inherit.contrasts is TRUE the contrasts will directly be copied from the ones used to define the model. Factors that have explicit contrasts defined in the model data frame will use those contrasts, regardless of the values defined for default.contrasts and inherit.contrasts. The analysis assumes that the contrasts are orthogonal to the intercept, which is the usual case if the default arguments are provided, and a warning will be issued if non-orthogonal contrasts are used; take special care of not using “treatment contrats” if inherit.contrasts is set to TRUE or default.contrasts is changed.
In generalized linear models, the adjusted means represent the expected values of the response by default, but the expected value of the link function may be shown by setting the argument link=FALSE. On the other hand, slope values and standard errors always refer to the link function.
For multivariate models, the arguments idata, and icontrasts may be used to define an intra-subjects model for multivariate repeated-measures data, as described for Anova or linearHypothesis in package car. Note, however, that the combinations of intra-subjects factor levels are defined in levels, and other arguments defined in those functions like idesign, imatrix or iterms will have no effect in testFactors.
The significance of adjusted values is tested by a call to linearHypothesis for each term, unless lht is set to FALSE. Extra arguments may be passed down to that function, for instance to specify the test statistic that will be evaluated.
Value
An object of class "testFactors", that contains the adjusted values and their standard errors for each term, and the otuput of the test, plus other variables used in the calculations. The summary method for this object will display those variables, unless they be omitted by setting the optional arguments predictors, matrices or covmat to FALSE. The argument predictors refers to the coefficients of specified combinations of factor levels, the values of covariates, and the contrast matrices used for terms that include factors; matrices refers to the “linear hypothesis matrix” used by linearHypothesis, and in multivariate linear models, to the “response transformation matrix” as well — if it exists; covmat refers to the variance-covariance matrix of the adjusted values.
Moreover, summary groups the results of the tests for all terms in one table. By default this table shows the test statistics, their degrees of freedom, and the p-values. If the model is of class "lm", it also shows the sums of squares; and if it is of class "mlm", only the first type of test statistic returned by linearHypothesis (by default “Pillai”) is shown. This variable shape of the ANOVA table is controlled by additional classes assigned to the object (either "testFactors.lm" or "testFactors.mlm", as suitable).
Note
The tests of mixed models are done under the assumption that the estimation of the random part of the model is exact.
Author(s)
Helios De Rosario-Martinez, helios.derosario@gmail.com
See Also
linearHypothesis in package car.
interactionMeans, and testInteractions as useful wrappers of testFactors.
Examples
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 | # Example with factors and covariates
# Analyse prestige of Canadian occupations depending on
# education, income and type of occupation
# See ?Prestige for a description of the data set
prestige.mod <- lm(prestige ~ (education+log2(income))*type, data=Prestige)
# Pairwise comparisons for factor "type", to see how it influences
# the mean value of prestige and interacts with log.income
# 1: "white collar" vs "blue collar"
wc.vs.bc <- list(type=c("wc", "bc"))
testFactors(prestige.mod, wc.vs.bc, terms.formula=~log2(income))
# 2: "professional" vs. "blue collar"
prof.vs.bc <- list(type=c("prof", "bc"))
testFactors(prestige.mod, prof.vs.bc, terms.formula=~log2(income))
# 3: "professional" vs. "white collar"
prof.vs.wc <- list(type=c("prof", "wc"))
testFactors(prestige.mod, prof.vs.wc, terms.formula=~log2(income))
# Interaction contrasts in a repeated-measures experiment
# See ?OBrienKaiser for a description of the data set
mod.ok <- lm(cbind(pre.1, pre.2, pre.3, pre.4, pre.5,
post.1, post.2, post.3, post.4, post.5,
fup.1, fup.2, fup.3, fup.4, fup.5) ~ treatment*gender,
data=OBrienKaiser)
# intra-subjects data:
phase <- factor(rep(c("pretest", "posttest", "followup"), each=5))
hour <- ordered(rep(1:5, 3))
idata <- data.frame(phase, hour)
Anova(mod.ok, idata=idata, idesign=~phase*hour)
# Contrasts across "phase", for different contrasts of "treatment"
# Using different definitions of the argument "levels"
# 1: Treatment "A" vs. treatment "B".
A.vs.B <- list(treatment=c("A", "B"))
# The following are equivalent:
# A.vs.B <- list(treatment=c(A=1, B=-1, control=0))
# A.vs.B <- list(treatment=c(A=1, B=-1))
# A.vs.B <- list(treatment ~ A - B)
# A.vs.B <- treatment ~ A - B
testFactors(mod.ok, A.vs.B, idata=idata, terms.formula=~0+phase)
# 2: Controls vs. treatments
control.vs.AB <- list(treatment=c(A=0.5, B=0.5, control=-1))
# The following is equivalent:
# control.vs.AB <- treatment ~ (A+B)/2 - control
testFactors(mod.ok, control.vs.AB, idata=idata, terms.formula=~0+phase)
# Shortcut to get only the adjusted values and simplified ANOVA tables
contr <- list(A.vs.B=A.vs.B, control.vs.AB=control.vs.AB)
anovaTables <- function(contrast) summary(testFactors(mod.ok, contrast,
idata=idata, terms.formula=~0+phase),
predictors=FALSE, matrices=FALSE)
lapply(contr,anovaTables)
|
Example output
Loading required package: car
Loading required package: carData
Call: testFactors(model = prestige.mod, levels = wc.vs.bc, terms.formula = ~log2(income))
Term (Intercept)
Adjusted mean:
-1.44478
Std. Error:
2.842545
Linear hypothesis test
Hypothesis:
typewc + 10.7951020408163 education:typewc + 12.5561779064803 log2(income):typewc = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 3666.007
2 89 3655.397 1 10.61044 0.25834 0.61252
------
Term log2(income)
Adjusted slope for log2(income):
-5.653038
Std. Error:
3.051886
Linear hypothesis test
Hypothesis:
log2(income):typewc = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 3796.316
2 89 3655.397 1 140.9197 3.43105 0.067296 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------
Call: testFactors(model = prestige.mod, levels = prof.vs.bc, terms.formula = ~log2(income))
Term (Intercept)
Adjusted mean:
10.62698
Std. Error:
3.784129
Linear hypothesis test
Hypothesis:
typeprof + 10.7951020408163 education:typeprof + 12.5561779064803 log2(income):typeprof = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 3979.313
2 89 3655.397 1 323.916 7.88657 0.0061205 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------
Term log2(income)
Adjusted slope for log2(income):
-6.535559
Std. Error:
2.616708
Linear hypothesis test
Hypothesis:
log2(income):typeprof = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 3911.609
2 89 3655.397 1 256.212 6.23814 0.014341 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------
Call: testFactors(model = prestige.mod, levels = prof.vs.wc, terms.formula = ~log2(income))
Term (Intercept)
Adjusted mean:
12.07176
Std. Error:
3.42693
Linear hypothesis test
Hypothesis:
typeprof - typewc + 10.7951020408163 education:typeprof - 10.7951020408163 education:typewc + 12.5561779064803 log2(income):typeprof - 12.5561779064803 log2(income):typewc = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 4165.051
2 89 3655.397 1 509.6539 12.40883 0.00067667 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------
Term log2(income)
Adjusted slope for log2(income):
-0.8825204
Std. Error:
3.104149
Linear hypothesis test
Hypothesis:
log2(income):typeprof - log2(income):typewc = 0
Model 1: restricted model
Model 2: prestige ~ (education + log2(income)) * type
Res.Df RSS Df Sum of Sq F Pr(>F)
1 90 3658.716
2 89 3655.397 1 3.319776 0.08083 0.77684
------
Type II Repeated Measures MANOVA Tests: Pillai test statistic
Df test stat approx F num Df den Df Pr(>F)
(Intercept) 1 0.96954 318.34 1 10 6.532e-09 ***
treatment 2 0.48092 4.63 2 10 0.0376868 *
gender 1 0.20356 2.56 1 10 0.1409735
treatment:gender 2 0.36350 2.86 2 10 0.1044692
phase 1 0.85052 25.61 2 9 0.0001930 ***
treatment:phase 2 0.68518 2.61 4 20 0.0667354 .
gender:phase 1 0.04314 0.20 2 9 0.8199968
treatment:gender:phase 2 0.31060 0.92 4 20 0.4721498
hour 1 0.93468 25.04 4 7 0.0003043 ***
treatment:hour 2 0.30144 0.35 8 16 0.9295212
gender:hour 1 0.29274 0.72 4 7 0.6023742
treatment:gender:hour 2 0.57022 0.80 8 16 0.6131884
phase:hour 1 0.54958 0.46 8 3 0.8324517
treatment:phase:hour 2 0.66367 0.25 16 8 0.9914415
gender:phase:hour 1 0.69505 0.85 8 3 0.6202076
treatment:gender:phase:hour 2 0.79277 0.33 16 8 0.9723693
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Call: testFactors(model = mod.ok, levels = A.vs.B, terms.formula = ~0 + phase, idata = idata)
Term phase
Adjusted mean at contrasts of phase:
phase1 phase2
-0.8750000 -0.9583333
Std. Error:
phase1 phase2
0.7010657 0.9763701
Response transformation matrix:
phase1 phase2
pre.1 -0.2 -0.2
pre.2 -0.2 -0.2
pre.3 -0.2 -0.2
pre.4 -0.2 -0.2
pre.5 -0.2 -0.2
post.1 0.0 0.2
post.2 0.0 0.2
post.3 0.0 0.2
post.4 0.0 0.2
post.5 0.0 0.2
fup.1 0.2 0.0
fup.2 0.2 0.0
fup.3 0.2 0.0
fup.4 0.2 0.0
fup.5 0.2 0.0
Sum of squares and products for the hypothesis:
phase1 phase2
phase1 1.934211 2.118421
phase2 2.118421 2.320175
Sum of squares and products for error:
phase1 phase2
phase1 12.41667 12.41667
phase2 12.41667 24.08333
Multivariate Tests:
Df test stat approx F num Df den Df Pr(>F)
Pillai 1 0.1359042 0.707756 2 9 0.51824
Wilks 1 0.8640958 0.707756 2 9 0.51824
Hotelling-Lawley 1 0.1572791 0.707756 2 9 0.51824
Roy 1 0.1572791 0.707756 2 9 0.51824
------
Call: testFactors(model = mod.ok, levels = control.vs.AB, terms.formula = ~0 + phase, idata = idata)
Term phase
Adjusted mean at contrasts of phase:
phase1 phase2
2.604167 2.145833
Std. Error:
phase1 phase2
0.6177004 0.8602678
Response transformation matrix:
phase1 phase2
pre.1 -0.2 -0.2
pre.2 -0.2 -0.2
pre.3 -0.2 -0.2
pre.4 -0.2 -0.2
pre.5 -0.2 -0.2
post.1 0.0 0.2
post.2 0.0 0.2
post.3 0.0 0.2
post.4 0.0 0.2
post.5 0.0 0.2
fup.1 0.2 0.0
fup.2 0.2 0.0
fup.3 0.2 0.0
fup.4 0.2 0.0
fup.5 0.2 0.0
Sum of squares and products for the hypothesis:
phase1 phase2
phase1 22.06921 18.18503
phase2 18.18503 14.98446
Sum of squares and products for error:
phase1 phase2
phase1 12.41667 12.41667
phase2 12.41667 24.08333
Multivariate Tests:
Df test stat approx F num Df den Df Pr(>F)
Pillai 1 0.6473884 8.261917 2 9 0.0091798 **
Wilks 1 0.3526116 8.261917 2 9 0.0091798 **
Hotelling-Lawley 1 1.8359816 8.261917 2 9 0.0091798 **
Roy 1 1.8359816 8.261917 2 9 0.0091798 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------
$A.vs.B
Adjusted values for factor combinations in model:
lm(formula = cbind(pre.1, pre.2, pre.3, pre.4, pre.5, post.1, post.2, post.3, post.4, post.5, fup.1, fup.2, fup.3, fup.4, fup.5) ~ treatment * gender, data = OBrienKaiser)
------
Adjusted values
Term phase
Adjusted mean at contrasts of phase:
phase1 phase2
-0.8750000 -0.9583333
Standard error:
phase1 phase2
0.7010657 0.9763701
---
------
Multivariate Test: Pillai test statistic
Df test stat approx F num Df den Df Pr(>F)
phase 1 0.1359 0.70776 2 9 0.5182
------------
$control.vs.AB
Adjusted values for factor combinations in model:
lm(formula = cbind(pre.1, pre.2, pre.3, pre.4, pre.5, post.1, post.2, post.3, post.4, post.5, fup.1, fup.2, fup.3, fup.4, fup.5) ~ treatment * gender, data = OBrienKaiser)
------
Adjusted values
Term phase
Adjusted mean at contrasts of phase:
phase1 phase2
2.604167 2.145833
Standard error:
phase1 phase2
0.6177004 0.8602678
---
------
Multivariate Test: Pillai test statistic
Df test stat approx F num Df den Df Pr(>F)
phase 1 0.64739 8.2619 2 9 0.00918 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
------------
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.