1

Perl（BioPerlではない）に2つの連続する各文字の数を見つける方法はありますか？

つまり、次のAA, AC, AG, AT, CC, CA, ...ようなシーケンスの数：

$sequence = 'AACGTACTGACGTACTGGTTGGTACGA'

PS：正規表現、つまりシーケンス内のGCの数を返す$ GC =（$ sequence =〜s / GC / GC / g）を使用して、手動で作成できます。

自動化された一般的な方法が必要です。

4

3 に答える 3

3

あなたは私をしばらく混乱させました、しかし私はあなたが与えられた文字列のジヌクレオチドを数えたいと思います。

コード：

my @dinucs = qw(AA AC AG CC CA CG);
my %count;
my $sequence = 'AACGTACTGACGTACTGGTTGGTACGA';

for my $dinuc (@dinucs) {
    $count{$dinuc} = ($sequence =~ s/\Q$dinuc\E/$dinuc/g);
}

データからの出力::ダンパー：

$VAR1 = {
          "AC" => 5,
          "CC" => "",
          "AG" => "",
          "AA" => 1,
          "CG" => 3,
          "CA" => ""
        };
于 2011-11-18T15:19:59.293 に答える
3

TLPの答えに近いが、置換なし：

my $sequence = 'AACGTACTGACGTACTGGTTGGTACGA';
my @dinucs = qw(AA AC AG AT CC CG);
my %count = map{$_ => 0}@dinucs;

for my $dinuc (@dinucs) {
    while($sequence=~/$dinuc/g) {
        $count{$dinuc}++;
    }
}

基準：

my $sequence = 'AACGTACTGACGTACTGGTTGGTACGA';
my @dinucs = qw(AA AC AG AT CC CG);
my %count = map{$_ => 0}@dinucs;

my $count = -3;
my $r = cmpthese($count, {
        'match' => sub {
            for my $dinuc (@dinucs) {
               while($sequence=~/$dinuc/g) {
                    $count{$dinuc}++;
               }
            }
        },
        'substitute' => sub {
            for my $dinuc (@dinucs) {
                $count{$dinuc} = ($sequence =~ s/\Q$dinuc\E/$dinuc/g);
            }
         }
});

出力：

              Rate substitute      Match
Substitute 13897/s         --       -11%
Match      15622/s        12%         --
于 2011-11-18T15:40:38.093 に答える
1

注意すれば正規表現は機能しますが、より高速で柔軟な substr を使用した簡単な解決策があります。

(この投稿の時点で、承認済みとしてマークされた正規表現ソリューションは、「AAAA...」のような繰り返される領域のジヌクレオチドを正しくカウントできません。これらの領域には、自然に発生するシーケンスが多数あります。

「AA」に一致すると、正規表現検索は 3 番目の文字で再開され、中央の「AA」ジヌクレオチドはスキップされます。ある位置に 'AC' がある場合、当然、次の塩基にそれがないことが保証されるため、これは他のジヌクレオチドには影響しません。質問で与えられた特定の配列は、塩基が 3 回続けて出現しないため、この問題に悩まされることはありません。)

私が提案する方法は、任意の長さの単語を数えることができるという点でより柔軟です。正規表現をより長い単語に拡張することは、正確な数を得るために正規表現でさらに体操を行う必要があるため、複雑です。

sub substrWise {
    my ($seq, $wordLength) = @_;

    my $cnt = {};

    my $w;
    for my $i (0 .. length($seq) - $wordLength) {
        $w = substr($seq, $i, $wordLength);
        $cnt->{$w}++;
    }

    return $cnt;
}

sub regexWise {
    my ($seq, $dinucs) = @_;

    my $cnt = {};
    for my $d (@$dinucs) {
        if (substr($d, 0,1) eq substr($d, 1,1) ) {
            my $n = substr($d, 0,1);
            $cnt->{$d} = ($seq =~ s/$n(?=$n)/$n/g); # use look-ahead
        } else {
            $cnt->{$d} = ($seq =~ s/$d/$d/g);
        }
    }

    return $cnt;
}


my @dinucs = qw(AA AC AG AT CA CC CG CT GA GC GG GT TA TC TG TT);

my $sequence = 'AACGTACTGACGTACTGGTTGGTACGA';

use Test::More tests => 1;
my $rWise = regexWise($sequence, \@dinucs);
my $sWise = substrWise($sequence, 2);
$sWise->{$_} //= '' for @dinucs; # substrWise will not create keys for words not found
# this seems like desirable behavior IMO,
# but i'm adding '' to show that the counts match
is_deeply($rWise, $sWise, 'verify equivalence');

use Benchmark qw(:all);
cmpthese(100000, {
    'regex' => sub {
        regexWise($sequence, \@dinucs);
    },
    'substr' => sub {
        substrWise($sequence, 2);
    }

出力：

1..1
ok 1 - verify equivalence
          Rate  regex substr
regex  11834/s     --   -85%
substr 76923/s   550%     --

より長い配列 (10 ～ 100 kbase) の場合、利点はそれほど顕著ではありませんが、それでも約 70% 勝っています。

于 2011-12-15T22:30:08.887 に答える