The TGF-β pathway genes were significantly differentially expressed in ccRCC and were closely related to prognosis
From KEGG TGF-β signaling pathway of GSEA website, we found 85 genes related to TGF-β signaling pathway (Table S1), and analyzed the RNA sequencing data of 539 patients with human renal clear cell carcinoma and 72 normal kidney samples. We found that 76 genes of TGF-β pathway were dysregulated in ccRCC, such as TGFB1, SMAD3, SMAD4, BMP2 and so on (Figure 1A, Table S2). Furthermore, we analyzed the effects of these differentially expressed genes related to TGF-β pathway on the overall survival of patients with clear cell renal cell carcinoma. The results showed that 55 genes could significantly affect the prognosis of patients with ccRCC, of which 19 genes were risk factors for the prognosis of patients with ccRCC (risk ratio > 1) and 36 genes were protective factors (risk ratio < 1) (Figure 1B, Table S3).
Unsupervised hierarchical clustering and Prognosis analysis
The unsupervised hierarchical clustering of the data revealed four different clusters of TGF-β pathways (Figure 2A). Cluster “1” showed that the expression level of 85 genes was similar to that of normal samples, indicating that the TGF-β pathway was in a normal state. The expression level of genes related to TGF- β pathway in cluster “2” increased, indicating that TGF-β scores are high, while cluster “3” and cluster “4” showed the low expression of TGF- β pathway genes, indicating that TGF-β scores are low. Figure 2B depicts the level of TGF-β scores in four clusters and further shows the differences in TGF-β scores among them. Table S4 showed the TGF-β score of each sample. However, subsequent survival analysis showed that the prognosis of cluster “3” was the worst and that of cluster “4” was the best (Figure 2C). Then we further classified all the samples into three groups, the “normal” group, the “KEGG-TGF-β high score group” and the “KEGG-TGF-β low score group”. That is, clusters “3” and “4” were merged into “KEGG-TGF-β low score group”. Then we analyzed the correlation between TGF-β score and clinicopathological features, and the results showed that the expression of TGF-β pathway genes was significantly correlated with fustat (Figure 2D). Further survival analysis showed that the prognosis of "KEGG-TGF-β low score group" group was the worst and the “KEGG-TGF-β high score group” group was the best (Figure 2E).
The disruption of TGF-β signaling pathway is closely related to the dysregulation of key genes and deacetylated genes in ccRCC
We explored the relationship between the expression of a variety of well-known key genes and TGF-β pathway in ccRCC (Figure 3A). VHL, TP53 and PTEN were down-regulated in the “KEGG-TGF-β low score group”. These results suggested that the destruction of TGF-β signal pathway was related to the promotion of tumor. EGFR, MYC and VEGFA and other oncogenes were highly expressed in “KEGG-TGF-β high score group”, and it may be more effective to target these genes in KEGG-TGF-β high score group. The analysis of transcriptome of ccRCC patients also showed that there was a strong correlation between the abnormal expression of sirtuin and HDAC and the abnormality of TGF-β pathway (Figure 3B).
Prediction of IC50 of different targeted drugs based on the TGF-β score
Considering that targeted drugs are commonly used in the treatment of metastatic renal cell carcinoma, we tried to evaluate the efficacy of different targeted drugs in TGF-β “KEGG-TGF-β high score group” and “KEGG-TGF-β low score group”. We have obtained a satisfactory prediction model by ridge regression on the GDSC cell line data set. Based on this model, we evaluated the IC50 of 12 targeted drugs. The results of the analysis showed that there was no significant difference in the IC50 of Pazopanib, Gefitinib, Bosutinib, and Lapatinib among the three groups (Figure 4A, 4D, 4F, 4H and 4M). Compared with KEGG-TGF-β high score group and KEGG-TGF-β low score group, the IC50 of Temsirolimus and Sunitinib in normal group was lower (Figure 4E and 4K). It is recommended that these two drugs be used in normal group. Compared with normal group and KEGG-TGF-β low score group, the IC50 of Imatinib, Nilotinib and Axitinib in KEGG-TGF-β high score group was lower (Figure 4B, 4C and 4L). These three drugs may have better results in KEGG-TGF-β high score group. Compared with normal group and KEGG-TGF-β high score group, the IC50 of Metformin, Tipifarnib and Sorafenib in KEGG-TGF-β low score group was lower (Figure 4G, 4I and 4J). These three drugs may have better results in KEGG-TGF-β low score group.
TGF-β pathway is related to immune regulation
Immune regulation plays an important role in the occurrence and development of tumors. Many genes related to TGF-β signaling pathway were associated with the infiltration of many kinds of immune cells, especially TNF, TGFB1 and IFNG genes (Figure 5A). Then we analyzed the correlation between TGF-β score and immune infiltration, and found that there was a close relationship between TGF-β score and immune infiltration(Figure 5B), especially Type_Ⅱ_IFN-responce and mast cells, which were positively correlated with TGF-β score(Figure 5C and 5D). Next, we analyzed the correlation between TGF-β score and immune checkpoints, and the results showed that TGF-β score was negatively correlated with CTLA4 and PDCD1(Figure 5E). We used the TIDE algorithm to predict the response of KEGG-TGF-β low score group and KEGG-TGF-β high score group to immune checkpoint inhibitors, but after Bonferroni correction, we did not find that there was a significant difference in the response of KEGG-TGF-β low score group and KEGG-TGF-β low score group to immune checkpoint inhibitors (Figure 5F).
Construction and verification of the new TGF-β-based survival model
Firstly, TGF-β genes related to survival were screened out according to survival analysis and P < 0.05. Then the lasso regression model was used to analyze and determine the most reliable prognostic markers. On this basis, five genes ZFYVE9, ACVR2A, IFNG, AMH and THBS3 were selected to establish a risk characteristic model based on the minimum criterion (Figure 6A and 6B). Then, according to the median risk score, we divided patients with ccRCC into low risk group and high risk group, and studied the prognostic performance of the new survival model composed of 5 genetic risk characteristics. Kaplan-Meier survival curve analysis showed that the survival rate of patients in the high-risk group was significantly lower than that in the low-risk group (Figure 6C). In addition, we also performed ROC curve analysis to analyze the prognostic performance of the new survival model in patients with ccRCC. The AUC value of 5-year survival rate was 0.728, and that of 7-year survival rate was 0.744. The AUC value of 10-year survival rate was 0.752. It is suggested that the prognostic model of ccRCC based on TGF-β pathway gene has high predictive value (Figure 6D, 6E and 6F). In order to better understand the relationship between TGF-β pathway genes and clinicopathological characteristics of ccRCC patients, we systematically analyzed the correlation between risk scores based on five TGF-β pathway genes and clinicopathological characteristics of ccRCC patients in TCGA database. We observed that the risk score was closely related to the clinicopathological features of patients with renal clear cell carcinoma, such as T (tumor size), M (tumor metastasis), tumor grade, tumor stage and fustat(Figure 6G). Univariate Cox regression analysis showed that age, grade, stage, tumor size, tumor metastasis and risk score were associated with Overall survival in patients with ccRCC (Figure 6H, Table S5). Multivariate Cox regression analysis showed that age, grade, stage and risk score were independent risk factors that affect the prognosis of ccRCC patients (Figure 6I, Table S6). We used nomogram to predict the risk of ccRCC patients. Five-year, seven-year and ten-year survival rates can be estimated based on the patient's age, grade, stage, and risk score (Figure 6J).
TGF-β pathway genes had a wide range of genetic changes across cancer types and affected the prognosis of many cancers
In order to further explore the genetic changes of TGF-β pathway genes in pan-tumor, we analyzed the single nucleotide variation and gene expression changes of TGF-β pathway genes across cancer types. The results showed that TGF-β pathway genes had a wide range of SNV (Figure 7A, Table S7) and gene expression differences (Figure 7B, Table S8) across cancer types. Then, we analyzed the effect of TGF-β pathway genes on the prognosis of cancer patients (Figure 7C, Table S9), and the results showed that most of TGF-β pathway genes were risk factors in other types of tumors, but they were very different in renal cell carcinoma. Most of TGF-β pathway genes were protective factors in KIRC (Kidney renal clear cell carcinoma), which was consistent with our previous analysis: the prognosis of KEGG-TGF-β low score group was worse.
Connectivity map analysis identified potential compounds / inhibitors targeting TGF-β
Considering that most TGF-β pathway genes are risk factors in tumors, we try to find compounds that can inhibit TGF-β pathway genes. We used a data-driven systematic method, Connectivity Map (CMAP), to explore the relationship between genes, compounds, and biological conditions to find candidate compounds that may target TGF-β pathway genes. We found that 54 compounds are enriched in different tumors, and they can inhibit TGF-β pathway genes (Figure 8A). At the same time, we explored the possible action mechanisms of 19 small molecular compounds, we found that the compounds involved 18 mechanisms, and two compounds had the same mechanism (Figure 8B). Therefore, we can select different compounds in different tumors and suppress TGF-β pathway genes according to different mechanisms, so as to achieve the purpose of tumor inhibition.