Brassica {pda} | R Documentation |
The scientist examined the effects of rate of
application of a certain herbicide (pyridate) on a range of closely
related crop
s. The experiment was conducted over two
year
s, with three block
s each year
.
The two year
s were conducted in roughly the same location,
but there was separate randomization in each year
. Three
herbicide rates (pyrrat
= 0, 1, 2) were randomly assigned
to separate rows within each crop
plot. Recorded values
consist of yield
per row for each crop*pyrrat
combination within each block
for each year
.
Since crop
2 was not used in the second year
, it
is dropped from further consideration.
data(Brassica)
Brassica data frame with 144 observations on 10 variables.
[,1] | crop | factor | crop identifier (1-7) |
[,2] | pyrrat | ordered | pyridate level (0-2) |
[,3] | block | factor | block (1-3) |
[,4] | inj7dat | numeric | injury after 7 days |
[,5] | in14dat | numeric | injury after 14 days |
[,6] | yldkga | numeric | yield in kilograms |
[,7] | pconyld | numeric | percent control yield |
[,8] | pcnblk | numeric | percent control in block yield |
[,9] | year | ordered | year - 1990 |
[,10] | trt | factor | treatment code |
Alan Miller, U WI Horticulture (http://www.hort.wisc.edu)
Miller A (1993) `Pyridate tolerance in Brassicaceae crops', PhD Dissertation, Department of Horticulture, University of Wisconsin-Madison.
# full data set--code not ready yet! # model fits for bill and leg data( Brassica ) # drop Crop 2 Brassica <- data.frame(Brassica[ Brassica$crop != 2, ]) Brassica$crop <- factor( as.character( Brassica$crop )) # Analysis based on lmer() library(lme4) Brassica1 <- data.frame( Brassica[ Brassica$year == 1, ] ) Brassica.lme <- lmer( yldkga ~ crop * pyrrat + (1|crop:block), data = Brassica1 ) anova( Brassica.lme ) VarCorr( Brassica.lme ) # Analysis based on aov() using Error() specification Brassica.aov <- aov ( yldkga ~ crop * pyrrat + Error( block/crop ), Brassica1 ) summary( Brassica.aov ) # Analysis of full data set: both years Brassica.full <- aov ( yldkga ~ crop * pyrrat * year + Error( year:block + crop:year:block ), Brassica ) Brassica.red <- aov ( yldkga ~ crop * year + pyrrat + crop:pyrrat + Error( year:block + crop:year:block ), Brassica) # Analysis based on lmer() Brassica1 <- data.frame( Brassica[ Brassica$year == 1, ] ) Brassica.lme <- lmer( yldkga ~ crop * pyrrat + (1|block:crop), data = Brassica1 ) Brassica$yrblk = interaction(Brassica[,c("year","block")]) Brassica$cropyrblk = interaction(Brassica[,c("crop","year","block")]) Brassica.fulle = lmer(yldkga ~ crop * pyrrat * year + (1|yrblk:cropyrblk), data = Brassica ) anova(Brassica.fulle) Brassica.rede = lmer(yldkga ~ crop * year + pyrrat + crop:pyrrat + (1|yrblk:cropyrblk), data = Brassica ) anova(Brassica.rede) anova(Brassica.rede,Brassica.fulle) lsd.plot(Brassica.rede,factor=c("pyrrat","crop"), xlab = "pyr rate by crop", ylab = "mean crop yield" ) ## Not run: # The following code is specialized to get plots for the book. # It is not recommended as standard practice, but presented to # show one torturous way to get LSMEANS for a nested design. Brassica.red <- aov (yldkga ~ crop * year + pyrrat + crop:pyrrat + Error(year:block + crop:year:block), Brassica, qr = TRUE) Brassica.proj <- proj ( Brassica.red ) Brassica.ref <- data.frame(Reference=c(Brassica.proj[[1]])) for (i in names(Brassica.proj)[-1]) Brassica.ref[[i]] <- apply( Brassica.proj[[i]], 1, sum ) # whole plot Brassica.wp <- nested( Brassica.proj, Brassica, c("block","year","crop"), "yldkga", c(1,3), Brassica.ref ) Brassica.wpaov <- aov( yldkga ~ year*block + crop + crop:year, Brassica.wp, weight = Brassica.wp$weight) # split plot Brassica.sp <- nested( Brassica.proj, Brassica, c("block","year","crop","pyrrat"), "yldkga", c(1,4), Brassica.ref ) Brassica.spaov <- aov( yldkga ~ year*block*crop + pyrrat + crop:pyrrat, Brassica.sp ) margin.plot( Brassica.wpaov, factors = c("crop","year"), xlab = "(a) crop by year", ylab = "yield in 10000 kg" ) #Figure H:24.1 Brassicaica split split plot anova lsd.plot( Brassica.spaov, factors = c("pyrrat","crop"), xpos = .6, ypos = 3.25, xlab = "(b) pyr rate by crop", ylab = "deviation from mean crop yield" ) ## End(Not run)