Эпигенетические достижения в педиатрическом многоразовом склерозе: ДНК-метилирование примечательно участие Эпштейна-Барра заражённых вирусом и принимающих иммунных сетей.

Introduction

Symptoms and disease activity in multiple sclerosis (MS) vary greatly from person to person. Despite recent advances in treatment, there is no cure for MS, and it remains challenging to disentangle the complex relationships between genetic susceptibility and environmental factors contributing to disease onset and activity over time. Pediatric-onset MS (PedMS), also known as early-onset MS or juvenile MS, is defined as MS with an onset before the age of 18 years. PedMS has distinctive features, and the disease course is different than in adults, 1,2 characterized by a higher burden of relapses and inflammation compared with adult MS 3,4 together with a higher capability to compensate brain damage, as suggested by the longer time to reach mild and severe disability. 5 Furthermore, the current understanding of PedMS suggests a narrow window of onset with unique cognitive decline. 6,7 In general, the relative contributions and timing of environmental exposures that promote MS onset are not clearly defined, as well as the knowledge about the mechanisms that link this exposure to genetic predisposition. Various early-life factors near gestation as well as exposures during pregnancy and breastfeeding 8 and exposure to various environmental factors in early phases of life including common viral antigens 9 could provide a unique opportunity to study their role in MS pathogenesis. 10 By focusing on PedMS, where the time lag between potential exposures and disease onset is shorter than in adult-onset MS (AOMS), with a disease duration less than 4 years as in our study, we can investigate these epigenetic mechanisms closer to their origin and before extensive confounding exerted by cumulative environmental exposures, prolonged disease progression, therapies, and aging in general.

Epigenome-wide association studies (EWASs) have consistently linked DNA methylation changes to MS, identifying among other signals a differential methylation at the human leukocyte antigen – DR beta 1 (HLA-DRB1)*15:01 locus as a hallmark of the disease. 11-13 Recent evidence suggests that these epigenetic signatures can occur independently of known genetic risk loci, offering superior diagnostic accuracy (area under the curve = 0.85) compared with genetic markers alone (area under the curve = 0.72). 13 These shifts seem to be cell specific, predominantly affecting B cells and monocytes, and extend beyond peripheral circulation; for instance, pathology-free brain regions in MS show distinct methylation patterns in genes regulating oligodendrocyte survival ( BC2L2 , NDRG1 ) and proteolytic processing ( LGMN , CTSZ ). 14 Furthermore, a complex genetic-epigenetic-transcriptional interplay at the 1q21.1 locus has been proposed to drive pathogenesis in primary progressive MS. 15

Despite these advances in adult cohorts, epigenetic data on PedMS remain scarce. Early investigations have begun to unravel this complexity, showing, for example, that childhood adversity may differentially shape the epigenome depending on the age of MS onset 16 or that PedMS patients exhibit accelerated epigenetic aging compared with pediatric healthy controls (HCs). 17 Crucially, recent evidence suggests that MS onset risk is determined by a concerted action of EBV infection, vitamin D levels, and genetic susceptibility, all converging through shared epigenetic pathways. 18 The aim of this work is to leverage the multicentric data and collaborative effort of the PEDIGREE consortium, 9 a collaboration of 29 Italian MS centers focused on identifying genetic and nongenetic risk factors as well as their interaction in PedMS, to identify critical epigenetic signatures associated with PedMS, with a focus on differentially methylated regions (DMRs). DMRs are clusters of CpG sites that exhibit coordinated and statistically significant differences in DNA methylation, which have been shown to play a role in many different diseases including MS because of their association with epigenetic control of gene regulation. 19

Methods

Study Population

This is a multicenter retrospective study, including overall 175 subjects (122 PedMS patients and 53 HC). The inclusion criteria for PedMS patients were (1) onset of MS before age 18 years, (2) diagnosis of MS according to Krupp's criteria, 20 and (3) disease duration of less than 4 years. Patients were excluded if they had an alternative diagnosis confirmed through appropriate testing (specifically, testing for anti-myelin oligodendrocyte glycoprotein (MOG) and anti-neuromyelitis optica antibodies using a live cell-based assay [CBA]). HCs were either healthy individuals or patients hospitalized for trauma, minor surgical procedures, or mild disorders who were not affected by systemic autoimmune disorders, including psoriasis, celiac disease, autoimmune thyroiditis, Systemic Lupus Erythematosus, Sjogren syndrome, or type 1 diabetes.

DNA Extraction and Bisulphite Treatment of DNA

Genomic DNA (gDNA) was extracted from peripheral blood collected by venipuncture from 175 subjects, drawn into tubes containing anticoagulant (EDTA), and stored at −80°C for further processing. Quantification was performed by the UV-Vis NanoDrop 2000 (Thermo Scientific, Waltham, MA) spectrophotometer or Quant-iT PicoGreen dsDNA Assay Kit (ThermoFisher Scientific Waltham, MA), Qubit dsDNA HS assay Kit, and Qubit dsDNA BR Assay Kit (Life Technologies, Grand Island, NY). In total, 500 ng of DNA from whole blood were converted by using bisulphite conversion technique (EZ DNA Methylation Kit) that converts unmethylated cytosines into uracils, while methylated cytosines remain unchanged.

Genome-Wide DNA Methylation Analysis

