Contando codones de ADN en un archivo de ADN

Contando codones de ADN en un archivo de ADN

Quiero crear un script bash que tome un archivo de ADN y verifique que no tenga caracteres de nueva línea o espacios en blanco, y luego genere los codones únicos junto con el recuento de la cantidad de veces que ocurren. He utilizado el siguiente código pero el codón sigue dándome un resultado de "bash-3.2$". Estoy muy confundido sobre si mi sintaxis es incorrecta y por qué no obtengo el resultado adecuado.

! /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

Por ejemplo, si un archivo llamado dnafile contiene el patrón aacacgaactttaacacg, el script tomará la siguiente entrada y salida

 $script dnafile              
 aac 3
 acg 2
 ttt 1

Respuesta1

Obtienes ese resultado porque la primera línea de tu script inicia un nuevo bashshell.

Esa línea debería decir

#!/bin/bash

(tenga en cuenta el #al principio).

Luego mezclas awkla sintaxis con el código shell de una manera que nunca funcionará.

En su lugar, manténgalo simple y divida su archivo en grupos de tres caracteres, ordénelos y cuente cuántos únicos obtiene:

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

Esto funcionaría siempre que la entrada siempre contenga un múltiplo de tres caracteres, sin espacios incrustados ni otros caracteres.

Respuesta2

(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
            }'

Respuesta3

awkUna versión un poco más larga

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*

Para esta entrada

archivo1: Aceptar

aacacgaactttaacacg

archivo2: espacio

aacacgaact ttaacacg

archivo3: salto de línea

aacacgaact
ttaacacg

archivo4: no es un múltiplo de 3 bases

aacacgaactttaacac

Usted obtiene

file1
aac 3
ttt 1
acg 2

file2
is broken

file3
is broken

file4
is broken

Si solo desea reparar los archivos y no tiene ninguno como file4entonces, catsus archivos pasarán trdesde un extremo awko desde el otro, como en su ejemplo.

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

información relacionada