what is reduced representation bisulfite sequencing

• Technology Overview Next-Generation Sequencing Digital Microfluidics Microarray Technology Technology Development
• Sequencing Method Explorer
• Publication Reviews
• Customer Stories Overview iCommunity Customer Interviews & Case Studies
• Genomics Education Overview Illumina Resources & Tools Genetic Analysis FFPE Sample Analysis Sequencing Workflow Accuracy Adventures in Genomics Genomics 101
• Genomics Podcast
• Genomics Consortia Overview Human Consortia
• Genomics Research Hub Overview Articles Publications Open-Source Bioinformatics Tools

RRBS-Seq/scRRBS

Reduced-representation bisulfite sequencing (RRBS-Seq ) or single-cell reduced-representation bisulfite sequencing (scRRBS ) is a protocol that uses one or multiple restriction enzymes on the genomic DNA to produce sequence-specific fragmentation. The fragmented genomic DNA is then treated with bisulfite and sequenced. This is the method of choice to study specific regions of interest. It is particularly effective where methylation is high, such as in promoters and repeat regions. 

• Genome-wide coverage of CpGs in islands at single-base resolution
• Areas dense in CpG methylation are covered
• Restriction enzymes cut at specific sites, providing biased sequence selection
• Method measures 10-15% of all CpGs in genome
• Cannot distinguish between 5mC and 5hmC
• Does not cover non-CpG areas, genome-wide CpGs, and CpGs in areas without the enzyme restriction site

Methods Links:

• RRBS: Meissner A., Gnirke A., Bell G. W., Ramsahoye B., Lander E. S., et al. (2005) Reduced representation bisulfite sequencing for comparative high-resolution DNA methylation analysis. Nucleic Acids Res 33: 5868-5877
• scRRBS: Guo H., Zhu P., Guo F., Li X., Wu X., et al. (2015) Profiling DNA methylome landscapes of mammalian cells with single-cell reduced-representation bisulfite sequencing. Nat Protoc 10: 645-659
• Infinium MethylationEPIC Kit
• TruSeq DNA PCR-Free Library Prep Kit
• Guo H., Zhu P., Guo F., Li X., Wu X., et al. (2015) Profiling DNA methylome landscapes of mammalian cells with single-cell reduced-representation bisulfite sequencing. Nat Protoc 10: 645-659
• Hannon E., Chand A. N., Evans M. D., Wong C. C. Y., Grubb M. S., et al. A role for CaV1 and calcineurin signalling to depolarization-induced changes in neuronal DNA methylation. Neuroepigenetics  
• Lin X., Stenvang J., Rasmussen M. H., Zhu S., Jensen N. F., et al. (2015) The potential role of Alu Y in the development of resistance to SN38 (Irinotecan) or oxaliplatin in colorectal cancer. BMC Genomics 16: 404
• Lee E.-J., Rath P., Liu J., Ryu D., Pei L., et al. Identification of Global DNA Methylation Signatures in Glioblastoma-Derived Cancer Stem Cells. Journal of Genetics and Genomics  
• Yang Y., Sebra R., Pullman B. S., Qiao W., Peter I., et al. (2015) Quantitative and multiplexed DNA methylation analysis using long-read single-molecule real-time bisulfite sequencing (SMRT-BS). BMC Genomics 16: 350
• Zhang G., Esteve P. O., Chin H. G., Terragni J., Dai N., et al. (2015) Small RNA-mediated DNA (cytosine-5) methyltransferase 1 inhibition leads to aberrant DNA methylation. Nucleic Acids Res 43: 6112-6124

• Complete Plasmid DNA Sequencing
• Mitochondrial DNA (mtDNA) Sequencing
• Viral Genome Sequencing
• Whole Genome Sequencing
• Whole Exome Sequencing
• Targeted Region Sequencing
• Amplicon Sequencing Services
• Chloroplast DNA (cpDNA) Sequencing
• Long Amplicon Analysis (LAA)
• Shallow Whole Genome Sequencing
• Human Mitochondrial DNA (mtDNA) Sequencing
• Long-Read Sequencing
• ctDNA Sequencing Services
• Gene Panel Sequencing Service
• Bacterial RNA Sequencing
• Small RNA Sequencing
• LncRNA Sequencing
• CircRNA Sequencing
• Ribosome Profiling (Ribo-seq)
• Total RNA Sequencing
• Targeted RNA Sequencing
• Degradome Sequencing
• Exosomal RNA Sequencing
• Ultra Low RNA Sequencing
• Digital RNA Sequencing (UMI-RNA-Seq)
• Dual RNA-seq
• microRNA Sequencing Service
• mRNA Sequencing Service
• EpiTYPER DNA Methylation Analysis
• Targeted Bisulfite Sequencing
• Reduced Representation Bisulfite Sequencing
• Whole Genome Bisulfite Sequencing (WGBS)
• MeDIP Sequencing (MeDIP-Seq)
• MethylRAD-Seq
• MeRIP Sequencing (m6A Analysis)
• DNA 6mA Sequencing
• DAP-Seq Service (DNA affinity purification sequencing)
• 5mC/5hmC Sequencing
• Nanopore RNA Methylation Sequencing Service
• RNA Methylation Sequencing Service
• 2'-O-RNA Methylation Sequencing Service
• Long-Read Metagenomic Sequencing
• Bacterial Whole Genome de novo Sequencing
• Fungal Whole Genome de novo Sequencing
• Full-Length Transcripts Sequencing (Iso-Seq)
• Human Whole Genome PacBio SMRT Sequencing
• Full-Length 16S/18S/ITS Amplicon Sequencing
• 16S/18S/ITS Amplicon Sequencing
• Metagenomic Shotgun Sequencing
• Viral Metagenomic Sequencing
• Metatranscriptomic Sequencing
• Microbial Whole Genome Sequencing
• Absolute Quantitative 16s/18s/ITS Amplicon Sequencing
• Microbial Identification
• Single-cell RNA Sequencing
• Single-cell DNA Methylation Sequencing
• Single-Cell DNA Sequencing
• 10x Genomics Single-Cell Sequencing
• Nanopore Full-Length Transcripts Sequencing
• Nanopore Direct RNA Sequencing
• CRISPR Sequencing
• CRISPR Screen Sequencing
• Antibody Screening Sequencing (Phage Display Library Screening)
• Immune Repertoire Sequencing
• Sanger Sequencing
• Pre-made Library Sequencing
• Lentiviral/Retroviral Integration Sites Analysis
• AAV Genome Sequencing
• Genotyping by Sequencing (GBS)
• SNP Microarray
• PCR-LDR SNP Genotyping Service
• RAD-seq Service
• MassARRAY SNP Genotyping
• SNaPshot Multiplex System for SNP Genotyping
• TaqMan SNP Genotyping
• CGH Microarray Service
• Microsatellite Genotyping Service
• Microsatellite Instability Analysis
• Microsatellite Development
• Genome-wide association study (GWAS)
• Variant Calling
• Population Evolution
• Genetic Linkage Map
• Bulk Segregant Analysis (BSA)
• Genomic Data Analysis
• Transcriptomic Data Analysis
• Epigenomics Data Analysis Service
• Long-Read Sequencing Data Analysis Service
• Single-cell RNA Sequencing Data Analysis Service
• MicroRNA Expression Profiling Microarray Service
• Gene Expression Profiling Microarray Service
• LncRNA Microarray Service
• Methylation Screening Array 270K Service
• Agriculture and Food Science
• Drug Development
• Inherited Disease Research
• Microbial Research
• Oncology Research
• Sample Submission Guidelines
• Order Online
• Genome Research

An Introduction to Reduced Representation Bisulfite Sequencing (RRBS)

What is DNA methylation?

DNA methylation refers to the addition of a methyl group (CH3) to the DNA strand which is catalyzed by DNA methyltransferases (DNMTs). DNA contains four kinds of nitrogenous bases, namely cytosine, guanine, thymine and adenine (Figure 1A). Researchers have found that cytosine and adenine can be methylated. Cytosine methylation often happens at the 5-carbon position of cytosine (5-methylcytosine or 5mC), which is found exclusively at symmetric CG sites on the DNA double helix across the entire genome, namely the CpG island (Figure 2B). Cytosine methylation is widespread in both eukaryotes and prokaryotes (Cloney, 2016; Cooper, 1983). Adenine methylation happens at the 6-nitrogen position of adenine (N6-methyladenine or N6mA). Adenine methylation has been observed in bacterial, plant and mammalian DNA (Ratel et al ., 2006; Wu et al ., 2016), but has received considerably less attention (Figure 2B).

DNA methylation is one of the most important epigenetic modifications that is correlated with gene repression (Figure 2) and is known to play an essential role in embryonic development, maintenance of pluripotency, genomic imprinting and X-chromosome inactivation through regulation of transcription, chromatin structure and chromosome stability (Robertson, 2005). DNA methylation can also affect health resulting in cancer, autoimmune disease, neurological disorders or other diseases. In many disease processes, gene promoter CpG islands acquire abnormal hypermethylation, which results in transcriptional silencing that can be inherited by daughter cells following cell division (Wang and Lei, 2018). Alterations of DNA methylation patterns have been recognized as an important component of cancer development (Figure 3).

DNA Methylation Detection and RRBS

Understanding the role of DNA methylation in development and disease requires knowledge of the distribution of these modifications in the genome (Harris et al ., 2010). Measuring the total amount of 5mC or 5hmC (5-hydroxymethylcytosine) allows researchers to gain insight into profound biological processes and identify biomarkers for disease. The detection of DNA methylation patterns is a rapidly advancing area of research and several methods have been available for the assessment of DNA methylation, with bisulfite treatment being a central procedure to a majority. A simple algorithm for choosing a suitable method for DNA methylation analysis is illustrated in Figure 4 (Kurdyukov and Bullock, 2016).

Bisulfite sequencing is the use of bisulfite treatment of DNA to determine its pattern of methylation (Frommer et al ., 1992). The bisulfite treatment of DNA mediates the deamination of cytosine into uracil, and these converted residues will be read as thymine, as determined by PCR-amplification and subsequent sequencing analysis (Figure 5).

Whole genome bisulfite sequencing (WGBS)  is the most comprehensive method of DNA methylation analysis, and now has been considered as the “gold standard” method in DNA methylation studies. The only limitations are the cost and difficulties in the analysis of NGS data.  Reduced representation bisulfite sequencing (RRBS) , which combines restriction enzymes and bisulfite sequencing to enrich for areas of the genome with a high CpG content, is an efficient and high-throughput technique for analyzing the genome-wide methylation profiles on a single nucleotide level. The method can reduce the number of nucleotides required to sequence to 1% of the genome (Meissner et al ., 2005). Due to its economical and productive manner,  RRBS has been used in quickly profiling of aberrant methylation in cancer (Smith et al ., 2009) and methylation states in development.

Workflow of RRBS

The protocol of RRBS is illustrated in Figure 6. The process consists of several steps:

• Enzyme digestion Genomic DNA is digested into DNA fragments of various sizes using a methylation-insensitive restriction enzyme. MspI is commonly used.
• Library preparation After enzyme digestion, library preparation, a process consisting of end repair, A-tailing and sequencing adapter ligation, is performed. The reaction of MspI digestion results in strands with sticky ends. End repair is necessary to fill in the 3’ terminal of the ends of the strands. A-tailing, in which an extra adenosine is added to both the plus and minus strands, is necessary for methylated sequence adapter ligation in the subsequent step.
• Size selection Afterward, a size selection of the resulting fragments is performed. Different sizes of the fragments are separated using gel electrophoresis and are purified using gel excising.
• Bisulfite conversion After the size selection, bisulfite conversion is executed, in which unmethylated cytosine is deaminated into uracil.
• PCR amplification and purification The bisulfite converted DNA is then amplified using PCR with primers that are complementary to the sequence adapters. A step for PCR purification is required before sequencing.
• Sequencing The sequencing is then performed using next generation sequencing technology. The Illumina platforms are most commonly performed.

Advantages and Limitations of RRBS

RRBS has some technical advantages over WGBS and other DNA methylation detection methods. But there are also some limitations in the specific protocol steps.

• Advantages:
• Taking advantage of using specific restriction enzyme, RRBS can enrich CpG regions of the genome, and decrease the amount of sequencing required as well as decrease the cost.
• Only a low sample concentration is required for accurate data analysis.
• Limitations
• The restriction enzyme digestion could not cover all the CG regions in the genome. Some CpG’s are missed and some regions thus have lower coverage.
• PCR amplification errors could occur.
• Complete bisulfite conversion of double stranded DNA (dsDNA) requires a denaturation step since bisulfite sequencing only converts single-stranded DNA (ssDNA).

