Как объединить последовательность FASTA в соответствии с идентификатором последовательности?

Question

У меня есть 9 файлов FASTA, представляющих последовательность ДНК 9 генов.

Каждый файл FASTA содержит 121 последовательность, представляющую 121 штамм.Название для каждой последовательности является идентификатором для каждого штамма.

Однако в каждом файле идентификатор не сортируется, например, в gene1.fasta:

>1
AAA
>16
TTT
>2
GGG
...

В gene2.fasta:

>2
CCC
>34
AAA
>1
GGG
...

Iхотите изменить эти файлы 9 генов FASTA на 121 штамм. Файлы FASTA в каждом файле просто объединяют 9 генов для одного штамма.Например, в Strain1.fasta:

AAAGGG

в Strain2.fasta:

GGGCCC

Как я могу сделать это в R?

Answer 1

Вот решение в R по запросу, использующее пакет Biostrings для чтения файлов fasta.

Это работает, но я должен сказать, что это один из самых уродливых кодов, которые я написал в длинномвремя.Я просто хотел посмотреть, смогу ли я как-нибудь это сделать - это 100% не лучшее решение.

library("Biostrings")
library("tidyverse")

convertStringSet = function(seq){
  return(df = data.frame("names" = names(seq), "seq" = paste(seq)))
}

# change the path accordingly
filenames = list.files("/home/x/y/z", pattern="gene*", full.names=TRUE)%>%
  lapply(readDNAStringSet)

fastaDF = filenames %>% lapply(convertStringSet) %>% 
  reduce(full_join, by = "names") %>% 
  unite("seq", -1,  sep="")

writeOutput = function(x){

  header = paste(">",x[1],sep="")
  fileName = paste("strain",x[1],".fasta",sep="")

  writeLines(c(header,x[2]), fileName)
}

apply(fastaDF, 1, writeOutput)

В качестве альтернативы, если вы работаете в системе Unix, эта строка awk должна выдавать то же самое.результаты:

awk '$0 ~ /^>/ {i=substr($0,2); next;} i != -1 {printf "%s", $0 >> "file"i; i=-1;}' gene*