Skip to main content


In-depth genomic data analyses revealed complex transcriptional and epigenetic dysregulations of BRAF V600E in melanoma

Article metrics



The recurrent BRAF driver mutation V600E (BRAF V600E) is currently one of the most clinically relevant mutations in melanoma. However, the genome-wide transcriptional and epigenetic dysregulations induced by BRAF V600E are still unclear. The investigation of this driver mutation’s functional consequences is critical to the understanding of tumorigenesis and the development of therapeutic strategies.

Methods and results

We performed an integrative analysis of transcriptomic and epigenomic changes disturbed by BRAF V600E by comparing the gene expression and methylation profiles of 34 primary cutaneous melanoma tumors harboring BRAF V600E with those of 27 BRAF WT samples available from The Cancer Genome Atlas (TCGA). A total of 711 significantly differentially expressed genes were identified as putative BRAF V600E target genes. Functional enrichment analyses revealed the transcription factor MITF (p < 3.6 × 10−16) and growth factor TGFB1 (p < 3.1 × 10−9) were the most significantly enriched up-regulators, with MITF being significantly up-regulated, whereas TGFB1 was significantly down-regulated in BRAF V600E, suggesting that they may mediate tumorigenesis driven by BRAF V600E. Further investigation using the MITF ChIP-Seq data confirmed that BRAF V600E led to an overall increased level of gene expression for the MITF targets. Furthermore, DNA methylation analysis revealed a global DNA methylation loss in BRAF V600E relative to BRAF WT. This might be due to BRAF dysregulation of DNMT3A, which was identified as a potential target with significant down-regulation in BRAF V600E. Finally, we demonstrated that BRAF V600E targets may play essential functional roles in cell growth and proliferation, measured by their effects on melanoma tumor growth using a short hairpin RNA silencing experimental dataset.


Our integrative analysis identified a set of BRAF V600E target genes. Further analyses suggested a complex mechanism driven by mutation BRAF V600E on melanoma tumorigenesis that disturbs specific cancer-related genes, pathways, and methylation modifications.


Next-generation sequencing has enabled us to identify numerous genetic alternations in melanoma genomes. These genetic alterations provide us with opportunities not only to investigate the novel insights into the molecular mechanisms of melanoma tumorigenesis but also to provide a new discovery basis for the identification of biomarkers for personalized targeted therapies [1-3]. So far, several driver genes including BRAF, NRAS, KIT, GNAQ, and GNA11 have been characterized and routinely used in clinical screenings for melanoma [4-6]. Other clinically relevant mutations or genes associated with those driver genes were systematically explored from 241 melanoma genomes [7]. Among these driver genes, the BRAF mutation at position 600 (BRAF V600) occurs in approximately 50% of melanoma patients, and among them, V600E accounts for approximately 79% [8]. The BRAF V600 in melanoma tumor genomes is currently one of the most clinically relevant mutation sites in melanoma [4,9]. Importantly, BRAF inhibitors such as vemurafenib and dabrafenib have been developed as targeted therapies for melanoma patients that harbor the BRAF V600E mutation. These compounds have provided tremendous clinical benefit to personalized cancer treatment; unfortunately, like other inhibitors, patients eventually develop resistance after treatment [10-12].

The exploration of the functional consequences of the transcriptional dysregulations of BRAF V600E is critical to the understanding of tumorigenesis and the potential discovery of targeted therapy. BRAF is part of the mitogen-activated protein kinase (MAPK) pathway that regulates cell growth and proliferation. The gain-of-function in BRAF V600E is well-known to highly activate the MAPK kinase pathway that promotes tumor cell growths in melanoma [13]. In recent years, several groups have explored the downstream genes promoted by BRAF V600E [14,15]. For example, Kannengiesser et al. [14] identified a few hundred genes associated with BRAF V600E through the differential analysis of microarray gene expression data in a survey of 69 human primary cutaneous melanoma tumors. Interestingly, they observed that most of those BRAF V600E regulated genes controlled by MITF were associated with over-expression. Through the investigation of transcriptome-wide changes using transduction BRAF V600E on primary human melanocytes, Flockhart et al. [15] have recently reported approximately one thousand mRNA transcripts that may be impacted by BRAF V600E.

Accumulating evidence has shown that aberrant methylation leads to the initiation and progression of tumorigenesis and this has been recognized as a hallmark of cancer [16,17]. Aberrant methylation in specific cancer genes has been reported to contribute to melanoma development [18-20]. Although the gain-of-function of BRAF V600E can promote specific target genes and pathways, to what extent the epigenetic modifications (i.e. DNA methylation) are involved and how to interplay within this process has been poorly understood.

The Cancer Genome Atlas (TCGA) project generated massive high-throughput genomic data, including mutation, DNA methylation, and transcription profiles for several hundred melanoma samples. These data provide us with an unprecedented opportunity for in-depth exploration of the functional consequences of a driver mutation (e.g., BRAF V600E) on tumors that integrate multiple types of genomic data. For this purpose, we performed an integrative analysis of the transcriptional and epigenetic alterations associated with a driver mutation (BRAF V600E) and applied it to the primary and metastatic tumor cutaneous melanoma samples available from TCGA.


