How to replace sequence identifiers in a fasta with OTU IDs from another file?

Question

I'm pretty new to Unix and bioinformatics and having a hard time accomplishing the following. I have one FASTA file with sequences and headers, and one OTU ID mapping file (.txt) with OTU IDs and their sequences. The sequence header in the .fasta file partially matches the sequence listed in the OTU ID file. What I want to do is replace the current headers in the FASTA file with just the OTU IDs from the .txt file. The files look like something like this:

OTU_ID.txt

denovo123 U1019.1010.2020_20384
denovo453 U1472.2992.2020_20443
denovo763 U1928.1898.2020_20343

.FASTA

>U1019.1010.2020_20384 0000-AHDJ-38347_223
GAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGA
>U1928.1898.2020_20343 0000-AHDJ-24247_242
GAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGA

and my desired output is that the .fasta file would look like this:

>denovo123
GAAATACCGAGAGAGAAATACCGAGAGAGAAAGAAATACCGAGAGAGAAATACCGAGAAATACCGAGAGAGAAATACCGAGAGA
>denovo763
GAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGA

Is it possible to do this in the command line?

Answer 1

You can do it with bash quiet easy. First create a backup of your fasta file with

cp otu.fasta otu_backup.fasta

Then you can run:

while IFS=" " read -r otu fastaId rest; 
do
sed -E -i "s/>$fastaId .+/>$otu/" otu.fasta; done < ids.txt;
done < OTU_ID.txt

What does the code?

With IFS=" " you declare a whitespace as field seperator and by reading each line of your OTU_ID.txt, you have the two variables otu and fastaId which hold the otu name and the id of the Fasta Header.

Within the while loop, sed is looking each loop for the corresponding id within a RegEx and substitute it with the corresponding otu name. By using the -i option, you overwrite the existing file. This is why you need to create a backup.

Answer 2

You can put you fasta file in a single line per sequence:

perl -pe 'if(/\>/){s/\n/\t/; s/>//; s/\n//}' file.fasta > file.tsv

This line will create the next output:

U1019.1010.2020_20384 0000-AHDJ-38347_223   GAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGA
U1928.1898.2020_20343 0000-AHDJ-24247_242   GAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGA

Then with awk you can match columns separated with tab:

awk -F' ' 'FNR==NR{a[$2]=$1;next} ($1 in a) {print ">"a[$1]"\n"$3}' otu_id.txt file.tsv

You will get this output:

>denovo123
GAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGA
>denovo763
AAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGAGAAATACCGAGAGA

You can add $1 to check if you are getting the right id.

awk -F' ' 'FNR==NR{a[$2]=$1;next} ($1 in a) {print ">"a[$1], $1, "\n"$3}' otu_id.txt file.tsv

If everything is alright you can redirect to a new file:

awk -F' ' 'FNR==NR{a[$2]=$1;next} ($1 in a) {print ">"a[$1]"\n"$3}' otu_id.txt file.tsv > file_ids.fasta

How it works, well awk implicitly loops through each file, one line at a time.

FNR==NR{a[$2]=$1;next}

NR is the number of lines that awk has read so far and FNR is the number of lines that awk has read from the current file. Thus, if FNR==NR, we are still reading the first-named file: otu_id.txt. For every line in file.tsv, we assign a[$2]=$1.

Here, a is an associative array, and a[$2]=$1 means saving file.tsv first column, denoted $1, as a value in array a using file.tsv second column, $2, as the key.

next tells awk to skip the rest of the commands and start over with the next line.

($1 in a) {print ">"a[$1]"\n"$3}

This means that we are reading the second file: otu_id.tsv. If match the first field of the line in file.tsv, as determined by the contents of array a, it will print out a line with the values of fields a[$1] and $3 separated with a new line character.