Reduced representation bisulfite sequencing for comparative high-resolution DNA methylation analysis

Abstract

We describe a large-scale random approach termed reduced representation bisulfite sequencing (RRBS) for analyzing and comparing genomic methylation patterns. BglII restriction fragments were size-selected to 500-600 bp, equipped with adapters, treated with bisulfite, PCR amplified, cloned and sequenced. We constructed RRBS libraries from murine ES cells and from ES cells lacking DNA methyltransferases Dnmt3a and 3b and with knocked-down (kd) levels of Dnmt1 (Dnmt[1(kd),3a-/-,3b-/-]). Sequencing of 960 RRBS clones from Dnmt[1(kd),3a-/-,3b-/-] cells generated 343 kb of non-redundant bisulfite sequence covering 66212 cytosines in the genome. All but 38 cytosines had been converted to uracil indicating a conversion rate of >99.9%. Of the remaining cytosines 35 were found in CpG and 3 in CpT dinucleotides. Non-CpG methylation was >250-fold reduced compared with wild-type ES cells, consistent with a role for Dnmt3a and/or Dnmt3b in CpA and CpT methylation. Closer inspection revealed neither a consensus sequence around the methylated sites nor evidence for clustering of residual methylation in the genome. Our findings indicate random loss rather than specific maintenance of methylation in Dnmt[1(kd),3a-/-,3b-/-] cells. Near-complete bisulfite conversion and largely unbiased representation of RRBS libraries suggest that random shotgun bisulfite sequencing can be scaled to a genome-wide approach.

Figure 1

Reduced representation bisulfite sequencing. Genomic DNA is digested to completion using a restriction enzyme (here BglII). After size-selection an adapter is added. The DNA is denatured, and unmethylated cytosines are bisulfite-converted to uracil. The two resulting C-poor strands are no longer complementary to each other. Primers specific for the converted adapter sequence are used to fill-in the second (G-poor) strand and for PCR amplification. PCR products are cloned and sequenced. Sequences generated from RRBS libraries are projected onto the genome by searching against a reduced representation database of BglII fragments that had been size-selected and bisulfite-converted in silico.

Figure 2

Generation of Dnmt1, Dnmt3a and Dnmt3b deficient ES cells. (A) Dnmt3a/3b homozygous double knockout ES cells have been described earlier (43). The knockdown virus is expressing a Dnmt1 shRNA, whereas the control is not. The infection was termed Passage 0. After the infection ES cells were passaged four times on feeders followed by two additional passages under feeder-free conditions (Passage 6). Number of viral integrations were determined by Southern blotting and clones with single integration were selected (data not shown). (B) Western blot analysis. The status of the different Dnmts is indicated above. The knockdown ES cells showed a significant reduction in Dnmt1 levels compared with their sister clone. c/c is a previously reported Dnmt1 null ES line (16).

Figure 3

Methylation status of the Dnmt-deficient ES cells. All knockdown and control ES cells were analyzed at Passage 6 after infection. (A) Minor satellite repeat methylation. HpaII digests of genomic DNA were hybridized to minor satellite probe pMR150. MspI is an isoschizomere of HpaII and cuts irrespective of the methylation status (i.e. appearance of a ladder in HpaII lane indicates loss of methylation). The status of the different Dnmts is shown above the Southern blot. All knockdown and control ES cell lines were generated as described in Figure 2. Each knockdown line contains a single lentiviral integration (data not shown). (B) IAP methylation. HpaII-digested genomic DNA was hybridized to an IAP probe. (C) COBRA analysis for imprinted genes. Genomic DNA was bisulfite treated and after PCR amplification of H19, Snrpn, Peg1 and Peg3 a restriction digest was performed to analyze the methylation status of the differentially methylated regions (U = unmethylated, M = methylated). The second (smaller) fragment of the methylated and digest product is not shown. (D) Total mCpG quantification by NNA. The spots corresponding to CpG and mCpG are indicated in the upper left panel. The per cent mCpG/(CpG+mCpG) are displayed in each panel (estimated error 5%).

Figure 4

Size distributions of the sequenced clones from each library. RRBS reads from wild-type ES cells (black) had a mean of 553 bp and an SD of 17 bp. Dnmt[1^kd,3a^−/−,3b^−/−] reads were (570 ± 20) bp in size (grey bars). The size distributions of the two libraries were overlapping but not identical.

Figure 5

Targeted bisulfite sequencing of specific loci. Ten loci for which RRBS sequencing indicated mCpGs in Dnmt-deficient cells and 10 loci that were devoid of methylation were bisulfite re-sequenced using specific primers in wild-type (top), 3a/b double knockout (middle) and Dnmt[1^kd,3a^−/−,3b^−/−] cells (bottom). Shown are two examples of each set. Each row represents a single sequenced molecule. Filled squares are methylated CpGs and empty ones indicate unmethylated sites. The asterisk indicates the original clone sequenced from the library.