Analysis of differential gene immune infiltration and clinical characteristics of skin cutaneous melanoma based on systems biology and drug repositioning methods to identify drug candidates for skin cutaneous melanoma.

The genes related to the pathogenesis of SKCM obtained in this study were obtained from three public databases (KEGG, OMIM, and Genotype). Fifteen datasets ( GSE100797 , GSE78220 , GSE133713 , GSE98394 , GSE190113 , GSE22153 , GSE54467 , GSE53118 , GSE22154 , GSE19234 , GSE59455 , GSE46517 , GSE99898 , GSE15605 , GSE6590293) were used to obtain SKCM expression profile data, obtain clinical information from the TCGA SKCM dataset, obtain normal tissue expression profiles from the GTEx database, and download 684 FDA-certified drugs and corresponding 7281 protein targets from Drugbank database ( https://go.drugbank.com/stats (version 5.1.3, released on April 2, 2019)). Protein interaction data was integrated from three databases (STRING database( https://string-db.org/ ), PNA database ( https://omics.bjcancer.org/pina/home.action ), and HuRI database ( http://www.interactome-atlas.org/ ), by integrating protein interaction data from three human protein interaction databases; 52,857,362 pairs of interactions with high confidence were obtained. In this study, we obtained SKCM pathogenesis-related pathway genes (map05218) from the KEGG database ( https://www.genome.jp/kegg/ ) and from the OMIM database ( https://www.omim .org/) search to obtain skin melanoma-related genes, and skin melanoma-related genes from the Genotype database ( https://www.ncbi.nlm.nih.gov/gap/phegeni ), and the genes from the above three databases were obtained. Union and mapping to HGNC human genes. Then, based on the STRING, PINA, and HuRI databases, a human protein-protein interaction (PPI) was constructed for the acquired genes related to the pathogenesis of SKCM. Finally, GO and KEGG analysis of the genes related to the pathogenesis of SKCM were performed based on the clusterProfiler R package (Yu et al. 2012 ). Taking “adj p 1” as the thresholds for differential gene screening, TCGA combined with GTEX database, GSE3189 , and GSE15605 were used to screen SKCM differential genes using R software Limma package, and then SKCM pathogenesis-related genes were screened. The intersection of SKCM and SKCM differential genes was obtained to obtain SKCM pathogenesis-related differential genes, and K-M and Cox regression analysis of SKCM pathogenesis-related differential genes was performed. Based on the results of multivariate Cox regression analysis, we established a prognostic nomogram model of differential genes related to the pathogenesis of SKCM, evaluated the prognostic value of the prognostic model, and analyzed the correlation between risk scores and immunity. Based on the prognostic analysis of differential genes related to the pathogenesis of SKCM, it was found that ZC3H12A affected the prognosis of SKCM patients in both K-M analysis (OS, DSS) and multivariate Cox regression analysis, revealing that ZC3H12A is an important factor affecting the progression of SKCM patients. Therefore, to verify the importance of ZC3H12A in SKCM, we analyzed the expression and prognosis of ZC3H12A in SKCM, and analyzed the potential function of ZC3H12A based on the clusterProfiler R package. Finally, the effect of ZC3H12A on TME was analyzed. To identify drugs for the treatment of SKCM, we used the protein interaction network method to identify drugs for the treatment of SKCM. First, based on the PPI network constructed by 294 SKCM pathogenesis-related genes, the network distance between the drug-corresponding targets and skin melanoma-related genes in the Drugbank database was analyzed, to quantify the relationship between each drug and SKCM. Then, we use formula 1 (Peng et al. 2020 ) to calculate the distance between any drug target and skin melanoma-related genes. For any drugs, record the corresponding target of S in the Drugbank database as T , w is the weight of the target gene, w  =  − ln ( D  + 1), D is the degree of the gene in the PPI network. 1 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$d\left(S,T\right)=\frac{1}{\left|T\right|}\sum\nolimits_{t\in T}{min}_{s\in S}(d\left(s,t\right)+w)$$\end{document} d S , T = 1 T ∑ t ∈ T min s ∈ S ( d s , t + w ) To assess the significance of each drug, we generated a distribution of reference distances for each drug. For any drug S , the target is T , we randomly select any gene R equal to T from the network, and then calculate the network distance d ( S , R ) between the random gene combination and SKCM-related genes. This process is repeated 10,000 times to calculate the mean \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\mu }_{d(S,R)}$$\end{document} μ d ( S , R ) and standard deviation \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\sigma }_{d(S,R)}$$\end{document} σ d ( S , R ) of the reference distance distribution generated by the random process, and the proximity value is calculated according to formula 2 (Peng et al. 2020 ). 2 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$z(S,T)=\frac{d\left(S,T\right)-{\mu }_{d(S,R)}}{{\sigma }_{d(S,R)}}$$\end{document} z ( S , T ) = d S , T - μ d ( S , R ) σ d ( S , R ) Finally, based on the reverse gene set enrichment analysis of drug features to identify potential drugs and drug feature analysis, according to the concept of pharmacogenomics, drugs that are effective for diseases should be able to affect the expression of disease genes. We downloaded drug treatment-related gene expression profiling data from the LINCS database, and each dataset was analyzed by GSEA (an enrichment score (ES) and p -value to evaluate the statistical significance of the drug. Then, the Benjamini and Hochberg algorithm was used for p -values. Values were corrected for multiple testing to obtain FDR. All drugs with FDR < 0.25 were considered to be significantly associated with cutaneous melanoma. The Wilcoxon rank sum test was used to compare variables between the two groups. The correlation between different genes was analyzed by Spearman’s or Pearson’s test, and the correlation coefficient was calculated. The K-M survival curve was used to determine the difference in survival rate, and the log-rank test was used to evaluate the significant difference in survival time between the two groups, determining independence using univariate and multivariate Cox analysis. p  < 0.05 was considered statistically significant. In this study, we used R v4.0.3 for data analysis and ggplot2 for visualization.

We obtained 72 pathway genes (map05218) related to the pathogenesis of skin melanoma from the KEGG database and obtained the pathogenesis of skin melanoma from the OMIM database search, There are 190 genes, and 43 genes related to skin melanoma were obtained from the Genotype database. Based on 3 databases (KEGG, OMIM, Genotype), we obtained 294 genes related to the pathogenesis of the SKCM gene (Fig.  1 A); based on the TCGA SKCM dataset, the expression analysis of 294 SKCM pathogenesis-related genes was carried out, and it was found that 294 genes were highly heterogeneous in TCGA SKCM patients, and there were significant differences in expression in patients with different stages (Fig.  1 B). Subsequently, a protein interaction network was constructed based on 294 SKCM pathogenesis-related genes, and it was found that the 294 SKCM pathogenesis-related genes were closely related to each other to form clusters (Fig.  1 C). It is suggested that the related genes that lead to the pathogenesis of skin melanoma interact with each other and participate in the related biological processes regulating the pathological changes of SKCM. To further study the biological functions involved in the regulation of these genes, we used the clusterProfiler package to perform GO and KEGG functional enrichment analysis and found that they were involved in the regulation of various cancer-related pathways and kinase activities (Fig.  1 D). Fig. 1 Acquisition and functional analysis of skin melanoma-related genes. A Based on 3 databases of KEGG, OMIM, and Genotype, 294 skin melanoma-related genes were obtained. B Heatmap of skin melanoma-related genes. C 294 SKCM pathogenesis-related gene-protein interaction networks. Red ovals represent skin melanoma-related genes in the OMIM database, green diamonds represent skin melanoma-related genes in the KEGG pathway, and purple triangles represent skin melanoma-related genes in the PheGenI database. D Shows the function bubble graph in which SKCM-related genes are significantly enriched. The color of the nodes represents significance, and the size of the nodes represents the number of genes enriched into functions. The top ten significantly enriched functions corresponding to BP, MF, CC, and KEGG, respectively Acquisition and functional analysis of skin melanoma-related genes. A Based on 3 databases of KEGG, OMIM, and Genotype, 294 skin melanoma-related genes were obtained. B Heatmap of skin melanoma-related genes. C 294 SKCM pathogenesis-related gene-protein interaction networks. Red ovals represent skin melanoma-related genes in the OMIM database, green diamonds represent skin melanoma-related genes in the KEGG pathway, and purple triangles represent skin melanoma-related genes in the PheGenI database. D Shows the function bubble graph in which SKCM-related genes are significantly enriched. The color of the nodes represents significance, and the size of the nodes represents the number of genes enriched into functions. The top ten significantly enriched functions corresponding to BP, MF, CC, and KEGG, respectively Based on TCGA combined with the GTEX database, 3187 SKCM differentially upregulated genes, 1933 SKCM differentially downregulated genes were found, and GSE3189 found 1395 SKCM differentially upregulated genes, 1884 SKCM differentially downregulated genes, GSE15605 differentially upregulated genes, and 1449 SKCM differentially downregulated genes. Subsequently, after the intersection of SKCM differential genes and related genes, 9 upregulated SKCM differentially related genes ( SLC45A2 , LEF1 , MITF , TYR , TNFRSF21 , PARP1 , STMN1 , E2F3 , BAX ) and 9 downregulated SKCM differentially related genes ( EGFR , PLAGL1 ), ZC3H12A , PERP , FILIP1L , EDN1 , FGFR3 , ACKR1 , PLA2G2A ) were found (Fig.  2 A, B). Correlation analysis found that 18 genes were more frequently positively correlated than negatively correlated in SKCM (Fig.  2 E). K-M survival analysis (OS, DSS) was performed on 18 SKCM differentially related genes, and 7 genes were found to have prognostic significance (Fig.  2 C, D), while multivariate COX regression analysis found that 3 genes ( EGFR , TNFRSF21 , ZC3H12A ) were for independent prognostic factors in SKCM patients (Fig.  2 F, G). Finally, the correlation between the expression of 18 differentially related genes and the level of immune cell infiltration was explored based on the TIMER2.0 online database (Li et al. 2020 ); it was found that most upregulated related genes were negatively correlated with the infiltration abundance of anti-tumor immune cells (B cells, CD8 + T cells, CD4 + T cells). while most downregulated genes were positively correlated with the infiltration abundance of anti-tumor immune cells (B cells, CD8 + T cells, CD4 + T cells) (Supple Fig. 1 (D)). Correlation analysis with T cell signature genes, immune cell marker genes, and immune checkpoint-related genes found that most upregulated genes were negatively correlated with the above genes, while most downregulated genes were positively correlated with the above genes (Supple Fig. 1  (A-C)). Fig. 2 Identification of SKCM-related differential genes and prognostic analysis ( A , B ). Identification of SKCM-related differential genes ( C , D ). KM prognostic analysis (OS, DSS) of 18 SKCM-related differential genes. E Correlation analysis of 18 SKCM-related differential genes. F , G COX regression analysis of 18 SKCM-related differential genes Identification of SKCM-related differential genes and prognostic analysis ( A , B ). Identification of SKCM-related differential genes ( C , D ). KM prognostic analysis (OS, DSS) of 18 SKCM-related differential genes. E Correlation analysis of 18 SKCM-related differential genes. F , G COX regression analysis of 18 SKCM-related differential genes A prognostic nomogram-based model was developed based on 3 independent prognostic genes (EGFR, TNFRSF21, ZC3H12A), independent prognostic indicators for T, M, and N stages (Fig.  3 A). The C-index was 0.732, and the ROC analysis showed that the AUC values for 2, 4, and 6 years were 0.823, 0.795, and 0.792, respectively (Fig.  3 C). Correspondingly, the C-index of the T, M, and N staging models was 0.703, and the AUC values of the 2-year, 4-year, and 6-year models were 0.795, 0.768, and 0.760, respectively (Fig.  3 E, F). This suggested that the new model had higher predictive accuracy. Patients with SKCM were categorized into low-risk and high-risk groups according to the model risk score thresholds. The Kaplan–Meier analysis showed significant differences in survival between the two risk groups, with patients in the low-risk group having a longer survival duration (Fig.  3 B). Subsequently, the associations of risk scores and pathological characteristics of SKCM patients were evaluated. Significant correlations between risk scores with pathological stage, T stage, N stage, melanoma ulceration, melanoma Clark level, and Breslow depth were observed (Fig.  3 G). In addition, according to the results of univariate and multivariate Cox regression analyses, the risk score was found to be an independent prognostic factor for SKCM (Supplementary Table 1 ). Correlational analysis suggested that risk score was negatively associated with the majority of antitumor immune cells (B cells, CD8 + T, and CD4 + T cells), while negatively correlated with AIM2 expression (Fig.  4 A). The ESTIMATE algorithm results suggested the risk score was negatively correlated with ESTIMATEScore, ImmuneScore, and StromalScore (Fig.  4 B). Eight immune checkpoint marker genes, type I antitumor response genes (IFNG, GZMB), and major immune activation cell marker genes were negatively correlated with risk scores (Fig.  4 E, F). Finally, the immune cell infiltrates in low-risk and high-risk groups were also examined. The results indicated a substantial increase in the level of infiltration of multiple immune cells in the patients belonging to the low-risk group (Fig.  4 C), At the same time, we also used the ESTIMATE algorithm to calculate the ESTIMATEScore and the immune score StromalScore between high- and low-risk groups. (Fig.  4 D). There is further strong support for a low-risk group having a role in antitumor immunity. Fig. 3 Prognostic value of the prognostic nomogram model. A Prognostic nomogram. B OS and DSS KM curves of high- and low-risk factor groups. C Prognostic model risk factor ROC curve. D Correlation between expression of 3 independent prognostic genes (EGFR, TNFRSF21, ZC3H12A) and risk factors. E , F T, M, N stage prognostic model risk factor ROC curve. G Association of risk factors with clinicopathological Fig. 4 Correlation of risk factors with immune infiltrating cells. A Bar graph showing the correlation of risk factors with the level of immune cell infiltration. B Correlation of risk factors with immune score, matrix score, and Estimate score. C Differences in the infiltration of 24 immune cells in the high- and low-risk factor groups. D Differences in immune score, matrix score, and Estimate score between high-and low-risk factor groups. E Association of risk factors with immune cell marker genes. F Association of risk factors with 8 immune checkpoints Prognostic value of the prognostic nomogram model. A Prognostic nomogram. B OS and DSS KM curves of high- and low-risk factor groups. C Prognostic model risk factor ROC curve. D Correlation between expression of 3 independent prognostic genes (EGFR, TNFRSF21, ZC3H12A) and risk factors. E , F T, M, N stage prognostic model risk factor ROC curve. G Association of risk factors with clinicopathological Correlation of risk factors with immune infiltrating cells. A Bar graph showing the correlation of risk factors with the level of immune cell infiltration. B Correlation of risk factors with immune score, matrix score, and Estimate score. C Differences in the infiltration of 24 immune cells in the high- and low-risk factor groups. D Differences in immune score, matrix score, and Estimate score between high-and low-risk factor groups. E Association of risk factors with immune cell marker genes. F Association of risk factors with 8 immune checkpoints In view of the prognostic analysis of 18 SKCM-related differential genes (K-M (OS, DSS), multivariate COX regression analysis), ZC3H12A was found to be an important factor affecting the prognosis of SKCM (Fig.  5 C). Therefore, first, based on TCGA combined with GTEX database, GSE3189 , GSE15605 , GSE46517 , and CCLE database to explore the expression of ZC3H12A mRNA found that ZC3H12A was in low expression (Fig.  5 A, B); based on GSE190113 and GSE46517 , it was found that ZC3H12A expression in primary tumor tissue was higher than that in metastatic tumors. However, the clinicopathological staging analysis found that patients with low ZC3H12A expression had tumor aggressive characteristics in SKCM patients, ZC3H12A showed high expression in M0, N0, and Stage I stages, and with tumor invasion, its expression gradually decreased (Fig.  5 E), revealing that ZC3H12A may be an important protective factor of SKCM. Fig. 5 Differential expression and clinicopathological characteristics of ZC3H12A in SKCM. A The expression of ZC3H12A in 54 SKCM cell lines in the CCLE database. B Differential expression of ZC3H12A in TCGA SKCM dataset, GSE3189 , GSE46517 , GSE15605 . C Identification of ZC3H12A as an important factor affecting the prognosis of SKCM. D Differential expression of ZC3H12A in primary and metastatic SKCM. E The relationship between ZC3H12A expression and M, N, and Stage stages in SKCM patients Differential expression and clinicopathological characteristics of ZC3H12A in SKCM. A The expression of ZC3H12A in 54 SKCM cell lines in the CCLE database. B Differential expression of ZC3H12A in TCGA SKCM dataset, GSE3189 , GSE46517 , GSE15605 . C Identification of ZC3H12A as an important factor affecting the prognosis of SKCM. D Differential expression of ZC3H12A in primary and metastatic SKCM. E The relationship between ZC3H12A expression and M, N, and Stage stages in SKCM patients To further test our hypothesis, the prognostic value of ZC3H12A was explored by K-M analysis on the TCGA SKCM dataset and GSE46517 , GSE54467 , GSE53118 , and GSE59455 ; it was found that with the increase of ZC3H12A expression, the prognosis of SKCM patients was significantly improved (Fig.  6 A–C), the prognostic value of ZC3H12A was subsequently explored by COX regression analysis on the TCGA SKCM dataset. It was found that ZC3H12A was an independent prognostic protective factor in SKCM patients (Fig.  6 D, E), which was consistent with that observed in SKCM patients. The association of lower ZC3H12A expression with more aggressive tumor characteristics was consistent. It was confirmed that ZC3H12A was an important protective factor of SKCM, but further experimental verification is needed. At the same time, the prognostic nomogram model constructed based on age and ZC3H12A had better predictive performance (Fig.  6 F, G). Taken together, these data collectively suggest that ZC3H12A was a potent and independent biomarker for early SKCM. Fig. 6 ZC3H12A is an important independent factor in SKCM patients and improves patient prognosis. A K-M method ZC3H12A is a good prognostic factor in GSE46517 , GSE54467 , GSE53118 , and GSE59455 datasets. B , C K-M method ZC3H12A is a good prognostic factor in the TCGA SKCM dataset. D , E Cox regression analysis ZC3H12A was an independent factor for improved prognosis in the TCGA SKCM dataset. F , G ZC3H12A prognostic nomogram was constructed based on the TCGA SKCM dataset ZC3H12A is an important independent factor in SKCM patients and improves patient prognosis. A K-M method ZC3H12A is a good prognostic factor in GSE46517 , GSE54467 , GSE53118 , and GSE59455 datasets. B , C K-M method ZC3H12A is a good prognostic factor in the TCGA SKCM dataset. D , E Cox regression analysis ZC3H12A was an independent factor for improved prognosis in the TCGA SKCM dataset. F , G ZC3H12A prognostic nomogram was constructed based on the TCGA SKCM dataset To characterize the biological characteristics of the ZC3H12A gene, GO, KEGG, and GSEA function enrichment analysis was performed on ZC3H12A based on the TCGA SKCM dataset using the clusterProfiler package of R language. The results showed that ZC3H12A was related to immune response, indicating that ZC3H12A may be involved in immune or inflammatory response (Supple Fig.  2 (A-D), Fig.  7 A). Then, in order to further verify the involvement of ZC3H12A in immune response, the correlation and differential expression between ZC3H12A and type I anti-tumor response genes ( IFNG , GZMB ) were first discussed based on the TCGA SKCM dataset; ZC3H12A was positively correlated ( p  < 0.01), and it was expressed in high- and low-expression groups of ZC3H12A high-expression group (Fig.  7 B, C), revealing that ZC3H12A was involved in the immune response of SKCM, and then based on TMER and CIBERSORT algorithms to analyze various immune responses between ZC3H12A high- and low-expression groups. Cell infiltration abundance was different, and it was found that most of the immune cell infiltration levels were significantly different between ZC3H12A high- and low-expression groups (Fig.  7 D, E). Among them, the infiltration level of CD4 memory-activated T cells and CD8 + T cells was higher in the ZC3H12A high-expression group. It is further revealed that ZC3H12A is involved in the immune response of SKCM patients and is sensitive to immunotherapy. Finally, the correlation between ZC3H12A expression and the infiltration abundance of various immune cells in the tumor microenvironment was analyzed based on 8 algorithms in 15 SKCM datasets in the GEO database and TCGA SKCM dataset. The results showed that with the increased expression of ZC3H12A, the infiltration level of anti-tumor immune cells (DC, CD4 + T cells, CD8 + Tcells, B cells) increased; especially, the infiltration of B cells and DC cells increased significantly, while the tumor-promoting immune cell ( MODS, Treg, CAF) infiltration levels did not change or even partially decreased (Fig.  7 A). ESTIMATE algorithm analysis found that ZC3H12A was positively correlated with ESTIMATEScore, immune score, and StromalScore except for GSE98394 , and the ESTIMATEScore, immune score, and StromalScore in the ZC3H12A high-expression group were higher than those in the ZC3H12A low-expression group (Fig.  7 A, F). The correlation and difference analysis of immune cell marker genes and T cell signature genes were carried out based on TCGA SKCM, and the relationship between ZC3H12A expression and chemokine and chemokine receptors was investigated based on 15 SKCM datasets in the GEO database and TCGA SKCM datasets. Based on the correlation of expression levels of the body, MHC molecules, immunosuppressive factors, and immune activators, it was found that ZC3H12A was positively correlated with almost all anti-tumor immune cell marker genes and T cell signature genes, and was significantly expressed in the ZC3H12A high-expression group (Fig.  7 B, C; Fig.  8 C, D), while in the GEO database, 15 except for a few MHC molecules in GSE133713 and GSE98314 datasets were positively correlated with ZC3H12A expression; the remaining 14 datasets found that MHC molecules were positively correlated with ZC3H12A expression, which revealed that the increase of ZC3H12A expression in SKCM patients increased antigen presentation and processing capacity. At the same time, it was found that the expressions of chemokines, chemokine receptors, and immune activators were positively correlated with ZC3H12A expression in most datasets, while immunosuppressive factors were positively correlated with ZC3H12A expression in a few datasets (Fig.  7 G). It was further confirmed that ZC3H12A is involved in the immune response in SKCM and can promote anti-tumor immunity. In addition, ZC3H12A expression was investigated based on the expression data of 6 SKCM single cells in the TISCH database (Sun et al. 2021 ). Interestingly, in the six SKCM single-cell datasets analyzed, ZC3H12A was expressed in antitumor immune cells (proliferating T cells, B cells, CD4 + T cells, DC cells) (F i g.  8 I), which further confirmed that ZC3H12A was involved in immune response and anti-tumor function. Fig. 7 Comprehensive analysis of ZC3H12A and immune cell infiltration in SKCM. A Heat map showing the effect of ZC3H12A on immune cell infiltration levels in 9 algorithms. B , C The correlation and differential expression of ZC3H12A with type I antitumor response genes (IFNG, GZMB) and immune cell marker genes. D The CIBERSORT algorithm analyzed the differences in the infiltration abundance of various immune cells between ZC3H12A high- and low-expression groups. E TMER algorithm to analyze the differences in the infiltration abundance of various immune cells between ZC3H12A high-and low-expression groups. F ESTIMATE algorithm ZC3H12A high- and low-expression groups ESTIMATEScore, immune score, StromalScore. G Heat map showing the correlation of ZC3H12A expression with the expression levels of chemokines, chemokine receptors, MHC molecules, immunosuppressive factors, and immune activators Fig. 8 Association of ZC3H12A with response to immunotherapy. A GSEA enrichment analysis showed that immune activation-related signals (B cell receptor signaling pathway, T cell receptor signaling pathway, chemokine signaling pathway,) were significantly enriched in patients with high ZC3H12A expression. B Differential expression of tumor stemness index between high- and low-ZC3H12A-expression groups and normal tissues. C , D The correlation and differential expression of ZC3H12A and T cell signature genes. E , F Correlation and differential expression of ZC3H12A with 8 predictors of immunotherapy response. G , H Differential expression of ZC3H12A in melanoma immunotherapy cohorts Homet cohort 2019 Anti-PD-1 and Riaz cohort 2018 Anti-PD-1/CTLA-4 ICB-responsive and non-responsive patients. I ZC3H12A expression in 6 SKCM single-cell expression datasets. J Analysis of the effect of ZC3H12A on prognosis in the Nathanson_2017_anti-CTLA-4, phs000452_anti-PD-1, and PRJEB23709_anti-PD-1 cohorts. K AUC values of ZC3H12A in 13 melanoma ICB-treated cohorts Comprehensive analysis of ZC3H12A and immune cell infiltration in SKCM. A Heat map showing the effect of ZC3H12A on immune cell infiltration levels in 9 algorithms. B , C The correlation and differential expression of ZC3H12A with type I antitumor response genes (IFNG, GZMB) and immune cell marker genes. D The CIBERSORT algorithm analyzed the differences in the infiltration abundance of various immune cells between ZC3H12A high- and low-expression groups. E TMER algorithm to analyze the differences in the infiltration abundance of various immune cells between ZC3H12A high-and low-expression groups. F ESTIMATE algorithm ZC3H12A high- and low-expression groups ESTIMATEScore, immune score, StromalScore. G Heat map showing the correlation of ZC3H12A expression with the expression levels of chemokines, chemokine receptors, MHC molecules, immunosuppressive factors, and immune activators Association of ZC3H12A with response to immunotherapy. A GSEA enrichment analysis showed that immune activation-related signals (B cell receptor signaling pathway, T cell receptor signaling pathway, chemokine signaling pathway,) were significantly enriched in patients with high ZC3H12A expression. B Differential expression of tumor stemness index between high- and low-ZC3H12A-expression groups and normal tissues. C , D The correlation and differential expression of ZC3H12A and T cell signature genes. E , F Correlation and differential expression of ZC3H12A with 8 predictors of immunotherapy response. G , H Differential expression of ZC3H12A in melanoma immunotherapy cohorts Homet cohort 2019 Anti-PD-1 and Riaz cohort 2018 Anti-PD-1/CTLA-4 ICB-responsive and non-responsive patients. I ZC3H12A expression in 6 SKCM single-cell expression datasets. J Analysis of the effect of ZC3H12A on prognosis in the Nathanson_2017_anti-CTLA-4, phs000452_anti-PD-1, and PRJEB23709_anti-PD-1 cohorts. K AUC values of ZC3H12A in 13 melanoma ICB-treated cohorts Based on the above analysis revealing that ZC3H12A was involved in immune response and was associated with anti-tumor immunity, it was speculated that ZC3H12A was involved in immune checkpoint blockade (ICB) response. To test our hypothesis, first, based on the hot tumor signature genes ( CCL5 , CD8A , PDCD1 , CD8B , CXCR3 , CXCL9 , CXCL10 , CD4 , CD3E , CXCL11 , CD274 , and CXCR4 ) (Dong et al. 2021 ), an unsupervised clustering method was used to group. In the TCGA database, SKCM patients were divided into two clusters of cold and hot tumors (Fig.  9 A, B). ZC3H12A showed high expression in hot tumor clusters, and ZC3H12A prognosis (OS, DSS) was significantly improved in hot tumor clusters (Fig.  9 C–E), indicating that it is associated with ICB response. Subsequently, we analyzed the association of ZC3H12A with 8 predictors of immunotherapy response ( SIGLEC15 , IDO1 , CD274 , HAVCR2 , PDCD1 , CTLA4 , LAG3 , PDCD1LG2 ) (Li et al. 2016 ); 8 immunotherapy response predictors were found to be positively correlated with ZC3H12A ( p  < 0.0001, Fig.  8 F) and were highly expressed in the ZC3H12A high-expression group (Fig.  8 E). The above analysis results strongly suggest that ZC3H12A may promote the therapeutic effect of tumor ICB. To verify its ICB efficacy in SKCM, based on the analysis of Homet cohort 2019 Anti-PD-1 and Riaz cohort 2018 Anti-PD-1/CTLA-4 in the melanoma immunotherapy cohort, it was found that ZC3H12A expression was significantly higher in ICB-treated patients than in ICB-treated patients. Non-responding patients, and in the Nathanson_2017_anti-CTLA-4, phs000452_anti-PD-1and, and PRJEB23709_anti-PD-1 cohorts, the prognosis of patients with high ZC3H12A was significantly improved (Fig.  8 G, H, J)), confirming that patients with high ZC3H12A expression responded significantly to ICB treatment. In addition, we also found that the AUC values of ZC3H12A in 13 melanoma ICB treatment cohorts were all greater than 0.5 (Fig.  8 K), revealing that ZC3H12A may be an indicator to predict the effect of SKCM immunotherapy. Fig. 9 Differential expression, prognosis, and drug identification of ZC3H12A in thermal tumor signature gene subtypes. A , B Classification of TCGA SKCMs into hot and cold 2 clusters based on hot tumor signature genes. C – E Differential expression and prognosis of ZC3H12A between cold and hot 2 clusters. F Drug and skin melanoma distance distribution. G The predicted 43 potential drugs, their corresponding targets, and the interaction network relationship between skin melanoma-related genes. All related genes in the network have direct interactions with drug targets Differential expression, prognosis, and drug identification of ZC3H12A in thermal tumor signature gene subtypes. A , B Classification of TCGA SKCMs into hot and cold 2 clusters based on hot tumor signature genes. C – E Differential expression and prognosis of ZC3H12A between cold and hot 2 clusters. F Drug and skin melanoma distance distribution. G The predicted 43 potential drugs, their corresponding targets, and the interaction network relationship between skin melanoma-related genes. All related genes in the network have direct interactions with drug targets In recent years, studies have found that patients with high tumor stemness index have poor efficacy of ICB therapy (Zhang et al. 2022 ). In this study, we used the OCLR algorithm constructed by Malta et al. (Malta et al. 2018 ) to explore the differential expression of tumor stemness in patients with high and low expression of ZC3H12A. Compared with patients with low ZC3H12A expression, the tumor stemness index was lower in patients with high ZC3H12A expression, but they were all higher than normal skin tissue (Fig.  8 B). These results reconfirmed that the invasive ability of patients with high ZC3H12A expression was weakened and the prognosis was improved. Significant improvement, inhibiting tumor immune escape and drug resistance, and ICB treatment have a good effect. We explored the distribution of distances between drug targets and skin melanoma-related genes (Fig.  9 F), The distance distribution between drug targets and SKCM-related genes and the reference distance distribution ranged from − 10.0 to 4.0, The density curves of the actual drug and the stochastic process simulated drug peaked between 1 and 2, followed by a second peak, indicating that the distances corresponding to most drugs were around this value. Drugs that fall in the second peak are unlikely candidates because most of the two distribution curves overlap in this range. On the other hand, the density curve of the drug drops sharply around point 1, but the density of the reference set is much smaller after that point, In other words, drugs with a distance of less than 1 were more likely to be effective in SKCM treatment. Therefore, we chose this value as the threshold for screening drug candidates for cutaneous melanoma. In this way, we excluded those drugs that were not associated with cutaneous melanoma and retained a sufficient number of reliable drugs. Network analysis of drug-cutaneous melanoma relationships helped us identify potential drug-SKCM associations, From the LINCS database, we searched gene expression datasets for screened drug effects at different concentrations and cell lines. Finally, we extracted a total of 165 gene expression datasets. Through an Inverted Gene Set Enrichment Analysis (IGSEA), we identified 43 drugs with FDR < 0.05 (Fig.  9 G). The full prediction results are presented in the supplementary material (Supplementary Table 2 ). This means they had a greater effect on the expression of genes associated with melanoma in the skin tested. All of these drugs target genes that overlap with genes associated with SKCM or with known drugs. In fact, there were more indirect interactions between the targets of these 43 drugs and SKCM-related genes, and the drug and disease proximity distances calculated based on the average shortest path by Eq.  1 were statistically significant.

