Quantitative Systems Pharmacology Analysis of KRAS G12C Covalent Inhibitors.

KRAS has proven difficult to target pharmacologically. Two strategies have recently been described for covalently targeting the most common KRAS mutant in lung cancer, KRAS G12C. Previously, we developed a computational model of the processes that regulate Ras activation. Here, we use this model to investigate KRAS G12C covalent inhibitors. We updated the model to include Ras protein turnover, and validation demonstrates that our model performs well in areas of G12C targeting where conventional wisdom struggles. We then used the model to investigate possible strategies to improve KRAS G12C inhibitors and identified GEF loading as a mechanism that could improve efficacy. Our simulations also found resistance-promoting mutations may reverse which class of KRAS G12C inhibitor inhibits the system better, suggesting that there may be value to pursuing both types of KRAS G12C inhibitors. Overall, this work demonstrates areas in which systems biology approaches can inform Ras drug development.

WHAT IS THE CURRENT KNOWLEDGE ON THE TOPIC?

☑ The individual mechanisms of Ras signal regulation and of covalent KRAS G12C inhibitors are known. How the combined actions of the different processes that regulate Ras signals influence the efficacy of covalent inhibitors is neither well known nor easily inferred.

WHAT QUESTION DID THIS STUDY ADDRESS?

☑ Can mathematical modeling account for the observed, nonintuitive, behaviors of Ras inhibiors? How can covalent inhibitor strategies be improved?

WHAT DOES THIS STUDY ADD TO OUR KNOWLEDGE?

☑ This study finds that observed, nonintuitive, behaviors are fully consistent with known Ras biology and could have been predicted a priori. The study also highlights the role of Ras GEFs in oncogenic Ras regulation and targeting.

HOW MIGHT THIS CHANGE DRUG DISCOVERY, DEVELOPMENT, AND/OR THERAPEUTICS?

☑ The identification of properties that could make NPIs efficacious suggests new assays for screening drug candidates. The identification that which of the two strategies is better can reverse upon resistance promoting mutations suggests benefit to continued pursuit of both strategies.

Human cancer cells, such as those from pancreatic cancer, colorectal cancer, and lung cancer, commonly include somatically acquired KRAS mutations.1, 2, 3 The Ras proteins ( KRAS, NRAS, and HRAS) bind to guanine nucleotides GDP and guanosine‐triphosphate (GTP) with high affinity.4 The Ras proteins can hydrolyze GTP to GDP, which is why they are commonly referred to as “GTPases”. The GTP‐bound form of Ras is considered the active form, and downstream signaling effectors specifically bind to the GTP‐bound form of Ras.5 The cancer‐promoting KRAS mutations most commonly occur at codon 12, 13, or 61, and result in increased levels of GTP‐bound KRAS, which in turn promotes downstream signaling.6 It has long been believed that drugs with the ability to block aberrant KRAS signaling would benefit patients with cancer.4

A major advance in the targeting of KRAS has been made with the development of molecules that covalently and irreversibly bind to the cysteine residue of the KRAS G12C mutant.7, 8, 9, 10 The KRAS G12C mutant is particularly common in lung cancers. Selectivity of a drug for the cysteine in this mutant protein should result in specific targeting of KRAS G12C mutant‐containing cancer cells. Compounds that covalently interact with the codon 12 cysteine can be grouped into two different classes based upon which pocket the molecule settles into after bonding with cysteine. The compounds described by the Gray laboratory rest in the nucleotide pocket (NP) in place of a guanine nucleotide.8 The G12C inhibitors that reside in the NP will be referred to as NP‐inhibitors (or NPIs). The compounds initially described by the Shokat laboratory rest in a distinct pocket, the “switch two pocket” that was previously unappreciated.7 Compounds using this strategy will be referred to as SIIP‐inhibitors (or SIIPIs). The originally described NPIs and SIIPIs were not able to reliably target KRAS G12C in cell culture. A new SIIPI, ARS‐853, has been described that has much improved binding properties and displays KRAS G12C inhibition in cell culture.11

