Fisrt pass on taxonomic rank fetcher

This commit is contained in:
Anton Nekrutenko
2008-03-03 19:40:36 +00:00
parent 135c9a0377
commit 9dca050a79
3 changed files with 143 additions and 0 deletions
+5
View File
@@ -0,0 +1,5 @@
33001686 9443 root Eukaryota Metazoa n n Chordata Craniata Gnathostomata Mammalia n Euarchontoglires Primates n n n n n n n n n n
23236241 9604 root Eukaryota Metazoa n n Chordata Craniata Gnathostomata Mammalia n Euarchontoglires Primates Haplorrhini Hominoidea Hominidae n n n n n n n
12583 9606 root Eukaryota Metazoa n n Chordata Craniata Gnathostomata Mammalia n Euarchontoglires Primates Haplorrhini Hominoidea Hominidae n n n Homo n Homo sapiens n
410771 40674 root Eukaryota Metazoa n n Chordata Craniata Gnathostomata Mammalia n n n n n n n n n n n n n
2286205 63221 root Eukaryota Metazoa n n Chordata Craniata Gnathostomata Mammalia n Euarchontoglires Primates Haplorrhini Hominoidea Hominidae n n n Homo n Homo sapiens Homo sapiens neanderthalensis
+74
View File
@@ -0,0 +1,74 @@
<tool id="Fetch Taxonomic Ranks" name="Fetch Taxonomic Ranks" version="1.0.0">
<description>from a list of GIs</description>
<command interpreter="python">tax.py $input $out_file1 -c $giField</command>
<inputs>
<param format="tabular" name="input" type="data" label="Fetch sequences corresponding to Query"></param>
<param name="giField" label="GIs column" type="data_column" data_ref="input" />
</inputs>
<outputs>
<data format="tabular" name="out_file1" />
</outputs>
<tests>
<test>
<param name="input" value="taxonomyGI.txt"/>
<param name="giField" value="1"/>
<output name="out_file1" file="taxonomyGI.dat"/>
</test>
</tests>
<help>
.. class:: infomark
Use *Filter and Sort->Filter* to restrict output of this tool to desired taxonomic ranks. You can also use *Text Manipulation->Cut* to remove unwanted columns from the output.
------
**What it does**
Fetches taxonomic information for a list of GI numbers (sequiences identifiers used by teh National Center for Biotechnology Information http://www.ncbi.nlm.nih.gov).
-------
**Example**
Suppose you have BLAST output that looks like this::
+-----------------------+----------+----------+-----------------+------------+------+--------+
| queryId | targetGI | identity | alignmentLength | mismatches | gaps | score |
+-----------------------+----------+----------+-----------------+------------+------+--------+
| 1L_EYKX4VC01BXWX1_265 | 1430919 | 90.09 | 212 | 15 | 6 | 252.00 |
| 1L_EYKX4VC01BXWX1_265 | 516142 | 90.09 | 212 | 15 | 6 | 252.00 |
+-----------------------+----------+----------+-----------------+------------+------+--------+
and you want to obtain full taxonomic representation for GIs listed in *targetGI* column. The output of this tool will add 21 columns to the input data shown above. The 21 columns will correspond to::
'root' :1,
'superkingdom':2,
'kingdom' :3,
'subkingdom' :4,
'superphylum' :5,
'phylum' :6,
'subphylum' :7,
'superclass' :8,
'class' :9,
'subclass' :10,
'superorder' :11,
'order' :12,
'suborder' :13,
'superfamily' :14,
'family' :15,
'subfamily' :16,
'tribe' :17,
'subtribe' :18,
'genus' :19,
'subgenus' :20,
'species' :21,
'subspecies' :22
</help>
</tool>
+64
View File
@@ -0,0 +1,64 @@
#!/usr/bin/env python
"""
Identify full taxonomic standing for sequences identified by gi number
usage: %prog gi_list_file out_file
-c, --cols=N: Number of column containing gi within the gi_list file
gi_list_file - user's input containing gi's in the column specified by option -c
taxonomy_db - database containing collapsed NCBI taxonomy generated by prepareTaxonomy.sh script
distributed with Galaxy. See prepareTaxonomy.readme for information on how to generate
necessary files
"""
import pkg_resources
pkg_resources.require( 'bx-python' )
pkg_resources.require( 'pysqlite' )
import traceback
import fileinput
from pysqlite2 import dbapi2 as sqlite
from warnings import warn
from bx.cookbook import doc_optparse
import string, sys
TAXONOMY = '/Users/anton/galaxy/static/taxonomy/taxonomy.db'
def main():
options, args = doc_optparse.parse( __doc__ )
if len(args) < 2:
sys.stderr.write('Not enough arguments\n')
sys.exit(0)
try:
gi_fname, out_fname = args
giCol = int( options.cols ) - 1
except:
doc_optparse.exception
con = sqlite.connect(TAXONOMY)
cur = con.cursor()
fg = open(gi_fname, 'r')
of = open( out_fname, "w" )
try:
for line in fg:
field = string.split(line.rstrip(), '\t')
sqlTemplate = string.Template('select gi2tax.gi, tax.* from gi2tax left join tax on gi2tax.taxId = tax.taxId where gi2tax.gi = $gi')
sql = sqlTemplate.substitute(gi = int(field[giCol]))
cur.execute(sql)
for item in cur.fetchall():
ranks = string.split(item[2], ",")
print >> of, str(item[0]) + "\t" + str(item[1]) + "\t" + "\t".join(ranks)
finally:
fg.close()
of.close()
if __name__ == "__main__":
main()