Development and validation of a necroptosis-related gene signature for predicting prognosis and immune infiltration in papillary thyroid cancer

Introduction

Thyroid cancer (THCA) is a malignancy originating from thyroid follicular epithelial cells or parafollicular C cells. According to joint data from the American Cancer Society and the International Agency for Research on Cancer (IARC), it ranked as the seventh most common malignancy globally and the fifth most prevalent cancer among females worldwide as of 2022 (1,2). In 2022, China’s World Age-Standardized Incidence Rate (WASIR) of THCA reached 24.64×10⁵, ranking as the third most prevalent malignancy (3). This represents a significant upward trend compared to 2016 statistics (4), which reported a THCA WASIR of 10.37×10⁵. Simultaneously, THCA mortality rates are also on the rise (5). Papillary thyroid cancer (PTC), which accounts for 80–85% of THCA cases, is the most common form (6). Although PTC typically portends good outcomes, a pronounced survival disparity exists between China (84% 5-year survival) and developed countries (99%) (7,8). Despite the continuous refinement of treatment strategies for PTC, some patients still face challenges, including local invasion, iodine refractoriness, metastasis, and recurrence, which contribute to poor prognosis. These difficulties are compounded by the heterogeneity of PTC, which exhibits diverse clinical, pathological, and molecular characteristics (9). The natural history of these aggressive variants remains incompletely characterized, frequently resulting in suboptimal therapeutic outcomes. This critical unmet need underscores the necessity to establish a robust risk stratification model for prognosis prediction and precision management across risk-differentiated cohorts. The prognosis of patients with THCA is influenced by multiple factors, including clinical features and genetic factors. Previous studies have identified age, tumor-nodes-metastasis (TNM) stage, extrathyroidal extension, and lymph node metastasis as independent prognostic factors for THCA (10,11). BRAF and TERT mutations serve as molecular drivers that independently correlate with diminished survival outcomes in THCA patients (12). Despite thyroglobulin’s (Tg) established role as a principal prognostic indicator in PTC, its clinical utility is significantly confounded by concurrent Tg antibody (TgAb) seropositivity (13). Therefore, novel biomarker discovery is critical to enhance PTC diagnostic reliability, therapeutic targeting, and prognostic modeling.

As a genetically programmed cell death pathway, apoptotic dysregulation constitutes a fundamental hallmark of oncogenesis (14). Before necroptosis was characterized as a programmed necrotic death pathway, necrosis was mechanistically distinguished from apoptosis (15). However, necroptosis shares mechanistic similarities with apoptosis and morphological characteristics with necrosis (16). Unlike the classical apoptotic pathway, which relies on caspase activation, necroptosis serves as an alternative when caspases are absent or inhibited (17). Multiple innate immune signaling cascades, particularly TNFR1 and Toll-like receptor (TLR)-dependent pathways, are capable of initiating necroptotic cell death (16). These pathways lead to the phosphorylation and activation of receptor-interacting protein kinase 1 (RIPK1) and RIPK3, which initiate a phosphorylation cascade converting mixed lineage kinase domain-like (MLKL) to its membrane-targeting conformation. Oligomerized MLKL compromises lipid bilayer integrity, facilitating damage-associated molecular patterns (DAMPs) release and effecting necroptotic demise (18).

Necroptotic signaling networks significantly impact cancer development, immunity, and metastasis. This death pathway demonstrates context-dependent duality: acting as a tumor constraint mechanism during apoptotic impairment (19), while its molecular executors (singly or synergistically) fuel malignant progression (20). Nevertheless, necroptosis-related genes (NRGs) in PTC represent an understudied area, and their mechanistic roles in disease evolution await systematic investigation.

Therefore, this study employed a bioinformatics approach to identify NRGs associated with the prognosis of patients with PTC, examining how these prognostic-related genes (PRGs) impact the PTC immune microenvironment. Additionally, leveraging NRGs to identify candidate therapeutic drugs, delivering a molecular roadmap to guide precision diagnosis and targeted therapy in PTC. We present this article in accordance with the TRIPOD reporting checklist (available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-2124/rc).

Methods

Data source

