File Coverage

blib/lib/Bio/Tools/Run/Phylo/SLR.pm
Criterion Covered Total %
statement 55 213 25.8
branch 14 104 13.4
condition 1 24 4.1
subroutine 15 24 62.5
pod 13 13 100.0
total 98 378 25.9


line stmt bran cond sub pod time code
1             #
2             # BioPerl module for Bio::Tools::Run::Phylo::SLR
3             #
4             # Please direct questions and support issues to
5             #
6             # Cared for by Albert Vilella
7             #
8             # Copyright Albert Vilella
9             #
10             # You may distribute this module under the same terms as perl itself
11              
12             # POD documentation - main docs before the code
13              
14             =head1 NAME
15              
16             Bio::Tools::Run::Phylo::SLR - Wrapper around the SLR program
17              
18             =head1 SYNOPSIS
19              
20             use Bio::Tools::Run::Phylo::SLR;
21             use Bio::AlignIO;
22             use Bio::TreeIO;
23             use Bio::SimpleAlign;
24              
25             my $alignio = Bio::AlignIO->new
26             (-format => 'fasta',
27             -file => 't/data/219877.cdna.fasta');
28              
29             my $aln = $alignio->next_aln;
30              
31             my $treeio = Bio::TreeIO->new
32             (-format => 'newick', -file => 't/data/219877.tree');
33              
34             my $tree = $treeio->next_tree;
35              
36             my $slr = Bio::Tools::Run::Phylo::SLR->new();
37             $slr->alignment($aln);
38             $slr->tree($tree);
39             # $rc = 1 for success, 0 for errors
40             my ($rc,$results) = $slr->run();
41              
42             my $positive_sites = $results->{'positive'};
43              
44             print "# Site\tNeutral\tOptimal\tOmega\t",
45             "lower\tupper\tLRT_Stat\tPval\tAdj.Pval\tResult\tNote\n";
46             foreach my $positive_site (@$positive_sites) {
47             print
48             $positive_site->[0], "\t",
49             $positive_site->[1], "\t",
50             $positive_site->[2], "\t",
51             $positive_site->[3], "\t",
52             $positive_site->[4], "\t",
53             $positive_site->[5], "\t",
54             $positive_site->[6], "\t",
55             $positive_site->[7], "\t",
56             $positive_site->[8], "\t",
57             "positive\n";
58             }
59              
60             =head1 DESCRIPTION
61              
62             This is a wrapper around the SLR program. See
63             http://www.ebi.ac.uk/goldman/SLR/ for more information.
64              
65             This module is more about generating the proper ctl file and
66             will run the program in a separate temporary directory to avoid
67             creating temp files all over the place.
68              
69             =head1 FEEDBACK
70              
71             =head2 Mailing Lists
72              
73             User feedback is an integral part of the evolution of this and other
74             Bioperl modules. Send your comments and suggestions preferably to
75             the Bioperl mailing list. Your participation is much appreciated.
76              
77             bioperl-l@bioperl.org - General discussion
78             http://bioperl.org/wiki/Mailing_lists - About the mailing lists
79              
80             =head2 Support
81              
82             Please direct usage questions or support issues to the mailing list:
83              
84             I
85              
86             rather than to the module maintainer directly. Many experienced and
87             reponsive experts will be able look at the problem and quickly
88             address it. Please include a thorough description of the problem
89             with code and data examples if at all possible.
90              
91             =head2 Reporting Bugs
92              
93             Report bugs to the Bioperl bug tracking system to help us keep track
94             of the bugs and their resolution. Bug reports can be submitted via the
95             web:
96              
97             http://redmine.open-bio.org/projects/bioperl/
98              
99             =head1 AUTHOR - Albert Vilella
100              
101             Email avilella-at-gmail-dot-com
102              
103             =head1 CONTRIBUTORS
104              
105             Additional contributors names and emails here
106              
107             =head1 APPENDIX
108              
109             The rest of the documentation details each of the object methods.
110             Internal methods are usually preceded with a _
111              
112             =cut
113              
114              
115             #' keep my emacs happy
116             # Let the code begin...
117              
118              
119             package Bio::Tools::Run::Phylo::SLR;
120 1     1   136585 use vars qw(@ISA %VALIDVALUES $MINNAMELEN $PROGRAMNAME $PROGRAM);
  1         1  
  1         62  
