Transcriptomic Alterations Induced By Vemurafenib after Treatment of Melanoma: A Comprehensive Bioinformatics Analysis

*Corresponding Author:

J. Shi
The Forth School of Nanjing Medical University,
Nanjing,
Jiangsu 210029,
China
E-mail: jpshi2021@163.com

Abstract

Melanoma is a highly aggressive kind of cancer with a very poor prognosis. v-Raf murine sarcoma viral oncogene homolog B1 inhibitor vemurafenib has indeed harvested substantial clinical benefits. Nevertheless, its drug resistance has also hampered scientists’ efforts towards successful melanoma treatment. In this study, we used data derived from the gene expression omnibus database to analyze the effect on vemurafenib sensitive cell lines after vemurafenib treatment. Gene expression omnibus datasets GSE42872 (cohort 1), GSE127988 (cohort 2), GSE110054 (cohort 3) were included in the analysis. We found 25 common differentially expressed genes in 3 datasets, including 10 upregulated genes and 15 downregulated genes after vemurafenib application. Analysis using web tool Timer showed a significant correlation of the upregulated genes with immune infiltration level in skin cell melanoma. Gene ontology enrichment analysis showed that after vemurafenib treatment, all datasets showed downregulation in deoxyribonucleic acid replication and cell cycle arrest. Meanwhile, genes related to neurogeneration, extracellular matrix and cell-cell adhesion were significantly enriched in all three datasets. Kyoto Encyclopedia of Genes and Genomes analysis showed that pathways like P53, phosphatidylinositol 3-kinase-protein kinase B and Ras-associated protein signaling pathways were enriched in differentially expressed genes after vemurafenib administration. The findings of the candidate differentially expressed genes and pathways may not only reveal the cellular sensitivity to vemurafenib treatment but also give rise to a better understanding of the mechanism of cancer cell cycle arrest and cellular resistance towards vemurafenib targeted therapy.

Keywords

Melanoma, vemurafenib, v-raf murine sarcoma viral oncogene homolog B1, drug resistance, bioinformatics, differentially expressed genes

Melanoma, a highly malignant tumor, has puzzled many scholars with its treatment [1]. Approximately 1.7 % of newly diagnosed primary malignancies worldwide are cases of cutaneous malignant melanoma which causes about 55 500 deaths annually [2]. Advanced melanoma is often associated with a life-threatening tumor [3]. Before the advent of immune checkpoint inhibitors and gene-targeted drugs, there were few effective plans to treat advanced melanoma and clinical trials did not yield satisfactory results [4]. Programmed Cell Death 1/Programmed Death Ligand 1 (PD1/PDL1) immune checkpoint inhibitors have achieved good results in the treatment of melanoma, but there is no decrease in enthusiasm for the research of gene-targeted therapy to improve the clinical effect [5]. v-Raf murine sarcoma viral oncogene homolog B1 (B-RAF) mutation, a genetic mutation found in 50 % to 60 % of melanomas, often occurs in the 600th codon, replacing valine with glutamic acid (so called B-RAF^V600E). This mutation leads to the activation of the BRAF/ERK kinase (MEK)/Extracellular-Signal-Regulated Kinase (ERK) pathway, causing continuous cell growth [6].

Vemurafenib, a novel drug first synthesized in 2005, inhibits BRAF primarily. This drug was first approved by the Food and Drug Administration (FDA) in 2011 to treat melanoma with BRAF mutations [7]. However, patients who receive targeted therapies initially respond well to treatment, but there is resistance in the long term [8]. At present, researchers speculate that drug resistance may be related to secondary mutations, bypass signaling pathways, tumor microenvironment changes and immune infiltration [9]. To provide supporting material for current drug resistance studies, we conducted a detailed bioinformatics analysis to explore what genetic changes were caused by vemurafenib in BRAF mutant melanoma cells.

Materials and Methods

Oncogenic properties and clinical manifestations of BRAF:

University of Alabama Cancer Database (UALCAN) is a user-friendly interactive web portal to perform indepth analyses of The Cancer Genome Atlas (TCGA) gene expression data [10]. We used the UALCAN databases to verify the differential Messenger Ribonucleic Acid (mRNA) expression of BRAF in humans. We investigated relative BRAF mRNA expression in variant types of tumors and different stages Skin Cutaneous Melanoma (SKCM) patients using UACLAN databases.

TCGA datasets:

TCGA, which is a joint effort between National Cancer Institute (NCI) and the National Human Genome Research Institute beginning in 2006, molecularly characterized over 20 000 primary cancer and matched normal samples spanning 33 cancer types, bringing together researchers from diverse disciplines and multiple institutions. We investigated the TCGA datasets in order to figure out the genetic mutations in BRAF genes within skin cutaneous melanoma patients and to analyze the survival curve between the mutated and not mutated cases.

