Introduction
Alzheimer’s disease (AD) is the most prevalent cause of dementia worldwide and is characterized by the combined presence of β-amyloid plaques and tau neurofibrillary tangle deposition1. AD pathogenesis is complex and may involve the disruption of many molecular pathways related to protein-protein interaction as well as protein-lipid interactions2,[3](#ref-CR3 “Hajjar, I., Liu, C., Jones, D. P. & Uppal, K. Untargeted metabolomics reveal dysregulations in sugar, methionine, and tyrosine pathways in the prodromal state of AD. Alzheimer’s. …
Introduction
Alzheimer’s disease (AD) is the most prevalent cause of dementia worldwide and is characterized by the combined presence of β-amyloid plaques and tau neurofibrillary tangle deposition1. AD pathogenesis is complex and may involve the disruption of many molecular pathways related to protein-protein interaction as well as protein-lipid interactions2,3,4,5,[6](#ref-CR6 “Guo, Q. et al. Heparin-enriched plasma proteome is significantly altered in Alzheimer’s Disease. Res. Sq. https://doi.org/10.21203/rs.3.rs-3933136/v1
(2024).“),7. The prevalence of cognitive impairment increases exponentially with advanced aging making this the most important risk factor for AD7. Given the complexity of aging and multifactorial pathology in AD, phenotypic heterogeneity in AD ranges from asymptomatic AD (AAD), to mild cognitive impairment (MCI), to dementia (specifically dementia due to AD, which we subsequently refer to as symptomatic AD, SAD)8; however, the ability to categorize and diagnose AD via observation of biochemical changes in the brain remains elusive. To gain insight into the molecular networks in AD brains, we aimed at conducting non-targeted lipidomics on the post-mortem brains from the Religious Orders Study (ROS) or Rush Memory and Aging Project (MAP) and integrated our lipidomics data with the matching non-targeted proteomic data obtained by our team9 on the individuals from the same cohort, with further integration of rich antemortem clinical and post-mortem pathological traits for the ROSMAP individuals.
The brain is the second-most lipid-rich organ after adipose tissue10 and lipid trafficking is essential for maintaining physiological functions such as membrane formation, neurotransmitter transport, and energy production11,12,13. Lipids constitute 50-60% of the brain dry weight14 with neural membranes and myelin comprised largely of cholesterol, sphingolipids, and glycerophospholipids15. Increased cholesterol levels in the brain directly regulate γ-secretase activity16 and sphingolipids, including ceramides, sphingomyelin, and glycosphingolipids (GSL or HexCer), are directly involved in amyloid precursor protein (APP) metabolism, leading to β-amyloid accumulation and aggregation17. β-amyloid oligomers induce the hydrolysis of glycerophospholipids such as phosphoinositol (PI) and phosphatidylcholine (PC) by phospholipase C (PLC) and PLD2/PLA2, respectively, and further elevate calcium release from the endoplasmic reticulum (ER) as well as inflammatory prostaglandin production, resulting in synaptotoxicity16,18. Recently, targeted and non-targeted lipidomic investigations of AD have demonstrated unbalanced lipid homeostasis associated with brain pathology in biological fluids and brain tissues8,19. Although abnormal lipid metabolism has been known as a critical contributor to AD pathogenesis, lipid profiles of the AD brain and how they relate to AD pathophysiology remain incompletely explored.
Recent advancements in liquid-chromatography coupled with high-resolution mass spectrometry permit deep interrogation, identification, and quantification of the lipidome. To investigate the lipidomic biomarkers associated with antemortem cognition and neuropathological changes in AD, we utilized lipidomic methodology to compare the lipid profiles of 316 post-mortem dorsolateral prefrontal cortex (DLPFC) samples from cases of control, asymptomatic AD (AAD), and symptomatic AD (SAD) within the Religious Orders Study and Rush Memory and Aging Project (ROSMAP). In this work, we identify differentially abundant LPE species (LPE 20:4 and LPE 22:6) between groups that are highly associated with AD case classifications. Using a multivariate, module-based approach, the lipid modules most strongly associated with AD clinical traits are enriched in LPE only or a mixture of LPE and LPC, and are highly correlated to previously reported AD-relevant protein modules including the following: MAPK/metabolism, Post-synaptic density, and cell-ECM interaction9. Our results highlight the importance of lipid-protein interaction networks in the AD pathogenesis and provide a robust framework for future multi-omics studies on AD brains.
Results
Overview of workflow
Our overall strategy was to discriminate significant differences in DLPFC-derived lipids between groups of AAD or SAD individuals relative to controls, and to identify the biological pathways associated with these significant lipids using integrated omics methods. Briefly, we processed 316 autopsied tissues from gray matter of DLPFC in 92 controls, 77 AAD subjects, and 147 SAD subjects. A non-targeted UPLC-MS/MS assay was applied to obtain a broad coverage of the brain lipidome through the combination of high-resolution accurate-mass and ddMS2 data acquisition. The raw data were processed, and all mass spectra were manually investigated through LipidSearch, yielding a total of 2200 annotated lipids (including deuterium-labeled lipids and internal standards). 343 among these annotated lipids met the data filtering criteria and were used for batch correction by SERRF (Supplementary Fig. 1; validated lipids are listed in Supplementary Data 3). To determine significant abundance differences in lipids among control and AD groups, the normalized lipid profiles were used for single-omic statistical analysis. Lipid modules were generated by WGCNA and further integrated with proteomic modules obtained from Johnson et al.9 using DIABLO. Figure 1a illustrates this study design, with detailed information given in methods.
Fig. 1: Overview of DLPFC dataset of ROSMAP.
a Graphic illustration of the pipeline to acquire DLPFC lipidome data as input for integrative multi-omics analysis. The illustration was created in BioRender. Chen, C. (https://BioRender.com/z77u571). b Principal component analysis (PCA) plots were obtained before (left) and after (right) SERRF normalization for human brain lipidomics data acquired in negative (-) and positive (+) electrospray mode (ESI). Pooled QC samples and human brain samples are represented as red and black colors, respectively. c The distribution of sample RSD (%) presented as histrograms. Raw- and SERRF corrected-data are highlighted as blue and red, respectively. d Pie chart showing the relative lipid class composition of all ROSMAP brain samples. AcCa Acyl carnitine, Cer Ceramide, ChE Cholesterol ester, GSL Glycosphingolipids, LPC Lysophosphatidylcholine, LPE Lysophosphatidylethanolamine, MAG Monoacylglyceride, PC Phosphatidylcholine, PC_e Plasmanyl PC, PE Phosphatidylethanolamine, PE_e Plasmanyl PE, PE_p Plasmenyl PE, PG Phosphatidylglycerol, PI Phosphatidylinositol, PS Phosphatidylserine, DAG Diacylglyceride, SM Sphingomyelins, TAG Triacylglyceride.
Demographics of participants
The basic characteristics of brain samples are summarized in Table 1. Study subjects from these three groups had similar education levels, and the majority of samples in each group were female. Relative to controls, both AD groups had higher global pathology, β-amyloid, tangles, and Braak, and lower scores for CERAD as well as Reagan (P < 0.05), where the last 2 metrics are inverted scores as recorded for the ROSMAP cohort compared to standard versions of these scores. AAD and SAD groups were older on average (P < 0.001) than controls, and importantly SAD had significantly lower global (all-domain) cognitive performance (GlobalCongFunc) (P < 1.0E-7), as expected.
Differential profiles of lipid classes in brain among the three groups
QC samples’ raw data showed systematic drifts in both positive and negative electrospray modes, and after SERRF normalization, QC distribution was largely tightened to one cluster (Fig. 1b). SERRF-normalized data showed consistently reduced RSD and a large increase in the number of lipids with less than 50% of RSD compared to raw data (Fig. 1c). The covariates: age, sex, and APOE risk, showed significant differences among the three classifications based on MANOVA analysis (P < 0.05). Therefore, all data analyses presented in the main manuscript are based on covariate-adjusted results to ensure robustness and comparability. For transparency, Supplementary Fig. 14 includes the main results of the study, shown both covariate-adjusted and unadjusted results for reference. A total of 343 lipid species across 15 lipid classes were measured, and the first and second in net abundance in brain were triacylglycerides (TAG) (15.2%) and phosphatidylcholine (PC) species (12.8%), respectively (Fig. 1d). As depicted in Fig. 2a, a heatmap displaying an overview of lipid classes reveals observable variations in the z-score patterns of acyl carnitines (AcCa), monoglycerides (MAG), diacylglycerides (DAG), TAG, and Cholesterol ester (ChE) across sample groups, in contrast to the more consistent lipid profiles in phospholipids. Of 15 lipid classes, the SAD diagnosis group had significantly lower levels of lyso-phospholipids including lysophosphatidylethanolamine (LPE) and lysophosphatidylcholine (LPC) species compared to the AAD and control groups (P < 0.05) (Fig. 2b). PC abundance was lower in the SAD group relative to the controls (P = 0.035) and AAD group (P = 0.054), and a marginally greater reduction was observed for phosphatidylethanolamine (PE) in SAD group, resulting in a slightly lower ratio of PC to PE in controls (Supplementary Fig. 2). Interestingly, the ratio of LPC to PC was significantly reduced in the SAD group compared to either control or AAD groups, indicating a decrease of potential for hydrolysis of phospholipids to lyso-phospholipids in SAD (Fig. 2c). Overall, the phospholipid classes were at similar levels in control and AAD groups, but trended lower in SAD (Supplementary Fig. 2). By dissecting lipid classes into subclasses by desaturation status and acyl chain length, the lowest mean abundances of phospholipids including PC, PE, PG, PS, PI, LPC, and LPE at either the level of desaturation or that of acyl chain length were mostly observed in the SAD group (Supplementary Fig. 3, 4a and 5a). The TAG species with 2 to 7 double bonds and acyl chain lengths between 50 and 58 carbons were generally higher in the SAD group than controls, with most species showing increased levels (Supplementary Figs. 3, 4b, 5b). In comparison to controls, phospholipids with acyl chain lengths of 22, 34, 35, 40, 44 net carbons (Supplementary Fig. 4a) and TAG with 52, 56, and 58 net carbons were significantly lower in AAD, while TAG with 47, 48, 51, 53 carbons was slightly higher in AAD (Supplementary Fig. 4b). Concerning desaturation status, AAD showed levels similar to SAD for TAG containing 0-4 double bonds, except those with 2 double bonds, while AAD had the lowest values compared to the other groups considering TAG containing 5–11 double bonds, except TAG containing 7 double bonds (Supplementary Fig. 5b). Lipid set enrichment analysis (LSEA), revealed LPE, LPC, and TAG classes were significantly different between SAD and AAD, and comparing SAD versus control groups. The PE class showed significant enrichment only in the comparison of SAD versus control groups (Fig. 2d). Leading-edge class members that drive significant enrichment of differentially abundant lipid classes are listed in Supplementary Data 4. Performing LSEA for differentially quantified lipids of specific net acyl chain length, lipid subsets with 20 and 22 carbons were significantly enriched in SAD comparing to either AAD or controls. Lipids with 56 carbons only showed significant enrichment comparing SAD to AAD (Fig. 2f). Lipids with a total acyl chain length of 22 carbons (LPE 22:6, LPC 22:6, LPE 22:5, and LPE 22:4) and of 20 carbons (LPC 20:4, LPE 20:4, and LPE 20:3) showed significant enrichment in SAD, and were mainly from LPE and LPC classes. TAG species containing oleic acid (C18:1) such as TAG 18:0_18:1_20:4, TAG 18:1_18:1_20:4, TAG 20:1_18:1_18:2, and TAG 16:0_18:1_22:6, are driving this net-56-carbon-lipid enrichment (Supplementary Data 5). The lipid subset with 12 double bonds was found to be significantly enriched in SAD versus either controls or AAD group (Fig. 2e). Intriguingly, DHA (22:6), containing 6 double bonds, is mainly attributed to the leading-edge of this subset, which also includes PE 22:6_22:6, TG 18:0_22:6_22:6, PG 22:6_22:6, PS 22:6_22:6, and TG 16:0_22:6_22:6 (Supplementary Data 6). At lipid class level, the most prevalent changes including in abundance and enrichment significance among lipids with acyl chain length of 20 and 22 carbons changing between SAD and control groups were found to be LPE or LPC, suggesting an important role of brain lyso-phospholipids in AD phenotypes.
Fig. 2: Characteristics of AD clinical traits and the non-targeted lipidome.
a The abundance of each lipid class was converted to a Z score and is represented in the heat map (blue for low expression and red for high expression). The heat map shows lipid classes differentially expressed across AD clinical traits and the three case classifications. b, c Box plots display median values (center line), all individual sample points, and the 25th and 75th percentiles (top and bottom edges of the boxes), with whiskers extending to 1.5× the interquartile range (IQR). Outliers beyond this range are plotted as gray points. Data are shown across three case classifications, SAD (n = 147), AAD (n = 77), and control (n = 92), for LPE, LPC, and the ratios of LPE/PE and LPC/PC. Statistical comparisons were performed using one-way ANOVA followed by Tukey’s multiple comparison test. The total lipid class levels are provided in a source data file. d–f Heatmap matrix of pairwise comparison for lipid class, unsaturation (double-bond) levels and total acyl chain length shown as the distribution of log2 fold change between pairs of case classes. Asterisk (*) indicates a significant enrichment (adjusted P < 0.05) based on the lipid set enrichment analysis (LSEA) method. Significance was determined using a one-sided Fisher’s exact test, with p-values adjusted using the Benjamini-Hochberg method.
LPE 22:6 [sn-1] is differentially abundant across the cohort
To understand the overall differences in lipid abundance among the three diagnosis groups, and to validate the results from Fig. 1, a circular heatmap was drawn to visualize averaged abundance of each lipid group after normalization. To highlight the most statistically significant differences between groups, the 50 lipids with the lowest P values were selected for visualization (Fig. 3a, abundance of all lipid species is shown in Supplementary Fig. 6). These were ranked from a total of 343 annotated lipids based on the results of the nonparametric Kruskal–Wallis test, all meeting an FDR threshold of less than 0.392 (raw P < 0.09). LPE 22:6 is the only significantly differentially abundant lipid across the three groups (FDR = 0.034, raw P value = 0.0001). LPE 20:4 reached statistical significance for the raw P value, 0.0014 (FDR = 0.092). LPE 22:6 regiospecificity (sn-1 vs sn-2) of the glycerol backbone has not been specified in previous studies and the physiological function of LPE might be different based on this structural distinction. In our study, lyso-phospholipid regiospecificity was determined by manual inspection of MS/MS fragment data. The fragmentation pattern of phospholipids is well established, resulting in peaks from the loss of headgroup fragments from the intact species, as well as peaks from each fatty acyl anion20,21. Lyso-phospholipids have one fatty acid and one hydroxyl group covalently bonded to the phosphoglycerol backbone. Polyunsaturated fatty acids (PUFA), are predominantly located at the sn-2 position of phospholipids; liberation of PUFA from membrane phospholipids by PLA2 results in lyso-phospholipids with a fatty acyl chain at the sn-1 position22. Fragmentation in the mass spectrometer at the sn-2 position of lyso-phospholipids with a PUFA, such as docosahexaenoic fatty acid (DHA; 22:6) is shown in Supplementary Fig. 7d, e. The mass spectrum shown are peaks resulting from the intact lipid as well as the fatty acyl anion, m/z 524 and m/z 327, respectively. Upon fragmentation, DHA rapidly loses carbon dioxide from the carboxy- terminus, resulting in a peak at m/z 283. When DHA is at the sn-1 position of lyso-phospholipids (Supplementary Fig. 7f-g), a peak at m/z 196 is also present21. This peak is reflective of the remaining hydroxy group on the newly fragmented lyso-phospholipid. A peak at m/z 196 is also present in commercially available internal standard sn-1 18:1D7 LPE (Supplementary Fig. 7b, c). For additional confirmation of positional isomers of LPE species, it has been previously reported that LPE [sn-2] isomer elutes prior to LPE [sn-1] on C18 chromatography23. This trend is appreciated in Supplementary Fig. 7d–g; here the predominant species is LPE 22:6 [sn-1].
Fig. 3: Identification of significant differentially abundant lipid species across the three groups.
a The abundance of 59 annotated lipid species were converted to Z scores and represented in a circos plot (blue for low expression and red for high expression above the mean, white). The innermost track shows color based on non-parametric Kruskal–Wallis tests across the three groups. The color gradient represents significance (P value), ranging from white (P = 0.115) to pink (P < 0.001). Specifically, orange indicates P < 0.05, green indicates P < 0.01, and pink represents the strongest significance (P < 0.001). b–d Volcano plots showing log2 fold change versus -log(P_value) for pairwise comparisons, as calculated by the R package, limma. Each lipid was tested using a moderated t-test. Lipid species demonstrating statistically significant differences in means between two case classes are highlighted in red for increased fold changes and in blue for decreased fold changes, with a significance level of P < 0.05. e Supervised PLS-DA (sPLS-DA) analysis of lipidomic data with the first two principle components (Comp 1 and Comp 2) presented. Confidence ellipses surrounding core Controls, AAD, and SAD group points (left) and the prediction background (right) as analyzed by the R package, mixOmics. f ROC analysis and AUC from sPLS-DA with 2 PCs. The AUC is presented with a 95% confidence interval (CI).
Using differential expression analysis, we identified magnitude of significant changes in lipid species between the two groups within these three clinical categories. The volcano plots highlight species with raw P value below 0.05 between groups, as presented in Fig. 3b-d. Relative to controls, a total of 25 lipids (for SAD) and 8 lipids (for AAD) were found to have apparent fold change. Intriguingly, lyso-phospholipids including LPE 22:6, LPE 22:5, LPE 22:4, LPE 20:4, LPC 22:6, LPC 20:4, LPC 18:1, and LPC 15:0, exhibited significant differences between SAD and control groups, and also between SAD versus AAD, but remained at a similar level comparing AAD earlier in the disease course to the control group. Between AAD and control groups, notable changes of annotated lipids were mostly found within phospholipids. To further investigate whether AD classification can be performed based on lipid profile, supervised partial least squares discriminant analysis (PLSDA) and feature selection was performed using the sparse PLSDA (sPLSDA) algorithm from the mixOmics R package. The first two principal components, Component 1 and Component 2 explained 42.54% and 6.13% of the total variance in lipid profile, respectively. These components were used for plotting sPLSDA, and dichotomizing the individual case samples for ROC analysis, which was performed descriptively rather than for predictive modeling. The sPLSDA result exhibited a poor classification across three groups (Fig. 3e, left panel), but interestingly, the clustering centers and background prediction areas between control and SAD groups were clearly separated (Fig. 3e, right panel). The ranking sPLSDA analysis revealed that LPE 22:6 has the highest eigenvalue of Component 1 and mainly contributes to the discrimination of controls from other groups (Supplementary Fig. 8a). The sPLS-DA algorithm, using Component 1 and Component 2, demonstrated classification accuracy (Fig. 3f) with an AUC of 0.761 (95% CI [0.708, 0.813]) under the ROC curve, specifically for the comparison of SAD vs. all other non-SAD cases (Wilcoxon test, P = 2.22E-15). However, the accuracy was slightly lower at 0.657 (95% CI [0.595, 0.720]) for distinguishing controls vs. all other cases (Wilcoxon test, P = 8.83E-06), and at an AUC of 0.673 (95% CI [0.604, 0.743]) for the classification of AAD vs. other cases (Wilcoxon test, P = 6.57E-06). The lipids with top 30 VIP scores were further used for ROC analysis using binominal regression model demonstrating the ability to distinguish between SAD vs. controls, AAD vs. controls, and SAD vs. AAD with AUC values of 0.784 (95%CI [0.725, 0.842]), 0.712 (95%CI [0.632, 0.792]), and 0.767 (95%CI [0.677,0.819]), respectively (Supplementary Fig. 8b). The sPLS-DA algorithm model yielded a substantially higher AUC than any individual annotated lipid, such as LPE 22:6 which had the highest individual species AUC value (AUC for SAD versus Control cases of 0.653, (95%CI [0.582, 0.725]), P = 3.0E-05) (Supplementary Fig. 8c). Our brain lipidomic results indicate that the most robust overall changes were found between control cases and those of the SAD group, and the most significant species changing in this cohort was LPE 22:6 [sn-1].
Construction of a lipid co-expression network
The log2-transformed data of 343 lipids were used to generate a lipid co-expression network using the WGCNA algorithm. As depicted in Fig. 4a, the resulting network consisted of 13 lipid modules ranging in size from 3 lipids to 73 lipids across 314 cases analyzed, after excluding outliers. Detailed information including module membership (module epigenlipid values, MEs) and P values are shown in Supplementary Data 7. Module membership (MM) of each lipid species within the module was highly relevant and greater than 0.60. t-distributed stochastic neighbor embedding (t-SNE) analysis was applied to lipid species within each AD lipid network by their kME values and no overlapping was observed, indicating that the lipid modules identified by WGCNA are robust and independent (Fig. 4c). To assess whether a given co-expression network module was related to AD clinical traits, we correlated the MEs to GlobalCogFunc, global AD pathology (gpath), β-amyloid, tangles, Reagan, CERAD, and Braak. We also correlated these MEs to AD case classification (Control, SAD, and AAD), determined by a combination of neuropathological and cognitive metrics as described in the method section. Two lipid modules (M2 and M13) were found to be significantly associated with all AD clinical traits (P < 0.05). These two modules have strong correlation with their intramodular LPE, and LPC class members as indicated by high intramodular kME for these members (Fig. 4b). Specifically, lipid species with MM greater than 0.9 are enriched in LPE and LPC classes; M13 module is enriched in LPC 22:6, LPC 15:0, and LPE 22:5. There was a negative correlation between the relevant module eigengenes (MEs) and the AD classification across the three groups, with M2 (cor = −0.15; P = 0.009), M5 (cor = −0.14; P = 0.012), M12 (cor = −0.14; P = 0.011), and M13 (cor = −0.17; P = 0.003) showing statistically significant associations (Fig. 4d–g). M5, M9, and M12, all enriched in different phospholipid species, were significantly associated with Braak, gpath, and tangles (|cor | (\ge ,)0.1) (Fig. 4a, b). We hypothesized that integration of these lipid modules of interest with proteomic modules across common cases of the same ROSMAP cohort generated by WGCNA could provide insight into biological processes, molecular functions, and cellular components altered in conjunction with the lipid profiles of these modules, either upstream, downstream, or coincidentally.
Fig. 4: Network analysis of the ROSMAP DLPFC lipidome.
Lipid levels in DLPFC from Controls, AAD, and SAD participants measured by UPLC-MS and analyzed by WGCNA and for differential abundance. a A lipid correlation network consisting of 13 lipid modules generated from 343 lipids in the ROSMAP cohort (n = 314 DLPFC case samples, 2 outliers were removed). Module eigenlipids, which represent PC1 of the lipid abundance within each module, are correlated with AD clinical traits including global cognitive function (GlobalCogFunc), β-amyloid, Braak, CERAD (inverted scale), Tangles, Reagan (inverted scale), and global pathology (gpath). A bubble heatmap shows the correlation and P value between lipid modules and clinical traits, in which the size of bubbles indicates P value (the larger the bubble size, the lower the P values are) and the color of bubbles indicates either positive (red) or negative correlation (blue). The numbers underneath bubbles are correlation coefficients. ME correlations and P values were performed using Pearson correlation and two-sided Student asymptotic P value for given correlations in the WGCNA R package. b The bubble heat map showing the lipid class enrichment in each of 17 lipid modules. The size of the bubble indicates the FDR (the larger the bubble size, the lower the P value.), and the color the bubble indicates the overlap (enrichment) of lipid class in each module. Significance was determined using a one-sided Fisher’s exact test, with P-values adjusted using the Benjamini-Hochberg method. c Dimensionality reduction and visualization by t-SNE analysis was applied to all the lipids within each lipid network module. d–g Module eigenlipid levels by AD case status for 4 lipid modules that had significant correlations with at least 4 AD clinical traits in (a) were assessed by Kruskal-Wallis test shown in boxplots. Box plots show the median (center line), first and third quartiles (top and bottom edges of boxes), and whiskers extending to the furthest data points within 1.5 × IQR. Outliers beyond this range are plotted as gray points; otherwise, individual samples are shown in their respective module colors. Number of samples in each AD classifications: SAD (n = 147), AAD (*n *= 76), and control (n = 91). The ME values from lipidomic co-expression modules are provided in a source data file.
Lipid-Proteome integration for AD-related DLPFC changes
A proteomic co-expression network was published by Johnson et al. in 2022 on the ROSMAP cohort DLPFC tissue samples9. Briefly, a total of 516 DLPFC tissue samples from ROSMAP24 and the Banner Sun Health Research Institute25, were analyzed by tandem mass tag mass spectrometry (TMT-MS), and a total of 8619 proteins were used to build a protein co-expression network using WGCNA. This network consisted of 44 proteomic modules and the biology represented by each module was determined using gene ontology enrichment for its constituent proteins. To identify a subset of lipid and protein modules that discriminate between the three case diagnoses in our study, we employed DIABLO26, a multi-block sPLS-DA approach. This method integrates multiple omics datasets using the input of modules (hereafter used interchangeably with MEs) from both the lipid and protein networks, along with the clinical traits common to both sets, thereby allowing for pairing of the networks. Using DIABLO’s model, although the averaged components (referred to X-variates) from proteomic, lipidomic, and clinical datasets did not show clear discrimination among the three groups, control samples were evident as a separate cluster and could already be discriminated from most SAD samples (Supplementary Fig. 9d). After including a Y-variate for the three case classes, clearer cluster separation was observed in a consensus plot of composite averaged variate axes (Fig. 5a). The clinical and proteomic data showed the highest discriminatory capacity and correlation (cor = 0.51). DIABLO analysis identified 5 lipid modules, 30 protein modules, and 5 clinical AD traits that distinguished SAD from the control group (Supplementary Fig. 9e-g). To visualize the multi-omics integration results, a network plot was used to show the inter-connection of nodes with edges kept representing absolute correlation of more than 0.7 (Fig. 5b). Lipid modules M6, and M8, consisting of TAG species, were highly correlated to protein modules for RNA splicing (lightgreen) and RNA binding (lightcyan); this 2-lipid module cluster was positively correlated to the post-synaptic density (darkgreen, cor=0.617) protein module and negatively correlated with the clinical trait for gpath (cor = -0.661). The M1 lipid module, primarily composed of PC species, was negatively correlated with clinical traits measuring Braak stage and gpath, and positively correlated with traits for Reagan and CERAD. This module was also linked to the MAPK/metabolism (black) protein module. Lipid modules M9 and M12, enriched in polyunsaturated PI, PC, and PE, were correlated to clinical trait for CERAD. Depicted in Fig. 5c, the M2 and M13 lipid modules were highly correlated to protein modules for MAPK/metabolism (black), a distinct post-synaptic density module (green), and cell-ECM interactions (green-yellow). We employed xMWAS network analysis to assess the individual correlation between top hub lipid species (MM > 0.9) extracted from lipid modules M2 and M13, and individual proteins. Subsequently, gene ontology enrichment analysis was utilized to elucidate the biology of the group of proteins strongly associated with individual lipid species. Consistently, LPE 22:6 [sn-1], which had a strong membership within the M2 lipid module, was highly associated with pathways involved in neurotransmitter-mediated neural activity and post-synaptic structural organization and synaptic plasticity (Fig. 5d). It also had significantly correlations (|cor | >0.5) with ~300 proteins, that were also broadly associated with multiple other lipid species (Supplementary Fig. 10). Utilizing a proteomic dataset from iPSC-derived neurons27 from matched ROSMAP individuals (Control = 10, SAD = 5, AAD = 12), we examined the correlations between LPE 22:6 [sn-1] and protein expression using Pearson correlation analysis. We found that LPE 22:6 [sn-1] is significantly associated with Echinoderm microtubule-associated protein-like 4 (EML4), which is essential for the stabilization of microtubules (cor = -0.605, P = 0.0008; Supplementary Fig. 13a). LPE 22:6 [sn-1] is also significant correlated to Protein phosphatase 1 and regulatory subunit 1 A (PPP1R1A), which is positively correlated with amyloid aggregates (cor = -0.592, P = 0.001; Supplementary Fig. 13b). Based on our results, we demonstrated that LPE, specifically LPE 22:6 [sn-1], is strongly associated with post-synaptic neuronal function, and could be a potential biomarker for distinguishing SAD cases from control ones.
Fig. 5: Lipid modules associated with protein modules relevant to AD neuropathology.
DIABLO-processed data for lipid modules, protein modules, and clinical traits as input was used for multi-omics integration analysis. a A DIABLO consensus component plot was generated from the identified multi-omics biomarker panel: test samples (dots) were overlaid with 95% confidence ellipses calculated from the training data (ellipses). X-variate indicates the average of the components from each dataset and Y-variate indicates different sample types, the three case classifications. b Relevance network showing the associations of the lipid network, protein network, and clinical traits at cutoff = 0.7. This plot was generated using Cytoscape. The correlation coefficients in the multi-omics integration are provided in a source data file. c A chordDiagram represents the associations of three lipid modules containing LPE and LPC to the protein network and AD clinical traits. The lipid M3 module containing LPE 20:4 and LPE 22:6 is highlighted in a black outline. d The hub lipids from lipid modules M2 (LPC 16:0, LPC 18:1, LPE 20:4, LPE 22:4, and LPE 22:6) and M13 (LPC 15:0, LPC 22:6, LPE 22:5) were further used for correlation analysis to protein expression by xMWAS. The proteins relevant to each hub lipid with correlation > 0.5 at P < 0.05 were applied to gene set enrichment analysis by the R package ClusterProfiler. The heat map of associations of each lipid species to Gene ontology (GO) terms and biological pathways. Significance was determined using a one-sided Fisher’s exact test, with p-values adjusted using the Benjamini-Hochberg method. Enrichment for a given ontology was adjusted by Bonferroni correction and is shown by -log10(P value). P values of 0.05 and 10-3 refer to -log(P values) of 1.3 and 3, respectively. ME values from proteomic co-expression modules, along with the corresponding protein expression profiles, are provided as a source data file. The pathway associated with the GO term is described in Supplementary Data 8.
Discussion
We measured non-targeted lipidomics on 316 post-mortem brain samples and employed a multi-omics network approach to integrate proteomic data available for the same brain region. Our key findings are lower abundance of total LPE and LPC in symptomatic stages of AD compared to control or asymptomatic cases. Although brain lipid composition differs between brain regions and shows dynamic changes during development and aging28, overall, the deficiency of brain PE and PE plasmalogens (PE_p) are associated with AD29,30,31,32, and a decreased ratio of PC/PE inhibits γ-secretase activity, resulting in less β-amyloid formation33. Consistent with this observation, the total abundance of PE was significantly decreased in SAD compared to controls (P = 0.0075). However, when analyzed individually, neither vinyl ether-linked PE (PE_p, P = 0.075) nor ether-linked PE (PE_e, *P *= 0.15) showed a statistically significant change. Additionally, significantly lower ether-linked PC (PC_e) was observed in SAD vs. controls (P = 0.021) (Supplementary Fig. 2). This resulted in a slightly lower PC/PE ratio in controls relative to the SAD group, in line with decreased BACE1 and APP (Supplementary Fig. 2 and 12). Lipid class variation may be due to membrane remodeling, consistent with phospholipid composition redistribution rather than a dramatic change in lipid metabolism, which is rarely seen in the brain34. Also, previous studies demonstrated that oleic acid-enriched TAG accumulates in ependymal cells in both postmortem AD brains and a transgenic AD mouse model is associated with the suppressed regeneration and homeostasis of neural stem cells35. In our LSEA results, the TAG lipid class showed significant enrichment across group comparisons in the DLPFC, with enrichment levels following the trend SAD > Control > AAD. Over half of the enriched TAG species contained oleic acid, suggesting a potential involvement of oleic acid-enriched TAG in AD pathology. Targeted lipidomics showed a decreased ratio of PC/PE was also observed in controls compared to SAD and AAD (Supplementary Fig. 11d). In general, controls tended to have slightly higher levels than other groups of phospholipids incorporating omega-3 and omega-6 PUFA, including DHA (omega-3), DPA (omega-3), AA (omega-6), and LA (omega-6) (Supplementary Fig. 11b, c). A decrease of total PUFA in phospholipids has been observed in the prefrontal cortex in other AD st