Word + bioinformática = CRISPR

Hola,
este será un artículo muy breve, porque tengo poco que escribir yo mismo, a parte de volver a comprobar una vez más que la ciencia básica, tan denostada a menudo, es a menudo fundamental para los avances más espectaculares de la ciencia aplicada.  En vez de contar yo aquí una historia hoy os invito a leer la increíble historia del descubrimiento de los sistemas CRISPR-Cas9, que ya están revolucionando la ciencia y pronto la medicina.

La fuente original, Cell, bajo subscripción:
http://www.cell.com/cell/fulltext/S0092-8674%2815%2901705-5

Un resumen en español en El Confidencial:
http://www.elconfidencial.com/tecnologia/2016-01-16/crispr-francis-mojica-charpentier-doudna-edicion-genomica_1136337/

Reconocimiento de una hebra sencilla de DNA por hibridación de un RNA guía de 32 bases (crRNA) por parte del complejo de proteínas Cascade. Tomado de [http://www.ncbi.nlm.nih.gov/pmc/articles/PMC4427192/]. Como se ve en nuestro repositorio 3D-footprint, la especificidad de reconocimiento es de las más altas observadas entre proteínas naturales: http://floresta.eead.csic.es/3dfootprint/complexes/4qyz_ABCDEFGHIJ.html

 Ah, y por cierto, esta historia una cura de humildad para los que nos dedicamos a inventar y probar complicados algoritmos. A veces, la herramienta "Buscar y reemplazar" de Word es todo lo que hace falta para hacer descubrimientos en Biología, como a veces me recuerda mi colega Ana Casas, y eso también es Biología computacional,
un saludo,
Bruno

Bioinformática en C (y perl6 async)

Hola, y feliz año!
tras la última entrada, y antes de que escriba algo más completo sobre perl6,
parece que la gestión de tareas asíncronas tiene mucho potencial...

Pero vamos a lo que vamos. Hoy quería compartir un documento sobre cómo programar en C con los compiladores actuales, porque algunos todavía seguimos usando sintaxis anticuada y esto es relevante en nuestro campo, dado que hay muchas aplicaciones bioinformáticas escritas en C, normalmente por cuestiones de eficiencia.

El documento en cuestión es: https://matt.sh/howto-c

y ofrece un montón de trucos y recomendaciones para modernizar y simplificar la escritura de código en C, con ejemplos sencillos como éste para vectores dinámicos:

 uintmax_t arrayLength = strtoumax(argv[1], NULL, 10);
    void *array[arrayLength];

 /* Hace innecesario liberar el puntero al final.

    Ojo: debe ser un tamaño razonable para el sistema donde se ejecute. */

Ha sido traducido al español aquí y comentado/criticado aquí y en menéame,
un saludo,
Bruno

Despedimos el 2015, viene Perl 6

Buenas,
quedan ya pocos días para terminar el año, o lo que es lo mismo, pocos días para que podamos probar la primera versión oficial de Perl 6, tras demasiados (15) años de desarrollo, de la mano del compilador Rakudo. Si no os suena nada de esto, podéis empezar por leer las preguntas frecuentes o directamente empezar tutoriales como perl6intro.

Ha habido mucho morbo con la muerte de Perl, y seguirá habiendo a lo largo de 2016 discusiones sobre la compatibilidad de perl 6 y 5, el que usamos ahora, sobre la eficiencia del nuevo lenguaje, o sobre la portabilidad de CPAN a perl 6. Todo esto habrá que ir viendo como se desarrolla muy pronto, pero desde luego es una noticia importante para los usuarios de Perl, esperemos que buena!

Un saludo y los mejores deseos para el próximo año, long live perl!
Bruno

Post-datas del 29 de Diciembre: 

1) Podéis ampliar esto en [http://perlweekly.com/archive/231.html]

2) Para instalar rakudo vía git debes tener abierto el puerto 9418; en caso contrario es posible hacerlo vía https añadiendo estas líneas a tu archivo ~/.gitconfig: 

[url "https://"]
    insteadOf = git://
[url "https://github.com"]
    insteadOf = git@github.com

Post-data del 1 de Enero: 

