Impact of metabolic syndrome on Overall Survival (OS) in Patients with PDAC
A total of 95 PDAC patients from the First Affiliated Hospital of Wenzhou Medical University were involved in our study. As of March of 2021, 57 patients died during follow-up, none lost follow-up. The details of clinical information were shown in Supplement materials 1.
Metabolic syndrome is conferred as central obesity, dyslipidemia, hyperglycemia, insulin resistance and hypertension. Metabolic disorders were proven to be associated with the increased tumor risk. A nomogram was constructed to predict 6-month, 1- and 2-year overall survival of the PDAC. Total scores were summations of each variable based on the intersection of the vertical line. By using this nomogram, we could convert each clinical index to the corresponding point, and then calculate the total point, which was used to evaluate the 6-month, 1- and 2-year survival rate. We found that metabolic syndromes and together with other clinical factors played an important role in the prognosis of PDAC patients (Figure 1A). Moreover, decision curve analysis showed the high accuracy of the predictive prognostic of MetS score for 6-month, 1- and 2-year OS possibility (Figure 1B-D).
Identification of Metabolism-related Signature in patients with PDAC.
In order to systematically characterized the distinct molecular mechanism of MetS genes in PDAC patients, patients from the TCGA database were set as a training dataset and applied to establish the Metabolism-related Signature. A total of 1466 metabolism-related genes were applied to univariate Cox regression analysis. A list of 20 genes associated most with the prognosis was selected for further study (P<0.0001). Through the Multivariate stepwise Cox regression analysis, a 5-gene prognosis model was successfully constructed.
Performance of the risk score in training dataset and validation dataset
Patients with PDAC were divided into a high-risk group (N=88) and a low-risk group (N=88) according to their median risk score. The performance of the risk score in predicting the prognosis of the PDAC patient was firstly verified in the training dataset. The survival analysis revealed that the PDAC patients with high-risk score had a significantly shorter overall survival time than patients in low-risk-score group (P<0.001, Figure 2A). Receiver operating characteristic (ROC) analysis was used to describe the discrimination accuracy of our model. The area under the ROC curve (AUC) of the model was 0.786 in the TCGA data set, which indicates that the reliability of our PDAC’s prognosis model. (Figure 2B). As Figure 2C showed, as the riskscore increased, death toll increased and the follow-up time decreased. In addition, we found that the overall survival (OS), disease-free survival (DFS) and progression-free survival (DFS) were shortened in the PDAC patients in TCGA database with an increasing risk score (P<0.00001, respectively) (Supplementary Figure 1).
Besides, progression-free survival (PFS) and disease-free survival (DFS) of two groups were compared with each other respectively. The results showed that patients with high riskscore had a significantly shorter PFS, so as DFS, than the patients with low riskscore (P<0.0001, respectively). The area under the ROC curve (AUC) for the model was 0.79 and 0.798, respectively (Figure 3A-B). These results suggest that the risk score act as a strong prognostic indicator for PDAC patients. Risk score distributions, survival status and expression profiles of the 5 genes of patients in ICGC dataset and GEO dataset was shown in Supplementary Figure 3.
Furthermore, the expression of the model's constituent genes between PDAC and normal tissue was investigated via CEPIA2, a web server. The results showed that the expressions of the CDA, GMPS, and PI4KB were significantly elevated in PDAC from the TCGA dataset compared with normal tissues from the GTEx dataset, while CA12 was downregulated in PDAC (Figure 3C). Survival analysis showed that CA12, GMPS, and PI4KB were associated with the prognosis of PDAC (Figure 3D)
Moreover, two independent validation cohorts were applied to evaluate the robustness of our model in PDAC patients. In the ICGC databases, based on the formula, PDAC patients were subdivided into high-risk and low-risk groups according to their median risk score. Consistent with training cohort results, we found that patients in the high-risk group had a significantly shorter OS than patients in the low-risk group (P=0.034, Figure 4A). Then, this 5-gene signature was further applied to the GEO dataset. In agreement with the TCGA dataset and ICGC dataset, the survival analysis showed that patients with high riskscore had a worse prognosis than those with low riskscore (P=0.005, Figure 4B).
The nomogram prediction model
Multivariate Cox hazard analysis was performed to compare the robustness of risk score with other clinical information in predicting the prognosis of PDAC patients. our analysis contains clinicopathologic indicators including Riskscore, Age, Gender, Histologic Grade, AJCC stage, T stage, and N stage. As shown in Figure 5A, Risk score maintained independence from other clinical indicators in predicting the OS of PDAC patients and was the only factor that was remarkably correlated with the prognosis of PDAC patients in the TCGA dataset (P<0.0001, HR (3.38(2.125-5.375))). According to the Cox regression, a nomogram was built to predict the prognosis of PDAC patients in clinical practice (Figure 5B).
Construction of a Weighted Correlation Network and the identification of target module.
The top 25% of most variant genes were selected to construct a co-expression network with Weighted gene co-expression network analysis (WGCNA) with WGCNA R package (Figure 6A). Clinical data including age, gender, DFS time, DFS status, Histologic grade, OS time, OS status, T stage, N stage, M stage, and risk score were included in our research (Figure 6B). The power of 8 (scale-free R2=0.88) was selected as the soft-thresholding. A total of 15 modules were identified among the patients with PDAC and assigned different colors. The association between each module and the clinical character was then analyzed. The results showed that the blue module (Correlation coefficient=-0.38, P=4e-04), lightcyan module (Correlation coefficient=0.53, P=3e-07), midnight blue module (Correlation coefficient=-0.49, P=2e-06) and purple module (Correlation coefficient=0.42, P=9e-05) were strongly correlated with the risk score of PDAC patients. (Figure 6C). Figure 6D showed a scatter plot of genes in blue modules.
Function enrichment analysis
The target module genes were further analyzed by GO and KEGG to obtain insights into the function of genes in the hub module. GO analysis showed that the genes in the blue module were significantly enriched in T cell activation, leukocyte cell-cell adhesion in biological process (BP); external side of plasma membrane, immunological synapse, alpha-beta T cell receptor complex in Cellular Component (CC); carbohydrate-binding, immune receptor activity, receptor-ligand activity in Molecular Function (MF) (Figure 6E). Moreover, hub genes were also highly represented in these top 5 KEGG pathways, including cytokine-cytokine receptor interaction, viral protein interaction with cytokine and cytokine receptor, hematopoietic cell lineage, chemokine signaling pathway, and T cell receptor signaling pathway (Figure 6F).
CIBERSORT analysis of tumor-infiltrating immune cells (TIICs) in pancreatic cancer
The function that chemokines and chemokine receptors can regulate the immune cells’ infiltration in tumors caught most of our attention. Furthermore, GO and KEGG analysis showed that genes in blue module enriched in the immune process, especially in T cells’ activities. Thus, we would like to explore whether the infiltration of immune cells was different between PDAC patients from high-risk and low-risk groups. CIBERSORT was applied to analyze the profile of tumor-infiltrating immune cells in PDAC. We divided the PDAC patients from TCGA into high-risk and low-risk group according to the median riskscore in our previous study. The differences in immune microenvironment between the high- and low-risk groups were investigated further. The TIICs composition of PDAC patients from two groups was analyzed by CIBERSORT. As shown in Figure 7A and 7C, the intergroup proportions of the 22 types of TIICs were similar while the intragroup proportions were varied. A visualization of the relative proportions of 22 TIICs between the high-risk and low-risk groups was shown by violin diagram (Figure 7C). The results showed that 2 TIICs (T cells CD8, T cells Regulatory) were in higher proportions in the low-risk group than those in the high-risk group, whereas 3 TIICs (Tregs NK cells Activated, Dendritic cells Activated, Mast cells Resting) were in higher proportions in the high-risk group (P<0.05, respectively). Moreover, the two most common TIICs in PDCA tissues were B and T lymphocytes, accounting for approximately 50% of all TIICs. Figure 7B showed the correlation between each TIIC.
Hub gene identification
Our study identified 693 genes in the blue module that correlated highly with the risk score of PDAC patients. Furthermore, we uploaded all genes in the blue module to the STRING database to construct a network of protein-protein interactions (PPI). In the PPI network, genes with a top 20 connectivity degrees were also defined as hub genes (Figure 7D). The expressions of these 20 genes between pancreatic tumor and normal tissue was investigated via CEPIA2, a web server. The results showed that expressions of CCR7, CXCR3, CXCR5, CCL5, CXCR4, CCL19, CCL21, CXCL13, CXCR6, SAA1, S1PR4, PONC were significantly elevated in pancreatic tumor tissue from the TCGA dataset compared with normal tissues from the GTEx dataset. (Figure 7E). In addition, the hub genes which expressed differently between pancreatic tumor tissue.
Gene set enrichment analysis of the immune status between high-risk and low‑risk group
GSEA analysis was applied to explore the potential Immune-related signaling pathways between high-risk and low-risk groups (Figure 8A-C). GSEA revealed that low riskscore was significantly associated with peptide injection OT2 thymocyte up (NES=1.57, P=0.002), 24H TLR1 TLR2 Ligand treated monocyte (NES=1.58, P=0.014) and PBMC CD4 T cell up (NES=1.63, P=0.017). The results elucidated that immune-related responses and processes played an essential role in the metabolism of PDAC patients.
Knockdown of GMPS can significantly represses the proliferation and migration ability of pancreatic cancer cells
To evaluate the potential function of Guanosine Monophosphate Synthetase (GMPS) in PDAC, we knocked down GMPS in PANC-1 and CFPAC-1 cells through the CRISPR-Cas9 system. The efficiency of GMPS knockdown was verified by western blot (Figure 9A). Moreover, the viability of cells was measured by Cell Counting Kit-8 (CCK-8) assay, and the result showed that GMPS knockdown significantly compromised the growth rate of PDAC cells (Figure 9B). Similarly, the colony formation assay revealed that the GMPS knockdown markedly reduced the number of clones after 14 days of cultivation, compared with the control (P<0.001) (Figure 9C). Together, these results suggest that GMPS plays a positive role in the proliferation of PDAC cells in vitro.
Furthermore, Transwell assays were used to assess the migration and invasion of PC cells. The results showed that the motilities of PANC-1 and CFPAC-1 were significantly affected by the expression level of GMPS. The migration rate of cells with GMPS depletion was much lower than the control cells (Figure 9D). These results demonstrated that GMPS depletion inhibited the migration of PDAC cells.
Overexpression of GMPS significantly enhances the proliferation and migration of pancreatic cancer cells
To further assess the oncogenic role of GMPS, GMPS was overexpressed in the pancreatic cancer cell lines PANC-1 and CFPAC-1. Western blot analysis confirmed a significant increase in GMPS expression in PANC-1- GMPS and CFPAC-1-GMPS cells relative to the expression of GMPS in control cells (Figure 9E). CCK-8 assay showed that the overexpression of GMPS increased viability of PANC-1 and CFPAC-1 cells (Figure 9F). Moreover, the colony-forming assay showed that PANC-1 and CFPAC-1 cells with GMPS overexpression demonstrated a considerable growth advantage relative to the respective control cells (Figure 9G). Meanwhile, a distinctly higher migration rate was observed in GMPS overexpression cell lines than in the control cells through Transwell assay (Figure 9H). Thus, GMPS overexpression significantly enhances the proliferation and migration ability of pancreatic cancer cells in vitro.