Renaming fasta files with their headers (gene name) ?

Question

[sunnykevin97]

Hi

I had around 10,0000 gene sequences in individual fasta files. I'd like to rename each file with their header name containing the gene name.

Original file

head 1.fasta

==> 1.fasta <==

> Gloriosasuperba; 8324-9004; -; atp6 

ATGACAGTAAGCCTTTTTGACCAATTTATGAGCCCCACACTACTAGGCATCCCCCTGCTC

Modified file

head atp6.fasta

==> atp6.fasta <==

> atp6 
ATGACAGTAAGCCTTTTTGACCAATTTATGAGCCCCACACTACTAGGCATCCCCCTGCTC

Answer 1

[Matthias Zepper]

I would not rename the files, but write the contents into a new file, which is trivial with awk:

cat *.fasta | tr "\n" ";" | tr ">" "\n>" | awk -F ";" '{a="../outdir/$4.fasta" ; print ">"$4"\n"$5 > a}'

This quick n' dirty solution assumes that there are no extra linebreaks in the sequences of your FASTA files and that there is a "outdir" folder in the parent directory to write the files into. You can also write your output to the current directory, but then your input/output might get mixed, so I would avoid that.

[sunnykevin97]

I executed the awk in the directory, containing all the files.

I created dir with "out"

cat *.fasta | tr "\n" ";" | tr ">" "\n>" | awk -F ";" '{print ">"$4"\n"$5 > "../out/$5.fasta"}' 

awk: cmd. line:1: (FILENAME=- FNR=1) fatal: cannot redirect to `../out/$5.fasta': No such file or directory

[Matthias Zepper]

I slightly edited my initial reply since the direct redirect to a mix of string and field variable doesn't seem to work. So it is now assigned to a variable a for being written out.

But this doesn't explain the "No such file or directory" error, because it would have created a file with the name $5.fasta in the out directory instead. Did you create the out folder in the parent directory or as a subfolder in the current directory? I am almost sure that this was the problem. Try again after running mkdir -p ../out.

[sunnykevin97]

Thanks Matthias,

It generates a concatenated file ($5.fasta), that is true.

with a fixed string length

head -n 6 $5.fasta

> trnW(tga)
AGAGACTTAGGCTAATATAAAACCAAGAGCCTTCAAAGCCCTAAATGAAAGTGAAAATCC
> trnA(gca)
AAAGTTTTAGCTTAATTAAAGTGTCTGTTTTGCGTACAGAAGATGTGGGTTAGTGTCCTG
> trnN(aac)
TAGATGGAGGCTCCTTGGTTTGAGCGTTTAGCTGTTAACTAAGAGTTTGTAGGATCGAAG
> trnC(tgc)
AGTCCCATGGTGTAACATATAAGATTGCAAATCTTAAGACGCAGATTAATATTTGCTGGG

> atp6
 GCTCTAGCTATTTCTCTTCCTTGATTAATATTCCCTGCCCCTTCAACTCGATGATTAAAT

> atp8
CTAATTATTCTTCCCCCTAAAGTGATTGCTCATACTTTCCCAAATGAACCAACCCTACAA

I'm looking for the whole sequence, how I do it ?

Further, the concatenated file needs to be splitted in to individual fasta files.

[Matthias Zepper]

That is weird, because tr should call the corresponding command to replace the respective letters:

echo "DNA" | tr "D" "R"
RNA 

But then let's go step by step and without using tr.

Does

cat *.fasta | paste - - | sed 's/;/\t/g'

give you

Gloriosasuperba  8324-9004   -   atp6   ATGACAGTAAGCCTTTTTGACCAATTTATGAGCCCCACACTACTAGGCATCCCCCTGCTC

and

cat *.fasta | paste - - | sed 's/;/\t/g' | awk '{fields[NF] = NF};END{for(i in fields) print fields[i]}'

returns only 5?

[sunnykevin97]

Yes, you're right.

[Matthias Zepper]

If you only have entries with 5 fields and the contents of each FASTA file are now on one line, this should then give you the desired output:

 cat *.fasta | paste - - | sed 's/;/\t/g' | awk -F ";" '{a="../out/$4.fasta" ; print ">"$4"\n"$5 > a}'

If you have multiple sequences for the same gene (e.g. transcripts), then use >> a such that they are concatenated.

[sunnykevin97]

In the cmd, what does the tr $4 and $5 represents for ?

'{print ">"$4"\n"$5 > "../out/$5.fasta"}'

[Matthias Zepper]

$4 and $5 are the respective columns.

echo "aaaa bbbb cccc dddd" | awk '{print $1 > $2; print $3 > $4}'

should give you the file bbbb with the content aaaa and the file dddd with the content cccc.

So what I am attempting is:

• get all of the FASTA header plus the actual DNA sequence in one line.
• combine the columns accordingly and reintroduce the newline "\n" characters at the appropriate spot.
• print that to a file, using the column with the gene name as output file name.