对数据框中的唯一组合进行统计测试

Question

对于下面的示例数据框，我想对所有独特的 DRUG - ADR 组合执行统计测试（例如 t 检验）。为此，我需要以下内容：

1) 每个独特的 DRUG - ADR 组合的 X 向量

2) 如果我感兴趣的 DRUG - ADR 是 D1 - A1，我想用向量测试 X 的向量（这里是 34）：

• D1 - 所有不是 A1 的 A（在示例 D1-A2 中，x = 37）
• A1 - 不是 D1 的所有 D（在示例中为 D4-A1，x = 65）

此过程应遍历数据框中的所有记录，并且应忽略 ID 变量，因为一个 ID 可以有多个 DRUG - ADR 组合。显然，我的数据集要大得多，从 X 得到的向量将包含超过 1 个值

dat <- data.frame(ID=c(1,2,3,4,4,4,5,6,6,7),
                  DRUG=c("D1","D2","D2","D3","D3","D3","D5","D1","D4","D2"),
                  ADR=c("A1","A2","A3","A6","A7","A8","A4","A2","A1","A2"),
                  X=c(34,76,34,45,2,41,56,37,65,12))


   ID DRUG ADR  X
1   1   D1  A1 34
2   2   D2  A2 76
3   3   D2  A3 34
4   4   D3  A6 45
5   4   D3  A7  2
6   4   D3  A8 41
7   5   D5  A4 56
8   6   D1  A2 37
9   6   D4  A1 65
10  7   D2  A2 12

期待您的建议！

其他信息： 我无法使用 boxplot() 解释这一点，所以我会口头尝试。在这个例子中，我想为每个独特的 DRUG - ADR 组合对 X 执行 2 次测试。因此，对于 D1-A1，我想测试与 D1-A1 和 D1 - A1 的对比。对于测试编号 1，我需要两个向量：一个包含 D1-A1 的 X 值（此处为 34），另一个包含 D1 - A1 的 X 值（因此在本例中为 D1-A2（值 = 37 ) ，但在我更大的数据集中，这将是一个更大的向量，其中包含例如 D1 - A23、D1-A573 等的值。对于测试号 2，我需要用于 D1-A1（我已经从测试 1 中获得）和 D1-A1（此处为 D4-A1，值 = 65）的 X 向量。然后应该对记录 2 (D2-A2)、记录 3 (D2-A3) 等重复相同的过程。因此，对于数据框中的每条记录，我需要来自 X 的 2 组 2 个向量。

接下来是记录 2。对于此记录，D2-A2 的向量将是 (76,12)，而 D2-A2 的向量将是 (37；来自记录 8)，而 D2 的向量-A2 将是（34；来自记录 3）。我希望这能解释它

附加注释：下面的语法描述了 D1-A1 的 2 个测试。但是，我想要一个通用函数，它遍历所有独特的 Dx-Ax 组合并执行相同的过程（因此在这种情况下，9 次 2 测试，因为记录 2 和 10 适用于相同的组合）。

x <- ifelse (dat$DRUG == "D1" & dat$ADR == "A1",dat$X, NA)
x <- x[!is.na(x)]

y <- ifelse (dat$DRUG != "D1" & dat$ADR == "A1",dat$X, NA)
y <- y[!is.na(y)]

z <- ifelse (dat$DRUG == "D1" & dat$ADR != "A1",dat$X, NA)
z <- z[!is.na(z)]

t.test(x,y)
t.test(x,z)

Answer 1

你可以试试这样的。 ifelse 非常适合分组

group1 <- ifelse(dat$DRUG == "D1", ifelse(dat$ADR  == "A1", 1, 2), NA)
group2 <- ifelse(dat$ADR == "A1", ifelse(dat$DRUG  == "D1", 1, 2), NA)
par(mfrow=c(1,2))
boxplot(dat$X ~ group1)
boxplot(dat$X ~ group2)

# for the test you can try following:
t.test(dat$X ~ group1)
t.test(dat$X ~ group2)

对于循环，你可以试试这个。如果小组观察太少，将被打断，并会为您提供所有结果的列表。否则你可以使用 t.test()p.value 仅报告 pvalues。

res <- list()
n <- 2 # change here the number of tests
for (x in 1:n){
  D <- paste0("D", x)
  A <- paste0("A", x)
  group1 <- factor(ifelse(dat$DRUG == D, ifelse(dat$ADR  == A, paste0(D,A), paste0(D,"Ax")), NA))
  group2 <- factor(ifelse(dat$ADR == A, ifelse(dat$DRUG  == D, paste0(D,A),  paste0("Dx",A)), NA))
  g1 <- paste0(levels(group1), collapse = "_vs_")
  g2 <- paste0(levels(group2), collapse = "_vs_")
  if(nlevels(group1) == 2){
    res1 <- t.test(dat$X ~ group1)$p.value
  }else{
    res1 <- NA
  }
  if(nlevels(group2) == 2){
    res2 <- t.test(dat$X ~ group2)$p.value
  }else{
    res2 <- NA
  }
  res[[x]] <- cbind.data.frame(g1, res1, g2, res2)
}
do.call("rbind", res)