在 R 中随机化裂区(和其他经典设计)

Randomizing a split-plot (and other classical designs) in R

我的问题是是否有任何 R 包提供的功能可以相当容易地随机化可能涉及交叉因素、嵌套、and/or 分块的标准实验设计。

为了具体起见,请具体说明如何对 nlme 包中作为数据集提供的 Oats 实验进行新的随机化。

> 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 

在那个实验中,有六个区块。每个块分为三个随机分配给品种的地块(每个块中每个品种一个地块,每个块单独随机分配)。每个地块被细分为 4 个子地块,并随机分配四种氮含量(0、0.2、0.4 和 0.6),每个地块中每个硝基水平分别随机分配一个子地块。在数据集中,这些图可识别为 BlockVariety 的组合。 (响应变量是 yield,因此这实际上不是治疗设计的一部分。)

第二个问题:鉴于你可以随机化 Oats,你能否使用相同的包轻松随机化其他经典设计,例如三因子 CRD、嵌套设计、三周期交叉设计,还是 5x5 希腊拉丁方块?

我实际上已经知道如何使用 R 语言中的基函数来完成此操作;所以我对看到程序化的答案不是特别感兴趣。我想知道是否有任何现有的 packages 使这变得容易。我可以确定一些可能有帮助的软件包,例如 randomizeRrandomizr,但快速阅读这些文档仍然无法实现非常明显(对我来说)如何做到这一点。

我已经准备好几年前为我的学生开发的通用随机化包,我正在尝试决定是否进一步开发它以在 CRAN 上发布。

以下是我使用 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))

关键是要意识到您需要对变体和块进行阻止,以将子地块随机化为硝基条件。 (这就是为什么我们需要调用 paste0)。

编辑

这是另一种方法(见评论)

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))