Theoretical analysis reveals a role for RAF conformational autoinhibition in paradoxical activation
Abstract
RAF kinase inhibitors can, under certain conditions, increase RAF kinase signaling. This process, which is commonly referred to as ‘paradoxical activation’ (PA), is incompletely understood. We use mathematical and computational modeling to investigate PA and derive rigorous analytical expressions that illuminate the underlying mechanism of this complex phenomenon. We find that conformational autoinhibition modulation by a RAF inhibitor could be sufficient to create PA. We find that experimental RAF inhibitor drug dose–response data that characterize PA across different types of RAF inhibitors are best explained by a model that includes RAF inhibitor modulation of three properties: conformational autoinhibition, dimer affinity, and drug binding within the dimer (i.e., negative cooperativity). Overall, this work establishes conformational autoinhibition as a robust mechanism for RAF inhibitordriven PA based solely on equilibrium dynamics of canonical interactions that comprise RAF signaling and inhibition.
Editor's evaluation
This important study uses mathematical modeling to demonstrate that conformational autoinhibition of the RAF kinase is an important feature of its paradoxical activation by pharmacological inhibitors. The provided theoretical analysis is highly compelling and increases the mechanistic understanding of how regulation by the Nterminal regulatory domains of RAF contributes to the paradoxical effect of the inhibitors. This article is poised to be of interest to biochemists, cell biologists, and cancer biologists, ultimately paving the way for improved RAF therapeutics.
https://doi.org/10.7554/eLife.82739.sa0Introduction
Cancer is a disease that is characterized by a collection of shared phenotypes (Hanahan and Weinberg, 2011). Among these phenotypes is selfsufficiency in the progrowth signals that aberrantly drive the continued proliferation of the cancerous cells. These excessive progrowth signals may result from somatic mutations within a proliferation signaling pathway, and they could also result from excessive stimulation of transmembrane receptor proteins. Such progrowth, proliferation, signals commonly pass through RAS GTPases (KRAS, NRAS, and HRAS) and RAF kinases (BRAF, CRAF, and ARAF) to initiate the RAF/MEK/ERK mitogenactivated protein kinase (MAPK) cascade.
Pharmacological agents that counteract these progrowth signals have proven to be therapeutically beneficial in cancer. FDAapproved agents for upstream receptor tyrosine kinases that can drive proliferation through RAS/RAF signaling have proven effective for a variety of cancers. For example, EGFR inhibitors (both small molecule and therapeutic antibodies) have established benefit in lung, colorectal, and head and neck adenocarcinomas (Ciardiello and Tortora, 2008). Downstream from receptor tyrosine kinases, the RAS GTPases have only recently become targetable, although for only one specific somatic mutation (KRAS G12C) (Skoulidis et al., 2021). Multiple drugs have been developed to inhibit the RAF kinases (Bollag et al., 2010; Nakamura et al., 2013; Peng et al., 2015; Zhang et al., 2015), and several of these agents have proven clinically valuable for RAF mutant melanoma and colorectal cancer.
The development of a drug that successfully inhibits its target protein is not sufficient to produce an effective drug. The drug must also avoid significant effects in nondiseased cells and tissues. This can be difficult, considering that the targeted protein could also be produced in the noncancerous cells of the treated patient. For example, it was originally hoped that RAF inhibitors would be able to block the transmission of RAS signals (Figure 1A). Unexpectedly, RAF inhibitors were found to amplify RAS signals within the RAF wildtype context through a process that is now commonly referred to as ‘paradoxical activation’ (PA) (Hatzivassiliou et al., 2010; Heidorn et al., 2010; Poulikakos et al., 2010). PA is also associated with the development of drugdependent squamous cell carcinomas in a subpopulation of those treated with RAF inhibitors (Su et al., 2012). Despite numerous studies, the mechanisms driving PA are still not fully understood (Karoulia et al., 2017; Karoulia et al., 2016; Peng et al., 2015).
A mathematical model of biochemical processes allows one to rigorously analyze what behaviors are possible for a given set of reaction mechanisms. Modeling can lead to mechanismdriven hypotheses for further experimental testing (Kholodenko, 2015; McFall et al., 2019; Rukhlenko et al., 2018; Stites et al., 2015; Stites et al., 2007; Wofsy et al., 1992). Some of these mathematical models have proven useful for understanding mechanisms of anticancer agents. For example, we recently used a mathematical model to identify the mechanism (McFall et al., 2019; McFall et al., 2020; McFall and Stites, 2021) that explains why colorectal cancer patients with a KRAS G13D mutation benefit from treatment with EGFR inhibitors despite KRAS mutations more commonly causing resistance to EGFR inhibitors (De Roock et al., 2010). These studies are significant because a mechanistic explanation for this highprofile problem evaded determination through traditional experimental methods for approximately a decade. Thus, mathematical models can complement experimental approaches to provide new understandings to clinically relevant biological systems.
Mechanistically, two behaviors are commonly implicated as contributing to the phenomenon of PA (Figure 1B). First, some RAF inhibitors have been shown to result in an increased level of RAF dimerization (Hatzivassiliou et al., 2010; Jin et al., 2017; Karoulia et al., 2016; Lavoie et al., 2013). This druginduced dimer potentiation (DP) is commonly thought of as manifesting in a higher affinity between RAF protomers when one (or both) are bound to a RAF inhibitor (Kholodenko, 2015). Second, many RAF inhibitors do not appear capable of binding to both protomers of a RAF dimer equally well (Karoulia et al., 2017; Karoulia et al., 2016; Peng et al., 2015). As only one protomer in a RAF dimer need be signaling competent for RAF signaling to propagate (Heidorn et al., 2010; Hu et al., 2013; Yuan et al., 2018), such negative cooperativity (NC) is posited to result in a reduced ability to fully inhibit RAF signaling. This NC is commonly thought of as manifesting in a reduced affinity for the binding of a second RAF inhibitor relative to the affinity with which the first RAF inhibitor bound (Kholodenko, 2015). Previous theoretical studies of PA have focused on these two mechanisms (Kholodenko, 2015; Rukhlenko et al., 2018).
Here, we report our mathematical analysis of RAF PA. We developed a mathematical model that describes key regulatory steps in RAF signaling. We developed the model to follow biochemical and thermodynamic principles, and we derive the behaviors that logically follow from the RAF regulatory mechanisms. Our model includes three different mechanisms that contribute to PA; two have received previous theoretical attention and one (RAF autoinhibitory conformational regulation) has not previously received a mathematical consideration. Our modeling reveals that, under certain conditions, RAF autoinhibitory conformational changes and their modulation by RAF inhibitor binding can be sufficient to drive PA. We fit our model to experimental data, and the fit model suggests links between the three types of RAF inhibitors and the three mechanisms that contribute to PA. Overall, this study suggests that all three mechanisms and their modulation by RAF inhibitor binding are necessary to explain PA.
Results
Analytical modeling of RAF autoinhibition
The regulation of RAF kinase activation is complex with multiple regulatory steps (Hu et al., 2013; Karoulia et al., 2016; Lavoie and Therrien, 2015; Wellbrock et al., 2004, Brummer and McInnes, 2020), and several of these processes have been described to play a role in PA (Cho et al., 2012; Heidorn et al., 2010; Holderfield et al., 2013; Kondo et al., 2019; McKay et al., 2011). Conformational changes of the RAF monomer contribute significantly to RAF kinase activation (Jin et al., 2017; Lavoie and Therrien, 2015; Martinez Fiesco et al., 2022; Wellbrock et al., 2004; Zhang et al., 2021). In the ‘autoinhibited’ form, associations between its Nterminus and its kinase domain maintain RAF in an inactive form that does not dimerize (Cutler et al., 1998; Lavoie and Therrien, 2015). In the ‘nonautoinhibited’ form, the kinase domain is no longer occluded, and other regulatory mechanisms contribute to full RAF kinase activation, such as kinase domain conformational changes and dimerization (Lavoie and Therrien, 2015; Lavoie et al., 2013). Recent experimental work reports that RAF inhibitors tend to promote a net transition to the nonautoinhibited conformation that is bound to RASGTP (Jin et al., 2017; Karoulia et al., 2017). It has previously been suggested that this biasing to the nonautoinhibited state may contribute to PA (Figure 1B; Jin et al., 2017). However, conformational autoinhibition (CA) has not previously received a theoretical analysis to evaluate whether, and how, it may contribute to PA.
We therefore developed a series of mathematical models that include CA of RAF, NC within a RAF dimer, and RAF DP. We systematically studied different subsets of the full mechanism to investigate how each contributes to PA. For each submodel, we derived analytic equations that provide the equilibrium solution for the system. One benefit of an analytic expression is that it can be utilized to determine the generality of a result and/or to identify specific criteria required for a given behavior to occur. Conclusions from most of our models are presented in theoremproof form in Appendix 1 and are thereby devoid of any sensitivities arising from specific parameter choices. We refer to these conclusions as arising from ‘analytical modeling.’ We complement our analysis with ‘numerical modeling’ by using specific parameter values to create representative figures that can illuminate some of the predictions for the modeled network motif.
Our first mathematical model focuses on RAF CA and RAF dimerization (Figure 2A). The model allows RAF to adopt two different conformations: one is autoinhibited and can neither dimerize nor bind drug, and the other is nonautoinhibited and can bind drug and/or dimerize (Lavoie et al., 2013). Drugbound RAF is assumed to only be able to transition back to an autoinhibited state only after the bound drug has dissociated. Within the model, wildtype RAF is implicitly assumed to be activated by RASGTP as binding to RASGTP is an essential step to wildtype RAF activation (Lavoie and Therrien, 2015). We define active RAF as the RAF protomers that are not bound to a drug and are part of a RAF dimer.
The mechanism shown in Figure 2A and its conversion to a set of equations that provide the system’s equilibrium solutions is detailed in Appendix 1. We use the principle of detailed balance (Kholodenko, 2015; Wofsy et al., 1992) and total protein and drug conservation equations to derive closedform, analytic expressions for the steadystate solution for this system (Table 1, with derivations in Appendix 1). We focus on steadystate solutions because PA is a steadystate phenomenon, as is reflected in the longterm outgrowth of secondary tumors with RAS mutations from patients treated with a RAF inhibitor (Su et al., 2012).
Paradoxical activation is a robust outcome of conformational autoinhibition
One benefit of an analytic expression is that it can be utilized to determine the generality of a result and/or to identify specific criteria required for a given behavior to occur. Along these lines, we investigated whether and under what conditions the mechanisms in our CA base model are sufficient to create PA. We can demonstrate analytically that the presence of both CA and stabilization of the active form by RAF inhibitors is sufficient to create PA for some, but not all, parameter values (Figure 2B, Table 1, Appendix 1). Additionally, our analysis yields algebraic expressions that define the conditions necessary for PA to occur by this mechanism (Table 1). One condition we find is that RAF states must favor the autoinhibited form of RAF at equilibrium and in the absence of a drug for PA to occur. In other words, PA does not always happen when CA is present. By substituting biophysical parameter values into these expressions, we can predict PA for a wide range of RAF concentration values and visualize the range of system parameters for which CA promotes PA (Figure 2C and D, Figure 2—figure supplement 1). In our investigations, we also analytically show that PA is not generated by simpler mechanisms such as one without the autoinhibitory conformational changes and one without dimerization (Appendix 1).
Paradoxical activation reflects a shifting balance of signaling complexes
To communicate how PA arises from this mechanism, we describe the proportion of RAF in each of its possible states: (i) autoinhibited RAF monomer, (ii) nonautoinhibited RAF monomer that is not bound to drug, (iii) nonautoinhibited monomer that is bound to drug, (iv) RAF dimer with no drug bound, (v) RAF dimer with one of two kinase domains bound to drug, and (vi) RAF dimer with both kinase domains bound to drug. We considered the total amount of kinase activity to be the number of RAF protomers within a dimer that are not bound to drug (Figure 2—figure supplement 2).
Before a drug is given, a significant fraction of RAF is autoinhibited and there are low levels of nonautoinhibited RAF and RAF dimers. As RAF inhibitor levels increase, the level of autoinhibited RAF progressively declines. Nonautoinhibited RAF distributes between drugbound monomeric and dimeric forms while the unbound monomeric form maintains equilibrium with the autoinhibited RAF. The increased quantity of RAF dimers reflects the increased pool of RAF proteins that are nonautoinhibited and therefore capable of dimerization. This results in a drugdependent increase in RAF dimers bound to drug in one site and thereby the increase in total RAF kinase activation that accounts for PA. The quantity of drugbound RAF monomer and doubly drugbound RAF dimer progressively increases to saturation as a function of the drug amount, resulting in the eventual reduction in RAF kinase activity that is associated with PA dose–responses.
Evaluation of distinct mechanisms that promote paradoxical activation
There is evidence that NC and RAF inhibitordriven DP also contribute to PA. Previous theoretical work on PA has focused on the contribution of DP and NC (Kholodenko, 2015; Rukhlenko et al., 2018). We were curious how CA might interact with these other processes that contribute to PA. We developed a mathematical model that combined the mechanisms studied in these previous studies with our model of CA; the complete mechanism that is modeled is shown in Figure 3A, and we adapt the same nomenclature for parameters (Kholodenko, 2015). Analysis of our model with CA, NC, and inhibitordriven DP finds that these three mechanisms are mathematically distinct, in that one cannot translate or rescale variables to interchange between the mechanisms (Table 1). Further, each mechanism brings a qualitatively distinct feature to the dose–response of active RAF protomers (Figure 3—figure supplement 1A and B) and widely diverging capacity to potentiate the fold change in PA (Figure 3—figure supplement 1A) and the width of PA (Figure 3—figure supplement 1B). Because of these differences, we reasoned that a mathematical model that includes these biochemically, mathematically, and qualitatively distinct mechanisms may enable the underlying mechanisms of PA to be inferred, as in an inverse problem, through computationally fitting our model to data that characterizes the PA response of RAF inhibitors.
We developed a computational process for obtaining sets of parameters that fit our model to published experimental data (Karoulia et al., 2016). Our approach finds multiple parameter sets for the same drug that are within 10% of the error for the parameter set with the least error. These published data characterize the PA response of SKMEL2 cells to nine different RAF inhibitors, each at seven different concentrations of inhibitor (Figure 3B). We reasoned that some parameters should be specific to the individual drugs (i.e., the dissociation constant for the drug binding to RAF, and the druginduced DP and NC parameters) and that some parameters should be the same no matter which drug was used in the experiment (i.e., the abundance of RAF, the equilibrium constant of RAF dimerization when no drug is bound, and the CA equilibrium constant). We set values for two key parameters, RAF abundance and RAF dimerization equilibrium constant, based on estimates and do not fit them. A table listing the bestfit parameters of each of the submodels is provided in Supplementary file 1. Of note, we have 1 global parameter (the autoinhibition equilibrium constant) and 27 drugspecific parameters. Our 28 total parameter estimates yield a model that matches the experimental data well (Figure 3B). Details of our approach for identifying model parameters are provided in the ‘Materials and methods’ section, and our code is provided in the supplement. We obtain similar quality overall fits if we do not prespecify the RAF abundance and the RAF dimerization constant and add these two parameters to the other 28 parameters and then perform the same parameter estimation procedure (Figure 3—figure supplement 1C–E). We observed that some parameters could be constrained to a narrow region through this procedure, while other parameters could vary widely and still match the same data (Figure 3—figure supplements 2 and 3).
Intriguingly, when we considered the collection of bestfit parameters we found the different RAF inhibitors appeared to separate in a manner reflective of their known identity as type I, I.5, and II inhibitors (Figure 3C, Supplementary file 1). This is notable as only the RAF inhibitor dose–response data was examined to completely characterize the drugdependent parameters, whereas the classification of kinase inhibitors is typically made by consideration of the structure of each kinase inhibitor complex (Dar and Shokat, 2011). We also plot alternative parameter estimates that are within 10% error of the considered parameter set with the least error, and this helps demonstrate the robustness of this observation. Thus, our model finds that the different mechanisms result in qualitatively distinct features and our fit parameters suggest that these differences enable RAF inhibitors to be grouped based on respective mechanisms of action.
Conformational autoinhibition and dimer potentiation together best explain PA
We next evaluated how essential all three processes are to explain the PA data. We considered the accuracy of model fits for our full model with all three mechanisms, as well as models that only included one or two of these processes. We found that the unified model (i.e., the version that includes all three mechanisms) fits the data much better than any of these submodels (Figure 3D). Only the simpler model that included both CA and DP (but not NC) could fit the data nearly as well.
With respect to NC, the addition of this process to the model did add an incrementally better fit. However, the data did not strongly constrain the value of ‘g’ that defines the extent of NC (Figure 3D, Figure 3—figure supplements 2 and 3). Additional experiments with a larger number of drug dose concentrations may in the future allow a better characterization of NC.
Altogether, we believe that our analysis provides strong support for CA being a critical factor to PA. Our submodel analysis suggests that both CA and DP are necessary for the mathematical model to reproduce the experimental data with the smallest degree of error. Even though the NC term ‘g’ was not strongly constrained, and even though it had the smallest impact on the ability of the model to fit the experimental data, the model fits consistently had values of g that were much larger than 1 (Supplementary file 1), which is consistent with NC. Thus, our analysis suggests that all three processes contribute to the overall PA phenomenon.
Discussion
From our studies, we conclude that CA, NC, and druginduced increases in dimer affinity are required to explain the PA observed with the three major types of RAF inhibitors. Our mathematical analysis suggests that the inhibitormodulated conformational regulation of RAF kinase activity combined with conformationdependent dimerization is a critical mechanism that drives PA. Although inhibitormodulated CA has been recognized as a possible contributor to the entire PA phenomenon (Heidorn et al., 2010; Jin et al., 2017), it has been underappreciated as a mechanism capable of driving PA. Supporting this assertion is that CA has neither been discussed as a motivation for the development of thirdgeneration RAF inhibitors (Nakamura et al., 2013; Peng et al., 2015; Zhang et al., 2015) nor has it been included in recent mathematical analyses of PA and of RAF signaling (Fröhlich et al., 2022; Kholodenko, 2015; Rukhlenko et al., 2018; Varga et al., 2017). This is notable as both drug development and mathematical analyses pay considerable attention to conformations within the kinase domain (i.e., whether the alphaC helix and DFG motif are in the ‘in’ or ‘out’ conformations) (Kholodenko, 2015; Nakamura et al., 2013; Peng et al., 2015; Rukhlenko et al., 2018; Zhang et al., 2015). That suggests that the concept of structural and conformational factors on PA is not foreign to drug developers and theoretical biologists, but that they may have limited their attention to the small changes in the kinase domain that have dominated the recent literature at the expense of also considering the large, autoinhibitory, conformational changes.
Additional processes have been claimed to also play a role in PA. For example, RAF protein phosphorylation changes (Holderfield et al., 2013), preferential binding of RAF inhibitor to the different RAF proteins and/or mutant RAF (Bollag et al., 2010), scaffold proteins (McKay et al., 2011), allosteric transactivation (Hu et al., 2013), and RAS nanoclusters (Cho et al., 2012) may all further tune the response to inhibitor, including in drugspecific manners. 1433 proteins can stabilize both RAF dimers (Kondo et al., 2019; Park et al., 2019; Wilker and Yaffe, 2004) and RAF CA (Park et al., 2019; Wilker and Yaffe, 2004), and it has been hypothesized that 1433 stabilization of dimers plays roles in PA (Kondo et al., 2019; Park et al., 2019). Intuitively, stabilization of RAF CA by itself would result in a stronger PA and could be approximated as an increase in the K_{A} term within our model. RAF stabilization of RAF dimers could be approximated as in increase in the K_{dim} term within our model. Although any of these other processes could potentially modulate PA, we find that our mathematical model that includes CA, NC, and druginduced increases in dimer affinity fits the experimental data well and that those fits reproduce orthogonal biochemical features of these inhibitors. Thus, we conclude that CA, NC, and druginduced increases in dimer affinity are the three principal mechanisms underlying PA, and that the other processes listed above may provide additional modulation of the overall behavior.
It has previously been difficult to reconcile PA for type I.5 inhibitors, which are sometimes thought of as dimer breakers because they position the alphaC helix in the ‘out’ position (in contrast to type I and II inhibitors). Studies with recombinant protein and analytic ultracentrifugation clearly found type I.5 inhibitors to predominantly be in the monomeric form (Lavoie et al., 2013). Withincell assays have similarly found type I.5 inhibitors to promote dimerization less than other type I and II RAF inhibitors (Hatzivassiliou et al., 2010; Peng et al., 2015; Thevakumaran et al., 2015); however, RAF inhibitors still appeared to promote some dimerization in those incell assays. 1433 binding proteins, which can help stabilize RAF dimers, may help explain this discrepancy (Kondo et al., 2019; Liau et al., 2020; Park et al., 2019). For example, by promoting the nonautoinhibited form, a type I.5 inhibitorbound RAF monomer is more dimerization capable than an autoinhibited (and noninhibitorbound) RAF monomer. Even if the dimerization affinity is reduced compared to a nonautoinhibited and noninhibitorbound RAF monomer, 1433 proteins may be able to bind and overcome the decrease in dimerization affinity. As our model does not explicitly include 1433 proteins, this effect may contribute to our parameter estimation process finding an elevated binding affinity for type I.5bound RAF monomers.
Although NC has been difficult to precisely measure experimentally, it has widely been assumed to be present to help explain the paradoxical activation caused by type I.5 inhibitors that do not promote dimerization as strongly as other RAF inhibitors. Our bestfit parameters did tend to have g values that were larger than 1, indicating that the model fit best when there was some NC. This could suggest that NC is more abundant than widely believed. Alternatively, the model without NC was able to fit the data nearly as well as the full model that included NC (i.e., Figure 3D). This may suggest that other processes not included in the model may be modulating paradoxical activation and the g parameter, as the only other term in the model, is contributing to the model’s ability to account for these otherwise not included effects.
We found parameter sets that reproduced available, published, data in order to test our model and investigate the potential for it to help illuminate aspects of PA. The bestfit parameter sets further support a role for CA and its modulation by RAF inhibitors in PA. However, it is also important not to read too deeply into the fits. For example, the data for the type II inhibitors AZ628, LY3009120, and TAK632 had small total fold change PA magnitudes, and our fits for them have even less PA. We anticipate that the modelfitting approach would converge to increasingly accurate estimates for the parameters as the set of data being fit to expands. Additionally, quantitative experimental measurements of the parameters being fit should also cascade to impact other parameters and result in better estimates (Gutenkunst et al., 2007).
Experimental testing of model predictions is an important next step. One possible approach for experimentally investigating the role of CA on PA would involve RAF proteins that are defective in CA. Of note, there are some RAF mutants in Noonan syndrome that are believed to be activating due to impaired regulation of CA. For example, stabilization of the autoinhibited form of CRAF involves the phosphorylation of serine at residue 259 followed by 1433 binding to the phosphorylated serine. The S257L, S259F, P261A, and N262K CRAF mutations are impaired at binding to 1433 and have severely impaired phosphorylation of serine 259 (Kobayashi et al., 2010; Pandit et al., 2007). Alternatively, mutants involving the autoinhibitory domain could be used (Hartsough et al., 2018; Poulikakos et al., 2011). In such experiments, comparisons of RAF inhibitor dose–responses between cells expressing one such mutant against control cells expressing the wildtype equivalent RAF protein would reveal whether there were changes in PA fold change and PA range (as defined in Figure 2). We would expect to see reduced PA fold change and reduced PA range for the CAimpaired mutants. We anticipate that PA range would be the more useful observable because nonlinearities in downstream signaling (including possible saturation of readouts like ERK phosphorylation) may make measurements of fold change less useful. As PA range is determined by when the signal passes back through the baseline level after the increase induced by a RAF inhibitor, it should be less impacted by nonlinear signal processing downstream from RAF. Whether the experiment utilizes Noonan syndrome RAF1 mutants or truncated RAF mutants, the presence of endogenous, wildtype BRAF, ARAF, and CRAF may make it difficult to observe differences in PA due to an introduced mutant, so such experiments may need to be performed in cells where the three endogenous RAF genes are knocked out or silenced.
Our analysis was motivated by RAF inhibitors and PA in RAS mutant cells treated with a RAF inhibitor. Our model, however, is generalizable to other systems that share the modeled features. We anticipate that PA will be observed for other proteins (i) that have a dynamicequilibrium of conformations, (ii) where not all conformations can dimerize, and (iii) where drug binding the protein stabilizes one or more of the conformations that can dimerize. As dimerization and CA are both common features for kinase regulation (Huse and Kuriyan, 2002; Lavoie et al., 2014), it seems reasonable to hypothesize that PA will be observed for more kinases through modulation of the conformation and dimerization dynamic equilibrium.
Materials and methods
Mathematical models and analysis
Request a detailed protocolWe focus on steadystate levels of the different states in which RAF can exist, as portrayed in the diagrams for each model. Between any two states, an equilibrium relationship can be expressed as the ratio of abundances in the two states. Conservation of total protein quantities and zero value of total Gibbs free energy change at equilibrium both provide mechanisms to algebraically combine these expressions. We thereby derive expressions that relate the relative abundance of the RAF within its different monomeric and dimeric states. We perform algebraic manipulations and derive analytic solutions which we crosscheck using Mathematica software v 12.0 (Wolfram Research). We perform numerical evaluations of these relationships and generate plots of these equations using Python packages, including numpy, scipy, and matplotlib.
Where applicable, parameters not varied in descriptive numerical plots are set as follows.
K_{A} = 10.0, K_{d} = 0.1 μM, K_{dim} = 0.1 μM, [RAF] = 0.04 μM. We utilize mathematical models to clarify specific mechanisms that contribute to RAF activation and drug response. The models were driven by specific hypotheses that cannot be generalized to include all the features of RAF activation simultaneously due to limitations of available experimental data that does not cover all the required conditions and variations thereof. Moreover, our analytic approach allows a parameterindependent prediction for components of canonical RAF activation cycle that would not be accessible if we attempted to include every potential feature independent of a focused hypothesis. This is similar to differences in complexity between in vitro and in vivo systems – the former often forms a more versatile and precise method of establishing mechanisms while not including most of the complexities present in vivo, thereby providing a precise testing ground for focused hypotheses.
Model fitting to dose–response data
Request a detailed protocolWe reasoned that some parameters (such as the abundance of the RAF kinases, their dimerization rate constant [K_{dim}], and the equilibrium constant for RAF conformational changes [K_{A}]) are intrinsic to a cell and should not vary between different RAF inhibitors. The values of RAF concentrations and dimerization constant are also effective parameters since all three types of RAF are implicitly included when we fit to in vitro data. In contrast, some properties will vary between RAF inhibitors. Specifically, these would be the affinity of a drug for binding to RAF (K_{d}), the druginduced change in RAF dimer affinity, or “f” to use the nomenclature of Kholodenko, 2015, and NC, or “g” to again use the nomenclature of Kholodenko, 2015.
Numerical fitting for unified model of PA mechanisms was performed using SLSQP algorithm from scipy library in Python language. Sixtythree experimental dose–response points from the PDF publication (Karoulia et al., 2016) were quantified using LICOR software (Image Studio V5.2) and used as raw data to fit 28 (Figure 3, Supplementary file 1) or all 30 parameters (Figure 3—figure supplement 1C–E, Supplementary file 1) in our model for the nine drugs included. A total of 1500 random initial states were chosen from a loguniform search domain for each of the parameters. The boundaries of the fitting search were identified in the units of μM (where applicable) as follows: ‘f,’ [10^{–5},100]; ‘g,’ [1,10^{4}]; ‘K_{A},’ [0.001,10^{2}]; ‘K_{d},’ [10^{–4},10^{4}]; ‘K_{dim},’ [10^{–4},10^{4}]; ‘RAF,’ [10^{–4},10^{3}]. Best fits were defined as within 10% of the lowest value of the fitting metric. Several fitting algorithms and metrics were evaluated to finally identify SLSQP with a chisquarelike fit metric leading to convergence with higher likelihood across different initial conditions. This metric identified deviation across all data points relative to the model prediction as follows:
To compare how well each of the submodels fit the data relative to one another, we optimized over the absolute error relative to the data as this metric has a straightforward interpretation of average proportionate deviation.
In 30 parameter fits that varied all the parameters in our models, the value for RAF concentration was best fit to be 0.033 ± 0.005 μM for SKMEL2 cells (Supplementary file 1). In these fits, we also found that the dimerization equilibrium constant of RAF is monotonically correlated to autoinhibition constant K_{A} creating a nonidentifiability problem (Figure 3—figure supplement 1E). Hence, other than Figure 3—figure supplement 1C–E, the model fits were performed with a representative choice of RAF concentration fixed at 0.04 μM and RAF intrinsic equilibrium dimerization dissociation rate fixed at 0.1 μM. These parameters vary with cellular context, and the chosen values are within the range of values observed in the literature (Fujioka et al., 2006; Lavoie et al., 2013; Sadaie et al., 2014).
Mathematica and Python files that allow for the reproduction of both our analytic and numerical analyses respectively are provided as supplementary data files and are made public on GitHub link (https://github.com/GMendiratta/RAFPA copy archived at Mendiratta, 2023).
Code availability
Request a detailed protocolAll codes needed to reproduce the work presented in the article are presented in the supplementary code. The codes are also available on GitHub (https://github.com/GMendiratta/RAFPA).
Appendix 1
Supplementary theory
The theory of ordinary differential equations (ODEs) was developed alongside calculus by Newton, Leibniz, Bernoulli, Riccati, and others in the 17th century. Nonlinear differential equations are extensively used in research including to model biological systems of neural signaling, protein interactions, and cardiac potentials; however, these are rife with singularities and numerical instabilities. Through this long history, differential equations that can be fully solved and are applicable to realworld conditions are quite rare. Hence, scientifically relevant, analytically accessible solutions are of considerable value to both progress of applied mathematics and quantitative sciences. In this work, we show such solutions under equilibrium conditions for a set of ODE models that describe the dimerization and drug binding of RAF protooncoproteins. Our work includes a series of mathematical models involving RAF dimerization, drug binding, RAF CA, and the interactions that affect these component processes. Each of these models is based on evidence from multiple sources in the literature that support the modeled RAF activation cycle and drug interactions. The value of such analytically solved models is both in their generalizability and in their global predictions. For example, in the first section, we will show that it is not possible to generate PA in a system where RAF monomer can dimerize and where drug binds equivalently well to monomeric and dimeric RAF proteins. This statement is derived independent of any parameter values that may be used – this is an ‘analytic’ result. Other examples in biological literature of similarly powerful analytic studies include Michaelis–Menten enzyme kinetics (Michaelis and Menten, 1913), negative feedbackbased oscillatory circuits (Glass and Mackey, 1988; Lewis, 2003), and the Hodgkin–Huxley model (Hodgkin and Huxley, 1952; Noble, 1962).
The molecular interaction rates for RAF activation and drug binding all occur on a time scale of protein–protein and protein–ligand interactions (Fujioka et al., 2006; Gillies et al., 2020), which are orders of magnitude smaller than the time scales of paradoxical activation and proliferation phenotypes associated with MAPK pathway activation (Hatzivassiliou et al., 2010; Heidorn et al., 2010) or the time scale of clinical impact (Gibney et al., 2013). We therefore focus on steadystate levels of the different states in which RAF can exist within our models. At equilibrium, we derive expressions for the relative abundance of the RAF within its different monomeric and dimeric states. We perform algebraic manipulations and derive analytic solutions that we reproduce using Mathematica software. We solve the equilibrium state numerically where it is instructive as visual aid to the analytical results or where analytic results cannot be obtained globally.
Within Appendix 1, we first describe a model of RAF dimerization and drug binding, followed by a model of autoinhibition without dimerization, and then a model that combines these two components. A numerical illustration of the mathematical model used that displays the drug dosedependent response of the different states of RAF is included as Figure 2—figure supplement 2, with narrative text in the main article under the section ‘ Paradoxical activation reflects a shifting balance of signaling complexes.’ The NC and DP mechanisms utilized in unified mechanisms model have been discussed in detail in the literature (Kholodenko, 2015). We show the resulting expressions in the unified model in Supplementary file 1 (and use these in Figure 3 plots and fits) and the derivations thereof are included in Mathematica code present in the included supplementary code files. We provide a numerical illustration of the unified model that displays the drug dosedependent response of the different states of RAF included as the last section of Appendix 1.
The promotion of RAF dimerization by a drug does not induce PA by itself
We begin by constructing a model of kinase activation to dispel a possible misunderstanding of the behavior possible through RAF dimerization and drug binding alone. In the mechanism that we model (Appendix 1—figure 1), a kinase exists in the competent/nonautoinhibited form (denoted by capital A) where it may dimerize (AA) and bind with the drug in the monomeric (Ad) or dimer form (AAd,AdAd). This model resembles a similar, older model of EGFR activation and dimerization (Wofsy et al., 1992). We can assign an equilibrium constant to each of the reactions that is equal to the ratio of reactants to the products:
Therefore, the process of drug binding has an equilibrium rate constant ${K}_{d}$:
The process of nonautoinhibited kinase dimerization corresponds to the equilibrium rate constant ${K}_{dim}$:
We derive closedform, analytic steadystate solutions relating the above rate constants with protein and drug concentrations using the principle of detailed balance (Kholodenko, 2015; Wofsy et al., 1992). At equilibrium, the principle of detailed balance or microscopic reversibility ensures that all cyclical processes have a zero change in Gibbs free energy. Free energy with a constant temperature, pressure, and pH is proportional to the logarithm of the equilibrium rate constant. Hence, corresponding to the network described in Appendix 1—figure 1, we obtain the following relations:
and
Since there are two binding locations for the drug in a dimer ($AA+d>AAd$), the reaction is twice as likely, and the dissociation constant is halved ${K}_{d}/2$. From the detailed balance equation, we calculate that ${K}_{Ad+A>AAd}={K}_{dim}/2$, hence connecting drug binding to the dimerization process via dynamical equilibrium relationships. In this first, simple, model, we assume that the drug does not induce any changes in the kinase conformation and therefore the dimerization of drugbound kinase proceeds at the same rate as unbound kinase, $K}_{Ad+Ad>AdAd}={K}_{dim$ . We can now calculate the last unknown rate, $K}_{AAd+d>AdAd}=2{K}_{d$. Thereby, only two rate constants determine the equilibrium state of this system.
To calculate the steadystate or equilibrium relationships, we begin with reactions that lead from a RAF protomer to a fully drugbound RAF dimer:
The following expressions correspond to the rate constants for the above reactions:
The total amount of drug and kinase is obtained by combining all states that are assumed unchanged over the short reaction times involved in the binding and unbinding molecular reactions, and are defined as below:
Replacing the equilibrium concentrations into the expression for total kinase concentration, we first solve for the total number of unbound RAF protomers:
Here, $RA{F}_{rel}=RAF/{K}_{dim}$ and ${d}_{rel}=\left[d\right]/{K}_{d}$. The total quantity of active protomers is a sum of the protomers in the kinase dimers and the dimers bound to one drug, ($2\left[AA\right]+\left[AAd\right]$). Substituting the equilibrium concentrations, the concentration of active kinase at equilibrium is obtained:
Utilizing the previously calculated equilibrium concentration of nonautoinhibited monomer kinase, we can reduce the above expression to a function of the dimensionless ratios $RA{F}_{rel}$ and ${d}_{rel}$:
It is clear that the above expression for active kinase protomers is a monotonically reducing function of the unbound drug concentration since that variable only occurs in the denominator. The second component required to analytically establish the absence of paradoxical activation in this model is to prove that the unbound drug is monotonically and proportionally related to total drug concentration. We do this by substituting the equilibrium relations, which provides the total drug amount:
Next, we substitute the equilibrium kinase monomer concentration [A] in the above expression to solve for equilibrium unbound drug concentration [d].
Here, ${\zeta}_{0}=DrugRAF{K}_{d}$. Therefore, the relationship between total and unbound drug concentrations is monotonic and directly related.
We therefore show that active RAF is an inverse function of drug concentration. PA is not possible while considering only simple dimerization and drug binding processes of RAF protein. This point may not have been sufficiently stressed in the literature, leading to confusion about the added dimers are induced by the drug causing PA by themselves. These new, partly drugbound, dimers will never create additional RAF signaling protomers unless there was more to the RAF activation processes than dimerization and drug binding. For example, DP, where the affinity for RAF dimerization is higher if one (or both) RAF protomers is bound to a RAF inhibitor, has previously been demonstrated as an addition to the above modeled mechanism that could result in a PA response to RAF inhibitor (Kholodenko, 2015).
Autoinhibited proteins active as monomers will not display PA from druginduced stabilization of the nonautoinhibited form alone
Extensive research has searched for explanations of PA that go beyond DP. These mechanisms include drug NC, RAS oligomerization, downstream negative feedback release, and more. One of these mechanisms, based on RAF CA, is particularly interesting because it does not involve any new, druginduced, structural conformations of the RAF protein. With a simple extension to our modeled network, we identify the conditions under which this phenomenon successfully provides a mechanism for PA in response to inhibition of dimerizing kinases.
We first model drugmediated stabilization of the active form of a conformation autoinhibitionregulated kinase for a kinase that does not dimerize and is fully active as a monomer. To study this, we create a model where the kinase may autoinhibit (denoted by ‘a’) to become reversibly inactive and where a drug will only bind the nonautoinhibited form (Appendix 1—figure 2). The following equilibria contribute to this network:
with corresponding equilibrium dissociation constants defined as follows:
The total drug and kinase amounts are calculated as below:
We solve the above four equations for the total nonautoinhibited protomers ([A]) as a function of total concentrations and equilibrium dissociation constants:
where ${C}_{1}=RAFDrug{K}_{d}\left(1+{K}_{A}\right)$. The function monotonically reduces to zero as the total drug increases to a large value. Therefore, a kinase that can signal downstream as a monomer does not show anomalous activation upon drug inhibition even though there is a reservoir available to shift the equilibrium to activated kinase.
The combination of conformational autoinhibition and stabilization of the active form by a drug can generate PA
Conformational changes of the RAF monomer contribute significantly to the regulation of RAF kinase activity (Lavoie and Therrien, 2015; Wellbrock et al., 2004). In the ‘autoinhibited’ form, associations between its Nterminus and its kinase domain maintain RAF in an inactive form that does not dimerize (Cutler et al., 1998; Lavoie and Therrien, 2015). In the ‘nonautoinhibited’ form, the kinase domain is no longer occluded, and other regulatory mechanisms contribute to full RAF kinase activation, such as kinase domain conformational changes and dimerization (Lavoie et al., 2013). Recent experimental work reports that RAF inhibitors tend to promote a net transition to the nonautoinhibited conformation that binds to RASGTP (Jin et al., 2017; Karoulia et al., 2017). It has previously been suggested that this biasing to the nonautoinhibited state may contribute to PA (Jin et al., 2017).
To study the conditions under which the stabilization of RAF in its nonautoinhibited state by RAF inhibitors may be sufficient to generate PA, we create a mathematical model that includes RAF CA, RAF dimerization, and RAF inhibitor (Appendix 1—figure 3). The model allows RAF to adopt two different conformations: one is autoinhibited and can neither dimerize nor bind drug, and the other is nonautoinhibited and can bind drug and/or dimerize (Lavoie et al., 2013). Drugbound RAF is assumed to only be able to transition back to an autoinhibited state after any bound drug has dissociated. Within the model, wildtype RAF is implicitly assumed to be activated by RASGTP as binding to RASGTP is an essential step to wildtype RAF activation (Lavoie and Therrien, 2015). We define active RAF as the RAF protomers that are not bound to a drug and are part of a RAF dimer. The following expressions identify the different interactions included in this network:
The processes in this model include autoinhibition, which, for simplicity, is modeled as a firstorder transition from autoinhibited state $a$ to nonautoinhibited state $A$ with a rate constant ${K}_{A}$:
The equilibrium equations defined in the previous section apply in the base model with the addition of the CA dynamical equilibrium represented in above equation. Over the reaction periods, we assume that the total RAF and drug concentrations remain unchanged, leading to the conservation equations below:
Following the analysis described in the first section, we solve the RAF conservation equation to find the concentration of unbound nonautoinhibited RAF:
Note that there is no explicit dependence on the drug binding rate ${K}_{d}$ and only on the ratio of the free drug concentration to ${K}_{d}$ appears explicitly (d_{rel} = [d]/K_{d}). Similarly, dimerization rate constant ${K}_{dim}$ can be combined with the total amount of RAF protein $RA{F}_{rel}=\left[RAF\right]/{K}_{dim}$ . This reduces the model from four varying parameters to two effective parameters, which govern all predictions. The total active kinase is the sum of drugunbound RAF protomers in a drugunbound RAF dimer and partly drugbound RAF dimer.
Substituting the above equation and the equilibrium conditions, we obtain the active RAF proportion shown in the first column of Supplementary file 1. This function is not monotonic, and the behavior can be defined by evaluating the zeroes of the derivatives of the function. The first derivative of the function representing active RAF protein simplifies to the following expression:
Here, $E1=1+{K}_{A}+{d}_{rel}$ and $E2=8{\left(1+{d}_{rel}\right)}^{2}RA{F}_{rel}$ . Maximal paradoxical activation, if it exists, occurs when the first derivative of active RAF relative to the drug concentration is zero. We calculate the values of ${d}_{rel}$, which satisfy this condition and find that only one solution can potentially be real and positive:
To additionally establish that this point of inflection is a maximum, we need to establish conditions under which the second derivative of the active RAF equation relative to ${d}_{rel}$ is negative. At the drug concentration in the above equation, the second derivative of active RAF simplifies to the following expression:
This function is negative definite as long as the parameters are all positive. The latter is a given since the parameters are either concentrations (nonnegative) or ratios of concentrations (equilibrium constants). The drug will show paradoxical activation via the suggested autoinhibition mechanism when the value of ${d}_{rel}$ is positive. We find that the maximum fold change occurs at a positive drug concentration only when the following condition is true to derive the necessary and sufficient condition for PA within this submodel:
As the concentration ratio $RA{F}_{rel}$ may not be smaller than 0, the righthand expression is always positive. Therefore, the condition that the autoinhibition equilibrium should favor the inactive RAF kinase naturally arises from the above expression (i.e., ${K}_{A}>1$). In other words, PA does not always happen when CA is present. It is also clear from the same expression that a model that does not show autoinhibition (${K}_{A}=0$) or a model where there are no dimers (${K}_{dim}\to \mathrm{\infty}:RA{F}_{rel}\to 0$) will not produce PA. Therefore, we demonstrate analytically that the presence of both CA and stabilization of the active form by RAF inhibitors is sufficient to create PA for a wide range of system parameters.
Next, we identify the predictions of this model for the baseline RAF signaling, defined as active RAF in the absence of the drug. The expression for baseline signaling in the basic model is shown in Supplementary file 1. This expression is positive definite and monotonic as a function of E3, where $E3=8\times RA{F}_{rel}/{\left(1+{K}_{A}\right)}^{2}$. This can be observed by taking the first derivative relative to E3, which reduces to the following positive definite function:
Therefore, baseline signaling monotonically increases as a function of RAF concentration and reduces as a function of equilibrium constant ${K}_{A}$ and equilibrium dissociation constant ${K}_{dim}$.
Moreover, expression for the maximum paradoxical activation can be calculated by substituting the corresponding value of drug concentration into the expression for active kinase concentration:
Note that in this mechanism the maximal kinase activity is a simple inverse function of ${K}_{A}$ . Since the minimum value of ${K}_{A}$ to produce PA limits to 1 from above, the maximum possible PA for a kinase is a function of $RAF/{K}_{dim}$ only. This is also because we have modeled the autoinhibitory mechanism in the simplest possible, firstorder approximation that represents the phenomena itself but may not encompass the detailed dynamics such as that modeled in the following sections. The kinase activity in the absence of the drug is shown in the following expression:
The maximal fold change (FC) in activity can therefore be calculated:
When the conditions of PA are satisfied, the fold change monotonically increases as a function of the autoinhibition equilibrium constant and reduces as a function of the RAF concentration to dimerization dissociation constant ratio. Both of these facts can be derived upon evaluating the first derivative of the fold change function above and applying the constraint relationship for PA existence, which results in a positive definite derivative relative to ${K}_{A}$ and a negative definite derivative relative to $RA{F}_{rel}$. The implication of the former is that when the autoinhibition equilibrium favors the inactive state, the fold change relative to baseline (drugfree) condition is enhanced. Further, when the dimerization rate of RAF increases (i.e., dissociation constant ${K}_{dim}$ reduces), the corresponding increase in $RA{F}_{rel}$ leads to a reduced fold change in RAF activity in response to the drug, relative to the control.
The relationship between unbound drug and total drug
We note here that the expressions characterizing the PA phenomenon thus far are done as a function of the freedrug concentration. This is to keep the expressions relatively simple. The total drug available in a cell also includes the drug concentration bound to RAF, including both monomeric and dimeric states of RAF. Our conclusions become final predictions of the model once we can show that unbound/free drug and the total drug concentrations are directly and monotonically related.
The relationship between total drug available (DTOT) and the unbound (free) drug (d) can be calculated by substituting the equilibrium conditions:
where the unbound kinase (A) is a function of free drug (d). After splitting this function into three additive terms, we substitute the expression for [A] and show that the derivative relative to unbound drug for each term is nonnegative. We do not show the explicit expression for these derivatives here and refer to the Mathematica notebook provided as Supplementary file 1. We note that these derivatives are positive definite for all positive values of the parameters in the model. Therefore, we obtain that the total drug quantity is a monotonic increasing function of the quantity of unbound drug (Appendix 1—figure 4). The monotonicity of this function implies that our conclusions for paradoxical activation follow through as a function of total drug as well. However, the nonlinear relationship implies that the map between drug IC50 values and the freedrug amount is directly related but also nonlinear. This is significant while comparing IC50 values predicted by a modeling framework with experimental results. When fitting our RAF activity functions to experimental data, we work with the total drug concentration.
Numerical model illustrations of three mechanisms of PA: NC, DP, and CA
Appendix 1—figure 5 shows the mechanisms of drug potentiation of RAF dimerization (DP, blue line), NC (red dotted line), CA (green dashed line), and the unified model with all three mechanisms (unified, black dashed line). We illustrate each of the states of RAF (i–vi) for each of the submodels, in a manner similar to what is shown in Figure 2—figure supplement 2, but now extended to include these additional models. The CA and DP mechanisms show distinct features in active RAF states. The DP model shows a rapid activation of RAF in response to a dimerizationinducing drug (blue, solid line) and then a rapid reduction of signaling. In contrast, the PA peak occurs at higher drug concentrations for the CA model, and the CA model also has a much wider tail of the distribution (green, dashed line). The NC model (red, dotted curves) also shows some increase of the partly active RAF dimer at low drug concentrations. However, this mechanism alone will not generate PA as the total abundance of active RAF protomers within RAF dimers will always be less than the levels obtained in the absence of drug. This point can also be shown as a parameterindependent, analytic feature of the NC model (Supplementary file 1) and is in agreement with prior in silico modeling of the NC mechanism (Kholodenko, 2015). We show in Figure 3—figure supplement 1A and B that in the presence of other mechanisms that induce PA, the NC model can modulate the dose–response curve to activate RAF at higher concentrations of the drug. The nature of their contributions remains quite distinct across the DP, NC, and CA models even as they combine to create more complex dose–response.
All of the mechanisms are combined into a unified model (black, dashed), which incorporates the unique features of each of the submodels. The curves in Appendix 1—figure 5 identify the differences within the models and also show the synergistic combination of the mechanisms. We also note the analytic expressions for each of the models in Supplementary file 1 are linearly independent, and therefore cannot be converted into each other by any simple translation or scaling transformation of any of the variables. The hypothesis that we could completely characterize these parameters using normalized data from signaling experiments arose from this insight into the analytic and numerical activation curves. These fits may be performed in the future for other cell lines, primary cells, mouse models, or any system of interest where such characterization may aid in drug development and precision therapy identification. The clear distinguishability of different mechanisms using fitting is an emergent feature due to the nonlinearities between how the models combine. The distinguishability is a rare feature in multiparameter models because of degeneracies that are common in multiparameter models. An example of a degeneracy in our model system involves the abundance of RAF and the dimerization equilibrium constant; these terms are difficult to isolate from each other in our fits. The coarsegrained data signaling data in our fits is much improved with narrower and more robust parameter regions predicted when the RAF concentration and dimerization rates are fixed to values based on the literature.
Finally, while NC and druginduced increases in dimerization are behaviors that can be tuned between RAF inhibitors, we highlight that CA is a feature that spans all RAF inhibitors. This means that this process and the secondary events that potentiate this process are important factors to consider in the development of any RAF inhibitor.
Data availability
All data needed to evaluate the conclusions in the paper are present in the paper or the Supplementary Materials. All materials are available upon request from the corresponding author.
References

Raf inhibitors target ras spatiotemporal dynamicsCurrent Biology 22:945–955.https://doi.org/10.1016/j.cub.2012.03.067

EGFR antagonists in cancer treatmentThe New England Journal of Medicine 358:1160–1174.https://doi.org/10.1056/NEJMra0707704

The evolution of protein kinase inhibitors from antagonists to agonists of cellular signalingAnnual Review of Biochemistry 80:769–795.https://doi.org/10.1146/annurevbiochem090308173656

Dynamics of the Ras/ERK MAPK cascade as monitored by fluorescent probesThe Journal of Biological Chemistry 281:8917–8926.https://doi.org/10.1074/jbc.M509344200

Paradoxical oncogenesis—the longterm effects of BRAF inhibition in melanomaNature Reviews Clinical Oncology 10:390–399.https://doi.org/10.1038/nrclinonc.2013.83

Oncogenic mutant RAS signaling activity is rescaled by the ERK/MAPK pathwayMolecular Systems Biology 16:e9518.https://doi.org/10.15252/msb.20209518

BookFrom Clocks to Chaos: The Rhythms of Life.Princeton University Press.https://doi.org/10.1515/9780691221793

Universally sloppy parameter sensitivities in systems biology modelsPLOS Computational Biology 3:0030189.https://doi.org/10.1371/journal.pcbi.0030189

Response and resistance to ParadoxBreaking BRAF inhibitor in Melanomas In Vivo and Ex VivoMolecular Cancer Therapeutics 17:84–95.https://doi.org/10.1158/15357163.MCT170705

A quantitative description of membrane current and its application to conduction and excitation in nerveThe Journal of Physiology 117:500–544.https://doi.org/10.1113/jphysiol.1952.sp004764

New perspectives for targeting RAF kinase in human cancerNature Reviews. Cancer 17:676–691.https://doi.org/10.1038/nrc.2017.79

Inhibitors that stabilize a closed RAF kinase domain conformation induce dimerizationNature Chemical Biology 9:428–436.https://doi.org/10.1038/nchembio.1257

Dimerizationinduced allostery in protein kinase regulationTrends in Biochemical Sciences 39:475–486.https://doi.org/10.1016/j.tibs.2014.08.004

Regulation of RAF protein kinases in ERK signallingNature Reviews. Molecular Cell Biology 16:281–298.https://doi.org/10.1038/nrm3979

Negative regulation of RAF kinase activity by ATP is overcome by 1433induced dimerizationNature Structural & Molecular Biology 27:134–141.https://doi.org/10.1038/s4159401903650

Discernment between candidate mechanisms for KRAS G13D colorectal cancer sensitivity to EGFR inhibitorsCell Communication and Signaling 18:179.https://doi.org/10.1186/s12964020006453

Sotorasib for Lung Cancers with KRAS p G12C MutationThe New England Journal of Medicine 384:2371–2381.https://doi.org/10.1056/NEJMoa2103695

RAS mutations in cutaneous squamouscell carcinomas in patients treated with BRAF inhibitorsThe New England Journal of Medicine 366:207–215.https://doi.org/10.1056/NEJMoa1105358

Crystal structure of a BRAF kinase domain monomer explains basis for allosteric regulationNature Structural & Molecular Biology 22:37–43.https://doi.org/10.1038/nsmb.2924

RAF1/BRAF dimerization integrates the signal from RAS to ERK and ROKαScience Signaling 10:eaai8482.https://doi.org/10.1126/scisignal.aai8482

The RAF proteins take centre stageNature Reviews. Molecular Cell Biology 5:875–885.https://doi.org/10.1038/nrm1498

1433 Proteinsa focus on cancer and human diseaseJournal of Molecular and Cellular Cardiology 37:633–642.https://doi.org/10.1016/j.yjmcc.2004.04.015

The AMPK inhibitor overcomes the paradoxical effect of RAF inhibitors through blocking phosphoSer621 in the C terminus of CRAFThe Journal of Biological Chemistry 293:14276–14284.https://doi.org/10.1074/jbc.RA118.004597
Article and author information
Author details
Funding
National Institutes of Health (K22CA216318)
 Edward Stites
National Institutes of Health (DP2AT011327)
 Edward Stites
Melanoma Research Alliance (Young Investigator Award)
 Edward Stites
Joe W. and Dorothy Dorsett Brown Foundation
 Edward Stites
Salk Institute for Biological Studies
 Gaurav Mendiratta
Conrad Prebys Foundation
 Edward Stites
The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.
Acknowledgements
We thank Tony Hunter, Geoff Wahl, Reuben Shaw, Bjorn Lillemeier, Dmitry Lyumkis, Joe Noel, Andrey Shaw, Rajasree Kalagiri, Leo Li, Michael Trogdon, Noah DeTal, Peter Carlip, and Melinda Tong for helpful conversations and feedback. We thank Amy Cao for help with illustrative graphics. This work was supported by a Salk Pioneer Fund Postdoctoral Scholar Award (GM); the National Institutes of Health NIH K22CA216318 (ECS), NIH DP2AT011327 (ECS), a Melanoma Research Alliance Young Investigator Award (ECS), the Joe W and Dorothy Dorsett Brown Foundation (ECS), and the Conrad Prebys Foundation (ECS).
Version history
 Preprint posted: November 20, 2019 (view preprint)
 Received: August 16, 2022
 Accepted: October 10, 2023
 Accepted Manuscript published: October 12, 2023 (version 1)
 Version of Record published: November 6, 2023 (version 2)
Copyright
© 2023, Mendiratta and Stites
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.
Metrics

 733
 views

 133
 downloads

 2
 citations
Views, downloads and citations are aggregated across all versions of this paper published by eLife.
Download links
Downloads (link to download the article as PDF)
Open citations (links to open the citations from this article in various online reference manager services)
Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)
Further reading

 Biochemistry and Chemical Biology
 Cell Biology
Numerous lipids are heterogeneously distributed among organelles. Most lipid trafficking between organelles is achieved by a group of lipid transfer proteins (LTPs) that carry lipids using their hydrophobic cavities. The human genome encodes many intracellular LTPs responsible for lipid trafficking and the function of many LTPs in defining cellular lipid levels and distributions is unclear. Here, we created a gene knockout library targeting 90 intracellular LTPs and performed wholecell lipidomics analysis. This analysis confirmed known lipid disturbances and identified new ones caused by the loss of LTPs. Among these, we found major sphingolipid imbalances in ORP9 and ORP11 knockout cells, two proteins of previously unknown function in sphingolipid metabolism. ORP9 and ORP11 form a heterodimer to localize at the ERtransGolgi membrane contact sites, where the dimer exchanges phosphatidylserine (PS) for phosphatidylinositol4phosphate (PI(4)P) between the two organelles. Consequently, loss of either protein causes phospholipid imbalances in the Golgi apparatus that result in lowered sphingomyelin synthesis at this organelle. Overall, our LTP knockout library toolbox identifies various proteins in control of cellular lipid levels, including the ORP9ORP11 heterodimer, which exchanges PS and PI(4)P at the ERGolgi membrane contact site as a critical step in sphingomyelin synthesis in the Golgi apparatus.

 Biochemistry and Chemical Biology
 Structural Biology and Molecular Biophysics
The mechanism underlying the preferential and cooperative binding of cofilin and the expansion of clusters toward the pointedend side of actin filaments remains poorly understood. To address this, we conducted a principal component analysis based on available filamentous actin (Factin) and Cactin (cofilins were excluded from cofilactin) structures and compared to monomeric Gactin. The results strongly suggest that Cactin, rather than FADPactin, represented the favourable structure for binding preference of cofilin. Highspeed atomic force microscopy explored that the shortened bare half helix adjacent to the cofilin clusters on the pointed end side included fewer actin protomers than normal helices. The mean axial distance (MAD) between two adjacent actin protomers along the same longpitch strand within shortened bare half helices was longer (5.0–6.3 nm) than the MAD within typical helices (4.3–5.6 nm). The inhibition of torsional motion during helical twisting, achieved through stronger attachment to the lipid membrane, led to more pronounced inhibition of cofilin binding and cluster formation than the presence of inorganic phosphate (Pi) in solution. FADPactin exhibited more naturally supertwisted half helices than FADP.Piactin, explaining how Pi inhibits cofilin binding to Factin with variable helical twists. We propose that protomers within the shorter bare helical twists, either influenced by thermal fluctuation or induced allosterically by cofilin clusters, exhibit characteristics of Cactinlike structures with an elongated MAD, leading to preferential and cooperative binding of cofilin.