Bisulphite-converted DNA was hybridized to the Infinium Methylation EPIC Array v2.0 (Illumina, San Diego), covering more than 930,000 methylation sites at the genome-wide level (reference genome hg38), according to the manufacturer's protocol. Signal intensity has been measured with the Illumina NextSeq 550 sequencer, and quality controls were supported by the GenomeStudio Methylation Module. In silico preprocessing of data and quality control was performed on raw methylation data applying standard Chip Analysis Methylation Pipeline (ChAMP; version 2.21.1). 21 Specifically, we removed probes with a detection p value >0.01 (n = 14,263), probes with a bead count <3 in at least 5% of samples (n = 5,711), non-CpG probes (n = 3,371), and probes overlapping with single nucleotide polymorphism sites (n = 119,867). Similar to other microarray experiments, methylation data are susceptible to various technical artifacts; thus, batch effect presence was verified within the exploratory analysis, checking for the presence of an association between the most important principal components and the Sentrix identifiers. After performing β-mixture quantile normalization (BMIQ) to correct for type II probe bias, the ComBat method 22 implemented in ChAMP was used to adjust for batch effect at the slide level. The obtained β values, reflecting CpG methylation levels as a ratio of normalized signal intensities, were used for further identification of differentially methylated probes (DMPs).

Cell Types Estimation

To account for the inherent cellular heterogeneity of whole blood samples, immune cell proportions (CD4 + T cells, CD8 + T cells, B cells, natural killer [NK] cells, monocytes, and granulocytes) were estimated using the DNAmAge calculator, a robust epigenetic deconvolution tool. 23 Sex was also bioinformatically inferred with the same tool based on the methylation levels of X-chromosomal markers. The estimated cell fractions were then statistically evaluated across different treatment subgroups using the Kruskal-Wallis test, followed by the Dunn post hoc test with Benjamini-Hochberg correction for multiple comparisons. This statistical step was crucial to identify whether specific disease-modifying treatments (DMTs) significantly altered the circulating immune landscape, potentially introducing systematic bias in the methylation signal.

Differential Methylation Analysis and Functional Annotation

Analyses were performed fitting multiple linear regression model on methylation M value matrix using the limma method 24 implemented in RnBeads package (version 2.21.2). 25 M values, which better satisfy the homoscedasticity assumption of statistical models commonly used in microarray analysis, were calculated as the log2 ratio of the β values. 26 Because gDNA samples originated from whole blood to contrast possible confounding effects caused by differences in blood cell composition between and within groups, including those caused by the concomitant use of DMT, we considered cellular heterogeneity using the method by Houseman et al. 27 already described. Thus, the linear model was fitted using as covariates the first 3 principal components of the matrix that includes immune cell estimated proportions, which captured most of the cellular variance (86%) and sex, to prevent multicollinearity.

DMPs were filtered using automatically selected rank cutoffs as implemented in RnBeads. The rank is calculated based on 3 different metrics: (1) the difference in mean methylation levels of the 2 groups, (2) the quotient in mean methylation, and (3) the statistical test performed with limma assessing if methylation values in the 2 groups originate from distinct distributions. A combined rank is computed as the maximum (i.e., worst) rank among the 3 ranks. The smaller the combined rank for a site, the more evidence for differential methylation it exhibits. Regional DNA methylation was analyzed using the DMRcate method, 28 updated for the EPICv2 platform. To ensure probe specificity, the workflow included the remapping of 4,130 cross-hybridizing probes to their most likely off-target locations and the filtering of replicate probes mapping to the same CpG site. DMRs were then detected by fitting the same linear model used in the prior analysis, setting the minimum number of consecutive DMPs constituting a DMR as 7 and a minimum gap between 2 different DMRs (defined as lambda value) of 1,000 base-pairs. DMRs were selected according to a false discovery rate (FDR)–adjusted p value of 0.05. A threshold in mean absolute delta-beta (mean Δβ) of 0.02 was used to discern whether a detected signal is caused by true biological effect or technical variation. 29 DMRs were ranked according to minimum FDR of the smoothed estimate (min FDR), where FDR smoothing results in increased power and cleaner spatial separation of signals from noise. Subgrouping strategies for sensitivity analyses were further defined based on treatment-induced cellular distribution, whereas an internal validation was performed based on a 1/3–2/3 partitioning criterion (Results). To further characterize the regulatory potential of the most relevant epigenetic signatures, top-ranked DMRs were integrated with functional genomic data from the ENCODE Project (Phase 4). This included H3K27ac histone acetylation tracks generated by the Broad Institute (Experiment Accession: ENCSR000AKP, File Accession: ENCFF381NDD) using the GM12878 lymphoblastoid cell line as a gold-standard reference for B-lymphocyte immune profiling. 30 In addition, we used the multicell registry of candidate cis-regulatory elements (cCREs, v3), accessed using the SCREEN portal. 30,31 This multilayered annotation, aligned to the GRCh38 genome assembly, allowed us to contextualize the identified differential methylation within established active enhancer or promoter frameworks, providing additional evidence for their functional profiling in the context of pediatric neuroinflammation.

Network Analysis

Using the NDEx tool, 32 we performed a subnetwork detection analysis, projecting genes overlapping with DMRs onto STRING v12.0 reference interactome. 33 We retained only those links with very high interaction evidence (score >0.95) in 1 of the 3 evidence domains: (1) protein-protein interaction derived from multiple interactomes such as IntAct, BioGrid, MINT, and others; (2) coexpression, which leverages gene expression data from multiple sources; and (3) databases, which collects evidence of interaction from curated pathway resources. From this filtered interactome (10,430 nodes and 58,559 edges), we explored the network using the 1-step neighborhood method to identify modules defined by the group of nodes connected to the query term(s) and all edges between these nodes. Centrality measures of the nodes were estimated with CentiScaPe tool, 34 and gProfiler 35 was used to perform pathway enrichment analysis on full set of network nodes.

Standard Protocol Approvals, Registrations, and Patient Consents

This study was conducted in accordance with the principles of the Declaration of Helsinki. The study protocol was approved by the Comitato Etico dell’Insubria (protocol: PEDIGREE study) and by the local Ethical Committees of all participating centers. Written informed consent was obtained from the patients or their legal guardians, and all data were anonymized.

Data Availability