Differential gene co-expression analyses identified putative targets of BRAF V600E

To identify the genes and related pathways perturbed by BRAF V600E, we developed a novel statistical approach, named Snowball, to identify differentially expressed genes based on their aggregated association between co-expression patterns and BRAF V600E mutation status. We identified the regulatory network modules that were significantly associated with BRAF V600E with a permutation p < 0.05, followed by a Weighted Gene Co-expression Network analysis (see Materials and methods, Figure 1A) [21]. As a result, a total of 711 putative target genes were identified including 330 down-regulated and 381 up-regulated genes (Additional file 1). Figure 2A shows a heat-map of expression patterns in the BRAF V600E and BRAF WT samples for those significantly associated genes identified by Snowball.

Figure 1

Workflow for identifying significantly altered transcriptional and epigenetic regulations associated with the BRAF V600E driver mutation in melanoma. Matched expression A) and methylation profiles B) are built for BRAF V600E and BRAF WT samples. Expression and methylation data in the column highlighted in the same color are derived from the same sample. We used the Snowball approach (top) to identify significantly and differentially expressed genes. A schematic demonstration of the Snowball approach for the identification of BRAF V600E putative targets is shown in the top right panel. Gene expression profiles for multiple cancer samples are measured in two groups, BRAF V600E and BRAF WT. All genes (g1-g7) can be powerfully detected based on their co-expression profiles from the BRAF mutation and wild-type groups (for details see Materials and methods). The LIMMA method was applied to detect differential methylation loci between BRAF V600E and BRAF WT samples (bottom).

Figure 2

Functional analysis of BRAF V600E target genes identified in primary tumor samples. A) Heat-map showing the differential signals for BRAF V600E target genes identified by Snowball approach. B) Enriched functional categories of the BRAF V600E target genes. * refers to genes in CGC catalogue.

Next, we used Ingenuity Pathway Analysis (IPA) to examine the functional categories and biological pathways of those putative BRAF V600E target genes. A significant portion of them were cancer-related (p < 1.6 × 10−8), including 20 genes from the Cancer Gene Census (CGC) catalogue (Figure 2B). In particular, cellular growth, proliferation, and development were found to be most overrepresented in molecular function, which was consistent with previous studies that found the BRAF mutation activates the MAPK pathway to facilitate various cellular processes [22-24] (Figure 2B). The most significantly enriched canonical pathway was the Aryl Hydrocarbon Receptor Signaling pathway (AHR pathway, p < 3.0 × 10−6), which belongs to the basic helix-loop-helix/Per-Arnt-Sim family of transcription factors (Figure 2B). This pathway has been reported to regulate xenobiotic metabolizing enzymes such as cytochrome P450 and has been demonstrated to cross-talk with the MAPK pathway [25,26]. Recent studies revealed that the AHR pathway is involved in various signaling pathways that are critical to cell proliferation and differentiation, gene regulation, cell motility and migration, and inflammation [27,28]. Another particular interesting pathway that was enriched was the IL-1 signaling pathway (p = 3.7 × 10−5), which had been reported to be dysregulated by BRAF V600E in a previous study [29]. This signaling pathway has also been shown to interact with the MAPK pathway [30,31] and contribute to multiple cancer progressions, including melanoma [32-34]. Taken together, our results indicate that BRAF V600E may regulate many genes and pathways that are crucial for melanoma development.

BRAF V600E target genes mediated by MITF and TGFB1

We next examined whether BRAF V600E target genes were regulated by specific up-regulators (i.e., transcription factors). The top two up-regulators identified using the IPA tool were oncogene MITF and tumor suppressor TGFB1; both were significantly enriched among BRAF V600E target genes (p < 3.6 × 10−16 and p < 3.1 × 10−9 for MITF and TGFB1, respectively; Figure 2B). Previous studies revealed that the BRAF mutation hyper-activated the MAPK signaling pathway and led to MITF promotion [35-37], whereas TGFB1 was reported to be down-regulated in multiple cancers, including melanoma [38-40]. Consistently, our results showed that MITF expression was also significantly higher in BRAF V600E than in BRAF WT samples, whereas TGFB1 showed significantly lower expression in BRAF V600E than in BRAF WT (Wilcoxon test, p < 0.05 and p < 0.01 for MITF and TGFB1, respectively; Figure 3).

Figure 3

MITF and TGFB1 dysregulated by BRAF V600E driver mutation in primary tumor samples. Boxplots show significantly higher expressions of MITF but lower expressions of TGFB1 in BRAF V600E, as compared to BRAF WT samples.

In addition, we repeated the Snowball analyses on TCGA melanoma metastatic samples. A total of 1010 putative BRAF V600Etarget genes were identified in the analysis, and we replicated a total of 213 (30%) BRAF V600E targets from the primary tumor samples (Additional file 2). Functional enrichment analysis of up-regulators using the IPA tool revealed that TGFB1 (p = 8.89 × 10−32) and MITF (p = 1.51 × 10−21) were again significantly and consistently enriched, suggesting that gain-of-function in BRAF V600E may generally influence the down-stream genes mediated by those two genes or pathways in different developmental stages of melanoma.

