2

私の質問は、文字列の大きなリスト (7mio と 3mio) を含む 2 つのファイルの比較に関するものです。私がする必要があるのは、1 つのリスト (「データベース」と呼ばれる) の各要素を他のリスト (「参照」と呼ばれる) のすべての要素と比較することです。

このための非常に基本的なスクリプトを作成しましたが、膨大な量の比較 (~3 兆) のため、検索を並列化したいと考えています。両方のファイルが非常に大きいため、私の考えでは、それらを 2 つの独立した配列 (@reference と @dbPep) に 1 回だけロードしてから、サブルーチンを使用して、別のスレッドでリストのチャンクを比較します。つまり、スレッド 1 で 0 から 1000000、1000001 を比較します。 - スレッド 2 では 2000000 など。チャンク サイズはスレッドの数によって決まりますが、最後のスレッドはモジュロを使用した後に残ったものも取得します。出力は、ファイルにパイプできるように STDOUT に出力する必要があります。結果の順序は重要ではありません。

今、私の問題は、$ARGV 入力に応じてスレッドを正しく動的に開始することにあります。

ここに私がこれまでに持っているもの、マルチスレッドがありますが、まだ実際には機能していないので、コメントアウトしました:

#!/usr/bin/perl
use Thread;

    $minLength      = $ARGV[0];
    $numThreads         = $ARGV[1];
    my $reference_File  = $ARGV[2];
    my $peptideDB_File  = $ARGV[3];

    if ($minLength == ""){$minLength = 6;}

    print STDERR ("Threads:\t".$numThreads."\n");
    print STDERR ("Reference:\t".$reference_File."\n");
    print STDERR ("peptide DB:\t".$peptideDB_File."\n----------------\n");
    print STDERR ("Loading reference...\n");

# Load reference peptides
    open (REF,$reference_File) || die ("File not found!");
    @reference = ();
    while (<REF>) {
        chomp;
        my @refPeptide = split("\t");
        if($refPeptide[1] eq "Sequence") {next;}
        push(@reference,$refPeptide[1]);
    }
    close (REF);
    print STDERR ("Done!\n");
    print STDERR ("Length:".@reference."\n");
    print STDERR ("Loading peptide database...\n");


# Open peptide DB to check
    open (PEP,$peptideDB_File) || die ("File not found!");
    @dbID  = ();
    @dbPep = ();
    while (<PEP>) {
        chomp;
        my @dbPeptide = split("\t");
        if($dbPeptide[1] eq "Sequence") {next;}
        push(@dbID,$dbPeptide[0]);
        push(@dbPep,$dbPeptide[1]);
    }
    close (PEP);
print STDERR ("Done!\n");
    print STDERR ("Length:".@dbPep."\n");
    print STDERR ("Starting analysis...\n");
    print STDERR ("(This can take a while)\n");
    print STDERR ("----------------\n");

#my %ref = map { $_ => 1 } @reference;

    my $dbLength        = @dbPep;
    my $rest            = $dbLength % $numThreads;
    my $chunkStart      = 0;
    my $chunkLength     = ($dbLength - $rest) / $numThreads;
    my $chunkEnd        = $chunkLength;

    #print STDERR($chunkEnd."\n");
    &peptideCheck($chunkStart,$chunkEnd);

    #for (my $i=0;$i<$numThreads;$i++) {
#       if ($i == $numThreads-1) {$chunkEnd += $rest;};
#       new Thread \&peptideCheck, $chunkStart, $chunkEnd;
#       $chunkStart += $chunkLength;
#       $chunkEnd   += $chunkLength;
#   }
    print STDERR ("----------------\nFinished!\n");
    sub peptideCheck {
        my @params = @_;
        for (my $j=@params[0]; $j<@params[1];$j++){
            if (length($dbPep[$j]) < $minLength) {next;}
            if ( grep { $_ eq $dbPep[$j]} @reference ) {
                #print STDERR ("RefCall:\n".$peptideCandidate[1]."\n");
                next;
            }
            else {
                #print STDERR ("Found novel peptide:\n");
                print (">".$dbID[$j]."\n".$dbPep[$j]."\n");
            }
        }
    }

