Clinical and genomic characteristics of metabolic syndrome in colorectal cancer

Association of metabolic syndrome with clinicopathological features

At first, a total of 716 CRC patients from Wenzhou Medical University met the inclusion criteria from 2014 January 1st to 2016 January 1st. As of March of 2020, 211 patients died during follow-up, none lost follow-up. About 49.4% of the patients only received surgery, while the others received both adjuvant postoperative therapy and surgery. Baseline clinicopathological parameters were summarized in Supplementary Table 1.

Metabolic syndrome is conferred as central obesity, dyslipidemia, hyperglycemia, insulin resistance, and hypertension. Metabolic disorders were proven to involve in increased tumor risk. In our study, the proportions of patients with overweight, hypertension, diabetes, and dyslipidemia were 27.5%, 12.2% and 59.5%, respectively. There were 108 patients diagnosed with metabolic syndrome (meeting over three requirements) and 608 patients failed to meet the diagnostic criteria. Correlation values between metabolic syndrome indexes and clinicopathological features were shown in Table 1.

Impact of metabolic syndrome on OS in patients with colorectal cancer

The following metabolic variables had significant differences in survival rates within their respective groups by using univariate analysis. There were statistically significant differences between the groups with respect to TG, HDL, HP, DG and BMI (each parameter with P< 0.05, Figure 1A–1E, Table 2). The univariate analysis demonstrated that compared with normal blood glucose level patients, CRC patients with diabetes showed a poor prognosis in OS (P = 0.00011). Low BMI can be considered as a predictor of poor survival in OS. In particular, the underweight group has the worst prognosis, while the normal and overweight groups have similar prognosis (P<0.0001). In the same way, an increasing in HDL could lead to the increase in the mortality of CRC (P = 0.0038). Overall, patients diagnosed with metabolic syndrome (score ≥ 3) had the worst survival outcome among the enrolled population (Figure 1F). The multivariate analysis showed that only BMI (P < 0.001), diabetes (P = 0.046), TG (P = 0.001), and metabolic syndrome (P < 0.001) were significantly independently associated with OS.

Figure 1. Kaplan-Meier curves for CRC patients stratified by metabolic syndrome. Kaplan-Meier analysis of overall survival (OS) of TG (A), HDL (B), hypertension (C), diabetes mellitus (D), BMI (E) and metabolic syndrome (F); (G). Nomogram developed by integrating metabolic syndrome and clinical pathological parameters for predicting 1-, 3-, 5-year survival of CRC patients; (H). Calibration curve for risk of 1-, 3-, 5-year survival of CRC metabolic syndrome.

An OS nomogram was constructed to predict 1-, 3- and 5-year overall survival of colorectal cancer (Figure 1G). Total scores were summations of each variable based on the intersection of the vertical line. As shown in Figure 1G, metabolic syndrome and age contributed the most risk points (ranged, 0-100), whereas the other clinical information contributed much less (ranged, 0-75). By using this nomogram, we could convert each clinical index to the corresponding point, and then calculated the total point, which was used to evaluate the 1-year, 3-year, and 5-year survival rate. Moreover, decision curve analysis showed the high accuracy of the predictive prognostic of MetS score for 1-, 3- and 5-year OS possibility (Figure 1H).

Summary and analysis of the genomic mechanism of metabolic syndrome

With the rapid development of oncotherapy, metabolism regulation, as a promising predictor, has received more and more attention. Aforementioned epidemiological data have proven the tight connection between metabolic syndrome and increased cancer risk. Nevertheless, the mechanisms that link metabolic disorders and cancer risk remain unknown. The possible association mechanism between these two risk factors will be firstly described in this paper, focusing on the shedding light on candidate signature genes and biological events occurring in CRC progression. We identified 11 key metabolic genes, including TGFB, IGF2BP1, IGFL1, IGF2BP3, CHS2, IGF1R, IL6, IL6R, CSF1, CSF3 and IGF2 from recognized literature [5, 10]. Subsequently, we constructed a gene metabolic score via Cox regression analyses. As shown in Figure 2A, 2B, CRC patients with high metabolic signature had significantly poorer OS than patients with low metabolic signature. Next, univariate and multivariate Cox regression analyses were constructed to evaluate the prognostic significance of the gene metabolic signature combined with various clinicopathologic factors (Figure 2C, 2D). Gene MetS score was also proven to be an independent prognostic factor and had strong predictive power.