Transcriptomic profiles of 510 PTC and 58 normal thyroid tissues were acquired from The Cancer Genome Atlas (TCGA) (accession: https://portal.gdc.cancer.gov/). Using the caret package (v6.0.93), 509 PTC cases with complete survival records were partitioned via createDataPartition into training (n=255) and testing sets (n=254) at 1:1 ratio. The Gene Expression Omnibus (GEO) dataset GSE29265 (http://www.ncbi.nlm.nih.gov/geo/) provided 20 matched normal-PTC samples. Additionally, 636 differentially expressed NRGs (DENRGs) with |logfold change (logFC)| >1 were curated from GeneCards (http://www.genecards.org/). The study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments.

DENRGs acquisition and functional enrichment

For the TCGA-PTC dataset, differential expression analysis of DENRGs between PTC and normal samples was performed using the “limma” package (version 3.46.0) (21), comparing tumor versus normal tissues. Differential expression criteria were set at |log₂FC| >0.5 and adjusted P<0.05. Volcano plots were visualized using the “ggplot2” package (version 3.3.6), and heat maps were created using “pheatmap” (version 1.0.12) (22). Functional annotation of significant DENRGs included Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment via clusterProfiler (v3.18.1) and Visualization of enrichment results as bubble plots using ggplot2 (v3.3.6). All analyses were implemented in the R statistical environment (v4.1.0).

PTC prognostic model construction and validation

Initially, Univariate Cox regression identified survival-associated DENRGs in the training set. Least absolute shrinkage and selection operator (LASSO) regularization via glmnet (v4.1-4) (23) refined gene selection, minimizing false-positive inclusions to enhance model robustness. LASSO-generated coefficients were applied to corresponding gene expression values to compute per-patient risk scores. Calculate the risk score according to the following formula:

RiskScore=exp(i1)×coef(i1)+exp(i2)×coef(i2)+…+exp(in)×coef(in)[1]

Among them, “exp” denoted the gene expression level, and “coef” stood for the coefficient.

Kaplan-Meier survival curves generated via the survminer package (v0.4.9) (24) delineated survival probability disparities between risk-stratified cohorts, with statistical significance assessed through log-rank testing. Model discrimination capacity was quantified using time-dependent receiver operating characteristic (ROC) analysis implemented in survivalROC (v1.0.3), evaluating area under the curve (AUC) metrics across both training and testing datasets to confirm predictive robustness. Diagnostic performance assessment integrated two independent cohorts (TCGA and GSE29625), where ROC curves constructed with pROC (v2.3.0) measured classification accuracy of the prognostic signature. Complementary differential expression analysis examined transcriptional profiles of model genes in malignant versus normal tissues within these datasets, visualized through box plots created with ggplot2 (v3.3.6) (25) to illustrate pathological dysregulation patterns. Furthermore, consolidated ROC evaluation across TCGA and GSE29625 using pROC specifically interrogated the model’s precision in distinguishing PTC from non-malignant states. Furthermore, the expression of PRGs in PTC and adjacent normal tissues in the TCGA and GSE29625 was detected using a nonparametric test. Meanwhile, in the disease samples of the TCGA and GSE29625 the expression differences of the genes associated with prognosis were compared between the high- and low-risk groups (P<0.05).

Evaluation of independent prognostic factors in patients with PTC

In the present study, clinicopathological variables, including gender, tumor stage, TNM categories, and risk score, were integrated into multivariate models using 137 training cohort cases with complete case records. Multivariable Cox proportional hazards regression identified independent prognostic determinants for PTC (significance threshold: P<0.05). These significant predictors were subsequently incorporated into a nomogram constructed with the rms package (v6.5-0), with model calibration assessed through observed-versus-predicted outcome curves.

Functional enrichment profiling of PRGs

Gene set enrichment analysis (GSEA) was performed on TCGA PTC samples stratified by prognostic gene expression levels (high versus low) to interrogate GO terms and KEGG pathways. Significant enrichments were defined by |normalized enrichment score (NES)| >1 and adjusted P<0.05.

Tumor mutational landscape and immune microenvironment characterization

Immune cell infiltration in the training cohort was quantified via xCell algorithm, with risk group disparities assessed by Wilcoxon signed-rank test. Somatic mutation patterns were decoded using maftools (v2.6.05) (26), including F-test evaluation of mutation co-occurrence/exclusivity. Tumor mutation burden (TMB) stratification informed Kaplan-Meier survival comparisons (survminer v0.4.9). Additionally, single-sample GSEA (ssGSEA) enumerated infiltration densities of 28 immune cell subsets in PTC specimens.

Immunotherapy responsiveness evaluation

Differential expression of immune checkpoint molecules was analyzed across risk strata. Immunophenoscores (IPS) were acquired from The Cancer Imaging Archive (https://tcia.at/) to model PD1/CTLA4 therapeutic efficacy disparities. Immune subtype distributions were resolved using ImmuneSubtypeClassifier.

Prospective therapeutic drug screening

Targeted pharmacological agents linked to PRGs were discovered by integrating data from the Genomics of Drug Sensitivity in Cancer (GDSC; https://www.cancerrxgene.org/) and the Drug-Gene Interaction Database (DGIdb; http://www.dgidb.org). The “oncoPredict” tool (version 0.2) was used to compute the half-maximal inhibitory concentration (IC₅₀) values for a variety of drugs across patients from both risk categories.

Experimental validation via quantitative reverse transcription polymerase chain reaction (qRT-PCR) and immunoblotting

Tissue specimens from tumor and adjacent normal regions of 10 classical PTC patients underwent molecular validation. Total RNA isolation from peripheral blood mononuclear cells (PBMCs) employed TRIzol reagent (Ambion, Austin, TX, USA) per manufacturer’s protocol, followed by complementary DNA (cDNA) synthesis using SureScript First Strand cDNA Synthesis Kit (Servicebio, Wuhan, China). Target-specific primers enabled qRT-PCR quantification on BIO-RAD systems with GAPDH normalization, applying the ∆Ct calculation method. Primer (Tsingke Biotechnology, Xi’an, China) sequences were designed as follows:

• FNDC4: forward, 5'-CTGACCGGCTACCTTCAAACA-3'; reverse, 5'-GCCTTCCCTGAGGACTCTGTT-3';
• TYRO3: forward, 5'-CAGCCACCAGCAGCAGTATC-3'; reverse, 5'-GGATGGCAAGGTCCAGAAGT-3’;
• CXCL5: forward, 5'-TGCGTTGCGTTTGTTTACAG-3'; reverse, 5'-TCTTCAGGGAGGCTACCACTT-3';
• GAPDH (internal control): forward, 5'-GAAGGTCGGAGTCAACGGATTT-3'; reverse, 5'-ATGGGTGGAATCATATTGGAAC-3'.

For protein immunodetection, radioimmunoprecipitation assay (RIPA) lysates (Solarbio, Beijin, China) were resolved on 10% sodium dodecyl sulfate-polyacrylamide gel electrophoresis (SDS-PAGE) gels and electrotransferred to polyvinylidene fluoride (PVDF) membranes (Millipore, Burlington, VT, USA). Immunoblotting utilized: primary antibodies: β-actin (1:25,000; Proteintech, Rosemont, IL, USA; 66009-1-Ig), CXCL5 (1:1,000; Abcam, Boston, MA, USA; ab126763), TYRO3 (1:1,000; CST, Danvers, MA, USA; #5585), FNDC4 (1:1,000; NOVUS, Littleton, CO, USA; NBP1-59690); secondary antibodies: goat anti-mouse immunoglobulin G (IgG)-horseradish peroxidase (HRP) (Servicebio; GB23301), goat anti-rabbit IgG-HRP (Invitrogen, Carlsbad, CA, USA; AB_228338). All procedures complied with institutional ethical guidelines (Ethics Committee approved).

Statistical analysis

R software (v4.2.2) was used to conduct all analyses. In the bioinformatics analysis, the Wilcoxon rank-sum test was used to examine the differences between the two groups, while the t-test was used in qRT-PCR to compare the differences between the two groups. A P value less than 0.05 was regarded as statistically significant.

Results

DENRGs profiling and functional annotation

Differential analysis of TCGA datasets discerned 32 DENRGs exhibiting significant dysregulation in PTC versus normal tissues (adjusted P<0.05; |log₂FC| >0.5), comprising 15 upregulated and 17 downregulated transcripts (Figure 1A,1B; Table S1). Functional annotation revealed 437 GO terms predominantly enriched in TNF-mediated signaling cascades, cellular responses to inflammatory stimuli, and extrinsic apoptotic pathways (Figure 1C). Parallel KEGG pathway analysis identified 58 significant entries, with necroptosis, nuclear factor-κB (NF-κB) signaling, and viral infection-associated pathogenesis [e.g., human immunodeficiency virus type 1 (HIV-1)/hepatitis C virus (HCV) infection, alcoholic liver disease] representing core mechanisms (Figure 1D).

Figure 1 Identification and enrichment analyses of DENRGs. (A) Volcano plot illustrating the DENRGs, with red dots indicating upregulated DENRGs, it reflected the change multiples of the expression levels of DENRGs, as well as the significance of the differences in gene expression. (B) Heatmap comparing the expression profiles of 32 DENRGs between normal and PTC tissues (red: high expression; blue: low expression). (C) GO analysis showing the top 10 enriched biological terms for DENRGs. (D) KEGG analysis highlighting the top 10 significantly associated pathways. BP, biological process; DENRG, differentially expressed necroptosis-related gene; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; MF, molecular function; PTC, papillary thyroid cancer.

Prognostic model development for patients with PTC

Univariate Cox regression screening in the training cohort nominated CXCL5, TYRO3, and FNDC4 as survival-predictive genes (P<0.05; Figure 2A). LASSO regularization refined these candidates into a tri-gene prognostic signature, with risk scores computed via expression-coefficient integration (Figure 2B). Stratification at the optimal cutoff (−7.06) demarcated distinct prognostic cohorts: high-risk patients demonstrated precipitously shortened survival durations in risk distribution curves (Figure 2C) and significantly inferior outcomes in Kaplan-Meier analysis (log-rank P=0.02; Figure 2D). Using ROC curve methodology, we objectively validated the model’s predictive capacity, with time-dependent AUC analyses yielding values >0.70 for all examined survival durations (1-, 3-, and 5-year), as visualized in Figure 2E.

Figure 2 Development of the prognostic signature. (A) Univariate Cox regression assessing OS associations for NRGs. (B) LASSO regression with 10-fold cross-validation to select three optimal PRGs. (C1) Patient stratification by risk score; (C2) survival status distribution (dotted line separates low- and high-risk groups); (C3) heat maps of the expression levels of CXCL5, TYRO3, and FNDC4. The gradient bars represent the expression levels of genes, with red indicating up-regulation and blue indicating down-regulation. (D) Kaplan-Meier survival analysis comparing OS between risk groups. (E) ROC curves evaluating the model’s predictive accuracy over time. AUC, area under the curve; CI, confidence interval; LASSO, least absolute shrinkage and selection operator; NRG, necroptosis-related gene; OS, overall survival; PRG, prognostic-related gene; ROC, receiver operating characteristic.

Verification of the PTC prognostic model

The prognostic risk model exhibited consistent validity across independent cohorts. In the testing set, high-risk patients demonstrated reduced survival times (Figure 3A) and significantly inferior outcomes via Kaplan-Meier analysis (log-rank P=0.048; Figure 3B), supported by time-dependent ROC curves showing robust predictive accuracy for 1-/3-/5-year survival (AUC >0.66; Figure 3C). Validation in the TCGA cohort corroborated these findings: high-risk patients manifested accelerated mortality (Figure 3D), striking survival disadvantage (Kaplan-Meier log-rank P=0.002; Figure 3E), and superior prognostic discrimination (AUC >0.70; Figure 3F). Diagnostic evaluation in both GSE29625 (AUC =0.70) and TCGA (AUC =0.72) datasets further confirmed strong discriminative capacity (Figure 3G,3H). Transcriptomic profiling revealed significant oncogenic dysregulation of PRGs (CXCL5/TYRO3/FNDC4), with marked overexpression in PTC versus normal tissues (P<0.001; Figure 3I,3J). In addition, CXCL5 showed markedly higher expression in the high-risk group, whereas TYRO3 was predominantly upregulated in the low-risk group (Figure S1).

Figure 3 Performance evaluation of the prognostic signature. (A1) Distribution of risk scores; (A2) corresponding survival outcomes; (A3) heatmap displaying expression patterns of the three signature genes. The gradient bars represent the expression levels of genes, with red indicating up-regulation and blue indicating down-regulation. (B) Kaplan-Meier survival curves for the test cohort. (C) Time-dependent ROC analysis (1/3/5 years) in TCGA-PTC. (D) Risk score, survival, and gene expression in the validation set. The gradient bars represent the expression levels of genes, with red indicating up-regulation and blue indicating down-regulation. (E) Kaplan-Meier curves for the validation cohort. (F) ROC analysis (1/3/5 years) in GSE29625. (G,H) ROC curves assessing model performance. (I,J) Expression profiles of PRGs in GSE29625 and TCGA-PTC datasets. **, P<0.01; ***, P<0.001. AUC, area under the curve; PRG, prognostic-related gene; PTC, papillary thyroid cancer; ROC, receiver operating characteristic; TCGA, The Cancer Genome Atlas.

Developing an independent prognostic risk score

Multivariable Cox regression incorporating clinicopathological variables (gender, tumor stage, TNM classification) and the genomic risk score demonstrated the latter’s capacity as an independent prognostic determinant for PTC {hazard ratio (HR) [95% confidence interval (CI)]: 6.805 (1.822–25.422); P<0.001; Figure 4A,4B}. This finding informed the construction of a clinically implementable nomogram integrating the risk score (Figure 4C). Calibration analysis revealed exceptional concordance between predicted and observed survival probabilities at 1-, 3-, and 5-year intervals (Figure 4D-4F), validating the model’s translational utility.

Figure 4 Prognostic evaluation and clinical prediction model. (A) HR assessment of individual variables through Cox proportional hazards regression. (B) Multivariable-adjusted survival analysis incorporating risk score and clinicopathological variables. (C) Clinically applicable nomogram integrating molecular risk stratification with conventional prognostic factors. Model accuracy was verified through (D) 1-year, (E) 3-year, and (F) 5-year survival probability calibration curves comparing predicted versus observed outcomes. CI, confidence interval; HR, hazard ratio; M, metastasis; N, node; OS, overall survival; T, tumor.

Functional landscape of PRGs

Comprehensive functional interrogation revealed that the PRGs orchestrate cell division and cytokine/chemokine-mediated immune responses in PTC. CXCL5 emerged as a multifunctional regulator, with GO annotation implicating its roles in cell cycle checkpoints (G1/S transition), MAPK cascade activation and ossification (adjusted P<0.05, q<0.25; Figure 5A). KEGG analysis further associated CXCL5 with oxidative phosphorylation and cytokine-receptor interactions (Figure 5B).

Figure 5 GSEA of the three PRGs in the TCGA cohort. (A,B) Enriched terms for CXCL5. (C,D) Enriched terms for TYRO3. (E,F) Enriched terms for FNDC4. GSEA, gene set enrichment analysis; PRG, prognostic-related gene; TCGA, The Cancer Genome Atlas.

TYRO3 demonstrated critical involvement in genomic integrity maintenance via GO-annotated cell cycle checkpoints (G1/S, G2/M transitions), DNA damage surveillance, and ubiquitin ligase complexes (Figure 5C). Its KEGG enrichment highlighted roles in endocytic trafficking, axon guidance tight junction dynamics (Figure 5D).

FNDC4 governed messenger RNA (mRNA) metabolism and cell cycle progression, evidenced by GO terms spanning nuclear-transcribed mRNA catabolism (nonsense-mediated decay) and mitotic regulation (Figure 5E). KEGG profiling connected FNDC4 to ribosomal functions, viral-cytokine interplay, chemokine signaling, and lysosomal pathways (Figure 5F). All enrichments met statistical thresholds (adjusted P<0.05, q<0.25).

Somatic mutation analysis and immune cell infiltration landscape in the high- and low-risk groups

Distinct immune microenvironment patterns emerged between risk-stratified cohorts. High-risk patients demonstrated significantly elevated immune scores and reduced stromal scores (P<0.05), though microenvironment scores remained comparable (Figure 6A). Mutational profiling of 123 training-set PTC specimens revealed somatic mutations in 99 samples (80.49%), with BRAF (68%), TTN (8%), and NRAS (5%) representing predominant alterations (Figure 6B). Low-risk counterparts (n=121) exhibited mutations in 87 samples (71.9%), featuring BRAF (47%), NRAS (9%), and MUC16 (7%) as top mutated genes (Figure 6C).

Figure 6 TMB and immune landscape. (A) Immune/stromal/microenvironment scores stratified by risk. (B,C) Mutation spectra of top 20 genes in (B) high- and (C) low-risk groups. Gene mutation co-occurrence/exclusivity patterns in the (D) high-risk and (E) low-risk group. (F) TMB distribution across risk groups. (G) Kaplan-Meier analysis by TMB status. (H) Immune cell infiltration (ssGSEA) differences. (I) Correlation between signature genes and immune cells. The gradient bar represents the degree of correlation, with red indicating positive correlation and blue indicating negative correlation. ns, not significant (P>0.05); *, P<0.05; **, P<0.01; ***, P<0.001. ssGSEA, single-sample gene set enrichment analysis; TMB, tumor mutation burden.

Mutational synergy analysis uncovered significant co-occurrence patterns. For example, TAFB and HMCN1 mutations showed significant co-occurrence in the high-risk group (see Figure 6D for details), while SMC3 and MACF1 along with several other gene mutations demonstrated significant co-occurrence in the low-risk group (see Figure 6E for details). Notably, statistical analysis revealed significant mutual exclusivity between BRAF and NRAS/HRAS mutations in both risk groups (Figure 6D,6E). While TMB showed no intergroup difference (P=0.75), high-TMB patients experienced significantly reduced survival (log-rank P=0.005; Figure 6F,6G).

Immune deconvolution identified enhanced infiltration of 27 immune subsets in high-risk patients, with multiple cell types such as activated/central memory/effector memory CD4⁺ T cells demonstrating significant differences (P<0.05; see Figure 6H for details). PRGs CXCL5/FNDC4 exhibited strong immunomodulatory correlations (|Cor| >0.3, P<0.05), contrasting with TYRO3’s minimal immune interactions (Figure 6I).

FNDC4 and CXCL5 are strongly correlated with immune cells

High-risk patients exhibited significant upregulation of 11 immune checkpoints (P<0.05), including inhibitory molecules (BTLA, CD274/PD-L1, CTLA4) and co-stimulatory receptors (ICOS, TNFRSF9) critical for T-cell exhaustion (Figure 7A). Immunotherapeutic vulnerability assessment revealed enhanced response to combined PD1/CTLA4 blockade in high-risk cohort (ΔIPS >2, P<0.05; Figure 7B), though monotherapy responses remained comparable between groups (Figure 7B). Immune subtyping identified C3 (inflammatory) as the dominant phenotype across risk strata (Figure 7C).

Figure 7 Immune checkpoint and phenotype analysis. (A) Differential expression of immune checkpoints. (B) IPS comparison for PD1/CTLA4 blockade response. (C) Immunophenotype distribution between risk group. ns, not significant (P>0.05); *, P<0.05; **, P<0.01; ***, P<0.001. IPS, immunophenoscore.

Potential therapeutic drug prediction for PTC

Pharmacogenomic screening of the GDSC database identified 198 candidate therapeutic agents for PTC. Among these, compounds AZD5991_1720 and IAP_5620_1428 demonstrated significantly elevated IC₅₀ values in high-risk patients compared to low-risk counterparts (P<0.05; table available at https://cdn.amegroups.cn/static/public/tcr-24-2124-1.xlsx, Figure 8A-8E). Concurrent interrogation of the DGIdb database revealed 18 TYRO3-targeting agents (Figure 8F, Table 1), including clinically relevant compounds such as DEZAPELISIB, TOZASERTIB, AST-487, SNS-314, and PD-0166285 (among others), whereas no compounds associated with CXCL5 or FNDC4 were discerned.

Figure 8 Drug sensitivity and biomarker interactions. (A-E) IC₅₀ values of selected drugs (AZD5991, AZD3759, etc.) in risk groups. (F) Drug-TYRO3 interaction network from DGIdb. DGIdb, Drug-Gene Interaction Database; IC₅₀, half-maximal inhibitory concentration.

PRGs were highly expressed in PTC samples

qRT-PCR (Figure 9A-9C) and western blot (Figure 9D-9F) analyses showed significantly higher mRNA expression levels of PRGs CXCL5, TYRO3, and FNDC4 in PTC tissues compared to normal tissues. Original western blot images are presented in Figure 10. Our protein-level findings corroborate the transcriptomic patterns previously observed in TCGA and GSE29652.

Figure 9 Experimental validation of gene expression. (A-C) RT-PCR and (D-F) WB quantification of CXCL5, TYRO3, and FNDC4 in 20 paired PTCt and ADt. *, P<0.05; **, P<0.01; ***, P<0.001. ADt, adjacent normal tissues; PTC, papillary thyroid cancer; PTCt, PTC tumors; RT-PCR, reverse transcription polymerase chain reaction; WB, western blot.

Figure 10 Protein expression of CXCL5, TYRO3, and FNDC4 detected by western blot analysis. ADt, adjacent normal tissues; PTC, papillary thyroid cancer; PTCt, PTC tumors.

Discussion

PTC accounts for over 80–85% of all thyroid malignancies and presents a rising global incidence (1,4). This trend underscores the necessity for more accurate prognostic tools to guide clinical management of PTC cases. Necroptosis, a regulated necrotic cell death pathway, exerts dual oncological functions: promoting tumor progression (27) while counteracting oncogenesis during apoptotic failure (28). Despite its pathological significance, the diagnostic and prognostic utility of NRGs in PTC remains largely uncharacterized. This study therefore establishes an NRG-derived molecular signature as a novel prognostic biomarker for PTC.

In this research, 32 DENRGs were identified through differential analysis in PTC. Functional enrichment analyses (GO/KEGG) delineated these genes’ pivotal involvement in programmed cell death modalities—particularly necroptosis, apoptosis, and immune-inflammatory-mediated cytolysis—alongside their regulation of TNF signaling, NF-κB transduction cascades, and associated inflammatory pathways. These molecular mechanisms are recognized as critical drivers of PTC pathogenesis, as substantiated by prior mechanistic studies, for example, Chen et al. found that SPRED3 promotes the proliferation of THCA cells through the NF-κB signaling pathway (29). Specifically, Wen et al. demonstrated that MEIS2 acts as a tumor suppressor in THCA, where its overexpression inhibits cellular proliferation and triggers intrinsic apoptosis through NF-κB pathway modulation (30). Complementarily, Luo et al. established that the FOXP4-AS1 constrains PTC tumorigenesis by impeding AKT-mediated proliferative and migratory signaling networks (31). A separate investigation elucidated the diagnostic utility of Galectin-3 and Tg in THCA while interrogating TNF-α’s pathogenetic involvement in disease progression (32). These findings suggest that the identified DENRGs may significantly impact THCA progression by regulating necroptosis. We established a three-gene prognostic signature (CXCL5, TYRO3, and FNDC4) for survival prediction in the TCGA-PTC cohort. Survival analysis based on the risk score revealed significantly better overall survival (OS) in the low-risk group. Additionally, AUC analysis confirmed the excellent predictive value of this model for the prognosis of patients with PTC. Univariate and multivariate Cox regression confirmed the risk score’s status as an independent prognostic determinant (P<0.001), with validation consistency across external cohorts TCGA and GSE29625. Initial univariate Cox regression analysis in our study demonstrated that CXCL5 acts as a risk factor in PTC. CXCL5, a member of the CXC-type chemokine family, plays a critical role in tumorigenesis and cancer progression (33). CXCL5 overexpression is a conserved oncogenic feature across multiple malignancies—notably nasopharyngeal, non-small cell lung cancer (NSCLC), hepatocellular carcinoma (HCC), breast and pancreatic carcinomas—where it drives advanced-stage progression through facilitating local invasion, neutrophil recruitment, and metastatic dissemination (34-38). The precise biological functions of CXCL5 in PTC are still under investigation. The study by Le et al. revealed that CXCL5 exhibits significant correlations with tumor immunology and the tumor microenvironment in PTC (39). Chang et al. showed that downregulating CXCL5 inhibited malignant behavior in THCA cells (40).

Our initial univariate Cox regression analysis revealed that TYRO3 and FNDC4 act as protective factors in PTC. The TAM receptor tyrosine kinase family (TYRO3, AXL, MER) mediates pivotal oncogenic processes through anti-apoptotic signaling and immunoregulation. Although AXL and MER functionalities are well-characterized, TYRO3’s mechanistic contributions remain poorly defined (41). Consistent with our findings, previous studies have reported the protective roles of TYRO3 in inflammatory conditions, brain injury, and renal diseases (42-44). However, in the oncological field, several recent studies have indicated that TYRO3 exerts tumor-promoting effects. For instance, the study by Hara et al. revealed that patients with TYRO3 overexpression in pancreatic cancer exhibited a poorer prognosis (45). Xiong et al. demonstrated that inhibition of TYRO3 suppresses the progression of liver cancer through the ERK signaling pathway (46). We hypothesize that the underlying reason may reside in the dual nature of TYRO3, a member of the TAM family of receptor tyrosine kinases, which exhibits context-dependent roles in tumorigenesis. Its function may manifest either protective or pro-tumorigenic effects depending on tumor type, immune microenvironment, and molecular regulatory networks (47). Furthermore, TYRO3 could exert bidirectional effects on cancer progression by acting as an NRG (17). FNDC4, a member of the fibroblast growth factor superfamily, has been less studied in the context of tumors. The study by Xiao et al. has confirmed the protective role of FNDC4 in chronic joint inflammation (48). Li et al. demonstrated that overexpression of FNDC4 suppresses the progression of ovarian cancer by promoting apoptosis and inhibiting cell proliferation (49). In addition, elevated transcriptional and translational expression of CXCL5, TYRO3, and FNDC4 in PTC specimens versus matched normal tissues (P<0.05) corroborates prior transcriptomic predictions from TCGA and GSE29625 cohorts.

Within the TCGA-PTC cohort, the three PRGs exhibited dichotomous expression patterns: CXCL5 was elevated in high-risk patients, whereas TYRO3 and FNDC4 showed increased expression in low-risk subgroups. Emerging evidence implicates necroptosis in oncogenesis, exemplified by Ando et al. (50) demonstrating CXCL5-mediated metastatic enhancement in pancreatic cancer via necroptotic signaling. Complementary work by Najafov et al. (51) identified TYRO3 as a regulator of MLKL oligomerization during necroptotic execution. Wang et al. (52) developed a prognostic risk model consisting of seven NRGs for patients with PTC. Their findings indicated that overexpression of genes such as IPMK, KLF9, and SPATA2 significantly inhibited the proliferation, invasion, and migration of PTC cells. Mechanistic understanding of necroptosis in PTC pathogenesis remains elusive. GSEA enrichment revealed CXCL5 potentiates PTC cell invasion/migration via mitogen-activated protein kinase (MAPK) cascades (c-Jun N-terminal kinase/extracellular signal-regulated kinase), and may drive osseous colonization through chemokine signaling and nuclear NF-κB activation-suggesting novel bone metastasis mechanisms. While only Ren et al. (53) report TYRO3’s protective role in PTC (aligning with present data), computational evidence implicates TYRO3 in JAK-STAT and PI3K/AKT pathway activation. Crucially, MAPK/NF-κB/PI3K-AKT/STAT pathways collectively modulate THCA progression (54). This study pioneers the association between elevated FNDC4 expression and PTC prognosis (validated at transcriptomic/proteomic levels), though mechanistic contributions of these DENRGs require further dissection.

Our molecular characterization revealed three differentially expressed necroptosis regulators showing significant enrichment in immunomodulatory pathways, indicating their potential role in shaping the immunological landscape of PTC. Notably, elevated immune scores and heightened infiltration of 27 immune cell subtypes were observed in the high-risk cohort, possibly attributable to the role of immune-inflammatory cells in necroptosis-driven tumor progression through angiogenesis stimulation, proliferation enhancement, and metastatic facilitation (16,55). Furthermore, the high-risk group displayed increased infiltration of immunosuppressive cell populations, such as γδT cells, myeloid-derived suppressor cells (MDSCs), regulatory T cells (Tregs), macrophages, monocytes, immature dendritic cells, plasmacytoid dendritic cells, and T follicular helper cells, indicative of an immunosuppressive microenvironment that facilitates tumor immune evasion and disease advancement (56). These findings align with the growing prominence of cancer immunotherapy as a promising therapeutic avenue (57). Differential expression analysis revealed upregulation of multiple immune checkpoints such as CTLA4, CD70/274, and TNFRSF9 in the high-risk group, suggesting potential exploitation of checkpoint pathways for immune escape. Although tumor mutational burden (TMB) did not differ significantly between risk groups, high-TMB cases exhibited poorer survival outcomes, underscoring its prognostic relevance in predicting immunotherapy response (58). The observed disparities in IPS following PD1 and CTLA4 blockade further support the therapeutic relevance of these checkpoints in high-risk PTC, given their established roles in T-cell immunomodulation and validated predictive utility in immunotherapy cohorts (59). While immune checkpoint inhibitors (ICIs) have gained clinical approval for multiple tumor types, their application in PTC remains investigational (60). Additionally, the DGIdb database identified 18 TYRO3-associated compounds as potential therapeutics, though clinical validation is warranted.

This study has several limitations, including reliance on computational methodologies susceptible to algorithmic biases and unresolved molecular mechanisms of DENRGs in PTC. Given the intricacy of tumor microenvironments and regulatory networks, further experimental validation is essential. Prospective studies are also needed to assess immunotherapy outcomes in PTC patients.

Conclusions

This study established and validated a three-DENRG prognostic signature for PTC, elucidating its molecular underpinnings, immune cell infiltration patterns, and implications for immunotherapy. The findings provide a foundation for exploring necroptosis in PTC pathogenesis and advance the potential for personalized therapeutic strategies.

Acknowledgments

None.

Footnote

Reporting Checklist: The authors have completed the TRIPOD reporting checklist. Available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-2124/rc

Peer Review File: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-2124/prf

Funding: None.

Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-2124/coif). The authors have no conflicts of interest to declare.

Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. The study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments.

Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.

References

1. Miller KD, Nogueira L, Devasia T, et al. Cancer treatment and survivorship statistics, 2022. CA Cancer J Clin 2022;72:409-36. [Crossref] [PubMed]
2. Bray F, Laversanne M, Sung H, et al. Global cancer statistics 2022: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin 2024;74:229-63. [Crossref] [PubMed]
3. Han B, Zheng R, Zeng H, et al. Cancer incidence and mortality in China, 2022. J Natl Cancer Cent 2024;4:47-53. [Crossref] [PubMed]
4. Zheng R, Zhang S, Zeng H, et al. Cancer incidence and mortality in China, 2016. J Natl Cancer Cent 2022;2:1-9. [Crossref] [PubMed]
5. Siegel RL, Kratzer TB, Giaquinto AN, et al. Cancer statistics, 2025. CA Cancer J Clin 2025;75:10-45. [Crossref] [PubMed]
6. Beheshtirouy S, Shayanfar A. A Review of Thyroid Cancer. JAMA 2024;331:1862-3. [Crossref] [PubMed]
7. National Cancer Center. Quality control index for standardized diagnosis and treatment of thyroid cancer in China (2022 edition). Zhonghua Zhong Liu Za Zhi 2022;44:902-7. [PubMed]
8. Siegel RL, Giaquinto AN, Jemal A. Cancer statistics, 2024. CA Cancer J Clin 2024;74:12-49. [Crossref] [PubMed]
9. Chen X, Lin S, Lin Y, et al. BRAF-activated WT1 contributes to cancer growth and regulates autophagy and apoptosis in papillary thyroid carcinoma. J Transl Med 2022;20:79. [Crossref] [PubMed]
10. Papaleontiou M, Norton EC, Reyes-Gastelum D, et al. Competing Causes of Death in Older Adults with Thyroid Cancer. Thyroid 2021;31:1359-65. [Crossref] [PubMed]
11. Yu L, Chen X, Liu J, et al. Associations between preoperative thyroid parameters, aggressive clinicopathological features and risk of recurrence in differentiated thyroid cancer. Ann Med 2025;57:2491153. [Crossref] [PubMed]
12. Tao Y, Wang F, Shen X, et al. BRAF V600E Status Sharply Differentiates Lymph Node Metastasis-associated Mortality Risk in Papillary Thyroid Cancer. J Clin Endocrinol Metab 2021;106:3228-38. [Crossref] [PubMed]
13. Cheng SH, Su YJ, Cheng RC. Research progress on thyroglobulin and anti-thyroglobulin antibodies in the diagnosis and treatment of differentiated thyroid cancer. Chinese Journal of Clinical Oncology 2020;47:1164-7.
14. Yuan J, Ofengeim D. A guide to cell death pathways. Nat Rev Mol Cell Biol 2024;25:379-95. [Crossref] [PubMed]
15. Linkermann A, Green DR. Necroptosis. N Engl J Med 2014;370:455-65. [Crossref] [PubMed]
16. Gong Y, Fan Z, Luo G, et al. The role of necroptosis in cancer biology and therapy. Mol Cancer 2019;18:100. [Crossref] [PubMed]
17. Li XY, Su LX, Chen WX, et al. Clinical Implications of Necroptosis Genes Expression for Cancer Immunity and Prognosis: A Pan-Cancer Analysis. Front Immunol 2022;13:882216. [Crossref] [PubMed]
18. Bedoui S, Herold MJ, Strasser A. Emerging connectivity of programmed cell death pathways and its physiological implications. Nat Rev Mol Cell Biol 2020;21:678-95. [Crossref] [PubMed]
19. Strilic B, Yang L, Albarrán-Juárez J, et al. Tumour-cell-induced endothelial cell necroptosis via death receptor 6 promotes metastasis. Nature 2016;536:215-8. [Crossref] [PubMed]
20. Höckendorf U, Yabal M, Herold T, et al. RIPK3 Restricts Myeloid Leukemogenesis by Promoting Cell Death and Differentiation of Leukemia Initiating Cells. Cancer Cell 2016;30:75-91. [Crossref] [PubMed]
21. Ritchie ME, Phipson B, Wu D, et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res 2015;43:e47. [Crossref] [PubMed]
22. Broom BM, Ryan MC, Brown RE, et al. A Galaxy Implementation of Next-Generation Clustered Heatmaps for Interactive Exploration of Molecular Profiling Data. Cancer Res 2017;77:e23-6. [Crossref] [PubMed]
23. Yu G, Wang LG, Han Y, et al. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS 2012;16:284-7. [Crossref] [PubMed]
24. In J, Lee DK. Survival analysis: part II - applied clinical data analysis. Korean J Anesthesiol 2019;72:441-57. [Crossref] [PubMed]
25. Robin X, Turck N, Hainard A, et al. pROC: an open-source package for R and S+ to analyze and compare ROC curves. BMC Bioinformatics 2011;12:77. [Crossref] [PubMed]
26. Wilkerson MD, Hayes DN. ConsensusClusterPlus: a class discovery tool with confidence assessments and item tracking. Bioinformatics 2010;26:1572-3. [Crossref] [PubMed]
27. Tong X, Tang R, Xiao M, et al. Targeting cell death pathways for cancer therapy: recent developments in necroptosis, pyroptosis, ferroptosis, and cuproptosis research. J Hematol Oncol 2022;15:174. [Crossref] [PubMed]
28. Gao W, Wang X, Zhou Y, et al. Autophagy, ferroptosis, pyroptosis, and necroptosis in tumor immunotherapy. Signal Transduct Target Ther 2022;7:196. [Crossref] [PubMed]
29. Chen Z, Wang C, Li M, et al. SPRED3 regulates the NF-κB signaling pathway in thyroid cancer and promotes the proliferation. Sci Rep 2024;14:20506. [Crossref] [PubMed]
30. Wen X, Liu M, Du J, et al. Meis homeobox 2 (MEIS2) inhibits the proliferation and promotes apoptosis of thyroid cancer cell and through the NF-κB signaling pathway. Bioengineered 2021;12:1766-72. [Crossref] [PubMed]
31. Luo X, Gao Q, Zhou T, et al. FOXP4-AS1 Inhibits Papillary Thyroid Carcinoma Proliferation and Migration Through the AKT Signaling Pathway. Front Oncol 2022;12:900836. [Crossref] [PubMed]
32. Okda TM, Atwa GMK, Eldehn AF, et al. A Novel Role of Galectin-3 and Thyroglobulin in Prognosis and Differentiation of Different Stages of Thyroid Cancer and Elucidation of the Potential Contribution of Bcl-2, IL-8 and TNF-α. Biomedicines 2022;10:352. [Crossref] [PubMed]
33. Hu B, Fan H, Lv X, et al. Prognostic significance of CXCL5 expression in cancer patients: a meta-analysis. Cancer Cell Int 2018;18:68. [Crossref] [PubMed]
34. Liu Q, Wang Y, Song M, et al. CCL20/CXCL5 Drives Crosstalk Between Anaplastic Thyroid Cancer Stem Cells and Tumor-Associated Macrophages to Promote Tumor Progression. Adv Sci (Weinh) 2025;12:e2405399. [Crossref] [PubMed]
35. Deng J, Ma X, Ni Y, et al. Identification of CXCL5 expression as a predictive biomarker associated with response and prognosis of immunotherapy in patients with non-small cell lung cancer. Cancer Med 2022;11:1787-95. [Crossref] [PubMed]
36. Nie Y, Jiang MC, Liu C, et al. CXCL5 Has Potential to Be a Marker for Hepatocellular Carcinoma Prognosis and Was Correlating With Immune Infiltrates. Front Oncol 2021;11:637023. [Crossref] [PubMed]
37. Wang T, Sheng J, Wang X, et al. CXCL5 Promotes the Malignant Phenotype of Pancreatic Cancer and Is Associated With Immune Infiltration. Clin Med Insights Oncol 2024;18:11795549241271691. [Crossref] [PubMed]
38. Liu Y, Qu HC, Huang Y. Renal clear cell carcinoma-derived CXCL5 drives tumor-associated fibroblast formation and facilitates cancer progression. Pathol Res Pract 2023;244:154319. [Crossref] [PubMed]
39. Le TN, Le MK, Dang MX, et al. CXCL5 expression is associated with active signals of macrophages in the microenvironment of papillary thyroid carcinoma. Pathol Res Pract 2024;261:155474. [Crossref] [PubMed]
40. Chang W, Chang Q, Lu H, et al. MicroRNA-873-5p suppresses cell malignant behaviors of thyroid cancer via targeting CXCL5 and regulating P53 pathway. Hum Vaccin Immunother 2022;18:2065837. [Crossref] [PubMed]
41. DeRyckere D, Huelse JM, Earp HS, et al. TAM family kinases as therapeutic targets at the interface of cancer and immunity. Nat Rev Clin Oncol 2023;20:755-79. [Crossref] [PubMed]
42. Kirimlioglu E, Katirci E, Simsek M. Distinct TYRO3 and PROS1 expression levels contribute to preeclampsia pathogenesis. Histochem Cell Biol 2025;163:29. [Crossref] [PubMed]
43. Zhang X, Li M, Xu Y, et al. Gal-3 activates Tyro3 to ameliorate ferroptosis of hippocampal neurons after traumatic brain injury. Mol Ther Nucleic Acids 2025;36:102433. [Crossref] [PubMed]
44. Zhang L, Jiang S, Shi J, et al. TYRO3 protects podocyte via JNK/c-jun-P53 pathway. Arch Biochem Biophys 2023;739:109578. [Crossref] [PubMed]
45. Hara K, Horikoshi Y, Morimoto M, et al. TYRO3 promotes chemoresistance via increased LC3 expression in pancreatic cancer. Transl Oncol 2023;28:101608. [Crossref] [PubMed]
46. Xiong Z, Sui X, Bai Y, et al. Hua Zheng San Ji Fang suppresses liver cancer progression by inhibiting TYRO3 expression via the ERK signaling pathway. Phytomedicine 2025;139:156497. [Crossref] [PubMed]
47. Miao YR, Rankin EB, Giaccia AJ. Therapeutic targeting of the functionally elusive TAM receptor family. Nat Rev Drug Discov 2024;23:201-17. [Crossref] [PubMed]
48. Xiao C, Xie S, Ouyang S. FNDC4 reduces inflammation, proliferation, invasion and migration of rheumatoid synovial cells by inhibiting CCL2/ERK signaling. Tissue Cell 2024;86:102294. [Crossref] [PubMed]
49. Li N, Li H, Zhou S, et al. Overexpression of FNDC4 constrains ovarian cancer progression by promoting cell apoptosis and inhibiting cell growth. J Cancer 2023;14:3416-28. [Crossref] [PubMed]
50. Ando Y, Ohuchida K, Otsubo Y, et al. Necroptosis in pancreatic cancer promotes cancer cell migration and invasion by release of CXCL5. PLoS One 2020;15:e0228015. [Crossref] [PubMed]
51. Najafov A, Mookhtiar AK, Luu HS, et al. TAM Kinases Promote Necroptosis by Regulating Oligomerization of MLKL. Mol Cell 2019;75:457-468.e4. [Crossref] [PubMed]
52. Wang Z, Wu P, Shi J, et al. A novel necroptosis-related gene signature associated with immune landscape for predicting the prognosis of papillary thyroid cancer. Front Genet 2022;13:947216. [Crossref] [PubMed]
53. Ren X, Du H, Cheng W, et al. Construction of a ferroptosis-related eight gene signature for predicting the prognosis and immune infiltration of thyroid cancer. Front Endocrinol (Lausanne) 2022;13:997873. [Crossref] [PubMed]
54. Jin S, Yang YT, Bao W. Signaling Pathways in Thyroid Cancer. Vitam Horm 2018;106:501-15. [Crossref] [PubMed]
55. Hanahan D, Weinberg RA. Hallmarks of cancer: the next generation. Cell 2011;144:646-74. [Crossref] [PubMed]
56. Hua X, Ge S, Zhang J, et al. A costimulatory molecule-related signature in regard to evaluation of prognosis and immune features for clear cell renal cell carcinoma. Cell Death Discov 2021;7:252. [Crossref] [PubMed]
57. Bergman PJ. Cancer Immunotherapy. Vet Clin North Am Small Anim Pract 2024;54:441-68. [Crossref] [PubMed]
58. Jardim DL, Goodman A, de Melo Gagliato D, et al. The Challenges of Tumor Mutational Burden as an Immunotherapy Biomarker. Cancer Cell 2021;39:154-73. [Crossref] [PubMed]
59. Charoentong P, Finotello F, Angelova M, et al. Pan-cancer Immunogenomic Analyses Reveal Genotype-Immunophenotype Relationships and Predictors of Response to Checkpoint Blockade. Cell Rep 2017;18:248-62. [Crossref] [PubMed]
60. Zhou W, Lin J, Liu J, et al. Thyroid cancer risk prediction model using m6A RNA methylation regulators: integrated bioinformatics analysis and histological validation. Aging (Albany NY) 2023;15:846-65. [Crossref] [PubMed]