At CD Genomics, we are dedicated to providing reliable epigenomics sequencing  services, including EpiTYPER DNA methylation analysis ,  targeted bisulfite sequencing ,  reduced representation bisulfite sequencing (RRBS) ,  whole genome bisulfite sequencing ,  MeDIP sequencing ,  ChIP-seq , and MethylRAD-seq .

References:

• Cloney, R. (2016). DNA methylation: A SMRT analysis of prokaryotic epigenomes. Nature reviews Genetics  17 , 195.
• Cooper, D.N. (1983). Eukaryotic DNA methylation. Human genetics  64 , 315-333.
• Frommer, M., McDonald, L.E., Millar, D.S., Collis, C.M., Watt, F., Grigg, G.W., Molloy, P.L., and Paul, C.L. (1992). A genomic sequencing protocol that yields a positive display of 5-methylcytosine residues in individual DNA strands. Proceedings of the National Academy of Sciences of the United States of America  89 , 1827-1831.
• Harris, R.A., Wang, T., Coarfa, C., Nagarajan, R.P., Hong, C., Downey, S.L., Johnson, B.E., Fouse, S.D., Delaney, A., Zhao, Y. , et al.  (2010). Comparison of sequencing-based methods to profile DNA methylation and identification of monoallelic epigenetic modifications. Nature biotechnology  28 , 1097-1105.
• Kurdyukov, S., and Bullock, M. (2016). DNA Methylation Analysis: Choosing the Right Method. Biology  5 .
• Meissner, A., Gnirke, A., Bell, G.W., Ramsahoye, B., Lander, E.S., and Jaenisch, R. (2005). Reduced representation bisulfite sequencing for comparative high-resolution DNA methylation analysis. Nucleic acids research  33 , 5868-5877.
• Mukherjee, K., Twyman, R.M., and Vilcinskas, A. (2015). Insects as models to study the epigenetic basis of disease. Progress in biophysics and molecular biology  118 , 69-78.
• Ratel, D., Ravanat, J.L., Berger, F., and Wion, D. (2006). N6-methyladenine: the other methylated base of DNA. BioEssays : news and reviews in molecular, cellular and developmental biology  28 , 309-315.
• Robertson, K.D. (2005). DNA methylation and human disease. Nature reviews Genetics  6 , 597-610.
• Smith, Z.D., Gu, H., Bock, C., Gnirke, A., and Meissner, A. (2009). High-throughput bisulfite sequencing in mammalian genomes. Methods  48 , 226-232.
• Wang, Y.P., and Lei, Q.Y. (2018). Metabolic recoding of epigenetics in cancer. Cancer communications  38 , 25.
• Wu, T.P., Wang, T., Seetin, M.G., Lai, Y., Zhu, S., Lin, K., Liu, Y., Byrum, S.D., Mackintosh, S.G., Zhong, M. , et al.  (2016). DNA methylation on N(6)-adenine in mammalian embryonic stem cells. Nature  532 , 329-333.

* is a required item.

• Next Generation Sequencing
• Bioinformatics Services
• Microarray Services

Reduced Representation Bisulfite Sequencing (RRBS)

Series: Methods In Molecular Biology > Book: Epigenomics

Protocol | DOI: 10.1007/978-1-0716-2724-2_3

• Department of Obstetrics and Gynecology, The Jikei University School of Medicine, Minato, Tokyo, Japan
• Department of Maternal-Fetal Biology, Research Institute, National Center for Child Health and Development, Setagaya, Tokyo, Japan
• Department of Obstetrics and Gynecology, Keio University School of Medicine, Shinjuku, Tokyo, Japan

Reduced representation bisulfite sequencing (RRBS) enriches CpG-rich genomic regions using the MspI restriction enzyme—which cuts DNA at all CCGG sites, regardless of their DNA methylation status at the CG site—and enables the measurement of DNA

Reduced representation bisulfite sequencing (RRBS) enriches CpG-rich genomic regions using the MspI restriction enzyme—which cuts DNA at all CCGG sites, regardless of their DNA methylation status at the CG site—and enables the measurement of DNA methylation levels at 5% ~ 10% of all CpG sites in the mammalian genome. RRBS has been utilized in a large number of studies as a cost-effective method to investigate DNA methylation patterns, mainly at gene promoters and CpG islands. Here, we describe protocols for gel-free preparation of RRBS libraries, quality control, sequencing, and data analysis. Our protocols typically require nine cycles of polymerase chain reaction (PCR) amplification to obtain a sufficient amount of library for sequencing when 100 ng of genomic DNA is used as a starting material; moreover, it takes 3 days to complete library preparation and quality control procedures for up to eight samples.

Figures ( 0 ) & Videos ( 0 )

Experimental specifications, other keywords.

Related articles

Studying dna methylation genome-wide by bisulfite sequencing from low amounts of dna in mammals, bisulfite sequencing using small dna amounts, epigenetic approaches in non-model plants, low-input whole-genome bisulfite sequencing, plant-rrbs: dna methylome profiling adjusted to plant genomes, utilizing efficient endonuclease combinations, for multi-sample studies, epigenomics: sequencing the methylome, reduced representation bisulfite sequencing to identify global alteration of dna methylation, from methylome to integrative analysis of tissue specificity, generating sequencing-based dna methylation maps from low dna input samples, next-generation bisulfite sequencing for targeted dna methylation analysis.

• Meissner A, Gnirke A, Bell GW et al (2005) Reduced representation bisulfite sequencing for comparative high-resolution DNA methylation analysis. Nucleic Acids Res 33:5868–5877
• Gu H, Smith ZD, Bock C et al (2011) Preparation of reduced representation bisulfite sequencing libraries for genome-scale DNA methylation profiling. Nat Protoc 6:468–481
• Boyle P, Clement K, Gu H et al (2012) Gel-free multiplexed reduced representation bisulfite sequencing for large-scale DNA methylation profiling. Genome Biol 13:R92
• Garrett-Bakelman FE, Sheridan CK, Kacmarczyk TJ et al (2015) Enhanced reduced representation bisulfite sequencing for assessment of DNA methylation at base-pair resolution. J Vis Exp 96:e52246
• Guo H, Zhu P, Guo F et al (2015) Profiling DNA methylome landscapes of mammalian cells with single-cell reduced-representation bisulfite sequencing. Nat Protoc 10:645–659
• Gu H, Raman AT, Wang X et al (2021) Smart-RRBS for single-cell methylome and transcriptome analysis. Nat Protoc 16:4004–4030
• Miura F, Enomoto Y, Dairiki R, Ito T (2012) Amplification-free whole-genome bisulfite sequencing by post-bisulfite adaptor tagging. Nucleic Acids Res 40:e136
• Katsuoka F, Yokozawa J, Tsuda K et al (2014) An efficient quantitation method of next-generation sequencing libraries by using MiSeq sequencer. Anal Biochem 466:27–29
• Toh H, Shirane K, Miura F et al (2017) Software updates in the Illumina HiSeq platform affect whole-genome bisulfite sequencing. BMC Genomics 18:31
• https://www.bioinformatics.babraham.ac.uk/projects/trim_galore/
• Krueger F, Andrews SR (2011) Bismark: a flexible aligner and methylation caller for Bisulfite-Seq applications. Bioinformatics 27:1571–1572
• https://rawgit.com/FelixKrueger/Bismark/master/Docs/Bismark_User_Guide.html
• Akalin A, Kormaksson M, Li S et al (2012) methylKit: a comprehensive R package for the analysis of genome-wide DNA methylation profiles. Genome Biol 13:R87
• https://bioconductor.org/packages/release/bioc/vignettes/methylKit/inst/doc/methylKit.html
• Thorvaldsdóttir H, Robinson JT, Mesirov JP (2013) Integrative Genomics Viewer (IGV): high-performance genomics data visualization and exploration. Brief Bioinform 14:178–192
• Nohara K, Nakabayashi K, Okamura K et al (2020) Gestational arsenic exposure induces site-specific DNA hypomethylation in active retrotransposon subfamilies in offspring sperm in mice. Epigenetics Chromatin 13:53
• MacConaill LE, Burns RT et al (2018) Unique, dual-indexed sequencing adapters with UMIs effectively eliminate index cross-talk and significantly improve sensitivity of massively parallel sequencing. BMC Genomics 19:30

Advertisement

• Sponsored Article

Brush Up: What Is Bisulfite Sequencing and How Do Researchers Use It to Study DNA Methylation?

Prior to dna methylation sequencing, researchers treat their samples with sodium bisulfite to distinguish methylated cytosine from unmethylated cytosine..

Deanna earned their PhD in cellular biology from McGill University in 2020 and has a professional background in medical writing. They joined The Scientist's Creative Services team as an assistant science editor in 2022.

View full profile.

Learn about our editorial policies.

What is DNA Methylation?

DNA methylation is a fundamental epigenetic mechanism in which a methyl group is added onto a nucleotide, commonly cytosine. DNA methylation affects many biological processes, including gene expression, embryonic development, inflammation, and cellular proliferation and differentiation. Several complex diseases have aberrant DNA methylation patterns, such as different types of cancer and neurodegenerative disorders, which are often associated with genomic instability and loss of DNA homeostasis. 1–3

Researchers need effective methods with high sensitivity and reliability to explore the importance of DNA methylation. The gold standard technology that scientists use to detect DNA methylation is bisulfite genomic sequencing. It is a qualitative, quantitative, and efficient approach to identify methylated cytosine at single base-pair resolution. 3

Bisulfite Conversion of Methylated Cytosine

Cytosine methylation occurs when a methyl group binds to the fifth carbon of a cytosine to form 5-methylcytosine (5mC). Researchers cannot detect 5mC with traditional molecular techniques such as PCR or cloning methods because methyl groups are not copied during PCR amplification. Instead, researchers must treat the methylated DNA samples with bisulfite prior to detection, to distinguish between methylated and unmethylated cytosines. 1,2

Bisulfite conversion is the chemical reaction that occurs when scientists treat DNA with sodium bisulfite. In 1992, researchers found that the amination reaction of sodium bisulfite with unmethylated cytosine is different than the reaction of sodium bisulfite with 5mC. Because of this difference, unmethylated cytosines in single-stranded DNA become uracil residues after exposure to sodium bisulfite, while 5mCs remain cytosines. Pretreatment with sodium bisulfite is the basis of many methylation detection and analysis techniques. 1 Following bisulfite conversion, researchers determine the methylation status in loci of interest because unmethylated cytosines are recognized as thymine after PCR amplification and sequencing. 3 

Challenges of Treating DNA with Bisulfite

Incomplete bisulfite treatment is a major challenge in methylation detection methods. Bisulfite conversion relies on the chemical modification of cytosine in single-stranded DNA, so scientists must ensure that their DNA sample is completely denatured prior to sodium bisulfite treatment. As such, the quality and quantity of purified DNA and the pH of the reaction are fundamental sample parameters in complete bisulfite conversion. 1 Incomplete conversion typically leads researchers to overestimate the amount of methylation in a sample. 2 An additional challenge to bisulfite-based methylation analyses is the sensitivity of the DNA sample to degradation during long incubation steps. 1

Types of Bisulfite Sequencing

Whole genome bisulfite next-generation sequencing.

Methylation detection methods are diverse in accuracy, sensitivity, speed, simplicity, and cost. Researchers investigate the methylation of an entire genome, also referred to as a methylome, with whole genome bisulfite sequencing (WGBS). WGBS consists of library preparation with bisulfite treatment, next-generation sequencing (NGS), and high throughput analysis in which researchers explore the biological processes connected to the observed methylome. 1,2,5   Challenges during WGBS arise from protocol changes that cause discrepancies in coverage depth, read quality, duplication rates, mapping efficiency, and methylation estimation. Additionally, the standard WGBS protocol requires a large amount of DNA, which may be a roadblock in many methylome studies. 2,5

Reduced Representation Bisulfite Sequencing

An alternative approach to WGBS is reduced representation bisulfite sequencing (RRBS). This method reduces the cost and complexity of methylome analysis. With RRBS, scientists enrich for GC-rich parts of the genome by digesting the DNA with an enzyme that generates fragments with CpG dinucleotides at both ends, independent of methylation status. This enrichment is informative because DNA methylation exists primarily in the context of symmetric CpG dinucleotides. After digestion, researchers use size selection to isolate short DNA fragments that correspond to GC-rich regions. Larger fragments that are not GC-rich and contain fewer potential methylation sites are discarded. This digestion and selection step essentially excludes uninformative DNA prior to bisulfite treatment, while capturing the majority of relevant genomic regions of the methylome. RRBS reduces the number of reads needed to obtain high coverage and requires much less DNA than WGBS, but researchers can only examine a fraction of the methylome with this technique due to the amount of excluded DNA. 4–7