Figure 2. The distribution of gene MetS Score in the TCGA cohort. (A). K-M survival curve of the low- and high- MetS Score for TCGA CRC patients; (B). The distributions of the MetS Score and survival status for each CRC patients; (C). Forest plot summary the univariable analyses of overall survival of gene MetS Score; (D). Forest plot summary the multivariable analyses of overall survival of gene MetS Score.

Identifying metabolic -related modules

To particularly describe molecular events relative to cancer metabolic progress, the WGCNA was used to identify the module eigengenes that were associated with MetS score. The power (β) of 6 was selected as the soft-power to ensure a scale-free network (scale R2 = 0.95), and ultimately identified 14 co-expressed gene modules (Figure 3A). The correlation between the module eigengenes and the clinical traits were shown in Figure 3B. Particularly, black and turquoise eigengenes were significantly associated with cancer MetS Score (Figure 3C).

Figure 3. Weight gene co-expression network analysis (WGCNA) identified metabolic -related modules eigengenes. (A). The clustering dendrograms for the CRC groups. (B). Heatmap of module-trait relationships. (C). The relationship of module membership and gene expression.

Metabolic signature was associated with drug metabolism

To better annotate the black module function, we singled out the 20 central genes in the co-expression network whose MM > 0.8. The known metabolic-related genes were strikingly marked in red (Figure 4A). Functional enrichment showed that these genes owe strong association with classic drug metabolism pathways, such as EGFR signaling, PI3K-mTOR signaling and JAK signaling (Figure 4B). In addition, we compared the gene expression profiles between high- and low- gene MetS score by GSEA. It suggested that the high-MetS score group was closely associated with key drug metabolism process, such as UDP glycosyltransferase, compound transmembrane transporter activity, drug metabolism cytochrome P450 and other enzymes (Supplementary Figure 1). Consistent with black module functional annotation, metabolic syndrome factors (specially for DG, BMI and HDL) also have an immense influence on chemotherapy drug metabolism in our above clinic retrospective investigation (Supplementary Figure 2).

Figure 4. Functional annotation of black module and drug response. (A). PPI network of black module genes. The known metabolic -related genes were strikingly marked in red; (B). Functional annotation for black module; (C). Correlation between drugs metabolism and metabolic –related genes, which affect major signaling pathways—IGF, EGFR, mTOR, ERK-MAPK, p53, JNK and p38-MAPK signaling. (D). The correlation of gene MetS Score and major drugs metabolism signaling pathways.

To further detect the effects of metabolic syndrome on drug metabolism, further in-depth study of the correlation of metabolic-related gene expression and chemotherapy drug sensitivity was performed by Pearson correlation, which clarified the pathological mechanisms that link its components with carcinogenesis. We focused on major targeted drug pathways—IGF, EGFR, mTOR, ERK-MAPK, p53, JNK and p38-MAPK signaling. Remarkably, IL-6 and IL6R remarkably improved the effectiveness of targeted therapy administration, yet IGF2 was closely related with drug tolerance (Figure 4C). Specially, low metabolism group was strongly associated with the mTOR pathway and VEGFR pathway without significantly association with AKT pathway (Figure 4D).

IL6 promoted oncogenic growth in CRC by stimulating mTOR signaling

We compared the differential expression of key metabolic genes and they all showed high diversity (P< 0.05, Supplementary Figure 3). IL6 was thought to promote tumor growth mainly by paracrine and autocrine methods, which explained the above experimental results to some extent. GSEA revealed that high IL6 expression groups in the TCGA CRC cohort were mainly enriched in KEGG pathways related to drug metabolism, especially the mTOR pathway (Figure 5A). We further specially delineated the protein expression patterns of IL6 and p-S6 by immunohistochemical staining, and they showed higher expression patterns in above MetS CRC patients (Figure 5B). We further detailed investigated the consequence of IL6 overexpression by stably expressing IL6 in SW480 and DLD1. Ectopic expression of IL6 was also sufficient to robustly promote CRC cell migration and transwell compared with control cells (Figure 5C, 5D). On the other hand, IL6 overexpression strongly potentiates the activity of mTOR pathway as shown by much enhanced S6K1 and mTOR phosphorylation (Figure 5E). Given the previous research, we speculated that whether IL6 participated in the regulation of the activity of PP242 activity, a mTOR kinase inhibition. As shown in Figure 5F, IL6 ectopic expression enhanced PP242 drug susceptibility in colorectal cancer cells, which completely abolished the activity of mTOR pathway. We also found that PP242 strongly inhibited the invasion of IL6-overexpressing CRC cells (Figure 5C, 5D). The above results suggested that overexpression of IL6 prompted colorectal cancer highly dependent on mTOR-S6K signaling and more sensitive to mTOR kinase inhibitors.

