• Publicidad

Mejorar script para que sea más rápido

Perl aplicado a la bioinformática

Mejorar script para que sea más rápido

Notapor damiancambiagno » 2015-08-01 06:18 @304

Estimados:

Soy novato en esto, y ando teniendo algunos problemitas.

Yo he adaptado un script en Perl para comparar dos archivos de unos análisis bioinformáticos y me está dando mucho dolor de cabeza porque es muy lento.

Mi pregunta es si alguien sabe cómo podría escribirlo para que funcione más rápido.


Básicamente lo que tengo que hacer es lo siguiente: Tengo dos archivo que tienen tres y cuatro columnas.
El primer archivo es de esta forma:
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
     Chr   posición   cuentas
      1      100       0.5
      1      300       0.1
      2      200       0.3
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4

y el segundo, esta:
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
     gen     Chr      inicio   fin
   AT1G01020  1        50      350
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4


Lo que yo quiero hacer es que cuando encuentre una fila del archivo 2º donde Chr coincida con el del archivo 1º, y la posición del archivo 1º esté entre el inicio y fin del archivo 2º, me imprima en un archivo de salida lo siguiente:
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
  gen    Chr    cuentas    posición
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4


Para eso adapté este script que me dieron pero resulta muy lento. ¿Alguien me podría ayudar con esto?

Desde ya, ¡muchas gracias!

Sintáxis: [ Descargar ] [ Ocultar ]
Using perl Syntax Highlighting
  1. my $input1 = "/home/damian/Escritorio/Arrays/archivo1";
  2. my $input2 = "/home/damian/Escritorio/Arrays/archivo2";
  3. my $output = "/home/damian/Escritorio/Arrays/salida";
  4. open( FILE1,      "<$input1.txt" )  or die "could not open $input1 file";
  5. open( OUTPUTFILE, ">>$output.txt" ) or die "could not $output.txt output file";
  6. while (<FILE1>) {    ##       1      100       0.5
  7.     my @split1    = split( /\s/, $_ );
  8.     my $Chr1      = $split1[0];
  9.     my $posicion1 = $split1[1];
  10.     my $cuentas   = $split1[2];
  11.     open( FILE2, "<$input2.txt" ) or die "could not open $input2 file";
  12.     while (<FILE2>) {    ##    AT1G01020  1            50      350
  13.         my @split2 = split( /\s/, $_ );
  14.         my $Gen    = $split2[0];
  15.         my $Chr2   = $split2[1];
  16.         my $incio  = $split2[2];
  17.         my $fin    = $split2[3];
  18.         if (    ( $Chr1 eq $Chr2 )
  19.             and ( $posicion1 > $incio )
  20.             and ( $posicion1 < $fin ) )
  21.         {
  22.             print OUTPUTFILE "$Gen\t$Chr2\t$cuentas\t$posicion1\n";
  23.             print "$Gen\t$Chr2\t$cuentas\t$posicion1\n";
  24.  
  25.         }
  26.     }
  27. }
  28. exit;
Coloreado en 0.017 segundos, usando GeSHi 1.0.8.4
damiancambiagno
Perlero nuevo
Perlero nuevo
 
Mensajes: 6
Registrado: 2015-07-31 16:57 @748

Publicidad

Re: Mejorar script para que sea más rápido

Notapor explorer » 2015-08-01 14:37 @651

Bienvenido a los foros de Perl en Español, damiancambiagno.

El problema de la lentitud es que estás leyendo todas las líneas del archivo 2º por cada línea del archivo 1º, así que estamos haciendo

<líneas archivo 1> + <líneas archivo 1> * <líneas archivo 2> lecturas, que pueden ser muchas, comparadas con las ideales:

<líneas archivo 1> + <líneas archivo 2>


La clave está en saber si alguno de los campos de estos archivos son claves. Por ejemplo, ¿la columna Chr del archivo 2º es única? ¿No se repite a lo largo del archivo, o puede repetirse?

Si no se repite, entonces podemos leer todas las líneas del archivo 2 y meterlas en un hash, y luego recorrer las líneas del archivo 1, procesando una por una.

Por este subforo de Bioinformática hay algunos casos parecidos. Dinos si ese campo es clave y te echaremos una mano.
JF^D Perl Programming Language
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14102
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Re: Mejorar script para que sea más rápido

Notapor damiancambiagno » 2015-08-01 17:31 @771

¡Muchas gracias por la pronta respuesta!

Claro, no especifiqué mucho de cómo son los archivos. Son así:

En el archivo 1, la columna "Chr" tiene 5 valores posibles (1 a 5) y de cada uno debe haber un millón de lineas.
De la columna "posición" son valores únicos para cada "Chr".

Para el archivo 2, los valores de Chr también pueden ser 5 (de 1 a 5) pero hay unas 6000 líneas de cada uno.
Y en el caso de los valores de las columnas inicio y fin van variando en cada fila.
Lo único que no se repite son los nombres de la columna gen.

No sé si sirve como dato, pero la combinación "Chr", "inicio" y "fin", es única.

¿Hay alguna forma entonces de acelerarlo?

Nuevamente, ¡muchas gracias!
damiancambiagno
Perlero nuevo
Perlero nuevo
 
Mensajes: 6
Registrado: 2015-07-31 16:57 @748

Re: Mejorar script para que sea más rápido

Notapor explorer » 2015-08-01 21:16 @928

¿Puede haber valores de Fin para las mismas combinaciones de Chr e Inicio?

Has dicho que los tres juntos forman combinaciones únicas pero, ¿podría darse este caso?
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
     gen     Chr      inicio   fin
   AT1G01020  1        50      350
   AT1G01021  1        50      500
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4
JF^D Perl Programming Language
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14102
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Re: Mejorar script para que sea más rápido

Notapor damiancambiagno » 2015-08-02 09:16 @427

Sí puede haber, o sea, puede ser que para una misma línea compartan el valor de Chr e Inicio pero no Fin, o el mismo valor de Chr y Fin pero no Inicio... Sería así.

Es decir, al menos una de las 3 va a variar entre líneas, pero cada combinación de las 3 es única.
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
gen  chr  inicio   fin
X     1     50     300
X1    1     50     290
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4

o también
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
X2    1     30     600
X3    1     50     600
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4

incluso
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
X4    1     200    400
X5    2     200    400
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4


¿Estoy en problemas?
damiancambiagno
Perlero nuevo
Perlero nuevo
 
Mensajes: 6
Registrado: 2015-07-31 16:57 @748

Re: Mejorar script para que sea más rápido

Notapor explorer » 2015-08-02 18:52 @828

Sí que veo problemas. Miro varias opciones, y no se me ocurre la idea feliz.

Hay una forma. Se trataría de leer el segundo archivo, y organizar la información en una estructura de varias dimensiones, que facilite la localización del gen, que es al final lo que nos interesa.

Usando el campo Chr como clave, se reduce el espacio de búsqueda a 1/5. Dentro de la clave que nos interesa, tenemos que hacer una búsqueda entre los valores Inicio que nos interesa (los que son mayores o iguales al valor Posición). Y dentro de ellos, buscar los valores finales que nos interesa (los que son menores o iguales a Posición).

En fin, se trata de una estructura tridimensional, compuesta de:

Primera dimensión: un hash cuyas claves son los valores posibles de Chr.
Segunda dimensión: un array que contiene los valores Inicio, asociadas al valor de la primera dimensión, ordenados numéricamente.
Tercera dimensión: un array que contiene los valores Final, asociados a los valores Inicio de la segunda dimensión, ordenados numéricamente.

Como vemos en los ejemplos, los valores de Chr son también numéricos, así que podemos pensar que la primera dimensión puede ser también un array, por simplificar.

Una vez construida la estructura, recorremos el primer archivo por todas sus líneas, y buscamos por la estructura qué valores son coincidentes.

Pero... puede haber mucha separación entre elementos. Los bucles pueden volverse muy lentos. Y el número de condicionales es casi el mismo, así que otra opción es la de obviar la ordenación de los elementos, y simplemente, buscar por todas las combinaciones.

Usaremos hash, porque así la construcción de la estructura será muy fácil.

Supongamos un archivo 1 así:
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
      1      100       0.5
      1      300       0.1
      2      200       0.3
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4

y un archivo 2 así:
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
AT1G01020  1    50      350
AT1G01021  1    50      500
AT1G010X0  1    50     300
AT1G010X1  1    50     290
AT1G010X2  1    30     600
AT1G010X3  1    50     600
AT1G010X4  1    200    400
AT1G010X5  2    200    400
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4

con este programa:
Sintáxis: [ Descargar ] [ Ocultar ]
Using perl Syntax Highlighting
  1. #!/usr/bin/perl
  2. use feature 'say';
  3. use autodie;
  4.  
  5. #my $dir    = '/home/damian/Escritorio/Arrays/';
  6. my $dir         = '.';
  7. my $input1      = "$dir/archivo1.txt";
  8. my $input2      = "$dir/archivo2.txt";
  9. my $output      = "$dir/salida.txt";
  10.  
  11. # Procesar el segundo archivo para crear la estructura
  12. my @matrix;
  13.  
  14. open my $INPUT2, '<', $input2;
  15.  
  16. while (<$INPUT2>) {
  17.     ##    AT1G01020  1       50         350
  18.     my(   $gen,      $chr,   $inicio,   $final) = split;
  19.    
  20.     $matrix[$chr]{$inicio}{$final} = $gen;                      # "Matrix te posee"
  21. }
  22.  
  23. close $INPUT2;
  24.  
  25. #use Data::Dumper;
  26. #say Dumper \@matrix;
  27.  
  28. # Procesar el primer archivo, por líneas
  29. open my $INPUT1, '<', $input1;
  30. open my $SALIDA, '>', $output;
  31.  
  32. while (<$INPUT1>) {
  33.     chomp;
  34.    
  35.     ##    1      100        0.5
  36.     my(   $chr,  $posicion, $cuentas) = split;
  37.    
  38.     #say "$chr, $posicion, $cuentas";
  39.    
  40.     for my $inicio ( keys %{ $matrix[$chr] } ) {                # recorremos todos los valores de inicio
  41.    
  42.         next if $inicio > $posicion;                            # saltar al siguiente si $posicion no está dentro
  43.        
  44.         for my $final ( keys %{ $matrix[$chr]{$inicio} } ) {    # recorremos los finales asociados a ese $inicio
  45.        
  46.             next if $final < $posicion;                         # saltar al siguiente si $posición no está dentro
  47.        
  48.             # gen    Chr    cuentas    posición
  49.             say $SALIDA join "\t", $matrix[$chr]{$inicio}{$final}, $chr, $cuentas, $posicion
  50.         }
  51.     }
  52. }
  53.  
  54. close $INPUT1;
  55. close $SALIDA;
  56.  
  57. __END__
Coloreado en 0.019 segundos, usando GeSHi 1.0.8.4
la salida es:
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
AT1G01020       1       0.5     100
AT1G01021       1       0.5     100
AT1G010X3       1       0.5     100
AT1G010X0       1       0.5     100
AT1G010X1       1       0.5     100
AT1G010X2       1       0.5     100
AT1G01020       1       0.1     300
AT1G01021       1       0.1     300
AT1G010X3       1       0.1     300
AT1G010X0       1       0.1     300
AT1G010X4       1       0.1     300
AT1G010X2       1       0.1     300
AT1G010X5       2       0.3     200
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4

Si descomentas las líneas de la salida del Data::Dumper, se muestra el formato que tiene la estructura:
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
$VAR1 = [
          undef,
          {
            '50' => {
                      '350' => 'AT1G01020',
                      '500' => 'AT1G01021',
                      '600' => 'AT1G010X3',
                      '300' => 'AT1G010X0',
                      '290' => 'AT1G010X1'
                    },
            '200' => {
                       '400' => 'AT1G010X4'
                     },
            '30' => {
                      '600' => 'AT1G010X2'
                    }
          },
          {
            '200' => {
                       '400' => 'AT1G010X5'
                     }
          }
        ];
Coloreado en 0.001 segundos, usando GeSHi 1.0.8.4

Prueba a ver. El tiempo consumido debería reducirse a un teórico 1/5 del anterior.
JF^D Perl Programming Language
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14102
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Re: Mejorar script para que sea más rápido

Notapor damiancambiagno » 2015-08-03 17:20 @763

Increíble... Esto me gusta cada vez más... ¡¡¡Lo pruebo y te aviso!!!

¡¡Mil gracias!!
damiancambiagno
Perlero nuevo
Perlero nuevo
 
Mensajes: 6
Registrado: 2015-07-31 16:57 @748

Re: Mejorar script para que sea más rápido

Notapor damiancambiagno » 2015-08-04 20:15 @885

¡Funcionó muy bien! Ahora es totalmente factible de correr en un tiempo razonable.

¡¡Muchísimas gracias, explorer!!
damiancambiagno
Perlero nuevo
Perlero nuevo
 
Mensajes: 6
Registrado: 2015-07-31 16:57 @748

Re: Mejorar script para que sea más rápido

Notapor explorer » 2015-08-05 08:01 @375

Pues... el caso es que se puede acelerar aún más, a costa de complicar un poco más el código, claro.

Consistiría en reconvertir la estructura para que la segunda dimensión estuviera ordenada numéricamente (pasar de hash a array), y así podríamos abortar el bucle con un last en cuanto viéramos que tenemos inicios superiores al valor que tenemos.

Sintáxis: [ Descargar ] [ Ocultar ]
Using perl Syntax Highlighting
  1. #!/usr/bin/perl
  2. use feature 'say';
  3. use autodie;
  4.  
  5. #my $dir    = '/home/damian/Escritorio/Arrays/';
  6. my $dir         = '.';
  7. my $input1      = "$dir/archivo1.txt";
  8. my $input2      = "$dir/archivo2.txt";
  9. my $output      = "$dir/salida.txt";
  10.  
  11. # Procesar el segundo archivo para crear la estructura
  12. my @matrix;
  13.  
  14. open my $INPUT2, '<', $input2;
  15.  
  16. while (<$INPUT2>) {
  17.     ##    AT1G01020  1            50      350
  18.     my($gen,$chr,$inicio,$final) = split;
  19.    
  20.     $matrix[$chr]{$inicio}{$final} = $gen;
  21. }
  22.  
  23. close $INPUT2;
  24.  
  25. #use Data::Dumper;
  26. #say Dumper \@matrix;
  27.  
  28. # Ordenación de los valores de inicio,
  29. # con la conversión a array de hashes
  30. for my $chr ( 1 .. $#matrix ) {
  31.  
  32.     my @inicios =                                               # creamos un array cuyos valores
  33.                 sort { $a->[0] <=> $b->[0] }                    # estarán ordenados numéricamente
  34.                 map { [ $_, $matrix[$chr]{$_} ] }
  35.                 keys %{ $matrix[$chr] }
  36.                 ;
  37.  
  38.     $matrix[$chr] = \@inicios;
  39. }
  40.  
  41. #say Dumper \@matrix;
  42.  
  43. # Procesar el primer archivo, por líneas
  44. open my $INPUT1, '<', $input1;
  45. open my $SALIDA, '>', $output;
  46.  
  47. while (<$INPUT1>) {
  48.     chomp;
  49.    
  50.     my($chr, $posicion, $cuentas) = split;
  51.    
  52.     #say "$chr, $posicion, $cuentas";
  53.  
  54.     for my $inicio_ref ( @{ $matrix[$chr] } ) {                 # recorremos todos los valores de inicio
  55.         my($inicio,$final_ref) = @{ $inicio_ref };
  56.  
  57.         last if $inicio > $posicion;                            # salimos si $posicion ya no está dentro
  58.  
  59.         for my $final ( keys %{ $final_ref } ) {                # recorremos los finales asociados a ese $inicio
  60.        
  61.             next if $final < $posicion;                         # saltar al siguiente si $posición no está dentro
  62.        
  63.             # gen    Chr    cuentas    posición
  64.             say $SALIDA join "\t", ${ $final_ref }{$final}, $chr, $cuentas, $posicion
  65.         }
  66.     }
  67. }
  68.  
  69. close $INPUT1;
  70. close $SALIDA;
Coloreado en 0.019 segundos, usando GeSHi 1.0.8.4
JF^D Perl Programming Language
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14102
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Re: Mejorar script para que sea más rápido

Notapor damiancambiagno » 2015-08-06 15:55 @705

¡¡Aplausos!! Veo que tengo que ponerme a estudiar más Perl.

La verdad que mejoró muchísimo el script.

Con mi computadora (que es una Notebook i3 muy básica) el primer script iba a andar durante 140 horas.
Con el segundo que tenía el hash y los arrays se redujo a solo 10 horas.
Y con el último a ¡¡6 horas!!

Eso me viene muy bien ya que lo tengo que correr unas 15 veces, así que ahora es totalmente factible.

¡¡¡Muchísimas gracias!!!
damiancambiagno
Perlero nuevo
Perlero nuevo
 
Mensajes: 6
Registrado: 2015-07-31 16:57 @748

Siguiente

Volver a Bioinformática

¿Quién está conectado?

Usuarios navegando por este Foro: No hay usuarios registrados visitando el Foro y 1 invitado