121 1     1   5 use strict;
  1         1  
  1         16  
122 1     1   4 use Bio::Root::Root;
  1         1  
  1         18  
123 1     1   454 use Bio::AlignIO;
  1         71776  
  1         29  
124 1     1   370 use Bio::TreeIO;
  1         13890  
  1         23  
125 1     1   6 use Bio::SimpleAlign;
  1         1  
  1         16  
126 1     1   395 use Bio::Tools::Run::WrapperBase;
  1         2  
  1         21  
127 1     1   4 use Cwd;
  1         1  
  1         71  
128 1     1   6 use File::Spec;
  1         2  
  1         209  
129              
130             @ISA = qw(Bio::Root::Root Bio::Tools::Run::WrapperBase);
131              
132             =head2 Default Values
133              
134             INCOMPLETE DOCUMENTATION OF ALL METHODS
135              
136             seqfile [incodon]
137             File from which to read alignment of codon sequences. The file
138             should be in PAML format.
139              
140             treefile [intree]
141             File from which tree should be read. The tree should be in Nexus
142             format
143              
144             outfile [slr.res]
145             File to which results are written. If the file already exists, it will
146             be overwritten.
147              
148             reoptimise [1]
149             Should the branch lengths, omega and kappa be reoptimized?
150             0 - no
151             1 - yes.
152              
153             kappa [2.0]
154             Value for kappa. If 'reoptimise' is specified, the value
155             given will be used as am initial estimate,
156              
157             omega [0.1]
158             Value for omega (dN/dS). If 'reoptimise' is specified, the value
159             given will be used as an initial estimate.
160              
161             codonf [0]
162             How codon frequencies are estimated:
163             0: F61/F60 Estimates used are the empirical frequencies from the
164             data.
165             1: F3x4 The frequencies of nucleotides at each codon position
166             are estimated from the data and then multiplied together to get the
167             frequency of observing a given codon. The frequency of stop codons is
168             set to zero, and all other frequencies scaled appropriately.
169             2: F1x4 Nucleotide frequencies are estimated from the data
170             (not taking into account at which position in the codon it occurs).
171             The nucleotide frequencies are multiplied together to get the frequency
172             of observing and then corrected for stop codons.
173              
174             freqtype [0]
175             How codon frequencies are incorporated into the substitution matrix.
176             0: q_{ij} = pi_{j} s_{ij}
177             1: q_{ij} = \sqrt(pi_j/pi_i) s_{ij}
178             2: q_{ij} = \pi_{n} s_{ij}, where n is the nucleotide that the
179             subsitution is to.
180             3: q_{ij} = s_{ij} / pi_i
181             Option 0 is the tradition method of incorporating equilibrium frequencies
182             into subsitution matrices (Felsenstein 1981; Goldman and Yang, 1994)
183             Option 1 is described by Goldman and Whelan (2002), in this case with the
184             additional parameter set to 0.5.
185             Option 2 was suggested by Muse and Gaut (1994).
186             Option 3 is included as an experiment, originally suggested by Bret Larget.
187             it does not appear to describe evolution very successfully and should not
188             be used for analyses.
189              
190             Kosakovsky-Pond has repeatedly stated that he finds incorporating codon
191             frequencies in the manner of option 2 to be superior to option 0. We find
192             that option 1 tends to perform better than either of these options.
193              
194             positive_only [0]
195             If only positively selected sites are of interest, set this to "1".
196             Calculation will be slightly faster, but information about sites under
197             purifying selection is lost.
198              
199             gencode [universal]
200             Which genetic code to use when determining whether a given mutation
201             is synonymous or nonsynonymous. Currently only "universal" and
202             "mammalian" mitochondrial are supported.
203              
204             nucleof [0]
205             Allow for empirical exchangabilities for nucleotide substitution.
206             0: No adjustment. All nucleotides treated the same, modulo
207             transition / transversion.
208             1: The rate at which a substitution caused a mutation from nucleotide
209             a to nucleotide b is adjust by a constant N_{ab}. This adjustment is
210             in addition to other adjustments (e.g. transition / transversion or
211             base frequencies).
212              
213             aminof [0]
214             Incorporate amino acid similarity parameters into substitution matrix,
215             adjusting omega for a change between amino acid i and amino acid j.
216             A_{ij} is a symmetric matrix of constants representing amino acid
217             similarities.
218             0: Constant omega for all amino acid changes
219             1: omega_{ij} = omega^{A_{ij}}
220             2: omega_{ij} = a_{ij} log(omega) / [ 1 - exp(-a_{ij} log(omega)) ]
221             Option 1 has the same form as the original codon subsitution model
222             proposed by Goldman and Yang (but with potentially different
223             constants).
224             Option 2 has a more population genetic derivtion, with omega being
225             interpreted as the ratio of fixation probabilities.
226              
227             nucfile [nuc.dat]
228             If nucleof is non-zero, read nucleotide substitution constants from
229             nucfile. If this file does not exist, hard coded constants are used.
230              
231             aminofile [amino.dat]
232             If aminof is non-zero, read amino acid similarity constants from
233             aminofile. If this file does not exist, hard coded constants are used.
234              
235             timemem [0]
236             Print summary of real time and CPU time used. Will eventually print
237             summary of memory use as well.
238              
239             ldiff [3.841459]
240             Twice log-likelihood difference used as a threshold for calculating
241             support (confidence) intervals for sitewise omega estimates. This
242             value should be the quantile from a chi-square distribution with one
243             degree of freedom corresponding to the support required.
244             E.g. qchisq(0.95,1) = 3.841459
245             0.4549364 = 50% support
246             1.323304 = 75% support
247             2.705543 = 90% support
248             3.841459 = 95% support
249             6.634897 = 99% support
250             7.879439 = 99.5% support
251             10.82757 = 99.9% support
252              
253             paramin []
254             If not blank, read in parameters from file given by the argument.
255              
256             paramout []
257             If not blank, write out parameter estimates to file given.
258              
259             skipsitewise [0]
260             Skip sitewise estimation of omega. Depending on other options given,
261             either calculate maximum likelihood or likelihood fixed at parameter
262             values given.
263              
264             seed [0]
265             Seed for random number generator. If seed is 0, then previously
266             produced seed file (~/.rng64) is used. If this does not exist, the
267             random number generator is initialised using the clock.
268              
269             saveseed [1]
270             If non-zero, save finial seed in file (~/.rng64) to be used as initial
271             seed in future runs of program.
272              
273             =head2 Results Format
274              
275             Results file (default: slr.res)
276             ------------
277             Results are presented in nine columns
278              
279             Site
280             Number of sites in alignment
281              
282             Neutral
283             (minus) Log-probability of observing site given that it was
284             evolving neutrally (omega=1)
285              
286             Optimal
287             (minus) Log-probability of observing site given that it was
288             evolving at the optimal value of omega.
289              
290             Omega
291             The value of omega which maximizes the log-probability of observing
292              
293             LRT_Stat
294             Log-likelihood ratio statistic for non-neutral selection (or
295             positive selection if the positive_only option is set to 1).
296             LRT_Stat = 2 * (Neutral-Optimal)
297              
298             Pval
299             P-value for non-neutral (or positive) selection at a site,
300             unadjusted for multiple comparisons.
301              
302             Adj. Pval
303             P-value for non-neutral (or positive) selection at a site, after
304             adjusting for multiple comparisons using the Hochberg procedure
305             (see the file "MultipleComparisons.txt" in the doc directory).
306              
307             Result
308             A simple visual guide to the result. Sites detected as having been
309             under positive selection are marked with a '+', sites under
310             purifying selection are marked with '-'. The number of symbols
311             Number symbols Threshold
312             1 95%
313             2 99%
314             3 95% after adjustment
315             4 99% after adjustment
316              
317             Occasionally the result may also contain an exclamation mark. This
318             indicates that the observation at a site is not significantly
319             different from random (equivalent to infinitely strong positive
320             selection). This may indicate that the alignment at that site is bad
321              
322             Note
323              
324             The following events are flagged:
325             Synonymous All codons at a site code for the same amino
326             acid.
327             Single character Only one sequence at the site is ungapped,
328             the result of a recent insertion for example.
329             All gaps All sequences at a site contain a gap
330             character.
331              
332             Sites marked "Single character" or "All gaps" are not counted
333             towards the number of sites for the purposes of correcting for
334             multiple comparisons since it is not possible to detect selection
335             from none or one observation under the assumptions made by the
336             sitewise likelihood ratio test.
337              
338             =cut
339              
340              
341             #' keep my emacs happy
342              
343             BEGIN {
344              
345 1     1   2 $MINNAMELEN = 25;
346 1         1 $PROGRAMNAME = 'Slr_Linux_static';
347 1 50       8 if ($^O =~ /darwin/i) {
    50          
348 0         0 $PROGRAMNAME = 'Slr_osx';
349             } elsif ($^O =~ /mswin/i) {
350 0         0 $PROGRAMNAME = 'Slr_windows.exe';
351             }
352 1 50       3 if( defined $ENV{'SLRDIR'} ) {
353 0 0       0 $PROGRAM = Bio::Root::IO->catfile($ENV{'SLRDIR'},$PROGRAMNAME). ($^O =~ /mswin/i ?'_windows.exe':'');;
354             }
355            
356             # valid values for parameters, the default one is always
357             # the first one in the array
358             # example file provided with the package
359             %VALIDVALUES = (
360 1         1594 'outfile' => 'slr.res',
361             'reoptimise' => [ 1,0],
362             'kappa' => '2.0',
363             'omega' => '0.1',
364             'codonf' => [ 0, 1,2],
365             'freqtype' => [ 0, 1,2,3],
366             'positive_only' => [ 0, 1],
367             'gencode' => [ "universal", "mammalian"],
368             'nucleof' => [ 0, 1],
369             'aminof' => [ 0, 1,2],
370             'nucfile' => '',
371             'aminofile' => '',
372             'timemem' => [ 0, 1],
373             'ldiff' => [ 3.841459, 0.4549364,1.323304,2.705543,6.634897,7.879439,10.82757],
374             'paramin' => '',
375             'paramout' => '',
376             'skipsitewise' => [ 0, 1],
377             'seed' => [0],
378             'saveseed' => [ 1, 0]
379             );
380             }
381              
382             =head2 program_name
383              
384             Title : program_name
385             Usage : $factory->program_name()
386             Function: holds the program name
387             Returns: string
388             Args : None
389              
390             =cut
391              
392             sub program_name {
393 6     6 1 26 return $PROGRAMNAME;
394             }
395              
396             =head2 program_dir
397              
398             Title : program_dir
399             Usage : ->program_dir()
400             Function: returns the program directory, obtained from ENV variable.
401             Returns: string
402             Args :
403              
404             =cut
405              
406             sub program_dir {
407 3 50   3 1 14 return Bio::Root::IO->catfile($ENV{SLRDIR}) if $ENV{SLRDIR};
408             }
409              
410              
411             =head2 new
412              
413             Title : new
414             Usage : my $obj = Bio::Tools::Run::Phylo::SLR->new();
415             Function: Builds a new Bio::Tools::Run::Phylo::SLR object
416             Returns : Bio::Tools::Run::Phylo::SLR
417             Args : -alignment => the Bio::Align::AlignI object
418             -save_tempfiles => boolean to save the generated tempfiles and
419             NOT cleanup after onesself (default FALSE)
420             -tree => the Bio::Tree::TreeI object
421             -params => a hashref of SLR parameters (all passed to set_parameter)
422             -executable => where the SLR executable resides
423              
424             See also: L, L
425              
426             =cut
427              
428             sub new {
429 1     1 1 98 my($class,@args) = @_;
430              
431 1         10 my $self = $class->SUPER::new(@args);
432 1         17 my ($aln, $tree, $st, $params, $exe,
433             $ubl) = $self->_rearrange([qw(ALIGNMENT TREE SAVE_TEMPFILES
434             PARAMS EXECUTABLE)],
435             @args);
436 1 50       9 defined $aln && $self->alignment($aln);
437 1 50       2 defined $tree && $self->tree($tree);
438 1 50       2 defined $st && $self->save_tempfiles($st);
439 1 50       2 defined $exe && $self->executable($exe);
440              
441 1         3 $self->set_default_parameters();
442 1 50       3 if( defined $params ) {
443 0 0       0 if( ref($params) !~ /HASH/i ) {
444 0         0 $self->warn("Must provide a valid hash ref for parameter -FLAGS");
445             } else {
446 0         0 map { $self->set_parameter($_, $$params{$_}) } keys %$params;
  0         0  
447             }
448             }
449 1         4 return $self;
450             }
451              
452              
453             =head2 prepare
454              
455             Title : prepare
456             Usage : my $rundir = $slr->prepare($aln);
457             Function: prepare the SLR analysis using the default or updated parameters
458             the alignment parameter must have been set
459             Returns : value of rundir
460             Args : L object,
461             L object
462              
463             =cut
464              
465             sub prepare{
466 0     0 1 0 my ($self,$aln,$tree) = @_;
467 0 0       0 unless ( $self->save_tempfiles ) {
468             # brush so we don't get plaque buildup ;)
469 0         0 $self->cleanup();
470             }
471 0 0       0 $tree = $self->tree unless $tree;
472 0 0       0 $aln = $self->alignment unless $aln;
473 0 0       0 if( ! $aln ) {
474 0         0 $self->warn("must have supplied a valid alignment file in order to run SLR");
475 0         0 return 0;
476             }
477 0 0       0 if( ! $tree ) {
478 0         0 $self->warn("must have supplied a valid tree file in order to run SLR");
479 0         0 return 0;
480             }
481 0         0 my ($tempdir) = $self->tempdir();
482 0         0 my ($tempseqFH,$tempseqfile);
483              
484             # Reorder the alignment according to the tree
485 0         0 my $ct = 1;
486 0         0 my %order;
487 0         0 foreach my $node ($tree->get_leaf_nodes) {
488 0         0 $order{$node->id_output} = $ct++;
489             }
490 0         0 my @seq; my @ids;
491 0         0 foreach my $seq ( $aln->each_seq() ) {
492 0         0 push @seq, $seq;
493 0         0 push @ids, $seq->display_id;
494             }
495             # use the map-sort-map idiom:
496 0         0 my @sorted = map { $_->[1] } sort { $a->[0] <=> $b->[0] } map { [$order{$_->id()}, $_] } @seq;
  0         0  
  0         0  
  0         0  
497 0         0 my $sorted_aln = Bio::SimpleAlign->new();
498 0         0 foreach (@sorted) {
499 0         0 $sorted_aln->add_seq($_);
500             }
501              
502             # Rename the leaf nodes in the tree from 1 to n
503 0         0 $ct = 1;
504 0         0 foreach my $node ($tree->get_leaf_nodes) {
505 0         0 $node->id($ct++);
506             }
507              
508 0 0       0 ($tempseqFH,$tempseqfile) = $self->io->tempfile
509             ('-dir' => $tempdir,
510             UNLINK => ($self->save_tempfiles ? 0 : 1));
511 0 0       0 my $alnout = Bio::AlignIO->new('-format' => 'phylip',
512             '-fh' => $tempseqFH,
513             '-interleaved' => 0,
514             '-idlinebreak' => 1,
515             '-idlength' => $MINNAMELEN > $aln->maxdisplayname_length() ? $MINNAMELEN : $aln->maxdisplayname_length() +1);
516              
517 0         0 $alnout->write_aln($sorted_aln);
518 0         0 $alnout->close();
519 0         0 undef $alnout;
520 0         0 close($tempseqFH);
521              
522 0         0 my ($temptreeFH,$temptreefile);
523 0 0       0 ($temptreeFH,$temptreefile) = $self->io->tempfile
524             ('-dir' => $tempdir,
525             UNLINK => ($self->save_tempfiles ? 0 : 1));
526              
527 0         0 my $treeout = Bio::TreeIO->new('-format' => 'newick',
528             '-fh' => $temptreeFH);
529              
530             # We need to add a line with the num of leaves ($ct-1) and the
531             # num of trees (1)
532 0         0 $treeout->_print(sprintf("%d 1\n",($ct-1)));
533 0         0 $treeout->write_tree($tree);
534 0         0 $treeout->close();
535 0         0 close($temptreeFH);
536              
537             # now let's print the ctl file.
538             # many of the these programs are finicky about what the filename is
539             # and won't even run without the properly named file.
540              
541 0         0 my ($treevolume,$treedirectories,$treefile) = File::Spec->splitpath( $temptreefile );
542 0         0 my ($alnvolume,$alndirectories,$alnfile) = File::Spec->splitpath( $tempseqfile );
543 0         0 my $slr_ctl = "$tempdir/slr.ctl";
544 0 0       0 open(SLR, ">$slr_ctl") or $self->throw("cannot open $slr_ctl for writing");
545 0         0 print SLR "seqfile\: $alnfile\n";
546 0         0 print SLR "treefile\: $treefile\n";
547 0         0 my $outfile = $self->outfile_name;
548 0         0 print SLR "outfile\: $outfile\n";
549              
550 0         0 my %params = $self->get_parameters;
551 0         0 while( my ($param,$val) = each %params ) {
552 0 0       0 next if $param eq 'outfile';
553 0         0 print SLR "$param\: $val\n";
554             }
555 0         0 close(SLR);
556 0         0 return $tempdir;
557             }
558              
559              
560              
561             =head2 run
562              
563             Title : run
564             Usage : my ($rc,$parser) = $slr->run($aln,$tree);
565             Function: run the SLR analysis using the default or updated parameters
566             the alignment parameter must have been set
567             Returns : Return code, L
568             Args : L object,
569             L object
570              
571              
572             =cut
573              
574             sub run {
575 0     0 1 0 my ($self) = shift;;
576 0         0 my $outfile = $self->outfile_name;
577 0         0 my $tmpdir = $self->prepare(@_);
578              
579             #my ($rc,$parser) = (1);
580 0         0 my ($rc,$results) = (1);
581             {
582 0         0 my $cwd = cwd();
  0         0  
583 0         0 my $exit_status;
584 0         0 chdir($tmpdir);
585 0         0 my $slrexe = $self->executable();
586 0 0 0     0 $self->throw("unable to find or run executable for SLR") unless $slrexe && -e $slrexe && -x _;
      0        
587 0         0 my $run;
588 0 0       0 open($run, "$slrexe |") or $self->throw("Cannot open exe $slrexe");
589 0         0 my @output = <$run>;
590 0         0 $exit_status = close($run);
591 0         0 $self->error_string(join('',@output));
592 0 0 0     0 if( (grep { /\berr(or)?: /io } @output) || !$exit_status) {
  0         0  
593 0         0 $self->warn("There was an error - see error_string for the program output");
594 0         0 $rc = 0;
595             }
596 0         0 eval {
597 0 0       0 open RESULTS, "$tmpdir/$outfile" or die "couldnt open results file: $!\n";
598 0         0 my $okay = 0;
599 0         0 my $sites;
600 0         0 my $type = 'default';
601 0         0 while () {
602 0         0 chomp $_;
603 0 0       0 if ( /^\#/ ) {next;}
  0         0  
604 0 0       0 if ( /\!/ ) {$type = 'random';} # random is last
  0 0       0  
    0          
    0          
    0          
    0          
    0          
605 0         0 elsif ( /\+/ ) {$type = 'positive';}
606 0         0 elsif ( /\-\s+/ ) {$type = 'negative';}
607 0         0 elsif ( /Constant/ ) {$type = 'constant';}
608 0         0 elsif ( /All gaps/ ) {$type = 'all_gaps';}
609 0         0 elsif ( /Single character/ ) {$type = 'single_character';}
610 0         0 elsif ( /Synonymous/ ) {$type = 'synonymous';}
611 0         0 else {$type = 'default'}
612 0 0       0 if ( /^\s+(\d+)\s+(\S+)\s+(\S+)\s+(\S+)\s+(\S+)\s+(\S+)\s+(\S+)\s+(\S+)\s+(\S+)/ ) {
613 0         0 push @{$sites->{$type}}, [$1,$2,$3,$4,$5,$6,$7,$8,$9];
  0         0  
614             } else {
615 0         0 $DB::single=1;1;
  0         0  
616             }
617             }
618 0         0 $results = $sites;
619 0         0 close RESULTS;
620             # TODO: we could have a proper parser object
621             # $parser = Bio::Tools::Phylo::SLR->new(-file => "$tmpdir/$outfile",
622             # -dir => "$tmpdir");
623             };
624 0 0       0 if( $@ ) {
625 0         0 $self->warn($self->error_string);
626             }
627 0         0 chdir($cwd);
628             }
629             # return ($rc,$parser);
630 0         0 return ($rc,$results);
631             }
632              
633             =head2 error_string
634              
635             Title : error_string
636             Usage : $obj->error_string($newval)
637             Function: Where the output from the last analysus run is stored.
638             Returns : value of error_string
639             Args : newvalue (optional)
640              
641              
642             =cut
643              
644             sub error_string{
645 0     0 1 0 my ($self,$value) = @_;
646 0 0       0 if( defined $value) {
647 0         0 $self->{'error_string'} = $value;
648             }
649 0         0 return $self->{'error_string'};
650              
651             }
652              
653             =head2 alignment
654              
655             Title : alignment
656             Usage : $slr->align($aln);
657             Function: Get/Set the L object
658             Returns : L object
659             Args : [optional] L
660             Comment : We could potentially add support for running directly on a file
661             but we shall keep it simple
662             See also: L
663              
664             =cut
665              
666             sub alignment{
667 0     0 1 0 my ($self,$aln) = @_;
668              
669 0 0       0 if( defined $aln ) {
670 0 0 0     0 if( -e $aln ) {
    0          
671 0         0 $self->{'_alignment'} = $aln;
672             } elsif( !ref($aln) || ! $aln->isa('Bio::Align::AlignI') ) {
673 0         0 $self->warn("Must specify a valid Bio::Align::AlignI object to the alignment function not $aln");
674 0         0 return undef;
675             } else {
676 0         0 $self->{'_alignment'} = $aln;
677             }
678             }
679 0         0 return $self->{'_alignment'};
680             }
681              
682             =head2 tree
683              
684             Title : tree
685             Usage : $slr->tree($tree, %params);
686             Function: Get/Set the L object
687             Returns : L
688             Args : [optional] $tree => L,
689              
690             Comment : We could potentially add support for running directly on a file
691             but we shall keep it simple
692             See also: L
693              
694             =cut
695              
696             sub tree {
697 0     0 1 0 my ($self, $tree, %params) = @_;
698 0 0       0 if( defined $tree ) {
699 0 0 0     0 if( ! ref($tree) || ! $tree->isa('Bio::Tree::TreeI') ) {
700 0         0 $self->warn("Must specify a valid Bio::Tree::TreeI object to the alignment function");
701             }
702 0         0 $self->{'_tree'} = $tree;
703             }
704 0         0 return $self->{'_tree'};
705             }
706              
707             =head2 get_parameters
708              
709             Title : get_parameters
710             Usage : my %params = $self->get_parameters();
711             Function: returns the list of parameters as a hash
712             Returns : associative array keyed on parameter names
713             Args : none
714              
715              
716             =cut
717              
718             sub get_parameters{
719 0     0 1 0 my ($self) = @_;
720             # we're returning a copy of this
721 0         0 return %{ $self->{'_slrparams'} };
  0         0  
722             }
723              
724              
725             =head2 set_parameter
726              
727             Title : set_parameter
728             Usage : $slr->set_parameter($param,$val);
729             Function: Sets a SLR parameter, will be validated against
730             the valid values as set in the %VALIDVALUES class variable.
731             The checks can be ignored if one turns off param checks like this:
732             $slr->no_param_checks(1)
733             Returns : boolean if set was success, if verbose is set to -1
734             then no warning will be reported
735             Args : $param => name of the parameter
736             $value => value to set the parameter to
737             See also: L
738              
739             =cut
740              
741             sub set_parameter{
742 0     0 1 0 my ($self,$param,$value) = @_;
743 0 0 0     0 unless (defined $self->{'no_param_checks'} && $self->{'no_param_checks'} == 1) {
744 0 0       0 if ( ! defined $VALIDVALUES{$param} ) {
745 0         0 $self->warn("unknown parameter $param will not be set unless you force by setting no_param_checks to true");
746 0         0 return 0;
747             }
748 0 0 0     0 if ( ref( $VALIDVALUES{$param}) =~ /ARRAY/i &&
749 0         0 scalar @{$VALIDVALUES{$param}} > 0 ) {
750            
751 0 0       0 unless ( grep { $value eq $_ } @{ $VALIDVALUES{$param} } ) {
  0         0  
  0         0  
752 0         0 $self->warn("parameter $param specified value $value is not recognized, please see the documentation and the code for this module or set the no_param_checks to a true value");
753 0         0 return 0;
754             }
755             }
756             }
757 0         0 $self->{'_slrparams'}->{$param} = $value;
758 0         0 return 1;
759             }
760              
761             =head2 set_default_parameters
762              
763             Title : set_default_parameters
764             Usage : $slr->set_default_parameters(0);
765             Function: (Re)set the default parameters from the defaults
766             (the first value in each array in the
767             %VALIDVALUES class variable)
768             Returns : none
769             Args : boolean: keep existing parameter values
770              
771              
772             =cut
773              
774             sub set_default_parameters{
775 1     1 1 2 my ($self,$keepold) = @_;
776 1 50       3 $keepold = 0 unless defined $keepold;
777            
778 1         6 while( my ($param,$val) = each %VALIDVALUES ) {
779             # skip if we want to keep old values and it is already set
780 19 50 33     48 next if( defined $self->{'_slrparams'}->{$param} && $keepold);
781 19 100       40 if(ref($val)=~/ARRAY/i ) {
782 12         31 $self->{'_slrparams'}->{$param} = $val->[0];
783             } else {
784 7         19 $self->{'_slrparams'}->{$param} = $val;
785             }
786             }
787             }
788              
789              
790             =head1 Bio::Tools::Run::WrapperBase methods
791              
792             =cut
793              
794             =head2 no_param_checks
795              
796             Title : no_param_checks
797             Usage : $obj->no_param_checks($newval)
798             Function: Boolean flag as to whether or not we should
799             trust the sanity checks for parameter values
800             Returns : value of no_param_checks
801             Args : newvalue (optional)
802              
803              
804             =cut
805              
806             sub no_param_checks{
807 0     0 1 0 my ($self,$value) = @_;
808 0 0       0 if( defined $value) {
809 0         0 $self->{'no_param_checks'} = $value;
810             }
811 0         0 return $self->{'no_param_checks'};
812             }
813              
814              
815             =head2 save_tempfiles
816              
817             Title : save_tempfiles
818             Usage : $obj->save_tempfiles($newval)
819             Function:
820             Returns : value of save_tempfiles
821             Args : newvalue (optional)
822              
823              
824             =cut
825              
826             =head2 outfile_name
827              
828             Title : outfile_name
829             Usage : my $outfile = $slr->outfile_name();
830             Function: Get/Set the name of the output file for this run
831             (if you wanted to do something special)
832             Returns : string
833             Args : [optional] string to set value to
834              
835              
836             =cut
837              
838             sub outfile_name {
839 0     0 1 0 my $self = shift;
840 0 0       0 if( @_ ) {
841 0         0 return $self->{'_slrparams'}->{'outfile'} = shift @_;
842             }
843 0 0       0 unless (defined $self->{'_slrparams'}->{'outfile'}) {
844 0         0 $self->{'_slrparams'}->{'outfile'} = 'out.res';
845             }
846 0         0 return $self->{'_slrparams'}->{'outfile'};
847             }
848              
849             =head2 tempdir
850              
851             Title : tempdir
852             Usage : my $tmpdir = $self->tempdir();
853             Function: Retrieve a temporary directory name (which is created)
854             Returns : string which is the name of the temporary directory
855             Args : none
856              
857              
858             =cut
859              
860             =head2 cleanup
861              
862             Title : cleanup
863             Usage : $slr->cleanup();
864             Function: Will cleanup the tempdir directory after an SLR run
865             Returns : none
866             Args : none
867              
868              
869             =cut
870              
871             =head2 io
872              
873             Title : io
874             Usage : $obj->io($newval)
875             Function: Gets a L object
876             Returns : L
877             Args : none
878              
879              
880             =cut
881              
882             sub DESTROY {
883 1     1   1768 my $self= shift;
884 1 50       8 unless ( $self->save_tempfiles ) {
885 1         10 $self->cleanup();
886             }
887 1         5 $self->SUPER::DESTROY();
888             }
889              
890             1;