Problem with extracting reverse-reverse read pairs
1
0
Entering edit mode
6 weeks ago
Wang Cong • 0

Hi, I am trying to extract reverse-reverse read pairs from a bam file (all reads are paired and mapped).

I used the following command:

samtools view --expr 'flag.paired && !flag.unmap && !flag.munmap && flag.reverse && flag.mreverse' Pooled_valid.bam -o Pooled_valid.RR.bam

However, I found that all the reads are extracted (the input and output have the exactly same number of reads). May I ask if there is a problem with my command? Thanks

bam read samtools pairs alignment • 263 views
ADD COMMENT
0
Entering edit mode

Please create a Minimal, Reproducible Example to enable others to exactly reproduce your issue on their machines.

ADD REPLY
1
Entering edit mode
6 weeks ago

The command doesnt' seem wrong to me but depending on the version of samtools, it might not work. It is arguably more robust to use the -f and -F options along with the numerical flags as it works with all samtools versions:

samtools view -f 49 -F 12 Pooled_valid.bam -o Pooled_valid.RR.bam

-f 49 means read paired (0x1), read reverse strand (0x10), mate reverse strand (0x20)

-F 12 mean NOT read unmapped (0x4) NOT mate unmapped (0x8)*

ADD COMMENT
0
Entering edit mode

Agreed the -f and -F options will be much faster, as the expression syntax is an interpretted string.

Note you can use -f and -F with strings, eg -f PAIRED,REVERSE,MREVERSE which may be more descriptive if you don't like the bit fields. Also the samtools flags command is a handy guide to the numeric and string names available.

ADD REPLY

Login before adding your answer.

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