A Phosphoproteomic Comparison of B-RAFV600E and MKK1/2 Inhibitors in Melanoma Cells*

Inhibitors of oncogenic B-RAFV600E and MKK1/2 have yielded remarkable responses in B-RAFV600E-positive melanoma patients. However, the efficacy of these inhibitors is limited by the inevitable onset of resistance. Despite the fact that these inhibitors target the same pathway, combination treatment with B-RAFV600E and MKK1/2 inhibitors has been shown to improve both response rates and progression-free survival in B-RAFV600E melanoma patients. To provide insight into the molecular nature of the combinatorial response, we used quantitative mass spectrometry to characterize the inhibitor-dependent phosphoproteome of human melanoma cells treated with the B-RAFV600E inhibitor PLX4032 (vemurafenib) or the MKK1/2 inhibitor AZD6244 (selumetinib). In three replicate experiments, we quantified changes at a total of 23,986 phosphosites on 4784 proteins. This included 1317 phosphosites that reproducibly decreased in response to at least one inhibitor. Phosphosites that responded to both inhibitors grouped into networks that included the nuclear pore complex, growth factor signaling, and transcriptional regulators. Although the majority of phosphosites were responsive to both inhibitors, we identified 16 sites that decreased only in response to PLX4032, suggesting rare instances where oncogenic B-RAF signaling occurs in an MKK1/2-independent manner. Only two phosphosites were identified that appeared to be uniquely responsive to AZD6244. When cells were treated with the combination of AZD6244 and PLX4032 at subsaturating concentrations (30 nm), responses at nearly all phosphosites were additive. We conclude that AZD6244 does not substantially widen the range of phosphosites inhibited by PLX4032 and that the benefit of the drug combination is best explained by their additive effects on suppressing ERK1/2 signaling. Comparison of our results to another recent ERK1/2 phosphoproteomics study revealed a surprising degree of variability in the sensitivity of phosphosites to MKK1/2 inhibitors in human cell lines, revealing unexpected cell specificity in the molecular responses to pathway activation.

Inhibitors of oncogenic B-RAF V600E and MKK1/2 have yielded remarkable responses in B-RAF V600E -positive melanoma patients. However, the efficacy of these inhibitors is limited by the inevitable onset of resistance. Despite the fact that these inhibitors target the same pathway, combination treatment with B-RAF V600E and MKK1/2 inhibitors has been shown to improve both response rates and progression-free survival in B-RAF V600E melanoma patients. To provide insight into the molecular nature of the combinatorial response, we used quantitative mass spectrometry to characterize the inhibitor-dependent phosphoproteome of human melanoma cells treated with the B-RAF V600E inhibitor PLX4032 (vemurafenib) or the MKK1/2 inhibitor AZD6244 (selumetinib). In three replicate experiments, we quantified changes at a total of 23,986 phosphosites on 4784 proteins. This included 1317 phosphosites that reproducibly decreased in response to at least one inhibitor. Phosphosites that responded to both inhibitors grouped into networks that included the nuclear pore complex, growth factor signaling, and transcriptional regulators. Although the majority of phosphosites were responsive to both inhibitors, we identified 16 sites that decreased only in response to PLX4032, suggesting rare instances where oncogenic B-RAF signaling occurs in an MKK1/2-independent manner. Only two phosphosites were identified that appeared to be uniquely responsive to AZD6244. When cells were treated with the combination of AZD6244 and PLX4032 at subsaturating concentrations (30 nM), responses at nearly all phosphosites were additive. We conclude that AZD6244 does not substantially widen the range of phosphosites inhibited by PLX4032 and that the benefit of the drug combination is best explained by their additive effects on suppressing ERK1/2 signaling. Comparison of our results to another recent ERK1/2 phosphoproteomics study revealed a surprising degree of variability in the sensitivity of phosphosites to Mitogen-activated protein kinase (MAPK) 1 pathways orchestrate key intracellular responses to a variety of extracellular signals including mitogenic stimuli and cellular stress. In the case of the RAF/MKK/ERK pathway, receptor tyrosine kinases activate the small GTPase Ras, which then binds members of the RAF family of kinases (RAF1, B-RAF, ARAF) leading to their activation. Activated RAF kinases phosphorylate and activate MAP kinase kinases 1 and 2 (MKK1, MKK2), which in turn phosphorylate and activate extracellular signalrelated kinases 1 and 2 (ERK1, ERK2). The specificity of this cascade is remarkable, as the only widely accepted targets of B-RAF are MKK1/2, and the only validated targets of MKK1/2 are ERK1/2 (1)(2)(3)(4). Once activated, ERK1/2 mediates the effects of pathway activation by phosphorylating scores of cytoplasmic and nuclear targets. However, the full scope of cellular substrates of ERK1/2 remains unknown.
It is not necessarily intuitive that two inhibitors that target the same pathway should lead to improved patient responses. It has been suggested that the combination of B-RAF and MKK1/2 inhibitors may be more effective because it provides a barrier to mechanisms of acquired resistance (MOR) that reactivate ERK1/2 signaling downstream of B-RAF V600E (16,18). Another explanation for the improved patient response is that the combination is more effective at inhibiting ERK signaling below the threshold required to achieve a positive clinical response (19,20). This could be because of either an additive or synergistic effect of the combination on ERK signaling. A third possibility is that, in addition to their shared targets, B-RAF V600E and MKK1/2 each have small number of unique targets and that inhibition of all MAPK pathway targets is therefore only possible with the combination. By profiling the changes in phosphorylation in response to B-RAF V600E and MKK1/2 inhibitors, phosphoproteomics can provide insight into which of these possibilities is most likely.
Several proteomics and phosphoproteomics strategies have been employed to identify ERK1/2 targets, including 2D-PAGE (21,22), analog sensitive ERK1/2 kinases (23, 24), negative ionization mass spectrometry (MS) (25), and shotgun phosphoproteomics using stable isotope labeling by amino acids in cell culture (SILAC)-based MS (24,26,27) or labelfree MS (28,29). Surprisingly, the degree of overlap between the ERK1/2 targets identified in these studies has been reported to be very low (28). This may be reflective of low sampling and indicate that only a fraction of ERK targets are identified in each study. Alternatively, it may be reflective of the variability in responses between different cell types and/or different treatment conditions. Importantly, the degree of variability in phosphorylation responses to a single kinase inhibitor in different cell types remains largely unexplored. Delineating the degree to which cells respond similarly to ERK pathway inhibition will be critical for identifying universal biomarkers that can report activation or inhibition of the pathway.
Here we use large-scale phosphoproteomics to identify targets of the B-RAF/MKK/ERK pathway in a human metastatic melanoma cell line that is sensitive to the B-RAF V600E inhibitor PLX4032 (PLX; vemurafenib) and the MKK1/2 inhibitor, AZD6244 (AZD; selumetinib). In total, changes at 23,986 phosphosites on 4784 proteins could be quantified from rep-licate experiments and 1317 phosphosites showed reproducible decreases with inhibitor. Several new candidate substrates of ERK1/2 with potential roles in melanoma progression or melanoma cell survival were identified. Although the vast majority of phosphosites responded similarly to both AZD and PLX, 16 sites decreased only in response to PLX, suggesting that there are MKK1/2-independent targets of B-RAF V600E . Only two sites decreased in an AZD-specific manner, indicating that unique targets of MKK1/2 inhibitors are exceedingly rare and unlikely to explain the benefits of combination treatment. Combinatorial experiments using subsaturating concentrations of AZD and PLX showed additive responses at nearly all phosphosites, consistent with a Bliss independence model. Therefore, a more complete inhibition of ERK1/2 signaling is the most likely explanation for the clinical efficacy of combination treatment with B-RAF V600E and MKK1/2 inhibitors. Finally, we compared phosphosites that responded to AZD in our dataset to those in a prior study of an MKK1/2 inhibitor in two other cell lines. This comparison revealed a surprisingly large degree of variability in the phosphosite responses across cell lines, even for well-established targets of the ERK1/2 pathway.

