if (FALSE) {
library(agridat)
data(student.barley)
dat <- student.barley
libs(lattice)
bwplot(yield ~ gen|district, dat, main="student.barley - yield")
dat$year <- factor(dat$year)
dat$income <- NULL
# convert to tons/ha
dat <- transform(dat, yield=yield*14 * 0.00112085116)
# Define 'loc' the way that Kempton does
dat$loc <- rep("",nrow(dat))
dat[is.element(dat$farmer, c("Allardyce","Roche","Quinn")),"loc"] <- "1"
dat[is.element(dat$farmer, c("Luttrell","Dooley")), "loc"] <- "2"
dat[is.element(dat$year, c("1904","1905","1906")) & dat$farmer=="Kearney","loc"] <- "2"
dat[dat$farmer=="Mulhall","loc"] <- "3"
dat <- transform(dat, loc=factor(paste(place,loc,sep="")))
libs(reshape2)
datm <- melt(dat, measure.var='yield')
# Kempton Table 9.5
round(acast(datm, loc+gen~year),2)
# Kempton Table 9.6
d2 <- dcast(datm, year+loc~gen)
mean(d2$Archer)
mean(d2$Goldthorpe)
mean(d2$Archer-d2$Goldthorpe)
sqrt(var(d2$Archer-d2$Goldthorpe)/51)
cor(d2$Archer,d2$Goldthorpe)
if(0){
# Kempton Table 9.6b
libs(lme4)
m2 <- lmer(yield~1 + (1|loc) + (1|year) +
(1|loc:year) + (1|gen:loc) + (1|gen:year), data=dat,
control=lmerControl(check.nobs.vs.rankZ="ignore"))
}
}
Run the code above in your browser using DataLab