Long-read single-molecule maps of the functional methylome

Abstract

We report on the development of a methylation analysis workflow for optical detection of fluorescent methylation profiles along chromosomal DNA molecules. In combination with Bionano Genomics genome mapping technology, these profiles provide a hybrid genetic/epigenetic genome-wide map composed of DNA molecules spanning hundreds of kilobase pairs. The method provides kilobase pair-scale genomic methylation patterns comparable to whole-genome bisulfite sequencing (WGBS) along genes and regulatory elements. These long single-molecule reads allow for methylation variation calling and analysis of large structural aberrations such as pathogenic macrosatellite arrays not accessible to single-cell second-generation sequencing. The method is applied here to study facioscapulohumeral muscular dystrophy (FSHD), simultaneously recording the haplotype, copy number, and methylation status of the disease-associated, highly repetitive locus on Chromosome 4q.

Figure 1.

ROM experimental scheme. (A, top) M.TaqI catalyzes the transfer of a TAMRA fluorophore from the cofactor AdoYnTAMRA onto the adenine residue that lies within its TCGA recognition site. (Bottom) If the CpG nested within the M.TaqI recognition site is methylated, the reaction is blocked. (B) Scatter plot comparing nonmethylation levels in two biological ROM replicates in 10-kbp windows along the human genome. (C) Representative dually labeled molecule (red dots indicate genetic barcode; green dots, methylation profile). The molecule's ROM fluorescence intensity profile is presented in green above the molecule image. (D) Representative field of view of DNA molecules (blue) fluorescently labeled in two colors and stretched in nanochannel arrays. Red dots indicate genetic labels; green dots, methylation profile. (E) Digitized representation of single molecules (yellow) aligned to an in silico generated reference (blue) according to their distinct genetic barcode (red dots). The positions of epigenetic labels (green dots) are inferred from alignment results.

Figure 2.

Global comparison of ROM and WGBS results in functional genomic elements. (A) Nonmethylation levels as a function of distance from histone modification peaks, for ROM and WGBS (blue indicates H3K4me3; red, H3K4ac). (B) Nonmethylation levels across gene bodies, in correlation with overall gene methylation level, determined by WGBS (red indicates low methylation level; blue, medium methylation level; gray, high methylation level). Gene lengths were normalized to 15 kbp, and 3 kbp was added upstream of the TSS and downstream from the TES.

Figure 3.

Genome-wide correlation between WGBS and ROM. (A) Comparison of coverage produced by both methods in a representative 500-kbp region from Chromosome 22 and the locus-specific wavelet decomposition correlation values in 10-kbp windows centered at each genomic position. (B) A zoom-in on a 40-kbp region from the same locus. In both cases, the WGBS data shown were inverted to represent nonmethylation levels and were smoothed using the wavelet transform. (C) Spearman's rank-order correlation values for genome-wide comparison of ROM and WGBS in different window sizes. (D) 2D density plot of nonmethylation levels in ROM in WGBS for regions identified as highly correlated by wavelet decomposition analysis.

Figure 4.

Comparison of long-range methylation profiles of single molecules. (A) Global view of Chromosome 1: 6,350,000–6,670,000 (bp). (B) Locations of CpG islands across the specified region. (C) Density of M.TaqI sites across the specified region in a 1.5-kb sliding window. (D) ROM fluorescence intensity profiles of three detected molecules aligned to the specified regions based on genetic labels (P-value <10⁻²⁰). (E) Gene body locations and corresponding HGNC gene symbols. Each gene is displayed as a purple arrow indicating gene orientation. (B–E) Black and red rectangles indicate methylated and nonmethylated gene promoters overlapping with CGIs, respectively.

Figure 5.

Copy number analysis. (A) The FSHD BAC model system contained a D4Z4 repeat array of unknown size (black triangles), genomic DNA upstream of the repeat array of unknown length (green), and the cloning vector (blue). (B) Box plots displaying the 25th, median, and 75th percentile (bottom, middle, and top of the box) of read coverage values along the repeat region (left box) and the nonrepetitive region (including the vector and the nonrepetitive genomic DNA; right box). The scale on the right is normalized to the median coverage along the nonrepetitive region. (C) Representative images of three intact model molecules labeled with Nb.BsmI (red dots) and stretched on modified glass surfaces. The labeling pattern can be aligned to the reference map presented below the images (expected labeling locations are shown in red). (D) Six hundred twenty-seven digital representations of labeled DNA molecules (yellow, blue dots represent detected labels) stretched and imaged in a nanochannel array chip. Consensus map of de novo assembly of the molecules is displayed in blue. The nonrepetitive regions are highlighted in gray (E) Comparative ROM profiles of nonmethylated and partially methylated BAC samples. Normalized averaged profiles of detected ROM labels are presented for the nonmethylated sample (green; 18,074 molecules) and partially methylated sample (blue; 9089 molecules). Each peak in the repetitive region corresponds to one repeat unit.

Figure 6.

Copy number and methylation analysis of the pathogenic contraction on Chromosome arm 4q of an FSHD patient compared with a healthy individual. (A) De novo assembly consensus optical maps of the two alleles, 4qA (red) and 4qB (gray) for both samples. Alignment to the in silico–generated map of this region (blue) is indicated by gray lines. Each repeat is represented by a black triangle. (B) Average profiles of detected ROM methylation labels for the two individuals (red indicates 4qA; gray, 4qB). Each peak in the plot corresponds to the average nonmethylation level of one repeat unit. (C) Table summarizing the number of detected repeat units in each allele. (D) Bar plot displaying the mean ROM nonmethylation value in the complete repeat array in both alleles for both samples.