Gene Expression Omnibus (GEO) datasets:

The GEO is a public repository at the National Center of Biotechnology Information (NCBI) for storing high-throughput gene expression datasets. We selected potential GEO datasets using search terms “melanoma" [MeSH Terms] and vemurafenib [All Fields] and "Homo sapiens" [Organism] “Expression profiling by high throughput sequencing” [Dataset Type] in the GEO datasets to identify potential datasets. Then, we further screened these datasets according to the inclusion criteria above, deleting invalid outcomes and samples that are too few for statistical analysis (with sample volume less than 3), as well as single-cell sequencing analysis. Finally, datasets GSE42872 (cohort 1), GSE127988 (cohort 2), GSE110054 (cohort 3) were included in our study, all of which were experiments conducted using melanoma cell lines M297.

R software and data processing:

R software (version 4.0.3) was used to download and process the raw data. After the expression data has been processed, we used log2 fold change to quantify the changes in expression between the control group and the vemurafenib group. R packages from Bioconductor were used in the analysis of the data; the volcano plots were drawn using the package “ggscatter”; the heatmaps were drawn using the packages “pheatmap”. p<0.01 and logFC>2 were considered statistically significant. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis were conducted using the package “ClusterProfiler” to investigate whether those significantly upregulated or downregulated genes were enriched in certain pathways.

Timer:

Timer is an open web tool that can be used to explore the association of genes in immune cell infiltration in different tumors. With Timer, B cells, CD4⁺ T cells, CD8⁺ T cells, macrophages, neutrophils and dendritic cells can all be studied. In this paper, we used Timer to explore the role of 10 highly expressed genes in immune infiltration in melanoma.

String database and Cytoscape software:

The Search Tool for the Retrieval of Interacting Genes/ Proteins (STRING) is a user-friendly, interactive database to analyze protein-protein interaction. In this study, we used the string database to investigate common differential genes in order to understand the protein interaction after using the Cytoscape, an opensource software platform for visualizing complex networks was used in our study to further analyze the genes that are analyzed using the STRING database.

The genes were ranked according to their degrees (meaning the number of a protein’s interaction with other proteins) and those genes with the highest degree scores (with scores ranging from 1 to 15) were recognized as hub genes.

Survival analysis:

Gene Expression Profiling Interactive Analysis (GEPIA) is an interactive web tool for analyzing gene expression within different types of cancers and normal samples, based on the data acquired from TCGA database [11].

We conducted the survival analysis and drew the survival curve of between the expression of the hub genes from the Protein-Protein Interaction (PPI) network and the patient’s survival outcome using the GEPIA web server. Therefore, the genes analyzed from GEO datasets are ultimately validated by TCGA data.

Results and Discussion

The process of our study was summarized in fig. 1. We enquired into the BRAF mutation and expression level within the TCGA database using TCGA, cBioportal and GEPIA. Firstly, a pan-cancer analysis showed that BRAF mRNA level in multiple melanomas was significantly higher than that in various other tumor tissues, demonstrating that BRAF may be a characteristic marker within multiple kinds of melanomas (fig. 2A). Secondly, using the TCGA database, we found that BRAF ranks 3^rd among the most frequently mutated genes in SKCM, with almost half of the total cases involved (fig. 2B). Further analysis using GEPIA exhibited that BRAF expression level was significantly up regulated as melanoma stage upgrades (fig. 2C). Finally using TCGA datasets, we plotted the survival curve of the patients with and without BRAF mutation (fig. 2D). The result was quite intriguing: the mutated group had a significantly better outcome than the not mutated group, implying that treatment targeting BRAF mutation had a significant effect, which is also the key point that we are going to address later on.

Fig. 1: Flow chart of our study

