A Prognostic Model of Differentiated Thyroid Cancer Based on Up-Regulated Glycolysis-Related Genes

Article type: Research Article

Keywords: thyroid cancer, prognosis, glycolysis-related genes, recurrence, decision tree.

Authors: Min Wu, Deng-jie Ou-yang, Bo Wei, Pei Chen, Qi-man Shi, Hai-long Tan, Bo-qiang Huang, Mian Liu, Zi-en Qin, Ning Li, Hui-yu Hu, Peng Huang, Shi Chang

Affiliations: Department of General Surgery, Xiangya Hospital Central South University, Changsha, China; National Clinical Research Center for Geriatric Disorders, Xiangya Hospital, Changsha, China; Clinical Research Center for Thyroid Disease in Hunan Province, Xiangya Hospital, Changsha, China; Hunan Provincial Engineering Research Center for Thyroid and Related Diseases Treatment Technology, Xiangya Hospital, Changsha, China

License: Copyright © 2022 Wu, Ou-yang, Wei, Chen, Shi, Tan, Huang, Liu, Qin, Li, Hu, Huang and Chang CC BY 4.0 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.

Article links: DOI: 10.3389/fendo.2022.775278  |  PubMed: 35528004  |  PMC: PMC9072639

Relevance: Moderate: mentioned 3+ times in text

Full text: PDF (21.2 MB)

Introduction

Thyroid cancer is one of the most common endocrine malignancies, with an incidence that has substantially increased worldwide over the last several decades (ref. 1–ref. 3). Differentiated thyroid carcinoma (DTC) is the most common pathological subtype, including papillary thyroid carcinoma (PTC), follicular thyroid carcinoma (FTC), and Hürthle cell carcinoma (ref. 4). Generally, DTC has an excellent prognosis due to its indolent features and a better survival rate than other carcinomas; however, disease ultimately recurs in approximately 5–21% of DTC patients (ref. 5). In DTC patients who experience recurrent disease, surgical treatment is generally required, and reoperation poses additional medical costs and significant morbidity compared to the initial surgery. Therefore, preventing recurrence in DTC patients can reduce the deteriorations in quality of life related to reoperation (ref. 6–ref. 8).

According to previous reports, age, lymph node (LN) metastasis, tumor size, extrathyroidal extension (ETE), tumor multiplicity, and extranodal extension (ENE) are known risk factors for recurrence of DTC (ref. 9). Thus, the American Thyroid Association (ATA) management guidelines have proposed a clinicopathological risk stratification system based on those risk factors, which subdivides patients into high-risk and low-risk cohorts, providing prognostic and predictive information to facilitate clinical decisions. Unfortunately, since most prognostic signatures lack reproducibility due to individual heterogeneity, the current guidelines are unable to effectively predict recurrence in patients with DTC.

Our previous studies confirmed that metabolic reprogramming and aberrant expression of glycolysis-related genes have been characterized as hallmarks of DTC (ref. 10). A high rate of aerobic glycolysis, known as the “Warburg effect”, has been gradually confirmed to occur in various cancers, and treatment with novel therapeutic targets can be used as an effective anticancer strategy (ref. 11). Our previous studies also found that aerobic glycolysis is active in DTC, promoting thyroid tumorigenesis, proliferation and invasion (ref. 12). Additionally, mitochondrial glycerophosphate dehydrogenase (MGPDH), the key enzyme connecting oxidative phosphorylation (OXPHOS) and glycolysis, promotes Warburg metabolism in thyroid cancer and can be effectively targeted by metformin (ref. 13). Studies have shown that high lactate production due to abnormal glycolytic activity promotes tumor angiogenesis and invasion, which may play a role in suppressing immune cells against cancer and promoting tumor recurrence (ref. 14, ref. 15). Therefore, targeting glycolysis is a promising strategy for cancer treatment. In addition, a glycolysis-related prognostic score for the prediction of prognosis and chemosensitivity of pancreatic ductal adenocarcinoma has shown better performance than conventional methods (ref. 16). However, limited studies have systematically investigated the metabolic status and its prognostic value in patients with DTC. The association between the genetic characteristics of glycolysis and the heterogeneity of DTC has never been reported.