Previously, we developed a computational model of the processes that together regulate Ras signaling.12, 13 The model has made multiple predictions about Ras signaling that went against conventional wisdom but have now been experimentally confirmed. For example, our model predicted that oncogenic Ras leads to wild‐type Ras activation,12 which has since been experimentally demonstrated by multiple laboratories.14, 15 More recently, we predicted that mutations to Ras GAP NF1 should frequently co‐occur with noncanonical cancer driver genes in the Ras pathway,16 which was subsequently detected in melanoma genomes.17 Here, we applied the model to the problem of targeting KRAS G12C. To do this, the model was updated to include protein turnover, the specific biochemical properties of the G12C mutant, and the described mechanisms of interaction for the two classes of inhibitors. Simulations of the updated model yield patterns of Ras activation that match well with previous experimental observations. Application of the model to NPIs finds that GEF loading is a property that could be optimized to make NPIs more effective. The model also suggests that mutations that promote resistance to SIIPIs may increase sensitivity to NPIs. Overall, this work demonstrates the potential of mechanistically accurate models of oncogenic signal regulation to contribute to quantitative systems pharmacology.

METHODS

The details of the Ras model have been previously published and described in detail.12, 13, 18, 19, 20 The model is specified as a set of coupled, nonlinear, ordinary differential equations. Simulations in MATLAB were used to find steady‐state levels of RasGTP and RasGTP‐Effector complex, both of which were considered as measures of Ras pathway signaling.12 Mutant Ras proteins are specified in terms of how each parameter describing a Ras biochemical reaction differs from the wild‐type Ras value. We specify these as multiplicative “scaling factors” that are equal to the ratio between the mutant parameter and the wild‐type parameter found within the same experimental study. Table 1 provides parameters unique to this model, or that are particularly relevant to understanding the present study.

Table 1

Parameters for the G12C model

Parameter	Value	Units	Description
kdeg	8 × 10−6	/s	Ras degradation rate
kproduction	3.2 × 10−12	/Ms	Total Ras production rate
kon,NPI	2 × 106	/Ms	On rate for NPI
kon,SIIPI,12	0.12	/Ms	On rate for SIIPI (compound 12)
kon,SIIPI,ARS	76	/Ms	On rate for SIIPI (ARS‐853)
Percent_KRAS	50	%	Percentage of total Ras that is KRAS
Percent_Mut	50	%	Percentage of KRAS that is mutant
ScalingfactorGTPase,G12C	72	%	Intrinsic GTPase rate of G12C relative to that of wild‐type Ras
ScalingfactorkdissG12C	120	%	Rate of effector dissociation of G12C relative to that of wild‐type Ras
GTP	1.8 × 10−4	M	Cellular abundance of GTP
GDP	1.8 × 10−5	M	Cellular abundance of GDP

GTPase, guanosine triphosphate; NPI, nucleotide pocket inhibitor; SIIPI, switch II pocket inhibitor.

Parameters for the expanded the Ras model that includes protein production and degradation, NPIs, SIIPIs, and the G12C mutant.

RESULTS

Model development and analysis

Model extension to include protein turnover and the G12C mutant

To study the G12C mutant targeted by a covalent inhibitor, the model was updated in two key ways. First, protein synthesis and degradation were added to the model. When the model was originally developed, we assumed that these processes would be slow relative to the processes that regulate steady‐state Ras nucleotide binding and could, therefore, be ignored when studying oncogenic mutant activation. Because covalent inhibitors target G12C irreversibly (and could potentially increase degradation), we added a degradative process to the model. A schematic with the updates to the original model is shown in Figure 1 a. Briefly, the model now allows for degradation of Ras in any form, and production of Ras results in nucleotide‐free Ras, which will rapidly bind to available cytoplasmic nucleotides upon modeled production. Protein degradation was modeled as a first order process with a half‐life of ∼24 hours (kdeg = 8 × 10−6/s).21, 22 Production was modeled to occur at a constant rate. This rate was set so the steady‐state level of total Ras remains at the level of our original estimate for total cellular Ras. Values for these parameters are provided in Table 1.

