Integrating scRNA-seq and bulk RNA-seq to characterize immune microenvironment and construct a prognostic model for HNSCC

Abstract

Head and neck squamous cell carcinoma (HNSCC) is among the leading cancers across the globe and continues to be related to unfavorable clinical outcomes. In many cases, limited survival and poor prognosis remain major challenges. Recent studies revealed that the tumor microenvironment (TME) is key to HNSCC progression and might partly account for the suboptimal response to immunotherapy observed in a large proportion of patients. To better characterize TME heterogeneity, the research combines single-cell RNA sequencing (scRNA-seq) with bulk RNA sequencing (bulk RNA-seq) data to develop a prognostic model for HNSCC using publicly available datasets. The scRNA-seq data were obtained from the Gene expression omnibus (GEO) database, while bulk RNA-seq data were retrieved from the UCSC Xena platform. Cell populations within HNSCC samples were annotated using R-based analytical workflows. Malignant cell populations were inferred through infer copy number variation (inferCNV) analysis. Differentially expressed genes (DEGs) were observed via bioinformatic approaches, and weighted gene co-expression network analysis (WGCNA) was applied to detect modules associated with tumor-related traits. Overlapping genes were subsequently selected for downstream analyses, including prognostic model construction, gene set enrichment analysis (GSEA), gene set variation analysis (GSVA), immune infiltration assessment, mutation profiling, along with drug sensitivity evaluation. 8 cell clusters were identified from scRNA-seq data. Among them, epithelial cells exhibited the strongest malignant features, as indicated by higher CNV levels compared with T cells. Integration of 106 epithelial marker genes, 279 DEGs, and 1089 module genes led to the identification of 2 key genes, KRT5 and TUBA1B. A prognostic model incorporating risk score, clinical stage, and age was subsequently established, showing a significant difference in overall survival between high- and low-risk groups. Enrichment analyses revealed that cancer-related pathways, including proteoglycans in cancer and HIF-1 signaling cascade, were prominently involved. GSVA indicated increased activity in telomere tethering at the nuclear periphery in the high-expression group, whereas pathways related to cilium movement were relatively suppressed. Immune infiltration analysis suggested a generally limited responsiveness of HNSCC to immunotherapy in the absence of specific clinical indicators. In addition, several compounds, including BRD.K37390332, NSC.74859, fluvastatin, and pifithrin-α, were identified as potential therapeutic candidates, although further validation is required. Through the combination of scRNA-seq and bulk RNA-seq data, this study establishes a prognostic model with moderate predictive performance for HNSCC. The genes KRT5 and TUBA1B emerge as potential biomarkers. Despite these findings, the predicted sensitivity to immunotherapy remains limited derived from computational analyses, highlighting the necessity for additional experimental and clinical confirmation.