Enriched Metabolic Pathway Analysis of Cytokine-Treated RPE Cells

Metabolomics data were first processed using an in-house R script. Processed data were first inspected for sample clustering and variance difference via Principal Component Analysis (PCA) and Partial Least Square Discriminant Analysis (PLSDA), in which PLSDA plots were generated using the R package mixOmics v6.18.1 [28 (link)]. Raw p-values were adjusted using the Benjamini–Hochberg correction and then tested for significance using the Welch method for metabolite profiling using the R package Omu v1.0.6 [29 (link)]. Differential metabolite analysis between the TNFα- and TGFβ2-treated samples compared to control (untreated) samples for both cell and media metabolites was generated by the R package Omu using the count_fold_change function, with the log FoldChange (FC) >1.5 or <−1.5 (TNFα-treated H-RPE vs. control H-RPE), logFC > 1.0, and logFC < −1.0 (TGFβ2-treated H-RPE vs. control H-RPE), with a Benjamini–Hochberg adjusted p-value cut off < 0.05 [29 (link)].
The Enriched Pathway Network Analysis was generated using the KEGG database of Homo sapiens [30 (link)], which was loaded locally using the functions buildGraphFromKEGGREST() and buildDataFromGraph() in the R package FELLA v1.14.0 [31 (link)]. The lists of metabolites that were significantly different from the TNFα-treated vs. control and TGFβ2-treated vs. control were extracted. The KEGG compound hierarchy was assigned to the extracted list of metabolites for the two comparisons using the function defineCompounds() in FELLA by mapping the metabolite compounds against the loaded database. The KEGG-assigned metabolomics data were then used as input for pathway enrichment analysis using the undirected heat diffusion model followed by statistical normalization using Z-scores for sub-network analysis in the R package FELLA v1.14.0 [31 (link)]. For the pathway enrichment analysis, metabolites with logFC > 1.5 or <−1.5 (TNFα-treated vs. control), logFC > 1.0, and logFC < −1.0 (TGFβ2-treated vs. control), with a Benjamini–Hochberg adjusted p-value < 0.05 were included. The purpose of two different fold-change cutoffs was due to differences in metabolite changes between TNFα-treated and TGFβ2-treated H-RPE, enabling the capture of sufficient differential metabolites in each group for generating a comprehensive network analysis. The enrichment analysis outputs were then mapped to the Homo sapiens (hsa) KEGG graphs and subsequently used for network analysis. Optimal visualization of the metabolic network graphs was generated with the number of nodes limit (nlimit) of 250 for TNFα-treated H-RPE vs. control and nlimit of 160 for TGFβ2-treated H-RPE vs. control using the generateResultsGraph() in FELLA. KEGG IDs unmapped to the KEGG graphs were retrieved and searched against the KEGG pathway database (https://www.genome.jp/kegg/pathway.html, accessed on 28 December 2022).