Figure 1

Modeled biochemical reactions of the Ras network, including production and degradation. (a) The original Ras model was based upon a dynamic equilibrium of Ras signaling states that followed from the biochemical reactions that influence Ras nucleotide binding state. This model was updated to include Ras protein production (green arrow) and degradation (red arrows). Parameters for these reactions are presented in Table 1. (b) Predicted Ras signaling levels for Ras mutant proteins. The updated model diagrammed in a was compared to the original model for predictions of Ras pathway activation for the networks with all Ras in the wild‐type form, or with both wild‐type Ras and a mutant form of Ras (G12V, G1D, G12C, or F28L). (c) Predicted Ras signaling levels for different proportions of basal GAP activity present. Reduced values could come for reasons such as deletion of tumor suppressor neurofibromin (NF1). b and c present computational predictions of the proportion of total Ras bound to GTP (left) and the proportion of total effector bound to RasGTP (right). Predictions with the original model without turnover are in blue, those with the model updated to include turnover are in maroon. GTPase, guanosine triphosphate.

We also updated the model to include the G12C mutant. Of note, oncogenic mutants at codons 12, 13, and 61 are biochemically similar in that they are insensitive to inactivation by Ras GAPs.6 However, there are subtle differences with other parameters between different oncogenic mutants that, in some cases, can have large effects.16, 23 For example, patients with colon cancer with the oncogenic G13D KRAS mutation seem to benefit from cetuximab, whereas patients with colon cancer with the oncogenic G12D or G12V KRAS mutations seem to not benefit from cetuximab.24 Therefore, we conclude that it is most prudent to model the G12C KRAS mutant by incorporating recently published data.25 Specifically, the intrinsic GTPase rate and the effector binding constant are two parameters of the G12C mutant that vary from wild type.25 The values used in our model are presented in Table 1.

Tests of the extended model with turnover and the G12C mutant

We first wanted to compare the behavior of the model extended to include protein turnover with the previous version of the model that did not include turnover. For test cases, we considered cells with all wild‐type Ras, that contain the G12C Ras mutant, and that include the G12D, G12V, or F28L Ras mutants studied in our original Ras model.12 We first evaluated the behavior of the updated model by examining the predicted level of RasGTP and RasGTP‐effector complex for these different cases. For these conditions, we modeled 25% of total RAS as mutant and 75% as wild type, such as might occur if 50% of the RAS protein in the cell was KRAS, the other 50% NRAS and HRAS, and one of the KRAS alleles was mutated. Model predictions for levels of RasGTP and RasGTP‐effector were essentially the same whether or not Ras turnover was included (Figure 1 b). This confirmed the validity of our original assumption that turnover could be ignored for predicting behavior of oncogenic mutants in the absence of covalent inhibitors.

We also investigated the behavior of the updated model by considering the loss of the tumor suppressor protein neurofibromin (NF1). The NF1 is a Ras GTPase activating protein (Ras‐GAP). Ras‐GAPs act as negative regulators of RasGTP by promoting the conversion of RasGTP to RasGDP. Germline absence of a single copy of NF1 results in increased RasGTP, and loss of both copies of NF1 results in further increases in RasGTP.26 NF1 is also one of the most commonly mutated driver genes in human cancer,27 in which one or both copies can be mutated, deleted, and/or silenced.

Modeling 100% and 50% loss of total Ras‐GAP in both the model with turnover and without turnover resulted in similar levels of RasGTP and RasGTP‐effector complex (Figure 1 c). We also considered fractions of total basal GAP loss other than 50% and 100% to consider conditions in which other GAPs contribute to RasGTP homeostasis.28 We considered 25% loss and 75% loss to span the range of levels of GAP activity, and found similar levels of RasGTP and RasGTP‐effector when the model with turnover was compared to the model without turnover. Overall, these simulations found that updating our model to include turnover results in essentially identical predictions for steady‐state levels of Ras signal, whether caused by oncogenes or by loss of tumor suppressor genes.

