☰

Zählen von DNA-Codons in einer DNA-Datei

tag-icon • tag-icon bash shell-script shell awk
Zählen von DNA-Codons in einer DNA-Datei

Ich möchte ein Bash-Skript erstellen, das eine DNA-Datei aufnimmt und überprüft, ob sie keine Zeilenumbruchzeichen oder Leerzeichen enthält, und dann die eindeutigen Codons zusammen mit der Anzahl ihrer Vorkommen ausgibt. Ich habe den folgenden Code verwendet, aber das Codon gibt mir immer die Ausgabe „bash-3.2$“ aus. Ich bin so verwirrt, ob meine Syntax falsch ist und warum ich nicht die richtige Ausgabe erhalte.

! /bin/bash

for (( pos=1; pos < length - 1; ++pos )); do
    codon = substr($1, $pos, 3)
    tr-d '\n' $1 | awk -f '{print $codon}' | sort | uniq -c
done

Wenn beispielsweise eine Datei mit dem Namen dnafile das Muster aacacgaactttaacacg enthält, dann übernimmt das Skript die folgende Eingabe und Ausgabe

 $script dnafile              
 aac 3
 acg 2
 ttt 1

Antwort1

Sie erhalten diese Ausgabe, weil die erste Zeile Ihres Skripts eine neue bashShell startet.

Diese Zeile sollte lauten

#!/bin/bash

(beachten Sie das #am Anfang).

Sie vermischen dann awkSyntax mit Shellcode auf eine Weise, die niemals funktionieren wird.

Machen Sie es sich stattdessen einfach, teilen Sie Ihre Datei in Gruppen von je drei Zeichen auf, sortieren Sie diese und zählen Sie, wie viele eindeutige Zeichen Sie erhalten:

$ fold -w 3 dnafile | sort | uniq -c
   3 aac
   2 acg
   1 ttt

Dies würde funktionieren, solange die Eingabe immer ein Vielfaches von drei Zeichen enthält, ohne eingebettete Leerzeichen oder andere Zeichen.

Antwort2

(echo aacacgaactttaacacg ;echo aacacgaactttaacacg ) |
  perl -ne '# Split input into triplets (A3)
            # use each triplet as key in the hash table count
            #   and increase the value for the key
            map { $count{$_}++ } unpack("(A3)*",$_);
            # When we are at the end of the file
            END{ 
                 # Remove the key "" (which is wrong)
                 delete $count{""};
                 # For each key: Print key, count
                 print map { "$_ $count{$_}\n" } keys %count
            }'

Antwort3

Eine etwas ausführlichere awkVersion

awk 'BEGINFILE{print FILENAME; delete codon}
     ENDFILE {
     if (NR!=1 || NF!=1 || length($0)%3!=0){
         print "is broken"}
     else{
         for (i=1; i<=length($0); i+=3) codon[substr($0,i,3)]++}; 
         for (c in codon) print c, codon[c]; 
         print ""}' file*

Für diese Eingabe

Datei1: OK

aacacgaactttaacacg

Datei2: Leerzeichen

aacacgaact ttaacacg

Datei3: Zeilenumbruch

aacacgaact
ttaacacg

Datei4: kein Vielfaches von 3 Basen

aacacgaactttaacac

Du erhältst

file1
aac 3
ttt 1
acg 2

file2
is broken

file3
is broken

file4
is broken

Wenn Sie nur die Dateien reparieren möchten und keine haben, wie file4dann catIhre Dateien durch trvon einem Ende awkoder dem anderen, genau wie Ihr Beispiel

<<< $(cat file[1..3] | tr -d "\n ")

verwandte Informationen