R многоядерный сверхвысокомерный одномерный анализ

Question

У меня есть фрейм данных, состоящий из 1М ковариат, где я хотел бы регрессировать каждый из них, независимо от конкретного столбца фрейма данных, используя R, но в многоядерном режиме. Под однофакторным анализом я подразумеваю либо биномиальную регрессию, либо критерий Вилкоксона.

мой текущий код это

library(MASS)
library(doParallel)
tt=dat.at.fil
nm.cores = detectCores() - 1
cl=makeCluster(nm.cores)
registerDoParallel(cl)
x <- foreach(cnt=1:nrow(tt),.combine=cbind) %dopar% {
whol.dat = data.frame(log10(t(tt)[,cnt]), y=factor(my.y))
deviance(glm(y~., data = whol.dat[-which(whol.dat[,1] == -Inf),], family = "binomial"))
}

или же

library(MASS)
library(doParallel)
tt=dat.at.fil
nm.cores = detectCores() - 1
cl=makeCluster(nm.cores)
registerDoParallel(cl)
x <- foreach(cnt=1:nrow(tt),.combine=cbind) %dopar% {
whol.dat = data.frame(t(tt)[,cnt], y=factor(my.y))
wilcoxon.test(y~., data = whol.dat))
     }

Интересно, как я могу улучшить это, чтобы быть еще более эффективным?

Answer 1

Это отличный пример того, как шаг разбивки вашей задачи на части может реально повлиять на преимущества параллельных вычислений. Я пропустил некоторые части вашего кода (например, преобразование логарифмов ковариат и устранение пропущенных значений), которые не являются существенными для проблемы. Я думаю, что вы хотите избежать транспонирования всей матрицы при каждом вызове - просто сделайте это один раз в верхней части вашего скрипта. AFAIK R хранит данные в главном порядке столбцов, поэтому простое избегание этого шага при работе со столбцами может сэкономить вам немного времени.

В первом испытании я сначала запустил версию в сериале, чтобы увидеть, насколько много улучшений. Это на четырехъядерном процессоре AMD Phenom 9850 с тактовой частотой 2,5 ГГц и 8 ГБ ОЗУ (такой старый).

library(doParallel)
library(iterators)
#make covariate data
N = 100
P = 100000 # number of predictors
tt = as.data.frame(matrix(rnorm(N*P),nrow=N,ncol=P))
my.y = rbinom(N,p=0.5,size=1)
y = factor(my.y)

# How fast to do it serially?
system.time(x1 <- foreach(cc = iter(tt, by='col'),.combine=c) %do% {
  deviance(glm(y~cc, family = "binomial"))
}) # elapsed 718 s

nm.cores = detectCores() - 1
cl=makeCluster(nm.cores)
registerDoParallel(cl)

# send entire dataframe to each worker, pull out the desired column
system.time(x2 <- foreach(cnt=1:ncol(tt),.combine=c) %dopar% {
  whol.dat = data.frame(tt[,cnt], y=factor(my.y))
  deviance(glm(y~., data = whol.dat, family = "binomial"))
}) # elapsed 276 s, so 3 x faster

all.equal(x1,x2) # TRUE, just checkin' ...

Моей первой мыслью было, что отправка всей матрицы каждому работнику каждый раз может повлечь за собой некоторые накладные расходы, поэтому я переписал foreach() использовать iter() просто отправить каждый столбец работнику:

system.time(x3 <- foreach(cc = iter(tt, by='col'),.combine=c) %dopar% {
  deviance(glm(y~cc, family = "binomial"))
}) # not much faster, 248s

И это немного ускоряет, но не сильно. Раньше я не использовал итераторы, поэтому, читая виньетку foreach, я натолкнулся на собственный итератор. iblkcol() это разбивает data.frame на фрагменты и отправляет каждый блок для экономии на накладных расходах по отправке данных и их получению от рабочих. Код для этого скрыт на Github (см. Строки 199-218).

## from vignette on foreach:
## use iblkcol() instead of iter in loop to send blocks of columns instead of one at a time
system.time(x4 <- foreach(cc = iblkcol(tt, chunks = nm.cores),.combine=c,.packages='foreach') %dopar% {
  foreach(x = 1:col(cc),.combine=c) %do% {
    deviance(glm(y~cc[,x], family = "binomial"))
  }
}) # 193 s! 

И это существенное улучшение по сравнению с отправкой каждого столбца по одному. Я думаю, что могут быть некоторые дополнительные ускорения путем настройки вызова glm(), чтобы воспользоваться тем фактом, что большая часть фрейма модели используется повторно от одного вызова к другому. То же самое должно работать с вызовом wilcoxon().