Extension of the model to covalent G12C inhibitors

We next extended the model to include the two types of KRAS G12C targeted covalent inhibitors (Figure 2 a,b). Our turnover model posits that NPIs bind to nucleotide‐free G12C mutant irreversibly, consistent with the manner in which NPIs have been described to date.8, 9 The SIIPIs are reported to bind to GDP bound G12C mutant (and presumably nucleotide‐free G12C mutant, as well), to bind irreversibly, and to prevent the inhibitor bound protein from exchanging GDP for GTP.11, 29 We modeled SIIPIs with these same activities. Table 1 includes the reaction parameters for NPIs and SIIPIs. Of note, the reported difference between the original SIIPI referred to as compound 127 and the new SIIPI ARS‐85311 are manifested exclusively in the on‐rate of drug binding.

Figure 2

Modeled biochemical reactions of KRAS G12C targeted covalent inhibitors. The Ras model was further updated to include the two classes of G12C inhibitors. (a) Nucleotide pocket inhibitors (NPIs) that bind to nucleotide‐free KRAS G12C are subject to degradation, and are incapable of binding other nucleotides once the drug is bound. (b) The switch II pocket inhibitors (SIIPIs) bind to GDP bound KRAS G12C and nucleotide‐free KRAS G12C, are subject to degradation, and GTP is not permitted to bind to drug bound G12C. All other interactions with GEFs, GAPs, and effectors are assumed to be impermissible, consistent with available data. (c) Simulated G12C inhibitor dose responses. The updated Ras model was used to simulate dose responses for NPI (blue) and SIIPI (red and black). For SIIPI, simulated dose responses were generated using both an on rate consistent with compound 12 described by Ostrem et al.7 (red) and an on‐rate consistent with ARS‐853 described by Patricelli et al.11 (black). For comparison, the dose response for covalent inhibitors in the case where protein turnover is not modeled is included (dashed line). GTPase, guanosine triphosphate.

Model validation: Testing the model's ability to reproduce the behaviors of G12C inhibitors

We then simulated dose responses for NPI and SIIPI (Figure 2 c). The NPI and SIIPI would seem similarly effective at steady‐state if turnover was not modeled. However, when turnover was included in the model, the SIIPI (ARS‐853 compound) was clearly superior to the other SIIPI and to the NPI. In this more physiological case, the maximal amount of drug that will bind to KRAS G12C is limited by the rate of each compound binding to KRAS G12C, which in the case of the NPI is restricted by the infrequency of nucleotide‐free mutated KRAS. This highlights the well‐recognized importance of considering turnover when evaluating covalent inhibitors.

Simulated behaviors of these drugs are consistent with their described abilities to inhibit KRAS G12C containing cancer cell lines, suggesting that our model includes the aspects of Ras biology needed to study KRAS G12C targeted covalent inhibitors. Experimentally, only the ARS‐853 compound is consistently effective at inhibiting KRAS G12C lung cancer cells. Lito et al.29 reported ARS‐853 has an half‐maximal inhibitory concentration of ∼2.5 microM; Patricelli et al.11 reported half‐maximal inhibitory concentration values ranging from 1 µM to 2 µM. In our simulations, the range of greatest change in Ras signal was in the low µmolar range, with ∼50% of Ras inhibition occurring at 1 µM, roughly mirroring the experimental observations. Lito et al.29 also reported an ∼95% reduction of Ras signal when 10 µM of ARS‐853 is applied, which is also comparable to the reductions in Ras signal found in our simulations. It is worth restating that these inhibitors are covalent and, therefore, do not have a Kd; they only have a kinetic reaction rate constant for the binding reaction. This highlights that the agreement between the model and experiment is not due to assigning a Kd for the drug in this range (as no such Kd exists to be assigned), but rather suggests that the agreement reflects an accurate portrayal of the system.