In this study, a glycolysis-related prognostic signature was developed based on expression profiling data in patients with DTC from several datasets. Histopathology and cellular experiments were used to investigate the utility of a simplified protein signature. In addition, the underlying molecular mechanisms in different risk subgroups were also explored. This study aimed to clarify the relationship between glycolysis and DTC, as well as to determine whether this prognostic signature can be used to detect a group of patients with DTC who have a high risk of recurrence.

Materials and Methods

Patients and Datasets

To identify the expression levels of glycolysis-related genes (GRGs) in DTCs, we retrospectively analyzed the gene expression profiles and clinical parameters of DTCs patients from 6 public cohorts, including 5 datasets (GSE33630, GSE60542, GSE58545, GSE35570, and GSE27155) from GEO and 1 dataset (TCGA-THCA) from TCGA. We downloaded the RNA-seq data of 510 cancer samples and 58 normal samples from TCGA, only patients meet the following criteria were included: 1) R0 surgical margins, 2) intact recurrence status, and 3) recurrence time > 30 days. Patients who were treated with RAI were also included in the analysis. Overall, a total of 379 patients were included in our study (Supplemental Table 1). Patients from those 5 GEO datasets (n=400, Supplemental Table 2) were enrolled to form the validation cohort of differentially expressed GRGs, while, patients from TCGA-THCA were randomly divided into training set (n=191), validation set (n=188) and test set (n=379) for model validation.

Candidate GRGs Selection and Signature Establishment

According to previous reports, 2528 protein-coding genes from GeneCards database were selected as GRGs (Supplemental Table 3) (ref. 17, ref. 18). Then, we performed gene set variation analysis (GSVA) in upregulated GRGs (log2FC ≥ 1 and FDR < 0.05), and GSVA score for each sample was defined as the glycolytic carcinogenicity score (ref. 19). The candidate GRGs were screened from those upregulated GRGs that co-expressed in TCGA-THCA, 5 GEO datasets, and functioned as prognosis-related genes. Then, prognosis-related DEGs were imported into LASSO regression to construct the most stable glycolysis-related risk score (GRS) as follows:

After the glycolysis-related gene signature (quantified by GRS) was established, univariate and multivariate Cox proportional hazards models were used to determine independent prognostic factors, thereafter, the independent prognostic factors, which were identified by multivariate Cox proportional hazards regression analysis, were used to create a composite risk score and displayed as a nomogram. Moreover, Schoenfeld residuals test were used to determine the proportional hazards assumption of the nomogram. Except for using ROC and calibration curves to measure the nomogram, we also performed decision curve analysis (DCA) to assess the clinical net benefit.

Mechanism Exploration

Differentially expressed genes (limma, |logFC|≥1 and FDR<0.05) between the GRS-high and GRS-low groups were screened and analyzed for GO and KEGG enrichment. We further analyzed the infiltration extent of 22 immune cells in different risk populations using the ESTIMATE and ssGSEA algorithm, as well as the correlation between GRS genes and immune cell infiltration via TIMER tool (ref. 20). In particular, the ESTIMATE algorithm estimates stromal and immune scores to predict the level of infiltrating stromal and immune cells based on their specific gene expression profile (ref. 21). For validation of each possible mechanism, we performed a pan-cancer analysis via GSCALite portal. Overall, 10 cancers including bladder urothelial carcinoma (BLCA), breast invasive carcinoma (BRCA), esophageal carcinoma (ESCA), head and neck squamous cell carcinoma (HNSC), kidney renal clear cell carcinoma (KIRC), kidney renal papillary cell carcinoma (KIRP), liver hepatocellular carcinoma (LIHC), lung adenocarcinoma (LUAD), lung squamous cell carcinoma (LUSC), and stomach adenocarcinoma (STAD) were analyzed.

Cell Culture and Plasmid Transfections

The DTC cancer cell lines TPC-1, BCPAP, and K1 and the normal thyroid Nthy-ori 3-1 cell line used in this study were provided by the Institute of Medical Sciences at Xiangya Hospital. All cell lines were cultured in a cell culture incubator at 37°C and 5% CO₂ according to our previous study (ref. 10). Blank plasmids or siRNAs targeting the four GRGs were manufactured by Jikai Company (Shanghai, China). Cells were transiently transfected with blank plasmid or siRNA using Lipofectamine 2000 (Invitrogen) according to the manufacturer’s instructions.