Figure 5. IL6 promotes oncogenic growth in CRC by stimulating mTOR signaling. (A). GSEA of IL6 in the TCGA CRC cohort; (B). Immunohistochemical staining of IL6 in MetS and non-MetS CRC patients; IL6 promotes growth-factor-induced migration (C) and invasion of CRC cells (D). Cancer cells with or without IL6 overexpression were treated with or without PP242; (E) IL6 promoted mTOR signaling in colorectal cancer cells. SW480 and DLD1 expressing ectopic IL6 or vector were analyzed for mTOR signaling by immunoblotting for P-mTOR and P-S6K. GAPDH was used as a loading control; (F) Pharmacological inhibition of mTOR signaling by PP242 abrogated IL6 overexpression-induced activation of mTOR signaling. SW480 and DLD1 cells were transfected with IL6 shRNAs or a control shRNA and then treated with PP242. Their effect on mTOR signaling was analyzed by immunoblotting.

Patients and clinical outcome assessment

This study utilized data from the First and the Second Affiliated Hospital of Wenzhou Medical University and public databases. We performed a retrospective study of colorectal cancer patients in these two institutions in Wenzhou, China (Wenzhou Medical University), from January 2010 to January 2016. The inclusion criteria were as follows: 1. histopathological diagnosis of colorectal cancer; 2. all data for patients diagnosed with colorectal cancer for the first time without any treatment; 3. complete pathology, laboratory, and follow-up data. Patients with unknown included variables were excluded. The following demographic, clinical, and pathology data were used: T stage, N stage, M stage, pathological stage, tumor history, laboratory test results (age, sex, body-mass index (BMI), TG, HDL-C, CEA, creatinine). Pathologists assessed the tumor stage according to the 7th edition of the AJCC TNM staging guidelines. Patients who received only neoadjuvant chemotherapy or surgery or those with unknown included variables were excluded. Metabolic syndrome was internationally defined as included more than three criteria: 1) BMI greater than or equal to 25.0 kg/m²; 2) diagnosed with diabetes; 3) diagnosed with hypertension, SBP/DBP >140/90 mmHg; 4) blood HDL-C < 0.9 mmol/L, 5) blood TG > 1.7 mmol/L. Totally, there were 716 eligible cases selected in our study, 15.1% were diagnosed with MetS. All of these patients were followed up, and recurrent and dead patients were recorded during the follow-up. The time was cut off in March 2020. The study protocols were approved by the Wenzhou Medical University Ethics Committee. All procedures adhere to the BRISQ Guidelines for reporting research on human biospecimens.

