Integrated omics profiling reveals systemic dysregulation and potential biomarkers in the blood of patients with neuromyelitis optica spectrum disorders

Abstract

Background

Neuromyelitis optica spectrum disorders (NMOSD) are autoimmune conditions that affect the central nervous system. The contribution of peripheral abnormalities to the disease's pathogenesis is not well understood.

Methods

To investigate this, we employed a multi-omics approach analyzing blood samples from 52 NMOSD patients and 46 healthy controls (HC). This included mass cytometry, cytokine arrays, and targeted metabolomics. We then analyzed the peripheral changes of NMOSD, and features related to NMOSD's disease severity. Furthermore, an integrative analysis was conducted to identify the distinguishing characteristics of NMOSD from HC. Additionally, we unveiled the variations in peripheral features among different clinical subgroups within NMOSD. An independent cohort of 40 individuals with NMOSD was utilized to assess the serum levels of fibroblast activation protein alpha (FAP).

Results

Our analysis revealed a distinct peripheral immune and metabolic signature in NMOSD patients. This signature is characterized by an increase in monocytes and a decrease in regulatory T cells, dendritic cells, natural killer cells, and various T cell subsets. Additionally, we found elevated levels of inflammatory cytokines and reduced levels of tissue-repair cytokines. Metabolic changes were also evident, with higher levels of bile acids, lactates, triglycerides, and lower levels of dehydroepiandrosterone sulfate, homoarginine, octadecadienoic acid (FA[18:2]), and sphingolipids. We identified distinctive biomarkers differentiating NMOSD from HC and found blood factors correlating with disease severity. Among these, fibroblast activation protein alpha (FAP) was a notable marker of disease progression.

Conclusions

Our comprehensive blood profile analysis offers new insights into NMOSD pathophysiology, revealing significant peripheral immune and metabolic alterations. This work lays the groundwork for future biomarker identification and mechanistic studies in NMOSD, highlighting the potential of FAP as a marker of disease progression.

Introduction

Neuromyelitis optica spectrum disorders (NMOSD) are inflammatory disorders targeting the central nervous system (CNS), primarily affecting the optic nerve and spinal cord. Previously categorized as a subtype of multiple sclerosis (MS), NMOSD exhibits distinct clinical characteristics, etiologies, and therapeutic approaches [1]. Prevalence is about 1 per 100,000 in Whites and 3.5 per 100,000 among East Asians [2]. Most affected are individuals aged 35–45, experiencing high relapses and disability rates; over 90% have multiple relapses within three years [3, 4]. Recent advancements have refined the concept of NMOSD through the identification of the aquaporin-4 (AQP4) autoantibody (NMO-IgG) in roughly 70% of patients [5] and the myelin oligodendrocyte glycoprotein (MOG) autoantibody in approximately 40% of AQP4-negative cases [6]. Though crucial for diagnosis, the role of NMO-IgG in managing disease progression and treatment remains debated [7,8,9]. The recurrent, severe nature of NMOSD highlights the pressing necessity for precise diagnostic tools and a deeper understanding of its pathogenesis, given its status as an incurable disease.

There is a well-established role for humoral immunity in the pathogenesis of NMOSD, while recent studies increasingly highlight the significant involvement of cellular immunity. Elevated levels of CD69⁺T cells and CD40L⁺ CD4T cells during acute phases, associated with greater disease severity, indicate a key role for cellular immunity [10, 11]. The proportion of IL-22-secreting CD4T cells and Th17 cells in the blood is also elevated in NMOSD [12, 13], and Th17 cells have a positive correlation with disease severity [14, 15]. These findings demonstrate that cellular immunity is involved in NMOSD pathogenesis, and their detailed understanding could contribute to the development of effective therapies. Moreover, many studies have explored the role of peripheral inflammatory factors in NMOSD pathogenesis, and IL-6 is among the most widely studied. This cytokine is significantly elevated in the blood and cerebrospinal fluid (CSF) of NMOSD patients and is believed to promote plasma cell survival, stimulate the production of AQP4 autoantibodies, disrupt the integrity and functionality of the blood–brain barrier, and enhance pro-inflammatory T lymphocyte differentiation and activation [16]. Besides, metabolic disturbances involving tryptophan metabolites also play a role in the development of NMOSD [17,18,19]. However, we still lack a comprehensive understanding of how cellular immunity and cytokines, and metabolic changes contribute to the onset and progression of NMOSD.

Considering NMOSD as an inflammatory autoimmune CNS disease increasingly implicates peripheral immune responses in exacerbating CNS pathology [20, 21]. The pronounced effects of an inflammatory milieu and the disruption of the blood–brain barrier on CNS integrity within NMOSD illuminate our limited understanding of the peripheral alterations in this condition. To address this gap and provide new potential biomarkers and therapeutic targets, we conducted a multi-omics analysis of blood samples collected from NMOSD patients and healthy controls (HC). We used this approach to comprehensively analyze the alteration of blood immune cell phenotypes, plasma proteins, and metabolites in NMOSD and to reveal disease-related alterations and severity-related systemic factors. We also established inter-correlation networks across multi-omics levels and discriminative signatures to differentiate NMOSD from HC. Furthermore, we compared peripheral differences between different clinical subgroups to identify potential biomarkers associated with disease severity. Besides, an independent cohort of NMOSD patients was used to confirm the findings. Our study provides an important basis for further investigation into the pathogenesis of NMOSD.

Materials and methods

Participants

A study cohort of 52 NMOSD patients and 46 HC were recruited from Ruijin Hospital and Huashan Hospital, Shanghai, China. An independent cohort of 40 NMOSD patients, including 20 with their first episode and 20 experiencing a relapse, was recruited from Ruijin Hospital, Shanghai, China. All the NMOSD patients fulfilled the international consensus diagnostic criteria for NMOSD [22]. The exclusion criteria were as follows: (1) cardiovascular and/or metabolic diseases; (2) psychiatric disorders and/or neurologic disease other than MS; (3) genetic disorders, cancer, or other autoimmune conditions; (4) body mass index (weight/height²) higher than 30; and (5) pregnancy. These criteria were based on in-depth medical history investigations and extensive laboratory tests, which included a complete blood count, liver function tests, kidney function tests, electrolyte levels, assessments for disseminated intravascular coagulation, tumor markers, autoimmune antibodies, and advanced imaging techniques such as head and spinal magnetic resonance imaging, chest computed tomography, abdominal ultrasound, and urinary system ultrasound. The study was approved by the Ethics Committee of Ruijin Hospital and Huashan Hospital and conducted in accordance with the principles of the Helsinki Declaration. All participants provided written informed consent.

The Expanded Disability Status Scale (EDSS) score is utilized to evaluate the disability level in all NMOSD patients. The Multiple Sclerosis Severity Score (MSSS) is used to determine the severity of NMOSD by adjusting the EDSS score to account for the duration of the disease based on a previously described method [23]. The level of plasma AQP4 and MOG antibodies of patients was determined by a cell-based assay. Briefly, this assay involves introducing the antigen gene into mammalian cells (HEK293T) to induce specific expression of the target antigen. Antibodies present in patient specimens specifically bind to the antigen, and then fluorescently labeled secondary antibodies bind to these primary antibodies. The results are interpreted based on the fluorescence observed under a microscope (DM IL LBD, Leica Microsystems CMS GmbH).

Mass cytometry analysis of blood immune cells

Peripheral blood samples were collected from both HC (n = 46) and NMOSD (n = 48) patients and collected into BD Vacutainer K2 EDTA tubes (BD, 366,643). The upper layer plasma was stored at -80 °C until further analysis after centrifuging at 350 g for 10 min. The lower layer cells were subjected to red blood cell lysis for 10 min at room temperature. The peripheral blood mononuclear cells were separated and analyzed by mass cytometry according to the previous study [24]. The cells were washed twice with phosphate-buffered saline (PBS) and stained with the Human Immune Monitoring Panel Kit (Fluidigm, Cat. No. #201,324). For viability staining, the cells were stained with Cisplatin (Fluidigm, Cat. No. #201,064) at a final concentration of 0.1 μM before surface staining for 4 min. Fc receptors were blocked by incubating the cells with Cell staining buffer for 10 min at room temperature. The surface antibodies cocktail (Supplementary Key Resource Table) was added to the cell suspension for 30 min on ice. The cells were then washed with staining buffer and fixed with 1.6% paraformaldehyde (Thermo Fisher, Cat.No. #28,908) for 10 min at room temperature. Cells were suspended in Ir-Interchelator (Fluidigm, Cat. No. #201192B) in Fix/Perm buffer (Fluidigm, Cat. No. #201,067) and incubated overnight at 2–8℃. Cells were then resuspended with Cell Acquisition Solution (Fluidigm, Cat.No. #201,237) with a 1:10 dilution of EQ Four Element Calibration beads (Fluidigm, Cat.No. #201,078) and filtered through a 35 μm nylon mesh filter cap (Corning, No. #352,235). The cells were acquired on a Helios Mass Cytometer (Fluidigm) at an event rate of 200–300 events/second. Mass cytometry data files were exported and analyzed using Cytobank analysis software (https://www.cytobank.org/), with CD45^highCD66⁻ gating to exclude granulocytes for their high proportion and instability.

Human cytokine antibody array analysis

Peripheral blood samples were obtained from both HC (n = 46) and NMOSD (n = 46) patients and placed into BD Vacutainer K2 EDTA tubes (BD, #366,643). Subsequently, the plasma samples were separated following the removal of blood cells and stored in a -80℃ until use. The detection of different cytokines in the plasma were conducted by RayBiotech (Guangzhou, China) with Quantibody® Human Cytokine Antibody Array 440 kit (RayBiotech, Inc., Cat# QAH-CAA-440). The experimental procedures were carried out following the manufacturer’s instructions, as previously described [25]. In brief, 100 μL of sample diluent was added to each well and incubated at room temperature for 30 min to block the slides. Then, 100 μL of samples or standard cytokines were added to the wells and incubated at room temperature for 2 h. After washing the wells 5 times with 150 μL of 1 × Wash Buffer I (5 min each) and 2 times with 150 μL of 1 × Wash Buffer II (5 min each) at room temperature, 80μL of detection antibody cocktail was added to each well and incubated at room temperature for 2 h. After another round of washing, 80μL of Cy3 equivalent dye-conjugated streptavidin was added and incubated for 1 h at room temperature. The signal was collected using a laser scanner (InnoScan 300 Microarray Scanner, Innopsys) after a final wash step.

MxP® Quant 500 kit metabolite measurements

Peripheral blood specimens were collected from HC (n = 46) and NMOSD (n = 47) patients and transferred into BD Vacutainer K2 EDTA tubes (BD, #366,643). Following this, the plasma was extracted by discarding the blood cells and stored in a -80℃ until use. Metabolites were measured with a targeted metabolomics approach using the MxP® Quant 500 kit (BIOCRATES Life Science AG, Innsbruck, Austria), with an ultra-performance liquid chromatography (UPLC)/MS/MS system [ExionLC UPLC (Sciex), QTRAP 6500⁺ triple quadrupole/linear ion trap MS/MS (Sciex)] which provides measurements of up to 630 endogenous and microbiome-derived metabolites quantitatively and/or semiquantitatively. The assay was conducted in accordance with the manufacturer’s instructions, as described in the previous study [24]. The MxP^® Quant 500 kit was developed and validated close to the guidelines for bioanalytical method validation by the US Food and Drug Administration (FDA) and European Medicines Agency (EMA) for its use with human plasma. It includes internal and calibration standards to ensure accurate quantification and result reproducibility. Plasma samples from 93 participants were analyzed, with each assay plate containing duplicates of a quality control pool from 10 human plasma samples. Raw MS data acquisition was performed using Analyst software version 1.6.3 (Sciex, Framingham, MA, USA). MS data processing and analyte concentration calculation were carried out by Biocrates MetIDQ™ software.

Enzyme-linked immunosorbent assay (ELISA)

In an independent cohort of 40 NMOSD patients, serum FAP levels were measured using the Human FAP ELISA Kit (cat. Code: ab193701, Abcam, USA) according to the manufacturer’s instructions. The FAP concentration (ng/mL) was calculated according to the standard curve.

Statistical analysis

The samples were from two cohorts: HC and NMOSD. To prepare the data for subsequent analysis, we excluded features of blood immune cells, plasma cytokines, and metabolites that were present in less than 30% of the participants, and then implemented K-nearest neighbor (KNN) imputation to address missing values. Coefficient of variation (CV) values of cytokines and metabolites that are less than 20% are acceptable for data analysis.

Linear regression models were employed to evaluate the differences in immune cell subsets, cytokines, and metabolites, while adjusting for confounding factors such as age and sex. This approach allows for the identification of variables that show significant changes across different groups. To account for the multiple testing problem and to control the false discovery rate, we applied the Benjamini–Hochberg procedure for adjusting p-values. Features exhibiting an adjusted p-value of less than 0.05 were identified as statistically significant. Box plots were used to visualize the significantly changed immune cell subsets or metabolites, with the central line representing the median, the box denoting the interquartile range (IQR), and the whiskers extending to the furthest data point within 1.5 times the IQR from the box to describe data spread. Bar plots were made to visualize the change of differential features between NMOSD and HC, EDSS-H (≥ 4) and EDSS-L (< 4), AQP4-H (≥ 100:1) and AQP4-L (< 100:1), T2 (> 3 years) and T1 (< 3 years), HR (≥ 1 relapse in the past year) and LR (no relapse in the past year), recurrence and first episode by ggpubr R package.

Principal component analysis (PCA) was an unconstrained ordination method for dimension reduction. PCA could determine principal components (PCs) which explain the most variance for data. Therefore, PCA projects high dimensional data on 2-dimensional scatterplot which enables the assessment of sample grouping.

Spearman correlation was applied between the selected features and MSSS. Spearman correlation analysis was performed to examine the correlation between immune cells, cytokines and metabolites, using the associate function of R package. The immune cells, cytokines and metabolites significantly differentially expressed between HC and NMOSD were identified by linear regression with adjusted for age and sex and then used in this analysis. The associations with correlation coefficient > 0.3 and p.adj < 0.05 were selected.

Integrative network analysis: identifying associations between blood-borne immune cells, plasma cytokines, and metabolites

To integrate blood immune cells, plasma cytokines, and metabolites data from HC and NMOSD groups, we used the Cross-Omics Multi-Way Association Study (xMWAS) package in R [26]. Integrative network analysis was performed using sparse partial least squares regression analysis, a multivariate approach for data integration to only include associations with |r|> 0.7 and p < 0.05, and network analysis techniques to evaluate the centrality (importance) of nodes in the integrative networks. The multilevel community detection method in xMWAS was used to identify communities of tightly connected immune cells, cytokines, metabolites. Community detection reveals topological modules comprised of functionally related features. The analysis also identifies nodes that undergo network changes, which are determined based on the Delta centrality (importance) measure (ECMnmosd—ECMcontrol), which a measure of the influence of a node in a network. A high eigenvector centrality means that a node is connected to any nodes who themselves have high scores. ECM, eigenvector centrality measure.

Integrative multi-omics analysis

To identify important signatures that highly correlated among multi-omics and potentially discriminative for HC and NMOSD, we performed Data Integration Analysis for Biomarker discovery using Latent Component (DIABLO) on omics data of the blood immune cells, the plasma cytokines and plasma metabolites (LC–MS/FIA). DIABLO is based on Partial Least Squares Discriminant Analysis (PLS-DA) and aims to integrate multi-omics data by maximizing covariance between all pairs of datasets. Prior to DIABLO, multi-omics data was log transformed. As discrimination is prioritized, the design matrix was set to 0.1. DIABLO performed 10 times repeat of tenfold cross validation by block.splsda and tune.block.splsda functions in mixomics, to tune successively model hyperparameters for a final model which minimizes classification error rate.

The first 2 components from the final model are considered and demonstrated in scatterplot using plotIndiv function for samples and plotVar function for features. Clustering of samples by group (HC/NMOSD) and features by omics dataset was assessed. The mixomics package also provides loadingplot function and cim function to reveal important features (selected by DIABLO) for each omics dataset, and to combine them in a heatmap, respectively. Heatmap that computed by cim function added a dendrogram with hierarchical clustering (Euclidean distance and complete linkage). Only the important features selected in component 1 of each omics dataset were assessed.

Results

Characteristics of participants and study design

This study enrolled 46 HC (13 males, 33 females) and 52 patients diagnosed with NMOSD (9 males, 43 females). The mean age was 42 ± 2.42 years in the HC group and 47.8 ± 2.11 years in the NMOSD group, with no significant age difference (P = 0.07) or sex difference (P = 0.23) between the groups. In the NMOSD cohort, the EDSS was 3.16 ± 0.34 and the MSSS was 5.38 ± 0.44. The mean duration of disease was 42.88 ± 9.74 months. Out of the NMOSD group, 31 patients exhibited an AQP4 antibody titer exceeding 100:1, whereas 21 patients had a titer below 100:1, including 13 who tested negative for the antibody. All NMOSD patients were negative for MOG antibodies. Detailed participant details information can be found in Table S1. For conducting a multi-omics analysis, blood samples were collected from both groups to investigate alterations in immune cell phenotypes, plasma cytokines, and metabolites in the NMOSD, analyzed using bioinformatics tools as depicted in Fig. 1A.

Fig. 1

Study design and changes of blood immune cell subsets in NMOSD. A Study design. B-M Box plots representing expression values of differential expressed blood immune cell subsets in HC (n = 46) and NMOSD (n = 48), with the central line representing the median, the box denoting the interquartile range (IQR), and the whiskers extending to the furthest data point within 1.5 times the IQR from the box. FDR corrected q-value: ^* P < 0.05, ^** P < 0.01, ^*** P < 0.001, ^**** P < 0.001. mDC, myeloid dendritic cell; pDC, plasmacytoid dendritic cell

Alterations of blood immune cell subtypes in NMOSD

To explore the influence of peripheral immune cells on the progression of NMOSD, we assessed thirty immune cell subtypes using mass cytometry (Fig. S1A). We identified many significantly changed immune cell subtypes in NMOSD patients as compared to HC group. Specifically, the NMOSD group exhibited an increased proportion of monocytes (Fig. 1B), along with decreased proportion of regulatory T cells (Treg) and its subtypes—activated Treg and secreting Treg, when compared to the HC group (Fig. 1C–E). Additionally, we observed reductions in dendritic cells (DC), encompassing both plasmacytoid DC (pDC) and myeloid DC (mDC), as well as natural killer (NK) cells, T cells, CD4⁺T cells, CD28⁺T cells, CD28⁺CD4T cells, CD28⁺CD8T cells, and Th17 cells (Fig. 1F–M, Fig. S1B-C). Additionally, we investigated potential correlations between blood immune cell subtypes and the MSSS; however, our findings indicated no significant correlations.

Alterations of plasma cytokines in NMOSD

Plasma cytokines play a pivotal role in modulating peripheral inflammation and may interact with peripheral immune cells, influencing NMOSD pathogenesis. To gain a comprehensive understanding of the role of plasma cytokines in NMOSD pathogenesis, we examined the expression profiles of 440 cytokines using a cytokine array. We found significant elevations in several cytokines in the NMOSD group compared to the HC group, including hepatocyte growth factor activator inhibitor type 2 (HAI-2), S100 calcium binding protein A8 (S100A8), growth differentiation factor 15 (GDF15), and Prostasin. In contrast, several cytokines were significantly decreased in NMOSD, such as stem cell factor (SCF), neurotrophin-3 (NT-3), bone morphogenetic protein 2 (BMP-2), galectin-2, and beta-nerve growth factor (b-NGF) (Fig. 2A).

Fig. 2

Changes of plasma cytokines in NMOSD and correlation with MSSS. A Bar plot of significantly differential expressed cytokines between HC (n = 46) and NMOSD (n = 46), determined by linear regression analysis controlling for age and sex. red: cytokines increased in NMOSD, blue: cytokines decreased in NMOSD. B Scatter plot showing the top 6 cytokines that have a negative correlation with the MSSS, as determined by spearman analysis, n = 45. C Scatter plot showing the top 6 cytokines that have a positive correlation with the MSSS, as determined by spearman analysis

Furthermore, we investigated correlations between plasma cytokine levels and the MSSS in NMOSD patients. Interestingly, our results indicated a higher number of cytokines negatively correlated with MSSS than those positively correlated (Table S2). Specifically, cytokines such as osteoactivin, alpha-fetoprotein (AFP), CD99, FAP, surfactant protein D (SP-D), and brain-derived neurotrophic factor (BDNF) displayed strong negative correlations with MSSS (Fig. 2B). On the contrary, angiogenin, insulin-like growth factor-binding protein 6 (IGFBP-6), C–C motif chemokine 5 (RANTES), thrombospondin-1 (TSP-1), chitinase-3-like protein 1 (CHI3L1), and dickkopf-3 (Dkk-3) demonstrated the most prominent positive correlations with MSSS (Fig. 2C).

Alterations of plasma metabolites in NMOSD

To gain insights into the alteration of plasma metabolites in NMOSD, we analyzed 630 plasma metabolites using targeted metabolomics. We employed liquid chromatography-mass spectrometry (LC–MS) for small molecular metabolites and flow injection analysis (FIA) for larger molecules such as cholesteryl esters, glycerophospholipids, glycerol esters, sphingolipids, and hexoses. The distinct profiles observed in NMOSD versus the HC group were notably separated in the LC–MS mode and in the PC2 dimension of FIA mode (Fig. S2A), indicating significant alterations in plasma metabolites in NMOSD.

Further analysis of differential metabolites revealed increased levels of small molecule metabolites, including bile acids (e.g., tauroursodeoxycholic acid [TDCA], taurocholic acid [TCA], glycodeoxycholic acid [GDCA], glycocholic acid [GCA]), hypoxanthine, lactate (Lac), sarcosine, p-Cresol sulfate (p-Cresol-SO4), α-aminobutyric acid (AABA), γ-aminobutyric acid (GABA), cysteine (Cys), and homocysteine (HCys) in NMOSD patients. Conversely, decreases were observed in dehydroepiandrosterone (DHEAS), homoarginine (HArg), octadecadienoic acid (FA[18:2]), serotonin, and docosahexaenoate (DHA) (Fig. 3A).

Fig. 3

Changes of plasma metabolites in NMOSD and correlation with disease severity. A Bar plot illustrates the log₂ fold change of LC–MS metabolites with significant differential expression between HC (n = 46) and NMOSD patients (n = 47), identified through statistical analysis involving linear regression models and FDR correction for multiple comparisons. B Box plots representing expression values of differential expressed sphingolipids in HC and NMOSD, with the median and interquartile range shown. FDR corrected q-value: ^** P < 0.01, ^**** P < 0.001. C Bar plot displayed enriched pathways of significantly differential expressed LC–MS and FIA metabolites using MetaboAnalyst enrichment analysis. D Scatter plots showing small molecular metabolites that have a negative correlation with the MSSS, as determined by spearman analysis

Lipid-related alterations included an increase in triglycerides (TGs), whereas levels of phosphatidylcholines (PCs), cholesterol esters (CEs), ceramides (Cer), and hexosylceramides (HexCers) decreased (Fig. S2B). Notably, various sphingomyelin (SM) subclasses consistently decreased in NMOSD, including SM.C16:0, SM.C16:1, SM.C18:0, SM.C18:1, SM.C20:2, SM.C24:0, SM.C26:0, SM(OH)C14:1, SM(OH)C16:1, SM(OH)C22:1, SM(OH)C22:2, and SM(OH)C24:1 (Fig. 3B), highlighting these metabolites' reduction in plasma as a significant NMOSD feature, potentially useful as diagnostic markers.

Pathway enrichment analysis of differential metabolites pointed to significant enrichment in numerous pathways, including aminoacyl-tRNA biosynthesis, arginine biosynthesis, taurine and hypotaurine metabolism, cysteine and methionine metabolism, primary bile acid biosynthesis, starch and sucrose metabolism, and sphingolipid metabolism (Fig. 3C). These results implicated the contribution of these altered metabolic pathways to NMOSD pathogenesis. Furthermore, correlation analysis between differential metabolites and the MSSS showed taurine, hippuric acid (HipAcid), creatinine, glutamate (Glu), ornithine (Orn), and Lac levels negatively correlated with MSSS (Fig. 3D, Table S3). Additionally, negative correlations with MSSS were also noted for lipid metabolites, including PCs, TGs, and various sphingolipid subclasses (Table S4), suggesting their potential as biomarkers for monitoring NMOSD progression.

Associations among immune cell subsets, cytokines, and metabolites

To elucidate the potential interactions among differential immune cell subsets, cytokines, and metabolites in relationship to the development of NMOSD, we conducted Spearman correlation analysis. This analysis identified significant correlations: cytokines such as HAI-2, S100A8, and GDF15 were positively correlated with monocytes, yet exhibited negative correlations with Tregs, DCs, NK cells, and various T cell subtypes. In contrast, cytokines including IL-3, DNAX accessory molecule 1 (DNAM-1), and Cadherin-11 displayed reverse correlation patterns (Fig. S3A). Furthermore, the analysis extended to correlations between immune cell subtypes and metabolites, uncovering that DHEAS, FA(18:2), and HArg negatively correlated with monocytes and positively with Tregs, DCs, NKs, and various T cell subtypes. Conversely, metabolites such as hypoxantine, xanthine, and methionine sulfoxide (Met-SO) showed positive correlations with monocytes, but negative correlations with Tregs, DCs, NKs, and various T cell subtypes (Fig. S3B). Sphingolipids (e.g., SM, Cer, HexCer) and PCs displayed negative correlations with monocytes and positive correlations with DCs, NKs, and T cell subsets (Th2, Th17, CD28⁺T, CD28⁺CD4T, CD28⁺CD8T), whereas TGs exhibited the opposite correlation pattern (Fig. S3C). These findings suggest intricate functional interactions between varying immune cell subsets and specific cytokines or metabolites.

Further, we used the xMWAS integrative analysis to contrast the multi-omics interactions between NMOSD and HC groups. In the HC group, the analysis selected 9 immune cell subtypes, 11 cytokines, and 14 metabolites, which were organized into seven communities with robust inter-correlations within the integrative network (γ ≥ 0.7, P < 0.05) (Fig. 4A and Table S5). In contrast, the NMOSD group analysis identified 11 immune cell subtypes, 54 cytokines, and 14 metabolites, forming ten communities with similarly high correlations (Fig. 4B and Table S6). Notably, this analysis underscored a significant increase in the number of plasma cytokines associated with NMOSD, highlighting the critical role of plasma cytokines in NMOSD pathogenesis. Moreover, using eigenvector delta centrality measure (DCM), we found that secreting Treg, Treg, TDCA, resting Treg, glycochenodeoxycholic acid (GCDCA), taurochenodeoxycholic acid (TCDCA), hypoxanthine, and GABA showed the largest change in eigenvector centrality between NMOSD and HC (Fig. 4C). This analysis underscores the impactful roles of Treg subtypes and secondary bile acids in NMOSD pathogenesis, demonstrating that immune cell subtypes, cytokines, and metabolite communities may cooperatively contribute to disease development.

Fig. 4

Integrative analysis of blood multi-omics data. A Integrative network plot of cytokines, immune cells, LC–MS and FIA metabolites in HC group (n = 46) by using xMWAS, correlation higher than 0.7, one color representing one community. B Integrative network plot of cytokines, immune cells, LC–MS and FIA metabolites in NMOSD group (n = 34) by using xMWAS, correlation higher than 0.7, one color representing one community. C The identification of nodes that undergo network changes, which is determined based on the Delta centrality (importance) measure (ECMnmo—ECMcontrol) in which EMC represents eigenvector centrality measure

Integrative features discriminating NMOSD from HC

The DIABLO analysis aimed to discern a multi-omics signature that could effectively distinguish NMOSD from HC. By integrating data on immune cells, metabolites, and cytokines, we identified distinct signatures within each category that distinguish NMOSD from HC (Fig. 5A–C). Notably, immune cell subtype analysis revealed higher proportions of monocytes and effector CD8T cells (TE CD8) in NMOSD, whereas HC exhibited higher proportions of CD28⁺T cells, CD28⁺CD4T cells, and activated Treg among others (Fig. 5D). In terms of metabolites, NMOSD was characterized by higher levels of GABA, sarcosine, and AABA in NMOSD, contrasting with higher levels of FA (18:2), HArg, DHEAS, SM(OH)C22:1, SM(OH)C24:1, SM(OH)C22:2, and SM.C24:0 in HC (Fig. 5E). For cytokine expression, NMOSD was marked by higher levels of HAI-2, GDF15, p53-associated parkin-like cytoplasmic protein (PARC), S100A8, and T-cell immunoglobulin mucin receptor 3 (TIM-3), while HC had higher levels of SCF, Syndecan-1, and retinol-binding protein 4 (RBP4), among others (Fig. 5F).

Fig. 5

Integrative analysis revealed discriminative features between NMOSD and HC. Omics integration analysis among immune cells, cytokines and metabolites. The highly correlated multi-omics signatures that discriminate HC (n = 46) and NMOSD (n = 34) are identified from supervised model using DIABLO. A–C scatterplot of samples on first 2 components for each block (Omics), colored by group (blue: HC, orange: NMOSD). D–F Loading plot of component 1 from supervised model for each block (omics), important signatures resulted from DIABLO are ordered by absolute importance (x-axis), color indicate the class for which the median expression value is the highest for each feature. G Omics integration analysis among immune cells, cytokines and metabolites. Heatmap of picked important signatures (of component 1 and 2 from each block) with dendrogram computed by Euclidean distance and complete linkage method, high correlations between features (γ > 0.5)

Further analysis involved heat-map clustering to visualize the correlations within the multi-omics signatures differentiating NMOSD from HC. This revealed two primary clusters of highly correlated communities (γ > 0.5) encompassing immune cell subtypes, metabolites, and cytokines (Fig. 5G). Notably, an elevated subcommunity in NMOSD, featuring elements like GABA, choline, Lac, AABA, S100A8, sarcosine, HAI-2, monocytes, TIM-3, PARC, and GDF15, demonstrated strong inter-correlations. Conversely, a subcommunity marked by decreased levels in NMOSD, essential for its differentiation from HC, included RBP4, Th17 cells, Th2 cells, CD28⁺T cells, CD28⁺CD4T cells, CD4T cells, Treg, and activated Treg, displaying close connections among its members. This intricate multi-omics DIABLO analysis provides a deep insight into NMOSD's distinguishing biomarker interactions, illuminating the intricate and dynamic interplay essential for understanding this disease.

Differences in peripheral characteristics among clinical subtypes of NMOSD

To elucidate the clinical heterogeneity of NMOSD, we stratified patients based on clinical parameters, including the EDSS, AQP4 antibody titer, disease duration, number of relapses in the past year, and episode history (Table S1). In our analysis, NMOSD patients were dichotomized into high and low EDSS score groups (EDSS-H, EDSS ≥ 4; EDSS-L, EDSS < 4), high and low AQP4 antibody titers (AQP4-H, antibody titer ≥ 100:1; AQP4-L, antibody titer < 100:1), by disease duration into T1 (< 3 years) and T2 (> 3 years), by relapse frequency into low relapse (LR, no relapses in the past year) and high relapse (HR, ≥ 1 relapse in the past year), and finally into first episode versus recurrent groups.

Firstly, the comparison between the EDSS-H and EDSS-L groups revealed that the former exhibited elevated levels of macrophage inflammatory protein 3a (MIP-3a), GRO (CXCL1), and matrix metalloproteinases (MMP-3, MMP-2), and alongside lower levels of IL-10, IFN-γ, and IL-7, etc. (Fig. 6A). Notably, there were no discernible differences in blood immune cell subtypes and metabolites between these subgroups.

Fig. 6

Difference in peripheral features between different clinical subgroups. Bar plot of significantly differential expressed cytokines (P < 0.05) between EDSS-H (n = 16) and EDSS-L (n = 24) (A), between AQP4-L (n = 17) and AQP4-H (n = 23) (B, between T2 (n = 16) and T1 (n = 24) (C), between HR (n = 24) and LR (n = 16) (D), between Recurrence (n = 18) and First episode (n = 22) (E). F The serum level of FAP in first episode (n = 20) and recurrence (n = 20) of NMOSD, Student’s t-test, ^*P < 0.05. G The serum level of FAP in T1 (< 3 years) (n = 27) and T2 (≥ 3 years) (n = 13) of NMOSD, Student’s t-test, ^*P < 0.05

Next, the analysis of the AQP4-L versus AQP4-H groups showed lower levels of secreting Treg in AQP4-L (Table S7), alongside significant lower levels in bone morphogenetic protein 4 (BMP-4), BDNF, and neuron-specific enolase (NSE), etc. (Fig. 6B). Metabolomic analysis showed lower levels of ceramides and higher levels of asparagine (Asn) and serine (Ser) in the AQP4-L group (Table S7).

Comparing T2 and T1 subgroups, the former displayed an increased presence of CD28⁺CD8T cells, osteopontin (OPN), and vascular cell adhesion molecule 1 (VCAM1), alongside reduced CA15-3, MIP-3b, neutrophil-activating protein 78 (ENA-78), FAP, and CCN family member 3 (NOV) (Fig. 6C, Table S8). Metabolic differences included elevated glutamine (Gln) with decreased phenylalanine (Phe) and Glu in T2 (Table S8).

In the HR versus LR comparison, the HR group had increased levels of resting Treg, endoglin, carcinoembryonic antigen-related cell adhesion molecule 1 (CEACAM-1), adiponectin, renin, and matrix extracellular phosphoglycoprotein (MEPE), but lower levels of NOV, FAP, Lac, and GLCAS (Fig. 6D, Table S9).

In the comparison between recurrent and first episodes, recurrent episode was characterized by higher levels of Th17 cells (Table S10), vascular endothelial growth factor receptor 2 (VEGF-R2), VCAM1. In contrast, markers such as CA15-3, kallikrein-5, ENA-78, and FAP were lower (Fig. 6E). Besides, Phe was notably lower in recurrent episodes (Table S10).

We examined the distinctions between the relapse and remission phases. Our findings indicate that the relapse phase is associated with higher levels of memory B cells and monocytes, and lower levels of CD28⁺CD8T cells and NK cells compared to the remission phase (Fig. S4A). Furthermore, during the relapse, there are elevated levels of SP-D and reduced levels of OPN (Fig. S4B), as well as increased AABA and decreased citrulline (Cit) levels (Fig. S4C). Besides, considering that medication might influence peripheral factors, we also assessed the effects of different treatments, including glucocorticoids, combinations of glucocorticoids with intravenous immunoglobulin (IVIG), and glucocorticoids with immunosuppressants. The analysis indicated no significant changes in the proportions of immune cell subtypes across the groups. Similarly, most cytokines and metabolites did not differ significantly among the various treatment groups. However, TSP-1 and B7 homolog 3 (B7-H3) exhibited changes with glucocorticoids combined with intravenous immunoglobulin (IVIG) or with immunosuppressants, while p-Cresol-SO4 levels varied between the use of glucocorticoids alone and glucocorticoids paired with immunosuppressants (Table S11). This suggests that these treatments have a minimal impact on peripheral factors.

Notably, FAP was differentially expressed in multiple clinical subgroups of NMOSD, such as lower expression in EDSS-H compared to EDSS-L (Fig. 6A), in T2 compared to T1 (Fig. 6C), in HR compared to LR (Fig. 6D), and in the recurrence group compared to the first episode (Fig. 6E), suggesting its potential as a key marker for NMOSD progression. To investigate the link between FAP and disease progression, we collected serum samples from an independent cohort of 40 NMOSD patients, equally divided between 20 with first episode and 20 with recurrent cases (Table S12). We used ELISA method to measure the serum FAP concentrations and compared the levels between the first episode and recurrent cases. Our findings showed that patients with recurrent cases had significantly lower serum FAP concentrations than those experiencing their first episode (Fig. 6F). Furthermore, we found that patients with a long disease duration (T2) had lower serum FAP levels compared to those with a short duration (T1) (Fig. 6G). These results collectively underscore the relevance of FAP as a crucial marker for the progression of NMOSD.

Discussion

NMOSD is a debilitating inflammatory disorder of the CNS that often leads to irreversible deficits and rapid disability progression. AQP4 autoantibodies serve as a sensitive and specific marker, facilitating the diagnosis of NMOSD. Previous studies showed that median serum AQP4-IgG titers are elevated in NMOSD during acute attack [7, 27]. However, the variability of AQP4-IgG levels both within individuals and across different cases complicates its utility in predicting relapses or the extent of disability [8, 9]. The temporal increase in serum AQP4-IgG levels preceding relapses necessitates frequent monitoring, which poses considerable challenge [28]. Current therapeutic strategies for NMOSD rely on broad-spectrum immunosuppressive treatments, making accurate diagnosis crucial to avoiding the inappropriate application of these therapies and underscoring the need for greater diagnostic precision and a deeper understanding of its pathogenesis. Prior investigations into NMOSD biomarkers have predominantly adopted singular approaches, which do not fully capture the complexity of the disorder. Our multi-omics analysis—encompassing immune cell subsets, cytokines, and metabolites in blood—aims to delineate the peripheral blood alterations in NMOSD, offering insights into its pathogenesis and identifying biomarkers.

In terms of immune cell subset alterations in NMOSD, our findings show an increase in monocytes, which might be connected to the inflammatory response involved in the development of NMOSD. Tregs play a critical role in maintaining immune system balance and preventing autoimmune diseases, including suppression of immune responses and inflammation, promotion of self-tolerance and preventing the development of autoimmune diseases [29, 30]. The observed reduction in Tregs and their subtypes—including activated, secreting, and resting Tregs, suggesting a compromised immune suppressive function of Tregs facilitated the development of NMOSD. Additionally, there was a concomitant decrease in DCs, including pDC and mDC subtypes, NK cells, and T cells, CD28⁺T and CD4T cells, implying impaired immune surveillance might contribute to NMOSD development. Notably, the reduction in NK cells aligns with findings reported in previous NMOSD study [31]. Although NMOSD is thought to exhibit show a predilection for the Th17 cell subtype [13, 14], our data reveal a significant decrease in Th17 cells, this discrepancy might be attributed to the exhaustion resulting from Th17 overactivation, but it requires further investigation.

Regarding plasma cytokines in NMOSD, our analysis identified notable changes. For instance, the simultaneous elevation of S100A8 and GDF15 suggests a pervasive peripheral inflammatory response. S100A8 activates MAP-kinase and NF-kappa-B signaling pathways through its interaction with Toll-like receptor 4 (TLR4) and the receptor for advanced glycation end products (AGER) [32,33,34,35]. Meanwhile, GDF15, a mitogenic factor induced by mitochondrial dysfunction and organ damage, mediates tolerance to inflammatory injury, underscoring a protective response to mitochondrial damage [36, 37]. HAI-2 plays crucial roles in the regulation of various biological processes involving regulation of protease activity and barrier integrity [38, 39], its increase may be a reaction to tissue injury. These findings collectively indicate that systemic inflammation plays a role in the pathogenesis of NMOSD. Conversely, reduced levels of cytokines related to neurotrophins and survival, such as NT-3, b-NGF and SCF, implying a decreased neural repair capacity facilitated the development of NMOSD [40,41,42,43]. Additionally, the decrease in Galectin-2, which binds beta-galactoside and possesses antimicrobial, antinematode, and anti-inflammatory properties [44,45,46,47] and plays an important role in the mucosal defense of the gastrointestinal tract [48,49,50], might signify a diminished anti-infective and physical barrier function are involved in the development of NMOSD.

There were also many changes in plasma metabolites in NMOSD. Notably, lactate was dramatically elevated in NMOSD, aligning with previous finding [18]. Dysregulations were also evident in bile acid metabolism, with elevations in both primary (TCA and GCA) and cytotoxic secondary bile acids (TDCA and GDCA)—metabolites modified by gut microbiota [51, 52]. Changes in gut microbiota and its metabolite bile acids have been closely linked to the pathogenesis of NMOSD [53, 54], and our findings further support the involvement of aberrant bile acid pathways in the development of NMOSD. Furthermore, the rise in p-Cresol-SO4, a neurotoxic and genotoxic substance produced by gut microbiota, suggests its role in the pathogenesis of NMOSD [55]. An imbalance was also noted between elevated homocysteine, associated with vascular injury, and decreased protective vascular markers such as HArg, FA(18:2), and DHA, suggesting that potential vascular dysfunction may be involved in the development of NMOSD [56]. Additionally, a decline in DHEAS, a neuroprotective steroid [57, 58], suggests the compromised neuroprotective function may have contributed to the development of NMOSD.

Our correlation analysis of peripheral blood alterations and MSSS scores, as well as comparison of immune cell subsets, cytokines, and metabolites across NMOSD subtypes, identified distinct peripheral markers. Notably, FAP emerged as a significant marker, with levels inversely related to MSSS and EDSS scores and lower in recurrent, longer-duration, and high-relapsed cases, which was validated in an independent cohort with lower levels in patients with recurrent episodes and those with a longer disease duration. FAP, a serine protease involved in tissue remodeling and immunosuppression [59,60,61], is typically associated with activated fibroblasts in tumor microenvironments [62] but also appears in the CNS under neuroinflammatory conditions [63,64,65,66]. The observed decrease in FAP expression, indicating reduced immunosuppressive and tissue repair capabilities, might contribute to the progression of NMOSD. Additionally, various other molecules showed significant correlations with MSSS. Osteoactivin, crucial for tissue repair, anti-inflammation, and neuroprotection [67, 68], and AFP, with known immunosuppressive effects during pregnancy, autoimmunity, and cancer [69], displayed negative correlations with MSSS, suggesting their involvement in the weakening of immunosuppressive and neuronal repair mechanisms that contribute to the development of NMOSD. Similarly, SP-D and BDNF, crucial for immune defenses and neuronal functions [70,71,72,73]. Conversely, proteins like angiogenin, IGFBP-6, RANTES, TSP-1, CHI3L1, and Dkk-3, which are implicated in neovascularization, inflammation, and tissue remodeling [73,74,75,76,77,78,79,80,81], showed positive correlations with MSSS, suggesting their involvement in the progression of NMOSD.

To identify a multi-omics signature distinguishing NMOSD from HC, DIABLO analysis was performed. This analysis revealed distinct omics-specific signatures that effectively separate NMOSD from HC, underscoring the multifaceted nature of this disease. The findings highlighted not only clear distinctions but also highly correlated communities among immune cell subtypes, metabolites, and cytokines, illustrating the intricate interplay within these biomolecular networks. Furthermore, Spearman correlation and xMWAS integrative analyses delineated the potential interactions among these components, emphasizing their collaborative impact on the pathogenesis of NMOSD. Particularly notable is the significant role of Treg subtypes and secondary bile acids.

However, our study has several limitations. First, it did not account for patients with other neurological disorders that require precise differential diagnoses, such as MS, Myelin Oligodendrocyte Glycoprotein Antibody-Associated Disease (MOGAD), and Acute Disseminated Encephalomyelitis (ADEM). Future research should explore these distinctions to identify unique peripheral characteristics and potential biomarkers specific to NMOSD compared to other diseases. Second, our analyses did not consider potential confounding factors, including dietary patterns, or genetic predispositions, which could influence the results. Third, given the rarity of NMOSD, our study's sample size was adequate but limited. Larger, multicenter and independent cohort studies are necessary to validate our findings and enhance the robustness of the results, Fourth, while our subgroup comparisons provided detailed insights into NMOSD, these findings must be validated through longitudinal studies to confirm their implications over time. Fifth, our sample comprised exclusively individuals from China; thus, findings might vary with a geographically diverse pool. The influence of ethnic and regional variations on NMOSD should be considered in future studies. Finally, our data is association study that does not imply a causal relationship, the functional implications of our findings require validation through animal studies and in vitro experiments.

Conclusions

In summary, our multi-omics approach successfully captures the nuanced alterations in peripheral blood immune cell phenotypes, cytokines, and metabolites associated with NMOSD. By establishing distinct discriminative signatures differentiating NMOSD from HC, this study has laid groundwork for future endeavors in biomarker discovery, therapeutic target development, and deeper understanding of disease pathogenesis. Notably, the identification of plasma FAP as a key biomarker highlights its potential importance in clinical applications and research. Our findings emphasize the need for further studies to elucidate the complex mechanisms by which the identified immune cell subsets, cytokines, and metabolites contribute to the pathogenesis of NMOSD, ultimately aiming to enhance diagnostic precision and therapeutic efficacy.