Deidentified individual participant data (including epigenetic methylation profiles and clinical characteristics) that underlie the results reported in this article will be shared for scientific purposes only. Related documents, including the study protocol and statistical analysis plan, will also be made available. Proposals should be directed to the corresponding author to gain access, and data requestors will need to sign a data access agreement. Data sharing is subject to compliance with the general data protection regulation and the specific restrictions of the informed consent forms signed by participants.

Results

Baseline characteristics of the study cohort are presented in Table 1 , including 122 PedMS and 53 HC. The study design aimed for balanced recruitment across age and sex; however, the exclusion of samples because of technical limitations (low DNA yield) and missing associated demographic and clinical data, together with the sex ratio disparity characteristic of PedMS, 36 led to an imperfect matching by age and sex between cases and controls. The median age was 16 years in cases and 14 in controls ( p = 0.008). The male:female ratio was 1:2.8 in cases and 1:1.2: in controls ( p = 0.02). The median disease duration was 24 months (IQR: 12–42). Seventy-eight percent of cases were under DMT, of which 28 treated with low-efficacy treatment agents (LETA) (23% of the total, of which 10 with glatiramer acetate, 10 with beta-interferons and 8 with dimethylfumarate), 68 patients were treated with high-efficacy treatment agents (HETA) (56% of the total, of which 44 with natalizumab, 23 with fingolimod, and 1 with rituximab) and 5 patients (4.1% of the total) with unspecified treatments. eFigure 1 presents the data specifically for B cells. Detailed results of the statistical evaluation for all estimated cell proportions are presented in eTable 1, showing that the LETA and other DMT-treated subgroup (n = 33), along with the treatment-naïve patients (n = 21), exhibited a cellular profile indistinguishable from HC (e.g., B cells p-adj = 0.61; CD4 + T cells p-adj = 0.233). By contrast, significant cellular redistribution was confined to the subgroup of HETA-treated patients (n = 68). Detailed clinical and demographic characteristics for each subgroup, including disease duration and treatment distribution, are summarized in eTable 2.

As an outcome of the filtering standard procedures, 175 samples and 771,881 annotated probes were retained for downstream data processing from a starting pool of 937,055 probes and 175 samples. Association analysis revealed 7,968 DMPs (eTable 3, overall association at single site level is shown in Figure 1 ) and 55 DMRs, which were annotated to 70 unique genes and 30 long noncoding RNAs based on genomic proximity and overlap with promoter or intragenic regions. Twenty-eight of 55 top ranked identified DMRs overlap with promoter regions of 1 or different protein coding genes, defined as the range that spans 1,000 bp upstream and 500 bp downstream of the transcription start site.

To validate the DMR detection results and control for potential confounders, 2 sensitivity analyses and an internal validation were conducted. The sensitivity analyses aimed to decouple disease-specific epigenetic signatures from treatment-induced effects, focusing on cohorts where the pharmacologic impact on the methylome was expected to be minimal. The first included only treatment-naïve PedMS patients (n = 21) compared with HC. The second, comprising naïve, LETA, and other DMT-treated patients (n = 54), was directly informed by cell-type estimation results: because patients under HETA exhibited significant shifts in immune cell composition, this broader subset allowed us to validate the DMRs in the absence of potent immunomodulatory interference. The DMR identification results for these 2 sensitivity analyses are presented in eTables 4 and 5. The internal validation was conducted by partitioning the total cohort into a discovery sample (n = 59) and a validation sample (n = 116) using a double-stratified sampling approach implemented in caret R package, while a more permissive detection threshold for individual CpG sites (FDR <0.1) in DMRcate was applied to compensate for the reduced statistical power inherent in the partitioned subgroups. Comprehensive data regarding the identified DMRs are presented in eTable 6. The most robust signal was identified on chromosome 1, mapping to CHI3L2 and DENND2D (min FDR = 8.005*10 −39 , mean Δβ = −0.024) ( Figure 2A ). This signal persisted throughout the sensitivity and validation analyses. Similarly, the DMR mapping to LPCAT1 (min FDR = 3.318*10 −30 , mean Δβ = 0.023) was independently confirmed in the naïve group (min FDR = 3.988*10 −7 ), in naïve + LETA and other DMT (min FDR = 2.461*10 −11 ), and in the internal validation (min FDR = 5.42*10 −8 and 1.38*10 −15 ) (eFigure 2). Other highly resilient signals include the locus overlapping with CD19 ( Figure 2B ) (min FDR = 2.072*10 −29 , mean Δβ = 0.025), which persisted in the naïve + LETA and other DMT sensitivity analysis (min FDR = 3.85*10 −5 , mean Δβ = 0.018). Furthermore, the DMR mapping to PARP1 (min FDR = 7.314*10 −28 , mean Δβ = 0.027) persisted in the second sensitivity analysis (min FDR = 1.21*10 −6 , Figure 3A ). Another highly significant and resilient epigenetic signal was identified on chromosome 22, mapping to the CPT1B/CHKB locus (min FDR = 1.67*10 −26 , mean Δβ = 0.029) (eFigure 3). This DMR, which encompasses a complex transcriptional unit including the CHKB-CPT1B read-through transcript, was replicated in the naïve + LETA and other DMT sensitivity analysis (min FDR = 9.93*10 −10 ) and confirmed in the internal validation cohort (min FDR = 6.9*10 −29 and 1.81*10 −5 ). The DMR that maps in the AIRN , IGF2R locus (min FDR = 1.5*10 −4 , mean Δβ = 0.026, Figure 3B ) also replicated in the 2 frameworks (min FDR = 1.9*10 −13 and 5.1*10 −12 and 1.3*10 −7 ). Complete list of the top 55 DMRs identified in their genomic context is presented in Table 2 ; min FDR values and mean Δβ are additionally reported for those signals that successfully replicated across the 3 analytical frameworks.

