In this study, we used a large bioinformatics approach to identify genes associated with exercise and survival prognosis in HNSCC. Transcriptome data for HNSCC were initially downloaded from The Cancer Genome Atlas (TCGA) database. Statistical analyses were conducted to identify genes exhibiting significant differential expression between HNSCC tissues and normal tissues. Subsequently, survival analysis was performed on these differently expressed genes to pinpoint those closely linked to patient prognosis. Gene sets related to exercise, encompassing various biological processes associated with exercise response, were obtained from the Gene Set Enrichment Analysis (GSEA) database. The intersection analysis on differentially expressed genes and exercise-related gene sets was conducted to obtain the key genes using public GSEA database. We also demonstrated that these genes were not only exercise-response related but also extremely correlated with the survival prognosis of HNSCC patients. The identification of these resplices provides new insights into the molecular biology underlying HNSCC and may yield novel.

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were conducted on the selected genes to elucidate their distribution in biological processes, molecular functions, and cellular components, as well as their involvement in specific metabolic and signaling pathways. Gene Set Enrichment Analysis (GSEA) was then employed to investigate the expression patterns of these genes under various biological conditions, aiming to identify gene sets that may influence the survival prognosis of HNSCC patients. And survival prognosis analysis was conducted to test relationship between these gene expression levels and overall survival (OS) rate with disease free survival rate. Single-gene expression analysis was also conducted to delve into the expression patterns of each gene and their associations with the clinical characteristics and prognosis of HNSCC patients. These identification results are helpful in deepening understanding of the molecular mechanism for HNSCC and lay a scientific foundation in exploring novel therapeutic targets and prognostic biomarkers.

In the study of HNSCC, tumor samples were divided into two groups based on high and low PIGF expression levels. Differential expression analysis was conducted using the limma package, identifying significantly differently expressed genes, which were visualized with volcano plots. Protein-protein interaction data were filtered through the ComPPI database to exclude biologically unreasonable interactions, and interaction scores were introduced to quantify data accuracy. ROC analysis was performed using the pROC package to calculate the 95% confidence interval and AUC, and ROC curves were plotted to evaluate the diagnostic efficiency of gene expression in distinguishing tumor from normal tissues. The data were sourced from TCGA-corrected RNA-seq data, generated through the Firehose pipeline and normalized. Z-score standardization was used to identify outliers, and the Wilcoxon rank-sum test assessed expression differences between tumor and normal tissues. Furthermore, gene expression was normalized using Z-score after pairing data from the GTEx database with TCGA data to eliminate outliers, followed by ROC analysis to evaluate diagnostic performance. The Wilcoxon rank-sum test was used to compare PIGF expression between tumor and normal tissues. Calibration curves and goodness-of-fit tests were applied to assess the predictive accuracy of the models. GEO datasets were processed by converting probe matrices to gene matrices and applying Z-score standardization with the Wilcoxon rank-sum test evaluating expression differences between tumor and normal tissues. Additionally, six molecular immune subtypes associated with tumor characteristics and prognosis were evaluated using median grouping and chi-square tests to assess the significance of subtype proportions. The Kruskal-Wallis rank-sum test compared PIGF expression across different molecular subtypes, while clinical variables were statistically analyzed in different expression groups using median grouping and chi-square tests.

In tumor tissues, the Pearson correlation between the target gene and both mRNA and miRNA was calculated, with scatter plots used to display these relationships. Results were reported only when the absolute value of the correlation coefficient exceeded 0.3. Gene expression levels were categorized based on their correlation strength with the target gene into four classes: strongly positive, moderately positive, weakly positive, and negative correlations. These were visualized using a heatmap of contingency tables, and Fisher’s exact test was employed for statistical analysis. Kaplan-Meier survival analysis was used to evaluate the correlation between gene expression levels and patient survival times. Detailed survival data analysis was performed using the survival package in R. The survminer package was used to identify optimal cut-off values for high and low expression groups, ensuring that the sample sizes for these groups met statistical requirements, typically not less than 30% of the total sample size. Moreover, a meta-analysis that based on univariate Cox proportional hazards model was performed by using the inverse variance method. We chose the hazard ratios (HR) as our primary measure of effect size, separating potential tumor-suppressive versus oncogenic actions. This simple classification approach does not consider the biological aspects of these genes. The statistical analysis and visualization were performed in R (version 4.3.2) using the meta package, which provides a range of functions for conducting meta-analyses and generating forest and funnel plots, visually presenting the combined effect sizes and assessing publication bias.

