N6-Methyladenosine RNA modification in cerebrospinal fluid as a novel potential diagnostic biomarker for progressive multiple sclerosis

Abstract

Background: Progressive multiple sclerosis (PMS) is an uncommon and severe subtype of MS that worsens gradually and leads to irreversible disabilities in young adults. Currently, there are no applicable or reliable biomarkers to distinguish PMS from relapsing-remitting multiple sclerosis (RRMS). Previous studies have demonstrated that dysfunction of N6-methyladenosine (m6A) RNA modification is relevant to many neurological disorders. Thus, the aim of this study was to explore the diagnostic biomarkers for PMS based on m6A regulatory genes in the cerebrospinal fluid (CSF).

Methods: Gene expression matrices were downloaded from the ArrayExpress database. Then, we identified differentially expressed m6A regulatory genes between MS and non-MS patients. MS clusters were identified by consensus clustering analysis. Next, we analyzed the correlation between clusters and clinical characteristics. The random forest (RF) algorithm was applied to select key m6A-related genes. The support vector machine (SVM) was then used to construct a diagnostic gene signature. Receiver operating characteristic (ROC) curves were plotted to evaluate the accuracy of the diagnostic model. In addition, CSF samples from MS and non-MS patients were collected and used for external validation, as evaluated by an m6A RNA Methylation Quantification Kit and by real-time quantitative polymerase chain reaction.

Results: The 13 central m6A RNA methylation regulators were all upregulated in MS patients when compared with non-MS patients. Consensus clustering analysis identified two clusters, both of which were significantly associated with MS subtypes. Next, we divided 61 MS patients into a training set (n = 41) and a test set (n = 20). The RF algorithm identified eight feature genes, and the SVM method was successfully applied to construct a diagnostic model. ROC curves revealed good performance. Finally, the analysis of 11 CSF samples demonstrated that RRMS samples exhibited significantly higher levels of m6A RNA methylation and higher gene expression levels of m6A-related genes than PMS samples.

Conclusions: The dynamic modification of m6A RNA methylation is involved in the progression of MS and could potentially represent a novel CSF biomarker for diagnosing MS and distinguishing PMS from RRMS in the early stages of the disease.

Keywords: Cerebrospinal fluid (CSF); Diagnostic biomarker; N6-methyladenosine (m6A); Progressive multiple sclerosis (PMS).

Fig. 1

The flow chart of this study in detail

Fig. 2

Identification of differentially expressed m6A-related genes. a The density plot of the two datasets before normalization. b The density plot of the two datasets after normalization. c The heat map of 13 differentially expressed m6A methylation regulators. d Differentially expressed m6A-related genes between MS and non-MS patients using the Mann–Whitney U test. (ALKBH5: MS vs. non-MS = 5.587 ± 0.807 vs. 5.003 ± 0.820, p < 0.001; FTO: MS vs. non-MS = 7.976 ± 2.128 vs. 6.471 ± 2.043, p = 0.030; HNRNPC: MS vs. non-MS = 6.735 ± 0.441 vs. 6.288 ± 0.410, p < 0.001; KIAA1429: MS vs. non-MS = 5.535 ± 0.970 vs. 4.664 ± 0.822, p < 0.001; METTL14: MS vs. non-MS = 5.657 ± 0.648 vs. 5.016 ± 0.476, p < 0.001; METTL3: MS vs. non-MS = 6.299 ± 1.097 vs. 5.512 ± 0.600, p < 0.001; RBM15: MS vs. non-MS = 6.023 ± 0.490 vs. 5.305 ± 0.554, p < 0.001; WTAP: MS vs. non-MS = 6.292 ± 0.644 vs. 5.742 ± 0.607, p < 0.001; YTHDC1: MS vs. non-MS = 6.068 ± 0.841 vs. 5.206 ± 0.516, p < 0.001; YTHDC2: MS vs. non-MS = 5.914 ± 0.645 vs. 5.341 ± 0.556, p < 0.001; YTHDF1: MS vs. non-MS = 6.791 ± 1.337 vs. 6.114 ± 1.393, p = 0.004; YTHDF2: MS vs. non-MS = 6.716 ± 1.016 vs. 6.045 ± 1.098, p < 0.001; ZC3H13: MS vs. non-MS = 6.645 ± 1.255 vs. 5.572 ± 1.154, p < 0.001). The significant levels were set at p < 0.05 (*), < 0.01 (**), and < 0.001 (***). e Correlation analysis of the relationships between different m6A-related genes

