Samtools View: Only Forward Or Reverse Strand
3
11
Entering edit mode
13.1 years ago
Gregor Rot ▴ 540

Hello,

i would like to only get reads that aligned to the forward / reverse strand separately. I am using this for now:

samtools view -F 20 ... : forward strand
samtools view -f 0x10 ... : reverse strand

Is this correct? With -F 20 i also ignore the not aligned reads.

Tnx!

Gregor

samtools strand • 35k views
ADD COMMENT
5
Entering edit mode
13.1 years ago
Pascal ★ 1.5k

Gregor,

For me it makes sense because I just verified that the the sum of:

samtools view -c -F 16 [BAM_FILE]
+
samtools view -c -f 16 [BAM_FILE]

is equal to:

samtools view -c [BAM_FILE]

;-)

Pascal

ADD COMMENT
0
Entering edit mode

Hey there,

I want to subtract the reads from forward and reverse strand, from a BAM file of Histone ChiP-seq data, to get the net reads mapping on either forward/reverse strand, how to do that. Any suggestions will be appreciated.

ADD REPLY
3
Entering edit mode
13.1 years ago
Chris Penkett ▴ 490

0x20 is if your data is paired end and will filter based on the strand of the mate. You want 0x10 for the query. For standard NGS data, mostly your pairs should be on the opposite strands anyway - except for weird events like DNA translocations, e.g., inversions.

[There is a slide (Using NGS Data to find chromosomal rearrangements) in this [PDF][1] to explain how the pairs will behave for these events.]

so do this for forward:

samtools view -F 0x10

and this for reverse:

samtools view -f 0x10

Use the -X flag in "samtools view" to convert the flag column from a numeric bit field to letters (as described here) like this:

% samtools view -Xf 0x10 accepted_hits.bam | head -1
24_727_529_F3   r   chr1    2   0   50M *   0   0   ATCCAGCATCCAGCATCCAGCATCCAGCATCCAGCATCCAGCATCCAGCA  @@!!",>4CS/!!;UPB!(QB>=EH5/!>T>9?S^ZUQZaZVY\_]^bcc  NM:i:0  NH:i:31 CC:Z:=  CP:i:9
% samtools view -XF 0x10 accepted_hits.bam | head -1
2118_596_811_F3     chr1    5367    3   50M *   0   0   CAGTGCAATCGTTTATAGACACAAACACCTTTTCCGCTGTTGATGGTGGT  ___ZZHFZ\]PPRUb_ZNQUNOQYYP!!!#ROLPZSRNO@/GH2=NL;..  NM:i:0  NH:i:2  CC:Z:=  CP:i:2868239
ADD COMMENT
0
Entering edit mode

Chris, I think that Gregor was right: observe he did say "20" and not "0x20" as you assumed in your response. 20 = 0x10 + 0x04.

ADD REPLY
0
Entering edit mode

True - it was confusing as he mixed hex and dec values.

ADD REPLY
7
Entering edit mode

I found for using single end reads it could be written also as

samtools view -F 20 ... : forward strand
samtools view -f 16 ... : reverse strand
samtools view -f 4 ... : unmapped

And since I don't see it mentioned on this particular page, I'll add that tool at https://broadinstitute.github.io/picard/explain-flags.html is particularly helpful for making and deciphering these commands. The fact that the flag 20 includes both unmapped and those mapped to the reverse strand was not obvious to me at first. Bitwise is a compact but not immediately obvious way to indicate these sets. Slides 28-31 of Pierre Lindenbaums's presentation, which he highlighted in this post SAM flags meaning helped a lot.

ADD REPLY
1
Entering edit mode

Thanks! Specially, the URL from broadInstitute was very useful.

ADD REPLY
0
Entering edit mode

The option -h is important to further use of bam files. Need to keep the SAM header which is mandatory:

samtools view -F 20 -h ...: forward strand

ADD REPLY
0
Entering edit mode

agree, it's not obvious at all.. wondering what would be the command for PE reads?

ADD REPLY
0
Entering edit mode

I also wonder the command for PE reads

ADD REPLY
1
Entering edit mode
4.6 years ago
brianpenghe ▴ 80

They forgot unmapped reads!

Please use these:

Forward

samtools view -c alignment.bam -h -F 20

Reverse

samtools view -c alignment.bam -h -f 16
ADD COMMENT

Login before adding your answer.

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