> winejudge.aov <- aov(rating ~ wine + judge) > anova(winejudge.aov) Analysis of Variance Table Response: rating Df Sum Sq Mean Sq F value Pr(>F) wine 3 184.00 61.333 57.5 1.854e-08 *** judge 5 173.33 34.667 32.5 1.549e-07 *** Residuals 15 16.00 1.067 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 > TukeyHSD(winejudge.aov,"wine") Tukey multiple comparisons of means 95% family-wise confidence level Fit: aov(formula = rating ~ wine + judge) $wine diff lwr upr p adj 2-1 2.0000000 0.2814186 3.718581 0.0202300 3-1 6.6666667 4.9480852 8.385248 0.0000001 4-1 6.0000000 4.2814186 7.718581 0.0000003 3-2 4.6666667 2.9480852 6.385248 0.0000061 4-2 4.0000000 2.2814186 5.718581 0.0000375 4-3 -0.6666667 -2.3852481 1.051915 0.6844296 > e <- residuals(winejudge.aov) > > > (winemean <- tapply(rating,wine,mean)) 1 2 3 4 20.00000 22.00000 26.66667 26.00000 > (judgemean <- tapply(rating,judge,mean)) 1 2 3 4 5 6 25 20 21 28 25 23 > > s <- length(judgemean) > r <- length(winemean) > > winemean_y <- rep(winemean,s) > judgemean_y <- rep(judgemean,each=r) > > winedev <- rating-winemean_y > judgedev <- rating-judgemean_y > > > ratedev.mat <- matrix(rep(0,(s*r)),ncol=r) > > for (i1 in 1:r) { + ratedev.mat[,i1] <- winedev[wine==i1] + } > > ratedev.mat [,1] [,2] [,3] [,4] [1,] 0 2 1.3333333 2 [2,] -5 -4 -3.6666667 -2 [3,] -2 -3 -2.6666667 -3 [4,] 6 4 3.3333333 4 [5,] 2 2 1.3333333 0 [6,] -1 -1 0.3333333 -1 > > > cov.mat <- matrix(rep(0,(r^2)),ncol=r) > > for (i1 in 1:r) { + for (i2 in 1:r) { + cov.mat[i1,i2] <- sum(ratedev.mat[,i1]*ratedev.mat[,i2])/(s-1) + }} > > cov.mat [,1] [,2] [,3] [,4] [1,] 14.0 11.0 9.200000 8.2 [2,] 11.0 10.0 8.200000 7.6 [3,] 9.2 8.2 7.066667 6.2 [4,] 8.2 7.6 6.200000 6.8 > > qqnorm(e); qqline(e) > > timeseq <- c(3,2,1,4,1,3,2,4,3,2,4,1,4,3,1,2,3,2,4,1,1,3,2,4) > > par(mfrow=c(2,3)) # Note the dataset doesn't have a variable for time > > plot(timeseq[1:4],e[1:4],main="Judge 1") > plot(timeseq[5:8],e[5:8],main="Judge 2") > plot(timeseq[9:12],e[9:12],main="Judge 3") > plot(timeseq[13:16],e[13:16],main="Judge 4") > plot(timeseq[17:20],e[17:20],main="Judge 5") > plot(timeseq[21:24],e[21:24],main="Judge 6") > > friedman.test(rating ~ wine|judge) Friedman rank sum test data: rating and wine and judge Friedman chi-squared = 16.6842, df = 3, p-value = 0.0008207