EXPERIMENTAL PROCEDURES
Cell Culture-The WM239A cell line was a gift from Meenhard Herlyn, Wistar Institute. Cells were grown in SILAC-RPMI1640 media (Thermo Fisher Scientific) supplemented with 10% FBS (Gemini) and were maintained at 37°C with 5% CO 2 .
Sample Preparation-Low passage (p10 -15) WM239A cells were cultured for a minimum of five passages in SILAC-RPMI1640 (Thermo) supplemented with light, medium, or heavy lysine (Cambridge Isotopes; 40 g/ml) and arginine (Cambridge Isotopes; 200 g/ml). Following incorporation of labeled amino acids, 35 ϫ 10 6 cells were seeded overnight in 15 cm dishes. Media was changed the following morning and cells were treated with AZD, PLX, or AZDϩPLX for 2 h. Cells were rinsed twice with PBS and processed for MS analysis according to a modified FASP protocol (30). Cells were lysed in 750 l 95°SDT buffer (0.1 M Tris, pH 7.5, 4% SDS, 0.1 M dithiothreitol) and returned to 95°for 10 min. Lysates were then sonicated with a microtip sonicator for 15 s. Light, medium, and heavy lysates were combined and brought to a final volume of 30 ml with buffer UA (8 M urea, 0.1 M Tris, pH 8.5), and then split evenly into two Amicon Ultra-15 10K filter units (Millipore) and centrifuged for 25 min at 4000 ϫ g. Samples were washed with 10 ml UA, carbamidomethylated with 5 ml 50 mM iodoacetamide in UA, then centrifuged and washed three times with 5 ml UA. After transfer to new filter units (filter units exposed to the high concentrations of SDS in SDT buffer leaked overnight), samples were washed once with 7 ml 0.5 M ammonium bicarbonate to bring the urea concentration to 1 M. Samples were digested at 37°C with 2% (w/w) sequence grade modified trypsin (Promega) in a total volume of 3 ml. After 16 -18 h, digests were centrifuged and each filter unit was washed once with 3 ml H 2 O. Eluates from the two filter units were combined, acidified to pH Յ 2 and desalted on Oasis HLB extraction cartridges (150 mg, Waters). Five micrograms of each sample was reserved for total protein measurements. The remainder was subject to TiO 2 enrichment.
MS/MS was performed using an LTQ Orbitrap Velos, scanning MS between 400 -1600 m/z (1 ϫ 10 6 ions, 60,000 resolution), and selecting the 10 most intense ions for MS/MS with 180 s dynamic exclusion, 10 ppm exclusion width (20 ppm for phosphopeptide samples), and a repeat count of one. Ions with unassigned charge states or a charge state of ϩ1 were excluded. For phosphopeptide samples, neutral loss masses of 19.5948, 24.4937, 32.6584, and 48.9879 were used for data dependent MS 3 scans. Maximal injection times were 500 ms for Orbitrap precursor scans (one microscan) and 150 ms for LTQ MS/MS (one microscan) with AGC 1 ϫ 10 4 . The normalized collision energy was 35%, with activation Q ϭ 0.25 for 10 ms.
Data Analysis-Raw MS files for both phosphopeptide and total peptide data sets were uploaded together and searched against the False discovery rates were 1% for both phosphopeptide and protein identifications. For protein quantification, the minimum number of total peptides was two (unique ϩ razor). Phosphorylated peptides and their unmodified counterparts were excluded for protein quantification. Raw files for phosphopeptide data and total protein data were defined as separate groups. The multiplicity was set to three with medium labels of Arg6 and Lys4 and heavy labels of Arg10 and Lys8 for both groups. For phosphopeptide files, Phospho(STY) was set as a variable modification. All files were searched with carbamidomethylation (Cys) as a fixed modification and acetylation (N-term) and oxidation (Met) as variable modifications. Default settings were used except the "match between runs" feature was enabled with default settings and the site quantification was changed to "use highest change." Localization probabilities were assessed with the MaxQuant PTM score. The enzyme specificity was trypsin/P and two missed cleavages were allowed. For phosphopeptide data, the phospho(STY) output file was uploaded into Perseus (v1.4.1.3). Reverse and contaminant rows were removed and the site table was expanded so that ratios for phosphosites with different multiplicities of phosphorylated residues [e.g. singly (one phosphate), doubly (two phosphates), and triply (three or more phosphates) phosphorylated peptides] were each treated as unique phosphosites. Rows not quantified in any of the experiments after expanding the site table were removed. Phosphosite positions mentioned in the text have been manually validated to match the reviewed human Uniprot entry for the gene of interest. Phosphosite positions in supplemental tables are directly from the MaxQuant output files. For total peptide samples, the protein groups output file was uploaded into the Perseus and rows designated as reverse, contaminant, and "only identified by site" were removed. Rows not quantified in any of the experiments were also removed. The mass spectrometry proteomics data have been deposited to the ProteomeXchange Consortium (32) via the PRIDE partner repository with the data set identifier PXD001560 (10 M SILAC experiment) and PXD001563 (30 nM SILAC experiment).
STRING and Gene Ontology Analyses-Genes containing Ser/Thr-Pro phosphosites observed to decrease with 10 M AZD or 10 M PLX (log 2 Ͻ Ϫ0.84) were uploaded into STRING using the multiple names option, selecting human as the organism. The required confidence score was set to 0.900 and disconnected nodes were removed from the image.
The same genes were also uploaded into the PANTHER database selecting human as the organism. Functional classification of genes was viewed as a pie chart of Gene Ontology molecular function. The 49 genes identified as having nucleic acid binding transcription factor activity were manually validated to identify genes with known activity as transcription factors, transcriptional coactivators, or transcriptional repressors. This resulted in the final list of 39 genes.
Identifying Unique Targets of AZD and PLX-The 61 phosphosites with median PLX/AZD log 2 ratios exceeding the significance threshold of log 2 value Ϯ0.840 were filtered to exclude sites where both the AZD/DMSO and PLX/DMSO ratios had log 2 values Ͻ Ϫ0.840 (indicating differential regulation at a common target site). All phosphosites with PLX/AZD ratios that showed large variability in replicate experiments (i.e. both positive and negative PLX/AZD ratios in replicate experiments) were also removed.
Matching WM239A Data to Galan et al.-All quantified phosphosites in this study were first matched to Galan et al. based on common Uniprot identifier and amino acid position. For sites that could not be matched based on these criteria, an attempt was made to match them based on gene name and phosphopeptide sequence, or phosphopeptide sequence alone. Phosphopeptide sequences were the sequences in the phospho(STY) probabilities column with phosphorylations assigned at residues with probabilities Ն0.5. Sites were then filtered to accept only those in which the localization score in both studies was Ͼ0.75. For the significant sites identified in all three cell lines (879 total), 691 were matched by Uniprot identifier and amino acid position, 178 were matched by gene and phosphopeptide sequence, and 10 were matched by phosphopeptide sequence alone. Phosphosites were matched to regulatory sites in the Phosphosite-Plus database using the sequence window from the phospho(STY) text file. SILAC ratios reported for AZD/DMSO (this study) are reported as a median value for sites with replicate measurements.
Immunoblotting-Cells were lysed in SDT or RIPA buffer, and protein concentrations determined using the DC protein assay (Bio-Rad). Equal amounts of protein were separated on 10% SDS-PAGE gels or 4 -20% protean-TGX precast gradient gels (Bio-Rad) and transferred to PVDF-SQ membranes (Millipore) for 1 h at 100 V. Membranes were incubated in primary antibody for 2 h at room temperature or overnight at 4°C, and in secondary antibody for 1 h at room temperature.

