DNA methylation alterations in neonatal blood associated with pre-B ALL development

Epigenome-wide analysis in the prospective MoBa (Norway) and the retrospective CCLS (USA) discovery datasets (Supplementary Table 1) identified significant (FDR<0.05) differentially methylated regions (DMRs) in the blood of newborns who later developed pediatric pre-B ALL, relative to controls (Supplementary Fig. 1; Supplementary Tables 2A-B). DMRs in both studies were significantly enriched in CpG Shores, Promoters, First Exons, Exon-Intron/Intron-Exon boundaries and imprinted genes while being depleted in Open Sea, Shelf and intronic regions (p<0.05) (Supplementary Fig. 2). Both studies significantly converged (p<0.01) on several CpGs (Supplementary Tables 2A-B), among which 7 passed all filters (Fig. 1A), including effect size ≥ 3% (considered sufficiently high for validation by targeted sequencing) (Supplementary Fig. 1). All 7 CpGs mapped to the imprinted VTRNA2-1, which showed hypermethylation in nested cases relative to controls (Fig. 1A). VTRNA2-1 encompassed 9 additional CpGs, which were significant (FDR<0.05) in both studies’ Crude Models (Supplementary Table 2A). VTRNA2-1 methylation showed sex-dependent but ethnicity-independent alterations in the European and Hispanic descents (Fig. 1B), with the effects being observed in females, though requiring validation in larger sample sizes given that stratification by sex reduces power. In the replication phase, findings were reproducible in independent samples from MEDC (Australia) using a different technology, EpiTyper (Fig. 1C).

Fig. 1

Discovery, validation and functional analysis of VTRNA2-1 methylation in association with pediatric pre-B ALL development. A Upper section: Prioritized differentially methylated genes with at least one CpG with effect size ≥3% after DMR analysis for blood samples taken from newborns of either MoBa or CCLS. 7 CpG sites were significantly enriched (p = 2.2 x 10-16) between MoBa and CCLS relative to the total number of array CpGs analyzed (470,963 CpGs); all CpGs mapped to the same gene, which was also significantly enriched (p = 4.4 x 10-3) relative to the total number of genes in the human genome (21,306 genes) (Fisher’s Exact Test). Lower section: The 7 significant CpGs within the DMR of VTRNA2-1 are symbolized in CCLS and MoBa by circles of varying sizes and colors, representing the effect sizes and directions of effect, respectively, as per the figure legend. The 7 CpGs are arranged in order according to their genomic position. The direction of effect is reported for the pre-B ALL nested cases relative to the controls: hypermethylation (Hypermeth) or hypomethylation (Hypometh). B VTRNA2-1 differential methylation in nested cases and controls was stratified by subject sex and ethnicity in MoBa and CCLS cohorts. Data points represent average methylation values at each CpG site, and the ribbons denote the 95% confidence intervals. CpG HM450 IDs are shown on the x-axes. In addition to the CpGs (in red) identified in the Adjusted Models in both MoBa and CCLS, we also show (in black) the additional CpGs identified in the Crude Models in both MoBa and CCLS (detailed in Supplementary Fig. 7 and Supplementary Table 2). C Validation, based on profiling of VTRNA2-1 methylation using EpiTyper, which is sequencing- rather than array-based, applied to an independent set of biological samples from MEDC. Data points represent average methylation values at each CpG site, and the error bars denote the 95% confidence intervals. The p-values indicate the statistical significance across the whole DMR region and were calculated by inverse variance based meta-analysis using METAL software. The DMR profiled by EpiTyper partially overlaps with that by HM450; specifically, CpGs 10 and 11 in (C) are identical to the last two CpGs in (B), cg16615357 and cg18797653, respectively. CpG1-2 and CpG3-4 each represents an average methylation value of two adjacent CpGs, as detected by EpiTyper. The genomic coordinates of the CpG ID numbers are detailed in Supplementary Fig. 8. D Box plots showing the methylation distribution of VTRNA2-1 across a panel of human tissue types using data extracted from the EWAS Open Platform [10]. The box plots encompass the first quartile (bottom border), the median (middle line), the fourth quartile (upper border) and the extreme values (dots). No statistically significant differences (p>0.05; Mann-Whitney test) were detected in VTRNA2-1 mean methylation between the target bone marrow and surrogate cord blood tissues. The sample sizes (N) are indicated for each tissue type. (n=5,023 tissues; 30 types) E Pearson correlation of VTRNA2-1 expression with the methylation of its CpGs in a panel of cancer tissues extracted from the MEXPRESS database [11] (n=2,273 tissues; 25 types). Cancer types and sample sizes are as follows: kidney renal papillary cell carcinoma (KIRP, N = 140), rectum adenocarcinoma (READ, N = 28), pheochromocytoma and paraganglioma (PCPG, N = 66), skin cutaneous carcinoma (SKCM, N = 74), testicular germ cell tumor (TGCT, N = 47), uveal melanoma (UVM, N = 28), thyroid carcinoma (THCA, N = 162) , kidney renal clear cell carcinoma (KIRC, N = 136), breast invasive carcinoma (BRCA, N = 106), pancreatic adenocarcinoma (PAAD, N = 67), colon adenocarcinoma (COAD, N = 95), prostate adenocarcinoma (PRAD, N = 75), liver hepatocellular carcinoma (LIHC, N = 90), bladder urothelial (BLCA, N = 106), uterine corpus endometrial carcinoma (UCEC, N = 121), head and neck squamous cell carcinoma (HNSC, N = 103), lung adenocarcinoma (LUAD, N = 74), mesothelium (MESO, N = 42), lung squamous cell carcinoma (LUSC, N = 44), glioblastoma multiforme (GBM, N = 44), sarcoma (SARC, N = 82), cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC, N = 127), brain lower grade glioma (LGG, N = 96), stomach adenocarcinoma (STAD, N = 217) and esophageal carcinoma (ESCA, N = 103). The asterisk (*) mark significant correlation after adjustment for multiple testing (FDR < 0.05). One CpG (cg11978884) was omitted from the analysis because it had no methylation values. F VTRNA2-1 methylation in MoBa paired samples over time. None of the VTRNA2-1 CpGs were significantly (p>0.05) differentially methylated in cord blood collected from the control subjects at age 0 (blue) versus paired peripheral blood collected from the same controls at age 3 (orange) years (Wilcoxon test). Methylation values at birth from nested unpaired controls (green) and cases (red) are shown as a reference. The y-axes represent the methylation (beta) values, and p values are reported for each CpG. In E-F, the orange rectangles represent CpGs common to Adjusted Models of both MoBa and CCLS. The remaining CpGs are those identified in the Crude Models of both datasets

