Multi-Omics Investigation into Acute Myocardial Infarction: An Integrative Method Revealing Interconnections amongst the Metabolome, Lipidome, Glycome, and Metallome

Acute myocardial infarction (AMI) is a leading cause of mortality and morbidity worldwide. This work aims to investigate the translational potential of a multi-omics study (comprising metabolomics, lipidomics, glycomics, and metallomics) in revealing biomechanistic insights into AMI. Following the N-glycomics and metallomics studies performed by our group previously, untargeted metabolomic and lipidomic profiles were generated and analysed in this work via the use of a simultaneous metabolite/lipid extraction and liquid chromatography–tandem mass spectrometry (LC–MS/MS) analysis workflow. The workflow was applied to blood plasma samples from AMI cases (n = 101) and age-matched healthy controls (n = 66). The annotated metabolomic (number of features, n = 27) and lipidomic (n = 48) profiles, along with the glycomic (n = 37) and metallomic (n = 30) profiles of the same set of AMI and healthy samples were integrated and analysed. The integration method used here works by identifying a linear combination of maximally correlated features across the four omics datasets, via utilising both block-partial least squares-discriminant analysis (block-PLS-DA) based on sparse generalised canonical correlation analysis. Based on the multi-omics mapping of biomolecular interconnections, several postulations were derived. These include the potential roles of glycerophospholipids in N-glycan-modulated immunoregulatory effects, as well as the augmentation of the importance of Ca–ATPases in cardiovascular conditions, while also suggesting contributions of phosphatidylethanolamine in their functions. Moreover, it was shown that combining the four omics datasets synergistically enhanced the classifier performance in discriminating between AMI and healthy subjects. Fresh and intriguing insights into AMI, otherwise undetected via single-omics analysis, were revealed in this multi-omics study. Taken together, we provide evidence that a multi-omics strategy may synergistically reinforce and enhance our understanding of diseases.


Introduction
Acute myocardial infarction (AMI), a condition classified under coronary heart disease (CHD), is one of the leading causes of mortality and morbidity worldwide [1][2][3]. Novel and reliable biomarkers are needed to assist in risk assessment, accurate diagnosis and prognosis, and may also act as mediators of disease. Thus, more in-depth research, especially using 'omics' tools, is needed to better understand the biology of these novel

Untargeted Metabolomics and Lipidomics Analysis and Data Acquisition 2.2.1. Reagents and Chemicals
The following comprises chemical/reagent information for the metabolomics and lipidomics analysis. As for the glycomics and metallomics analysis, details regarding the materials used can be found in their respective publications [10,11]. All reagents used were of analytical grade unless otherwise stated. Ultrapure water (18.2 Ω) used for all chemical and sample preparation was obtained from an Ultra Clear TM water purification system by Siemenssie (Munich, Germany). For the mobile phases and extraction solvents used in the liquid chromatography-mass spectrometry (LC-MS)-based metabolomics/lipidomics analysis, LC-MS grade isopropanol, acetonitrile, methanol and formic acid were purchased from Fisher Chemicals (Waltham, MA, USA), and HPLC Plus grade tert-butyl methyl ether (MTBE) was purchased from Sigma Aldrich (St. Louis, MO, USA). For the internal standards used, SPLASH II Lipidomix Mass Spec Standard was purchased from Avanti Polar Lipids (Birmingham, AL, USA), and BOC-Leucine (99% for HPLC) standard was purchased from Sigma Aldrich.

Sample Extraction and Preparation
The polar metabolite and non-polar lipid extracts used for the untargeted analysis were obtained simultaneously using the Matyash extraction method [12]. Prior to the extraction, 40 uL of blood plasma sample was thawed on ice. This was followed by the addition of extraction solvents: 300 µL of methanol containing SPLASH II Lipidomix as well as BOC-Leucine as internal standards, and 1 mL of MTBE. Samples were then incubated at room temperature in a Vortemp Shaking Incubator from UniEquip (Munich, Germany) (Munchen, DE) for 1 h. Next, 250 µL MS grade water was then added to induce phase separation, before centrifugation at 1000× g for 10 min. After phase separation, the upper organic phase was transferred to a new tube while the bottom layer was re-extracted with extraction solvents. After centrifugation, the two phases were carefully separated and dried at 4 • C in a SpeedVac Concentrator until the solvents were fully evaporated. A process blank was also prepared with the same preparation steps, with the use of ultrapure water in place of blood plasma.
Before the sample run, dried non-polar extracts used for lipidomics analysis were reconstituted in 300 µL of isopropanol, while dried polar extracts used for the metabolomics analysis were reconstituted in 50 µL 50% acetonitrile. Quality control (QC) samples were prepared by mixing 10 µL of each plasma sample and aliquots of the pooled QC sample were subjected to the same sample preparation procedure. This pooled QC sample also provided a representation of all the analytes present in the samples.

