我正在使用以下脚本:
import csv
import time
import sys
from ete3 import NCBITaxa
ncbi = NCBITaxa()
def get_desired_ranks(taxid, desired_ranks):
lineage = ncbi.get_lineage(taxid)
names = ncbi.get_taxid_translator(lineage)
lineage2ranks = ncbi.get_rank(names)
ranks2lineage = dict((rank,taxid) for (taxid, rank) in lineage2ranks.items())
return{'{}_id'.format(rank): ranks2lineage.get(rank, '<not present>') for rank in desired_ranks}
if __name__ == '__main__':
file = open(sys.argv[1], "r")
taxids = []
contigs = []
for line in file:
line = line.split("\n")[0]
taxids.append(line.split(",")[0])
contigs.append(line.split(",")[1])
desired_ranks = ['superkingdom', 'phylum']
results = list()
for taxid in taxids:
results.append(list())
results[-1].append(str(taxid))
ranks = get_desired_ranks(taxid, desired_ranks)
for key, rank in ranks.items():
if rank != '<not present>':
results[-1].append(list(ncbi.get_taxid_translator([rank]).values())[0])
else:
results[-1].append(rank)
i = 0
for result in results:
print(contigs[i] + ','),
print(','.join(result))
i += 1
file.close()
该脚本从文件中获取taxid,并从NCBI的分类数据库的本地副本中获取它们各自的血统。奇怪的是,当我在一小组taxid(~70,~100)上运行这个脚本时,它运行得很好,但是我的大多数数据集都是280k以上的taxid,这些都破坏了脚本。你知道吗
我得到一个完全错误:
Traceback (most recent call last):
File "/data1/lstout/blast/scripts/getLineageByETE3.py", line 31, in <module>
ranks = get_desired_ranks(taxid, desired_ranks)
File "/data1/lstout/blast/scripts/getLineageByETE3.py", line 11, in get_desired_ranks
lineage = ncbi.get_lineage(taxid)
File "/data1/lstout/.local/lib/python2.7/site-packages/ete3/ncbi_taxonomy/ncbiquery.py", line 227, in get_lineage
result = self.db.execute('SELECT track FROM species WHERE taxid=%s' %taxid)
sqlite3.Warning: You can only execute one statement at a time.
回溯的前两个文件只是我上面引用的脚本,第三个文件是ete3的一个文件。正如我所说的,这个脚本可以很好地处理小数据集。你知道吗
我尝试过的:
转到ete3代码的第227行。。。你知道吗
result = self.db.execute('SELECT track FROM species WHERE taxid=%s' %merged_conversion[taxid])
并将“execute”函数更改为“executescript”,以便能够同时处理多个查询(这似乎是问题所在)。这产生了一个新的错误,并导致了一个兔子洞,我改变了他们的脚本中的小东西,试图捏造这个工作。没有结果。这是完整的冒犯功能:
def get_lineage(self, taxid):
"""Given a valid taxid number, return its corresponding lineage track as a
hierarchically sorted list of parent taxids.
"""
if not taxid:
return None
result = self.db.execute('SELECT track FROM species WHERE taxid=%s' %taxid)
raw_track = result.fetchone()
if not raw_track:
#perhaps is an obsolete taxid
_, merged_conversion = self._translate_merged([taxid])
if taxid in merged_conversion:
result = self.db.execute('SELECT track FROM species WHERE taxid=%s' %merged_conversion[taxid])
raw_track = result.fetchone()
# if not raise error
if not raw_track:
#raw_track = ["1"]
raise ValueError("%s taxid not found" %taxid)
else:
warnings.warn("taxid %s was translated into %s" %(taxid, merged_conversion[taxid]))
track = list(map(int, raw_track[0].split(",")))
return list(reversed(track))
让我如此困扰的是,它只处理少量数据!我在学校的高性能计算机上运行这些脚本,并尝试在它们的头部节点和交互式moab调度程序中运行。什么都没用。你知道吗
目前没有回答
相关问题 更多 >
编程相关推荐