Searching for paired reads around a location

[Ash]

Hi folks,

I'm looking at the 1000 genomes data using samtools. So far the flag filter has been great in helping me locate regions of interest, but for what I want to do next, I'm not sure this is the easiest/best option... Is there anyway to search for paired end reads that span a specific location, ie, one read would be to the left of the location in question, and the other read would be to the right of the location in question? Both reads would have to map properly. I suspect that at certain locations individuals are polymorphic for an insert, and the presence/absence of such reads as described above would confirm/deny this.

Apologies if this has already been addressed - I suspect I may not be searching using the correct terminology as I find it hard to believe I'm the only person interested in this!

All help gratefully received!

A

[swbarnes2]

Let's say your center point is 1000.

samtools view -f 34 -r chr:500-1000 my.bam > forward.sam
samtools view -f 18 -r chr:1001-1500 my.bam > reverse.sam

a flag of 34 = properly paired, and the mate is reversed, meaning the read itself is forward, since its properly paired. So it's left of your center point. A flag of 18 means properly paired, and reversed, and those fall after the center point. You want the reads whose names appear on both lists.

[Ash]

Many thanks for your prompt response, swbarnes2. I'd gotten as far as having to do 2 searches using :

samtools view -f 3 -F 1564 my.bam chr:500-1000 > forward.sam
samtools view -f 19 -F 1548 my.bam chr:1001-1500 > reverse.sam

which I believe does the same as what you've described above(?). I thought there may be some more sophisticated way of doing this that doesn't require you having to compare 2 files for each position, but I guess not.

All the best,

A