1) Tutorial exprés en [https://learnxinyminutes.com/docs/perl6]
2) Si ya usas perl5 [http://doc.perl6.org/language/5to6-nutshell]

Post-data de Junio de 2016: 

Mientras no sea pueda instalar perl6 directamente en Ubuntu prefiero usar rakudobrew para instalar perl6 como un ecosistema aislado, como se explica en: [http://rakudo.org/how-to-get-rakudo/#Installing-Rakudo-Star-Source-Rakudobrew]

Si los médicos tuvieran que trabajar de bioinformáticos...

Hoy os publico una viñeta que me ha hecho gracia y se adapta muy bien a la bioinformática... desgraciadamente se cumple totalmente, sobretodo lo del conserje...

El original aquí: http://sinergiasincontrol.blogspot.com/2008/10/33-suposiciones.html

Anotando datos del NCBI de forma automática con un script de Perl

Esta semana me di cuenta que un antiguo script de Perl que leía datos bibliográficos de Pubmed automáticamente había dejado de funcionar.

La explicación es que el NCBI ha dejado de dar soporte a su servicio web con el protocolo SOAP desde el 1 de julio de 2015. No soy un experto informático, pero hace unos años ese protocolo estaba de moda y ahora parece que ya no (imagino que por razones de seguridad).

He aquí el código OBSOLETO que da el error:

 use SOAP::Lite;  
   
 my $pubmed_id = '24234003';  
   
 my $WSDL = 'http://www.ncbi.nlm.nih.gov/soap/v2.0/efetch_pubmed.wsdl';  
 my $soap = SOAP::Lite->service($WSDL)  
                ->on_fault( sub { my($self, $res) = @_;  
                     die  
                     "faultcode:", ref $res ? $res->faultcode : $self->transport->status, "\n" ,  
                     "faultstring:", ref $res ? $res->faultstring : $self->transport->status, "\n";  
                });  
 my $response = $soap->run_eFetch(     SOAP::Data->name(id => $pubmed_id),  
                          SOAP::Data->name(db => "pubmed")  
                     );  
   
 my %results;  
   
 if (ref($response) eq "HASH") {  
      $results{'title'} = $response->{'PubmedArticle'}{'MedlineCitation'}{'Article'}{'ArticleTitle'};  
      $results{'journal'} = $response->{'PubmedArticle'}{'MedlineCitation'}{'Article'}{'Journal'}{'Title'};  
      $results{'pubmed'} = $response->{'PubmedArticle'}{'MedlineCitation'}{'PMID'};  
      $results{'volume'} = $response->{'PubmedArticle'}{'MedlineCitation'}{'Article'}{'Journal'}{'JournalIssue'}{'Volume'};  
      $results{'issue'} = $response->{'PubmedArticle'}{'MedlineCitation'}{'Article'}{'Journal'}{'JournalIssue'}{'Issue'};  
      $results{'year'} = $response->{'PubmedArticle'}{'MedlineCitation'}{'Article'}{'Journal'}{'JournalIssue'}{'PubDate'}{'Year'};  
      $results{'pages'} = $response->{'PubmedArticle'}{'MedlineCitation'}{'Article'}{'Pagination'}{'MedlinePgn'};  
      $results{'url'} = 'http://www.ncbi.nlm.nih.gov/pubmed/'.$results{'pubmed'};  
      if (ref($response->{'PubmedArticle'}{'MedlineCitation'}{'Article'}{'AuthorList'}{'Author'}) eq "ARRAY") {  
           my @authors;  
           foreach my $author (@{$response->{'PubmedArticle'}{'MedlineCitation'}{'Article'}{'AuthorList'}{'Author'}}){  
                my $initials = join('.',split('',$author->{'Initials'}));  
                push(@authors, $author->{'LastName'}." ".$initials);  
           }  
           $results{'authors'} = join(", ",@authors);  
      } elsif (ref($response->{'PubmedArticle'}{'MedlineCitation'}{'Article'}{'AuthorList'}{'Author'}) eq "HASH") {  
           my $author = $response->{'PubmedArticle'}{'MedlineCitation'}{'Article'}{'AuthorList'}{'Author'};  
           my $initials = join('.',split('',$author->{'Initials'}));  
           $results{'authors'} = $author->{'LastName'}." ".$initials;  
      }  
      foreach my $key (keys %{$results}){  
           if (!defined($results{$key})) {  
                $results{$key} = '';  
           }  
      }  
 }  
   
 print %results;  
   

Ahora las consultas automáticas al NCBI se han de realizar insertando nuestra consulta en una URL e indicando el formato deseado de respuesta. Por ejemplo, la consulta equivalente al código SOAP anterior se puede realizar con la siguiente URL:

• http://eutils.ncbi.nlm.nih.gov/entrez/eutils/efetch.fcgi?db=pubmed&id=24234003 

 O si preferimos el resultado en XML:

• http://eutils.ncbi.nlm.nih.gov/entrez/eutils/efetch.fcgi?db=pubmed&id=24234003&retmode=xml

De esta forma, el código anterior queda reescrito de la siguiente forma:

 use LWP::Simple;  
   
 my $pubmed_id = '24234003';  
   
 my $eutils_url = 'http://eutils.ncbi.nlm.nih.gov/entrez/eutils/';  
 my $eutil = 'efetch';  
 my $db = 'pubmed';  
 my $type = 'docsum';  
 my $mode = 'text';  
   
 my $url = sprintf('%s%s.fcgi?db=%s&id=%s&rettype=%s&retmode=%s', $eutils_url, $eutil, $db, $pubmed_id, $type, $mode);  
   
 my $data = get($url);  
   
 my %results;  
   
 if (defined($data)){  
      $data =~ s/\012\015?|\015\012?/ /g;  
      if ($data =~ /\d+\:\s*(.+?)\.(.+?)\.(.+?)\.(.+?)\s.+?\;(\d+)\((\d+)\)\:([\d\-]+)./){  
           @{$results{'authors'}} = split(/,\s?/,$1);  
           $results{'title'} = $2;  
           $results{'journal'} = $3;  
           $results{'year'} = $4;  
           $results{'volume'} = $5;  
           $results{'issue'} = $6;  
           $results{'pages'} = $7;  
           $results{'pubmed'} = $pubmed_id;  
           $results{'url'} = 'http://www.ncbi.nlm.nih.gov/pubmed/'.$pubmed_id;  
      } elsif ($data =~ /\d+\:\s*(.+?)\.(.+?)\.(.+?)\.(.+?)\s/){  
           @{$results{'authors'}} = split(/,\s?/,$1);  
           $results{'title'} = $2;  
           $results{'journal'} = $3;  
           $results{'year'} = $4;  
           $results{'pubmed'} = $pubmed_id;  
           $results{'url'} = 'http://www.ncbi.nlm.nih.gov/pubmed/'.$pubmed_id;  
      }  
 }  
   
 while (my ($key, $value) = each %results) {  
      if (ref($value) ne ARRAY){  
           printf("%s: %s\n", $key, $value);  
      } else {  
           printf("%s: %s\n", $key, join(', ',@$value));  
      }  
 }  
   

Para otro tipo de consultas automáticas (secuencias, artículos científicos...) recomiendo leer la documentación oficial del NCBI, donde se explica la utilidad de las diferentes Entrez Programming Utilities (E-utilities): ESearch, EFetch, ESummary, EPost y ELink.

Seguro que alguno de vosotros está pensando porqué no usar BioPerl o BioPython para estas consultas. La respuesta es porque quería explicar el fin de la era SOAP y la estructura de las nuevas consultas vía URL.

Aquí el código en BioPerl:

 use Bio::DB::EUtilities;  
   
 my $pubmed_id = '24234003';  
   
 my $eutil = 'efetch';  
 my $db = 'pubmed';  
 my $type = 'docsum';  
 my $mode = 'text';  
 my $email = 'mymail@foo.bar';  
    
 my $factory = Bio::DB::EUtilities->new(     -eutil  => $eutil,  
                          -db   => $db,  
                          -id   => $pubmed_id,  
                          -email  => $email,  
                          -rettype => $type,  
                          -retmode => $mode  
                     );  
    
 my $data = $factory->get_Response->content;  
   
 my %results;  
   
 if (defined($data)){  
      $data =~ s/\012\015?|\015\012?/ /g;  
      if ($data =~ /\d+\:\s*(.+?)\.(.+?)\.(.+?)\.(.+?)\s.+?\;(\d+)\((\d+)\)\:([\d\-]+)./){  
           @{$results{'authors'}} = split(/,\s?/,$1);  
           $results{'title'} = $2;  
           $results{'journal'} = $3;  
           $results{'year'} = $4;  
           $results{'volume'} = $5;  
           $results{'issue'} = $6;  
           $results{'pages'} = $7;  
           $results{'pubmed'} = $pubmed_id;  
           $results{'url'} = 'http://www.ncbi.nlm.nih.gov/pubmed/'.$pubmed_id;  
      } elsif ($data =~ /\d+\:\s*(.+?)\.(.+?)\.(.+?)\.(.+?)\s/){  
           @{$results{'authors'}} = split(/,\s?/,$1);  
           $results{'title'} = $2;  
           $results{'journal'} = $3;  
           $results{'year'} = $4;  
           $results{'pubmed'} = $pubmed_id;  
           $results{'url'} = 'http://www.ncbi.nlm.nih.gov/pubmed/'.$pubmed_id;  
      }  
 }  
   
 while (my ($key, $value) = each %results) {  
      if (ref($value) ne ARRAY){  
           printf("%s: %s\n", $key, $value);  
      } else {  
           printf("%s: %s\n", $key, join(', ',@$value));  
      }  
 }  
   

Y el código en BioPython:

 from Bio import Entrez  
   
 pubmed_id = '24234003';  
   
 handle = Entrez.efetch(db='pubmed', id=pubmed_id, rettype='docsum', retmode="xml", email='mymail@foo.bar')  
   
 results = Entrez.read(handle)  
   
 for key, value in results[0].iteritems():  
      print "%s: %s" % (key,value)