[求助]在5万个长3k的字符串中快速定位定长字符串

Lonki

UID: 21498
帖子: 1
积分: 2
在线时间: 10 分钟

11^# Lonki 发表于 2007-09-24 12:26

QUOTE:

原帖由 galaxy001 于 2007-9-23 23:16 发表
for (@lines) {

$name = $_;

}
光这个循环就309s
开来瓶颈在文件IO，……

上面的瓶颈主要在于频繁的读取磁盘, 试一试这个, 唯一改变的code是一开始就一次性读入整个文件到内存...

open(FIN, $ARGV[0]) or die "Where the hell is $ARGV[0]?";
    my @lines = (<FIN>, ">");

    my ($pattern, $name, $str) = ('CAAT', '', '');

    for (@lines) {
        if (substr($_, 0, 1) eq '>') {
            while ($str =~ /$pattern/g) {
                print "Found pattern in pos ",  $-[0] + 1, " of $name";
            }
            $name = $_;
            $str  = '';
        } else {
            $str .= $_;
        }
    }

修正: 读完文件后补一个'>'

Lonki

UID: 21498
帖子: 1
积分: 2
在线时间: 10 分钟

12^# Lonki 发表于 2007-09-24 23:12

简单做了个对比测试, 一次性读入大文件后处理 VS Tie::File后处理
100M的文本文件, 机器Pentium M 1.5G + 768M

测试用例 1:
open FIN, $ARGV[0];
my @lines = <FIN> ### 一次性读入该文件到内存
for (@lines) {
### 每次采用正则式对若干行处理...
}

测试用例 2:
tie my @lines, 'Tie::File', $ARGV[0];
for (@lines) {
### 每次采用正则式对若干行处理...
}

对于耗时, 内存占用(内存+虚拟内存) 均在连续运行后取均值.

测试结果:
测试用例 1:
耗时: 25秒
内存峰值: 600M, 运行过程中增长较快.

测试用例 2:
耗时: 5分钟
内存峰值: 120M, 运行过程中增长较慢, 硬盘灯闪的较频繁.

本来还想多测几种(sysread...), 时间不够了:<

个人觉得应根据实际文件大小和硬件配置, 来寻找空间和时间的平衡点.
从脚本语言的角度来说, 效率低一点无所谓, 重要的是快速简洁的完成代码.

仅供参考.

galaxy001

UID: 40457
帖子: 169
积分: 388
在线时间: 1 天 15 小时

13^# galaxy001 发表于 2007-09-25 18:32

今早有些困,没看到第二页,自己采取了先读入内存的办法:

[Copy to clipboard] [ - ]

CODE:

tie my @lines, 'Tie::File', $ARGV[0] or die "Where the hell is $ARGV[0]?";
my ($pattern, $str) = ('TTCAATCTCCCTTCGTTCCCTGTCC', '');
my @FastaCache;
# Begin to read fasta seq.
my $name=$lines[0]; # for the 1st line
for (@lines) {
      if (substr($_, 0, 1) eq '>') {
         my $tmppos=index $name,' ';
         my $tmpID=substr $name,1,$tmppos-1;
         my $tmpDesc=substr $name,$tmppos+1;
         push @FastaCache,[$tmpID,$tmpDesc,$str] if ! $str eq '';
         $name = $_;
         $str  = '';
      } else {
         $str .= $_;
      }
}
my $tmppos=index $name,' ';
my $tmpID=substr $name,1,$tmppos-1;
my $tmpDesc=substr $name,$tmppos+1;
push @FastaCache,[$tmpID,$tmpDesc,$str] if ! $str eq ''; # for the last seq.
# Finishing reading fasta file.
untie @lines;
for my $row (@FastaCache) {
while (@$row[2] =~ /$pattern/g) {
   # print "Found pattern in pos ",  $-[0] + 1, " of $name\n";
   print @$row[0],"\t",$-[0] + 1,"\n";
}
}

把字符串连接一起做了,应该更快.
对源文件的处理我只能想到目前这个方法,Lonki兄有更美观的代码吗?
目前一次搜索在0.5s内完成,由于我有631191次检索,正在考虑用index试试.

Lonki

UID: 21498
帖子: 1
积分: 2
在线时间: 10 分钟

14^# Lonki 发表于 2007-09-25 19:18

用Tie肯定是比较慢的办法...

index的话, 我测下来基本上没什么改观.

open FIN, $ARGV[0] or die "Where the hell is $ARGV[0]?";
    my @lines = (<FIN>, ">");

    my ($pattern, $name, $str) = ('ACCAGGCCGCCAGCCACAA', '', '');
    my $ps = 0;

    for (@lines) {
        if (substr($_, 0, 1) eq '>') {
            $ps = 0;
            while ((($ps = index($str, $pattern, $ps)) > -1 )) {
                print "Found pattern in pos ",  $ps + 1, " of $name\n";
                $ps += length($pattern);
            }
            $name = $_;
            $str  = '';
        } else {
            $str .= $_;
        }
    }
    close FIN;

如果只考虑速度的话, 用sysread吧
还不够的话, 就用c++写吧, 反正代码也简单的.

galaxy001

UID: 40457
帖子: 169
积分: 388
在线时间: 1 天 15 小时

15^# galaxy001 发表于 2007-09-25 21:36

为啥我用index的还要花6.7倍的时间？

[Copy to clipboard] [ - ]

CODE:

