Quantifying ERK activity in response to inhibition of the BRAFV600E-MEK-ERK cascade using mathematical modelling

Mitogen-activated protein kinase (MAPK) cascades regulate key cellular processes. In BRAFV600E-mutant melanoma, constitutively active BRAF hyperactivates the BRAF-MEK-ERK cascade, driving proliferation. Vertical inhibition with BRAF and MEK inhibitors improves outcomes compared with monotherapy but resistance emerges rapidly. The study addresses the mechanistic question of how combined inhibition by dabrafenib (ATP-competitive BRAFV600E inhibitor) and trametinib (allosteric MEK1/2 inhibitor) modulates ERK activity over time, and how intracellular factors such as BRAF and ATP concentrations influence drug sensitivity and resistance. The purpose is to build a mechanistic, quantitatively parameterized model to predict ERK activity under mono- and combination-therapy to guide dosing strategies.

Prior work established frequent BRAF mutations in melanoma and clinical efficacy of BRAF and MEK inhibitors, with combination therapy delaying resistance compared to BRAF inhibitor monotherapy. Multiple resistance mechanisms exist, including BRAFV600E amplification. Xue et al. showed BRAF amplification confers a drug-dependent growth advantage and introduced a fitness threshold concept; Sale et al. reported MEK inhibitor addiction in BRAF-amplified cells. The foundational Huang and Ferrell MAPK cascade model demonstrated ultrasensitivity arising from dual phosphorylation. Subsequent studies explored ultrasensitivity and bistability, receptor-mediated inputs, and quantitative signaling behaviors. These works motivate a pathway-specific, mechanistic model incorporating ATP-dependent phosphorylation and explicit drug actions for BRAFV600E melanoma.

• Cascade and drug mechanism: The model specializes the MAPK cascade to BRAFV600E-MEK-ERK with constitutively active BRAF. Double phosphorylation activates MEK and ERK. ATP-dependent phosphorylation is explicitly modeled, producing ADP. Dabrafenib (DBF) competitively binds the ATP site on BRAF, blocking ATP binding and downstream phosphorylation of MEK/pMEK. Trametinib (TMT) allosterically binds MEK, preventing phosphorylation of ERK/pERK despite ATP and substrate binding.
• Reaction system: 36 chemical reactions (R.1–R.36) capture binding, phosphorylation, dephosphorylation, and drug interactions, including phosphatases for MEK and ERK. Order of binding events is considered.
• ODE formulation: Using mass action kinetics, reactions are converted to a system of ODEs (O.1–O.36) with state variables as concentrations of all species over time. Conservation laws (C.1–C.7) enforce total BRAF, MEK, ERK, phosphatases, and drug conservation across free, bound, and phosphorylated states, yielding a differential-algebraic system.
• Parameters and initialization: Forward (a_j), reverse (d_j), and catalytic (k_j) rate constants are set using literature data. Forward rates a_j are taken equal (a_j=a1) following Michaelis constant considerations K_m=(d_j+k_j)/a_j. Initial total concentrations (e.g., ERK_tot) and baseline intracellular BRAF (e.g., 3 nM) and ATP (e.g., 1 mM) are from literature sources. Full parameter tables and initial conditions are in the Supplementary Material.
• Output metric: Activated ERK(t) is defined as ppERK/ERK_tot.
• Implementation: MATLAB code solves the system numerically; code and usage instructions are available on GitHub (per Supplementary Material).

• Monotherapy dynamics: Dabrafenib responses are highly time-dependent and sensitive to intracellular BRAF and ATP; higher BRAF or ATP reduces sensitivity to DBF. Activated ERK under DBF can converge to similar steady states across different DBF-BRAF-ATP combinations despite differing transients. Trametinib monotherapy drives distinct steady states; beyond ~8 h, TMT efficacy is largely insensitive to ATP variation (1–5 mM) and only sensitive to low BRAF concentrations (<3 nM) in the tested range.
• Combination therapy and synergy: Heatmaps across DBF-TMT doses show that combinations below the isobole connecting monotherapy doses achieving activated ERK=0.5 are synergistic. Increasing BRAF or ATP from baseline shifts responses toward higher activated ERK, indicating reduced combination efficacy.
• Resistance determinants: Elevated BRAF concentration increases the minimum drug exposure required to suppress activated ERK below 0.5 at 24 h; similarly, elevated ATP increases resistance. Maps of minimum BRAF (with ATP fixed at 1 mM) and minimum ATP (with BRAF fixed at 3 nM) required to exceed activated ERK=0.5 quantify these thresholds.
• Clinical concordance: Simulated reduction in phosphorylated ERK (pERK+ppERK) versus dabrafenib concentration at BRAF=3 nM aligns with clinical biopsy data from Falchook et al. Under increased BRAF (10 nM), the model predicts marked attenuation of ERK suppression by DBF, consistent with BRAF amplification–mediated resistance observed in patient-derived xenografts.
• Mechanistic insight: ATP-competitive inhibition at BRAF makes DBF sensitivity dependent on ATP levels, suggesting elevated ATP as a contributor to DBF resistance; TMT’s allosteric mechanism renders its efficacy comparatively ATP-insensitive within the tested range.

The model quantifies how vertical inhibition modulates ERK signaling and clarifies determinants of resistance. It supports the clinical rationale for BRAF+MEK combination therapy and explains why resistance persists: modulation of pathway component levels (e.g., BRAF amplification) or metabolic state (ATP) sustains ERK activity under treatment. The ATP dependence of DBF efficacy introduces a metabolic dimension to resistance. Intratumoral heterogeneity in BRAF and ATP likely shapes spatial and temporal patterns of drug response, complicating regimen design. The framework can guide preclinical prioritization of dose combinations achieving sustained ERK suppression across G1. Extensions to include ERK inhibitors, schedule optimization, and multi-scale models (e.g., agent-based or age-structured) could connect intracellular ERK dynamics to cell-cycle progression and population-level resistance evolution. Comparison with clinical data demonstrates plausibility, while highlighting the need for accurate kinetic and concentration parameters for translational application.

A mechanistic, parameterized model of the BRAFV600E-MEK-ERK cascade with explicit ATP dependence and drug actions of dabrafenib and trametinib quantifies ERK activity under mono- and combination-therapy. It predicts synergistic regions for DBF-TMT dosing and identifies elevated BRAF and ATP as drivers of reduced sensitivity and potential resistance. The in silico approach can prioritize dose combinations for empirical testing and is extensible to include ERK inhibitors and treatment scheduling within multi-scale tumor models. Future work should refine kinetic parameters, incorporate patient-specific pathway component measurements, and link ERK dynamics to proliferation and cell death to enable personalized therapy design.

Model predictions depend on the accuracy of kinetic rate constants and intracellular concentrations of pathway components and phosphatases; these were drawn from literature and may vary across tumors and patients. The ODE framework assumes well-mixed compartments and does not capture spatial heterogeneity, stochasticity, or dynamic changes in protein expression and ATP beyond preset values. Validation is limited to literature comparisons and a clinical dataset for DBF; broader experimental validation and patient-specific parameterization are needed for clinical translation.