Phosphoproteomics Responses to B-RAF
ERLIC has been reported to separate phosphopeptides relatively evenly across fractions, permitting higher numbers of phosphopeptide identifications (34,35). Consistent with these studies, we observed that phosphopeptides separated quite evenly across the 24 ERLIC fractions ( Fig. 2A). The resolution was also excellent, with the vast majority of phosphopeptides (85%) eluting in only one or two fractions (Fig.  2B). The majority of phosphopeptides identified were singly phosphorylated, with multiply phosphorylated peptides abundant only in fractions 22-24 (Fig. 2C). In three biological replicates, we were able to quantify 16,890, 18,833, and 18,103 phosphosites (Fig. 2D). In all, we quantified 23,986 phosphosites, 20,553 (86%) of which were class I sites (localization probability Ͼ0.75 and delta score Ͼ5) (36), showing that that most phosphosites could be localized with high confidence (supplemental Tables S1 and S2).
We performed three biological replicates, alternating SILAC labeling between experiments (Fig. 3A). At concentrations of 10 M, both AZD and PLX reduced phosphorylation of ERK1/2 at its regulatory phosphosites to levels undetectable by Western blotting (Fig. 3B). In addition, the SILAC peptide containing these phosphosites was undetectable under drug-treated conditions confirming that both inhibitors achieve nearly complete pathway inhibition at these concentrations (Fig. 3C). Changes in phosphorylation after inhibitor treatment were considered significant at a given phosphosite when the SILAC ratio reflected more than a 1.8-fold change (log 2 value Ϯ0.840). This threshold is the log 2 value exceeded by 1% of phosphosite ratio measurements in a control experiment in which light-, medium-, and heavy-labeled cells with were each treated with DMSO (supplemental Fig. S1). These values correspond to z-scores of Ϫ3.09 and ϩ3. 13.
In each biological replicate, 7-8% of phosphosites decreased significantly upon treatment with either AZD or PLX (Table I). As would be expected for treatment with kinase inhibitors, a smaller percentage of phosphosites were observed to increase in response to inhibitor treatment (1.6 -3.6%). There was a high degree of overlap between experimental replicates, and we focused our analysis on the 17,790 phosphosites that were quantified in at least two experiments (Fig. 2D). Responses at these sites were considered significant if the median value of the SILAC ratios from individual experiments exceeded the significance threshold of a 1.8-fold change (for sites quantified in only two replicates, the median value is the average of the two ratios) ( Table I).
Histograms of the median SILAC ratios (AZD/DMSO, PLX/ DMSO, and PLX/AZD) for these phosphosites centered at log 2 ϭ 0, and AZD/DMSO and PLX/DMSO histograms showed a clear enrichment for phosphosites that decreased following drug treatment (Fig. 4A). In total, 1317 phosphosites on 763 proteins decreased in response to either AZD or PLX. These included 664 phosphosites that contained the minimal ERK1/2 consensus phosphorylation motif (Ser/Thr-Pro), of which 277 sites contained the full consensus motif (Pro-Xxx-Ser/Thr-Pro). Nearly all drug-responsive phosphosites reflected true changes in phosphorylation. Only three phosphosites that decreased in response to inhibitor changed in a manner that could easily be explained by altered protein abundance (supplemental Fig. S2, supplemental Table S3). These were pSer74 , pThr189, and pSer191 on FOSL1, which is destabilized by the inhibitor-dependent dephosphorylation of pSer265 (37). Overall, we observed very high overlap between the phosphosites that changed in response to AZD and PLX ( Fig. 4B and 4C, blue symbols). Many phosphosites that scored as significant with only one inhibitor actually responded to both, but fell just below the 1% FDR cutoff with the other inhibitor (Fig. 4C, green symbols). Thus, the responses to inhibition of B-RAF and MKK1/2, representing different tiers of the B-RAF/MKK/ERK pathway, were remarkably uniform.
Cellular Networks Responding to AZD and PLX-Proteins with phosphosites that decrease in response to AZD and PLX were examined for high confidence interactions using the STRING database (38,39). We focused on proteins containing phosphosites with the minimal ERK1/2 consensus motif (Ser/ Thr-Pro), a subset that should be enriched in direct ERK1/2 substrates. Analysis of these proteins identified subnetworks representing the nuclear pore complex, regulators of transcription, growth factor signaling, small GTPase guanine nu- cleotide exchange factors (GEFs) and hydrolysis-activating proteins (GAPs), and proteins associated with the centrosome and mitotic regulation (Fig. 5).
The most interconnected network was the nuclear pore complex, of which NUP50, NUP98, NUP153, NUP214, and TPR have previously been identified as ERK1/2 substrates (22,28). In addition to these components, we identified sites that decreased on NUP35 and NUP188. Interestingly, the phosphosites that decreased in response to ERK pathway inhibition clustered to the cytoplasmic and nuclear FG nucleoporins (NUPs) as well as NUPs of the inner ring (Fig. 6A). Sites within the transmembrane ring, the linker, and outer ring NUPs were unaffected by inhibitor. Thus, there appears to be spatial clustering of the phosphosites regulated by AZD and PLX within the nuclear pore complex.
The largest interaction network identified by STRING analysis was comprised of transcriptional activators, coactivators, and repressors (Fig. 5). This network included many well characterized ERK1/2 substrates, such as ELK1 (40), JUN (41), FOSL1 (42), MYC (43), and STAT3 (4), and prompted us to interrogate the data set for additional transcriptional regulators that are responsive to MAPK pathway inhibition in melanoma. To do this, proteins harboring Ser/Thr-Pro phosphosites that decreased upon inhibitor treatment were sorted by Gene Ontology molecular function using the PANTHER classification system (44). Proteins classified as having nucleic acid transcription factor activity were examined further. This analysis identified 39 proteins with bona fide roles in transcription, 19 of which have previously been characterized as downstream targets of MAPK signaling (Fig 6B, supplemental Table S4). Importantly, of the 30 phosphosites that decreased in response to inhibitor on these 19 known targets, 21 have previously been shown to affect transcriptional activity, protein localization, or protein stability, showing that many of the sites identified have regulatory functions.
The remaining 20 transcriptional regulators are new candidates for ERK1/2 substrates. Fourteen of these proteins have at least one consensus docking motif for ERK1/2, which strongly predicts direct substrates (supplemental Table S4). Additionally, several of the phosphosites that decrease on these candidate targets occur in regions of the protein likely to affect protein function (Fig. 6B). The pSer73 site in PTOV1 is immediately adjacent to the PTOV-A domain, a domain that can affect transcriptional activity, although its role in PTOV1 function is currently unclear (45). The drug-responsive phosphosites in ATF6 (pSer13/pSer16) and TCF12 (pSer67) lie  within transactivation domains (46,47) suggesting B-RAF V600E may control their transcriptional activity. Similarly, pSer353 on NFIL3 occurs in a previously characterized repression domain in which phosphorylation might affect protein activity (48). Finally, the inhibitor-responsive phosphosite in LIMD1 (pSer424) lies within a region necessary for its interaction with Rb, suggesting that phosphorylation may affect LIMD1:Rb binding and repression of E2F target genes (49).
Mapping Inhibitor Responses onto the Human Kinome-We next examined the 17,790 phosphosites for which we had replicate measurements with the goal of addressing which protein kinases are controlled by B-RAF V600E signaling. In total, there were 904 phosphosites that matched 217 enzymes in the human kinome (50) (supplemental Table S5). Of these, 120 sites on 72 kinases showed significant responses to AZD or PLX. Nineteen of these sites have validated functions in regulating enzyme activity, localization, molecular interactions, or protein stability (supplemental Table S5). The majority of these sites occurred on kinases within the RAF/ MKK/ERK pathway, including A-RAF, B-RAF, RAF1, MKK1/2, ERK1/2, and p90RSK1/2/3. These included activating phosphorylation sites on MKK1 (pSer218, pSer222), MKK2 (pSer222, pSer226), ERK1 (pThr202, pTyr204), ERK2 (pThr185, pTyr187), and RSK2 (pSer369, pSer577). Each of these was inhibited in response to drug, as expected. A negative feedback site on RAF1 (pThr642) that is targeted by ERK1/2 (51) also decreased with inhibitor, reflecting loss of feedback inhibition following drug treatment.
Some responsive phosphosites occurred within protein kinases associated with other signaling pathways. Those with known regulatory function were linked to the activation of CHK1 (pSer280) and MKK4 (pSer257). These decreased with inhibitor treatment, indicating that B-RAF V600E enhances signaling through stress-response pathways. Likewise, activating sites in EPHA2 (pSer897) and p70S6K (pSer427) were inhibited by drug treatment, revealing the role for B-RAF V600E in promoting growth-regulatory mechanisms downstream of cell adhesion receptors and mTOR pathways. Drug-responsive phosphosite inhibition of GRK2/ADRBK1 (pSer670) suggests suppression of beta-adrenergic receptor signaling by B-RAF V600E . These results show crosstalk between oncogenic B-RAF V600E and other signaling pathways in melanoma. as MKK-independent targets of RAF1, the only known targets of B-RAF are MKK1/2 (1,2). Similarly, the only known targets of MKK1/2 are ERK1 and ERK2 (3). Nevertheless, the existence of unique targets of B-RAF V600E and MKK1/2 could potentially explain the improved clinical responses seen with the combination of B-RAF V600E and MKK1/2 inhibitors. Therefore, we examined phosphosites that were uniquely responsive to PLX or AZD. Of the 1317 phosphosites that decreased in response to either drug, 61 had PLX/AZD ratios that surpassed the significance threshold of a Ϯ1.8-fold change. After removing sites that decreased with both AZD and PLX but with different magnitudes, and sites with large variability between replicate experiments ("Experimental Procedures"), we identified 16 phosphosites FIG. 5. STRING network analysis of proteins harboring Ser-Thr/Pro phosphosites that decrease in response to inhibitor. STRING network of high-confidence interactions (minimum confidence score of 0.900) among proteins containing Ser/Thr-Pro phosphosites that decrease in response to AZD or PLX. Major subnetworks are highlighted in blue. Disconnected nodes are not shown. Color of the nodes is arbitrary. Large nodes indicate that structural information is available for that protein.

