Расчет корреляционной матрицы с помощью пакета pspearman с функцией apply ()

Я пытаюсь вычислить корреляцию Спирмена и p-значение для кадра данных. Для лучшего приближения значения p я должен придерживаться пакета pspearman. Я ожидаю аналогичного результата с функцией rcorr(). Но у меня проблема при выполнении pspearman:spearman.test() построчно.

Мой фрейм данных содержит 5000 строк (гены) и 200 столбцов (пятна). И я хочу получить матрицу корреляции и матрицу значений p для этих 5000 * 5000 пар ген-ген. Корреляция рассчитывается только в том случае, если оба гена не являются НА более чем в двух точках.

Я могу добиться этого с помощью циклов, но это слишком медленно для моего большого набора данных. У меня возникают проблемы, когда я пытаюсь использовать apply(),sapply(),mapply() для увеличения скорости.

Вот что я пробовал:

data = data.frame(matrix(rbinom(10*100000, 50, .5), ncol=200))
dim(data) #5000, 200
rownames(data) <- paste("gene", 1:5000, sep="") 
colnames(data) <- paste("spot",1:200,sep='')

library(pspearman)
spearFunc = function(x,y=data) {
  df = rbind(x,y)
  # Check the number of complete spots.There are no NAs in this set.
  complete = sum(!(is.na(x)) & !(is.na(y)))
  if (complete >=2 ) {
    pspearman::spearman.test(as.numeric(x),as.numeric(y))
    # This function returns a list containing 8 values, like pvalue,correlation
    }}

pair.all1 = mapply(spearFunc,data,data)
dim(pair.all1)
# 8 200, 200 is the number of columns 
pair.all2 = apply(data,1,spearFunc)

Что приводит к ошибке:

Ошибка в pspearman :: spearman.test (as.numeric (x), as.numeric (y)): объект (list) не может быть приведен к типу 'double'

Я надеюсь использовать spearman.test для каждой пары генов с помощью apply (), чтобы сделать

spearman.test(data[gene1],data[gene1]) 
spearman.test(data[gene1],data[gene2])
....
spearman.test(data[gene1],data[gene5000])
...
spearman.test(data[gene5000],data[gene5000])

Он должен вернуть фрейм данных из 8 строк и 25 000 000 столбцов (5000 * 5000 пар генов).

Можно ли использовать apply () внутри apply () для достижения моей цели?

Спасибо!


r apply correlation sapply mapply
person FewKey    schedule 24.08.2018    source источник

Ответы (1)


arrow_upward
0
arrow_downward

Рассмотрите возможность создания попарных комбинаций генов из row.names с combn, а затем итерации по списку пар с помощью определенной функции. Обязательно возвращайте структуру NA из логики if, чтобы избежать NULL в выводе матрицы.

Однако имейте в виду, что парные перестановки 5000 генов (choose(5000, 2)) приводят к очень высоким результатам - 12 497 500 элементов! Следовательно, sapply (сам цикл) может не сильно отличаться по производительности от for. Рассмотрите возможность распараллеливания итерации.

gene_combns <- combn(row.names(data), 2, simplify = FALSE)

spear_func <- function(x) {
  # EXTRACT ROWS BY ROW NAMES  
  row1 <- as.numeric(data[x[1],])
  row2 <- as.numeric(data[x[2],]) 

  # Check the number of complete spots.There are no NAs in this set.
  complete = sum(!(is.na(x)) & !(is.na(y)))

  if (complete >=2 ) {
    pspearman::spearman.test(row1, row2)        
  } else {
    c(statistic=NA, parameter=NA, p.value=NA, estimate=NA, 
      null.value=NA, alternative=NA,   method=NA, data.name=NA)
  }
}

pair.all2 <- sapply(gene_combns, spear_func)

Тестирование

Выше было протестировано с cor.test (точно такие же входные аргументы и выходной список, что и spearman.test, но более точный p-value) с использованием небольшой выборки данных (50 obs, 20 vars):

set.seed(82418)
data <- data.frame(matrix(rbinom(10*100000, 50, .5), ncol=200))[1:50, 1:20]
rownames(data) <- paste0("gene", 1:50) 
colnames(data) <- paste0("spot", 1:20)

gene_combns <- combn(row.names(data), 2, simplify = FALSE)
# [[1]]
# [1] "gene1" "gene2"    
# [[2]]
# [1] "gene1" "gene3"    
# [[3]]
# [1] "gene1" "gene4"    
# [[4]]
# [1] "gene1" "gene5"    
# [[5]]
# [1] "gene1" "gene6"    
# [[6]]
# [1] "gene1" "gene7"

test <- sapply(gene_combns, spear_func)  # SAME FUNC BUT WITH cor.test
test[,1:5]

#             [,1]                              [,2]                             
# statistic   885.1386                          1659.598                         
# parameter   NULL                              NULL                             
# p.value     0.1494607                         0.2921304                        
# estimate    0.3344823                         -0.2478179                       
# null.value  0                                 0                                
# alternative "two.sided"                       "two.sided"                      
# method      "Spearman's rank correlation rho" "Spearman's rank correlation rho"
# data.name   "row1 and row2"                   "row1 and row2"                  
#             [,3]                              [,4]                             
# statistic   1554.533                          1212.988                         
# parameter   NULL                              NULL                             
# p.value     0.4767667                         0.7122505                        
# estimate    -0.1688217                        0.08797877                       
# null.value  0                                 0                                
# alternative "two.sided"                       "two.sided"                      
# method      "Spearman's rank correlation rho" "Spearman's rank correlation rho"
# data.name   "row1 and row2"                   "row1 and row2"                  
#             [,5]                             
# statistic   1421.707                         
# parameter   NULL                             
# p.value     0.7726922                        
# estimate    -0.06895299                      
# null.value  0                                
# alternative "two.sided"                      
# method      "Spearman's rank correlation rho"
# data.name   "row1 and row2"
person Parfait    schedule 24.08.2018
comment
Большое спасибо за ваше терпение. Функция combn () и попытка избежать NULL очень полезны. Я использовал spearman.test для вычисления p-значения, потому что точность p-value cor.test () зависит от нулевого распределения данных и статистических связей. Хотя это называется точным, но невозможно вычислить реальное значение p, когда n ›22. При n ›22 лучше использовать приближения, как в функции spearman.test. Это значительно снижает количество ложных срабатываний. (ссылка: Быстрое вычисление точного нулевого распределения L-статистики Спирмена и Пейджа для выборок со связями и без них) - person FewKey; 25.08.2018
comment
Здорово! Рад помочь. Спасибо за информацию - TIL. Удачного кодирования! - person Parfait; 25.08.2018