Использование multicore в R для анализа данных GWAS

Я использую R для анализа данных широкогеномного ассоциативного исследования. У меня есть около 500 000 потенциальных предикторных переменных (однонуклеотидных полиморфизмов, или SNPs), и я хочу проверить ассоциацию между каждой из них и непрерывным результатом (в данном случае концентрацией липопротеина низкой плотности в крови).

Я уже написал сценарий, который делает это без проблем. Чтобы кратко объяснить, у меня есть объект данных, называемый "Data". Каждая строка соответствует определенному пациенту в исследовании. Есть столбцы для возраста, пола, индекса массы тела (ИМТ) и концентрации ЛПНП в крови. Есть также полмиллиона других столбцов с данными о SNP.

В настоящее время я использую цикл for для запуска линейной модели полмиллиона раз, как показано на рисунке:

# Repeat loop half a million times
for(i in 1:500000) {

# Select the appropriate SNP
SNP <- Data[i]

# For each iteration, perform linear regression adjusted for age, gender, and BMI and save the result in an object called "GenoMod"
GenoMod  <- lm(bloodLDLlevel ~ SNP + Age + Gender + BMI, data = Data)

# For each model, save the p value and error for each SNP. I save these two data points in columns 1 and 2 of a matrix called "results"
results[i,1] <- summary(GenoMod)$coefficients["Geno","Pr(>|t|)"]
results[i,2] <- summary(GenoMod)$coefficients["Geno","Estimate"]
}

Все это работает хорошо. Однако мне бы очень хотелось ускорить анализ. Поэтому я экспериментировал с пакетами multicore, DoMC и foreach.

Мой вопрос в том, не мог бы кто-нибудь помочь мне адаптировать этот код, используя схему foreach?

Я запускаю скрипт на сервере Linux, который, очевидно, имеет 16 ядер. Я пробовал экспериментировать с пакетом foreach, и мои результаты с его использованием были сравнительно хуже, то есть на выполнение анализа с использованием foreach уходит больше времени.

Например, я попробовал сохранить объекты линейной модели, как показано на рисунке:

library(doMC)
registerDoMC()
results <- foreach(i=1:500000) %dopar% { lm(bloodLDLlevel ~ SNP + Age + Gender + BMI, data = Data) }

Это занимает более чем в два раза больше времени, чем использование обычного цикла for. Любой совет о том, как сделать это лучше или быстрее, будет принят с благодарностью! Я понимаю, что использование параллельной версии lapply может быть вариантом, но не знаю, как это сделать.

Всего наилучшего,

Alex

7
задан Alexander 13 December 2011 в 11:35
поделиться