Phosphorylation Events Uniquely Responsive to PLX or AZD-Although a handful of proteins have been implicated
that reproducibly decreased only in response to PLX, and two phosphosites that decreased only in response to AZD (Table II) FIG. 6. Selected targets of B-RAF V600E signaling. A, The nuclear pore network, showing nuclear pore proteins harboring Ser/Thr-Pro phosphosites identified in this study (circles) and their positions within the nuclear pore complex (center). Phosphosites indicated in yellow were significantly inhibited by AZD or PLX, and phosphosites in black were unchanged. B, Transcriptional regulators harboring Ser/Thr-Pro phosphosites significantly inhibited by AZD or PLX. Top: Phosphosites on known substrates for ERK1/2. Bottom: Phosphosites on candidate substrates for ERK1/2. Color indicates the log 2 value of the median SILAC ratio for the phosphosite. For proteins in which multiple phosphosites decrease, the color reports the site with the largest change. Asterisks indicate validated regulatory sites where phosphorylation affects transcriptional activity, protein localization, or protein stability. Double asterisks indicate phosphosites located within or immediately adjacent to functional domains.

Targets of ERK1/2 with known roles in transcriptional regulation
The one phosphopeptide that would have been expected to respond differently to AZD and PLX is the MKK1/2 phosphopeptide harboring its activating sites, which is directly phosphorylated by B-RAF V600E , and would be expected to decrease only in response to PLX. Indeed, Western blotting showed that phosphorylation at this site decreased in response to PLX but increased slightly in response to AZD (Table II, Fig 7A). The other 15 sites identified as unique targets of PLX include four protein kinases (RIPK2, LIMK1, MAP2K4, STK10), a lipid kinase (PI4KB), a transcription factor (ATF2), and two membrane trafficking proteins (SNX13, VPS51) (Table II). A particularly intriguing target was RIPK2, a dual-specificity kinase best characterized for its role in innate immune signaling downstream of NOD1 and NOD2 (52). Previous work has shown that RAF1 activates RIPK2, which in turn directly phosphorylates and activates ERK1/2 (53). Identification of RIPK2 as a unique target of PLX suggests that B-RAF V600E may be capable of activating ERK1/2 in an MKK1/ 2-independent manner by regulating RIPK2. Sequence alignment of the 16 phosphosites identified as unique targets of PLX failed to identify a consensus motif around the phosphorylated residue (Fig. 7B). Thus, to the extent that these sites represent direct targets of B-RAF V600E , no easily identifiable consensus sequence was indicated.
The two phosphosites identified as unique targets of AZD occurred in the follistatin-related protein, FSTL1, and the RAF family member, A-RAF. Interestingly, the ARAF phosphopeptide exhibits a differential response to AZD and PLX only in its diphosphorylated form. Although both monophosphorylated peptides decrease in response to AZD and PLX, the diphosphorylated form decreases eightfold with AZD but increases as much as 16-fold with PLX. This behavior suggests that PLX  either facilitates A-RAF phosphorylation or protects it from dephosphorylation, perhaps by promoting A-RAF:B-RAF dimerization.
Responses to Combination Treatment with AZD and PLX-We next examined the effects of combination treatment with B-RAF V600E and MKK1/2 inhibitors by designing a triplelabeled SILAC experiment to compare the responses of melanoma cells to AZD, PLX, or the combination of AZD ϩ PLX (Fig. 8A). At 10 M concentrations, we observed no difference in phosphosite responses to the combination of AZD ϩ PLX versus AZD or PLX alone (supplemental Fig. S3). Thus, at 10 M concentrations, both AZD and PLX completely inhibit ERK signaling. To identify concentrations at which combinatorial effects might be observed, we treated cells with different concentrations of either inhibitor and monitored the decrease in phospho-ERK1/2 levels. At concentrations of 30 nM, each inhibitor reduced ERK1/2 phosphorylation well below the maximum level and within an intermediate range in which additive or synergistic interactions could be observed (Fig.  8B). SILAC experiments were then conducted at this concentration, measuring phosphosite responses to AZD, PLX, and the combination in duplicate (Fig. 8A). Because most phosphosite responses were small at 30 nM inhibitor (supplemental Tables S6 and S7), we evaluated responses at those sites most likely to represent real decreases, based on their overlap with the 1317 sites that clearly decreased with 10 M AZD or PLX (Fig. 8C). This yielded 319 sites where responses could be quantified in duplicate, and where inhibition was measurable with each drug individually. We used a Bliss independence model (54) to predict the additive response to the combination for each of these sites, based on the response of that site to each single drug. We calculated the predicted response using the equation: where R E is the predicted response and R A and R P are the individual responses measured for AZD and PLX, respectively. A histogram of the differences between the predicted responses to the combination and the measured responses showed a normal distribution centered at zero (Fig. 8C). Only four phosphosites exhibited differences that exceeded a 3 threshold. Manual inspection of the measured ratios for each of these four sites showed large variability between duplicate measurements suggesting these instances of significance may be spurious. Thus, at 30 nM concentrations, the molecular responses to the combination of AZD and PLX are additive, showing independence between the two inhibitors.
Variability in the Phosphorylation of RAF/MKK/ERK Signaling Targets-Finally, we examined the extent to which phosphorylation events targeted by inhibitors vary between different cell systems. To address this, we compared the phosphosite responses to 10 M AZD in our study against responses to 10 M PD184532, a MKK1/2 inhibitor used in large-scale phosphoproteomics experiments performed in human A375 cells expressing B-RAF V600E and HEK293 cells stimulated with PMA (27). In total, there were 4119 phosphosites that were quantified in all three cell lines from these two studies (supplemental Table S8). From these, we identified 879 class I phosphosites on 616 proteins that were inhibited by drug in at least one cell line (Ͼ1.8-fold decrease). These included 213 sites inhibited by AZD in WM239A cells, 360 sites inhibited by PD184352 in A375 cells, and 760 sites inhibited by PD184352 in HEK293 cells (Fig. 9A, supplemental Table S8). Surprisingly, only 142 of these 879 phosphosites passed the significance threshold in all three cell lines (Fig.  9A). This reveals that the majority of inhibitor-responsive phosphorylation events vary between cells. Because some of these differences might reflect phosphosites that narrowly miss the significance threshold in one or two cell types, we identified sites that were significantly responsive in one cell line (greater than 1.8-fold change; log 2 Ͻ Ϫ0.84) but unresponsive in others (less than 1.2-fold change; log 2 Ͼ Ϫ0.28). Using these criteria, we identified 217, 27, and 15 phosphosites that were uniquely responsive in HEK293 cells, A375 cells, and WM239A cells, respectively (Fig. 9B). The large number of phosphosite responses unique to HEK293 cells suggest that the spectrum of phosphorylation events regulated by the ERK1/2 pathway may be largely dependent on the cell type and/or mechanism of pathway activation. These unique responses do not appear to be because of different off-target effects of the two inhibitors given that both HEK293 and A375 cells were treated with PD184352.
Specific examples of phosphosite variation were instructive. NUP153 and NUP214 are both direct substrates for phosphorylation by ERK1/2 (22). Of the 10 phosphosites with inhibitor/control ratios that could be quantified in all three cell lines, seven showed significant responses in at least one cell line, but only two sites showed a significant response in all three cell lines (Fig. 9C). This reveals unexpected variability in nucleoporin phosphorylation, and implies that there is variation in the mechanism by which ERK1/2 regulates nuclear import in different cell systems. Other examples of proteins displaying variability in phosphorylation included the wellcharacterized ERK1/2 substrates, stathmin and cortactin. Here, three of five sites in stathmin, and three of four sites in cortactin showed varying responses to drug across the three cell lines (Fig. 9C). Two of these phosphorylation sites (pSer38 in STMN1 and pSer418 in CTTN) are direct substrates for ERK1/2 and also are key regulators of protein function. Taken together, our findings indicate that targets of a single pathway can show significantly different responses in different cell systems.

