How can I select reads Based on the best MAPQ and exclude reads with the same name and MAPQ value?

Question

[jleandroj]

I have a sam file. I want to select the reads with the best MAPQ (mapping quality).

Firstly, for example, I can have three reads with the same name but these have different MAPQ value, therefore I want to exclude the two reads with less MAPQ value and select the read with the best MAPQ.

Second, if there are three reads with the same name and the same number of MAPQ I want to exclude these three reads.

Thanks

Answer 1

[Devon Ryan]

For the first one you can more simply exclude secondary alignments (bit 256 in the flag, see -f and -F).

For the second requirement, you'll either need to know the range of possible MAPQ values where this can occur for the aligner you're using and exclude them or write a script to do this (e.g., in python with pysam).

[jleandroj]

Thanks Devon Can I do it in shell (awk, perl)? Are there any way to do it using samtools?

[Devon Ryan]

For the second requirement you can undoubtedly do that with perl. If the file is name sorted then you might be able to put something together with awk, but it'd be more trouble than it's worth. Stick to perl if that's what you know (or play with the code that Pierre posted if you're OK with javascript).

Answer 2

[Pierre Lindenbaum]

For fun, using nashorn (the java-based javascript engine) and htsjdk:

(not tested, It needs real data... )