30 de noviembre de 2012

comparando secuencias de igual longitud

Hola,
hoy quería compartir recetas en Perl para la comparación eficiente de secuencias de igual longitud. Las dos funciones que vamos a ver devuelven como resultado listas con las posiciones (contando desde 0) en las que difieren dos secuencias.

Para qué sirve esto? En general, para calcular distancias de edición (edit distances); en bioinformática es una manera, por ejemplo, de comparar oligos de igual longitud o secuencias extraídas de un alineamiento múltilple, incluyendo gaps. El siguiente código fuente incluye una subrutina en Perl (con el operador binario ^ y la función pos) y otra en C embebido para hacer la misma operación por dos caminos:

 use strict;  
 use warnings;  
 my $seq1 = 'ACTGGA';  
 my $seq2 = 'AGTG-A';  
   
 my @Pdiffs = find_diffs($seq1,$seq2);  
 my @Cdiffs = find_diffs_C($seq1,$seq2);  
   
 printf("# version Perl\ndiferencias = %d\n%s\n",scalar(@Pdiffs),join(',',@Pdiffs));  
 printf("# version C\ndiferencias = %d\n%s\n",scalar(@Cdiffs),join(',',@Cdiffs));  
   
 sub find_diffs  
 {  
    my ($seqA,$seqB) = @_;  
    my $XOR = $seqA ^ $seqB; 
    my @diff;   
    while($XOR =~ /([^\0])/g)  
    {   
       push(@diff,pos($XOR)-1);  
    }  
    return @diff;  
 }  
   
 use Inline C => <<'END_C';  
 void find_diffs_C(char* x, char* y)   
 {                        
   int i;   
    Inline_Stack_Vars;                               
   Inline_Stack_Reset;                                                                
   for(i=0; x[i] && y[i]; ++i) {                          
    if(x[i] != y[i]) {                              
     Inline_Stack_Push(sv_2mortal(newSViv(i)));                 
    }                                       
   }   
    Inline_Stack_Done;                                                                  
 }     
 END_C  
Podéis ver en stackoverflow otras versiones y su comparación en base a su tiempo de cálculo,
un saludo,
Bruno




21 de noviembre de 2012

Comprimiendo BLAST

Hola,
ya hemos hablado antes en este blog de las crecientes aplicaciones de los algoritmos de compresión en la bioinformática. Hoy precisamente quería enlazar a dos artículos recientes que los explotan para acelerar algo que a priori parece imposible, la herramienta más universal de la biología computacional, BLAST.

En un paper en Nature Biotech, Loh, Maym y Berger nos presentan el prototipo Compression-accelerated BLAST (fuente C++ aquí), que en pruebas empíricas acelera varios órdenes de magnitud las búsquedas de BLAST simplemente haciendo las operaciones en un espacio comprimido. Como ventaja añadida se generan archivos de resultados de tamaños significativamente menores. Obviamente pagamos una pequeña pérdida de sensibilidad, pero puede ser perfectamente asumible para tareas de mapeo de lecturas (reads) de secuenciación de última generación:

Original de http://www.nature.com/nbt/journal/v30/n7/full/nbt.2241.html.

En otro paper recientemente publicado en Bioinformatics, Koskinen y Holm, aplican el marco de los vectores de sufijos (ya discutidos aquí) a la búsqueda de proteínas homólogas con una identidad superior al 50%, acelerando de nuevo varias órdenes de magnitud por encima de BLAST.

Familias de algoritmos para la búsqueda inexacta de secuencias, según Koskinen y Hol (http://bioinformatics.oxfordjournals.org/content/28/18/i438.full).
El programa que implementa estas ideas se llama SANS, escrito en FORTRAN 90, está disponible en este enlace,
un saludo,
Bruno