opensource/galaxy
3
0
Fork 0
mirror of https://github.com/galaxyproject/galaxy.git synced 2026-09-24 16:30:27 +08:00
Code Issues Packages Projects Releases Wiki Activity

megablast wrapper tool with functional test data.

Browse Source
This commit is contained in:
Wen-Yu Chung
2008-02-22 16:37:23 +00:00
parent 6f8ba9a02c
commit e43d84131a
2 changed files with 18 additions and 67 deletions
Download Patch File Download Diff File Expand all files Collapse all files
tools/metag_tools/megablast_wrapper.py
+14 -53
View File
@@ -1,53 +1,15 @@
#! /usr/bin/python
"""
run megablast for metagenomics data
database builds file
usage: %prog database/reference query output
Wen-Yu Chung
"""
import sys, os, tempfile
DB_LOC = "/depot/data2/galaxy/blastdb.loc"
def print_bed_format(fields, file_handle):
#print " ".join(fields)
print >> file_handle, "\t".join(fields)
return 0
def print_maf_format(fields, file_handle):
score = fields.pop(0)
print >> file_handle, "a score=%s" %score
print >> file_handle, "s %s" %" ".join(fields)
print >> file_handle
return 0
def parse_megablast_output(filename, output_format, file_handle):
fields = []
file_to_be_parsed = open(filename, 'r')
for i, line in enumerate(file_to_be_parsed):
line = line.strip('\r\n')
if (not line.startswith("#")):
[query_id, subject_id, iden, align_length, mismatches, gaps, q_start, q_end, s_start, s_end, evalue, bit_score] = line.split()
if (int(s_start) > int(s_end)):
strand = "-"
temp = s_start
s_start = s_end
s_end = temp
else:
strand = "+"
if (output_format == "bed"):
fields = [subject_id, s_start, s_end, query_id, str(0), strand]
print_bed_format(fields, file_handle)
else:
fields = [bit_score, subject_id, s_start, str(int(align_length)-int(gaps)), strand, "srcSize", "text"]
print_maf_format(fields, file_handle)
return 0
def __main__():
# file I/O
db_build = sys.argv[1]
query_filename = sys.argv[2]
#output_format = sys.argv[3]
output_filename = sys.argv[3]
# megablast parameters
@@ -63,44 +25,43 @@ def __main__():
db = {}
db_file = open(DB_LOC, "r")
for i, line in enumerate(db_file):
line = line.strip('\r\n')
line = line.rstrip('\r\n')
fields = line.split()
db[(fields[0])] = []
for j in xrange(1, len(fields)):
db[(fields[0])].append(fields[j])
#print "\n".join(db[(db_build)])
# prepare to run megablast
if os.system('which megablast 2>&1'): print >> sys.stderr, "Cannot locate megablast."
for chunk in db[(db_build)]:
#chunk = db[(db_build)][0] # test
#if 1: # test
#if os.path.exists('/usr/bin/megablast'):
megablast_output_file = tempfile.NamedTemporaryFile('w')
megablast_output_filename = megablast_output_file.name
megablast_arguments = ["megablast", "-d", chunk, "-i", query_filename, "-o", megablast_output_filename]
megablast_parameters = ["-m", "8", "-D", "3", "-a", "1"]
megablast_parameters = ["-m", "8", "-D", "3", "-a", "8"]
megablast_user_inputs = ["-W", mega_word_size, "-p", mega_iden_cutoff, "-t", mega_disc_word, "-N", mega_disc_type, "-F", mega_filter]
megablast_command = " ".join(megablast_arguments) + " " + " ".join(megablast_parameters) + " " + " ".join(megablast_user_inputs) + " 2>&1"
#print megablast_command
os.system(megablast_command)
#parse_megablast_output(megablast_output_filename, output_format, output_file)
megablast_output_handle = open(megablast_output_filename, 'r')
for i, line in enumerate(megablast_output_handle):
line = line.strip('\r\n')
for i, line in enumerate( file(megablast_output_filename) ):
line = line.rstrip('\r\n')
fields = line.split()
if (not line.startswith("#")):
# replace subject id with gi number
# remove this after re-build blastdb
subject_id_fields = fields[1].split('|')
gi = subject_id_fields[1]
if len(subject_id_fields) > 1: gi = subject_id_fields[1]
else: gi = subject_id_fields[0]
fields[1] = gi
print >> output_file, "\t".join(fields)
megablast_output_file.close()
output_file.close()
# megablast generates a file called error.log, if empty, delete it, if not, show the contents
if os.path.exists('./error.log'):
for i, line in enumerate( file('./error.log') ):
line = line.rstrip('\r\n')
print >> sys.stdout, line
os.remove('./error.log')
if __name__ == "__main__" : __main__()
tools/metag_tools/megablast_wrapper.xml
+4 -14
View File
@@ -1,19 +1,12 @@
<tool id="megablast_wrapper" name="Run Megablast">
<description>on genomic data</description>
<description>for Metagenomics Projects</description>
<command interpreter="python">megablast_wrapper.py $source_select $input_query $output1 $word_size $iden_cutoff $disc_word $disc_type $filter_query</command>
<inputs>
<param name="input_query" type="data" format="fasta" label="Query Sequence"/>
<param name="source_select" type="select" display="radio" label="Target database">
<!-- <options from_file="/depot/data2/galaxy/blastdb.loc" /> -->
<option value="nt">nt (for nucleotides)</option>
</param>
<!-- <param name="something" type="blastdb" label="database"/> -->
<!--
<param name="output_format" type="select" display="radio" label="Output Type">
<option value="bed">BED</option>
<option value="maf">MAF</option>
</param>
-->
<option value="test">Ecoli_K12</option>
</param>
<param name="word_size" type="select" label="Word size (-W, length of best perfect match)">
<option value="28">28</option>
<option value="16">16</option>
@@ -33,12 +26,10 @@
<outputs>
<data name="output1" format="tabular"/>
</outputs>
<!--
<tests>
<test>
<param name="input_query" value="megablast_query.fa" ftype="fasta"/>
<param name="source_select" value="nt" />
<param name="output_format" value="bed" />
<param name="source_select" value="test" />
<param name="word_size" value="28" />
<param name="iden_cutoff" value="99.0" />
<param name="disc_word" value="0" />
@@ -47,7 +38,6 @@
<output name="output1" file="megablast_test.out"/>
</test>
</tests>
-->
<help>
.. class:: warningmark