Introduction
Thyroid carcinoma (TC) is the most common endocrine malignant tumor [1], categorized into four main subtypes based on histology [2]. Anaplastic thyroid carcinoma (ATC) is the most aggressive form, accounting for only 2% of TC cases but responsible for over 40% of thyroid cancer deaths [3]. Clinically, ATC is characterized by localized aggressive disease, a high metastasis rate, high mortality, and a low preoperative diagnosis rate [4]. Despite advances in treatments like radioiodine therapy, chemotherapy, and surgery, survival outcomes of ATC remain poor [5]. Therefore, novel and curative treatment strategies are urgently needed.
Recent research identified multiple critical factors in the development and progression of ATC. For instance, it was indicated that upregulation of abnormal spindle-like microcephaly-associated protein (ASPM) promoted ATC tumorigenesis, suggesting ASPM as a potential therapeutic target [6]. Furthermore, it was also demonstrated that 5-hydroxymethylcytosine (5hmC) had significant prognostic and therapeutic values among ATC patients [7]. However, while these studies enhanced our understanding of ATC, they primarily focus on single gene regulation and lack comprehensive genomic analysis from larger clinical datasets. Therefore, further systematic studies are essential to elucidate the molecular mechanisms underlying ATC and identify new diagnostic biomarkers and therapeutic targets.
Though several studies attempted to investigate the molecular mechanisms of ATC using high-throughput microarrays [8], they primarily focused on identifying differentially expressed genes (DEGs) from limited gene sets, largely overlooking the relationships between these genes [9–11]. Moreover, reliance on DEGs alone may miss potential therapeutic targets [12]. In this regard, weighted gene co-expression network analysis (WGCNA) offers a more effective approach to elucidate gene association patterns across different samples [13] and has been widely applied in cancer research [14, 15]. Nevertheless, although a study using WGCNA showed that TMEM158 could serve as a diagnostic biomarker for ATC, comprehensive associated survival analysis was still absent [16]. Therefore, there is still a need to analyze and identify new biomarkers for ATC from larger datasets combined with the WGCNA method.
Given this background, we systematically analyzed gene expression profiles from 4 Gene Expression Omnibus (GEO) microarray datasets between ATC and normal samples based on 2 different network analyses (protein–protein interaction [PPI] and WGCNA). Key genes at the intersection of these 2 methods were employed to identify the real hub genes that were essential to the development and progression of ATC. Additionally, we established an independent prognostic signature and a nomogram with clinical application. This study identified reliable biomarkers for the development, progression, and prognosis of ATC and provided a novel scientific basis for future therapy in ATC.
Materials and methods
Data collection and data processing
The microarray datasets (GSE33630, GSE65144, GSE29265, GSE53072) [17–20] and corresponding clinical features were downloaded from the Gene Expression Omnibus (GEO, https://www.ncbi.nlm.nih.gov/geo/) database. The datasets GSE33630, GSE65144, and GSE29265 were utilized to screen for differentially expressed genes (DEGs) by constructing PPI networks. The “sva” R package was applied to eliminate the batch effect of the 3 datasets. The dataset GSE53072 was utilized to perform WGCNA analysis. If multiple probes corresponded to the same gene, the maximum expression value was considered as the gene expression level.
Meanwhile, the RNA sequencing (RNA-seq) data of (thyroid cancer) THCA was downloaded from The Cancer Genome Atlas (TCGA, https://portal.gdc.cancer.gov/) database. The dataset included 502 THCA tumor samples and 58 normal samples. After excluding 52 samples without survival time, 450 tumor samples and corresponding clinical features were ultimately included in the study. The RNA-seq data of raw count normalized from the above databases was utilized for differential expression analysis. Then, the count values were converted to transcripts per kilobase million (TPM) values, which were further transformed to log2 (TPM + 1) for subsequent analysis.
In both databases, the p-values of DEGs were identified using a Wilcox test, and the cutoff value was set to |log2FC |> 1 and adjusted p < 0.05. The screening process was conducted using the “limma” R package.
PPI network construction and analysis of modules
PPI network was used to reveal the specific and unspecific interactions of protein and then to identify key genes and pivotal modules involved in the development of ATC. The online analysis tool STRING (https://www.string-db.org/) was used to construct PPI networks and the cutoff confidence score was set as > 0.7. Afterwards, Cytoscape software (version 3.9.0) was used to visualize the PPI network. Five methods in the plug-in CytoHubba were utilized to screen the key genes, namely MCC (Maximal clique centrality), MNC (Maximal neighborhood component), Degree (Node connect degree), EPC (Edge percolated component), and Closeness (Node connect degree). The top 20 genes in each method were selected, and then intersections were identified as key genes in PPI networks. In addition, the plug-in Molecular Complex Detection (MCODE) was employed to screen the significant module.
Weighted gene co-expression network analysis
Weighed Gene Co-expression Network Analysis (WGCNA) was an algorithm capable of transforming expression data into co-expression gene modules and exploring the relationship between modules and phenotypes. The top 5000 genes in GSE53072 were used to construct a co-expression network via the “WGCNA” R package. Firstly, we used the function “goodSamplesGenes” to check if the genes and samples were outliers. Subsequently, Pearson’s correlation coefficient of any two genes was calculated, and a correlation matrix was established. Then, a weighted value β conforming to the scale-free network law was set (scale-free R2 = 0.85). The correlation index αij of any gene pair was calculated based on the β square of the correlation coefficient Sij to construct an adjacency matrix. After determining the optimal β value for the soft threshold parameter, the adjacency matrix was turned into a topological overlap matrix (TOM). Eventually, the degree of dissimilarity dij between nodes was calculated based on topological coefficients ωij, and a dynamic hybrid tree cutting algorithm was further utilized to build a hierarchical clustering tree to divide modules. The specific calculation formulas were as follows:where cor represents Pearson’s correlation test. i, j denote gene i and gene j, respectively. , and represented node connectivity.
Identify significant relevant module
To identify the modules that were significantly associated with ATC, the module eigengene (ME) was utilized to characterize the expression profiles of each module and the correlation between ATC. Subsequently, the values of Gene Significance (GS) were utilized to calculate the association of the individual gene with ATC. In addition, the Module Membership (MM) was defined as the ME correlation and the gene expression profile for each module. If a gene within a module had both large MM and GS, then the gene was considered to be a key gene in the module.
Hub gene identification and verification
In this study, |GS |> 0.9 and |MM|> 0.9 were set as the threshold for screening candidate hub genes in the module. The hub genes in key modules were compared with key genes derived from PPI to obtain decisive genes in the development of ATC. For further validation of the hub genes, the immunohistochemistry (IHC) database (The Human Protein Atlas, https://www.proteinatlas.org/) was also employed.
Functional and pathway enrichment analysis
Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways analysis were performed using the “clusterProfiler” R package [21]. GO enrichment analysis included biological processes (BP), cellular components (CC), and molecular functions (MF). KEGG was a public database that integrates genomic, chemical, and systemic functional information. A P-value cutoff of < 0.05 was set to determine enriched GO terms and KEGG pathways. In addition, to increase the reliability and credibility of the results, the Gene Set Enrichment Analysis (GSEA) was employed to explore the potential molecular mechanisms.
Survival analysis and establishment of predictive nomogram
The clinical information was extracted from the TCGA database for survival analysis. After removing samples without survival time, we randomly divided the samples into a training cohort (225 samples) and a validation cohort (225 samples). We used the training cohort to establish a model and verified it in the validation cohort. Univariate Cox regression was performed to screen for genes significantly associated with prognosis (p < 0.05). Subsequently, LASSO Cox regression using the “glmnet” R package was applied to identify the optimal candidate genes combination to construct the prognostic gene signature [22]. The optimal value of the penalty parameter λ was determined by tenfold cross-validation based on minimum criteria. According to the coefficients calculated by LASSO Cox regression and expression level of each mRNA, the individual risk score of each patient could be calculated through the following formula:n was the number of prognostic genes, Expi was the expression value of gene i, and β was the regression coefficient. Simultaneously, all THCA patients were separated into high- and low-risk categories according to the median risk score. Kaplan–Meier survival curves were carried out to assess the survival difference of the two categories via the “survminer” R package. The “survival” and “timeROC” R package were utilized to perform the time-dependent receiver operating characteristic (ROC) curve analysis, which was applied to evaluate the prognostic gene signature’s predictive performance. Moreover, the independent prognostic significance of this identified gene signatures as well as other clinicopathological features were further investigated according to both univariate and multivariate Cox regression analysis. The expression levels of the prognostic gene signature in different risk groups were visualized by boxplots. Then, a three-dimensional principal component analysis (PCA) was carried out to explore the difference of the distribution between the low-risk group and the high-risk group. In addition, to make our signature to be more effectively applied to the clinical process, we took the clinical features (age, gender, stage, T, N, M) into the prognostic model and constructed a nomogram by the “rms” and “survival” R package. Moreover, the same methods were employed to further validate the prognostic capacity of the gene signature in the validation cohort.
Cell culture and quantitative RT-PCR (RT-qPCR)
Human ATC cell lines (CAL62 cell) and human normal thyroid follicular epithelial cell line (Nthy-ori 3–1 cell) were purchased from Ybio (Shanghai, China). Cells were cultured in Dulbecco’s modified Eagle’s medium (DMEM) containing 10% fetal bovine serum and 1% streptomycin at 37 °C with 5% CO2. In addition, total RNA was respectively extracted from CAL-62 cells and Nthy-ori 3–1 cells using Super Total RNA Extraction Kit (promega, China). Next, quantitative real-time polymerase chain reaction (RT-qPCR) was performed with the qPCR Master Mix Kit (promega, China) in a total volume of 20 μL. Relative mRNA expression levels of target genes were calculated using the 2−ΔΔCt method and GAPDH was used as the internal control. Primer information is shown in Online Resource Table S1.
Statistical analysis
All statistical analyses and visualizations were conducted using R software (version 4.2.0). Spearman correlation coefficients were calculated to assess associations. The Wilcoxon rank sum test was employed to evaluate differences in the expression levels of prognostic genes across different groups. Student's t-test was utilized to compare significant differences between two groups in the in vitro studies. Kaplan−Meier survival curves and the Log-rank test were applied to analyze patient survival rates. A p-value of less than 0.05 was considered statistically significant.
Results
Identification of differentially expressed genes
This study was carried out based on the flowchart (Fig. 1). We first collected 37 tumor samples and 82 normal samples of ATC from the GEO database (the detailed information was shown in Online Resource Table S2). The batch correction was performed to eliminate the batch effect of 3 datasets: GSE33630, GSE65144, and GSE29265 (Fig. 2A, B). Then, 1274 significantly up-regulated DEGs and 1431 significantly down-regulated DEGs in merged GEO microarray datasets were identified (Fig. 2C). The heatmap of DEGs was shown in Fig. 2D. These DEGs were utilized to perform subsequent PPI analysis. The intersection of 4 GEO datasets was shown in Fig. 2E, F, including 434 significantly up-regulated and 563 significantly down-regulated genes.
Fig. 1 [Images not available. See PDF.]
Data sources and analysis flow chart. DEGs, differentially expressed genes; GEO, The Gene Expression Omnibus; GO, Gene Ontology; HPA, Human Protein Atlas; KEGG, Kyoto Encyclopedia of Genes and Genomes; PCA, principal component analysis; PPI, protein–protein interaction; ROC, receiver operating characteristic; RT-qPCR, reverse transcription-quantitative polymerase chain reaction; THCA, thyroid cancer; WGCNA, weighed gene co-expression network analysis
Fig. 2 [Images not available. See PDF.]
Data preprocessing and identification of DEGs. A PCA analysis for 3 GEO datasets. B PCA analysis for 3 GEO datasets after elimination of batch effects. C Volcano plots of DEGs in merged GEO datasets. D Heatmap of DEGs expression from 3 GEO datasets. E Venn diagrams of up-regulated DEGs. F Venn diagrams of down-regulated DEGs. DEGs, differentially expressed genes; PCA, principal component analysis; GEO, Gene Expression Omnibus
Construction of PPI network and functional annotation
The PPI network was constructed by Cytoscape software based on the STRING database, comprising 416 nodes and 2514 edges (Fig. 3A). The top 20 genes screened by the 5 methods in CytoHubba were selected, and then intersections were identified as key genes for ATC in PPI analysis. These genes were: KIF2C, PBK, TOP2A, CDK1, KIF20A, and ASPM, which might play essential roles in ATC progression (Online Resource Table S3). Subsequently, MCODE within Cytoscape was employed to conduct module examination, and 5 modules were obtained from the PPI network (Online Resource Table S4). Among these results, module 1 stood out as the most important module (MCODE score = 41.238). This module included 43 vertices and 1732 links (Fig. 3B). Notably, the 43 genes within this module were all upregulated and the majority of the top 20 genes identified by five methodologies were also encompassed in this module. Moreover, to further explore the potential biological roles and pathways of these genes, the DEGs within this module were subjected to GO and KEGG analyses. GO enrichment analysis revealed that the genes within this module were significantly enriched in biological processes related to mitosis nuclear division and nuclear division (Fig. 3C). KEGG analysis proposed that genes were chiefly implicated in the cell cycle (Fig. 3D).
Fig. 3 [Images not available. See PDF.]
PPI network construction and modules analysis. A PPI networks diagram. Red represents upregulated genes and blue represents downregulated genes. B The most significant module of PPI network. The gradation of color positively associated with the expression level of this gene. C GO analysis and D KEGG pathway analysis of the most significant module in PPI analysis. PPI, protein–protein interaction; GO, gene ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes
Construction of weighted gene co-expression network and identification of key modules
We used the “WGCNA” R package to construct the weighted co-expression network. No outlier samples were identified by cluster analysis (Fig. 4A). A power of β = 20 (scale-free R2 = 0.85) was selected as a soft threshold to establish a gene regulatory network (Fig. 4B). Subsequently, 12 modules were identified using a dynamic pruning method (Fig. 4C). Among these, blue, brown and turquoise modules showed the strongest correlation with cancer traits (Fig. 4D). Additionally, we calculated the module eigengene and the Pearson correlation coefficients between modules based on clustering. The results indicated that the blue and brown modules exhibited the highest consistency (Fig. 4E). Ultimately, MM and GS values above the threshold of all genes in modules were screened as hub genes. The gene distribution of the blue, brown and turquoise modules was shown in Fig. 5A–C. Interestingly, we found that most of the top 20 genes in the 5 PPI methods were in the blue module. Therefore, we chose the blue module including 210 genes for further analysis.
Fig. 4 [Images not available. See PDF.]
WGCNA analysis of the GSE53072 dataset. A Cluster analysis of samples. B Analysis of scale-free index and average connectivity for various soft-threshold powers. C Dendrogram of the top 5000 genes of a variance cluster based on a dissimilarity measure (1-TOM). Each black vertical line represents a gene, and all genes are clustered into 12 different color modules. D Heat map of the correlation between ME and clinical traits. E Hierarchical clustering dendrogram and heatmap of MEs. WGCNA, weighed gene co-expression network analysis; TOM, topological overlap matrix; ME, module eigengene
Fig. 5 [Images not available. See PDF.]
Identification of key modules. A–C Scatterplots of GS versus MM in the blue module, brown module and turquoise module. D GO analysis of the blue module. GS, gene significance; MM, module membership; GO, gene ontology
Blue module functional annotation
The GO analysis results, as shown in Fig. 5D and Online Resource Fig. S1A, indicated a strong association of the blue module with glycosaminoglycan binding, chromosome segregation, mitotic nuclear division, and chromosomal region. Additionally, KEGG analysis revealed that genes within the blue module were predominantly enriched in the cell cycle and DNA replication pathways (Online Resource Fig. S1B). Furthermore, GSEA analysis indicated that genes were principally involved in cellular senescence, oocyte meiosis and progesterone-mediated oocyte maturation (Online Resource Fig. S1C).
Identification and verification of real hub gene
To further screen for the most significant hub genes, we combined two methods (PPI and WGCNA). Our findings revealed that the 6 key genes (KIF2C, PBK, TOP2A, CDK1, KIF20A, and ASPM) identified in the PPI analysis were all encompassed within the blue module in the WGCNA analysis. Additionally, KEGG analysis indicated that the blue module was related to the cell cycle, which was consistent with PPI submodule analysis. Furthermore, GO analysis highlighted a close correlation between the blue module, PPI submodule, and mitotic nuclear division. These results strongly suggested that the 6 genes identified played a critical role in the development of ATC. Consequently, we defined these 6 genes as the real hub genes relevant to ATC. Immunohistochemistry (IHC) results from the Human Protein Atlas (HPA) database were examined, revealing a significant up-regulated expression of 5 genes (KIF2C, PBK, TOP2A, CDK1, and KIF20A) in ATC (Fig. 6A–E and Online Resource Table S5). Furthermore, RT-qPCR analysis revealed significantly higher expression levels of 4 genes (KIF2C, CDK1, KIF20A, and ASPM) in CAL62 cells compared to Nthy-ori 3–1 cells (Fig. 6F and I–K), which aligns with mRNA expression data from 4 GEO datasets (Online Resource Fig. S2A–D). However, validation through RT-qPCR showed no significant difference in the relative expression levels of PBK and TOP2A between ATC and normal samples (Fig. 6G–H).
Fig. 6 [Images not available. See PDF.]
IHC and RT-qPCR validation about significant hub genes. Immunohistochemistry images of KIF2C (A), PBK (B), TOP2A (C), CDK1 (D), and KIF20 (E) in normal and ATC tissues obtained from HPA database. The RT-qPCR was used to compare mRNA levels of KIF2C (F), PBK (G), TOP2A (H), CDK1 (I), KIF20A (J), and ASPM (K) in CAL-62 cells and Nthy-ori 3–1 cells. HPA, Human Protein Atlas; IHC, immunohistochemistry; RT-qPCR, reverse transcription-quantitative polymerase chain reaction; ATC, anaplastic thyroid carcinoma. (P values: ns, nonsignificant; *P < 0.05; **P < 0.01; ***P < 0.001)
Survival analysis
To establish an effective model for predicting prognostic status, univariate Cox regression and LASSO Cox regression analysis were employed to screen genes. Firstly, 32 genes of the blue module were significantly associated with prognosis according to univariate Cox regression analysis in the training cohort (Online Resource Table S6). Subsequently, based on 1000 times tenfold cross-validation in LASSO Cox regression analysis, the minimum of λ value was selected as the optimum λ value (0.0121). Then, 6 of 32 genes with non-zero coefficients, namely GPSM2, FGF5, ASXL3, CYP4B1, CLMP, and DUXAP9 were screened to construct the gene signature according to the optimum λ (Fig. 7A, B).
Fig. 7 [Images not available. See PDF.]
Construction of a prognostic signature in training cohort. A Different colors represent different genes. B Cross-validation for turning the parameter selection in the LASSO Cox regression. C Risk scores, survival time, survival status, and 6-gene signature expression in the training set. D Kaplan–Meier OS curves for high-risk and low-risk groups. E The 1-, 3-, and 5-year ROC curve to predict the survival status. LASSO, the least absolute shrinkage and selection operator; OS, overall survival; ROC, receiver operating characteristic
The risk score calculated using the following formula: Risk score = (1.0505* expression of GPSM2) + (0.8805* expression of FGF5) + (0.1120* expression of ASXL3) + (− 0.1338* expression of CYP4B1) + (0.2679* expression of CLMP) + (0.3317* expression of DUXAP9). Subsequently, the 225 patients were divided into high-risk (n = 112) and low-risk (n = 113) groups based on the median risk score. Scatter plots and the corresponding heatmap were plotted based on the survival time and risk score in the training cohort (Fig. 7C) and validation cohort (Online Resource Fig. S3A). We found that the death cases were focused on the high-risk group, while the surviving cases were centralized in the low-risk group. These results demonstrated that the definition of “risk score” was effective. To evaluate the impact of risk scores on prognosis, we conducted a Kaplan–Meier curve analysis. The results showed that patients with low-risk scores had significantly better overall survival than those with high-risk scores in both the training cohort (p = 0.0059, Fig. 7D) and validation cohort (p = 0.0073, Online Resource Fig. S3B). Moreover, time-sensitive ROC analysis uncovered that this genetic signature displayed a beneficial predictive performance in the instructional group with AUCs of 0.989, 0.921, and 0.959 at 1-, 3-, and 5-year intervals, respectively (Fig. 7E). Correspondingly, in the confirmation group, the AUCs were NA, 0.680, and 0.678, respectively (Online Resource Fig S3C). Among the six genes, GPSM2, FGF5, ASXL3, CLMP, and DUXAP9 were notably upregulated in the high-risk cohort, correlating with poor prognosis, while CYP4B1 was significantly upregulated in the low-risk cohort and associated with favorable prognosis (Fig. 8A–C and Online Resource Fig. S3D). To further confirm these results, the mRNA levels of these 6 genes were determined in CAL-62 cells and normal Nthy-ori 3–1 cells with RT-qPCR. The findings showed that the mRNA expression levels of GPSM2, FGF5, ASXL3, CLMP and DUXAP9 were higher in cancer tissues, whereas CYP4B1 expression was higher in normal tissue (Fig. 8D–I), consistent with our bioinformatics analysis results.
Fig. 8 [Images not available. See PDF.]
Prognostic study and the RT-qPCR results of the training cohort via the signature model. A The bar charts display the variables and corresponding coefficients in the prognostic signature. B The forest plot shows the univariate regression results of 6 genes. C Box plot visualization of the expression levels of 6 prognostic gene in different risk groups. The RT-qPCR was used to compare mRNA levels of GPSM2 (D), FGF5 (E), ASXL3 (F), CYP4B1 (G), CLMP (H), and DUXAP9 (I) in CAL-62 cells and Nthy-ori 3–1 cells. RT-qPCR, reverse transcription-quantitative polymerase chain reaction. (P values: ns, nonsignificant; *P < 0.05; **P < 0.01; ***P < 0.001)
Independent prognostic role of the six-gene signature and establishment of a nomogram
To evaluate the risk score as an independent prognostic factor, we performed univariate and multivariate Cox regression analyses. Results of both analyses indicated demonstrated that the risk score was a significant independent predictor of poor survival (univariate: p < 0.001, HR = 3.176, 95% CI 3.101–19.148; multivariate: p = 0.015, HR = 44.801, 95% CI 2.090–960.279; Fig. 9A, B). We then explored the correlation between the risk score and clinicopathological parameters in the training cohort using the Wilcoxon signed-rank test (Online Resource Fig. S4A–F). Additionally, we developed a nomogram that integrates the risk score and clinicopathological parameters to predict overall survival (OS) rates at 1, 3, and 5 years (Fig. 9C, D). PCA analyses further showed identifiable dimensions between the high-risk and low-risk categories (Fig. 9E).
Fig. 9 [Images not available. See PDF.]
Analysis of the prognostic signature. A Univariate and B multivariate Cox regression of the risk score and other clinical characteristics in the training cohort. C A nomogram was constructed to predict 1-, 3-, and 5-year survival probabilities. D Calibration curves for 1-, 3-, and 5- year prediction was utilized to depict the consistency between predicted and actual survival. E 3-dimensional PCA plot showing distribution in the 2 groups. PCA, principal component analysis
Discussion
ATC is a rare but highly aggressive malignancy with a poor prognosis, and its rapid progression and limited therapeutic effect cause much concern [23, 24]. However, although several studies attempted to investigate molecular mechanisms of ATC by high-throughput microarrays, these studies were restricted to the identification of DEGs by analyzing fewer gene sets and failed to perform comprehensively associated survival analysis [9, 10, 16]. In this study, we identified 6 genes as hub genes for the occurrence and development of ATC and constructed a reliable and valid independent prognostic signature for ATC consisting of 6 survival-related genes that exhibited robust capabilities in predicting survival outcomes in ATC patients in both training and validation cohorts. Moreover, these 12 genes were also verified in CAL62 cells and Nthy-ori 3–1 cells by RT-qPCR. These findings offer further insight into tumorigenesis and development of ATC and provide novel prognostic biomarkers for ATC patients.
The WGCNA algorithm constructs a network based on the correlation between genes whereas the PPI network is based on protein interactions [25], both of which constitute an efficient and comprehensive strategy frequently employed for the exploration of disease biomarkers. For instance, previous studies identified several biomarkers that were closely related to the treatment and prognosis of lung squamous cancer using this method [8]. Furthermore, it was also indicated that, based on this approach, CSTA might play important roles in the progression of head and neck cancer and serve as a potential biomarker for future diagnosis and treatment [26]. In this respect, the present study for the first time combined WGCNA and PPI network analysis to excavate hub genes associated with ATC occurrence, progression, and prognosis, which will contribute to further improving the individualized prevention and treatment for ATC patients.
In this study, through a series of rigorous screening, 6 genes were finally certified as diagnostic biomarkers (namely KIF2C, PBK, TOP2A, CDK1, KIF20A, and ASPM). Among these diagnostic biomarkers, kinesin family member 2C (KIF2C), a mitotic centromere associated with kinesin, is an essential regulator of the cell cycle and its aberrant expression can cause abnormal spindle formation and failed chromosome segregation, thereby promoting tumor progression [27]. In the present study, we further revealed that this gene was significantly overexpressed in ATC, which was in accordance with previous studies indicating that KIF2C was closely associated with the occurrence and development of ATC [10]. Moreover, the other 2 genes (CDK1, and KIF20A) were also previously suggested to be closely correlated with the genesis and development of ATC and to exhibit tumor-suppressive effects [28, 29], which was also confirmed in the present study. Therefore, these 2 genes might also be key genes in the diagnosis of ATC though further studies are required. Additionally, PBK and TOP2A were both found to promote tumor cell proliferation and play positive regulatory roles in chromosome separation [30–32]. However, although our bioinformatics analysis was consistent with previous work, current RT-qPCR results showed no significant difference in the expression of PBK and TOP2A between ATC and normal samples. These results may be attributed to the relatively low incidence of ATC as well as the influence of ethnicity. Therefore, further study is still needed to elucidate the pivotal role of PBK and TOP2A during ATC progression.
Additionally, as a positive regulator of the canonical Wnt signaling pathway [33, 34], ASPM was found to contribute to proliferation and invasive properties of cancer cells and to be overexpressed in many cancers, such as lung carcinoma, gastric carcinoma, and hepatocellular carcinoma [35–37]. Notably, in our study, ASPM also showed a significant expression in ATC compared with normal tissues and was closely associated with ATC occurrence and development. Furthermore, a recent study was highly concordant with our results that ASPM promoted the progression of ATC [38]. Therefore, ASPM is expected to be a prospective target for the diagnosis and treatment of ATC. Regrettably, we were unable to extract ASPM related information between ATC and normal tissues in the HPA database mainly due to the rarity of ATC and a lack of more pathological tissue, which should be validated in future studies.
To effectively improve the outlook and therapy of ATC patients, the current research presented an innovative mark showcasing 6 predictive genes (GPSM2, FGF5, ASXL3, CYP4B1, CLMP, and DUXAP9), which demonstrated excellent prognostic ability in ATC patients and was also confirmed as an autonomous prognosticator of ATC outlook. Additionally, a diagram that merged this mark with clinical features was developed to precisely foretell 1-, 3- and 5-year survival probabilities of ATC patients, which could offer clinicians a handy and hopeful approach to evaluate the outlook of ATC patients and carry out individualized therapy.
Among these prognostic genes, the expressions of GPSM2, FGF5, ASXL3, CLMP, and DUXAP9 were significantly higher in high-risk group compared in low-risk group and were related to poor survival outcome, whereas CYP4B1 was overexpressed in low-risk group and might be a protective factor for ATC. GPSM2, a key member of the G-protein-signaling modulator family, has been identified in previous studies as a potential drug target and prognostic biomarker in hepatocellular carcinoma [39], lung adenocarcinoma [40], and liver cancer [41]. Our study further demonstrated that GPSM2 exerts a significant oncogenic effect and acts as an unfavorable factor in ATC. Furthermore, FGF5, part of the fibroblast growth factor family, promotes angiogenesis, tumor cell transformation, and proliferation [42]. Recent findings have highlighted the oncological significance of FGF5, noting its overexpression in various solid tumors and its association with poor survival outcomes [43]. This aligns with our results, indicating that FGF5 may be developed as a prognostic biomarker for ATC. Additionally, ASXL3, CLMP, DUXAP9, and CYP4B1 have previously been reported to be involved in various cancers [44–46] and were included in our prognostic signature, further validating their reliability as biomarkers.
Although we excavated the primary genes connected with the development, progression, and prediction of ATC as thoroughly as possible, some limitations need to be mentioned. Firstly, in the absence of suitable survival data for ATC, we used the clinical data from TCGA thyroid cancer to characterize the impact of identified key genes on survival. Secondly, the potential mechanisms of the current 12 primary genes needed to be further investigated to gain a better understanding of the tumorigenesis and growth of ATC. Thirdly, even though the suggested gene signature model was authenticated by cell experiments, additional clinical trials were still needed to validate its clinical usefulness. Finally, this study lacks validation in an independent ATC cohort to enhance the generalizability of the results as well as comparison of ATC with other types of thyroid cancer (papillary, follicular, and medullary). These constraints should be surmounted in forthcoming fundamental and clinical investigations.
In summary, we comprehensively identified 6 hub genes intricately linked to the onset and progression of ATC. We subsequently formulated a prognostic signature comprising these 6 genes, which was validated as an independent predictor of OS among ATC patients. This pivotal study offers profound insights into the mechanisms underlying ATC development and prognosis and presents novel and robust biomarkers for its precise prevention and treatment. It is expected to provide molecular targets for clinical diagnosis and drug design.
Acknowledgements
All authors sincerely acknowledge the contributions from the TCGA, GEO and HPA projects.
Author contributions
GGH conceived and designed the study and obtained the funding. The investigation, data collection, analysis, and writing were performed by KQ, JRW, ZZ, and JDZ. QF, PW, YG, LC, and QWZ participated in the interpretation of data and verification. KQ and TZ contributed to visualization. GGH and JDZ revised the manuscript. All gave final approval and agreed to be accountable for all aspects of work ensuring integrity and accuracy.
Funding
This work was supported by grants from the Applied Basic Research Program Yunnan Province of China (Joint Special Project of Kunming Medical University) (No. 202101AY070001-300) and the Young and Middle-aged Academic and Technical Leaders Reserve Talent Project of Yunnan Province (No. 202405AC350037).
Data availability
The original contributions presented in the study are included in the article/Supplementary Material; the raw data that support the findings of this study are available at the following URL: https://www.ncbi.nlm.nih.gov/geo/https://portal.gdc.cancer.gov/
Declarations
Ethics approval and consent to participate
Our study is based on open-source data, thus no ethical issues or other conflicts of interest exist.
Competing interests
The authors declare no competing interests.
Publisher's Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
References
1. Kitahara, CM; Sosa, JA. The changing incidence of thyroid cancer. Nat Rev Endocrinol; 2016; 12,
2. Li, X; Zhang, S; Zhang, Q et al. Diagnosis of thyroid cancer using deep convolutional neural network models applied to sonographic images: a retrospective, multicohort, diagnostic study. Lancet Oncol; 2019; 20,
3. Pan, Z; Fang, Q; Li, L et al. HN1 promotes tumor growth and metastasis of anaplastic thyroid carcinoma by interacting with STMN1. Cancer Lett; 2021; 31,
4. Li, L; Zhu, M; Huang, H; Wu, J; Meng, D. Identification of hub genes in anaplastic thyroid carcinoma: evidence from bioinformatics analysis. Technol Cancer Res Treat; 2022; 19,
5. Pan, Z; Li, L; Fang, Q et al. Integrated bioinformatics analysis of master regulators in anaplastic thyroid carcinoma. Biomed Res Int; 2019; 2, pp. 1-13.[COI: 1:CAS:528:DC%2BC1MXhvVaqt7rI] [DOI: https://dx.doi.org/10.1155/2019/9734576]
6. Jiang, L; Zhang, S; An, N; Chai, G; Ye, C. ASPM promotes the progression of anaplastic thyroid carcinomas by regulating the Wnt/β-catenin signaling pathway. Int J Endocrinol; 2022; 2022, 5316102.[COI: 1:CAS:528:DC%2BB38XhsF2js7vF] [DOI: https://dx.doi.org/10.1155/2022/5316102] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/35387319][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC8977346]
7. Seok, JY; Astvatsaturyan, K; Peralta-Venturina, M; Lai, J; Fan, X. TROP-2, 5hmC, and IDH1 expression in anaplastic thyroid carcinoma. Int J Surg Pathol; 2021; 29,
8. Gao, M; Kong, W; Huang, Z; Xie, Z. Identification of key genes related to lung squamous cell carcinoma using bioinformatics analysis. Int J Mol Sci; 2020; 21,
9. Ma, Y; Cang, S; Li, G et al. Integrated analysis of transcriptome data revealed MMP3 and MMP13 as critical genes in anaplastic thyroid cancer progression. J Cell Physiol; 2019; 234,
10. Wang, S; Wu, J; Guo, C et al. Identification and validation of novel genes in anaplastic thyroid carcinoma via bioinformatics analysis. Cancer Manag Res; 2020; 12, pp. 9787-9799.[COI: 1:CAS:528:DC%2BB3MXjs1Wgs78%3D] [DOI: https://dx.doi.org/10.2147/CMAR.S250792] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/33116838][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7550107]
11. Pan, Z; Xu, T; Bao, L et al. CREB3L1 promotes tumor growth and metastasis of anaplastic thyroid carcinoma by remodeling the tumor microenvironment. Mol Cancer; 2022; 21,
12. Wang, Y; Chen, L; Ju, L et al. Novel biomarkers associated with progression and prognosis of bladder cancer identified by co-expression analysis. Front Oncol; 2019; 9, 1030. [DOI: https://dx.doi.org/10.3389/fonc.2019.01030] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/31681575][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6799077]
13. Tian, Z; He, W; Tang, J et al. Identification of important modules and biomarkers in breast cancer based on WGCNA. Onco Targets Ther; 2020; 13, pp. 6805-6817.[COI: 1:CAS:528:DC%2BB3cXisVOqt7zO] [DOI: https://dx.doi.org/10.2147/OTT.S258439] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32764968][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7367932]
14. Liu, J; Zhou, S; Li, S et al. Eleven genes associated with progression and prognosis of endometrial cancer (EC) identified by comprehensive bioinformatics analysis. Cancer Cell Int; 2019; 19, pp. 1-17.[COI: 1:CAS:528:DC%2BC1MXhtFKktbzO] [DOI: https://dx.doi.org/10.1186/s12935-019-0859-1]
15. Wu, J; Fang, X; Xia, X. Identification of key genes and pathways associated with endometriosis by weighted gene co-expression network analysis. Int J Med Sci; 2021; 18,
16. Li, HN; Du, YY; Xu, T et al. TMEM158 may serve as a diagnostic biomarker for anaplastic thyroid carcinoma: an integrated bioinformatic analysis. Curr Med Sci; 2020; 40,
17. Wang, DP; Tang, XZ; Liang, QK; Zeng, XJ; Yang, JB; Xu, J. Overexpression of long noncoding RNA SLC26A4-AS1 inhibits the epithelial-mesenchymal transition via the MAPK pathway in papillary thyroid carcinoma. J Cell Physiol; 2020; 235,
18. Tu, Y; Fan, G; Xi, H et al. Identification of candidate aberrantly methylated and differentially expressed genes in thyroid cancer. J Cell Biochem; 2018; 119,
19. Liu, L; He, C; Zhou, Q; Wang, G; Lv, Z; Liu, J. Identification of key genes and pathways of thyroid cancer by integrated bioinformatics analysis. J Cell Physiol; 2019; 234,
20. Du, L; Zhao, Q; Li, J; Wang, M; Qiao, H. Expression of colorectal neoplasia differentially expressed in anaplastic thyroid carcinoma and its effect on cancer cell proliferation. Ann Transl Med.; 2022; 10,
21. Wu, T; Hu, E; Xu, S et al. clusterProfiler 4.0: a universal enrichment tool for interpreting omics data. Innovation (Camb).; 2021; 2,
22. Tibshirani, R. The lasso method for variable selection in the Cox model. Stat Med; 1997; 16,
23. Jannin, A; Escande, A; Al Ghuzlan, A et al. Anaplastic thyroid carcinoma: an update. Cancers (Basel); 2022; 14,
24. Pan, Z; Li, L; Qian, Y et al. The differences of regulatory networks between papillary and anaplastic thyroid carcinoma: an integrative transcriptomics study. Cancer Biol Ther; 2020; 21,
25. Kong, L; Chen, J; Ji, X et al. Alcoholic fatty liver disease inhibited the co-expression of Fmo5 and PPARα to activate the NF-κB signaling pathway, thereby reducing liver injury via inducing gut microbiota disturbance. J Exp Clin Cancer Res; 2021; 40,
26. Li, CY; Cai, JH; Tsai, JJP; Wang, CCN. Identification of hub genes associated with development of head and neck squamous cell carcinoma by integrated bioinformatics analysis. Front Oncol.; 2020; 10, 681. [DOI: https://dx.doi.org/10.3389/fonc.2020.00681] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32528874][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7258718]
27. Zhang, X; Li, Y; Hu, P; Xu, L; Qiu, H. KIF2C is a biomarker correlated with prognosis and immunosuppressive microenvironment in human tumors. Front Genet; 2022; 13,[COI: 1:CAS:528:DC%2BB38XhsVyrtrzO] [DOI: https://dx.doi.org/10.3389/fgene.2022.891408] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/35685442][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC9171145]
28. Gong, Y; Xu, F; Deng, L; Peng, L. Recognition of key genes in human anaplastic thyroid cancer via the weighing gene coexpression network. Biomed Res Int; 2022; 2022, 2244228.[COI: 1:CAS:528:DC%2BB38XhvFCrur7L] [DOI: https://dx.doi.org/10.1155/2022/2244228] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/35782055][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC9247818]
29. Liang, B; Zhou, Y; Jiao, J et al. Integrated analysis of transcriptome data revealed AURKA and KIF20A as critical genes in medulloblastoma progression. Front Oncol; 2022; 12,[COI: 1:CAS:528:DC%2BB3sXhtlGqsL3K] [DOI: https://dx.doi.org/10.3389/fonc.2022.875521] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/35574421][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC9092218]
30. Wen, H; Chen, Z; Li, M et al. An integrative pan-cancer analysis of PBK in human tumors. Front Mol Biosci; 2021; 8,[COI: 1:CAS:528:DC%2BB38XhtVCktLs%3D] [DOI: https://dx.doi.org/10.3389/fmolb.2021.755911] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34859049][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC8631476]
31. Uusküla-Reimand, L; Wilson, MD. Untangling the roles of TOP2A and TOP2B in transcription and cancer. Sci Adv.; 2022; 8,
32. Dong, C; Fan, W; Fang, S. PBK as a potential biomarker associated with prognosis of glioblastoma. J Mol Neurosci; 2020; 70,
33. Yi, X; Wan, Y; Cao, W; Peng, K; Li, X; Liao, W. Identification of four novel prognostic biomarkers and construction of two nomograms in adrenocortical carcinoma: a multi-omics data study via bioinformatics and machine learning methods. Front Mol Biosci; 2022; 9,[COI: 1:CAS:528:DC%2BB38XhsFyqurjJ] [DOI: https://dx.doi.org/10.3389/fmolb.2022.878073] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/35693556][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC9174903]
34. Hou, S; Xu, H; Liu, S et al. Integrated bioinformatics analysis identifies a new stemness index-related survival model for prognostic prediction in lung adenocarcinoma. Front Genet; 2022; 13,[COI: 1:CAS:528:DC%2BB38Xht12qs7fP] [DOI: https://dx.doi.org/10.3389/fgene.2022.860268] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/35464867][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC9026767]
35. Hsu, CC; Liao, WY; Chang, KY. Correction to: a multi-mode Wnt- and stemness-regulatory module dictated by FOXM1 and ASPM isoform I in gastric cancer. Gastric Cancer; 2021; 24,
36. Deng, T; Liu, Y; Zhuang, J; Tang, Y; Huo, Q. ASPM is a prognostic biomarker and correlates with immune infiltration in kidney renal clear cell carcinoma and liver hepatocellular carcinoma. Front Oncol; 2022; 12,[COI: 1:CAS:528:DC%2BB3sXhtlGqs7fO] [DOI: https://dx.doi.org/10.3389/fonc.2022.632042] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/35515103][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC9065448]
37. Wu, B; Hu, C; Kong, L. ASPM combined with KIF11 promotes the malignant progression of hepatocellular carcinoma via the Wnt/β-catenin signaling pathway. Exp Ther Med; 2021; 22,
38. Zhang, J; Lou, W. A key mRNA-miRNA-lncRNA competing endogenous RNA triple sub-network linked to diagnosis and prognosis of hepatocellular carcinoma. Front Oncol; 2020; 10, 340. [DOI: https://dx.doi.org/10.3389/fonc.2020.00340] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32257949][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7092636]
39. Deng, M; Liu, B; Zhang, Z et al. Loss of G-protein-signaling modulator 2 accelerates proliferation of lung adenocarcinoma via EGFR signaling pathway. Int J Biochem Cell Biol; 2020; 122,[COI: 1:CAS:528:DC%2BB3cXjslynsLc%3D] [DOI: https://dx.doi.org/10.1016/j.biocel.2020.105716] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32058048]
40. Yang, D; Ji, F; Li, Y; Jiao, Y; Fang, X. GPSM2 serves as an independent prognostic biomarker for liver cancer survival. Technol Cancer Res Treat; 2020; 19, 1533033820945817.[COI: 1:CAS:528:DC%2BB3MXks1Kiurw%3D] [DOI: https://dx.doi.org/10.1177/1533033820945817] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32812493][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7440740]
41. Catenacci, DVT; Rasco, D; Lee, J et al. Phase I escalation and expansion study of bemarituzumab (FPA144) in patients with advanced solid tumors and FGFR2b-selected gastroesophageal adenocarcinoma. J Clin Oncol; 2020; 38,
42. Zhao, T; Qian, K; Zhang, Y. High expression of FGF5 is an independent prognostic factor for poor overall survival and relapse-free survival in lung adenocarcinoma. J Comput Biol; 2020; 27,
43. Szczepanski, AP; Zhao, Z; Sosnowski, T; Goo, YA; Bartom, ET; Wang, L. ASXL3 bridges BRD4 to BAP1 complex and governs enhancer activity in small cell lung cancer. Genome Med; 2020; 12,
44. Wang, X; Duanmu, J; Fu, X; Li, T; Jiang, Q. Analyzing and validating the prognostic value and mechanism of colon cancer immune microenvironment. J Transl Med; 2020; 18,
45. Tan, L; Tang, Y; Li, H. N6-methyladenosine modification of LncRNA DUXAP9 promotes renal cancer cells proliferation and motility by activating the PI3K/AKT signaling pathway. Front Oncol; 2021; 11,[COI: 1:CAS:528:DC%2BB3sXht1elur7O] [DOI: https://dx.doi.org/10.3389/fonc.2021.641833] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34168980][PubMedCentral: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC8217835]
46. Liu, X; Jia, Y; Shi, C. CYP4B1 is a prognostic biomarker and potential therapeutic target in lung adenocarcinoma. PLoS ONE; 2021; 16,
You have requested "on-the-fly" machine translation of selected content from our databases. This functionality is provided solely for your convenience and is in no way intended to replace human translation. Show full disclaimer
Neither ProQuest nor its licensors make any representations or warranties with respect to the translations. The translations are automatically generated "AS IS" and "AS AVAILABLE" and are not retained in our systems. PROQUEST AND ITS LICENSORS SPECIFICALLY DISCLAIM ANY AND ALL EXPRESS OR IMPLIED WARRANTIES, INCLUDING WITHOUT LIMITATION, ANY WARRANTIES FOR AVAILABILITY, ACCURACY, TIMELINESS, COMPLETENESS, NON-INFRINGMENT, MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. Your use of the translations is subject to all use restrictions contained in your Electronic Products License Agreement and by using the translation functionality you agree to forgo any and all claims against ProQuest or its licensors for your use of the translation functionality and any output derived there from. Hide full disclaimer
© The Author(s) 2024. This work is published under http://creativecommons.org/licenses/by-nc-nd/4.0/ (the “License”). Notwithstanding the ProQuest Terms and Conditions, you may use this content in accordance with the terms of the License.
Abstract
Background
Anaplastic thyroid carcinoma (ATC) is a rare but the most aggressive type of thyroid carcinoma. Nevertheless, limited advances were made to reduce mortality and improve survival over the last decades. Therefore, identifying novel diagnostic biomarkers and therapeutic targets for ATC patients is still needed.
Materials and methods
RNA sequencing data and corresponding clinical features were available from GEO and TCGA databases. We integrated WGCNA and PPI network analysis to identify hub genes associated with ATC development, and RT-qPCR was employed for data verification. Univariate and LASSO Cox regression analyses were used to generate prognostic signatures.
Results
Based on PPI and WGCNA, 6 hub genes were identified, namely KIF2C, PBK, TOP2A, CDK1, KIF20A, and ASPM, which play vital roles in ATC development. Subsequently, RT-qPCR experiments showed that most of these genes were significantly upregulated in CAL-62 cells compared to Nthy-ori 3–1 cells. Moreover, a prognostic signature featuring GPSM2, FGF5, ASXL3, CYP4B1, CLMP, and DUXAP9 was generated, which was also verified by RT-qPCR results and proved as an independent predictor of poorer prognosis of ATC. Additionally, a nomogram incorporating the risk score and clinicopathological parameters was further constructed for accurate prediction of 1-, 3- and 5-year survival probabilities of ATC.
Conclusions
Our study identified 6 key genes critical to ATC development and constructed a prognostic signature. These findings provide reliable biomarkers and a relatively comprehensive tumorigenesis profile of ATC, which may inform future strategies for clinical diagnosis and pharmaceutical design.
You have requested "on-the-fly" machine translation of selected content from our databases. This functionality is provided solely for your convenience and is in no way intended to replace human translation. Show full disclaimer
Neither ProQuest nor its licensors make any representations or warranties with respect to the translations. The translations are automatically generated "AS IS" and "AS AVAILABLE" and are not retained in our systems. PROQUEST AND ITS LICENSORS SPECIFICALLY DISCLAIM ANY AND ALL EXPRESS OR IMPLIED WARRANTIES, INCLUDING WITHOUT LIMITATION, ANY WARRANTIES FOR AVAILABILITY, ACCURACY, TIMELINESS, COMPLETENESS, NON-INFRINGMENT, MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. Your use of the translations is subject to all use restrictions contained in your Electronic Products License Agreement and by using the translation functionality you agree to forgo any and all claims against ProQuest or its licensors for your use of the translation functionality and any output derived there from. Hide full disclaimer
Details
1 920th Hospitalof Joint Logistics Support Force of People’s Liberation Army, Department of Clinical Pharmacy, Kunming, China; Dali University, College of Pharmacy, Dali, China (GRID:grid.440682.c) (ISNI:0000 0001 1866 919X)
2 920th Hospital of JointLogistics Support Force of People’s Liberation Army, Department of Pathology, Kunming, China (GRID:grid.440682.c)
3 Kunming Tongren Hospital, Oncology Department, Kunming, China (GRID:grid.490275.d)
4 920th Hospitalof Joint Logistics Support Force of People’s Liberation Army, Department of Clinical Pharmacy, Kunming, China (GRID:grid.490275.d); Dali University, College of Pharmacy, Dali, China (GRID:grid.440682.c) (ISNI:0000 0001 1866 919X)
5 920th Hospitalof Joint Logistics Support Force of People’s Liberation Army, Department of Clinical Pharmacy, Kunming, China (GRID:grid.440682.c)
6 920th Hospitalof Joint Logistics Support Force of People’s Liberation Army, Department of Clinical Pharmacy, Kunming, China (GRID:grid.440682.c); Dali University, College of Pharmacy, Dali, China (GRID:grid.440682.c) (ISNI:0000 0001 1866 919X)
7 Dali University, College of Pharmacy, Dali, China (GRID:grid.440682.c) (ISNI:0000 0001 1866 919X)
8 The 306th Hospital of People’s LiberationArmy (PLA), Medical Engineering Section, Beijing, China (GRID:grid.413138.c)
9 920th Hospitalof Joint Logistics Support Force of People’s Liberation Army, Department of Clinical Pharmacy, Kunming, China (GRID:grid.413138.c)