henhaobio
Sign in

medRxivGenetic & Genomic MedicinePreprint

Population-scale detection of methylation outliers from long-read genome sequencing

Tanner D. JenseniD, Rhina Kaur, Devon E. BonneriD, Jon Nguyen, Chloe M. ReuteriD, Euan Angus AshleyiD, Jonathan Adam BernsteiniD, Matthew Thomas WheeleriD, Stephen B. MontgomeryiD

PDFDOICodeData
In-depth full-text digest available · Jump to digest ↓

Abstract

BackgroundAberrant DNA methylation can mediate the functional effects of rare genetic variation and contribute to imprinting disorders, repeat expansion diseases, and other pathogenic regulatory mechanisms. Long-read sequencing technologies now enable genome-wide detection of CpG methylation alongside genetic variation from a single assay. However, methods for systematic identification and interpretation of methylation outliers from long-read sequencing data remain limited. MethodsWe developed METAFORA, a computational workflow for detecting methylation outlier regions from PacBio and Oxford Nanopore long-read sequencing data. METAFORA constructs population-level methylation references, segments the genome into correlated CpG blocks, infers technical and biological sources of variation through hidden factor estimation, models uncertainty due to variable depth sequencing, and computes covariate-adjusted methylation outlier scores for individual samples. We applied METAFORA across large long-read sequencing cohorts and integrated methylation outliers with multi-omic data. METAFORA is implemented as a snakemake workflow available at https://github.com/tjense25/METAFORA. ResultsMETAFORA identified methylation outlier regions associated with rare structural variants, tandem repeat expansions, and imprinting abnormalities. We found outlier regions were enriched for molecular outliers across transcriptomic and chromatin accessibility datasets, supporting their functional relevance in gene regulation. In a representative case, METAFORA identified an imprinting defect affecting the GNAS locus associated with an STX16 deletion. ConclusionsMETAFORA enables scalable detection and interpretation of methylation outliers from long-read sequencing data and provides a framework for integrating epigenetic outliers with genomic and multi-omic analyses. These approaches may improve interpretation of rare regulatory variation and support discovery of clinically relevant epigenetic abnormalities in genomic medicine.

In-depth summary

Claude Opus 5.5 · xhigh effort · 8 Nov 2026 · from the full text of version 1 (11 Jun 2026) · CC BY 4.0

METAFORA calls depth-aware, covariate-adjusted methylation outlier regions from PacBio and Nanopore long reads; in 551 rare-disease genomes, outliers tracked structural variants, repeat expansions, and imprinting defects.

Array-based outlier calls profile only ∼80,000 CpGs and are made at single sites or promoters, other methods rely on fixed tiles or predefined annotations, episignature classifiers need pre-trained case-control models, and prior long-read methylation studies examined single individuals or families. METAFORA learns correlated CpG blocks from data, models depth uncertainty, and removes hidden technical and biological factors so PacBio and Oxford Nanopore genomes can be pooled into larger references.

Key results

  1. ShownAcross the multi-platform cohort, METAFORA called tens of thousands of outlier regions, and the number per genome varied by platform. Outliers were often haplotype-specific and enriched near rare structural variants, for example a TMEM97 promoter deletion on the same haplotype.

    32,808 outlier regions (26,816 blood; 5,992 fibroblast), 58.8% hypomethylated; median 17 large (>30 CpGs) autosomal outliers per sample (PacBio 21, ONT R9 13, ONT R10 6); 1,221 sex-chromosome outliers; 16.4% of PacBio outliers allelic; rare SVs 2.47 fold enriched near outlier genes (p < 2.2e-16). · Outlier burden per genome; Phased long-read genomes reveal…

  2. ShownOutliers sit preferentially in CpG islands, promoters, and constrained noncoding sequence. Nearby genes are more often outliers in chromatin accessibility, expression, and protein data, as seen in a FAM193B repeat-expansion carrier.

    Vs background: CpG islands OR=4.27, promoters OR=1.89, ABC elements OR=1.43; Gnocchi t value=57.1. Genes within 10Kb: ATAC outliers OR=6.41 (p=1.03e-28), RNA OR = 1.93 (p=3.7e-19), protein OR = 1.37 (p=0.011). FAM193B 675-copy CGG expansion: hypermethylation z-score = 6.68, ATAC log2 fold-change = −0.78, expression z-score = −4.04. · Outlier characterization reveals their…; Outlier characterization and enrichment analysis; Multi-omics enrichments

  3. ShownIntersecting outliers with 12 curated imprinting DMRs and filtering left three candidate cases. In one, a maternally phased STX16 deletion co-occurred with GNAS loss of methylation and GNAS overexpression, supporting pseudohypoparathyroidism type 1B.

    43 candidate outliers across 40 individuals before filtering. Proband (ONT R9, fibroblasts): 3 kb STX16 deletion of exons 4-6 plus GNAS outlier delta = −0.49, z-score = −3.84 on the same maternal haplotype; blood RNA-seq GNAS gene z-score=2.05, isoform ENST00000485673 z-score=2.93; affected sibling carries the deletion. · Case study: METAFORA outliers aid in rare…; STX16 deletion leads to loss of imprinting of…

  4. ShownSegment-level CpG methylation agreed closely across ONT R9, ONT R10, PacBio, and sequencing sites. The remaining platform differences tracked effective depth and error rate, which supports pooled references that include depth and platform covariates.

    Median Spearman r = 0.976 across blood samples; R10-R10 mean r = 0.984 vs R10-PacBio 0.974; same- vs cross-site difference ≈ 0.0024 (R10) and ≈ 0.0004 (PacBio); chrM error proxy R10 0.6%, PacBio ∼2%, R9 13.8%; 2.3% of segments differ by at least 0.1 median beta (ONT vs PacBio); depth-PC2 correlation 0.74. · Comparison of segment methylation profiles; Pairwise correlation across platforms and…; Sequencing depth is major source of variation in…

  5. ShownIn simulations, detection power rose with sequencing depth, effect size, and cohort size. False positives came mainly from noise, and size and effect filters removed nearly all of them.

    100 outliers per parameter set; at 5× depth only d ≥ 0.7 detectable; at 20× and above near-complete sensitivity for d ≥ 0.3; only 30× exceeded 50% power at d = 0.2; dropping calls < 10 CpGs or delta < 0.2 removed 99.95% of false positives; at d = 0.2, 84.4% detected with 100 samples vs 81% with 15. · Benchmarking METAFORA on simulated data; Sensitivity and power; False positives; Cohort Size

