1 de julio de 2010

A qué grupo taxonómico pertenece una especie?

Desde hace ya tiempo hay más secuencias genómicas disponibles que las que podemos recordar y a menudo nos encontramos con nombres de organismos que no podemos ubicar, no sabemos si son bacterias aguas termales o artrópodos, por ejemplo. El siguiente programa puede ayudarnos precisamente a colocar dentro de la clasificación taxonómica vigente del NCBI un nombre de especie, o en general, cualquier taxón de nombre aceptado. Antes de ejecutarlo debes acceder al servidor FTP en ftp://ftp.ncbi.nih.gov/pub/taxonomy y descomprimir los archivos de nombres y nodos, como se explica en la cabecera del código Perl, haciendo apuntar la variable global $TAXONOMYPATH al directorio donde se encuentren. Una vez hechos estos ajustes el código ya está listo para usarse desde el terminal. Por ejemplo, si consultamos la taxonomía de Oryza sativa obtendremos esta respuesta, que podemos comparar con la figura, tomada de www.biomedcentral.com/1471-2164/8/283:


$ ./prog6.pl Oryza sativa 

# input taxon: Oryza sativa

# selected taxonID: 4530 \\ 
   URL: http://www.ncbi.nlm.nih.gov/Taxonomy/Browser/wwwtax.cgi?id=4530

# variants:
> Oryza sativa
> Oryza sativa L.

# lineage : Oryza sativa, Oryza, Oryzeae, Ehrhartoideae, BEP clade, Poaceae, 
Poales,commelinids, Liliopsida, Magnoliophyta, Spermatophyta, Euphyllophyta, 
Tracheophyta, Embryophyta, Streptophytina, Streptophyta, Viridiplantae, 
Eukaryota, cellular organisms

# rank : species, genus, tribe, subfamily, no_rank, family, order, subclass, 
class, no_rank, no_rank, no_rank, no_rank, no_rank, no_rank, phylum, kingdom, 
superkingdom, no_rank,

