Research Paper Volume 16, Issue 8 pp 7405—7425
Identification and validation of an anoikis-related genes signature for prognostic implication in papillary thyroid cancer
- 1 Postgraduate Training Base at Shanghai Gongli Hospital, Ningxia Medical University, Shanghai 200135, China
- 2 School of Medicine, Shanghai University, Shanghai 200444, China
- 3 Department of Endocrinology, Gongli Hospital of Shanghai Pudong New Area, Shanghai 200135, China
- 4 Department of Otolaryngology Head and Neck Surgery, Gongli Hospital of Shanghai Pudong New Area, Shanghai 200135, China
Received: November 22, 2023 Accepted: March 3, 2024 Published: April 24, 2024
https://doi.org/10.18632/aging.205766How to Cite
Copyright: © 2024 Zhao et al. This is an open access article distributed under the terms of the Creative Commons Attribution License (CC BY 4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Abstract
Thyroid cancer, notably papillary thyroid cancer (PTC), is a global health concern with increasing incidence. Anoikis, a regulator of programmed cell death, is pivotal in normal physiology and, when dysregulated, can drive cancer progression and metastasis. This study explored the impact of anoikis on PTC prognosis. Analyzing data from GEO, TCGA, and GeneCards, we identified a prognostic signature consisting of six anoikis-related genes (ARGs): EZH2, PRKCQ, CD36, INHBB, TDGF1, and MMP9. This signature independently predicted patient outcomes, with high-risk scores associated with worse prognoses. A robust predictive ability was confirmed via ROC analysis, and a nomogram achieved a C-index of 0.712. Differences in immune infiltration levels were observed between high- and low-risk groups. Importantly, the high-risk group displayed reduced drug sensitivity and poor responses to immunotherapy. This research provides insights into anoikis in PTC, offering a novel ARG signature for predicting patient prognosis and guiding personalized treatment strategies.
Introduction
Thyroid cancer (THCA) is the most prevalent endocrine malignancy worldwide, with its incidence and morbidity steadily increasing [1]. Papillary thyroid cancer (PTC) is the predominant subtype of THCA, accounting for approximately 80% of all THCA cases [2]. The prognosis for PTC is typically favorable, with a 5-year survival rate exceeding 95% and a 10-year survival rate over 90% when treated with radioiodine ablation and/or revision surgery [3]. But it is important to note that approximately 20% of patients experience a reduced survival rate due to factors such as recurrence, metastasis, and other complications [4, 5]. Therefore, exploring the key molecules that promote the progression of PTC is necessary to improve the prognosis and quality of life of PTC patients [6].
Anoikis is a term used to describe a form of programmed cell death triggered by the detachment of cells from the extracellular matrix (ECM) or neighboring cells [7]. It is a crucial mechanism that regulates tissue homeostasis, development, and the prevention of cancer cell metastasis. When cells lose contact with their surrounding ECM, they undergo a series of biochemical and morphological changes that ultimately lead to their death [7]. During anoikis, various signaling pathways are activated, including integrin-mediated signaling, growth factor receptor signaling, and apoptotic pathways. These pathways converge to induce apoptosis, characterized by DNA fragmentation, mitochondrial dysfunction, cytoskeletal rearrangement, and caspase activation [7, 8]. Anoikis plays a significant role in normal physiological processes such as tissue remodeling and the maintenance of epithelial cell layers [9, 10]. Dysregulation of anoikis can have detrimental effects, contributing to cancer progression and metastasis [8]. Cancer cells often acquire resistance to anoikis, allowing them to survive and invade distant tissues [11].
Recently, there is a growing number of research focusing on the significance of anoikis resistance in tumorigenesis. The endocytic degradation of epidermal growth factor receptor (EGFR) can induce cancer cells to detach from the extracellular matrix, ultimately triggering apoptosis, known as anoikis. Targeting EGFR could be effective for anoikis suppression [12]. Upon stimulation by specific extracellular matrix components, such as cancer-related fibroblasts (CAFs), tumors acquire anoikis resistance and develop invasive and metastatic capabilities [13]. Emerging studies have demonstrated that the anoikis-related signature holds the potential to predict treatment response and prognosis in a wide range of cancer types. Until now, the establishment of clinical prognostic models based on anoikis-related genes has explored in head and neck squamous cell carcinoma, clear cell renal cell carcinoma, hepatocellular carcinoma, and others [14–16], but remains unexplored in PTC. Gaining insight into the mechanisms of action of anoikis-related genes may enhance our understanding of the development process of PTC.
In this study, we conducted an analysis of differentially expressed genes (DEGs) between PTC and adjacent normal tissues using publicly available databases. By integrating this data with clinical information, we developed a six-gene signature related to anoikis to predict the prognosis of PTC patients and elucidate immune cell infiltration patterns. This approach aims to enhance treatment options and improve patient outcomes in PTC.
Materials and Methods
Data acquisition
The messenger RNA (mRNA) expression data of PTC and normal tissues were downloaded from the GEO (including the GSE29265, GSE33630, and GSE60542 datasets; https://www.ncbi.nlm.nih.gov/geo). Additionally, the mRNA sequencing data (FPKM value) and corresponding clinical-pathologic information of PTC patients were obtained from the TCGA (https://cancergenome.nih.gov/). A total of 498 PTC patients were obtained and randomly assigned to either the training cohort (n = 349) or the test cohort (n = 149).
Identification of anoikis-related DEGs and functional enrichment analysis
Anoikis-related genes (ARGs) were extracted from GeneCards [17]. In the GeneCards database, each gene has a corresponding relevance score value, which is used to evaluate the correlation between the gene and the elements (chemical substances and diseases). The higher the score, the stronger the statistical correlation between genes and related elements [17]. Using a relevance score >0.4, 551 ARGs were selected. The “limma” R package was employed to screen DEGs (|log2(fold change)| > 1 and adjusted P-value < 0.05) between tumor tissue and normal tissues in GEO datasets. Then, the anoikis-related DEGs (ADGs) were obtained by overlapping the intersection of DEGs and ARGs using Venny (version 2.1) [18]. The ADGs were submitted to Metascape Online (https://metascape.org/), which incorporates a core set of default ontologies such as Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways to conduct functional analysis and construct a PPI network [19].
Establishment of the prognostic genes signature of the ADGs
To establish a prognostic model for ADGs in PTC patients from TCGA database, a univariate Cox analysis of progression-free interval (PFI) was first conducted using the R packages “survival” and “survminer” to identify the survival-related genes among the ADGs with a prognosis value (P < 0.05). These survival-related genes were then included in the subsequent Least Absolute Shrinkage and Selection Operator (LASSO) Cox regression analysis by using “glmnet” R package. The 10-fold cross validation and 1,000 iterations were performed to minimize the potential risk of overfitting and select optimal prognostic genes. The prognostic genes signature was calculated based on the normalized gene expression level and corresponding regression coefficients. The risk score for prognosis was determined by using a linear combination of the regression coefficient in the LASSO regression and the gene’s expression level.
The assessment of the prognosis model
PTC patients in both the training cohort and test cohort were respectively divided into high- and low-risk groups by the risk scores. The Kaplan-Meier (K-M) survival curve was performed to assess prognostic significance. The receiver operating characteristic curve (ROC curve) was created by using the R package “timeROC.” The area under the curve (AUC) in 1-, 3-, and 5-year progress-free survival was used to evaluate the performance of the prognostic model. The GEPIA database (http://gepia.cancer-pku.cn/) was used to perform survival analysis for the prognostic genes [20].
Nomogram construction and evaluation
Univariate and multivariate Cox analyses of patient risk score of the prognostic signature and clinical factors, including gender, pathological stage, and focal type, were performed to determine the significance of each factor in predicting PFI in patients with PTC. Based on the results of multivariate Cox regression, the R package “rms” was used to construct a nomogram to provide the 1-, 3- and 5-year survival probabilities. The calibration curve evaluated the predictive ability of the nomogram.
Functional analysis related to risk score
The “limma” R package was employed to screen DEGs (|log2(fold change)| > 1 and adjusted P-value < 0.05) between high- and low-risk groups in TCGA database. Heatmap was performed using the “pheatmap” and “ggplot2” R packages. GO and KEGG enrichment analyses were utilized to explore the biological functions of DEGs.
Evaluation of tumor immune microenvironment
The ESTIMATE algorithm was used to calculate the stromal score, immune score, and tumor purity between high-risk and low-risk groups [21]. The CIBERSORT algorithm was used to estimate the abundance of 22 immune cell types. The infiltration levels of 22 immune cell types between the high-risk and low-risk subgroups were compared [22].
Drug sensitivity and immunotherapeutic response analyses
The “pRRophetic” R package was employed to calculate the half-maximal inhibitory concentration (IC50) values of chemotherapeutic and targeted drugs for each PTC patients in TCGA cohort [23]. The tumor immune dysfunction and exclusion (TIDE; http://tide.dfci.harvard.edu/) was used to calculate the TIDE score of each patient according to myeloid-derived suppressor cell (MDSC), macrophage M2, T cell Dysfunction and Exclusion [24].
Statistical analysis
R software (version 4.3.0) was used to perform data analysis. The statistical value P < 0.05 indicates that the difference is statistically significant.
Data accessibility
The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request.
Results
Figure 1 illustrates the cohort design and analytical concepts for the entire study.
Figure 1. Flow chart of the current study.
ADGs in PTC patients and functional enrichment analysis
The data of PTC was obtained from three GEO cohorts and TCGA database, and the details are shown in Table 1. By comparing tumor tissues and adjacent non-neoplastic tissues, 561 DEGs were obtained from GEO: GSE29265, with 288 genes showing upregulation and 273 showing downregulation in tumor tissues (Figure 2A). 806 DEGs were obtained from GEO: GSE33630, with 459 genes showing upregulation and 347 showing downregulation in tumor tissues (Figure 2B). 745 DEGs were obtained from GEO: GSE60542, with 381 genes showing upregulation and 364 showing downregulation in tumor tissues (Figure 2C). The DEGs obtained from the GEO datasets were intersected with the 551 anoikis-related genes to obtain ADGs. 64 ADGs overlapped between the four datasets (Figure 2D). The modulation-specific biological processes and pathways features of the 64 ADGs were then analyzed. Results were visualized by bar graphs and PPI networks. The GO and KEGG analyses showed that enrichments were mainly focused on positive regulation of cell migration, pathways in cancer, response to wounding, proteoglycans in cancer, and positive regulation of programmed cell death et al. (Figure 2E, 2F).
Table 1. Overview of each dataset associated with PTC from GEO and TCGA.
| Platform | Cases of normal | Cases of tumor | Scanned items | Clinical files | |
| GSE29265 | GPL570 | 20 | 20 | mRNA | No |
| GSE33630 | GPL570 | 45 | 49 | mRNA | No |
| GSE60542 | GPL570 | 30 | 33 | mRNA | No |
| TCGA | NA | − | 498 | mRNA | Yes |
Figure 2. Overview of the differentially expressed anoikis-related genes in PTC. (A–C) Volcano plots of differentially expressed gens (DEGs) between PTC and normal tissues in GSE29265, GSE33630, and GSE60542. (D) Venn diagram showing the dysregulated anoikis-related genes common to the four datasets. (E) Bar graph showing the GO and KEGG analysis. (F) PPI network showing the distribution and relationship of the different enriched functions.
Construction of the ADGs-related prognostic model
Gene expression data (FPKM) with corresponding patient information were obtained from TCGA database of 498 PTC patients. The clinical characteristics of these two cohorts are summarized in Table 2, and there are no significant clinical differences between these two cohorts. Univariate Cox regression analysis in training cohort showed that six ADGs were significantly associated with tumor progression (Supplementary Table 1). Both six genes were subsequently identified as prognostic signature in a LASSO Cox regression analysis (Figure 3A, 3B). Thus, the anoikis-related genes signature was established based on the coefficient in the LASSO regression analysis. A risk score for each patient was calculated as follows: 1.4829 × (expression of EZH2) + (−0.2218) × (expression of PRKCQ) + (−0.0188) × (expression of CD36) + (−0.2417) × (expression of INHBB) + (−1.1662) × (expression of TDGF1) + (0.1198) × (expression of MMP9).
Table 2. The clinical characteristics of 498 PTC patients in TCGA.
| Categories | Total (n = 498) | Training (n = 349) | Test (n = 149) | P-value |
| Age at diagnosis | ||||
| <55 | 333 (66.9) | 239 (68.5) | 94 (63.1) | 0.242 |
| ≥55 | 165 (33.1) | 110 (31.5) | 55 (36.9) | |
| Gender | ||||
| Male | 134 (26.6) | 86 (24.6) | 48 (32.2) | 0.081 |
| Female | 364 (73.1) | 263 (75.4) | 101 (67.8) | |
| T stage | ||||
| T1 | 143 (28.7) | 97 (27.8) | 46 (30.9) | 0.174 |
| T2 | 162 (32.5) | 118 (33.8) | 44 (29.5) | |
| T3 | 168 (33.7) | 112 (32.1) | 56 (37.6) | |
| T4 | 23 (4.6) | 20 (5.7) | 3 (2.0) | |
| Tx | 2 (0.4) | 2 (0.6) | 0 (0.0) | |
| Pathological stage | ||||
| Stage I | 281 (56.4) | 208 (59.6) | 73 (49.0) | 0.051 |
| Stage II | 51 (10.2) | 34 (9.7) | 17 (11.4) | |
| Stage III | 110 (22.1) | 66 (18.9) | 44 (29.5) | |
| Stage IV | 54 (10.8) | 39 (11.2) | 15 (10.1) | |
| Unknown | 2 (0.4) | 2 (0.6) | 0 (0.0) | |
| Status | ||||
| Progress free | 448 (90.0) | 314 (90.0) | 134 (89.9) | 0.990 |
| Progress | 50 (10.0) | 35 (10.0) | 15 (10.1) | |
Figure 3. Prognostic analysis of the six anoikis-related genes signature in the TCGA cohort. (A, B) LASSO regression was performed with the minimum criteria. (C–E) Kaplan-Meier survival analyses comparing the progress-free survival of patients in the high- and low-risk groups were conducted in both the TCGA training cohort (P = 0.0022), test cohort (P = 0.043), and overall cohort (P = 0.00049). (F–H) The ROC curves validated the prognostic performance of the risk score in both the TCGA training, test, and overall cohort.
Estimations of the prognosis signature in the TCGA database
We employed K-M and ROC curves to evaluate the prognostic significance of the anoikis-related prognostic signature across training, test, and overall cohorts. The K-M analysis showed that patients with high-risk scores had significantly worse progression-free survival (FPS) rates compared to patients with low-risk score in both the training (P = 0.0022), test (P = 0.043), and overall (P = 0.00049) cohorts (Figure 3C–3E). The AUC values for 1-, 3-, and 5-year FPS rates obtained from the prognostic signature in the training cohort were 0.805, 0.705, and 0.703, respectively (Figure 3F), while those in the test cohort were 0.632, 0.649, and 0.813, respectively (Figure 3G), and in the overall cohort were 0.703, 0.683, and 0.69, respectively (Figure 3H). The distribution of the risk score is shown in Figure 4. The progression rate of PTC patients increases proportionally with their risk scores (Figure 4A, 4D). Furthermore, this pattern remains consistent across both the test (Figure 4B, 4E) and overall (Figure 4C, 4F) cohorts. Significantly, patients with high-risk scores exhibited elevated expression levels of EZH2 and MMP9, while PRKCQ, CD36, INHBB, and TDGF1 demonstrated decreased expression levels across all three cohorts (Figure 4G–4I). Utilizing the GEPIA database, we conducted survival analysis on the six genes featured in the prognostic signature. Our findings revealed that PRKCQ, CD36, INHBB, and TDGF1 were correlated with better prognosis and disease-free survival (DFS) in THCA patients, whereas EZH2 and MMP9 were linked to an unfavorable prognosis and reduced DFS among individuals with THCA (Figure 5).
Figure 4. Verification of the anoikis-related genes signature was conducted in the training, test, and overall cohorts. (A–F) Dot plots illustrating the survival and risk score for the training, test, and overall cohorts. (G–I) The heatmaps of six anoikis-related genes in the training, test, and overall cohorts.
Figure 5. Kaplan-Meier curves of genes associated with the six gene prognostic risk signature.
Establishment and evaluation of the nomogram
The univariate and multivariate Cox regression analyses showed that the risk score and T stage were the independent significant prognostic factors in predicting PFS in patients with PTC (Supplementary Table 2 and Figure 6A). We subsequently constructed a predictive nomogram with these factors (Figure 6B). Besides, the C-index calculated by R for the nomogram was 0.712, suggesting that the nomogram had a superior predictive performance. The calibration curves were used to evaluate the accuracy of the nomogram, in which a standard curve represented the best prediction. The predicted outcomes of 1-, 3- and 5-year progress-free survival rates showed excellent consistency (Figure 6C–6E).
Figure 6. Forest plot of the multivariate Cox regression analysis and the construction and evaluation of the nomogram. (A) Forest plot showing the risk score and T stage were the significant prognostic factors in predicting progress-free survival in patients with PTC. (B) Nomogram based on the risk score of the model and clinical information of PTC patients in the TCGA cohort. (C–E) Calibration curves of the nomogram for the probability of 1-, 3- and 5-years.
Functional enrichment analysis of the DEGs between the high- and low-risk patients
567 DEGs were obtained between high- and low-risk groups (Supplementary Figure 1). In order to gain a deeper comprehension of the potential biological mechanisms underlying the DEGs associated with high- and low-risk patients, we performed an analysis utilizing the GO and KEGG databases. The investigation revealed that the biological processes primarily encompassed the production of molecular mediator of immune response, immunoglobulin production, and immune response-regulating and activating cell surface receptor signaling pathway. Additionally, the molecular functions were found to be associated with antigen binding, receptor ligand activity, cytokine activity, and cytokine receptor binding. The cellular components analyzed in this study encompassed the external side of plasma membrane, plasma membrane signaling receptor complex, T cell receptor complex, and immunoglobulin complex. (Figure 7A, Supplementary Table 3). Through KEGG analysis, it was determined that these DEGs were primarily associated with cytokine-cytokine receptor interaction, viral protein interaction with cytokine and cytokine receptor, chemokine signaling pathway, and hematopoietic cell lineage (Figure 7B, Supplementary Table 4).
Figure 7. Enrichment analysis of differentially expressed genes between patients in high- and low-risk groups in the TCGA cohort. (A, B) GO (A) and KEGG (B) analysis results of differentially expressed genes.
Assessment of the immune microenvironment in PTC
The ESTIMATE algorithm was used to calculate the stromal score, tumor purity, immune score, and ESTIMATE score. The high-risk group had a higher stromal score, immune score, and ESTIMATE score, and a lower tumor purity than the low-risk group (Figure 8A). To explore the correlation between anoikis-related genes signature and immune landscape in PTC, we employed the CIBERSORT algorithm to assess the relative proportions of 22 immune cell types in a cohort of 498 PTC patients from the TCGA database (Figure 8B). Immune cells, including B cells naïve, plasma cells, T cells CD8+, T cells CD4+ memory activated, T cells follicular helper, T cells regulatory (Tregs), T cells gamma delta, NK cells resting, monocytes, macrophages M0, macrophages M1, macrophages M2, mast cells resting, eosinophils, and neutrophils were statistically different between the high- and low-risk groups (Figure 8C).
Figure 8. Correlation between immune cell infiltration and different risk scores. (A) Differences of TME between different risk groups by ESTIMATE algorithm. (B) Immune cell proportions for each tumor patient. (C) Analyzing the immune cell infiltration levels of PTC samples between different risk groups by CIBERSORT algorithm. *P < 0.05; **P < 0.01; ***P < 0.001.
Prediction of drug sensitivity and immunotherapy efficacy
Using the R package “pRRophetic”, we investigated the IC50 values of six common chemotherapeutic and targeted drugs (Cisplatin, Doxorubicin, Paclitaxel, Sorafenib, Axitinib, and Sunitinib) in both the high- and low-risk groups of PTC patients. PTC patients in the low-risk group demonstrated greater sensitivity to Sorafenib and Axitinib. Conversely, PTC patients in the high-risk group exhibited increased sensitivity to Cisplatin, Doxorubicin, Paclitaxel, and Sunitinib (Figure 9A–9F). To evaluate the risk model’s value in immunotherapy, we compared the TIDE and other immune-related scores between high- and low-risk groups. The high-risk group showed higher TIDE score and Dysfunction score, and lower Exclusion score (Figure 9G–9I), indicating a limited immunotherapy benefit for PTC patients in the high-risk group.
Figure 9. Evaluation the value of the anoikis-related prognostic model in drug sensitivity and immunotherapy. (A–F) Correlation between risk score and drug sensitivity. (G–I) Differences in TIDE, Exclusion, and Dysfunction between low- and high-risk groups. *P < 0.05; **P < 0.01; ***P < 0.001; ****P < 0.0001.
Discussion
An important aspect of cancer progression is the ability of tumor cells to detach from their primary site and invade surrounding tissues or spread to distant organs [25, 26]. In PTC, Lymph node metastasis has been associated with a lower overall survival rate [27], furthermore, distant metastasis serves as the primary cause of death in PTC patients [28]. Therefore, construction of a metastasis-related signature to predict the prognosis may provide important insights into disease management and early intervention [29]. Anoikis, a form of programmed cell death, is initiated when cells lose adhesion to the extracellular matrix [30]. A pivotal process in metastasis involves the capability of cancer cells to survive without adhesion, evading anchorage-dependent cell death [31].
Under normal physiological conditions, anoikis serves as a safeguard, thwarting the survival of detached cells and the formation of secondary tumors [32]. Nonetheless, cancer cells can evolve strategies to elude anoikis, thus bolstering their metastatic potential [33]. Anoikis resistance in PTC arises from intricate pathways, underscoring the significance of targeting anoikis-related genes to counteract PTC progression and metastasis. Multifaceted analyses reveal the intricate interplay among diverse factors influencing anoikis resistance in PTC [34–38]. Understanding the molecular mechanisms underlying anoikis resistance in PTC is crucial for developing targeted therapies to prevent metastasis and improve patient outcomes [39].
This study was the first attempt at establishing an anoikis-related prognostic genes signature in PTC. In our study, we analyzed DEGs in PTC patients from three GEO datasets. 64 ADGs were subsequently selected by intersecting DEGs and anoikis-related genes. Functional analysis showed positive regulation of cell migration, pathways in cancer, response to wounding, and et al were major enriched, which highlighted the role of ADGs in tumor progression. By using univariate and LASSO Cox regression analyses, we proposed 6-gene signature, including EZH2, PRKCQ, CD36, INHBB, TDGF1, and MMP9. Their role in the process of anoikis resistance has been demonstrated.
EZH2 could reduce ITGα2 transcription, leading to decreased focal adhesions between the extracellular matrix and the cytoskeleton, enhancing cell mobility and increase anoikis resistance [40]. PRKCQ could promote anoikis resistance via kinase-activity-dependent stimulation of Erk/MAPK in breast cancer [41]. Blocking CD36 function or reducing its transcription could limit demethyl fruticulin A intake and integrin sequestration, restoring cell division and preventing anoikis [42]. INHBB could suppress anoikis resistance and migration of nasopharyngeal carcinoma cells by the TGF-β signaling pathway [43]. TDGF1 could enhance the anoikis resistance and the invasion ability of breast cancer cells [44]. MMP9 exerts its effects on the epithelium by cleaving one or more components of cell-cell junctions and triggering anoikis [45].
The risk score was subsequently calculated, and the ROC curves and AUC values confirmed the predictive power of our signature. In addition, our K-M analyses showed that the expression of each of these six genes was associated with the prognosis of PTC. Furthermore, we developed a nomogram incorporating the risk score and T stage. The nomogram demonstrated excellent predictive ability, as indicated by a high C-index of 0.712. This may have the potential to introduce novel insights into clinical decision-making processes. Additionally, we divided patients into high- and low-risk groups by the median risk score, and performed enrichment analysis, immune infiltration analysis, and drug sensitivity evaluation between two groups of PTC patients in TCGA.
The GO enrichment analysis showed that the biological processes involved the production of molecular mediator of immune response, immunoglobulin production, and immune response signaling pathway, etc., which indicated that our prognostic model may be closely related to the biological immune process. Thus, these results motivated us to study the relationship between our prognostic model and tumor immune microenvironment (TME). Additionally, the KEGG enrichment analysis results showed the significance of cytokine-cytokine receptor interaction, which could regulate the immune mechanism of PTC patients and then affect the progression, metastasis and migration of the disease [46]. Besides, the KEGG enrichment analysis results highlighted the significance of pathways such as chemokine signaling pathway. The chemokine could be involved in recruitment and maintenance of immune cells within the THCA microenvironment [47], and this process may involve thyroid cancer cell growth, aggressiveness and metastasis [48].
The immune microenvironment refers to the complex network of immune cells, cytokines, chemokines, and other factors present within the tumor microenvironment [49]. A balanced TME can recognize and eliminate cancer cells, a process known as immunosurveillance [50]. However, tumors often develop various strategies to evade immune recognition and destruction, leading to immune evasion and tumor growth [51]. Therefore, we explored the TME between our two groups. The high-risk group exhibited higher stromal score, immune scores, and ESTIMATE score, and lower tumor purity compared to the low-risk group. Studies have suggested that an increase in stromal cells is associated with a poorer prognosis in THCA, a finding consistent with our results [52]. Furthermore, the augmentation of stromal cells can enhance the expression of diverse immune checkpoints, facilitating tumor cell immune evasion [52].
The differences in immune scores suggested different infiltration of immune cells between the two groups. We found that the infiltration of immune cells, including B cells naïve, plasma cells, T cells CD4+ memory activated, T cells follicular helper, T cells regulatory (Tregs), T cells gamma delta, macrophages M0, and macrophages M1, were higher in high-risk group, while T cells CD8+, NK cells resting, monocytes, macrophages M2, mast cells resting, eosinophils, and neutrophils were lower in high-risk group. These findings showed that aberrant immune cell infiltration may facilitate the progression of PTC, and may provide guidance for us to further analyze the correlation between TME and anoikis-related genes in PTC.
In this study, we analyzed the drug sensitivity of chemotherapy drugs and targeted drugs commonly used for PTC patients. The high-risk group showed increased sensitivity to Sorafenib and Axitinib, but reduced sensitivity to Cisplatin, Doxorubicin, Paclitaxel, and Sunitinib. Furthermore, tumor immunotherapy has gained significant popularity, but the immunotherapy of THCA is still in the initial stage of exploration [53]. TIDE scores serve as a predictive tool for patient response to immunotherapy, as they predict the tumor’s potential capacity for immune evasion [24]. In our study, the high-risk group exhibited higher TIDE scores and Dysfunction score, and lower Exclusion score. These findings indicate that PTC patients in the high-risk group may experience high immune escape potential due to immune cell dysfunction, and the immunotherapy efficacy may be limited.
However, our study has certain limitations. Firstly, our research primarily relies on data from public databases, necessitating the need for further cohort studies with larger sample sizes to validate our findings. Secondly, the underlying mechanisms of the association between ARGs and PTC warrant further experimental verification.
Conclusion
In conclusion, we developed a six anoikis-related genes signature capable of predicting progression in PTC patients. Further investigations into the molecular mechanisms underlying these features will enhance our understanding of their impact on PTC progression, potentially providing valuable insights for precision medicine approaches.
Author Contributions
R-Y Z, Z-H W, and X-P C conceived and designed the study. R-Y Z and Z-H W conducted the database mining. P-P Q and L-Y Y analyzed data and compiled charts. R-Y Z and Y-Y L wrote the paper. X-P C, S-X H, and Y Z reviewed and edited the manuscript.
Conflicts of Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
Funding
This study was supported by the Research Grant for Health Science and Technology of Pudong Health Bureau of Shanghai (Grant No. PW2019D-4), the Subject Construction Project of Pudong Health Committee of Shanghai (Grant No. PWZy2020-06), the Pudong New Area Health System Leading Talent Training Program (Grant No. PWRI2021-07), and the Pudong New Area Clinical Characteristic Discipline (Grant No. PWYts2021-15).
References
- 1. Kitahara CM, Sosa JA. The changing incidence of thyroid cancer. Nat Rev Endocrinol. 2016; 12:646–53. https://doi.org/10.1038/nrendo.2016.110 [PubMed]
- 2. Lee YC, Na SY, Park GC, Han JH, Kim SW, Eun YG. Occult lymph node metastasis and risk of regional recurrence in papillary thyroid cancer after bilateral prophylactic central neck dissection: A multi-institutional study. Surgery. 2017; 161:465–71. https://doi.org/10.1016/j.surg.2016.07.031 [PubMed]
- 3. Margolick J, Chen W, Wiseman SM. Systematic Review and Meta-Analysis of Unplanned Reoperations, Emergency Department Visits and Hospital Readmission After Thyroidectomy. Thyroid. 2018; 28:624–38. https://doi.org/10.1089/thy.2017.0543 [PubMed]
- 4. Carling T, Udelsman R. Thyroid cancer. Annu Rev Med. 2014; 65:125–37. https://doi.org/10.1146/annurev-med-061512-105739 [PubMed]
- 5. Mazzaferri EL, Massoll N. Management of papillary and follicular (differentiated) thyroid cancer: new paradigms using recombinant human thyrotropin. Endocr Relat Cancer. 2002; 9:227–47. https://doi.org/10.1677/erc.0.0090227 [PubMed]
- 6. Li MY, Tang XH, Fu Y, Wang TJ, Zhu JM. Regulatory Mechanisms and Clinical Applications of the Long Non-coding RNA PVT1 in Cancer Treatment. Front Oncol. 2019; 9:787. https://doi.org/10.3389/fonc.2019.00787 [PubMed]
- 7. Frisch SM, Francis H. Disruption of epithelial cell-matrix interactions induces apoptosis. J Cell Biol. 1994; 124:619–26. https://doi.org/10.1083/jcb.124.4.619 [PubMed]
- 8. Gilmore AP. Anoikis. Cell Death Differ. 2005 (Suppl 2); 12:1473–7. https://doi.org/10.1038/sj.cdd.4401723 [PubMed]
- 9. Santagostino SF, Assenmacher CA, Tarrant JC, Adedeji AO, Radaelli E. Mechanisms of Regulated Cell Death: Current Perspectives. Vet Pathol. 2021; 58:596–623. https://doi.org/10.1177/03009858211005537 [PubMed]
- 10. Hindupur SK, Balaji SA, Saxena M, Pandey S, Sravan GS, Heda N, Kumar MV, Mukherjee G, Dey D, Rangarajan A. Identification of a novel AMPK-PEA15 axis in the anoikis-resistant growth of mammary cells. Breast Cancer Res. 2014; 16:420. https://doi.org/10.1186/s13058-014-0420-z [PubMed]
- 11. Paoli P, Giannoni E, Chiarugi P. Anoikis molecular pathways and its role in cancer progression. Biochim Biophys Acta. 2013; 1833:3481–98. https://doi.org/10.1016/j.bbamcr.2013.06.026 [PubMed]
- 12. Iradyan M, Iradyan N, Hulin P, Hambardzumyan A, Gyulkhandanyan A, Alves de Sousa R, Hessani A, Roussakis C, Bollot G, Bauvais C, Sakanyan V. Targeting Degradation of EGFR through the Allosteric Site Leads to Cancer Cell Detachment-Promoted Death. Cancers (Basel). 2019; 11:1094. https://doi.org/10.3390/cancers11081094 [PubMed]
- 13. Yuan Z, Li Y, Zhang S, Wang X, Dou H, Yu X, Zhang Z, Yang S, Xiao M. Extracellular matrix remodeling in tumor progression and immune escape: from mechanisms to treatments. Mol Cancer. 2023; 22:48. https://doi.org/10.1186/s12943-023-01744-8 [PubMed]
- 14. Chi H, Jiang P, Xu K, Zhao Y, Song B, Peng G, He B, Liu X, Xia Z, Tian G. A novel anoikis-related gene signature predicts prognosis in patients with head and neck squamous cell carcinoma and reveals immune infiltration. Front Genet. 2022; 13:984273. https://doi.org/10.3389/fgene.2022.984273 [PubMed]
- 15. He Z, Gu Y, Yang H, Fu Q, Zhao M, Xie Y, Liu Y, Du W. Identification and verification of a novel anoikis-related gene signature with prognostic significance in clear cell renal cell carcinoma. J Cancer Res Clin Oncol. 2023; 149:11661–78. https://doi.org/10.1007/s00432-023-05012-6 [PubMed]
- 16. Chen Y, Lin QX, Xu YT, Qian FJ, Lin CJ, Zhao WY, Huang JR, Tian L, Gu DN. An anoikis-related gene signature predicts prognosis and reveals immune infiltration in hepatocellular carcinoma. Front Oncol. 2023; 13:1158605. https://doi.org/10.3389/fonc.2023.1158605 [PubMed]
- 17. Rebhan M, Chalifa-Caspi V, Prilusky J, Lancet D. GeneCards: integrating information about genes, proteins and diseases. Trends Genet. 1997; 13:163. https://doi.org/10.1016/s0168-9525(97)01103-7 [PubMed]
- 18. Oliveros JC Venny. An interactive tool for comparing lists with Venn's diagrams. 2007-2015.
- 19. Zhou Y, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, Benner C, Chanda SK. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun. 2019; 10:1523. https://doi.org/10.1038/s41467-019-09234-6 [PubMed]
- 20. Tang Z, Li C, Kang B, Gao G, Li C, Zhang Z. GEPIA: a web server for cancer and normal gene expression profiling and interactive analyses. Nucleic Acids Res. 2017; 45:W98–102. https://doi.org/10.1093/nar/gkx247 [PubMed]
- 21. Yoshihara K, Shahmoradgoli M, Martínez E, Vegesna R, Kim H, Torres-Garcia W, Treviño V, Shen H, Laird PW, Levine DA, Carter SL, Getz G, Stemke-Hale K, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun. 2013; 4:2612. https://doi.org/10.1038/ncomms3612 [PubMed]
- 22. Newman AM, Liu CL, Green MR, Gentles AJ, Feng W, Xu Y, Hoang CD, Diehn M, Alizadeh AA. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods. 2015; 12:453–7. https://doi.org/10.1038/nmeth.3337 [PubMed]
- 23. Geeleher P, Cox N, Huang RS. pRRophetic: an R package for prediction of clinical chemotherapeutic response from tumor gene expression levels. PLoS One. 2014; 9:e107468. https://doi.org/10.1371/journal.pone.0107468 [PubMed]
- 24. Jiang P, Gu S, Pan D, Fu J, Sahu A, Hu X, Li Z, Traugh N, Bu X, Li B, Liu J, Freeman GJ, Brown MA, et al. Signatures of T cell dysfunction and exclusion predict cancer immunotherapy response. Nat Med. 2018; 24:1550–8. https://doi.org/10.1038/s41591-018-0136-1 [PubMed]
- 25. Nguyen DX, Bos PD, Massagué J. Metastasis: from dissemination to organ-specific colonization. Nat Rev Cancer. 2009; 9:274–84. https://doi.org/10.1038/nrc2622 [PubMed]
- 26. Jaguri A, Ahmad A. Breaching the Fortress of Tumor Microenvironment to Control Cancer Metastasis. Cancers (Basel). 2023; 15:4562. https://doi.org/10.3390/cancers15184562 [PubMed]
- 27. Agarwal S, Chand G, Jaiswal S, Mishra A, Agarwal G, Agarwal A, Verma AK, Mishra SK. Pattern and risk factors of central compartment lymph node metastasis in papillary thyroid cancer: a prospective study from an endocrine surgery centre. J Thyroid Res. 2012; 2012:436243. https://doi.org/10.1155/2012/436243 [PubMed]
- 28. Sipos JA, Mazzaferri EL. Thyroid cancer epidemiology and prognostic variables. Clin Oncol (R Coll Radiol). 2010; 22:395–404. https://doi.org/10.1016/j.clon.2010.05.004 [PubMed]
- 29. Cruchaga C, Del-Aguila JL, Saef B, Black K, Fernandez MV, Budde J, Ibanez L, Deming Y, Kapoor M, Tosto G, Mayeux RP, Holtzman DM, Fagan AM, et al, and Dominantly Inherited Alzheimer Network (DIAN), and Disease Neuroimaging Initiative (ADNI), and NIA-LOAD family study. Polygenic risk score of sporadic late-onset Alzheimer's disease reveals a shared architecture with the familial and early-onset forms. Alzheimers Dement. 2018; 14:205–14. https://doi.org/10.1016/j.jalz.2017.08.013 [PubMed]
- 30. Supruniuk K, Radziejewska I. MUC1 is an oncoprotein with a significant role in apoptosis (Review). Int J Oncol. 2021; 59:68. https://doi.org/10.3892/ijo.2021.5248 [PubMed]
- 31. Kim YS, Tang PW, Welles JE, Pan W, Javed Z, Elhaw AT, Mythreye K, Kimball SR, Hempel N. HuR-dependent SOD2 protein synthesis is an early adaptation to anchorage-independence. Redox Biol. 2022; 53:102329. https://doi.org/10.1016/j.redox.2022.102329 [PubMed]
- 32. Shen M, Kang Y. Stresses in the metastatic cascade: molecular mechanisms and therapeutic opportunities. Genes Dev. 2020; 34:1577–98. https://doi.org/10.1101/gad.343251.120 [PubMed]
- 33. Tang W, Feng X, Zhang S, Ren Z, Liu Y, Yang B, lv B, Cai Y, Xia J, Ge N. Caveolin-1 Confers Resistance of Hepatoma Cells to Anoikis by Activating IGF-1 Pathway. Cell Physiol Biochem. 2015; 36:1223–36. https://doi.org/10.1159/000430292 [PubMed]
- 34. Zhao W, Chen S, Hou X, Chen G, Zhao Y. CHK2 Promotes Anoikis and is Associated with the Progression of Papillary Thyroid Cancer. Cell Physiol Biochem. 2018; 45:1590–602. https://doi.org/10.1159/000487724 [PubMed]
- 35. Pan Y, Zhu X, Wang K, Chen Y. MicroRNA-363-3p suppresses anoikis resistance in human papillary thyroid carcinoma via targeting integrin alpha 6. Acta Biochim Biophys Sin (Shanghai). 2019; 51:807–13. https://doi.org/10.1093/abbs/gmz066 [PubMed]
- 36. Li J, Luo M, Ou H, Liu X, Kang X, Yin W. Integrin β4 promotes invasion and anoikis resistance of papillary thyroid carcinoma and is consistently overexpressed in lymphovascular tumor thrombus. J Cancer. 2019; 10:6635–48. https://doi.org/10.7150/jca.36125 [PubMed]
- 37. Kuo CY, Hsu YC, Liu CL, Li YS, Chang SC, Cheng SP. SOX4 is a pivotal regulator of tumorigenesis in differentiated thyroid cancer. Mol Cell Endocrinol. 2023; 578:112062. https://doi.org/10.1016/j.mce.2023.112062 [PubMed]
- 38. Tang M, Luo W, Zhou Y, Zhang Z, Jiang Z. Anoikis-related gene CDKN2A predicts prognosis and immune response and mediates proliferation and migration in thyroid carcinoma. Transl Oncol. 2024; 40:101873. https://doi.org/10.1016/j.tranon.2023.101873 [PubMed]
- 39. Kim H, Sung JY, Park EK, Kho S, Koo KH, Park SY, Goh SH, Jeon YK, Oh S, Park BK, Jung YK, Kim YN. Regulation of anoikis resistance by NADPH oxidase 4 and epidermal growth factor receptor. Br J Cancer. 2017; 116:370–81. https://doi.org/10.1038/bjc.2016.440 [PubMed]
- 40. Ferraro A, Mourtzoukou D, Kosmidou V, Avlonitis S, Kontogeorgos G, Zografos G, Pintzas A. EZH2 is regulated by ERK/AKT and targets integrin alpha2 gene to control Epithelial-Mesenchymal Transition and anoikis in colon cancer cells. Int J Biochem Cell Biol. 2013; 45:243–54. https://doi.org/10.1016/j.biocel.2012.10.009 [PubMed]
- 41. Byerly J, Halstead-Nussloch G, Ito K, Katsyv I, Irie HY. PRKCQ promotes oncogenic growth and anoikis resistance of a subset of triple-negative breast cancer cells. Breast Cancer Res. 2016; 18:95. https://doi.org/10.1186/s13058-016-0749-6 [PubMed]
- 42. Giannoni P, Narcisi R, De Totero D, Romussi G, Quarto R, Bisio A. The administration of demethyl fruticulin A from Salvia corrugata to mammalian cells lines induces "anoikis", a special form of apoptosis. Phytomedicine. 2010; 17:449–56. https://doi.org/10.1016/j.phymed.2009.07.007 [PubMed]
- 43. Zou G, Ren B, Liu Y, Fu Y, Chen P, Li X, Luo S, He J, Gao G, Zeng Z, Xiong W, Li G, Huang Y, et al. Inhibin B suppresses anoikis resistance and migration through the transforming growth factor-β signaling pathway in nasopharyngeal carcinoma. Cancer Sci. 2018; 109:3416–27. https://doi.org/10.1111/cas.13780 [PubMed]
- 44. Normanno N, De Luca A, Bianco C, Maiello MR, Carriero MV, Rehman A, Wechselberger C, Arra C, Strizzi L, Sanicola M, Salomon DS. Cripto-1 overexpression leads to enhanced invasiveness and resistance to anoikis in human MCF-7 breast cancer cells. J Cell Physiol. 2004; 198:31–9. https://doi.org/10.1002/jcp.10375 [PubMed]
- 45. Vermeer PD, Denker J, Estin M, Moninger TO, Keshavjee S, Karp P, Kline JN, Zabner J. MMP9 modulates tight junction integrity and cell viability in human airway epithelia. Am J Physiol Lung Cell Mol Physiol. 2009; 296:L751–62. https://doi.org/10.1152/ajplung.90578.2008 [PubMed]
- 46. Wang B, Shen W, Yan L, Li X, Zhang L, Zhao S, Jin X. Reveal the potential molecular mechanism of circRNA regulating immune-related mRNA through sponge miRNA in the occurrence and immune regulation of papillary thyroid cancer. Ann Med. 2023; 55:2244515. https://doi.org/10.1080/07853890.2023.2244515 [PubMed]
- 47. Ostan R, Lanzarini C, Pini E, Scurti M, Vianello D, Bertarelli C, Fabbri C, Izzi M, Palmas G, Biondi F, Martucci M, Bellavista E, Salvioli S, et al. Inflammaging and cancer: a challenge for the Mediterranean diet. Nutrients. 2015; 7:2589–621. https://doi.org/10.3390/nu7042589 [PubMed]
- 48. Coperchini F, Croce L, Marinò M, Chiovato L, Rotondi M. Role of chemokine receptors in thyroid cancer and immunotherapy. Endocr Relat Cancer. 2019; 26:R465–78. https://doi.org/10.1530/ERC-19-0163 [PubMed]
- 49. Berraondo P, Umansky V, Melero I. Changing the tumor microenvironment: new strategies for immunotherapy. Cancer Res. 2012; 72:5159–64. https://doi.org/10.1158/0008-5472.CAN-12-1952 [PubMed]
- 50. Sun J, Shi R, Zhang X, Fang D, Rauch J, Lu S, Wang X, Käsmann L, Ma J, Belka C, Su C, Li M. Characterization of immune landscape in papillary thyroid cancer reveals distinct tumor immunogenicity and implications for immunotherapy. Oncoimmunology. 2021; 10:e1964189. https://doi.org/10.1080/2162402X.2021.1964189 [PubMed]
- 51. Muenst S, Läubli H, Soysal SD, Zippelius A, Tzankov A, Hoeller S. The immune system and cancer evasion strategies: therapeutic concepts. J Intern Med. 2016; 279:541–62. https://doi.org/10.1111/joim.12470 [PubMed]
- 52. Yang Z, Wei X, Pan Y, Xu J, Si Y, Min Z, Yu B. A new risk factor indicator for papillary thyroid cancer based on immune infiltration. Cell Death Dis. 2021; 12:51. https://doi.org/10.1038/s41419-020-03294-z [PubMed]
- 53. Pani F, Yasuda Y, Rousseau ST, Bermea KC, Roshanmehr S, Wang R, Yegnasubramanian S, Caturegli P, Adamo L. Preconditioning of the immune system modulates the response of papillary thyroid cancer to immune checkpoint inhibitors. J Immunother Cancer. 2022; 10:e005538. https://doi.org/10.1136/jitc-2022-005538 [PubMed]