To further evaluate whether BRAF V600E target genes are mostly mediated by MITF, we collected 5,579 MITF target genes that were reported in a ChIP-Seq experiment and 732 MITF-induced targets inferred from a small interfering RNA (siRNA)-mediated MITF knockdown (siMITF) experiment in a melanoma cell line [41]. We found that genes targeted by MITF ChIP-Seq binding and siMITF-induced genes were more highly expressed overall than randomly selected background genes, regardless of BRAF V600E mutation status (Figure 4A, Wilcoxon test, p < 5.0 × 10−30 for all comparisons). Furthermore, a random subset of non-target genes with the same range of expression levels was selected as a background to compare to MITF ChIP-Seq binding and siMITF induced target genes, and the result showed that both MITF ChIP-Seq binding and siMITF induced targets showed a significantly higher expression change in BRAF V600E versus BRAF WT than the randomly selected background genes (Figure 4B; Wilcoxon test, p < 3.0 × 10−11 for all comparisons). These results support the conclusion that BRAF V600E leads to an increased in the level of the MITF gene, which likely subsequently results in the overall activation of many MITF target genes.

Figure 4

MITF and its target genes dysregulated by BRAF V600E driver mutation in primary tumor samples. A) and B) Boxplots show the relative gene expression level (median value) in tumor samples and fold changes (absolute value) between BRAF V600E and BRAF WT samples for:1) non-MITF target genes, 2) MITF ChIP-Seq binding targets), 3) MITF induced genes, and 4) MITF ChIP-Seq binding and induced targets.

Snowball identified BRAF V600E targets in response to BRAF inhibition

To further evaluate the identified BRAF V600E regulated genes, we analyzed a publicly available gene expression dataset of A375 melanoma cells that harbor the BRAF V600E mutation. This dataset contains the gene expression profiles before and after treatment with BRAF inhibitor vemurafenib (RAFi) [42]. Interestingly, we found that BRAF V600E regulated genes identified by Snowball from both the TCGA primary and metastasis tumor samples as well as the MTIF ChIP-Seq targets [41] showed a significant response when compared to randomly selected control genes (Figure 5). This suggested that BRAF V600E regulated genes identified by Snowball are highly reliable and BRAF V600E may regulate MITF targets likely mediated via MITF.

Figure 5

Boxplot of gene expression fold changes (absolute value) on A375 melanoma cells before and after the treatment of vemurafenib, for Snowball identified BRAF V600E targets from both TCGA primary and metastatic tumour samples, MITF targets from literature and control genes. P value for each comparison was derived from Wilcoxon test.

In particular, we also found that TGFB1 exhibited significantly elevated gene expression levels in cells with BRAF inhibitor induction, supporting that low TGFB1 expression level is associated with BRAF V600E. However, MITF itself exhibited an opposite trend for gene expression. Recent work by Konieczkowski et al. has suggested that most drug-sensitive cell lines exhibit high MITF expression and activity, but this was not observed in A375 cells based on the analysis of 29 BRAF V600E-mutant melanoma cell lines [21]. This discrepancy might indicate that MITF plays a complex role in melanoma drug response.

Global loss of DNA methylation associated with BRAF V600E

We next compared DNA methylation profiles between BRAF V600E and BRAF WT samples (see Materials and methods). A genome-wide DNA methylation loss was observed in BRAF V600E samples based on the comparison of DNA methylation profiles between BRAF V600E and BRAF WT (Figure 6A). After carrying out the differential DNA methylation analysis, we identified 523 aberrant methylation loci (i.e. CpG sites) using criteria of raw p < 1×10−3 and absolute intercept ≥ 0.2 (see Materials and methods). Surprisingly, 97.9% (512 of 523) showed hypomethylation in BRAF V600E relative to BRAF WT samples. This indicates a consistent, dominant loss of DNA methylation associated with BRAF V600E (Figure 6B). We further repeated the same differential DNA methylation analysis on 56 BRAF V600E versus 37 BRAF WT TCGA metastatic samples (see Materials and methods). The same trend of genome-wide DNA hypomethylation as well as significant aberrant methylation sites was found to be associated with BRAF V600E (Additional file 3). These results were consistent with the previous findings [43,44].

Figure 6

Methylation alterations associated with BRAF V600E driver mutations in primary tumor samples. A) Density plots of the median methylation intensity of each CpG site in BRAF V600E samples and BRAF WT samples. Methylation loci with Δβ > 0.1 were labeled with black dots. The plot indicates a global DNA methylation loss associated with BRAFV600. B) Heat-map showing differential methylation signals between BRAF V600E and BRAF WT samples, indicating a dormant methylation loss in BRAF V600E samples. C) Boxplots showing significantly lower expressions of DNMT3A in TCGA BRAF V600E versus BRAF WT samples. D) Boxplot of DNMT3A expression level on A375 melanoma cells before and after the treatment of vemurafenib.