Éste es el código:

 #!/usr/bin/perl -w  
   
 # prog6  
 # programa que devuelve informacion taxonomica para un taxon de entrada,   
 # que se apoya en la version texto de la clasificacion Taxonomy del NCBI,   
 # que se pueden descargar de ftp://ftp.ncbi.nih.gov/pub/taxonomy , en   
 # particular el fichero taxdump.tar.gz, que debes descomprimir para obtener   
 # los dos ficheros que se utilizan aqui  
 # Bruno Contreras Moreira, Junio de 2010  
   
 use strict;  
 $|=1;  
   
 if(!$ARGV[0])  
 {   
    die "# usage: $0 <taxon name from superkingdom to species, such as:\n".  
    "       'Arabidopsis thaliana' or 'Viridiplantae' but not 'thaliana'>\n";   
 }  
   
 my $VERBOSE   = 0 ;  
 my $TAXONOMYPATH = '/path/taxonomy/'; # directorio de taxdump.tar.gz  
 my $NAMESFILE  = $TAXONOMYPATH . 'names.dmp';  
 my $NODESFILE  = $TAXONOMYPATH . 'nodes.dmp';  
 my $NCBIURL   = 'http://www.ncbi.nlm.nih.gov/Taxonomy/Browser/wwwtax.cgi?id=';  
   
 ###################################################################  
   
 my ($n_of_taxa,$taxonID,$taxon,$ID,$name,$rank,$parentID) = (0,'','');  
 my (@names,%node,@lineage,@sorted_lineage,%lineage_name);  
 my ($lineage_names,$lineage_ranks) = ('','');  
   
 ## 1) captura el nombre del taxon de entrada, como Arabidopsis thaliana, y  
 ## sustituye caracteres '_' por espacios si los hubiera   
 if(scalar(@ARGV) > 1){ @ARGV = @ARGV[0,1]; } # corta tercera palabra < especie  
 else{ @ARGV = split(/_/,$ARGV[0]); @ARGV = @ARGV[0,1]; }  
   
 $taxon = join(' ',@ARGV);  
 $taxon =~ s/^\s+//g; # elimina espacios iniciales y finales  
 $taxon =~ s/\s+$//g;   
   
 print "# input taxon: $taxon\n\n";  
   
 ## 2) busca el taxonID del taxon de entrada  
 open(NAMES,$NAMESFILE) || die "# $0 : cannot open $NAMESFILE\n";  
 while(<NAMES>)  
 {  
    next if(!/$taxon/i);  
      
    #3702   |   Arabidopsis thaliana   |      |   scientific name  
    if(/^(\d+)\s+\|\s+(.*?)\|/)  
    {  
       ($ID,$name) = ($1,$2);   
       $name =~ s/^\s+//g;   
       $name =~ s/\s+$//g;   
            
       # comprueba que el taxon buscado coincide con el principio de $name,   
       # evita aceptar cadenas que en realidad coinciden con la especie,   
       # que a veces se pueden encontrar en generos muy distintos  
       next if(index($name,$taxon) != 0);  
         
       # si el taxon buscado contiene una sola palabra evita  
       # aceptar cadenas de especies, que tendran dos palabras  
       next if($taxon !~ /\s+/ && $name =~ /\s+/);  
            
       if(!$n_of_taxa || $ID eq $taxonID) # guarda primer taxonID   
       {  
          $taxonID = $ID;  
          push(@names,$name);  
          $n_of_taxa++;  
               
          if($VERBOSE){ print }  
       }  
       else  
       {  
          if($VERBOSE)  
          {   
             # muestra otras posibles apariciones del taxon  
             print "# skip taxon: $ID\n";   
          }     
          else{ last } # termina la lectura del fichero   
       }     
    }        
 }  
 close(NAMES);  
   
 if(!$n_of_taxa){ die "# cannot find taxon $taxon in $NAMESFILE , exit\n"; }  
 else  
 {   
    print "# selected taxonID: $taxonID  URL: $NCBIURL$taxonID\n";  
    print "\n# variants:\n";   
    foreach $name (@names){ print "> $name\n"; }  
 }  
   
 ## 3) busca taxon padre del taxonID encontrado  
 open(NODES,$NODESFILE) || die "# $0 : cannot open $NODESFILE\n";  
 while(<NODES>)  
 {  
    #132596   |   79512   |   genus...  
    #tax_id         -- node id in GenBank taxonomy database  
     #parent tax_id      -- parent node id in GenBank taxonomy   
     #rank         -- rank of this node (superkingdom, kingdom, ...)   
    next if(/forma/ || /varietas/); # evita rangos < especie  
    if(/^$taxonID\t/)  
    {  
       my @data = split(/\t|\t/,$_);  
       $parentID = $data[2];  
       $rank   = $data[4];   
       if($VERBOSE){ print }  
       last;  
    }  
 }  
 close(NODES);   
   
 if(!$parentID){ die "\n# cannot find taxon '$taxon' in $NODESFILE , exit\n"; }  
   
 ## 4) lee arbol taxonomico del nivel de especie hacia arriba   
 ## (son en 2010 35MB, caben en cualquier RAM)  
 open(NODES,$NODESFILE) || die "# $0 : cannot open $NODESFILE\n";  
 while(<NODES>)  
 {  
    next if(/forma/ || /varietas/);  
    if(/^(\d+)\t\|\t(\d+)\t\|\t(\S+)/)  
    {  
       $node{$1}{'p'} = $2; # p = parent node  
       $node{$1}{'r'} = $3; # r = rank  
    }  
 }  
 close(NODES);   
   
 ## 5) reconstruye linaje de $taxonID en terminos de IDs  
 # 5.1) primero guarda datos de $taxonID  
 push(@lineage,$taxonID);  
 $lineage_name{$taxonID} = $rank;  
 $node{$taxonID}{'r'} = $rank;  
   
 # 5.2) ahora guarda los de los taxones superiores  
 while($parentID ne $taxonID && $parentID > 1)  
 {  
    push(@lineage,$parentID);  
    $parentID = $node{$parentID}{'p'};  
 }  
 @sorted_lineage = sort{$a<=>$b}(@lineage);  
 $taxon = shift(@sorted_lineage);   
   
 ## 6) recupera los nombres del linaje de $taxonID leyendo NAMES de nuevo,  
 ## que contiene los taxones ordenados de menor a mayor  
 open(NAMES,$NAMESFILE) || die "# $0 : cannot open $NAMESFILE\n";  
 while(<NAMES>)  
 {  
    next if(!/scientific name/);  
    if(/^(\d+)\s+\|\s+(.*?)\|/)  
    {  
       ($ID,$name) = ($1,$2);   
       next if($ID ne $taxon);  
         
       $name =~ s/^\s+//g;   
       $name =~ s/\s+$//g;   
       $lineage_name{$ID} = $name;  
       if(@sorted_lineage){ $taxon = shift(@sorted_lineage) }  
       else{ last } # termina si no quedan mas IDs  
    }        
 }  
 close(NAMES);  
   
 ## 7) enumera el linaje de $taxoID  
 foreach $ID (@lineage)  
 {  
    $rank = $node{$ID}{'r'};  
    if($rank eq 'no'){ $rank = 'no_rank' }  
      
    $lineage_names .= "$lineage_name{$ID}, ";  
    $lineage_ranks .= "$rank, ";  
 }  
   
 print "\n# lineage : $lineage_names\n# rank  : $lineage_ranks\n";  
   

No hay comentarios:

Publicar un comentario