Skip to main content

ORIGINAL RESEARCH article

Front. Immunol., 09 June 2023
Sec. Cancer Immunity and Immunotherapy
This article is part of the Research Topic Application of Advances in Molecular and Genetic Researches to Tackle Pancreatic Cancer: Screening, Diagnosis, Prognostication and Drug Response Prediction. View all 5 articles

Establishment and validation of a ubiquitination-related gene signature associated with prognosis in pancreatic duct adenocarcinoma

Yangyang Guo,&#x;Yangyang Guo1,4†Zhixuan Wu&#x;Zhixuan Wu2†Kenan Cen&#x;Kenan Cen1†Yongheng BaiYongheng Bai3Ying DaiYing Dai1Yifeng Mai*Yifeng Mai1*Kai Hong*Kai Hong1*Liangchen Qu*Liangchen Qu4*
  • 1Department of General Surgery, The First Affiliated Hospital of Ningbo University, Ningbo, China
  • 2Key Laboratory of Diagnosis and Treatment of Severe Hepato-Pancreatic Diseases of Zhejiang Province, The First Affiliated Hospital of Wenzhou Medical University, Wenzhou, China
  • 3National Key Clinical Specialty (General Surgery), The First Affiliated Hospital of Wenzhou Medical University, Wenzhou, China
  • 4Department of Emergency, Taizhou Hospital of Zhejiang Province Affiliated to Wenzhou Medical University, Taizhou, Zhejiang, China

Background: Patients with pancreatic duct adenocarcinoma (PDAC) have varied prognoses that depend on numerous variables. However, additional research is required to uncover the latent impact of ubiquitination-related genes (URGs) on determining PDAC patients’ prognoses.

Methods: The URGs clusters were discovered via consensus clustering, and the prognostic differentially expressed genes (DEGs) across clusters were utilized to develop a signature using a least absolute shrinkage and selection operator (LASSO) regression analysis of data from TCGA-PAAD. Verification analyses were conducted across TCGA-PAAD, GSE57495 and ICGC-PACA-AU to show the robustness of the signature. RT-qPCR was used to verify the expression of risk genes. Lastly, we formulated a nomogram to improve the clinical efficacy of our predictive tool.

Results: The URGs signature, comprised of three genes, was developed and was shown to be highly correlated with the prognoses of PAAD patients. The nomogram was established by combining the URGs signature with clinicopathological characteristics. We discovered that the URGs signature was remarkably superior than other individual predictors (age, grade, T stage, et al). Also, the immune microenvironment analysis indicated that ESTIMATEscore, ImmuneScores, and StromalScores were elevated in the low-risk group. The immune cells that infiltrated the tissues were different between the two groups, as did the expression of immune-related genes.

Conclusion: The URGs signature could act as the biomarker of prognosis and selecting appropriate therapeutic drugs for PDAC patients.

Introduction

Pancreatic duct adenocarcinoma (PDAC)is often called the “king of cancer” (1, 2). Its prognosis is extremely poor making it the fourth leading contributor to cancer-associated death globally (3). Due to the lack of specific clinical manifestations in the first stages of PDAC, only a small proportion of individuals are identified with certainty at such early stages, and the vast majority are diagnosed at a more advanced level (4). Despite recent breakthroughs in the systematic treatment of PDAC, the prognoses for those with the advanced disease remain dismal owing to the disease’s rapid local progression and frequent distant metastasis (5). Therefore, it is very important to find some pivotal genes that may regulate the onset and advancement of PDAC and serve as novel therapeutic targets for PDAC.

Recently, with the considerable progress of cutting-edge high-throughput sequencing technologies and the growing improvement of public databases, two authoritative databases, the Cancer Genome Atlas (TCGA) and the International Cancer Genome Consortium (ICGC), have collected considerable clinical, pathological, and biological data of cancer patients (6, 7). Researchers can use these data and a variety of bioinformatics analysis methods to screen and predict new diagnostic and prognostic markers for various cancers. Like most cancers, PDAC is a complex malignant disease involving multiple molecules. At present, researchers have successfully established a variety of effective polygene prognostic risk models using bioinformatics technology (8, 9). A multigene prognosis model is helpful to evaluate the total survival period and recurrence risk of patients, and identify high-risk patients with poor prognosis and timely and systematic treatment, while for low-risk patients, unnecessary treatment burden can be appropriately avoided. For example, based on DNA methylation, autophagy, and immune-related genes, the prognosis prediction models for PDAC have strong prediction ability, which can be used for early diagnosis, prognosis evaluation, and treatment (1012). Ubiquitination-related genes (URGs) have been reported as regulators of tumors, affecting tumor cell cycle regulation, gene expression, and progression (13). However, there is much less understanding of ubiquitination in the PDAC microenvironment and prognosis.

