python - Bio.Entrez EFetch 的替代方案，用于从 NCBI 下载完整基因组序列

Question

我的目标是从 NCBI 下载完整的后生动物基因组序列。我有一个我需要的基因组序列的唯一 ID 号列表。我计划使用 Bio.Entrez 模块 EFetch 下载数据，但今天通过 2011 年 11 月 2 日的发行说明 (http://1.usa.gov/1TA5osg) 了解到 EFetch 不支持“基因组”数据库。任何人都可以建议替代包/模块或其他方法来解决这个问题吗？先感谢您!

Answer 1

这是一个适合您的脚本——尽管您可能需要修改它才能使其工作。将脚本命名为您喜欢的任何名称，但是当您调用脚本时，请按如下所示操作:python name_of_script[带有 .py 扩展名] your_email_address。

您需要将电子邮件添加到通话末尾，否则它将无法正常工作。如果您有一个包含登录号(1/行)的文本文件，则选择选项 2。如果您选择选项 1，它将要求您输入生物体名称、菌株名称和关键字等项目。使用尽可能多的关键字——只需确保用逗号分隔它们即可。如果您选择第一个选项，将搜索 NCBI 并返回 GI 编号 [注意:NCBI 将在 9.2016 中逐步淘汰 GI 编号，因此此脚本在此之后可能无法工作]，然后将其用于获取入藏号。一旦出现所有入藏号，就会创建一个文件夹，并为每个入藏号(命名为入藏号)创建一个子文件夹。在每个子文件夹中，将下载相应的 fasta AND genbank 文件。这些文件将以登录号作为文件名(例如 accession_number.fa、accession_number.gb)。根据您的目的编辑脚本。

另外...请注意脚本的警告 (ACHTUNG) 部分。有时规则可以被改变……但如果你足够过分，你的 IP 可能会被 NCBI 屏蔽。您已收到警告。

import os
import os.path
import sys
import re #regular expressions
from Bio import Entrez
import datetime
import time
import glob

arguments = sys.argv
Entrez.email = arguments[1] #email

accession_ids = []
print('Select method for obtaining the accession numbers?\n')
action = input('1 -- Input Search Terms\n2 -- Use text file\n')

if action == '1':
    print('\nYou will be asked to enter an organism name, a strain name, and keywords.')
    print('It is not necessary to provide a value to each item (you may just hit [ENTER]), but you must provide at least one item.\n')

    organism = input('Enter the organism you wish to search for (e.g. Escherichia coli [ENTER])\n')
    strain = input('Enter the strain you wish to search for. (e.g., HUSEC2011 [ENTER])\n')
    keywords = input('Enter the keywords separated by a comma (e.g., complete genome, contigs, partial [ENTER])\n')
    search_phrase = ''

    if ',' in keywords:
        keywords = keywords.split(',')

    ncbi_terms = ['organism', 'strain', 'keyword']
    ncbi_values = [organism, strain, keywords]

    for index, n in enumerate(ncbi_values):
        if index == 0 and n != '':
            search_phrase = '(' + n + '[' + ncbi_terms[index] + '])'
        else:
            if n != '' and index != len(ncbi_values)-1:
                search_phrase = search_phrase + ' AND (' + n + '[' + ncbi_terms[index] + '])'
        if index == len(ncbi_values)-1 and n != '' and type(n) is not list:
            search_phrase = search_phrase + ' AND (' + n + '[' + ncbi_terms[index] + '])'
        if index == len(ncbi_values)-1 and n != '' and type(n) is list:
            for name in n:
                name = name.lstrip()
                search_phrase = search_phrase + ' AND (' + name + '[' + ncbi_terms[index] + '])'


    print('Here is the complete search line that will be used: \n\n', search_phrase)

    handle = Entrez.esearch(db='nuccore', term=search_phrase, retmax=1000, rettype='acc', retmode='text')
    result = Entrez.read(handle)
    handle.close()
    #print(result['Count'])
    gi_numbers = result['IdList']

    fetch_handle = Entrez.efetch(db='nucleotide', id=result['IdList'], rettype='acc', retmode='text')
    accession_ids = [id.strip() for id in fetch_handle]
    fetch_handle.close()

if action == '2': #use this option if you have a file of accession #s
    file_name = input('Enter the name of the file\n')

    with open(file_name, 'r') as input_file:
        lines = input_file.readlines()

        for line in lines:
            line = line.replace('\n', '')
            accession_ids.append(line)
#--------------------------------------------------------------------------------------------------------------

#----------------------------------- Make directory to store files --------------------------------------------
new_path = 'Genbank_Files/'
if not os.path.exists(new_path):
    os.makedirs(new_path)

print('You have ' + str(len(accession_ids)) + ' file(s) to download.') #print(accession_ids)

ending='.gb'

files = []
##CHECK IF FILE HAS BEEN DOWNLOADED
for dirpath, dirnames, filenames in os.walk(new_path):
    for filename in [f for f in filenames if f.endswith(ending)]: #for zipped files
        files.append(os.path.join(dirpath,filename))


for f in files:
    f = f.rsplit('/')[-1]
    f = f.replace('.gb', '')

    if f in accession_ids:
        ind = accession_ids.index(f)
        accession_ids.pop(ind)

print('')
print('You have ' + str(len(accession_ids)) + ' file(s) to download.')
#--------------------------------------------------------------------------

###############################################################################
#---ACHTUNG--ACHTUNG--ACHTUNG--ACHTUNG--ACHTUNG--ACHTUNG--ACHTUNG--ACHTUNG----#
###############################################################################
# Call Entrez to download files
# If downloading more than 100 files...
# Run this script only between 9pm-5am Monday - Friday EST
# Send E-utilities requests to http://eutils.ncbi.nlm.nih.gov
# Make no more than 3 requests every 1 second (Biopython takes care of this).
# Use URL parameter email & tool for distributed software
# NCBI's Disclaimer and Copyright notice must be evident to users of your service. 
#
# Use this script at your own risk. 
# Neither the script author nor author's employers are responsible for consequences arising from improper usage 
###############################################################################


# CALL ENTREZ: Call Entrez to download genbank AND fasta (nucleotide) files using accession numbers.
###############################################################################
start_day = datetime.date.today().weekday() # 0 is Monday, 6 is Sunday
start_time = datetime.datetime.now().time()
print(str(start_day), str(start_time))
print('')

if ((start_day < 5 and start_time > datetime.time(hour=21)) or (start_day < 5 and start_time < datetime.time(hour=5)) or start_day > 5 or len(accession_ids) <= 100 ):
    print('Calling Entrez...')

    for a in accession_ids:
        if ((datetime.date.today().weekday() < 5 and datetime.datetime.now().time() > datetime.time(hour=21)) or 
            (datetime.date.today().weekday() < 5 and datetime.datetime.now().time() < datetime.time(hour=5)) or 
            (datetime.date.today().weekday() == start_day + 1 and datetime.datetime.now().time() < datetime.time(hour=5)) or 
            (datetime.date.today().weekday() > 5) or len(accession_ids) <= 100 ):

            print('Downloading ' + a)

            new_path = 'Genbank_Files/' + a  + '/'
            if not os.path.exists(new_path):
                os.makedirs(new_path)

            handle=Entrez.efetch(db='nucleotide', id=a, rettype='gb', retmode='text', seq_start=0)
            FILENAME = new_path + a + '.gb'
            local_file=open(FILENAME,'w')
            local_file.write(handle.read())
            handle.close()
            local_file.close()

            handle=Entrez.efetch(db='nucleotide', id=a, rettype='fasta', retmode='text') 
            FILENAME =  new_path + a + '.fna'
            local_file=open(FILENAME,'w')
            local_file.write(handle.read())
            handle.close()
            local_file.close()

else:
    print('You have too many files to download at the time. Try again later.')
#-------

<小时/>