使用 Biopython 替换文件之间的序列

Replace sequences between files using Biopython

我有两个蛋白质序列 FASTA 文件:

nsp.fasta --> 原文件

wsp.fasta --> 来自信号肽预测工具的输出文件,该工具返回 nsp.fasta 中的蛋白质并去除了信号。

例如:

记录在nsp.fasta:

>gi|564250271|ref|XP_006264203.1| PREDICTED: apolipoprotein D [Alligator mississippiensis]
MRGMLALLAALLGLLGLVEGQTFHMGQCPNPPVQEDFDPSKYLGKWYEIEKLPSGFEQER
CVQANYSLKANGKIKVLTKMVRSAQHLTCLQHRMMLLVSSPVMPASPYWVVATDYENYAL
VYSCTSFFWLFHVDYAWIRSRTPQLHPETVEHLKSVLRSYRIQTGMMLPTDQMNCPSDM

记录在wsp.fasta:

>gi|564250271|ref|XP|006264203.1|  PREDICTED: apolipoprotein D [Alligator mississippiensis]; MatureChain: 21-179
QTFHMGQCPNPPVQEDFDPSKYLGKWYEIEKLPSGFEQERCVQANYSLKANGKIKVLTKM
VRSAQHLTCLQHRMMLLVSSPVMPASPYWVVATDYENYALVYSCTSFFWLFHVDYAWIRS
RTPQLHPETVEHLKSVLRSYRIQTGMMLPTDQMNCPSDM

然而,并非 nsp.fasta 中的所有蛋白质都包含信号肽,因此 wsp.fasta 是 nsp.fasta 中包含信号的蛋白质子集。我需要的是一个包含所有蛋白质记录的独特文件,包括未发现信号肽的蛋白质和去除信号肽的成熟链。

我试过以下方法:

from Bio import SeqIO

file1 = SeqIO.parse(r"c:\Users\Sergio\Desktop\nsp.fasta", "fasta")

file2 = SeqIO.parse(r"c:\Users\Sergio\Desktop\wsp.fasta", "fasta")

for seq1 in file1:
    for seq2 in file2:
        if seq2.id == seq1.id:
            seq1.seq = seq2.seq
            SeqIO.write(seq1, r"c:\Users\Sergio\Desktop\nuevsp.fasta", "fasta")

但是根本没有输出。我试过将 SeqIO.write 放在循环之外,它 returns 是一个空白文件。我究竟做错了什么?是否已经存在合并两个文件或用另一个文件中的序列替换一个文件中的序列的方法?

提前致谢!!

塞尔吉奥

编辑代码,我添加了一个 elif 子句,试图在 nsp.fasta 中添加与 wsp.fasta 不匹配的记录,但它不起作用:

to_write = []

for seq1 in SeqIO.parse(r"c:\Users\Sergio\Desktop\nsp.txt", "fasta"):
    for seq2 in SeqIO.parse(r"c:\Users\Sergio\Desktop\wsp.txt", "fasta"):
        if seq1.id == seq2.id:
            seq1.seq = seq2.seq
            to_write.append(seq1)
        elif seq1.id != seq2.id:
            to_write.append(seq1)

SeqIO.write(to_write, r"c:\Users\Sergio\Desktop\nuevsp.txt", "fasta")

就像你写的那样,每次你写一个新的序列,你都会覆盖前一个序列。尝试将记录存储在列表中,然后在循环完成时写出列表。

to_write = []
for seq1 in SeqIO.parse(r"c:\Users\Sergio\Desktop\nsp.fasta", "fasta"):
    for seq2 in SeqIO.parse(r"c:\Users\Sergio\Desktop\wsp.fasta", "fasta"):
        if seq2.id == seq1.id:
            seq1.seq = seq2.seq
            to_write.append(seq1)
SeqIO.write(to_write, r"c:\Users\Sergio\Desktop\nuevsp.fasta", "fasta")

编辑建议另一种使用列表理解的方法:

ids_to_save = [x.id for x in SeqIO.parse(r"c:\Users\Sergio\Desktop\nsp.fasta", "fasta")]
records_to_save = [x for x in SeqIO.parse(r"c:\Users\Sergio\Desktop\wsp.fasta", "fasta") if (x.id in ids_to_save)]
SeqIO.write(records_to_save, r"c:\Users\Sergio\Desktop\nuevsp.fasta", "fasta")

编辑以解决 "add the records in nsp.fasta that doesn't match wsp.fasta" 需要 - 一般方法,不一定是确切代码:

ids_not_wanted = [x.id for x in SeqIO.parse(r"c:\Users\Sergio\Desktop\wsp.fasta", "fasta")]
records_to_save_2 = [x for x in SeqIO.parse(r"c:\Users\Sergio\Desktop\wsp.fasta", "fasta") if (x.id not in ids_not_wanted)]

records_to_save.append(records_to_save_2)
# If duplicate records are a problem, eliminate them using "set"
records_to_save = list(set(records_to_save))
SeqIO.write(records_to_save, r"c:\Users\Sergio\Desktop\nuevsp.fasta", "fasta")