[Bioperl-l] PAML nssites model result object

[Edward Chuong]

Hi Jason,

Thanks for the help.

The code seems to get stuck at 

 if( ! $node->ancestor || ! $node->has_tag('t') ) {
(this condition turns out true for every node, not just root, so it
always hits "next")

I used Data::Dumper to check on the node and I've pasted the
results--it seems like those tags aren't being sent in?


Thanks!
-Ed

'_root_cleanup_methods' => [
		  sub { "DUMMY" }
		],
'_creation_id' => 0,
'_branch_length' => '0.613722',
'_desc' => {},
'_id' => 'NP_033437.2_mus',
'_ancestor' => bless( {
	 '_root_cleanup_methods' => [
					  $VAR1->{'_root_cleanup_methods'}[0]
					],
	 '_creation_id' => 3,
	 '_desc' => {
				  '2' => bless( {
					  '_root_cleanup_methods' => [
						   $VAR1->{'_root_cleanup_methods'}[0]
						 ],
					  '_creation_id' => 2,
					  '_branch_length' => '0.768322',
					  '_desc' => {},
					  '_id' => 'PM_BWp0001H02f',
					  '_ancestor' => $VAR1->{'_ancestor'},
					  '_root_verbose' => 0
					}, 'Bio::Tree::Node' ),
				  '0' => $VAR1,
				  '1' => bless( {
					  '_root_cleanup_methods' => [
												   $VAR1->{'_root_cleanup_methods'}[0]
												 ],
					  '_creation_id' => 1,
					  '_branch_length' => '0.366319',
					  '_desc' => {},
					  '_id' => 'NP_742070.1_rat',
					  '_ancestor' => $VAR1->{'_ancestor'},
					  '_root_verbose' => 0
					}, 'Bio::Tree::Node' )
				},
	 '_id' => '',
	 '_height' => undef,
	 '_root_verbose' => 0
   }, 'Bio::Tree::Node' ),
'_root_verbose' => 0
}, 'Bio::Tree::Node' );


On Wed, 9 Mar 2005 18:01:34 -0500, Jason Stajich <jason.stajich at duke.edu> wrote:
> Resend with code pasted....
> 
> #!/usr/bin/perl -w
> use strict;
> use Bio::Tools::Phylo::PAML;
> 
> my $outcodeml = shift(@ARGV);
> 
> my $paml_parser = new Bio::Tools::Phylo::PAML(-file => "./$outcodeml",
>                                               -dir => "./");
> my $result = $paml_parser->next_result();
> my $MLmatrix = $result->get_MLmatrix(); # get MaxLikelihood Matrix
> my @otus = $result->get_seqs;
> if( $#{$MLmatrix} < 0 ) {
>      for my $tree ($result->next_tree ) {
>         for my $node ( $tree->get_nodes ) {
>             my $id;
>             if( $node->is_Leaf() ) {
>                 $id = $node->id;
>             } else {
>                 $id = "(".join(",", map { $_->id } grep { $_->is_Leaf }
>                            $node->get_all_Descendents) .")";
>             }
>             if( ! $node->ancestor || ! $node->has_tag('t') ) {
>                 # skip when no values have been associated with this node
>                 # (like the root node)
>                 next;
>             }
>              # I know this looks complicated
>             # but we use the get_tag_values method to pull out the annotations
>             # for each branch
>             # The ()[0] around the call is because get_tag_values returns a
> list
>             # if we want to just get the 1st item in the list we have
>             # to tell Perl we are treating it like an array.
>             # in the future get_tag_values needs to be smart and just
>             # return the 1st item in the array if called in scalar
>             # context
> 
>             printf
> "%s\tt=%.3f\tS=%.1f\tN=%.1f\tdN/
> dS=%.4f\tdN=%.4f\tdS=%.4f\tS*dS=%.1f\tN*dN=%.1f\n",
>             $id,
>             map { ($node->get_tag_values($_))[0] }
>             qw(t S N dN/dS dN dS), 'S*dS', 'N*dN';
>         }
>      }
> } else {
>      my $i =0;
>      my @seqs = $result->get_seqs;
>      for my $row ( @$MLmatrix ) {
>         print $seqs[$i++]->display_id, join("\t",@$row), "\n";
>      }
> }
> 
> On Mar 9, 2005, at 5:41 PM, Jason Stajich wrote:
> 
> > I just updated things last week so this is brand-spanking-new.  I
> > don't know if I connected everything up for NSsites stuff quite yet
> > as that is handled in - the branch-specific parsing should work now.
> > I don't know if the synopsis code is really up to snuff either.  When
> > I get around to it I will try and see what still needs to be connected
> > in NSsites parsing.
> >
> > I don't think $node->param() is going to work -
> > $node->get_tag_values() is the way I've implemented it.
> >
> > <00parse_codeml.pl>
> >
> > -jason
> > --
> > Jason Stajich
> > jason.stajich at duke.edu
> > http://www.duke.edu/~jes12/
> >
> > On Mar 9, 2005, at 5:23 PM, Edward Chuong wrote:
> >
> >> Hi all,
> >>
> >> I'm trying to parse PAML results, and running into some trouble. I'm
> >> using branch specific omega model, and I want to get the branch
> >> specific ka/ks values out.
> >> http://cvs.bioperl.org/cgi-bin/viewcvs/viewcvs.cgi/bioperl-live/Bio/
> >> Tools/Phylo/PAML.pm?rev=HEAD&cvsroot=bioperl&content-type=text/
> >> vnd.viewcvs-markup
> >> says that $node->param('omega') should work, but Data::Dumper shows
> >> that this value isn't stored in the node (only branch lengths and seq
> >> IDs appear to be stored).
> >>
> >> I'm assuming that I can get these values out of the
> >> get_NSSite_result() Bio::Tools::Phylo::PAML::ModelResult object, but
> >> I'm not sure how to call it. The current synopsis uses
> >> "get_model_params" but it seems to be out of date because it's not in
> >> the current souce. The docs at
> >> http://cvs.bioperl.org/cgi-bin/viewcvs/viewcvs.cgi/bioperl-live/Bio/
> >> Tools/Phylo/PAML/Result.pm?rev=HEAD&cvsroot=bioperl&content-
> >> type=text/vnd.viewcvs-markup
> >> say to use my
> >> @results = @{$self->get_NSSite_results};
> >> --that looks like a mistake, and I've tried
> >> @result = $result->get_NSSite_results but that doesn't work either
> >> (just get undefined objs).
> >>
> >> Am I doing something wrong, or is this functionality still being
> >> worked on? I've tried using both 1.4 and the LIVE versions. Any help
> >> is appreciated, thanks!
> >>
> >> -Ed
> >> _______________________________________________
> >> Bioperl-l mailing list
> >> Bioperl-l at portal.open-bio.org
> >> http://portal.open-bio.org/mailman/listinfo/bioperl-l
> >>
> 
> 


-- 
Edward Chuong
(949) 939-2732
AIM: edawad85