连接来自同一物种名称下两个文件的 DNA 序列

Question

我有两个 FASTA 文件，其中包含编码两种不同蛋白质的 DNA 序列。我想将不同蛋白质和同一物种的序列连接成一个长序列。

例如，我有：

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()

但是由于 DNA 序列很长并且占文件的很多行，我不知道如何 select 它们。

请帮忙！

Answer 1

通过使用字典，您可以将 fasta 序列附加到每个 ID。然后，将它们打印到输出文件。

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 打开。

Answer 2

既然你标记了 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是你合并后的输出文件。

此外，请注意，从生物学角度来看，我认为这是一件奇怪的事情，跨多个文件连接序列可能会导致创建 'fake' 个连续序列，并且连接的顺序肯定很重要，所以小心

连接来自同一物种名称下两个文件的 DNA 序列

Joining DNA sequences from two files under the same species name

python

bioinformatics

biopython

fasta