Construction of a ceRNA Network and Analysis of Tumor Immune Infiltration in Pancreatic Adenocarcinoma

Abstract

Pancreatic adenocarcinoma (PAAD) is characterized by high malignancy, frequent metastasis, and recurrence with an unfavorable prognosis. This study is aimed at constructing a prognostic model for tumor-infiltrating immune cells and a competing endogenous RNA (ceRNA) network in PAAD and analyzing susceptibilities of chemotherapy and immunotherapy of PAAD. Gene expression profiles and clinical information of PAAD were downloaded from The Cancer Genome Atlas (TCGA) database and divided into the tumor group and the normal group. A total of five PAAD survival-related key genes in the ceRNA network and three survival-related immune infiltrating cells were uncovered, and two survival risk models and nomograms were constructed. The efficiency and performance of the two models were verified using multi-index area under the curve analysis at different time points, decision curve analysis, and calibration curves. Co-expression analysis showed that LRRC1, MIR600HG, and RNF166 in the ceRNA network and tumor-infiltrating immune cells including CD8 T cells and M1 macrophages were likely related to the PAAD prognosis, and the expression of key ceRNA-related genes was experimently validated in tissues and cell lines by RT-qPCR. Patients with low risk scores for key genes in the ceRNA network displayed a positive response to anti-programmed death-1 (PD-1) treatment and greater sensitivity to chemotherapeutic drugs such as docetaxel, lapatinib, and paclitaxel. More importantly, our results suggested that the IC50 values of gemcitabine in PAAD were not significantly different between the high and low risk groups. The expression levels of immune checkpoints were significantly different in the high-risk and low-risk groups. The prognostic model, nomogram, and drug analysis may provide an essential reference for PAAD patient management in the clinic.

Keywords: Pancreatic adenocarcinoma; TCGA; competing endogenous RNA network; prognostic model; tumor-infiltrating immune cell.

FIGURE 1

Flow chart of the analytical process.

FIGURE 2

Differentially expressed genes in PAAD. Heatmap (A) and volcano plot (B) of all differentially expressed mRNAs; heatmap (C) and volcano plot (D) of all differentially expressed lncRNAs; heatmap (E) and volcano plot (F) of all differentially expressed miRNAs. PAAD, pancreatic adenocarcinoma; miRNAs, microRNAs; lncRNAs, long non-coding RNAs; mRNAs, messenger RNAs.

FIGURE 3

Constructed ceRNA networks via Cytoscape (A); GO (B) and KEGG (C) analyses for ceRNA-related differentially expressed genes; visual network (D).

FIGURE 4

Kaplan–Meier survival curves of the former 12 genes in the ceRNA networks: FAM83A (A), EVL (B), TMC7 (C), B3GNT3 (D), hsa-miR-196b-5p (E), PCDH1 (F), SDR16C5 (G), hsa-miR-424-5p (H), LY6D (I), CDC6 (J), PADI1 (K), and ANLN (L).

FIGURE 5

Construction of the nomogram for predicting the survival probability of PAAD based on the identified hub genes in ceRNA networks. Results of Lasso regression (A, B); Kaplan–Meier survival curve in high- and low-risk groups based on the multivariate Cox regression analysis (C); forest plot of the multivariate Cox regression analysis (D); ROC curves of the multivariate Cox model (E); area under the multi-indicator ROC curve (F) and decision-making curve (G); nomogram (H) and its calibration curve (I).

FIGURE 6

The model acted as a potential predictor of chemosensitivity as high risk scores were related to a lower IC50 for chemotherapeutics such as doxorubicin, mitomycin, and cisplatin, whereas they were related to a higher IC50 for gemcitabine Box plots depicted the differences in the estimated IC50 levels of (A) Gemcitabine; (B) Dasatinib; (C) Lapatinib; (D) Docetaxel; (E) Cyclopamine; (F) Paclitaxel; (G) A.443654; (H) BI.2536 between the high and low risk score groups. The results illustrated that the group with low risk score was more likely to respond to immune checkpoint inhibitors (PD-1 and CTLA-4 inhibitors) than the group with high risk score (I); expression of immune checkpoints among high and low PAAD risk groups (J).

FIGURE 7

Composition (A) and heatmap (B) of 22 subsets of immune cells in PAAD; violin plot (C) of immune cell infiltration in tumor and normal groups.

FIGURE 8

(A): Co-expression patterns among fractions of immune cells; Kaplan–Meier survival curves of fractions of memory B cells (B) and activated dendritic cells (C).

FIGURE 9

Construction of the nomogram for predicting the survival probability of PAAD based on the prognostic-related immune cells. Results of Lasso regression (A, B); Kaplan–Meier survival curve in high- and low-risk groups based on the multivariate Cox regression analysis (C); forest plot of the multivariate Cox regression analysis (D); ROC curves of the multivariate Cox model (E); heatmap of the three immune cells in the Cox regression model (F); area under the multi-indicator ROC curve (G) and decision-making curve (H); nomogram (I); nomogram-predicted probability of three-year overall survival (J).

FIGURE 10

Co-expression patterns among fractions of three immune cells and five hub genes in the ceRNA network (A); MIR600HG and M1 macrophages (B); MIR600HG and CD8 T cells (C); RNF166 and CD8 T cells (D); LRRC1 and M1 macrophages (E).

FIGURE 11

Validation of the expression levels of the four hub genes between normal pancreatic samples and PAAD samples based on TCGA and GTEx data in GEPIA. Expression levels of LRRC1, RNF166, LY6D, and MIR600HG in normal pancreatic samples and PAAD samples (A–D). RT-qPCR was applied to test LRRC1 (E), RNF166 (F), LY6D (G), and MIR600HG (H) expressions in PC cell lines.Validation of the expression levels of the four hub genes between normal pancreatic samples (n = 12) and PAAD samples (n = 12) by PCR analysis (I). All of the data are presented as mean ± SD. Significant differences are defined by a p value < 0.01.