Model validation: Testing the model's ability to reproduce epidermal growth factor receptor modulation of Ras inhibition

Two recent studies demonstrate that the NPI ARS853 does not bind to GTP‐bound KRAS G12C, but rather binds to GDP‐bound KRAS.11, 29 That oncogenic Ras can exist in a GDP‐bound state was considered a surprising result by many Ras biologists.11, 29, 30, 31 However, multiple groups had previously appreciated that oncogenic Ras exists in a dynamic equilibrium with GTP and GDP‐bound fractions.12, 32, 33, 34 Of note, our original Ras model is based upon such a dynamic equilibrium.12 These new experiments may, therefore, serve as good tests for validating the application of our model to G12C inhibitors.

The recent SIIPI studies investigated how modulating the dynamic equilibrium would influence drug binding. Specifically, these studies probed the effects of epidermal growth factor receptor (EGFR) activation (or inhibition) on the inhibition of RAS signaling by SIIPI. The activation of EGFR was found to delay the kinetics of inhibition, and the inhibition of EGFR was found to expedite the kinetics of inhibition.29 This suggests that the dynamics of RAS activation can modify the effect of the drugs. To simulate these experiments, we modeled EGFR inhibition by decreasing the concentration of basally active GEFs in our model and we modeled EGFR activation by increasing the concentration of basally active GEFs. We found that our model reproduces the kinetic changes that follow from EGFR activation and inhibition (Figure 3 a). Of note, the changes in total RasGTP that our model predicts were also observable in the previous experiments. We also considered how EGFR/RTK activation would influence the ARS‐853 dose response. We found that the SIIPI dose response would shift to require higher levels of drug in situations of EGFR activation and to require less drug for situations of EGFR inactivation (Figure 3 b), consistent with the measured dose responses.11, 29

Figure 3

The Ras model predicts how receptor tyrosine kinase activation and inhibition influence ARS‐853 inhibition. Simulations of Ras G12C network inhibition by the switch II pocket inhibitor (SIIPI) ARS‐853 were performed, including conditions of upstream activation and upstream inhibition. (a) Simulations of the kinetics of the inhibition found increased GEF activity (upstream activation) delays the kinetics of inhibition (green), whereas decreased GEF activity results in faster inhibition (red). Kinetics for the basal model (no change in GEF activity) is shown in black. (b) The steady‐state dose responses of ARS‐853 in these same conditions of increased (dotted line) or decreased (dashed lined) GEF activity were also found through model simulations. Basal GEF activity is shown with the solid line. For both a and b, the plot on the left is normalized to the total amount of Ras signal in each case, whereas the plot on the right shows total RasGTP in each case.

Model validation: Testing the model's ability to simulate resistance promoting mutations

In their characterization of the behaviors of the SIIPI, Lito et al.29 engineered specific KRAS mutations that had both the G12C mutation as well as another mutation at a separate KRAS residue to introduce additional biochemical defects.29 They considered secondary mutations that would be anticipated to potentiate nucleotide exchange (e.g., Y40A, N116H, or A146V) and that would be anticipated to impair GTPase activity (e.g., A59G, Q61L, or Y64A). Lito et al.29 hypothesized and demonstrated that these secondary mutations would have reduced sensitivity to ARS‐853.

We evaluated our model's ability to reproduce these experiments. To simulate combinations of mutations within the same Ras protein, such as combinations of G12C with secondary mutations that are anticipated to reduce the intrinsic GTPase rate (e.g., KRAS G12C/A59G) or secondary mutations that are anticipated to increase the spontaneous exchange of nucleotide (e.g., KRAS G12C/A146V), we maintained the G12C parameters and adjusted only the rate constant(s) that were anticipated to change as a consequence of the secondary mutation. Our model simulations showed that increasing impairment of GTPase activity should make the mutant increasingly less sensitive to ARS853 (Figure 4 a) and that faster cycling secondary mutants should also be less sensitive to ARS853 (Figure 4 b); both model behaviors are qualitatively consistent with the previous experimental observations.29