Targeted Bisulfite Sequencing

In contrast to the broad enrichment of GC-rich regions with RRBS, targeted bisulfite sequencing allows researchers to investigate the methylation status of each cytosine in a specific genomic region. Scientists apply targeted methods to validate differentially-methylated regions or analyze candidate regions rather than the entire methylome. Targeted bisulfite sequencing relies on two PCR amplification steps after bisulfite conversion, followed by deep sequencing. The PCR steps increase copy number and introduce adaptors required for sequencing as well as sample-specific identifiers. These adaptors and identifiers allow researchers to process a large number of samples simultaneously on a high-throughput instrument. 7,8

• R. Halabian et al., “Laboratory methods to decipher epigenetic signatures: a comparative review,” Cell Mol Biol Lett, 26:1-30, 2021.
• T. Gong et al., “Analysis and performance assessment of the whole genome bisulfite sequencing data workflow: currently available tools and a practical guide to advance DNA methylation studies,” Small Methods , 6:e2101251, 2022.
• Y. Li, T.O. Tollefsbol, “DNA methylation detection: bisulfite genomic sequencing analysis,” Methods Mol Biol , 791:11-21, 2011.
• Y. Zhang, A. Jeltsch, “The application of next generation sequencing in DNA methylation analysis,” Genes (Basel) , 1:85-101, 2010.
• Q. Wang et al., “Tagmentation-based whole-genome bisulfite sequencing,” Nat Protoc , 8:2022-32, 2013.
• H. Gu et al., “Preparation of reduced representation bisulfite sequencing libraries for genome-scale DNA methylation profiling,” Nat Protoc , 6:468-81, 2011.
• D.A. Moser et al., “Targeted bisulfite sequencing: A novel tool for the assessment of DNA methylation with high sensitivity and increased coverage,” Psychoneuroendocrinology , 120:1-8, 2020.
• E. Leitão et al., “Locus-specific DNA methylation analysis by targeted deep bisulfite sequencing,” Methods Mol Biol , 1767:351-66, 2018.

Related genetics Research Resources

The Six-Base Genome Reveals Multimodal Data from a Single DNA Sample

Tapping into Immunotherapy’s Potential to Help More People

Epigenetics in a Dish

Cancer Genomics and Proteomics pp 23–31 Cite as

Reduced Representation Bisulfite Sequencing to Identify Global Alteration of DNA Methylation

• Arvindhan Nagarajan Ph.D. 3 ,
• Christine Roden B.S. 3 &
• Narendra Wajapeyee Ph.D. 3  
• First Online: 01 January 2014

2955 Accesses

7 Citations

1 Altmetric

Part of the book series: Methods in Molecular Biology ((MIMB,volume 1176))

Reduced representation bisulfite sequencing is a cost-effective high-throughput sequencing-based method to obtain DNA methylation status at a single-nucleotide level. DNA methylation status is determined by utilizing DNA methylation-specific restriction enzymes to selectively amplify for genomic regions that are rich in methylated DNA. Although the method is genome-wide, DNA methyl sequencing does not require the sequencing of the whole genome, hence the name “reduced representation.” However, a large majority of CpG islands are covered by reduced representation bisulfite sequencing allowing for the acquisition of comprehensive information of the methylation landscape in diseases like cancer. Data generated by this approach is typically reproducible and often covers between 65 and 75 % of the whole genome.

This is a preview of subscription content, log in via an institution .

Buying options

• Available as PDF
• Read on any device
• Instant download
• Own it forever
• Available as EPUB and PDF
• Compact, lightweight edition
• Dispatched in 3 to 5 business days
• Free shipping worldwide - see info
• Durable hardcover edition

Tax calculation will be finalised at checkout

Purchases are for personal use only

Springer Nature is developing a new tool to find and evaluate Protocols.  Learn more

Deaton AM, Bird A (2011) CpG islands and the regulation of transcription. Genes Dev 25(10):1010–1022, PubMed PMID: 21576262. Pubmed Central PMCID: 3093116. Epub 2011/05/18. eng

Article   CAS   PubMed Central   PubMed   Google Scholar  

Jones PA, Baylin SB (2007) The epigenomics of cancer. Cell 128(4):683–692, PubMed PMID: 17320506. Epub 2007/02/27. eng

Wu SC, Zhang Y (2010) Active DNA demethylation: many roads lead to Rome. Nat Rev Mol Cell Biol 11(9):607–620, PubMed PMID: 20683471. Pubmed Central PMCID: 3711520. Epub 2010/08/05. eng

Gnyszka A, Jastrzebski Z, Flis S (2013) DNA methyltransferase inhibitors and their emerging role in epigenetic therapy of cancer. Anticancer Res 33(8):2989–2996, PubMed PMID: 23898051. Epub 2013/07/31. eng

CAS   PubMed   Google Scholar  

Bock C, Tomazou EM, Brinkman AB et al (2010) Quantitative comparison of genome-wide DNA methylation mapping technologies. Nat Biotechnol 28(10):1106–1114, PubMed PMID: 20852634. Pubmed Central PMCID: 3066564. Epub 2010/09/21. eng

Boyle P, Clement K, Gu H et al (2012) Gel-free multiplexed reduced representation bisulfite sequencing for large-scale DNA methylation profiling. Genome Biol 13(10):R92, PubMed PMID: 23034176. Pubmed Central PMCID: 3491420. Epub 2012/10/05. Eng

Meissner A, Gnirke A, Bell GW et al (2005) Reduced representation bisulfite sequencing for comparative high-resolution DNA methylation analysis. Nucleic Acids Res 33(18):5868–5877, PubMed PMID: 16224102. Pubmed Central PMCID: 1258174. Epub 2005/10/15. eng

Gu H, Bock C, Mikkelsen TS et al (2010) Genome-scale DNA methylation mapping of clinical samples at single-nucleotide resolution. Nat Methods 7(2):133–136, PubMed PMID: 20062050. Pubmed Central PMCID: 2860480. Epub 2010/01/12. eng

Gu H, Smith ZD, Bock C et al (2011) Preparation of reduced representation bisulfite sequencing libraries for genome-scale DNA methylation profiling. Nat Protoc 6(4):468–481, PubMed PMID: 21412275. Epub 2011/03/18. eng

Article   CAS   PubMed   Google Scholar  

Download references

Acknowledgements

N.W. is a translational scholar of the Sidney Kimmel Foundation for Cancer Research and is supported by the Young Investigator Awards from National Lung Cancer Partnership, Uniting Against Lung Cancer, and International Association for the Study of Lung Cancer, Melanoma Research Alliance, and Melanoma Research Foundation.

Author information

Authors and affiliations.

Department of Pathology, Yale University School of Medicine, LH-214A, 310 Cedar Street, New Haven, CT, 06510, USA

Arvindhan Nagarajan Ph.D., Christine Roden B.S. & Narendra Wajapeyee Ph.D.

You can also search for this author in PubMed   Google Scholar

Corresponding author

Correspondence to Narendra Wajapeyee Ph.D. .

Editor information

Editors and affiliations.

Department of Pathology& Yale Cancer Center, Yale University School of Medicine, New Haven, Connecticut, USA

Narendra Wajapeyee

Rights and permissions

Reprints and permissions

Copyright information

© 2014 Springer Science+Business Media New York

About this protocol

Cite this protocol.

Nagarajan, A., Roden, C., Wajapeyee, N. (2014). Reduced Representation Bisulfite Sequencing to Identify Global Alteration of DNA Methylation. In: Wajapeyee, N. (eds) Cancer Genomics and Proteomics. Methods in Molecular Biology, vol 1176. Humana Press, New York, NY. https://doi.org/10.1007/978-1-4939-0992-6_3

Download citation

DOI : https://doi.org/10.1007/978-1-4939-0992-6_3

Published : 20 June 2014

Publisher Name : Humana Press, New York, NY

Print ISBN : 978-1-4939-0991-9

Online ISBN : 978-1-4939-0992-6

eBook Packages : Springer Protocols

• Publish with us

Policies and ethics

• Find a journal
• Track your research

Thank you for visiting nature.com. You are using a browser version with limited support for CSS. To obtain the best experience, we recommend you use a more up to date browser (or turn off compatibility mode in Internet Explorer). In the meantime, to ensure continued support, we are displaying the site without styles and JavaScript.

• View all journals
• Explore content
• About the journal
• Publish with us
• Sign up for alerts
• Published: 06 May 2021

Extended-representation bisulfite sequencing of gene regulatory elements in multiplexed samples and single cells

• Sarah J. Shareef 1 , 2 ,
• Samantha M. Bevill 1 , 2 ,
• Ayush T. Raman 1 , 2 ,
• Martin J. Aryee   ORCID: orcid.org/0000-0002-6848-1344 1 , 2 , 3 ,
• Peter van Galen   ORCID: orcid.org/0000-0002-0735-1570 1 , 2 , 4 ,
• Volker Hovestadt   ORCID: orcid.org/0000-0002-3480-6649 1 , 2 , 5 , 6 &
• Bradley E. Bernstein   ORCID: orcid.org/0000-0002-5726-6278 1 , 2  

Nature Biotechnology volume  39 ,  pages 1086–1094 ( 2021 ) Cite this article

9242 Accesses

24 Citations

55 Altmetric

Metrics details

• DNA methylation
• Epigenomics
• Methylation analysis
• Next-generation sequencing

The biological roles of DNA methylation have been elucidated by profiling methods based on whole-genome or reduced-representation bisulfite sequencing, but these approaches do not efficiently survey the vast numbers of non-coding regulatory elements in mammalian genomes. Here we present an extended-representation bisulfite sequencing (XRBS) method for targeted profiling of DNA methylation. Our design strikes a balance between expanding coverage of regulatory elements and reproducibly enriching informative CpG dinucleotides in promoters, enhancers and CTCF binding sites. Barcoded DNA fragments are pooled before bisulfite conversion, allowing multiplex processing and technical consistency in low-input samples. Application of XRBS to single leukemia cells enabled us to evaluate genetic copy number variations and methylation variability across individual cells. Our analysis highlights heterochromatic H3K9me3 regions as having the highest cell-to-cell variability in their methylation, likely reflecting inherent epigenetic instability of these late-replicating regions, compounded by differences in cell cycle stages among sampled cells.

This is a preview of subscription content, access via your institution

Access options

Access Nature and 54 other Nature Portfolio journals

Get Nature+, our best-value online-access subscription

24,99 € / 30 days

cancel any time

Subscribe to this journal

Receive 12 print issues and online access

195,33 € per year

only 16,28 € per issue

Rent or buy this article

Prices vary by article type

Prices may be subject to local taxes which are calculated during checkout

Similar content being viewed by others

Smart-RRBS for single-cell methylome and transcriptome analysis

Hongcang Gu, Ayush T. Raman, … Andreas Gnirke

scTEM-seq: Single-cell analysis of transposable element methylation to link global epigenetic heterogeneity with transcriptional programs

Kooper V. Hunt, Sean M. Burnard, … Heather J. Lee

De novo detection of somatic mutations in high-throughput single-cell profiling data sets

Francesc Muyas, Carolin M. Sauer, … Isidro Cortés-Ciriano

Data availability

The datasets generated during this study have been deposited in the Gene Expression Omnibus with primary accession code GSE149954 .

Luo, C., Hajkova, P. & Ecker, J. R. Dynamic DNA methylation: in the right place at the right time. Science 361 , 1336–1340 (2018).

Article   CAS   PubMed   PubMed Central   Google Scholar  

Yoder, J. A., Walsh, C. P. & Bestor, T. H. Cytosine methylation and the ecology of intragenomic parasites. Trends Genet. 13 , 335–340 (1997).

Article   CAS   PubMed   Google Scholar  

Greenberg, M. V. C. & Bourc’his, D. The diverse roles of DNA methylation in mammalian development and disease. Nat. Rev. Mol. Cell Biol. 20 , 590–607 (2019).

Deaton, A. M. & Bird, A. CpG islands and the regulation of transcription. Genes Dev. 25 , 1010–1022 (2011).