Interactome analysis mapped 33 of the 70 DMR-associated coding genes, resulting in a subnetwork comprising 165 nodes and 680 edges ( Figure 4 ). Of the initial seed genes, 19 were integrated into this main subnetwork, while the remaining 14 mapped to small, disconnected components and were excluded. Pathway enrichment analysis performed on this subnetwork identified a strong dominance of immune-related and stress-response pathways, revealing 43 high significant terms (adjusted p < 1*10 −18 , Table 3 ; eTable 7 for a complete list of significant terms), which were subsequently evaluated within the sensitivity framework of the naïve + LETA and other DMT-treated cohort to ensure a sufficient density of DMRs. The network analysis performed on this subset of the entire cohort revealed a subnetwork of 87 nodes and 320 edges (eFigure 4), and a pathway analysis was performed on this network of genes (eTable 8). The pathway of Epstein-Barr virus infection emerged as the top-ranked term in the primary analysis (adj-p = 3.23*10 −40 ), and it was confirmed in the sensitivity analysis (adj-p = 5.22*10 −5 ). Other pathways which were replicated refer to immune system process, the positive regulation of biological process which shows the highest number of intersecting genes (n = 111), and the cell surface receptor signaling pathway.

Discussion

Standard EWAS often lack the statistical power to detect the small effect sizes characteristic of epigenetic research, particularly in complex traits studies. Regional DNA methylation profiles have shown greater predictive power for disease status compared with the methylation status of individual CpG sites. 37 Therefore, the identification of DMRs has emerged as a strategy to increase statistical power and to reduce the multiple testing burden, suggesting a more complex epigenetic regulation and potentially a greater biological relevance than changes in individual sites. DMRs may indeed occur throughout the genome but have been identified particularly around the promoter regions of genes, within the body of genes, and at intergenic regulatory regions. Indeed, 55 top-ranked DMRs identified in this study overlapped with a total of 90 genes, 28 of which mapped to gene promoters and for most of them an association with MS has already been described. Although 78.8% of our cohort was under DMT, cell-type deconvolution confirmed that immune profiles in patients treated with LETA and other agents did not deviate significantly from controls. This lack of cellular redistribution, particularly when contrasted with HETA effects, justified focusing our sensitivity analyses on these subgroups to isolate intrinsic disease signals. A third step of internal validation was performed in the absence of an external data replication, which was difficult to find because of the rarity of epigenetics data on PedMS patients. We chose to study epigenetic modulation via CpG methylation because it is unique in reflecting environmental exposures and lifestyle factors. The top-ranked DMR, consistently identified across both sensitivity analyses and internal validation, overlaps with the first exon and the promoter region of CHI3L2 . The gene is structurally similar to CHI3L1 , which is a validated biomarker for MS activity but with a different function. 38 Elevated CHI3L1/2 expression levels are typically observed in the CSF of AOMS, 39 also being a CSF marker of long-term disability progression. 40 Our study identified a significant hypermethylation of the CHI3L2 promoter in the peripheral blood of PedMS patients. This epigenetic silencing, localized within a region encompassing a cluster of multiple active candidate cis-regulatory elements (cCREs) and marked by a prominent H3K27ac peak ( Figure 3A ), suggests a systemic downregulation of the protein expression. This divergence may reflect tissue-specific regulatory mechanisms or a compartmentalized immune response. Another highly consistent epigenetic signal maps to LPCAT1 (Lysophosphatidylcholine Acyltransferase 1) (eFigure 2). LPCAT1 governs the Land's cycle conversion of lysophosphatidylcholine into phosphatidylcholine, a process critical for myelin integrity and immune membrane modulation, while also providing a regulatory link to the de novo Kennedy pathway. 41 Parallel dysregulation at the CPT1B/CHKB locus (eFigure 3), involved in mitochondrial fatty acid transport and phospholipid biosynthesis, supports this connection, suggesting that altered lipid homeostasis and metabolic flux may contribute to the characteristic immune signature of PedMS. Another resilient signal maps to the PARP1 gene, involved in DNA repair, neuroinflammation, and oxidative stress. This hypomethylated DMR is located within an intronic region characterized by a prominent H3K27ac enrichment peak, overlapping with multiple candidate cis-regulatory elements (cCREs) ( Figure 3A ) suggesting an increased transcriptional accessibility of this regulatory locus in PedMS patients. Although the therapeutic potential of PARP1 inhibition in MS remains a subject of debate, 42 our findings highlight its epigenetic dysregulation as consistent with the PARP1 up-regulation previously reported in MS immune cells. 43 Other signals of DMR point to the role of B cells in PedMS including CD19 , BLK , and BLNK , but only CD19 ( Figure 2B ) replicates in sensitivity analyses, suggesting that they likely reflect treatment-induced shifts in B-lymphocyte proportions rather than a primary disease-specific epigenetic marker. The CD19 gene emerged as upregulated in a group of PedMS vs AOMS patients, 44 compatibly with the fact that promoter hypomethylation detected in PedMS subjects affects gene expression positively, as well as BLNK which was shown to be upregulated in the same study. An additional finding of the study refers to the involvement of specific loci associated with genomic imprinting, a mechanism whose role in MS is increasingly recognized. 19 In particular, IGF2R ( Figure 3B ) emerged among the most resilient signals, maintaining its statistical significance in the naïve + LETA subset and in the internal validation. The identified DMR overlaps with 2 imprinting control regions (ICRs) previously characterized. 45 Specifically, IGF2R has been implicated in modulating the internalization of cytotoxic granzyme B in MS models, suggesting its role in protecting neurons from inflammatory-mediated damage. 46 Similarly, the stability of ZFP57 , a master regulator essential for the maintenance of imprinting at multiple loci, 47 suggests a systemic disturbance in the epigenetic machinery rather than a localized stochastic event.

