Counting DNA codons in DNA file

Question

I want to create a bash script that takes in a dna file and checks that it has no newline characters or white space characters, and then outputs the unique codons along with their count of the number of times they occur. I have used the following code but the codon keeps giving me an output of "bash-3.2$". I am so confused as to whether my syntax is wrong and why I'm not getting the proper output.

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

For example if a file named dnafile contains the pattern aacacgaactttaacacg then the script will take the following input and output

 $script dnafile              
 aac 3
 acg 2
 ttt 1

Answer 1

You get that output because the first line of your script starts a new bash shell.

That line should read

#!/bin/bash

(note the # at the start).

You then mix awk syntax with shell code in a way that will never work.

Instead, keep it simple and chop up your file in groups of three characters, sort these and count how many unique ones you get:

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

This would work as long as the input always contains a multiple of three characters, with no embedded spaces or other characters.

Answer 2

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

Answer 3

A slightly more long-winded awk version

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*

For this input

file1 : OK

aacacgaactttaacacg

file2 : space

aacacgaact ttaacacg

file3 : linebreak

aacacgaact
ttaacacg

file4 : not a multiple of 3 bases

aacacgaactttaacac

You get

file1
aac 3
ttt 1
acg 2

file2
is broken

file3
is broken

file4
is broken

If you just want to repair the files and have none like file4 then cat your files through tr from one end of awk or the other, just like your example

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