pdf("E:\\sta4211\\KNNL_TA21_01.pdf") risk <- read.table("http://www.stat.ufl.edu/~winner/sta4211/datasets/CH21TA01.txt",header=F, col.names=c("conf","block", "method")) attach(risk) block <- factor(block) method <- factor(method); levels(method) <- c("Utility","Worry","Comparison") mean(conf) tapply(conf,method,mean) tapply(conf,block,mean) interaction.plot(method,block,conf,xlab="Method",ylim=c(0,20)) risk.aov <- aov(conf ~ method + block) # Obtain ANOVA & F-test anova(risk.aov) TukeyHSD(risk.aov,"method") yhat <- predict(risk.aov) e <- residuals(risk.aov) plot(yhat,e,xlab="Yhat",ylab="Residual",main="Residuals versus Yhat") qqnorm(e); qqline(e) stripchart(e ~ block,xlab="Residuals",ylab="Block") friedman.test(conf ~ method | block) # Friedman's Test ### Tukey's Additivity Test (mu_hat <- mean(conf)) (methodmean <- as.vector(tapply(conf,method,mean))) (blockmean <- as.vector(tapply(conf,block,mean))) (r <- length(methodmean)) (nb <- length(blockmean)) (tau_hat <- methodmean-mu_hat) (rho_hat <- blockmean-mu_hat) (tau_hat_y <- rep(tau_hat,nb)) (rho_hat_y <- rep(rho_hat,each=r)) (SSTO <- sum((conf-mean(conf))^2)) (SSTR <- sum(tau_hat_y^2)) (SSBL <- sum(rho_hat_y^2)) (SSBL.TR <- SSTO-SSTR-SSBL) (SSBL.TR_Tukey <- ((sum(tau_hat_y*rho_hat_y*conf))^2)/ ((sum(tau_hat^2))*(sum(rho_hat^2)))) (SSRem <- SSBL.TR-SSBL.TR_Tukey) (F_BL.TR_Tukey <- (SSBL.TR_Tukey/1)/(SSRem/((r-1)*(nb-1)-1))) (F_05 <- qf(0.95,1,(r-1)*(nb-1)-1)) (P_F_BL.TR <- 1-pf(F_BL.TR_Tukey,1,(r-1)*(nb-1)-1)) MSBL <- SSBL/(nb-1) MSBL.TR <- SSBL.TR/((nb-1)*(r-1)) (eff_RBD <- ((nb-1)*MSBL + nb*(r-1)*MSBL.TR) /((nb*r-1)*MSBL.TR)) df1 <- nb*(r-1); df2 <- (nb-1)*(r-1) (eff_RBD_mod <- (df2+1)*(df1+3)*eff_RBD/((df2+3)*(df1+1))) dev.off()