Development and implementation of a prognostic model for clear cell renal cell carcinoma based on heterogeneous TLR4 expression

Abstract

Objective: Clear cell renal cell carcinoma (ccRCC) is the most common subtype among renal cell carcinomas and has the worst prognosis, originating from renal tubular epithelial cells. Toll-like receptor 4 (TLR4) plays a crucial role in ccRCC proliferation, infiltration, and metastasis. The aim of this study was to construct a prognostic scoring model for ccRCC based on TLR4 expression heterogeneity and to explore its association with immune infiltration, thereby providing insights for the treatment and prognostic evaluation of ccRCC.

Methods: Using R software, a differential analysis was conducted on normal samples and ccRCC samples, and in conjunction with the KEGG database, a correlation analysis for the clear cell renal cell carcinoma pathway (hsa05211) was carried out. We observed the expression heterogeneity of TLR4 in the TCGA-KIRC cohort and identified its related differential genes (TRGs). Based on the expression levels of TRGs, consensus clustering was employed to identify TLR4-related subtypes, and further clustering heatmaps, principal component, and single-sample gene set enrichment analyses were conducted. Overlapping differential genes (ODEGs) between subtypes were analysed, and combined with survival data, univariate Cox regression, LASSO, and multivariate Cox regression were used to establish a prognostic risk model for ccRCC. This model was subsequently evaluated through ROC analysis, risk factor correlation analysis, independent prognostic factor analysis, and intergroup differential analysis. The ssGSEA model was employed to explore immune heterogeneity in ccRCC, and the performance of the model in predicting patient prognosis was evaluated using box plots and the oncoPredict software package.

Results: In the TCGA-KIRC cohort, TLR4 expression was notably elevated in ccRCC samples compared to normal samples, correlating with improved survival in the high-expression group. The study identified distinct TLR4-related differential genes and categorized ccRCC into three subtypes with varied survival outcomes. A risk prognosis model based on overlapping differential genes was established, showing significant associations with immune cell infiltration and key immune checkpoints (PD-1, PD-L1, CTLA4). Additionally, drug sensitivity differences were observed between risk groups.

Conclusion: In the TCGA-KIRC cohort, the expression of TLR4 in ccRCC samples exhibited significant heterogeneity. Through clustering analysis, we identified that the primary immune cells across subtypes are myeloid-derived suppressor cells, central memory CD4 T cells, and regulatory T cells. Furthermore, we successfully constructed a prognostic risk model for ccRCC composed of 17 genes. This model provides valuable references for the prognosis prediction and treatment of ccRCC patients.

Keywords: Consensus clustering; Prognostic model; Renal clear cell carcinoma; TLR4.

Fig. 1

Heterogeneity of TLR4 in TCGA-KIRC. A: Inter-group comparison between normal and ccRCC TLR4 expression, with light green representing the normal group, and brown representing the ccRCC group. B: Kaplan-Meier survival analysis diagram comparing high and low TLR4 expression levels in ccRCC. Blue represents the low TLR4 expression group, while yellow signifies the high TLR4 expression group. C: Single-gene GSEA chart, with the peak on the left indicating that the number of positively correlated genes in the differentially expressed genes between TLR4 high and low expression groups exceeds the number of negatively correlated genes. D: Volcano plot, with red representing 3557 upregulated genes, and blue representing 1119 downregulated genes. E: Venn diagram. Blue denotes high and low TLR4 differentially expressed genes (TLR4_Deg), and red denotes genes recorded in renal cell carcinoma (hsa05211; RCC_gene). The overlapping nine genes are therapeutic response genes (TRGs).

Fig. 2

Heterogeneity Analysis of TRGs. A: Dumbbell plot, where red dots represent copy number increase (GAIN) and blue dots represent copy number reduction (LOSS). The y-axis represents the proportion (%) of samples with gene alterations. B: Boxplot for intergroup differential analysis. Red represents ccRCC samples, and blue represents normal tissue samples. The x-axis indicates the names of TRGs genes, while the y-axis represents the expression levels of these genes.C–K: Kaplan-Meier survival plots for individual TRGs. Blue represents low expression group for a given TRG, and red indicates the high expression group for that TRG.

Fig. 3

