Extract reads by type +/-

Question

[ksadree]

Dear all,

I am trying to extract reads from a bam file that align to genes in one direction only; using samtools. See image:

https://postimg.org/image/inh0uxq9z/

The reads to extract are circled.

thanks for your help!

[Brian Bushnell]

How was the mapping done? The process will be different for genome versus transcriptome mapping.

[ksadree]

Reads were aligned against transcriptome with Bowtie2

Answer 1

[mbk0asis]

check this page

[ksadree]

I tried something like that, but it extracts only + or - reads. I am interested in reads that coincide with the direction for the gene. thanks.

[Brian Bushnell]

Since you mapped to the transcriptome, all + reads will be in the same direction as the gene, and minus reads will be in the opposite direction. It only gets complicated when you map to the genome. You can also use BBMap's SplitSam tool for this purpose:

splitsam.sh mapped.sam plus.sam minus.sam unmapped.sam

This will produce 3 output files from the input file "mapped.sam", such that "plus.sam" contains reads in the forward direction, and "minus.sam" contains reverse-direction reads.

P.S. It works on bam files as well, if samtools is installed.

Answer 2

[swbarnes2]

In .bam format, reads which are aligned in the reverse direction have the '16' flag set, reads in the forward do not. samtools view can pull out reads based on whether or not they have the 16 flagged.