使用Biopython通过BLAST获取未知序列的详细信息

2024-10-01 17:40:06 发布

您现在位置:Python中文网/ 问答频道 /正文

我是第一次用生物圈。我有来自未知生物的序列数据,并试图用BLAST来判断它们最有可能来自哪个生物体。为此,我编写了以下函数:

def find_organism(file):
    """
    Receives a fasta file with a single seq, and uses BLAST to find
    from which organism it was taken.
    """
    # get seq from fasta file
    seqRecord = SeqIO.read(file,"fasta")
    # run BLAST
    blastResult = NCBIWWW.qblast("blastn", "nt", seqRecord.seq)
    # get first hit
    blastRecord = NCBIXML.read(blastResult)
    firstHit = blastRecord.alignments[0]
    # get hit's gi number
    title = firstHit.title
    gi = title.split("|")[1]
    # search NCBI for the gi number
    ncbiResult = Entrez.efetch(db="nucleotide", id=gi, rettype="gb", retmode="text")
    ncbiResultSeqRec = SeqIO.read(ncbiResult,"gb")
    # get organism
    annotatDict = ncbiResultSeqRec.annotations
    return(annotatDict['organism'])

它工作得很好,但是每种物种都需要大约2分钟的时间来恢复,这在我看来很慢。我只是想知道我能不能做得更好。我知道我可能会创建NCBI的本地副本来提高性能,我可能会这样做。但是,我怀疑首先查询BLAST,然后使用id查询Entrez不是一种方法。你还有其他的改进建议吗?
谢谢!在


Tags: fromreadgettitlefindseqrecordseqfasta
1条回答
网友
1楼 · 发布于 2024-10-01 17:40:06

你可以通过以下方法获得有机体:

[...]
blastResult = NCBIWWW.qblast("blastn", "nt", seqRecord.seq)
blastRecord = NCBIXML.read(blastResult)

first_organism = blastRecord.descriptions[0]

这至少可以保存efetch查询。无论如何,“blastn”可能需要太长时间,如果你打算大规模查询NCBI,你将被禁止。在

相关问题 更多 >

    热门问题