# # Exercise 2.4 # # Statistical Methods in Biology: Design and Analysis of Experiments and Regression # by S.J. Welham, S.A. Gezan, S.J. Clark & A. Mead (2014) # Chapman & Hall/CRC Press, Boca Raton, Florida. ISBN: 978-1-4398-0878-8 # Version 1, 16/09/2014 # Set working directory - use setwd() function or from Session menu in RStudio # e.g. setwd("d:/stats4biol/data) # Read data & assign Field as factor carbon <- read.table("carbon.dat",sep="",header=T) carbon$Field <- as.factor(carbon$Field) summary(carbon) # Perform two-sample t-test directly t.test(Carbon~Field, var.equal=TRUE, data=carbon) # Rest of file verifies t-test calculations # Get number of observations, mean and variance for each field N.df <- aggregate(Carbon~Field, FUN="length", data=carbon); N.df mean.df <- aggregate(Carbon~Field, FUN="mean", data=carbon); mean.df var.df <- aggregate(Carbon~Field, FUN="var", data=carbon); var.df # Calculate pooled estimate of variance df <- N.df$Carbon[1]+N.df$Carbon[2]-2; df pooled.var <- ( (N.df$Carbon[1]-1)*var.df$Carbon[1] + (N.df$Carbon[2]-1)*var.df$Carbon[2] )/df pooled.var # Calculate t statistic SED <- sqrt(pooled.var*((1/N.df$Carbon[1])+(1/N.df$Carbon[2]))); SED t <- (mean.df$Carbon[1]-mean.df$Carbon[2])/SED; t # Get critical value # = 97.5th quantile of t distribution with 10 df critical <- qt(0.975,df); critical # Get observed significance level # = proportion of t distribution with 10 df > abs(t) or < -abs(t) P <- 2*(1-pt(abs(t),df)); P # End of file