Извлечь имена чтения и связанные нуклеотиды в определенных позициях из файла BAM (в R)

Предположим, у меня есть файл BAM и несколько позиций, которые я хотел бы более внимательно изучить в этом выравнивании. Моя цель - выяснить, находятся ли эти позиции в одних и тех же считываниях и какие нуклеотиды находятся в этих двух положениях. Таким образом, мой вывод должен представлять собой таблицу всех операций чтения, которые происходят в позициях, и список связанных нуклеотидов.

Результат должен выглядеть следующим образом (с 980, 1350 и 1400 первые три позиции; я составил названия для чтения):

В конце я хотел бы иметь таблицу, которая связывает нуклеотиды посредством чтения.

До сих пор я пытался прочитать файл Bam в желаемых местах:

      #Positions of interest
pos <- c(980, 1350, 1400)

#path to BAM file
path <- pathToBamFile

#name of reference used for mapping
refname <- "myreference"

#loading the bam file
bamFile <- BamFile(path)
gr <- GRanges(seqnames = refname, ranges = IRanges(start = pos, end = pos))
params <- ScanBamParam(which = gr, what = scanBamWhat())
aln <- scanBam(bamFile, param = params)

С помощью следующего кода я могу получить доступ к чтению, но как мне получить соответствующие нуклеотиды в отдельных позициях?

      pos.no <- c(1:length(pos))

i <- 1
aln[[pos.no[i]]]$qname

Я уже пробовал работать со струнами СИГАР с последовательностями чтобы извлечь нуклеотиды вручную, но это приводило к все большим и большим проблемам.

Есть ли простой способ извлечь нуклеотиды?

Я долго искал в Интернете простое решение, но безуспешно. Надеюсь, здесь мне кто-нибудь поможет! Большое спасибо за ваше время заранее.

0 ответов

Другие вопросы по тегам