Schübeler, D. Function and information content of DNA methylation. Nature 517 , 321–326 (2015).

Article   PubMed   CAS   Google Scholar  

Baylin, S. B. & Jones, P. A. A decade of exploring the cancer epigenome - biological and translational implications. Nat. Rev. Cancer 11 , 726–734 (2011).

Zemach, A., McDaniel, I. E., Silva, P. & Zilberman, D. Genome-wide evolutionary analysis of eukaryotic DNA methylation. Science 328 , 916–919 (2010).

Yin, Y. et al. Impact of cytosine methylation on DNA binding specificities of human transcription factors. Science 356 , eaaj2239 (2017).

Article   PubMed   PubMed Central   CAS   Google Scholar  

Song, Y. et al. Dynamic enhancer DNA methylation as basis for transcriptional and cellular heterogeneity of ESCs. Mol. Cell 75 , 905–920 (2019).

Flavahan, W. A. et al. Insulator dysfunction and oncogene activation in IDH mutant gliomas. Nature 529 , 110–114 (2016).

Gu, H. et al. Preparation of reduced representation bisulfite sequencing libraries for genome-scale DNA methylation profiling. Nat. Protoc. 6 , 468–481 (2011).

Meissner, A. et al. Genome-scale DNA methylation maps of pluripotent and differentiated cells. Nature 454 , 766–770 (2008).

Boyle, P. et al. Gel-free multiplexed reduced representation bisulfite sequencing for large-scale DNA methylation profiling. Genome Biol. 13 , 1–10 (2012).

Article   CAS   Google Scholar  

Akalin, A. et al. Base-pair resolution DNA methylation sequencing reveals profoundly divergent epigenetic landscapes in acute myeloid leukemia. PLoS Genet. 8 , e1002781 (2012).

Garrett-Bakelman, F. E. et al. Enhanced reduced representation bisulfite sequencing for assessment of DNA methylation at base pair resolution. J. Vis. Exp . e52246 (2015).

Li, G. et al. Joint profiling of DNA methylation and chromatin architecture in single cells. Nat. Methods 16 , 991–993 (2019).

Tanaka, K. & Okamoto, A. Degradation of DNA by bisulfite treatment. Bioorg. Med. Chem. Lett. 17 , 1912–1915 (2007).

Kint, S., De Spiegelaere, W., De Kesel, J., Vandekerckhove, L. & Van Criekinge, W. Evaluation of bisulfite kits for DNA methylation profiling in terms of DNA fragmentation and DNA recovery using digital PCR. PLoS ONE 13 , e0199091 (2018).

Ben-Hattar, J. & Jiricny, J. Methylation of single CpG dinucleotides within a promoter element of the Herpes simplex virus tk gene reduces its transcription in vivo. Gene 65 , 219–227 (1988).

Rao, S. S. P. et al. A 3D map of the human genome at kilobase resolution reveals principles of chromatin looping. Cell 159 , 1665–1680 (2014).

Liu, X. S. et al. Editing DNA methylation in the mammalian genome. Cell 167 , 233–247 (2016).

Fu, Y., Sinha, M., Peterson, C. L. & Weng, Z. The insulator binding protein CTCF positions 20 nucleosomes around its binding sites across the human genome. PLoS Genet. 4 , e1000138 (2008).

Schuyler, R. P. et al. Distinct trends of DNA methylation patterning in the innate and adaptive immune systems. Cell Rep. 17 , 2101–2111 (2016).

Wiehle, L. et al. DNA (de)methylation in embryonic stem cells controls CTCF-dependent chromatin boundaries. Genome Res. 29 , 750–761 (2019).

Li, Y. et al. Alterations of specific chromatin conformation affect ATRA-induced leukemia cell differentiation. Cell Death Dis. 9 , 200 (2018).

Varley, K. E. et al. Dynamic DNA methylation across diverse human cell lines and tissues. Genome Res. 23 , 555–567 (2013).

Zhou, W. et al. DNA methylation loss in late-replicating domains is linked to mitotic cell division. Nat. Genet. 50 , 591–602 (2018).

Michalowsky, L. A. & Jones, P. A. Differential nuclear protein binding to 5-azacytosine-containing DNA as a potential mechanism for 5-aza-2′-deoxycytidine resistance. Mol. Cell Biol. 7 , 3076–3083 (1987).

CAS   PubMed   PubMed Central   Google Scholar  

Aran, D. & Hellman, A. DNA methylation of transcriptional enhancers and cancer predisposition. Cell 154 , 11–13 (2013).

Aran, D., Sabato, S. & Hellman, A. DNA methylation of distal regulatory sites characterizes dysregulation of cancer genes. Genome Biol. 14 , R21 (2013).

Creyghton, M. P. et al. Histone H3K27ac separates active from poised enhancers and predicts developmental state. Proc. Natl Acad. Sci. USA 107 , 21931–21936 (2010).

Bell, A. C. & Felsenfeld, G. Methylation of a CTCF-dependent boundary controls imprinted expression of the Igf2 gene. Nature 405 , 482–485 (2000).

Hark, A. T. et al. CTCF mediates methylation-sensitive enhancer-blocking activity at the H19/Igf2 locus. Nature 405 , 486–489 (2000).

Moore, J. E. et al. Expanded encyclopaedias of DNA elements in the human and mouse genomes. Nature 583 , 699–710 (2020).

Corces, M. R. et al. Lineage-specific and single-cell chromatin accessibility charts human hematopoiesis and leukemia evolution. Nat. Genet. 48 , 1193–1203 (2016).

Barretina, J. et al. The Cancer Cell Line Encyclopedia enables predictive modelling of anticancer drug sensitivity. Nature 483 , 603–607 (2012).

McGahon, A. J. et al. Downregulation of Bcr-Abl in K562 cells restores susceptibility to apoptosis: characterization of the apoptotic death. Cell Death Differ. 4 , 95–104 (1997).

Charlton, J. et al. Global delay in nascent strand DNA methylation. Nat. Struct. Mol. Biol. 25 , 327–332 (2018).

Hansen, R. S. et al. Sequencing newly replicated DNA reveals widespread plasticity in human replication timing. Proc. Natl Acad. Sci. USA 107 , 139–144 (2010).

Angeloni, A. & Bogdanovic, O. Enhancer DNA methylation: implications for gene regulation. Essays Biochem. 63 , 707–715 (2019).

Mulqueen, R. M. et al. Highly scalable generation of DNA methylation profiles in single cells. Nat. Biotechnol. 36 , 428–431 (2018).

Luo, C. et al. Robust single-cell DNA methylome profiling with snmC-seq2. Nat. Commun. 9 , 3824 (2018).

Luo, C. et al. Single-cell methylomes identify neuronal subtypes and regulatory elements in mammalian cortex. Science 357 , 600–604 (2017).

Hou, Y. et al. Single-cell triple omics sequencing reveals genetic, epigenetic, and transcriptomic heterogeneity in hepatocellular carcinomas. Cell Res. 26 , 304–319 (2016).

ENCODE Project Consortium. An integrated encyclopedia of DNA elements in the human genome. Nature 489 , 57–74 (2012).

Kacmarczyk, T. J. et al. ‘Same difference’: comprehensive evaluation of four DNA methylation measurement platforms. Epigenetics Chromatin 11 , 21 (2018).

Chen, C. et al. Single-cell whole-genome analyses by Linear Amplification via Transposon Insertion (LIANTI). Science 356 , 189–194 (2017).

Hovestadt, V. et al. Decoding the regulatory landscape of medulloblastoma using DNA methylation sequencing. Nature 510 , 537–541 (2014).

Guo, F. et al. Active and passive demethylation of male and female pronuclear DNA in the mammalian zygote. Cell Stem Cell 15 , 447–459 (2014).

Gaiti, F. et al. Epigenetic evolution and lineage histories of chronic lymphocytic leukaemia. Nature 569 , 576–580 (2019).

Ghandi, M. et al. Next-generation characterization of the Cancer Cell Line Encyclopedia. Nature 569 , 503–508 (2019).

Download references

Acknowledgements

B.E.B. is the Bernard and Mildred Kayden Endowed MGH Research Institute Chair and an American Cancer Society Research Professor. This research was supported by a Pioneer Award from the National Institutes of Health Common Fund and National Cancer Institute (DP1CA216873) and by the Gene Regulation Observatory at the Broad Institute. S.J.S. was supported by a Medical Scientist Training Award from the National Institute of General Medical Sciences (T32GM007753). V.H. was supported by a Human Frontier Science Program long-term fellowship (LT000596/2016-L). We thank R. Boursiquot for technical assistance and L. Gaskell, W. Flavahan and other Bernstein lab members for helpful discussions.

Author information

Authors and affiliations.

Department of Pathology and Center for Cancer Research, Massachusetts General Hospital and Harvard Medical School, Boston, MA, USA

Sarah J. Shareef, Samantha M. Bevill, Ayush T. Raman, Martin J. Aryee, Peter van Galen, Volker Hovestadt & Bradley E. Bernstein

Broad Institute of Harvard and MIT, Cambridge, MA, USA

Department of Biostatistics, Harvard T.H. Chan School of Public Health, Boston, MA, USA

Martin J. Aryee

Division of Hematology, Brigham and Women’s Hospital and Harvard Medical School, Boston, MA, USA

Peter van Galen

Department of Pediatric Oncology, Dana-Farber Cancer Institute, Boston, MA, USA

Volker Hovestadt

Division of Hematology/Oncology, Boston Children’s Hospital, Boston, MA, USA

You can also search for this author in PubMed   Google Scholar

Contributions

S.J.S., P.v.G, V.H. and B.E.B. conceptualized and designed experiments. S.J.S. optimized and executed XRBS experiments. S.J.S. and S.M.B. performed decitabine experiments. V.H. performed computational analyses. A.T.R. contributed to computational analyses of single-cell data. M.J.A., V.H. and B.E.B. provided senior guidance. S.J.S., V.H. and B.E.B. wrote the manuscript with assistance from other authors, all of whom approved the final submission.

Corresponding authors

Correspondence to Volker Hovestadt or Bradley E. Bernstein .

Ethics declarations

Competing interests.

B.E.B. discloses financial interests in Fulcrum Therapeutics, HiFiBio, Arsenal Biosciences, Cell Signaling Technologies and BioMillenia. The remaining authors declare no competing financial interests.

Additional information

Peer review information Nature Biotechnology thanks Jonas Demeulemeester and the other, anonymous, reviewer(s) for their contribution to the peer review of this work.

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Extended data

Extended data fig. 1 evaluation of single mspi anchor design for methyl-cpg profiling..

a , Plots show results from an in silico MspI restriction digest analysis of the human genome. The cumulative number of MspI fragments (total of 2.3 million, left), of basepairs (total of 3.1 billion, middle), and of CpGs (total of 29.4 million, right) is shown relative to increasing MspI fragment length. Vertical dotted lines show the size range of fragments captured in typical RRBS experiments. This analysis shows that RRBS of MspI fragments 40–120 bases in length covers only 0.9% of the genome, but enriches for 5.6% of genomic CpGs. Recent implementations of RRBS (for example enhanced RRBS 14 , 15 ) that consider fragments up to 320 bases in length cover an additional 9.7% of CpGs. Approximately 35.0% of CpGs that are located within 300 bases of a single MspI site are not captured by these techniques. b , Histogram shows coverage depth of MspI restriction sites for individual replicates of a 10 ng XRBS library (left, middle), and both replicates combined (right). c , Heatmap shows coverage depth of CpGs between replicates of a 10 ng XRBS library (Pearson’s r = 0.90). d , Histogram shows coverage depth of CpGs in the combined dataset of both replicates (n = 2 independently generated libraries). e , Plot shows unique reads as a function of aligned reads in millions. With greater sequencing depth the fraction of unique reads decreases, as the chance of sampling a non-unique read (that is PCR duplicate) increases.

Extended Data Fig. 2 Comparison of MspI fragment detection and CpG coverage over regulatory elements.

a , Plot shows number of detected fragments plotted as a function of calculated MspI fragment length from XRBS 10 ng library replicates and from public RRBS and enhanced RRBS (ERRBS) datasets. Because of the random hexamer-primed second strand elongation step, XRBS efficiently detects fragments that exceed the selected fragment size range in RRBS (40–220 bp) and ERRBS (70–320 bp). XRBS less efficiently captures short fragments (<70 bp) compared to ERRBS and RRBS. Peaks in the graph correspond to three fragments commonly generated from Alu repetitive elements. b , Plot compares CpG coverage as a function of sequencing depth (x-axis) for XRBS (red), WGBS (blue, ENCODE), ERRBS (orange, obtained from ref. 47 ) and RRBS (green, ENCODE). c , Downsampling analysis plot as in panel B but restricted to CpGs within CpG islands (top) and gene promoters (bottom). d , Downsampling analysis plot as in panel B but restricted to CpGs within H3K27ac peaks (top) and CTCF binding sites (bottom).

