My question is whether any R packages out there provide function(s) that make it reasonably easy to randomize a standard experimental design that may involve crossed factors, nesting, and/or blocking.
For concreteness, show me specifically how to produce a new randomization of the Oats experiment provided as a dataset in the nlme package.
> data(Oats, package="nlme")
> summary(Oats)
Block Variety nitro yield
VI :12 Golden Rain:24 Min. :0.00 Min. : 53.0
V :12 Marvellous :24 1st Qu.:0.15 1st Qu.: 86.0
III:12 Victory :24 Median :0.30 Median :102.5
IV :12 Mean :0.30 Mean :104.0
II :12 3rd Qu.:0.45 3rd Qu.:121.2
I :12 Max. :0.60 Max. :174.0
In that experiment, there are six blocks. Each block is divided into three plots which are randomly assigned to varieties (one plot per variety in each block, each block separately randomized). Each plot is subdivided into 4 subplots, and randomly assigned to four amounts of nitrogen (0, 0.2, 0.4, and 0.6), one subplot per nitro level randomized separately in each plot. In the dataset, the plots are identifiable as the combinations of Block and Variety. (The response variable is yield, so that is not really part of the treatment design.)
Second question: Given you can randomize Oats, can you use the same package to easily randomize other classical designs, e.g., a three-factor CRD, a nested design, a three-period crossover design, or a 5x5 Graeco-Latin square?
I actually already know how to do this using base functions in the R language; so I'm not particularly interested in seeing a programmatic answer. I want to know if any existing packages make this easy. I can identify a few packages that may help, e.g., randomizeR and randomizr, but a quick reading of the documentation of these still does not make it very obvious (to me) how to do this.
I have the makings of a general-purpose randomization package I developed a few years ago for my students, and I am trying to decide whether to develop it further for release on CRAN.
Here's how I'd do it with randomizr:
data(Oats, package="nlme")
# get the latest version from github
install.packages("devtools")
devtools::install_github("DeclareDesign/randomizr")
library(randomizr)
Oats <- within(Oats,{
Variety_new <- block_ra(block_var = Block,
condition_names = c("Golden Rain", "Marvellous", "Victory"))
nitro_new <- block_ra(block_var = paste0(Block, Variety_new),
condition_names = c(0, 0.2, 0.4, 0.6))
})
# Original Random Assignment
with(Oats, table(Block, Variety))
with(Oats, table(Block, nitro))
with(Oats, table(Block, nitro, Variety))
# New Random Assignment
with(Oats, table(Block, Variety_new))
with(Oats, table(Block, nitro_new))
with(Oats, table(Block, nitro_new, Variety_new))
The key is to realize that you need to block on both variety and Block to randomize the subplots into nitro conditions. (that's why we need the call to paste0).
EDIT
here is another way to do this (see comments)
library(randomizr)
des <- rev(expand.grid(subplot=1:4, wholeplot=1:3, block=1:6))
des <- within(des,{
plot_id <- paste0(block, "_", wholeplot)
Variety <- block_and_cluster_ra(
block_var = block,
clust_var = plot_id,
condition_names = c("Golden Rain", "Marvellous", "Victory"))
nitro <- block_ra(block_var = plot_id,
condition_names = c(0, 0.2, 0.4, 0.6))
})
with(des, table(Variety, block))
with(des, table(Variety, nitro))
with(des, table(Variety, nitro, block))