Bioinformatic Profiling Identifies a Fatty Acid Metabolism-Related Gene Risk Signature for Malignancy, Prognosis, and Immune Phenotype of Glioma

Cancer cells commonly have metabolic abnormalities. Aside from altered glucose and amino acid metabolism, cancers cells often share the attribute of fatty acid metabolic alterations. However, fatty acid metabolism related-gene set has not been systematically investigated in gliomas. Here, we provide a bioinformatic profiling of the fatty acid catabolic metabolism-related gene risk signature for the malignancy, prognosis and immune phenotype of glioma. In this study, a cohort of 325 patients with whole genome RNA-seq expression data from the Chinese Glioma Genome Atlas (CGGA) dataset was used as training set, while another cohort of 667 patients from The Cancer Genome Atlas (TCGA) dataset was used as validating set. After confirmed that fatty acid catabolic metabolism-related gene set could distinguish clinicopathological features of gliomas, we used LASSO regression analysis to develop a fatty-acid metabolism-related gene risk signature for glioma. This 8-gene risk signature was found to be a good predictor of clinical and molecular features involved in the malignancy of gliomas. We also identified that this 8-gene risk signature had high prognostic values in patients with gliomas. Correlation analysis showed that our risk signature was closely associated with the immune cells involved in the microenvironment of glioma. Furthermore, the fatty acid catabolic metabolism-related gene risk signature was also found to be significantly correlated with immune checkpoint members B7-H3 and Tim-3. In summary, we have identified a fatty acid metabolism-related gene risk signature for malignancy, prognosis, and immune phenotype of glioma; and our study might contribute to better understanding of metabolic pathways and further developing of novel therapeutic approaches for gliomas.


Introduction
Abnormality of metabolism is the hallmark of cancer, and identification of the metabolic weaknesses of cancer cells has prompted new therapeutic approaches toward tumor treatments [1]. For example, increased glucose metabolism has been frequently seen as a characteristic of cancer cells [2]. Except for glucose metabolism, altered fatty acid metabolism in cancer cells has received increasing attention recently [3]. Fatty acid is the cornerstone of cell membrane formation, energy storage, and signaling molecule production in carcinogenesis; thus targeting at the pathway of fatty acid metabolism might inhibit rapid proliferation of the cancer cells [4].
In this study, we focus on one of the fatty acid metabolism-related gene setsthe catabolic metabolism of fatty acid gene set in gliomas. Gliomas are the most prevalent and malignant primary brain tumors in adults, and glioblastoma (GBM) is the most common and devastating type among all grades of glioma. Despite of the novel therapy of GBM, patients with GBM only have a median overall survival time of 14.6-16.7 months in clinical trials [5]. Metabolic profiling analysis of gliomas might contribute to better understanding of molecular pathways and further developing of novel therapies in gliomas [6]. For example, bioinformatic profiling had demonstrated a glucose metabolism-related risk signature and an amino acid metabolism-related risk signature were closely associated with malignancy and prognosis of glioma [7,8]. Through mass spectrometry, lipidomic signatures were also found to be approaches for classification of gliomas [9]. However, the role of the fatty acid metabolism-related gene set in glioma still remains unclear.
In our study, we firstly identified that the fatty acid catabolic metabolism-related gene set had the ability to distinguish clinicopathological features of gliomas. Then, we generated a fatty-acid catabolic metabolism-related gene risk signature in CGGA dataset, and further validated in TCGA dataset. We observed that our risk signature was associated with molecular features of gliomas and could serve as an independent prognostic factor for both all grade gliomas and GBM. Lastly, we also found that this risk signature was closely related to tumor infiltrating lymphocytes, which indicated an potential association between fatty acid metabolism and immune phenotype of gliomas. We believe that our results might provide a new insight for understanding the metabolic mechanism of gliomas.