Fig. 3

Functional annotation of the DEGs between MS and non-MS patients. a The integrated analysis of gene-protein interactions of top 94 DEGs and m6A-related genes. b, c The enriched GO terms of these m6A-related genes and DEGs. d The significant KEGG pathways of these m6A-related genes and DEGs. e The PPI network of these m6A-related genes with high confidence (> 0.7)

Fig. 4

Non-supervision consensus clustering analysis of the 61 MS patients on the expression similarity of m6A-related genes. a The cumulative distribution function (CDF) of consensus clustering for k from 2 to 9. b Relative change in area under the CDF curve for k from 2 to 9. c The consensus clustering matrix for k = 2. d The tracking plot was presented to verify the principal component for k from 2 to 9

Fig. 5

The differential clinical characteristics and gene expression between cluster 1 and cluster 2. a The PCA was used to verify the two distinct subgroups divided by non-supervision consensus clustering analysis of m6A-related genes. b The total m6A-related gene expression differences of individual patient between clusters. c The correlation heatmap showed a significant association between clusters and MS subtypes. d The gene expression differences of each m6A-related gene between PMS patients and RRMS patients using the Mann–Whitney U test. (ALKBH5: RRMS vs. PMS = 6.378 ± 0.995 vs. 4.259 ± 1.000, p < 0.001; FTO: RRMS vs. PMS = 8.411 ± 1.153 vs. 6.172 ± 3.121, p < 0.001; HNRNPC: RRMS vs. PMS = 8.006 ± 1.277 vs. 4.897 ± 0.494, p < 0.001; KIAA1429: RRMS vs. PMS = 6.132 ± 1.322 vs. 4.494 ± 1.000, p < 0.001; METTL14: RRMS vs. PMS = 6.212 ± 1.026 vs. 4.805 ± 0.787, p < 0.001; METTL3: RRMS vs. PMS = 7.688 ± 2.088 vs. 4.267 ± 1.530, p < 0.001; RBM15: RRMS vs. PMS = 6.963 ± 1.186 vs. 4.477 ± 0.671, p < 0.001; WTAP: RRMS vs. PMS = 7.177 ± 1.119 vs. 4.810 ± 0.874, p < 0.001; YTHDC1: RRMS vs. PMS = 6.563 ± 1.085 vs. 5.167 ± 1.185, p < 0.001; YTHDC2: RRMS vs. PMS = 6.749 ± 0.944 vs. 4.533 ± 0.951, p < 0.001; YTHDF1: RRMS vs. PMS = 9.171 ± 2.812 vs. 3.223 ± 1.457, p < 0.001; YTHDF2: RRMS vs. PMS = 8.856 ± 2.558 vs. 3.445 ± 1.191, p < 0.001; ZC3H13: RRMS vs. PMS = 7.239 ± 1.191 vs. 5.413 ± 1.938, p < 0.001). e The gene expression differences of each m6A-related gene between SPMS patients and PPMS patients using the Mann–Whitney U test. (ALKBH5: SPMS vs. PPMS = 4.189 ± 0.737 vs. 4.346 ± 1.362, p = 0.687; FTO: SPMS vs. PPMS = 6.136 ± 3.202 vs. 6.217 ± 2.911, p = 0.947; HNRNPC: SPMS vs. PPMS = 6.136 ± 0.507 vs. 6.217 ± 0.467, p = 0.588; KIAA1429: SPMS vs. PPMS = 4.643 ± 0.823 vs. 4.311 ± 1.268, p = 0.391; METTL14: SPMS vs. PPMS = 4.815 ± 0.737 vs. 4.792 ± 0.852, p = 0.941; METTL3: SPMS vs. PPMS = 4.272 ± 1.501 vs. 4.260 ± 1.564, p = 0.984; RBM15: SPMS vs. PPMS = 4.551 ± 0.633 vs. 4.387 ± 0.700, p = 0.531; WTAP: SPMS vs. PPMS = 4.767 ± 0.821 vs. 4.862 ± 0.965, p = 0.779; YTHDC1: SPMS vs. PPMS = 5.189 ± 1.183 vs. 5.139 ± 1.187, p = 0.915; YTHDC2: SPMS vs. PPMS = 4.625 ± 0.974 vs. 4.421 ± 0.863, p = 0.583; YTHDF1: SPMS vs. PPMS = 3.062 ± 0.504 vs. 3.421 ± 2.311, p = 0.526; YTHDF2: SPMS vs. PPMS = 3.220 ± 0.890 vs. 3.721 ± 1.568, p = 0.276; ZC3H13: SPMS vs. PPMS = 5.228 ± 1.899 vs. 5.641 ± 1.993, p = 0.584). The significant levels were set at p < 0.05 (*), < 0.01 (**), and < 0.001 (***)

