[PerlChina] Re: 在我遇到的这种情况下散列和数组哪个快？

更新了一下代码，使用OligoNuc.ctl保存相关配置内容。其格式如下：
# 寡核苷酸串长度
len of oligonuc string = 12

#来源文件目录
directory of genome sequences = /home/ribozyme/data/genomes/galga/

#是否遍历该目录
traversal the directory = N

遍历目录时所用的正则信息
key words of the files = _temp_\d+\.fsa

#不遍历目录时，给定的文件名
name of the files = Gallus_gallus.WASHUC2.55.dna_rm.chromosome.28.fa

#结果文件名称
result file name = /home/ribozyme/bin/ChrY/hash_12bp.dat

# Finish the ctrl file.

Perl程序文件中也做了修改，如下：

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

my (@ctl,$len,$dir_a,$i,$file,@order,$m,%hash,@array,$temp,$line,@seq,
$seq,$str,@str,$index,$char,$n,$key,$val,$num);
open FA,q(/home/ribozyme/bin/ChrY/OligoNuc.ctl) or die "Can't find the
control file\n";
$i = 0;
while ( $line = <FA> ) {
next if $line=~m{^#};
$line =~ m{\s\=\s(.+)};
$ctl[$i]=$1;
$i++;
}
close FA;

$len = $ctl[0];
$dir_a = $ctl[1];
opendir DIR, $dir_a or die "Directory can't be opened!\n";
if ( $ctl[2] =~ m{^Y$}i ) {
$i = 0;
while ( $file = readdir(DIR) ) {
if ( $file =~ m{$ctl[3]} ) { $order[$i] = $dir_a.$file; $i++; }
}
} elsif ( $ctl[2] !~ m{^Y$}i ) {
$order[0] = $dir_a.$ctl[4];
}
close DIR;

my %code_of = (
A => 0,
G => 1,
C => 2,
T => 3
);
%hash = ();
$m = 0;
$i = 0;

foreach $file (@order) {
open FA,qq($file); # 打开基因组序列文件
$temp = '';
while ( $line = <FA> ) {
if ( $line =~ m{^>} ) { #如果是以大于号开头，则表示是fasta格式的序列标题行，用N代替之。一个文件中可能有多个
fasta格式的序列
$temp .= 'N';
} elsif ( $line !~ m{^>} ) { # 如果不是以大于号开头，则表示是真正的序列，其中包含大写的ATGCN五种字母。
chomp $line; # N表示测序不清楚，只能确定那里有一个核苷酸的位置。
$temp .= $line;
}
else {}
}
close FA;
print length($temp)." raw data obtained for $file.\n";

$i=0;
# 把所有连续的非N的长度在所要测的寡核苷酸串以上的序列都存入序列数组。
while ($temp=~m{([^N]{$len,})}g) { $seq[$i]=$1; $i++; }
print "$i sequences for the file...\n";

foreach $seq (@seq) {
$i = 0;
do {
$str[0] = uc( substr($seq,$i,$len) ); #从序列中截取一个短串
$str[1] = lc( reverse($str[0]) );
$str[1] =~ s{a}{T}g;
$str[1] =~ s{t}{A}g;
$str[1] =~ s{g}{C}g;
$str[1] =~ s{c}{G}g; #得到该短串的互补串，对于序列来说，其互补的结果也是要考虑的。
foreach $str (@str) {
# 计算一个15bp核苷酸串对应的下标，每个字符占用两个二进制位
my $index = 0;
for $char (split q[], $str) { # 按字符分开
$index = $index << 2; # 左移两位
$index += $code_of{$char}; # 将当前位加入末两位
}
if ( exists( $hash{$index} ) ) {
$hash{$index}++; # 如hash中有键值则计数
} else {
$hash{$index} = 1; # 如hash中没有则添加键值
$m++;
}
}
$i++;
} until ( length($seq) < $i + $len );
}
print qq($i 序列已处理，散列中有\t $m 个记录。\n);
}

#遍历散列并转换数字下标为核苷酸串
open FB,qq(>$ctl[5]);
$n=0;
while ( ($index,$val) = each(%hash) ) {
$str='';
$i=0;
do {
$num = $index & 2;
foreach $key ( keys %code_of ) {
if ( $code_of{$key} == $num ) {
$str .= $key;
}
}
$index = $index >> 2;
$i++;
} until ( $i == $len );
print FB qq($str\t$val\n);
$n++;
}
close FB;
print "There are totally $m kinds of oligonucs.\n";
print qq(Normal quit!\n);

On 9月25日, 下午1时56分, 空格 <ribozyme2...@gmail.com> wrote:
--~--~---------~--~----~------------~-------~--~----~
您收到此信息是由于您订阅了 Google 论坛"PerlChina Mongers 讨论组"论坛。
要在此论坛发帖，请发电子邮件到 perlchina@googlegroups.com
要退订此论坛，请发邮件至 perlchina+unsubscribe@googlegroups.com
更多选项，请通过 http://groups.google.com/group/perlchina?hl=zh-CN 访问该论坛
-~----------~----~----~----~------~----~------~--~---