DISCUSSION
This study is the first of its scale to compare phosphoproteome responses to BRAF V600E and MKK1/2 inhibitors in melanoma cells. Our first objective was to evaluate the overlap between the molecular responses to B-RAF V600E and MKK1/2 inhibitors. Our results show that nearly all phosphosites that respond to AZD also respond to PLX. This high degree of overlap reinforces the idea that signaling from B-RAF V600E to MKK1/2 is remarkably linear. Our second objective was to use phosphoproteomics to understand the interactions between B-RAF V600E and MKK1/2 inhibitors and to provide molecular insight into why combination treatments lead to improved outcomes in clinical trials. Our results show that at subsaturating concentrations of inhibitors, phosphosite responses to combination treatment are additive, such that the magnitude of the responses to the combination can be predicted from the responses to each inhibitor individually in nearly all cases. This shows independence between BRAF V600E and MKK1/2 inhibitors and additive interactions in suppressing ERK signaling.
Collectively, these results provide insight into the clinical benefit of adding a MKK1/2 inhibitor to the B-RAF V600E inhibitor regimen. The fact that nearly all phosphosites that respond to AZD also respond to PLX make it unlikely that unique targets of MKK1/2 inhibitors account for the enhanced benefit of the combination. Instead, the interactions between the two drugs suggest that the benefit of the combination is best explained by the additive inhibition of ERK signaling. It has been reported that Ͼ80% inhibition of phospho-ERK1/2 is needed for an effective clinical response (20), and that this level of inhibition is difficult to achieve in many patients receiving single agent therapy with B-RAF V600E or MKK1/2 inhibitors (19). Our data suggest that the combination is likely more effective at reaching this level of inhibition than either drug individually.
Another common hypothesis used to explain the enhanced benefit of combination therapy is that the combination of inhibitors prevents the emergence of resistance mechanisms that still require MKK1/2 signaling (16,18). This would include mechanisms of resistance such as N-RAS mutations, B-RAF V600E amplification, B-RAF splice variants, or receptor tyrosine kinase up-regulation. However, it was recently reported that two of five patients developing acquired resistance to combination treatment harbored resistance mechanisms that reactivated ERK signaling upstream of MKK1/2 (55). One resistant tumor showed B-RAF V600E amplification, whereas the other expressed a B-RAF V600E splice variant (56)  (73) showing the overlap between phosphosites with significant responses to MKK1/2 inhibitor in WM239A cells (this study), A375 cells (27), or HEK293 cells (27). B, Area-proportional Euler diagram showing the overlap between phosphosites with a significant response in at least one cell line (log 2 Ͻ Ϫ0.84) and unchanged (log 2 Ͼ Ϫ0.28) in all lines where the response is not significant. C, Examples of proteins that display variability in phosphosite responses to MKK1/2 inhibitor. Charts above each example indicate the response in each cell line (displayed as log 2 values of SILAC ratios). The horizontal red line within each chart indicates the significance threshold of log 2 ϭ Ϫ0.84. known to confer resistance to B-RAF V600E inhibition. This suggests that the MKK1/2 inhibition achieved at the doses used in combination treatments is not sufficient for complete inhibition of MKK1/2 signaling and reinforces the notion that the enhanced benefit of the combination is likely because of the additive inhibition of ERK signaling.
Our high coverage of phosphosites was made possible by the use of the TiO 2 -ERLIC enrichment protocol, where batch TiO 2 enrichment prior to ERLIC fractionation significantly decreased sample processing time with no loss of phosphopeptide identifications. In addition, the formate buffer system used for ERLIC fractionation yielded phosphopeptides that were predominantly singly phosphorylated (Ͼ96%). This contrasts with triethylammonium phosphate (TEAP)-based ERLIC separations, in which as many as 75% of phosphopeptides are multiply phosphorylated (57)(58)(59)(60). Identification of singly phosphorylated peptides is advantageous for phosphoproteomics studies using MaxQuant because phosphosite ratios for these peptides reflect phosphorylation only at the site of interest and are easily interpreted. In contrast, phosphosite ratios that arise from multiply phosphorylated peptides can be derived from the SILAC ratios of more than one multiply phosphorylated peptide in which the site of interest is phosphorylated. Thus, quantifying changes in phosphosites present in multiply phosphorylated peptides requires manual analysis, which can be difficult with large data sets.
Our analysis of phosphosites that decreased in response to AZD and PLX led to the identification of potentially novel ERK sites on 20 transcription factors. Several of the 20 candidate target genes have previously been linked to melanoma survival or associated with melanoma progression (supplemental Table S4). PTOV1 (61) and the short isoform of PHF19 (62) are both overexpressed in advanced melanoma indicating that their regulation by B-RAF V600E might be important for disease progression. ATF6 has been shown to contribute to the intrinsic resistance of many melanoma cell lines to ER-stress induced apoptosis by upregulating antiapoptotic genes. Thus, phosphorylation of ATF6 within its transactivation domain by B-RAF V600E might increase transcription of antiapoptotic genes and promote cell survival (63,64). CRTC2 is a coactivator of ATF6 (65) and its regulation by B-RAF V600E may also modulate ATF6 activity to affect transcription and promote survival. Finally, the transcriptional repressor, TSC22D4, promotes the bypass of B-RAF V600E -mediated oncogene-induced senescence when overexpressed in fibroblasts (66). It is tempting to speculate that regulation of TSC22D4 by B-RAF V600E might affect a similar process in melanocytes to allow melanomas to bypass senescence.
Our identification of phosphosites that respond to PLX, but not to AZD, suggests that oncogenic B-RAF V600E has cellular targets that are independent of MKK1/2 and that signaling branch points between B-RAF V600E and MKK1/2 may be more prevalent than currently appreciated. Some of these targets may be therapeutically relevant. RIPK2, for example, has been shown to bypass MKK1/2 and directly phosphorylate ERK1/2 downstream of RAF1 (53). Thus, RIPK2 may play a role similar to COT1, a kinase that when overexpressed also directly phosphorylates and activates ERK1/2 (67). Further studies are needed to address the extent to which the unique responses to PLX or AZD represent direct substrates of B-RAF V600E and MKK1/2, and whether these phosphorylation events affect melanoma survival or progression.
In our analysis of phosphosite responses, all available ratios for phosphosites on multiply phosphorylated peptides were treated as unique. This increased the number of phosphosite ratios by 16.0% (from 20,760 to 23,986) and made it possible to compare phosphosite ratios for singly and multiply phosphorylated peptides to assess potential hierarchal dependences between neighboring phosphosites. This strategy revealed that the unique sensitivity of ARAF pSer186 to AZD was specific to the diphosphorylated form of the peptide, suggesting regulation at this site is influenced by phosphorylation of ARAF at pThr181.
The degree to which phosphorylation events targeted by specific pathways vary between different cell systems, or in response to different treatment conditions, is undetermined. Courcelles et al. recently reported a striking lack of overlap in responses to MKK1/2 inhibitors among five published studies (28). However, these comparisons were made at the level of proteins rather than phosphosites, and the five studies compared used cell lines from different species and differed widely in their methods for quantifying changes. Therefore, it was difficult to conclusively determine the extent of variation in phosphosite regulation. To better address this question, we compared our study against that of Galan et al. (27), who examined human HEK293 cells and A375 melanoma cells treated with 10 M PD184352. Both AZD and PD184352 bind to the same allosteric pocket in MKK1/2, and show similar effects in preclinical cancer models (68 -70). Despite this, we observed a surprising degree of heterogeneity in the inhibitor responses between these three cell lines. Of the phosphosites that could be quantified across all three cell lines, only 16% (142) passed our threshold for significance in all cells. In contrast, nearly 37% of sites (321) showed clear variation between cell lines. In addition, 25% (217) were responsive to MKK1/2 inhibition only in PMA-stimulated HEK293 cells suggesting that activated B-RAF V600E in A375 and WM239A cells regulates only a small subset of potential ERK targets. These comparisons reveal that even within a single species, the molecular responses to activation of a pathway depend heavily on both the cell line and the means of pathway activation.
Recent analyses of phosphoproteomics data sets have shown that orthologous proteins show large variations in phosphosite localization between species, even after accounting for sequence conservation (71,72). These analyses also showed that conserved phosphosites have regulatory functions more often than nonconserved sites. Therefore, we asked whether phosphosites that consistently respond to MKK1/2 inhibition across cell lines were also more likely to have regulatory function. To do this, we examined the 4119 phosphosites that could be quantified in all three cell lines, and found 237 that matched known human regulatory sites in the PhosphositePlus database (50). Sites with regulatory function accounted for 13 of the 142 sites (9.1%) that responded in all three cell lines, and 20 of the 333 sites (6.2%) that varied in at least one cell line. Thus, phosphosites that show consistent responses in all cell lines are not enriched for sites with regulatory function. Even with well-characterized ERK1/2 substrates such as NUP153/NUP214, STMN1, and CTTN, phosphorylation events controlling nuclear translocation, microtubule binding, and cell motility, respectively, show variability in their responses across cell systems. Determining the mechanisms by which phosphorylation responses vary is key to predicting responses to activation of a specific pathway. This is needed for the long-term goal of defining phosphorylation signatures that can serve as biomarkers for pathway activation.