Figure 4

Modeled secondary mutations to G12C and ARS‐853 dose response. Simulations were performed to evaluate secondary mutations in KRAS G12C that impair guanosine triphosphate (GTPase) activity and that altered the rate of nucleotide exchange. (a) Decreased rates of intrinsic GTPase activity are predicted to result in a decreased response to ARS‐853. Solid red: 10× slower GTPase activity than G12C, dotted red: 100× slower GTPase activity than G12C. (b) Increased intrinsic nucleotide dissociation is predicted to result in a decreased response to ARS‐853. Solid green: 10× faster, dashed green: 100× faster, both relative to G12C.

Model validation: Conclusion

We considered several different experiments that tested the inhibition of Ras signaling with the newly developed G12C inhibitors. We found that the model could readily reproduce experimentally observed behaviors, suggesting that our model is valid for problems involving G12C inhibitors. Additionally, it is worth noting that multiple experimental behaviors that were interpreted as contrary to conventional wisdom could have been predicted a priori had our readily available model first been applied to these problems.

Analysis of NPIs

We next used the model to investigate NPIs. We considered why NPIs are less effective than SIIPIs (Figure 5). A key step for the drug, as described, is binding to a nonoccupied nucleotide binding pocket.9 A nonoccupied nucleotide binding pocket is likely to be rarely encountered within the cell due both to the high affinity of the pocket for guanine nucleotides and also the high concentration of guanine nucleotides within the cell.

Figure 5

Simulations identify methods to improve nucleotide pocket inhibitor (NPI) effectiveness. The NPI described by Lim et al.8 and Hunter et al.9 was studied to determine how optimizing different biochemical rate constants might augment the effectiveness of NPIs. (a) An increased and decreased on rate of the compound was modeled: red: 10× slower, solid green: 10× faster, dashed green: 100× and 1,000× faster. B. The ability of GEFs to load the NPI, comparable to GEF loading of nucleotides, was considered. Solid blue: GEF loading of NPI proportional to the abundance of NPI and nucleotides. Dashed blue: GEF loading of NPI favoring NPI over nucleotides by a factor of 10 over their proportional abundance. Dotted blue: GEF loading of NPI favoring NPI over nucleotides by a factor of 100 over their proportional abundance.

We initially modeled NPI to bind the Ras NP with an on‐rate similar to that for nucleotides (kon = 2 × 105/Ms). We next considered on‐rates that were one, two, and three orders of magnitude faster to evaluate how much better inhibition would be if the on‐rate could be increased. Computational simulations suggested that a significant enhancement of the on‐rate would result in a much improved dose response, even when considering that nucleotide‐free Ras protein is very limited within the cell (Figure 5 a). However these levels are approaching the theoretical diffusion limit,35 suggesting that it would not be possible to engineer a better NPI by only optimizing the forward reaction rate constant.

We considered whether GEFs could facilitate the loading of an NPI. The authors of the NPI studies evaluated the ability of their compound to compete with GTP and GDP.9 Their assay included recombinant Ras protein, but did not include RasGEFs, Ras‐GAPs, nor Ras‐effectors. Thus, it did not assess how the network of reactions that together influence nucleotide binding would impact drug binding. As this compound cannot pass through the cell membrane, a cell‐based assay is not a simple option for investigating this problem. Mathematical modeling provides one mechanism by which we can extrapolate to the relevant, physiological, network and investigate this problem.

We were particularly interested in whether or not Ras GEFs could facilitate NPI loading. Our previous modeling studies of Ras suggest that a small level of basal GEF activity is needed to explain the experimental data of basal nucleotide exchange and that GEF inhibition should result in less oncogenic Ras signaling.12 We, therefore, considered the potential effects of GEF loading on NPIs.