Despite advances in diagnosis and treatment, SKCM is currently a disease with high morbidity and mortality. In recent years, the repurposing of drugs for other diseases as potential therapeutic agents for the treatment of SKCM has led to scientific efforts (Khosravi et al. 2019 ) and widespread attention. At the same time, the prognosis was strongly dependent on the diagnostic stage, and early detection of SKCM has great potential in improving the rate of radical resection and reducing the disease burden. In this study, for the first time, we screened out the drugs related to the treatment of SKCM by using the protein interaction network method. And for the first time, ZC3H12A was identified and validated as a promising and robust biomarker in SKCM stage I patients. We also constructed a clinical prognostic nomogram model based on three SKCM pathogenesis-related differentially independent prognostic genes and T, N, and M staging, and found that the combined application of their expression and clinical nomogram could improve the prognostic ability of SKCM. ZC3H12A is a key molecule regulating immune response (Akira 2013 ), and targeting ZC3H12A can lead to fatal inflammatory diseases (Miao et al. 2013 ). Based on the analysis of ZC3H12A function and immune infiltration, it was found that ZC3H12A was involved in regulating the immune response of SKCM, remodeling the tumor immune activation microenvironment, and promoting the infiltration of anti-tumor immune cells (CD8 + T cells, B cells, CD4 + T cells). This was basically consistent with previous research reports. Furthermore, our findings also suggest that patients with lower ZC3H12A gene expression in SKCM patients have more aggressive tumor characteristics and that lower ZC3H12A expression was significantly associated with worse prognosis. The possible mechanism is that the highly expressed ZC3H12A in SKCM patients promotes the increased infiltration of anti-tumor immune cells CD8 + T cells, thereby remodeling the tumor immune activation microenvironment. However, the causal relationship between ZC3H12A and remodeling of the tumor immune-activating microenvironment remains elusive and needs to be further investigated in follow-up studies. In the pathogenesis of SKCM, MAPK, PI3K, AKT, and other signaling pathways mediate the proliferation, infiltration, and metastasis of SKCM. Oncogenic signaling plays an important role in the occurrence and metastasis of SKCM. At the same time, many drug targets for the treatment of SKCM are highly close to SKCM gene modules, and drugs usually selectively act on specific targets to exert their biological functions. In this study, through the functional enrichment analysis of the screened 294 SKCM risk genes, it was found that a variety of cancer-related signaling pathways such as MAPK, PI3K, and AKT were involved in the occurrence and development of SKCM, which was consistent with the aforementioned literature on the pathogenesis of SKCM (Huang et al. 2020 ), and further confirmed that SKCM was caused by the overactivation of multiple pathogenic signaling pathways caused by the interaction of multiple genes. This is consistent with our understanding of the molecular mechanisms of SKCM and lays the foundation for the development of new drugs for the treatment of SKCM. SKCM was a complex disease caused by abnormal gene activity caused by polygenic hypermutation. Earlier studies have shown that the genes associated with SKCM are not randomly distributed but are functional modules clustered together. By measuring the relationship between drug targets and SKCM-related genes on the PPI network, we could exclude drugs with low proximity to SKCM-related genes from the drug list extracted from the drug library and retain potential drug candidates. In this study, we chose d ( s , t ) < 1 as the threshold for measuring whether a drug is close to SKCM based on distance distribution. Although it was a computational threshold, it can more accurately screen potential drugs that may treat SKCM. Passing this threshold, 1496 potential drugs for SKCM were discovered, and in the screening of FDA-approved SKCM-targeted inhibitors and immune checkpoint inhibitors, three of them (vemurafenib (targets: CYP1A2, ABCB1, ABCG2, ABCC1, CYP2C8, CYP2B6, ORM1, BRAF), cobimetinib (targets: ABCB1, ABCG2, SLCO1B1, SLCO1B3, MAP2K1), binimetinib (targets: TNF, CYP1A2, CYP2C19, UGT1A1, IL1B, IL6, MAP3K1, MAP2K2) were included in the 1496 SKCM potential drug list, indicating that this method is effective in screening SKCM drugs. Given that SKCM was a malignancy caused by aberrant activities of multiple genes, we selected SKCM cell lines with different concentrations as a gene expression dataset for screening drugs for IGSEA analysis. Finally, we obtained 165 SKCM gene expression datasets, and we screened 43 drugs as potential drugs for the treatment of SKCM through IGSEA analysis. By reviewing the literature, we found that the screened 21 drugs have strong clinical trials and animal experimental evidence to support their potential as drugs for the treatment of SKCM (Table 1 ). Some of these drugs were originally approved for the treatment of metastatic breast cancer (e.g., palbociclib, tucatinib), but have now been found to inhibit the proliferation and apoptosis of SKCM by inhibiting the CDK4/6 pathway and inhibiting the HER2 tyrosine kinase with high selectivity (AbuHammad et al. 2019 ; Liu et al. 2018b ; Yoshida et al. 2016 ). FDA-approved drugs for the treatment of lymphatic blood system diseases such as vorinostat, prednisolone acetate, arsenic trioxide, and tazemetostat are now found to have a certain effect on SKCM. Among them, vorinostat can inhibit the proliferation of melanoma by modifying the transcription level of selenoprotein and inhibiting the activity of BRAF + MEK and histone deacetylase (Ecker et al. 2020 ; Haas et al. 2014 ; Wang et al. 2018 ); vorinostat has entered the SKCM clinical phase II trial (Haas et al. 2014 ). Arsenic trioxide has now been found that it can downregulate the tumor procoagulant activity of MCF-7 and WM-115 cell lines and induce growth inhibition and apoptosis of melanoma cell lines (Hoffman et al. 2015 ). Table 1 Molecular and supporting information of anti-melanoma drug candidates based on network group analysis sigdrug Indications Proximity FDR Targets Possible mechanism of action References Palbociclib Hormone-receptorpositive breast cancer –4.529 0 ABCB1, ABCG2, SLC22A1, SULT2A1, CDK4, CDK6 Inhibition of CDK4/6 (AbuHammad et al. 2019 ; Liu et al. 2018b ; Yoshida et al. 2016 ). Tucatinib (HER2)–positive metastatic breast cancer –4.737 0 ABCB1, ABCG2, SLC22A2, CYP2C8, SLC47A1, SLC47A2, ERBB3 Highly selective inhibition of HER2 tyrosine kinase (Murthy et al. 2020 ) Tazemetostat 1. Epithelioid sarcoma 2. Follicular lymphoma –4.656 0 ABCB1, ABCG2, SLC47A1, SLC47A2, EZH2, EZH1 1. Strong selective inhibition of histone methyltransferase EZH2 2. Inhibit some EZH2 gain-of-function mutations (including Y646X and A687V), as well as EZH1 1. Histone Methyltransferase (HMT) Inhibitor 2. Histone deacetylase (HDAC) inhibitor NA Prednisolone acetate Arthritis, blood problems, immune system disorders, cancer, and severe allergies. –4.989 0 ABCB1, SLCO1A2, NR3C1, SERPINA6 NA NA Arsenic trioxide Chronic Myelogenous Leukemia (CML), Chronic Myelogenous Leukemia, Multiple Myeloma, Acute Myelogenous Leukemia –4.67 0 CYP3A4, TXNRD1, AKT1, JUN, IKBKB, MAPK3, CCND1, MAPK1, CDKN1A, HDAC1, PML 1. Induces growth inhibition and apoptosis; 2. Downregulates cancer procoagulant activity in MCF-7 and WM-115 cell lines (Hoffman et al. 2015 ) Vorinostat Mycosis fungoides, Sézary syndrome, non-Hodgkin's lymphomas, AML\MDS –4.589 0 HDAC2, HDAC1, HDAC8, HDAC3, HDAC6, acuC1 1. Modify selenoprotein transcript levels 2. BRAF+MEK inhibitor 3. The histone deacetylase inhibitor (Ecker et al. 2020 ) (Wang et al. 2018 ) (Haas et al. 2014 ) Tixocortol Rhinitis, pharyngitis, ulcerative colitis and oral, inflammatory conditions –4.989 0 HDAC2, NR3C1 NA NA Norethynodrel Contraceptives, endometriosis and polyhedron –4.704 0 SHBG, ESR1, AR NA NA Isosorbide Coronary artery disease, angina pectoris –4.989 0 BCL2, BCL2L1, MCL1 1. Release of NO, mediates apoptosis of tumor cells by generating physiological messenger circular GMP; 2. Increase the efficacy of chemoimmunotherapies in vivo, promoting the survival and increasing the function of injected cells by targeting a key pathway incisplatin-induced cytotoxicity. (Perrotta et al. 2007 ) Afatinib Non-small-cell lung cancer –4.66 0.01 EGFR, ABCB1, ABCG2, ERBB4 1. Inhibition of tyrosine kinase receptors; 2. Downregulated pRPS6KB1 and pRPS6, downstream of mTOR signaling, and IRS-1 in the insulin signaling pathway, specifically ablating IRS-1 nuclear signal. (Das et al. 2020 ) Salinomycin Plasmodium, lung cancer –4.989 0.01 ABCB1 1. By activating endogenous K+ and exogenous H+ exchange, it causes oxidative stress, induces mitochondrial permeability transformation (MPT), and produces cytotoxic effects on cancer cells 2. Induce autophagy in melanoma cell lines (M21, M29, SK-MEL-1, SK-MEL-19, SK-MEL-103 and A375) that are sensitive to SAL. (Ohkubo et al. 2018 ) (Liu et al. 2020 ) Tezacaftor Cystic fibrosis –4.489 0.01 ABCB1, ABCG2, CFTR NA NA Baricitinib Rheumatoid arthritis, inflammatory bowel disease, psoriasis or ankylosing spondylitis; CoviD-19 treatment trial; –4.308 0.01 ABCB1, SLC22A6, ABCG2, SLCO1B3, SLC22A2, SLC22A8, SLC47A2, PTK2B, JAK1, JAK2, JAK3 Inhibits Janus kinases 1 and 2 NA Ripretinib Gastrointestinal tumors –4.217 0.01 ABCB1, ABCG2, PDGFRB, PDGFRA, CYP2C8, SLC47A1, BRAF, KDR, KIT, TEK 1. Abrogates KIT/PDGFRa-mediated tumor cell signaling and prevents proliferation in KIT/PDGFRa-driven cancers. 2. Inhibits several other kinases: VEGFR2, KDR, TIE2, TEK, PDGFR-beta, FMS, CSF1R, inhibite tumor cell growth. National Library of Medicine Mometasone furoate Nasal allergy –4.641 0.01 CYP2C8, NR3C1, PGR NA NA Ractopamine Animal feed additive, clenbuterol –4.949 0.01 ESR1 NA NA Turpentine Paint and varnish solvents, infections of the stomach and intestines, autism, and pain –5.293 0.01 MAPK1 NA NA Estriol Multiple sclerosis –4.642 0.02 ABCB1, SLCO1A2, SHBG, ESR1, ESR2 Inhibit the production of t-PA; The expression of tyrosinase in B16 melanoma cells induced by imidazole was blocked NA Panobinostat Multiple myeloma –4.989 0.02 ABCB1 1. A histone deacetylase inhibitor; 2. Activate MEK1/2 and ERK1/2 , key melanoma-genic kinases blocked; 3. Interfere with the appropriate folding of HSP90-client proteins (including AKT and RAF) that are critical to cancer cell growth; 4. Increased apoptosis as well as a G1 cell cycle arrest; 5. Increased ability to activate antigen-specific T cells; 6. Restore BRAF-inhibitor sensitivity in melanoma. (Xia et al. 2014 ; Gallagher et al. 2018 ) Fidaxomicin Antibiotic –4.989 0.02 ABCB1, sigA1 Inhibits RNA polymerase sigma subunit resulting in inhibition of protein synthesis and cell death in susceptible organisms NA Ponatinib Leukemia –4.131 0.02 ABCB1, FGFR2, FGFR1, FGFR4, FGFR3, ABCG2, PDGFRA, CYP2C8, ABL1, KDR, FLT3, KIT, RET, TEK, BCR, SRC, LCK, LYN 1. Adhesive to KITD816V mutant protein and inhibit tumor growth 2. FGFR inhibitors (Han et al. 2019 ) Technetium Tc-99m tetrofosmin Diagnostic Agent, radiopharmaceutical imaging –4.989 0.02 ABCB1 NA NA Duvelisib Chronic lymphocytic leukemia, follicular lymphoma or small lymphocytic leukemia –4.506 0.02 ABCB1, ABCG2, PIK3CD, PIK3CG NA NA Grapiprant Immunomodulating, antineoplastic activities –4.989 0.02 ABCB1, PTGER4 NA NA Drospirenone Prevent pregnancy –4.344 0.02 PTGS2, NR3C1, AR, PGR, NR3C2, PON1 NA NA Segesterone acetate Prevent ovulation and pregnancy –4.581 0.02 NR3C1, AR, PGR NA NA Homosalate Sunscreen, skin care products with sun protection –4.567 0.02 ESR1, AR, PGR Src activation in melanoma and Src inhibitors (Homsi et al. 2009 ) Enzacamene Sunscreen –4.423 0.02 ESR1, ESR2, AR, PGR NA NA Bosutinib Leukemia –4.136 0.03 ABCB1, CYP2C8, ABL1, BCR, SRC, LYN, HCK, CDK2, MAP2K1, MAP2K2, MAP3K2, CAMK2G An SRC inhibitor, induces caspase-independent cell death associated with permeabilization of lysosomal membranes (Noguchi et al. 2018 ) Ledipasvir Chronic (long-lasting) hepatitis C, a viral infection of the liver –4.989 0.03 ABCB1, ABCG2, NS5A NA NA Omadacycline Community-Acquired Bacterial Pneumonia (CABP), Acute Bacterial Skin and Skin Structure Infections (ABSSSI) –4.989 0.03 ABCB1, rpsN, rpsS, rpsG, rpsC, rpsH NA NA Berberine Diabetes, obesity, and inflammation –4.989 0.03 BIRC5, qacR 1. via Affecting the FAK, uPA, and NF-κB Signaling Pathways and Inhibits PLX4032 Resistant A375.S2 Cell Migration; 2. Inhibits melanoma cancer cell migration by reducing the expressions of cyclooxygenase-2, prostaglandin E2and prostaglandin E2 Receptors suppressed epithelial mesenchymal transition through cross-talk regulation of PI3K/AKT and RARα/RARβ in melanoma cells (Kou et al. 2016 ; Liu et al. 2018a ; Ren et al. 2020 ) Oxybenzone Sunscreen –4.423 0.03 ESR1, ESR2, AR, PGR Inhibited cell growth, DNA synthesis and retarded cycle progression from G1 (Xu and Parsons 1999 ) Prucalopride Chronic idiopathic constipation –4.989 0.04 ABCB1, HTR4 NA NA Boceprevir Hepatitis C –4.989 0.04 ABCB1, NS3/4A NA NA Favipiravir Influenza A, B, and C virus infections. –4.989 0.04 ABCB1, SLC22A6, SLC22A8, AOX1, CYP2C8, SLC22A12, PB1 NA NA Pemigatinib Bile duct cancer –4.166 0.04 ABCB1, FGFR2, FGFR1, FGFR4, FGFR3, ABCG2, SLC22A2, SLC47A1 Acts on FGFR2 gene to prevent tumor growth and proliferation NA Ingenol mebutate Actinic keratosis. –4.61 0.04 PRKCA, PRKCD Upregulation of Erdr1 and a significant decrease in p53 and bcl-2 activates a Broad range of protein kinase C (PKC) (a, b, c, d, e, g and h) isoenzymes or Direct proapoptotic effects (Mansuy et al. 2014 ; Woo et al. 2017 ) Cabozantinib Renal cell carcinoma (RCC) Medullary Thyroid Cancer (MTC) –4.321 0.04 CYP2C8, KDR, RET, MET Inhibit both VEGFR2 and c-Met simultaneously (Roy et al. 2015 ) Fluoroestradiol F-18 An imaging agent –4.469 0.04 SHBG, ESR1, ESR2 NA NA Diflorasone Skin irritation or skin rashes –4.989 0.04 NR3C1 NA NA Arzoxifene –4.469 0.04 ESR1, ESR2 NA NA Molecular and supporting information of anti-melanoma drug candidates based on network group analysis 1. Epithelioid sarcoma 2. Follicular lymphoma 1. Strong selective inhibition of histone methyltransferase EZH2 2. Inhibit some EZH2 gain-of-function mutations (including Y646X and A687V), as well as EZH1 1. Histone Methyltransferase (HMT) Inhibitor 2. Histone deacetylase (HDAC) inhibitor 1. Induces growth inhibition and apoptosis; 2. Downregulates cancer procoagulant activity in MCF-7 and WM-115 cell lines 1. Modify selenoprotein transcript levels 2. BRAF+MEK inhibitor 3. The histone deacetylase inhibitor (Ecker et al. 2020 ) (Wang et al. 2018 ) (Haas et al. 2014 ) 1. Release of NO, mediates apoptosis of tumor cells by generating physiological messenger circular GMP; 2. Increase the efficacy of chemoimmunotherapies in vivo, promoting the survival and increasing the function of injected cells by targeting a key pathway incisplatin-induced cytotoxicity. 1. Inhibition of tyrosine kinase receptors; 2. Downregulated pRPS6KB1 and pRPS6, downstream of mTOR signaling, and IRS-1 in the insulin signaling pathway, specifically ablating IRS-1 nuclear signal. 1. By activating endogenous K+ and exogenous H+ exchange, it causes oxidative stress, induces mitochondrial permeability transformation (MPT), and produces cytotoxic effects on cancer cells 2. Induce autophagy in melanoma cell lines (M21, M29, SK-MEL-1, SK-MEL-19, SK-MEL-103 and A375) that are sensitive to SAL. (Ohkubo et al. 2018 ) (Liu et al. 2020 ) 1. Abrogates KIT/PDGFRa-mediated tumor cell signaling and prevents proliferation in KIT/PDGFRa-driven cancers. 2. Inhibits several other kinases: VEGFR2, KDR, TIE2, TEK, PDGFR-beta, FMS, CSF1R, inhibite tumor cell growth. Inhibit the production of t-PA; The expression of tyrosinase in B16 melanoma cells induced by imidazole was blocked 1. A histone deacetylase inhibitor; 2. Activate MEK1/2 and ERK1/2 , key melanoma-genic kinases blocked; 3. Interfere with the appropriate folding of HSP90-client proteins (including AKT and RAF) that are critical to cancer cell growth; 4. Increased apoptosis as well as a G1 cell cycle arrest; 5. Increased ability to activate antigen-specific T cells; 6. Restore BRAF-inhibitor sensitivity in melanoma. 1. Adhesive to KITD816V mutant protein and inhibit tumor growth 2. FGFR inhibitors 1. via Affecting the FAK, uPA, and NF-κB Signaling Pathways and Inhibits PLX4032 Resistant A375.S2 Cell Migration; 2. Inhibits melanoma cancer cell migration by reducing the expressions of cyclooxygenase-2, prostaglandin E2and prostaglandin E2 Receptors suppressed epithelial mesenchymal transition through cross-talk regulation of PI3K/AKT and RARα/RARβ in melanoma cells Renal cell carcinoma (RCC) Medullary Thyroid Cancer (MTC) The following shortcomings remain in this study: First, identifying risk genes and drug targets for SKCM is a long-term challenge in the development of computational drug repositioning. Our research protocol relies on known risk genes associated with SKCM. There are many unreported SKCM-related risk genes and drug targets that have not been studied in depth. Second, the human PPI network is extremely complex and huge. Due to the limitations of experimental techniques and algorithms, the currently available information may be biased, covering only a small part of the proteome, which needs further experimental verification. Finally, the drugs collected in the Drugbank database are still increasing, and there may be many drugs that are effective in the treatment of SKCM that have not yet been discovered. Furthermore, the association of ZC3H12A expression with SKCM tumor stage should be validated in large prospective clinical trials before translation into clinical practice. Meanwhile, functional studies are needed to fully understand and exploit the role of ZC3H12A in the activation of early SKCM immune responses. In short, through the combination of molecular network methods and drug epidemic diseases, we have determined several candidates that may be used to treat SKCM. Although the result is only initial, it can provide some clues for further SKCM drug development. At the same time, it was found that ZC3H12A was rising in phase I SKCM and was related to the improvement of prognosis. It may be due to its participation in the immune process to reshape the tumor immune activation micro-environment. These findings may have potential clinical significance.

Skin cutaneous melanoma (SKCM) is a common skin tumor with extremely high invasiveness and mortality(Qing Lv et al. 2017 ). In the past decade, targeted therapy and immunotherapy have made great progress in the treatment of SKCM, and the survival and quality of life of patients with SKCM have improved significantly (Baroudjian et al. 2016 ; Singh et al. 2011 ). However, there are still some defects such as secondary drug resistance and serious adverse reactions. In addition, this treatment approach has limitations in clinical practice, such as the lack of effective therapeutic targets in some patients, which leads to the failure of the corresponding targeted therapy, or the ineffective treatment of immune checkpoint inhibitors due to the lack of sufficient immune infiltrating cells in tumor tissues. These limitations lead to the failure of these two highly specific treatments to benefit patients with SKCM. Given the limitations of SKCM treatment, the development of new drugs and the identification of novel biomarkers for risk stratification are needed to improve clinical outcomes in SKCM. The tumor microenvironment (TME) refers to the complex tumor cells including peripheral blood vessels, immune cells, fibroblasts, bone marrow-derived inflammatory cells, various signaling molecules, and extracellular matrix (ECM). The dynamic ecosystem on which survival and development depend (Binnewies et al. 2018 ; Hanahan and Weinberg 2011 ). The immune cells and their regulation in this microenvironment play a crucial role in the occurrence and development of tumors (Pandey et al. 2022 ). At present, many studies have confirmed that the deterioration of malignant melanoma is closely related to the immune microenvironment, but no more studies have confirmed which immune cells infiltrate the metastasis of malignant melanoma (Kalaora et al. 2022 ). As for SKCM pathogenesis-related differential genes, microenvironment interactions are even less reported. The research and development of new SKCM drugs was a time-consuming process with a large investment and a low success rate. Therefore, drug repositioning has attracted extensive attention of medical researchers. Drug repositioning, also known as “reuse of old drug” and “re-examination of old drug,” refers to the discovery of new indications or new uses of listed drugs, including the re-use, re-evaluation, and repositioning of the treatment direction (Pushpakom et al. 2019 ). Because the drugs for repositioning studies have usually passed several stages of clinical trials or have been marketed, and the uncertainties in safety and pharmacokinetics are significantly reduced, the development cycle is shortened, and the development costs and risks are obviously reduced. In addition, drug repositioning was currently known because of the advantages of shorter time-to-market and greater likelihood of discovering differences in drug action compared to the strategy of reformulation after patent licensing. One of the strategies is the (Sonawane et al. 2019 ) development strategy (Pushpakom et al. 2019 ; Sonawane et al. 2019 ). With the development of network biology research, network pharmacology technology has provided a new means for new drug research and development, and has been continuously applied in drug repositioning research, becoming one of the important technologies in tumor drug repositioning research (Fiscon et al. 2022 ). Given the available findings, we know that SKCM pathogenesis-related genes play an important role in the development of SKCM, however, its impact on the TME has not been investigated. Therefore, we conducted a systematic study to identify differential genes associated with SKCM risk and to explore the prognostic value of these genes and their association with the tumor immune microenvironment; subsequently, network analysis was performed on genes and drug targets related to the pathogenesis of skin melanoma to screen potential drugs. The candidate drugs were further screened according to the gene expression profiles of the drugs, and finally, some drugs that may be effective in the treatment of skin melanoma were obtained.

ESM 1 (DOCX 10 kb) ESM 2 (PNG 386 kb) High Resolution Image (TIF 4119 kb) ESM 3 (PNG 771 kb) High Resolution Image (TIF 5237 kb) ESM 3 (DOCX 19 kb) ESM 4 (TXT 94 kb) (DOCX 10 kb) (PNG 386 kb) High Resolution Image (TIF 4119 kb) (PNG 771 kb) High Resolution Image (TIF 5237 kb) (DOCX 19 kb) (TXT 94 kb)