以下は、私のファイルの (非常に) 短い例です。最初はデータベース、2 番目はリファレンスです。

Protein_Name    Sequence    Unique_ID   Monoisotopic_Mass   Predicted_NET   Tryptic_Name
IPI:80000000.1  MGSTSEK 1   738,321789  0,0756  t1.1
IPI:80000001.1  MFLCSSFIFHHDCEASPAMLNCG 2   2559,0512972    0,5426  t1.1
IPI:80000002.1  MIVR    3   517,3046228 0,0924  t1.1
IPI:80000002.1  MIVRPPQPC   4   1039,5306652    0,2739  t1.2
IPI:80000002.1  PPQPC   5   540,2366066 0,1163  t2.1
IPI:80000003.1  MPVLMEK 6   846,4343096 0,2742  t1.1
IPI:80000003.1  MPVLMEKGGGR 7   1173,5998042    0,2599  t1.2
IPI:80000003.1  MPVLMEKGGGREPECSSPLGEQTR    8   2587,2192178    0,2881  t1.3
IPI:80000003.1  MPVLMEKGGGREPECSSPLGEQTRVASVGTCWIEVR    9   3887,878977 0,3859  t1.4


Protein_Name    Sequence    Unique_ID   Monoisotopic_Mass   Predicted_NET   Tryptic_Name
sp|P31946|1433B_HUMAN   MTMDK   1   624,2611054 0,1124  t1.1
sp|P31946|1433B_HUMAN   MTMDKSELVQK 2   1308,6417266    0,2304  t1.2
sp|P31946|1433B_HUMAN   MTMDKSELVQKAK   3   1507,7737968    0,2278  t1.3
sp|P31946|1433B_HUMAN   MTMDKSELVQKAKLAEQAER    4   2305,1769438    0,3107  t1.4
sp|P31946|1433B_HUMAN   SELVQK  5   702,3911854 0,1286  t2.1
sp|P31946|1433B_HUMAN   SELVQKAK    6   901,5232556 0,1294  t2.2
sp|P31946|1433B_HUMAN   SELVQKAKLAEQAER 7   1698,9264026    0,2544  t2.3
sp|P31946|1433B_HUMAN   SELVQKAKLAEQAERYDDMAAAMK    8   2695,330871 0,3435  t2.4
sp|P31946|1433B_HUMAN   AKLAEQAER   9   1014,5457814    0,1198  t3.2
IPI:80000000.1  MGSTSEK 1   738,321789  0,0756  t1.1
IPI:80000001.1  MFLCSSFIFHHDCEASPAMLNCG 2   2559,0512972    0,5426  t1.1

どんな助けでも大歓迎です:)

感謝と乾杯

4

1 に答える 1

0

わかりました、ピルクロウが推奨する「スレッド」モジュールに切り替えた後、コードに取り組んでいました(ありがとうございます:))、スレッドの生成は正常に機能しているようです。ただし、正常に動作していないように見えるのは、生成されたスレッドが終了するのを待っていることです (以下を参照)。現在、threads->detach() を使用しています。これは、スレッドからの戻り値を必要とせず、STDOUT への出力が唯一の目的であるためです。プログラムを実行すると、次の出力が得られます。

Threads:    3
Reference:  testRef.txt
peptide DB: testDB.txt
----------------
Loading reference...
Done!
Length:1999
Loading peptide database...
Done!
Length:1999
Starting analysis...
(This can take a while)
----------------
1999    1   0   666 666
main: 0 666
createThr: 0    666
main: 666   1332
createThr: 666  1332
peptideCheck: 0 666
MTMDKSELVQK
main: 1332  1999
createThr: 1332 1999
peptideCheck: 666   1332
MAVMAPRTLLLLLSGALALTQTWAGSHSMR
----------------
Finished!
peptideCheck: 1332  1999
THVIWGTSKDFGISGFR

