I have a region on the genome that was PCR amplified and bisulfite sequenced using ion-proton sequencer to detect methylation site. Can some one suggest how to align these reads back to the genome.
I came across BSSeeker however, their website do not work.... The deal here is it is just two regions on the genome that we are interested in, while most of the softwares are meant for the whole genome sequencing?
Any directions to detect methylation sites for such short stretches would be really great!
They're totally overkill for the amount of data that you have, but you can use bismark or bison (of which I'm the author and it'd be even more overkill for you than bismark) to do the alignments. Since you have targeted sequencing, you might just create a miniature fasta file containing just the targeted regions and align to that. It could be that someone has created a better solution for your use-case, but I haven't seen one.
BSmap will be by far the slowest, followed by bismark and then bison. The effective speed difference between bismark and bison will generally depend on the resources you have. I made bison with clusters in mind, so each instance can use a few nodes. If, instead, you just have a single server that's being used, the speed benefit with bison will generally no be that big (it depends on how many threads you give bowtie, fewer threads=smaller difference, at least in my hands). Regarding accuracy, I don't recall seeing a comparison against BSmap, though the bismark paper suggests that BSmap would likely perform the worst. In all of my tests, bison is more accurate than bismark, but the difference is very small. Bison's error rate ends up being that of bowtie2, which occasionally returns alignments with a high MAPQ that are actually completely wrong and should have a low MAPQ. I should note that in the next bison release, all MAPQs will be recalculated (the current version just passed the MAPQ through, which was a bad idea, though I guess it might still be better than not providing them, as done by bismark) and the methylation extractor has an option to only include reads with a minimum MAPQ. Bismark's error rate is only ever so slightly higher because of how it chooses between possible alignments in cases where there are valid alignments to multiple real and virtual strands. Bismark is quite flexible, so it ends up fitting more use cases than bison (I can't really say how BSmap compares in this regard). For someone that is only doing a one-off experiment, then the easier setup of bismark would make it a better choice. If you have more samples to process, then bison can make more sense.
For locus specific targetted bisulfite sequencing in HT (as in my case), BiQAnalyser HT really works well. It infact supports both genome wide and locus specific and is written in Java and so compatible with multiple operating systems. This works with needleman-wunsch algorithm which could be a potential limitation when working with bizzare gene sequence regions. There are several parameters that can be tweeked to increase the accuracy of the methylation call. The program is relatively fair in speed. More detailed information is available here in the paper
For regular sanger sequencing kind of approach, a similar tool called BiQAnalyser is avalaible
We have carried out targeted sequencing using ION TORRENT PGM of the PCR products amplified from bisulfite converted whole blood DNA. We have aligned the reads and now analysing the reads against reference sequences using BiQ Analyzer HiMod. But the problem is, can we just rely on the mean modification levels for every potential modification site of very sample for one specified amplicon generated by the software. We tried calculating the methylation percentage from the [0,1] we get from result files, however less reads per CpG site may generate higher methylation percentage. Therefore is there any way to deal with the problem. What cut off I should give for min number of reads per site?
BS Seeker is right here; you have the wrong URL.
Hi Devon,
Yes, we are comparing groups of samples. Thank You for your suggestions. We will try to use beta-binomial based method.