BLAT wrapper for short read analysis.

include functional test data.
This commit is contained in:
Wen-Yu Chung
2008-02-22 16:10:22 +00:00
parent 0d79c5ee56
commit 6f8ba9a02c
2 changed files with 152 additions and 0 deletions
+83
View File
@@ -0,0 +1,83 @@
#! /usr/bin/python
import os, sys, tempfile
nib_file = "/depot/data2/galaxy/alignseq.loc"
twobit_file = "/depot/data2/galaxy/twobit.loc"
def check_nib_file( dbkey ):
nib_path = ''
nibs = {}
for line in open( nib_file ):
line = line.rstrip( '\r\n' )
if line and not line.startswith( "#" ):
fields = line.split()
if len( fields ) < 3:
continue
if ( fields[0] == 'seq' ):
nibs[( fields[1] )] = fields[2]
if nibs.has_key( dbkey ):
nib_path = nibs[( dbkey )]
return nib_path
def check_twobit_file( dbkey ):
twobit_path = ''
twobits = {}
for line in open( twobit_file ):
line = line.rstrip( '\r\n' )
if line and not line.startswith( "#" ):
fields = line.split()
if len( fields ) < 2:
continue
twobits[( fields[0] )] = fields[1]
if twobits.has_key( dbkey ):
twobit_path = twobits[( dbkey )]
return twobit_path
def __main__():
# I/O
source_format = sys.argv[1] # 0: dbkey; 1: upload file
target_file = sys.argv[2]
query_file = sys.argv[3]
output_file = sys.argv[4]
min_iden = sys.argv[5]
tile_size = sys.argv[6]
one_off = sys.argv[7]
all_files = []
if (source_format == '0'):
# check target genome
dbkey = target_file
if dbkey == '?':
print >> sys.stdout, "No genome build specified. please check your dataset."
sys.exit()
nib_path = check_nib_file( dbkey )
twobit_path = check_twobit_file( dbkey )
if not os.path.exists( nib_path ) and not os.path.exists( twobit_path ):
print >> sys.stdout, "No sequences are available for %s, request them by reporting this error." % dbkey
sys.exit()
# check the query file, see whether all of them are legitimate sequence
if (nib_path):
compress_files = os.listdir(nib_path)
target_path = nib_path
elif (twobit_path):
compress_files = twobit_path
target_path = ""
else:
print >> sys.stdout, "Requested genome build has no available sequence."
sys.exit()
for file in compress_files:
file = target_path + '/' + file
file = os.path.normpath(file)
all_files.append(file)
else:
all_files = [target_file]
for detail_file_path in all_files:
output_tempfile = tempfile.NamedTemporaryFile().name
command = 'blat ' + detail_file_path + ' ' + query_file + ' ' + output_tempfile + ' -oneOff=' + one_off + ' -tileSize=' + tile_size + ' -minIdentity=' + min_iden + ' -mask=lower -noHead -out=pslx 2>&1'
os.system(command)
os.system('cat %s >> %s' %(output_tempfile, output_file))
if __name__ == '__main__': __main__()
+69
View File
@@ -0,0 +1,69 @@
<tool id="blat_wrapper" name="Run Sequence Alignment" version="1.0.0">
<description>on short reads against genome builds</description>
<command interpreter="python">
#if $source.source_select=="database":#blat_wrapper.py 0 $source.dbkey $input_query $output1 $iden $tile_size $one_off
#else:#blat_wrapper.py 1 $source.input_target $input_query $output1 $iden $tile_size $one_off
#end if
</command>
<inputs>
<conditional name="source">
<param name="source_select" type="select" label="Target Source">
<option value="database">Genome Build</option>
<option value="input_ref">Your Upload File</option>
</param>
<when value="database">
<param name="dbkey" type="genomebuild" label="Genome" />
</when>
<when value="input_ref">
<param name="input_target" type="data" format="fasta" label="Reference Sequence" />
</when>
</conditional>
<param name="input_query" type="data" format="fasta" label="Query Sequence"/>
<param name="iden" type="float" size="15" value="90.0" label="Minimal Identity (-minIdentity)" />
<param name="tile_size" type="integer" size="15" value="11" label="Minimal Size of Match (-tileSize)" help="for shorter reads, use size = 8"/>
<param name="one_off" type="integer" size="15" value="0" label="Number of Mismatch in the Word (-oneOff)" />
</inputs>
<outputs>
<data name="output1" format="tabular"/>
</outputs>
<tests>
<test>
<param name="source_select" value="database" />
<param name="dbkey" value="eschColi_K12_1" />
<param name="input_query" value="blat_wrapper_test1.fa" ftype="fasta"/>
<param name="iden" value="90.0" />
<param name="tile_size" value="11" />
<param name="one_off" value="0" />
<output name="output1" file="blat_wrapper_test1.out" />
</test>
</tests>
<help>
**What it does**
This tool runs sequence alignment program on your short reads dataset against a genome build.
-----
**Example**
Input as multiple fasta file::
&gt;seq1
TGGTAATGGTGGTTTTTTTTTTTTTTTTTTATTTTT
Search against ce2 (C. elegans March 2004)
(This is output type = pslx)::
25 1 0 0 0 0 0 0 + seq1 36 10 36 chrI 15080483 9704438 9704464 1 26, 10, 9704438, ggttttttttttttttttttattttt, ggtttttttttttttttttttttttt,
27 0 0 0 0 0 1 32 + seq1 36 9 36 chrI 15080483 1302536 1302595 2 21,6, 9,30, 1302536,1302589, tggtttttttttttttttttt,attttt, tggtttttttttttttttttt,attttt,
-----
.. class:: infomark
BLAT: Kent, W James, BLAT--the BLAST-like alignment tool. (2002) Genome Research:12(4) 656-664.
</help>
</tool>