LC-QTOF Analysis
For both the polar metabolite and non-polar lipid fractions, the samples were analysed on an Agilent 6540 UHD Accurate-Mass Q-TOF LC/MS coupled with an Agilent 1290 Infinity LC system and operated using the B.08.00 Agilent Mass Hunter Software from Agilent Technologies.
For the chromatographic separation, the columns and LC conditions that were used were different for the polar metabolite fraction versus the non-polar lipid fraction. For the polar metabolite fraction, 5 µL injections on a Kinetex Polar C18 column (100 mm × 2.1 mm, 2.6 µm) from Phenomenex (Torrance, CA, USA) at a flow rate of 0.3 mL/min at 40 • C were done. Mobile phase A consisted of 0.1% formic acid in water, and mobile phase B was 0.1% formic acid in acetonitrile. The gradient was as follows: 5% B at 0-2 min, 77% B at 8 min, 95% B at 12-14 min, and 5% B at 14.2-18 min. As for the non-polar lipid fractions, chromatographic separations of 5 µL injections were performed on an XSelect Acquity CSH C18 column (100 mm × 2.1 mm, 3.5 µm) from Waters Corporation (Milford, MA, USA) at a flow rate of 0.35 mL/min at 40 • C. Mobile phase A consisted of 60:40 acetonitrile/water with 10mM ammonium formate and 0.1% formic acid, and mobile phase B consisted of 90:10 isopropanol/acetonitrile with 10mM ammonium formate and 0.1% formic acid. The gradient was as follows: 10% B at 0-2 min, 40% B at 4.5 min, 79% B at 9.5 min, 95% B at 18.5-21.5 min, and 5% B at 21.6-25.1 min.
For the MS data acquisition after chromatographic separation, compounds were ionised with both polarities (i.e., positive and negative) using Dual Agilent Jet Stream Electrospray Ionisation with the following source conditions for the Q-TOF mass spectrometer system: a positive capillary voltage of 4 kV in positive ion mode, negative capillary voltage of 3.5 kV in negative ion mode, drying gas flow of 10 L/min, and the gas temperature of 320 • C. The nebuliser pressure was set at 35 psi. The Fragmentor voltage of the method was set at 75 V. Data was acquired over a mass range of 50-1000 m/z. Agilent Masshunter Auto MS/MS mode were used for MS/MS data acquisition. A collision energy of 10, 20, and 40 eV was applied for Data Dependent Acquisition (DDA) for both positive and negative modes, with an acquisition rate of 4 spectra/s. A maximum of 10 precursors were selected for each MS cycle for MS/MS acquisition. 'Iterative' mode was applied, in which MS precursors selected for a particular run were not selected for subsequent duplicate MS/MS runs.
After data acquisition, raw data files were converted to the mzXML format with the open-source ProteoWizard software [20], before processing with the XCMS Online program (https://xcmsonline.scripps.edu (accessed on 30 June 2022)). In XCMS, blanks, QC and plasma samples were uploaded for pre-treatment, including procedures such as peak picking and grouping, retention time correction, peak alignment, and annotation of isotopes and adducts. Specifically, we followed and adapted the recommended parameters for the UPLC/QTOF mode for untargeted metabolomics works available in XCMS to ensure accuracy and consistency of the integration of peaks. For feature detection, 15 ppm was set as the maximal tolerated m/z deviation in consecutive scans, with minimum and maximum peak widths of chromatographic peaks being 5 and 20 s respectively. Retention time was corrected by using 1 m/z step size for profile generation from the raw data files using the 'obiwarp' method. Alignment of the peaks across samples was done with the following parameters: (a) 5 s allowance for retention time deviations and (b) half of the samples in at least one sample group must have the peak for it to be valid. The adducts included during the annotation step include [M+H] + and [M+Na] + for the positive mode, and [M-H] − and [M+FA-H] − , which are potentially key ESI adducts formed based on the mobile phases used. After checking for the accurate and consistent integration of peaks across samples, information containing the sample retention times, m/z, peak abundance and isotope ions was downloaded from XCMS to be used for further data clean-up. Isotope ions were filtered and only the main isotope was retained for each feature.

Identification of Highly Contributing Features and Pathway Analysis
Identification of highly contributing features for the untargeted metabolomics and lipidomics datasets was based on the criteria considering both univariate and multivariate results obtained from MetaboAnalyst (http://www.metaboanalyst.ca (accessed on 30 June 2022)). Univariate tests for comparing mean peak areas between subject groups were done using the Mann-Whitney test. p-values obtained were false-discovery rate (FDR)-adjusted. Principal component analysis (PCA) and projection to latent structuresdiscriminant analysis (PLS-DA) were then performed. The PLS-DA model was validated with 10-fold cross-validation and a 100-iteration permutation test. The criteria for selecting important features for identification were as follows: a variable importance in projection (VIP) score > 2.0 (includes metabolites/lipids identified in both principal components (PC) 1 and 2) and a univariate test p-value of <0.05 with fold-change >1.5 or <0.67.
For those significant features based on the aforementioned criteria, a further structural identification step was required. For that purpose, a 3-step process was employed for putative identification. Firstly, matching of accurate m/z (accuracy threshold = 10 ppm) to online databases was done. The Human Metabolome Database (HMDB; http://hmdb.ca (accessed on 30 June 2022)), was used for the majority of the matches and class identifications. This was supplemented with matches using the METLIN database (https://metlin.scripps.edu (accessed on 30 June 2022)). Secondly, MS/MS spectra were obtained either from DDA- based results or additional MS/MS runs, and final identities were determined by comparing our MS/MS spectra with the in silico fragmentation spectra of potential molecular identities shortlisted from the first step. Thirdly, a final validation was done by checking that the retention time was reasonable based on the polarity of a proposed structure. As such, a level 2 metabolite identification confidence level may be achieved, according to the Metabolite Standards Initiative [21]. For features with no matching MS/MS spectra, they were grouped under the level 4 identification level, and are left out from network analysis due to insufficient molecular information.

Glycomics and Metallomics Analysis and Data Acquisition
N-glycan extracts and elemental digestates from blood plasma samples were obtained as described in our previous works [10,11], and also presented in the "Supplementary Methods" section found within the Supplementary Material file available online. Briefly, 2 µL blood plasma was treated with PNGase F and digestion buffer for the release of N-glycans. Before instrumental analysis, the N-glycans were reduced with 2-picolane borane and labelled with a fluorescence tag (8-aminopyrene-1,3,6-trisulfonic acid (APTS)). The APTS-labelled glycans were then purified via the use of a magnetic stand-enabled solid phase extraction step. As for the elemental digestates, briefly, equi-volumes of ultrapure nitric acid and ultrapure water were added to 80 µL blood plasma. Acid digestion proceeded on a hotplate at 98 • C for 2 h. The glycan extracts and elemental digestates were then subsequently diluted to suitable levels and analysed by capillary electrophoresislaser-induced fluorescence (CE-LIF) and inductively coupled plasma-mass spectrometry (ICP-MS) respectively, again as described in our previous works [10,11]. For the glycomics dataset, peak areas were obtained for each N-glycan variable, while concentrations were obtained for the metallomic variables via the building of external calibration curves with appropriate standards.

Data Processing
Data pre-processing was done to filter out low-quality metabolite/lipid/glycan/element peaks/concentrations via the following criteria: (1) peaks with less than 3-fold average intensity in the samples compared to the blanks; (2) peaks present in less than 75% of the samples in either patient/control group; (3) relative standard deviations (RSD) of intensities more than 20% in pooled QC samples (except for RSD threshold of 30% for the untargeted metabolomics portion). After filtering, any missing values were replaced with half the minimum value of each respective feature.
For the glycomics dataset, normalisation of peak areas was done with respect to the sum of peak areas as per the glycomics convention. On the other hand, for the metabolomics and lipidomics datasets, normalisation was done with respect to the QC samples by the internal standard [22]. No normalisation was done for the metallomics dataset as accurate concentrations were measured.
To prepare the four omics datasets for cross-omics multivariate analysis, log-transformation and scaling via Pareto scaling were done as appropriate.

Statistical Analysis and Multi-Omics Integrative Analysis
All statistical analyses were conducted in R (version 4.04, R Core Team, Vienna, Austria) [23]. The normality of all metabolites, lipids, and glycan features was first assessed via the Kolmorov-Smirnov test. Since non-normal distributions (p < 0.05 based on the Kolmorov-Smirnov test) were found with some of the features, non-parametric tests were done subsequently. Baseline characteristics among patients and controls were compared using the Mann-Whitney test for continuous variables, and Fisher's Exact test for categorical variables. Differences in plasma levels of the various omics features between patient groups were also assessed using the Mann-Whitney test (with a Benjamini-Hochberg multiplicity correction if comparisons between multiple patient groups were endeavoured).
For the multi-omics analysis, the 'mixOmics' package (version 6.0.0) was used and implemented in R [24]. The DIABLO framework from the 'mixOmics' package was used to integrate the metabolomics and glycomics datasets via block-PLS-DA based on sparse generalised canonical correlation analysis (CCA). The PLS-DA model was tuned and performed with 10-fold validation and using centroid distance for estimating the error rate. The visualisation of the relationships between variables from all four omics layers (based on results from the DIABLO framework's block-PLS-DA analysis) was done by plotting relevant correlation plots and relevance network maps using the 'mixOmics' package.
Classification modelling was done on MetaboAnalyst (http://www.metaboanalyst.ca (accessed on 30 June 2022)). The receiving operating characteristic (ROC) curve for biomarker identification and performance evaluation was generated based on random forest classification modelling with repeated random sub-sampling cross-validation (30 iterations; each iteration uses 2/3 samples for feature selection and model training, and the remaining 1/3 for testing). Table 1 summarises the baseline characteristics such as demographic information (i.e., age, gender, race), medical history related to the disease, as well as traditional cardiovascular risk factors and predictors of AMI and adverse events of the study populations. The following discussion regarding the demographic information and baseline characteristics of the samples had been mentioned in our previous works [10,11]. They are briefly reiterated in this work. Firstly, we observed no significant differences in the distributions of patient ages between the AMI and healthy groups, as they were age-matched in the study. On the other hand, other demographic characteristics such as gender, smoking status, and race were found to be significantly different between the two groups (p < 0.0001). Secondly, with reference to Table 1, we confirmed that the medical history for CHD was also closely matched between AMI and healthy groups, with the only exception being dyslipidaemia (p < 0.001). Differences in high-sensitive Troponin T (hsTnT) and N-terminal-pro hormone brain natriuretic peptides (NTproBNP) levels were found to be highly significant between patient groups (p < 0.0001). Nonetheless, this observation was expected as they are biomarkers of AMI, and raised troponin was also an inclusion criterion for the AMI cases in this study. As for other traditional cardiovascular risk factors, several of them (diastolic blood pressure (DBP), systolic blood pressure (SBP), high-density lipoprotein-cholesterol (HDL-C), white blood cell count (WBC), and creatinine levels) are also statistically significant variables (p < 0.0001), which again was expected due to their associations with CHD.

Analytical Validation for Untargeted Metabolomics and Lipidomics Analysis
Prior to the formal analysis of study samples, the combination of the extraction method with our LC-QTOF analytical workflow for metabolomics and lipidomics analysis was validated by assessing instrumental analytical precision, and inter-day method repeatability. As for the glycomics and metallomics analytical workflows, again, they were validated and relevant figures-of-merit had been presented in our previous works [10,11].
Firstly, for checking the instrument's analytical precision, especially necessitated by LC-MS's known signal drift and batch issues, a single extract was injected 10 times in succession, in both positive and negative modes. RSDs of 7.9% (positive mode) and 5.5% (negative mode) were obtained for the 10 consecutive injections, indicating satisfactory instrument precision. As such, for the actual analysis of study samples, the following run conditions were set: (1) 10 consecutive injections of a QC sample prior to sample injections for conditioning, and (2) one QC injection between every 10 sample injections for quality control.
For assessing inter-day repeatability/precision, three independent sets of extractions were done using the same pooled plasma sample, and they were analysed over different days to evaluate any significant batch effects and uncertainties in the overall workflow. Extraction and analysis of both the metabolite and lipid fractions were done in both positive and negative modes (four sets in total). A threshold of up to 30% RSD was used for reference, as it is commonly accepted as a standard when filtering samples based on QC in metabolomics works [25]. For the metabolite fraction, 58.0% of peaks from the positive mode and 88.0% of total peaks from the negative mode had an RSD < 20%. As for the lipid fraction, median inter-day RSDs of 14.0% (positive mode) and 12.9% (negative mode) were achieved. Additionally, 79.0% of the total peaks from the positive mode and 86.2% of total peaks from the negative mode had RSDs of <20%, indicating great precision, which is beyond comparable with the standard. While the results for the polar metabolite fraction in the positive mode are less desirable, the overall results for all four sets gave satisfactory intra-day precisions. Based on these results, we thus also set the RSD threshold during QC filtering of the actual sample run to be at 30% for the metabolomics analysis and 20% for the lipidomics portion, i.e., features with more than 20% RSD in QC samples for the lipidomics portion were removed.

Generation of Untargeted Lipidomics and Metabolomics Datasets and Identification of Significant Features
Firstly, the human plasma untargeted metabolomics and lipidomics profiles were generated as described, based on the Matyash extraction method, and according to the run/QC conditions ascertained in the previous section. Metabolite/lipid features extracted using XCMS were first defined by their LC retention time and accurate mass-to-charge ratios (m/z), as their structural identities were unknown at this point. For the metabolite fraction, 6361 and 4228 features were extracted from the positive and negative modes' results, respectively. For the lipid fraction, 15,968 and 3366 features were extracted from XCMS from the positive and negative modes' results, respectively.
Secondly, data clean-up such as the removal of features with RSDs of >20 or 30% in the QC, as well as filtering after data normalisation, was performed as described in the methods Section 2.4. Based on these data processing steps, the following number of features were left for further statistical filtering: 319 (positive mode) and 1420 (negative mode) for the metabolite fraction, and 3292 (positive mode) and 973 (negative mode) for the lipid fraction.
Thirdly, these validated and "cleaned-up" datasets were then further filtered based on statistical criteria as described in the methods Section 2.2.4. This process served to generate a focused list of highly contributing features in differentiating AMI and healthy patients, for further structural identification. A total of 40 metabolomic and 99 lipidomic significant (p < 0.05) features were filtered. Since one of the criteria for the filtering of significant features is based on the VIP score generated from PLS-DA modelling, the models' performances were checked. Figure 1A,B show the constructed three-dimensional PLS-DA models (with 10-fold cross-validation) for the metabolomics and lipidomics portions, respectively. Figure 1C,D illustrate the various performance measures (accuracy, R 2 , and Q 2 ) against the number of components used in the PLS-DA models. For the metabolomics analysis, the cumulated R 2 and Q 2 values at the optimal number of components were 0.79 and 0.22 for metabolomics, and 0.61 and 0.43 for the lipidomics analysis. The moderately high R 2 and closeness of Q 2 to R 2 values for the lipidomics analysis demonstrate good predictive relevance and validity, based on PLS model performance standards [26]. On the other hand, for the metabolomics portion, there may be possible overfitting, and discretion may be advised when interpreting results generated from PLS-DA-based models for our metabolomics dataset.
Finally, a three-step structural identification process as described in Section 2.2.4 was applied to annotate the identities of the 40 metabolomic and 98 lipidomic features. In this untargeted metabolomics/lipidomics study, a Level 2 metabolite identification confidence level was achieved for most of the significant features (according to the Metabolite Standards Initiative) [21]. Unfortunately, there was difficulty identifying potential structures based on an accurate m/z and in the matching of the MS/MS spectral patterns for multiple features. Subsequent analyses and interpretations were thus limited to what could be surmised from the identified features.
the other hand, for the metabolomics portion, there may be possible overfitting, and discretion may be advised when interpreting results generated from PLS-DA-based models for our metabolomics dataset. Finally, a three-step structural identification process as described in Section 2.2.4 was applied to annotate the identities of the 40 metabolomic and 98 lipidomic features. In this untargeted metabolomics/lipidomics study, a Level 2 metabolite identification confidence level was achieved for most of the significant features (according to the Metabolite Standards Initiative) [21]. Unfortunately, there was difficulty identifying potential structures based on an accurate m/z and in the matching of the MS/MS spectral patterns for multiple features. Subsequent analyses and interpretations were thus limited to what could be surmised from the identified features.
As such, a total of 76 (27 metabolomic and 48 lipidomic) features could be successfully annotated, and their most pertinent details are given in Table 2 (metabolomics) and Table 3 (lipidomics). As for those unannotated features, they continued to be included in classification modelling. LC-MS results and statistical information on these unannotated significant features are given in Table S1 (metabolomics) and Table S2 (lipidomics) in the Supplementary File, together with the additional structural information and the VIP scores of the identified features. However, they were excluded from multi-omics integration due to insufficient information. As such, a total of 76 (27 metabolomic and 48 lipidomic) features could be successfully annotated, and their most pertinent details are given in Table 2 (metabolomics) and Table 3 (lipidomics). As for those unannotated features, they continued to be included in classification modelling. LC-MS results and statistical information on these unannotated significant features are given in Table S1 (metabolomics) and Table S2 (lipidomics) in the Supplementary File, together with the additional structural information and the VIP scores of the identified features. However, they were excluded from multi-omics integration due to insufficient information.  a "UP" denotes presence of higher amounts of the analyte (i.e., up-regulated) in AMI patients as compared to healthy volunteers, while "DOWN" denotes the presence of lower amounts of the analyte (i.e., down-regulated) in AMI patients. b Fold-change was calculated by dividing the mean peak area of analyte in AMI patients over the mean peak area of analyte in healthy volunteers. A fold-change > 1 also implies up-regulation in AMI patients, while a fold change < 1 implies down-regulation. c p-value was FDR-adjusted.  a "UP" denotes presence of higher amounts of the analyte (i.e., up-regulated) in AMI patients as compared to healthy volunteers, vice versa. b Fold-change was calculated by dividing the mean peak area of analyte in AMI patients over the mean peak area of analyte in healthy volunteers. A fold-change > 1 also implies up-regulation in AMI patients, while a fold change < 1 implies down-regulation. c p-value was FDR-adjusted.
First, the metabolite/lipid classes that were significantly altered in AMI were assessed. With reference to Table 2, it is shown that there is a spread of metabolite classes involved with AMI risk. They include semi-polar compounds such as purine nucleosides and indole/derivatives, to non-polar compounds such as steroids/derivatives (overlap with lipidomics' coverage). We note that a few of the metabolites found at higher levels in AMI patients do not seem to be endogenous/non-naturally occurring in the human body. They include (1) perfluorododecanoic acid, a polyfluoroalkyl chemical found in stainresistant furniture, grease-resistant paper, kitchen wares, etc., whereby its exposure and accumulation in human serum has been reported to be associated with cardiovascular risk and specifically angina pectoris (chest pain) [44], and (2) carbosulfan, a pesticide, and a derivative of carbofuran, the exposure to which has several reported associations with CVD risk and AMI [45]. This untargeted metabolomics study thus augments the evidence pointing at such chemical exposures as potential CVD risk factors.
Second, unlike the metabolomics portion, great coverage for lipidomics was observed (Table 3), with significant lipid features constituting five out of the eight main lipid categories defined by LIPID MAPS [46]. We thus further investigated the distribution of lipid classes that were significantly altered (p < 0.05) in AMI. The distribution of lipid classes was plotted as a pie chart as shown in Figure 2. Based on Figure 2, glycerophospholipids (27%), fatty acids (16%), and triacylglycerols (14%) are the top three abundant lipid classes that were altered in AMI. As shown in Table 3, these lipid classes display distinct regulation patterns: glycerophospholipids were generally downregulated in AMI patients, all detected fatty acids were downregulated, while all triacylglycerols (as well as most diacylglycerols) were unanimously upregulated in AMI.   As for the down-regulation trend found in glycerophospholipid metabolism, this has already been reported to play a prominent role in CHD progression recently [47]. According to Chen et al., plasmalogens (including various phosphatidylcholines (PC) and phosphatidylethanolamines (PE) identified in our work and classified under glycerophospholipids in Table 3) were proposed to be protective against atherosclerosis [47]. Moreover, since plasmalogens are also notably more susceptible to oxidation under oxidative stress, low levels of such glycerophospholipids can be considered as a biomarker of oxidative stress and the negative actions of reactive oxygen species, which drive or can be driven by disease progression [48,49]. Additionally, glycerophospholipids have been postulated to be possible inflammatory mediators as glycerophospholipid metabolic pathway was associated with low-grade inflammatory states and general systemic-immune inflammatory states [50]. This link suggests the potential and value of integrating this lipidomics dataset with glycomics, to identify interactive metabolic pathways commonly associated with inflammation.
As for the simultaneous down-regulation of fatty acids and up-regulation of di-and tri-acylglycerols, this can be explained by the imbalance between fatty acid uptake and oxidation in AMI patients. This is since these two categories of lipids share a dynamic balance, as di-and triacylglycerols are lipid intermediates (alongside phospholipids and sphingolipids) that will accumulate within the myocardium when fatty acid oxidation (FAO) is unable to match the fatty acid delivered to the heart [51,52]. Overall, this simultaneous dysregulation of fatty acids and di-/triacylglycerols (as well as up-regulation of other lipid intermediates attributed to the dysfunction of FAO) served as additional evidence for the ongoing discussion on "lipotoxicity" in the human heart. As summarised by a key review by Schulze et al., there have been various supporting reports on the role of intracellular lipid accumulation leading to chronic states with ATP production dysfunction and energy depletion of the failing myocardium [53]. Such reports range from studies of animal models, studies of patients with inborn errors of FAO who develop cardiac abnormalities (e.g., sudden cardiac death, cardiac/skeletal myopathies, insulin resistance), and observational studies reporting the association of cardiac lipid accumulation with obesity and/or metabolic cardiovascular complications (e.g., diabetes mellitus), as well as more advanced imaging studies showing increased intramyocardial lipid content in patients with heart failure [53]. However, while it seems apparent that lipid accumulation is linked to the failing myocardium and thereby also AMI onset and progression, what initiates or drives it is still unclear (i.e., whether it is due to elevated fatty acid uptake, elevated di-/triacylglycerol synthesis, impaired degradation of lipids, or combinations of the above) [53].

Multi-Omics Integration and Analysis of Cross-Omics Relationships
Next, the annotated metabolomic (number of features, n = 27) and lipidomic (n = 48) profiles, along with the glycomic (n = 37; [10]) and metallomic (n = 30; [11]) profiles of the same set of AMI and healthy samples were integrated and analysed. The integration method used here works by identifying a linear combination of maximally correlated features across the four omics datasets, via utilising both block-PLS-DAs based on sparse generalised canonical correlation analysis [54]. A series of correlation plots (Figure 3) display the overall correlations amongst the omics blocks.
Primarily, Figure 3A showed that the latent components of each omics block were highly correlated between metabolomics and lipidomics (r = 0.67), and were moderately correlated between metallomics and lipidomics (r = 0.37), between glycomics and lipidomics (r = 0.33), and between glycomics and lipidomics (r = 0.33). This highlights the ability in modelling satisfactory agreement between the datasets. The clustering of the AMI and healthy samples, while not distinct, are inhomogeneous and thus preliminarily demonstrates the ability of the integrated omics models to discriminate the outcome of interest as well. However, low correlations (r < 0.33) were observed between glycomics and metallomics, and between glycomics and metabolomics data. Nevertheless, beyond general agreement among the various omics blocks in disease classification, further assessment for specific components highlights the correlations between individual features from each omics block ( Figure 3B). Additionally, the circos plot ( Figure 3C) displays the significantly contributing features across all four omics blocks in a circle, with links between each omics feature indicating at least a moderate correlation (r cut-off = 0.33; p-value cut-off at 0.05). The density of both red (positive correlation) and green (negative correlation) links in Figure 3C demonstrates the high-volume interactions amongst biomolecules from all four omics layers.  Next, through the generation of a series of relevance network plots amongst various omics combinations, various cross-omics clusters of interest have been identified. Figures S1-S9 in the Supplementary File display the individual relevance network plots (including one tetraomics network, two tri-omics networks, and all six combinations of the bi-omics networks; a high correlation cut-off at r ≥ 0.5 and p-value < 0.05). On the other hand, Figures 4 and 5 display the main relevance network plots, which reveal pertinent clusters.    Firstly, with reference to Figure 4A which shows the key cluster (r ≥ 0.67) from all four omics blocks, it mainly revolves around the FA2G1 [6] glycan and A2G2S2(2,6) glycan. FA2G1 [6] is a core-fucosylated and unsialylated glycan, while A2S2S2(2,6) is a fully-sialylated biantennary and noncore-fucosylated glycan. These glycans, which are abundantly found on the Fc region of human IgG, are both also correlated with the Cu/Se ratio, and lipids of various classes (acyl carnitine, glycerophospholipid, diacylglycerol, cholesterol and cholesterol ester). However, it is interesting to note that the directions of correlations are opposite for the two glycans-FA2G1 [6] is positively correlated with various lipids and negatively correlated with the Cu/Se ratio, while A2G2S2(2,6) is negatively correlated with various lipids but positively correlated with the Cu/Se ratio. The other two glycans found in the key cluster (FA2G2 and A2G2S1(2,6)) follow the same trends in terms of glycosylation traits versus the correlation direction with both lipids and the Cu/Se ratio. In fact, when looking at Figure 4B and specific bi-omics (glycomics-lipidomics and glycomics-metallomics) relevance networks at r = 0.4 cut-offs in Figure 5's Clusters 1, 2, and 5, the aforementioned correlation trend is consistent and even more apparent. This reveals intriguing insights regarding the immunoregulatory capacity of IgG in AMI and its possible interplay with lipid membrane domains.
We also looked at how the various detected glycans are linked to IgG-FcγR functions, as well as how IgG-FcγR interactions could be inter-regulated by lipid metabolism. Fc region glycans can directly influence the affinity of IgGs to FcγRs and thereby change their ability to recruit immune effector cells to activate distinct immunomodulatory pathways, either by changing the conformation of the Fc region or via glycan-glycan interactions [55,56]. For example, N-glycans with core fucose or terminal sialic acid present reduce the affinity of IgGs to Fcγ receptors in general [55], but specifically, unsialylated IgGs were reported to primarily interact with type I Fc receptors, which include all FcγRs. On the other hand, sialylated IgGs interacted primarily with type II FcγRs, which are mainly expressed in dendritic cells (DC) [55,57]. Next, lipid metabolism, and thereby lipid content/composition in immune cells, may be inter-regulated by IgG-FcγR interactions [58]. Now, as seen from the correlation trends with the various lipids seen in Figure 4, most of the lipids (11 lipids) linked with the glycans associated with the IgG-FcγR immune complex are glycerophospholipids (one of the three main lipid classes composing plasma membranes). Furthermore, it has been previously reported that pro-inflammatory activated lipid-rich human monocytederived DCs, with their remodelled plasma phospholipid and cholesterol content, might be involved in atherosclerosis [59]. All of these collaboratively suggest that our previous postulation regarding AMI atherosclerosis involving a switch from a pro-inflammatory to an anti-inflammatory state could have also been modulated by lipid metabolism. On the other hand, we postulate that it could also have simultaneously triggered a change in plasma membrane composition to elicit desired immunoregulatory responses as a response to atherosclerotic onset or progression. The involvement of Cu/Se ratio in this equation remains unclear, although Se has been shown to regulate murine lipid levels via restoring lipid peroxidation and ameliorating disruptions of membrane dynamics [60]. It is likely that Cu and Se are indirectly involved in this relevance network due to their strong but general role in inflammation and immunity [61]. This is further shown as we look at Figure 5's bi-omics (metallomics-lipidomics) Cluster 3, where Se is only positively correlated with a fatty acid and its acylcarnitine (oleic acid and o-linoleoylcarnitine) with r ≥ 0.4. This also highlights the need to investigate the relationship between copper and lipids, which has 15 links r ≥ 0.4, yet lacks evidence of their interaction in CVD pathophysiology.
Beyond the key cluster identified in the multi-omics network in Figure 4A, a side cluster was also revealed when the correlation cut-off was lowered to r = 0.5 ( Figure 4B). This side cluster revolves around Ca and its interactions with two unsialylated and non-fucosylated glycans (A2BG1 [3] and M6) and three phosphatidylethanolamine (PE) lipids under the glycerophospholipid class. It is noted that Ca's negative correlation with the PEs seemed to be specific, as correlations with other types of lipids, glycerophospholipids or otherwise, were absent at this correlation cut-off. Thus, while Ca ions have been reported to interact with lipid membranes as a whole (as evidenced by the presence of Ca binding sites in lipid bilayers and changes in membrane bilayer structures with Ca binding, etc.) [62], the lack of Ca-PC and Ca-phosphatidylserine (PS) correlations suggest a different and more specific mode of interaction with the PEs. Here, our evidence (along with the correlations with Na and Ca×P products) points towards the role of the sarcoplasmic reticulum Ca-ATPase (SRCA) in heart disease. The PE headgroups had been reported to modulate SRCA function by facilitating "dynamic structural changes involving the phosphorylation domain, which enhances the catalytic function of Ca-ATPase" [63]. This explains the link between the PEs with Ca, and thereby also possibly Na, since Na together with Na/K-ATPase play interactive roles in the regulation of intracellular Ca as well as Ca stores at the SRCA [64]. All of these have direct consequences on cardiac myocyte functions and are linked with heart contractions [64]. Since SRCA is a potent therapeutic target for CVD [65,66], our work suggests to also scrutinise specific intermolecular interactions that result from the phosphoethanolamine moiety of PEs. Moreover, it is also interesting to see from Figure 5 Cluster 4, that PEs not just have strong correlations with Ca and Na, but also moderate correlations (r ≥ 0.4) with many other elements as well (e.g., essential elements such as Al, Cu, Mg, and Zn). It thus highlights the multi-faceted interactions that PEs have with such essential elements, and their potential importance in phospholipid homeostasis in the cardiovascular system [67,68]. Lastly, further studies should also investigate how a terminal GlcNAc glycan and a high-mannose glycan as identified in this work may be involved, e.g., whether inhibition of the ATPases impairs N-glycan expression/function [69], or vice versa [70].

Classification Modelling and Performance of the Multi-Omics Model
As a culmination of all four omics studies, the performance of a multi-omics classifier in discriminating between AMI and healthy was compared with that of the individual omics. First, for quick visualisation of the classification performances, PLS-DA modelling was applied to depict the classifiers' ability to cluster/separate AMI versus controls. For the multi-omics classifier, it was built based on the block-PLS-DA analysis described earlier.
As shown in Figure 6A, the ability to cluster/separate AMI and healthy samples via the consensus multi-omics classifier is seemingly better than the individual PLS-DA models.
Secondly, besides the quick visualisation based on supervised PLS-DA models, an unsupervised method was also applied to observe if the discrimination between AMI and controls was indeed the main source of variability measured in omics variables. For that purpose, hierarchical cluster analysis was done. As shown in Figure 6B, the use of all four omics datasets was able to render the clustering of AMI from healthy relatively well. This can be seen from the grouping of the blue rows (AMI) mostly at the bottom, while the orange rows (controls) are mostly on top. Beyond that, the intermixing of features from the four omics layers could also be observed from the top dendrogram in Figure 6B. This was shown from the interspersed spread of the colored columns (showing features from different omics layers). The significance of this observation is that it possibly denotes an interactive relationship among the four omics in AMI discrimination.
Thirdly, to more objectively assess the differences in performance between individual and multi-omics classifiers, AUC ROC s were generated. In this case, random forest models were chosen as the standard modelling technique for comparison purposes, to account for any nonlinear interactions across the omics layers. As a note, random forest modelling was also chosen instead of PLS-DA modelling as PLS-DA is a linear-based method and may not be sensitive to the nonlinear cross-omics relationships. Again, the random forest modelling was done with repeated random sub-sampling cross-validation (30 iterations; each iteration uses two-thirds of the samples for feature selection and model training, and the remaining one-third for testing). No external validation was done. The results are summarised in Table 4 below. The performance of the multi-omics classifier (AUC ROC = 0.953) was higher than all the individual omics classifiers. The multi-omics classifier's lower confidence limit (CI 0.911-0.987) was also even higher than all the upper confidence limits of the individual classifiers, except for metabolomics. This exemplifies the synergistic effect when these four omics datasets were cohesively used to build a disease classification model. Although, there is one caveat here: the observation that the multi-omics classifier achieving a higher AUC ROC value is not unexpected mathematically. The most important question remains to be whether the multi-omics classifier's performance can be reproducible when tested in a large and independent cohort of AMI patients. This is the litmus test.   As a whole, despite the limitation that an external validation was not achieved in this portion, the potential of the multi-omics classifier, having synergistically enhanced performance, cannot be denied.

Limitations and Further Work
One of the biggest analytical challenges still pertains to the coverage in metabolomics analyses. As such, for future large-scale multi-omics studies, researchers should consider performing targeted metabolomics with a defined list of pathways and metabolites of interest [71], or to attempt, as comprehensively as possible, an untargeted metabolomics study by perhaps performing both GC-MS and LC-MS analyses for covering both volatile and non-volatiles [71], as well as running LC-MS with both normal and reverse-phased columns (e.g., HILIC + T3 + C8 combination) [72], after testing at different mobile phase pH conditions (acidic, neutral, basic) for optimal/complementary selectivity [73].
Secondly, in this study, only age was matched between the case and control groups. However, other demographic characteristics (gender, race, smoking status), as well as cardiovascular risk factors are also important factors that may be translated to measurable differences in omics profiles as well. Unfortunately, our study size was severely limited by patients' recruitment in local hospitals. Thus, we could only match the cases and controls based on age. Moreover, the inequity in distributions was further compounded by the fact that AMI patients are usually males. Additionally, adjustment for confounding factors by statistical analysis cannot completely eliminate the effect of mismatching between the two groups on the results. Therefore, in future work, sensitivity analyses (e.g., by removing some female controls and checking if the results would be affected) should be attempted. Otherwise, studying larger cohorts with appropriate matching would be ideal.
Lastly, we also note the challenges pertaining to visualisation capacities in our work, which were limited by the functions achievable using the available open-sourced software/packages. It would be greatly beneficial for any such omics studies in the future to have further developments in bioinformatics tools to streamline and hasten the data processing, analysis, and visualisation works.

Conclusions
Overall, an integrative multi-omics study of omics beyond the central dogma was attempted on a cohort of AMI and healthy patients to reveal intriguing biological insights pertaining to the cross-omics mechanisms associated with AMI pathophysiology. They included elucidating the potential roles of glycerophospholipids in N-glycan-modulated IgG-FcγR's immunoregulatory functions, the importance of SRCA in CVD, and the contributions of PEs for SRCA functions. Moreover, combining the four omics datasets synergistically enhanced the classifier performance in discriminating between AMI and healthy patients. All these discoveries were otherwise not attainable when these omics were singularly and independently analysed. The transitioning of a single omics to a multi-omics strategy is therefore supported by the work presented in this study.
It should also be noted that the significant contributions of the respective state-of-theart analytical workflows employed in this work enabled efficient multi-omics research. We highlight that the CE-LIF-based glycomics workflow selected in this multi-omics study used only 2 µL of blood plasma, while the untargeted metabolomics and lipidomics datasets generated via simultaneous extraction used only 50 µL plasma. Furthermore, even though the metallomics workflow required acid digestion, which makes it hard to integrate with other omics extraction workflows, the glycomics and metabolomics/lipidomics extractions could potentially be unified. One of our recent works [74], which demonstrates the validity of a simultaneous polar metabolite and N-glycan extraction bi-omics workflow, supports such an outlook.

Data Availability Statement:
The data presented in this study are available in our previous glycomics and metallomics works [10,11] and in the supplementary material for the metabolomics/lipidomics data.