2

我的问题围绕着比较两个包含大字符串列表(7mio 和 3mio)的文件。我需要做的是将一个列表(称为“数据库”)中的每个元素与另一个列表(称为“参考”)中的每个元素进行比较。

我为此编写了一个非常基本的脚本,但是由于大量的比较(约 3 万亿),我想并行化搜索。由于这两个文件都很大,我的想法是只将它们加载到两个独立的数组(@reference 和@dbPep)中一次,然后使用子例程比较单独线程中的列表块,即线程 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");
            }
        }
    }

这是我的文件的(非常)简短示例,第一个是数据库,第二个是参考。

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

好的,在切换到 pilcrow 推荐的“线程”模块后,我一直在处理代码(谢谢 :)),现在生成线程似乎工作正常。但是,似乎无法正常工作的是等待生成的线程完成(见下文)。目前,我使用 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 回答