如何生成行列频率均衡的40名参与者饮品随机分配序列？

问题 #

需要给40名参与者随机分配6种饮品的领取顺序（每日1种），需满足两个核心条件：

• 每位参与者每种饮品各领取一次
• 每日每种饮品在参与者中的出现次数尽可能均匀

最初直接给每个参与者生成随机排列的代码如下：

library(ggplot2)
set.seed(123)
random_order <- as.data.frame(replicate(40, sample(1:6, 6, replace = F)))
random_order$trial <- c(1:6)

但检查每日各饮品的出现次数时，发现分布差异明显：

tidyr::pivot_longer(random_order, cols = c(1:40), 
                    names_to = "participant", values_to = "drink_order") |>
  dplyr::group_by(trial, drink_order) |>
  dplyr::summarise(count = dplyr::n())

输出显示部分天数的饮品出现次数波动大（比如第1天饮品4仅出现4次，饮品1出现9次），密度图也显示频率分布分散。需要调整生成逻辑，让每日各饮品的出现次数更集中，缩小count的数值范围。

解决方案 #

独立生成每个参与者的随机排列无法保证每日的均匀性，我们可以通过构造平衡排列或使用专业包生成区块设计来实现需求：

方法1：手动构造平衡排列 #

总人数40，6种饮品，每日每种饮品的理想次数为40%/%6=6，剩余4个名额分配给4种饮品（即这4种每日出现7次）。先固定每日的饮品分布，再调整每个参与者的序列确保是完整排列：

library(dplyr)
library(tidyr)
library(ggplot2)

set.seed(123)

# 构造每日的饮品分配（4种7次，2种6次）
daily_drinks <- lapply(1:6, function(t) {
  selected <- sample(1:6, 4)
  drinks <- c(rep(selected, each=7), rep(setdiff(1:6, selected), each=6))
  sample(drinks)
})

# 转换格式并调整每个参与者的序列为完整排列
daily_df <- do.call(cbind, daily_drinks) %>%
  as.data.frame() %>%
  mutate(participant = 1:40) %>%
  pivot_longer(-participant, names_to = "trial", values_to = "drink") %>%
  mutate(trial = as.integer(gsub("V", "", trial)))

balanced_order <- daily_df %>%
  group_by(participant) %>%
  mutate(
    drink = ifelse(any(duplicated(drink)),
                   sample(setdiff(1:6, drink), sum(duplicated(drink))),
                   drink)
  ) %>%
  ungroup() %>%
  pivot_wider(names_from = participant, values_from = drink) %>%
  mutate(trial = 1:6)

# 验证每日计数
balanced_order %>%
  pivot_longer(-trial, names_to = "participant", values_to = "drink_order") %>%
  group_by(trial, drink_order) %>%
  summarise(count = n())

方法2：使用blockrand包生成平衡区块设计 #

blockrand可直接生成多周期的平衡随机化分配，更高效严谨：

library(blockrand)

set.seed(123)
# 生成40个参与者的6周期平衡分配
balanced_design <- blockrand(n = 40, num.levels = 6, block.sizes = rep(6, 40), 
                            id.prefix = "P", block.prefix = "T", 
                            randomize = TRUE)

# 整理为行是trial、列是participant的格式
balanced_order <- balanced_design %>%
  mutate(trial = as.integer(gsub("T", "", block)),
         participant = as.integer(gsub("P", "", id))) %>%
  pivot_wider(names_from = participant, values_from = treatment) %>%
  arrange(trial)

效果验证 #

运行以下代码查看新的计数分布和密度图：

balanced_order %>%
  pivot_longer(-trial, names_to = "participant", values_to = "drink_order") %>%
  group_by(trial, drink_order) %>%
  summarise(count = n()) %>%
  ggplot(aes(count)) +
  geom_density()

此时count只会在6和7之间波动，分布高度集中，完全满足需求。

内容的提问来源于stack exchange，提问作者Yaaman