使用dplyr / tidyverse同时对多个变量进行多重配对t检验

tmf*_*mnk 7 r dplyr

假设这样的数据结构:

   ID testA_wave1 testA_wave2 testA_wave3 testB_wave1 testB_wave2 testB_wave3
1   1           3           2           3           6           5           3
2   2           4           4           4           3           6           6
3   3          10           2           1           4           4           4
4   4           5           3          12           2           7           4
5   5           5           3           9           2           4           2
6   6          10           0           2           6           6           5
7   7           6           8           4           6           8           3
8   8           1           5           4           5           6           0
9   9           3           2           7           8           4           4
10 10           4           9           5          11           8           8
Run Code Online (Sandbox Code Playgroud)

我想要实现的是分别为每个测试计算一个配对的t检验(在这种情况下,这意味着testA和testB,但在现实生活中,我有更多的测试)。我想这样做,以便将给定测试的第一波与同一测试的所有其他后续波进行比较(在testA的情况下,意味着testA_wave1与testA_wave2以及testA_wave1与testA_wave3)。

这样,我能够实现:

df %>%
 gather(variable, value, -ID) %>%
 mutate(wave_ID = paste0("wave", parse_number(variable)),
        variable = ifelse(grepl("testA", variable), "testA",
                     ifelse(grepl("testB", variable), "testB", NA_character_))) %>%
 group_by(wave_ID, variable) %>% 
 summarise(value = list(value)) %>% 
 spread(wave_ID, value) %>% 
 group_by(variable) %>% 
 mutate(p_value_w1w2 = t.test(unlist(wave1), unlist(wave2), paired = TRUE)$p.value,
        p_value_w1w3 = t.test(unlist(wave1), unlist(wave3), paired = TRUE)$p.value) %>%
 select(variable, matches("(p_value)"))

  variable p_value_w1w2 p_value_w1w3
  <chr>           <dbl>        <dbl>
1 testA           0.664        0.921
2 testB           0.146        0.418
Run Code Online (Sandbox Code Playgroud)

但是,我希望看到不同/更优雅的解决方案,它们给出相似的结果。我主要是在寻找dplyr/ tidyverse解决方案,但是如果有一种完全不同的方法来实现它,那么我并不反对。

样本数据:

set.seed(123)
df <- data.frame(ID = 1:20,
testA_wave1 = round(rnorm(20, 5, 3), 0),
testA_wave2 = round(rnorm(20, 5, 3), 0),
testA_wave3 = round(rnorm(20, 5, 3), 0),
testB_wave1 = round(rnorm(20, 5, 3), 0),
testB_wave2 = round(rnorm(20, 5, 3), 0),
testB_wave3 = round(rnorm(20, 5, 3), 0))
Run Code Online (Sandbox Code Playgroud)

Ron*_*hah 10

dplyr0.8.0开始,我们可以group_split用来将数据帧拆分为数据帧列表。

我们gather将数据框转换为长格式,然后separate将列的名称(key)转换为不同的列(testwave)。然后group_split,我们用于将数据框基于test列拆分为列表。对于列表中的每个数据帧,我们将spread其转换为宽格式,然后计算t.test值,然后使用将其重新绑定到一个数据帧中map_dfr

library(tidyverse)

df %>%
  gather(key, value, -ID) %>%
  separate(key, c("test", "wave")) %>%
  group_split(test) %>% #Previously we had to do split(.$test) here
  map_dfr(. %>%
          spread(wave, value) %>%
          summarise(test = first(test),
                    p_value_w1w2 = t.test(wave1, wave2, paired = TRUE)$p.value, 
                    p_value_w1w3 = t.test(wave1, wave3, paired = TRUE)$p.value))


# A tibble: 2 x 3
#  test  p_value_w1w2 p_value_w1w3
#  <chr>        <dbl>        <dbl>
#1 testA        0.664        0.921
#2 testB        0.146        0.418
Run Code Online (Sandbox Code Playgroud)