Regarding network-integrated pathway analysis, EBV infection emerged as the top-ranked term in the primary analysis and it was confirmed in the sensitivity analysis. Notably, among the subjects for whom serologic data were available, 96.1% of patients (49/51) and 83.3% of controls (25/30) tested positive for EBV (VCA IgG and EBNA-1 IgG). The EBV pathway signal arises from network enrichment of DMR-associated genes independently of serologic status, and the 2 EBV-seronegative PedMS subjects displayed methylation profiles indistinguishable from their diagnostic group (eFigures 1 and 4), consistent with a multihit model of MS pathogenesis. In contrast to specific viral infection terms, broader regulatory categories and response to stimulus pathways showed greater stability across different cohort subsets. This functional robustness is epitomized by the central roles of PARP1 and IGF2R within our interaction network, which is maintained also in the sensitivity framework ( Figure 4 and eFigure 5). PARP1 , identified as the most central hub according to the betweenness centrality metric (0.56), links the EBV-related signaling clusters (comprising GRB2 , CBL , and CD79B ) to broader mechanisms of DNA repair and genotoxic stress management (interacting with XRCC5/6 and PRKDC ). The centrality of IGF2R underscores the involvement of logistical pathways, such as ER-Phagosome ( p = 3.05 * 10 −19 ) and vesicular trafficking (through GGA2 and PLIN3 ), as stable regulatory hubs in PedMS. By contrast, the isolation of top signals such as CHI3L2 and the LPCAT1/CHKB axis, which remained unconnected in the high-confidence network, suggests that they represent specialized phenotypic outputs rather than central drivers. Notably, although differential methylation at the HLA-DRB1*15:01 locus is a hallmark consistently reported in AOMS studies, 11-13 our analysis did not identify significant DMRs mapping to this region in the pediatric cohort. Nevertheless, our findings align with recent evidence highlighting how environmental triggers and genetic risk factors converge on shared epigenetic pathways leading to MS onset. 18 This supports the hypothesis that DNA methylation serves as a critical regulatory interface between early exposures and the PedMS immune profile. This is further corroborated by the high EBV seroprevalence observed in our cohort, although a contribution of the age gap between cases and controls (16 vs 14 years) to this difference cannot be excluded.

This study is subject to limitations common in EWAS involving rare clinical populations. First, the cross-sectional design of this study precludes the establishment of causality, leaving open whether the detected signatures are a cause or consequence of the disease. Second, the use of whole blood DNA, despite rigorous in silico correction for cell heterogeneity, might still mask subtle tissue-specific or cell-type–specific methylation changes. Nevertheless, it has been demonstrated that even small methylation differences in whole blood can reflect a nonlinear correlation with gene expression in specific tissues, particularly for imprinted genes 48 and pediatric subjects. 49 Moreover, given the rarity of the pediatric form, which accounts for only 5% of all MS cases, sample size and statistical power are inherently constrained. We acknowledge the demographic imbalance between PedMS cases and controls. The rigorous statistical framework adjusting for sex and cellular heterogeneity and the resilience of our findings across both sensitivity analyses and internal validation subsets, despite demographic skewing, suggests these signals represent core biological hallmarks rather than artifacts, although larger prospectively matched cohorts and external validation remain essential to further refine and confirm these observations.

The present multicentric Italian study provides novel epigenetic evidence linking DNA methylation alterations to the pathogenesis of PedMS. We successfully identified 55 DMRs, many of which maps to promoter regions of genes critical for immune regulation and previously implicated in MS, such as CHI3L2 , LPCAT1 , PARP-1 , IGF2R , and CD19 . Our network-based enrichment analysis, which identified PARP1 and IGF2R as high-centrality regulatory hubs, supports the involvement of EBV infection in the molecular pathways dysregulated by these epigenetic changes, suggesting a mechanism where specific DNA methylation patterns may mediate the interaction between genetic predisposition and environmental triggers in the disease's early onset. Furthermore, the identification of DMRs overlapping with imprinting control regions and master regulators such as ZFP57 suggests a potential role for parent-of-origin effects (POEs) in shaping the epigenetic vulnerability to PedMS. These findings highlight the importance of future studies integrating targeted sequencing to move beyond suggestive genomic overlaps.

Although the cross-sectional nature of this study limits causal inference, these findings underscore the potential of regional DNA methylation analysis to uncover reproducible epigenetic patterns even in the earliest stages of the disease. Further study prospects concern the validation of our results at the transcriptomic level, as well as the leveraging of exposure data available for our cohort to establish a functional link between the identified epigenetic signature and the subject's environmental history.

Acknowledgment

The authors thank all individuals, patients, families, contributors, and institutions who participated in this study, particularly FISM (Fondazione Italiana Sclerosi Multipla).

Author Contributions

