{'fraction_non_alphanumeric': 0.03128120872394575, 'fraction_numerical': 0.014496955639315743, 'mean_word_length': 5.3016646366220055, 'pattern_counts': {'":': 0, '<': 6, '<?xml version=': 0, '>': 4, 'https://': 2, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Identification of Metabolic-Associated Genes for the Prediction of Colon and Rectal Adenocarcinoma

Background and Aim Uncontrolled proliferation is the most prominent biological feature of tumors. In order to rapidly proliferate, tumor cells regulate their metabolic behavior by controlling the expression of metabolism-related genes (MRGs) to maximize the utilization of available nutrients. In this study, we aimed to construct prognosis models for colorectal adenocarcinoma (COAD) and rectum adenocarcinoma (READ) using MRGs to predict the prognoses of patients. Methods We first acquired the gene expression profiles of COAD and READ from the TCGA database, and then utilized univariate Cox analysis, Lasso regression, and multivariable Cox analysis to identify the MRGs for risk models. Results Eight genes (CPT1C, PLCB2, PLA2G2D, GAMT, ENPP2, PIP4K2B, GPX3, and GSR) in the colon cancer risk model and six genes (TDO2, PKLR, GAMT, EARS2, ACO1, and WAS) in the rectal cancer risk model were identified successfully. Multivariate Cox analysis indicated that these two models could accurately and independently predict overall survival (OS) for patients with COAD or READ. Furthermore, functional enrichment analysis was used to identify the metabolism pathway of MRGs in the risk models and analyzed these genes comprehensively. Then, we verified the prognosis model in independent COAD cohorts (GSE17538) and detected the correlations of the protein expression levels of GSR and ENPP2 with prognosis for COAD or READ. Conclusion In this study, 14 MRGs were identified as potential prognostic biomarkers and therapeutic targets for colorectal cancer.


Introduction
Colorectal cancer is a common malignant tumor of the digestive tract, which is a serious threat to human health. 1 Given that the early symptoms of colorectal cancer are atypical, diagnosed cases are often in the middle and late stages of the disease. The development of surgical techniques and adjuvant chemotherapy and radiotherapy chemoradiotherapy has improved the survival time of colorectal cancer patients to a certain extent, but the 5-year survival rate of patients with colorectal cancer remains extremely low because of chemical drug resistance and other reasons. 2 Therefore, biomarkers that can provide accurate prognosis and predict therapeutic effects are beneficial for patients with colorectal cancer. The prognoses of patients with colorectal cancer considerably vary, and thus establishing an accurate and effective prognosis model is necessary. Such model can be useful in assessing patients' risk and in treatment selection. 3 Current studies have shown that the blood supply, function, and environmental carcinogens of colon and rectum are different. Therefore, clinical characteristics, pathological stage, and the prognosis of colorectal cancer are always different among primary sites. [4][5][6] The analysis of the differences in the molecular biology of colorectal cancer among different primary sites is necessary in the selection of individualized treatment for patients with colorectal cancer. At present, with the deepening understanding of tumor biology and the complexity of tumor metabolism, metabolic reprogramming has been found to be a marker of malignancy. 7 Compared with normal tissues, cancer cells have a high degree of metabolic heterogeneity. 8 Metabolic phenotypes evolve at the different stages of tumorigenesis, early tumor growth requires nutrient uptake and biosynthesis, and other subtypes of selective metabolic requirements occur during local invasion and distant metastasis. 9 Two studies proposed prognostic models of prognosis-related metabolic genes for colorectal cancer prognosis. 10,11 In the present study, we modeled the prognosis genes of colorectal adenocarcinoma (COAD) and rectum adenocarcinoma (READ) separately for the individualized detection of patients with colorectal cancer.
We constructed the prognostic models of prognosisrelated metabolic genes for COAD and READ, utilizing the expression profiles of TCGA databases. The prognostic models were further optimized through Lasso regression and Cox regression analyses, and the specificity and sensitivity of the models were evaluated through ROC curve analysis. Meanwhile, metabolism-related genes (MRGs) related to OS were screened out. Our data established two models specific to COAD and READ, which are useful in guiding individualized prognosis and treatment.

Data Collection
The transcriptome data and the relevant clinical data of COAD (including 398 cancer and 39 normal tissues) and READ (including 84 cancer and 2 normal tissues) were downloaded from TCGA database (https://portal.gdc.can cer.gov/) to construct risk prediction models, and GSE17538 of GEO database (https://www.ncbi.nlm.nih. gov/geo/) including 238 COAD patients was used as the validation dataset. Meanwhile, the gene mutation data were obtained from GSCALite (http://bioinfo.life.hust.edu.cn/ web/GSCALite/). The clinical and pathological characteristics of patients in TCGA are summarized in Table 1. We obtained MRGs by searching "metabolism" from the Molecular Signatures Database v7.0 (http://software.broadin stitute.org/gsea/msigdb) and selected the C2 sub-collection (c2.cp.kegg.v7.0.symbols.gmt) as the needed gene set, then 944 MRGs were obtained. Then, MRGs expression was extracted from GSE17538 dataset through R language merge package. Finally, a gene expression matrix consisting of 772 gene expression values was obtained.

Differential Expression Analysis of MRGs and Functional Enrichment Analyses
For COAD, the differences of MRGs between cancer tissues and non-tumor samples were analyzed by Wilcoxon signed-rank test. An adjusted P < 0.01 and a | log2 (FC) |> 2 were considered the cutoffs for identifying differential expression MRGs (DEMRGs). KEGG analysis was then performed to discover the primary biological characteristics of these genes.

Construction and Validation of Risk Prediction Models Using MRGs
We first used univariate Cox analysis to initially identify potential MRGs, and least absolute shrinkage and selection operator (Lasso) regression analysis was subsequently applied by using "glmnet" and "survival" packages to eliminate false parameters caused by overfitting. Multivariable regression analysis was used to screen to MRGs with independent prognostic significance, and these MRGs were used to construct prognostic models with the following formula: Risk score = (Coefficient mRNA1 ×mRNA1 expression) + (Coefficient mRNA2 ×mRNA2 expression) + ⋯ + (Coefficient mRNAn × mRNAn expression). Finally, Cox regression was used to establish OS prognostic risk model. In addition, univariate and multivariate regressions were performed to determine whether the model was an independent prognostic factor.

Gene Set Enrichment Analysis (GSEA)
GSEA was performed to explore the characteristics and pathways that were enriched in the high-risk group and low-risk group. Using normalized enrichment score (NES) and normalized P-value calculation predefined tags and KEGG pathway enrichment. Terms with |NES|>1 and P<0.05 were considered significantly enriched.

Patient Samples
A total of 117 COAD and 114 READ were obtained from Tianjin cancer hospital. All tissue sections were confirmed by specialists to make a final diagnosis. Histopathological diagnoses were made using the World Health Organization criteria.

Immunohistochemistry Staining and Evaluation
Paraffin-embedded tissue sections (4µm) were dewaxed and rehydrated with xylene and fractional alcohol solution, the endogenous peroxidase activity was hardened with 3% hydrogen peroxide, and the sections were exposed to antigen by boiling in 10 mM citric acid buffer (pH 6.0). The sections were incubated with primary antibody at 4°C for 18 h, incubated with PV6001at 37°C for 30 min, and stained with DAB (Zhongshan Goldbridge Biotechnology Company) for 1 ~ 2 min. The control group was incubated with PBS instead of primary antibody.
To score staining of tissue sections, we selected 5 consecutive high-power fields on each section to estimate the average percentage of stained cells. For the expression of GSR in the nucleus, the mean value of DovePress expression in all cancer tissues was truncated, expression level ≤50% in tissues was classified as low expression group, and 50% of > was classified as high expression group ( Figure 1A). For the expression of GSR in cytoplasm, ENPP2 and GAMT, both percentage and intensity were considered in the semi-quantitative evaluation. The percentage of positive cells was divided into 0(0% positive cells), 1(1-25% positive cells), 2 (26-50% positive cells), 3(50-75% positive cells), 4 (>75% positive cells), and the strength was divided into 0(negative), 1(weak), 2(medium), and 3(strong). The intensity score (0-3) is multiplied by the percentage score (0-4), 0-4 being low expression and 5-12 being high expression.

Statistical Analysis
All statistical analyses were performed using R software (version 3.6.0), and the related source codes were provided in Supplementary

Flow Chart of This Study
The detailed workflow for constructing the risk models and downstream analysis was shown in Figure 2. We first used the TCGA database to extract MRGs from the COAD and READ data sets and then observed the differential expression of MRGs at the gene expression level. The specific prognostic risk models of COAD and READ were constructed using the MRG set data, and the predictive ability of each risk model was tested through ROC analysis. The risk model of COAD was validated in the GEO database. Finally, Kaplan-Meier analysis and expression level detection of the genes in the risk models were performed, and the mutation information of each gene in COAD and READ in the models was examined.

Differential Expression and Functional Annotation of MRGs in COAD
From the TCGA database, we first downloaded mRNA expression data and the clinical information of 398 COAD  (Table 1). A total of 772 MRGs were detected in this database, and 196 differentially expressed MRGs were obtained. The expression pattern of the differentially expressed MRGs in COAD and non-tumor tissues was shown in the volcanic and heat maps ( Figure 3A and B). Among the 196 differentially expressed genes, 98 genes were down-regulated, and 98 genes were upregulated in the tumor tissues. In order to observe metabolic differences between colon cancer and normal tissues, we used the KEGG pathway to conduct a functional enrichment analysis of the differentially expressed MRGs. We found that the up-regulated genes are mainly involved in the following metabolic pathways: purine metabolism, amino acid biosynthesis, carbon metabolism, and cysteine and methionine metabolism ( Figure 3C). Meanwhile, the down-regulated genes are mainly involved in fatty acid degradation, starch and sucrose metabolism, retinol metabolism, and galactose metabolism ( Figure 3D). We also downloaded the mRNA expression data and corresponding clinical information of 84 READ tissue samples and 2 non-tumor tissues (Table 1).
Owing to the small number of normal rectal tissue samples, we cannot obtain valid differential gene data.

Construction and Verification of the COAD Risk Model
To investigate the relationship between MRGs and COAD prognosis, we constructed a prognostic risk model. Initially, univariate regression analysis was performed on 13 genes that were significantly associated with prognosis, including 2 low-risk MRGs (GSR and SUCLG2P2) and 11 high-risk MRGs (CPT1C, PLCB2,  PLCG2, PLA2G2D, PTGDS, GAMT, ENPP2, SU LG2P2, PIK3CD, PIP4K2B, GPX3, and MAT1A; Figure 4A). Then, Lasso regression and multivariate regression were used to generate a final prognosis model with eight MRGs, including CPT1C, PLCB2, PLA2G2D, GAMT, ENPP2, PIP4K2B, GPX3, and GSR ( Figure 4B and C). The specific information of these genes, including full name, coefficient, and relevant pathway, was shown in Table 3. Then, the patients with COAD in the TCGA database were divided into low-and high-risk groups, and Kaplan-Meier survival analysis was performed. The results showed that patients with high-risk scores had significantly worse OS than patients with low-risk scores ( Figure 4D). Figure 4E shows the association between survival time and risk score, and survival time decreased significantly with increasing risk score. Time-dependent ROC analysis indicated that this prognosis model can accurately predict the OS for patients with COAD ( Figure 4F). The heatmap developed to show the gene expression profiles indicated that CPT1C, PLCB2, PLA2G2D, GAMT, ENPP2, PIP4K2B, and GPX3 were identified as positive risk-correlated genes, whereas GSR was identified as negative risk-correlated gene ( Figure 4G).
To verify the validity of the prediction model, we used the COAD data of GSE17538 from the Gene Expression Omnibus (GEO) as validation data. Consistent with the results obtained from TCGA database, our finding showed that the prediction model could provide reliable prognoses for the patients in GSE17538, and the patients with highrisk scores had considerably poor OS ( Figure 5A). ROC curve analysis showed the sensitivity and specificity of the risk score for the OS when combined with clinical characteristics ( Figure 5B). Furthermore, the association between survival time and risk score indicated that the survival time decreased with the significant increase in risk score ( Figure 5C). Meanwhile, the results of the heat map were consistent with those in the TCGA database ( Figure 5D).

Construction of the READ Risk Model
By the approach used in COAD analysis, we performed univariate regression analysis to construct a prognostic model for READ, and seven MRGs were used ( Figure  6A). Then, a final prognosis model including six MRGs (TDO2, PKLR, GAMT, EARS2, ACO1, and WAS) was constructed using Lasso regression and multivariate Cox regression ( Figure 6B and C). Specific information on the genes is shown in Table 4, including full name, coefficient, and relevant pathways. Meanwhile, the results of Kaplan-Meier survival analysis showed that patients with low-risk scores had significantly better OS than those with high-risk scores, suggesting that the established prognostic model was available ( Figure 6D). Survival time increased when the risk score decreased significantly ( Figure 6E), and ROC curve analysis showed the sensitivity and specificity of the risk score for OS when combined with clinical characteristics ( Figure 6F). Furthermore, the heatmap showed that TDO2, PKLR, and GAMT were identified as positive risk-correlated genes, while EARS2, ACO1, and WAS were identified as negative risk-correlated genes ( Figure  6G). Given the small number of READ cases in the GEO database, we were unable to use these data to validate the risk model.

Prognostic Risk Models of COAD and READ Were Independently Related to OS
To further validate this prediction model, univariate and multivariate regression analyses were performed, and TCGA and GSE17538 were used to investigate the correlations of clinical parameters and risk score with OS in patients with COAD. The TCGA results of univariate Cox regression showed that pathological stage; T, N, and M stages; and risk score were all significantly correlated with OS (all P<0.001; Figure 7A), and multivariate Cox regression showed that T stage and risk score were correlated with OS in patients with COAD (P<0.05; Figure 7B). The GSE17538 results showed that not only univariate Cox regression analysis but also multivariate regression analysis indicated that the prognostic risk model of COAD was independently related to OS (P<0.05; Figure 7C and D). As to READ, the results of univariate regression and multivariate Cox regression analysis indicated that the obtained risk model was independently related to OS (P<0.05; Figure 7E and F). Comprehensive Analysis of All the Genes in the Prognostic Models of COAD and READ According to the previous results, the COAD risk model has eight genes, and the READ risk model has six genes. To further evaluate the prognostic value of the above genes, we used the GEPIA database in performing Kaplan-Meier analysis. We found that high expression of CPT1C, PLCB2, or GPX3 indicated significantly poor prognosis, whereas the high expression of GSR indicated excellent prognosis in patients with COAD in the GEPIA database ( Figure 8A). Although PLA2G2D, GAMT, ENPP2, and PIP4K2B cannot significantly predict prognosis, patients with high expression level of PLA2G2D, GAMT, ENPP2, or PIP4K2B always had poor prognosis ( Figure 8A). In READ, high TDO2 and GAMT expression levels were correlated with poor prognosis, whereas high EARS2 and ACO1 expression levels predicted good prognosis ( Figure 8B). Overall, the results of Kaplan-Meier analysis were generally consistent with the above univariable Cox analysis results.
Meanwhile, we analyzed the mRNA expression levels of the genes in relative tumor and normal tissues with the GEPIA database. As shown in Figure 9A, in the COAD risk model, the gene expression levels of PLA2G2D, ENPP2, and GPX3 in COAD tissues were significantly lower than corresponding normal tissues, but changes in the other genes were not obvious. As to the READ risk model, only GAMT had significantly lower expression in READ tissues than in corresponding normal tissues ( Figure 9B).
In addition, we investigated genetic alterations of genes in the prognosis models with the data of Gene Set Cancer Analysis Lite (http://bioinfo.life.hust.edu.cn/web/ GSCALite/). The results showed that the genes in the COAD model were changed in all the 38 queried samples (100%; Figure 10A), and the genes in the READ model were changed in all the 5 queried samples (100%; Figure  10B). In COAD, genetic alterations included single nucleotide polymorphism (SNP), delete mutation, and insert mutation, whereas the genetic alterations in READ just included SNP ( Figure 10). Meanwhile, on Copy Number Variation (CNV) module, the statistics of heterozygous and homozygous CNV of each cancer type are displayed in Figure 10C.

Involved Metabolic Pathways of the Genes in the Risk Models Analyzed by GSEA
To understand the metabolic characteristics of COAD and READ tissues, we used the GSEA approach to analyze the major differences in metabolic pathways between high-and low-risk patients. In COAD, inositol phosphate metabolism, tyrosine metabolism, arachidonic acid metabolism, glycerophospholipid metabolism, and ether lipid metabolism were activated in high-risk tumors, whereas fatty acid metabolism, butanoate metabolism, amino sugar and nucleotide sugar metabolism, propanoate metabolism, and glycolysis gluconeogenesis were inhibited in low-risk tumors ( Figure 11A). In READ, glycerophospholipid metabolism, arachidonic acid metabolism, alpha linoleic acid metabolism, linoleic acid metabolism, and glutathione metabolism were activated in high-risk tumors, whereas the TCA cycle, propanoate metabolism, one carbon pool by folate, lysine, valine leucine and isoleucine degradation were inhibited in low-risk tumors ( Figure 11B). All the number of genes and full details for the GSEA analysis of COAD and READ were provided in Supplementary Files 1 and 2.

Expression of GSR, ENPP2 and GAMT in COAD and READ Detected Through Immunohistochemistry Staining
To test the validity of the prognostic models, we detected the protein expression levels of GSR, ENPP2 and GAMT  Immunohistochemical studies had demonstrated that GSR was expressed in the cytoplasm and nuclei of cancer cells ( Figure 1A), ENPP2 and GAMT were mainly expressed in the cytoplasm of cancer cells ( Figure 1H and K). Then, all the patients were divided into two groups (low expression group and high expression group) based on these protein expression levels in the tumor tissue. Kaplan-Meier survival analysis showed that the high expression levels of GSR in the cytoplasm was significantly correlated with the poor prognosis of patients with COAD ( Figure 1B) but not with the prognosis of patients with READ ( Figure 1E). No significant correlation was found between the expression levels of GSR in nucleus and prognoses of patients with COAD ( Figure 1C) and READ ( Figure 1F). Meanwhile, the high expression of GSR in both cytoplasm and nucleus could not indicate the prognoses of patients of COAD ( Figure 1D) or READ ( Figure 1G). In addition, high ENPP2 expression levels in cancer tissues was significantly correlated with the poor prognosis for COAD ( Figure 1I), and no significant correlation between ENPP2 expression and patient prognosis was found in READ ( Figure 1J). Furthermore, the significant correlations were found between the expression levels of GAMT and the prognoses of patients with COAD ( Figure 1L) and READ ( Figure 1M).

Discussion
In this study, the prognosis model of COAD with eight MRGs (CPT1C, PLCB2, PLA2G2D, GAMT, ENPP2, PIP4K2B, GPX3, and GSR), and the prognosis mode of READ with six MRGs (TDO2, PKLR, GAMT, EARS2, ACO1, and WAS) were identified successfully, and multivariable Cox regression analysis showed that the risk score calculated by the model could independently predicted the prognosis of patient with COAD or READ.
Furthermore, some protein expression levels of these MRGs in prognosis models were also significantly correlated with the prognosis of patients with COAD or READ.
The Warburg effect is one of the important characteristics of tumor cells, which shows hyperactive glycolysis even when oxygen is sufficient. However, this effect is only a part of metabolic reprogramming, and metabolic abnormalities in tumor cells are far more complex than previously thought. Tumor cells reprogram many metabolic pathways to the meet nutrient, energy, and REDOX requirements of rapid proliferation. 12 Metabolic reprogramming causes changes in the levels and types of specific metabolites inside and outside tumor cells. Such changes can promote tumor growth and metastasis by influencing gene expression, tumor microenvironment, and cell status. 13 Currently, drugs targeting these metabolic reprogramming can inhibit the malignant progression of tumors and promote the death of tumor cells. [14][15][16] Therefore, using metabolic gene changes in constructing a prognosis model for cancer is feasible. More and more studies have focused on the relationship between metabolic reprogramming and the progression of colorectal cancer. Compared with normal colorectal tissues, glucose, lipid, amino acid and nucleotide metabolism have been found significantly perturbed in cancer tissues. 17,18 Currently, metabolomics' tests have identified a variety of dysregulated metabolites associated with metabolic pathways in colorectal cancer samples, and these metabolites may serve as useful biomarkers for indicating the prognosis of patients with colorectal cancer, but lots of researches are still needed before clinical application. 19,20 Recently, two study teams analyzed DEMRGs between colorectal cancer and normal tissues using TCGA database, then used DEMRGs in constructing a prognosis MRG model. One study constructed a model comprising six genes (NAT2, XDH, GPX3, AKR1C4, SPHK1, and ADCY5) for colorectal cancer 10 and another model comprising AOC2, ENPP2, ADA, ACADL, GPD1L, and CPT2. 11 Although they used the same data, they obtained completely different prognostic models. We considered that the difference might be due to differences in methods used in data processing or different P values defined. In these studies, the authors modeled COAD and READ together for prognosis. As the blood supply, function, and environmental carcinogens vary among different sites, the prognosis of colorectal cancer varies among different primary sites. [21][22][23] Thus, in our study, we DovePress constructed a prognostic model of prognosis-related metabolic genes for COAD and another for READs, using the expression profiles of TCGA databases. Differentially expressed genes in the tumor and relative normal tissues have been used in constructing prognostic models in many studies. We think that by doing so, numerous genes that actually have prognostic significance would be overlooked. The reason is that genes that predict prognosis in cancer patients are in fact not necessarily expressed differently in tumors unlike in normal tissues. We did not directly use differential genes in constructing the prognostic models. Instead, we analyzed all the genes involved in metabolism. Finally, we constructed models using Cox regression and Lasso regression in COAD and READ, respectively. The final prognosis model of COAD with eight MRGs (CPT1C, PLCB2, PLA2G2D, GAMT, ENPP2, PIP4K2B, GPX3, and GSR), and the mode READ with six MRGs (TDO2, PKLR, GAMT, EARS2, ACO1, and WAS) were identified, and multivariable Cox regression analysis showed that the risk score calculated by the model independently predicted the prognosis of patient with COAD or READ. Thus, prognostic models for COAD and READ are quite different. Only one gene, GAMT, appeared in the risk models for COAD and READ. This finding further showed that modeling COAD and READ separately for the prediction of the prognoses of patients is necessary and reliable. Among the eight genes in the prognostic model of COAD, three genes (ENPP2, CPT1C, and PLA2G2D) are involved in lipid metabolism, two genes (PIP4K2B and PLCB2) are involved in phosphoinositol metabolism, and two genes (GPX3 and GSR) are involved in glutathione metabolism, indicating that the metabolism of lipid, phosphoinositol, and glutathione is significantly important to COAD. GPX3, as a risk gene, is responsible for the peroxide metabolism of glutathione, and GSR, as a protect gene, is responsible for the reductive metabolism of glutathione. Thus, glutathione peroxidation may be associated with the malignancy of COAD. In our study, we used immunohistochemistry to verify the relationship between GSR and the prognosis of colorectal cancer. We found that high GSR expression level of the cytoplasm was significantly correlated with the poor COAD prognosis but not with READ prognosis. Among the six genes in the model of READ, two genes (TDO2 and GAMT) are involved in amino acid metabolism, suggesting that amino acid metabolism is active, and the prognosis of READ patient is poor. Meanwhile, two genes (EARS2 and WARS) are involved in tRNA aminoacylation of tryptophany and glutamine, which is associated with good READ prognosis. In addition, PKLR is involved in pyruvate metabolism, and ACO1 is involved in the TCA cycle. GAMT, which is mainly involved in creatine synthesis, was present in the two tumor prognostic models, but the creatine system in relation to malignancy requires extensive investigation. 24 In addition to supporting the biosynthesis of creatine, GAMT supports the metabolism of polyamine and methionine synthesis in tumor cells, both of which are in high demand in proliferating sarcoma cells. 25 The relationship between GAMT protein expression level and the prognosis of colorectal cancer was also detected in this study, and we found that high GAMT expression was significantly correlated with the poor COAD prognosis and READ prognosis.
Notably, two genes, GPX3 and ENPP2, were also identified in the other two studies for MRGs in colorectal cancer. 10,11 GPX3, the only extracellular GPX of the family of oxidoreductases, can catalyze the detoxification of lipid hydroperoxides by reduced glutathione. GPX3 has dual properties in tumors, acting as a tumor suppressor in some tumors and as a tumor promoter in others. 26 In colorectal cancer, tumor cells with low GPX3 expression were more sensitive to chemotherapy drugs, such as oxaliplatin and cisplatin. 27 ENPP2 can encode a secreted glycoprotein autotaxin, which converts lysophosphatidylcholine into lysophosphatidic acid, and the autotaxin-LPA axis is essential for cell survival, migration, proliferation, and differentiation in many cancers. [28][29][30][31] Moreover, ATX expression is positively correlated with microvascular density and tumor angiogenesis in colorectal cancer. 32 In our study, patients with COAD and high ENPP2 expression in cancer tissues always had poor survival, but no significant correlation between ENPP2 expression and patient prognosis was found in READ.
Gene mutations were significantly varied between the two cancer types. In COAD, genetic alterations included SNP, delete mutation, and insert mutation, whereas the genetic alterations in READ included SNP only. Moreover, the major differences in metabolic pathways between the low-and high-risk patients with COAD or READ also differed significantly. In COAD, arachidonic acid metabolism, tyrosine metabolism, glycerophospholipid metabolism, inositol phosphate metabolism, and ether lipid metabolism were activated in the high-risk OncoTargets and Therapy 2021:14 submit your manuscript | www.dovepress.com DovePress tumors, whereas in READ, glycerophospholipid metabolism, arachidonic acid metabolism, linoleic acid metabolism, alpha linoleic acid metabolism, and glutathione metabolism were activated in the low-risk tumors. Basing on this result, we found that arachidonic acid metabolism and glycerophospholipid metabolism were all activated in high-risk groups for COAD and READ. However, no changes in the same metabolic pathway were found in low-risk patients with COAD or READ. In addition, we screened DEMRGs between COAD and non-tumor tissues and performed KEGG analysis on the DEMRGs to understand metabolism in COAD. We found that the genes involved in amino acids biosynthesis and carbon metabolism were up-regulated, whereas the genes involved in fatty acid degradation were downregulated. These results suggested that COAD cells are more prone to anabolism than normal tissues and inhibit catabolism.
Our study has some limitations. Although our prognostic models were conducted based on massive data mining, additional functional and mechanism experiments are still needed for exploring the effect of each gene in the prognostic model on the metabolism of colorectal cancer. Moreover, some important clinical factors and lifestyle habits, such as dietary habits, obesity, smoking, and diabetes, were not included in our analysis. However, these factors influenced the prognoses of patients with colorectal cancer. [33][34][35][36]

Conclusions
In conclusion, applying TCGA and GEO data sets, we successfully identified and validated the risk prediction models with metabolism-related genes for colon cancer and rectal cancer, respectively. Meanwhile, through immunohistochemistry, we verified the relationship between the expression levels of some proteins in the models and colorectal cancer prognosis. The results of our study may provide novel insights for further research on the metabolic characteristics and provide potential metabolic targets for treatment of colorectal cancer. However, there are still many survival analyses and functional experiments that need to be validated.