☰

DNA配列の計算

tag-icon • tag-icon text-processing bioinformatics
DNA配列の計算

スペースで区切られたDNA配列があります。スペースを削除し、スペース文字なしでシーケンス数を返す必要があります。これを行う方法について助けていますか？ bashで端末を使用しています。

たとえば、シーケンスは次のようになります。

GTCGATTGCAAGGATCCGCATGGGATAAAGGAATCGCAGTTCGAACAGGCAATGCCGCAG
CTATGATAGGACATCTCTTGGAGACACCTATTAATGTTTCAGAAACGGATACCTTGGTTG
TCCAGTACGAAATTAAGTTGGACAATTCTTTGACGTGCGGC
CTATATTAAAATTGTGGGTACATCACTCTCTTACCTGAGAATTCCAACAGAGCAGGACGC
TAACCCAGTGTCTATACCAGTCTGTGGCTTTGAAAGATTAGACACATTTCTGGATGAATT
TTCAAATTCTAAATTGATCGTTCAGTCTACACTAAGACATTCGTACGTTAGTCTTGAGAA

スペースを削除し、正確にいくつの塩基があるかどうかを計算したいと思います。あるいは、A、C、G、またはTが何個あるかを計算し、スペースを計算せずに追加することもできます。

答え1

GNU awkを使用したマルチ文字RSとRT：

$ awk -v RS='[^\n]' 'RT{cnt[RT]++} END{for (base in cnt) print base, cnt[base]}' file
A 101
C 68
T 98
G 74

あなたの説明では、「基本」はあなたの例では改行文字ではなくすべての文字であると仮定します。

答え2

空白行や末尾の空白などがないと仮定すると、fold個々の文字ストリームを作成してsort結合し、uniq -c次を使用して各文字の数を計算できます。

$ fold -w 1 file | sort | uniq -c
 101 A
  68 C
  74 G
  98 T

入力にジャンク空白文字がある場合は、最初のtr手順を使用してこれらの文字を削除します。

$ tr -d -c 'ACGT' <file | fold -w 1 | sort | uniq -c
 101 A
  68 C
  74 G
  98 T

ここのコマンドは、またはを除くtr入力のすべての文字を削除します。ACGT

sort | uniq -cパイプの終わりにあるビットは、入力内のawk各文字の発生回数を数え、次を報告する単一のコマンドで置き換えることができます。

$ tr -d -c 'ACGT' <file | fold -w 1 | awk '{ count[$0]++ } END { for (ch in count) printf "%4d %s\n", count[ch], ch }'
 101 A
  68 C
  74 G
  98 T

しかし、私たちがそれを導入しようとしているならばawk、それを取り除く方が良いかもしれませんfold：

$ tr -d -c 'ACGT' <file | awk '{ for (i = 1; i <= length; ++i) count[substr($0,i,1)]++ } END { for (ch in count) printf "%4d %s\n", count[ch], ch }'
 101 A
  68 C
  74 G
  98 T

...次のような場合もありますtr。

$ awk '{ gsub("[^ACGT]", ""); for (i = 1; i <= length; ++i) count[substr($0,i,1)]++ } END { for (ch in count) printf "%4d %s\n", count[ch], ch }' file
 101 A
  68 C
  74 G
  98 T

awk美しく印刷されたコード：

{
    gsub("[^ACGT]", "")  # removes anything not A, C, G, or T
    for (i = 1; i <= length; ++i)
        count[substr($0, i, 1)]++
}
END {
    for (ch in count) {
        printf "%4d %s\n", count[ch], ch
    }
}

gsub()代わりに、最初のブロック（入力の各行を解析）を再構築して使用できますsubstr()。

{
    count["A"] += gsub("A", "A")
    count["C"] += gsub("C", "C")
    count["G"] += gsub("G", "G")
    count["T"] += gsub("T", "T")
}
END {
    for (ch in count) {
        printf "%4d %s\n", count[ch], ch
    }
}

...しかし、入れ子になったコードが少し減った以外は、以前のコードと比較して大幅に改善されません（一部のユーザーの読みやすさに役立たない限り）。

答え3

Perl 1つのライナーを使用してください。

perl -F'' -e '
    BEGIN{my %h}
    map { /\S/ and $h{$_}++ } @F;
    END{print map { "$_ $h{$_}\n" } keys %h}
' file

出力

C 68
A 101
G 74
T 98

答え4

使用幸せ（以前のPerl_6）

raku -e '.say for slurp.comb(/\S/).Bag.pairs;'

出力例：

G => 74
T => 98
A => 101
C => 68

またはタブ区切りの出力（.sayに変更.put）：

~$ raku -e '.put for slurp.comb(/\S/).Bag.pairs;' file
G   74
A   101
T   98
C   68

出力をソートする必要がある場合は、.sort最後に以下を追加してください。

~$ raku -e '.put for slurp.comb(/\S/).Bag.pairs.sort;' file
A   101
C   68
G   74
T   98

または、最高のヌクレオチド数に基づいてソートします。

~$ raku -e '.put for slurp.comb(/\S/).Bag.pairs.sort: -*.value;' file
A   101
T   98
G   74
C   68

またはそれを数えてください（フォントなしスペース）：

~$ raku -e '.put for slurp.comb(/\S/).elems;' file
341

最後に、非常に大きなファイルで作業している場合は、より良いメモリ管理のために代わりに試してみるlines.joinことができます。slurp

入力例：

GTCGATTGCAAGGATCCGCATGGGATAAAGGAATCGCAGTTCGAACAGGCAATGCCGCAG
CTATGATAGGACATCTCTTGGAGACACCTATTAATGTTTCAGAAACGGATACCTTGGTTG
TCCAGTACGAAATTAAGTTGGACAATTCTTTGACGTGCGGC
CTATATTAAAATTGTGGGTACATCACTCTCTTACCTGAGAATTCCAACAGAGCAGGACGC
TAACCCAGTGTCTATACCAGTCTGTGGCTTTGAAAGATTAGACACATTTCTGGATGAATT
TTCAAATTCTAAATTGATCGTTCAGTCTACACTAAGACATTCGTACGTTAGTCTTGAGAA

https://raku.org

関連情報