Extended Data Fig. 3 XRBS efficiently covers CpGs in regulatory elements and repetitive regions.

a , Plots show the number of proximal enhancer-like, distal enhancer-like, and CTCF-only elements (as defined in the ENCODE SCREEN database 35 ) with at least 25-fold combined coverage as a function of sequencing depth for XRBS (red), WGBS (blue), ERRBS (orange), and RRBS (green). Enrichment for functional elements at a uniform sequencing depth of 10 billion base pairs is indicated. Vertical grey line indicates break in x-axis scale. b , Heatmaps depict genomic regions (rows, n = 3,725,365 LTRs, SINEs, and LINEs) containing different repeat element families (as defined by RepeatMasker). Individual repeat elements are divided into 50 equally sized windows (5’ and 3’ position indicated). Upstream and downstream regions (±200 bp) are divided into 25 equally sized windows. Panels from left to right show DNA methylation calls from 450k methylation array, RRBS, ERRBS, XRBS, and WGBS. c , Plots compare CpG coverage within different repeat element families as a function of sequencing depth for XRBS (red), WGBS (blue), ERRBS (orange), and RRBS (green). CpG enrichment relative to WGBS is indicated. In comparison to RRBS, XRBS enriches for most repeat families, with the exception of Alu and ERV1 elements that frequently contain MspI restriction sites and are also efficiently captured (see also Extended Data Fig. 2a ).

Extended Data Fig. 4 Correlation of DNA methylation with histone marks and compartment calls.

a , Plot shows unique reads as a function of aligned reads in low-coverage XRBS libraries from K562, HL-60, OCI-AML3, and Kasumi-1 cells. b , Plot shows unique reads as a function of aligned reads in low-coverage libraries from K562, Kasumi-1, HL-60, OCI-AML3 cells. Libraries were generated from 1,000 (green) and 100 (orange) cells sorted directly into lysis buffer. Libraries generated from 1,000 cells are comparable to libraries generated from 10 ng of purified DNA (panel A), whereas 100 cell libraries show reduced complexity. c , Heatmap shows Pearson correlation of XRBS methylation profiles of 100 kb windows generated from 10 ng gDNA, 1,000 or 100 sorted cells across four cell lines. Dendrogram derived from unsupervised clustering is indicated to the left. Sample grouping by DNA methylation is consistent with cell identity, indicating low technical variability between input material. d , Heatmaps show correlation between average DNA methylation values and signal for H3K9me3 (left), H3K27me3 (center), and H3K36me3 (right) in 100kb-windows for K562 cells. e , Heatmap shows correlation between DNA methylation and the Hi-C-derived first eigenvector indicating compartment A (positive values) and compartment B (negative values) in 100kb-windows for K562 cells. f , Heatmap shows correlation between average DNA methylation values and ChIP seq signal for H3K9me3 (top), H3K27me3 (middle), and H3K36me3 (bottom) in 100kb-windows for human H1 embryonic stem cells, primary T cells and mammary epithelial cells, and cultured IMR90, GM12878 and K562 cells. g , Heatmap as in panel F, but shows correlation between average DNA methylation values and the Hi-C-derived first eigenvector (x-axis). Positive values correspond to compartment A and negative values correspond to compartment B. While hypomethylation of compartment B is most pronounced in K562 cells, a similar trend is also observed in other cultured cell lines and in primary mammary epithelial cells.

Extended Data Fig. 5 Characterization of decitabine treatment of cancer cell lines.

a , Plot shows dose response curve for decitabine treatment of three cell lines Kasumi, HL-60, and OCI-AML3. Viability was measured using cell titer glo and is reported as luminescence relative to control DMSO treated cells (n = 3 independently treated replicates, error bars represent standard deviation). b , Images show HL60 and OCI-AML3 cells treated with 300 nM decitabine and a DMSO vehicle control. Morphology of decitabine treated cells similar to control, repeated three times. Scale bar is indicated and applies to all images. c , Plot shows unique reads as a function of aligned reads in XRBS libraries from DMSO- and decitabine-treated HL-60 and OCI-AML3 cells. d , Barplot shows average DNA methylation values across island (dark grey) and non-island (light grey) CpGs in DMSO- and decitabine-treated HL-60 and OCI-AML3 cells. For example, average methylation of non-island CpGs in HL-60 cells is reduced from 0.68 to 0.54 by decitabine treatment (20.2% reduction, n = 1 library per treatment). e , Heatmap shows correlation between Hi-C-derived first eigenvectors from K562 and HL-60 cell lines in 100kb-windows, indicating high agreement in compartment structure between both cell lines. f , Heatmaps show correlation between average DNA methylation values and Hi-C-derived eigenvector in 100kb-windows for DMSO- (left) and decitabine-treated HL-60 cells (center). Heatmap on the right shows relative DNA methylation values of decitabine- and DMSO-treated cells. Despite compartment B showing lower methylation compared to compartment A at baseline, induced DNA hypomethylation with decitabine treatment affects compartment A and B equally.

Extended Data Fig. 6 Differential DNA methylation of gene promoters.

a , Plot shows unique reads as a function of aligned reads in 1,000 cell high-coverage libraries of four cell lines. b , Heatmap depicts 8 kb regions (rows, n = 3,972 promoters) centered at transcription start sites that show cell line-specific hyper- or hypomethylation (as in Fig. 4a ) and divided into 100 equally sized windows. Panels from left to right show methylation calls from 450k methylation array, RRBS, XRBS, and WGBS in K562 cells. All datasets except XRBS were retrieved from ENCODE 46 . c , Plot shows expression levels for genes that were associated with cell line-specific promoter hyper- and hypo-methylation. Genes with an expression level larger than 0.5 are considered as expressed. Average gene expression levels are indicated by horizontal lines. P -values were generated using a two-sided Mann-Whitney U test. In K562, the majority (74.0%) of hypomethylated promoters are associated with non-expressed genes, which is unique to this cell line, consistent with global hypo-methylation in K562. d , Scatterplot compares gene expression level and H3K4me3 ChIP-seq signal for gene promoters that were identified as differentially methylated between all four cell lines. Individual promoters (dots) are colored if specifically hypermethylated (red) and hypermethylated (blue) in K562 cells. This analysis shows that promoters which are not expressed and specifically hypomethylated in K562 (n = 1,624 promoters) are generally negative for the H3K4me3 histone mark (98.7%), whereas promoters that are hypomethylated and expressed (n = 570) are more frequently positive for H3K4me3 (45.0%).

Extended Data Fig. 7 Evaluating the use of XRBS DNA methylation profiling to predict H3K27 acetylation and CTCF binding.

a , Heatmap depicts 8 kb regions (rows, n = 15,202 peaks) centered on H3K27ac peaks, grouped into regions that are hypomethylated specifically in K562 or OCI-AML3 cells (as in Fig. 4b ). Peaks that are not specifically hypomethylated (‘Others’) are downsampled for visualization. Regions are divided into 100 equally sized windows. Panels from left to right show: methylation calls from 450k methylation array, RRBS, XRBS, and WGBS in K562 cells. All datasets except XRBS were retrieved from ENCODE 46 . b , Scatterplot shows merged H3K27ac peaks from OCI-AML3 and K562 ChIP-seq datasets. Individual peaks (dots) are colored if specifically hypomethylated in K562 (blue) or OCI-AML3 (red) cells. c , Line plot (bottom) shows difference in methylation between K562 and OCI-AML3 cells over merged H3K27ac peaks (n = 15,202 peaks). Of these peaks, 8.3% and 2.3% are specifically hypomethylated in K562 (methylation difference ≤ −0.5) and OCI-AML3 ( ≥ 0.5) cells, respectively. Bar plot (top) shows the fraction of cell line-specific H3K27ac peaks within 100 equally sized bins grouped by difference in methylation. Shared peaks are indicated in gray. d , Receiver operating characteristic (ROC) curve shows performance of predicting cell line-specific H3K27ac peaks based on difference in DNA methylation over peaks that are covered by XRBS. Sensitivity and specificity are indicated at different thresholds (±0.2 and ±0.5, as in panel C). e , Heatmap depicts 2 kb regions (rows, n = 7,629 peaks) centered at merged CTCF peaks from K562 and HL-60 Chip-seq datasets. Individual peaks (dots) are colored if specifically hypermethylated in K562 or HL-60 cells (as in Fig. 4c ). Peaks not specifically hypermethylated (‘Others’) are downsampled for visualization. Panels from left to right show methylation calls from 450k methylation arrays, RRBS, XRBS, and WGBS in K562 cells. All datasets except XRBS were retrieved from ENCODE 46 . f , Scatterplot shows merged CTCF peaks from K562 and HL-60 ChIP-seq datasets. Individual CTCF binding sites (dots) are colored if specifically hypermethylated in K562 (red) or HL-60 (blue) cells. g , Line plot (bottom) shows difference in methylation between K562 and HL-60 cells over merged CTCF peaks (n = 7,629 peaks). Bar plot (top) shows the fraction of cell line-specific CTCF peaks within 100 equally sized bins grouped by difference in methylation. Shared peaks are indicated in gray. h , ROC curve shows performance of predicting cell line-specific CTCF peaks based on difference in DNA methylation over peaks that are covered by XRBS. Sensitivity and specificity are indicated at different thresholds (±0.2 and ±0.5, as in panel G).

Extended Data Fig. 8 XRBS profiling of limited human bone marrow cell types.

a , Plots show the gating strategy for fluorescence assisted cell sorting (FACS) of human bone marrow of CD34 + HSPCs, CD3 + T cells, and CD14 + monocytes. Singlets (FSC-W vs. -H) and viable cells (PI vs. FSC-A) were sorted based on cell surface marker signal. b , Plot shows unique reads as a function of aligned reads in libraries from unsorted human bone marrow, HSPCs, monocytes, and T cells. Libraries were generated from 100 sorted cells. 1000 cells were used for the unsorted bone marrow library. c , Heatmap depicts 4 kb regions (rows, n = 2,170 regions) centered over elements defined in the ENCODE SCREEN database. Only differentially methylated elements between monocytes and T cells are shown. Elements were stratified by their methylation status in HSPCs (hypomethylated: top; hypermethylated: bottom). Methylation levels of the unsorted bone marrow are shown for comparison (left). ATAC-seq signal for sorted hematopoietic stem cells (HSCs), monocyte, and CD4 + T cells (obtained from 36 ).

Extended Data Fig. 9 Evaluation of single cell XRBS profiles.

a , Plot shows unique reads as a function of aligned reads in single cell XRBS profiles (n = 96 cells). With greater sequencing depth the fraction of unique reads decreases, as the chance of sampling a non-unique read (that is PCR duplicate) increases. b , Boxplots compare DNA methylation profiles from human scXRBS (n = 59 cells) and three published scRRBS datasets generated from human cells: Chronic lymphocytic leukemia (n = 282 cells) 51 , hepatocellular carcinoma and HepG2 cells (n = 34 cells) 45 , and oocytes, sperm and pronuclei (n = 35 cells) 50 . Single cells from Hou et al. were generated using the scTrio-seq protocol that in part resembles scRRBS. Only CpGs within 75 bases of an MspI cut site were considered for scRRBS libraries to adjust for differences in read lengths. Libraries from Gaiti et al. were sequenced at 2×51 bases. Left plot shows the number of paired-end reads sequenced for each cell. Other plots show the number of CpGs covered (≥1-fold) across all CpGs in the genome, CpGs within distal enhancer-like regions, and CpGs within ‘CTCF-only’ regions (SCREEN database 35 ). Both strands of a CpG dinucleotide are assessed individually. Although sequenced at the lowest depth, scXRBS libraries on average capture the most CpGs, particularly in CpG-sparse regions. Boxplots were generated in R using default settings: Bounds of box and thick horizontal line represent 25 th , 75 th , and 50 th percentile of observations, whiskers represent minimum and maximum observations, and outliers are indicated as dots. c , Barplot shows the fraction of unique reads (that is reads not representing PCR duplicates) per single cell library. Within the same PCR reaction, the duplicate rate was very similar, irrespective of the total number of aligned reads per single cell. Each bar plot represents a single cell XRBS library. Twenty four barcoded cells were in each of 4 independent libraries. d , Heatmap compares alternate allele frequencies from SNP array data for K562 and HL-60 cell lines. Cell line-specific homozygous alleles are indicated by white boxes boxes and were used for single cell SNP analysis in Fig. 5d . e , Plots show copy number variation calls from combined single cell XRBS profiles (top) and whole exome sequencing data (middle) for K562 cells. A number of chromosomes show differences in copy number between XRBS and whole exome sequencing (bottom). However, these differences likely represent true copy number variations between K562 cells used for these experiments. f , Heatmap shows pairwise correlation coefficients of single cell methylation profiles. Dendrogram shows unsupervised clustering. Single cell XRBS profiles cluster by cell type. g , Barplot shows K562 single cell average DNA methylation values within various early and late replicating regions. Each bar represents an individual K562 single cell library. There are 32 single cell libraries plotted for each cell cycle phase. h , Heatmap shows pairwise correlation of average DNA methylation values within various early and late replicating regions. Late replicating regions (G2 phase) cluster separately. These results suggest that one source of single cell DNA methylation heterogeneity is decreased fidelity of maintenance DNA methylation in late replicating domains.

