MAPQ filtering in BS-Seq pipeline

Question

[int11ap1]

I have aligned my BS-Seq data with Bismark and I removed as well the duplicated fragments (I am following the commands from this link). My question is: do you normally filter those reads with a mapping quality, MAPQ, below a threshold in BS-Seq? Because the link provided does not say anything about it.

Thank you very much.

Answer 1

[Devon Ryan]

Yes, I personally recommend excluding alignments with MAPQ < 10. You can see why in the Bison paper and since bismark uses the same algorithm to calculate MAPQ scores (it just ported my C code to perl) you would be best off doing the same. I also exclude base calls with Phred scores < 5 or so. Presumably bismark's tools allow doing this, though if not you can always use PileOMeth, which can plot methylation bias graphs (and exclude biased regions) as well.

[int11ap1]

Thank you Devon. However, I read in the Bison's paper that Bismark produces slightly different MAPQ values compared to Bison, since the latter recalculates the MAPQ score (besides, all my reads have a score of 255). So, how can they use the same algorithm?

[Devon Ryan]

Bismark was modified after the paper was published (or at least submitted, I don't recall which). If you still have alignments with 255 MAPQ scores then you must be using an older version of bismark.

[int11ap1]

Well, I think it's because I used bowtie1 instead of bowtie2... Do you recommend me to align again using bowtie2?

[Devon Ryan]

Ah, yeah, that would do it. Yes, use bowtie2.

[int11ap1]

Do you recommend me filtering those methylated positions with very low coverage (e.g. supported by 1 read)? Thank you Devon.