Western-Blot

Western blotting was performed as previously described (ref. 10) with Abs targeting ADM (ab190819, Abcam), CD44 (ab189524, Abcam), MKI67 (ab92742, Abcam) and TYMS (15047‐1‐AP, Proteintech).

Cell Proliferation

For the cell proliferation assay, after transfection for 72 hours, the BrdU incorporation efficiency of cells was measured using the Cell Proliferation Assay kit (Roche, 11669915001) according to the manufacturer’s instructions.

Statistical Analyses

Statistical analyses were performed using R software (version 3.6.1). The ESTIMATE package was used to calculate the immunoreactivity of the samples. The LASSO algorithm was performed via glmnet package. Survival analysis, and optimal cut-off values were performed and obtained through the survival and survminer packages. The rms, pROC and ggDCA package were used to plot calibration, ROC curves and DCA curves. The nomogram was plotted using the regplot package. The chi-square test was used to assess clinical correlation. All statistical tests with p-values less than 0.05 were considered significant if not otherwise stated.

Results

Schematic Diagram of the Study Design

First, glycolysis was identified as an important risk factor affecting relapse in DTC patients. Then, a combination of univariate Cox regression analysis and the LASSO algorithm was used to screen for robust glycolysis-related genes and create a glycolysis-related gene signature to predict recurrence and survival in DTC patients. Subsequently, the prognostic value of this gene signature was validated using the training and validation cohorts. In addition, model assessments were used to further evaluate its prognostic ability, and a decision tree was constructed to improve risk stratification for recurrence based on GRSs and other clinicopathological variables. Finally, histopathological examination and cellular experiments were used to explore potential molecular mechanisms in the high-risk group (f1).

Glycolysis Levels in Tumor Tissues Are Significantly Associated With Poor Prognosis in DTC Patients

To understand the correlation between aberrant expression of GRGs and the prognosis of DTC, we downloaded the RNA-seq data and matched clinical information of 379 DTC patients from TCGA. In our results, we identified 3258 differentially expressed genes (DEGs), including 1709 upregulated genes and 1549 downregulated genes (f2). Meanwhile, we screened 2353 glycolysis-related genes (GRGs) from the GeneCards website, 175 of which were selected as candidate genes by taking the intersection between DEGs and GRGs (f2 and Supplemental Figure 1). Glycolysis scores based on 175 GRGs were used to assess the glycolysis level in tumor tissue. We found that the glycolysis scores in RFS (+) DTC were significantly higher than those in RFS (-) DTC (p < 0.05, f2). Moreover, using ROC curve analysis, we set the optimal cut-off value for glycolysis scores to 0.582 and found that DTC patients with high glycolysis scores exhibited worse prognosis (f2). Furthermore, glycolysis scores were found to be significantly associated with the clinicopathological features of DTC, including diagnosis age (p < 0.001), clinical stage (p < 0.001), T stage (p = 0.021), N stage (p < 0.001), and histology grading (p < 0.001) (f2). Besides, when we compared the levels of glycolysis scores between different clinical subgroups using ANOVA tests (Supplemental Figure 2). We found that overall glycolysis scores trended upwards with disease progression (AJCC stage). And patients with positive lymph node metastases had significantly higher glycolysis scores than those without lymph node metastases. Among the different histological subgroup types, the glycolysis score was significantly higher in the Tall-cell subgroup than in the other subgroups. However, the difference in glycolysis scores between the elderly and the young was not significant. This may be related to statistical methods.

A GRS Model Based on Four Targeting Genes Has Potential for Predicting Recurrence in DTC Patients

