- 1Department of Obstetrics and Gynecology, The Second Affiliated Hospital, Chongqing Medical University, Chongqing, China
- 2Department of Hepatobiliary Surgery, The First Affiliated Hospital, Chongqing Medical University, Chongqing, China
- 3Department of General Surgery, The Second Affiliated Hospital of Tianjin, University of Traditional Chinese Medicine, Tianjin, China
Background: Cervical cancer (CC) is a major health threat to females, and distal metastasis is common in patients with advanced CC. Anoikis is necessary for the development of distal metastases. Understanding the mechanisms associated with anoikis in CC is essential to improve its survival rate.
Methods: The expression matrix of long non-coding RNAs (lncRNAs) from cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC) patients was extracted from The Cancer Genome Atlas (TCGA), and highly relevant anoikis-related lncRNAs (ARLs) were identified by the single sample gene set enrichment analysis (ssGSEA) method. ARLs-related molecular subtypes were discerned based on prognosis-related ARLs. ARLs-related prognostic risk score (APR_Score) was calculated and risk model was constructed using LASSO COX and COX models. In addition, we also assessed immune cell activity in the immune microenvironment (TME) for both subtypes and APR_Score groups. A nomogram was utilized for predicting improved clinical outcome. Finally, this study also discussed the potential of ARLs-related signatures in predicting response to immunotherapy and small molecular drugs.
Results: Three ARLs-related subtypes were identified from TCGA-CESC (AC1, AC2, and AC3), with AC3 patients having the highest ARG scores, higher angiogenesis scores, and the worst prognosis. AC3 had lower immune cell scores in TME but higher immune checkpoint gene expression and higher potential for immune escape. Next, we constructed a prognostic risk model consisting of 7-ARLs. The APR_Score exhibited a greater robustness as an independent prognostic indicator in predicting prognosis, and the nomogram was a valuable tool for survival prediction. ARLs-related signatures emerged as a potential novel indicator for immunotherapy and small molecular drug selection.
Conclusion: We firstly constructed novel ARLs-related signatures capable of predicting prognosis and offered novel ideas for therapy response in CC patients.
Introduction
Cervical cancer (CC) is a serious health threat to females, and according to WHO data, CC is the fourth leading cause of cancer deaths among females, accounting for 342,000 deaths in 2020 (Sung et al., 2021). The most common histological type of CC is squamous cell carcinoma that accounts for about 75% of all cases, followed by the second most common type of adenocarcinoma that accounts for 20% of all cases (Meng et al., 2021). Unlike other cancers among females, effective measures including effective early screening and HPV vaccination could significantly reduce the incidence of CC and early targeted treatment (Canfell et al., 2020). Current mainstream treatment options for CC patients including systemic surgical treatment and radiotherapy could improve the 5-year survival rate of early CC patients to 91.5% (Bhatla et al., 2018). However, survival outcomes for patients with advanced CC due to postoperative recurrence and distant metastases are unsatisfactory (Marth et al., 2017). Therefore, there was an immediate demand for monitoring biomarkers to evaluate the CC metastasis risk.
Anoikis, a defense mechanism to prevent abnormal cell migration, is mainly manifested as the loss of cell adhesion to ECM and cell apoptosis and death in the process of migration (Paoli et al., 2013). Anoikis plays an essential component in inhibiting tumor cell migration and is important in delaying cancer progression (Adeshakin et al., 2021). Studies have found that tumor cells secrete growth factors and activate EMT signaling pathway to resist Anoikis, which also became a prerequisite for tumor cell metastasis (Kim et al., 2012). Anti-anoikis has become a landmark event in the occurrence of remote cancer metastasis (Kim et al., 2012). However, few studies focused on the correlation between anoikis and distal metastasis of CC. In view of the inability to accurately evaluate the prognosis of advanced CC patients, it was indispensable to explore the mechanism of anoikis in CC to improve its treatment.
Therefore, we developed an assessment system for prognostic risk stratification and constructed an anoikis-related long non-coding RNAs (ARLs)-related prognostic model in CC. We further investigated the relationship among prognostic indicators and immune microenvironment (TME), immunotherapeutic response, and chemotherapeutic drug suitability. The purpose of this study was to formulate a novel prognostic scoring system for CC, which was designed to accurately guide the prognosis of CC and improve its treatment options.
Materials and methods
Dataset download and processing
RNA-sequencing data and corresponding clinical information of cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC) samples were available in The Cancer Genome Atlas (TCGA, https://portal.gdc.cancer.gov/) database using TCGA GDC API. We also downloaded mutect2-processed single-nucleotide variants (SNVs) data from TCGA. The expression matrix with TPM format was transformed into log 2 (TPM+1). Samples were filtered using the Sangerbox (http://sangerbox.com/) (Shen et al., 2022) database as follows: 1) removing samples without follow-up information; 2) retaining samples with survival time greater than 0; 3) remove samples without Status. A total of 291 samples were included after screening. RNA annotation file in GENCODE (https://www.gencodegenes.org/) was used to obtain mRNAs and lncRNAs expression matrix. Anoikis-related genes (ARGs) were available in GeneCards (https://www.genecards.org/).
Recognition of anoikis-related lncRNAs (ARLs) and molecular subtypes
Based on the expression profiles of ARGs, the ARG scores of samples were attained via the single sample gene set enrichment analysis (ssGSEA) algorithm (Barbie et al., 2009). ARLs were identified using the rcorr function in the Hmisc package under screening thresholds of |cor|>0.3 and p < 0.01 (Hmisc.pdf). Univariate COX models were performed on the ARLs to screen for prognostically relevant ARLs under the threshold of p < 0.01. A consistency clustering analysis was carried out on samples in TCGA-CESC to identify molecular subtypes, according to the methods of Wilkerson et al. (Wilkerson and Hayes, 2010).
Gene set variation analysis (GSVA) and gene mutation landscape in ARLs-related subtypes
For patients in different subtype groups, hallmark gene sets were captured from the Molecular Signatures Database (MSigDB, https://www.gsea-msigdb.org/gsea/msigdb/) and GSVA was conducted using the GSVA package to explore differences in biological pathway variants across groups (Hanzelmann et al., 2013). Next, oncogenic pathway signatures were obtained from Sanchez-Vega et al. (Sanchez-Vega et al., 2018) and differences in oncogenic pathway scores were assessed in subtypes using the ssGSEA method. SNV data from CC were processed in The Genome Analysis Toolkit (GATK) software in the mutect2 plugin. Genes with mutation frequencies ≥3 were screened (fisher test, p < 0.05), and the mutation landscape was mapped using the maftools package (Mayakonda et al., 2018).
Construction of ARLs-related prognostic risk model (APR_Score)
The TCGA-CESC samples were clustered into training and test sets at 1:1 ratio for analysis. In the training set, LASSO and multivariate COX models were executed on prognosis-related ARLs to discriminate ARLs significantly affecting CC prognosis. ARLs-related prognostic risk model was constructed based on the following formula (Simon et al., 2011).
In the formula, expression of ARLs represented the expression data of prognosis-related ARLs, and βARLs indicated the COX regression coefficients after normalization. We calculated the ARLs-related prognostic risk score (APR_Score) of each patient in TCGA. Patients were assigned to the high-APR_Score and low-APR_Score groups according to the optimal p-value based on survminer package.
Clinical meaning of APR_Score and association with prognosis
We analyzed the variation in APR_Score and prognosis amongst the molecular subtypes. In this study, for different molecular subtypes and APR_Score groups, K-M curves were evaluated in the training set and the validation set to assess prognostic differences, and ROC curves were used to assess the accuracy of APR_Score. Furthermore, univariate and multivariate COX model analyses were performed on the TCGA-CESC cohort.
Gene set enrichment analysis (GSEA)
For patients in different APR_Score groups, GSEA was performed to explore differences in the biological pathways involved (p < 0.05, FDR<0.25).
Relationship between CC patients and tumor microenvironment (TME)
For ARLs-related molecular subtypes and APR_Score groups, we quantified the relative abundance of 22 immune cell species in the TME of CC patients using the CIBERSORT algorithm (Chen et al., 2018). Next, we further assessed the immune scores of 10 immune cell species in the TME using the MCP-Count method (Becht et al., 2016). Finally, 47 classes of immune checkpoint genes were obtained from the study by Danilova et al. (Danilova et al., 2019) to assess their expression levels in subtypes of CC patients.
Construction of nomogram
In this study, we constructed a nomogram for predicting 1-year, 3-year and 5-year survival rates of patients using the APR_Score (rms.pdf), and calibration curves were applied to evaluate the predictive accuracy of the nomogram. Finally, the ROC curves were exploited to validate the clinical use of APR_Score and nomogram.
Immunotherapy and small molecular drug sensitivity analysis
We assessed TIDE scores in the Tumor Immune Dysfunction and Exclusion (TIDE, http://tide.dfci.harvard.edu/) website between the APR_Score groups. In addition, we computed the half maximal inhibitory concentration (IC50) values of commonly used small molecular drugs in the pRRophetic package (Geeleher et al., 2014).
Statistical analysis
The packages included in this study were downloaded from R (version 4.1.1, https://www.r-project.org/) and analyzed using R Studio, an integrated development environment (IDE) for the R language. The packages included Hmisc, ConsensusClusterPlus, GSVA, survival, survminer, glmnet, maftools, timeROC, rms, and pRRophetic. Sangerbox was deployed for sample screening and data processing. In this study, p < 0.05 was considered statistically significant.
Results
Identification of ARLs in TCGA-CESC
Based on annotation file in GENCODE, a total of 14,176 lncRNAs were obtained. Firstly, we computed the ARG score for samples of TCGA-CESC via ssGSEA, and 574 ARLs (|cor|>0.3 and p < 0.01) were screened based on correlation analysis. This was followed by univariate COX model analysis, filtering a total of 37 prognosis-related ARLs (26 risk factors and 11 protective factors, p < 0.01) (Figure 1A). Furthermore, we analyzed the correlation among 37-ARLs, and the correlation heat map was illustrated in Figure 1B. Finally, the expression of 37-ARLs in 291 TCGA-CESC samples was counted, and we observed that 37-ARLs were differentially expressed, showing changes of ARG score in TCGA-CESC patients (Figure 1C).
FIGURE 1. Identification of ARLs TCGA-CESC (A) Forest map of ARLs affected OS (B) Heatmap of ARLs correlation affected OS (C) 37 Expression heatmap OS related ARLs.
ARLs-related molecular subtypes in TCGA-CESC
To determine the potential connections among ARG scores and clinical information, this study counted the expression of clinical information and 37-ARLs in TCGA-CESC patients. We found that risk factors were highly expressed in AC3 and protective factors were highly expressed in AC1 (Figure 2). Classification of TCGA-CESC patients was conducted based on consistent clustering of the 37-ARLs expression matrix. We found that the TCGA-CESC patients could be significantly clustered into three clusters (Figures 3A–C), therefore we derived three ARLs-related molecular subtypes, namely, ARLs-Cluster 1 (AC1), ARLs-Cluster 2 (AC2), and ARLs-Cluster 3 (AC3). Furthermore, the K-M curves showed significant prognostic differences between patients with AC1-3, with AC1 patients having the optimal prognosis and AC3 having the poorest prognosis (Figure 3D). Survival status statistics also showed more deaths among AC3 patients (p < 0.05, Figure 3E). We examined the ARG scores of patients in AC1-3 and found that the highest ARG scores were in AC3 patients (Figure 3F).
FIGURE 2. Heat map shows the expression of 37 ARLs in ARLs-related molecular subtypes in CC patients and clinical information.
FIGURE 3. ARLs-related molecular subtypes in TCGA-CESC. (A) Cumulative distribution function (CDF) of consensus clustering (B) CDF Delta area curve for k = 2–10 (C) Heatmap of sample clusters shows consensus matrix k=3 (D) K-M curves of AC1-3 patients in TCGA-CESC (E) Statistical graph of patient survival status for patients in AC1-3 (F) Distribution of ARG score amongst clusters AC1-3. *p < 0.05.
Features of TME in AC1-3 subtypes
To determine the relationship between ARG score and TME in CC, we assessed the relative infiltration abundance of 22 immune cell species in TME using the CIBERSORT method. We determined that the infiltration abundance of B cells naïve, T cells CD8, T cells regulatory (Tregs), and Mast cells resting was significantly lower in AC3 than in AC1 and AC2. In contrast, the infiltration risk of Macrophages M0, Mast cells activated, and Neutrophils was remarkably higher in AC3 than in AC1 and AC2 (Figure 4A). Accumulating evidence indicated that angiogenesis could influence tumor development and metastasis (Detmar, 2000; Parmar and Apte, 2021). Therefore, we computed the Angiogenesis score in AC1-3 and clearly observed that the poor prognosis AC3 subtype had the highest Angiogenesis score (p < 0.05, Figure 4B). In addition, we found that 17 immune checkpoint genes were highly expressed in clusters AC3 (Figure 4C). These results suggested that a higher ARG score was associated with immunosuppressive activity in TME of CC, which could lead to higher Angiogenesis scores, and that the development and metastasis of CC might be inextricably linked to immune escape.
FIGURE 4. Features of TME in AC1-3 subtypes. (A) Distribution of immune cell scores estimated by CIBERSORT amongst three subtypes (B) Differences of Angiogenesis scores amongst three subtypes (C) Differences of immune checkpoint gene expression amongst three subtypes. ns p > 0.05; *p < 0.05; **p < 0.01; *** p < 0.001; ****p < 0.0001.
Genetic landscapes and biological pathways in AC1-3 subtypes
Based on the genomic mutational landscape, we found variation among AC1-3 subtypes. AC3 showed markedly higher mutation frequencies than in AC1 and AC2, with the highest mutation frequencies for PCLO, HDWE1 and CDK12 in AC3 (Figures 5A–C). Based on the results of ssGSEA analysis, we noted that AC3 was substantially enriched in tumorigenesis-related pathways containing RAS, TGF-Beta, NRF1, NOTCH, MYC, HIPPO, CellCyle (Figure 5D). Similarly, GSVA results further validated the result because TGF BETA SIGNALING, PI3K AKT MTOR SIGNALING, TNFA SIGNALING VIA NFKB, WNT BETA CATENIN SIGNALING and NOTCH SIGNALING were remarkably activated in AC3 but remarkably inhibited in AC1 (Figure 5E). These results further supported a higher risk of metastasis in patients with CC in AC3, which in turn could lead to the negative prognosis.
FIGURE 5. Genetic Landscapes and Biological Pathways in AC1-3 subtypes (A–C) The mutation frequencies of top 15 genes in clusters AC1-3 (D) Distribution of tumorigenesis-related pathways score in AC1-3 (E) Heatmap shows pathway characteristic based on GSVA results. Red represents ssGSEA score is greater than 0 while blue represents ssGSEA score is lower than 0. ns p > 0.05; *p < 0.05; **p < 0.01; *** p < 0.001; ****p < 0.0001.
Construction and validation of APR_Score
The prognostic risk model was constructed based on ARLs by randomly dividing patients in TCGA-CESC into training and validation set at a ratio of 1:1, with the clinical information of the samples in each group shown in Table 1. In the training set, the LASSO COX model was applied to optimize the model, and 11 ARLs were identified based on the penalty parameter lambda and the model trajectory change curve (Supplementary Figure S1). Seven ARLs affecting prognosis were selected based on the multivariate COX model, namely AC092614.1, AL158071.3, AC016394.2, LINC02749, MIR100HG, AC079313.1, and ATP2A1-AS1. According to APR_Score=1.847*AC092614.1+(−0.88*AL158071.3)+(−0.675*AC016394.2)+(−1.099*LINC02749)+0.514*MIR100HG+0.766*AC079313.1+(-0.323* ATP2A1-AS1), the prognostic risk model was assessed. Patients in the training and validation sets were classified into high-APR_Score and low-APR_Score groups by the optimal p-value in the survminer. Based on the K-M curves, we noted that high APR_Score scores predicted poorer OS in the training set, validation set and total TCGA-CESC cohort (Figures 6A–C). In the training set, the AUCs were 0.73, 0.89, 0.91, 0.91 and 0.9 for 1, 2, 3, 4 and 5 year(s) of survival, respectively (Figure 6D). While in the validation set and the total TCGA-CESC cohort, APR_Score still showed higher AUC values in predicting 1-, 2-, 3-, 4-, and 5-year OS in CC patients (Figures 6E, F). These results suggested that APR_Score had a high accuracy in predicting OS risk in CC patients and might be a novel prognostic indicator.
FIGURE 6. Development and validation of APR_Score (A–C) K-M curves for the 7 ARLs in the training set, validation set and TCGA-CESC cohort (D–F) ROC curves of the training set, validation set and TCGA-CESC cohort.
Correlation between APR_Score and clinical characteristics
To determine the relationship between APR_Score and clinical characteristics, this study discussed the interaction between APR_Score and different clinical features. First, we found that the APR_Score was the highest in the poor prognosis AC3 subtype, with more death cases in the high APR_Score group (Figures 7A, B). Expression differences of AC092614.1, AL158071., AC016394.2, LINC02749, MIR100HG, AC079313.1, ATP2A1- AS1 in different types of patients are shown in Figure 7B. In addition, to compare the meaning of APR_Score in clinicopathological subgroups for CC prognosis, K-M survival analysis was performed on patients with high APR_Score and low APR_Score in clinicopathological subgroups. The results showed that patients with high APR_Score had markedly poorer OS than those with low APR_Score in Stage, Grade and Age subgroups (p < 0.0001, Figures 7C–E). Additionally, we compared the differences of clinical characteristics (N satge and M stage) between risk groups in TCGA-CESC cohort. The patients with metastasis (including lymph node metastasis and distal metastasis) in the high-APR_Score group were more than those in the low-APR_Score group (Supplementary Figure S2).
FIGURE 7. Correlation between APR_Score and clinical characteristics of CC (A) APR_Score in AC1-3 subtypes (B) Heatmap showing the expression of 7 OS-related ARLs in ARLs-related molecular subtypes and APR_Score groupings (C–E) K-M curves for patients with high APR_Score and low APR_Score in Stage, Grade, and Age subgroups.
Nomogram for predicting survival rate of CC
To further discuss the clinical value of APR_Score in CC patients, univariate and multivariate COX model analyses were performed, and we determined that APR_Score and Stage were clinically significant for CC prognosis (p < 0.05, Figures 8A, B). Due to a close correlation between APR_Score and Stage and prognosis, we created a nomogram based on APR_Score and Stage to assess OS of CC patients (Figure 8C). The calibration curves showed a great overlap between Nomogram predicted survival at 1-, 3- and 5-year OS with actual observations, indicating that nomogram was a reliable tool for predicting OS (Figure 8D). We also found that the nomogram and APR_Score showed higher accuracy in predicting 1-year, 3-year and 5-year OS when compared to Age, TNM Stage, Stage and Grade (Figures 8E–G).
FIGURE 8. Nomogram for predicting survival rate of CC (A, B) Forest plots for univariate and multivariate COX models for APR_Score (C) Nomogram of APR_Score and Stage (D) Calibration curves for 1-year, 3-years and 5-years survival (E–G) The ROC curves of a variety of clinical features for overall OS at 1-year, 3-year and 5-year.
TME activity assessment in APR_Score groups
The CIBERSORT and MCP-count algorithms were applied to assess the abundance of immune cell infiltration and immune score in TME of APR_Score groups. CIBERSORT results, as presented in Figure 9A, showed that T cells CD8, Tregs infiltration abundance was lower in high-APR_Score group. MCP-count results also showed lower immune score of CD8 T cells in high-APR_Score group (Figure 9B). Then, to discuss the potential connection between APR_Score and immunotherapy response, we assessed TIDE scores and Exclusion scores in the APR_Score groups. TIDE scores and Exclusion scores were slightly higher in the high APR_Score group than in the low APR_Score group, suggesting that patients with high APR_Score were more likely to experience immune escape and less responsive to immunotherapy (Figures 9C, D).
FIGURE 9. TME activity assessment in APR_Score groups (A) Results of immune cell scores estimated by CIBERSORT between APR_Score groups in TCGA-CESC cohort (B) Results of immune cell scores using MCP-count between APR_Score groups in TCGA-CESC cohort (C) Distribution of TIDE score between APR_Score groups in TCGA-CESC cohort (D) Exclusion scores in high-APR_Score group and low-APR_Score group. ns p > 0.05; *p < 0.05; **p < 0.01; *** p < 0.001; ****p < 0.0001.
Association between APR_Score and chemotherapy drug sensitivity
Until immunotherapy was proposed as an alternative treatment for CC, conventional resection and radiotherapy were the dominant treatments (Serkies and Jassem, 2018). In this study, to discuss the potential of APR_Score as therapeutic response marker for predicting chemotherapeutic agents, we evaluated the IC50 values of 20 chemotherapeutic agents in TCGA-CESC. Initially, correlations between APR_Score and drug IC50 values were computed, and highly correlated drugs were selected for comparison in the APR_Score groups (Supplementary Figure S3). We found positive responses to small molecular drugs including Rapamycin, KIN001-135, Roscovitine, Phenformin treatments in the low APR_Score group, while the high APR_Score group responded positively to Sunitinib, MG-132, Paclitaxel, AZ628, Sorafenib, Saracatinib, Dasatinib, CGP-60474, A-770041, WH-4-023, WZ-1 -84, CMK, Bortezomib, Lapatinib, Midostaurin, Embelin treatments (Figure 10). Overall, the APR_Score was correlated with the sensitivity to small molecular drugs.
FIGURE 10. Differences of IC50 values for 20 chemotherapy drugs between APR_Score groups in TCGA-CESC cohort. ns p > 0.05; *p < 0.05; **p < 0.01; *** p < 0.001; ****p < 0.0001.
Biological pathway characterization of APR_Score groups
Moreover, to discuss the biological pathway variation in APR_Score groups, we evaluated the markedly enriched pathways in distinct APR_Score groups using GSEA methods via h.all.v7.5.1.symbols.gmt signatures. We noted that TGF BETA SIGNALING, APICAL JUNCTION, HYPOXIA, EPITHELIAL MESENCHYMAL TRANSITION, APOPTOSIS, UV RESPONSE DN, KRAS SIGNALING UP, TNFA SIGNALING VIA NFKB, PROTEIN SECRETION, ANGIOGENESIS were markedly enriched in high-APR_Score group (p < 0.05, FDR<0.25, Figure 11).
Discussion
Anoikis is a protective mechanism via which the organism self-corrects abnormal disorders in the presence of disorders and damage, and is essential for the normal growth and development process of living organisms (Zhong and Rescorla, 2012). Advanced cancer patients are affected by tumor metastasis, which occurs in large tumors and could result in death (Steeg, 2016). The prerequisite for tumor cells to migrate would be their own resistance to Anoikis, otherwise the self-correcting capacity of the organism could prevent the formation of metastatic lesions (Simpson et al., 2008; Kim et al., 2012). In fact, the existing studies related to the principles of the anoikis mechanism in CC remained limitations.
In this study, we defined three potential subtypes based on anoikis-related lncRNAs in CC that exhibited remarkable prognostic variation as well as immunoreactivity. Notably, AC3 patients exhibited higher ARG scores and angiogenesis scores. Angiogenesis is essential and important in the development of malignant tumor (Viallard and Larrivee, 2017). In contrast, tumor cell anoikis resistance was sufficiently necessary for the formation of new blood vessels and for migration (Kim et al., 2012). The interaction between anoikis and tumor angiogenesis was reported. Gao and colleagues showed that osteosarcoma cells resist through anoikis by activating Src kinase, which in turn activates JNK/ERK/VEGF-A to promote the formation of tumor metastases (Gao et al., 2019). In the present study, we also found marked activation of angiogenesis and NOTCH signaling in the AC3 subtype (Figure 5E). The survival rate of AC3 patients was not satisfactory. This allowed us to speculate that AC3 patients might have a higher probability of postoperative recurrence and higher prognostic risk. However, distal metastases occurred in the majority of cases at the time of diagnosis in CC patients at present (Dereje et al., 2020; Sengayi-Muchengeti et al., 2020; Ryzhov et al., 2021). Accurate determination of disease staging would be crucial for patient treatment modalities and clinical outcomes. Therefore, the ARLs-related subtypes defined in this study might assist M Stage and more accurately calculate the risk of metastasis in CC patients.
Here, we have constructed prognostic risk stratification models for CC based on the expression of ARLs. We explored the potential utility of APR_Score in predicting patient survival and prognosis. In addition, this study also explored the potential function of APR_Score in guiding immunotherapy and chemotherapy. Immunotherapy offered promising opportunities for patients as an emerging option for the rehabilitation of advanced CC patients (Duska et al., 2020; Colombo et al., 2021). Pembrolizumab was approved by the Food and Drug Administration (FDA) as an immunotherapeutic agent for CC (De Felice et al., 2021), however, the current treatment with Pembrolizumab for CC remained limited due to the absence of efficacy biomarkers to determine treatment benefit (Marret et al., 2019). In this study, there were discrepancies in immune cell scores between high- and low-APR_Score patients and in the TIDE scores. We observed higher T cell CD8, T cell CD4 memory resting and TIDE scores in patients with high APR_Score. Tumor cells blocked T cell killing by expression of PD-1/PD-L1, causing immune escape of tumor cells (Pardoll, 2012; Chabanon et al., 2016). Another study concluded that CC cells had a greater tendency of immune escape and were accompanied by increased T cells (Mortezaee, 2020). Moreover, acquisition of anoikis-resistance enhances the abilities of invasiveness, escaping from immune surveillance and therapeutic agents in cancer cells. Fanfone and colleagues have demonstrated that mechanically stressed and anoikis-resistant cancer cells had increased cell motility and escape from killing by natural killer cells (Fanfone et al., 2022), suggesting that anoikis was strongly associated with immune escape of tumor cells. We constructed a prognostic model based on ARLs, and patients with high APR_Score had suboptimal immunotherapy benefit and were more suitable for taking conventional chemotherapy, suggesting that ARLs-based prognostic risk model was a promising predictor in immunotherapy and chemotherapy. Furthermore, the COX model demonstrated that APR_Score was an independent prognostic indicator, and that the nomogram developed in conjunction with APR_Score was a valid prognostic tool. Thus, this finding confirmed the positive utility of ARLs-related signatures in guiding the prognosis and treatment selection, which might provide a theoretical basis for the development of precise and personalized treatment guidelines for CC.
ARLs-related prognostic risk model consisted of AC092614.1, AL158071.3, AC016394.2, LINC02749, MIR100HG, AC079313.1, ATP2A1-AS1. We found that MIR100HG and ATP2A1-AS1 were correlated with multiple cancers and the remaining ARLs were the first identified tumor prognostic markers. Several studies suggested that aberrant expression of MIR100HG was associated with poor clinical outcome and pathological features, and that it was involved in multiple pathways related to tumorigenesis (Wu et al., 2022). aTP2A1-AS1 is a potential prognostic biomarker for CC (Feng et al., 2021). The mechanism of anoikis in CC was not completely elucidated, and the ARLs identified in this study were important for the elucidation of the molecular mechanisms of CC. However, limitations of this study should be equally noted. In other datasets, there are few cervical cancer-related data sets and they lack clinical information (survival time). Due to the difficulty in annotating lncRNAs in other databases, only the TCGA database was used and randomly divided into training set and validation set to construct and validate the risk model in this study. The accuracy of APR_Score should be validated in the future using clinical samples or sequencing data from multiple centers. Finally, the mechanism of 7-ARLs in CC still required further study, which was our follow-up research plan.
Conclusion
In this study, we defined three molecular subtypes of anoikis-related lncRNAs and generated ARLs-related signatures for assessing the prognosis of patients with CC. The molecular subtypes contributed to a better understanding of the mechanism of CC metastasis and the signatures held potential clinical value in predicting response to therapy.
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 author.
Author contributions
All authors contributed to this present work: HL and LX designed the study and acquired the data. HW drafted the manuscript, YL and JZ revised the manuscript. All authors read and approved the manuscript.
Funding
The subject was funded by the General Program of Chongqing Natural Science Foundation: Study on the Role of CircRNA-0008193 Targeting miR-1226-5p in Regulating PD-L1 in Immune Escape of Cervical Cancer(cstc2021jcyj-msxmX0219). We are grateful for the support of the Tianjin Municipal Education Commission Scientific Research Program Project (2018KJ036).
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/fphar.2023.1135626/full#supplementary-material
Abbreviations
ARLs, Anoikis-related lncRNAs; APR_Score, Anoikis prognostic risk model; ARGs, Anoikis-related genes; AUC, area under ROC curve; CC, Cervical cancer; FC, fold change; GSVA, Gene set variation analysis; GSEA, Gene set enrichment analysis; HR, hazard ratio; IC50, the half maximal inhibitory concentration; IDE, integrated development environment; K-M, Kaplan-Meier; LASSO, least absolute shrinkage and selection operator; MSigDB, Molecular Signatures Database; OS, overall survival; ROC, receiver operating characteristic analysis; ssGSEA, single-sample gene set enrichment analysis; TCGA, The Cancer Genome Atlas; TIDE, Tumor Immune Dysfunction and Exclusion; TME, Tumor microenvironment.
References
Adeshakin, F. O., Adeshakin, A. O., Afolabi, L. O., Yan, D., Zhang, G., and Wan, X. (2021). Mechanisms for modulating anoikis resistance in cancer and the relevance of metabolic reprogramming. Front. Oncol. 11, 626577. doi:10.3389/fonc.2021.626577
Barbie, D. A., Tamayo, P., Boehm, J. S., Kim, S. Y., Moody, S. E., Dunn, I. F., et al. (2009). Systematic RNA interference reveals that oncogenic KRAS-driven cancers require TBK1. Nature 462 (7269), 108–112. doi:10.1038/nature08460
Becht, E., Giraldo, N. A., Lacroix, L., Buttard, B., Elarouci, N., Petitprez, F., et al. (2016). Estimating the population abundance of tissue-infiltrating immune and stromal cell populations using gene expression. Genome Biol. 17 (1), 218. doi:10.1186/s13059-016-1070-5
Bhatla, N., Aoki, D., Sharma, D. N., and Sankaranarayanan, R. (2018). Cancer of the cervix uteri. Int. J. Gynaecol. Obstet. 143 (2), 22–36. doi:10.1002/ijgo.12611
Canfell, K., Kim, J. J., Brisson, M., Keane, A., Simms, K. T., Caruana, M., et al. (2020). Mortality impact of achieving WHO cervical cancer elimination targets: A comparative modelling analysis in 78 low-income and lower-middle-income countries. Lancet 395 (10224), 591–603. doi:10.1016/S0140-6736(20)30157-4
Chabanon, R. M., Pedrero, M., Lefebvre, C., Marabelle, A., Soria, J. C., and Postel-Vinay, S. (2016). Mutational landscape and sensitivity to immune checkpoint blockers. Clin. Cancer Res. 22 (17), 4309–4321. doi:10.1158/1078-0432.CCR-16-0903
Chen, B., Khodadoust, M. S., Liu, C. L., Newman, A. M., and Alizadeh, A. A. (2018). Profiling tumor infiltrating immune cells with CIBERSORT. Methods Mol. Biol. 1711, 243–259. doi:10.1007/978-1-4939-7493-1_12
Colombo, N., Dubot, C., Lorusso, D., Caceres, M. V., Hasegawa, K., Shapira-Frommer, R., et al. (2021). Pembrolizumab for persistent, recurrent, or metastatic cervical cancer. N. Engl. J. Med. 385 (20), 1856–1867. doi:10.1056/NEJMoa2112435
Danilova, L., Ho, W. J., Zhu, Q., Vithayathil, T., De Jesus-Acosta, A., Azad, N. S., et al. (2019). Programmed cell death ligand-1 (PD-L1) and CD8 expression profiling identify an immunologic subtype of pancreatic ductal adenocarcinomas with favorable survival. Cancer Immunol. Res. 7 (6), 886–895. doi:10.1158/2326-6066.CIR-18-0822
De Felice, F., Giudice, E., Bolomini, G., Distefano, M. G., Scambia, G., Fagotti, A., et al. (2021). Pembrolizumab for advanced cervical cancer: Safety and efficacy. Expert Rev. Anticancer Ther. 21 (2), 221–228. doi:10.1080/14737140.2021.1850279
Dereje, N., Gebremariam, A., Addissie, A., Worku, A., Assefa, M., Abraha, A., et al. (2020). Factors associated with advanced stage at diagnosis of cervical cancer in addis ababa, Ethiopia: A population-based study. BMJ Open 10 (10), e040645. doi:10.1136/bmjopen-2020-040645
Detmar, M. (2000). Tumor angiogenesis. J. Investig. Dermatol Symp. Proc. 5 (1), 20–23. doi:10.1046/j.1087-0024.2000.00003.x
Duska, L. R., Scalici, J. M., Temkin, S. M., Schwarz, J. K., Crane, E. K., Moxley, K. M., et al. (2020). Results of an early safety analysis of a study of the combination of pembrolizumab and pelvic chemoradiation in locally advanced cervical cancer. Cancer 126 (22), 4948–4956. doi:10.1002/cncr.33136
Fanfone, D., Wu, Z., Mammi, J., Berthenet, K., Neves, D., Weber, K., et al. (2022). Confined migration promotes cancer metastasis through resistance to anoikis and increased invasiveness. Elife 11, e73150. doi:10.7554/eLife.73150
Feng, Q., Wang, J., Cui, N., Liu, X., and Wang, H. (2021). Autophagy-related long non-coding RNA signature for potential prognostic biomarkers of patients with cervical cancer: A study based on public databases. Ann. Transl. Med. 9 (22), 1668. doi:10.21037/atm-21-5156
Gao, Z., Zhao, G. S., Lv, Y., Peng, D., Tang, X., Song, H., et al. (2019). Anoikis-resistant human osteosarcoma cells display significant angiogenesis by activating the Src kinase-mediated MAPK pathway. Oncol. Rep. 41 (1), 235–245. doi:10.3892/or.2018.6827
Geeleher, P., Cox, N., and Huang, R. S. (2014). pRRophetic: an R package for prediction of clinical chemotherapeutic response from tumor gene expression levels. PLoS One 9 (9), e107468. doi:10.1371/journal.pone.0107468
Hanzelmann, S., Castelo, R., and Guinney, J. (2013). Gsva: Gene set variation analysis for microarray and RNA-seq data. BMC Bioinforma. 14, 7. doi:10.1186/1471-2105-14-7
Harrell Jr, F. E., Harrell Jr, M. F. E., and Hmisc, D. (2017). Package ‘rms’[J]. Vanderbilt University 17 (4), 229.
Kim, Y. N., Koo, K. H., Sung, J. Y., Yun, U. J., and Kim, H. (2012). Anoikis resistance: An essential prerequisite for tumor metastasis. Int. J. Cell Biol. 2012, 306879. doi:10.1155/2012/306879
Marret, G., Borcoman, E., and Le Tourneau, C. (2019). Pembrolizumab for the treatment of cervical cancer. Expert Opin. Biol. Ther. 19 (9), 871–877. doi:10.1080/14712598.2019.1646721
Marth, C., Landoni, F., Mahner, S., McCormack, M., Gonzalez-Martin, A., Colombo, N., et al. (2017). Cervical cancer: ESMO clinical practice guidelines for diagnosis, treatment and follow-up. Ann. Oncol. 28 (4), iv72–iv83. doi:10.1093/annonc/mdx220
Mayakonda, A., Lin, D. C., Assenov, Y., Plass, C., and Koeffler, H. P. (2018). Maftools: Efficient and comprehensive analysis of somatic variants in cancer. Genome Res. 28 (11), 1747–1756. doi:10.1101/gr.239244.118
Meng, Y., Chu, T., Lin, S., Wu, P., Zhi, W., Peng, T., et al. (2021). Clinicopathological characteristics and prognosis of cervical cancer with different histological types: A population-based cohort study. Gynecol. Oncol. 163 (3), 545–551. doi:10.1016/j.ygyno.2021.10.007
Mortezaee, K. (2020). Immune escape: A critical hallmark in solid tumors. Life Sci. 258, 118110. doi:10.1016/j.lfs.2020.118110
Paoli, P., Giannoni, E., and Chiarugi, P. (2013). Anoikis molecular pathways and its role in cancer progression. Biochim. Biophys. Acta 1833 (12), 3481–3498. doi:10.1016/j.bbamcr.2013.06.026
Pardoll, D. M. (2012). The blockade of immune checkpoints in cancer immunotherapy. Nat. Rev. Cancer 12 (4), 252–264. doi:10.1038/nrc3239
Parmar, D., and Apte, M. (2021). Angiopoietin inhibitors: A review on targeting tumor angiogenesis. Eur. J. Pharmacol. 899, 174021. doi:10.1016/j.ejphar.2021.174021
Ryzhov, A., Corbex, M., Pineros, M., Barchuk, A., Andreasyan, D., Djanklich, S., et al. (2021). Comparison of breast cancer and cervical cancer stage distributions in ten newly independent states of the former soviet union: A population-based study. Lancet Oncol. 22 (3), 361–369. doi:10.1016/S1470-2045(20)30674-4
Sanchez-Vega, F., Mina, M., Armenia, J., Chatila, W. K., Luna, A., La, K. C., et al. (2018). Oncogenic signaling pathways in the cancer Genome Atlas. Cell 173 (2), 321–337. e10. doi:10.1016/j.cell.2018.03.035
Sengayi-Muchengeti, M., Joko-Fru, W. Y., Miranda-Filho, A., Egue, M., Akele-Akpo, M. T., N'da, G., et al. (2020). Cervical cancer survival in sub-saharan africa by age, stage at diagnosis and human development index: A population-based registry study. Int. J. Cancer 147 (11), 3037–3048. doi:10.1002/ijc.33120
Serkies, K., and Jassem, J. (2018). Systemic therapy for cervical carcinoma - current status. Chin. J. Cancer Res. 30 (2), 209–221. doi:10.21147/j.issn.1000-9604.2018.02.04
Shen, W., Song, Z., Zhong, X., Huang, M., Shen, D., Gao, P., et al. (2022). Sangerbox: A comprehensive, interaction-friendly clinical bioinformatics analysis platform. iMeta 1 (3), e36. doi:10.1002/imt2.36
Simon, N., Friedman, J., Hastie, T., and Tibshirani, R. (2011). Regularization paths for cox's proportional hazards model via coordinate descent. J. Stat. Softw. 39 (5), 1–13. doi:10.18637/jss.v039.i05
Simpson, C. D., Anyiwe, K., and Schimmer, A. D. (2008). Anoikis resistance and tumor metastasis. Cancer Lett. 272 (2), 177–185. doi:10.1016/j.canlet.2008.05.029
Steeg, P. S. (2016). Targeting metastasis. Nat. Rev. Cancer 16 (4), 201–218. doi:10.1038/nrc.2016.25
Sung, H., Ferlay, J., Siegel, R. L., Laversanne, M., Soerjomataram, I., Jemal, A., et al. (2021). Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J. Clin. 71 (3), 209–249. doi:10.3322/caac.21660
Viallard, C., and Larrivee, B. (2017). Tumor angiogenesis and vascular normalization: Alternative therapeutic targets. Angiogenesis 20 (4), 409–426. doi:10.1007/s10456-017-9562-9
Wilkerson, M. D., and Hayes, D. N. (2010). ConsensusClusterPlus: A class discovery tool with confidence assessments and item tracking. Bioinformatics 26 (12), 1572–1573. doi:10.1093/bioinformatics/btq170
Wu, Y., Wang, Z., Yu, S., Liu, D., and Sun, L. (2022). LncmiRHG-MIR100HG: A new budding star in cancer. Front. Oncol. 12, 997532. doi:10.3389/fonc.2022.997532
Keywords: TCGA, cervical cancer, anoikis, LncRNAs, prognosis, immunotherapy, small molecular drug
Citation: Liang H, Xiang L, Wu H, Liu Y, Tian W and Zeng J (2023) Anoikis-related long non-coding RNA signatures to predict prognosis and small molecular drug response in cervical cancer. Front. Pharmacol. 14:1135626. doi: 10.3389/fphar.2023.1135626
Received: 01 January 2023; Accepted: 07 March 2023;
Published: 20 March 2023.
Edited by:
Zhi-Qian Zhang, Southern University of Science and Technology, ChinaReviewed by:
Zhenyao Chen, Fudan University, ChinaYan Wang, Beijing Obstetrics and Gynecology Hospital, Capital Medical University, China
Zigang Zhao, Hainan Hospital of PLA General Hospital, China
Copyright © 2023 Liang, Xiang, Wu, Liu Tian and Zeng. 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: Jianhua Zeng, 300541@hospital.cqmu.edu.cn