 # Read data

> data <- read.table("~/Documents/subversion/concurrency/ITICSE 2015/stats_R/anova/theory/data.csv", sep=",", header=T)

 # Calculate means

> means<- round(tapply(data$nmark, data$row, mean), digits=2)
> means
  r1   r2   r3   r4 
1.16 1.05 0.88 0.99 

 # Two ways of plotting raw data

> library(gplots)  #I load the "gplots" package to plot means

> plotmeans(data$nmark~data$row, digits=2, ccol="red", ylab = "normalized mark",xlab = "row", main="Plot of normalized mark by row", mean.labels=TRUE)

> boxplot(data$nmark ~ data$row, main="Normalized mark by row (mean is black dot)", xlab="row", ylab="normalized mark", col=rainbow(7))
> points(means, col="black", pch=18)

 # Volumes

> volumes <- c(length(data[data[, "row"]=="r1", ]$row), length(data[data[, "row"]=="r2", ]$row), length(data[data[, "row"]=="r3", ]$row), length(data[data[, "row"]=="r4", ]$row))
> volumes
[1] 127 163 218 135
> round(100*(volumes/643), digits=2)
[1] 19.75 25.35 33.90 21.00
> sum(volumes)
[1] 643

 # ANOVA

> aov_cont<- aov(data$nmark ~ data$row)
> summary(aov_cont)
             Df Sum Sq Mean Sq F value   Pr(>F)    
data$row      3   6.78  2.2611   8.348 1.89e-05 ***
Residuals   639 173.07  0.2708                     
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1



If the overall P value is small, then it is unlikely that the differences you observed are due to random sampling. You can reject the idea that all the populations have identical means. This doesn't mean that every mean differs from every other mean, only that at least one differs from the rest. Look at the results of post tests to identify where the differences are.

  # TukeyHSD: To see the differences between each row
  # r5-r3 are practically the same (p = 0.9937223)
  # r3-r1 there is a significance difference (p = 0.0000000)

> tuk<- TukeyHSD(aov_cont)
> tuk
  Tukey multiple comparisons of means
    95% family-wise confidence level

Fit: aov(formula = data$nmark ~ data$row)

$`data$row`
             diff         lwr          upr     p adj
r2-r1 -0.11365080 -0.27231408  0.045012475 0.2532406
r3-r1 -0.28048935 -0.43013101 -0.130847683 0.0000103
r4-r1 -0.17122186 -0.33693426 -0.005509469 0.0397484
r3-r2 -0.16683854 -0.30564623 -0.028030853 0.0109949
r4-r2 -0.05757106 -0.21356960  0.098427482 0.7774458
r4-r3  0.10926748 -0.03754578  0.256080754 0.2218324

# In following plot: Significant differences are the ones which not cross the zero value.

> plot (tuk)
psig=as.numeric(apply(tuk$`data$row`[,2:3],1,prod)>=0)+1
op=par(mar=c(2,3,3.8,2))
plot(tuk,col=psig,yaxt="n")
for (j in 1:length(psig)){
axis(2,at=j,labels=rownames(tuk$`data$row`)[length(psig)-j+1],
     las=1,cex.axis=.8,col.axis=psig[length(psig)-j+1])
}
axis
par(op)



Para hacer mantel test creamos un nuevo conjunto de datos donde las row dejan de ser r1, r2,etc para ser 1,2,etc.

> data <- read.table("~/Documents/subversion/concurrency/ITICSE 2015/stats_R/anova/theory/data_mantel.csv", sep=",", header=T)

> diffrow.dist <- dist(data$row)
> diffmark.dist <- dist(data$nmark)

> mantel.rtest(diffmark.dist, diffrow.dist, nrepet = 999)
Monte-Carlo test
Observation: 0.01634641 
Call: mantel.rtest(m1 = diffmark.dist, m2 = diffrow.dist, nrepet = 999)
Based on 999 replicates
Simulated p-value: 0.069 