We next validated the expression of 175 GRGs in five GEO cohorts (GSE33630, GSE60542, GSE58545, GSE35570 and GSE27155), and 51 GRGs were confirmed to be upregulated in all cohorts (f3). However, only 27 GRGs were further confirmed to be significantly associated with RFS in DTC patients (f3). Next, the LASSO Cox algorithm was conducted to screen for robust glycolysis-related prognostic biomarkers. We finally constructed glycolysis risk score (GRS) predictive models based on four targeting genes (ADM, CD44, MKI67 and TYMS), which were highly upregulated in tumor samples, including papillary and follicular thyroid cancers (Supplemental Figure 3). The GRS was calculated via the combination GRG panel, GRS = (0.176×expression value of ADM) + (0.452×expression value of MKI67) + (0.73×expression value of CD44) + (0.833×expression value of TYMS) (f3). ROC curve analysis indicated that our GRS model has good predictive value for recurrence in DTC patients (AUC=0.739) (f3). Meanwhile, validation of the model showed that GRS was significantly elevated in RFS (+) patients compared to RFS (-) patients in the training, validation and test sets (f3).

Increased Expression of the Four Targeting Genes Promotes Cell Proliferation and Is Significantly Related to Lymph Node Metastasis in DTC Patients

According to previous reports, the functions of CD44, ADM, TYMS, and MKI67 are primarily in promoting the Warburg effect and tumor growth in various cancers, while the function of those four genes in DTC has rarely been reported (Supplemental Table 4). Based on the available data from the HPA database (https://www.proteinatlas.org/), we observed high expression of CD44, ADM, TYMS, and MKI67 in tumor tissues compared to normal thyroid tissue (f4). The same results were confirmed using the data from TCGA (f4). Moreover, expression of TYMS in tumors with lymph node metastasis (N1) was higher than in tumors without (N0) (f4).

To elucidate the cellular function of the four targeting genes, we cultured four DTC cell lines and assessed changes in cell proliferation after alterations in their expression. We chose TPC-1 and BCPAP as cellular models since all four targeting genes are overexpressed in these lines (f4). Cell proliferation was determined by performing BrdU incorporation assays. In this study, silencing the four GRGs with siRNA resulted in a 40–70% decrease in cell proliferation in the two DTC cell lines (f4).

The Four GRGs Involved in Remodeling of the Tumor Immune Microenvironment

To clearly clarify the molecular mechanisms of the four GRGs that may be involved, we divided the DTC patients into a low-GRS group and a high-GRS group. Using the limma algorithm, we identified 159 downregulated genes and 494 upregulated genes in the high-GRS group (f5 and Supplemental Table 5), and those genes were primarily related to the activation of the immune signaling pathway (f5 and Supplemental Figure 4). Next, we compared the distribution of immune scores in different GRS groups. The results revealed that the high-GRS group exhibited higher immunoreactivity and worse prognosis (f5 and Supplemental Figure 5). In addition, we found that immunoreactivity scores were significantly correlated with glycolysis scores (f5). Next, we assessed the level of immune cell infiltration in different immune score groups, glycolysis score groups and risk groups using the ssGSEA algorithm. The results showed that the high-risk group population was also clustered into both the high glycolysis score group and the high immune score group, and the high-risk group displayed a higher level of immune cell infiltration (f5). Finally, we analyzed the correlation between the GRS genes and the degree of infiltration of the six immune cells using the TIMER algorithm. We found that CD44, MKI67 and TYMS were significantly associated with all six types of cells. In particular, the expression of GRS genes was positively correlated with the degree of B-cell infiltration in all of them. In addition, we explored the mutations of four GRS genes in thyroid cancer, their regulatory roles in ten classical tumor pathways, and their interactions (Supplemental Figures 6–8). We found that all GRS genes showed different degrees of mutation in thyroid cancer and were closely related to the cell proliferation cycle and epithelial-mesenchymal transition (EMT). Overall, our results suggest that the level of glycolysis may affect the patients’ immune status and thus the prognosis of thyroid cancer patients (f5).

The Four GRGs Model Exhibited a Generalizability Signature in Pan-Cancer

In the pan-cancer analysis, we found that expression levels of the GRS model were significantly upregulated in all 10 cancer tissues (f6 and Supplemental Table 6). Similarly, we observed that gene expression levels in the GRS model were strongly associated with recurrence-free survival for most tumors (f6 and Supplemental Table 7).

Single nucleotide variant analysis revealed that MKI67 was the most frequently mutated gene (61%) among the four genes, missense mutation was the most common type of mutation, and C-base mutation to T-base was the most frequent base mutation (f6). We also explored the copy number variation of GRS genes across cancers. The results showed that GRS genes exhibited amplification and deletion in most cancers, with BLCA, ESCA and LUSC ranking highest in the frequency of GRS gene variants (f6). In the analysis of the degree of immune cell infiltration, we found that the mRNA expression level of the GRS genes was highly correlated with the degree of infiltration of six immune cells in BLCA, KIRC, and LIHC. In particular, LIHC was the tumor with the highest degree of correlation (f6). Besides, based on the genetic model, we calculated the GRS score for pan-cancer. The effect of GRS score on the prognosis and recurrence of pancytopenia was also examined using Kaplan-Meier analysis (Supplemental Figure 9). The results showed that GRS score for most tumors such as, GRS in BLCA, ESCA, KIRC, KIRP, LIHC and LUAD still showed as a risk factor for tumor outcome. High GRS scores were associated with poorer prognosis (p<0.05). Overall, our results suggest that GRS model may also be prognostic and therapeutic targets for other cancers.

The Excellent Prognostic Predictive Efficacy of the GRS Model in DTC Patients Compared to Traditional Features

According to previous reports, age, lymph node metastasis (LNM), tumor size, extrathyroidal extension (ETE), tumor multiplicity, and extranodal extension (ENE) are known risk factors for the recurrence of DTC (ref. 22). Consistent with the results of successive studies, we found that GRS also GRG also plays a role in most clinical subgroups. Besides, in this study (f7), stage, T stage, N stage and GRS were significantly correlated with the recurrence of DTC (f7). However, when comparing the predictive performance of GRS to traditional prognostic indicators, we observed better prognostic predictive efficacy of the GRS model for DTC patients than for traditional features. ROC curve analysis showed that the AUC of GRS was 0.739, while the AUC of sex, age, stage and grade was 0.431, 0.565, 0.515 and 0.634, respectively (f7). Furthermore, we compared the sensitivity and accuracy of the discriminatory efficacy for the GRS model in the first, third, and fifth years. As shown in Supplemental Figure 10, the AUCs of the GRS model in the first, third and fifth years were 0.774, 0.751 and 0.744, respectively.

Integrated Model Constructed Using GRS and N-Stage Optimizes the Risk Stratification and RFS in DTC Patients

Since we confirmed that GRS and N-stage were independent risk factors for RFS in DTC (f7), we constructed an integrated model based on these two factors. To quantify the risk assessment for individual DTC patients, a nomogram was generated using GRS together with N-stage, and the red arrow shows an example (f8). According to the nomogram, when a DTC patient with high GRS and N-stage was N0, the risk of recurrence at 1 year, 3 years, and 5 years was 0.127, 0.403, and 0.426, respectively. The underlying proportional hazards assumptions of the nomogram model were verified using Schoenfeld residual tests (f8). Furthermore, we compared the value of the nomogram model to the GRS model and traditional clinicopathological features and found that the nomogram exhibited the most powerful capacity for RFS prediction (f8). Similarly, both calibration tests and DCA graphically illustrated that the nomogram model exhibited more accuracy in estimating patient recurrence than other parameters (f8).

Finally, to aid clinical decision-making, we constructed a recurrence decision tree with N-stage and GRS to optimize risk stratification. As shown in f8, three different risk subgroups were defined based on two primary components: patients with low GRS were defined as the “low risk” group, while patients with high GRS & N0 stage and high GRS & N1 stage were defined as “intermediate risk” and “high risk” groups, respectively. As shown in the K-M curve, RFS differed markedly between the three risk subgroups (f8).

Discussion

The increasing prevalence and incidence of DTC each year makes it a worldwide health concern. Although DTC is an indolent disease associated with excellent overall survival, its high recurrence rate remains a challenge for current postoperative management. Since recurrent disease often necessitates remedial surgery, additional administration of radioactive iodine (RAI) and more intensive long-term surveillance (ref. 23), the ability to determine the risk of disease recurrence may be a more meaningful outcome for both patients and clinicians. In the present study, we are the first to construct an integrated model incorporating clinicopathological features and glycolysis-related gene signatures, which performed well in differentiating between low‐ and high‐risk cohorts and showed good performance for predicting RFS in patients with DTC.

Although several follow-up data suggest minimal differences in disease-specific survival (DSS) between total thyroidectomy and lobectomy, the extent of surgery has a definite impact on recurrence in DTC patients (ref. 24). In this study, all patients with DTC who were used for analysis underwent total thyroidectomy; however, approximately 10% of patients relapsed after surgery. Therefore, a larger surgical scope does not prevent recurrence and increases the likelihood of surgical complications. According to previous reports, age, ETE, LNM, tumor size and tumor multifocality were found to be risk factors for the recurrence of DTC (ref. 22, ref. 25). Consistent with previous studies, we also found that stage, T stage and N stage were significantly associated with RFS in DTC patients; however, only N stage was an independent risk factor for RFS.

Emerging evidence has implied that predictive models based on gene signatures can quantify the risk of tumor recurrence, but few studies have attempted to apply them in DTC. Our previous studies confirmed that abnormal glycolysis plays an important role in promoting the invasion and metastasis of DTC (ref. 10, ref. 12). Similarly, in the present study, we found that the glycolysis level in tumor tissues was significantly associated with poor prognosis in DTC patients. Through univariate/multivariate Cox, LASSO and ROC analyses, we demonstrated that the GRS could be an independent prognostic factor in DTC patients, as calculated using the gene expression signatures of four GRGs (ADM, CD44, MKI67 and TYMS). It has been reported that CD44, MKI67 and TYMS are overexpressed in many cancer types, including DTC, and regulate glucose metabolism by targeting different genes (ref. 26–ref. 29). CD44, a marker of cancer stem-like cells and epithelial-mesenchymal transition, was found to increase CREB phosphorylation and sustain the proliferation of thyroid cancer cells (ref. 27). However, the function of the remaining 3 genes in DTC remains unclear. In the mechanistic validation of this study, we found that all four GRGs participate in promoting cell proliferation and metastasis, which may be attributed to their function in remodeling the tumor immune microenvironment of DTC. Interestingly, similar results were found in another ten cancers, which indirectly illustrates that the prognostic assessment function of the four-GRG signature may be generalizable across cancer types.

Multivariate Cox regression analysis demonstrated that GRS was a fairly strong predictor of risk for overall recurrence and survival, even stronger than age and clinical stage. To quantify the risk assessment for individual patients, we generated a nomogram including GRS and N-stage. Calibration analysis revealed that the nomogram showed very close to accurate prediction of actual survival. In addition, ROC analysis revealed that the nomogram showed more stable and stronger predictive power at different time points during follow-up compared to any other individual variable. We combined clinicopathological characteristics to construct a decision tree to improve risk stratification and optimize existing clinical guidelines. In the decision tree, GRS is the primary determinant. However, in the second node, the decision tree demonstrated that risk stratification was improved if GRS was replaced with N stage. Together, the decision tree and multivariate Cox results suggest that the glycolytic gene signature is indeed a strong risk factor for e recurrence and survival in DTC patients.

Some limitations of our study should be acknowledged. First, this was a retrospective study, and therefore, the prognostic robustness and clinical utility of the glycolysis-related gene signature need to be further validated in larger prospective trials. Second, further experimental studies are needed to elucidate the biological functions underlying the GRS genes related to tumor glycolysis.

In summary, we established a novel glycolysis-related gene signature to discriminate patients with DTC who are at high risk of recurrence. Integrating this with clinicopathological features, we constructed a decision tree to optimize risk stratification for recurrence-free survival and a nomogram to quantify risk assessment for individual patients. The GRG signature-based model may help clinicians provide a prognosis for DTC patients with a high risk of recurrence after surgery and provide further personalized treatment to decrease the chance of relapse.

Data Availability Statement

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

Ethics Statement

This study is based on published or public datasets and does not include new data that require ethical approval or consent.

Author Contributions

MW, PH, and SC conceived and planned the study design. MW, BW, PC, Q-mS, H-lT, D-jO-y, B-qH, ML, Z-eQ, NL, and H-yH performed formal analysis and data interpretation. MW and BW wrote the original draft. PH and SC provided critical revisions and contributed to the editing of the paper. All authors read and approved the final manuscript.

Funding

This work was supported by the National Natural Science Foundation of China (grant numbers 81974423, 81902729), the Key Research and Development Program of Hunan Province (grant number 2019SK2031), the China Postdoctoral Science Foundation (grant number 2020M672517, 2021T140749), the Natural Science Foundation of Hunan Province (grant number 2020JJ5904), and the Xiangya Hospital Foundation for Young Scholars (grant number 2018Q01).

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.

References

1. C La Vecchia, M Malvezzi, C Bosetti, W Garavello, P Bertuccio, F Levi. Thyroid Cancer Mortality and Incidence: A Global Overview.. Int J Cancer (, 2015. [DOI]
2. P Lin, Y-N Guo, L Shi, XJ Li, H Yang, Y He. Development of a Prognostic Index Based on an Immunogenomic Landscape Analysis of Papillary Thyroid Cancer.. Aging (Albany NY (, 2019. [DOI | PubMed]
3. T Zha, H Wu. Expression of Serum AMPD1 in Thyroid Carcinoma and its Clinical Significance.. Exp Ther Med (, 2018. [DOI]
4. J Park, K Kim, DJ Lim, JS Bae, JS Kim. Male Sex Is Not an Independent Risk Factor for Recurrence of Differentiated Thyroid Cancer: A Propensity Score-Matching Study.. Sci Rep (, 2021. [DOI | PubMed]
5. FH Liu, SF Kuo, C Hsueh, TC Chao, JD Lin. Postoperative Recurrence of Papillary Thyroid Carcinoma With Lymph Node Metastasis.. J Surg Oncol (, 2015. [DOI]
6. L Zhang, Y Wang, X Li, Y Wang, K Wu, J Wu. Identification of a Recurrence Signature and Validation of Cell Infiltration Level of Thyroid Cancer Microenvironment.. Front Endocrinol (Lausanne) (, 2020. [DOI]
7. H Zhao, S Zhang, S Shao, H Fang. Identification of a Prognostic 3-Gene Risk Prediction Model for Thyroid Cancer.. Front Endocrinol (Lausanne) (, 2020. [DOI]
8. M Wu, H Yuan, X Li, Q Liao, Z Liu. Identification of a Five-Gene Signature and Establishment of a Prognostic Nomogram to Predict Progression-Free Interval of Papillary Thyroid Carcinoma.. Front Endocrinol (Lausanne) (, 2019. [DOI]
9. CH Park, CM Song, YB Ji, JY Pyo, KJ Yi, YS Song. Significance of the Extracapsular Spread of Metastatic Lymph Nodes in Papillary Thyroid Carcinoma.. Clin Exp Otorhinolaryngol (, 2015. [DOI]
10. P Huang, LF Mao, ZP Zhang, WW Lv, XP Feng, HJ Liao. Down-Regulated miR-125a-5p Promotes the Reprogramming of Glucose Metabolism and Cell Malignancy by Increasing Levels of CD147 in Thyroid Cancer.. Thyroid (, 2018. [DOI]
11. JK Byun, YK Choi, YN Kang, BK Jang, KJ Kang, YH Jeon. Retinoic Acid-Related Orphan Receptor Alpha Reprograms Glucose Metabolism in Glutamine-Deficient Hepatoma Cells.. Hepatology (, 2015. [DOI]
12. P Huang, S Chang, X Jiang, J Su, C Dong, X Liu. RNA Interference Targeting CD147 Inhibits the Proliferation, Invasiveness, and Metastatic Activity of Thyroid Carcinoma Cells by Down-Regulating Glycolysis.. Int J Clin Exp Pathol (, 2015
13. S Thakur, B Daley, K Gaskins, VV Vasko, M Boufraqech, D Patel. Metformin Targets Mitochondrial Glycerophosphate Dehydrogenase to Control Rate of Oxidative Phosphorylation and Growth of Thyroid Cancer.. In Vitro In Vivo. Clin Cancer Res (, 2018. [DOI]
14. N Lee, W-J Jang, JH Seo, S Lee, C-H Jeong. 2-Deoxy-D-Glucose-Induced Metabolic Alteration in Human Oral Squamous SCC15 Cells: Involvement of N-Glycosylation of Axl and Met.. Metabolites (, 2019. [DOI]
15. Z Husain, Y Huang, P Seth, VP Sukhatme. Tumor-Derived Lactate Modifies Antitumor Immune Response: Effect on Myeloid-Derived Suppressor Cells and NK Cells.. J Immunol (, 2013. [DOI]
16. XP Zhang, Q Chen, Q Liu, Y Wang, F Wang, ZM Zhao. Development and Validation of Glycolysis-Related Prognostic Score for Prediction of Prognosis and Chemosensitivity of Pancreatic Ductal Adenocarcinoma.. J Cell Mol Med (, 2021. [DOI]
17. X Zhao, Z Liu, Z Ren, H Wang, Z Wang, J Zhai. Triptolide Inhibits Pancreatic Cancer Cell Proliferation and Migration via Down-Regulating PLAU Based on Network Pharmacology of Tripterygium Wilfordii Hook F.. Eur J Pharmacol (, 2020. [DOI | PubMed]
18. M Safran, I Dalah, J Alexander, N Rosen, TI Stein, M Shmoish. GeneCards Version 3: The Human Gene Integrator.. Database (, 2010. [DOI]
19. C Zou, C Yuan, J Ye, Z Liu, X Gao, X Piao. Identification and Validation of a Ten-Gene Set Variation Score as a Diagnostic and Prognostic Stratification Tool in Hepatocellular Carcinoma.. Am J Transl Res (, 2020
20. T Li, J Fan, B Wang, N Traugh, Q Chen, JS Liu. TIMER: A Web Server for Comprehensive Analysis of Tumor-Infiltrating Immune Cells.. Cancer Res (, 2017. [DOI]
21. H Wang, X Wu, Y Chen. Stromal-Immune Score-Based Gene Signature: A Prognosis Stratification Tool in Gastric Cancer.. Front Oncol (, 2019. [DOI | PubMed]
22. H Yan, X Zhou, H Jin, X Li, M Zheng, X Ming. A Study on Central Lymph Node Metastasis in 543 Cn0 Papillary Thyroid Carcinoma Patients.. Int J Endocrinol (, 2016. [DOI | PubMed]
23. S Filetti, C Durante, D Hartl, S Leboulleux, LD Locati, K Newbold. Thyroid Cancer: ESMO Clinical Practice Guidelines for Diagnosis, Treatment and Follow-Updagger.. Ann Oncol (, 2019. [DOI]
24. RM Gartland, CC Lubitz. Impact of Extent of Surgery on Tumor Recurrence and Survival for Papillary Thyroid Cancer Patients.. Ann Surg Oncol (, 2018. [DOI]
25. X Lan, W Sun, H Zhang, W Dong, Z Wang. Zhang T. A Meta-Analysis of Central Lymph Node Metastasis for Predicting Lateral Involvement in Papillary Thyroid Carcinoma.. ORG Head Neck Surg (, 2015. [DOI]
26. Z Zhang, T Xu, W Qin, B Huang, W Chen, S Li. Upregulated PTPN2 Induced by Inflammatory Response or Oxidative Stress Stimulates the Progression of Thyroid Cancer.. Biochem Biophys Res Commun (, 2020. [DOI]
27. V De Falco, A Tamburrino, S Ventre, MD Castellone, M Malek, SN Manié. CD44 Proteolysis Increases CREB Phosphorylation and Sustains Proliferation of Thyroid Cancer Cells.. Cancer Res (, 2012. [DOI]
28. MA Hossain, TA Asa, MM Rahman, S Uddin, AA Moustafa, JMW Quinn. Network-Based Genetic Profiling Reveals Cellular Pathway Differences Between Follicular Thyroid Carcinoma and Follicular Thyroid Adenoma.. Int J Environ Res Public Health (, 2020. [DOI]
29. R Gao, D Li, J Xun, W Zhou, J Li, J Wang. CD44ICD Promotes Breast Cancer Stemness via PFKFB4-Mediated Glucose Metabolism.. Theranostics (, 2018. [DOI]