Excessive clipping when mapping to a genomic fragment

[limepoutine]

I have four genomic fragments from a tetraploid potato and want to identify splicing junctions computationally using NM_001247643.3 as a query. MSA suggests the sequences match well except for a few indels. However, magicblast seems to be clipping the query sequence excessively.

Running magicblast -query NM_001247643.fasta -subject Fragment1.fasta -out Fragment1.sam -limit_lookup F -score 20 generates CIGAR 510S531M108N51M. The expected result is that the start codon matches.

The sequences are attached here as a text file: sequences.txt.

[boratyng]

Hi @limepoutine, I am sorry you ran into problems. Yes, magicblast clips excessively here. Unfortunately the fix requires code changes. What system are you running magicblast on? How did you install it? I will see if I can get you a pre-release binary with this issue fixed.

[limepoutine]

I am running magicblast on Linux with conda. I do have the privilege to run any binaries on that cluster, so any channel would be OK for me. Thanks for the quick follow-up!

[boratyng]

Here is the Linux binary: https://ftp.ncbi.nlm.nih.gov/blast/temp/ncbi-magicblast-1.7.2-x64-linux/magicblast. Please, let us know if it works for you.

[limepoutine]

Hi, I can confirm the issue described has been fixed. However, I ran into another issue during testing. Magicblast appears to generate different alignments based on the direction of subject sequence.

When I run ./magicblast -query NM_001247643.fasta -subject Fragment2.fasta -limit_lookup F -score 20, the SAM file is now fully consistent with my expectation:

NM_001247643.3	0	Fragment2	103	60	28M6D1013M106N51M	*	0	0	 ATGTATCCTCCAAGCAACAATAACTGCAGCCCAATTTTGTCTTCTTTGATATGCCAAAATATTCCATCTTCTCCTTGTATGCAATATGAACACGAACTATACTTTCAAAGCTTTAATCATGATAACCAATATTATTTTCAACAACAGCAACTAGTTCCCTCGATAGATGATTTGAGTCCTCACATCTTAGCTGACAGCTGCACCGAGATTATTACTAAGCCTTCGAATTGCAACCACGAACTACAAGGAATGGAAGAAGGCCGAGGCGAAAAGAAAGGAGATGATGATGTTATGAGTAGCAGAATTAGTGGACGGATCTCAAAAAATAATAAGAGATCTTCCAATAAAGATCGACACAGCAAGATCAACACCGCTCGTGGTCCAAGAGATCGAAGGATGAGACTTTCACTTGATGCTGCTCGCAAGTTTTTCCGTTTGCAGGACTTATTGGGATTCGATAAGGCCAGCAAAACTGTTGAATGGTTGCTTACTCAATCGGATTCTGCAATTGAAGAGCTTGTTGCCGCTAAAGGCAATGATGCACAGGTTGCTCAGCAAACTAGCTGCAATACCCCCACTACTACTACTGGAATTGGTGCAATTTGTGCATCTAATTCTATTTCTGAGTCGTGTGAAGTTATATCAGGAACTGATGAAACTTCCTCTAATGACAAAAACAAGGAAACCGCTCAAGATGAGGAGAAGAAGAAAAGGAAGAAGGTGGTTAACACAGCTCGTAGAGCTGTGTTAGAACCTCTTACGAAGGAATCGAGGAATCAAGCAAGAGCCAGGGCTAGAGAGAGAACAAAATCAAAGAAAATGAGCCAAACTGGAAAATCCAAATCCCTAGCTAATGATTTGAACCCTTCAGGATCTCGGAGGCCGGCTAATAAAACTTGTGAAGAACCTGGAACACATGAAGAACTCAACTTCCATCAAGAGAAGAACACTGTCGATGACTGTAATTTTATGGTAAATGGAAATTGGAATCCATTTACAATCTTTAGCTATCATGAGCAATACGCTGGAATTTCCAACGAGCATCAATTGGTTACAGACTTGCAATTTTGTGGAAAGCTATGGGAAGGCTAG	*	NH:i:1	AS:i:986	NM:i:24	XS:A:+

When genomic fragment is from the template strand, magicblast tends to overextend and miss the splice signal. This is the SAM output when I replace Fragment2 with its reverse complement.

NM_001247643.3	16	Fragment2	11	60	53M2I108D1016M6D21M	*	0	0	CTAGCCTTCCCATAGCTTTCCACAAAATTGCAAGTCTGTAACCAATTGATGCTCGTTGGAAATTCCAGCGTATTGCTCATGATAGCTAAAGATTGTAAATGGATTCCAATTTCCATTTACCATAAAATTACAGTCATCGACAGTGTTCTTCTCTTGATGGAAGTTGAGTTCTTCATGTGTTCCAGGTTCTTCACAAGTTTTATTAGCCGGCCTCCGAGATCCTGAAGGGTTCAAATCATTAGCTAGGGATTTGGATTTTCCAGTTTGGCTCATTTTCTTTGATTTTGTTCTCTCTCTAGCCCTGGCTCTTGCTTGATTCCTCGATTCCTTCGTAAGAGGTTCTAACACAGCTCTACGAGCTGTGTTAACCACCTTCTTCCTTTTCTTCTTCTCCTCATCTTGAGCGGTTTCCTTGTTTTTGTCATTAGAGGAAGTTTCATCAGTTCCTGATATAACTTCACACGACTCAGAAATAGAATTAGATGCACAAATTGCACCAATTCCAGTAGTAGTAGTGGGGGTATTGCAGCTAGTTTGCTGAGCAACCTGTGCATCATTGCCTTTAGCGGCAACAAGCTCTTCAATTGCAGAATCCGATTGAGTAAGCAACCATTCAACAGTTTTGCTGGCCTTATCGAATCCCAATAAGTCCTGCAAACGGAAAAACTTGCGAGCAGCATCAAGTGAAAGTCTCATCCTTCGATCTCTTGGACCACGAGCGGTGTTGATCTTGCTGTGTCGATCTTTATTGGAAGATCTCTTATTATTTTTTGAGATCCGTCCACTAATTCTGCTACTCATAACATCATCATCTCCTTTCTTTTCGCCTCGGCCTTCTTCCATTCCTTGTAGTTCGTGGTTGCAATTCGAAGGCTTAGTAATAATCTCGGTGCAGCTGTCAGCTAAGATGTGAGGACTCAAATCATCTATCGAGGGAACTAGTTGCTGTTGTTGAAAATAATATTGGTTATCATGATTAAAGCTTTGAAAGTATAGTTCGTGTTCATATTGCATACAAGGAGAAGATGGAATATTTTGGCATATCAAAGAAGACAAAATTGGGCTGCAGTTATTGTTGCTTGGAGGATACAT	*	NH:i:1	AS:i:960	NM:i:134

This occurs with multiple previous releases and the pre-release version. I am able to work around this problem by adjusting the direction manually before mapping. Thank you for paying attention to this!

[boratyng]

Thank you for this report I am glad that the binary fixes the previous problem. Unfortunately we do not have a quick fix for the new issue. We will try to fix it in the next magicblast release.