A. Corona: drafting/revision of the manuscript for content, including medical writing for content; study concept or design; analysis or interpretation of data. M. Tosi: drafting/revision of the manuscript for content, including medical writing for content; study concept or design. A. Zollo: drafting/revision of the manuscript for content, including medical writing for content. N. Barizzone: drafting/revision of the manuscript for content, including medical writing for content. N. Pomella: drafting/revision of the manuscript for content, including medical writing for content. M. Simone: major role in the acquisition of data. A. Protti: major role in the acquisition of data. A.L. Berardinelli: major role in the acquisition of data. Antonio Gallo: major role in the acquisition of data. C. Canavese: major role in the acquisition of data. D. Vecchio: major role in the acquisition of data. E. Cocco: major role in the acquisition of data. L. Moiola: major role in the acquisition of data. M.Z. Conti: major role in the acquisition of data. S. Malucchi: major role in the acquisition of data. Amanda Papa: major role in the acquisition of data. A. Mingione: drafting/revision of the manuscript for content, including medical writing for content. C. Monzani: drafting/revision of the manuscript for content, including medical writing for content; analysis or interpretation of data. P. Annovazzi: major role in the acquisition of data. L.M.E. Grimaldi: major role in the acquisition of data. R. Lanzillo: major role in the acquisition of data. S. Rasia: major role in the acquisition of data. S. Bova: major role in the acquisition of data. P. Ravanini: major role in the acquisition of data. F. Caushi: major role in the acquisition of data. V.L.A.M. Torri Clerici: major role in the acquisition of data. S. Sotgiu: major role in the acquisition of data. A. Priori: drafting/revision of the manuscript for content, including medical writing for content. M. Trojano: major role in the acquisition of data. M.P. Amato: major role in the acquisition of data. R. Bergamaschi: major role in the acquisition of data. S. Pilotto: major role in the acquisition of data. C. Pozzilli: major role in the acquisition of data. S. Cottone: major role in the acquisition of data. G. Santangelo: major role in the acquisition of data. G. De Luca: major role in the acquisition of data. M. Pugliatti: major role in the acquisition of data. A. Ghezzi: drafting/revision of the manuscript for content, including medical writing for content; major role in the acquisition of data; study concept or design. S. D'Alfonso: drafting/revision of the manuscript for content, including medical writing for content; major role in the acquisition of data; study concept or design. F. Martinelli Boneschi: drafting/revision of the manuscript for content, including medical writing for content; major role in the acquisition of data; study concept or design.

Study Funding

This study was funded by the Fondazione Italiana Sclerosi Multipla (grant N. 12/18/F14, Project 2022/R-Multi/013 FISM, Research Special Project FISM-2022, P.I. Angelo Ghezzi).

Disclosure

The authors report no relevant disclosures. Go to Neurology.org/NN for full disclosures.

Coinvestigators

Coinvestigators are listed at Neurology.org/NN .

Glossary

adult-onset multiple sclerosis

beta-mixture quantile normalization

cell-based assay

disease-modifying treatment

differentially methylated probe

differentially methylated region

Epstein-Barr virus

epigenome-wide association study

genomic DNA

pediatric healthy control

high-efficacy treatment agent

imprinting control region

low-efficacy treatment agent

pediatric MS

parent of origin effect

References

1. Waldman A , Ghezzi A , Bar-Or A , Mikaeloff Y , Tardieu M , Banwell B . Multiple sclerosis in children: an update on clinical diagnosis, therapeutic strategies, and research . Lancet Neurol. 2014 ; 13 ( 9 ): 936 - 948 . doi: 10.1016/S1474-4422(14)70093-6 PMC4443918 25142460

2. Ghezzi A , Baroncini D , Zaffaroni M , Comi G . Pediatric versus adult MS: similar or different? Mult Scler Demyelinating Disord. 2017 ; 2 ( 1 ): 5 . doi: 10.1186/s40893-017-0022-6

3. Waubant E , Chabas D , Okuda DT , et al . Difference in disease burden and activity in pediatric patients on brain magnetic resonance imaging at time of multiple sclerosis onset vs adults . Arch Neurol. 2009 ; 66 ( 8 ): 967 - 971 . doi: 10.1001/archneurol.2009.135 19667217

4. Chitnis T , Aaen G , Belman A , et al .; US Network of Paediatric Multiple Sclerosis Centers . Improved relapse recovery in paediatric compared to adult multiple sclerosis . Brain 2020 ; 143 ( 9 ): 2733 - 2741 . doi: 10.1093/brain/awaa199 32810215

5. Renoux C , Vukusic S , Mikaeloff Y , et al .; Adult Neurology Departments KIDMUS Study Group . Natural history of multiple sclerosis with childhood onset . N Engl J Med 2007 ; 356 ( 25 ): 2603 - 2613 . doi: 10.1056/NEJMoa067597 17582070

6. McKay KA , Manouchehrinia A , Berrigan L , Fisk JD , Olsson T , Hillert J . Long-term cognitive outcomes in patients with pediatric-onset vs adult-onset multiple sclerosis . JAMA Neurol. 2019 ; 76 ( 9 ): 1028 - 1034 . doi: 10.1001/jamaneurol.2019.1546 PMC6580443 31206130

7. Amato MP , Krupp LB , Charvet LE , Penner I , Till C . Pediatric multiple sclerosis . Neurology. 2016 ; 87 ( 9_Supplement_2 ). s82 - 7 doi: 10.1212/WNL.0000000000002883 27572867

8. Graves JS , Chitnis T , Weinstock-Guttman B , et al . Maternal and perinatal exposures are associated with risk for pediatric-onset multiple sclerosis . Pediatrics. 2017 ; 139 ( 4 ): e20162838 . doi: 10.1542/peds.2016-2838 PMC5369674 28562303

9. Ghezzi A , D'Alfonso S , Pugliatti M , et al .; PEDIGREE Study Group . A narrative review of genetic and environmental factors and risk for multiple sclerosis. The design of the Italian multicentric PEDIGREE study . Neurol Sci. 2025 ; 46 ( 10 ): 4913 - 4923 . doi: 10.1007/s10072-025-08385-3 40815455

10. Pohl D , Krone B , Rostasy K , et al . High seroprevalence of epstein–Barr virus in children with multiple sclerosis . Neurology. 2006 ; 67 ( 11 ): 2063 - 2065 . doi: 10.1212/01.wnl.0000247665.94088.8d 17159123

11. Maltby VE , Lea RA , Sanders KA , et al . Differential methylation at MHC in CD4+ T cells is associated with multiple sclerosis independently of HLA-DRB1 . Clin Epigenetics. 2017 ; 9 : 71 . doi: 10.1186/s13148-017-0371-1 PMC5516341 28729889

12. Kular L , Liu Y , Ruhrmann S , et al . DNA methylation as a mediator of HLA-DRB1*15:01 and a protective variant in multiple sclerosis . Nat Commun. 2018 ; 9 ( 1 ): 2397 . doi: 10.1038/s41467-018-04732-5 PMC6008330 29921915

