Сделайте так, чтобы функции, использующие apply, stringr, stringi и rbind, работали быстрее

Question

Фон: Я собираюсь предоставить фон для применения этого кода и программный фон.Надеюсь, оба помогут.Я занимаюсь вычислительной работой в области геномики.Да, просто еще один биолог, выдававший себя за компьютерщика.Я работаю над сценарием, который позволит мне интегрировать наборы данных по каждой позиции в геноме человека.Это переводит в фрейм данных, который составляет более 3 миллиардов строк на 12 столбцов.В качестве тестового набора данных я строю свой конвейер анализа с использованием дрожжевого генома, который сгенерирует кадр данных с примерно 25 миллионами строк и 12 столбцами.

Проблема: Мой текущий код работаетхорошо, но жестоко медленно.Например, я запустил свой трубопровод 45 минут назад, и он составляет около 1/3 пути через дрожжевой геном.Это означает, что, вероятно, потребуется 135 минут, чтобы закончить один образец дрожжей, или 270 часов для 1 образца человека ... теперь умножьте это на 90 образцов человека, которые я готовлю к анализу, и вы можете надеяться увидеть мою проблему.Мне нужно ускорить это.Я буду распараллеливать это, но даже тогда я думаю, что сам код слишком неуклюжий.Я хочу помочь моей существующей функции намного, намного быстрее.Пожалуйста, не говорите мне, что мне нужно распараллелить это (это получит отрицательный голос).

Пример данных:

chrom <- c("chr1", "chr1", "chr1", "chr1")
start <- c("0","1","2","6")
stop <- c("1","2","6","7")
sequence <- c("a", "t", "tcag", "a")
seqData <- data.frame(chrom, start, stop, sequence)

Пример вывода:

chrom_out <- c("chr1", "chr1", "chr1", "chr1", "chr1", "chr1", "chr1")
start_out <- c("0", "1", "2", "3", "4", "5", "6")
stop_out <- c("1", "2", "3", "4", "5", "6", "7")
sequence_out <- c("a", "t", "t", "c", "a", "g", "a")
out_seqdata <- data.frame(chrom_out, start_out, stop_out, sequence_out)

Текущий код:

library(dplyr)
library(stringi)
library(stringr) 


wl = function(x){

  length<- stri_length(x["sequence"])
  if(length ==1){
    tmpseq<- x["sequence"]
    tmpstart <- as.numeric(x["start"])
    tmpstop <- as.numeric(x["stop"])
    tmpchrom <- x["chrom"]
    tmpdf <- data.frame(tmpseq, tmpstart, tmpstop, tmpchrom)
    colnames(tmpdf)<- c("tmpseq", "tmpstart", "tmpstop", "tmpchrom")
    print(tmpdf)
  }else{
    tmpseq<- strsplit(x["sequence"], "(?<=.{1})", perl = TRUE)
    tmpstart <- as.numeric(x["start"])+(1:length-1)
    tmpstop<- as.numeric(x["start"])+(1:length)
    tmpdf <- data.frame(tmpseq, tmpstart, tmpstop)
    tmpdf$tmpchrom <- x["chrom"]
    colnames(tmpdf)<- c("tmpseq", "tmpstart", "tmpstop", "tmpchrom")
    print(tmpdf)
  }
}

Объяснение кода: Я использую apply для итерации по каждой строке кадра данных.Фрейм данных представляет собой список координат и геномную последовательность для этих координат.Chrom = хромосома, start = начальная позиция в хромосоме, stop = позиция остановки, а последовательность - это фактическая последовательность.Данные в настоящее время представлены в сжатом формате, примером которого является третий ряд данных.Я хочу расширить эти данные таким образом, чтобы каждая геномная буква стала отдельной строкой, а затем соответствующим образом скорректировать диапазон координат.Функция wl (от широкого до длинного) выполняет это.Сначала он определяет длину строки последовательности.Если длина равна 1, она возвращает эту строку в качестве кадра данных без дальнейших манипуляций;иначе он разбивает строку на отдельные символы, определяет координаты для каждого символа и возвращает этот фрейм данных.В результате получается список фреймов данных, которые затем объединяются, производя выходные данные примера.

Что мне нужно: Я собираюсь разбить геном на части, создав список, что позволит мнераспараллелить этот список.Чанки приведут к серии кадров длиной около 25 миллионов строк.Я собираюсь распараллелить несколько выборок тоже.Распараллеливание внутри распараллеливания ... звучит как отличный способ разбить кластер.Я знаю, как это сделать (и написать этот код, и разбить кластер).Мне нужна помощь, чтобы ускорить фактическую функцию.25 миллионов строк все еще занимает много времени для обработки с использованием моей текущей функции.Любые идеи очень приветствуются.Пожалуйста, измените мою функцию или порекомендуйте новый подход - все идеи приветствуются.Я не знаю о более быстрых методах, кроме добавления большего количества лошадиных сил.

Answer 1

Вы можете векторизовать все свои операции:

# Generate vector of start positions
# Goes from 0 (minimal position in given data) to maximum base position in chromosome
foo <- 0:max(as.numeric(as.character(seqData$start)))
# Split sequence into a character vector
bar <- unlist(strsplit(as.character(seqData$sequence), ""))
# Generate final data frame
data.frame(start = foo, end = foo + 1, seq = bar)
#   start end seq
# 1     0   1   a
# 2     1   2   t
# 3     2   3   t
# 4     3   4   c
# 5     4   5   a
# 6     5   6   g
# 7     6   7   a

Вы можете использовать этот код для перебора хромосом по одной за раз.

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

wl <- function(data, chr) {
    startPos <- 0:max(as.numeric(as.character(data$start)))
    nucs     <- unlist(strsplit(as.character(data$sequence), ""))
    data.frame(chr, start = startPos, end = startPos + 1, seq = nucs)
}
library(foreach)
# use dopar for parallel computations 
foreach(i = unique(seqData$chr), .combine = rbind) %do% {
    wl(subset(seqData, chrom == i), i)
}

PS: я бы никогда не использовал бы координаты генома как символьный вектор.Кроме того, создание столбца end - это пустая трата пространства, поскольку вы знаете, что он расположен на 1 из start.