Functional analysis of VTRNA2-1

As expected for an imprinted gene, VTRNA2-1 methylation levels were centered around 50±10% and were similar across various tissue types, including bone marrow (pre-B ALL target tissue) and cord blood (surrogate tissue) (p>0.05), with the exception for placenta (being extra-embryonal) and sperm (being hypomethylated given that VTRNA2-1 is maternally imprinted) (Fig. 1D). In various tissue types, VTRNA2-1 methylation negatively correlated with its gene transcription (FDR<0.05), including in pediatric pre-B ALL tissues (Supplementary Fig. 3 and Fig. 1E).

Longitudinal analysis in controls and cases pre-diagnosis

Assessment of VTRNA2-1 methylation stability over time showed that methylation at all 16 CpGs was similar across birth years in cord blood samples collected from controls and pre-diagnostic cases from 2000-2008 in MoBa (Supplementary Fig. 4A), and this was replicated in neonatal blood spots collected from 1984-1994 in an independent population, UKCS (UK) (Supplementary Fig. 4B). Moreover, there was no significant difference (p>0.05) in VTRNA2-1 methylation levels between cord blood (i.e. at birth) and peripheral blood collected from the same individuals (controls) at age three years (Fig. 1F), which is the peak incidence age of pre-B ALL [12]. This highlights the stability of VTRNA2-1 methylation over critical time windows and across matched cord and peripheral blood tissues, reinforcing its observed methylation stability in a panel of human tissues (Fig. 1D).

Longitudinal analysis post-diagnosis

In line with the VTRNA2-1 hypermethylation observed at birth in nested cases versus controls, VTRNA2-1 was significantly hypermethylated in pediatric pre-B ALL tissues at diagnosis compared to controls, regardless of diagnostic matrix (surrogate=blood or target=bone marrow, Fig. 2A). This was based on the NOPHO (Nordic countries) cohort and validated in the QcALL cohort (Canada) (Supplementary Fig. 5). At remission, VTRNA2-1 methylation levels were reset to normal, and, at relapse, they re-increased to above control levels (Fig. 2A). This trend was validated in subjects matched at diagnosis and remission in QcALL (Fig. 2B). Overall, these results suggest VTRNA2-1 methylation as a marker of pre-B ALL prognosis, including leukemic state.

Focusing on individual-level data, we observed two QcALL patient clusters: C1 exhibited low methylation levels (10±10%), stable from diagnosis to remission (p>0.05), and C2 exhibited at diagnosis higher methylation levels (~20-100%), which converged to 50±10% at remission (p<0.05) (Fig. 2B and Supplementary Fig. 6A-B). The C1 and C2 patterns were also observed in blood at birth (Supplementary Fig. 6C) and in NOPHO (Supplementary Fig. 6D) and are in line with recent observations showing that VTRNA2-1 is imprinted in ~75% of individuals [13], as in C2, and non-imprinted in the remaining portion (C1).

