Объединение последовательностей ДНК из двух файлов под одним названием вида

У меня есть два файла FASTA с последовательностями ДНК, кодирующими два разных белка. Я хочу объединить последовательности для разных белков и одного вида в одну длинную последовательность.

например, у меня есть:

Protein 1
>sce
AGTAGATGACAGCT
>act
GCTAGCTAGCT
Protein 2
>sce
GCTACGATCGACT
>act
TACGATCAGCTA
Protein 1+2
>sce
AGTAGATGACAGCTGCTACGATCGACT
>act
GCTAGCTAGCTTACGATCAGCTA

Что-то, что может быть проблемой, заключается в том, что виды не появляются в одном и том же порядке в обоих файлах, и есть несколько последовательностей, которые находятся в одном, но не в другом (файлы имеют длину около 110 видов, с несоответствие 4 или 5).

Моя первая попытка написать код для него была:

gamma = open('gamma.fas', 'w')
spc = open("spc98.fas", 'w')
outfile = open("joined.fas", 'w')
for line in gamma:
    if line.startswith(">"):
        for line2 in spc:
             if line2.startswith(">"):
                if line == line2:
                    outfile.write(line)
    else:
        outfile.write(line)
fh.close()

но поскольку последовательности ДНК очень длинные и занимают много строк в файле, я не знаю, как их выбрать.

Пожалуйста помоги!


person Ana Catarina Vitorino    schedule 14.08.2019    source источник


Ответы (2)


Используя словарь, вы можете добавлять последовательности fasta к каждому идентификатору. А затем распечатайте их в выходной файл.

outfile = open("joined.fas", 'w')

d = dict()

for file in ('gamma.fas', 'spc98.fas'):
    with open(file, 'r') as f:
        for line in f:
            line = line.rstrip()
            if line.startswith('>'):
                key = line
            else:
                d.setdefault(key, '')
                d[key] += line

for key, seq in d.items():
    outfile.write(key + "\n" + seq + "\n")

outfile.close()

РЕДАКТИРОВАТЬ: Кстати, вы открываете два файла для чтения как открытые для записи, что приведет к стиранию двух входных файлов.

gamma = open('gamma.fas', 'w') spc = open("spc98.fas", 'w')

Их следует открывать с помощью r вместо w.

person Chris Charley    schedule 14.08.2019

Поскольку вы отметили Biopython, вот компактное решение. Обратите внимание, что он помещает весь файл в память (как и большинство простых подходов):

from Bio.Seq import Seq
from Bio import SeqIO

d = SeqIO.to_dict(SeqIO.parse('1.fasta', 'fasta'))

for r in SeqIO.parse('2.fasta', 'fasta'):
    d[r.id] = d.setdefault(r.id, Seq('')) + r.seq

SeqIO.write(d.values(), 'output.fasta', 'fasta')

Здесь 1.fasta и 2.fasta — два входных файла fasta, а output.fasta — объединенный выходной файл.

Кроме того, обратите внимание, что биологически я думаю, что это странная вещь, объединение последовательностей в нескольких файлах может привести к созданию «поддельных» смежных последовательностей, и порядок объединения, безусловно, важен, поэтому будьте осторожны.

person Chris_Rands    schedule 15.08.2019