2
$\begingroup$

when mapping short reads with bwa mem, if a read has alternative mapping positions they are reported by bwa mem in the X0 and XA tags. Now, let's say I want to compute the coverage of my bam file. How would I do for the alternative hits to be taken into accounts. What I mean is that, if a read QUERY has a selected best HIT position and 3 alternate HIT positions, I would like that QUERY read to be counted as 1 unit of coverage for the best HIT position, but also as 1 unit of coverage for each of the alternate HIT position.

My understanding is that samtools depth will only consider the best hit (which makes sense in 99.99% of the cases).

I guess I could edit the bam file with custom scripts so that the alternative hits are also reported as primary hits but I was wondering if there was an already existing tool (as I am not that confident with my ability to mess around with sam files).

Thanks a lot.

$\endgroup$
1
$\begingroup$

You would want to instead use the -a option with bwa mem to output all alignments. You would be best off filtering these with a custom script to only keep those with a close enough score to the primary hit. samtools depth can't be told to include secondary alignments (it's hard-coded), so you'll need to use something else like bamCoverage or samtools mpileup.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.