In this study, to examine the link between URGs and the prognosis of PDAC individuals, we used multivariate Cox and LASSO regression analyses to identify three ubiquitination genes that have the most impact on the prognosis of PDAC individuals and constructed a three genes prognosis model. In addition, we integrated clinicopathological parameters and risk scores to develop a novel nomogram for clinical application, which can more directly assess PDAC patients’ prognoses, and help achieve personalized therapy. We then found the different immune statuses between various URGs risk groups. Additionally, we examined the prognostic model’s biological roles and signaling pathways to further evaluate the probable molecular processes that influence PDAC patients’ survival and prognoses.

Materials and methods

Data retrieval

The TCGA database (TCGA-PAAD, https://portal.gdc.cancer.gov, 2022.11.21) was searched to obtain the gene expression data (FPKM) of 178 PDAC tumor samples, 4 normal tissue samples, and the related clinical data. Genotype-Tissue Expression (GTEx) database (https://commonfund.nih.gov/GTEx) was used to obtain the gene expression data of pancreatic normal tissue. TCGA-PAAD was categorized into the train and test groups according to the 1:1 ratio with R software (Supplementary Table S1). The ICGC database (ICGC-PACA-AU, https://dcc.icgc.org/) was searched, and 88 pancreatic cancer samples along with their associated prognostic data were retrieved for external verification after normalization. The GEO database (GSE57495, https://www.ncbi.nlm.nih.gov/geo, 2022.11.21) was also obtained for external verification. A search of the MSigDB database (http://www.broad.mit.edu/gsea/msigdb/, 2022.11.21) yielded 79 URGs, and the genes are listed in Supplementary Table S2. The research flowchart was shown in Figure 1.

FIGURE 1
www.frontiersin.org

Figure 1 The flowchart of the current study. *P < 0.05.

Consensus clustering analysis of URGs

The R package “limma” and “ConsensusClusterPlus” were used for consistent cluster classification of PDAC (14). The association between clusters and overall survival (OS) was analyzed by R packet “survival”. The results were analyzed by R packages “pheatmap”, “survival” and “survminer” as heat maps and Kaplan-Meier (KM) curves. The “limma” program was employed to determine DEGs between two clusters with the criteria of |log fold change (FC)| >1 and FDR < 0.05. Scores of infiltrating immune cells were derived via the CIBERSORT method, and the difference in infiltration of immune cells between the two subtypes was assessed by “limma” package.

Development and verification of the URGs prognostic signature

A univariate Cox regression analysis was conducted to determine DEGs associated with prognosis. Then, using multivariate Cox and least absolute shrinkage and selection operator (LASSO) regression analyses by “caret”, “glmnet” and “survival” packages, thress genes were selected and integrated into the prognostic signature (15). The median risk score was used to classify individuals with PDAC into two categories (low- and high-risk categories). Subsequently, the OS was compared by KM analysis. OS and Receiver Operating Characteristics (ROC) of subgroups were analyzed with the “survival”, “survminer” and “timeROC” R packages for 1, 3, and 5 years. Specifically, the “ggplot2” R program was employed to conduct a principal component analysis (PCA). By incorporating risk assessment with clinical data, a nomogram was developed. Next, multifactor ROC was implemented to verify the predictive accuracy of the nomogram.

Comparative analysis of the tumor microenvironment between high- and low-risk categories

Immune cell abundance (ImmuneScores) and stromal cell abundance (StromalScores) were evaluated by the ESTIMATE (16). To examine the variation in immune cell infiltration between high-risk and low-risk categories, we used the TIMER, CIBERSORT-ABS, QUANTISEQ, EPIC, MCPCOUNTER, and CIBERSORT, XCELL, algorithms. The correlation of risk score and immune cell was evaluated by “corrplot” packages. Differential immune cell infiltration and immune function were probed via single-sample gene set enrichment analysis (ssGSEA) using “GSEABase” package. The expression patterns of immune-related genes were also determined. Tumor immune dysfunction and exclusion (TIDE) acted as a vital biomarker for immunotherapy response. Additionally, we assessed whether there was any link between TIDE scores and risk scores.

Pathway analysis of the URGs signature

We examined the DEGs in the high- and low-risk categories. The underlying pathway analysis associated with DEGs was enriched through the Gene Ontology (GO), Disease Ontology (DO), and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis using “DOSE”, “org.Hs.eg.db” R packages. To assess the probable biological functioning differences between high- and low-risk categories, a gene set variation analysis (GSVA) was carried out.

Drug sensitivity analysis

We investigated the potential for URGs signature to serve as a predictor for medications used in chemotherapy and targeted treatment. Subsequently, the half-maximal inhibitory concentration (IC50) was computed with the “pRRophetic” method (1, 8).

Reverse transcription quantitative polymerase chain reaction

Pancreatic tissue samples were collected from the Ningbo First Hospital, including eight normal pancreatic tissue samples and eight pancreatic cancer tissue samples. The study was approved by the Ethics Committee of the Ningbo First Hospital. All research was performed in accordance with relevant guidelines/regulations. Trizol was employed to isolate total RNA, after which it was reverse-transcribed into the cDNA template. Next, RT-qPCR was conducted with the aid of SYBR Green Real-Time PCR Master Mix Plus (Toyobo). Analyses were conducted according to MIQE guidelines. The internal reference gene utilized was β-Actin. Supplementary Table S3 outlines the amplification primer sequences.

Cell culture and S100A2-knockdown by siRNA

PANC-1 cells were cultured in DMEM with 10% fetal bovine serum under standard culture conditions. Using Lipofectamine 2000 (Invitrogen), siRNA (100 nm) was transfected into cells 48 hours after transfection according to manufacturer’s instructions. Colony formation assay was evaluated by crystal violet staining methods. PANC-1 cells were seeded at 1000 cells per well to six‐well plates then cultivated for 14 days. Sequentially, cell number at each well was counted after staining.

Transwell assay

Transfected PANC-1 cells were seeded at 2× 105 cells per upper transwell chamber, cultivated with or without 100 mL of reconstituted Matrigel-coated membrane for 36-48 hours. Then stained the cells with crystal violet. Thereafter, number of migration or invasion cells was count.

Statistical analysis

Data were presented as the mean ± SEM. Significant differences were evaluated by performing Student’s t-test using Prism software v6.02. Moreover, the Kruskal-Wallis test was used for variables with more than two groups. The Kruskal-Wallis test and Wilcoxon rank sum test were applied to analyze correlations. Correlation analysis between two groups of variables was used spearman correlation coefficient. Statistical significance was set at P < 0.05.

Results

Identification of URGs clusters in PDAC

The link between URGs expression and PDAC subtypes was first analyzed using a consensus clustering method. As depicted in Figures 2A–C, the CDF curve was applied to categorize patients with PDAC into two clusters (C1 and C2). In contrast with C2, C1 individuals diagnosed with PDAC had remarkably lower survival duration (Figure 2D). The correlation between URG clusters, clinical characteristics, and URGs expression in PDAC patients is depicted in Figure 2E. Most URGs were expressed higher in C1 than in C2, and the Grade of the samples in C1 was higher, while T, N, stage, age and gender had no significant difference between C1 and C2.

FIGURE 2
www.frontiersin.org

Figure 2 URG clusters and clinical characteristics between PDAC samples in two clusters. (A, B) The cumulative distribution function curve illustrates the most effective way of URG clustering. (C) The consensus matrix of the clustering analysis via k-means clustering (k = 2). (D) Kaplan–Meier (KM) curves for the overall survival (OS) of PDAC patients among different URG groups. (E) Heatmap of URG expression in PDAC patients with different clinical characteristics and URG clusters.

Since immune cells perform an instrumental function in the onset and advancement of PDAC, we next evaluated the variations in infiltrating immune cells between the two clusters. In cluster 2, the level of monocytes, resting mast cells, naive B cells, and CD8 T cells were higher than in cluster 1, while Tregs, Eosinophils, Macrophages M0 and Mast cells activated were lower than in cluster 1 (Supplementary Figure S1).

Development and validation of the ubiquitination-related prognostic signature

Using the “limma” program, 996 DEGs were found between two clusters with the criteria of |log fold change (FC)| >1 and FDR < 0.05. Thereafter, 43 ubiquitination-related DEGs whose expression levels were remarkably different between PDAC and normal tissues were identified (Figures 3A, B). Next, 12 prognosis-related DEGs were found by the univariate Cox analysis. Subsequently, we completed a LASSO analysis to remove the overfitting genes and the URGs signature of 3 genes (SLC22A17, UCHL1 and S100A2) was created (Figures 3C, D). The equation applied to derive the risk score is as indicated: risk score= (SLC22A17 × (-0.260926538020362) + (UCHL1 × (-0.286371148071792) + (S100A2 × (0.157355046660652).

FIGURE 3
www.frontiersin.org

Figure 3 Construction of the prognostic signature. (A, B) The different expression of DEGs between PDAC and normal tissue. (C) LASSO coefficient profiles (y-axis) of the gene sets and the optimal penalization coefficient (l) via 3-fold cross-validation based on partial likelihood deviance. (D) The dotted vertical lines represent the optimal values of l. The top x-axis has the numbers of gene sets, whereas the lower x-axis revealed the log (λ). (E) Risk score and survival outcome of each case. (F) Heatmap showed the expression of 2 genes in two risk groups. (G) The KM curve showed that patients in the high-risk group had a worse prognosis. (H) The AUC for 1-, 2- and 3-years survival.

Patients with PDAC were classified into low- and high-risk categories as per the median risk score value (Figure 3E). The variations in the expression of these two genes between the two risk categories are illustrated in Figure 3F. Also, patients having elevated risk scores had a greater fatality rate (Figure 3G). Moreover, the ROC curve was performed to assess the URGs signature, which manifested that the AUC values for 1-, 2- and 3-year periods were 0.710, 0.706, and 0.748, respectively (Figure 3H).

Furthermore, we verified the aforementioned findings in test datasets. All patients with PADC in the test datasets were also classified into low- and high-risk categories. The KM analysis disclosed that the low-risk individuals exhibited a more favorable prognosis in contrast to those at high risk in TCGA-test, TCGA-all, and ICGC-PACA-AU (Figures 4A–C). The AUC values of the ROC curve of 1-, 2-, and 3-year periods were 0.738, 0.654, and 0.723, correspondingly, in TCGA-test (Figure 4D), 0.723, 0.707, and 0.704 in TCGA-all (Figure 4E), and 0.726, 0.738, and 0.759 in ICGC-PACA-AU (Figure 4F). In addition, another test dataset GSE57495 was also used for validation. Patients having high risk scores had a greater fatality rate, the AUC values for 1-, 2- and 3-year periods were 0.706, 0.733, and 0.871, respectively (Supplementary Figure S2).

FIGURE 4
www.frontiersin.org

Figure 4 Validation of the prognostic signature. KM curve showed that patients in the high-risk group had a worse prognosis in TCGA-test (A), TCGA-all (B), and ICGC-PACA-AU (C). The AUC for 1-, 2- and 3-years survival in TCGA-test (D), TCGA-all (E), and ICGC-PACA-AU (F).

Construction of a nomogram for PDAC

Multivariate and univariate Cox regression analyses proved that risk score independently acted as a robust prognostic marker (P < 0.05) (Figures 5A, B). An innovative nomogram was developed using the URG signature and clinical variables from the training dataset to further exploit the URG signature’s prognostic potential (Figure 5C). The predictive 1, 2, 3-survival rate was close to the actual observation (Figure 5D). In addition, the ROC analysis was conducted to evaluate the nomogram’s prognosis-predicting value. For 1-year survival times, the AUC value was 0.745 (nomogram), 0.731(risk score) (Figure 5E). For 2-year survival times, the AUC value was 0.776 (nomogram), 0.680(risk score) (Figure 5F). For 3-year survival times, the AUC value was 0.802 (nomogram), 0.708 (risk score) (Figure 5G). These results revealed that this novel nomogram could act as an admirable prognosis prediction model.

FIGURE 5
www.frontiersin.org

Figure 5 Construction and assessment of nomogram. (A) Univariate Cox regression (B) multivariate Cox regression analyses. (C) The prediction of nomogram in the TCGA-train dataset. (D) Calibration plots for the nomogram. The multifactor AUC for 1- (E), 2- (F), and 3-years (G) survival.

The tumor microenvironment analysis in high- and low-risk groups

The TME serves as a crucial indicator of the biological behavior of the tumor. ESTIMATE analysis revealed that the ImuneScores, StromalScores, and ESTIMATEScores were all lower in the high-risk category in contrast with the low-risk category (Figure 6A). ssGSEA analysis found less infiltration of the B cells, CD8+ T cells, DCs, iDCs, Neutrophils, Mast cells, T helper cells, Tumor-infiltrating cell (TIL), and T cells regulatory (Treg) in the high-risk patients in contrast with the low-risk patients (Figure 6B). Some immunologic functions, including T cells co-stimulation, CCR, Type II IFN response, and T cell co-inhibition were also improved in the low-risk patients (Figure 6C). Additionally, the distinctions of immune cell levels between the two risk groups were also investigated through CIBERSORT, MCPCOUNTER, QUANTISEQ, EPIC, TIMER, CIBERSORT-ABS, and XCELL. As per the findings, the low-risk category had remarkably higher levels in most immune cells, including naive CD4 T cells, CD8 T cells, DCs, Cancer associated fibroblast, NK cell, B cell and Monocyte (Figure 6D). This may explain why the low-risk category has a superior prognosis. Additionally, Figure 6E depicted the distribution of low- and high-risk individuals across multiple immune subtypes.

FIGURE 6
www.frontiersin.org

Figure 6 Analysis of immune conditions of high- and low-risk groups. (A) Differences in immune microenvironment scores between the two groups. (B) The analysis of differences in immune cell infiltration between the two groups with ssGSEA. (C) The analysis of differences in immune functions between the two groups with ssGSEA. (D) The analysis of differences in immune cell infiltration between the two groups with Multiple algorithms. (E) The distribution of patients with high- and low-risk in different immune subtypes. *P <0.05; **P <0.01; ***P <0.001.

We next examined the low- and high-risk patients in terms of the expression patterns of immune-related genes. A majority of immune-related genes were discovered to be expressed at low levels in the high-risk category (Figures 7A—D). TIDE scores acted as a vital biomarker for immunotherapy response. The link between the TIDE score and risk score was also investigated. TIDE scores were found to be lower in the high-risk category in contrast with the low-risk category (Figure 7E). Furthermore, high-risk patients respond better to immunotherapy in contrast with those at low risk (Figure 7F).

FIGURE 7
www.frontiersin.org

Figure 7 Assessment of Immunotherapy response of high- and low-risk groups. (A–D) The immune-related gene expression levels in different groups. (E) The TIDE scores in different groups. (F) Prediction immunotherapy response in IMVigor210. *P <0.05; **P <0.01; ***P <0.001.

Comparative analysis of mutations and functional evaluation of the URG signature

The somatic mutations of PDAC with high and low URG scores were examined to adequately describe the pathobiological significance of the URG score. The high-URG score group had a greater somatic mutation frequency (94.05%) in contrast with the low-URG score group (70.51%) (Figures 8A, B). The GO, KEGG, and DO enrichment analyses were conducted to investigate the latent biological roles of the DEGs. Accordingly, the GO result revealed that the DEGs were primarily enriched in signal release, external side of plasma membrane, and antigen binding (Figures 9A, C). The KEGG result suggested that the DEGs were primarily enriched in insulin secretion, cAMP signaling pathway, and Chemokine signaling pathway (Figure 9B). The DO result suggested that the DEGs were primarily enriched in cell type benign neoplasm, adenocarcinoma and pancreas disease (Figure 9C). In addition, the GSVA indicated that many pathways were substantially altered between the high- and low-risk PDAC patients (Figure 9D).

FIGURE 8
www.frontiersin.org

Figure 8 The somatic gene mutations in the high-risk group (A) and low-risk group (B).

FIGURE 9
www.frontiersin.org

Figure 9 Function analysis. (A) GO analysis of differential genes between high and low-risk groups. (B) KEGG analysis of differential genes between high- and low-risk groups. (C) DO analysis of differential expression genes. (D) GSVA enrichment analysis in high- and low-risk groups.

The correlation analysis of drug sensitivity and risk score

We correlated the PDAC patients’ risk scores with the IC50 values of chemotherapy and targeted treatment medications to learn more about the possible variations in drug sensitivity between low- and high-risk categories. The IC50 values of 17-AAG and PD-0325901 were significantly higher in low-risk group, whereas the IC50 values of Phenformin, Axitinib, AZD8055 and TAK-715 were lower in low-risk group (Figure 10).

FIGURE 10
www.frontiersin.org

Figure 10 Drug sensitivity analysis in high and low-risk groups.

Knockdown of S100A2 inhibits the malignant biological behavior of pancreatic cancer

First, we detected the expression of three risk genes by RT-qPCR. The results showed that the expressions of three risk genes in pancreatic cancer tissue were higher than normal pancreatic tissue (Figure 11A). Next, we knocked down of S100A2 in PANC-1 cells to explore the role of S100A2 in pancreatic cancer. The results showed that S100A2 knockdown reduced the cloning, migration, and invasion ability of PANC-1 cells (Figures 11B–D). In addition, it was found that S100A2 knockdown increased the ubiquitination of β-catenin, thereby reducing its protein expression (Figure 11E). These results suggested that S100A2 knockdown might inhibit the malignant biological behavior of pancreatic cancer cells by increasing the ubiquitination of β-catenin.

FIGURE 11
www.frontiersin.org

Figure 11 S100A2 knockdown inhibits the malignant biological behavior of PDAC cell. (A) Detection of risk gene expression in pancreatic cancer tissues and normal tissues by RT-qPCR. (B) S100A2 knockdown inhibits the cloning of PANC-1 cell. (C) S100A2 knockdown inhibits the migration ability of PANC-1 cell. (D) S100A2 knockdown inhibits the invasion ability of PANC-1 cell. (E) S100A2 knockdown increased the ubiquitination of β-catenin. *P <0.05.

Discussion

Pancreatic cancer is a particularly deadly malignancy of the digestive system (18). Although there is a new understanding of the molecular mechanism of PDAC, and new progress in surgery, adjuvant therapy, and chemotherapy, patients diagnosed with PDAC have not seen a substantial improvement in their prognoses (19). For a long time, the prognosis of PDAC has been judged mainly based on clinical manifestations, tumor stage, pathological grade, lymph node metastasis, neurovascular invasion, and other pathological characteristics (19, 20). Consequently, predicting a patient’s prognosis in PDAC using these markers is challenging. Additionally, in the early diagnosis of PDAC, the detection of single tumor markers, such as CEA, CA199, and CA125, is difficult to reach the ideal level of sensitivity and specificity (2126), which can no longer meet the current clinical needs. Recently, Wu et al. found a nine-gene signature, which could anticipate the OS duration of PDAC patients (27). As per the expression of immune-related genes such as MMP14 and INHBA, Xu et al. developed a predictive risk model to assess PADC patients’ prognoses and identify PDAC therapy opportunities (28). However, there is much less understanding of ubiquitination in the PDAC prognosis and microenvironment. Here, we built a risk model using three SLC22A17, UCHL1 and S100A2, to indicate the outcomes of PDAC. The risk model also presented as a potential biomarker to reflect the sensitivity of targeted therapy and immune status in tumor tissues. Our study analyzed the mRNA expression data of these three risk model genes.

First, we identified 996 DEGs between the two clusters of ubiquitinated subtypes, and further screened 12 genes correlated with the prognosis from these DEGs utilizing univariate Cox regression analysis. Thereafter, we utilized LASSO regression analysis to obtain the optimal genes for predicting outcomes. Lastly, a prognostic risk model of PDAC containing three ubiquitination-related genes (SLC22A17, UCHL1 and S100A2) was constructed, and KM and ROC curve analyses were employed in TCGA, ICGC and GEO data sets to confirm the effectiveness of the model for prognostic evaluation of PDAC. Wei et al. found that high expression of SLC22A17 indicates poorer prognosis in gastric cancer (29). UCHL1 might play a role in the malignant progression of triple-negative breast cancer by maintaining dryness and promoting cell invasion (30). Li et al. revealed S100A2 promoted glycolysis and proliferation of colorectal cancer through GLUT1 regulation (31). In this study, we found S100A2 knockdown could inhibit the malignant biological behavior of pancreatic cancer cells by increasing the ubiquitination of β-catenin. However, the significance of other two genes in PDAC remains uncertain and this is a direction of pancreatic cancer research in the future.

The tumor microenvironment (TME) is constituted of diverse immune cells, interstitial cells, extracellular matrix, and tumor blood vessels, which stimulate the onset and advancement of cancer (31). The infiltrating immune cell levels in TME usually change with tumorigenesis and progression (32). Our analysis illustrated that PDAC patients having high risk scores recorded lower ImuneScores, StromalScores and ESTIMATEScores. We found that most of the immune cells (B cells, CD8+ T cells, Treg, T helper cells, Neutrophils, TIL, and Mast cells) were substantially reduced in the high-risk patients in contrast with the low-risk patients. Additionally, the majority of immune-related genes tended to be downregulated in the high-risk population, whereas the low-risk category illustrated considerable improvement in immunologic functioning. Research suggests that immune cells are important components of anti-tumor immunity (33). One reason high-risk individuals have such a dismal prognosis is that they have fewer immune cells and attenuated immunological functioning. TIDE scores acted as vital biomarkers for immunotherapy response (34). Our results highlighted that high-risk individuals with PDAC responded more positively to immunotherapy compared to those in the low-risk category and that the TIDE scores were lower in the high-risk patients in contrast with those in the low-risk category. The findings of this research shed light on the involvement of ubiquitination in PDAC and may be utilized to direct immunotherapeutic and chemotherapeutic interventions for PDAC patients.

Nevertheless, our investigation does have a few drawbacks. Case selection bias could be present since the vast majority of analyses use data from publicly available data sets and all samples are retrieved retroactively. Furthermore, additional in vitro and in vivo tests are warranted to corroborate our findings.

In summary, we designed a molecular cluster and prognostic signature based on URGs, which aid in anticipating survival, directing immunotherapy, and determining clinical outcomes. This research potentially provides deeper insights into the function of ubiquitination in PDAC and facilitates the development of more effective therapies for this disease.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding authors.

Author contributions

YG designed the study. LQ, YB and ZW wrote the manuscript. KC, YD, YM and KH performed data collection and curation. ZW performed RT-qPCR experiments. All authors have read and agreed to the published version of the manuscript.

Funding

This study was sponsored by the Ningbo University Institute of Geriatrics (LNBYJS-2021), Key Laboratory of Diagnosis and Treatment of Severe Hepato-Pancreatic Diseases of Zhejiang Province (2018E10008) and the Wenzhou Municipal Science and Technology Plan Project (Y20220185).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2023.1171811/full#supplementary-material

Supplementary Figure 1 | The differences in immune cell infiltration between two clusters.

Supplementary Figure 2 | Validation of the Ubiquitination-related prognostic signature in GSE57495.

References

1. Guo X, Li Y, Che X, Hou K, Qu X, Li C. microRNA-569 inhibits tumor metastasis in pancreatic cancer by directly targeting NUSAP1. Aging (Albany NY) (2022) 14(8):3652–65. doi: 10.18632/aging.204035

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Wang L, Wang L, Wang S, Zhou Z, Liu Z, Xu P, et al. N2E4, a monoclonal antibody targeting neuropilin-2, inhibits tumor growth and metastasis in pancreatic ductal adenocarcinoma via suppressing FAK/Erk/HIF-1α signaling. Front Oncol (2021) 11:657008. doi: 10.3389/fonc.2021.657008

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Chao MW, Chang LH, Tu HJ, Chang CD, Lai MJ, Chen YY, et al. Combination treatment strategy for pancreatic cancer involving the novel HDAC inhibitor MPT0E028 with a MEK inhibitor beyond K-ras status. Clin Epigenetics (2019) 11(1):85. doi: 10.1186/s13148-019-0681-6

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Lee MKC, Grimmond SM, McArthur GA, Sheppard KE. PRMT5: an emerging target for pancreatic adenocarcinoma. Cancers (Basel) (2021) 13(20):5136. doi: 10.3390/cancers13205136

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Zhan HX, Xu JW, Wu D, Wu ZY, Wang L, Hu SY, et al. Neoadjuvant therapy in pancreatic cancer: a systematic review and meta-analysis of prospective studies. Cancer Med (2017) 6(6):1201–19. doi: 10.1002/cam4.1071

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Deng T, Gong Y, Liao X, Wang X, Zhou X, Zhu G, et al. Integrative analysis of a novel eleven-small nucleolar RNA prognostic signature in patients with lower grade glioma. Front Oncol (2021) 11:650828. doi: 10.3389/fonc.2021.650828

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Kong L, Liu P, Fei X, Wu T, Wang Z, Zhang B, et al. A prognostic prediction model developed based on four CpG sites and weighted correlation network analysis identified DNAJB1 as a novel biomarker for pancreatic cancer. Front Oncol (2020) 10:1716. doi: 10.3389/fonc.2020.01716

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Tian X, Zheng J, Mou W, Lu G, Chen S, Du J, et al. Development and validation of a hypoxia-stemness-based prognostic signature in pancreatic adenocarcinoma. Front Pharmacol (2022) 13:939542. doi: 10.3389/fphar.2022.939542

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Zhu J, Huang Q, Peng X, Luo C, Liu S, Liu Z, et al. Identification of LncRNA prognostic signature associated with genomic instability in pancreatic adenocarcinoma. Front Oncol (2022) 12:799475. doi: 10.3389/fonc.2022.799475

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Vrba L, Futscher BW, Oshiro M, Luo C, Liu S, Liu Z, et al. Liquid biopsy, using a novel DNA methylation signature, distinguishes pancreatic adenocarcinoma from benign pancreatic disease. Clin Epigenetics (2022) 14(1):28. doi: 10.1186/s13148-022-01246-2

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Zhang Q, Lv L, Ma P, Zhang Y, Deng J, Zhang Y. Identification of an autophagy-related pair signature for predicting prognoses and immune activity in pancreatic adenocarcinoma. Front Immunol (2021) 12:743938. doi: 10.3389/fimmu.2021.743938

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Zhang X, Li X, Xie J, Zhu Q, Yuan Y. A novel immune-related prognostic signature predicting survival in patients with pancreatic adenocarcinoma. J Oncol (2022) 2022:8909631. doi: 10.1155/2022/8909631

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Wang G, Bie F, Qu X, Yang X, Liu S, Wang Y, et al. Expression profiling of ubiquitin-related genes in LKB1 mutant lung adenocarcinoma. Sci Rep (2018) 8(1):13221. doi: 10.1038/s41598-018-31592-2

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Huang G, Zhou J, Chen J, Liu G. Identification of pyroptosis related subtypes and tumor microenvironment infiltration characteristics in breast cancer. Sci Rep (2022) 12(1):10640. doi: 10.1038/s41598-022-14897-1

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Zhao Y, Zhang Y, Dai C, Hong K, Guo Y. A signature constructed with mitophagy-related genes to predict the prognosis and therapy response for breast cancer. Aging (Albany NY) (2022) 14(15):6169–86. doi: 10.18632/aging.204209

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Liang L, Liu Y, Jiang S, Huang J, He H, Shen L, et al. Novel circulating tumour cell-related risk model indicates prognosis and immune infiltration in lung adenocarcinoma. J Immunol Res (2022) 2022:6521290. doi: 10.1155/2022/6521290

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Ye Y, Li L, Dai Q, Liu Y, Shen L. Comprehensive analysis of histone methylation modification regulators for predicting prognosis and drug sensitivity in lung adenocarcinoma. Front Cell Dev Biol (2022) 10:991980. doi: 10.3389/fcell.2022.991980

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Morgan RG, Wormsley KG. Progress report. cancer of the pancreas. Gut (1977) 18(7):580–96. doi: 10.1136/gut.18.7.580

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Mao W, Mai J, Peng H, Wan J, Sun T. YAP in pancreatic cancer: oncogenic role and therapeutic strategy. Theranostics (2021) 11(4):1753–62. doi: 10.7150/thno.53438

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Tong W, Zhu L, Bai Y, Yang L, Liu Z, Zhang Y. Enhancer RNA LINC00242-induced expression of PHF10 drives a better prognosis in pancreatic adenocarcinoma. Front Oncol (2022) 11:795090. doi: 10.3389/fonc.2021.795090

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Glaß M, Dorn A, Hüttelmaier S, Haemmerle M, Gutschner T. Comprehensive analysis of LincRNAs in classical and basal-like subtypes of pancreatic cancer. Cancers (Basel) (2020) 12(8):2077. doi: 10.3390/cancers12082077

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Chang JC, Kundranda M. Novel diagnostic and predictive biomarkers in pancreatic adenocarcinoma. Int J Mol Sci (2017) 18(3):667. doi: 10.3390/ijms18030667

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Xu HX, Li S, Wu CT, Qi ZH, Wang WQ, Jin W, et al. Postoperative serum CA19-9, CEA and CA125 predicts the response to adjuvant chemoradiotherapy following radical resection in pancreatic adenocarcinoma. Pancreatology (2018) 18(6):671–7. doi: 10.1016/j.pan.2018.05.479

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Grady WM, Yu M, Markowitz SD, Chak A. Barrett's esophagus and esophageal adenocarcinoma biomarkers. Cancer Epidemiol Biomarkers Prev (2020) 29(12):2486–94. doi: 10.1158/1055-9965.EPI-20-0223

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Zhao Y, Zhang Y, Dai C, Hong K, Guo Y. A signature constructed with mitophagy-related genes to predict the prognosis and therapy response for breast cancer. Aging (Albany NY) (2022) 14(15):6169–86. doi: 10.18632/aging.204209

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Zheng Y, Weng X, Hu D, He J. Identification of a signature based on non-apoptotic regulatory cell death to improve prognosis prediction in acute myeloid leukaemia [published online ahead of print, 2022 Dec 9]. Br J Haematol (2022), 10.1111/bjh.18601. doi: 10.1111/bjh.18601

CrossRef Full Text | Google Scholar

27. Wu M, Li X, Zhang T, Liu Z, Zhao Y. Identification of a nine-gene signature and establishment of a prognostic nomogram predicting overall survival of pancreatic cancer. Front Oncol (2019) 9:996. doi: 10.3389/fonc.2019.00996

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Xu H, Sun J, Zhou L, Du QC, Zhu HY, Chen Y, et al. Development of a lipid metabolism-related gene model to predict prognosis in patients with pancreatic cancer. World J Clin Cases (2021) 9(35):10884–98. doi: 10.12998/wjcc.v9.i35.10884

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Wei J, Gao X, Qin Y, Liu T, Kang Y. An iron metabolism-related SLC22A17 for the prognostic value of gastric cancer. Onco Targets Ther (2020) 13:12763–75. doi: 10.2147/OTT.S287811

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Tian C, Liu Y, Liu Y, Hu P, Xie S, Guo Y, et al. UCHL1 promotes cancer stemness in triple-negative breast cancer. Pathol Res Pract (2022) 240:154235. doi: 10.1016/j.prp.2022.154235

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Li C, Chen Q, Zhou Y, Niu Y, Wang X, Li X, et al. S100A2 promotes glycolysis and proliferation via GLUT1 regulation in colorectal cancer. FASEB J (2020) 34(10):13333–44. doi: 10.1096/fj.202000555R

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Dong S, Hou D, Peng Y, Chen X, Li H, Wang H. Pan-cancer analysis of the prognostic and immunotherapeutic value of MITD1. Cells (2022) 11(20):3308. doi: 10.3390/cells11203308

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Wang L, Wu W, Zhu X, Ng W, Gong C, Yao C, et al. The ancient Chinese decoction yu-Ping-Feng suppresses orthotopic Lewis lung cancer tumor growth through increasing M1 macrophage polarization and CD4+ T cell cytotoxicity. Front Pharmacol (2019) 10:1333. doi: 10.3389/fphar.2019.01333

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Sun Y, Wu J, Yuan Y, Lu Y, Luo M, Lin L, et al. Construction of a promising tumor-infiltrating CD8+ T cells gene signature to improve prediction of the prognosis and immune response of uveal melanoma. Front Cell Dev Biol (2021) 9:673838. doi: 10.3389/fcell.2021.673838

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: ubiquitination, PDAC, immune microenvironment, immunotherapy, prognosis

Citation: Guo Y, Wu Z, Cen K, Bai Y, Dai Y, Mai Y, Hong K and Qu L (2023) Establishment and validation of a ubiquitination-related gene signature associated with prognosis in pancreatic duct adenocarcinoma. Front. Immunol. 14:1171811. doi: 10.3389/fimmu.2023.1171811

Received: 22 February 2023; Accepted: 26 May 2023;
Published: 09 June 2023.

Edited by:

Saleh AlGhamdi, King Fahad Medical City, Saudi Arabia

Reviewed by:

Theodore Welling, New York University, United States
Haotian Fu, Zhejiang University, China
Edmund Ui-Hang Sim, University of Malaysia Sarawak, Malaysia

Copyright © 2023 Guo, Wu, Cen, Bai, Dai, Mai, Hong and Qu. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Liangchen Qu, 1106083035@qq.com; Kai Hong, hongkai0629@163.com; Yifeng Mai, fymaiyifeng@nbu.edu.cn

These authors have contributed equally to this work

Disclaimer: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.