Зацикливание переменной и выполнение перестановочных тестов на разность означает переменную

Question

Итак, у меня есть набор данных, который является гораздо большей версией этого:

value <- c(4, 6, 7, 2, 3, 4)
category <- rep(c("good", "bad"), c(6, 6))
gene <- rep(c("gene_1", "gene_2"), 3)
df <- as.data.frame(cbind(category, gene, value))

   category   gene value
1      good gene_1     4
2      good gene_2     6
3      good gene_1     7
4      good gene_2     2
5      good gene_1     3
6      good gene_2     4
7       bad gene_1     4
8       bad gene_2     6
9       bad gene_1     7
10      bad gene_2     2
11      bad gene_1     3
12      bad gene_2     4

Мне нужно сделать тест на перестановку, чтобы найти средние различия в «значении» между «категориями» для каждого гена .

Мой вопрос состоит из двух частей:

1) Как мне создать внешний цикл, который будет циклически проходить через каждый ген и выполнять тест перестановки,

и

2) является ли этот код для правильного выполнения теста перестановки?

y <- df$value
x <- df$category

obs <- mean(y[x == "good"]) - mean(y[x == "bad"])
## Permutation test
perm <- NULL
for(i in 1:1000){
  y.perm <- sample(y, replace=FALSE, size=100) 
  perm[i] <- mean(y.perm[x == "good"]) - mean(y.perm[x == "bad"])
}

## Two-sided p-value:
twosided <- sum(abs(perm) >= abs(obs))/1000

Вот пример хвоста моего набора данных:

# A tibble: 6 x 3
# Groups:   category, gene [2]
  category gene     value
  <fct>    <chr>    <dbl>
1 bad      gene_145  8.54
2 good     gene_145  8.40
3 good     gene_145  8.46
4 bad      gene_145  8.56
5 good     gene_145  8.25
6 bad      gene_145  8.43

Итак, вы можете видеть, что он в основном идентичен образцу df.Кроме того, вот строка вложенных фактических данных и вложенного образца df:

> str(nest(df, -gene)) #sample df
'data.frame':   2 obs. of  2 variables:
 $ gene: chr  "gene_1" "gene_2"
 $ data:List of 2
  ..$ :'data.frame':    6 obs. of  2 variables:
  .. ..$ category: Factor w/ 2 levels "bad","good": 2 2 2 1 1 1

> str(nest(merged_df, -gene)) #actual data
Classes ‘tbl_df’, ‘tbl’ and 'data.frame':   200 obs. of  3 variables:
 $ category: Factor w/ 2 levels "bad","good": 1 2 1 2 1 2 1 2 1 2 ...
 $ gene    : chr  "gene_1" "gene_1" "gene_3" "gene_3" ...
 $ data    :List of 200
  ..$ :Classes ‘tbl_df’, ‘tbl’ and 'data.frame':    19 obs. of  1 variable:
  .. ..$ value: num  4.84 5.38 5.3 5.47 4.88 ...

Answer 1

Давайте сначала возьмем второй вопрос.Возможно, это лучший вопрос для перекрестной проверки, но, поскольку он тесно связан с вашим вопросом программирования, мы разберем его здесь.Когда дело доходит до перестановочной выборки, следует иметь в виду две основные вещи:

1. Каждая перестановка должна быть такого же размера, что и наблюдаемая выборка.Это гарантирует, что переставленные образцы будут иметь (или, по крайней мере, будут иметь потенциал) ту же дисперсию, что и наблюдаемые.

2. Когда это возможно, следует перемешать метки ,не значения (см. этот перекрестный проверенный вопрос ).Это соответствует первоначальной цели выборки перестановок: создание альтернатив, в которых собранные данные остались прежними, но категории данных оказались произвольно разными.В таких ситуациях, как ваша, обычно следует выбирать метки без замены, чтобы в итоге вы получили такое же количество наблюдений на группу, как и у вас изначально (дисбалансы и все).

Предположим, у нас есть данные этого примера (обратите внимание, что я сделал несколько небольших изменений для исправления ошибок в вашем оригинале):

value <- round(runif(12) * 10)
category <- rep(c("good", "bad"), c(6, 6))
gene <- rep(c("gene_1", "gene_2"), 3)
df <- data.frame(category, gene, value)

   category   gene value
1      good gene_1     8
2      good gene_2     0
3      good gene_1     4
4      good gene_2     8
5      good gene_1     3
6      good gene_2     9
7       bad gene_1     0
8       bad gene_2     3
9       bad gene_1     7
10      bad gene_2     0
11      bad gene_1     5
12      bad gene_2     2

Код теста перестановки может выглядеть следующим образом:

perm <- rep(NA, 1000)
  for (i in 1:1000) {
    labels <- sample(df$category, nrow(df), replace=FALSE)
    perm[i] <- mean(df$value[labels == "good"]) - mean(df$value[labels == "bad"])
  }

Для другихвопрос, мы можем выполнить эту операцию для каждого гена во фрейме данных с помощью tidyverse .Полное объяснение того, как работают эти функции, выходит далеко за рамки вопроса, но вкратце, мы используем nest, чтобы «свернуть» фрейм данных для каждого уровня «гена», во вложенный фрейм данных, обычно называемый «данные".Затем мы используем функции mutate и map для работы с вложенными кадрами, реализуя код перестановки сверху.

library(tidyverse)

df.nest <- nest(df, -gene) %>% 
  dplyr::mutate(
    obs = map_dbl(data, function(x) mean(x$value[x$category == 'good']) - mean(x$value[x$category == 'bad']) ), 
    permutes = map(data, function(x) {
      perm <- rep(NA, 1000)
      for (i in 1:1000) {
        labels <- sample(x$category, nrow(x), replace=FALSE)
        perm[i] <- mean(x$value[labels == "good"]) - mean(x$value[labels == "bad"])
      }
      return(perm)
    }),
    p.val = map2_dbl(obs, permutes, function(obs, permutes) {
      mean(abs(permutes) >= abs(obs))
    })
  )

  gene   data                   obs permutes      p.val
  <fct>  <list>               <dbl> <list>        <dbl>
1 gene_1 <data.frame [6 × 2]>     1 <dbl [1,000]> 0.9  
2 gene_2 <data.frame [6 × 2]>     4 <dbl [1,000]> 0.388