Fig. 6

The m6A-related feature gene selection and the diagnostic gene signature construction. a The random forest algorithm revealed that the error is small and stable after 400 nTree in the training set. b Eight feature genes were selected according to the cutoff value of 0.4, including KIAA1429, WTAP, YTHDF1, ALKBH5, YTHDF2, HNRNPC, METTL3, and YTHDC2. c The ROC curve for assessing the performance of this diagnostic gene signature in the training set. d The ROC curve for assessing the performance of this diagnostic gene signature in the test set

Fig. 7

The total m6A level and qRT-PCR validation of the m6A-related feature genes in patients with MS. a The total m6A RNA methylation context of total RNA between PMS patients and RRMS patients (PMS vs. RRMS = 0.515 ± 0.154% vs. 1.488 ± 0.611%, p = 0.036). b–i The gene expression of the feature genes ALKBH5, HNRNPC, KIAA1429, METTL3, WTAP, YTHDF1, YTHDF2, and YTHDC2 between PMS patients and RRMS patients (ALKBH5: PMS vs. RRMS = 1.004 ± 0.091 vs. 1.702 ± 0.110, p = 0.002; HNRNPC: PMS vs. RRMS = 1.002 ± 0.070 vs. 1.549 ± 0.262, p = 0.046; KIAA1429: PMS vs. RRMS = 1.003 ± 0.068 vs. 1.760 ± 0.242, p = 0.013; METTL3: PMS vs. RRMS = 1.015 ± 0.178 vs. 1.680 ± 0.202, p = 0.025; WTAP: PMS vs. RRMS = 1.002 ± 0.072 vs. 1.399 ± 0.123, p = 0.017; YTHDF1: PMS vs. RRMS = 1.003 ± 0.083 vs. 2.432 ± 0.407, p = 0.008; YTHDF2: PMS vs. RRMS = 1.012 ± 0.153 vs. 1.364 ± 0.057, p = 0.038; YTHDC2: PMS vs. RRMS = 1.078 ± 0.372 vs. 1.478 ± 0.180, p = 0.243). j–q The gene expression of the feature genes ALKBH5, HNRNPC, KIAA1429, METTL3, WTAP, YTHDF1, YTHDF2, and YTHDC2 between PMS patients and non-MS patients (ALKBH5: PMS vs. non-MS = 1.000 ± 0.026 vs. 1.034 ± 0.030, p = 0.295; HNRNPC: PMS vs. non-MS = 1.075 ± 0.429 vs. 0.830 ± 0.024, p = 0.465; KIAA1429: PMS vs. non-MS = 1.000 ± 0.025 vs. 1.085 ± 0.048, p = 0.092; METTL3: PMS vs. non-MS = 1.000 ± 0.024 vs. 0.920 ± 0.034, p = 0.005; WTAP: PMS vs. non-MS = 1.052 ± 0.351 vs. 0.866 ± 0.010, p = 0.495; YTHDC2: PMS vs. non-MS = 1.051 ± 0.301 vs. 0.753 ± 0.149, p = 0.277; YTHDF1: PMS vs. non-MS = 1.000 ± 0.007 vs. 0.972 ± 0.047, p = 0.457; YTHDF2: PMS vs. non-MS = 1.000 ± 0.006 vs. 0.973 ± 0.135, p = 0.790). The significant levels were set at p < 0.05 (*), < 0.01 (**), and < 0.001 (***)