The study had used a stratified approach to classify samples as high or low expression group basing on the top 30% of most expressed samples versus bottom 30%. This second classification aimed to find the most extreme gene expression changes associated with disease states. This was followed up by a differential expression analysis using the limma package, producing log2 fold changes (log2FC) and ranking genes that were statistically altered. The additional analysis utilizes gene sets from the KEGG database and was performed with fgsea function in the R package fgsea. Enrichment scores (ES) of gene sets with significant P values were calculated by using GSEA analysis, and tested for significance/multiple hypothesis correction was applied. Gene Sets with uncorrected p-value < 0.05 and corrected p-value < 0.25 were assumed to have biological significance and visualized as before based on species partitioning. To explore the states of tumor cells better in single-cell level, we performed CancerSEA analysis. This platform creates a merged view of the datasets and uncovers 14 new functional states that are able to articulate significant aspects of tumor cell function, making it usable as an effective resource for conducting meaningful experiments. Using the z-score algorithm proposed by Lee et al., gene set values were calculated and converted to z-scores with the GSVA algorithm in the R package GSVA. Pearson correlation analysis was then employed to explore the relationships between gene expression and functional states, calculating the correlation between gene expression and gene set z-scores. Finally, the gsva function in the GSVA package was used to score 73 metabolic gene sets from the KEGG database. Based on these GSVA scores, the limma package was again used to compare metabolic pathway activities between the high and low expression groups, revealing the role of metabolic pathways in disease progression.

Gene expression data from multiple publicly available datasets of cancer patients undergoing immunotherapy were utilized. To assess the diagnostic performance of PIGF expression in distinguishing between responders and non-responders to immunotherapy, ROC curve analysis was conducted using the pROC package in R. The area under the curve (AUC) and 95% confidence intervals (CI) were calculated, and smoothing techniques were applied to the ROC curves for improved visualization. Patients were categorized into high and low PIGF expression groups based on the median expression level of PIGF. A Chi-square test was performed to examine the differences in the proportions of responders and non-responders between these two groups. Furthermore, the Wilcoxon rank-sum test was employed to compare gene expression differences between responders and non-responders, and gene expression levels were standardized into Z-scores for statistical comparison.

The TIMER 2.0 database was employed to collect and analyze immune infiltration data from TCThe TIMER 2.0 database was used in this study to obtain immune infiltration data from the TCGA tumor samples5. We developed a comprehensive immune-related long non-coding RNAs database and discovered that multiple algorithms can evaluate the abundance of different types of immune cells in tumor tissues, and analyze their relationship with gene expression. This approach served data quality and provenance, giving a holistic overview on the interaction of immune cells to gene expression. The bar of scatter plots shows us correlation coefficient, which clearly describes the relationship between immune cell type and gene expression. Samples were bimodalized based on median of gene expression as an robust estimator to differentiate between low and hight expressions people — the healthy group in blue line while ICH is presented by red solid line. To find whether the immune cell contents in high and low expression groups had statistical differences, Wilcoxon rank sum test was performed as a nonparametric method applied to comparisons among different kinds of data distributions. Heatmaps were generated to show key immune cell types in detail.

Whole-genome CRISPR screening data were obtained from the DepMap portal, and dependency scores for approximately 17,000 candidate genes were analyzed using the CERES algorithm. The pan-cancer mutation landscape of core genes was visualized using the plotmafSummary function from the maftools package. To assess the independence between gene expression levels and specific gene mutation types the independence_test function from the R coin package, based on permutation tests, was employed. Genes with a mutation rate exceeding 10% and a p-value less than 0.01 were identified and visualized to highlight significant associations between gene expression and mutation types. In the TCGA-HNSC project, copy number variation (CNV) analysis was performed using the GISTIC score method to identify genomics CNVs. The CNV profile of 451 samples was visualized using bar plots, which reflected chromosomal copy number changes. Quantitative measures of genomics alterations, such as FGA, FGG, and FGL, were defined based on the genomics distance of clonal regions. When analyzing differences between specific gene expression subgroups, ANOVA was used, followed by Tukey’s Honest Significant Difference (TukeyHSD) test for multiple comparisons if ANOVA indicated significance. This approach was used to identify specific group differences. The correlation between CNV scores and gene expression levels was analyzed using scatter plots combined with Spearman’s rank correlation coefficient, which measures the monotonic relationship between two variables. CNV data were obtained from the TCGA Genome Characterization Center andmeasured through whole-genome arrays. Gene-level copy number estimates were derived using the TCGA FIREHOSE pipeline and GISTIC2 method. The Kruskal-Wallis test, a non-parametric method for comparing multiple samples, was employed to compare gene expression differences among different CNV types (ranging from −2 to 2).