そのことから、スレッドは正しいパラメーターで適切に生成されていると結論付けていますが、計算に関しては、メイン スクリプトが終了するため中断されます。detach() を使用する場合、これは発生しませんが、何かを見逃したのでしょうか、間違った方法で使用していますか?

新しいコードは次のとおりです。

#!/usr/bin/perl

    use threads;
    use threads::shared;
    use strict;
    use warnings;

    my $minLength :shared   = $ARGV[0];
    my $numThreads          = $ARGV[1];
    my $reference_File      = $ARGV[2];
    my $peptideDB_File      = $ARGV[3];

    if ($minLength eq ""){$minLength = 6;}

    print STDERR ("Threads:\t".$numThreads."\n");
    print STDERR ("Reference:\t".$reference_File."\n");
    print STDERR ("peptide DB:\t".$peptideDB_File."\n----------------\n");
    print STDERR ("Loading reference...\n");

# Load reference peptides
    open (REF,$reference_File) || die ("File not found!");
    my @reference :shared = ();
    while (<REF>) {
        chomp;
        my @refPeptide = split("\t");
        if($refPeptide[1] eq "Sequence") {next;}
        push(@reference,$refPeptide[1]);
    }
    close (REF);
    print STDERR ("Done!\n");
    print STDERR ("Length:".@reference."\n");
    print STDERR ("Loading peptide database...\n");

# Open peptide DB to check
    open (PEP,$peptideDB_File) || die ("File not found!");
    my @dbID :shared  = ();
    my @dbPep :shared = ();
    while (<PEP>) {
        chomp;
        my @dbPeptide = split("\t");
        if($dbPeptide[1] eq "Sequence") {next;}
        push(@dbID,$dbPeptide[0]);
        push(@dbPep,$dbPeptide[1]);
    }
    close (PEP);
    print STDERR ("Done!\n");
    print STDERR ("Length:".@dbPep."\n");
    print STDERR ("Starting analysis...\n");
    print STDERR ("(This can take a while)\n");
    print STDERR ("----------------\n");

#my %ref = map { $_ => 1 } @reference;

    my $dbLength        = @dbPep;
    my $rest            = $dbLength % $numThreads;
    my $chunkStart      = 0;
    my $chunkLength     = ($dbLength - $rest) / $numThreads;
    my $chunkEnd        = $chunkLength;

    print STDERR ($dbLength."\t".$rest."\t".$chunkStart."\t".$chunkLength."\t".$chunkEnd."\n");

    #&peptideCheck($chunkStart,$chunkEnd);

    for (my $i=0;$i<$numThreads;$i++) {
        if ($i == $numThreads-1) {$chunkEnd += $rest;};
        print STDERR ("main: ".$chunkStart."\t".$chunkEnd."\n");
        createThr($chunkStart, $chunkEnd);
        $chunkStart += $chunkLength;
        $chunkEnd   += $chunkLength;
    }
    print STDERR ("----------------\nFinished!\n");

    sub peptideCheck {
        my @params = @_;
        print STDERR ("peptideCheck: ".$params[0]."\t".$params[1]."\n");
        for (my $j=$params[0]; $j<$params[1];$j++){
            if (length($dbPep[$j]) < $minLength) {next;}
            print ($dbPep[$j]."\n");
            if ( grep { $_ eq $dbPep[$j]} @reference ) {
                next;
            }
            else {
                print STDOUT (">".$dbID[$j]."\n".$dbPep[$j]."\n");
            }
        }
    }

    sub createThr {
        my @thrParams = @_;
        print STDERR ("createThr: ".$thrParams[0]."\t".$thrParams[1]."\n");
        my $thr = threads->create(\&peptideCheck, $thrParams[0], $thrParams[1]);
        $thr->detach();
    }
于 2012-08-09T07:41:59.490 に答える