Supplementary information

Reporting summary, supplementary tables 1–5, rights and permissions.

Reprints and permissions

About this article

Cite this article.

Shareef, S.J., Bevill, S.M., Raman, A.T. et al. Extended-representation bisulfite sequencing of gene regulatory elements in multiplexed samples and single cells. Nat Biotechnol 39 , 1086–1094 (2021). https://doi.org/10.1038/s41587-021-00910-x

Download citation

Received : 23 August 2020

Accepted : 31 March 2021

Published : 06 May 2021

Issue Date : September 2021

DOI : https://doi.org/10.1038/s41587-021-00910-x

Share this article

Anyone you share the following link with will be able to read this content:

Sorry, a shareable link is not currently available for this article.

Provided by the Springer Nature SharedIt content-sharing initiative

This article is cited by

Reduced representative methylome profiling of cell-free dna for breast cancer detection.

• Qingmo Yang
• Xingqiang Zhu
• Hongliang Chen

Clinical Epigenetics (2024)

Methylomics and cancer: the current state of methylation profiling and marker development for clinical care

• Chengyin Liu

Cancer Cell International (2023)

Characterizing cis-regulatory elements using single-cell epigenomics

• Sebastian Preissl
• Kyle J. Gaulton

Nature Reviews Genetics (2023)

Joint single-cell profiling resolves 5mC and 5hmC and reveals their distinct gene regulatory effects

• Emily B. Fabyanic

Nature Biotechnology (2023)

Biomarkers of aging

• Guang-Hui Liu

Science China Life Sciences (2023)

Quick links

• Explore articles by subject
• Guide to authors
• Editorial policies

Sign up for the Nature Briefing newsletter — what matters in science, free to your inbox daily.

• Order Login
• At-Home Testing

Omics Services

• Research Services
• Clinical Services
• Whole Genome Sequencing
• Whole Exome Sequencing
• Targeted Gene Sequencing
• Microarray Analysis
• Bacterial Genome Sequencing
• Sanger Sequencing
• Whole Plasmid Sequencing
• Long-Read Sequencing
• Whole Genome Bisulfite Sequencing
• Reduced Representation Bisulfite Sequencing
• 16S/ITS Amplicon Sequencing
• Shotgun Metagenomic Sequencing
• mRNA Sequencing
• Total RNA Sequencing

Small RNA Sequencing

• Low Input RNA Sequencing
• Single-Cell Analysis
• Spatial Biology
• Biomarker Profiling
• Sanger Sequencing (SNP validation)
• Our Accreditations
• Our Newest Instruments
• Omics Insights Blog
• Customer Publications
• Our Newsletter
• Our Network

Omics Research // Epigenomics

Reduced representation bisulfite sequencing (rrbs).

Sequence the important targets of methylation of CpG sites throughout the genome.

What Is Reduced Representation Bisulfite Sequencing (RRBS)?

RRBS is used to get a “reduced representation" of the genome, with a focus on CpG islands. The restriction enzyme is used to digest the DNA during the fragmentation step. It enables us to profile genome-wide methylation on a single nucleotide scale using a MspI restriction enzyme digestion at CpG sites.

RRBS is an alternative to whole genome bisulfite sequencing (WGBS) that requires less sequencing coverage from a single sample.

Why Is Whole-Genome Bisulfite Sequencing Important?

RRBS enriches regions of high methylation, resulting in a focus on the regions of most interest. RRBS aids research in epigenomic-wide studies by studying the methylation status at CpG locations. It can be used to identify gene expression regulation and transcription modification that is used to identify different cell types and disease states.

Using RRBS, the methylation status of cysteine residues can be identified. The genomic region of the CpG and its function can also be identified and reported.

What Can You Do?

Focus your methylation studies.

Study areas of high methylation at CpG dinucleotides by using a restriction enzyme digestion

Reduce Costs, Lower Sequencing Coverage

Use less raw sequencing throughput than WGBS to study the methylation that plays an important role in regulation of transcription

Genome-Wide Studies Across Multiple Applications

Study regions of high methylation across the entire genome to inform on epigenomic changes in fields of immunology, cancer, and neuroscience

On-Going Promotions

We're always adding new deals — check in often for updated promotions.

RRBS Best Practices

Sample preparation & dna extraction.

DNA extraction from cells, tissue, FFPE, blood, and swab samples

Get your high-quality genomic DNA isolated from most sample types at Psomagen, or submit your extracted genomic DNA.

Sample Requirements

250 ng genomic DNA and DIN >7.0

Library Preparation

Zymo-seq RRBS library

• MspI restriction enzyme digestion
• Library QC >10 nM at 150-750 bp

Sequencing Instrumentation & Read Length

Illumina NovaSeq 6000

• 2 x 151 bp read length
• Human - at least 30 million reads raw sequencing depth (approximately 10X average coverage of the detected CpG sites)

Bioinformatic Analysis

Basic analysis:

• Data quality control

Advanced analysis:

• Methylation calling
• Functional location mapping
• Differentially methylated CpGs

Methylation-specific analyses are available.

Why Our Clients Choose Psomagen

Superior data quality.

Exceeds manufacturer benchmarks

Fast Turnaround

Starting at one week for sequencing projects

US-Based Processing

Sequencing and customer support in Rockville, MD

Ph.D. Level Support

Expert assistance every step of the way, with real-time project updates through our online system

Automated Workflows

Increase scalability and improve reproducibility

CLIA-Certified, CAP-Accredited Labs

Clinical research options available (inquire for more information)

Complete Sequencing Solutions

From extraction to data analysis, with many sample types accepted

Quotes From Our Clients

I chose Psomagen because of the price and US-based labs. We've used Psomagen before and were very satisfied. The people were responsive and helpful and the work done was very good.

University Research Laboratory

We're Ready for Your Order!

Here's what to expect in your project's life cycle. If you have a question during sequencing, it's easy to get in touch!

Request a Consultation

Create an account through our client portal and get in touch with your local field team to start designing a project.

We’ll help drive your research forward with the right tools for your research needs.

Place an Order

After finalizing the scope of work, the official quote will be issued for signing. Then, your sales representative will send specific shipping instructions.

Once all the samples and necessary documents are ready, we can start your project. Need more analysis? Reach out to our bioinformatics experts to dig deeper.

QC and Data Delivery

All Psomagen sequencing projects come with a built-in QC step to ensure you get the information you need. An account manager will update you on your project's progress in a timely manner.

You can also see some project information and securely view all valid reports on the client portal site .

Customer Publications and Other Resources

A chromosome-level reference genome and pangenome for barn swallow population genomics, whole genome sequencing identifies structural variants contributing to hematologic traits in the nhlbi topmed program, first complete mitochondrial genome of the saharan striped polecat (ictonyx libycus), go beyond rrbs.

There's other options for your epigenomics research! Find the right fit for you with other Psomagen services.

Whole-Genome Bisulfite Sequencing (WGBS)

Sequence genome-wide methylation status quickly and accurately via US-based, massively parallel, next-generation sequencing (NGS).

Understand the impact of small RNA on transcription and gene silencing with small RNA sequencing.

PSO Comprehensive!

What forms of whole-genome sequencing will get you the insights you need, subscribe to our newsletter.

• Our Certifications
• News & Events
• © 2023 Psomagen
• Cookie Settings
• HIPAA Statement
• Privacy Policy
• Website Terms & Conditions

Reduced Representation Bisulfite Sequencing (RRBS)

Affiliations.

• 1 Department of Maternal-Fetal Biology, Research Institute, National Center for Child Health and Development, Setagaya, Tokyo, Japan. [email protected].
• 2 Department of Obstetrics and Gynecology, The Jikei University School of Medicine, Minato, Tokyo, Japan.
• 3 Department of Obstetrics and Gynecology, Keio University School of Medicine, Shinjuku, Tokyo, Japan.
• 4 Department of Maternal-Fetal Biology, Research Institute, National Center for Child Health and Development, Setagaya, Tokyo, Japan.
• PMID: 36173564
• DOI: 10.1007/978-1-0716-2724-2_3

Reduced representation bisulfite sequencing (RRBS) enriches CpG-rich genomic regions using the MspI restriction enzyme-which cuts DNA at all CCGG sites, regardless of their DNA methylation status at the CG site-and enables the measurement of DNA methylation levels at 5% ~ 10% of all CpG sites in the mammalian genome. RRBS has been utilized in a large number of studies as a cost-effective method to investigate DNA methylation patterns, mainly at gene promoters and CpG islands. Here, we describe protocols for gel-free preparation of RRBS libraries, quality control, sequencing, and data analysis. Our protocols typically require nine cycles of polymerase chain reaction (PCR) amplification to obtain a sufficient amount of library for sequencing when 100 ng of genomic DNA is used as a starting material; moreover, it takes 3 days to complete library preparation and quality control procedures for up to eight samples.

Keywords: DNA methylation; Methylome; MspI; Next-generation sequencing; Reduced representation bisulfite sequencing.

© 2023. The Author(s), under exclusive license to Springer Science+Business Media, LLC, part of Springer Nature.

Publication types

• Research Support, Non-U.S. Gov't
• CpG Islands
• DNA Methylation*
• High-Throughput Nucleotide Sequencing
• Mammals / genetics
• Sequence Analysis, DNA / methods
• hydrogen sulfite

An official website of the United States government

The .gov means it’s official. Federal government websites often end in .gov or .mil. Before sharing sensitive information, make sure you’re on a federal government site.

The site is secure. The https:// ensures that you are connecting to the official website and that any information you provide is encrypted and transmitted securely.

• Publications
• Account settings

Preview improvements coming to the PMC website in October 2024. Learn More or Try it out now .

• Advanced Search
• Journal List
• Biol Proced Online

Performances of Different Fragment Sizes for Reduced Representation Bisulfite Sequencing in Pigs

Xiao-long yuan.

1 Guangdong Provincial Key Lab of Agro-Animal Genomics and Molecular Breeding, National Engineering Research Centre for Breeding Swine Industry, College of Animal Science, South China Agricultural University, Guangzhou, Guangdong China

2 Section of Comparative Pediatrics and Nutrition, Department of Veterinary Clinical and Animal Sciences, University of Copenhagen, Frederiksberg, Denmark

Rong-Yang Pan

3 Department of Animal Sciences, Georg-August University, Albrecht Thaer-Weg 3, Göttingen, Germany

Per Torp Sangild

Associated data.

The datasets generated during and/or analysed during the current study are available in the European Nucleotide Archive (accession number: PRJEB14111).

Reduced representation bisulfite sequencing (RRBS) has been widely used to profile genome-scale DNA methylation in mammalian genomes. However, the applications and technical performances of RRBS with different fragment sizes have not been systematically reported in pigs, which serve as one of the important biomedical models for humans. The aims of this study were to evaluate capacities of RRBS libraries with different fragment sizes to characterize the porcine genome.

We found that the Msp I-digested segments between 40 and 220 bp harbored a high distribution peak at 74 bp, which were highly overlapped with the repetitive elements and might reduce the unique mapping alignment. The RRBS library of 110–220 bp fragment size had the highest unique mapping alignment and the lowest multiple alignment. The cost-effectiveness of the 40–110 bp, 110–220 bp and 40–220 bp fragment sizes might decrease when the dataset size was more than 70, 50 and 110 million reads for these three fragment sizes, respectively. Given a 50-million dataset size, the average sequencing depth of the detected CpG sites in the 110–220 bp fragment size appeared to be deeper than in the 40–110 bp and 40–220 bp fragment sizes, and these detected CpG sties differently located in gene- and CpG island-related regions.

