Decodificando la Gene Ontology

Hola,
a estas alturas es difícil no haber oído hablar de la Gene Ontology (GO), el vocabulario controlado más utilizado en Biología para describir secuencias y sus funciones. Diría yo que la principal aplicación de la GO es facilitar la anotación, manipulación y comparación computacional de grandes volúmenes de secuencias. De hecho, los repositorios de referencia, como por ejemplo UniProt, llevan tiempo anotando sus secuencias usando la jerarquía de la GO. Por lo tanto, actualmente es fácil encontrarse con un montón de secuencias, por ejemplo procedentes de un experimento o generadas por un colaborador, que en vez de descripciones legibles tienen asociados códigos de la GO. Por ejemplo, puedes encontrarte con una proteína de A.thaliana,

> AT3G04220  GO:0005524  

ATGGATTCTT CTTTTTTACT CGAAACTGTT GCTGCTGCAA CAGGCTTCTT
CACACTTTTG GGTACAATAC TTTTTATGGT TTACAGAAAA TTCAAAGACC
ATCGAGAAAA CAAAGAAAAT GATTCTTCTT CTTCAACACA ...

y preguntarte, qué será? Cuál es su función? En el Taller de (bio)Perl ya apuntamos una manera de averiguar la genealogía del término GO:0005524, pero dadas las limitaciones del módulo GO::Parser aquí os muestro una manera más eficiente de extraer el significado y la genealogía de cualquier identificador de la GO, que puede pertenecer a una de las 3 ramas fundamentales de la jerarquía  (proceso biológico, función molecular y localización celular). Para ello es necesario descargar la versión más reciente de la GO en un archivo de texto en formato OBO, y guardar el siguiente módulo escrito en Perl:

 package myGOparser;  
 require Exporter;  
   
 @ISA = qw(Exporter);  
 @EXPORT = qw( $DEFGOFLATFILE parse_GO_flat_file get_full_annotation_for_id );   
 use strict;  
   
 our $DEFGOFLATFILE = "/path/to/GO/gene_ontology.1_2.obo";  
   
 sub parse_GO_flat_file  
 {  
    my $flatfile = $_[0] || $DEFGOFLATFILE;  
    my ($id,$alt_id,$name,$namespace,%synonim,%GO);  
      
    open(GOFILE,$flatfile) ||   
       die "# parse_GO_flat_file : cannot read $flatfile\n";  
    while(<GOFILE>)  
    {  
       #[Term]  
       #id: GO:0006355  
       #name: regulation of transcription, DNA-dependent  
       #namespace: biological_process  
       #alt_id: GO:0032583  
       #is_a: GO:0010468 ! regulation of gene expression  
       if(/^id: (GO:\d+)/){ $id = $1; }  
       elsif(/^name: (.+?)\n/){ $GO{$id}{'name'} = $1; }  
       elsif(/^alt_id: (GO:\d+)/){ push(@{$synonim{$id}},$1); }  
       elsif(/^is_a: (GO:\d+)/){ push(@{$GO{$id}{'parents'}},$1); }  
    }  
    close(GOFILE);  
      
    # link synonims  
    foreach $id (keys(%synonim))  
    {  
       foreach my $syn (@{$synonim{$id}}){ $GO{$syn} = $GO{$id} }  
    }  
      
    return \%GO;  
 }  
   
 sub get_full_annotation_for_id  
 {  
    my ($id,$ref_parsed_GO) = @_;  
    my $annot;  
      
    if(!$ref_parsed_GO->{$id})  
    {  
       return "[cannot find GO term $id]";  
    }  
    else  
    {  
       if(!$ref_parsed_GO->{$id}{'parents'})  
       {  
          $annot .= $ref_parsed_GO->{$id}{'name'}."($id) | ";  
       }  
       else  
       {  
          foreach my $pid (@{$ref_parsed_GO->{$id}{'parents'}})  
          {  
             $annot .= $ref_parsed_GO->{$id}{'name'}."($id)".', '.  
               get_full_annotation_for_id($pid,$ref_parsed_GO);  
          }  
       }     
    }  
      
    return $annot;  
 }  
 1;  

Ahora podemos crear un programa como éste, que invocado desde el terminal nos devuelve la descripción jerárquica completa de un termino GO:

 use strict;  
 use myGOparser;  
   
 my $input_id = $ARGV[0] || die "# usage: ./go_parser.pl \n\n";  
   
 print "# parsing GO flat file ... ";  
 my $ref_GO = parse_GO_flat_file();  
 print "done\n\n";  
   
 print get_full_annotation_for_id($input_id,$ref_GO)."\n";  

Para el identificador GO:0046983 obtendremos la siguiente descripción:

protein dimerization activity(GO:0046983), protein binding(GO:0005515), binding(GO:0005488), molecular_function(GO:0003674) |