lib/Bio/Tools/Run/Phylo/QuickTree.pm

   1 # $Id$
   2 #
   3 # BioPerl module for Bio::Tools::Run::Phylo::QuickTree
   4 #
   5 # Please direct questions and support issues to <bioperl-l@bioperl.org>
   6 #
   7 # Cared for by Sendu Bala <bix@sendu.me.uk>
   8 #
   9 # Copyright Sendu Bala
  10 #
  11 # You may distribute this module under the same terms as perl itself
  12
  13 # POD documentation - main docs before the code
  14
  15 =head1 NAME
  16
  17 Bio::Tools::Run::Phylo::QuickTree - Wrapper for rapid reconstruction of
  18                                     phylogenies using QuickTree
  19
  20 =head1 SYNOPSIS
  21
  22   use Bio::Tools::Run::Phylo::QuickTree;
  23
  24   #  Make a QuickTree factory
  25   @params = ();
  26   $factory = Bio::Tools::Run::Phylo::QuickTree->new(@params);
  27
  28   #  Pass the factory an alignment
  29   $inputfilename = 't/data/cysprot.stockholm';
  30   $tree = $factory->run($inputfilename); # $tree is a Bio::Tree::Tree object.
  31   # or get a Bio::Align::AlignI (SimpleAlign) object from somewhere
  32   $tree = $factory->run($aln);
  33
  34 =head1 DESCRIPTION
  35
  36 This is a wrapper for running the QuickTree application by Kevin Howe. You
  37 can download it here: http://www.sanger.ac.uk/Software/analysis/quicktree/
  38
  39 Currently only input with alignments and output of trees is supported. (Ie.
  40 no support for distance matrix in/out.)
  41
  42 You will need to enable this QuickTree wrapper to find the quicktree program.
  43 This can be done in (at least) three ways:
  44
  45  1. Make sure the QuickTree executable is in your path.
  46  2. Define an environmental variable QUICKTREEDIR which is a
  47     directory which contains the 'quicktree' application:
  48     In bash:
  49
  50     export QUICKTREEDIR=/home/username/quicktree_1.1/bin
  51
  52     In csh/tcsh:
  53
  54     setenv QUICKTREEDIR /home/username/quicktree_1.1/bin
  55
  56  3. Include a definition of an environmental variable QUICKTREEDIR in
  57     every script that will use this QuickTree wrapper module, e.g.:
  58
  59     BEGIN { $ENV{QUICKTREEDIR} = '/home/username/quicktree_1.1/bin' }
  60     use Bio::Tools::Run::Phylo::QuickTree;
  61
  62 =head1 FEEDBACK
  63
  64 =head2 Mailing Lists
  65
  66 User feedback is an integral part of the evolution of this and other
  67 Bioperl modules. Send your comments and suggestions preferably to
  68 the Bioperl mailing list.  Your participation is much appreciated.
  69
  70   bioperl-l@bioperl.org                  - General discussion
  71   http://bioperl.org/wiki/Mailing_lists  - About the mailing lists
  72
  73 =head2 Support
  74
  75 Please direct usage questions or support issues to the mailing list:
  76
  77 I<bioperl-l@bioperl.org>
  78
  79 rather than to the module maintainer directly. Many experienced and
  80 reponsive experts will be able look at the problem and quickly
  81 address it. Please include a thorough description of the problem
  82 with code and data examples if at all possible.
  83
  84 =head2 Reporting Bugs
  85
  86 Report bugs to the Bioperl bug tracking system to help us keep track
  87 of the bugs and their resolution. Bug reports can be submitted via
  88 the web:
  89
  90   http://redmine.open-bio.org/projects/bioperl/
  91
  92 =head1 AUTHOR - Sendu Bala
  93
  94 Email bix@sendu.me.uk
  95
  96 =head1 APPENDIX
  97
  98 The rest of the documentation details each of the object methods.
  99 Internal methods are usually preceded with a _
 100
 101 =cut
 102
 103 package Bio::Tools::Run::Phylo::QuickTree;
 104 use strict;
 105
 106 use Bio::AlignIO;
 107 use Bio::TreeIO;
 108
 109 use base qw(Bio::Tools::Run::WrapperBase);
 110
 111 our $PROGRAM_NAME = 'quicktree';
 112 our $PROGRAM_DIR = $ENV{'QUICKTREEDIR'};
 113
 114 =head2 program_name
 115
 116  Title   : program_name
 117  Usage   : $factory>program_name()
 118  Function: holds the program name
 119  Returns : string
 120  Args    : None
 121
 122 =cut
 123
 124 sub program_name {
 125     return $PROGRAM_NAME;
 126 }
 127
 128 =head2 program_dir
 129
 130  Title   : program_dir
 131  Usage   : $factory->program_dir(@params)
 132  Function: returns the program directory, obtained from ENV variable.
 133  Returns : string
 134  Args    : None
 135
 136 =cut
 137
 138 sub program_dir {
 139     return $PROGRAM_DIR;
 140 }
 141
 142 =head2 new
 143
 144  Title   : new
 145  Usage   : $factory = Bio::Tools::Run::Phylo::QuickTree->new(@params)
 146  Function: creates a new QuickTree factory
 147  Returns : Bio::Tools::Run::Phylo::QuickTree
 148  Args    : Optionally, provide any of the following (default in []):
 149            -upgma  => boolean # Use the UPGMA method to construct the tree [0]
 150            -kimura => boolean # Use the kimura translation for pairwise
 151                               # distances [0]
 152            -boot   => int     # Calculate bootstrap values with n iterations [0]
 153
 154 =cut
 155
 156 sub new {
 157     my ($class, @args) = @_;
 158     my $self = $class->SUPER::new(@args);
 159
 160     # for consistency with other run modules, allow params to be dashless
 161     my %args = @args;
 162     while (my ($key, $val) = each %args) {
 163         if ($key !~ /^-/) {
 164             delete $args{$key};
 165             $args{'-'.$key} = $val;
 166         }
 167     }
 168
 169     my ($upgma, $kimura, $boot) = $self->_rearrange([qw(UPGMA
 170                                                         KIMURA
 171                                                         BOOT)], %args);
 172
 173     $self->upgma(1) if $upgma;
 174     $self->kimura(1) if $kimura;
 175     $self->boot($boot) if $boot;
 176
 177     return $self;
 178 }
 179
 180 =head2 upgma
 181
 182  Title   : upgma
 183  Usage   : $factory->upgma(1);
 184  Function: Choose to use the UPGMA method to construct the tree.
 185  Returns : boolean (default 0)
 186  Args    : None to get, boolean to set.
 187
 188 =cut
 189
 190 sub upgma {
 191     my ($self, $bool) = @_;
 192     if (defined ($bool)) {
 193         $self->{upgma} = $bool;
 194     }
 195     return $self->{upgma} || 0;
 196 }
 197
 198 =head2 kimura
 199
 200  Title   : kimura
 201  Usage   : $factory->kimura(1);
 202  Function: Choose to use the kimura translation for pairwise distances.
 203  Returns : boolean (default 0)
 204  Args    : None to get, boolean to set.
 205
 206 =cut
 207
 208 sub kimura {
 209     my ($self, $bool) = @_;
 210     if (defined ($bool)) {
 211         $self->{kimura} = $bool;
 212     }
 213     return $self->{kimura} || 0;
 214 }
 215
 216 =head2 boot
 217
 218  Title   : boot
 219  Usage   : $factory->boot(100);
 220  Function: Choose to calculate bootstrap values with the supplied number of
 221            iterations.
 222  Returns : int (default 0)
 223  Args    : None to get, int to set.
 224
 225 =cut
 226
 227 sub boot {
 228     my ($self, $int) = @_;
 229     if (defined ($int)) {
 230         $self->{boot} = $int;
 231     }
 232     return $self->{boot} || 0;
 233 }
 234
 235 =head2 run
 236
 237  Title   : run
 238  Usage   : $factory->run($stockholm_file);
 239            $factory->run($align_object);
 240  Function: Runs QuickTree to generate a tree
 241  Returns : Bio::Tree::Tree object
 242  Args    : file name for your input alignment in stockholm format, OR
 243            Bio::Align::AlignI compliant object (eg. Bio::SimpleAlign).
 244
 245 =cut
 246
 247 sub run {
 248     my ($self, $in) = @_;
 249
 250     if (ref $in && $in->isa("Bio::Align::AlignI")) {
 251         $in = $self->_writeAlignFile($in);
 252     }
 253     elsif (! -e $in) {
 254         $self->throw("When not supplying a Bio::Align::AlignI object, you must supply a readable filename");
 255     }
 256
 257     return $self->_run($in);
 258 }
 259
 260 sub _run {
 261     my ($self, $file)= @_;
 262
 263     my $exe = $self->executable || return;
 264     my $param_str = $self->arguments." ".$self->_setparams;
 265     my $command = $exe." $param_str ".$file;
 266
 267     $self->debug("QuickTree command = $command");
 268
 269     open(my $result, "$command |") || $self->throw("QuickTree call ($command) crashed: $?");
 270     my $treeio = Bio::TreeIO->new(-format => 'nhx', -fh => $result);
 271     my $tree = $treeio->next_tree;
 272     close($result);
 273
 274     # if bootstraps were enabled, the bootstraps are the ids; convert to
 275     # bootstrap and no id
 276     if ($self->boot) {
 277         my @nodes = $tree->get_nodes;
 278         my %non_internal = map { $_ => 1 } ($tree->get_leaf_nodes, $tree->get_root_node);
 279         foreach my $node (@nodes) {
 280             next if exists $non_internal{$node};
 281             $node->bootstrap && next; # protect ourselves incase the parser improves
 282             $node->bootstrap($node->id);
 283             $node->id('');
 284         }
 285     }
 286
 287     return $tree;
 288 }
 289
 290 =head2 _setparams
 291
 292  Title   : _setparams
 293  Usage   : Internal function, not to be called directly
 294  Function: Creates a string of params to be used in the command string
 295  Returns : string of params
 296  Args    : none
 297
 298 =cut
 299
 300 sub _setparams {
 301     my $self = shift;
 302     my $param_string = '-in a -out t';
 303
 304     $param_string .= ' -upgma' if $self->upgma;
 305     $param_string .= ' -kimura' if $self->kimura;
 306     $param_string .= ' -boot '.$self->boot if $self->boot;
 307
 308     return $param_string;
 309 }
 310
 311 =head2 _writeAlignFile
 312
 313  Title   : _writeAlignFile
 314  Usage   : obj->_writeAlignFile($seq)
 315  Function: Internal(not to be used directly)
 316  Returns : filename
 317  Args    : Bio::Align::AlignI
 318
 319 =cut
 320
 321 sub _writeAlignFile{
 322     my ($self, $align) = @_;
 323
 324     my ($tfh, $tempfile) = $self->io->tempfile(-dir=>$self->tempdir);
 325
 326     my $out = Bio::AlignIO->new('-fh'     => $tfh,
 327                                 '-format' => 'stockholm');
 328     $out->write_aln($align);
 329     $out->close();
 330     $out = undef;
 331     close($tfh);
 332     undef $tfh;
 333     return $tempfile;
 334 }
 335
 336 1;