Conclusions

In this study, our results demonstrated that selections of fragment sizes could affect the numbers and sequencing depth of detected CpG sites as well as the cost-efficiency. No single solution of RRBS is optimal in all circumstances for investigating genome-scale DNA methylation. This work provides the useful knowledge on designing and executing RRBS for investigating the genome-wide DNA methylation in tissues from pigs.

In mammals, DNA methylation preferably occurs at CpG dinucleotides, and it modifies many key biological processes, including gene transcription [ 1 ], genomic imprinting [ 2 ], tissue differentiation [ 3 ] and phenotypic variation [ 4 ]. By profiling the DNA methylome, it is possible to mirror the epigenetic patterns that may regulate gene expression. However, the asymmetric distribution of CpG sites and the short length of sequencing reads make whole-genome bisulfite sequencing (WGBS) relatively costly [ 5 ]. Therefore, reduced representation bisulfite sequencing (RRBS) is developed, which uses Msp I, a restriction enzyme that cuts C|CGG sites, to select the CG-rich regions and reduces the required amount of sequencing to study the genome-wide DNA methylation [ 6 , 7 ].

RRBS has a low cost for per detected CpG site, and it is highly sensitive to low DNA input, while providing single-nucleotide resolution to quantify the DNA methylation level distribution [ 7 , 8 ]. The cost-effective RRBS processes allow for large-scale mapping of DNA methylation in a large number of samples (e.g. >100 per week) [ 9 ]. Recently, RRBS has been performed on samples from humans [ 9 ], pigs [ 10 ], sheep [ 5 ] and many model organisms [ 11 – 13 ] to generate the genome-scale DNA methylation and screen the dynamic changes in the methylomes. Additionally, the DNA methylation information obtained by RRBS could be scaled up to WGBS data [ 14 ] and used in studies on methylation quantitative trait loci [ 15 ]. New approaches and developments based on RRBS have been developed, such as the laser capture microdissection-RRBS [ 16 ], single-cell RRBS [ 17 ], double-enzyme RRBS [ 18 ], and high-throughput targeted repeat element RRBS [ 19 ], to profile the genome-wide DNA methylation of the representative genome or of the specific genomic features. These processes and developments have greatly advanced investigations of DNA methylomes in mammals.

In the vertebrate genome, the fragment size of 40–220 bp has been suggested for RRBS [ 20 ]. The applications and technical analyses of RRBS with this fragment size, such as the genomic coverage and coverage depth, have been systematically investigated in mice [ 20 ] and humans [ 21 ]. Furthermore, performances and technical assessments of RRBS with different fragment sizes but not the 40–220 bp fragment size have also been discussed to resolve the fragment size selection and sequence depth in livestock, e.g. sheep [ 5 ]. However, the applications and technical performances of RRBS with different fragment sizes have not been systematically reported in pigs, which serve as one of the important biomedical models for humans [ 10 , 22 ].

In this study, we attempted to report and discuss technical applications of RRBS with different fragment sizes on the sample from pigs. We first bioinformatically predicted distribution characteristics of the different Msp I-digested fragment sizes of 40–110 bp, 110–220 bp and 40–220 bp in the porcine genome, respectively. Then, RRBS libraries were built with these fragment sizes from the same porcine DNA sample, and these RRBS libraries were sequenced to show the mapping efficiencies, optimal sequencing quantities, and distributions of the detected CpG sites across the locations of genes and CpG islands (CGIs) at the genome scale. This work would provide the methodological information about RRBS for use in epigenomic investigations of pigs, and it sheds additional light on how to design RRBS with the appropriate fragment size for comprehensively representing the methylome of pigs.

Ethic Statement

