Checking that two multiline FASTAs are identical (allowing for different order and lower/upper case)

Question

Is there a good CLI tool that can diff two multiline fasta files, ignoring order, and comparing sequences with identical names with each other pairwise.

Bonus points if lower/upper case is tolerated (as in, it's considered identical) and if also indels can be tolerated, maybe through a flag?

I could script this in Python myself, but I was wondering whether there was a custom tool for this purpose already out there. Maybe something more performant than a Biopython script.

If there's no custom tool, I'd also be happy with a short script that uses a number of tools, e.g. seqkit samtools diff etc.

Answer 1

There are a lot of ambiguities in the question, but I'll provide one approach by means of a fasta/fastq sorting program (see here), using process substitution to avoid intermediary files. This will only work for small-ish fasta files, which can fit into computer memory:

diff -q <(fastx-sort.pl a.fasta) <(fastx-sort.pl b.fasta)

Note that this will compare the entire file. If there are any sequences in one file but not the other, then the result will be that the files are different.

The -q option for diff will just output whether or not the two files differ; you didn't specify in your question whether you wanted any more detail. This option can be removed to show differences, or replaced with -u to add context. Context and/or line changes are not necessarily going to be helpful for multi-line fasta sequences, because the compared files may have a different order from the input, and the context lines may not include the name of the sequence. This could be fixed by changing the files to single-line fasta format (e.g. by piping through fasta_formatter from the FASTX-Toolkit):

diff -u <(fastx-sort.pl a.fasta | fasta_formatter) <(fastx-sort.pl b.fasta | fasta_formatter)

Note that this approach will output the entire sequence for any different sequences, together with adjacent context (including the sequence name). There are ways round that as well... but it would be far better if you could make your question more specific so that we don't need to explore the universe when answering.

If you want to ignore case, convert everything to upper (or lower) case:

diff -q <(fastx-sort.pl a.fasta | perl -pe 'tr/a-z/A-Z/') <(fastx-sort.pl b.fasta | perl -pe 'tr/a-z/A-Z/')

Answer 2

Try md5fa from samtools. It computes two MD5s with and without ordering, respectively. Alternatively, there is seqkit sum, which is perhaps more powerful.