Methods, data and reproducibility

Snakemake workflow: CpG beta/depth profiles from long-read BAMs are summarized over fastseg-derived correlated CpG blocks; per-CpG depth-aware beta-distribution deviance scores are segmented into candidate regions, whose M-values are residualized on hidden-factor PCs, sex, and depth to call |z| > 3 outliers, then annotated with regulatory and haplotype features.

  • Stanford UDN and GREGoR R02 long-read genomes (ONT R9.4, ONT R10.4, PacBio Revio) from blood (n=503) and fibroblasts (n=48) across several sequencing sites · human · n = 551 genomes (77 ONT R9, 98 ONT R10, 376 PacBio) · dbGaP phs003047 (GREGoR); dbGaP phs001232 (UDN)
  • Simulated outliers spiked into chr1 population mean betas, varying median beta, delta, cohort size, CpG count, depth, and noise · n = 100 outliers per parameter combination
  • GREGoR ATAC-seq BAMs (ENCODE pipeline), jointly peak-called with MACS2, for chromatin accessibility outliers via EpiOut · human · n = 184 GREGoR participants (155 shared with methylation) · dbGaP phs003047
  • Whole-blood RNA-seq BAMs (STAR, hg38) used for expression outlier z-scores · human · n = 455 RNA-seq samples (108 shared with methylation) · dbGaP phs003047; phs001232
  • Olink HT plasma protein abundance (PCNormalizedNPX) run on 3 plates, used for protein outlier z-scores · human · n = 172 GREGoR participants (120 shared with methylation) · dbGaP phs003047
  • Curated imprinting-disorder DMRs from a recent review, with boundaries refined using allele-specific methylation in blood · human · n = 12 DMRs
  • Annotation tracks: ABC enhancer predictions (Nasser2021, lifted hg19 to hg38), gnomAD Gnocchi noncoding constraint, UCSC CpG islands, ENCODE cCREs, GENCODE v32 genes · human
Methods & tools
METAFORA (Snakemake workflow), modBam2Bed (ONT) and pb-CpG-tools (PacBio), fastseg change-point segmentation, edgeR correctBatchEffects, PCA with Gavish-Donoho threshold; PCAtools, limma lmFit, plyranges, Sniffles2 and jasmine (SV calling and merging), EpiOut with MACS2 peaks, STAR, RNASeQC, RSEM, Olink HT proteomics, IGV HTML reports
Software & environment
Snakemake workflow; GRCh38/hg38 with GENCODE v32 gene models; software package versions and containers are not stated.
Reproducibility gaps
Long-read and multi-omic data are accessible only through dbGaP (phs003047, phs001232); tool versions, runtime, and hardware requirements are not reported.
Limitations (6)
  • No genome-wide methylation truth set exists. Error was proxied by within-platform correlation and chrM CpG methylation, and the chrM estimate reflects error only where methylation is extremely low and may not transfer to nuclear CpGs. · stated by the authors
  • Power is limited for small effects (delta ≤0.2) and low depth. High-depth, high-noise data inflate false positives, which may partly explain the excess outlier burden in PacBio genomes. · stated by the authors
  • PacBio genomes were generated with the older SMRTbell prep kit 3.0 chemistry, so the finding that ONT R10 outperformed PacBio may not hold for newer chemistries such as SPRQ, which were not benchmarked. · stated by the authors
  • Haplotype-specific methylation is not explicitly modeled during calling. Haplotype coverage imbalance can mimic outliers at imprinted loci, and the vast majority of imprinting-DMR candidates showed such imbalance, which manual pileup inspection was used to filter out. · stated by the authors
  • No head-to-head comparison against array-based, fixed-tile, or other long-read outlier callers is reported. The simulations use a beta-plus-Gaussian generative model close to METAFORA's own beta test, which may make the sensitivity estimates look better than they are. · noted in this summary
  • Sample counts are inconsistent: the 77/98/376 platform split is described both as all 551 genomes (including 48 fibroblasts) and as whole-blood only, and the simulations are described as 100 outliers in one place and 10 replicates per setting in another. · noted in this summary

Read this if you call methylation from PacBio or Oxford Nanopore genomes in rare-disease cohorts and need a depth-aware, covariate-adjusted way to flag outlier regions, imprinting defects, and cis-acting variants.

Citations
0

Reproducibility ledger

Nobody has reported trying to reproduce this paper yet. Did you replicate an experiment or run the code? Your report saves the next lab a week.

Sign in to file a reproducibility report.

Related papers

Similarity Map
SelectedSame fieldCross-field
THISGigante '19Tian '25Costeira '26Hansen '25O’Neill '24Oliver '22
63% semantic matchGenomics
Scott Gigante, Quentin A. Gouil, Alexis Lucattini+8 more·Nucleic Acids Research·Feb 17, 2019· 137 citations