文字列内の文字の長さnの組み合わせの数をカウントする方法

Question

ATCGの組み合わせの出現を次の1つのライナーで示し、長さ6の文字列を出力します。 0マッチの正規表現やその他の部分を "0 ATTTAG"のような場所に変更する方法はありますか？

#!/bin/bash 
for file in e_coli.fa 
do 
    base=$(basename $file .fa) 
    cat $file | perl -nE 'say for /(?<=([ATCG]{6}))/g' \ 
     | sort | uniq -c >> ${base}_hexhits_6mer.txt 
done 

stdout: 
    465 AAAAAA 
    607 AAAAAC 
    661 AAAAAG 
    581 AAAAAT 
    563 AAAACA 
    807 AAAACC 
    770 AAAACG 
    373 AAAACT 
    663 AAAAGA 
    1213 AAAAGC 

Answer 1

uniq -cはラインが発生した回数をカウントするので要求された変更は、完全に書き直さを必要と0、それはおそらく返すことはできません。

perl -e' 
    while (<>) { 
     ++$counts{$_} for /(?=([ATCG]{6}))/g; 
    } 

    for my $seq (glob("{A,C,G,T}" x 6)) { 
     printf("%7d %s\n", $counts{$seq}, $seq); 
    } 
' "$file" >"${base}_hexhits_6mer.txt" 

Answer 2

あなたがしたいことははるかに複雑です。 が表示されていないことを確認するには、が表示されている必要があります。可能なすべての組み合わせの文字をフィルタリングできます。ここで

、私はハッシュにタリー全て（__DATA__から読んで）「み」ATCGAsの文字列の文字、Ts、CsとGsを見つけるために、Perlでsubstrを使用してスライディングウィンドウアプローチを使用します。これらは並べ替えられ、最も一般的に観察される6量体が最初に示され、印刷される。

use strict; 
use warnings; 

my @bases = qw/ A G C T /; 

my %data; 

for my $a1(@bases){ 
    for my $a2(@bases){ 
     for my $a3(@bases){ 
      for my $a4(@bases){ 
       for my $a5(@bases){ 
        for my $a6(@bases){ 
         $data{"$a1$a2$a3$a4$a5$a6"} = 0; 
        } 
       } 
      } 
     } 
    } 
} 

my $nucs = <DATA>; 

my $len = length($nucs); 

for (my $i = 0; $i <= $len - 6; $i++) { 
    my $kmer = substr($nucs, $i, 6); 
    next if $kmer =~ tr/ACGT//c; 
    $data{$kmer}++; # populate hash with "seen" 6-mers 
} 

# print out sorted hash 
foreach my $seq (sort { $data{$b} <=> $data{$a} } keys %data){ 
    print "$seq,$data{$seq}\n"; 
} 

__DATA__ 
ATGCCCGTCGTAGTCATGCATGCATCGATCGATGCATGCTACGTGTTGT 

は明らかに私が何をやったかよりも、文字列内のすべての文字の順列を計算するより良い/きれいな方法があるように起こっているが、それは動作します。

ボロディン氏によれば、これは主に "見えない"文字列のバリエーションを表示します。

Answer 3

最も簡単な方法は、各パターンについて出現回数のハッシュを構築し、その後、このプログラムは形成された全ての可能な6文字の文字列のリストを生成するglobトリックを使用して、すべての可能なパターン

のカウントを印刷しますA、T、CとG

から

use strict; 
use warnings 'all'; 

my @files = qw/ e_coli.fa /; 

my %counts; 

for my $file (@files) { 

    open my $fh, '<', $file or die qq{Unable to open "$file" for input: $!}; 

    while (<$fh>) { 
     ++$counts{$1} while /(?= ([ATCG]{6})) /gx; 
    } 
} 

for my $pattern (glob '{A,T,C,G}' x 6) { 
    printf "%4d %s\n", $counts{$pattern} // 0, $pattern; 
} 

Answer 4

あなたは多くのデータを持っていて、少し速く何かを必要とする場合は、ここではCソリューションです：

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

void reader(FILE* in, unsigned long hist[4096]) { 
    for (unsigned long key=0, count=0;;) { 
    switch(getc(in)) { 
     case EOF:      return; 
     case 'A': key <<= 2;   break; 
     case 'C': key <<= 2; key += 1; break; 
     case 'G': key <<= 2; key += 2; break; 
     case 'T': key <<= 2; key += 3; break; 
     default: count=0;    continue; 
    } 
    if (count == 5) ++hist[key & 0xFFF]; 
    else ++count; 
    } 
} 

int putkey(FILE* out, unsigned long key) { 
    char s[6]; 
    for (int j=6; j--; key >>= 2) s[j] = "ACGT"[key&3]; 
    return fprintf(out, "%.6s", s); 
} 

void writer(FILE* out, unsigned long hist[4096]) { 
    for (unsigned long key = 0; key < 4096; ++key) { 
    fprintf(stdout, "%7lu ", hist[key]); 
    putkey(out, key); 
    putchar('\n'); 
    } 
} 

int main(int argc, char** argv) { 
    FILE* in = stdin; 
    if (argc > 1) in = fopen(argv[1], "r"); 
    if (!in) { perror(argv[1]); exit(1); } 
    unsigned long hist[4096] = {0}; 
    reader(in, hist); 
    writer(stdout, hist); 
    return 0; 
} 

31MBのfastqサンプルを処理するのに半分以下の時間がかかりました（これには、4096個の可能な6文字シーケンスがすべて含まれています）。 Perlの解はそれぞれ12秒（fugu）と18秒（ikegami/borodin）であった。