Fig. 2: Analysis of genetic features of BRAF within melanoma, (A) pan-cancer analysis of relative genetic expression within variant types of human cancers; (B) Distribution of most frequently mutated genes within SKCM; (C) Expression of BRAF in different stages of SKCM; (D) Survival analysis of different mutation situations (BRAF-V600Emutated and (BRAF-V600E not mutated)

To understand the changes in RNA expression after the application of vemurafenib cell lines, we used the GEO datasets GSE42872 (cohort 1), GSE127988 (cohort 2), GSE110054 (cohort 3) to conduct the differential analysis using R. Among them, GSE127988 contains expression profile data of 3 different vemurafenib concentrations: 0.1 M, 0.316 M and 1 M. Therefore, 5 parallel experiments were included in the study. Using the Venn diagram, we identified 25 common differential genes that co-existed in all those 5 parallel experiments (we deleted one invalid output from the 26 genes that served as an outcome (fig. 3A). Fig. 4A is the volcano map of gene expression after vemurafenib treatment in Cohort 1. Fig. 4C is a volcanic map of three different drug concentrations in Cohort 2 and fig. 4E is a volcanic map of Cohort 3. Next to the volcano map, we also made a heat map of gene expression (fig. 4B, fig. 4D and fig. 4F). These genes are potentially helpful in identifying a melanoma patient’s response to vemurafenib treatment. The relationship between 10 highly expressed genes and immune infiltration was studied in Timer, as shown in fig. 5.

Fig. 3: The venn gram

Fig. 4: Expression profile of datasets, (A) Volcano plot of gene expression in cohort 1, the significantly upregulated and downregulated genes are labeled in red and blue; (B) the genes with the greatest fold changes are and heatmap of 25 common differentially expressed genes (DEGs) in cohort 1; (C) Volcano plot of gene expression in cohort 2 with different vemurafenib concentration of 0.1 M, 0.316 M and 1 M; (D) heatmap of 25 common DEGs in cohort 2 with different vemurafenib concentration of 0.1 M, 0.316 M and 1 M; (E) Volcano plot of 25 common DEGs in cohort 3; (F) heatmap of 25 common DEGs in cohort 3

Fig. 5: The relationship between 10 highly expressed genes and immune infiltration

We then used the R package ‘ClusterProfiler’ to further investigate where those differential genes are significantly enriched. Since Cohort 2 contains three different concentrations, GO and KEGG enrichment analyses were conducted in five parallel experiments. Not surprisingly, GO enrichment of the downregulated genes (fig. 6B, fig. 6D, fig. 6F, fig. 6H, fig.7B) showed that the down-regulated genes were predominantly enriched within the biological processes such as “Deoxyribonucleic Acid (DNA) replication”, “Cell Cycle”, “DNA packaging”, implying that after vemurafenib administration, the cellular replication activity was significantly reduced. Moreover, cellular component enrichment analysis (fig. 6B, fig. 6D, fig. 6F, fig. 6H, fig. 7B) showed that those downregulated genes were most frequently distributed within the chromosomal regions, as well as those regions that are extremely vital to cellular replication and metabolism, such as replication fork, spindle and kinetochore. Besides, molecular function GO analysis of down-regulated genes showed that enzymes that are essential in DNA replication and formation were also significantly enriched (fig. 6B, fig. 6D, fig. 6F, fig. 6H, fig. 7B). Finally, KEGG analysis of down-regulated genes was enriched in pathways such as cell cycle, DNA replication (fig. 7D, fig. 7F, fig. 7H, fig. 8B, fig. 8D).

Fig. 6: (A) Results of GO analysis of upregulated DEGs in cohort 1; (B) Results of GO analysis of downregulated DEGs in cohort 1; (C, E, G) results of GO analysis of upregulated DEGs in cohort 2 with different vemurafenib concentration of 0.1 M, 0.316 M and 1 M; (D, F, H) results of GO analysis of downregulated DEGs in cohort 2 with different vemurafenib concentration of 0.1 M, 0.316 M and 1 M

Fig. 7: (A) Results of GO analysis of upregulated DEGs in cohort 3; (B) Results of GO analysis of downregulated DEGs in cohort 3; (C) Results of KEGG analysis of upregulated DEGs in cohort 1; (D) Results of KEGG analysis of downregulated DEGs in cohort 1; (E,G) Results of KEGG analysis of upregulated DEGs in cohort 2 with different vemurafenib concentration of 0.1 M, 0.316 M; (F, H) Results of KEGG analysis of downregulated DEGs in cohort 2 with different vemurafenib concentration of 0.1 M, 0.316 M

Fig. 8: (A) Results of KEGG analysis of upregulated DEGs in cohort 2 with vemurafenib concentration of 1 M; (B) Results of KEGG analysis of downregulated DEGs in cohort 2 with vemurafenib concentration of 1 M; (C) Results of KEGG analysis of upregulated DEGs in cohort 3; (D) Results of KEGG analysis of downregulated DEGs in cohort 3

Interestingly, GO enrichment analysis of up-regulated genes in all five GEO experiments showed that mRNAs concerned with extracellular matrix or cell-cell adhesion were significantly enriched (fig. 6A, fig. 6C, fig. 6E, fig. 6G and fig. 7A). This finding is important in that it implies vemurafenib may enhance cellular adhesion via acting on extracellular matrix, thus reducing migration in melanoma cells, which also requires further investigation. Besides, we also found that molecular functions related to neurogenesis were significantly enriched in cohort 1 and cohort 2. KEGG analysis of upregulated genes (fig. 7C, fig. 7E, fig. 7G, fig. 8A and fig. 8C) showed that multiple cellular signaling pathways were enriched after using the BRAF inhibitor, such as the P53 signaling pathway was deactivated after using vemurafenib, while phosphatidylinositol 3-kinaseprotein kinase B (PI3K-AKT) signaling pathways, Role of Ras-Associated Protein (Rap) signaling pathways, as well as pathways related to focal adhesion, were somehow activated after vemurafenib application.

We first drew a Venn diagram of the differentially expressed genes within the GEO datasets GSE42872, GSE 127988 (which contains 3 different vemurafenib concentrations 0.1 M, 0.316 M and 1 M) and GSE234231 (fig. 3A). Not surprisingly, the differentially expressed genes have multiple genes in common. However, due to a lack of sufficient samples, the datasets of GSE110054 only had 98 differential genes that met our previous criteria. Therefore, we selected the genes that are included in at least 4 parallel experiments, which are illustrated as samples encircled in red lines (fig. 3B).

By using the STRING database and the Cytoscape software, we conducted a protein-protein interaction analysis of those common genes (fig. 9A). The result of our PPI analysis showed multiple hub genes that are at the core of the PPI network, which was further confirmed by survival analysis. Using the GEPIA web tool, we examined the hub genes (meaning genes that have the maximal degree of freedom-interaction with other proteins) in our PPI network. We eventually found that the overexpression of most of the hub genes was related to poor survival outcomes. And most of those genes were genes that are essential in cell cycle and replication, for instance, gene Cyclin Dependent Kinase 1 (CDK1) (fig. 9B). Most of the hub genes were significantly downregulated (Table 1), showing the drug’s core effect relies on melanoma cell cycle arrest.

Fig. 9: (A) Protein-Protein interaction of 605 Genes that tare enriched within at least 4 experiments, the genes that have the biggest degree are situated in the center of the diagram; (B) The Kaplan-Meier survival plot of the 15 hub genes that have the highest degree of interaction

Table 1: Significant Values of Genes

Targeted therapies for melanoma are evolving rapidly with the advance of basic research [12]. In the past, advanced melanoma therapy has baffled many researchers, with an estimated 5 y survival rate of only 16 % for patients with advanced malignant melanoma [13]. With the development of precision medicine, B-RAF mutations and immune checkpoint changes have been found in the pathogenesis of melanoma, which is an exciting discovery and has greatly promoted the treatment of melanoma [14]. This was followed by the emergence of B-RAF inhibitors (vemurafenib, dabrafenib) and immune checkpoint inhibitors (Ipilimumab, Pembrolizumab) that gave patients with advanced melanoma a ray of life [15]. A growing number of clinical trials have shown that these two types of drugs can significantly improve survival in patients with advanced melanoma [16].

In 2002, researchers found that melanoma cells often harbor B-RAF mutations [14]. This opens the door to targeted therapies. The activation process of this pathway is: growth factor binding Receptor Tyrosine Kinases (RTKs) leads to automatic phosphorylation of the receptor, thus phosphorylation of RAS-Guanosine Diphosphate (GDP) into RAS-Guanosine Triphosphate (GTP) (active state). Activated RAS binds to RAF to activate downstream sites such as MEK/ERK [17-19]. This cascade is highly conserved, with many protein-protein interaction networks that mediate cell proliferation, differentiation and apoptosis [20]. Normally, B-RAF proteins are highly regulated by the body and are often bound to Heat Shock Proteins (HSP) in a low active state [21]. Once the B-RAF mutation occurs, it leads to the continuous activation of functional domains, causing a cascade of downstream sites, thereby promoting cell growth and inhibiting apoptosis [22,23].

Vemurafenib, which was primarily synthesized in 2005, targets B-RAF effectively [24,25]. It was found to be selective for melanoma cells with B-RAF mutation compared to non B-RAF mutant melanoma cells and has the good clinical effect and strong compliance [24]. Vemurafenib had a combined response rate of about 50 % and significantly improved progression-free survival compared with the previous chemotherapy regimen-dacarbazine [26]. Because B-RAF mutations occur in 50 % to 60 % of melanomas, the drug could be a promising treatment for B-RAF mutated melanomas[25]. In addition, vemurafenib has shown promising results in the treatment of lung, colorectal, breast and ovarian cancer with B-RAF mutations[27,28]. It mainly acts on the "The Asp-Phe-Gly (DFG) motif", which is continuously activated in B-RAF mutations, to achieve inhibitory pathway activity [29,30]. Because of these benefits, it was approved by the FDA in 2011 as the first drug to treat cancers harboring B-RAF mutations [31,32]. Despite all its virtues, it did not escape the fate of drug resistance [33]. Like other kinase inhibitors, vemurafenib shows striking clinical effects at first, but resistance develops after 6-8 mo [34]. The resurrection of mutations, contradictory activation of drugs, tumor microenvironment and immune infiltration are speculated to be the reasons for the development of drug resistance and tumor progression [35,36]. Studies have found that the Mitogen-Activated Protein Kinases (MAPK) pathway is highly expressed in drug-resistant tumors, so scholars speculated that B-RAF inhibitor and MEK inhibitor could be used in combination to reduce drug resistance [37]. In order to prove the view points above, many clinical studies have been carried out, which indeed received good clinical effects, showing its potential practicability [38].

In this study, using the datasets derived from the vemurafenib cell lines, we investigated the genes which are significantly upregulated and downregulated and identified 25 genes that have significant changes in expression in all 3 variant GEO datasets: Olfactomedin like 2A (OLFML2A), Aurora kinase B (AURKB), RAS Related (RRAS), ADAM Metallopeptidase With Thrombospondin Type 1 Motif 15 (ADAMTS15), CDK2, Polymerase Epsilon (POLE), Musculoaponeurotic Fibrosarcoma (MAF), Sprouty RTK Signaling Antagonist 4 (SPRY4), Matrix Metallopeptidase 19 (MMP19), SRY-Box Transcription Factor 4 (SOX4), Replication Factor C Subunit 3 (RFC3), Forkhead Box D3 (FOXD3), Methylenetetrahydrofolate Dehydrogenase (MTHFD2), Bloom syndrome helicase (BLM), ETS translocation variant 4 (ETV4), Cell Division Cycle Associated 7 (CDCA7), DNA Primase Subunit 1 (PRIM1), ATPase Family AAA Domain Containing 2 (ATAD2), Tumor Protein P53 Inducible Nuclear Protein 1 (TP53INP1), SH2B adapter protein 3 (SH2B3), Zinc Finger E-Box Binding Homeobox 2 (ZEB2), Integrin Subunit Beta 3 (ITGB3), Jun Proto-Oncogene (JUN), Non-SMC Condensin I Complex Subunit G (NCAPG). These genes are potentially useful in identifying a patient’s responsiveness to vemurafenib treatment, although more studies are needed to verify this issue. Moreover, Enrichment results showed a significant upregulation in genes that are concerned with extracellular matrix and cellular adhesion. Whether vemurafenib is directly correlated with reduced cancer cell migration and invasion requires more validation. Despite all these findings, we have seen enrichment in focal adhesion, which has been proved to participate in drug resistance by activating B-Cell Lymphoma 2 (BCL-2) through ERK signaling pathways. Moreover, we also found multiple tumorigenic genes downregulated that can significantly reduce patient survival, including CDK1, CDK2, CDC45, Protein Regulator Of Cytokinesis 1 (PRC1), Polo Like Kinase 1 (PLK1), Cyclin A2 (CCNA2)—most of them being key players in cell cycles, were at the core in the PPI network after vemurafenib application.

The function of downregulated genes after vemurafenib treatment was enriched in cell replication (particularly in the cell cycle), which was closely related to the expression of cyclin and CDK, revealing its significant role in inhibiting the growth of tumor cells. Tumorigenesis is the process of loss of normal regulation resulting in continuous growth. Cell cycle refers to the time that a cell goes through mitosis, including the Gap 1 (G1) phase, Synthesis (S) phase, Gap 2 (G2) phase and Mitosis (M) phase [39]. The G1, S and G2 phases are collectively referred to as the intermitotic phase, which is the preparation period for the M phase [40]. G1/S checkpoint is an important site of cell regulation, which is completed by dephosphorylation and phosphorylation of cyclin, CDK and casein kinase I (CKI) [41]. Cyclin and CDK positively promote the cell cycle, while CKI is involved in negative regulation. High expression of cyclin and CDK was detected in many tumors [42]. Cyclin D1 may be related to the activation of the MAPK signaling pathway and vemurafenib can inhibit this pathway. In our analysis, we also found that vemurafenib can reduce the expression of cyclin and CDK, suggesting the rationality of vemurafenib treatment. In conclusion, vemurafenib inhibits genes highly related to the cyclin and CDK expression, thereby reducing tumor cell replication and slowing tumor growth.

The tumor microenvironment is one of the hotspots of current research, among which the Extracellular Matrix (ECM) is widely mentioned in the resistance of melanoma.

Fibronectin, a class of extracellular matrix, is often found to be expressed in drug-resistant cells with B-RAF mutations. Its mechanism of resistance is the activation of the alpha-5 beta-1 (α5β1)/Akt pathway [43]. In addition, Hirata et al. found that the abnormal activation of vemurafenib on melanoma-associated fibroblasts is associated with the stromal generation and induction of Integrin β1 (ITGB1)/Focal adhesion kinase (FAK)/steroid receptor coactivator (Src) signaling in melanoma cells, which provides a mechanism of drug resistance [44]. Our analysis found that the up-regulated genes after vemurafenib administration are associated with ECM function, which may be a potential mechanism for drug resistance after targeted therapy.

There are many adhesion molecules on the surface of melanoma cells, such as cadherin, integrin and immunoglobulin superfamily [45]. The cell adhesion molecule melanoma cell adhesion molecule (MCAM)/ molecule MCAM (MUC18) enables melanoma cells to have metastatic potential and enhanced tumorigenicity. MCAM/MUC18 mediated homomorphic and heteromorphic adhesion between melanoma cells and endothelial cells, respectively [46]. These two interactions may promote transfer at different stages of the cascade [47]. The loss of E-cadherin and the increase of N-cadherin in melanoma cells often disrupt the normal homeostasis, promote the occurrence and development of tumor and may mediate the occurrence of Epithelial-Mesenchymal Transformation (EMT) [48]. Therefore, some scholars speculate that the increase of intercellular adhesion may be the mechanism of drug resistance in targeted therapy. In our analysis, there was also an increase in cellular adhesion after vemurafenib treatment, which may be related to the development of drug resistance.

Neurogenesis often occurs after repairing nerve damage, such as a stroke. It is often accompanied by changes in the local microenvironment, causing the regeneration of blood vessels and nerves. Such local microenvironmental alterations are also present in melanoma cells. Prakash et al. found that neurogenesis could increase the brain metastasis of melanoma [49]. According to their study, neurogenesis can promote the action of melanoma on components such as ECM, potentially leading to the activation of some signaling pathways and increasing the metastatic potential of melanoma. In our study, vemurafenib led to the overexpression of genes whose functions were enriched in neurogenesis. This may also be related to the development of drug resistance after vemurafenib treatment, resulting in adverse outcomes.

Consistent with current clinical observations, vemurafenib has a dual effect in the treatment of melanoma. On the one hand, it targets B-RAF mutations, causing downstream cascades (such as the MAPK pathway) to be suppressed, resulting in low expression of the associated proto-oncogenes. The results of our analysis showed that the suppressed genes were involved in promoting cell replication, suggesting that vemurafenib inhibited tumor cell proliferation. However, on the other hand, vemurafenib can also cause some adverse events, such as the promotion of cell adhesion, changes in ECM and neurogenesis, which may be the result of the contradictory activation of vemurafenib and may be related to the acquisition of drug resistance, but it still needs to be confirmed by further studies.

Currently, the combination of B-RAF inhibitor and MEK inhibitor has been found to reduce contradictory activation, thus avoiding resistance to vemurafenib, improving clinical outcomes and promoting patient survival. Our findings can provide a direction for future exploration of the mechanism of vemurafenib drug resistance. But there are limitations to our study. First, we lack in vivo and in vitro trials to verify the results. And then, there was heterogeneity in the different experimental data we included.

Vemurafenib is one of the drugs commonly used in targeted therapies that inhibit B-RAF mutations. In clinical trials, vemurafenib is found to have excellent clinical results at first, but resistance emerges after a few months. Researches on drug resistance continue to mount. Our data analysis also found that vemurafenib had a dual effect: inhibiting the expression of harmful genes and promoting patient survival, but also promoting adverse events leading to drug resistance and tumor progression. Therefore, it is necessary to further study the drug resistance mechanism of vemurafenib and how to reduce its drug resistance, so as to shed light on the treatment of melanoma.

Authors' contributions:

Acquisition of data, analyzing and interpretation of data, drafting the article was done by Jiaheng Xie, Yuan Cao, Zhechen Zhu; Designing, revising and guiding the study was performed by Jingping Shi, Shujie Ruan and Ming Wang. The authors read and approved the final manuscript.

Acknowledgements:

We are very grateful for data provided by databases such as TCGA, GEO and GEPIA.

Conflict of interests:

The authors declared no conflicts of interest.

References

1. Schadendorf D, van Akkooi ACJ, Berking C, Griewank KG, Gutzmer R, Hauschild A, et al. Melanoma. Lancet 2018;392(10151):971-84.
2. International Agency for Research on Cancer, World Health Organization. GLOBOCAN 2012: estimated cancer incidence, mortality and prevalence worldwide in 2012; 2013.
3. Hagen B. Managing side effects of vemurafenib therapy for advanced melanoma. J Adv Pract Oncol 2014;5(6):400-10.
4. Niezgoda A, Niezgoda P, Czajkowski R. Novel approaches to treatment of advanced melanoma: a review on targeted therapy and immunotherapy. Biomed Res Int 2015;2015:851387.
5. Schadendorf D, Fisher DE, Garbe C, Gershenwald JE, Grob JJ, Allan Halpern A, et al. Melanoma. Nat Rev Dis Primers 2015;1:15003.     
6. Wan PT, Garnett MJ, Roe SM, Lee S, Niculescu-Duvaz D, Good VM, et al. Mechanism of activation of the RAF-ERK signaling pathway by oncogenic mutations of B-RAF. Cell 2004;116(6):855-67.
7. Shelledy L, Roman D. Vemurafenib: first-in-class BRAF-mutated inhibitor for the treatment of unresectable or metastatic melanoma. J Adv Pract Oncol 2015;6(4):361-5.
8. Manzano JL, Layos L, Bugés C, de los Llanos Gil M, Vila L, Martinez-Balibrea E, et al. Resistant mechanisms to BRAF inhibitors in melanoma. Ann Transl Med 2016;4(12):237.
9. Poulikakos PI, Persaud Y, Janakiraman M, Kong X, Ng C, Moriceau G, et al. RAF inhibitor resistance is mediated by dimerization of aberrantly spliced BRAF (V600E). Nature 2011;480(7377):387-90.
10. Chandrashekar DS, Bashel B, Balasubramanya SA, Creighton CJ, Ponce-Rodriguez I, Chakravarthi BV, et al. UALCAN: a portal for facilitating tumor subgroup gene expression and survival analyses. Neoplasia 2017;19(8):649-58.
11. Tang Z, Li C, Kang B, Gao G, Li C, Zhang Z. GEPIA: a web server for cancer and normal gene expression profiling and interactive analyses. Nucleic Acids Res 2017;45(W1):W98-102.
12. Garber K. Cancer research. Melanoma drug vindicates targeted approach. Science 2009;326(5960):1619.
13. DeSantis CE, Lin CC, Mariotto AB, Siegel RL, Stein KD, Kramer JL, et al. Cancer treatment and survivorship statistics, 2014. CA Cancer J Clin 2014;64(4):252-71.
14. Davies H, Bignell GR, Cox C, Stephens P, Edkins S, Clegg S, et al. Mutations of the BRAF gene in human cancer. Nature 2002;417(6892):949-54.
15. Menzies AM, Long GV. Systemic treatment for BRAF-mutant melanoma: where do we go next?. Lancet Oncol 2014;15(9):e371-81.
16. Lo JA, Fisher DE. The melanoma revolution: from UV carcinogenesis to a new era in therapeutics. Science 2014;346(6212):945-9.
17. Howe LR, Leevers SJ, Gómez N, Nakielny S, Cohen P, Marshall CJ. Activation of the MAP kinase pathway by the protein kinase raf. Cell 1992;71(2):335-42.
18. Jamal S, Ziff E. Transactivation of c-fos and beta-actin genes by raf as a step in early response to transmembrane signals. Nature 1990;344(6265):463-466.
19. Luo Z, Tzivion G, Belshaw PJ, Vavvas D, Marshall M, Avruch J. Oligomerization activates c-Raf-1 through a Ras-dependent mechanism. Nature 1996;383(6596):181-5.
20. Halaban R, Zhang W, Bacchiocchi A, Cheng E, Parisi F, Ariyan S, et al. PLX4032, a selective BRAF(V600E) kinase inhibitor, activates the ERK pathway and enhances cell migration and proliferation of BRAF melanoma cells. Pigment Cell Melanoma Res 2010;23(2):190-200.
21. Roskoski R Jr. RAF protein-serine/threonine kinases: structure and regulation. Biochem Biophys Res Commun 2010;399(3):313-7.
22. Pollock PM, Harper UL, Hansen KS, Yudt LM, Stark M, Robbins CM, et al. High frequency of BRAF mutations in nevi. Nat Genet 2003;33(1):19-20.
23. Dong J, Phelps RG, Qiao R, Yao S, Benard O, Ronai Z, et al. BRAF oncogenic mutations correlate with progression rather than initiation of human melanoma. Cancer Res 2003;63(14):3883-5.
24. Tsai J, Lee JT, Wang W, Zhang J, Cho H, Mamo S, et al. Discovery of a selective inhibitor of oncogenic B-Raf kinase with potent antimelanoma activity. Proc Natl Acad Sci U S A 2008;105(8):3041-6.
25. Bollag G, Hirth P, Tsai J, Zhang J, Ibrahim PN, Cho H, et al. Clinical efficacy of a RAF inhibitor needs broad target blockade in BRAF-mutant melanoma. Nature 2010;467(7315):596-9.
26. Chapman PB, Hauschild A, Robert C, Haanen JB, Ascierto P, Larkin J, et al. Improved survival with vemurafenib in melanoma with BRAF V600E mutation. N Engl J Med 2011;364(26):2507-16.
27. Prahallad A, Sun C, Huang S, Nicolantonio FD, Salazar R, Zecchin D, et al. Unresponsiveness of colon cancer to BRAF(V600E) inhibition through feedback activation of EGFR. Nature 2012;483(7387):100-3.
28. Corcoran RB, Ebi H, Turke AB, Coffee EM, Nishino M, Cogdill AP, et al. EGFR-mediated re-activation of MAPK signaling contributes to insensitivity of BRAF mutant colorectal cancers to RAF inhibition with vemurafenib. Cancer Discov 2012;2(3):227-35.
29. Zhang C, Bollag G. Scaffold-based design of kinase inhibitors for cancer therapy. Curr Opin Genet Dev 2010;20(1):79-86.
30. Pratilas CA, Solit DB. Targeting the mitogen-activated protein kinase pathway: physiological feedback and drug response. Clin Cancer Res 2010;16(13):3329-34.
31. Larkin J, Del Vecchio M, Ascierto PA, Krajsova I, Schachter J, Neyns B, et al. Vemurafenib in patients with BRAFV600 mutated metastatic melanoma: an open-label, multicentre, safety study. Lancet Oncol 2014;15(4):436-44.
32. Kim G, McKee AE, Ning YM, Hazarika M, Theoret M, Johnson JR, et al. FDA approval summary: vemurafenib for treatment of unresectable or metastatic melanoma with the BRAFV600E mutation. Clin Cancer Res 2014;20(19):4994-5000.
33. Nazarian R, Shi H, Wang Q, Kong X, Koya RC, Lee H, et al. Melanomas acquire resistance to B-RAF (V600E) inhibition by RTK or N-RAS upregulation. Nature 2010;468:973-7.
34. Deng W, Vashisht Gopal YN, Scott A, Chen G, Woodman SE, Davies MA. Role and therapeutic potential of PI3K‐mTOR signaling in de novo resistance to BRAF inhibition. Pigment Cell Melanoma Res 2012;25(2):248-58.
35. Kaplan FM, Shao Y, Mayberry MM, Aplin AE. Hyperactivation of MEK–ERK1/2 signaling and resistance to apoptosis induced by the oncogenic b-RAF inhibitor, PLX4720, in mutant n-RAS melanoma cells. Oncogene 2011;30(3):366-71.
36. Basile KJ, Abel EV, Aplin AE. Adaptive upregulation of FOXD3 and resistance to PLX4032/4720-induced cell death in mutant B-RAF melanoma cells. Oncogene 2012;31(19):2471-9.
37. Larkin J, Ascierto PA, Dréno B, Atkinson V, Liszkay G, Maio M, et al. Combined vemurafenib and cobimetinib in BRAF-mutated melanoma. N Engl J Med 2014;371(20):1867-76.
38. Niehr F, von Euw E, Attar N, Guo D, Matsunaga D, Sazegar H, et al. Combination therapy with vemurafenib (PLX4032/RG7204) and metformin in melanoma cell lines with distinct driver mutations. J Transl Med 2011;9(1):1-3.
39. Malumbres M, Barbacid M. Cell cycle, CDKs and cancer: a changing paradigm. Nat Rev Cancer 2009;9(3):153-66.
40. Malumbres M. Cyclins and related kinases in cancer cells. J BUON 2007;12:S45-S52.
41. Malumbres M, Barbacid M. Mammalian cyclin-dependent kinases. Trends Biochem Sci 2005;30(11):630-41.
42. Haneke K, Schott J, Lindner D, Hollensen AK, Damgaard CK, Mongis C, et al. CDK1 couples proliferation with protein synthesis. J Cell Biol 2020;219(3):e201906147.
43. Fedorenko IV, Abel EV, Koomen JM, Fang B, Wood ER, Chen YA, et al. Fibronectin induction abrogates the BRAF inhibitor response of BRAF V600E/PTEN-null melanoma cells. Oncogene 2016;35(10):1225-35.
44. Hirata E, Girotti MR, Viros A, Hooper S, Spencer-Dene B, Matsuda M, et al. Intravital imaging reveals how BRAF inhibition generates drug-tolerant microenvironments with high integrin β1/FAK signaling. Cancer cell 2015;27(4):574-88.
45. Koukoulis GK, Patriarca C, Gould VE. Adhesion molecules and tumor metastasis. Hum Pathol 1998;29(9):889-892.
46. Rapanotti MC, Campione E, Spallone G, Orlandi A, Bernardini S, Bianchi L. Minimal residual disease in melanoma: circulating melanoma cells and predictive role of MCAM/MUC18/MelCAM/CD146. Cell Death Discov 2017;3:17005.
47. Melnikova VO, Bar-Eli M. Bioimmunotherapy for melanoma using fully human antibodies targeting MCAM/MUC18 and IL-8. Pigment Cell Res 2006;19(5):395-405.
48. Kuphal S, Bosserhoff AK. E-cadherin cell–cell communication in melanogenesis and during development of malignant melanoma. Arch Biochem Biophys 2012;524(1):43-7.
49. Prakash R, Izraely S, Thareja NS, Lee RH, Rappaport M, Kawaguchi R, et al. Regeneration enhances metastasis: a novel role for neurovascular signaling in promoting melanoma brain metastasis. Front Neurosci 2019;13:297.