This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Get the coverage after read alignment with bwa-mem2

Hi,

I have aligned genomic reads to a reference genome using bwa-mem2 (using bwa-mem2 mem). From what I understood, if a read map to 2 different locations in the genome with the same score, bwa-mem2 will choose at random one of the location and this will have a MAPQ score of 0 (I did not run mem with the -a option).

What I want to do then is to get the mean read coverage at different part of the genome using bedtools (Thats why I absolutely need reads to be mapped to only 1 genomic location) and compare these different coverages.

However, I have several questions :

  • Is my first asumption about bwa-mem correct ?
  • Will bedtools coverage take into account reads flagged with a MAPQ of 0 ?

Thanks for any help and tips !

Maxime

coverage alignment bwa-mem2 mapq

You may check out deeptools.

1 answer

1) yes

2)

$ samtools depth -a -r "chr1:123-456" in.bam -q 0  | awk '{D+=int($3)} END{print D/(NR*1.0);}'

Log in to answer this question.