esta entrada es para compartir un script Perl para enraizar árboles, muchos árboles, de manera automática, y ordenar los nodos por distancia. Pongo el título en inglés para los buscadores.
Rubén Sancho y yo estamos probando el software GRAMPA, que requiere una colección de árboles de genes precalculados en formato Newick, pero además los quiere enraizados y ordenados de manera ascendente. Como son más de mil árboles no lo queríamos hacer uno a uno con FigTree, por poner un ejemplo; queríamos automatizarlo y para ello aprovechamos los módulos Bio::TreeIO y Bio::Phylo. El primero es de Bioperl y ya lo tenía instalado, el segundo lo instalé con: $ sudo cpanm Bio::Phylo
El árbol de ejemplo, contenido en el fichero unrooted.ph, es:
(1_Sbic:0.047,2_Osat:0.068,(((((((((((3_B422:0.007,(((4_Barb:0.000,((5_Bret:0.003,(6_BsyE:0.000,7_BsyE:0.000):0.000):0.000,8_BsyC:0.000):0.000):0.000,9_Bpho:0.002):0.000,(10_BsyC:0.000,11_BsyC:0.000):0.001):0.000):0.008,((12_B422:0.005,13_Bpho:0.005):0.002,14_Barb:0.000):0.005):0.000,((((15_Bdis:0.000,16_Bdis:0.000):0.000,17_Bmex:0.033):0.013,18_Bhyb:0.001):0.014,(19_Bboi:0.021,((20_Bmex:0.017,21_Bsta:0.034):0.012,22_Bsta:0.000):0.003):0.020):0.012):0.002,23_Barb:0.000):0.000,24_Brup:0.012):0.000,25_BsyG:0.000):0.000,(26_Bpin:0.000,27_BsyG:0.000):0.000):0.011,28_Bsta:0.000):0.053,29_BsyG:0.000):0.036,(30_Barb:0.005,(31_Bpho:0.024,(32_BsyC:0.000,33_BsyE:0.000):0.027):0.029):0.000):0.005,34_Bboi:0.030):0.035);
Invocando el script lo enraizamos con el taxón elegido, el outgroup '_Sbic':
$ perl reroot_tree.pl unrooted.ph > rooted.ph
Obtenemos el fichero rooted.ph, que podemos visualizar con FigTree:
Árbol enraizado y con orden creciente de nodos. |
#!/usr/bin/perl -w
# Re-roots an input Newick tree with a user-defined outgroup and
# prints the resulting tree in ascending or descending node order.
# Based on https://github.com/phac-nml/snvphyl-tools/blob/master/rearrange_snv_matrix.pl
use strict;
use Bio::TreeIO;
use Bio::Phylo::IO;
use Bio::Phylo::Forest::Tree;
my $OUTGROUPSTRING = '_Sbic'; # change as needed
my $NODEORDER = 1; # 1:increasing, 0:decreasing
my ($outfound,$outnode,$sorted_newick) = (0);
die "# usage: $0 <tree.newick>\n" if(!$ARGV[0]);
# read input tree
my $input = new Bio::TreeIO(-file=>$ARGV[0],-format=>'newick');
my $intree = $input->next_tree();
# find outgroup taxon
for my $node ( $intree->get_nodes() ) {
if(defined($node->id()) && $node->id() =~ m/$OUTGROUPSTRING/) {
$outnode = $node;
$outfound = 1;
last;
}
}
if($outfound == 0) {
die "# cannot find outgroup $OUTGROUPSTRING in input tree $ARGV[0]\n";
}
# root in outgroup and sort in increasing order
$intree->reroot_at_midpoint($outnode);
# sort nodes in defined order and print
my $unsorted_tree = Bio::Phylo::IO->parse(
'-string' => $intree->as_text('newick'),
'-format' => 'newick'
)->first();
$unsorted_tree->ladderize($NODEORDER);
$sorted_newick = $unsorted_tree->to_newick();
$sorted_newick =~ s/'//g;
print $sorted_newick;
Un saludo,
Bruno