Single-cell gene expression data for HNSCC were obtained from the TISCH database. Heatmaps, generated using the pheatmap package, effectively revealed gene expression patterns at the single-cell level across different cancer types. Hierarchical clustering analysis, performed using Euclidean distance and Ward’s minimum variance method, uncovered intrinsic patterns of gene expression and their conservation across various cancers. Additionally, UMAP (Uniform Manifold Approximation and Projection) was utilized to explore expression patterns in high-dimensional data, preserving the original data topology during dimensionality reduction. By using UMAP analysis of gene expression data for CENPF, we could render a clearer depiction of the patterns underlying adding new observations to key biological discovery. To determine significant differences in specific gene expression between cell types, the Kruskal-Wallis rank sum (KW) test was applied as a non-parametric method suitable for un-normally distributed samples. Furthermore, UMAP visualization of the AUC cell scoring capturing heterogeneity of pathway activity in individual cells was performed. Methodology for spatial transcriptomics can be found in the supplementary methods.

Cell proliferation was determined by the cell counting Kit-8 (CCK-8, Beyotime Biotechnology Co., Ltd., Shanghai, China) at 0 h, 24 h and 48 through various treatments. The cells were then cultured in 96-well plates (Thermo Fisher Scientific. MA, United States) and exposed to respective interventions for another 24 h essentially as described above briefly After that, 10 μL of CCK-8 solution was added to each well and incubated for another 2 h. Absorbance (450 nm) was determined by a Microplate Reader (Bio-Rad, Hercules, CA). Cell viability was evaluated by calculating the ratio of the average absorbance of the treated groups to that of the control group, expressed as a percentage ([Absorbance of treated group/Absorbance of control group] × 100%).

Cells in the logarithmic growth phase were collected and diluted to a concentration of 500 cells/mL. Each well of a 6-well plate was pre-wetted with 1 mL of culture medium before adding 1 mL of the cell suspension. Three replicate wells were prepared for each group. The cells were incubated overnight in a 37°C, 5% CO₂ incubator to allow for attachment. Subsequently, cells were collected and 5 × 10^4 cells per well were added to the corresponding wells, with the medium being changed every 2 days. After 12 days, the medium was discarded, and the wells were washed twice with PBS. Cells were then fixed by adding 1 mL of methanol to each well and incubating at room temperature for 20 min. After removing the methanol, 1 mL of 0.1% crystal violet was added for staining at room temperature for 20 min. The wells were then washed with PBS until the background was clear, followed by photographing and counting the colonies.

SCC4 cells were seeded into 6-well plates at a density of 5 × 10⁵ cells per well. After 12 h of incubation, a sterile pipette tip was used to create scratches along predefined tracks. Detached cells were washed away with PBS, and photographs were taken. The plates were then returned to the incubator for an additional 24 h, followed by another round of photography to capture cell migration. For the Transwell assay, treated SCC4 cells were seeded into the upper chamber of Transwell inserts pre-coated with Matrigel at a density of 1 × 10⁵ cells per well. The lower chamber was filled with 600 μL of culture medium containing 5% fetal bovine serum. The cells were incubated for 24 h, after which they were fixed with 4% paraformaldehyde for 20 min and stained with 0.1% crystal violet for 10 min. The number of invading cells was then observed under a microscope and photographed.

The statistical analysis was conducted using SPSS software version 26.0 (SPSS Inc., Chicago, United States). The results are expressed as the mean ± standard deviation. For comparisons between two groups, Student’s t-test was utilized, while one-way ANOVA was applied for comparisons across multiple groups. Statistical significance was defined as P-value < 0.05.