由于只有2个值需要计算,因此我们手动执行了t检验。如果wave...列数更多,那么这将变得很麻烦。在这种情况下,我们可以

df %>%
   gather(key, value, -ID) %>%
   separate(key, c("test", "wave")) %>%
   group_split(test) %>% 
   map_dfr(function(data) 
              data %>%
                   spread(wave, value) %>%
                   summarise_at(vars(setdiff(unique(data$wave), "wave1")), 
                   function(x) t.test(.$wave1, x, paired = TRUE)$p.value) %>%
                   mutate(test = first(data$test)))

#  wave2 wave3 test 
#  <dbl> <dbl> <chr>
#1 0.664 0.921 testA
#2 0.146 0.418 testB
Run Code Online (Sandbox Code Playgroud)

在这里,它将对每个带有“ wave1”列的“ wave ..”列执行t检验。


由于您也愿意接受其他解决方案,因此请尝试使用纯基础R解决方案

sapply(split.default(df[-1], sub("_.*", "", names(df[-1]))), function(x) 
 c(p_value_w1w2 = t.test(x[[1]], x[[2]],paired = TRUE)$p.value, 
   p_value_w1w3 = t.test(x[[1]], x[[3]],paired = TRUE)$p.value))


#                 testA     testB
#p_value_w1w2 0.6642769 0.1456059
#p_value_w1w3 0.9209554 0.4184603
Run Code Online (Sandbox Code Playgroud)

我们根据划分列,test*并创建数据框列表,然后t.test为每个数据框应用不同的列组合。


dip*_*kov 5

这是一种做到这一点的方法,需要花费purrr很多时间。

library("tidyverse")

set.seed(123)
df <- tibble(
  ID = 1:20,
  testA_wave1 = round(rnorm(20, 5, 3), 0),
  testA_wave2 = round(rnorm(20, 5, 3), 0),
  testA_wave3 = round(rnorm(20, 5, 3), 0),
  testB_wave1 = round(rnorm(20, 5, 3), 0),
  testB_wave2 = round(rnorm(20, 5, 3), 0),
  testB_wave3 = round(rnorm(20, 5, 3), 0)
)

pvalues <- df %>%
  # From wide tibble to long tibble
  gather(test, value, -ID) %>%
  separate(test, c("test", "wave")) %>%
  # Not stricly necessary; will order the waves alphabetically instead
  mutate(wave = parse_number(wave)) %>%
  inner_join(., ., by = c("ID", "test")) %>%
  # If there are two waves w1 and w2,
  # we end up with pairs (w1, w1), (w1, w2), (w2, w1) and (w2, w2),
  # so filter out to keep the pairing (w1, w2) only
  filter(wave.x == 1, wave.x < wave.y) %>%
  nest(ID, value.x, value.y) %>%
  mutate(pvalue = data %>%
           # Perform the test
           map(~t.test(.$value.x, .$value.y, paired = TRUE)) %>%
           map(broom::tidy) %>%
           # Also not strictly necessary; you might want to keep all
           # information about the test: estimate, statistic, etc.
           map_dbl(pluck, "p.value"))
pvalues
#> # A tibble: 4 x 5
#>   test  wave.x wave.y data              pvalue
#>   <chr>  <dbl>  <dbl> <list>             <dbl>
#> 1 testA      1      2 <tibble [20 x 3]>  0.664
#> 2 testA      1      3 <tibble [20 x 3]>  0.921
#> 3 testB      1      2 <tibble [20 x 3]>  0.146
#> 4 testB      1      3 <tibble [20 x 3]>  0.418

pvalues %>%
  # Drop the data in order to pivot the table
  select(- data) %>%
  unite("waves", wave.x, wave.y, sep = ":") %>%
  spread(waves, pvalue)
#> # A tibble: 2 x 3
#>   test  `1:2` `1:3`
#>   <chr> <dbl> <dbl>
#> 1 testA 0.664 0.921
#> 2 testB 0.146 0.418
Run Code Online (Sandbox Code Playgroud)

reprex软件包(v0.2.1)创建于2019-03-08