This document was created by Jihong Zhang and Jonathan Templin, based on an example developed by Lesa Hoffman.
The models for this example come from Hoffman (2015) chapter 2. We will be examining the extent to which cognition (as measured by an information test outcome) can be predicted from age (centered at 85 years) grip strength (centered at 9 pounds), sex (with men as the reference group) and dementia status (none, future, current, with none as the reference) in a sample of 550 older adults. In example 2a we examined main effects only; in example 2b we examined interactions involving age, grip strength, and sex (a binary predictor that was treated as continuous). Now we examine the interaction of sex by dementia group when treated as “continuous” via manual contrasts or as “categorical” via R’s factor type.
First, renouncing ANOVA dogma: You can ask for Type I error correction to pairwise follow-ups if desired, although I don’t believe they should have any special status relative to any other default-provided p-values because all group differences and their SEs are already predicted by the model.
In addition, it is possible to have significant pairwise group comparisons even if the “omnibus” test is not significant because it is evaluated per df. So if one group is very different than all the others, then the average effect per df may not be different than 0 (but you’d be missing the real story by only examining the “omnibus” result). So I do not think the omnibus F-test should be the gate-keeper to examining group comparisons. But just remember, for every 20 p-values, one is significant by chance!
#Step 1: Find path of data file and copy path (if using windows)
#Step 2: In the console below, type readClipboard()
#Step 3: Copy and paste R's path to the line below in quotes
#CHANGE BETWEEN QUOTES IN THIS LINE TO REFLECT DIRECTORY OF DATA:
#IMPORT DATA AND PUT INTO DATASET
data01 = read.csv("mv18epsy905_lecture03.csv", header = TRUE)
#AUTOMATING PACKAGES NEEDED FOR ANALYSES
haspackage = require("multcomp")Loading required package: multcomp
Loading required package: mvtnorm
Loading required package: survival
Loading required package: TH.data
Loading required package: MASS
Attaching package: ‘TH.data’
The following object is masked from ‘package:MASS’:
geyserif (haspackage==FALSE){
install.packages("multcomp")
}
library("multcomp")
### GLANCE AT DATA
summary(data01) PersonID cognition age grip sexMW demgroup
Min. : 1.0 Min. : 0.00 Min. :80.02 Min. : 0.000 Min. :0.0000 Min. :1.000
1st Qu.:138.2 1st Qu.:17.25 1st Qu.:82.33 1st Qu.: 7.000 1st Qu.:0.0000 1st Qu.:1.000
Median :275.5 Median :25.00 Median :84.33 Median : 9.000 Median :1.0000 Median :1.000
Mean :275.5 Mean :24.82 Mean :84.93 Mean : 9.113 Mean :0.5873 Mean :1.351
3rd Qu.:412.8 3rd Qu.:33.00 3rd Qu.:86.19 3rd Qu.:11.000 3rd Qu.:1.0000 3rd Qu.:2.000
Max. :550.0 Max. :44.00 Max. :96.97 Max. :19.000 Max. :1.0000 Max. :3.000 ### CENTERING CONTINOUSE PREDICTORS;
data01$age85 = data01$age - 85
data01$grip9 = data01$grip - 9
### CHECK AGAIN, age85/grip9 have been added
summary(data01) PersonID cognition age grip sexMW demgroup age85 grip9
Min. : 1.0 Min. : 0.00 Min. :80.02 Min. : 0.000 Min. :0.0000 Min. :1.000 Min. :-4.98351 Min. :-9.0000
1st Qu.:138.2 1st Qu.:17.25 1st Qu.:82.33 1st Qu.: 7.000 1st Qu.:0.0000 1st Qu.:1.000 1st Qu.:-2.67264 1st Qu.:-2.0000
Median :275.5 Median :25.00 Median :84.33 Median : 9.000 Median :1.0000 Median :1.000 Median :-0.67323 Median : 0.0000
Mean :275.5 Mean :24.82 Mean :84.93 Mean : 9.113 Mean :0.5873 Mean :1.351 Mean :-0.07321 Mean : 0.1127
3rd Qu.:412.8 3rd Qu.:33.00 3rd Qu.:86.19 3rd Qu.:11.000 3rd Qu.:1.0000 3rd Qu.:2.000 3rd Qu.: 1.18658 3rd Qu.: 2.0000
Max. :550.0 Max. :44.00 Max. :96.97 Max. :19.000 Max. :1.0000 Max. :3.000 Max. :11.96728 Max. :10.0000 #Create FACTOR variables for demgroup and Sex
data01$demgroupF = as.factor(data01$demgroup)
data01$sexMWF = as.factor(data01$sexMW)Here, all tests from anova() are df=1 because all fixed effects are single df (treated as separate).
#MODEL #1 -- Using 0/1 coding instead of factors
model1 = lm(cognition~ age85 +grip9 + age85*grip9 + sexMWF + demgroupF + sexMWF*demgroupF, data=data01)
model11 <- update(model1, .~.-sexMWF:demgroupF, data= data01)
summary(model1)
Call:
lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMWF +
demgroupF + sexMWF * demgroupF, data = data01)
Residuals:
Min 1Q Median 3Q Max
-27.9650 -5.7937 -0.1291 6.6792 19.7428
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 29.07015 0.74850 38.838 < 2e-16 ***
age85 -0.33480 0.11989 -2.793 0.005414 **
grip9 0.61789 0.14808 4.173 3.51e-05 ***
sexMWF1 -2.87559 1.01124 -2.844 0.004629 **
demgroupF2 -6.05590 1.63513 -3.704 0.000234 ***
demgroupF3 -11.97073 2.24495 -5.332 1.43e-07 ***
age85:grip9 0.12215 0.04035 3.027 0.002587 **
sexMWF1:demgroupF2 0.16427 2.07048 0.079 0.936792
sexMWF1:demgroupF3 -7.87510 3.02454 -2.604 0.009475 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 9.272 on 541 degrees of freedom
Multiple R-squared: 0.2984, Adjusted R-squared: 0.2881
F-statistic: 28.77 on 8 and 541 DF, p-value: < 2.2e-16summary(model11)
Call:
lm(formula = cognition ~ age85 + grip9 + sexMWF + demgroupF +
age85:grip9, data = data01)
Residuals:
Min 1Q Median 3Q Max
-27.7181 -6.0173 -0.1319 6.9618 19.9896
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 29.40780 0.69491 42.319 < 2e-16 ***
age85 -0.33396 0.12036 -2.775 0.005715 **
grip9 0.61942 0.14874 4.164 3.63e-05 ***
sexMWF1 -3.45564 0.88727 -3.895 0.000111 ***
demgroupF2 -5.92254 1.01363 -5.843 8.85e-09 ***
demgroupF3 -16.30040 1.51255 -10.777 < 2e-16 ***
age85:grip9 0.12302 0.04054 3.035 0.002522 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 9.315 on 543 degrees of freedom
Multiple R-squared: 0.2894, Adjusted R-squared: 0.2815
F-statistic: 36.85 on 6 and 543 DF, p-value: < 2.2e-16anova(model1)Analysis of Variance Table
Response: cognition
Df Sum Sq Mean Sq F value Pr(>F)
age85 1 1926 1926.2 22.4046 2.823e-06 ***
grip9 1 3039 3039.2 35.3506 4.942e-09 ***
sexMWF 1 1609 1609.3 18.7191 1.804e-05 ***
demgroupF 2 11811 5905.7 68.6924 < 2.2e-16 ***
age85:grip9 1 799 799.1 9.2944 0.002411 **
sexMWF:demgroupF 2 600 300.2 3.4919 0.031130 *
Residuals 541 46511 86.0
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1anova(model1, model11)Analysis of Variance Table
Model 1: cognition ~ age85 + grip9 + age85 * grip9 + sexMWF + demgroupF +
sexMWF * demgroupF
Model 2: cognition ~ age85 + grip9 + sexMWF + demgroupF + age85:grip9
Res.Df RSS Df Sum of Sq F Pr(>F)
1 541 46511
2 543 47111 -2 -600.42 3.4919 0.03113 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1# Contrasts
# Conditional Main Effects of Age85
conAgeeffect = matrix(c(0,1,0,0,0,0,0,0,0),1)
condeffects1= glht(model1,linfct=conAgeeffect,test=Ftest())
print("Conditional Effect for Age85")[1] "Conditional Effect for Age85"summary(condeffects1,test=Ftest())
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -0.3348
Global Test:
F DF1 DF2 Pr(>F)
1 7.799 1 541 0.005414# R function used to calculate Conditional Main effects that match SPSS and SAS for "Regression" style output
condEffects <- function(model = model1) {
vars <- names(coef(model)) # name vector for predictors
k <- 2:length(vars) # place rank for each predictor
diagmatrix <- diag(length(vars)) # make a diag matrix
for (i in k) {
print(paste("conditional effect of", vars[i]))
effectmatrix <- matrix(diagmatrix[i,],1)
condeffects= glht(model,linfct=effectmatrix,test=Ftest())
print(summary(condeffects,test=Ftest()))
cat("----------------------------------\n")
}
}
condEffects(model = model1)[1] "conditional effect of age85"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -0.3348
Global Test:
F DF1 DF2 Pr(>F)
1 7.799 1 541 0.005414
----------------------------------
[1] "conditional effect of grip9"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 0.6179
Global Test:
F DF1 DF2 Pr(>F)
1 17.41 1 541 3.506e-05
----------------------------------
[1] "conditional effect of sexMWF1"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -2.876
Global Test:
F DF1 DF2 Pr(>F)
1 8.086 1 541 0.004629
----------------------------------
[1] "conditional effect of demgroupF2"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -6.056
Global Test:
F DF1 DF2 Pr(>F)
1 13.72 1 541 0.0002344
----------------------------------
[1] "conditional effect of demgroupF3"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -11.97
Global Test:
F DF1 DF2 Pr(>F)
1 28.43 1 541 1.427e-07
----------------------------------
[1] "conditional effect of age85:grip9"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 0.1222
Global Test:
F DF1 DF2 Pr(>F)
1 9.163 1 541 0.002587
----------------------------------
[1] "conditional effect of sexMWF1:demgroupF2"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 0.1643
Global Test:
F DF1 DF2 Pr(>F)
1 0.006295 1 541 0.9368
----------------------------------
[1] "conditional effect of sexMWF1:demgroupF3"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -7.875
Global Test:
F DF1 DF2 Pr(>F)
1 6.779 1 541 0.009475
----------------------------------# Cell means (predicted intercepts)
mean1 = matrix(c(1,0,0,0,0,0,0,0,0),1); rownames(mean1)="Men-None"
mean2 = matrix(c(1,0,0,1,0,0,0,0,0),1); rownames(mean2)="Women-None"
mean3 = matrix(c(1,0,0,0,1,0,0,0,0),1); rownames(mean3)="Men-Future"
mean4 = matrix(c(1,0,0,1,1,0,0,1,0),1); rownames(mean4)="Women-Future"
mean5 = matrix(c(1,0,0,0,0,1,0,0,0),1); rownames(mean5)="Men-Current"
mean6 = matrix(c(1,0,0,1,0,1,0,0,1),1); rownames(mean6)="Women-Current"
meansvec = rbind(mean1,mean2,mean3,mean4,mean5,mean6)
means = glht(model1,linfct=meansvec)
summary(means)
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMWF +
demgroupF + sexMWF * demgroupF, data = data01)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
Men-None == 0 29.0701 0.7485 38.838 <1e-05 ***
Women-None == 0 26.1946 0.6388 41.004 <1e-05 ***
Men-Future == 0 23.0142 1.4928 15.417 <1e-05 ***
Women-Future == 0 20.3029 1.1186 18.150 <1e-05 ***
Men-Current == 0 17.0994 2.1402 7.990 <1e-05 ***
Women-Current == 0 6.3487 1.9479 3.259 0.0071 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)# Simple effects of sex and dementia groups;
## Simple effect of sex
effect1 = matrix(c(0,0,0,1,0,0,0,0,0),1); rownames(effect1) = "Sex Diff for None"
effect2 = matrix(c(0,0,0,1,0,0,0,1,0),1); rownames(effect2) = "Sex Diff for Future"
effect3 = matrix(c(0,0,0,1,0,0,0,0,1),1); rownames(effect3) = "Sex Diff for Current"
## Simple effect of dementia group
effect4 = matrix(c(0,0,0,0,1,0,0,0,0),1); rownames(effect4) = "None-Fut Diff in Men"
effect5 = matrix(c(0,0,0,0,1,0,0,1,0),1); rownames(effect5) = "None-Fut Diff in Women"
effect6 = matrix(c(0,0,0,0,0,1,0,0,0),1); rownames(effect6) = "None-Cur Diff in Men"
effect7 = matrix(c(0,0,0,0,0,1,0,0,1),1); rownames(effect7) = "None-Cur Diff in Women"
effect8 = matrix(c(0,0,0,0,-1,1,0,0,0),1); rownames(effect8) = "Fut-Current Diff in Men"
effect9 = matrix(c(0,0,0,0,-1,1,0,-1,1),1); rownames(effect9) = "Fut-Current Diff in Women"
effectsvec1 = rbind(effect1,effect2,effect3,effect4,effect5,effect6,effect7,effect8,effect9)
effects1 = glht(model1,linfct=effectsvec1)
summary(effects1)
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMWF +
demgroupF + sexMWF * demgroupF, data = data01)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
Sex Diff for None == 0 -2.876 1.011 -2.844 0.03406 *
Sex Diff for Future == 0 -2.711 1.874 -1.447 0.61742
Sex Diff for Current == 0 -10.751 2.899 -3.708 0.00183 **
None-Fut Diff in Men == 0 -6.056 1.635 -3.704 0.00200 **
None-Fut Diff in Women == 0 -5.892 1.278 -4.611 < 0.001 ***
None-Cur Diff in Men == 0 -11.971 2.245 -5.332 < 0.001 ***
None-Cur Diff in Women == 0 -19.846 2.029 -9.783 < 0.001 ***
Fut-Current Diff in Men == 0 -5.915 2.587 -2.287 0.14540
Fut-Current Diff in Women == 0 -13.954 2.239 -6.233 < 0.001 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)# Diffs in simple effects = interactions;
effect10 = matrix(c(0,0,0,0,0,0,0,1,0),1); rownames(effect10) = "A: Sex Effect differ between None and Future?"
effect11 = matrix(c(0,0,0,0,0,0,0,1,0),1); rownames(effect11) = "A: None-Future Effect differ by Sex?"
effect12 = matrix(c(0,0,0,0,0,0,0,0,1),1); rownames(effect12) = "B: Sex Effect differ between None and Current?"
effect13 = matrix(c(0,0,0,0,0,0,0,0,1),1); rownames(effect13) = "B: None-Current Effect differ by Sex?"
effect14 = matrix(c(0,0,0,0,0,0,0,-1,1),1); rownames(effect14) = "C: Sex Effect differ between Future and Current?"
effect15 = matrix(c(0,0,0,0,0,0,0,-1,1),1); rownames(effect15) = "C: Future-Current Effect differ by Sex?"
effectsvec2 = rbind(effect10,effect11,effect12,effect13,effect14,effect15)
effects2 = glht(model1,linfct=effectsvec2)
summary(effects2)
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMWF +
demgroupF + sexMWF * demgroupF, data = data01)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
A: Sex Effect differ between None and Future? == 0 0.1643 2.0705 0.079 0.9964
A: None-Future Effect differ by Sex? == 0 0.1643 2.0705 0.079 0.9964
B: Sex Effect differ between None and Current? == 0 -7.8751 3.0245 -2.604 0.0242 *
B: None-Current Effect differ by Sex? == 0 -7.8751 3.0245 -2.604 0.0242 *
C: Sex Effect differ between Future and Current? == 0 -8.0394 3.4152 -2.354 0.0470 *
C: Future-Current Effect differ by Sex? == 0 -8.0394 3.4152 -2.354 0.0470 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)R-created dummy code for ref group is omitted – and the output only shows the effects included in the model. You have to know which one is left out to understand the reference group.
#Continuous Sex (0=Men) and Categorical Dementia (Ref=Current)
data01$demgroupF <- relevel(data01$demgroupF, ref = "3")
model2 = lm(cognition~age85 +grip9 + age85*grip9 + sexMW + demgroupF +
sexMW*demgroupF, data=data01)
summary(model2)
Call:
lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMW +
demgroupF + sexMW * demgroupF, data = data01)
Residuals:
Min 1Q Median 3Q Max
-27.9650 -5.7937 -0.1291 6.6792 19.7428
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 17.09942 2.14022 7.990 8.21e-15 ***
age85 -0.33480 0.11989 -2.793 0.00541 **
grip9 0.61789 0.14808 4.173 3.51e-05 ***
sexMW -10.75069 2.89932 -3.708 0.00023 ***
demgroupF1 11.97073 2.24495 5.332 1.43e-07 ***
demgroupF2 5.91483 2.58676 2.287 0.02261 *
age85:grip9 0.12215 0.04035 3.027 0.00259 **
sexMW:demgroupF1 7.87510 3.02454 2.604 0.00947 **
sexMW:demgroupF2 8.03937 3.41516 2.354 0.01893 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 9.272 on 541 degrees of freedom
Multiple R-squared: 0.2984, Adjusted R-squared: 0.2881
F-statistic: 28.77 on 8 and 541 DF, p-value: < 2.2e-16anova(model2)Analysis of Variance Table
Response: cognition
Df Sum Sq Mean Sq F value Pr(>F)
age85 1 1926 1926.2 22.4046 2.823e-06 ***
grip9 1 3039 3039.2 35.3506 4.942e-09 ***
sexMW 1 1609 1609.3 18.7191 1.804e-05 ***
demgroupF 2 11811 5905.7 68.6924 < 2.2e-16 ***
age85:grip9 1 799 799.1 9.2944 0.002411 **
sexMW:demgroupF 2 600 300.2 3.4919 0.031130 *
Residuals 541 46511 86.0
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1levels(data01$demgroupF) # Ref=Current(3), demgroupF1 = None, demgroupF2=Future[1] "3" "1" "2"# F tests
condEffects(model = model2) # Marginal Effect of sexMW and demgroup are not correct[1] "conditional effect of age85"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -0.3348
Global Test:
F DF1 DF2 Pr(>F)
1 7.799 1 541 0.005414
----------------------------------
[1] "conditional effect of grip9"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 0.6179
Global Test:
F DF1 DF2 Pr(>F)
1 17.41 1 541 3.506e-05
----------------------------------
[1] "conditional effect of sexMW"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -10.75
Global Test:
F DF1 DF2 Pr(>F)
1 13.75 1 541 0.0002304
----------------------------------
[1] "conditional effect of demgroupF1"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 11.97
Global Test:
F DF1 DF2 Pr(>F)
1 28.43 1 541 1.427e-07
----------------------------------
[1] "conditional effect of demgroupF2"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 5.915
Global Test:
F DF1 DF2 Pr(>F)
1 5.228 1 541 0.02261
----------------------------------
[1] "conditional effect of age85:grip9"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 0.1222
Global Test:
F DF1 DF2 Pr(>F)
1 9.163 1 541 0.002587
----------------------------------
[1] "conditional effect of sexMW:demgroupF1"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 7.875
Global Test:
F DF1 DF2 Pr(>F)
1 6.779 1 541 0.009475
----------------------------------
[1] "conditional effect of sexMW:demgroupF2"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 8.039
Global Test:
F DF1 DF2 Pr(>F)
1 5.541 1 541 0.01893
----------------------------------# Marginal dementia group mean differences averaged across sex means (what F-tests give);
model2_marGroupmean1 <- matrix(c(1,0,0,0.5,1,0,0,0.5,0),1);
rownames(model2_marGroupmean1) = "Marginal group mean for None"
model2_marGroupmean2 <- matrix(c(1,0,0,0.5,0,1,0,0,0.5),1);
rownames(model2_marGroupmean2) = "Marginal group mean for Future"
model2_marGroupmean3 <- matrix(c(1,0,0,0.5,0,0,0,0,0),1);
rownames(model2_marGroupmean3) = "Marginal group mean for Current"
model2_effectsvec0 = rbind(model2_marGroupmean1, model2_marGroupmean2,
model2_marGroupmean3)
model2_groupmean = glht(model2,linfct=model2_effectsvec0)
summary(model2_groupmean)
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMW +
demgroupF + sexMW * demgroupF, data = data01)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
Marginal group mean for None == 0 27.6323 0.4780 57.802 <1e-10 ***
Marginal group mean for Future == 0 21.6586 0.9283 23.331 <1e-10 ***
Marginal group mean for Current == 0 11.7241 1.4443 8.118 <1e-10 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)
# Simple effects of sex per dementia group;
levels(data01$demgroupF)[1] "3" "1" "2"model2_effectSexDiff1 = matrix(c(0,0,0,1,0,0,0,1,0),1); rownames(model2_effectSexDiff1) = "Sex Diff for None"
model2_effectSexDiff2 = matrix(c(0,0,0,1,0,0,0,0,1),1); rownames(model2_effectSexDiff2) = "Sex Diff for Future"
model2_effectSexDiff3 = matrix(c(0,0,0,1,0,0,0,0,0),1); rownames(model2_effectSexDiff3) = "Sex Diff for Current"
model2_effectsvec1 = rbind(model2_effectSexDiff1, model2_effectSexDiff2,
model2_effectSexDiff3)
model2_effects = glht(model2,linfct=model2_effectsvec1)
summary(model2_effects)
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMW +
demgroupF + sexMW * demgroupF, data = data01)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
Sex Diff for None == 0 -2.876 1.011 -2.844 0.013814 *
Sex Diff for Future == 0 -2.711 1.874 -1.447 0.381839
Sex Diff for Current == 0 -10.751 2.899 -3.708 0.000691 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)# Simple effects of dementia group per sex;
#summary(model2)
model2_effectDemDiff1 = matrix(c(0,0,0,0,-1,1,0,0,0),1); rownames(model2_effectDemDiff1) = "None-Fut Diff in Men"
model2_effectDemDiff2 = matrix(c(0,0,0,0,-1,1,0,-1,1),1); rownames(model2_effectDemDiff2) = "None-Fut Diff in Women"
model2_effectDemDiff3 = matrix(c(0,0,0,0,-1,0,0,0,0),1); rownames(model2_effectDemDiff3) = "None-Cur Diff in Men"
model2_effectDemDiff4 = matrix(c(0,0,0,0,-1,0,0,-1,0),1); rownames(model2_effectDemDiff4) = "None-Cur Diff in Women"
model2_effectDemDiff5 = matrix(c(0,0,0,0,0,-1,0,0,0),1); rownames(model2_effectDemDiff5) = "Fut-Current Diff in Men"
model2_effectDemDiff6 = matrix(c(0,0,0,0,0,-1,0,0,-1),1); rownames(model2_effectDemDiff6) = "Fut-Current Diff in Women"
model2_effectsvec2 = rbind(model2_effectDemDiff1, model2_effectDemDiff2,
model2_effectDemDiff3, model2_effectDemDiff4,
model2_effectDemDiff5, model2_effectDemDiff6)
model2_effects2 = glht(model2,linfct=model2_effectsvec2)
summary(model2_effects2)
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMW +
demgroupF + sexMW * demgroupF, data = data01)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
None-Fut Diff in Men == 0 -6.056 1.635 -3.704 0.00131 **
None-Fut Diff in Women == 0 -5.892 1.278 -4.611 < 1e-04 ***
None-Cur Diff in Men == 0 -11.971 2.245 -5.332 < 1e-04 ***
None-Cur Diff in Women == 0 -19.846 2.029 -9.783 < 1e-04 ***
Fut-Current Diff in Men == 0 -5.915 2.587 -2.287 0.10788
Fut-Current Diff in Women == 0 -13.954 2.239 -6.233 < 1e-04 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)# Differences in simple effects = interactions;
model2_effectIntDiff1 = matrix(c(0,0,0,0,0,0,0,-1,1),1); rownames(model2_effectIntDiff1) = "A: Sex Effect differ between None and Future?"
model2_effectIntDiff2 = matrix(c(0,0,0,0,0,0,0,-1,1),1); rownames(model2_effectIntDiff2) = "A: None-Future Effect differ by Sex?"
model2_effectIntDiff3 = matrix(c(0,0,0,0,0,0,0,-1,0),1); rownames(model2_effectIntDiff3) = "B: Sex Effect differ between None and Current?"
model2_effectIntDiff4 = matrix(c(0,0,0,0,0,0,0,-1,0),1); rownames(model2_effectIntDiff4) = "B: None-Current Effect differ by Sex?"
model2_effectIntDiff5 = matrix(c(0,0,0,0,0,0,0,0,-1),1); rownames(model2_effectIntDiff5) = "C: Sex Effect differ between Future and Current?"
model2_effectIntDiff6 = matrix(c(0,0,0,0,0,0,0,0,-1),1); rownames(model2_effectIntDiff6) = "C: Future-Current Effect differ by Sex?"
model2_effectsvec3 = rbind(model2_effectIntDiff1, model2_effectIntDiff2,
model2_effectIntDiff3, model2_effectIntDiff4,
model2_effectIntDiff5, model2_effectIntDiff6)
model2_effects3 = glht(model2,linfct=model2_effectsvec3)
summary(model2_effects3)
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMW +
demgroupF + sexMW * demgroupF, data = data01)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
A: Sex Effect differ between None and Future? == 0 0.1643 2.0705 0.079 0.9964
A: None-Future Effect differ by Sex? == 0 0.1643 2.0705 0.079 0.9964
B: Sex Effect differ between None and Current? == 0 -7.8751 3.0245 -2.604 0.0242 *
B: None-Current Effect differ by Sex? == 0 -7.8751 3.0245 -2.604 0.0242 *
C: Sex Effect differ between Future and Current? == 0 -8.0394 3.4152 -2.354 0.0470 *
C: Future-Current Effect differ by Sex? == 0 -8.0394 3.4152 -2.354 0.0470 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)In the ANOVA table, Demgroup effects are df=2 now – a simultaneous test of the marginal main effect of Demgroup.
# Categorical Sex (Ref=Women) and Categorical Dementia (Ref=Current)
data01$demgroupF <- relevel(data01$demgroupF, ref = "3")
data01$sexMWF <- relevel(data01$sexMWF, ref = "1")
model3 = lm(cognition~age85 +grip9 + age85*grip9 + sexMWF + demgroupF +
sexMWF*demgroupF, data=data01)
levels(data01$sexMWF)[1] "1" "0"summary(model3)
Call:
lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMWF +
demgroupF + sexMWF * demgroupF, data = data01)
Residuals:
Min 1Q Median 3Q Max
-27.9650 -5.7937 -0.1291 6.6792 19.7428
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.34872 1.94788 3.259 0.00119 **
age85 -0.33480 0.11989 -2.793 0.00541 **
grip9 0.61789 0.14808 4.173 3.51e-05 ***
sexMWF0 10.75069 2.89932 3.708 0.00023 ***
demgroupF1 19.84583 2.02858 9.783 < 2e-16 ***
demgroupF2 13.95420 2.23892 6.233 9.24e-10 ***
age85:grip9 0.12215 0.04035 3.027 0.00259 **
sexMWF0:demgroupF1 -7.87510 3.02454 -2.604 0.00947 **
sexMWF0:demgroupF2 -8.03937 3.41516 -2.354 0.01893 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 9.272 on 541 degrees of freedom
Multiple R-squared: 0.2984, Adjusted R-squared: 0.2881
F-statistic: 28.77 on 8 and 541 DF, p-value: < 2.2e-16anova(model3)Analysis of Variance Table
Response: cognition
Df Sum Sq Mean Sq F value Pr(>F)
age85 1 1926 1926.2 22.4046 2.823e-06 ***
grip9 1 3039 3039.2 35.3506 4.942e-09 ***
sexMWF 1 1609 1609.3 18.7191 1.804e-05 ***
demgroupF 2 11811 5905.7 68.6924 < 2.2e-16 ***
age85:grip9 1 799 799.1 9.2944 0.002411 **
sexMWF:demgroupF 2 600 300.2 3.4919 0.031130 *
Residuals 541 46511 86.0
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1# F-tests
condEffects(model = model3)[1] "conditional effect of age85"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -0.3348
Global Test:
F DF1 DF2 Pr(>F)
1 7.799 1 541 0.005414
----------------------------------
[1] "conditional effect of grip9"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 0.6179
Global Test:
F DF1 DF2 Pr(>F)
1 17.41 1 541 3.506e-05
----------------------------------
[1] "conditional effect of sexMWF0"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 10.75
Global Test:
F DF1 DF2 Pr(>F)
1 13.75 1 541 0.0002304
----------------------------------
[1] "conditional effect of demgroupF1"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 19.85
Global Test:
F DF1 DF2 Pr(>F)
1 95.71 1 541 6.398e-21
----------------------------------
[1] "conditional effect of demgroupF2"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 13.95
Global Test:
F DF1 DF2 Pr(>F)
1 38.84 1 541 9.239e-10
----------------------------------
[1] "conditional effect of age85:grip9"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 0.1222
Global Test:
F DF1 DF2 Pr(>F)
1 9.163 1 541 0.002587
----------------------------------
[1] "conditional effect of sexMWF0:demgroupF1"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -7.875
Global Test:
F DF1 DF2 Pr(>F)
1 6.779 1 541 0.009475
----------------------------------
[1] "conditional effect of sexMWF0:demgroupF2"
General Linear Hypotheses
Linear Hypotheses:
Estimate
1 == 0 -8.039
Global Test:
F DF1 DF2 Pr(>F)
1 5.541 1 541 0.01893
----------------------------------# Simple effects of sex per dementia group;
levels(data01$sexMWF) # Ref=Women(1); sexMWF0=Men[1] "1" "0"levels(data01$demgroupF) # Ref=Current(3); demgroupF1=None, demgroupF2=Future[1] "3" "1" "2"model3_effectSexDiff1 = matrix(c(0,0,0,-1,0,0,0,-1,0),1); rownames(model3_effectSexDiff1) = "Sex Diff for None"
model3_effectSexDiff2 = matrix(c(0,0,0,-1,0,0,0,0,-1),1); rownames(model3_effectSexDiff2) = "Sex Diff for Future"
model3_effectSexDiff3 = matrix(c(0,0,0,-1,0,0,0,0,0),1); rownames(model3_effectSexDiff3) = "Sex Diff for Current"
model3_effectsvec1 = rbind(model3_effectSexDiff1, model3_effectSexDiff2,
model3_effectSexDiff3)
model3_effects = glht(model3,linfct=model3_effectsvec1)
summary(model3_effects)
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMWF +
demgroupF + sexMWF * demgroupF, data = data01)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
Sex Diff for None == 0 -2.876 1.011 -2.844 0.013817 *
Sex Diff for Future == 0 -2.711 1.874 -1.447 0.381836
Sex Diff for Current == 0 -10.751 2.899 -3.708 0.000691 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)# Simple effects of dementia group per sex;
#summary(model3)
levels(data01$sexMWF) # Ref=Women(1); sexMWF0=Men[1] "1" "0"levels(data01$demgroupF) # Ref=Current(3); demgroupF1=None, demgroupF2=Future[1] "3" "1" "2"model3_effectDemDiff1 = matrix(c(0,0,0,0,-1,1,0,-1,1),1); rownames(model3_effectDemDiff1) = "None-Fut Diff in Men"
model3_effectDemDiff2 = matrix(c(0,0,0,0,-1,1,0,0,0),1); rownames(model3_effectDemDiff2) = "None-Fut Diff in Women"
model3_effectDemDiff3 = matrix(c(0,0,0,0,-1,0,0,-1,0),1); rownames(model3_effectDemDiff3) = "None-Cur Diff in Men"
model3_effectDemDiff4 = matrix(c(0,0,0,0,-1,0,0,0,0),1); rownames(model3_effectDemDiff4) = "None-Cur Diff in Women"
model3_effectDemDiff5 = matrix(c(0,0,0,0,0,-1,0,0,-1),1); rownames(model3_effectDemDiff5) = "Fut-Current Diff in Men"
model3_effectDemDiff6 = matrix(c(0,0,0,0,0,-1,0,0,0),1); rownames(model3_effectDemDiff6) = "Fut-Current Diff in Women"
model3_effectsvec2 = rbind(model3_effectDemDiff1, model3_effectDemDiff2,
model3_effectDemDiff3, model3_effectDemDiff4,
model3_effectDemDiff5, model3_effectDemDiff6)
model3_effects2 = glht(model3,linfct=model3_effectsvec2)
summary(model3_effects2)
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMWF +
demgroupF + sexMWF * demgroupF, data = data01)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
None-Fut Diff in Men == 0 -6.056 1.635 -3.704 0.00123 **
None-Fut Diff in Women == 0 -5.892 1.278 -4.611 < 1e-04 ***
None-Cur Diff in Men == 0 -11.971 2.245 -5.332 < 1e-04 ***
None-Cur Diff in Women == 0 -19.846 2.029 -9.783 < 1e-04 ***
Fut-Current Diff in Men == 0 -5.915 2.587 -2.287 0.10786
Fut-Current Diff in Women == 0 -13.954 2.239 -6.233 < 1e-04 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)# Differences in simple effects = interactions;
model3_effectIntDiff1 = matrix(c(0,0,0,0,0,0,0,-1,1),1); rownames(model3_effectIntDiff1) = "A: Sex Effect differ between None and Future?"
model3_effectIntDiff2 = matrix(c(0,0,0,0,0,0,0,-1,1),1); rownames(model3_effectIntDiff2) = "A: None-Future Effect differ by Sex?"
model3_effectIntDiff3 = matrix(c(0,0,0,0,0,0,0,-1,0),1); rownames(model3_effectIntDiff3) = "B: Sex Effect differ between None and Current?"
model3_effectIntDiff4 = matrix(c(0,0,0,0,0,0,0,-1,0),1); rownames(model3_effectIntDiff4) = "B: None-Current Effect differ by Sex?"
model3_effectIntDiff5 = matrix(c(0,0,0,0,0,0,0,0,-1),1); rownames(model3_effectIntDiff5) = "C: Sex Effect differ between Future and Current?"
model3_effectIntDiff6 = matrix(c(0,0,0,0,0,0,0,0,-1),1); rownames(model3_effectIntDiff6) = "C: Future-Current Effect differ by Sex?"
model3_effectsvec3 = rbind(model3_effectIntDiff1, model3_effectIntDiff2,
model3_effectIntDiff3, model3_effectIntDiff4,
model3_effectIntDiff5, model3_effectIntDiff6)
model3_effects3 = glht(model3,linfct=model3_effectsvec3)
summary(model3_effects3)
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = cognition ~ age85 + grip9 + age85 * grip9 + sexMWF +
demgroupF + sexMWF * demgroupF, data = data01)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
A: Sex Effect differ between None and Future? == 0 -0.1643 2.0705 -0.079 0.9964
A: None-Future Effect differ by Sex? == 0 -0.1643 2.0705 -0.079 0.9964
B: Sex Effect differ between None and Current? == 0 7.8751 3.0245 2.604 0.0242 *
B: None-Current Effect differ by Sex? == 0 7.8751 3.0245 2.604 0.0242 *
C: Sex Effect differ between Future and Current? == 0 8.0394 3.4152 2.354 0.0470 *
C: Future-Current Effect differ by Sex? == 0 8.0394 3.4152 2.354 0.0470 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)Moral of the story: Given a significant interaction, it is very likely that the marginal tests provided by default will be misleading and/or will not help you test any meaningful hypotheses. Get simple effects!