Go Back   SEQanswers > Bioinformatics > Bioinformatics

Similar Threads
Thread Thread Starter Forum Replies Last Post
Bismark, bisulfite and chloroplast DNA colindaven Bioinformatics 5 08-12-2014 03:39 AM
Parallelization for Bismark using Bowtie 1 KlChen Bioinformatics 2 10-21-2013 03:14 PM
Questions about Bismark and bisulfite-seq sero Bioinformatics 4 07-11-2013 01:09 PM
Bismark and unique matches mixter Bioinformatics 3 12-12-2012 05:06 AM
Bismark paired-end positions mixter Epigenetics 1 01-19-2012 11:45 AM

Thread Tools
Old 07-22-2014, 08:59 PM   #1
Junior Member
Location: California

Join Date: Sep 2012
Posts: 6
Default Help with Bismark methylation extractor


Could I please get some advice on the following issue. I am working on targted capture genome-wide DNA methylation assay. The aligner I use is Bismark/bowtie2. After obtaining output from bismark, I isolated targeted regions in the bismark BAM output file by using bedtools 'intersect' between the bismark BAM output file and a targeted region bed file. I then ran the bismark methylation extractor. I obtained this error:

The IDs of Read 1 (HWI-D00356:33:C3D9RACXX:1:1101:1482:2234_1:N:0:TGACCA) and Read 2 (HWI-D00356:33:C3D9RACXX:1:1101:1551:2226_1:N:0:TGACCA) are not the same. This might be a result of sorting the paired-end SAM/BAM files by chromosomal position which is not compatible with correct methylation extraction. Please use an unsorted file instead

I have not done any sorting prior to methylation extractor so I am unsure about the cause of this issue. If I ommit bedtools 'intersect' and used bismark methylation extractor to process the bismark BAM output, the methylation extractor works well. It therefore suggests that something is going on with the bedtools 'intersect' which I am not aware of. Could someone please advice?

Thanks very much.
Dipro is offline   Reply With Quote
Old 07-23-2014, 01:39 AM   #2
Senior Member
Location: Cambridge, UK

Join Date: Sep 2009
Posts: 625

Hi Dipro,

If you run the methyation extractor in paired-end mode it expects read 1 and read 2 of a pair to follow each other line by line in the BAM file. If you disturb this R1/R2 order you would team up two reads that have nothing to do with each other, and thus features such as the --no_overlap detection would not work. Also, this will interfere with the strand origin of the reads. Because of this we included a check at the start that would not allow the methylation extraction of such files since it almost never is what you are intending to do. What appears to have happened with bedtools intersect is that it accepted some hits where only R1 or R2 of a pair overlapped with a desired region, but its partner read didn't, thus getting the read pairs out of sync.

As solutions to this you could:

a) run the methylation extractor first and simply use all data for downstream analysis

b) run the methylation extractor first and then filter on the position of the methylation calls afterwards (even though this might not work out of the box with bedtools intersect if it requires BAM files)

c) write a quick script to go through the filtered BAM file that reads in two lines at a time and sorts R1 and R2 into a paired-end file if they follow each other, and sort singleton reads into another file. Then you should be able to process both files in paired-end (-p) and single-end mode (-s), respectively.

I personally would probably stick with a) since it is the most straight forward method and you don't exclude data that is actually there. Let me know if I can be of any help with whatever you choose to do.
fkrueger is offline   Reply With Quote
Old 07-23-2014, 05:18 AM   #3
Junior Member
Location: California

Join Date: Sep 2012
Posts: 6

Hi fkrueger,

Thanks very much for your reply. I'll work on a) as adviced.

Much appreciated.
Dipro is offline   Reply With Quote

bedtools intersect, bismark, methylation extractor

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off

All times are GMT -8. The time now is 12:04 AM.

Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2021, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO