Abstract Background and Objectives People with multiple sclerosis (MS) have a dysregulated circulating metabolome, but the metabolome of MS brain lesions has not been studied. The aims of this study were to identify differences in the brain tissue metabolome in MS compared with controls and to assess its association with the cellular profile of corresponding tissue. Methods MS tissues included samples from the edge and core of chronic active or inactive lesions and periplaque white matter (WM). Control specimens were obtained from normal WM. Metabolomic analysis was performed using mass-spectrometry coupled with liquid/gas chromatography and subsequently integrated with single-nucleus RNA-sequencing data by correlating metabolite abundances with relative cell counts, as well as individual genes using Multiomics Factor Analysis (MOFA). Results Seventeen samples from 5 people with secondary progressive MS and 8 samples from 6 controls underwent metabolomic profiling identifying 783 metabolites. MS lesions had higher levels of sphingosines (false discovery rate-adjusted p-value[q] = 2.88E-05) and sphingomyelins and ceramides (q = 2.15E-07), but lower nucleotide (q = 0.05), energy (q = 0.001), lysophospholipid (q = 1.86E-07), and monoacylglycerol (q = 0.04) metabolite levels compared with control WM. Periplaque WM had elevated sphingomyelins and ceramides (q = 0.05) and decreased energy metabolites (q = 0.01) and lysophospholipids (q = 0.05) compared with control WM. Sphingolipids and membrane lipid metabolites were positively correlated with astrocyte and immune cell abundances and negatively correlated with oligodendrocytes. On the other hand, long-chain fatty acid, endocannabinoid, and monoacylglycerol pathways were negatively correlated with astrocyte and immune cell populations and positively correlated with oligodendrocytes. MOFA demonstrated associations between differentially expressed metabolites and genes involved in myelination and lipid biosynthesis. Discussion MS lesions and perilesional WM demonstrated a significantly altered metabolome compared with control WM. Many of the altered metabolites were associated with altered cellular composition and gene expression, indicating an important role of lipid metabolism in chronic neuroinflammation in MS. Introduction Multiple sclerosis (MS) is a chronic inflammatory disease of the CNS characterized by demyelination and neurodegeneration.^[42]1 MS lesions, which can be identified radiologically and histologically, are defined by myelin loss, inflammatory cell infiltrates, and glial reactions.^[43]2 Following the acute phase, lesions can evolve into chronic active (CA) or chronic inactive (CIna) types, based on their inflammatory activity, with the CA lesions having a hypercellular rim.^[44]3 Both lesions are characterized by hypocellular centers or cores.^[45]2[46]-[47]4 In a recent single-nucleus transcriptomic study, there were robust differences in cell composition and gene expression between brain tissues obtained from different stages of white matter (WM) MS lesions, periplaque WM, and control WM.^[48]4 We have previously demonstrated that the circulating metabolome is altered in people with MS.^[49]5,[50]6 Other studies suggest the lipid profile of MS brain lesions and normal appearing WM to be abnormal,^[51]7[52]-[53]9 but the metabolome of MS brain tissue has not been studied. Metabolism plays a critical role in determining the phenotype of immune and glial cells,^[54]10 and altering metabolic function has been shown to be a promising therapeutic strategy in neurologic diseases, including MS.^[55]11 Thus, a better understanding of the metabolites that characterize proinflammatory cell phenotypes and those that are associated with oligodendrocyte differentiation and remyelination could be the first step in identifying novel therapeutic targets. In this study, we processed brain tissues obtained from the edge of CA vs CIna lesions, demyelinated lesion core, and periplaque WM from people with secondary progressive MS, as well as normal WM from neurologically healthy brains, to identify differences in the metabolomic profile of MS lesions. Integrating single-nucleus RNA-sequencing (snRNA-seq) data obtained from the same tissue samples, we identified metabolites associated with specific cell phenotypes, which could provide novel insights into the pathophysiology of MS. Methods Brain Tissue and Sampling Strategy Frozen and paraffin-embedded mirror brain tissue blocks were provided by the Netherlands Brain Bank (5 people with secondary progressive MS and 3 age-matched and sex-matched nonaffected, nondementia controls) and the Rocky Mountain Brain Bank (3 nonaffected, nondementia controls). For tissue selection in MS samples, initial pathologic staging and characterization of MS lesions were provided by the Netherlands Brain Bank. Using a 3-mm punch-biopsy tool, frozen samples (∼50 mg) were isolated from the edge or the core of chronic active and inactive lesions, respectively, as well as from the periplaque WM; control specimens were similarly obtained from the normal WM, as previously described.^[56]4 At each location, adjacent samples were processed for both snRNA-seq and metabolomic analyses. Detailed information about the brain tissue sampling strategy along with the main clinical characteristics of MS and control cases is listed in eTable 1. Metabolomic Analysis Brain tissue obtained from 25 samples (eTable 1) was sent to Metabolon Inc (Durham, NC) for metabolomic analysis. A fraction of the mass (either 24 mg or 16 mg) was extracted and processed as follows. Samples were prepared using the automated MicroLab STAR system. Proteins were precipitated with methanol under vigorous shaking for 2 minutes (Glen Mills GenoGrinder 2000) followed by centrifugation to remove proteins and isolate protein-bound small metabolites. The extracts underwent chromatography and spectrometry analyses, as described in detail elsewhere.^[57]12 Briefly, liquid chromatography with tandem mass spectrometry and gas chromatography with mass spectrometry were performed, and the subsequent mass spectra were matched to a reference library for metabolite identification. Subsequently, the area under the curve was used for quantification of each metabolite level. Spatially Corresponding snRNA-Sequencing Data Nuclei extraction from the same brain tissue blocks was performed according to a standardized protocol, and snRNA-seq was performed using the 10X Genomics platform, as published in detail elsewhere.^[58]4 Statistical Analysis A total of 783 metabolites were identified through a global metabolomics panel (25 brain samples). Since there was 1 sample with different amounts of mass extracted for the analysis (16 mg vs 24 mg for the rest), we normalized the metabolites by the extracted mass. Metabolites that were missing in ≥30% of the samples were removed (n = 107 metabolites). Missing values were imputed using k-nearest neighbors (3 neighbors used for each imputation). Metabolite values were then log transformed and subsequently scaled. To identify any potential outliers, we used an Euclidean distance-based sample network (and a standardized network connectivity measure [Zk]); samples with Zk values <−3 would classify as outliers. Using this definition, we did not detect outliers. We compared individual log-transformed metabolites between control and MS samples using generalized estimating equation (GEE) models adjusted for age and sex (unstructured correlation structure) to account for multiple tissue samples from the same participants. Sensitivity analyses additionally adjusting for postmortem interval were also performed. To identify meaningful relationships and metabolic signatures across the samples, we took a pathway approach using weighted gene correlation network analysis (WGCNA), which can also be applied to other settings besides gene expression.^[59]13 Briefly, the goal of WGCNA is to identify clusters of highly interconnected nodes (metabolites here) using a correlation network. Using the module scores from WGCNA, we assessed whether there was any differential association between each module and MS lesions or periplaque WM compared with controls (using GEE adjusting for sex and age, as above). In a subsequent analysis, we compared differences in these scores between different types of lesions, as well as periplaque and control WM (5 group analysis), using ANOVA with Tukey post hoc tests for demonstrating evidence of differences in overall models. Integration of Metabolomic Data With the Corresponding Single-Cell Transcriptomic Data In an approach to identify which metabolic pathways were associated with specific cell populations, a metabolite set enrichment analysis (MSEA) was performed based on predefined biological classification of the metabolites provided by Metabolon, as well as on the WGCNA-derived modules. MSEA performs the same steps as gene set enrichment analysis (GSEA).^[60]14 In brief, metabolites are ranked based on the Spearman correlations with the specific proportional cell-type abundance (provided by snRNA-seq analysis), and a normalized enrichment score (NES) is calculated for each predetermined pathway, based on the size of the pathway and the degree to which the pathway's metabolites are overrepresented at the extremes (positives or negatives) of the entire ranked metabolite list. Pathways with 3 or fewer metabolites were excluded from the analysis. To further integrate metabolomics with snRNA-seq data, we used Multiomics Factor Analysis (MOFA).^[61]15 MOFA reduces the dimensionality of the data by inferring a set of factors that capture biological and technical sources of variability in multiomics data sets, in a way similar to Principal Component Analysis. These factors are derived by finding linear combinations of the original variables that maximize the explained variance, and each factor captures a different source of variability in the data. Factors of interest can then be used to discern any possible relationship with the phenotype. The features (i.e., metabolites and genes) are assigned a weight for each factor, which provides a score for how strongly they relate to each factor. Features with strong association with a factor have large absolute values, whereas features with little to no association have values close to zero. The sign of the weight indicates the directionality of the effect. In this way, highly correlated multiomics features that explain differences across sample phenotypes can be recognized. Pathway enrichment analysis (GSEA) of the top genes (using gene ontology [GO], KEGG, and Reactome databases) and metabolites (using a priori defined pathways from Metabolon) was then performed to identify the most important pathways. GSEA was performed separately for the positively and negatively weighted features. Both metabolomics and snRNA-seq data were scaled before running MOFA, and a number of 4 factors were chosen. Statistical analysis was performed using the R software, version 4.2.0. False discovery rates (q-values) were calculated using the Benjamini-Hochberg adjustment method. Statistical significance was defined as q < 0.05. Standard Protocol Approvals, Registrations, and Patient Consents Regarding tissue samples obtained from the Netherlands Brain Bank, written informed consent for the use of brain tissue and clinical data for research purposes were obtained from all donors or their next of kin in accordance with Europe's Code of conduct for Brain Banking and the International Declaration of Helsinki. For the samples obtained from the Rocky Mountain Brain Bank, tissue donation was approved by the Colorado Multiple Institutional Review Board, and informed consent was obtained from all donors or their next of kin. Data Availability Deidentified data and the code to reproduce the main findings of this study can be found at the GitHub repository ([62]github.com/dladakis/ms_tissue_metabolomics). All snRNA-seq raw data as well as expression matrix and annotations had previously been made available for reanalysis in GEO ([63]ncbi.nlm.nih.gov/geo/ under accession number [64]GSE180759). Results Individual Metabolites The results from individual metabolite comparisons across groups are included in eTable 2. It is worth mentioning that while most of the sphingolipids (sphingosines, ceramides, and sphingomyelins) were uniformly elevated in MS tissue, sphingosine-1-phosphate was not significantly different (eTable 2). WGCNA Modules Our metabolomic analysis of the brain tissue identified a total of 16 metabolite modules by WGCNA. As detailed in the Methods, since the number of metabolites was significantly larger than our sample size, we used an agnostic approach for dimension reduction and derived clusters (modules) of related metabolites using WGCNA. Modules were summarized into single scores (module scores) for each sample, which reflect the levels of the individual metabolites that they contain, with more important metabolites (high-module membership) having a larger influence on the score. Therefore, high-module scores correspond to increased levels of their respective metabolites and low-module scores to decreased levels. The 3 metabolites with the highest module membership for each metabolite module are listed in [65]Table 1. All the metabolites with their respective module membership are summarized in eTable 3. Table 1. Top 3 Metabolites for Each Module Based on Module Membership Module Module composition Metabolite Metabolite subpathway (from Metabolon) Black Nucleotide and histidine metabolism pseudouridine Pyrimidine metabolism, uracil containing 5,6-dihydrouridine Pyrimidine metabolism, uracil containing dimethylarginine (ADMA + SDMA) Urea cycle; arginine and proline metabolism Blue Dipeptides glycylleucine Dipeptides 6-phosphogluconate Pentose phosphate pathway gluconate Food component/plant Brown Unsaturated fatty acids, endocannabinoids, lysophospholipids 10-nonadecenoate (19:1n9) Long-chain monounsaturated fatty acid 10-heptadecenoate (17:1n7) Long-chain monounsaturated fatty acid eicosenoate (20:1n9 or 1n11) Long-chain monounsaturated fatty acid Cyan Benzoate metabolism & food component trimethylamine N-oxide Phospholipid metabolism guaiacol sulfate Benzoate metabolism hippurate Benzoate metabolism Green Phosphatidylethanolamines & phosphatidylcholines 1-stearoyl-2-arachidonoyl-GPE (18:0/20:4) Phosphatidylethanolamines 1-stearoyl-2-arachidonoyl-GPC (18:0/20:4) Phosphatidylcholines 1,2-dipalmitoyl-GPC (16:0/16:0) Phosphatidylcholines Green-yellow Hexosylceramides & plasmalogens glycosyl-N-palmitoyl-sphingosine (d18:1/16:0) Hexosylceramides 1-(1-enyl-palmitoyl)-2-oleoyl-GPC (P-16:0/18:1)* Plasmalogens glycosyl-N-stearoyl-sphingadienine (d18:2/18:0)* Hexosylceramides Light cyan Nucleotide metabolism 5′- GMP Purine metabolism, guanine containing AMP Purine metabolism, adenine containing UMP Pyrimidine metabolism, uracil containing Magenta Energy metabolites guanosine Purine metabolism, guanine containing pyridoxamine phosphate Vitamin B6 metabolism glutamate, gamma-methyl ester Glutamate metabolism Midnight blue Nonspecific erythronate* Aminosugar metabolism ribonate Pentose metabolism erythritol Food component/plant Pink Lysophospholipids 1-arachidonoyl-GPI* (20:4)* Lysophospholipids 3′-dephospho-CoA-glutathione* Glutathione metabolism 1-linoleoyl-GPC (18:2) Lysophospholipids Purple Hexosylceramides glycosyl ceramide (d18:1/23:1, d17:1/24:1)* Hexosylceramides 1-stearoyl-2-oleoyl-GPS (18:0/18:1) Phosphatidylserines 1-(1-enyl-oleoyl)-2-oleoyl-GPE (P-18:1/18:1)* Lysoplasmalogens Red Acyl carnitines palmitoleoylcarnitine (C16:1)* Fatty acid metabolism (acyl carnitine, monounsaturated) 3-hydroxypalmitoylcarnitine Fatty acid metabolism (acyl carnitine, hydroxy) oleoylcarnitine (C18:1) Fatty acid metabolism (acyl carnitine, monounsaturated) Salmon Sphingosines sphingosine Sphingosines hexadecasphingosine (d16:1)* Sphingosines heptadecasphingosine (d17:1) Sphingosines Tan Monoacylglycerols 2-oleoylglycerol (18:1) Monoacylglycerols 1-arachidonylglycerol (20:4) Monoacylglycerols 2-linoleoylglycerol (18:2) Monoacylglycerols Turquoise Amino acids metabolism serine Glycine, serine, and threonine metabolism lysine Lysine metabolism gamma-glutamylleucine Gamma-glutamyl amino acids Yellow Sphingomyelins and ceramides palmitoyl dihydrosphingomyelin (d18:0/16:0)* Dihydrosphingomyelins sphingomyelin (d18:1/17:0, d17:1/18:0, d19:1/16:0) Sphingomyelins palmitoyl sphingomyelin (d18:1/16:0) Sphingomyelins [66]Open in a new tab Identification of Extensive Changes in Lipid Metabolism in MS vs Control Brain Tissue Differences in metabolite modules (as derived from GEE models adjusted for sex and age) between MS lesions (including CA edge, CIna edge, and lesion core samples) and periplaque WM in comparison with control WM are listed in [67]Table 2. Table 2. Metabolomic Differences Between MS Lesions and Periplaque WM Compared With Control WM Using WGCNA Module Module composition Module size MS lesions vs control WM Periplaque WM vs control WM Standardized mean difference (95% CI) p Value q Value Standardized mean difference (95% CI) p Value q Value Black Nucleotide and histidine metabolism 34 0.04 (−0.27 to 0.36) 0.79 0.85 0.05 (−0.25 to 0.36) 0.74 0.84 Blue Dipeptides 81 0.19 (0.01 to 0.36) 0.04 0.1 0.08 (−0.16 to 0.33) 0.49 0.7 Brown Unsaturated fatty acids, endocannabinoids, lysophospholipids 69 −0.14 (−0.26 to −0.03) 0.02 0.06 0.04 (−0.11 to 0.2) 0.58 0.73 Cyan Benzoate metabolism & food component 14 −0.08 (−0.3 to 0.14) 0.49 0.7 −0.1 (−0.33 to 0.12) 0.35 0.6 Green Phosphatidylethanolamines & phosphatidylcholines 47 0.05 (−0.13 to 0.22) 0.62 0.75 −0.09 (−0.31 to 0.13) 0.44 0.69 Green-yellow Hexosylceramides & plasmalogens 18 −0.07 (−0.28 to 0.15) 0.54 0.71 −0.06 (−0.31 to 0.2) 0.66 0.78 Gray N/A 68 0.02 (−0.23 to 0.28) 0.85 0.88 0.01 (−0.25 to 0.27) 0.92 0.92 Light cyan Nucleotide metabolism 13 −0.28 (−0.5 to −0.07) 0.01 0.05 −0.22 (−0.46 to 0.03) 0.08 0.18 Magenta Energy metabolites 30 −0.31 (−0.48 to −0.15) 1.40E-04 0.001 −0.29 (−0.46 to −0.11) 0.002 0.01 Midnight blue Nonspecific 14 −0.08 (−0.26 to 0.09) 0.34 0.6 −0.29 (−0.69 to 0.11) 0.15 0.3 Pink Lysophospholipids 31 −0.37 (−0.49 to −0.24) 5.48E-09 1.86E-07 −0.13 (−0.23 to −0.03) 0.01 0.05 Purple Hexosylceramides 24 −0.21 (−0.38 to −0.04) 0.02 0.06 −0.02 (−0.18 to 0.13) 0.8 0.85 Red Acyl carnitines 40 −0.12 (−0.41 to 0.17) 0.41 0.66 −0.21 (−0.42 to 0.01) 0.06 0.14 Salmon Sphingosines 15 0.24 (0.14 to 0.35) 2.54E-06 2.88E-05 0.03 (−0.07 to 0.13) 0.54 0.71 Tan Monoacylglycerols 15 −0.23 (−0.4 to −0.06) 0.008 0.04 −0.14 (−0.3 to 0.02) 0.09 0.19 Turquoise Amino acids metabolism 115 0.17 (0.02 to 0.31) 0.03 0.08 0.19 (−0.13 to 0.5) 0.24 0.45 Yellow Sphingomyelins and ceramides 48 0.26 (0.17 to 0.35) 1.26E-08 2.15E-07 0.09 (0.02 to 0.17) 0.01 0.05 [68]Open in a new tab Abbreviations: WGCNA = weighted gene correlation network analysis; WM = white matter. Bold p- and q-values indicate statistical significance. Compared with control WM, MS lesions demonstrated increased sphingosines (standardized mean difference [SMD] = 0.24, 95% confidence intervals [CI] 0.14 to 0.35; q = 2.88E-05), sphingomyelins, and ceramides (SMD = 0.26, 95% CI 0.17 to 0.35; q = 2.15E-07), and decreased nucleotide (SMD = −0.28, 95% CI −0.5 to −0.07; q = 0.05), energy (SMD = −0.31, 95% CI −0.48 to −0.15; q = 0.001), lysophospholipid (SMD = −0.37, 95% CI −0.49 to −0.24; q = 1.86E-07), and monoacylglycerol metabolites (SMD = −0.23, 95% CI −0.4 to −0.06; q = 0.04). Compared with control WM, the perilesional WM had increased sphingomyelins and ceramides (SMD = 0.09, 95% CI 0.02 to 0.17; q = 0.05) and decreased energy (SMD = −0.29, 95% CI −0.46 to −0.11; q = 0.01) and lysophospholipid (SMD = −0.13, 95% CI −0.23 to −0.03; q = 0.05) metabolites. Results were similar in sensitivity analyses additionally adjusting for postmortem intervals (results not shown). Divergent Inner-Outer Tissue Gradients of Specific Lipid Modules In lesion cores, sphingomyelins and ceramides were higher compared with control WM, periplaque WM, and both CA and CIna lesion edge. Core lesions also demonstrated higher sphingosines and dipeptides as well as lower unsaturated fatty acids (UFAs), endocannabinoids, lysophospholipids, and nucleotide metabolites compared with control WM. CA and CIna edge lesions showed increased sphingomyelins and ceramides, while CIna edge lesions had additionally increased sphingosines compared with control WM. All tissues obtained from brains with MS had decreased lysophospholipids compared with control WM. The significant results by comparing all 5 analytical groups are shown in [69]Figure 1 and listed in eTable 4. Figure 1. Differences in Metabolic Pathways Between Various Pathologic Groups Using WGCNA. Figure 1 [70]Open in a new tab Jitter plots of the significantly different modules in the different tissue samples. Y axis denotes the module score derived from the WGCNA analysis (high-module scores represent increased metabolites of the respective modules), and X axis denotes the different types of samples (control white matter [WM], periplaque WM, chronic active edge of MS WM lesion, chronic inactive edge of MS WM lesion, and core of MS WM lesion). Horizontal bars represent the median values for each group. Comparisons were performed using the ANOVA and Tukey post hoc test. p-values were adjusted for multiple comparisons using Benjamini-Hochberg (q-values). *q < 0.05, **q < 0.01, ***q < 0.001. WGCNA = weighted gene correlation network analysis. Integrating Tissue Metabolome to Single-Cell Transcriptomic Cellular Composition snRNA-seq and cell count data of the same tissue samples were available in 19 samples.^[71]4 Spearman correlations were calculated between each individual metabolite and each cell subpopulation relative abundance and were used as inputs for MSEA ([72]Figure 2A). Metabolites were included in specific metabolic pathways based on a priori biological knowledge. A score (NES) and a p-value for each pathway and each cell type were derived. High-positive NES indicated strong positive correlations between the metabolites of the respective pathway and the cell type, while low-negative NES indicated strong negative correlations. Pathways with the strongest associations with specific cell populations are shown in [73]Figure 2B. The list of all results is presented in eTable 5 and eFigure 1. Figure 2. Integration of Metabolomics and snRNA-seq Data. [74]Figure 2 [75]Open in a new tab (A) Tissue samples were obtained from MS lesions (shown), periplaque WM, or control WM. snRNA-seq and metabolomics analyses were performed. Integration of the results was performed in 2 ways: (1) correlating individual metabolites with cell populations and performing subsequent metabolite pathway enrichment analysis and (2) using MOFA. Created with [76]BioRender.com. (B) Bubble heatmap of normalized enrichment scores derived from pathway enrichment analysis of metabolite and relative cell abundance correlation. Metabolic pathways are shown on the Y axis, while the different types of cells are shown on the X axis. The color of the bubble denotes the directionality of the relationship and the significance; red hues denote positive correlations, whereas blue hues denote negative correlations; darker colors denote significant relationships. The size of the bubbles corresponds to the absolute normalized enrichment score. Pathways positively correlated with an inflammatory cell phenotype are denoted with red font color, while pathways negatively correlated with an inflammatory cell phenotype are denoted with blue font color. MIMS = microglia inflamed in multiple sclerosis, mono/moDC = monocytes and monocyte-derived dendritic cells, Oligo = oligodendrocytes, OPC = oligodendrocyte progenitor cells. Sphingolipids, diacylglycerols, and various cell membrane lipid metabolites were positively correlated with the number of astrocyte and immune cell subpopulations and negatively correlated with oligodendrocyte and premyelinating oligodendrocyte progenitor cell (OPC) populations ([77]Figure 2B). On the other hand, endocannabinoids, long-chain fatty acids, monoacylglycerols, lysophospholipids, and pantothenate/CoA metabolites were negatively correlated with astrocytic and immune cell subpopulations and positively correlated with oligodendrocyte and premyelinating OPC populations. Similar findings were observed using the WGCNA-derived pathways for MSEA (data not shown). Multiomics Factor Analysis We used MOFA to identify metabolites and genes that explain the variability across the different samples. Of the 4 identified factors that explained most of the variance, factor 1 accounted for 40% and 12% of the gene and metabolite variability, respectively ([78]Figure 3A). Factor 1 showed a strong association with the type of sample with low factor 1 values in samples from control and periplaque WM and incrementally increasing values in CIna edge, CA edge, and finally core samples ([79]Figure 3B). The top factor 1-weighted genes and metabolites and their sample expression are demonstrated in [80]Figure 3C. Features with large positive weights were more abundant in samples with high factor 1 (MS samples) compared with samples with low factor 1 (control and periplaque WM samples), while features with large negative weights were decreased in samples with high and increased in samples with low factor 1. The top 20 factor 1-weighted genes and metabolites and their expression in each sample are shown in [81]Figure 3C. The mRNA expression pattern of top genes correlated strongly with the levels of the top metabolites. When performing pathway enrichment analysis (GO enrichment analysis for genes and MSEA for metabolites) in the top positive factor 1 features, circulatory system development (q = 2.7e-06), heart system development (q = 1.3e-05), cell-to-cell signaling (q = 1.7e-05), and synaptic signaling (q = 3.5e-04) were some of the notable gene pathways that were enriched in MS lesions (and particularly in core lesions), while phosphatidylethanolamines (q = 9e-07), dipeptides (q = 1.5e-06), and sphingomyelins (q = 1.6e-04) were the metabolic pathways enriched in MS ([82]Figure 3E). On the other hand, ensheathment of neurons (q = 6e-07), oligodendrocyte differentiation (q = 2.1e-06), cellular lipid metabolic process (q = 5.5e-04), and phospholipid biosynthetic process (q = 0.003) were the gene pathways, while hexosylceramides (q = 9.6e-12), endocannabinoids (q = 5.2e-09), lysophospholipids (q = 4.3e-05), and monoacylglycerols (q = 0.014) were the metabolic pathways downregulated in MS lesions (and particularly in the core lesions) ([83]Figure 3D). KEGG and Reactome enrichment analyses demonstrated similar results (data not shown). Figure 3. MOFA Analysis. [84]Figure 3 [85]Open in a new tab (A) Proportion of total variance explained (R^2) by individual factors derived from MOFA analysis for mRNA and metabolic data and cumulative proportion of total variance explained. (B) Violin and jitter plot of Factor 1 values by pathologic sample. (C) Top 20 weighted genes and top 20 weighted metabolites in factor 1 and their relative expression in the pathologic samples. (D) Pathway analysis based on factor 1 weights. GO biological process pathway analysis for genes is shown in the top and metabolite set enrichment analysis (MSEA) is shown in the bottom. Discussion In this study, we identified significant differences in the tissue metabolome between MS lesions and periplaque WM compared with control WM. MS lesions demonstrated increased sphingolipids, in addition to decreased nucleotide and energy metabolites, lysophospholipids, and monoacylglycerols. Increased sphingolipids and decreased energy metabolites and lysophospholipids were also found in periplaque WM compared with control WM. These differences in the metabolite levels were more pronounced in the core of the lesions, where they showed significant difference compared with both the edge of the lesion and periplaque WM samples. Based on previous findings, demonstrating a shift of cellularity and gene expression toward an inflammatory phenotype in MS lesions,^[86]4 a difference in the metabolomic profile with control WM was expected. Periplaque WM had more similarities in cell composition with control WM, but the presence of immune cells, ‘stressed’ oligodendrocytes, and astrocytes, notably increased in periplaque WM,^[87]4 could potentially explain the observed shift in the tissue metabolome. Sphingolipids were increased in MS lesions and periplaque WM compared with control WM. In addition, they were positively correlated with the relative abundance of stressed and inflammatory astrocytes, immune cells, and immune-like OPCs (which have been shown to be involved in antigen presentation and cytotoxicity in EAE models^[88]16), but negatively correlated with oligodendrocytes. Besides their structural function, ceramides participate in inflammatory processes and can induce proinflammatory astrocytic activation, reactive oxygen species production, and oligodendrocyte damage.^[89]17,[90]18 Long-chain ceramides (C16 and C18 ceramides)—which make up the majority of the ceramides in our analyses (eTable 2)—and sphingosines have been associated with cell death, necroptosis, autophagy, and mitophagy.^[91]19 Furthermore, increased levels of de novo synthesized ceramides have been reported in human brain tissue^[92]8,[93]18,[94]20 and serum^[95]6 from people with MS as well as in mouse models of EAE.^[96]18 This supports that the increased ceramides found in lesions in our study might be contributing to neuroinflammation and may not just be the result of myelin breakdown. Furthermore, the well-studied balance between ceramides and S1P—which promotes cell survival^[97]21—is disturbed in our study with increased ceramide levels but no difference in S1P levels in MS compared with control tissue. These results further reinforce the notion of potentially targeting enzymes involved in ceramide production as a therapeutic strategy in MS.^[98]22[99]-[100]24 Fingolimod seems to have a direct effect on ceramide production in the CNS in mouse models of Alzheimer disease,^[101]25 which might contribute to its therapeutic effects in addition to the lymphocyte sequestration effects. Furthermore, natural sphingomyelinase 2 (nSMase2) inhibition demonstrated restoration of the healthy balance between ceramides and sphingosines and improved behavioral outcomes in mouse models of HIV.^[102]26 Therefore, nSMase2 inhibitors might be a promising approach to test in models of MS. Polyunsaturated UFAs (PUFAs), endocannabinoids (eCBs), and lysophospholipids were diminished in the core of the lesions compared with control and periplaque WM. PUFAs have anti-inflammatory and neuroprotective functions and are mainly derived from blood.^[103]27,[104]28 Higher supply and levels of PUFAs have been associated with beneficial effects in neuronal development during infancy and childhood.^[105]29 Higher PUFAs in serum have been consistently associated with decreased risk of cognitive impairment and development of dementia.^[106]30[107]-[108]32 eCBs modulate many aspects of the immune system, both in the periphery and in CNS, and generally exert immunosuppressive effects by preventing glutamate excitotoxicity^[109]33 and attenuating proinflammatory cytokine production, leukocyte migration at inflammation sites, reactive oxygen species production, and overall immune cell activation.^[110]34,[111]35 Dysregulation of eCB metabolism and pathways has been found in both serum and CSF of people with MS as well as in brain tissue from mice with EAE.^[112]36 Numerous studies in murine models of MS have consistently demonstrated the beneficial actions of eCBs, including suppressing inflammation, ameliorating disease symptoms, and promoting neuroprotection.^[113]37[114]-[115]40 In our study, not only were their concentrations decreased in the core of MS lesions, but they were also negatively associated with inflammatory astrocytes and immune-like OPCs and positively correlated with oligodendrocytes. This suggests that paucity of these metabolites could promote proinflammatory responses and decreased myelination capacity. There have been 3 studies with ω-3 PUFA supplementation with inconclusive results or lack of sustained benefit.^[116]41[117]-[118]43 Cannabinoids have been studied for symptom management with negative results or small clinical—and mainly subjective—improvement.^[119]44 The underwhelming results of the aforementioned studies could suggest that these metabolites may not be suitable as therapeutic targets. Lysophospholipids were reduced in all types of MS samples compared with control WM. They are closely related to PUFAs and have multiple functions in the CNS mediated by different receptors.^[120]45 They are believed to promote neuronal development and attenuate EAE severity in mouse models demonstrating an anti-inflammatory and neuroprotective effect.^[121]45 Another finding in our study was that MS lesions and periplaque WM had decreased nucleotide and energy metabolites, which were profoundly reduced in the core lesions. We speculate that inflammatory cells may have altered nucleotide metabolism, and the increased numbers of apoptotic or necroptotic cells, especially in the core, could explain the decreased nucleotide and energy metabolites observed in the MS lesions compared with control WM. In addition, the mitochondrial dysfunction that exists in MS lesions^[122]46 could be a contributor to these findings. Of interest, energy metabolites were reduced in periplaque WM as well, indicating an existing pathology that could contribute to lesion expansion. An additional point of view would be that the observed metabolic differences are a result primarily of inflammation. Ceramides and sphingomyelins, along with diacylglycerols and phospholipids (including phosphatidylcholines, phosphatidylethanolamines, and plasmalogens), were correlated with the proinflammatory cellular profile. These metabolites are the major components of myelin^[123]47 and could be derived from myelin breakdown. Myelin has been identified as an inhibitor of OPC differentiation and maturation,^[124]48 and thus, this effect could be exerted through its metabolites. On the other hand, endocannabinoids, various long-chain fatty acids, lysophospholipids, monoacylglycerols, and pantothenate and coenzyme-A metabolites were correlated with promyelinating and anti-inflammatory cellular profile. These products could indicate the resolution of inflammation and the transition into the remyelinating phase. In line with these findings, a study of tissue lipidomics in a mouse model of demyelination (lysolecithin-induced) showed that ceramides and phospholipids were increased only in the acute phase of demyelination and inflammation, while free fatty acids increased in the later phase of remyelination and resolution of inflammation.^[125]49 This could potentially be attributed to the effect of various lipases that can produce either proinflammatory or anti-inflammatory lipid mediators. Using MOFA as a different unsupervised approach to integrate metabolomics and snRNA-seq data, we identified a factor that was associated with the pathologic stage. The top metabolic pathways of that factor overlapped with the modules derived from WGCNA that differed between pathologic groups and also with many of the pathways associated with cell populations in our enrichment analyses. The top genes of the same factor—which were strongly correlated with these metabolites—were involved in axon myelination, oligodendrocyte differentiation, and lipid and phospholipid biosynthesis. This suggests that transcriptomic changes may underlie the changes in cell populations, which were correlated with observed metabolic changes. An additional interesting finding from the MOFA analysis was the upregulation of circulatory and heart system development in the lesions. In a recent study, vessel density in MS lesions was deemed to be higher compared with normal-appearing WM,^[126]50 consistent with our finding on a molecular level. This could indicate that neoangiogenesis or increased vasculature could play a role in neuroinflammation. There are certain limitations of this study that need to be addressed. Our sample size was small, and larger cohorts are required to validate these findings. Although the patients were untreated for MS within the past 3 months before death and the majority of them were never treated with a disease-modifying treatment, interactions between other medications and the metabolomic profile cannot be completely excluded. In addition, the nature of metabolomic analysis, with homogenization of the brain tissue, limits our ability to distinguish between intracellular and extracellular metabolite concentrations or their cellular origin, information that could help validate some of our hypotheses. Furthermore, it is not known whether changes in the metabolomic profile of the different lesions are the main drivers of a subsequent inflammatory reaction that leads to demyelination or rather are the byproduct of demyelination and myelin breakdown. Finally, correlating metabolite levels with cell populations does not imply causation, and definite conclusions cannot be drawn about the origin of these metabolites. In conclusion, we found that metabolomic profiles differ in MS lesions and periplaque WM compared with control WM, with lipids being the most dysregulated class of metabolites. These changes might derive from a proinflammatory state (as evidenced by the cell type correlations) and may therefore play an important role in the pathophysiology of MS. Lesion cores demonstrated increased sphingolipids and decreased UFAs, eCBs, and energy compared with other MS tissues, highlighting the severe pathology in these areas. These metabolites could act as potential therapeutic targets for MS treatment, due to their interplay with the immune system. Further larger studies are needed to verify these findings, and mouse models of MS could be helpful for identifying potential therapeutic targets. Acknowledgment