From f64f71e7ff2b4523bfbd8ad7257a977028d3f99d Mon Sep 17 00:00:00 2001 From: Jeremy Goecks Date: Mon, 30 Apr 2012 11:41:04 -0400 Subject: [PATCH] Revert Tophat wrapper due to parameter incompatibility issues. --- tools/ngs_rna/tophat_wrapper.py | 25 +- tools/ngs_rna/tophat_wrapper.xml | 564 ++++++++++++++++++++----------- 2 files changed, 378 insertions(+), 211 deletions(-) diff --git a/tools/ngs_rna/tophat_wrapper.py b/tools/ngs_rna/tophat_wrapper.py index 1ac5179479c..4bc91f923ef 100644 --- a/tools/ngs_rna/tophat_wrapper.py +++ b/tools/ngs_rna/tophat_wrapper.py @@ -52,16 +52,16 @@ def __main__(): supplied GFF file. (ignored without -G)") parser.add_option( '', '--no-novel-indels', action="store_true", dest='no_novel_indels', help="Skip indel search. Indel search is enabled by default.") # Types of search. + parser.add_option( '', '--microexon-search', action="store_true", dest='microexon_search', help='With this option, the pipeline will attempt to find alignments incident to microexons. Works only for reads 50bp or longer.') parser.add_option( '', '--closure-search', action="store_true", dest='closure_search', help='Enables the mate pair closure-based search for junctions. Closure-based search should only be used when the expected inner distance between mates is small (<= 50bp)') parser.add_option( '', '--no-closure-search', action="store_false", dest='closure_search' ) - parser.add_option( '', '--min-closure-exon', dest='min_closure_exon', help='Minimum length for exonic hops in potential splice graph' ) - parser.add_option( '', '--min-closure-intron', dest='min_closure_intron', help='Minimum intron length that may be found during closure search' ) - parser.add_option( '', '--max-closure-intron', dest='max_closure_intron', help='Maximum intron length that may be found during closure search' ) - parser.add_option( '', '--microexon-search', action="store_true", dest='microexon_search', help='With this option, the pipeline will attempt to find alignments incident to microexons. Works only for reads 50bp or longer.') parser.add_option( '', '--coverage-search', action="store_true", dest='coverage_search', help='Enables the coverage based search for junctions. Use when coverage search is disabled by default (such as for reads 75bp or longer), for maximum sensitivity.') parser.add_option( '', '--no-coverage-search', action="store_false", dest='coverage_search' ) parser.add_option( '', '--min-segment-intron', dest='min_segment_intron', help='Minimum intron length that may be found during split-segment search' ) parser.add_option( '', '--max-segment-intron', dest='max_segment_intron', help='Maximum intron length that may be found during split-segment search' ) + parser.add_option( '', '--min-closure-exon', dest='min_closure_exon', help='Minimum length for exonic hops in potential splice graph' ) + parser.add_option( '', '--min-closure-intron', dest='min_closure_intron', help='Minimum intron length that may be found during closure search' ) + parser.add_option( '', '--max-closure-intron', dest='max_closure_intron', help='Maximum intron length that may be found during closure search' ) parser.add_option( '', '--min-coverage-intron', dest='min_coverage_intron', help='Minimum intron length that may be found during coverage search' ) parser.add_option( '', '--max-coverage-intron', dest='max_coverage_intron', help='Maximum intron length that may be found during coverage search' ) @@ -73,6 +73,21 @@ def __main__(): (options, args) = parser.parse_args() + # output version # of tool + try: + tmp = tempfile.NamedTemporaryFile().name + tmp_stdout = open( tmp, 'wb' ) + proc = subprocess.Popen( args='tophat -v', shell=True, stdout=tmp_stdout ) + tmp_stdout.close() + returncode = proc.wait() + stdout = open( tmp_stdout.name, 'rb' ).readline().strip() + if stdout: + sys.stdout.write( '%s\n' % stdout ) + else: + raise Exception + except: + sys.stdout.write( 'Could not determine Tophat version\n' ) + # Color or base space space = '' if options.color_space: @@ -166,7 +181,7 @@ def __main__(): opts += ' --no-closure-search' if options.microexon_search: opts += ' --microexon-search' - if options.single_paired == 'paired' and options.mate_std_dev: + if options.single_paired == 'paired': opts += ' --mate-std-dev %s' % options.mate_std_dev if options.initial_read_mismatches: opts += ' --initial-read-mismatches %d' % int( options.initial_read_mismatches ) diff --git a/tools/ngs_rna/tophat_wrapper.xml b/tools/ngs_rna/tophat_wrapper.xml index bc7a494cc0f..d28c2518b76 100644 --- a/tools/ngs_rna/tophat_wrapper.xml +++ b/tools/ngs_rna/tophat_wrapper.xml @@ -7,111 +7,150 @@ tophat_wrapper.py - - ## Change this to accommodate the number of threads you have available. - --num-threads="4" + ## Change this to accommodate the number of threads you have available. + --num-threads="4" - ## Provide outputs. - --junctions-output=$junctions - --hits-output=$accepted_hits + ## Provide outputs. + --junctions-output=$junctions + --hits-output=$accepted_hits - ## Handle reference file. - #if $refGenomeSource.genomeSource == "history": - --own-file=$refGenomeSource.ownFile - #else: - --indexes-path="${refGenomeSource.index.fields.path}" - #end if - - ## Are reads single-end or paired? - --single-paired=$singlePaired.sPaired - - ## First input file always required. - --input1=$input1 - - ## Second input only if input is paired-end. - #if $singlePaired.sPaired == "paired" - --input2=$singlePaired.input2 - -r $singlePaired.mate_inner_distance - --mate-std-dev=$singlePaired.mate_std_dev - #end if - - ## Set params. - --settings=$params.settingsType - #if $params.settingsType == "full": - -a $params.anchor_length - -m $params.splice_mismatches - -i $params.min_intron_length - -I $params.max_intron_length - -g $params.max_multihits - --min-segment-intron $params.min_segment_intron - --max-segment-intron $params.max_segment_intron - --initial-read-mismatches=$params.initial_read_mismatches - --seg-mismatches=$params.seg_mismatches - --seg-length=$params.seg_length - --library-type=$params.library_type - - ## Closure search. - #if $params.closure_search.use_search == "Yes": - --closure-search - --min-closure-exon $params.closure_search.min_closure_exon - --min-closure-intron $params.closure_search.min_closure_intron - --max-closure-intron $params.closure_search.max_closure_intron + ## Handle reference file. + #if $refGenomeSource.genomeSource == "history": + --own-file=$refGenomeSource.ownFile #else: - --no-closure-search - #end if - - ## Indel search. - #if $params.indel_search.allow_indel_search == "Yes": - ## --allow-indels - --max-insertion-length $params.indel_search.max_insertion_length - --max-deletion-length $params.indel_search.max_deletion_length - #else: - --no-novel-indels + --indexes-path="${refGenomeSource.index.fields.path}" #end if - ## Supplying junctions parameters. - #if $params.own_junctions.use_junctions == "Yes": - #if $params.own_junctions.gene_model_ann.use_annotations == "Yes": - -G $params.own_junctions.gene_model_ann.gene_annotation_model + ## Are reads single-end or paired? + --single-paired=$singlePaired.sPaired + + ## First input file always required. + --input1=$input1 + + ## Set params based on whether reads are single-end or paired. + #if $singlePaired.sPaired == "single": + --settings=$singlePaired.sParams.sSettingsType + #if $singlePaired.sParams.sSettingsType == "full": + -a $singlePaired.sParams.anchor_length + -m $singlePaired.sParams.splice_mismatches + -i $singlePaired.sParams.min_intron_length + -I $singlePaired.sParams.max_intron_length + -g $singlePaired.sParams.max_multihits + --min-segment-intron $singlePaired.sParams.min_segment_intron + --max-segment-intron $singlePaired.sParams.max_segment_intron + --initial-read-mismatches=$singlePaired.sParams.initial_read_mismatches + --seg-mismatches=$singlePaired.sParams.seg_mismatches + --seg-length=$singlePaired.sParams.seg_length + --library-type=$singlePaired.sParams.library_type + + ## Indel search. + #if $singlePaired.sParams.indel_search.allow_indel_search == "Yes": + ## --allow-indels + --max-insertion-length $singlePaired.sParams.indel_search.max_insertion_length + --max-deletion-length $singlePaired.sParams.indel_search.max_deletion_length + #else: + --no-novel-indels + #end if + + ## Supplying junctions parameters. + #if $singlePaired.sParams.own_junctions.use_junctions == "Yes": + #if $singlePaired.sParams.own_junctions.gene_model_ann.use_annotations == "Yes": + -G $singlePaired.sParams.own_junctions.gene_model_ann.gene_annotation_model + #end if + #if $singlePaired.sParams.own_junctions.raw_juncs.use_juncs == "Yes": + -j $singlePaired.sParams.own_junctions.raw_juncs.raw_juncs + #end if + ## TODO: No idea why a string cast is necessary, but it is: + #if str($singlePaired.sParams.own_junctions.no_novel_juncs) == "Yes": + --no-novel-juncs + #end if + #end if + + #if $singlePaired.sParams.closure_search.use_search == "Yes": + --closure-search + --min-closure-exon $singlePaired.sParams.closure_search.min_closure_exon + --min-closure-intron $singlePaired.sParams.closure_search.min_closure_intron + --max-closure-intron $singlePaired.sParams.closure_search.max_closure_intron + #else: + --no-closure-search + #end if + #if $singlePaired.sParams.coverage_search.use_search == "Yes": + --coverage-search + --min-coverage-intron $singlePaired.sParams.coverage_search.min_coverage_intron + --max-coverage-intron $singlePaired.sParams.coverage_search.max_coverage_intron + #else: + --no-coverage-search + #end if + ## TODO: No idea why the type conversion is necessary, but it seems to be. + #if str($singlePaired.sParams.microexon_search) == "Yes": + --microexon-search + #end if #end if - #if $params.own_junctions.raw_juncs.use_juncs == "Yes": - -j $params.own_junctions.raw_juncs.raw_juncs - #end if - ## TODO: No idea why a string cast is necessary, but it is: - #if str($params.own_junctions.no_novel_juncs) == "Yes": - --no-novel-juncs + #else: + --input2=$singlePaired.input2 + -r $singlePaired.mate_inner_distance + --settings=$singlePaired.pParams.pSettingsType + #if $singlePaired.pParams.pSettingsType == "full": + --mate-std-dev=$singlePaired.pParams.mate_std_dev + -a $singlePaired.pParams.anchor_length + -m $singlePaired.pParams.splice_mismatches + -i $singlePaired.pParams.min_intron_length + -I $singlePaired.pParams.max_intron_length + -g $singlePaired.pParams.max_multihits + --min-segment-intron $singlePaired.pParams.min_segment_intron + --max-segment-intron $singlePaired.pParams.max_segment_intron + --initial-read-mismatches=$singlePaired.pParams.initial_read_mismatches + --seg-mismatches=$singlePaired.pParams.seg_mismatches + --seg-length=$singlePaired.pParams.seg_length + --library-type=$singlePaired.pParams.library_type + + ## Indel search. + #if $singlePaired.pParams.indel_search.allow_indel_search == "Yes": + ## --allow-indels + --max-insertion-length $singlePaired.pParams.indel_search.max_insertion_length + --max-deletion-length $singlePaired.pParams.indel_search.max_deletion_length + #else: + --no-novel-indels + #end if + + ## Supplying junctions parameters. + #if $singlePaired.pParams.own_junctions.use_junctions == "Yes": + #if $singlePaired.pParams.own_junctions.gene_model_ann.use_annotations == "Yes": + -G $singlePaired.pParams.own_junctions.gene_model_ann.gene_annotation_model + #end if + #if $singlePaired.pParams.own_junctions.raw_juncs.use_juncs == "Yes": + -j $singlePaired.pParams.own_junctions.raw_juncs.raw_juncs + #end if + ## TODO: No idea why type cast is necessary, but it is: + #if str($singlePaired.pParams.own_junctions.no_novel_juncs) == "Yes": + --no-novel-juncs + #end if + #end if + + #if $singlePaired.pParams.closure_search.use_search == "Yes": + --closure-search + --min-closure-exon $singlePaired.pParams.closure_search.min_closure_exon + --min-closure-intron $singlePaired.pParams.closure_search.min_closure_intron + --max-closure-intron $singlePaired.pParams.closure_search.max_closure_intron + #else: + --no-closure-search + #end if + #if $singlePaired.pParams.coverage_search.use_search == "Yes": + --coverage-search + --min-coverage-intron $singlePaired.pParams.coverage_search.min_coverage_intron + --max-coverage-intron $singlePaired.pParams.coverage_search.max_coverage_intron + #else: + --no-coverage-search + #end if + ## TODO: No idea why the type conversion is necessary, but it seems to be. + #if str ($singlePaired.pParams.microexon_search) == "Yes": + --microexon-search + #end if #end if #end if - - #if $params.coverage_search.use_search == "Yes": - --coverage-search - --min-coverage-intron $params.coverage_search.min_coverage_intron - --max-coverage-intron $params.coverage_search.max_coverage_intron - #else: - --no-coverage-search - #end if - ## TODO: No idea why the type conversion is necessary, but it seems to be. - #if str($params.microexon_search) == "Yes": - --microexon-search - #end if - #end if - - - - - - - - - - - - - - - + @@ -129,109 +168,221 @@ - - - - + + + + - - - - - - - - - - - - - - - - + + + + + + + + + + + + + - - - - - - - alignments (number of reads divided by average depth of coverage)" help="0.0 to 1.0 (0 to turn off)" /> - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + +alignments (number of reads divided by average depth of coverage)" help="0.0 to 1.0 (0 to turn off)" /> + + + + + + + + + + - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + @@ -320,11 +471,11 @@ tophat -o tmp_dir -p 1 tophat_in1 test-data/tophat_in2.fastqsanger Rename the files in tmp_dir appropriately --> - - + + @@ -335,13 +486,13 @@ tophat -o tmp_dir -p 1 -r 20 tophat_in1 test-data/tophat_in2.fastqsanger test-data/tophat_in3.fastqsanger Rename the files in tmp_dir appropriately --> - - + + - + @@ -349,15 +500,15 @@ - - + + @@ -395,13 +546,13 @@ Replace the + with double-dash Rename the files in tmp_dir appropriately --> - - + + - + @@ -493,15 +644,16 @@ This is a list of implemented Tophat options:: -j/--raw-juncs [juncs file] Supply TopHat with a list of raw junctions. Junctions are specified one per line, in a tab-delimited format. Records look like: [chrom] [left] [right] [+/-], left and right are zero-based coordinates, and specify the last character of the left sequenced to be spliced to the first character of the right sequence, inclusive. -no-novel-juncs Only look for junctions indicated in the supplied GFF file. (ignored without -G) --no-closure-search Disables the mate pair closure-based search for junctions. Currently, has no effect - closure search is off by default. - --closure-search Enables the mate pair closure-based search for junctions. Closure-based search should only be used when the - --min-closure-exon During closure search for paired end reads, exonic hops in the potential splice graph must be at least this long. The default is 50. - --min-closure-intron The minimum intron length that may be found during closure search. The default is 50. - --max-closure-intron The maximum intron length that may be found during closure search. The default is 5000. expected inner distance between mates is small (about or less than 50bp) + --closure-search Enables the mate pair closure-based search for junctions. Closure-based search should only be used when the expected inner distance between mates is small (about or less than 50bp) --no-coverage-search Disables the coverage based search for junctions. --coverage-search Enables the coverage based search for junctions. Use when coverage search is disabled by default (such as for reads 75bp or longer), for maximum sensitivity. --microexon-search With this option, the pipeline will attempt to find alignments incident to microexons. Works only for reads 50bp or longer. + --butterfly-search TopHat will use a slower but potentially more sensitive algorithm to find junctions in addition to its standard search. Consider using this if you expect that your experiment produced a lot of reads from pre-mRNA, that fall within the introns of your transcripts. --segment-mismatches Read segments are mapped independently, allowing up to this many mismatches in each segment alignment. The default is 2. --segment-length Each read is cut up into segments, each at least this long. These segments are mapped independently. The default is 25. + --min-closure-exon During closure search for paired end reads, exonic hops in the potential splice graph must be at least this long. The default is 50. + --min-closure-intron The minimum intron length that may be found during closure search. The default is 50. + --max-closure-intron The maximum intron length that may be found during closure search. The default is 5000. --min-coverage-intron The minimum intron length that may be found during coverage search. The default is 50. --max-coverage-intron The maximum intron length that may be found during coverage search. The default is 20000. --min-segment-intron The minimum intron length that may be found during split-segment search. The default is 50.