How To Get The Correct Alignment Statistics After Using Bwa Mem With -M Option
1
0
Entering edit mode
7.3 years ago
skm770 ▴ 150

I want to get uniquely mapped reads from bam files aligned with bwa-mem using the -M option (to make it compatible with picard-tools). If I try to use the following command with samtools to get uniquely mapping reads( or reads with mapping quality greater than 10. I get a very high number of around 90% uniquely mapping. I know I can decrease this number using a more strong threshold of mapping quality of 30 or so but I want to know what would be the ideal way to do so and if the bwa-mem needs to be run with different set of parameters to get more reliable alignments. Following are the bwa and samtools commands that I am using.

bwa mem -t 4  -M ref_genome.fasta  fastq[0] fastq[1] > samfile
samtools view -c -q10 -f 1 -F 12 -S samfile

Thanks

picard samtools • 4.0k views
ADD COMMENT
0
Entering edit mode

You can get non-secondary alignments with samtools view -F 0x100 (see: http://samtools.sourceforge.net/SAMv1.pdf)

ADD REPLY
0
Entering edit mode
7.3 years ago
Pavel Senin ★ 1.9k

I am not sure if I understand the purpose of the procedure aiming on the reliability of an alignment right, but you can make bwa to search for seeds in more rigorous way - exhaustively - by lowering the r paremeters which is set by default to 1.5. Alternatively or in addition, you could adjust A and B parameters to penalize mismatches more aggressively - if your reads are of good quality, then it makes sense not to align a read if there is mismatch, IMO. Then, you can scrutinize what is going on with gaps and tighten that too as well as clipping and unpaired read penalties.

ADD COMMENT

Login before adding your answer.

Traffic: 2111 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6