Extracting Reads by position
1
1
Entering edit mode
8.7 years ago

I have an alignment file for a protein encoded by about ~3000 bases (I have both bam and sam files), and I am looking to extract all the reads that from bases 877-981. Is there any way to do this with samtools?

I tried using this command:

samtools view 8.2alnCorrected.bam "877-981" > filteredreads.bam

However I got this message with an empty file:

region "877-981" specifies an unknown reference name. Continue anyway.
bam sam alignment • 4.5k views
ADD COMMENT
1
Entering edit mode

Take a look at the manual to correctly specify the region (RNAME[:STARTPOS[-ENDPOS]]): http://www.htslib.org/doc/samtools.html

ADD REPLY
2
Entering edit mode
8.7 years ago

You need to specify the chromosome or contig as well. Your bam will have lines that look like:

chr1     12345    ...

So for reads are mapped to that chromosome or contig, you want:

samtools view 8.2alnCorrected.bam chr1:877-981
ADD COMMENT
0
Entering edit mode

It worked out. Thanks Chris!

ADD REPLY
0
Entering edit mode

Glad to hear it! Consider hitting the "Accept" button next to the answer - that lets anyone else who stumbles across this question find the solution easily.

ADD REPLY

Login before adding your answer.

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