Consensus Clustering Analysis and Evaluation. A–C: Consensus matrix plots when clustering into 2, 3, and 4 subclasses, respectively. Deeper shades of blue indicate higher similarity among base clusters, lighter shades of blue or colorless squares indicate lower similarity. The size of each square represents the weight of the cluster in the final classification. D: Consensus Cumulative Distribution Function (CDF) plot. The x-axis represents the similarity threshold of clustering results, and the y-axis represents the frequency of sample assignments to different clusters. A flatter curve in the plot suggests closer frequencies of sample assignments to clusters, indicating better clustering performance. E: Delta area plot in consensus clustering. When the Delta area value tends to stabilize, and the line in the plot becomes flat, the optimal number of clusters is considered to be reached. F: Kaplan-Meier survival plot. Different colors represent different cluster subtypes, where red signifies subtype A, yellow represents subtype B, and gray denotes subtype C. G: In the Principal Component Analysis (PCA) plot, each point represents the position of a sample in the principal component space, which in turn represents the concentration of gene expression for that sample. H: Multidimensional heatmap: red represents high gene expression, while blue indicates low gene expression. I: Single-sample Gene Set Enrichment Analysis (ssGSEA) plot. The x-axis represents the names of infiltrating immune cells in the tissue, while the y-axis represents the relative abundance of immune cell infiltration. “ns” denotes no statistical significance.

Fig. 4

Construction and Validation of the Prognostic Model. A: Venn diagram comparing the overlap of three groups of differentially expressed genes obtained from pairwise differential analysis. B: LASSO regression coefficient shrinkage path plot: shows how gene coefficients are reduced with the increase in regularization strength. C: Cross-validation error plot for LASSO regression: depicts the model's prediction error changing with the change in regularization strength. The prediction error is minimal when λ = λmin. D: Forest plot for the prognostic risk model: The x-axis represents the HR value. Each gene displays a 95% confidence interval. Genes are generally considered tumor suppressors when HR < 1, and oncogenes when HR > 1. E: Kaplan-Meier survival plot for the training set: Blue represents the low-risk group, yellow represents the high-risk group, and the shadow indicates the 95% confidence interval. F: Kaplan-Meier survival plot for the test set: The color scheme is the same as in Fig. 5E. Because the cutoff value is based on the training set, the number of samples in the high and low-risk groups differs.

Fig. 5

Series of Plots Evaluating the Performance of the Prognostic Risk Model. A: 1, 3, and 5-year ROC curves for riskScore. B: Risk factor correlation plot. The top diagram shows each ccRCC sample ordered from smallest to largest by riskScore. The middle diagram plots the prognostic outcomes of samples with riskScore on the x-axis and survival time on the y-axis. Red indicates death, blue indicates survival/loss to follow-up, which corresponds to the order of samples in the top diagram. The bottom diagram is a heatmap matrix of model genes in the samples, where red represents high expression, and blue indicates low expression, which also corresponds to the top and middle diagrams. C: Multi-ROC curve plot that combines age, sex, and TNM stage. D: Univariate prognostic analysis plot. E: Multivariate prognostic analysis plot. F: Sankey plot, from left to right: subtype classification, high and low-risk groups, and outcome classification. G: Box plot of riskScore among subtypes, with the y-axis representing riskScore and the x-axis representing subtype names. H: Box plot of prognostic grouping, with the y-axis representing riskScore and the x-axis representing survival outcomes. I: Stacked histogram for riskScore grouping. The y-axis represents percentage values, and the x-axis represents high and low-riskScore groups, showing the proportion of dead and living samples within each group. J: Box plot of TMB versus riskScore, showing a comparison of TMB values between high and low-riskScore groups. K: Trend plot among TMB, riskScore, and subtypes. The x-axis represents the specific values of riskScore, the y-axis represents TMB values, different colors represent samples from different subtypes, straight lines represent linear regression, and shadows represent the 95% confidence interval. L: KM survival plot for TMB. M: Joint KM survival analysis plot for TMB and riskScore.

Fig. 6

Study on the Correlation between the Prognostic Risk Model and Immune Infiltrating Cells. A: ssGSVA analysis plot between high and low-risk groups. The x-axis represents immune cells infiltrating ccRCC tissues, and the y-axis represents the relative abundance of infiltration. B: The left diagram is for the matrix score, and the right diagram is for the immune score. The x-axis represents the high and low-risk groups, and the y-axis represents the relative scores. C: Association plot of the model construction genes and immune infiltrating cells. The right-side y-axis represents the names of the genes used in model construction, and the x-axis represents the names of immune infiltrating cells. D: Correlation analysis plot of model construction genes and immune checkpoint genes: the x-axis represents the names of immune checkpoint genes. In both C and D, red indicates positive correlation, and blue indicates negative correlation.

Fig. 7

The role of riskScore in immunotherapy and anti-cancer drug sensitivity. A: Lollipop chart where the left y-axis represents immune cell types, and the right y-axis represents the p-value of correlation analysis, and the x-axis represents the correlation coefficient between the type of immune cell and the riskScore of the sample. B: Chord plot, where the size and depth of color of the red band respectively represent the strength of correlation between the modules. C–F: Box plots of risk model for immune checkpoint genes. Blue and red represent the low-risk and high-risk groups, respectively, and the y-axis represents the gene expression level. G–J: Box plots of drug sensitivity analysis risk model. Green and red represent the low-risk and high-risk groups, respectively, and the y-axis represents the drug sensitivity coefficient.