for my $row (@FastaCache) {
my $tmppos=-1;
while ($tmppos) {
$tmppos=index(@$row[2],$pattern2,$tmppos)+1;
print @$row[0],"\t",$tmppos,"\n" if $tmppos;
}
}

[Copy to clipboard] [ - ]

CODE:

my $ps;
for my $row (@FastaCache) {
while ((($ps = index(@$row[2], $pattern2, $ps)) > -1 )) {
print @$row[0],"\t",$ps+1,"\n";
$ps += length($pattern);
}
}

2个都是6.7倍时间。用my $tmp=@$row[2]给index作参数更慢。

文件只读一次，用tie不会慢很多。我只是懒得换了。反正重点在搜索循环。

不会说用数组的数组会慢吧？

Lonki

UID: 21498
帖子: 1
积分: 2
在线时间: 10 分钟

16^# Lonki 发表于 2007-09-25 22:20

QUOTE:

原帖由 galaxy001 于 2007-9-25 21:36 发表
为啥我用index的还要花6.7倍的时间？

for my $row (@FastaCache) {
my $tmppos=-1;
while ($tmppos) {
$tmppos=index(@$row[2],$pattern2,$tmppos)+1;
print @$row[0],"\t",$tmppos," ...

tie不是读一次哦~
它是一次将文件的"索引"映射到内存中的数组, 该数组的一个元素对应文件的一行.
在每一次循环访问该数组时, 就会发生对文件的读取, 所以慢.
如果你在自己机器上测试, 会看到磁盘灯狂闪的

你要搜60W多个pattern的话, 得要点时间....

刚才用C写了一个, 我的机器上4S之内跑完100M(之前发在11楼的代码跑25S), Dev-C++(gcc,g++) for Win.

#include <stdio.h>
#include <string.h>

int main(int argc, char *argv[]) {
    FILE  *fp;
    char  line[256] = ""; /* keeps one line of input in a loop */
    char  name[256] = ""; /* keeps the name of each DNA group */
    char  str[4096] = ""; /* keeps each DNA group */
    char  *p = NULL;          /* position of the pattern in string */


    if (argc != 3) {
        printf("Usage: %s pattern filename\n", argv[0]);
        return 1;
    }

    if (NULL == (fp = fopen(argv[2], "r"))) {
        printf("Cannot open %s\n", argv[2]);
        return 99;
    }

    /* parse previous DNA group when read '>' */
    while (NULL != (fgets(line, sizeof(line), fp))) {
        if ('>' == line[0]) {
            p = str;
            while (p = strstr(p, argv[1])) {
                printf("Found pattern in pos %d of %s", p - str + 1, name);
                p += strlen(argv[1]);
            }
            strcpy(name, line);
            memset(str, 0, sizeof(str));
        }
        else {
            strcat(str, line);
        }
    }

    fclose(fp);

    /* parse the last DNA group */
    p = str;
    while (p = strstr(p, argv[1])) {
        printf("Found pattern in pos %d of %s", p - str + 1, name);
        p += strlen(argv[1]);
    }


    return 0;
}

Updated: 09/26/2007 3:11 PM

galaxy001

UID: 40457
帖子: 169
积分: 388
在线时间: 1 天 15 小时

17^# galaxy001 发表于 2007-09-25 23:36

看了，还是那个问题，第一条序列丢失（name未定义），最后一条也丢失（无>则不处理上一条）。
关键是我目前能看懂C却不会语法，下笔老错。
由于pattern[64] 是来源于另一tsv文件，如果用perl的``调用就要写服务daemon来保存序列信息了，看来还是没法避开C。
我决定试试C，同时先把perl版跑着，但估计是perl先干完。

我的环境debian 4.

Lonki

UID: 21498
帖子: 1
积分: 2
在线时间: 10 分钟

18^# Lonki 发表于 2007-09-26 02:20

QUOTE:

原帖由 galaxy001 于 2007-9-25 23:36 发表
看了，还是那个问题，第一条序列丢失（name未定义），最后一条也丢失（无>则不处理上一条）。
关键是我目前能看懂C却不会语法，下笔老错。
由于pattern[64] 是来源于另一tsv文件，如果用perl的``调用就要写服 ...

这个我测试过程中到没有遇到...
你把出现问题的input贴几行呢? 以及相应的pattern.

galaxy001

UID: 40457
帖子: 169
积分: 388
在线时间: 1 天 15 小时

19^# galaxy001 发表于 2007-09-26 09:20

这个,人工解释执行一下不就看到了.

[Copy to clipboard] [ - ]

CODE:

>s111
a111---111NNN
a222---222NNN
>s222
a111---111NNN
a222---222NNN
>s333
a111---111NNN
a222---222NNN

假设是找“---”。
第一行，>会直接去打印，此时name尚未赋值。
第三条，扫完后，没有更多>，不打印。
-----------------------------------
想错了了，第一行没问题，刚循环时str=NULL，绝不会有匹配的。

Lonki

UID: 21498
帖子: 1
积分: 2
在线时间: 10 分钟

20^# Lonki 发表于 2007-09-26 10:21

QUOTE:

原帖由 galaxy001 于 2007-9-26 09:20 发表
这个,人工解释执行一下不就看到了.

>s111
a111---111NNN
a222---222NNN
>s222
a111---111NNN
a222---222NNN
>s333
a111---111NNN
a222---222NNN

假设是找“---”。
第一行，>会直接去打印，此时 ...

刚死机了.... 最后应该再push一个>吧

Updated: 09/26/2007 11:23 AM
11楼 Perl 代码
16楼 C 代码