Ignore Reads In Bwa
1
1
Entering edit mode
9.6 years ago
Ric ▴ 360

Hello, Is there a way to ignore reads which can align in multiple position equally good using bwa?

Thank you in advance.

bwa paired galaxy samtools • 7.2k views
ADD COMMENT
0
Entering edit mode

Optional flags in BWA, that is, XT:A:U flag in the sam file denotes unique read and XT:A:R denotes multiple mappings for that read. For paired-end reads, you might also want to consider the flag XT:A:M (one-mate recovered) which means that one of the pairs is uniquely mapped and the other isn't.

ADD REPLY
0
Entering edit mode

Thank you. Where did you get the flags from? I am interested on paired-end where both reads mapped unique. Which command have I to use?

ADD REPLY
0
Entering edit mode

I usually considered XT:A:U flag alone and discarded any other. But on seqanswers.com I came to know that people retain both XT:A:U and XT:A:M and removing any other read. So, I'd suggest that. If you're using perl, you can just use a regular expression just like that: if( $_ =~ /XT:A:U/ || $_ =~ /XT:A:M/ ) { print FILE_PTR $_ } else { next; }

ADD REPLY
0
Entering edit mode

Because you want both reads to be mapped uniquely (essentially all reads are mapped uniquely then), you can just use awk and regular expression as follows: awk '$0 ~ /XT:A:U/ {print}' in.sam > out.sam

ADD REPLY
0
Entering edit mode
9.6 years ago
Wen.Huang ★ 1.2k

bwa assigns non-unique reads a "0" mapping quality, if you filter by mapping quality by "samtools view -q 1", this gives you uniquely mapped reads.

ADD COMMENT
1
Entering edit mode

I think thats not true anymore. I see many reads which map to multiple location and have mapping quality in 10-23 range.

ADD REPLY
0
Entering edit mode

Thank you. Does non-unique reads mean that they match on different position equally good? But "1" mapping quality is different to what MarkDuplicates does? How to include in the above samtools command only paired-end (where both are unique)?

ADD REPLY

Login before adding your answer.

Traffic: 2488 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