Direct detection of DNA methylation during single-molecule, real-time sequencing

Abstract

We describe the direct detection of DNA methylation, without bisulfite conversion, through single-molecule, real-time (SMRT) sequencing. In SMRT sequencing, DNA polymerases catalyze the incorporation of fluorescently labeled nucleotides into complementary nucleic acid strands. The arrival times and durations of the resulting fluorescence pulses yield information about polymerase kinetics and allow direct detection of modified nucleotides in the DNA template, including N6-methyladenine, 5-methylcytosine and 5-hydroxymethylcytosine. Measurement of polymerase kinetics is an intrinsic part of SMRT sequencing and does not adversely affect determination of primary DNA sequence. The various modifications affect polymerase kinetics differently, allowing discrimination between them. We used these kinetic signatures to identify adenine methylation in genomic samples and found that, in combination with circular consensus sequencing, they can enable single-molecule identification of epigenetic modifications with base-pair resolution. This method is amenable to long read lengths and will likely enable mapping of methylation patterns in even highly repetitive genomic regions.

Figure 1

Principle and corresponding example of detecting DNA methylation during SMRT sequencing. (a) Schematics of polymerase synthesis of DNA strands containing a methylated (top) or unmethylated (bottom) adenosine. (b) Typical SMRT sequencing fluorescence traces from these templates. Letters above the fluorescence trace pulses indicate the identity of the nucleotide incorporated into the growing complementary strand. The dashed arrows indicate the IPD before incorporation of the cognate T, and, for this typical example, the IPD is ∼5× larger for mA in the template compared to A.

Figure 2

SMRT sequencing-mediated detection of methylated DNA bases. All three panels show the ratio of the average IPD in the methylated template to the average IPD in the control template, plotted versus DNA template position. In the region shown, the two templates are identical except at the two positions marked by triangles. Polymerase synthesis runs in the direction of increasing position number. While in all cases the two templates have a circular topology and are 199 bases in length, only 90 base segments surrounding the methylated regions are shown for clarity. Error bars indicate the s.e.m. IPD ratio at each template position (average n = 346 measurements for each position in (a), average n = 504 for each position in (b), and average n = 393 for each position in (c), computed using bootstrapping techniques (Online Methods).

Figure 3

Principal component analysis of C, mC, and hmC IPD and PW signatures. Each principal component is a linear combination of the mean IPD and PW at positions 71-79, which surround the variably modified template position 73. The weightings of IPD and PW at each position are shown in the Supplementary Table 1. Data points on the plot were computed by projecting a random 20% subsample of the IPD and PW values onto the first two principal components (PC1 and PC2) and then converting to a z-score (Online Methods).

Figure 4

IPD distributions for A and mA in synthetic DNA templates. (a) Schematic of the DNA templates with a total length of 199 bases. (b) IPD distributions at the indicated positions for both templates. For each row, the histograms depict the distributions of mean IPD (averaged over the indicated number of circular subreads). (c) Receiver operating characteristic (ROC) curves, based on the IPD distributions from the differentially methylated position in (b) and parameterized by IPD threshold, for assigning a methylation status to an adenosine nucleotide after one (gray), three (red), or five (blue) circular consensus sequencing subreads. The black dashed line depicts the ROC curve for randomly guessing the methylation status. Note that because the templates have a length of 199 bases, five full circular subreads correspond to read lengths of nearly 1000 bases.

Figure 5

Comparison of SMRT sequencing kinetics for DNA samples propagated within dam+ E. coli and for the same samples after whole-genome amplification (WGA). The sample comprises a 3.7-kb subregion of a C. elegans fosmid cloned into an E. coli vector. (a) 50-bp window GC-content of the sample, plotted versus template position. (b) Average IPD at each template position within the dam+ sample. (c) Average IPD at each template position within the WGA sample. (d) Ratio of the average IPDs (dam+ in (b) divided by WGA in (c)), plotted versus template position. Positions with a GATC context, where methylation of adenine at the sequence motif GATC is expected, are denoted by black squares, and all other positions are denoted by open blue circles. Error bars at the GATC positions denote the s.e.m. IPD ratio at those positions (average n = 106 measurements at each position). For comparison, the mean ± s.d. of all IPD ratios at non-GATC positions (open blue circles) is 1.00 ± 0.24 (n is ∼389,000 measurements). Average sequencing coverage across this fosmid region was 121-fold for the dam+ sample and 91-fold for the WGA sample.