Is there a way to remove split alignments of single-end long reads from BWA-SW output? BWA-MEM has an option to include or flag split alignments in the result, but it seems BWA-SW method include them in the result unless the split occur in paired-end reads.
$\begingroup$
$\endgroup$
Add a comment
|
1 Answer
$\begingroup$
$\endgroup$
From BWA manual, it is advised to use BWA-MEM in your longs read case (unless you have a good reason not to):
BWA-MEM and BWA-SW share similar features such as long-read support and split alignment, but BWA-MEM, which is the latest, is generally recommended for high-quality queries as it is faster and more accurate.
I guess that BWA-SW returns several lines per transcript (based on local alignments). In that case, the issue has already been reported. The user asks:
How can I differentiate between a read that has multiple hits in the reference and a read that has a split alignment?
The author replies:
In principle, it is hard to tell. Bwa-mem makes an arbitrary choice. Mappings flagged with 0x100 are multiple hits. Mappings flagged with 0x800 are from a chimeric alignment.
So you could actually filter by flag. But I am afraid the output could include something you do not necessarily want to (since there is no specific flag for the split reads).
Maybe a better solution: you could switch to some mapper that will report the proper split in the hit of the read map, such as minimap2 from the same author. It is suitable for long reads and comes with a splice -x feature.
Then, the laziest solution I see here would be to remove the split alignments exploiting the sam format and in particular the CIGAR string. Meaning, to remove the sam fields containing a N in their CIGAR (aka, skipped region from the reference). Something like this would work:
cat alignment.sam | awk '$6 !~ /N/ {print $0}'
Hope this can be somewhat helpful.
