Bio/Tools/Phylo/PAML/Result.pm

   1 #
   2 # BioPerl module for Bio::Tools::Phylo::PAML::Result
   3 #
   4 # Please direct questions and support issues to <bioperl-l@bioperl.org>
   5 #
   6 # Cared for by Jason Stajich <jason-at-bioperl.org>
   7 #
   8 # Copyright Jason Stajich, Aaron Mackey
   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::Phylo::PAML::Result - A PAML result set object
  17
  18 =head1 SYNOPSIS
  19
  20   # see Bio::Tools::Phylo::PAML for example usage
  21   use Bio::Tools::Phylo::PAML;
  22   my $parser = Bio::Tools::Phylo::PAML->new
  23     (-file => "./results/mlc", -dir => "./results/");
  24
  25   # get the first/next result; a Bio::Tools::Phylo::PAML::Result object,
  26   # which isa Bio::SeqAnalysisResultI object.
  27   my $result = $parser->next_result();
  28
  29   my @seqs         = $result->get_seqs;
  30   my %input_params = $result->get_input_parameters;
  31   my @basfreq      = $result->get_codon_pos_basefreq;
  32   my $MLmatrix     = $result->get_MLmatrix; # get MaxLikelihood Matrix
  33   my $NGmatrix     = $result->get_NGmatrix; # get Nei-Gojoburi Matrix
  34
  35
  36   # for AAML runs
  37   my $AAmatrix   = $result->get_AADistMatrix;
  38   my $AAMLmatrix   = $result->get_AAMLDistMatrix;
  39
  40   # if -dir contains an rst file get list of
  41   # Bio::PrimarySeq ancestral state reconstructions of the sequences
  42   my @rsts          = $result->get_rst_seqs;
  43
  44
  45   # if you want to print the changes on the tree
  46   # this will print out the
  47   # anc_aa       => ANCESTRAL AMINO ACID
  48   # anc_prob     => ANCESTRAL AA PROBABILITY
  49   # derived_aa   => DERIVED AA
  50   # derived_prob => DERIVE AA PROBABILITY (where appropriate - NA for extant/tip taxas)
  51   # site         => which codon site this in the alignment
  52     @trees = $result->get_rst_trees;
  53     for my $t ( @trees ) {
  54         for my $node ( $t->get_nodes ) {
  55             next unless $node->ancestor; # skip root node
  56             my @changes = $node->get_tag_values('changes');
  57             my $chgstr = '';
  58             for my $c ( @changes ) {
  59                 for my $k ( sort keys %$c ) {
  60                     $chgstr .= "$k => $c->{$k} ";
  61                 }
  62                 $chgstr .= "\n\t";
  63             }
  64
  65             printf "node:%s n=%s s=%s\n\t%s\n",
  66             $node->id,
  67             $node->get_tag_values('n'),
  68             $node->get_tag_values('s'),
  69             $chgstr;
  70         }
  71     }
  72
  73   # Persite probabilities
  74   my $persite = $result->get_rst_persite;
  75   # let's score site 1
  76   $site = $persite->[2];
  77   # so site 2, node 2 (extant node, node 2)
  78   print $site->[2]->{'codon'}, ' ',$site->[2]->{'aa'},"\n";
  79   # site 2, node 3
  80   print $site->[3]->{'codon'}, ' ',$site->[3]->{'aa'}, "\n";
  81
  82   # ancestral node 9, codon, aa, marginal probabilities; Yang95 is listed as
  83   #  (eqn. 4 in Yang et al. 1995 Genetics 141:1641-1650) in PAML rst file.
  84   print $site->[9]->{'codon'}, ' ',$site->[9]->{'aa'}, ' ', $site->[9]->{'prob'}, ' ',
  85         $site->[9]->{'Yang95_aa'},' ', $site->[9]->{'Yang95_aa_prob'},"\n";
  86
  87
  88 =head1 DESCRIPTION
  89
  90 This is a container object for PAML Results.
  91
  92 =head1 FEEDBACK
  93
  94 =head2 Mailing Lists
  95
  96 User feedback is an integral part of the evolution of this and other
  97 Bioperl modules. Send your comments and suggestions preferably to
  98 the Bioperl mailing list.  Your participation is much appreciated.
  99
 100   bioperl-l@bioperl.org                  - General discussion
 101   http://bioperl.org/wiki/Mailing_lists  - About the mailing lists
 102
 103 =head2 Support
 104
 105 Please direct usage questions or support issues to the mailing list:
 106
 107 I<bioperl-l@bioperl.org>
 108
 109 rather than to the module maintainer directly. Many experienced and
 110 reponsive experts will be able look at the problem and quickly
 111 address it. Please include a thorough description of the problem
 112 with code and data examples if at all possible.
 113
 114 =head2 Reporting Bugs
 115
 116 Report bugs to the Bioperl bug tracking system to help us keep track
 117 of the bugs and their resolution. Bug reports can be submitted via
 118 email or the web:
 119
 120   https://github.com/bioperl/bioperl-live/issues
 121
 122 =head1 AUTHOR - Jason Stajich, Aaron Mackey
 123
 124  Email jason-at-bioperl-dot-org
 125  Email amackey-at-virginia-dot-edu
 126
 127 =head1 CONTRIBUTORS
 128
 129 Albert Vilella avilella-AT-gmail-DOT-com
 130
 131 =head1 APPENDIX
 132
 133 The rest of the documentation details each of the object methods.
 134 Internal methods are usually preceded with a _
 135
 136 =cut
 137
 138
 139 # Let the code begin...
 140
 141
 142 package Bio::Tools::Phylo::PAML::Result;
 143 use strict;
 144
 145
 146 use base qw(Bio::Root::Root Bio::AnalysisResultI);
 147
 148 =head2 new
 149
 150  Title   : new
 151  Usage   : my $obj = Bio::Tools::Phylo::PAML::Result->new(%data);
 152  Function: Builds a new Bio::Tools::Phylo::PAML::Result object
 153  Returns : Bio::Tools::Phylo::PAML::Result
 154  Args    : -trees     => array reference of Bio::Tree::TreeI objects
 155            -MLmatrix  => ML matrix
 156            -seqs      => array reference of Bio::PrimarySeqI objects
 157            -codonpos  => array reference of codon positions
 158            -codonfreq => array reference of codon frequencies
 159            -version   => version string
 160            -model     => model string
 161            -patterns  => hashref with the fields '-patterns', '-ns', '-ls'
 162            -stats     => array ref of misc stats   (optional)
 163            -aafreq    => Hashref of AA frequencies (only for AAML)
 164            -aadistmat => Bio::Matrix::PhylipDist   (only for AAML)
 165            -aamldistmat => Bio::Matrix::PhylipDist   (only for pairwise AAML)
 166            -ntfreq    => array ref of NT frequencies (only for BASEML)
 167            -seqfile    => seqfile used
 168            -kappa_mat => Bio::Matrix::PhylipDist of kappa values (only for BASEML)
 169            -alpha_mat => Bio::Matrix::PhylipDist of alpha values (only for BASEML)
 170            -NSSitesresult => arrayref of PAML::ModelResult
 171            -input_params  => input params from .ctl file
 172            -rst       => array reference of Bio::PrimarySeqI objects
 173                          of ancestral state reconstruction
 174            -rst_persite=> arrayref of persite data, this is a complicated set of AoH
 175            -rst_trees  => rst trees with changes coded on the tree
 176
 177 See Also: L<Bio::Tree::TreeI>, L<Bio::PrimarySeqI>, L<Bio::Matrix::PhylipDist>, L<Bio::Tools::Phylo::PAML>
 178
 179
 180 =cut
 181
 182 sub new {
 183   my($class,@args) = @_;
 184
 185   my $self = $class->SUPER::new(@args);
 186   my ($trees,$mlmat,$seqs,$ngmatrix,
 187       $codonpos,$codonfreq,$version,
 188       $model,$patterns, $stats,
 189       $aafreq, $aadistmat,
 190       $aamldistmat,
 191       $ntfreqs, $seqfile, $kappa_mat, $alpha_mat,
 192       $NSSitesresults,$input_params,$rst,$rst_persite,$rst_trees ) =
 193           $self->_rearrange([qw
 194                              (TREES MLMATRIX
 195                               SEQS NGMATRIX
 196                               CODONPOS CODONFREQ
 197                               VERSION MODEL PATTERNS
 198                               STATS AAFREQ AADISTMAT
 199                               AAMLDISTMAT
 200                               NTFREQ SEQFILE
 201                               KAPPA_DISTMAT
 202                               ALPHA_DISTMAT
 203                               NSSITESRESULTS
 204                               INPUT_PARAMS
 205                               RST RST_PERSITE RST_TREES)],
 206                             @args);
 207   $self->reset_seqs;
 208   if( $trees ) {
 209       if(ref($trees) !~ /ARRAY/i ) {
 210           $self->warn("Must provide a valid array reference to initialize trees");
 211       } else {
 212           foreach my $t ( @$trees ) {
 213               $self->add_tree($t);
 214           }
 215       }
 216   }
 217   $self->{'_treeiterator'} = 0;
 218
 219   if( $mlmat ) {
 220       if( ref($mlmat) !~ /ARRAY/i ) {
 221           $self->warn("Must provide a valid array reference to initialize MLmatrix");
 222       } else {
 223           $self->set_MLmatrix($mlmat);
 224       }
 225   }
 226   if( $seqs ) {
 227       if( ref($seqs) !~ /ARRAY/i ) {
 228           $self->warn("Must provide a valid array reference to initialize seqs");
 229       } else {
 230           foreach my $s ( @$seqs ) {
 231               $self->add_seq($s);
 232           }
 233       }
 234   }
 235   if( $ngmatrix ) {
 236       if( ref($ngmatrix) !~ /ARRAY/i ) {
 237           $self->warn("Must provide a valid array reference to initialize NGmatrix");
 238       } else {
 239           $self->set_NGmatrix($ngmatrix);
 240       }
 241   }
 242   if( $codonfreq ) {
 243       if( ref($codonfreq) =~ /ARRAY/i ) {
 244           $self->set_CodonFreqs($codonfreq);
 245       } else {
 246           $self->warn("Must provide a valid array reference to initialize codonfreq");
 247       }
 248   }
 249
 250   if( $codonpos ) {
 251       if( ref($codonpos) !~ /ARRAY/i ) {
 252           $self->warn("Must provide a valid array reference to initialize codonpos");
 253       } else {
 254           $self->set_codon_pos_basefreq(@$codonpos);
 255       }
 256   }
 257
 258   $self->version($version)   if defined $version;
 259   $self->seqfile($seqfile)   if defined $seqfile;
 260   $self->model($model)       if defined $model;
 261   if( defined $patterns ) {
 262       if( ref($patterns) =~ /HASH/i ) {
 263           $self->patterns($patterns);
 264       } else {
 265           $self->warn("Must provide a valid array reference to initialize patterns");
 266       }
 267   }
 268
 269   $self->{'_aafreqs'} = {};
 270   if( $aafreq ) {
 271       if( ref($aafreq) =~ /HASH/i ) {
 272           $self->set_AAFreqs($aafreq);
 273       } else {
 274           $self->warn("Must provide a valid hash reference to initialize aafreq");
 275       }
 276   }
 277   if( $stats ) {
 278       if( ref($stats) =~ /HASH/i ) {
 279           while( my ($stat,$val) = each %$stats) {
 280               $self->add_stat($stat,$val);
 281           }
 282       } else {
 283           $self->warn("Must provide a valid hash reference initialize stats");
 284       }
 285   }
 286   $self->set_AADistMatrix($aadistmat) if defined $aadistmat;
 287   $self->set_AAMLDistMatrix($aamldistmat) if defined $aamldistmat;
 288
 289   if( defined $NSSitesresults ) {
 290       if( ref($NSSitesresults) !~ /ARRAY/i ) {
 291           $self->warn("expected an arrayref for -NSSitesresults");
 292       } else {
 293           foreach my $m ( @$NSSitesresults ) {
 294               $self->add_NSSite_result($m);
 295           }
 296       }
 297   }
 298
 299   $self->{'_ntfreqs'} = {};
 300   if( $ntfreqs ) {
 301       if( ref($ntfreqs) =~ /HASH/i ) {
 302           $self->set_NTFreqs($ntfreqs);
 303       } else {
 304           $self->warn("Must provide a valid hash reference to initialize ntfreq");
 305       }
 306   }
 307
 308   if( $kappa_mat ) {
 309       $self->set_KappaMatrix($kappa_mat);
 310   }
 311   if( $alpha_mat ) {
 312       $self->set_AlphaMatrix($alpha_mat);
 313   }
 314
 315   if( $input_params ) {
 316       if(  ref($input_params) !~ /HASH/i ) {
 317           $self->warn("Must provide a valid hash object for input_params\n");
 318       } else {
 319           while( my ($p,$v) = each %$input_params ) {
 320               $self->set_input_parameter($p,$v);
 321           }
 322       }
 323
 324   }
 325   $self->reset_rst_seqs;
 326   if( $rst ) {
 327       if( ref($rst) =~ /ARRAY/i ) {
 328           for ( @$rst ) {
 329               $self->add_rst_seq($_);
 330           }
 331       } else {
 332           $self->warn("Need a valid array ref for -rst option\n");
 333       }
 334   }
 335   if( defined $rst_persite ) {
 336       $self->set_rst_persite($rst_persite);
 337   }
 338   $self->reset_rst_trees;
 339   if( $rst_trees ) {
 340       if( ref($rst_trees) =~ /ARRAY/i ) {
 341           for ( @$rst_trees ) {
 342               $self->add_rst_tree($_);
 343           }
 344       } else {
 345           $self->warn("Need a valid array ref for -rst_trees option\n");
 346       }
 347   }
 348
 349   return $self;
 350 }
 351
 352 =head2 next_tree
 353
 354  Title   : next_tree
 355  Usage   : my $tree = $factory->next_tree;
 356  Function: Get the next tree from the factory
 357  Returns : L<Bio::Tree::TreeI>
 358  Args    : none
 359
 360 =cut
 361
 362 sub next_tree{
 363    my ($self,@args) = @_;
 364    return $self->{'_trees'}->[$self->{'_treeiterator'}++] || undef;
 365 }
 366
 367 =head2 get_trees
 368
 369  Title   : get_trees
 370  Usage   : my @trees = $result->get_trees;
 371  Function: Get all the parsed trees as an array
 372  Returns : Array of trees
 373  Args    : none
 374
 375
 376 =cut
 377
 378 sub get_trees{
 379    my ($self) = @_;
 380    return @{$self->{'_trees'} || []};
 381 }
 382
 383 =head2 rewind_tree_iterator
 384
 385  Title   : rewind_tree_iterator
 386  Usage   : $result->rewind_tree_iterator()
 387  Function: Rewinds the tree iterator so that next_tree can be
 388            called again from the beginning
 389  Returns : none
 390  Args    : none
 391
 392 =cut
 393
 394 sub rewind_tree_iterator {
 395     shift->{'_treeiterator'} = 0;
 396 }
 397
 398 =head2 add_tree
 399
 400  Title   : add_tree
 401  Usage   : $result->add_tree($tree);
 402  Function: Adds a tree
 403  Returns : integer which is the number of trees stored
 404  Args    : L<Bio::Tree::TreeI>
 405
 406 =cut
 407
 408 sub add_tree{
 409    my ($self,$tree) = @_;
 410    if( $tree && ref($tree) && $tree->isa('Bio::Tree::TreeI') ) {
 411        push @{$self->{'_trees'}},$tree;
 412    }
 413    return scalar @{$self->{'_trees'}};
 414 }
 415
 416
 417 =head2 set_MLmatrix
 418
 419  Title   : set_MLmatrix
 420  Usage   : $result->set_MLmatrix($mat)
 421  Function: Set the ML Matrix
 422  Returns : none
 423  Args    : Arrayref to MLmatrix (must be arrayref to 2D matrix whic is
 424            lower triangle pairwise)
 425
 426
 427 =cut
 428
 429 sub set_MLmatrix{
 430    my ($self,$mat) = @_;
 431    return unless ( defined $mat );
 432    if( ref($mat) !~ /ARRAY/i ) {
 433        $self->warn("Did not provide a valid 2D Array reference for set_MLmatrix");
 434        return;
 435    }
 436    $self->{'_mlmatrix'} = $mat;
 437 }
 438
 439 =head2 get_MLmatrix
 440
 441  Title   : get_MLmatrix
 442  Usage   : my $mat = $result->get_MLmatrix()
 443  Function: Get the ML matrix
 444  Returns : 2D Array reference
 445  Args    : none
 446
 447
 448 =cut
 449
 450 sub get_MLmatrix{
 451    my ($self,@args) = @_;
 452    return $self->{'_mlmatrix'};
 453 }
 454
 455 =head2 set_NGmatrix
 456
 457  Title   : set_NGmatrix
 458  Usage   : $result->set_NGmatrix($mat)
 459  Function: Set the Nei & Gojobori Matrix
 460  Returns : none
 461  Args    : Arrayref to NGmatrix (must be arrayref to 2D matrix whic is
 462            lower triangle pairwise)
 463
 464
 465 =cut
 466
 467 sub set_NGmatrix{
 468    my ($self,$mat) = @_;
 469    return unless ( defined $mat );
 470    if( ref($mat) !~ /ARRAY/i ) {
 471        $self->warn("Did not provide a valid 2D Array reference for set_NGmatrix");
 472        return;
 473    }
 474    $self->{'_ngmatrix'} = $mat;
 475 }
 476
 477 =head2 get_NGmatrix
 478
 479  Title   : get_NGmatrix
 480  Usage   : my $mat = $result->get_NGmatrix()
 481  Function: Get the Nei & Gojobori matrix
 482  Returns : 2D Array reference
 483  Args    : none
 484
 485
 486 =cut
 487
 488 sub get_NGmatrix{
 489    my ($self,@args) = @_;
 490    return $self->{'_ngmatrix'};
 491 }
 492
 493
 494 =head2 add_seq
 495
 496  Title   : add_seq
 497  Usage   : $obj->add_seq($seq)
 498  Function: Add a Bio::PrimarySeq to the Result
 499  Returns : none
 500  Args    : Bio::PrimarySeqI
 501 See also : L<Bio::PrimarySeqI>
 502
 503 =cut
 504
 505 sub add_seq{
 506    my ($self,$seq) = @_;
 507    if( $seq ) {
 508        unless( $seq->isa("Bio::PrimarySeqI") ) {
 509            $self->warn("Must provide a valid Bio::PrimarySeqI to add_seq");
 510            return;
 511        }
 512        push @{$self->{'_seqs'}},$seq;
 513    }
 514
 515 }
 516
 517 =head2 reset_seqs
 518
 519  Title   : reset_seqs
 520  Usage   : $result->reset_seqs
 521  Function: Reset the OTU seqs stored
 522  Returns : none
 523  Args    : none
 524
 525
 526 =cut
 527
 528 sub reset_seqs{
 529    my ($self) = @_;
 530    $self->{'_seqs'} = [];
 531 }
 532
 533 =head2 get_seqs
 534
 535  Title   : get_seqs
 536  Usage   : my @otus = $result->get_seqs
 537  Function: Get the seqs Bio::PrimarySeq (OTU = Operational Taxonomic Unit)
 538  Returns : Array of Bio::PrimarySeq
 539  Args    : None
 540 See also : L<Bio::PrimarySeq>
 541
 542 =cut
 543
 544 sub get_seqs{
 545    my ($self) = @_;
 546    return @{$self->{'_seqs'}};
 547 }
 548
 549 =head2 set_codon_pos_basefreq
 550
 551  Title   : set_codon_pos_basefreq
 552  Usage   : $result->set_codon_pos_basefreq(@freqs)
 553  Function: Set the codon position base frequencies
 554  Returns : none
 555  Args    : Array of length 3 where each slot has a hashref
 556            keyed on DNA base
 557
 558
 559 =cut
 560
 561 sub set_codon_pos_basefreq {
 562     my ($self,@codonpos) = @_;
 563     if( scalar @codonpos != 3 ) {
 564         $self->warn("invalid array to set_codon_pos_basefreq, must be an array of length 3");
 565         return;
 566     }
 567     foreach my $pos ( @codonpos ) {
 568         if( ref($pos) !~ /HASH/i ||
 569             ! exists $pos->{'A'} ) {
 570             $self->warn("invalid array to set_codon_pos_basefreq, must be an array with hashreferences keyed on DNA bases, C,A,G,T");
 571         }
 572     }
 573     $self->{'_codonposbasefreq'} = [@codonpos];
 574 }
 575
 576 =head2 get_codon_pos_basefreq
 577
 578  Title   : get_codon_pos_basefreq
 579  Usage   : my @basepos = $result->get_codon_pos_basefreq;
 580  Function: Get the codon position base frequencies
 581  Returns : Array of length 3 (each codon position), each
 582            slot is a hashref keyed on DNA bases, the values are
 583            the frequency of the base at that position for all sequences
 584  Args    : none
 585  Note    : The array starts at 0 so position '1' is in position '0'
 586            of the array
 587
 588 =cut
 589
 590 sub get_codon_pos_basefreq{
 591    my ($self) = @_;
 592    return @{$self->{'_codonposbasefreq'}};
 593 }
 594
 595 =head2 version
 596
 597  Title   : version
 598  Usage   : $obj->version($newval)
 599  Function: Get/Set version
 600  Returns : value of version
 601  Args    : newvalue (optional)
 602
 603
 604 =cut
 605
 606 sub version{
 607    my $self = shift;
 608    $self->{'_version'} = shift if @_;
 609    return $self->{'_version'};
 610 }
 611
 612 =head2 seqfile
 613
 614  Title   : seqfile
 615  Usage   : $obj->seqfile($newval)
 616  Function: Get/Set seqfile
 617  Returns : value of seqfile
 618  Args    : newvalue (optional)
 619
 620
 621 =cut
 622
 623 sub seqfile{
 624    my $self = shift;
 625    $self->{'_seqfile'} = shift if @_;
 626    return $self->{'_seqfile'};
 627 }
 628
 629 =head2 model
 630
 631  Title   : model
 632  Usage   : $obj->model($newval)
 633  Function: Get/Set model
 634  Returns : value of model
 635  Args    : on set, new value (a scalar or undef, optional)
 636
 637
 638 =cut
 639
 640 sub model{
 641     my $self = shift;
 642
 643     return $self->{'_model'} = shift if @_;
 644     return $self->{'_model'};
 645 }
 646
 647
 648 =head2 patterns
 649
 650  Title   : patterns
 651  Usage   : $obj->patterns($newval)
 652  Function: Get/Set Patterns hash
 653  Returns : Hashref of pattern data
 654  Args    : [optional] Hashref of patterns
 655          : The hashref is typically
 656          : { -patterns => \@arrayref
 657          :   -ns       => $ns
 658          :   -ls       => $ls
 659          : }
 660
 661 =cut
 662
 663 sub patterns{
 664     my $self = shift;
 665     return $self->{'_patterns'} = shift if @_;
 666     return $self->{'_patterns'};
 667 }
 668
 669 =head2 set_AAFreqs
 670
 671  Title   : set_AAFreqs
 672  Usage   : $result->set_AAFreqs(\%aafreqs);
 673  Function: Get/Set AA freqs
 674  Returns : none
 675  Args    : Hashref, keys are the sequence names, each points to a hashref
 676            which in turn has keys which are the amino acids
 677
 678
 679 =cut
 680
 681 sub set_AAFreqs{
 682    my ($self,$aafreqs) = @_;
 683
 684    if( $aafreqs && ref($aafreqs) =~ /HASH/i ) {
 685        foreach my $seqname ( keys %{$aafreqs} ) {
 686            $self->{'_aafreqs'}->{$seqname} = $aafreqs->{$seqname};
 687        }
 688    }
 689 }
 690
 691 =head2 get_AAFreqs
 692
 693  Title   : get_AAFreqs
 694  Usage   : my %all_aa_freqs = $result->get_AAFreqs()
 695             OR
 696            my %seq_aa_freqs = $result->get_AAFreqs($seqname)
 697  Function: Get the AA freqs, either for every sequence or just
 698            for a specific sequence
 699            The average aa freqs for the entire set are also available
 700            for the sequence named 'Average'
 701  Returns : Hashref
 702  Args    : (optional) sequence name to retrieve aa freqs for
 703
 704
 705 =cut
 706
 707 sub get_AAFreqs{
 708    my ($self,$seqname) = @_;
 709    if( $seqname ) {
 710        return $self->{'_aafreqs'}->{$seqname} || {};
 711    } else {
 712        return $self->{'_aafreqs'};
 713    }
 714 }
 715
 716 =head2 set_NTFreqs
 717
 718  Title   : set_NTFreqs
 719  Usage   : $result->set_NTFreqs(\%aafreqs);
 720  Function: Get/Set NT freqs
 721  Returns : none
 722  Args    : Hashref, keys are the sequence names, each points to a hashref
 723            which in turn has keys which are the amino acids
 724
 725
 726 =cut
 727
 728 sub set_NTFreqs{
 729    my ($self,$freqs) = @_;
 730
 731    if( $freqs && ref($freqs) =~ /HASH/i ) {
 732        foreach my $seqname ( keys %{$freqs} ) {
 733            $self->{'_ntfreqs'}->{$seqname} = $freqs->{$seqname};
 734        }
 735    }
 736 }
 737
 738 =head2 get_NTFreqs
 739
 740  Title   : get_NTFreqs
 741  Usage   : my %all_nt_freqs = $result->get_NTFreqs()
 742             OR
 743            my %seq_nt_freqs = $result->get_NTFreqs($seqname)
 744  Function: Get the NT freqs, either for every sequence or just
 745            for a specific sequence
 746            The average nt freqs for the entire set are also available
 747            for the sequence named 'Average'
 748  Returns : Hashref
 749  Args    : (optional) sequence name to retrieve nt freqs for
 750
 751
 752 =cut
 753
 754 sub get_NTFreqs{
 755    my ($self,$seqname) = @_;
 756    if( $seqname ) {
 757        return $self->{'_ntfreqs'}->{$seqname} || {};
 758    } else {
 759        return $self->{'_ntfreqs'};
 760    }
 761 }
 762
 763 =head2 add_stat
 764
 765  Title   : add_stat
 766  Usage   : $result->add_stat($stat,$value);
 767  Function: Add some misc stat valuess (key/value pairs)
 768  Returns : none
 769  Args    : $stat  stat name
 770            $value stat value
 771
 772
 773 =cut
 774
 775 sub add_stat{
 776    my ($self,$stat,$value) = @_;
 777    return if( ! defined $stat || !defined $value );
 778    $self->{'_stats'}->{$stat} = $value;
 779    return;
 780 }
 781
 782 =head2 get_stat
 783
 784  Title   : get_stat
 785  Usage   : my $value = $result->get_stat($name);
 786  Function: Get the value for a stat of a given name
 787  Returns : scalar value
 788  Args    : name of the stat
 789
 790
 791 =cut
 792
 793 sub get_stat{
 794    my ($self,$statname) = @_;
 795    return $self->{'_stats'}->{$statname};
 796 }
 797
 798 =head2 get_stat_names
 799
 800  Title   : get_stat_names
 801  Usage   : my @names = $result->get_stat_names;
 802  Function: Get the stat names stored for the result
 803  Returns : array of names
 804  Args    : none
 805
 806
 807 =cut
 808
 809 sub get_stat_names{
 810    my ($self) = @_;
 811    return keys %{$self->{'_stats'} || {}};
 812 }
 813
 814 =head2 get_AADistMatrix
 815
 816  Title   : get_AADistMatrix
 817  Usage   : my $mat = $obj->get_AADistMatrix()
 818  Function: Get AADistance Matrix
 819  Returns : value of AADistMatrix (Bio::Matrix::PhylipDist)
 820  Args    : none
 821
 822
 823 =cut
 824
 825 sub get_AADistMatrix{
 826     my $self = shift;
 827     return $self->{'_AADistMatix'};
 828 }
 829
 830 =head2 set_AADistMatrix
 831
 832  Title   : set_AADistMatrix
 833  Usage   : $obj->set_AADistMatrix($mat);
 834  Function: Set the AADistrance Matrix (Bio::Matrix::PhylipDist)
 835  Returns : none
 836  Args    : AADistrance Matrix (Bio::Matrix::PhylipDist)
 837
 838
 839 =cut
 840
 841 sub set_AADistMatrix{
 842    my ($self,$d) = @_;
 843    if( ! $d ||
 844        ! ref($d) ||
 845        ! $d->isa('Bio::Matrix::PhylipDist') ) {
 846        $self->warn("Must provide a valid Bio::Matrix::MatrixI for set_AADistMatrix");
 847    }
 848    $self->{'_AADistMatix'} = $d;
 849    return;
 850 }
 851
 852 =head2 get_AAMLDistMatrix
 853
 854  Title   : get_AAMLDistMatrix
 855  Usage   : my $mat = $obj->get_AAMLDistMatrix()
 856  Function: Get AAMLDistance Matrix
 857  Returns : value of AAMLDistMatrix (Bio::Matrix::PhylipDist)
 858  Args    : none
 859
 860
 861 =cut
 862
 863 sub get_AAMLDistMatrix{
 864     my $self = shift;
 865     return $self->{'_AAMLDistMatix'};
 866 }
 867
 868 =head2 set_AAMLDistMatrix
 869
 870  Title   : set_AAMLDistMatrix
 871  Usage   : $obj->set_AAMLDistMatrix($mat);
 872  Function: Set the AA ML Distrance Matrix (Bio::Matrix::PhylipDist)
 873  Returns : none
 874  Args    : AAMLDistrance Matrix (Bio::Matrix::PhylipDist)
 875
 876
 877 =cut
 878
 879 sub set_AAMLDistMatrix{
 880    my ($self,$d) = @_;
 881    if( ! $d ||
 882        ! ref($d) ||
 883        ! $d->isa('Bio::Matrix::PhylipDist') ) {
 884        $self->warn("Must provide a valid Bio::Matrix::MatrixI for set_AAMLDistMatrix");
 885    }
 886    $self->{'_AAMLDistMatix'} = $d;
 887    return;
 888 }
 889
 890 =head2 add_NSSite_result
 891
 892  Title   : add_NSSite_result
 893  Usage   : $result->add_NSSite_result($model)
 894  Function: Add a NSsite result (PAML::ModelResult)
 895  Returns : none
 896  Args    : Bio::Tools::Phylo::PAML::ModelResult
 897
 898
 899 =cut
 900
 901 sub add_NSSite_result{
 902    my ($self,$model) = @_;
 903    if( defined $model ) {
 904        push @{$self->{'_nssiteresult'}}, $model;
 905    }
 906    return scalar @{$self->{'_nssiteresult'}};
 907 }
 908
 909 =head2 get_NSSite_results
 910
 911  Title   : get_NSSite_results
 912  Usage   : my @results = @{$self->get_NSSite_results};
 913  Function: Get the reference to the array of NSSite_results
 914  Returns : Array of PAML::ModelResult results
 915  Args    : none
 916
 917
 918 =cut
 919
 920 sub get_NSSite_results{
 921    my ($self) = @_;
 922    return @{$self->{'_nssiteresult'} || []};
 923 }
 924
 925 =head2 set_CodonFreqs
 926
 927  Title   : set_CodonFreqs
 928  Usage   : $obj->set_CodonFreqs($newval)
 929  Function: Get/Set the Codon Frequence table
 930  Returns : value of set_CodonFreqs (a scalar)
 931  Args    : on set, new value (a scalar or undef, optional)
 932
 933
 934 =cut
 935
 936 sub set_CodonFreqs{
 937     my $self = shift;
 938
 939     return $self->{'_codonfreqs'} = shift if @_;
 940     return $self->{'_codonfreqs'};
 941 }
 942
 943 =head2 get_CodonFreqs
 944
 945  Title   : get_CodonFreqs
 946  Usage   : my @codon_freqs = $result->get_CodonFreqs()
 947  Function: Get the Codon freqs
 948  Returns : Array
 949  Args    : none
 950
 951
 952 =cut
 953
 954 sub get_CodonFreqs{
 955    my ($self) = @_;
 956    return @{$self->{'_codonfreqs'} || []};
 957 }
 958
 959
 960 =head2 BASEML Relevant values
 961
 962 =cut
 963
 964 =head2 get_KappaMatrix
 965
 966  Title   : get_KappaMatrix
 967  Usage   : my $mat = $obj->get_KappaMatrix()
 968  Function: Get KappaDistance Matrix
 969  Returns : value of KappaMatrix (Bio::Matrix::PhylipDist)
 970  Args    : none
 971
 972
 973 =cut
 974
 975 sub get_KappaMatrix{
 976     my $self = shift;
 977     return $self->{'_KappaMatix'};
 978 }
 979
 980 =head2 set_KappaMatrix
 981
 982  Title   : set_KappaMatrix
 983  Usage   : $obj->set_KappaMatrix($mat);
 984  Function: Set the KappaDistrance Matrix (Bio::Matrix::PhylipDist)
 985  Returns : none
 986  Args    : KappaDistrance Matrix (Bio::Matrix::PhylipDist)
 987
 988
 989 =cut
 990
 991 sub set_KappaMatrix{
 992    my ($self,$d) = @_;
 993    if( ! $d ||
 994        ! ref($d) ||
 995        ! $d->isa('Bio::Matrix::PhylipDist') ) {
 996        $self->warn("Must provide a valid Bio::Matrix::MatrixI for set_NTDistMatrix");
 997    }
 998    $self->{'_KappaMatix'} = $d;
 999    return;
1000 }
1001
1002
1003 =head2 get_AlphaMatrix
1004
1005  Title   : get_AlphaMatrix
1006  Usage   : my $mat = $obj->get_AlphaMatrix()
1007  Function: Get AlphaDistance Matrix
1008  Returns : value of AlphaMatrix (Bio::Matrix::PhylipDist)
1009  Args    : none
1010
1011
1012 =cut
1013
1014 sub get_AlphaMatrix{
1015     my $self = shift;
1016     return $self->{'_AlphaMatix'};
1017 }
1018
1019 =head2 set_AlphaMatrix
1020
1021  Title   : set_AlphaMatrix
1022  Usage   : $obj->set_AlphaMatrix($mat);
1023  Function: Set the AlphaDistrance Matrix (Bio::Matrix::PhylipDist)
1024  Returns : none
1025  Args    : AlphaDistrance Matrix (Bio::Matrix::PhylipDist)
1026
1027
1028 =cut
1029
1030 sub set_AlphaMatrix{
1031    my ($self,$d) = @_;
1032    if( ! $d ||
1033        ! ref($d) ||
1034        ! $d->isa('Bio::Matrix::PhylipDist') ) {
1035        $self->warn("Must provide a valid Bio::Matrix::MatrixI for set_NTDistMatrix");
1036    }
1037    $self->{'_AlphaMatix'} = $d;
1038    return;
1039 }
1040
1041 =head2 set_input_parameter
1042
1043  Title   : set_input_parameter
1044  Usage   : $obj->set_input_parameter($p,$vl);
1045  Function: Set an Input Parameter
1046  Returns : none
1047  Args    : $parameter and $value
1048
1049
1050 =cut
1051
1052 sub set_input_parameter{
1053    my ($self,$p,$v) = @_;
1054    return unless defined $p;
1055    $self->{'_input_parameters'}->{$p} = $v;
1056 }
1057
1058 =head2 get_input_parameters
1059
1060  Title   : get_input_parameters
1061  Usage   : $obj->get_input_parameters;
1062  Function: Get Input Parameters
1063  Returns : Hash of key/value pairs
1064  Args    : none
1065
1066
1067 =cut
1068
1069 sub get_input_parameters{
1070    my ($self) = @_;
1071    return %{$self->{'_input_parameters'} || {}};
1072 }
1073
1074 =head2 reset_input_parameters
1075
1076  Title   : reset_input_parameters
1077  Usage   : $obj->reset_input_parameters;
1078  Function: Reset the Input Parameters hash
1079  Returns : none
1080  Args    : none
1081
1082
1083 =cut
1084
1085 sub reset_input_parameters{
1086    my ($self) = @_;
1087    $self->{'_input_parameters'} = {};
1088 }
1089
1090 =head1 Reconstructed Ancestral State relevant options
1091
1092 =head2 add_rst_seq
1093
1094  Title   : add_rst_seq
1095  Usage   : $obj->add_rst_seq($seq)
1096  Function: Add a Bio::PrimarySeq to the RST Result
1097  Returns : none
1098  Args    : Bio::PrimarySeqI
1099 See also : L<Bio::PrimarySeqI>
1100
1101 =cut
1102
1103 sub add_rst_seq{
1104    my ($self,$seq) = @_;
1105    if( $seq ) {
1106        unless( $seq->isa("Bio::PrimarySeqI") ) {
1107            $self->warn("Must provide a valid Bio::PrimarySeqI to add_rst_seq");
1108            return;
1109        }
1110        push @{$self->{'_rstseqs'}},$seq;
1111    }
1112
1113 }
1114
1115 =head2 reset_rst_seqs
1116
1117  Title   : reset_rst_seqs
1118  Usage   : $result->reset_rst_seqs
1119  Function: Reset the RST seqs stored
1120  Returns : none
1121  Args    : none
1122
1123
1124 =cut
1125
1126 sub reset_rst_seqs{
1127    my ($self) = @_;
1128    $self->{'_rstseqs'} = [];
1129 }
1130
1131 =head2 get_rst_seqs
1132
1133  Title   : get_rst_seqs
1134  Usage   : my @otus = $result->get_rst_seqs
1135  Function: Get the seqs Bio::PrimarySeq
1136  Returns : Array of Bio::PrimarySeqI objects
1137  Args    : None
1138 See also : L<Bio::PrimarySeq>
1139
1140 =cut
1141
1142 sub get_rst_seqs{
1143    my ($self) = @_;
1144    return @{$self->{'_rstseqs'} || []};
1145 }
1146
1147
1148 =head2 add_rst_tree
1149
1150  Title   : add_rst_tree
1151  Usage   : $obj->add_rst_tree($tree)
1152  Function: Add a Bio::Tree::TreeI to the RST Result
1153  Returns : none
1154  Args    : Bio::Tree::TreeI
1155 See also : L<Bio::Tree::TreeI>
1156
1157 =cut
1158
1159 sub add_rst_tree{
1160    my ($self,$tree) = @_;
1161    if( $tree ) {
1162        unless( $tree->isa("Bio::Tree::TreeI") ) {
1163            $self->warn("Must provide a valid Bio::Tree::TreeI to add_rst_tree not $tree");
1164            return;
1165        }
1166        push @{$self->{'_rsttrees'}},$tree;
1167    }
1168 }
1169
1170 =head2 reset_rst_trees
1171
1172  Title   : reset_rst_trees
1173  Usage   : $result->reset_rst_trees
1174  Function: Reset the RST trees stored
1175  Returns : none
1176  Args    : none
1177
1178
1179 =cut
1180
1181 sub reset_rst_trees{
1182    my ($self) = @_;
1183    $self->{'_rsttrees'} = [];
1184 }
1185
1186 =head2 get_rst_trees
1187
1188  Title   : get_rst_trees
1189  Usage   : my @otus = $result->get_rst_trees
1190  Function: Get the trees Bio::Tree::TreeI
1191  Returns : Array of Bio::Tree::TreeI objects
1192  Args    : None
1193 See also : L<Bio::Tree::TreeI>
1194
1195 =cut
1196
1197 sub get_rst_trees{
1198    my ($self) = @_;
1199    return @{$self->{'_rsttrees'} || []};
1200 }
1201
1202 =head2 set_rst_persite
1203
1204  Title   : set_rst_persite
1205  Usage   : $obj->set_rst_persite($newval)
1206  Function: Get/Set the per-site RST values
1207  Returns : value of set_rst_persite (a scalar)
1208  Args    : on set, new value (a scalar or undef, optional)
1209
1210
1211 =cut
1212
1213 sub set_rst_persite{
1214     my $self = shift;
1215
1216     return $self->{'_rstpersite'} = shift if @_;
1217     return $self->{'_rstpersite'};
1218 }
1219
1220 =head2 get_rst_persite
1221
1222  Title   : get_rst_persite
1223  Usage   : my @rst_persite = @{$result->get_rst_persite()}
1224  Function: Get the per-site RST values
1225  Returns : Array
1226  Args    : none
1227
1228
1229 =cut
1230
1231 sub get_rst_persite{
1232    my ($self) = @_;
1233    return $self->{'_rstpersite'} || [];
1234 }
1235
1236
1237
1238 1;