library(psych) sat.act$education <- factor(sat.act$education, labels = paste0("L", 0:5)) sat.act$gender <- factor(sat.act$gender, labels = c("male", "female")) genderDiff <- aov(ACT ~ gender, data = sat.act) objanova <- summary(genderDiff) objttest <- t.test(ACT ~ gender, data = sat.act, var.equal = TRUE) # Notice: equal p, equal df, t^2 = F objanova[[1]][1,4] # Get F value objttest[["statistic"]] # Get t value educDiff <- aov(ACT ~ education, data = sat.act) aggregate(ACT ~ education, data = sat.act, FUN = mean) plot(aggregate(ACT ~ education, data = sat.act, FUN = mean), type = "l") # Post hoc comparison TukeyHSD(educDiff, "education") # Contrast library(multcomp) pairwise <- glht(educDiff, linfct = mcp(education = "Tukey")) summary(pairwise, test = adjusted(type = "bonferroni")) confint(pairwise) confint(pairwise, test = adjusted(type = "bonferroni")) plot(confint(pairwise)) ?p.adjust # x - mean(x) 0:5 - mean(0:5) ctr <- c(-2.5, -1.5, -0.5, 0.5, 1.5, 2.5) ctr <- matrix(ctr, nrow = 1) ctr <- rbind(ctr) linear <- glht(educDiff, linfct = mcp(education = ctr)) summary(linear) linear <- c(-5, -3, -1, 1, 3, 5) quadratic <- c(5, -1, -4, -4, -1, 5) cubic <- c(-5, 7, 4, -4, -7, 5) quartic <- c(1, -3, 2, 2, -3, 1) quintic <- c(-1, 5, -10, 10, -5, 1) ctr2 <- rbind(linear, quadratic, cubic, quartic, quintic) polynomial <- glht(educDiff, linfct = mcp(education = ctr2)) summary(polynomial, test = adjusted(type = "bonferroni")) ctr3 <- t(poly(1:6, 4)) polynomial2 <- glht(educDiff, linfct = mcp(education = ctr3)) summary(polynomial2, test = adjusted(type = "bonferroni")) # Two-way ANOVA twoway <- aov(ACT ~ education + gender + education:gender, data = sat.act) twoway <- aov(ACT ~ education*gender, data = sat.act) summary(twoway) # Equal to type I in SAS library(car) Anova(twoway, type = 2) Anova(twoway, type = 3) aggregate(ACT ~ education*gender, data = sat.act, FUN = mean) library(phia) interactionMeans(twoway) plot(interactionMeans(twoway)) testInteractions(twoway, fixed = "education", across = "gender") testInteractions(twoway, fixed = "gender", across = "education") testInteractions(twoway, fixed = "gender", pairwise = "education") testInteractions(twoway, fixed = "gender", across = "education", adjustment = "bonferroni") ctr <- cbind(linear, quadratic) total.ctr <- list(education = ctr) testInteractions(twoway, custom = total.ctr, fixed = "gender")