Our simulations suggest that GEF loading could be an important variable for NPIs. If GEF loads the G12C inhibitor as well as it loads the nucleotide (e.g., it loads the G12C inhibitor like GTP and GDP, and at a rate that is proportional to the abundances of drug and nucleotides), our simulations suggest that there would be a modest increase in the ability of NPIs to target the G12C mutant (Figure 5 b, blue solid line). If GEF loading could favor drug loading by a factor of 10 to 100 over nucleotide, the effectiveness of this class of drug is predicted to be much higher (Figure 5 b, blue dotted and dashed lines). This computational analysis, therefore, suggests that efforts to optimize NPIs should evaluate candidates for whether or not GEFs facilitate their loading. Of note, the predicted levels of optimized (GEF‐facilitated) NPI (oNPI) needed to inhibit 50% of the Ras signal were on the order of 1 µM, much smaller than the near mM levels of nucleotide that are found within the cell.

Analysis of oNPI on resistance‐promoting mutations

We next considered how oNPI would respond to the mechanisms described to promote resistance to SIIPI. Our simulations found that increased intrinsic nucleotide cycling, such as what occurs with Y40A, N116H, and A146V compound mutations with G12C, could actually result in an increased sensitivity to NPI (Figure 6 a). This suggests that there may be value in the continued development of oNPIs. In contrast, our simulations found that secondary mutations to KRAS that impair GTPase activity (Figure 6 b) or mutations that increase RTK activation (Figure 6 c) would both result in less sensitivity to NPI. This suggests that an oNPI would not be able to combat all potential forms of G12C‐acquired resistance. It is worth noting that although the oNPI was less effective than the SIIPI for a G12C mutant, the oNPI was more effective than the SIIPI for all three resistance‐promoting mutations.

Figure 6

Simulations find switch II pocket inhibitor (SIIPI) resistant promoting mutations can make nucleotide pocket inhibitors (NPIs) the more effective KRAS G12C inhibitor. Potential SIIPI resistance promoting mutations were studied for both an optimized NPI (oNPI) and an SIIPI. The oNPI was assumed to have its loading facilitated by GEFs by a factor of four over the proportional abundance of NPI to total guanine nucleotides. (a) For the G12C mutant, the SIIPI ARS‐853 (black solid) was predicted to be superior to the oNPI (red solid), but once a secondary mutation that caused increased nucleotide exchange was modeled (100× over G12C) simulations found a decreased effectiveness of the SIIPI (black dashed), but a higher effectiveness for the oNPI (red dashed). (b) When a secondary mutation that causes decreased intrinsic GTPase activity was modeled (100× relative to G12C), both the SIIPI (black dashed) and oNPI (red dashed) were less effective. In the conditions of the secondary mutation, the oNPI became relatively more effective than the SIIPI. (c) When a secondary mutation increased RAS activation through increased GEF activation (e.g., as might happen with an RTK mutation), the oNPI (red dashed line) became more effective than the SIIPI (black dashed), although both classes of inhibitor are predicted to be overall less effective against mutations that result in increased GEF activation.

DISCUSSION

There is great interest in the covalent targeting of KRAS G12C,4, 10, 30 and multiple inhibitors exist at different stages of development.7, 8, 9, 11, 36, 37, 38, 39 Few of these compounds seem capable of passing through the cell membrane. Within the cell, Ras does not function in isolation but is rather a part of a complex signaling network. Modeling provides an approach for extrapolating to cellular conditions to evaluate how candidates in development may function in more physiological conditions. Modeling may also provide an approach for identifying which properties could be optimized within these physiological conditions, even before a derivative capable of being reliably delivered into the cell is developed. For investigating how a mutant will influence a signaling network, it is important that the model includes the proteins with which Ras directly interacts. Many models have been developed of the signaling networks that include Ras, and these typically include at least one representative GEF, GAP, and effector.40, 41, 42, 43, 44 It is much less common for pathway models that include Ras to dive into the details of Ras nucleotide exchange at the level needed to study oncogenic Ras activation,20 although a few models at this level of detail have been developed.12, 32, 43, 45, 46, 47 Systems pharmacology may be an area in which this additional detail becomes critical.48

