+1 for Python/ETE - excellent package!
I'm trying to remove a node from a newick tree using bioperl. The tree file contains this data :
(((A:5,B:5)90:2,C:4)25:3,D:10);
Below is the code:
use Bio::TreeIO;
use Bio::Tree::TreeFunctionsI;
use strict;
use warnings;
my $input = new Bio::TreeIO(-file => "tree.newick", -format => "newick");
my $tree = $input->next_tree;
my $bool = $tree-> remove_Node('B');
print $bool;
output: 1;
So the remove_Node is working, but the problem is when I open the tree file I find the node is still present!... What is wrong with the code ? how to remove the node from the tree??
Thanks in advance.
4 answers
Cool, Sphinx'd version of the BioPython tutorial phylogenetics chapter here http://www.bio-cloud.info/Biopython/en/ch12.html
You need to save the tree to a file:
my $treeout = Bio::TreeIO->new(-format => 'newick', -file => ">outputree.nw");
$treeout->write_tree($tree)
That's pretty rubbish! Can you not view the updated tree without writing to a new file?
I don't think so because $tree is a Bio::Tree::TreeI object. I know in Bio::Phylo you can do something like print $tree->as_string but I'm not sure you can do that with BioPerl tree objects. EDIT: This is only partially true, see my latest answer for a more thorough description of these methods.
This question has an accepted answer, but it is not complete. Bio::TreeIO has an as_text method that is available if you upgrade your BioPerl to a recent version, and you may use this method instead of saving to a file.
#!/usr/bin/env perl
use strict;
use warnings;
use Bio::TreeIO;
my $treeio = Bio::TreeIO->new(-fh => \*DATA, -format => 'newick');
my $tree = $treeio->next_tree;
$tree->remove_Node('B');
print $tree->as_text('newick');
__DATA__
(((A:5,B:5)90:2,C:4)25:3,D:10);
Explanation: If we had tried to print $treeio above, instead of calling as_text on the tree, we will see that $treeio is Bio::TreeIO::newick object.
print ref $treeio; # prints Bio::TreeIO::newick
Printing the bare object handle will technically print the object name and a reference showing the underlying data structure, and while you can inspect this, it is usually not very informative unless you know what you are looking at.
Similar to the above, Bio::Phylo has the ability to write the tree string to STDOUT. Here is the equivalent code using Bio::Phylo, which is slightly different.
#!/usr/bin/env perl
use strict;
use warnings;
use Bio::Phylo::IO;
my $tree = Bio::Phylo::IO->parse( -handle => \*DATA,
-format => 'newick',
-keep => ['A','C','D'])->first;
print $tree->to_newick, "\n";
__DATA__
(((A:5,B:5)90:2,C:4)25:3,D:10);
The same applies here about printing with a method versus printing the object itself.
print ref $treeio; # prints Bio::Phylo::Forest::Tree
The actual answer as to why your print statement is printing '1' is that you are trying to assign the modified tree to a variable, but really you are assigning the return value of the method, which is '1' (meaning success). Just call the method on the tree object itself (as in the code above), without assignment, and print or save the tree.
How would you flip this example to specifically remove a node? Something like this?
#!/usr/bin/env perl
use strict;
use warnings;
use Bio::Phylo::IO;
my $tree = Bio::Phylo::IO->parse( -handle => \*DATA,
-format => 'newick',
-remove => ['A'])->first; #####
print $tree->to_newick, "\n";
__DATA__
(((A:5,B:5)90:2,C:4)25:3,D:10);
I guess that would be more expressive and easier to write (also a little redundant), but '-remove' is not an argument that is available to the Bio::Phylo::Parsers::Newick class. So, I think using '-keep' is the way to go.
I would recommend you also consider Bio::Phylo for doing phylogenetic analyses in Perl. For simplicity, I like to do these things in R because you can really see what you are doing. Something like this should work:
library(ape)
phylo <- read.tree("mytree.newick")
drop.tip(phylo, taxon) ## taxon can be a single taxon or a vector of taxon names
+1 for Bio::Phylo mention - much better package for working with trees in Perl!
I tried the bove and R displayed the msg with the taxon removed from the tree.. But I am unable to find the new tree file. do i need to print the new tree again. M using R for the first time. please help
Log in to answer this question.