Skip to content

Instantly share code, notes, and snippets.

Show Gist options
  • Save nickloman/2828015 to your computer and use it in GitHub Desktop.
Save nickloman/2828015 to your computer and use it in GitHub Desktop.
lib/qiime/assign_taxonomy.py.patch
--- lib/qiime/assign_taxonomy.py 2011-09-11 20:53:06.000000000 +0100
+++ ../16s/Qiime-1.2.1/build/lib/qiime/assign_taxonomy.py 2011-02-22 21:00:06.000000000 +0000
@@ -89,7 +89,7 @@
""" Initialize the object
"""
_params = {
- 'Min percent identity': 90,
+ 'Min percent identity': 0.90,
'Max E value': 1e-30,
'Application': 'blastn/megablast'
}
@@ -174,12 +174,9 @@
# if the user provided a result_path, write the
# results to file
of = open(result_path,'w')
- for seq_id, stuff in result.items():
- import sys
- print >>sys.stderr, stuff
- lineage, confidence, blast_hit_id, identity, aln_length = stuff
- of.write('%s\t%s\t%s\t%s\t%s\t%s\n' %
- (seq_id, lineage, confidence, blast_hit_id, identity, aln_length))
+ for seq_id, (lineage, confidence, blast_hit_id) in result.items():
+ of.write('%s\t%s\t%s\t%s\n' %
+ (seq_id, lineage, confidence, blast_hit_id))
of.close()
result = None
logger.info('Result path: %s' % result_path)
@@ -232,11 +229,11 @@
for query_id, hit in hits.items():
query_id=query_id.split()[0]
try:
- hit_id, e_value, idy, aln_length = hit
+ hit_id, e_value = hit
hits[query_id] = \
- (id_to_taxonomy_map.get(hit_id, None),e_value,hit_id,idy, aln_length)
+ (id_to_taxonomy_map.get(hit_id, None),e_value,hit_id)
except TypeError:
- hits[query_id] = ('No blast hit', None, None, None, None)
+ hits[query_id] = ('No blast hit', None, None)
return hits
@@ -250,7 +247,7 @@
blast_result = blast_seqs(\
seqs,Blastall,blast_db=blast_db,\
- params={'-p':'blastn','-n':'T','-a':'8'},\
+ params={'-p':'blastn','-n':'T'},\
add_seq_names=False)
if blast_result['StdOut']:
@@ -262,7 +259,7 @@
for seq_id in seq_ids:
blast_result_id = seq_id.split()[0]
try:
- result[seq_id] = [(e['SUBJECT ID'],float(e['E-VALUE']),e['% IDENTITY'],e['ALIGNMENT LENGTH']) \
+ result[seq_id] = [(e['SUBJECT ID'],float(e['E-VALUE'])) \
for e in blast_result[blast_result_id][0]
if (float(e['E-VALUE']) <= max_evalue and \
float(e['% IDENTITY']) >= min_percent_identity)]
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment