Aumenta a velocidade do script Bash que usou grep em um loop while

1

fez com que esse script funcionasse contra um arquivo, composto por muita linha (> 500Mb) com este esquema:

odd lines: >BLA_BLA lenght_XX cov.XX
even lines: AGCAGCAGACTCAGACTACAGAT  # on even lines there's a DNA sequence

Sua função é recalcular valor após "cov." usando parâmetros passados por argumentos e substitua o mais antigo e calcule a quantidade percentual de "G" e "C" no DNA seq, em linhas pares.

Então, a saída é semelhante a:

> BLA_BLA lenght_XX
> nucleotidic_cov XX
> DNA seq (the same of even lines)
> GC_CONT: XX

Aqui está o código (apenas o loop):

K=$(($READLENGHT - $KMER + 1))
Y=$(echo "scale=4; $K / $READLENGHT" | bc)

while read odd; do
    echo -n "${odd##}" | cut -d "_" -f 1,2,3,4 && printf "nucleotide_cov: " 
    echo "scale=4;${odd##*_} / $Y" | bc 
    read even
    echo "${even##}" &&
    ACOUNT=$(echo "${even##}" |  sed -e "s/./&\n /g" | grep -c "A")  
    GCOUNT=$(echo "${even##}" |  sed -e "s/./&\n /g" | grep -c "G")
    CCOUNT=$(echo "${even##}" |  sed -e "s/./&\n /g" | grep -c "C")
    TCOUNT=$(echo "${even##}" |  sed -e "s/./&\n /g" | grep -c "T")
    TOTALBASES=$(($ACOUNT+$GCOUNT+$CCOUNT+$TCOUNT))
    GCCONT=$(($GCOUNT+$CCOUNT))
    printf "GC_CONT: " 
    echo "scale=2;$GCCONT / $TOTALBASES *100" | bc  
done < "$1"

É incrivelmente lento quando é executado em um arquivo de texto enorme (maior que 500Mb) em um servidor de 16 núcleos. Alguma idéia de como aumentar a velocidade desse script?

EDITAR

Conforme solicitado, a E / S desejada é fornecida via pastebin: link

    
por Shred 27.02.2018 / 10:13

4 respostas

2

O cálculo da porcentagem pode ser reduzido a uma única operação como esta

 echo "${even##}" | awk '{x=gsub(/[ACT]/,""); y=gsub(/G/,""); printf "GC_CONT : %.2f%%\b", (y*100)/(x+y) }'

O gsub substitui um padrão e retorna a contagem de substituições que ele fez. Então, isso pode ser usado para calcular rapidamente a porcentagem.

Você também pode processar as linhas pares e ímpares no awk. Não está claro o que você está fazendo com linhas estranhas, mas sua função completa pode ser colocada em um único awk -

awk -F '_' -v Y="$Y" '{ if(NR%2==1) {
    printf "%s %s %s %s %s\nnucleotidic_cov : %.4f\n",$1,$2,$3,$4,$5, ($6 / Y)
} else {
    x=gsub(/[AT]/,""); 
    y=gsub(/[GC]/,""); 
    printf "GC_CONT : %.2f%%\n", (y*100)/(x+y)
    }
 }' large_file

EDIT : Com base no requisito do OP, o bloco if foi alterado para linhas ímpares. O gsub removeria o "cov". do número. Depois de passar a variável de shell $ Y para awk, agora podemos dividir e imprimir no formato necessário.

Usar um único script awk em vez de várias operações acelerará significativamente a operação.

    
por 27.02.2018 / 11:09
3

Você chegou (para dizer o mínimo) aos limites do que pode ser razoavelmente feito no shell - você deve reescrever seu script em algo como AWK, ou Perl ou Python. Usando uma linguagem mais avançada como essas evitará ter que executar vários processos para todo o seu processamento de texto; você poderá fazer isso usando funções internas.

    
por 27.02.2018 / 10:22
1

O número de núcleos pouco importa se o programa não estiver em paralelo (muito).

Você pode usar wc e tr em vez de sed e grep, o que pode acelerar um pouco as coisas:

ACOUNT=$(echo "${even##}" | tr -d [^A] | wc -m)

Mas, na verdade, acho que o maior problema é que o shell, embora seja fácil de programar para trabalhos rápidos e sujos, não é a ferramenta certa para o trabalho quando se trata de poder de processamento bruto. Eu sugeriria uma linguagem de programação mais sofisticada, como Perl ou Python, que também possui habilidades de threading (permitindo que você use todos os seus núcleos).

Você poderia fazer isso em perl um pouco assim:

#!/usr/bin/perl -w
use strict;
use warnings;

my $y = ...;                              # calculate your Y value here
while(my $odd = <ARGV>) {                 # Read a line from the file(s) passed
                                          # on the command line
    chomp $odd;                           # lose the newline
    my @split = split /_/, $odd;          # split the read line on a "_" boundary
                                          # into an array
    print join("_", @split[0..3]) . "\n"; # print the first four elements of the
                                          # array, separated by "_"
    print $split[$#split] / $y . "\n";    # Treat the final element of the
                                          # @split array as a number, divide it
                                          # by $y, and output the result
    my %charcount = (                     # Initialize a hash table
        A => 0,
        G => 0,
        C => 0,
        T => 0
    );
    my $even = <ARGV>;                    # read the even line
    chomp $even;
    foreach my $char(split //,$even) {    # split the string into separate
                                          # characters, and loop over them
        $charcount{$char}++;              # Count the correct character
    }
    my $total = $charcount{A} + $charcount{G} + $charcount{C} + $charcount{T};
    my $gc = $charcount{G} + $charcount{C};
    my $perc = $gc / $total;
    print "GC_CONT: $perc\n";             # Do our final calculations and
                                          # output the result
}

Nota: não testado (além de "o perl aceita este código")

Se você quiser saber mais sobre o perl, execute perldoc perlintro e comece; -)

    
por 27.02.2018 / 10:34
0

Você está lendo um longo arquivo linha por linha e executando vários comandos em cada iteração. O principal problema que você está enfrentando é a latência de executar esses cálculos e ler pedaços muito pequenos do arquivo por vez.

A resposta de Stephen Kitt é boa, você quer reescrever isso em uma linguagem de nível mais alto na qual você pode armazenar em cache o conteúdo do arquivo e executar suas operações de string com muito mais eficiência.

Se você quiser excluir o desempenho do armazenamento e do sistema de arquivos, é possível carregar o arquivo da RAM usando:

# mkdir /mnt/tmpfs
# mount -t tmpfs -o size=1024M tmpfs /mnt/tmpfs
# cp <input_file> /tmp/tmpfs
# <script> /tmp/tmpfs/<input_file>

Isso deve tornar o processo mais rápido, tanto quanto você está restrito a E / S. Mas nunca será tão bom quanto poderia ser se reescrito em C ou ruby ou python.

    
por 27.02.2018 / 10:31