intersectbed between a bam and a bed file gives different results with or without -sorted mode

[wang...@gmail.com]

Hi, I have a sorted bam file and a bed file.  I am trying to find all the reads in the bam file that overlap with the entries in the bed file.

I tried two ways:

1) intersectBed -abam sample_sorted.bam -b interval.bed -f 1 > intersect.bam

2) Now, I sorted the bed file based on the order of chromosomes appearing in the bam header using bedtools sort:

      bedtools sort -faidx chromo_order_sam_header.txt -i interval.bed > interval_sorted.bed

 And now I do intersect in the sorted mode (chromosome_and_size.txt that is fed to -g is retrieved from UCSC using the fetchChromSizes binary):

      intersectBed -abam sample_sorted.bam -b interval_sorted.bed -sorted -g chromosome_and_size.txt -f 1 > intersect_sorted.bam

But the results are different.  From a quick check on genome browser I found 1) gives the right answer while 2) is missing quite a few reads.  Have you observed this before and any idea why this happens?  How can I make the two match?

Qingqing

[Braulio Valdebenito Maturana]

Hi there!

I’m about to write about something like this soon. I have observed something like this before. I have gotten some unexpected results when using some BAM files (i.e., less results). The workaround I used was to convert them to BED format, since I was interesting in counting overlaps. 

Could you check the number of overlaps you get by converting the BAM files into BED files (with bedtools bamtobed)?

Which version of BEDtools are you using?

Best regards,

Braulio.

-- 
You received this message because you are subscribed to the Google Groups "bedtools-discuss" group.
To unsubscribe from this group and stop receiving emails from it, send an email to bedtools-discu...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Qingqing Wang]

Whoops, sorry, I forgot I should add a -u parameter when I compare bed intersect bed, in comparison to bam intersect bed, since bed intersect bed will print multiple original entries in bed if intersected, while bam intersect bed will just print out one read one time if intersected.

Yes the workaround works!!!  Thanks a lot :)

Qingqing

On Fri, Jul 6, 2018 at 10:21 AM, Qingqing Wang <wang...@gmail.com> wrote:

Yes I did try the workaround.

bedtools bamtobed -i sample_sorted.bam > samplebam.bed

sort -k1,1 -k2,2n samplebam.bed > samplebam_sorted.bed

sort -k1,1 -k2,2n interval.bed > interval_sorted.bed

intersectBed -a samplebam_sorted.bed -b interval_sorted.bed -sorted -wa -f 1 > intersect.txt

However, intersect.txt reported three times more reads than intersect with bam and interval.bed directly without the sorted mode...  Please help...

I am using v2.26.0.

Qingqing

To unsubscribe from this group and stop receiving emails from it, send an email to bedtools-discuss+unsubscribe...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.

--
You received this message because you are subscribed to the Google Groups "bedtools-discuss" group.

To unsubscribe from this group and stop receiving emails from it, send an email to bedtools-discuss+unsubscribe@googlegroups.com.

[bvaldebenitom]

Hi there!

I’m about to write about something like this soon. I have observed something like this before. I have gotten some unexpected results when using some BAM files (i.e., less results). The workaround I used was to convert them to BED format, since I was interesting in counting overlaps. 

Could you check the number of overlaps you get by converting the BAM files into BED files (with bedtools bamtobed)?

Which version of BEDtools are you using?

Best regards,

Braulio.

[bvaldebenitom]

Glad you were able to make it work.

Best regards,

Braulio.

To unsubscribe from this group and stop receiving emails from it, send an email to bedtools-discuss+unsub...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.

--
You received this message because you are subscribed to the Google Groups "bedtools-discuss" group.

To unsubscribe from this group and stop receiving emails from it, send an email to bedtools-discu...@googlegroups.com.