To test whether VTRNA2-1 methylation affects clinical outcomes, ten-year follow-up data was used in NOPHO (n=598, including 134 relapse and 79 death events; Supplementary Table 3). VTRNA2-1 methylation at two CpG sites was significantly inversely associated with overall survival (hazard ratio>1; p<0.05), after adjusting for sex, age and prognosis risk groups (Fig. 2C), hence, reinforcing the prognostic potential of this gene. No significant difference (p>0.05) on patient overall or relapse-free survival was detected between C1 and C2 (Supplementary Fig. 6E-F). Also, no significant associations (p>0.05) were observed between VTRNA2-1 methylation and relapse-free survival (Fig. 2D), although relapse events were more frequent (exhibiting higher statistical power) than death events, suggesting that VTRNA2-1 methylation may likely associate with patient survival more than relapse.

Fig. 2

Longitudinal analysis of VTRNA2-1 methylation post-diagnosis and its hypothesized role in pre-B ALL development. A Methylation of VTRNA2-1 CpGs in peripheral blood and bone marrow of cases versus controls at diagnosis, remission and relapse in NOPHO. In purple: Methylation of VTRNA2-1 CpGs in peripheral blood samples from sorted B-cells of normal subjects (N=26) and from pediatric pre-B ALL patients collected at diagnosis (n=74) from NOPHO. In green: Methylation of VTRNA2-1 CpGs in sorted B-cells (N=26) from bone marrow of fetuses (N=8) and in bone marrow samples from cases of pediatric pre-B ALL collected at diagnosis (n=535), remission (n=82) and relapse (n=32) from NOPHO. Whiskers represent the minimum and the maximum, while the top, the bottom, and the band in the box represent the first and third quartile and the median respectively. Significant differences between methylation of normal and tumor samples are marked for each CpG with an asterisk (Wilcoxon test). B Methylation of VTRNA2-1 CpGs in 46 pediatric pre-B ALL samples collected at diagnosis (red) and remission (blue) from the same patients in QcALL. Significant differences between methylation at diagnosis and remission are marked for each CpG with an asterisk (Wilcoxon test). The data are represented in the form of a dot plot to better visualize the paired samples (a line links each pair). Red and blue box plots are also shown for each time point (diagnosis and remission, respectively), representing the first quartile (bottom border), the median (middle line) and the fourth quartile (upper border) for each condition. C and D Methylation of VTRNA2-1 CpGs in relation to overall and relapse-free survival, respectively, represented by hazard ratios. In NOPHO, 598 pre-B ALL patients were followed up for ten years or more. VTRNA2-1 methylation at two CpG sites significantly affected overall survival (denoted by *, Wald test), after adjusting for patient sex, age and risk groups using a Multivariate Cox Regression model. Risk group variables also affected overall and relapse-free survival (denoted by ** or ***, Wald test). HR: high risk, IR: intermediate risk and SR: standard risk. * denotes p<0.05, ** denotes p<0.01, *** denotes p<0.001. In A-D, the orange rectangles represent CpGs common to Adjusted Models of both MoBa and CCLS. The remaining CpGs are those identified in the Crude Models of both datasets. E The tumor surveillance model offering a biologically plausible mechanism of VTRNA2-1 in pediatric pre-B ALL development. The basal methylation and expression level of VTRNA2-1 determines the degree of gradients (narrow: RIGHT versus wide: LEFT), which is important to shift the balance from cell survival (RIGHT) to cell death (LEFT) via PKR activation. Graphic icons used to construct the figure were retrieved from thenounproject.com. F Summary of the study’s time points, sample types, and VTRNA2-1 results

Hypothetical model of VTRNA2-1 role pre- and post-diagnosis

VTRNA2-1 is a 100-nucleotide non-coding RNA with central roles in multiple cancer types based on cell, animal and clinical models and is a known regulator of Protein Kinase R (PKR)-mediated cell death (Supplementary Table 4). A tumor surveillance model for eliminating pre-cancerous cells has been proposed requiring VTRNA2-1 hypermethylation and the ‘drop’ in its expression as a critical event for PKR activation [14] (Fig. 2E). The effectiveness of this mechanism is speculated to be weaker in individuals in which VTRNA2-1 expression is already low (i.e. hypermethylated) at birth, prohibiting any possible drop of VTRNA2-1 levels and subsequent cell death. This can lead since birth to the accumulation of precancerous cells that can malignantly transform over time (Fig. 2E). VTRNA2-1 methylation could also affect PKR-mediated immune regulation [15], which can serve as an additional hit to activate pre-leukemic clones to progress into malignancy and/or could enable existing tumor cells to evade immune attack (hence, worsening prognosis).