Der schnellste Weg, um einen String mit zwei Unstimmigkeiten in großen Dateien in Perl zu suchen?

Question

Ich habe den Code geändert, um mit zwei Dateien zu arbeiten. to_search.txt muss durchsucht werden. big_file.fastq hat Zeilen, die durchsucht werden sollen, und wenn eine Zeichenfolge gefunden wird (2 Nichtübereinstimmungen mit exakter Länge zulässig, die von 8-10 reichen, kein Hinzufügen und Löschen), in den entsprechenden Namen einfügen. So wird jede Zeichenfolge in allen Zeilen (2. Zeile) in big_file.fastq gesucht.

# to_search.txt: (length can be from 8-20 characters) 

1  TCCCTTGT 
2  ACGAGACT 
3  GCTGTACG 
4  ATCACCAG 
5  TGGTCAAC 
6  ATCGCACA 
7  GTCGTGTA 
8  AGCGGAGG 
9  ATCCTTTG 
10  TACAGCGC 
#2000 search needed 

# big_file.fastq: 2 billions lines (each 4 lines are associated: string search is in second line of each 4 lines). 
# Second line can have 100-200 characters 
@M04398:19:000000000-APDK3:1:1101:21860:1000 1:N:0:1 
TCttTTGTGATCGATCGATCGATCGATCGGTCGTGTAGCCTCCAACCAAGCACCCCATCTGTTCCAAATCTTCTCCCACTGCTACTTGAAGACGCTGAAGTTGAAGGGCCACCTTCATCATTCTGG 
+ 
#8ACCDGGGGGGGGGGGGGEGGGGGGGGDFFEGGG@FFGGGGGGGGGGGGGGGGGCF@<FFGGGGGFE9FGGGFEACE8EFFFGGGGGF9F?CECEFEG,CFGF,CF@CCC?BFFG<,9<9AEFG,, 
@M04398:19:000000000-APDK3:1:1101:13382:1000 1:N:0:1 
NTCGATCGATCGATCGATCGATCGATCGTTCTGAGAGGTACCAACCAAGCACACCACGGGCGACACAGACAGCTCCGTGTTGAACGGGTTGTTCTTCTTCTTGCCTTCATCATCCCCATCCTCAGTGGACGCAGCTTGCTCATCCTTCCTC 
+ 
#8BCCGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG 
@M04398:19:000000000-APDK3:1:1101:18888:1000 1:N:0:1 
NCAGAATGAGGAAGGATGAGCCCCGTCGTGTCGAAGCTATTGACACAGCGCTATTCCGTCTTTATGTTCACTTTAAGCGGTACAAGGAGCTGCTTGTTCTGATTCAGGAACCGAACCCTGGTGGTGTGCTTGGTTGGCAAGTTTACGGCTC 
+ 
#8BCCGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGCGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGE 

Hier ist der Code für zwei Fehlpaarungen. Ich habe versucht mit genauem Match, Geschwindigkeit ist nicht schlecht: dauert ungefähr einen Tag. Ich habe das Time :: Progress-Modul verwendet. Wenn ich 2 Mismatch verwende: zeigt 115 Tage bis zum Ende. Wie kann die Geschwindigkeit hier verbessert werden?

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

$| = 1; 

open(IN_P1, "big_file.fastq") or die "File not found"; 

my (@sample_file_names, @barcode1); 
open(BC_FILE, "to_search.txt") or die "No barcode file"; 
my @barcode_file_content = <BC_FILE>; 

foreach (@barcode_file_content) { 
    chomp $_; 
    $_ =~ s/\r//; 
    $_ =~ s/\n//; 

    #print $_; 
    my @elements = split("(\t|,|)", $_); 
    push @sample_file_names, $elements[0]; 
    push @barcode1, $elements[2]; 
} 

# open FH 
my @fh_array_R1; 
foreach (@sample_file_names) { 
    chomp $_; 
    local *OUT_R1; 
    open(OUT_R1, ">", "$_\.fq") or die "cannot write file"; 
    push @fh_array_R1, *OUT_R1; 
} 

# unknown barcode file 
open(UNKNOWN_R1, ">unknown-barcode_SE.fq") or die "cannot create unknown-r1 file"; 

while (defined(my $firstp1 = <IN_P1>)) { 

    my $p1_first_line = $firstp1; 
    my $p1_second_line = <IN_P1>; 
    my $p1_third_line = <IN_P1>; 
    my $p1_fourth_line = <IN_P1>; 

    chomp($p1_first_line, $p1_second_line, $p1_third_line, $p1_fourth_line,); 

    my $matched_R1 = "$p1_first_line\n$p1_second_line\n$p1_third_line\n$p1_fourth_line\n"; 

    for (my $j = 0 ; $j < scalar @barcode1 ; $j++) { 
     chomp $barcode1[$j]; 

     my $barcode1_regex = make_barcode_fragments($barcode1[$j]); 

     if ($p1_second_line =~ /$barcode1_regex/i) { 
      # keep if matched 
      print { $fh_array_R1[$j] } $matched_R1; 
      last; 
     } 
     else { 
      #print to unknown; 
      print UNKNOWN_R1 $matched_R1; 
     } 
    } 
} 

