R commands for Tukey test Data at Table 5-10 of text y <- c(5, 4, 6, 3, 5, 3, 1, 4, 2, 3, 1, 1, 3, 1, 2) temp <- as.factor(rep(c(100, 125, 150), each=5)) press <- as.factor(rep(c(25, 30, 35, 40, 45), times=3)) data <- data.frame(y, temp, press) data # Use the function tukey.1df() to carry out Tukey's 1 df test # Copy and paste everything between the lines into R ============================================================ tukey.1df <- function(y, factorA, factorB) { factorA <- as.factor(factorA) factorB <- as.factor(factorB) a <- length(unique(factorA)) b <- length(unique(factorB)) Y <- matrix(nrow=a, ncol=b) for(i in 1:a) {for (j in 1:b) Y[i,j] <- y[factorA == unique(factorA)[i] & factorB == unique(factorB)[j]]} Ameans <- apply(Y, 1, mean) Bmeans <- apply(Y, 2, mean) term1 <- t(Ameans)%*%Y%*%Bmeans y.. <- mean(Y) SSA <- b*sum((Ameans-y..)^2) SSB <- a*sum((Bmeans-y..)^2) term2 <- SSA + SSB + a*b*y..^2 SSN <- a*b*((term1 - y..*term2)^2)/(SSA*SSB) SST <- sum((Y-y..)^2) SSE <- SST-SSN-SSA-SSB SS <- c(SSA, SSB, SSN, SSE, SST) df <- round(c(a-1, b-1, 1, (a-1)*(b-1)-1, a*b-1)) MS <- SS[1:4]/df[1:4] F0 <- c(MS[1]/MS[4], MS[2]/MS[4], MS[3]/MS[4]) p <- vector(length=3) for(i in 1:3) p[i] <-1-pf(F0[i], df[i], df[4]) SS <- round(SS,3) MS <- round(MS,3) F0 <- round(F0,3) p <- round(p,4) MS <- c(MS, " ") F0 <- c(F0, " ", " ") p <- c(p, " ", " ") out <- cbind(SS, df, MS, F0, p) dimnames(out) <- list(c("A","B","N","Err","Tot"),c("SS","df","MS", "F0", "p")) print.matrix(out, quote=F) list(out = out) } =============================================================== > h <- tukey.1df(y,temp,press) SS df MS F0 p A 23.333 2 11.667 42.949 1e-04 B 11.6 4 2.9 10.676 0.0042 N 0.099 1 0.099 0.363 0.566 Err 1.901 7 0.272 Tot 36.933 14 # Compare with Table 5-11 of the text.