KLB is lowly expressed in tumor tissueFigure 1 illustrates the study flow. To determine whether KLB plays a role in human cancers, we first examined the expression of the KLB gene in different cancers in the Cancer Genome Atlas (TCGA) database through the database TIMER 2.0. As shown in Fig. 2A, KLB is lowly expressed in 15 epithelial-type tumors compared to standard samples, including Breast invasive carcinoma, Cholangio carcinoma, Colon adenocarcinoma, Glioblastoma multiforme, Thyroid carcinoma, Kidney chromophobe, Kidney renal clear cell carcinoma, Kidney renal papillary cell carcinoma, Liver hepatocellular carcinoma, Lung adenocarcinoma, Lung squamous cell carcinoma, Prostate adenocarcinoma, Rectum adenocarcinoma, Stomach adenocarcinoma, and Head and Neck Squamous Cell Carcinoma. We compared the expression of KLB in HNSC with that of normal samples, and the results in Fig. 2B,C show that KLB is lower expressed in HNSC. In addition, Fig. 2E–H shows the immunohistochemical staining results of the HNSC tumor and paraneoplastic tissues, with KLB expressed at low levels in tumor tissue.Figure 1The whole analytical process of the study.Figure2Differences in KLB expression between tumor and paraneoplastic tissues. (A) KLB expression in pan-cancer. (B) KLB expression was lower in tumor tissues than in paraneoplastic tissues in TCGA-HNSC data. (C) Lower KLB expression in tumor tissues than paraneoplastic tissues in paired samples of TCGA-HNSC data. (D) In collected clinical samples, there was a higher KLB expression in paraneoplastic tissues (N = 8) than in tumor tissues (N = 8). (E–H) Immunohistochemical staining of head and neck tissues against KLB expression levels in the HPA database. (F) and (H) are from different parts of sections from the same patient, and (E) and (G) are from different patients’ normal mucosal epithelial tissues of the head and neck. (E) Specimen number HPA021136, site: head and neck normal tissue, patient id: 3362, (F) Specimen number HPA021136, site: head and neck, squamous carcinoma tissue, patient id: 2547, (G) Specimen number HPA021136, site: head and neck, normal nasopharyngeal tissue, patient id: 3402, (H) Specimen number HPA021136, site: head and neck, head and neck, squamous carcinoma tissue, patient id: 2547).To validate the results of bioinformatics analysis, we selected pathological specimens of cancerous (n = 8) and paraneoplastic (n = 8) tissues from patients with squamous cell carcinoma of the head and neck and extracted total tissue RNA. We analyzed the relative expression of KLB in cancerous and paraneoplastic tissues using the RT-qPCR method. The outcomes showed that the word of KLB was significantly lower in cancer tissues than in the paraneoplastic tissues, consistent with the bioinformatics analysis (Fig. 2D).KLB expression correlates with immune infiltrationNext, we assessed the relationship between KLB expression and immune infiltration using the CIBERSORT algorithm. As shown in Fig. 3A–H, the transcriptional expression of KLB was negatively correlated with Mast cells activated, Macrophages M2, Macrophages M0, and NK cells resting, and positively correlated with T cells regulatory (Tregs), T cells follicular helper, T cells CD8, Plasma cells and B cells naïve. On the other hand, KLB highly expressing HNSC patients with B cells naive, Plasma cells, T cells follicular helper and T cells regulatory (Tregs) high infiltration and NK cells resting, Macrophages M0 and Macrophages M2 low infiltration (Fig. 3I).Figure 3Relationship between KLB expression level and immune infiltration. (A–H) The transcriptional expression of KLB was negatively correlated with Mast cells activated, Macrophages M2, Macrophages M0, and NK cells resting, and positively correlated with T cells regulatory (Tregs), T cells follicular helper, T cells CD8, Plasma cells and B cells naïve. (I) Differences in immune cell infiltration between KLB high and low expression groups.Correlation of KLB expression with the prognosis of HNSC patientsTo investigate the prognostic value of KLB in HNSC, we used the Kaplan–Meier Plotter to predict the relationship between KLB expression and the overall survival of HNSC patients. Figure 4A showed a significant difference between the high and low KLB expression groups. Median survival was higher in the high KLB expression group than in the low KLB expression group, which indicated that KLB expression might be associated with a better prognosis for HNSC patients. In addition, we compared the impact of KLB on overall survival in male and female patients separately. The overall survival rate was higher in male patients with higher KLB expression than in patients with lower KLB expression. In contrast, there was no statistically significant difference in the effect of KLB on overall survival in female patients (Fig. 4B,C).Figure 4KLB expression levels correlate with survival in HNSC patients. (A–C). K–M plotter plotted survival curves for HNSC patients, with more prolonged survival in the high KLB expression group in overall and male patients and no statistically significant difference in survival between high and low KLB expression in female patients. (A). Survival of different KLB expression subgroups in the overall situation. (B) survival of different KLB expression subgroups when sex is female. survival of different KLB expression subgroups when sex is male). (D) Differential genes between high and low KLB expression groups. (E) Survival differences between the three subtypes. (F) and (G) Enrichment analysis between the three subtypes. (F) Comparative results of enrichment analysis between subtype B and subtype C. Twenty differentially expressed pathways between the two isoforms are illustrated in the figure. (G) Comparative results of enrichment analysis between subtype A and subtype B. Twenty differentially expressed pathways between the two isoforms are illustrated in the figure).To further investigate the prognostic value of KLB, we merged HNSC transcriptomic and clinical data from TCGA (normal = 44, tumor = 500) and a GEO dataset containing HNSC transcriptomic expression data and survival time and survival status (GSE65858, n = 270). We divided the 772 HNSC samples into a KLB high-expression group and a KLB low-expression group. We found the differential genes between the two groups, with the screening criteria of P < 0.05 and Log FC > 1 (Fig. 4D). We found 2798 genes that differed between groups and used the R package “Consensus Cluster Plus” to cluster all HNSC samples based on the expression patterns of the differential genes. The clustering was best when all samples were divided into three subtypes with good intra-subtype stability. We named the KLB subtypes A–C, with KLB cluster A containing 281 illustrations, KLB cluster B containing 409 samples, and KLB cluster C collecting 82 samples. Survival analysis showed that samples in cluster C had a better prognosis than those in clusters A and B (Fig. 4E). In addition, we used GSEA analysis to explore the internal characteristics of the three subtypes. We obtained enrichment analysis results showing that C subtypes are associated with immune activation, such as the Fc epsilon RI signaling pathway, B cell receptor signaling pathway, T cell receptor signaling pathway, the intestinal immune network for Ig A production, and natural killer cell-mediated cytotoxicity. In contrast, the B subtype is associated with pathways that promote tumor progression, such as the Wnt signaling pathway (Fig. 4F,G). In addition, immune cell infiltration-related pathways, such as the B cell receptor signaling pathway, T cell receptor signaling pathway, and natural killer cell-mediated cytotoxicity, were up-regulated in the A subtype compared to the B subtype.Prognostic risk assessment of HNSC patients using differential genes grouped by KLB expressionTo assess the prognostic risk of HNSC patients, we performed unifactorial Cox analysis on differential genes grouped by KLB expression to find prognosis-related genes; all samples were randomly divided into a training set and a test set and then performed lasso regression analysis and multifactorial Cox regression analysis on these prognosis-related genes to finally obtain 16 genes for assessing the KLB-Related Prognostic Risk Score (KRPRS) of HNSC patients. Eight of these positive factors were: ODF4, FAM187B, RAB37, CALN1, SRY, PLA2G2D, PIWIL2, CALML5 and eight negative factors were: TRIML2, KIF1A, SRPX, DPY19L2P1, CDKL2, RLN2, TMEFF2, C4BPB. The KLB-Related Prognostic Risk Score = expression level of ODF4*(− 1.135722681) + FAM187B*(-1.12696817) + RAB37*(− 0.494219462) + CALN1*(− 0.364682077) + SRY*(-0.324556358) + PLA2G2D*(− 0.220799715) + PIWIL2*(− 0.175252875) + CALML5*(− 0.060812431) + TRIML2*0.141813192 + KIF1A*0.16220514 + SRPX*0.178929297 + DPY19L2P1*0.428554182 + CDKL2*0.436786829 + RLN2*0.703418091 + TMEFF2*0.779739064 + C4BPB*1.317576808.Taking the median KRPRS as the cutoff value, in both the training and validation sets, the high KLB-Related Prognostic Risk Score group had a poorer prognosis than the low KRPRS group. We also analyzed the difference in KLB expression between the high-risk and low-risk groups. We showed that KLB expression was lower in the high-risk group (Fig. 5A), suggesting that lower KLB expression predicted a higher risk score and a worse prognosis. The results of the ROC analysis show that the area under the curve at 1, 3, and 5 years is above 0.65 (Fig. 5B). The conclusions we obtained in the training and test sets are consistent.Figure 5Construction of a prognostic risk score for HNSC patients based on differential genes between the three subtypes. (A) Survival differences between high and low-risk groups. (B) Accuracy of AUC curves to assess prognostic modeling. (C) and (D) Difference in expression of KLB, PD-L1 between high and low-risk groups. (E) Difference in immune infiltration between high and low-risk groups. (F) and (G) Mutation differences in transcription factors associated with tumor progression between high and low-risk groups. (F. Waterfall plot of gene mutations in the low-risk group. G. Waterfall plot of gene mutations in high-risk group).The low-risk score group had higher PD-L1 expression and immune cell infiltrationTo explore the differences between the high and low-risk score groups, we compared the differences in immune cell infiltration between the high and low-risk groups. Previous accumulating evidence suggests that high tumor mutation load implies durable anti-PD-1/PD-L1 immunotherapy; we compared the differences in PD-1 and PD-L1 expression between the high-risk and low-risk groups and showed that PD-L1 expression was higher in the low-risk group compared to the high-risk group (Fig. 5C,D). Therefore, the above results indirectly demonstrate that risk scores predict clinical response to anti-PD-1/PD-L1 immunotherapy. In addition, the results showed that the high-risk score group had higher Macrophages M0, Mast cells activated, and T cells CD4 memory resting infiltration. The low-risk score group had higher Macrophages M1, T cells CD4 memory started, T cells CD8, T cells follicular helper, T cells gamma delta, T cells regulatory (Tregs), and B cells naïve infiltration (Fig. 5E). Next, we analyzed high- and low-risk scores for tumor mutations in the TCGA-HNSC cohort using the map tools package. Figure 5F,G shows that the high-risk score group predicted a higher tumor mutation load than the low-risk score.Evaluation of response to drug candidates in patients at high and low risk of HNSCTo calculate drug sensitivity for each sample of patients in the high-risk and low-risk groups, we used the R package pRRophetic to identify 73 compounds that differed between the high-risk and low-risk scoring groups. The five drugs with the strongest positive correlation between risk score and drug sensitivity were ATRA(R = 0.49), THZ-2-49(R = 0.47), KIN001-102(R = 0.45), GSK690693(R = 0.36), MK-2206(R = 0.36) (Fig. 6A–J). Of these, ATRA (trans-retinoic acid), THZ-2-49 (selective CDK7 inhibitor), KIN001 (a patented drug combination of Pamapimod and pioglitazone, two active ingredients of an oral drug), GSK690693 (Pan-Akt inhibitor), and MK-2206, a selective Akt1/2/3 inhibitor. To predict the responsiveness of different risk groups to immune checkpoint inhibitor therapy, TIDE scores and the proportion of samples expected to be “responders” in various risk groups were compared between the high-risk and low-risk groups. The TIDE scores of the low-risk group were significantly lower than those of the high-risk group, and the proportion of “responder” samples in the low-risk group was higher than that in the high-risk group (Fig. 6K–L). This suggests that the samples in the low-risk group have a better response to immune checkpoint inhibitors.Figure 6Drug sensitivity analysis between high and low-risk groups. ( A–E) Risk Scores were positively correlated with IC50 for ATRA (R = 0.49), THZ-2-49 (R = 0.47), KIN001-102 (R = 0.45), GSK690693 (R = 0.36), and MK-2206 (R = 0.36). F–J) HNSC patients in the low-risk group had higher drug sensitivity to ATRA, THZ-2-49, KIN001-102, GSK690693, and MK-2206. K–L. TIDE score for assessing immune checkpoint inhibitor therapy responsiveness in individual tumor patients. (K) Percentage of Responder and Non-responder in high and low Risk groups. (L) Difference in TIDE scores between high- and low-risk groups).