We also retrospectively selected colorectal cancer gene expression and relative clinicopathological data from the TCGA database (https://portal.gdc.cancer.gov/). Raw microarray data Affymetrix were downloaded and normalized using the limma package. The detailed working algorithm was demonstrated in the Figure 6.

Figure 6. Flow chart of the experimental design and main process.

Construction nomogram models

An OS nomogram was constructed based on the prognostic factors derived from multivariate Cox regression analysis to predict 1-, 3- and 5-year survival rate. Each patient could sum up corresponding variable score and finally establish predictive measures of OS. The nomogram was generated using ggplot packages together with R software. The calibration curve for predicting 1-, 3- and 5-year OS indicated that the nomogram-predicted survival was closely corresponded with actual survival outcomes. The survival analysis was conducted using rms, survivalROC, survcomp and survival package. Hazard ratios (HRs) and 95% confidence intervals (CIs) were recorded.

Evaluation of TCGA gene metabolic score

R language survival package was used to perform the Cox regression model analysis. The key metabolic gene lists were obtained from literature investigation the molecular mechanism of MetS from Pubmed [10, 35], and then applied to multivariate Cox regression analysis. The MetS score of each patient was calculated with the expression level and its relative coefficient. On the basis of the median score as the cut-off setline, these patients were divided into high- and low-MetS subgroups. Log-rank test was performed to calculate the corresponding hazard ratios (HRs) and 95% confidence interval (CI).

Weighted gene co-expression network analysis

The weighted gene co-expression network analysis (WGCNA) hierarchically clustered module eigengenes of the clinical features based on co-expression relationships, thoroughly explored the biological processes and molecular mechanisms behind cancer metabolic disorders [36]. In this study, we constructed WGCNA to explore the colorectal cancer dataset to identify gene modules associated with expression patterns of MetS score and pathology factors. The co-expression network was constructed by the R package WGCNA. The connectivity degree of each node of the network was calculated by STRING database and reconstructed via Cytoscape software.

Functional annotation

The signaling pathway and molecular function underlying MetS score were explored with GSEA (Version 4.0.1). The number of permutations was set at 1000, and P < 0.05 and an FDR < 0.25 were considered statistically significant. Gene ontology (GO) enrichment analysis was performed with the DAVID platform.

GDSC (Genomics of drug sensitivity in cancer)

GDSC, a large-scale drug screening data screened on a panel of 990 human cancer cell lines, concentrates on providing publicly available tumor treatment genome data and identifying potential tumor treatment targets. (http://www.cancerrxgene.org/gdsc1000/) [37–39]. Gene mutations of the cancer genomes greatly affect clinical outcomes and drug targets respond. To analyze the correlation of expression and drug sensitivity, the Pearson correlation coefficient of transcript levels was calculated.

Immunohistochemistry

Above-mentioned 20 CRC specimens were collected, including 10 MetS and 10 non-MetS CRC tissues. Two researchers evaluated the staining results independently and scored the staining intensity of immunostaining as: 0 (negative), 1 (weakly positive), 2 (moderately positive) and 3 (strongly positive). The primary antibody against P-S6 and IL6 was used at a concentration of 1:200. For quantitative analysis, a Histo (H) score was calculated based on the staining intensity and percentage of stained cells using Aperio ScanScope system (Vista, CA, USA).

Colony formation and transwell migration assay

A number of 1 × 10³ DLD1 and SW480 were inoculated in six-well plates and incubated at 37° C for 5-7 days. Then, we treated the colonies with PP242 or vector for 5 days. Colonies were fixed with 4% paraformaldehyde formaldehyde (Solarbio, Beijing, China) followed by staining with crystal violet (Sigma-Aldrich). The number of colonies with more than 50 cells was calculated. Transwell migration experiments were used to confirm the migration ability of DLD1 and SW480, and 5 × 10⁴ cells were added to the upper chamber placed in a 24-well plate, with serum-free medium. Meanwhile, medium containing 15% serum was added to the lower chamber. Taking cell images at 100× magnification 5 random microscopic fields of view were selected to count the number of migrated cells.

Antibodies and western blot analysis

An equal amount of protein was subjected to SDS-PAGE. Proteins were transferred into PVDF membranes, and the blots were incubated with the following different primary antibodies: Rabbit mTOR, p-mTOR (S2448), S6K, p-S6K (T389), p-S6 from Cell Signaling Technology, and IL6 from Proteintech. HRP-labeled GAPDH, anti-mouse and anti-rabbit antibodies were purchased from Santa-Cruz Biotechnology. All primary antibodies were confirmed to be reactive only to the targets by the manufacturer and used at 1:1000, and secondary antibodies were used at 1:5000.

Statistical analysis

R software and Stats were used for statistical analyses. Continuous variables were exhibited for means, medians, range, and standard deviation (SD) and compared using an independent t test or Wilcoxon test; Spearmen correlation coefficient was used for variable correlation; Chi-square test was used to analyze categorical variables. The log-rank survival test and the results were shown in the forest plot. All statistical tests were two-sided and P <0.05 was considered statistically significant.

Date accessibility

The data are available in the TCGA datasets.