13. Xavier A , Maltby VE , Ewing E , et al . DNA methylation signatures of multiple sclerosis occur independently of known genetic risk and are primarily attributed to B cells and monocytes . Int J Mol Sci. 2023 ; 24 ( 16 ): 12576 . doi: 10.3390/ijms241612576 PMC10454485 37628757

14. Huynh JL , Garg P , Thin TH , et al . Epigenome-wide differences in pathology-free regions of multiple sclerosis–affected brains . Nat Neurosci. 2014 ; 17 ( 1 ): 121 - 130 . doi: 10.1038/nn.3588 PMC3934491 24270187

15. Pahlevan Kakhki M , Giordano A , Starvaggi Cucuzza C , et al . A genetic-epigenetic interplay at 1q21.1 locus underlies CHD1L-mediated vulnerability to primary progressive multiple sclerosis . Nat Commun. 2024 ; 15 ( 1 ): 6419 . doi: 10.1038/s41467-024-50794-z PMC11289459 39079955

16. O'Neill KA , van der Veer BK , Charvet L , et al . Childhood adversity, allostatic load and epigenetic signatures in paediatric and adult-onset multiple sclerosis . Brain Commun. 2026 ; 8 ( 1 ): fcaf512 . doi: 10.1093/braincomms/fcaf512 PMC12917236 41728265

17. Yang JH , Race J , Sumowski P , et al . Epigenetic age and telomere length correlations in pediatric-onset multiple sclerosis . Mult Scler. 2025 ; 31 ( 14 ): 1641 - 1650 . doi: 10.1177/13524585251386217 PMC12614491 41195728

18. Simpson‐Yap S , Morwitch E , Tanner SA , et al .; Ausimmune Investigator Group . Epstein–Barr virus, lower vitamin D, low sun exposure, and HLA ‐ DRB1 *1501 risk variant share common epigenetic pathways leading to multiple sclerosis onset . Ann Neurol. 2026 ; 99 ( 2 ): 341 - 355 . doi: 10.1002/ana.78043 PMC12894492 41070760

19. Baulina N , Kiselev I , Favorova O . Imprinted genes and multiple sclerosis: what do we know? Int J Mol Sci. 2021 ; 22 ( 3 ): 1346 . doi: 10.3390/ijms22031346 PMC7866243 33572862

20. Krupp LB , Tardieu M , Amato MP , et al .; International Pediatric Multiple Sclerosis Study Group . International pediatric multiple sclerosis study group criteria for pediatric multiple sclerosis and immune-mediated central nervous system demyelinating disorders: revisions to the 2007 definitions . Mult Scler. 2013 ; 19 ( 10 ): 1261 - 1267 . doi: 10.1177/1352458513484547 23572237

21. Tian Y , Morris TJ , Webster AP , et al . ChAMP: updated methylation analysis pipeline for illumina BeadChips . Bioinformatics. 2017 ; 33 ( 24 ): 3982 - 3984 . doi: 10.1093/bioinformatics/btx513 PMC5860089 28961746

22. Johnson WE , Li C , Rabinovic A . Adjusting batch effects in microarray expression data using empirical bayes methods . Biostatistics. 2007 ; 8 ( 1 ): 118 - 127 . doi: 10.1093/biostatistics/kxj037 16632515

23. Horvath S . DNA methylation age of human tissues and cell types . Genome Biol. 2013 ; 14 ( 10 ): 3156 . doi: 10.1186/gb-2013-14-10-r115 PMC4015143 24138928

24. Ritchie ME , Phipson B , Wu D , et al . Limma powers differential expression analyses for RNA-sequencing and microarray studies . Nucleic Acids Res. 2015 ; 43 ( 7 ): e47 . doi: 10.1093/nar/gkv007 PMC4402510 25605792

25. Müller F , Scherer M , Assenov Y , et al . RnBeads 2.0: comprehensive analysis of DNA methylation data . Genome Biol. 2019 ; 20 ( 1 ): 55 . doi: 10.1186/s13059-019-1664-9 PMC6419383 30871603

26. Du P , Zhang X , Huang C-C , et al . Comparison of beta-value and M-value methods for quantifying methylation levels by microarray analysis . BMC Bioinformatics. 2010 ; 11 : 587 . doi: 10.1186/1471-2105-11-587 PMC3012676 21118553

27. Houseman EA , Accomando WP , Koestler DC , et al . DNA methylation arrays as surrogate measures of cell mixture distribution . BMC Bioinformatics. 2012 ; 13 : 86 . doi: 10.1186/1471-2105-13-86 PMC3532182 22568884

28. Peters TJ , Buckley MJ , Statham AL , et al . De novo identification of differentially methylated regions in the human genome . Epigenetics Chromatin. 2015 ; 8 : 6 . doi: 10.1186/1756-8935-8-6 PMC4429355 25972926

29. Campagna MP , Xavier A , Lechner-Scott J , et al . Epigenome-wide association studies: current knowledge, strategies and recommendations . Clin Epigenetics. 2021 ; 13 ( 1 ): 214 . doi: 10.1186/s13148-021-01200-8 PMC8645110 34863305

30. ENCODE Project Consortium , Moore JE , Purcaro MJ , Pratt HE , et al . Expanded encyclopaedias of DNA elements in the human and mouse genomes . Nature. 2020 ; 583 ( 7818 ): 699 - 710 . doi: 10.1038/s41586-020-2493-4 PMC7410828 32728249

31. Luo Y , Hitz BC , Gabdank I , et al . New developments on the encyclopedia of DNA elements (ENCODE) data portal . Nucleic Acids Res. 2020 ; 48 ( D1 ): D882 - D889 . doi: 10.1093/nar/gkz1062 PMC7061942 31713622