To explore the possible mechanism associated with the genome-wide hypomethylation, we systematically examined whether BRAF V600E dysregulated any previously reported chromatin regulatory factors, such as DNA methyltransferases. Interestingly, only DNMT3A, functioning as a de novo DNA methyltransferase, showed significantly lower expression in BRAF V600E (Figure 6C) and was identified as a putative BRAF V600E target. In the preceding analysis, DNMT3A also exhibited significantly elevated expression levels in the A375 melanoma cells after being induced with the BRAF inhibitor (Figure 6D). It should be noted that no significantly differential gene expression patterns were observed based on the analysis of metastatic samples. One possible explanation is that DNMT3A might play a critical role in the initiation of tumorigenesis but may not be necessary in the later metastatic stage to maintain global hypomethylation. Taken together, these results suggest that BRAF V600E might initiate genome-wide epigenetic modifications through the regulation of DNMT3A, facilitating the initiation of melanoma tumorigenesis [35,37].

BRAF V600E targets associated with melanoma proliferation

We have shown that putative BRAF V600E target genes may play essential roles in melanoma tumorigenesis. To further verify the effects of those genes on cancer cell proliferation, we used a publicly available, large-scale gene silencing dataset from the short hairpin RNA (shRNA) screens of three melanoma cell lines (see Materials and methods, [45]). Among the 711 putative BRAF V600E target genes, we found that down-regulated genes significantly increased cell growth and proliferation, whereas up-regulated genes slightly decreased both (Figure 7). According to the effects of putative BRAF V600E targets on melanoma cell proliferation, we identified top tumor suppressor genes including TGFB1, TGFB1I1, PRODH, NAT6, ZNF205, ZNF142, FRS3, RUNX3, IGFBP5, HPGD, MAPK11, and NFIC, which significantly increased cell growth and proliferations. In contrast, top oncogenes including MET, BFSP1, CDH19, and ST6GALNAC3 were found to be associated with decreased cell growth and proliferations. In summary, our results suggest that BRAF V600E may play essential functional roles in cell growth and proliferation.

Figure 7

BRAF V600E target genes in melanoma cell proliferation. Boxplots of shRNA values (measuring melanoma relative cell proliferation) for BRAF V600E target genes in down-regulation (D), up-regulation (U), and control genes (C). BRAF V600E target genes in down-regulations showed a statistically significant increase of melanoma cell proliferation relative to the control, whereas up-regulation genes exhibited a slight decrease in melanoma cell proliferation.


To our knowledge, this study is among the first attempts at an in-depth exploration of the functional consequences of a single driver mutation using an integrative genomic data analysis strategy. We applied our recently developed Snowball approach and identified 711 putative BRAF V600E target genes, many of which are known to be involved in tumorigenesis. We further demonstrated that BRAF V600E might dysregulate specific cancer-related pathways and epigenetic modifications in melanoma tumorigenesis.

Although previous findings based on the differential gene expression analysis of BRAF V600E and BRAF WT samples provided novel insights into the understanding of BRAF-driven biology in melanoma [14,15], applying a sensitive detection method towards the functional consequences of the driver mutation BRAF V600E from a cohort of clinical samples is challenging. Traditional analyses used to detect significantly and differentially expressed genes are based on statistical parametric models, typically in a regression framework [46,47]. Those approaches assume expression independence among genes and apply a gene-by-gene strategy. However, analyses based on those approaches may result in ineffective detection of important gene or pathway targets due to the fact that a driver mutation is typically part of a small sample size situation and is expected to alter the expression of its cognate genes and genes in the same downstream pathways. The Snowball approach was implemented to meet the specific challenges in identifying the functional consequences of a driver mutation on clinical samples. Based on gene transcriptions and their interactions in the context of regulatory networks and the fact that a driver mutation may disturb the gene regulatory networks and generate differential co-expression profiles, the Snowball approach specifically utilizes a multivariate, distance-based regression to provide a more sensitive detection based on the co-expression profile of a given set of genes and its association with the mutation status. It is known that each tumor genome may have about 2–8 driver mutations, and their frequencies in the population are typically not high [48]. It is a real challenge to detect them using traditional differential signals. The Snowball approach is specifically developed to amplify the detection power by aggregating with resampling. Our simulation and real data analyses demonstrate that it is a more powerful approach for identifying large numbers of potential targets for downstream analyses. Additionally, the genetic background of patients and their tumor samples exhibit high heterogeneity with patient-specific and sample-specific variation. This heterogeneous predisposition to the driver mutation perturbation may lead to different gene expression patterns per gene from sample to sample as well as from patient to patient. Snowball utilizes a distance-based regression based on the gene co-expression profiles and assigns a robust ranking index to genes even when they have different predispositions [49].

This study revealed a key regulatory mechanism in melanoma, where BRAF V600E may play dual roles as a positive regulator of the MITF pathway and as a negative regulator of the TGFB1 pathway in the initiation of melanoma development. Recently, several studies have pinpointed an alteration of MITF in patients that may fail to eradicate tumors due to chemoresistance, which reactivates the MAPK signaling pathway [50-52]. Our work, together with those findings, highlights the potentially important role of MITF in melanomagenesis. In contrast to MITF, BRAF V600E represses the TGFB1 pathway, which may lead to the deactivation of the apoptosis process and the consequent cause of uncontrolled cell proliferation [53]. Moreover, DNMT3A, which also acts as a potential TGFB1 target [54], was found to be likely to mediate BRAF V600E epigenetic modifications in this study, thus facilitating melanoma development. While these findings are insightful, future studies using in vitro and in vivo assays are warranted to verify these results.