# make two mismatch patterm of barcode 
sub make_barcode_fragments { 
    my ($in1) = @_; 
    my @subpats; 
    for my $i (0 .. length($in1) - 1) { 
     for my $j ($i + 1 .. length($in1) - 1) { 
      my $subpat = join('', 
       substr($in1, 0, $i), 
       '\\w', substr($in1, $i + 1, $j - $i - 1), 
       '\\w', substr($in1, $j + 1), 
      ); 
      push @subpats, $subpat; 
     } 
    } 
    my $pat = join('|', @subpats); 

    #print $pat; 
    return "$pat"; 
} 
exit; 

Answer 1

Wenn Ihr Algorithmus nicht geändert werden kann/sich in Perl verbessert, können Sie immer noch Speedup erhalten, indem Sie die Zeit mit dem Schreiben Teile in C raubend Hier ist ein Beispiel unter Verwendung von Inline-C:

use strict; 
use warnings; 
use Benchmark qw(timethese); 
use Inline C => './check_line_c.c'; 

my $find = "MATCH1"; 
my $search = "saasdadadadadasd"; 

my %sub_info = (
    c => sub { check_line_c($find, $search) }, 
    perl => sub { check_line_perl($find, $search) }, 
); 

timethese(4_000_000, \%sub_info); 

sub check_line_perl { 
    my ($find, $search) = @_; 

    my $max_distance = 2; 

    for my $offset (0 .. length($search) - length($find)) { 
     my $substr = substr($search, $offset, length($find)); 
     my $hd = hd($find, $substr); 
     if ($hd <= $max_distance) { 
      return ($hd, $substr); 
     } 
    } 
    return (undef, undef); 
} 

sub hd { 
    return ($_[0]^$_[1]) =~ tr/\001-\377//; 
} 

wo check_line_c.c ist:

void check_line_c(char* find, char * search) { 
    int max_distance = 2; 
    int flen = strlen(find); 
    int last_ind = strlen(search) - flen; 

    SV *dis = &PL_sv_undef; 
    SV *match = &PL_sv_undef; 
    for (int ind = 0; ind <= last_ind; ind++) 
    { 
     int count = 0; 
     for (int j = 0; j < flen; j++) 
     { 
      if (find[j]^search[ind+j]) count++; 
     } 
     if (count < max_distance) 
     { 
      match = newSV(flen); 
      sv_catpvn(match, search+ind, flen); 
      dis = newSViv(count); 
      break; 
     } 
    } 
    Inline_Stack_Vars; 
    Inline_Stack_Reset; 
    Inline_Stack_Push(sv_2mortal(dis)); 
    Inline_Stack_Push(sv_2mortal(match)); 
    Inline_Stack_Done; 
} 

Der Ausgang ist (Ubuntu-Laptop mit Intel Core i7-4702MQ CPU @ 2.20GHz):

Benchmark: timing 4000000 iterations of c, perl... 
     c: 2 wallclock secs (0.76 usr + 0.00 sys = 0.76 CPU) @ 5263157.89/s (n=4000000) 
     perl: 19 wallclock secs (18.30 usr + 0.00 sys = 18.30 CPU) @ 218579.23/s (n=4000000) 

Dies ergibt also eine 24-fache Beschleunigung für diesen Fall.

Answer 2

Ich würde vorschlagen, einen wirklich schlechten Hash-Algorithmus zu erstellen. Etwas nettes und reversibel und ineffizient, wie die Summe der Charaktere. Oder vielleicht die Summe eindeutiger Werte (1-4), die durch die Zeichen repräsentiert werden.

Berechnen Sie die Zielsummen und berechnen Sie auch die maximal zulässige Varianz. Das heißt, wenn das Ziel eine Übereinstimmung mit zwei Substitutionen ist, was ist der maximal mögliche Unterschied? (4-1 + 4-1 = 6).

Dann berechnen Sie für jedes "Fenster" von Text der entsprechenden Länge in der Zieldatendatei eine laufende Punktzahl. (Fügen Sie am Ende ein Zeichen hinzu, lassen Sie ein Zeichen von Anfang an los und aktualisieren Sie den Hash-Wert.) Wenn der Wert für ein Fenster innerhalb des zulässigen Bereichs liegt, können Sie weitere Untersuchungen durchführen.

Möglicherweise möchten Sie dies als verschiedene Durchgänge implementieren. Möglicherweise sogar als verschiedene Stufen in einer Shell-Pipeline oder einem Skript. Die Idee ist, dass Sie möglicherweise Teile der Suche parallelisieren können. (Zum Beispiel könnten alle Übereinstimmungszeichenfolgen mit der gleichen Länge von einem Prozess durchsucht werden, da die Hash-Fenster identisch sind.)

Außerdem ist es natürlich vorteilhaft, dass Sie Ihre frühe Arbeit behalten können, wenn Ihr Programm stürzt in den späteren Phasen ab. Und Sie können sogar die ersten Teile des Prozesses laufen lassen, während Sie noch die Endstufen entwickeln.