The results of this study present the identification and multigene analysis of exercise-survival prognosis-related genes in HNSCC, as illustrated in Figure 1. The Venn diagram (Figure 1A) highlights the intersection of exercise-related genes, tumor-related genes, and survival prognosis-related genes in HNSCC, pinpointing core genes crucial for survival prognosis influenced by exercise. The distribution of these core genes across various biological processes and pathways is depicted in the GO analysis pie chart (Figure 1B), emphasizing their involvement in diverse cellular functions. Further insights are provided by the KEGG pathway enrichment analysis (Figure 1C), which identifies significant pathways involving these core genes. The GSEA plot (Figure 1D) demonstrates the enrichment of core genes in specific biological pathways or processes, indicating their functional relevance in the context of HNSCC and exercise. Survival analyses using the ssGSEA score are presented for different HNSCC patient cohorts (Figures 1E, F, datasets GSE126 and GSE525), revealing the association between core gene expression levels and overall survival. Additional Kaplan-Meier survival analyses (Figures 1G, H, datasets GSE407 and GSE123) consistently show the correlation between core gene expression and patient survival, reinforcing their prognostic value. A meta-analysis of survival data (Figure 1I) confirms the significant prognostic impact of these genes across multiple datasets. Multivariate Cox regression analyses (Figures 1J–M) further illustrate the independent prognostic value of these genes, highlighting their potential as biomarkers for HNSCC prognosis. These comprehensive analyses demonstrate the critical role of exercise-survival prognosis-related genes in HNSCC and their potential as prognostic biomarkers, providing valuable insights into the molecular mechanisms of HNSCC influenced by exercise.

The study focused on analyzing the expression pattern of the PIGF gene in HNSCC using data from The Cancer Genome Atlas (TCGA). Figure 2A displays a volcano plot generated through differential expression analysis. The diagnostic potential of PIGF in distinguishing tumor tissues from normal tissues is evaluated by the ROC curve in Figure 2B, which presents an AUC value indicating PIGF’s predictive ability. Further analysis in Figure 2C shows a violin plot comparing PIGF expression levels between normal and tumor tissues, approaching statistical significance (p = 0.06), while Figure 2D presents a paired sample plot indicating a non-significant difference (p = 0.318). The expression pattern across different tumor stages is shown in Figure 2E, with significant differences (p = 0.021) between early-stage (Stage I-II) and advanced-stage (Stage III-IV) HNSCC. Figure 2F’s ROC curve assesses the diagnostic performance of PIGF in differentiating early from advanced stages, demonstrating moderate efficacy. Figure 2G visualizes the gene interaction network with PIGF as the central node, indicating its interactions with other genes. Finally, Figure 2H shows the median expression levels of PIGF at different tumor stages, emphasizing the dynamic changes in expression during tumor progression. The expression level of PIGF was significantly upregulated in stage III and IV HNSCC, showing a linear relationship, indicating a correlation with disease severity. However, the expression level of PIGF in stage II HNSCC is lower than that in stage I. These results collectively underscore the significant upregulation of PIGF in HNSCC, its diagnostic potential, and its variable expression across different tumor stages, enhancing our understanding of its role in cancer progression.

A comprehensive analysis of the PIGF gene in HNSCC using TCGA-GTEX data provides significant insights into its diagnostic and clinical relevance (Figures 3A–F). Figure 3A shows the ROC curve assessing the diagnostic efficiency of PIGF expression between tumor and normal groups, including the Hosmer-Lemeshow test with a p-value of 0.096, indicating a good model fit. Figure 3B further highlights the diagnostic performance of PIGF expression in distinguishing tumor from normal tissue in HNSCC, with an area under the curve (AUC) of 0.660 and a 95% confidence interval of 0.597–0.724, demonstrating high diagnostic accuracy. The differential expression of the PIGF gene across various tumor subtypes of HNSCC is depicted in Figure 3C. The violin plot illustrates the distribution of PIGF mRNA levels in atypical, basal, classical, and mesenchymal subtypes, with a significant p-value of less than 0.001, indicating notable differences among subtypes. Figure 3D presents a heatmap showing expression differences of the PIGF gene in HNSCC patients from TCGA, categorizing expression levels into high and low groups and correlating these with patient subtypes. Figure 3E shows a density map comparing the estimated expression of PIGF, showing the obvious differences between normal and tumor tissues. The results showed that the expression level of PIGF in tumor tissues was significantly higher than that in normal tissues. Chi square test results (Figure 3F) highlighted the association between PIGF expression level and various clinical characteristics of hNSC patients, including alcohol consumption, HPV status, lymph node involvement, tumor stage, gender, targeted therapy, radiotherapy, tumor grade and patient age, discriminated between high expression group and low expression group and provided the relevant p value. The results showed that patients with high PIGF expression had later tumor stage, more lymph node involvement, and were more sensitive to targeted therapy and radiotherapy. These findings collectively reveal the multidimensional characteristics of the PIGF gene in HNSCC, underscoring its potential as a diagnostic biomarker and its association with clinical traits.

