如果我正确解释它,这应该有效。
我们将进行两步采样:
- 首先,对处理组本身进行采样,从而更容易确定块中的特定行是否与前一个块的同一行位于同一处理组中。
- 其次,从每个经过验证的安全组中抽取一个。
我将在这里使用随机种子以获得可重复性,不要set.seed(.)
在生产中使用。
set.seed(42)
nBlocks <- 7
treatments <- list(1:7, 8:14, 15:21, 22:28, 29:35, 36:42, 43:49, 50:56, 57:63)
blocks <- Reduce(function(prev, ign) {
while (TRUE) {
this <- sample(length(treatments))
if (!any(this == prev)) break
}
this
}, seq.int(nBlocks)[-1], init = sample(length(treatments)), accumulate = TRUE)
blocks <- do.call(cbind, blocks)
blocks
# [,1] [,2] [,3] [,4] [,5] [,6] [,7]
# [1,] 1 3 4 2 8 2 1
# [2,] 5 1 2 4 5 7 9
# [3,] 9 8 9 3 1 3 5
# [4,] 7 9 3 6 7 9 3
# [5,] 2 4 8 5 4 1 4
# [6,] 4 7 1 9 6 4 2
# [7,] 8 6 5 7 2 6 8
# [8,] 3 5 6 8 9 5 6
# [9,] 6 2 7 1 3 8 7
这里每一列都是一个“块”,每个数字代表分配给每一行的治疗组。您可以看到后续列中没有任何行包含相同的组。
例如,第一列(“块 1”)将在第一行包含来自治疗 1 组的内容,在第二行包含治疗 5 组,等等。此外,检查将显示所有治疗都包含在每个块列中,一个推断实验设计的要求。
(仅供参考,根据随机条件,理论上这可能需要一段时间。因为它每列重复,但它应该是相对有效的。我在这里没有针对太长执行的保护措施,但我没有t 认为这是必需的:这里的条件不会导致需要大量重复的“失败”的高可能性。)
下一步是将这些组编号中的每一个转换为来自相应治疗组的编号。
apply(blocks, 1:2, function(ind) sample(treatments[[ind]], 1))
# [,1] [,2] [,3] [,4] [,5] [,6] [,7]
# [1,] 6 17 22 11 54 14 3
# [2,] 30 3 13 22 33 48 58
# [3,] 63 55 61 15 4 21 33
# [4,] 49 60 21 36 43 58 21
# [5,] 12 25 55 32 27 7 25
# [6,] 24 46 4 58 38 28 11
# [7,] 53 38 35 49 11 36 56
# [8,] 16 29 36 56 63 29 40
# [9,] 36 8 47 3 19 50 43
为了验证,在第一个矩阵中,我们的前三行(块 1)是 1、5 和 9,它们应该分别转换为 1-7、29-35 和 57-63。“6”在 1-7 以内,“30”在 29-35 以内,“63”在 59-63 以内。检查将显示其余部分是正确的。
由于首先确定治疗组的步骤,验证/保证您不会在两个相邻块之间连续重复治疗组要简单得多。
编辑
规则:
- 同一治疗组可能不在相邻列的同一行;和
- 相同的处理(不是组)可能不在相邻列的任何行中。
我们可以使用与以前相同的方法。请注意,随着任何组变得更小,迭代时间可能会增加,但我不认为它可能会进入无限循环。(但是,如果你不小心有一组长度为 1,那么......这将永远不会结束。)
nBlocks <- 7
treatments <- list(1:7, 8:14, 15:21, 22:28, 29:35, 36:42, 43:49, 50:56, 57:63)
# helper function for randomized selection of treatments given groups
func <- function(grp) cbind(grp, sapply(treatments[grp], sample, size = 1))
set.seed(42)
func(c(1,3,5))
# grp
# [1,] 1 1
# [2,] 3 19
# [3,] 5 29
然后是同样的Reduce
心态:
set.seed(42)
blocks <- Reduce(function(prev, ign) {
while (TRUE) {
this1 <- sample(length(treatments))
if (!any(this1 == prev[,1])) break
}
while (TRUE) {
this2 <- func(this1)
if (!any(this2[,2] %in% prev[,2])) break
}
this2
}, seq.int(nBlocks-1), init = func(sample(length(treatments))), accumulate = TRUE)
blocks <- do.call(cbind, blocks)
groups <- blocks[, seq(1, by = 2, length.out = nBlocks)]
treats <- blocks[, seq(2, by = 2, length.out = nBlocks)]
由此,我们有两种产品(尽管您可能只关心第二种):
治疗组,很好地验证了上面的规则 1:相邻列中没有组可以在同一行中:
groups
# grp grp grp grp grp grp grp
# [1,] 1 3 1 7 8 5 1
# [2,] 5 1 2 8 2 7 3
# [3,] 9 8 5 2 1 4 6
# [4,] 7 9 6 3 4 8 5
# [5,] 2 4 7 9 3 9 4
# [6,] 4 7 4 5 7 1 2
# [7,] 8 6 9 1 9 6 7
# [8,] 3 5 8 6 5 2 9
# [9,] 6 2 3 4 6 3 8
处理本身,对于上面的规则 2,相邻列中可能没有处理:
treats
#
# [1,] 7 19 2 47 51 33 3
# [2,] 35 4 12 50 8 44 15
# [3,] 60 51 35 10 1 22 41
# [4,] 43 58 41 21 26 55 31
# [5,] 12 24 43 57 17 57 26
# [6,] 27 49 26 34 48 6 11
# [7,] 53 36 62 6 62 36 47
# [8,] 16 33 54 42 32 10 62
# [9,] 37 9 15 27 37 18 56
编辑 2:
另一个规则:
- 每个治疗组必须在每行和每列中恰好出现一次(需要方形实验设计)。
我认为这有效地生成了一个类似数独的治疗组矩阵,一旦满足,回填规则#2(在相邻列中不重复处理)。https://gamedev.stackexchange.com/a/138228建议了一种方法(虽然它很仓促):
set.seed(42)
vec <- sample(9)
ind <- sapply(cumsum(c(0, 3, 3, 1, 3, 3, 1, 3, 3)), rot, x = vec)
apply(ind, 1, function(z) all(1:9 %in% z)) # all rows have all 1-9, no repeats
# [1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
apply(ind, 1, function(z) all(1:9 %in% z)) # ... columns ...
# [1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
ind
# [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
# [1,] 1 7 8 3 5 2 4 6 9
# [2,] 5 2 3 6 9 4 8 1 7
# [3,] 9 4 6 1 7 8 3 5 2
# [4,] 7 8 1 5 2 3 6 9 4
# [5,] 2 3 5 9 4 6 1 7 8
# [6,] 4 6 9 7 8 1 5 2 3
# [7,] 8 1 7 2 3 5 9 4 6
# [8,] 3 5 2 4 6 9 7 8 1
# [9,] 6 9 4 8 1 7 2 3 5
考虑到组的限制,这使得随机组安排的风格相当固定。由于这是一种实验设计,如果您要使用这种方法(并且块之间的接近性完全是一个问题),那么您应该在对治疗本身进行采样之前ind
随机化矩阵的列和/或行。(你可以做列和行,只是分段做,它应该保留约束。)
ind <- ind[sample(9),][,sample(9)]
ind
# [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
# [1,] 2 3 8 1 4 7 9 6 5
# [2,] 7 8 4 6 2 9 5 3 1
# [3,] 1 7 9 4 5 6 3 2 8
# [4,] 8 1 6 9 3 4 2 5 7
# [5,] 5 2 7 8 9 1 6 4 3
# [6,] 3 5 1 7 6 8 4 9 2
# [7,] 4 6 3 5 8 2 7 1 9
# [8,] 6 9 5 2 1 3 8 7 4
# [9,] 9 4 2 3 7 5 1 8 6
从这里,我们可以制定规则 2:
treatments <- list(1:7, 8:14, 15:21, 22:28, 29:35, 36:42, 43:49, 50:56, 57:63)
mtx <- do.call(rbind, Reduce(function(prev, ind) {
while (TRUE) {
this <- sapply(treatments[ind], sample, size = 1)
if (!any(prev %in% this)) break
}
this
}, asplit(ind, 2)[-1],
init = sapply(treatments[ind[,1]], sample, size = 1),
accumulate = TRUE))
mtx
# [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
# [1,] 11 44 4 52 30 15 23 41 59
# [2,] 16 56 49 3 12 33 39 57 27
# [3,] 52 24 60 40 46 2 20 29 13
# [4,] 1 37 23 63 56 48 32 12 17
# [5,] 24 10 30 16 58 39 50 2 47
# [6,] 49 57 41 25 6 52 11 17 34
# [7,] 59 31 19 14 38 23 47 51 7
# [8,] 41 17 11 33 24 61 5 43 54
# [9,] 29 4 51 45 20 8 58 28 40