Added SOLiD (colorspace) support to Bowtie wrapper tool

This commit is contained in:
Kelly Vincent
2010-01-11 17:04:45 -05:00
parent 093a9700ed
commit b5fa090092
4 changed files with 1365 additions and 522 deletions
+31
View File
@@ -0,0 +1,31 @@
#This is a sample file distributed with Galaxy that enables tools
#to use a directory of colorspace Bowtie indexed sequences data files.
#You will need to create these data files and then create a
#bowtie_indices_color.loc file similar to this one (store it in this
#directory) that points to the directories in which those files are
#stored. The bowtie_indices_color.loc file has this format (white
#space characters are TAB characters):
#
#<build> <file_base>
#
#So, for example, if you had hg18 indexed stored in
#/depot/data2/galaxy/bowtie/hg18/,
#then the bowtie_indices_color.loc entry would look like this:
#
#hg18 /depot/data2/galaxy/bowtie/hg18/hg18
#
#and your /depot/data2/galaxy/bowtie/hg18/ directory
#would contain hg18.*.ebwt files:
#
#-rw-r--r-- 1 james universe 830134 2005-09-13 10:12 hg18.1.ebwt
#-rw-r--r-- 1 james universe 527388 2005-09-13 10:12 hg18.2.ebwt
#-rw-r--r-- 1 james universe 269808 2005-09-13 10:12 gh18.3.ebwt
#...etc...
#
#Your bowtie_indices_color.loc file should include an entry per line
#for each index set you have stored. The "file" in the path does not
#actually exist, but it is the prefix for the actual index files. For
#example:
#
#hg18 /depot/data2/galaxy/bowtie/hg18/hg18
#hg19 /depot/data2/galaxy/bowtie/hg19/hg19
+311 -136
View File
@@ -2,6 +2,61 @@
"""
Runs Bowtie on single-end or paired-end data.
For use with Bowtie v. 0.12.1
usage: bowtie_wrapper.py [options]
-t, --threads=t: The number of threads to run
-i, --input1=i: The (forward or single-end) reads file in Sanger FASTQ format
-I, --input2=I: The reverse reads file in Sanger FASTQ format
-o, --output=o: The output file
-4, --dataType=4: The type of data (SOLiD or Solexa)
-2, --paired=2: Whether the data is single- or paired-end
-g, --genomeSource=g: The type of reference provided
-r, --ref=r: The reference genome to use or index
-s, --skip=s: Skip the first n reads
-a, --alignLimit=a: Only align the first n reads
-T, --trimH=T: Trim n bases from high-quality (left) end of each read before alignment
-L, --trimL=L: Trim n bases from low-quality (right) end of each read before alignment
-m, --mismatchSeed=m: Maximum number of mismatches permitted in the seed
-M, --mismatchQual=M: Maximum permitted total of quality values at mismatched read positions
-l, --seedLen=l: Seed length
-n, --rounding=n: Whether or not to round to the nearest 10 and saturating at 30
-P, --maqSoapAlign=P: Choose MAQ- or SOAP-like alignment policy
-w, --tryHard=: Whether or not to try as hard as possible to find valid alignments when they exist
-v, --valAlign=v: Report up to n valid arguments per read
-V, --allValAligns=V: Whether or not to report all valid alignments per read
-G, --suppressAlign=G: Suppress all alignments for a read if more than n reportable alignments exist
-b, --best=b: Whether or not to make Bowtie guarantee that reported singleton alignments are 'best' in terms of stratum and in terms of the quality values at the mismatched positions
-B, --maxBacktracks=B: Maximum number of backtracks permitted when aligning a read
-R, --strata=R: Whether or not to report only those alignments that fall in the best stratum if many valid alignments exist and are reportable
-j, --minInsert=j: Minimum insert size for valid paired-end alignments
-J, --maxInsert=J: Maximum insert size for valid paired-end alignments
-O, --mateOrient=O: The upstream/downstream mate orientation for valid paired-end alignment against the forward reference strand
-A, --maxAlignAttempt=A: Maximum number of attempts Bowtie will make to match an alignment for one mate with an alignment for the opposite mate
-f, --forwardAlign=f: Whether or not to attempt to align the forward reference strand
-E, --reverseAlign=E: Whether or not to attempt to align the reverse-complement reference strand
-F, --offrate=F: Override the offrate of the index to n
-8, --snpphred=8: SNP penalty on Phred scale
-6, --snpfrac=6: Fraction of sites expected to be SNP sites
-7, --keepends=7: Keep extreme-end nucleotides and qualities
-S, --seed=S: Seed for pseudo-random number generator
-d, --dbkey=d: Dbkey of reference genome
-C, --params=C: Whether to use default or specified parameters
-u, --iauto_b=u: Automatic or specified behavior
-K, --ipacked=K: Whether or not to use a packed representation for DNA strings
-Q, --ibmax=Q: Maximum number of suffixes allowed in a block
-Y, --ibmaxdivn=Y: Maximum number of suffixes allowed in a block as a fraction of the length of the reference
-D, --idcv=D: The period for the difference-cover sample
-U, --inodc=U: Whether or not to disable the use of the difference-cover sample
-y, --inoref=y: Whether or not to build the part of the reference index used only in paired-end alignment
-z, --ioffrate=z: How many rows get marked during annotation of some or all of the Burrows-Wheeler rows
-W, --iftab=W: The size of the lookup table used to calculate an initial Burrows-Wheeler range with respect to the first n characters of the query
-X, --intoa=X: Whether or not to convert Ns in the reference sequence to As
-N, --iendian=N: Endianness to use when serializing integers to the index file
-Z, --iseed=Z: Seed for the pseudorandom number generator
-c, --icutoff=c: Number of first bases of the reference sequence to index
-x, --indexSettings=x: Whether or not indexing options are to be set
-H, --suppressHeader=H: Suppress header
"""
import optparse, os, shutil, sys, tempfile
@@ -13,161 +68,281 @@ def stop_err( msg ):
def __main__():
#Parse Command Line
parser = optparse.OptionParser()
parser.add_option('', '--threads', dest='threads', help='The number of threads to run')
parser.add_option('', '--input1', dest='input1', help='The (forward or single-end) reads file in Sanger FASTQ format')
parser.add_option('', '--input2', dest='input2', help='The reverse reads file in Sanger FASTQ format')
parser.add_option('', '--output', dest='output', help='The output file')
parser.add_option('', '--paired', dest='paired', help='Whether the data is single- or paired-end')
parser.add_option('', '--genomeSource', dest='genomeSource', help='The type of reference provided')
parser.add_option('', '--ref', dest='ref', help='The reference genome to use or index')
parser.add_option('', '--skip', dest='skip', help='Skip the first n reads')
parser.add_option('', '--alignLimit', dest='alignLimit', help='Only align the first n reads')
parser.add_option('', '--trimH', dest='trimH', help='Trim n bases from high-quality (left) end of each read before alignment')
parser.add_option('', '--trimL', dest='trimL', help='Trim n bases from low-quality (right) end of each read before alignment')
parser.add_option('', '--mismatchSeed', dest='mismatchSeed', help='Maximum number of mismatches permitted in the seed')
parser.add_option('', '--mismatchQual', dest='mismatchQual', help='Maximum permitted total of quality values at mismatched read positions')
parser.add_option('', '--seedLen', dest='seedLen', help='Seed length')
parser.add_option('', '--rounding', dest='rounding', help='Whether or not to round to the nearest 10 and saturating at 30')
parser.add_option('', '--maqSoapAlign', dest='maqSoapAlign', help='Choose MAQ- or SOAP-like alignment policy')
parser.add_option('', '--tryHard', dest='tryHard', help='Whether or not to try as hard as possible to find valid alignments when they exist')
parser.add_option('', '--valAlign', dest='valAlign', help='Report up to n valid arguments per read')
parser.add_option('', '--allValAligns', dest='allValAligns', help='Whether or not to report all valid alignments per read')
parser.add_option('', '--suppressAlign', dest='suppressAlign', help='Suppress all alignments for a read if more than n reportable alignments exist')
parser.add_option('', '--best', dest='best', help="Whether or not to make Bowtie guarantee that reported singleton alignments are 'best' in terms of stratum and in terms of the quality values at the mismatched positions")
parser.add_option('', '--maxBacktracks', dest='maxBacktracks', help='Maximum number of backtracks permitted when aligning a read')
parser.add_option('', '--strata', dest='strata', help='Whether or not to report only those alignments that fall in the best stratum if many valid alignments exist and are reportable')
parser.add_option('', '--minInsert', dest='minInsert', help='Minimum insert size for valid paired-end alignments')
parser.add_option('', '--maxInsert', dest='maxInsert', help='Maximum insert size for valid paired-end alignments')
parser.add_option('', '--mateOrient', dest='mateOrient', help='The upstream/downstream mate orientation for valid paired-end alignment against the forward reference strand')
parser.add_option('', '--maxAlignAttempt', dest='maxAlignAttempt', help='Maximum number of attempts Bowtie will make to match an alignment for one mate with an alignment for the opposite mate')
parser.add_option('', '--forwardAlign', dest='forwardAlign', help='Whether or not to attempt to align the forward reference strand')
parser.add_option('', '--reverseAlign', dest='reverseAlign', help='Whether or not to attempt to align the reverse-complement reference strand')
parser.add_option('', '--offrate', dest='offrate', help='Override the offrate of the index to n')
parser.add_option('', '--seed', dest='seed', help='Seed for pseudo-random number generator')
parser.add_option('', '--dbkey', dest='dbkey', help='')
parser.add_option('', '--params', dest='params', help='Whether to use default or specified parameters')
parser.add_option('', '--iauto_b', dest='iauto_b', help='Automatic or specified behavior')
parser.add_option('', '--ipacked', dest='ipacked', help='Whether or not to use a packed representation for DNA strings')
parser.add_option('', '--ibmax', dest='ibmax', help='Maximum number of suffixes allowed in a block')
parser.add_option('', '--ibmaxdivn', dest='ibmaxdivn', help='Maximum number of suffixes allowed in a block as a fraction of the length of the reference')
parser.add_option('', '--idcv', dest='idcv', help='The period for the difference-cover sample')
parser.add_option('', '--inodc', dest='inodc', help='Whether or not to disable the use of the difference-cover sample')
parser.add_option('', '--inoref', dest='inoref', help='Whether or not to build the part of the reference index used only in paried-end alignment')
parser.add_option('', '--ioffrate', dest='ioffrate', help='How many rows get marked during annotation of some or all of the Burrows-Wheeler rows')
parser.add_option('', '--iftab', dest='iftab', help='The size of the lookup table used to calculate an initial Burrows-Wheeler range with respect to the first n characters of the query')
parser.add_option('', '--intoa', dest='intoa', help='Whether or not to convert Ns in the reference sequence to As')
parser.add_option('', '--iendian', dest='iendian', help='Endianness to use when serializing integers to the index file')
parser.add_option('', '--iseed', dest='iseed', help='Seed for the pseudorandom number generator')
parser.add_option('', '--icutoff', dest='icutoff', help='Number of first bases of the reference sequence to index')
parser.add_option('', '--indexSettings', dest='index_settings', help='Whether or not indexing options are to be set')
parser.add_option('', '--suppressHeader', dest='suppressHeader', help='Suppress header')
parser.add_option( '-t', '--threads', dest='threads', help='The number of threads to run' )
parser.add_option( '-4', '--dataType', dest='dataType', help='The type of data (SOLiD or Solexa)' )
parser.add_option( '-i', '--input1', dest='input1', help='The (forward or single-end) reads file in Sanger FASTQ format' )
parser.add_option( '-I', '--input2', dest='input2', help='The reverse reads file in Sanger FASTQ format' )
parser.add_option( '-o', '--output', dest='output', help='The output file' )
parser.add_option( '-2', '--paired', dest='paired', help='Whether the data is single- or paired-end' )
parser.add_option( '-g', '--genomeSource', dest='genomeSource', help='The type of reference provided' )
parser.add_option( '-r', '--ref', dest='ref', help='The reference genome to use or index' )
parser.add_option( '-s', '--skip', dest='skip', help='Skip the first n reads' )
parser.add_option( '-a', '--alignLimit', dest='alignLimit', help='Only align the first n reads' )
parser.add_option( '-T', '--trimH', dest='trimH', help='Trim n bases from high-quality (left) end of each read before alignment' )
parser.add_option( '-L', '--trimL', dest='trimL', help='Trim n bases from low-quality (right) end of each read before alignment' )
parser.add_option( '-m', '--mismatchSeed', dest='mismatchSeed', help='Maximum number of mismatches permitted in the seed' )
parser.add_option( '-M', '--mismatchQual', dest='mismatchQual', help='Maximum permitted total of quality values at mismatched read positions' )
parser.add_option( '-l', '--seedLen', dest='seedLen', help='Seed length' )
parser.add_option( '-n', '--rounding', dest='rounding', help='Whether or not to round to the nearest 10 and saturating at 30' )
parser.add_option( '-P', '--maqSoapAlign', dest='maqSoapAlign', help='Choose MAQ- or SOAP-like alignment policy' )
parser.add_option( '-w', '--tryHard', dest='tryHard', help='Whether or not to try as hard as possible to find valid alignments when they exist' )
parser.add_option( '-v', '--valAlign', dest='valAlign', help='Report up to n valid arguments per read' )
parser.add_option( '-V', '--allValAligns', dest='allValAligns', help='Whether or not to report all valid alignments per read' )
parser.add_option( '-G', '--suppressAlign', dest='suppressAlign', help='Suppress all alignments for a read if more than n reportable alignments exist' )
parser.add_option( '-b', '--best', dest='best', help="Whether or not to make Bowtie guarantee that reported singleton alignments are 'best' in terms of stratum and in terms of the quality values at the mismatched positions" )
parser.add_option( '-B', '--maxBacktracks', dest='maxBacktracks', help='Maximum number of backtracks permitted when aligning a read' )
parser.add_option( '-R', '--strata', dest='strata', help='Whether or not to report only those alignments that fall in the best stratum if many valid alignments exist and are reportable' )
parser.add_option( '-j', '--minInsert', dest='minInsert', help='Minimum insert size for valid paired-end alignments' )
parser.add_option( '-J', '--maxInsert', dest='maxInsert', help='Maximum insert size for valid paired-end alignments' )
parser.add_option( '-O', '--mateOrient', dest='mateOrient', help='The upstream/downstream mate orientation for valid paired-end alignment against the forward reference strand' )
parser.add_option( '-A', '--maxAlignAttempt', dest='maxAlignAttempt', help='Maximum number of attempts Bowtie will make to match an alignment for one mate with an alignment for the opposite mate' )
parser.add_option( '-f', '--forwardAlign', dest='forwardAlign', help='Whether or not to attempt to align the forward reference strand' )
parser.add_option( '-E', '--reverseAlign', dest='reverseAlign', help='Whether or not to attempt to align the reverse-complement reference strand' )
parser.add_option( '-F', '--offrate', dest='offrate', help='Override the offrate of the index to n' )
parser.add_option( '-S', '--seed', dest='seed', help='Seed for pseudo-random number generator' )
parser.add_option( '-8', '--snpphred', dest='snpphred', help='SNP penalty on Phred scale' )
parser.add_option( '-6', '--snpfrac', dest='snpfrac', help='Fraction of sites expected to be SNP sites' )
parser.add_option( '-7', '--keepends', dest='keepends', help='Keep extreme-end nucleotides and qualities' )
parser.add_option( '-d', '--dbkey', dest='dbkey', help='Dbkey of reference genome' )
parser.add_option( '-C', '--params', dest='params', help='Whether to use default or specified parameters' )
parser.add_option( '-u', '--iauto_b', dest='iauto_b', help='Automatic or specified behavior' )
parser.add_option( '-K', '--ipacked', dest='ipacked', help='Whether or not to use a packed representation for DNA strings' )
parser.add_option( '-Q', '--ibmax', dest='ibmax', help='Maximum number of suffixes allowed in a block' )
parser.add_option( '-Y', '--ibmaxdivn', dest='ibmaxdivn', help='Maximum number of suffixes allowed in a block as a fraction of the length of the reference' )
parser.add_option( '-D', '--idcv', dest='idcv', help='The period for the difference-cover sample' )
parser.add_option( '-U', '--inodc', dest='inodc', help='Whether or not to disable the use of the difference-cover sample' )
parser.add_option( '-y', '--inoref', dest='inoref', help='Whether or not to build the part of the reference index used only in paired-end alignment' )
parser.add_option( '-z', '--ioffrate', dest='ioffrate', help='How many rows get marked during annotation of some or all of the Burrows-Wheeler rows' )
parser.add_option( '-W', '--iftab', dest='iftab', help='The size of the lookup table used to calculate an initial Burrows-Wheeler range with respect to the first n characters of the query' )
parser.add_option( '-X', '--intoa', dest='intoa', help='Whether or not to convert Ns in the reference sequence to As' )
parser.add_option( '-N', '--iendian', dest='iendian', help='Endianness to use when serializing integers to the index file' )
parser.add_option( '-Z', '--iseed', dest='iseed', help='Seed for the pseudorandom number generator' )
parser.add_option( '-c', '--icutoff', dest='icutoff', help='Number of first bases of the reference sequence to index' )
parser.add_option( '-x', '--indexSettings', dest='index_settings', help='Whether or not indexing options are to be set' )
parser.add_option( '-H', '--suppressHeader', dest='suppressHeader', help='Suppress header' )
(options, args) = parser.parse_args()
# make temp directory for placement of indices and copy reference file there if necessary
tmp_index_dir = tempfile.mkdtemp()
# get type of data (solid or solexa)
if options.dataType == 'solid':
colorspace = '-C'
else:
colorspace = ''
# index if necessary
if options.genomeSource == 'history':
if options.genomeSource == 'cHistory' or options.genomeSource == 'xHistory':
# set up commands
if options.index_settings =='index_pre_set':
indexing_cmds = ''
if options.index_settings =='cIndexPreSet' or options.index_settings == 'xIndexPreSet':
indexing_cmds = '%s' % colorspace
else:
try:
indexing_cmds = '%s %s %s %s %s %s %s --offrate %s %s %s %s %s %s' % \
(('','--noauto')[options.iauto_b=='set'],
('','--packed')[options.ipacked=='packed'],
('','--bmax %s'%options.ibmax)[options.ibmax!='None' and options.ibmax>=1],
('','--bmaxdivn %s'%options.ibmaxdivn)[options.ibmaxdivn!='None'],
('','--dcv %s'%options.idcv)[options.idcv!='None'],
('','--nodc')[options.inodc=='nodc'],
('','--noref')[options.inoref=='noref'], options.ioffrate,
('','--ftabchars %s'%options.iftab)[int(options.iftab)>=0],
('','--ntoa')[options.intoa=='yes'],
('--little','--big')[options.iendian=='big'],
('','--seed %s'%options.iseed)[int(options.iseed)>0],
('','--cutoff %s'%options.icutoff)[int(options.icutoff)>0])
if options.iauto_b == 'set':
iauto_b = '--noauto'
else:
iauto_b = ''
if options.ipacked == 'packed':
ipacked = '--packed'
else:
ipacked = ''
if options.ibmax != 'None' and int( options.ibmax ) >= 1:
ibmax = '--bmax %s' % options.ibmax
else:
ibmax = ''
if options.ibmaxdivn != 'None' and int( options.ibmaxdivn ) >= 0:
ibmaxdivn = '--bmaxdivn %s' % options.ibmaxdivn
else:
ibmaxdivn = ''
if options.idcv != 'None' and int( options.idcv ) > 0:
idcv = '--dcv %s' % options.idcv
else:
idcv = ''
if options.inodc == 'nodc':
inodc = '--nodc'
else:
inodc = ''
if options.inoref == 'noref':
inoref = '--noref'
else:
inoref = ''
if options.iftab != 'None' and int( options.iftab ) >= 0:
iftab = '--ftabchars %s' % options.iftab
else:
iftab = ''
if options.intoa == 'yes':
intoa = '--ntoa'
else:
intoa = ''
if options.iendian == 'big':
iendian = '--big'
else:
iendian = '--little'
if int( options.iseed ) > 0:
iseed = '--seed %s' % options.iseed
else:
iseed = ''
if int( options.icutoff ) > 0:
icutoff = '--cutoff %s' % options.icutoff
else:
icutoff = ''
indexing_cmds = '%s %s %s %s %s %s %s --offrate %s %s %s %s %s %s %s' % \
( iauto_b, ipacked, ibmax, ibmaxdivn, idcv, inodc,
inoref, options.ioffrate, iftab, intoa, iendian,
iseed, icutoff, colorspace )
except ValueError:
indexing_cmds = ''
try:
shutil.copy(options.ref, tmp_index_dir)
except Exception, erf:
stop_err('Error creating temp directory for indexing purposes\n' + str(erf))
options.ref = os.path.join(tmp_index_dir,os.path.split(options.ref)[1])
cmd1 = 'bowtie-build %s -f %s %s 2> /dev/null' % (indexing_cmds, options.ref, options.ref)
shutil.copy( options.ref, tmp_index_dir )
except Exception, e:
stop_err( 'Error creating temp directory for indexing purposes\n' + str( e ) )
options.ref = os.path.join( tmp_index_dir, os.path.split( options.ref )[1] )
cmd1 = 'bowtie-build %s -f %s %s 2> /dev/null' % ( indexing_cmds, options.ref, options.ref )
try:
os.chdir(tmp_index_dir)
os.system(cmd1)
except Exception, erf:
stop_err('Error indexing reference sequence\n' + str(erf))
os.chdir( tmp_index_dir )
os.system( cmd1 )
except Exception, e:
stop_err( 'Error indexing reference sequence\n' + str( e ) )
# set up aligning and generate aligning command options
# automatically set threads in both cases
if options.params == 'pre_set':
aligning_cmds = '-p %s -S' % options.threads
if options.suppressHeader == 'true':
suppressHeader = '--sam-nohead'
else:
suppressHeader = ''
if options.params == 'csPreSet' or options.params == 'cpPreSet' or \
options.params == 'xsPreSet' or options.params == 'xpPreSet':
aligning_cmds = '-p %s -S %s -q %s ' % ( options.threads, suppressHeader, colorspace )
else:
try:
aligning_cmds = '%s %s %s %s %s %s %s %s %s %s %s %s %s %s ' \
'%s %s %s %s %s %s %s %s %s %s -p %s -S' % \
(('','-s %s'%options.skip)[options.skip!='None'],
('','-u %s'%options.alignLimit)[int(options.alignLimit)>0],
('','-5 %s'%options.trimH)[int(options.trimH)>=0],
('','-3 %s'%options.trimL)[int(options.trimL)>=0],
('','-n %s'%options.mismatchSeed)[options.mismatchSeed=='0' or options.mismatchSeed=='1' or options.mismatchSeed=='2' or options.mismatchSeed=='3'],
('','-e %s'%options.mismatchQual)[int(options.mismatchQual)>=0],
('','-l %s'%options.seedLen)[int(options.seedLen)>=5],
('','--nomaqround')[options.rounding=='noRound'],
('','-v %s'%options.maqSoapAlign)[options.maqSoapAlign!='-1'],
('','-I %s'%options.minInsert)[options.minInsert!='None'],
('','-X %s'%options.maxInsert)[options.maxInsert!='None'],
('','--%s'%options.mateOrient)[options.mateOrient!='None'],
('','--pairtries %s'%options.maxAlignAttempt)[options.maxAlignAttempt!='None' and int(options.maxAlignAttempt)>=0],
('','--nofw')[options.forwardAlign=='noForward'],
('','--norc')[options.reverseAlign=='noReverse'],
('','--maxbts %s'%options.maxBacktracks)[options.maxBacktracks!='None' and (options.mismatchSeed=='2' or options.mismatchSeed=='3')],
('','-y')[options.tryHard=='doTryHard'],
('','-k %s'%options.valAlign)[options.valAlign!='None' and int(options.valAlign)>=0],
('','-a')[options.allValAligns=='doAllValAligns' and int(options.allValAligns)>=0],
('','-m %s'%options.suppressAlign)[int(options.suppressAlign)>=0],
('','--best')[options.best=='doBest'],
('','--strata')[options.strata=='doStrata'],
('','-o %s'%options.offrate)[int(options.offrate)>=0],
('','--seed %s'%options.seed)[int(options.seed)>=0],
options.threads)
except ValueError, erf:
stop_err('Something is wrong with the alignment parameters and the alignment could not be run\n' + str(erf))
if options.skip != 'None' and int( options.skip ) > 0:
skip = '-s %s' % options.skip
else:
skip = ''
if int( options.alignLimit ) >= 0:
alignLimit = '-u %s' % options.alignLimit
else:
alignLimit = ''
if int( options.trimH ) > 0:
trimH = '-5 %s' % options.trimH
else:
trimH = ''
if int( options.trimL ) > 0:
trimL = '-3 %s' % options.trimL
else:
trimL = ''
if options.mismatchSeed == '0' or options.mismatchSeed == '1' or options.mismatchSeed == '2' or options.mismatchSeed == '3':
mismatchSeed = '-n %s' % options.mismatchSeed
else:
mismatchSeed = ''
if int( options.mismatchQual ) >= 0:
mismatchQual = '-e %s' % options.mismatchQual
else:
mismatchQual = ''
if int( options.seedLen ) >= 5:
seedLen = '-l %s' % options.seedLen
else:
seedLen = ''
if options.rounding == 'noRound':
rounding = '--nomaqround'
else:
rounding = ''
if options.maqSoapAlign != '-1':
maqSoapAlign = '-v %s' % options.maqSoapAlign
else:
maqSoapAlign = ''
if options.minInsert != 'None' and int( options.minInsert ) > 0:
minInsert = '-I %s' % options.minInsert
else:
minInsert = ''
if options.maxInsert != 'None' and int( options.maxInsert ) > 0:
maxInsert = '-X %s' % options.maxInsert
else:
maxInsert = ''
if options.mateOrient != 'None':
mateOrient = '--%s' % options.mateOrient
else:
mateOrient = ''
if options.maxAlignAttempt != 'None' and int( options.maxAlignAttempt ) >= 0:
maxAlignAttempt = '--pairtries %s' % options.maxAlignAttempt
else:
maxAlignAttempt = ''
if options.forwardAlign == 'noForward':
forwardAlign = '--nofw'
else:
forwardAlign = ''
if options.reverseAlign == 'noReverse':
reverseAlign = '--norc'
else:
reverseAlign = ''
if options.maxBacktracks != 'None' and int( options.maxBacktracks ) > 0 and \
( options.mismatchSeed == '2' or options.mismatchSeed == '3' ):
maxBacktracks = '--maxbts %s' % options.maxBacktracks
else:
maxBacktracks = ''
if options.tryHard == 'doTryHard':
tryHard = '-y'
else:
tryHard = ''
if options.valAlign != 'None' and int( options.valAlign ) >= 0:
valAlign = '-k %s' % options.valAlign
else:
valAlign = ''
if options.allValAligns == 'doAllValAligns':
allValAligns = '-a'
else:
allValAligns = ''
if options.suppressAlign != 'None' and int( options.suppressAlign ) >= 0:
suppressAlign = '-m %s' % options.suppressAlign
else:
suppressAlign = ''
if options.best == 'csDoBest' or options.best == 'cpDoBest' or \
options.best == 'xsDoBest' or options.best == 'xpDoBest':
best = '--best'
else:
best = ''
if options.strata == 'doStrata':
strata = '--strata'
else:
strata = ''
if options.offrate != 'None' and int( options.offrate ) >= 0:
offrate = '-o %s' % options.offrate
else:
offrate = ''
if options.seed != 'None' and int( options.seed ) >= 0:
seed = '--seed %s' % options.seed
else:
seed = ''
if options.snpphred != 'None' and int( options.snpphred ) >= 0:
snpphred = '--snpphred %s' % options.snpphred
else:
snpphred = ''
if options.snpfrac != 'None' and float( options.snpfrac ) >= 0:
snpfrac = '--snpfrac %s' % options.snpfrac
else:
snpfrac = ''
if options.keepends != 'None' and options.keepends == 'doKeepends':
keepends = '--col-keepends'
else:
keepends = ''
aligning_cmds = '%s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s ' \
'%s %s %s %s %s %s %s %s %s %s %s -p %s -S %s -q' % \
( skip, alignLimit, trimH, trimL, mismatchSeed, mismatchQual,
seedLen, rounding, maqSoapAlign, minInsert, maxInsert,
mateOrient, maxAlignAttempt, forwardAlign, reverseAlign,
maxBacktracks, tryHard, valAlign, allValAligns, suppressAlign,
best, strata, offrate, seed, colorspace, snpphred, snpfrac,
keepends, options.threads, suppressHeader )
except ValueError, e:
stop_err( 'Something is wrong with the alignment parameters and the alignment could not be run\n' + str( e ) )
# prepare actual aligning commands
if options.paired == 'paired':
cmd2 = 'bowtie %s %s -1 %s -2 %s > %s 2> /dev/null' % (aligning_cmds, options.ref, options.input1, options.input2, options.output)
if options.paired == 'cPaired' or options.paired == 'xPaired':
cmd2 = 'bowtie %s %s -1 %s -2 %s > %s 2> /dev/null' % ( aligning_cmds, options.ref, options.input1, options.input2, options.output )
else:
cmd2 = 'bowtie %s %s %s > %s 2> /dev/null' % (aligning_cmds, options.ref, options.input1, options.output)
cmd2 = 'bowtie %s %s %s > %s 2> /dev/null' % ( aligning_cmds, options.ref, options.input1, options.output )
# align
try:
os.system(cmd2)
except Exception, erf:
stop_err("Error aligning sequence\n" + str(erf))
# remove header if necessary
if options.suppressHeader == 'true':
tmp_out = tempfile.NamedTemporaryFile()
cmd3 = 'cp %s %s' % (options.output, tmp_out.name)
try:
os.system(cmd3)
except Exception, erf:
stop_err("Error copying output file before removing headers\n" + str(erf))
output = file(tmp_out.name, 'r')
fout = file(options.output, 'w')
header = True
line = output.readline()
while line.strip() != '':
if header:
if line.startswith('@HD') or line.startswith('@SQ') or line.startswith('@RG') or line.startswith('@PG') or line.startswith('@CO'):
pass
else:
header = False
fout.write(line)
else:
fout.write(line)
line = output.readline()
fout.close()
tmp_out.close()
os.system( cmd2 )
except Exception, e:
stop_err( 'Error aligning sequence\n' + str( e ) )
# clean up temp dir
if os.path.exists(tmp_index_dir):
shutil.rmtree(tmp_index_dir)
if os.path.exists( tmp_index_dir ):
shutil.rmtree( tmp_index_dir )
if __name__=="__main__": __main__()
File diff suppressed because it is too large Load Diff
+12 -6
View File
@@ -2,15 +2,21 @@ import os
def exec_before_job(app, inp_data, out_data, param_dict, tool):
try:
refFile = param_dict['refGenomeSource']['indices'].value
try:
refFile = param_dict[ 'solidOrSolexa' ][ 'cRefGenomeSource' ][ 'cIndex' ].value
except:
refFile = param_dict[ 'solidOrSolexa' ][ 'xRefGenomeSource' ][ 'xIndex' ].value
except:
try:
refFile = param_dict['refGenomeSource']['ownFile'].dbkey
try:
refFile = param_dict[ 'solidOrSolexa' ][ 'cRefGenomeSource' ][ 'cOwnFile' ].dbkey
except:
refFile = param_dict[ 'solidOrSolexa' ][ 'xRefGenomeSource' ][ 'xOwnFile' ].dbkey
except:
out_data['output'].set_dbkey('?')
out_data[ 'output' ].set_dbkey( '?' )
return
dbkey = os.path.split(refFile)[1].split('.')[0]
dbkey = os.path.split( refFile )[1].split( '.' )[0]
# deal with the one odd case
if dbkey.find('chrM') >= 0:
if dbkey.find( 'chrM' ) >= 0 or dbkey.find( 'chr_m' ) >= 0:
dbkey = 'equCab2'
out_data['output'].set_dbkey(dbkey)
out_data[ 'output' ].set_dbkey(dbkey)