# # Exercise 12.2 # # 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 # Data from R Harrington & C Shortall, Rothamsted Research # Version 1, 14/08/2015 # Set working directory - use setwd() function or from Session menu in RStudio # e.g. setwd("d:/stats4biol/data) # Set up packages to be used later - available from CRAN library(ggplot2) # for plotting # Read data hereford <- read.table("hereford.dat",sep="",header=T) summary(hereford) # Get log transformation of wet weights hereford$logWt <- log10(hereford$WetWeight+0.5) # Plot data qplot(y=logWt, x=Year, data=hereford) # Fit SLR model hereford.slr <- lm(logWt~Year, data=hereford) # Use ANOVA to check for strength of relationship anova(hereford.slr) # Summarise using adjusted R-squared summary(hereford.slr)$adj.r.squared # Plot residuals plot(hereford.slr, ask=FALSE) # Regression model parameter estimates summary(hereford.slr)$coefficients # Prediction for 2010 - extrapolation! hereford.pred.slr <- predict(hereford.slr, new=data.frame(Year=c(2010)), interval="confidence", se.fit=TRUE) hereford.pred.slr$fit hereford.pred.slr$se.fit # Back-transform predictions hereford.pred.back <- (10^hereford.pred.slr$fit)-0.5 hereford.pred.back # Get fitted model plot # Set up basic plot using model object plot.obs <- ggplot(data=hereford.slr, aes(y=logWt, x=Year) ) # Plot observations and fitted line plot.obs + geom_point() + geom_line(aes(y=fitted(hereford.slr))) + ggtitle("Observed data with fitted SLR") # Look for temporal correlation res <- residuals(hereford.slr) nv <- length(res) nvm1 <- nv-1 # correlation between residual and residual from previous year cor(res[1:nvm1],res[2:nv]) # plot residual vs that from previous year plot(y=res[1:nvm1],x=res[2:nv]) # plot residuals across time plot(y=res,x=hereford$Year,type='b') # End of file - but see also Exercises 13.4 and 15.1