• Publicidad

Búsqueda de valores en rangos

Perl aplicado a la bioinformática

Búsqueda de valores en rangos

Notapor Dixan_Megaperl » 2011-09-06 03:44 @197

Hola,

necesito ayuda con una tarea que en principio parece sencilla, pero que dado el tamaño de los ficheros de datos que se manejan se me está volviendo compleja.

El caso es que tengo dos ficheros, uno de ellos con 22.500 líneas y en cada línea varias columnas de las que me interesan las dos primeras (cromosoma y posición). En otro fichero, de 180.000 líneas, tengo 3 columnas (cromosoma, inicio y fin). La cosa es: necesito por cada posición del primer fichero ver si está o no en alguno de los rangos del segundo.

He conseguido una solución "a lo bestia" que consiste en cargar los dos ficheros en sendos arrays y recorrer el primero. Por cada posición del primero me recorro el segundo para ver si está en algún rango (con alguna medida para no recorrerme el fichero entero, claro).

Ejemplo fichero de valores:
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
  1. #CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT 
  2. chr10   76210   .       A       G       88.8    .       DP=2;AF1=1;CI95=0.5,1;DP4=0,0,1,1;MQ=60;FQ=-33  GT:PL:GQ        1/1:120,6,0:10
  3. chr10   94238   .       A       T       28      .       DP=3;AF1=1;CI95=0.5,1;DP4=0,0,0,3;MQ=45;FQ=-36  GT:PL:GQ        1/1:60,9,0:16
  4. chr10   131577  .       A       T       4.77    .       DP=14;AF1=0.4999;CI95=0.5,0.5;DP4=6,6,1,1;MQ=60;FQ=6.99;PV4=1,4.9e-07,1,0.053   GT:PL:GQ        0/1:33,0,255:33
  5. chr10   146528  .       G       A       19.8    .       DP=2;AF1=1;CI95=0.5,1;DP4=0,0,1,1;MQ=60;FQ=-33  GT:PL:GQ        1/1:51,6,0:10
  6.  
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4


Ejemplo fichero de rangos:
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
  1. track:exome and plus captures          
  2. chr1    29775   30931
  3. chr1    68569   70529
  4. chr1    227733  228854
  5. chr1    227971  229211
  6. chr1    367147  369108
  7. chr1    470471  471830
  8. chr1    620584  622545
  9. chr1    740665  741785
  10.  
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4


Mi "solución":
Sintáxis: [ Descargar ] [ Ocultar ]
Using perl Syntax Highlighting
  1. #! /usr/bin/perl
  2.  
  3. my @ranges;
  4. my @values;
  5.  
  6. open( RANGES, $ARGV[0] ) or die "Can't open file: $!";
  7. while ( my $line = <RANGES> ) {
  8.     chomp $line;
  9.     my @val = split( "\t", $line );
  10.     push( @ranges, $val[0] . "_" . $val[1] . "_" . $val[2] );
  11. }
  12. close(RANGES);
  13.  
  14. open( VALUES, $ARGV[1] ) or die "Can't open the file $!";
  15. while ( my $line = <VALUES> ) {
  16.     chomp $line;
  17.     my @val = split( "\t", $line );
  18.     if ( $val[0] !~ "#" ) { $leer = 1; }
  19.     if ( $leer == 1 ) {
  20.         push( @values, $val[0] . "_" . $val[1] );
  21.     }
  22. }
  23. close(VALUES);
  24.  
  25. foreach $Val (@values) {
  26.     my @aux = split( "_", $Val );
  27.     foreach $range (@ranges) {
  28.         my @aux2 = split( "_", $range );
  29.         if ( $aux[1] >= $aux2[1] && $aux[1] <= $aux2[2] ) {
  30.             chomp $aux2[2];
  31.             print "$aux[1] esta en el rango [$aux2[1]-$aux2[2]] \n";
  32.             last;
  33.         }
  34.         elsif ( $aux[1] < $aux2[1] ) {
  35.             print "$aux[1] no esta en ningun rango\n";
  36.             last;
  37.         }
  38.     }
  39. }
  40.  
Coloreado en 0.003 segundos, usando GeSHi 1.0.8.4


Viendo el problema, ¿existe alguna solución más optimizada que la de recorrerse ambos arrays a piñón?

¡Muchas gracias!
Dixan_Megaperl
Perlero nuevo
Perlero nuevo
 
Mensajes: 2
Registrado: 2011-09-06 01:40 @111

Publicidad

Re: Búsqueda de valores en rangos

Notapor explorer » 2011-09-06 05:33 @273

Bienvenido a los foros de Perl en español, Dixan_Megaperl.

En este hilo comentamos otra solución.

Primero, leemos el fichero con los cromosomas y sus posiciones. Guardamos la información en una estructura de hash de array de array (por cada cromosoma tenemos un array, que guarda arrays con la posición e ID.

Segundo, hacemos una reordenación numérica (de menor a mayor) de las posiciones de los cromosomas.

Tercero, abrimos el fichero con los rangos. Por cada línea, primero comprobamos si tenemos información de ese cromosoma, y si es así, hacemos un bucle por todas las posiciones de ese cromosoma. Como están ordenadas numéricamente, es fácil saltar por las posiciones que son menores al inicio del rango. Una vez que encontramos una posición mayor que el inicio del rango, comprobamos si está antes del final del rango, y si es así, pintamos la información de la coincidencia, y terminamos el bucle.

Y es casi lo mismo a lo que tienes hecho, salvo que tu no miras por la coincidencia del nombre del cromosoma, sino que estás mirando todas las posiciones para todos los rangos.

Tu solución quizás se pueda acelerar, si eliminas los split(), y cambias la estructura de datos, de una dimensión a dos (array de array) (habría que lanzar un Benchmark, para comprobarlo).
JF^D Perl programming & Raku programming. Grupo en Telegram: https://t.me/Perl_ES
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14476
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Re: Búsqueda de valores en rangos

Notapor Dixan_Megaperl » 2011-09-06 07:40 @361

¡Muchas gracias!

Con la información del enlace he conseguido mejorar el tiempo de ejecución del programa. He bajado de 6 minutos a 1 y poco.
Dixan_Megaperl
Perlero nuevo
Perlero nuevo
 
Mensajes: 2
Registrado: 2011-09-06 01:40 @111


Volver a Bioinformática

¿Quién está conectado?

Usuarios navegando por este Foro: No hay usuarios registrados visitando el Foro y 4 invitados