Pearson correlation analysis is used to identify significant interactions between PIGF expression and various related mRNAs and miRNAs. Figure 4A presents a heatmap displaying the correlation matrix between PIGF mRNA and other mRNAs, with color intensity indicating the correlation strength—blue for negative and red for positive correlations. Figure 4B features a scatter plot with a regression line, illustrating the Pearson correlation between PIGF mRNA expression and a specific related mRNA. The x-axis represents PIGF mRNA levels, and the y-axis represents the related mRNA levels, highlighting the direction and strength of the correlation. Furthermore, Figures 4C–L depict scatter plots showing the correlation between PIGF mRNA and various miRNAs. Each plot, such as Figure 4C, includes the correlation coefficient and p-value, indicating a statistically significant positive correlation. The analysis reveals moderate to strong correlations between PIGF mRNA and these mRNAs and miRNAs, with significant p-values. These findings suggest that PIGF may play a crucial role in the regulatory network involving these genes, providing a foundation for future studies on the functional implications of PIGF interactions in the studied condition.

A comprehensive analysis is conducted to evaluate the prognostic significance of PIGF gene expression across various survival outcomes using internal and external datasets. Kaplan-Meier survival analysis is performed for four key survival outcomes: Overall Survival (OS), Disease-Specific Survival (DSS), Progression-Free Interval (PFI), and Disease-Free Interval (DFI). Survival curves are stratified by quartiles (Q1-Q4) of PIGF expression levels (Figures 5A–D). The results indicate no significant difference in OS (log-rank test P = 0.13; Figure 5A), a significant difference in PFI (log-rank test P = 0.002; Figure 5B), no significant difference in DSS (log-rank test P = 0.232; Figure 5C), and no significant association with DFI (log-rank test P = 0.363; Figure 5D). A meta-analysis of univariate Cox regression survival analysis across multiple datasets determines that higher PIGF expression correlates with an increased risk of poor survival outcomes (Figures 5E, F). External GEO datasets analysis validates these findings in head and neck squamous cell carcinoma (HNSCC), with significant associations observed in datasets GSE10406 (P = 0.004; Figure 5G), GSE84318 (P = 0.016; Figure 5H), and GSE53161 (P = 0.001; Figure 5I), though no significant association is found in GSE24362 (P = 0.348; Figure 5J). Further independent prognostic analysis confirmed the significance of PIGF expression independent of clinical variables through both univariate and multivariate Cox regression analyses (Supplementary Figures 1A–D). Expression levels are correlated with overall survival status (Supplementary Figure 1E), and a Chi-square test indicated no significant distribution difference across expression quartiles (P = 0.706; Supplementary Figure 1F). Restricted cubic spline models explored non-linear risk associations adjusted for relevant covariates, showing estimated log hazard ratios with 95% confidence bands (Supplementary Figures 1G–J). Kaplan-Meier survival curves stratified by expression levels reveal significant differences, particularly for OS and DSS (Supplementary Figures 1K–N). This analysis underscores the significant prognostic value of PIGF gene expression for various survival outcomes, suggesting its potential as a biomarker for cancer prognosis.

The analysis of core genes in HNSCC provides comprehensive insights into their roles in tumor progression and patient survival outcomes through various Kaplan-Meier (KM) survival analyses, GSEA, and GSVA. Figures 6A–E demonstrate KM survival analyses for the four subgroups of dual-gene molecular subtypes based on CD274 expression. These figures show that high CD274 expression correlates with poorer survival outcomes. Similarly, Figures 6F–J demonstrate that elevated PDCD1 expression is associated with reduced survival rates. Figure 6K highlights significant pathways identified through GSEA for hallmark gene sets, revealing differently regulated biological processes. The KEGG gene set enrichment analysis in Figure 6L identifies critical signaling pathways involved in HNSCC, including those related to immune response and apoptosis. Figure 6M uses the clusterProfiler package to compare high and low expression groups, showcasing normalized enrichment scores (NES) for significantly enriched pathways, thereby emphasizing the diverse biological functions influenced by core genes. Differential GSVA scores for metabolic pathways in Figure 6N indicate altered metabolic activities linked to core gene expression. The heatmap in Figure 6O illustrates the correlation between immune response signatures and genome state, highlighting the interplay between immune activity and genetic alterations. Figure 6P represents PIGF expression across various immune stimulators, underscoring its role in modulating immune responses. Finally, Figure 6Q presents Pearson correlation analyses between z-scores of core gene expression and tumor state parameters, offering insights into the relationships between core gene expression and tumor-related parameters. This comprehensive analysis provides a detailed perspective on the functional roles of core genes in HNSCC, emphasizing their potential as prognostic biomarkers and therapeutic targets.

