Introduction

Library construction of standard directional BS-Seq samples is a multi-step process typically consisting of several consecutive steps:

  1. Sonication
  2. End-repair
  3. A-tailing
  4. Adapter ligation
  5. Bisulfite conversion

The end-repair step typically uses a mix of nucleotides to fill in fragments with a 5′ overhang which are frequently generated by sonication, leaving blunt ended fragments ready for the A-tailing procedure. This would not be apparent in other applications such as ChIP-Seq or RNA-Seq because the genomic sequence is simply regenerated in these fragments. The situation changes in BS-Seq applications though because it does not really interrogate the sequence itself but the methylation state of cytosine residues. The nucleotide-mix used in the lab typically uses a standard mix of the four nucleotides C, A, T and G for the fill-in reaction, whereby C is typically unmethylated for cost reason.

The Symptoms

This schematic drawing shows that cytosines located towards the 3′ end of reads may be filled in with unmethylated cytosines during the end-repair step. Since the fragments undergo bisulfite treatment after the adapters have been ligated on, these cytosines (marked in red) will be converted to T. Since Read 2 starts right from this end carrying the artificial bases they will be (incorrectly) read as unmethylated cytosines. If the sequenced fragments happen to be very short you might also see a read these filled-in positions at the 3′ end of Read 1, but this will be arguably quite rare, whereas the Read 2 bias definitely affects every single molecule with filled in cytosines.

endrepair

It is actually pretty much impossible to spot the symptoms of this bias unless you really look for them as a post-alignment QC step (see below).

Diagnosis

The trouble with this kind of bias is that it only becomes apparent when you really look for it. To spot the bias one has to look at the average methylation levels across the entire length of the reads and plot out a methylation bias (or M-bias) graph:

R1

Read 1 looks fairly normal with a little bit of a ‘wobble’ at the 5-end but generally fairly even methylation levels across the entire read.

 

Read 2 however tells an entirely different story:

R2 allIn

The methylation state of the first couple of bases or Read 2 drops from the average ~70% down to ~3%, and this steep drop can certainly be explained by the filled-in unmethylated cytosines. Since there is no reason to believe that there would be a biological role for lower methylation at the start of reads we have to assume that this low-methylation is artefactual, and thus leaving this untreated will introduce, in this case several hundred thousand, incorrect methylation calls (and thus noise) into the dataset.

Mitigation

Fortunately it is relatively easy to fix this problem by simply disregarding the biased positions when running the methylation extractor. Since this problem only affects Read 2 of paired-end libraries a typical command like this should fix the problem:

bismark_methylation_extractor --ignore_r2 2 --gzip sample1_bismark_bt2_pe.bam

 

The same M-bias plot for Read 2 should then look like this:

R2 remedied

Prevention

Since this is a common phenomenon with any standard directional library the option –ignore_r2 2 of the methylation extractor should be added to the default data processing pipeline for paired-end applications.

Lessons Learnt

Read 2 of paired-end BS-Seq libraries suffers from artificial hypomethylation at the start of Read 2 of paired-end reads. This will add spurious hypomethylated calls and should be removed from the dataset before proceeding to downstream analysis.

References

The M-bias plot was first described in this paper: Hansen et al., Genome Biology, 2012, 13:R83

Datasets

Any paired-end directional BS-Seq dataset will show this problem.

Software

To get to the stage where you can see the M-bias plot data was adapter and quality trimmed with Trim Galore and aligned to the mouse genome using Bismark. Methylation information was extracted using the bismark_methylation_extractor and reports shown in this article were generated using bismark2report.

February 12, 2016

2 thoughts on "Library end-repair reaction introduces methylation biases in paired-end (PE) Bisulfite-Seq applications"

  • Xinlei Hu

    In Introduction part , it says that fragments with a 3′ overhang are filled in with nucleotides. But in The Symptoms part, it shows that the 3′ end of reads may be filled in with unmethylated cytosines—-I mean it’s 5′ overhang. According to Illumina’s support file of WGBS, the 3ʹ to 5ʹ exonuclease activity of the kit removes the 3ʹ overhangs and the polymerase activity fills in the 5ʹ overhangs. I get a little confused that which kind of overhang is filled in, or both overhangs are filled in .

  • Felix Krueger

    You are absolutely right that it should read 5′ overhang in the introduction section, I have now updated the main article. Thanks for spotting this!

Leave a Reply

Your email address will not be published.