улучшить скорость разбора fastq

Question

@solved C# с тем же кодом в два раза быстрее

я разборе файла fastq phred33 в Perl и занимает значительное количество времени (порядка 15 минут). Файл fastq составляет около 3 концертов. Есть ли разумные способы сделать это быстрее?

$file=shift; 
open(FILE,$file); 
open(FILEFA,">".$file.".fa"); 
open(FILEQA,">".$file.".qual"); 
while($line=<FILE>) 
{ 
    chomp($line); 
    if($line=~m/^@/) 
    { 


    $header=$line; 
    $header =~ s/@/>/g; 
    $seq=<FILE>; 
    chomp($seq); 
    $nothing=<FILE>; 
    $nothing=""; 
    $fastq=<FILE>; 

    print FILEFA $header."\n"; 
    print FILEFA $seq."\n"; 
    $seq=""; 
    print FILEQA $header."\n"; 

     @elm=split("",$fastq); 
     $i=0; 
     while(defined($elm[$i])) 
     { 
      $Q = ord($elm[$i]) - 33; 
      if($Q!="-23") 
      { 
      print FILEQA $Q." "; 
      } 
      $i=$i+1; 
     } 
     print FILEQA "\n"; 
    } 
} 
print $file.".fa\n"; 
print $file.".qual\n"; 

Answer 1

Здесь не используется процессор. он привязан к IO, поэтому в основном это время для чтения через 3 ГБ. Существуют микрооптимизации (и другие очистки), которые можно выполнить.

Во-первых, всегда используйте use strict; use warnings;.

Основной код

my @elm = split(//, $fastq); 
my $i=0; 
while(defined($elm[$i])) { 
    my $Q = ord($elm[$i]) - 33; 
    if($Q!="-23") { 
     print FILEQA $Q." "; 
    } 
    $i=$i+1; 
} 

Цель if($Q!="-23"), чтобы проверить, если символ новой строки, которые вы не должны делать, если вы сделали chomp($fastq);. (Что с кавычками -23?!)

chomp($fastq); 
my @elm = split(//, $fastq); 
my $i=0; 
while(defined($elm[$i])) { 
    my $Q = ord($elm[$i]) - 33; 
    print FILEQA $Q." "; 
    $i=$i+1; 
} 
print FILEQA "\n"; 

Использование while цикл только усложняет ситуацию. Используйте цикл for, когда у вас есть известное количество итераций.

chomp($fastq); 
for (split(//, $fastq)) { 
    print FILEQA (ord($_)-33)." "; 
} 
print FILEQA "\n"; 

Это может помочь немного повернуть это наизнанку.

chomp($fastq); 
print FILEQA join(' ', map ord($_)-33, split //, $fastq), "\n"; 

На второй мысли, а не наизнанку достаточно :)

$fastq =~ s/(.)/(ord($1)-33) . " "/eg; 
print FILEQA $fastq; 

Но что мы предварительно вычислены переводы? Тогда нам не пришлось бы повторно называть sub (код /e).

my %map = map { chr($_) => ($_-33)." " } 0x00..0xFF; 

$fastq =~ s/(.)/$map{$1}/g; 
print FILEQA $fastq; 

---

После немного больше очистки, мы получаем:

use strict; 
use warnings; 

my %map = map { chr($_) => ($_-33)." " } 0x00..0xFF; 

my $file = shift; 

my $fa_file = "$file.fa"; 
my $qual_file = "$file.qual"; 

open(my $FILE, '<', $file ) or die $!; 
open(my $FILEFA, '>', $fa_file ) or die $!; 
open(my $FILEQA, '>', $qual_file) or die $!; 

while (my $header = <$FILE>) { 
    next if $header !~ /^@/; 

    my $seq = <$FILE>; 
    <$FILE>; 
    my $fastq = <$FILE>; 

    $header =~ s/@/>/g; 
    $fastq =~ s/(.)/$map{$1}/g; 

    print $FILEFA $header; 
    print $FILEFA $seq; 

    print $FILEQA $header; 
    print $FILEQA $fastq; 
} 

print "$fa_file\n"; 
print "$qual_file\n";