Our study provides an in-depth analysis of PIGF gene expression in HNSCC using single-cell sequencing data. Figures 7A–C present UMAP plots displaying major cell lineages (A), single-cell PIGF gene expression (B), and density contour lines (C) in HNSCC, highlighting the distribution and expression patterns of PIGF across different cell populations within the tumor microenvironment. Figures 7D–F show UMAP plots depicting differential expression of core genes across various cell types in HNSCC, illustrating the variability in gene expression among different cell populations. Figure 7G represents pathway differences between core gene-positive and core gene-negative groups across various cell types, with a dot plot displaying enriched pathways in each cell type, indicating functional pathways associated with PIGF expression. Figure 7H utilizes spatial transcriptomic deconvolution to show cellular composition with the maximum value for each spot, providing a spatial map of cell distribution within the tumor. Figure 7I illustrates the Spearman correlation between gene expression and microenvironment components at single-cell resolution, revealing interactions between gene expression and the tumor microenvironment. Figure 7J demonstrates the differential expression of PIGF in malignant, mixed malignant, and normal regions, suggesting its potential role in tumor progression. Figure 7K shows PIGF expression across different tumor stages allowing comparison across various stages of tumor development. Figure 7L presents the active landscape of core gene set scores in microzones, identifying regions within the tumor with high or low activity of specific gene sets. Figure 7M shows differences in AUC scores of gene sets between malignant, mixed malignant, and normal microzones, highlighting significant differences in gene set activity among the different regions. This detailed single-cell sequencing analysis reveals the complex role of PIGF in the tumor microenvironment of HNSCC. Additionally, in vitro cell experiments demonstrates that silencing PIGF affects cell proliferation, apoptosis, and the expression of related factors, suggesting potential molecular targets for the HNSCC treatment.

This study evaluated the role of PIGF in the proliferation and colony formation of HNSCC cells, as well as its potential as a drug target. By knocking down (sh1-PIGF) or overexpressing (PIGF-OE) PIGF in different experimental groups, real-time quantitative PCR data indicated a significant reduction in PIGF mRNA expression in the sh1-PIGF group, while the PIGF-OE group showed a notable increase, confirming the efficiency of both knockdown and overexpression constructs (Figure 8A). Functional analysis revealed that colony formation in the sh1-PIGF group was significantly reduced, suggesting that PIGF knockdown inhibited cell proliferation. Conversely, PIGF-OE reversed this effect, further demonstrating the critical role of PIGF in promoting HNSCC cell proliferation (Figure 8B). Similarly, the CCK-8 assay showed that PIGF knockdown significantly reduced cell viability, while PIGF overexpression restored the proliferative ability of the sh1-PIGF group (Figure 8C). Immunofluorescence staining demonstrated that the expression of proliferation markers MYC and Ki67 was markedly reduced in the sh1-PIGF group, whereas both markers were upregulated in the PIGF-OE group, supporting the hypothesis that PIGF promotes tumor cell growth by regulating proliferation-related pathways (Figure 8D). To further explore the potential of PIGF as a drug target, we analyzed its predictive performance across multiple cancer datasets (Figure 8E). In the GBM dataset, PIGF exhibited strong predictive power, particularly in the PIGF_GBM-PRJNA482620 dataset. Moreover, analysis of several cancer immunotherapy-related metrics (e.g., TIDE, MSI score, TMB, CD274 expression, CD8 infiltration, IFNγ expression) revealed significant variability in PIGF’s predictive performance across different datasets, with high predictive accuracy observed in glioblastoma and certain melanoma datasets (Figure 8F). These findings suggest that PIGF may play a role in modulating the immune microenvironment, making it a potential target for cancer immunotherapy.