One interesting finding from the present systems pharmacology study is that it suggests NPIs could be more effective if GEFs can promote their loading into the NP. This possibility seems consistent with known GEF biology and needs to be considered experimentally. However, the ability of Ras GEFs to load NPIs has not, to the best our knowledge, been addressed to date. The reason it has not been addressed previously may be because it is widely believed that the increased activation of oncogenic Ras is independent of GEFs. Our work suggests that efforts to screen derivatives of the known NPIs and new NPI candidates should include the evaluation of how well Ras GEFs can load the compound. This could be done experimentally by adapting the recombinant protein based assay used by Hunter et al.9 to also include a recombinant GEF and determining how NPI binding to RAS is enhanced by GEFs. Our work further suggests that compounds that are loaded by GEFs preferentially over nucleotides would be most beneficial.

Another interesting finding was that, theoretically, the better inhibitor could change as secondary resistance mutations are acquired. For example, the modeled GEF‐facilitated NPI is inferior to the SIIPI with respect to inhibiting KRAS G12C signaling. However, the model suggests that the GEF‐facilitated NPI should be superior to the SIIPI if a G12C mutant picks up a second mutation to the same allele that result in faster nucleotide dissociation.

A major challenge in developing targeted therapies for cancer is that it can be difficult to anticipate the response of a biological network to an inhibitor. Indeed, multiple experiments using KRAS G12C inhibitors have been interpreted as unexpected and contradictory to conventional wisdom. Computational systems biology models can provide an alternative viewpoint to expert opinion. It is possible to use a computational model to find the systems level behaviors that naturally emerge from the constituent reactions when the model is based upon the fundamental reactions of the network. We have done that here, and we demonstrate that multiple experimental observations that were interpreted as unexpected were foreseeable by computational modeling. As there are ongoing efforts to advance both strategies,36, 37, 38, 39 it may be valuable to integrate modeling into the drug development process.

Mathematical models have previously contributed to the study of G12C inhibitors. Modeling has been used to help analyze and interpret kinetic data.11 Modeling has also been used to extrapolate from data on reversible inhibitors to irreversible inhibitors. That study also included an argument that Ras GEFs may accelerate binding kinetics.49 Our model differs from these other models in that it includes a more complete set of Ras regulatory reactions and also protein turnover. Although this detail may be more than is needed for the questions addressed by the other models, this additional information is needed to evaluate signal inhibition in the cellular context (Figure 2). We anticipate that models like these will be increasingly utilized to inform efforts to develop and use targeted therapies.

Our consideration of methods to improve covalent Ras inhibitors made the assumption that increased Ras inhibition would be beneficial. Currently, most of these inhibitors display limited activity and limited selectivity, so increasing the level of Ras inhibition seems like a valid goal. However, it is possible that partial inhibition is sufficient to kill cancer cells and/or that modest changes smaller in magnitude than those studied here may have large clinical effects. At this early stage of development, we believe that our assumption that “more inhibition is better” is reasonable. As this class of drugs matures, it will be worthwhile to revisit this assumption.

Altogether, we have applied a computational systems biology approach to the analysis of targeted covalent inhibitors. Specifically, we used a model of the biochemical reaction network that regulates Ras signals to study how two specific KRAS G12C oncogenic mutant inhibition strategies alter Ras signaling. Our model finds the behaviors that logically follow from what is known about Ras signaling and about these inhibitors. We find that our model naturally reproduces the experimentally observed behaviors of both NPIs and SIIPIs, including results that were widely perceived as unexpected and inconsistent with known Ras biology. Our model also suggests strategies to improve the effectiveness of one of these classes of inhibitors.

Source of Funding

E.S. was supported by the Training Program in Immunology and Immunogenetics under award 5T32AI7163.

Conflict of Interest

A.S. is an employee of Genentech. E.S. declares that he has no competing interests.

Author Contributions

E.S. and A.S. wrote the manuscript. E.S. and A.S. designed the research. E.S. and A.S. performed the research. E.S. and A.S. analyzed the data.

Model Code

Click here for additional data file.