In conclusion, we performed an integrative analysis to exhaustively interrogate mutation, expression, and methylation datasets in an attempt to detect putative target genes and their regulations that are associated with the BRAF V600E driver mutation in melanoma. Our analyses identified not only known genes that contribute to melanoma pathogenesis but also many novel genes with potential clinical relevance. Importantly, our analysis indicated that a substantial proportion of the putative BRAF V600E target genes were significantly regulated by the transcription factor MITF and tumor suppressor TGFB1, suggesting that BRAF V600E may control specific cancer-related pathways via MITF and TGFB1 in order to initiate tumorigenesis. In particular, DNMT3A, one of the putative BRAF V600E targets, may reprogram epigenetic modifications to facilitate cancer development. These target genes were further shown to be essential in melanoma cell proliferation. Our analysis strategy provides a novel way to explore the functional consequences of a driver mutation and can be similarly applied to other driver mutations in complex diseases.

Materials and methods


We retrieved genomic data, including somatic mutations, from whole exome sequencing (total 385 samples), DNA methylation (n = 413), and RNA-Seq (n = 371) of cutaneous melanoma samples from the TCGA data portal ( We analyzed only those samples derived from primary tumors with matched mutation, expression, and methylation profiles in each sample. The DNA somatic mutation data (TCGA level 2) was retrieved from the TCGA somatic mutation annotation file (summarized as “maf” file), from which the BRAF V600E mutation status of each sample was examined and determined. We systemically examined the mutational profiles across all TCGA melanoma samples (N = 385). Since multiple driver mutations may co-exist on the same sample, the samples with known driver mutations including NRAS, CDKN2A, GNAQ, KIT and GNA11 have been removed from both the case and negative control groups to reduce the confounding effects. We finally included a dataset of 34 BRAF V600E samples and 27 BRAF WT samples; here, wide-type (WT) denotes pan-negative samples (those without any mutations in the above driver genes) (Figure 1). The gene-level expression data (TCGA level 3) was generated using Illumina HiSeq 2000 and measured by normalized RSEM (RNA-Seq by Expectation-Maximization) read counts. The DNA methylation data (TCGA level 3) was generated using the Illumina HumanMethylation450 BeadChip Array. Each methylation CpG locus was measured by a β value representing a ratio of M/(U + M), where M is the methylated probe intensity and U is the unmethylated probe intensity. The β value ranged from 0 to 1 (0: unmethylated; 1: fully methylated).

The transcription factor MITF’s binding targets from the ChIP-Seq data and its induced genes detected by the small interfering RNA (siRNA)-mediated MITF knockdown (siMITF) experiment were collected from a previous work [41].

DNA methylation analysis

We started methylation analysis from methylation profiles (β value, TCGA level 3) and then converted methylation β value to M value, which is compatible with the typical assumptions of linear models. We applied the R package Minfi [55] to detect differential methylation loci between BRAF V600E and BRAF WT samples. Significantly differentially methylated sites were detected used an F-test implemented in the function ‘dmpFinder’. The significantly aberrant methylation loci were identified by applying raw p value < 1 × 10−3 and absolute intercept ≥ 0.2.

Gene expression analysis

We developed the Snowball algorithm to identify a set of genes or gene modules that are likely regulated by a driver mutation [21]. This approach takes into account gene-gene interactions by evaluating each gene in a group of other genes. It is a more effective learning approach for the identification of functionally relevant genes or gene modules medicated by driver mutations that spread their genetic turbulence in the gene regulatory network to penetrate its functional impact. By applying the Snowball approach to the BRAF V600E and BRAF WT sample sets, we identified 1072 genes with significantly aggregated association with the mutation BRAF V600E. We further applied the Weighted Gene Co-expression Network Analysis [48] and identified 9 gene modules, each significantly associated with the BRAF mutation status when assessed using a generalized distance-based regression [21] at a permutation P < 0.05. The fold change of each gene’s expression in BRAF V600E relative to BRAF WT was calculated based on log2-transformed RSEM measurement.

To evaluate how Snowball identified BRAF regulated genes in response to the BRAF inhibitor vemurafenib, we also analyzed gene expression data on A375 melanoma cells harboring the BRAF V600E mutation from recent literature. This dataset contained the gene expression profiles before and after treatment with the BRAF inhibitor (GEO: GDS5085) [6]. Briefly, using LIMMA, we compared the gene expression profiles of A375 melanoma cells before and after treatment with the BRAF inhibitor, and the fold change of each gene was computed and reported. A total of 5000 genes were randomly selected from the genome as control genes for the comparative analysis.

Functional analysis

For the abovementioned 711 putative BRAF V600E regulated target genes, we examined their functional enrichment in gene networks and biological pathways, using the Ingenuity Pathway Analysis (IPA) tool ( The top 5 ranked gene networks and biological pathways were present.

Effect of gene silencing on cell proliferation using RNA interference data

To estimate the effect of an individual gene on cancer cell proliferation, we downloaded a comprehensive dataset from a genome-wide shRNA analysis of 10,941 genes (comprising of 52,209 probes) for three melanoma cell lines: A2058, HS944, and IGR39 (from the previous study) [45]. The effect of an individual gene’s silence for each of the three melanoma cell proliferations (measured by shRNA value) was computed using the log2 ratio of cell abundance in the pool generated by shRNA sequences at the endpoint, relative to the initial reference pool (details described in [45]).


  1. 1.

    Berger MF, Hodis E, Heffernan TP, Deribe YL, Lawrence MS, Protopopov A, et al. Melanoma genome sequencing reveals frequent PREX2 mutations. Nature. 2012;485:502–6.

  2. 2.

    Hodis E, Watson IR, Kryukov GV, Arold ST, Imielinski M, Theurillat JP, et al. A landscape of driver mutations in melanoma. Cell. 2012;150:251–63.

  3. 3.

    Gartner JJ, Parker SC, Prickett TD, Dutton-Regester K, Stitzel ML, Lin JC, et al. Whole-genome sequencing identifies a recurrent functional synonymous mutation in melanoma. Proc Natl Acad Sci U S A. 2013;110:13481–6.

  4. 4.

    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:949–54.

  5. 5.

    Curtin JA, Busam K, Pinkel D, Bastian BC. Somatic activation of KIT in distinct subtypes of melanoma. J Clin Oncol. 2006;24:4340–6.

  6. 6.

    Curtin JA, Fridlyand J, Kageshita T, Patel HN, Busam KJ, Kutzner H, et al. Distinct sets of genetic alterations in melanoma. N Engl J Med. 2005;353:2135–47.

  7. 7.

    Xia J, Jia P, Hutchinson KE, Dahlman KB, Johnson D, Sosman J, et al. A meta-analysis of somatic mutations from next generation sequencing of 241 melanomas: a road map for the study of genes with potential clinical relevance. Mol Cancer Ther. 2014;13:1918–28.

  8. 8.

    Lovly CM, Dahlman KB, Fohn LE, Su Z, Dias-Santagata D, Hicks DJ, et al. Routine multiplex mutational profiling of melanomas enables enrollment in genotype-driven therapeutic trials. PLoS One. 2012;7:e35309.

  9. 9.

    Rubinstein JC, Sznol M, Pavlick AC, Ariyan S, Cheng E, Bacchiocchi A, et al. Incidence of the V600K mutation among melanoma patients with BRAF mutations, and potential therapeutic response to the specific BRAF inhibitor PLX4032. J Transl Med. 2010;8:67.

  10. 10.

    Johannessen CM, Boehm JS, Kim SY, Thomas SR, Wardwell L, Johnson LA, et al. COT drives resistance to RAF inhibition through MAP kinase pathway reactivation. Nature. 2010;468:968–72.

  11. 11.

    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.

  12. 12.

    Villanueva J, Vultur A, Lee JT, Somasundaram R, Fukunaga-Kalabis M, Cipolla AK, et al. Acquired resistance to BRAF inhibitors mediated by a RAF kinase switch in melanoma can be overcome by cotargeting MEK and IGF-1R/PI3K. Cancer Cell. 2010;18:683–95.

  13. 13.

    Johannessen CM, Johnson LA, Piccioni F, Townes A, Frederick DT, Donahue MK, et al. A melanocyte lineage program confers resistance to MAP kinase pathway inhibition. Nature. 2013;504:138–42.

  14. 14.

    Kannengiesser C, Spatz A, Michiels S, Eychene A, Dessen P, Lazar V, et al. Gene expression signature associated with BRAF mutations in human primary cutaneous melanomas. Mol Oncol. 2008;1:425–30.

  15. 15.

    Flockhart RJ, Webster DE, Qu K, Mascarenhas N, Kovalski J, Kretz M, et al. BRAFV600E remodels the melanocyte transcriptome and induces BANCR to regulate melanoma cell migration. Genome Res. 2012;22:1006–14.

  16. 16.

    Jones PA, Baylin SB. The epigenomics of cancer. Cell. 2007;128:683–92.

  17. 17.

    Sandoval J, Esteller M. Cancer epigenomics: beyond genomics. Curr Opin Genet Dev. 2012;22:50–5.

  18. 18.

    Dawson MA, Kouzarides T. Cancer epigenetics: from mechanism to therapy. Cell. 2012;150:12–27.

  19. 19.

    Rodriguez-Cerdeira C, Molares-Vila A. New perspectives of “omics” applications in Melanoma research. Open Biochem J. 2011;5:60–6.

  20. 20.

    van den Hurk K, Niessen HE, Veeck J, van den Oord JJ, van Steensel MA, Zur Hausen A, et al. Genetics and epigenetics of cutaneous malignant melanoma: a concert out of tune. Biochim Biophys Acta. 1826;2012:89–102.

  21. 21.

    Konieczkowski DJ, Johannessen CM, Abudayyeh O, Kim JW, Cooper ZA, Piris A, et al. A melanoma cell state distinction influences sensitivity to MAPK pathway inhibitors. Cancer Discov. 2014;4:816–27.

  22. 22.

    Flaherty KT, Robert C, Hersey P, Nathan P, Garbe C, Milhem M, et al. Improved survival with MEK inhibition in BRAF-mutated melanoma. N Engl J Med. 2012;367:107–14.

  23. 23.

    Flaherty KT, Infante JR, Daud A, Gonzalez R, Kefford RF, Sosman J, et al. Combined BRAF and MEK inhibition in melanoma with BRAF V600 mutations. N Engl J Med. 2012;367:1694–703.

  24. 24.

    Falchook GS, Lewis KD, Infante JR, Gordon MS, Vogelzang NJ, DeMarini DJ, et al. Activity of the oral MEK inhibitor trametinib in patients with advanced melanoma: a phase 1 dose-escalation trial. Lancet Oncol. 2012;13:782–9.

  25. 25.

    Chiaro CR, Patel RD, Marcus CB, Perdew GH. Evidence for an aryl hydrocarbon receptor-mediated cytochrome p450 autoregulatory pathway. Mol Pharmacol. 2007;72:1369–79.

  26. 26.

    Borlak J, Jenke HS. Cross-talk between aryl hydrocarbon receptor and mitogen-activated protein kinase signaling pathway in liver cancer through c-raf transcriptional regulation. Mol Cancer Res. 2008;6:1326–36.

  27. 27.

    Feng S, Cao Z, Wang X. Role of aryl hydrocarbon receptor in cancer. Biochim Biophys Acta. 1836;2013:197–210.

  28. 28.

    Fan Y, Boivin GP, Knudsen ES, Nebert DW, Xia Y, Puga A. The aryl hydrocarbon receptor functions as a tumor suppressor of liver carcinogenesis. Cancer Res. 2010;70:212–20.

  29. 29.

    Khalili JS, Liu S, Rodriguez-Cruz TG, Whittington M, Wardell S, Liu C, et al. Oncogenic BRAF(V600E) promotes stromal cell-mediated immunosuppression via induction of interleukin-1 in melanoma. Clin Cancer Res. 2012;18:5329–40.

  30. 30.

    Adefuye AO, Sales KJ, Katz AA. Seminal plasma induces the expression of IL-1alpha in normal and neoplastic cervical cells via EP2/EGFR/PI3K/AKT pathway. J Mol Signal. 2014;9:8.

  31. 31.

    Liu X, Ye F, Xiong H, Hu D, Limb GA, Xie T, et al. IL-1beta Upregulates IL-8 Production in Human Muller Cells Through Activation of the p38 MAPK and ERK1/2 Signaling Pathways. Inflammation. 2014;37:1486–95.

  32. 32.

    Ma L, Lan F, Zheng Z, Xie F, Wang L, Liu W, et al. Epidermal growth factor (EGF) and interleukin (IL)-1beta synergistically promote ERK1/2-mediated invasive breast ductal cancer cell migration and invasion. Mol Cancer. 2012;11:79.

  33. 33.

    Saito M, Fan D, Lachman LB. Antitumor effects of liposomal IL1 alpha and TNF alpha against the pulmonary metastases of the B16F10 murine melanoma in syngeneic mice. Clin Exp Metastasis. 1995;13:249–59.

  34. 34.

    Tartour E, Adams D, Besancenot JF, Fridman WH, Schlumberger M. Poems syndrome with high interleukin (IL)6 and IL1 beta serum levels, in a patient with thyroid carcinoma and melanoma. Eur J Cancer. 1994;30A:893–4.

  35. 35.

    Garraway LA, Widlund HR, Rubin MA, Getz G, Berger AJ, Ramaswamy S, et al. Integrative genomic analyses identify MITF as a lineage survival oncogene amplified in malignant melanoma. Nature. 2005;436:117–22.

  36. 36.

    Wellbrock C, Marais R. Elevated expression of MITF counteracts B-RAF-stimulated melanocyte and melanoma cell proliferation. J Cell Biol. 2005;170:703–8.

  37. 37.

    Wellbrock C, Rana S, Paterson H, Pickersgill H, Brummelkamp T, Marais R. Oncogenic BRAF regulates melanoma proliferation through the lineage specific factor MITF. PLoS One. 2008;3:e2734.

  38. 38.

    Matsushita M, Matsuzaki K, Date M, Watanabe T, Shibano K, Nakagawa T, et al. Down-regulation of TGF-beta receptors in human colorectal cancer: implications for cancer development. Br J Cancer. 1999;80:194–205.

  39. 39.

    Lasfar A, Cohen-Solal KA. Resistance to transforming growth factor beta-mediated tumor suppression in melanoma: are multiple mechanisms in place? Carcinogenesis. 2010;31:1710–7.

  40. 40.

    Elliott RL, Blobe GC. Role of transforming growth factor Beta in human cancer. J Clin Oncol. 2005;23:2078–93.

  41. 41.

    Strub T, Giuliano S, Ye T, Bonet C, Keime C, Kobi D, et al. Essential role of microphthalmia transcription factor for DNA replication, mitosis and genomic stability in melanoma. Oncogene. 2011;30:2319–32.

  42. 42.

    Parmenter TJ, Kleinschmidt M, Kinross KM, Bond ST, Li J, Kaadige MR, et al. Response of BRAF-mutant melanoma to BRAF inhibition is mediated by a network of transcriptional regulators of glycolysis. Cancer Discov. 2014;4:423–33.

  43. 43.

    Hou P, Liu D, Dong J, Xing M. The BRAF(V600E) causes widespread alterations in gene methylation in the genome of melanoma cells. Cell Cycle. 2012;11:286–95.

  44. 44.

    Maddodi N, Bhat KM, Devi S, Zhang SC, Setaluri V. Oncogenic BRAFV600E induces expression of neuronal differentiation marker MAP2 in melanoma cells by promoter demethylation and down-regulation of transcription repressor HES1. J Biol Chem. 2010;285:242–54.

  45. 45.

    Cheung HW, Cowley GS, Weir BA, Boehm JS, Rusin S, Scott JA, et al. Systematic investigation of genetic vulnerabilities across cancer cell lines reveals lineage-specific dependencies in ovarian cancer. Proc Natl Acad Sci U S A. 2011;108:12372–7.

  46. 46.

    Smyth GK. Linear models and empirical bayes methods for assessing differential expression in microarray experiments. Stat Appl Genet Mol Biol. 2004;3(1):1–25.

  47. 47.

    Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010;26:139–40.

  48. 48.

    Cancer Genome Atlas Research N. Integrated genomic analyses of ovarian carcinoma. Nature. 2011;474:609–15.

  49. 49.

    Xu Y, Guo X, Sun J, Zhao Z. Snowball: resampling combined with distance-based regression to discover transcriptional consequences of a driver mutation. Bioinformatics. 2015;31:84–93.

  50. 50.

    Haq R, Yokoyama S, Hawryluk EB, Jonsson GB, Frederick DT, McHenry K, et al. BCL2A1 is a lineage-specific antiapoptotic melanoma oncogene that confers resistance to BRAF inhibition. Proc Natl Acad Sci U S A. 2013;110:4321–6.

  51. 51.

    Van Allen EM, Wagle N, Sucker A, Treacy DJ, Johannessen CM, Goetz EM, et al. The genetic landscape of clinical resistance to RAF inhibition in metastatic Melanoma. Cancer Discov. 2014;4:94–109.

  52. 52.

    Saez-Ayala M, Montenegro MF, Sanchez-Del-Campo L, Fernandez-Perez MP, Chazarra S, Freter R, et al. Directed phenotype switching as an effective antimelanoma strategy. Cancer Cell. 2013;24:105–19.

  53. 53.

    Sharpless E, Chin L. The INK4a/ARF locus and melanoma. Oncogene. 2003;22:3092–8.

  54. 54.

    Zhang Q, Chen L, Helfand BT, Jang TL, Sharma V, Kozlowski J, et al. TGF-beta regulates DNA methyltransferase expression in prostate cancer, correlates with aggressive capabilities, and predicts disease recurrence. PLoS One. 2011;6:e25168.

  55. 55.

    Aryee MJ, Jaffe AE, Corrada-Bravo H, Ladd-Acosta C, Feinberg AP, Hansen KD, et al. Minfi: a flexible and comprehensive Bioconductor package for the analysis of Infinium DNA methylation microarrays. Bioinformatics. 2014;30:1363–9.

Download references

Grant support

This work was partially supported by Ingram Professorship Funds (to Z.Z.) and National Institutes of Health grants (R01LM011177 and P30CA68485). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Author information

Correspondence to Zhongming Zhao.

Additional information

Competing interests

The authors declare that they have no competing interests.

Authors’ contributions

ZZ and XG conceived and designed the project. XG collected the data. XG and YX performed the data analysis. XG, ZZ, and YX drafted the manuscript. All authors read and approved the final manuscript.

Additional files

Additional file 1:

A total of 711 BRAF V600E target genes were listed based on the analysis of primary samples.

Additional file 2:

A total of 1010 BRAF V600E target genes were listed based on the analysis of primary samples.

Additional file 3:

Methylation alterations associated with BRAF V600E driver mutations based on TCGA metastatic samples. A) Density plots of the median methylation intensity of each CpG site in BRAF V600E samples and BRAF WT samples. Methylation loci with Δβ > 0.1 were labeled with black dots. B) Heat-map showing differential methylation signals between BRAF V600E and BRAF WT samples, indicating a dominant methylation loss in BRAF V600E samples.

Rights and permissions

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Guo, X., Xu, Y. & Zhao, Z. In-depth genomic data analyses revealed complex transcriptional and epigenetic dysregulations of BRAF V600E in melanoma. Mol Cancer 14, 60 (2015) doi:10.1186/s12943-015-0328-y

Download citation


  • Melanoma
  • Expression
  • DNA methylation
  • Driver mutation
  • BRAF
  • MITF
  • TGFB1
  • DNMT3A