The landscape of genetic variation of TRP-related genes in pancreatic adenocarcinoma
At first, the incidence frequency of somatic mutations of 28 TRP family genes in PAAD was summarized based on cBioPortal database (http://cbioportal.org) [16]. TRP-related gene mutations were found in 52 of the 149 cases, with a frequency of 35%. According to the findings, the highest mutation frequency was found in TRPA1, followed by TRPM4 and TRPV6 (Fig. 1 A). The incidence of copy number variation (CNV) mutations was shown to be prevalent in 28 TRP family genes, with more number of copy number deletion genes as compared to amplification genes (Fig. 1 B). The site of CNV alteration of TRP family genes on chromosomes was illustrated in Fig. 1 C. Univariate Cox regression evaluated the prognostic significance of 28 TRP family genes in PAAD patients. A network was used to describe the interaction of TRP family genes and their effect on the prognosis of PAAD patients. We discovered that TRP family genes’ expression was significantly correlated (Fig. 1 D).
The differential expression and prognostic value of TRP-related genes
We first explored the differential expression of TRP-related genes in pancreatic adenocarcinoma using GEPIA [17] and CCLE databases [18]. As shown in Supplementary Fig. 1, we found that TRPC1, TRPM2, TRPM4, MCOLN1, PDK2, TRPV2 and TRPV4 were highly expressed in PAAD tissues, while MCOLN3 and TRPV6 were lowly expressed. In addition, the expression of TRP-related genes in pancreatic cancer cell line was shown in Supplementary Fig. 2. More importantly, we also evaluated the prognostic value of TRP-related genes in pancreatic adenocarcinoma using the KM-plotter database (https://kmplot.com/analysis/). The result showed that TRPC1, TRPC4, TRPC5, TRPC6, TRPC7, TRPM6, MCOLN1, MCOLN3, PDK1 and TRPV1 were favorable to the prognosis of patients, while TRPV1, TRPM1, TRPM4, TRPM8 and TRPV6 were unfavorable (Supplementary Fig. 3). In addition, we constructed the PPI network based TRP channel proteins by STRING website. And we found that TRPC1 was associated with almost all TRP channel proteins, suggesting that TRPC1 may be a key molecule in TRP ion channels (Supplementary Fig. 4A-B).
Unsupervised clustering based on 28 TRP-related genes
Three different subtypes were subsequently identified on the basis of the expression of 28 TRP-related genes with 122 patients in cluster A, 60 patients in cluster B, and 60 patients in cluster C (Figs. 2 A-B). Figure 2 C illustrated remarkable differences in OS among these three clusters (p < 0.05). GSVA enrichment analysis helped us understand the biological mechanisms involved in different clusters. As demonstrated in Fig. 2 D, cluster A was highly enriched in carcinogenic signaling pathways including the TGF beta signaling pathway, MAPK signaling pathways, focal adhesion, and pathways in cancer. Cluster C was enriched in immune-related pathways including T cell receptor signaling pathway, Toll-like receptor signaling pathways, JAK/STAT signaling pathway, chemokine signaling pathway, and mTOR signaling pathway (Fig. 2 E). Furthermore, PCA revealed remarkable differences among the three clusters (Fig. 3 A). Following investigations of immune cell infiltration, it was discovered that the abundance level of tumor-infiltrating immune cells in cluster B was considerably lower than that in other clusters (Fig. 3 B), which was linked to poor survival in cluster B. In addition, we also compared the differential expression of TRP-related genes among the three clusters to screen out the relatively important TRP-related genes in each cluster. The result showed that TRPC4, MCOLN3, TRPV5 and TRPV6 were highly expressed in cluster A. And TRPV1, TRPC1, TRPC3, MCOLN1, PKD1, PKD2L1 and TRPV2 were highly expressed in cluster B. Only TRPM4 was highly expressed in cluster C (Supplementary Fig. 5) .
Identification of TRP-associated lncRNAs with significant prognostic value in pancreatic adenocarcinoma
The Limma R package was performed for the identification of TRP-associated lncRNAs by Pearson correlation analysis. Finally, 345 TRP-associated lncRNAs were used for the following analysis (Fig. 3 C). And the detailed lncRNA-mRNA pairs information was presented in Supplementary Table 1. lncRNAs with significant differences between adjacent normal pancreatic tissues and pancreatic tumor tissues were shown in the heatmap (|logFC| > 2 & FDR < 0.05, Fig. 3 D). Univariate Cox proportional hazards analysis was used to screen 12 TRP-related lncRNAs which were prominently related to patients’ prognosis (p < 0.05). As shown in Fig. 4 A, LINC01091, AC012213.4, LINC02251, AP000802.1 and TRPC7 − AS1 were low risk with hazard ration (HR) < 1, whereas LINC01133, AC068580.2, AP005233.2, LINC02323, LINC00973, LINC02041 and AC009065.5 were high risk with hazard ration (HR) > 1. Afterward, a prognostic signature comprising 4 critical lncRNAs, including LINC01091, LINC01133, TRPC7-AS1, and LINC00973, was selected by multivariate Cox analysis in the training set (Fig. 4 B).
Validation of the 4-TRlncRNAs prognostic signature
On the basis of their median risk level, patients in the training cohort were categorized into high and low-risk groups. The survival curve highlighted, in particular, that PAAD patients in the low-risk group have considerably longer survival times as compared to the high-risk group (Fig. 4 C). Moreover, the AUC of risk score at 1, 3, and 5-years was 0.678, 0.780, and 0.837 respectively, indicating a high predictive value for the prognosis of PAAD patients (Fig. 4 E). Additionally, we also looked into the prognostic model’s predictive values in the validation cohort. The survival curve in the validation cohort demonstrated that PAAD patients in the low-risk group had better survival outcomes than those in the high-risk group, which was consistent with the findings obtained in the training cohort (Fig. 4 D). The respective AUC values in the validation cohort at 1, 3, and 5 years were 0.736, 0.650, and 0.857 (Fig. 4 F).
Evaluation of the diagnostic efficiency of 4-TRlncRNAs prognostic signature
In addition, we did univariate and multivariate Cox regression analyses for investigating if the 4-TRlncRNAs signature is not dependent on other clinicopathological factors in predicting the survival of PAAD patients. The risk score’s hazard ratio (HR) and 95% CI in the training cohort were 1.402 and 1.247–1.577 (p < 0.001) in univariate Cox regression analysis (Fig. 5 A), and 1.423 and 1.251–1.618 (p < 0.001) in multivariate Cox regression analysis (Fig. 5 B), showing that 4-TRlncRNAs signature is an independent prognostic factor for PAAD patients. Moreover, findings in the validation cohort verified the signature’s independent prognostic value (Figs. 5 C-D). After dimensionality reduction, we performed PCA on the model genes in the training and validation cohorts based on the high and low-risk groups, and the scatter plot revealed that high and low-risk patients were considerably distinguished, demonstrating that the model genes were the main factor influencing the distribution of high-low risk groups in pancreatic adenocarcinoma (Figs. 5 E-F).
Establishment and verification of nomogram model for predicting the survival rate of PAAD patients
ROC curves for 4-TRlncRNAs signature and clinicopathological features were created to compare diagnostic effectiveness between the risk score and clinicopathological features, as shown in Fig. 6 A. In PAAD patients, the AUC of 4-TRlncRNAs signature was higher as compared to those of the clinical indexes (AUC = 0.711). The risk score’s C-index was also considerably higher as compared to other clinical features (Fig. 6 B). In clinical decision making, DCA further showed that risk score was a viable prognostic indicator as compared to other variables (Fig. 6 C). After that, we effectively constructed a prognostic nomogram for predicting the 1/3/5-year survival rate. The final OS prediction model incorporated age, gender, grade, stage, and 4-TRlncRNAs (Fig. 6 D). To examine the congruence between actual and expected survival, calibration curves for 1–3-, and 5-year survival were generated, which almost fell on the diagonal of 45°, indicating that the model has good accuracy (Fig. 6 E). Then, using GSEA, it was discovered that cancer-related pathways, such as the cell cycle and the p53 signaling pathway, were considerably enriched in the high-risk group (Figs. 7 A-B).
The association between the risk signature with TMB
In the high-risk group, LINC01091 and TRPC7 − AS1 had low expression but LINC01133 and LINC00973 had high expression according to the clustering heatmap of the four TRP-associated IncRNAs. Additionally, after comparing high and low-risk groups, there were considerable differences in the immune score and cluster grouping (Fig. 7 C). Moreover, data on somatic mutation were further studied to investigate mutation differences in our risk model. TMB and risk score had a substantial positive link as displayed in Fig. 7 D, with patients in the high-risk group having a higher TMB. Figures 7 E-F depicted the top 20 mutated genes in two risk groups. With an overall mutation rate of 88.24% versus 61.4%, the high-risk group showed a more extensive tumor mutation burden as compared to the low-risk group.
Relationship between risk score model and immune infiltration
Furthermore, the MCP-counter algorithm revealed that reduced CD8 T cells, NK cells, cytotoxic lymphocytes, and endothelial cells infiltration were associated with a high risk (Fig. 8 A). The expression levels of POLE2, MCM6 CD274, and LOXL2 were found to be higher in the high-risk group (Fig. 8 B), indicating that immune checkpoint inhibitors (ICIs) therapy may benefit high-risk patients more. Moreover, as illustrated in Fig. 8 C, we found considerable changes in the immune subtypes between the high and low-risk groups (p < 0.05).
Validation of TRP-associated lncRNAs in prognostic signature
Using the ENCORI Pan-Cancer Analysis Platform (https://starbase.sysu.edu.cn/panCancer.php), we looked at the expression of four TRP-associated lncRNAs in normal pancreatic tissues and pancreatic tumor tissues. LINC01091 and LINC01133 expression levels in PAAD tissues were significantly different, as demonstrated in Figs. 9 A-D, whereas TRPC7-AS1 and LINC00973 expression levels were not. However, the ENCORI database contained only four normal pancreatic samples, so we further analyzed the expression of four lncRNAs using the GEPIA website, which incorporated normal pancreas data from the GTEx database. The results showed that only LINC01133 was significantly overexpressed in PAAD tissues, although other lncrnas had no statistical significance, the up-regulation trend was obvious (Figs. 9 E-H). In pancreatic cancer cells and tissues, qRT-PCR indicated the relative expression of LINC01091, LINC01133, TRPC7-AS1, and LINC00973. Figure 10 A-D illustrated the expression levels of four lncRNAs in four cell lines: HPDE6-C7, BXPC-3, SW1990, and PANC-1. In PAAD tissues, LINC01091 and LINC01133 expression levels were significantly different, whereas TRPC7-AS1 and LINC00973 expression levels were not (Fig. 10 E-H). These findings showed that prognostic TRP-associated lncRNAs may play a role in PC tumorigenesis.