The ovary was collected from one female Landrace × Yorkshire crossed gilt aged 180 days. Pig cares and experiments were approved by the Animal Care and Use Committer of the South China Agricultural University, Guangzhou, China (approval number: SCAU#2013–10). All experiments and conductions were performed in accordance with the guidelines and regulations of the Administration of Affairs Concerning Experimental Animals (Ministry of Science and Technology, China, revised in June 2004).

Simulation of Msp I Digestion and Size Selection

The porcine reference genome (Sscrofa 10.2), which was downloaded from the Ensembl Genome Browser ( http://www.ensembl.org/Sus_scrofa/Info/Index ), was digested using Msp I in the simulation. The restriction enzyme cutting site for Msp I was C|CGG. The single sequences between two consecutive restriction sites were extracted as a Msp I-digested segment. The fragment size of 40–220 bp was recommended for RRBS in mammalian genomes, and based on this, fragment sizes of 40–110 bp, 110–220 bp and 40–220 bp were selected to predict the number and distributions of Msp I-digested fragments. The gene locations were downloaded from the Ensembl Genome Browser ( http://www.ensembl.org/Sus_scrofa/Info/Index ). The upstream regions were 5 kb upstream regions to the transcription start sites, and the downstream regions were 5 kb downstream regions to the transcription end sites. Additionally, 5′ untranslated (5’UTR), coding sequence (CDS), intronic and 3′ untranslated (3’UTR) and non-coding regions were denoted the same to the Ensembl Genome Browser. The outside regions of the upstream, 5’UTR, CDS, intron, 3’UTR, downstream, and non-coding regions were defined as the intergenic regions. CGI locations were downloaded from UCSC ( http://hgdownload.soe.ucsc.edu/goldenPath/susScr3/database/ ). The 2 kb upstream and downstream regions of the CGIs were defined as the CGI shores. The 2 kb upstream and downstream regions of the CGI shores were defined as the CGI shelves. The outside regions of the CGIs, shores and shelves were defined as the inter CGI regions. The locations of repetitive elements were also downloaded from UCSC ( http://hgdownload.soe.ucsc.edu/goldenPath/susScr3/database/ ). The simulations and calculations were completed by Perl and R scripts.

RRBS Library Preparation and Sequencing

The library constructions and sequencing services were provided by RiboBio Co., Ltd. (Guangzhou, China). The processes and procedures for building RRBS libraries were based on the technical processes described by previously published RRBS studies [ 7 , 20 ]. Briefly, the porcine ovarian genomic DNA was first extracted using a DNeasy Blood & Tissue Kit (Qiagen, Beijing). After checking the quality of the extracted DNA, the ovarian genomic DNA was digested overnight with Msp I (New England Biolabs, USA). The sticky ends were filled with CG nucleotides and 3′ A overhangs were added to the Msp I-digested segments. Second, the methylated Illumina sequencing adapters with 3′ T overhangs were ligated to the digested segments; then, the products were purified. Afterwards, the 40–110 bp, 110–220 bp and 40–220 bp fragments were separately selected and converted by bisulfite using an EZ DNA Methylation Gold Kit (Zymo Research, USA). Finally, libraries of 40–110 bp, 110–220 bp and 40–220 bp fragments were PCR amplified, and each library was sequenced with one lane of an Illumina HiSeq 2500 as well as 100-bp paired-end reads (PE100). All reads were trimmed using Trim Galore (v0.4.0) software (Babraham Bioinformatics, http://www.bioinformatics.babraham.ac.uk/projects/trim_galore/ ) and a Phred quality score of 20 as the minimum. The adaptor pollution reads and multiple N reads (where N  > 10% of one read) were removed off to generate the clean reads. The quality control checks were performed by FastQC (v0.11.3) software (Babraham Bioinformatics, http://www.bioinformatics.babraham.ac.uk/projects/fastqc/ ). The clean RRBS data were mapped to the porcine reference genome (Sscrofa 10.2, http://www.ensembl.org/Sus_scrofa/Info/Index ) and were called the DNA methylation by Bismark software (v0.14.5) [ 23 ]. The first two nucleotides were trimmed from all the second read sequences to blunt-end the Msp I site. For the overlapped reads, only the methylation calls of read 1 were used for in the process by Bismark with the option “-- no_overlap”, in order to avoid scoring the overlapping methylation calls twice. The bisulfite conversion rates were calculated as the number of covered cytosines in the non-CpG context, which were converted, was divided by the total number of covered cytosines in the non-CpG context [ 20 ]. The conversion efficiencies of the fragment sizes of 40–110 bp, 110–220 bp and 40–220 bp were 99.22, 99.60 and 99.31%, respectively. The RRBS data of these three fragment sizes were submitted to the European Nucleotide Archive (accession number: PRJEB14111).

Sub-sampling of RRBS Data and Bioinformatic Analysis

To assess the cost performances of different fragment sizes with different dataset sizes, we randomly sampled subsets of the paired reads from the whole RRBS data of 40–110 bp, 110–220 bp and 40–220 bp in triplicate and then investigated the features of these sub-sampled data. We respectively triplicated the generations of 10, 7, and 15 increasingly sub-sampled data sets for the fragment sizes of 40–110 bp, 110–220 bp and 40–220 bp. During the sampling process, the minimum number of reads was set as 10 million, which were paired reads (5 million from read 1 and 5 million from read 2). Then the dataset size was increasingly paired by 10 million reads (5 million from read 1 and 5 million from read 2). All sub-sampled data were aligned to the porcine genome, and the detected CpG sites were extracted by Bismark [ 23 ]. The mapping efficiencies for the different sub-samples of data were the same as those of the total data, suggesting that the sampling process was successful. Uniquely mapped reads were retained for further calculations. For RRBS data, the detected CpG sites with more than three covered reads were remained for further analyses. The average values of triplications were presented in this study.

Length Distribution of Msp I-digested Segments in Three Fragment Sizes

Based on the simulating processes of one previous study [ 24 ], the porcine reference genome (Sscrofa10.2) was digested by Msp I in the simulation. Fragment sizes of 40–110 bp, 110–220 bp and 40–220 bp were selected to evaluate the performance of RRBS on pigs. The fragment sizes of 40–110 bp, 110–220 bp and 40–220 bp contained 385,352, 281,798 and 664,080 Msp I-digested segments (Table ​ (Table1) 1 ) and covered 3,550,514, 3,718,483 and 7,234,567 CpG sites, respectively (Table ​ (Table1). 1 ). The average counts of CpG sites per segment was 9.21, 13.20 and 10.89 in the fragment sizes of 40–110 bp, 110–220 bp and 40–220 bp (Table ​ (Table1), 1 ), respectively.

The contents of the three fragment sizes

The length distribution of the Msp I-digested segments between 40 and 220 bp is shown in Fig. ​ Fig.1. 1 . Compared with that of 110–220 bp, the specific feature of the 40–110 bp fragment size was that it harbored a high distribution peak at 74 bp, which contained 66,733 Msp I-digested segments (Fig. ​ (Fig.1). 1 ). Moreover, we found that 95.40% of 74-bp length digested segments overlapped with the repetitive elements. Then, the single base sequences belonging to the 74-bp length digested segments were extracted, and these sequences were aligned with the porcine reference genome by bowtie2 (v2.2.5) [ 25 ]. We found that there were only 10.61% uniquely aligned segments; moreover, 89.39% of the 74-bp length digested segments had multiple alignments.

The length distribution of the Msp I-digested segments between 40 and 220 bp

Mapping Efficiencies of These Three Fragment Sizes

To characterize the RRBS performances of these three fragment sizes, we generated more than 100, 70 and 150 million RRBS reads for the 40–110 bp, 110–220 bp and 40–220 bp libraries from the same DNA sample, respectively. Bismark [ 23 ] was used to map these RRBS data to the porcine genome. We found that the unique mapping efficiencies were 40.70, 59.18 and 48.06% for the 40–110 bp, 110–220 bp and 40–220 bp libraries, respectively. The multiple mapping efficiencies were 19.70, 9.99 and 13.24% for the 40–110 bp, 110–220 bp and 40–220 bp libraries, respectively (Fig. ​ (Fig.2). 2 ). These results indicated that the fragment size of 110–220 bp had the highest unique mapping alignment with the lowest multiple alignment (Fig. ​ (Fig.2). 2 ). Moreover, compared with the mapping efficiency of 110–220 bp, the relatively higher multiple alignments and the relatively lower unique mapping alignments of 40–110 bp and 40–220 bp libraries suggested that the 74-bp length digested segments, which were highly overlapped with the repetitive elements, might reduce the unique mapping efficiency.

Mapping efficiencies of the 40–110 bp, 110–220 bp and 40–220 bp fragment sizes

Optimal Sequencing Quantities for These Three Fragment Sizes

To investigate the optimal sequencing quantities of these fragment sizes, we randomly sampled different subsets of reads from the whole RRBS data of these three fragment sizes and generated 10, 7, and 15 increasingly sub-sampled data sets for the fragment sizes of 40–110 bp, 110–220 bp and 40–220 bp in triplicate, respectively. Uniquely mapped reads and the detected CpG sites with ≥3 covered reads (3×) were retained for further calculations.

The number distributions and sequencing saturations of the detected CpG sites with ≥5, 10 and 15 covered reads (5X, 10X and 15X) are shown in Fig. ​ Fig.3. 3 . As expected, the numbers of 5X, 10X and 15X detected CpG sites increased with increasing data sizes for these three fragment sizes (Fig. ​ (Fig.3a 3a – e ). However, the increasing speed of 10X detected CpG sites decreased when the data size was more than 70, 50 and 110 million reads for the 40–110 bp, 110–220 bp and 40–220 bp fragment sizes, respectively (Fig. ​ (Fig.3a 3a – e ). Moreover, the saturations and percentages of the 10X detected CpG sites over the 3× detected CpG sites also descended when the data size was more than 70, 50 and 110 million reads for the 40–110 bp, 110–220 bp and 40–220 bp fragment sizes, respectively (Fig. ​ (Fig.3d 3d – f ), suggesting that the cost effectiveness of these fragment sizes might decrease with a higher dataset size.

Distributions of detected CpG sites in the differently sub-sampled RRBS data. The number distributions of detected CpG sites with ≥5, 10 and 15 covered reads (5X, 10X and 15X) in the differently sub-sampled RRBS data for 40–110 bp ( a ), 110–220 bp ( c ) and 40–220 bp ( e ) fragment sizes in triplications. The percentages of 5X, 10X and 15X detected CpG sites over the 3× detected CpG sites in the differently sub-sampled RRBS data for 40–110 bp ( b ), 110–220 bp ( d ) and 40–220 bp ( f ) fragment sizes in triplications

Sequencing Depth of These Three Fragment Sizes

Considering the cost and acquired number of detected CpG sites, we selected the data size of 50 million reads to evaluate the sequencing depth of these three fragment sizes. Given a 50-million read dataset size, the fragment size of 40–220 bp detected more CpG sites with 5X than 40–110 bp and 110–220 bp (Fig. ​ (Fig.4). 4 ). The 40–220 bp fragment size detected almost the same number of CpG sites with 10X as for 110–220 bp and detected more sites than for 40–110 bp. However, the 40–220 bp fragment size detected fewer CpG sites with 15X than for 110–220 bp, while it was almost the same as that for 40–110 bp (Fig. ​ (Fig.4). 4 ). These results suggested that the average sequencing depth of the detected CpG sites in the 110–220 bp fragment size appeared to be deeper than in the 40–110 bp and 40–220 bp fragment sizes.

Distributions of detected CpG sites versus the differently covered depth for the three fragment sizes in 50 million reads in triplications

Distribution of Detected CpG Sites in These Three Fragment Sizes

Given a 50-million dataset size, the fragment size of 40–220 bp detected the highest number of CpG sites with 5X within gene- and CGI-related regions compared with 40–110 bp and 110–220 bp fragment sizes (Fig. ​ (Fig.5a, 5a , ​ ,b). b ). For the CpG sites with 10X, the 40–220 bp fragment size detected the highest number of CpG sites within gene-related regions and CGIs (Fig. ​ (Fig.5c, 5c , ​ ,d), d ), but the 110–220 bp fragment size detected the highest number of CpG sites within the CGI shores and CGI shelves compared with the other two fragment sizes (Fig. ​ (Fig.5d). 5d ). Furthermore, considering the CpG sites with 15X, the fragment size of 110–220 bp detected the most CpG sites within gene- and CGI-related regions compared with the 40–110 bp and 40–220 bp fragment sizes (Fig. ​ (Fig.5e, 5e , ​ ,f). f ). Interestingly, considering the CpG sites with 5X, 10X and 15X, the 40–110 bp fragment size always detected more CpG sites within the 5’UTR regions than those with a fragment size of 110–220 bp (Fig. 5a–e ).

Coverage of the detected CpG sites of these three fragment sizes across the gene-related and CGI-related regions for the whole porcine genome. The coverages of 5X, 10X and 15X detected CpG sites across the gene-related regions for 40–110 bp ( a ), 110–220 bp ( c ) and 40–220 bp ( e ) fragment sizes in triplications. The coverages of 5X, 10X and 15X detected CpG sites across the CGI-related regions for 40–110 bp ( b ), 110–220 bp ( d ) and 40–220 bp ( f ) fragment sizes in triplications

DNA methylation is an important epigenetic modification that plays a critical function in many biological processes. Profiling of the genome-wide DNA methylation allows for investigations of DNA methylation dynamics and epigenetic mechanisms of many key biological processes [ 10 , 26 ]. Compared with the other sequencing strategies, RRBS is a cost-effective pipeline to generate genome-wide DNA methylation at the single-nucleotide resolution [ 8 , 20 ]. By enriching in the CpG-rich regions and relying on bisulfite sequencing methods, RRBS reduces the sequencing requirement and enhances the sequencing depth and accuracy of DNA methylation information for the targeted regions of the genome-scale DNA methylation [ 27 ]. In this study, the 40–110 bp, 110–220 bp and 40–220 bp fragment sizes were selected to assess the capabilities, data utilization efficiencies and cost performances of RRBS in pigs. We found that there was a length distribution peak at 74 bp, which highly overlapped with the repetitive elements and might reduce the unique mapping alignment. The 110–220 bp library displayed the highest unique mapping alignment and the lowest multiple alignment. When the data sizes are more than 70, 50 and 110 million reads for the fragment sizes of 40–110 bp, 110–220 bp and 40–220 bp, respectively, the cost effectiveness of these fragment sizes might decrease. Given a 50-million dataset size, the average sequencing depth of the detected CpG sites in the 110–220 bp fragment size appeared to be deeper than in the 40–110 bp and 40–220 bp fragment sizes.

In vertebrate genomes, the 40–220 bp fragment size is commonly used to resolve the mammalian methylome with RRBS. With a fragment size of 40–220 bp, the alignment efficiencies were approximately 30–40% for different mouse cells by 36-bp single-end sequencing [ 6 ]. However, the percentages of uniquely mapped reads only ranged from 27.0 to 32.7% for zebrafish with 100-bp single-end sequencing [ 11 ]. The mapping ratio was approximately 40–50% for human embryos and sperm cells with PE100 [ 28 ]. The uniquely mapped reads were approximately 48% for pigs with 50-bp paired-end sequencing (PE50) [ 10 ]. In this study, we also found that the unique mapping efficiency was approximately 48% for the 40–220 bp fragment size in pigs (Fig. ​ (Fig.2). 2 ). However, the unique mapping efficiency increased to approximately 60% for the 110–220 bp fragment size in pigs (Fig. ​ (Fig.2 2 ).

Multiple factors contributed to the relatively low utilization efficiency for the 40–220 bp fragment size. First, there were many short segments in the range of 40–220 bp (Fig. ​ (Fig.1). 1 ). The short segments were easily aligned to multiple locations with the present mapping model, which reduced the unique mapping ratio. Second, the selection of the sequencing strategy might be not appropriate. For example, considering the case of sequencing the 40–220 bp library with PE50, there was always a region (i.e., 50 bp for 150-bp segments and 120 bp for 220-bp segments in the center) that could not be covered by any read, resulting in loss of methylation information harbored by the uncovered regions. Third, large repetitive sequences might be located in the 40–220 bp fragment size, aggravating the ratio of multiple alignments [ 29 ].

One previous study suggested that short reads and a large number of repetitive elements might decrease the unique mapping efficiency of RRBS data [ 29 ]. In this study, we found that 95.40% of the 74-bp Msp I-digested segments overlapped with repetitive elements, and 89.39% of these segments were aligned to multiple locations. Furthermore, the unique mapping efficiency of the 40–110 bp fragment size was 40.10% (Fig. ​ (Fig.2), 2 ), which was the lowest compared with the other two fragment sizes (Fig. ​ (Fig.2). 2 ). In addition, one previous study recommended that there are redundant microsatellites, one of the repetitive elements, located in the Msp I-digested segments for 40–220 bp in the mouse genome [ 30 ], and these microsatellites might result in the low alignment efficiencies of the RRBS data from mouse cells by the 36-bp single-end sequencing [ 6 ]. These results showed that repetitive elements might decrease the unique mapping efficiency.

Given the 50, 70, 90 and 110 million reads of RRBS data, the 40–220 bp fragment size detected 2,233,833; 3,122,473.33; 3,699,966 and 4,058,334.33 CpG sites with 10X, respectively (Fig. ​ (Fig.3e). 3e ). Given the 130 and 150 million reads of RRBS data, the 40–220 bp fragment size detected 4,285,318.67 and 4,438,622.67 CpG sites with 10X, respectively (Fig. ​ (Fig.3e), 3e ), respectively. Compared with the 50 million reads, the dataset sizes of 70, 90 and 110 million reads increased by 40, 80 and 120%, and the numbers of detected CpG sites with 10X increased by 39.78, 65.63 and 81.68%, respectively (Fig. ​ (Fig.3e). 3e ). However, compared with the 50 million reads, the dataset sizes of the 130 and 150 million reads increased by 160 and 200%, but the numbers of detected CpG sites with 10X only increased by 91.84 and 98.70%, respectively (Fig. ​ (Fig.3e). 3e ). As a result, when the data set was more than 110 million reads, the cost-efficiency decreased for the fragment size of 40–220 bp.

The selections of fragment sizes could affect the numbers and sequencing depth of detected CpG sites as well as the mapping efficiency. For the 50 million reads, the fragment size of 40–220 bp detected 2,233,833 CpG sites with 10X, which was almost the same number of CpG sites with 10X for 110–220 bp, while it was more than for 40–110 bp (Fig. ​ (Fig.4). 4 ). However, the average sequencing depth per detected CpG site was 18.82 for the fragment size of 40–220 bp, which was lower than for 110–220 bp (26.56) and 40–110 bp (21.60). Moreover, the unique mapping efficiency of 40–110 bp fragment size was 40.70% with PE100, but it approximately increased to 60% for the 110–220 bp fragment size (Fig. ​ (Fig.2). 2 ). The unique mapping efficiency of 50–150 bp fragment size was 38.3% with PE100 in sheep, but it increased to 61.4% for the 150–220 bp fragment size [ 5 ]. Therefore, the selections of fragment sizes could affect the numbers and sequencing depth of detected CpG sites as well as the cost-efficiency. Taken together, when designing methylome studies using RRBS, researchers should consider the fragment size, mapping efficiency, sequencing depth, covered CpG sites, and dataset size. No single solution of RRBS is optimal in all circumstances for investigating genome-scale DNA methylation.

Acknowledgements

The authors thank RiboBio Co., Ltd. (Guangzhou, China) for providing the library constructions and sequencing services.

This work was supported by the Special Program for Applied Research on Super Computation of the NSFC-Guangdong Joint Fund (the second phase), Basic Work of Science and Technology Project (2014FY120800), Guangdong Sailing Program (2014YT02H042) and the earmarked fund for China Agriculture Research System (CARS-36).

Availability of Data and Materials

Authors’ contributions.

Conceived and designed the experiments: XLY, JQL and ZZ. Prepared biological samples: XLY, XD. Bioinformatic analyses: XLY, RYP and NG. Wrote the paper: XLY, ZZ and BL. Revised the paper: JQL, HZ and PTS. All authors reviewed the manuscript. All authors read and approved the final manuscript.

Authors’ Information

Not applicable.

Competing Interests

The authors declare that they have no competing interests.

Consent for Publication

Ethics approval, publisher’s note.

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Abbreviations

Contributor information.

Xiao-Long Yuan, Email: moc.361@230gnolxy .

Zhe Zhang, Email: nc.ude.uacs@gnahzehz .

Rong-Yang Pan, Email: nc.ude.uacs.uts@napyr .

Ning Gao, Email: ten.haey@liam_oaggnin .

Xi Deng, Email: moc.361@4090gnedx .

Bin Li, Email: nc.ude.uacs.uts@ilnib .

Hao Zhang, Email: nc.ude.uacs@oahgnahz .

Per Torp Sangild, Email: kd.uk.dnus@stp .

Jia-Qi Li, Email: nc.ude.uacs@ilqj .