32. Pillich RT , Chen J , Churas C , et al . NDEx IQuery: a multi-method network gene set analysis leveraging the network data exchange . Bioinformatics. 2023 ; 39 ( 3 ): btad118 . doi: 10.1093/bioinformatics/btad118 PMC10023220 36882166

33. Szklarczyk D , Kirsch R , Koutrouli M , et al . The STRING database in 2023: protein-protein association networks and functional enrichment analyses for any sequenced genome of interest . Nucleic Acids Res. 2023 ; 51 ( D1 ): D638 - D646 . doi: 10.1093/nar/gkac1000 PMC9825434 36370105

34. Scardoni G , Petterlini M , Laudanna C . Analyzing biological network parameters with CentiScaPe . Bioinformatics. 2009 ; 25 ( 21 ): 2857 - 2859 . doi: 10.1093/bioinformatics/btp517 PMC2781755 19729372

35. Kolberg L , Raudvere U , Kuzmin I , Adler P , Vilo J , Peterson H . g:Profiler-interoperable web service for functional enrichment analysis and gene identifier mapping (2023 update) . Nucleic Acids Res. 2023 ; 51 ( W1 ): W207 - W212 . doi: 10.1093/nar/gkad347 PMC10320099 37144459

36. Jeong A , Oleske DM , Holman J . Epidemiology of pediatric-onset multiple sclerosis: a systematic review of the literature . J Child Neurol. 2019 ; 34 : 705 - 712 . doi: 10.1177/0883073819845827 31185780

37. Lehmann-Werman R , Neiman D , Zemmour H , et al . Identification of tissue-specific cell death using methylation patterns of circulating DNA . Proc Natl Acad Sci U S A. 2016 ; 113 ( 13 ): E1826 - E1834 . doi: 10.1073/pnas.1519286113 PMC4822610 26976580

38. Areshkov PO , Avdieiev SS , Balynska OV , Leroith D , Kavsan VM . Two closely related human members of chitinase-like family, CHI3L1 and CHI3L2, activate ERK1/2 in 293 and U373 cells but have the different influence on cell proliferation . Int J Biol Sci. 2012 ; 8 ( 1 ): 39 - 48 . doi: 10.7150/ijbs.8.39 PMC3226031 22211103

39. Hinsinger G , Galéotti N , Nabholz N , et al . Chitinase 3-like proteins as diagnostic and prognostic biomarkers of multiple sclerosis . Mult Scler. 2015 ; 21 ( 10 ): 1251 - 1261 . doi: 10.1177/1352458514561906 25698171

40. Comabella M , Sastre-Garriga J , Borras E , et al . CSF chitinase 3-Like 2 is associated with long-term disability progression in patients with progressive multiple sclerosis . Neurol Neuroimmunol Neuroinflamm. 2021 ; 8 ( 6 ): e1082 . doi: 10.1212/NXI.0000000000001082 PMC8428018 34497102

41. Moessinger C , Klizaite K , Steinhagen A , et al . Two different pathways of phosphatidylcholine synthesis, the kennedy pathway and the lands cycle, differentially regulate cellular triacylglycerol storage . BMC Cell Biol. 2014 ; 15 ( 1 ): 43 . doi: 10.1186/s12860-014-0043-3 PMC4293825 25491198

42. Farez MF , Quintana FJ , Gandhi R , Izquierdo G , Lucas M , Weiner HL . Toll-like receptor 2 and poly(ADP-ribose) polymerase 1 promote central nervous system neuroinflammation in progressive EAE . Nat Immunol. 2009 ; 10 ( 9 ): 958 - 964 . doi: 10.1038/ni.1775 PMC2746562 19684606

43. Meira M , Sievers C , Hoffmann F , et al . PARP-1 deregulation in multiple sclerosis . Mult Scler J Exp Transl Clin. 2019 ; 5 ( 4 ): 2055217319894604 . doi: 10.1177/2055217319894604 PMC6918498 31897308

44. Menascu S , Khavkin Y , Zilkha-Falb R , et al . Clinical and transcriptional recovery profiles in pediatric and adult multiple sclerosis patients . Ann Clin Transl Neurol. 2021 ; 8 ( 1 ): 81 - 94 . doi: 10.1002/acn3.51244 PMC7818128 33197148

45. Jima DD , Skaar DA , Planchart A , et al . Genomic map of candidate human imprint control regions: the imprintome . Epigenetics. 2022 ; 17 ( 13 ): 1920 - 1943 . doi: 10.1080/15592294.2022.2091815 PMC9665137 35786392

46. Wang Y , MacDonald RG , Thinakaran G , Kar S . Insulin-like growth Factor-II/Cation-Independent mannose 6-Phosphate receptor in neurodegenerative diseases . Mol Neurobiol. 2017 ; 54 ( 4 ): 2636 - 2658 . doi: 10.1007/s12035-016-9849-7 PMC5901910 26993302

47. Jiang W , Shi J , Zhao J , et al . ZFP57 dictates allelic expression switch of target imprinted genes . Proc Natl Acad Sci U S A. 2021 ; 118 ( 5 ): e2005377118 . doi: 10.1073/pnas.2005377118 PMC7865185 33500348

48. Murphy SK , Huang Z , Hoyo C . Differentially methylated regions of imprinted genes in prenatal, perinatal and postnatal human tissues . PLoS One. 2012 ; 7 : e40924 . doi: 10.1371/journal.pone.0040924 PMC3396645 22808284

49. Breton CV , Marsit CJ , Faustman E , et al . Small-magnitude effect sizes in epigenetic end points are important in children's environmental health studies: the children's environmental health and disease prevention research center's epigenetics working group . Environ Health Perspect. 2017 ; 125 ( 4 ): 511 - 526 . doi: 10.1289/EHP595 PMC5382002 28362264