python - 使用 Python 检索丢失的序列 -'split' 命令不起作用

Question

我有一组使用软件找到​​的(蛋白质)序列，但它们的长度比数据库中原始序列的长度短。我下载了整个数据库，现在我有这些不完整的序列已找到的序列以及从中找到序列的原始数据库。

来自软件的示例结果:

>tr|E7EWP2|E7EWP2_HUMAN  Uncharacterized protein OS=Homo sapiens GN=TRIO PE=4 SV=2
KEFIMAELIQTEKAYVRDLRECMDTYLWEMTSGVE

数据库中的序列:

>tr|E7EWP2|E7EWP2_HUMAN  Uncharacterized protein OS=Homo sapiens GN=TRIO PE=4 SV=2
ARRKEFIMAELIQTEKAYVRDLRECMDTYLWEMTSGVEEIP

所以缺失的残基是“ARR”，最后是“EIP”，我有大约 70 个这样的不完整序列？我想编写一个 Python 程序，可以自动从数据库中检索完整的序列。我真的是 python 的新手，当然我会尝试编写自己的代码，我想知道是否有任何库或类似 biopython 模块的东西可以做到这一点。 我的计划是从我的结果中提取间隔，扩展它们并在原始数据库中选择它，但我不知道如何进一步进行。

我想得到 list_seq = [ARR,KEFIMAELIQTEKAYVRDLRECMDTYLWEMTSGVE,EIP] 以便我可以进一步使用 list_seq[0] r.strip(3) 和 list_seq[1] l.strip[3] 这样我就能得到完整的序列。但 list_seq 不起作用。

提前致谢

Answer 1

您可以使用 BioPython SeqIO 中的 index 方法。 index 方法通过蛋白质 id 对数据库记录进行索引，并且不会将整个数据库加载到内存中，从而允许使用完整的大型数据库进行高效搜索(或者您可以使用传统的数据库例如 sqlite 首先存储您的记录并对其执行搜索):

from Bio import SeqIO

db1 = "dbase.fasta"
db2 = "my_collection.fasta"

dbase_dict = SeqIO.index(db1, "fasta")
my_record_dict = SeqIO.index(db2, "fasta")

for record in my_record_dict:
    if record in dbase_dict:
        rec_dbase = dbase_dict[record]
        rec_mine = my_record_dict[record]
        if rec_dbase.seq != rec_mine.seq:
            print rec_dbase

这个程序只是打印有差异的记录。从这一点开始，您可以将它们保存在列表中或写入 fasta 文件