Methods and Materials
2.1. Data Collection. The whole genome RNA-seq expression data and clinical information of 325 glioma patients from CGGA dataset (http://www.cgga.org.cn) were used as the training set [7,8]. RNA-seq data and clinical information from TCGA dataset (http://cancergenome.nih.gov) were used as validation set [10,11]. After dropping the samples with severely incomplete data (e.g. lack of critical clinical information such as overall survival time and IDH status), the eventual size of validation set was 667.

Bioinformatics Analysis.
Fatty acid catabolic metabolismrelated gene set (GO_FATTY_ACID_CATABOLIC_PRO-CESS), consisted of 73 genes in total, was extracted from Molecular Signatures Database v6.2 (http://www.broad.mit .edu/gsea/msigdb/) [7]. Most variable genes of the gene set, determined by their median absolute deviation (MAD), were selected for further consensus clustering [12]. Consensus clustering was carried out in R programming language (http://cran.r-project.org) for detecting the fatty acid catabolic metabolism-related glioma subgroups of the training set. The optimal number of the glioma clusters was determined by quantitative stability evidence in an unsupervised analysis. For evaluating the correlation between risk signature and the immune phenotype of glioma, the ESTIMATE package of R programming language was conducted to calculate the immune score which represented the infiltration of immune cell in the microenvironment of glioma [13]. The association between immune cells and glioma risk signature was analyzed by Gene Set Variation Analysis (GSVA) in R programming language as described by Zhang and colleagues [14,15].

Statistical Analysis.
Screened by univariate Cox regression analysis in the training set, 46 genes with prognostic significance (p < 0.05) in fatty acid catabolic metabolism-related gene set were selected for Least Absolute Shrinkage and Selection Operator (LASSO) regression analysis [16,17]. The generalized linear model produced by LASSO regression analysis was further analyzed with 10-fold cross validation in order to generate the minimum cross validated error. Based on the cross validation, 8 genes with their regression coefficients (Coef) were eventually achieved. Then the risk score for each patient in the training set and validation set was calculated by the following formula: All patients in the training set and validation set were then separated into high or low risk group according to the median risk score. Survival analysis based on the risk score was evaluated by Kaplan-Meier survival curve by using R programming language. Univariate and multivariate survival analysis was performed by using Cox proportional hazards model in R programming language. Other main statistical analysis including Student's t-test, chi-square test, and Pearson's test were all performed in R programming language. Statistical significance was considered at the level of p < 0:05.

Identification of an 8-Gene Risk Signature Associated with
Fatty Acid Catabolic Metabolism. Through univariate Cox regression analysis, 46 fatty acid catabolic metabolismrelated genes with prognostic significance (p < 0:05) were selected for further analysis in the training set. To identify the gene risk signature associated with fatty acid catabolic metabolism, these 46 genes were undergone the LASSO regression analysis. After cross validation, LASSO regression analysis generated 8 genes (ABCD1, ACADSB, CEL, CPT2, GCDH, NUDT19, PCCA, PEX13) in total as active covariates to calculate the risk score ( Figure 2 and Table 1). The signature risk score of each patient in the training set and validating set was then calculated with the LASSO regression coefficients and expression value of these 8 genes through formula (1) mentioned above.

Prognostic Value of 8-Gene Risk Signature in All Grade
Gliomas and Glioblastoma. In the CGGA dataset, Kaplan-Meier survival analysis revealed that patients in high risk group (n = 162) had a significantly poorer prognosis compared with patients in low risk group (n = 163; median OS: 10.5 vs 37.1 months; p < 0:001; Figure 4(a)). In TCGA dataset, patients in high risk group (n = 333) were also found to have much shorter overall survival times than patients in low risk group (n = 334, median OS: 15.5 vs 24.3 months; p < 0:001; Figure 4(d)). After taking important clinical and molecular factors (including age, gender, WHO grade, IDH status, chemotherapy and radiotherapy) into account, univariate and multivariate Cox analysis further demonstrated that this risk score was an independent prognostic factor of prognosis in CGGA dataset (Table 2). Cox proportional hazard model also found risk score could serve as an independent prognostic factor in TCGA dataset (Table 2). When focusing on the GBM phenotype, we also observed that patients in high risk group (n = 72) had a shorter OS than patients in low risk group (n = 72) of GBM phenotype in CGGA dataset(median OS: 8.5 vs 11.5 months; p < 0:001; Figure 4(b)). Results in TCGA dataset further validated the prognostic value of the risk signature in GBM phenotype (Figure 4(e)). In addition, the progression-free survival time of the high risk group was much shorter than low risk group in both CGGA and TCGA datasets (Figures 4(c)-4(f)). These

Discussion
Altered cancer metabolic processes such as glucose metabolism and amino acid metabolism are the hallmarks of cancers [18]. Metabolomic signatures can provide a better understanding of the molecular pathways of gliomas and offer great potentials for developing novel therapeutic approaches in glioma treatments [6]. For example, metabolic pathways including cysteine metabolism, nucleotides metabolism and 2-hydroxyglutarate have been demonstrated to be helpful for classification of gliomas [19,20]. Myo-inositol, an important osmolyte and substrate in phosphatidylinositol lipid      family, was also found to be associated with glioma grade [21]. Previous studies have also identified an amino acid metabolism-related gene risk signature and a glucose metabolism-related gene risk signature both have high prognostic value in glioma through bioinformatic analysis [7,8]. Nevertheless, the role of the fatty acid metabolism-related gene set in glioma still remains unclear.
To date, this is the first study introducing a fatty acid catabolic metabolism-related gene risk signature for the malignancy of gliomas and the survival of patients with gliomas. After confirmed that 73 fatty acid catabolic metabolism-related genes had the ability to distinguish the key clinicopathological features of gliomas in both CGGA and TCGA datasets, we built a fatty acid catabolic metabolism-related gene risk signature through LASSO regression analysis. Lower grade gliomas (LGG, WHO grade II and III) have a preferentially better clinical and prognostic characteristics compared with higher grade gliomas (HGG, WHO grade IV) [22]. IDH status and TCGA molecular subtypes(classical, mesenchymal proneural and neural) were key features for classification and prognosis of glioma [11,23]. Using our risk signature, patients in higher risk group tend to be associated with the higher grade, the more invasive TCGA molecular subtypes (classical and mesenchymal) and IDH wide type, which represented worse prognosis. By contrast, the lower grade, the less invasive TCGA subtypes (proneural and neural subtype) and IDH mutation were preferentially associated with patients in lower risk group. Furthermore, several oncometabolites have been confirmed to be the accumulated metabolic products of IDH mutation. For example, the abnormal glucose metabolite 2-hydroxyglutarate was a metabolic biomarker in gliomas and was useful in classification of gliomas [24]. However, the association between fatty acid metabolite and IDH mutation still remains unclear and our risk signature might provide some clues in this researching field. In summary, our results indicated that fatty acid catabolic metabolism might be involved in the progression of gliomas.
To further explored the potential role of fatty acid catabolic metabolism in gliomas, we evaluated the correlations between our risk signature and immune cell populations. Previous studies have demonstrated that tumor infiltrating lymphocytes(TILs), especially the CD4 + T cells and CD8 + T cells, are correlated with clinical prognosis in gliomas [25,26]. Our study found fatty acid metabolic gene risk signature was highly associated with CD4 + T cells and moderately correlated with CD8 + T cells, which indicated that patients of higher risk tend to have an unfavorable TILs pattern as previously demonstrated: lower level of CD8 + T cells combined with higher level CD4 + T cells [26]. In addition, our risk signature was also strongly correlated with innate immune cells including monocytes, macrophages, and NK cells, which might be consistent with previous studies that tryptophan metabolic adaptations in GBM were associated with evasion of innate immune system by tumors cells [27]. Our findings indicated that fatty acid catabolic metabolism-related gene risk signature might, to some extent, involve in the altered immune microenvironment of the gliomas. Furthermore, we evaluated the correlation of our gene risk signature and common immune checkpoint members. B7-H3 and Tim-3, the novel targets of immunotherapy against solid tumors [28,29], both showed close association with risk signature. Clinical trials of anti-B7-H3 (NCT02475213) and anti-Tim-3 (NCT02817633) are  carrying on, and our study showed the fatty acid catabolic metabolism-related gene risk signature might be a possible metabolic marker of the immunotherapy for gliomas. In our fatty acid catabolic metabolism-related signature, the protein encoded by ABCD1 is one of the superfamily of ATP-binding cassette transporters and is involved in the catabolic metabolism of very long chain fatty acid. ABCD1 is associated with altered white matter microvascular perfusion [30] and may contribute to the cell differentiation with parallel to tumorigenesis [31]. In previous study, ABCD1 transcript levels were overexpressed in breast cancer [32]. The protein encoded by CPT2 is a nuclear protein transported to the mitochondrial membrane. CPT2 plays a critical role in regulation of fatty acid oxidation [33] and might promote carcinogenesis in liver cancer by leading hepatocellular carcinoma to lipid-rich environment [34]. PCCA, encoding the mitochondrial enzyme Propionyl-CoA carboxylase, was also found to be altered in gastric and colorectal cancer [35]. Relationship of gliomas and other proteins encoded by the genes of our risk signature remains unclear and needs further researches. In summary, our fatty acid metabolic gene risk signature model may provide new insights into the carcinogenesis and therapeutic approaches of gliomas.

Conclusion
In conclusion, we identified a fatty acid catabolic metabolismrelated gene risk signature for the malignancy, prognosis and immune phenotype of gliomas, and our study might provide better understanding of fatty acid metabolic role in glioma carcinogenesis and in glioma immune phenotype.

Data Availability
All the data websites were confirmed to be available online.

Conflicts of Interest
The authors declare that they have no conflicts of interest.

Authors' Contributions
Ying Qi and Di Chen contributed equally to this work.

Supplementary Materials
Figure S1 Consensus matrixes for k = 4 to k = 10 of the 325 patients in the CGGA datasets by clustering the gene expression profile of the 73 fatty acid metabolism genes. Figure S2: Correlation of 8-gene risk signature and immune checkpoints related molecules in CGGA and TCGA datasets. Table S1: Clinicopathological features of two clusters classified by consensus clustering based on fatty acid catabolic metabolism-related gene set in CGGA dataset.