Growth rate-driven modelling reveals how phenotypic adaptation drives drug resistance in BRAFV600E-mutant melanoma

Summary Phenotypic adaptation, the ability of cells to change phenotype in response to external pressures, has been identified as a driver of drug resistance in cancer. To quantify phenotypic adaptation in BRAFV600E-mutant melanoma, we develop a theoretical model that emerges from data analysis of WM239A-BRAFV600E growth rates in response to drug challenge with the BRAF-inhibitor encorafenib. Our model constitutes a cell population model in which each cell is individually described by one of multiple discrete and plastic phenotype states that are directly linked to drug-dependent net growth rates and, by extension, drug resistance. Data-matched simulations reveal that phenotypic adaptation in the cells is directed towards states of high net growth rates, which enables evasion of drug-effects. The model subsequently provides an explanation for when and why intermittent treatments outperform continuous treatments in vitro, and demonstrates the benefits of not only targeting, but also leveraging, phenotypic adaptation in treatment protocols.


Introduction
Over 50% of melanoma cases present with BRAFV600E mutations that result in constituent activation of the Mitogen Activated Protein Kinase (MAPK) pathway and, by extension, uncontrolled cell growth and division. 1,2Targeted small molecule therapies that selectively inhibit mutant BRAF signalling have consequently been developed over the last decades and are now part of the standard of care for BRAFmutant melanoma.While these targeted therapies have lead to a remarkable increase in progression free survival, treatment resistance inevitably develops and limits long-term survival. 3This drug resistance in BRAF-mutant melanoma, and other solid tumours, is increasingly understood as a result of intratumoural heterogeneity, which has long-been viewed through the lens of inter-clonal differences. 4,5,6owever, a growing body of work implicates cell plasticity, i.e., phenotypic adaptation without associated genetic mutations, as a key component of therapy resistance in melanoma. 1,7,8dvances in transcriptomics and multi-omics have led to the identification of reversible and adaptive phenotype changes that confer treatment resistance and are mediated by a number of complex physiological factors on the single-cell level. 9,10,11As systems-level experiments of the complex interactions driving these adaptations are currently intractable, computational models have become increasingly important to our understanding of how therapeutic selection pressure shape the epigenetic evolution of malignant tumours. 12,13,14Many existing computational approaches leverage precise mechanistic representations of the MAPK signalling pathway that are informed by extensive perturbation experiments. 13,15hile these models permit a detailed description of the molecular signalling culminating in phosphorylated ERK, they are not directly linked with the cellular response.Further, due to the large number of chemical reactions being modelled, these models are computationally expensive to train and simulate.Here, we present a novel approach to understanding cell plasticity that is computationally inexpensive and directly bridges cellular phenotype with net growth rates.
Specifically, we develop a computational model that is directly informed by growth rate analysis of BRAFV600E-WM239A melanoma cells exposed to the BRAF-inhibitor encorafenib (LGX818) in a recent study by Kavran et al. 9 In our model, each individual cell is described by a plastic phenotype state that can take one of  discrete values, where  is an arbitrary integer larger than 0. Further, each phenotype state is directly related to drug dose-dependent growth rates.Modelled phenotype states are thus mathematical abstractions that link cells to drug resistance.To investigate phenotypic adaptation in response to drug exposure and removal, we propose four candidate cell-level phenotype adaptation strategies.Subsequent comparisons of model-generated cell populations and in vitro cell count data reveal that phenotypic adaptation is directed towards phenotype states with high net growth rates.As such, our individual-based model demonstrates the emergence of drug resistance at the population-level as a result of directed phenotypic adaptation on the individual cell-level.Further, our modelling captures the experimental observation that intermittent treatment strategies can outperform their continuous counterparts in vitro, despite halved cumulative drug doses. 9Accordingly, our framework can identify when intermittent treatments may be more effective than continuous treatments directly from multiple-dose growth rate data from untreated and treated cells.Building on the above results, we perform further simulation experiments to quantitatively assess the potential of not only targeting, but also leveraging, cell plasticity in melanoma towards improved clinical outcomes.

Empirical growth rates are mapped to phenotype states and give rise to a theoretical model
Our theoretical model formulation emerges from data analysis of drug dose-dependent growth rates of BRAFV600E-mutant melanoma cells exposed to the BRAF-inhibior encorafenib.The growth rates are shown in Figure 1a and are extracted from a recent experimental study by Kavran et al. 9 in which cell count numbers of WM239A-BRAFV600E melanoma cells were measured at 0 and 72 hours in response to drug challenge following 0-4 weeks of 500 nM encorafenib incubation.Cell populations that have been continuously incubated with drug for 1-4 weeks present dose-dependent net growth rates that are not significantly different from each other.We therefore assume that cells reach a drug-adapted state within 1 week of encorafenib incubation and stratify the observed growth rates as belonging to one of two extreme phenotypes: drug-naive cells and drug-adapted cells, i.e, cells that have been incubated in 500 nM encorafenib for 0 and 1-4 weeks, respectively (Figure 1b).Further, upon drug removal, drugadapted cells in the regarded in vitro experiment resumed growth rates similar to those of drug-naive cells 9 (results not shown) hence we allow for bi-directional phenotype updates in the model.
Building on our growth rate stratification, we formulate a discrete adaptive phenotype model in which each individual cell in a cell population is described by a discrete and plastic phenotype state that can be updated in response to external pressures, here encorafenib presence or absence.We use the variable   to describe the phenotype state of cell , where a cell's phenotype state corresponds to its sensitivity to the BRAF-inhibitor such that, arbitrarily, cells in state  = 0 are the most drug-sensitive and cells in state  = 1 are the most drug-resistant.In our model, growth rates for these phenotypes correspond to those for drug-naive and drug-adapted cells in Figure 1b, respectively.We further assume that the cells can exist in  (here  = 9) intermittent phenotype states between the extrema  = 0 and  = 1, where each state is associated with drug dose-dependent net growth rates that are estimated through linear interpolation of growth rate data and are mapped onto a fitness matrix (Figure 1c-d).In our discrete adaptive phenotype model, positive values in the fitness matrix correspond to daily probabilities that cells divide, whereas negative values correspond to daily probabilities that cells die.Our presented growth rate-to-phenotype state mapping enables growth rate-driven phenotype characterisation of cells, which we use to elucidate the directionality of phenotypic adaptation, predict treatment responses, and investigate the possibility of both targeting and leveraging cell plasticity in treatments.

Phenotypic adaptation in BRAFV600E-mutant melanoma is directed towards high fitness
Curiously, it was observed that four weeks of intermittent (one week on/one week off) treatments with 500 nM encorafenib are more effective at suppressing WM239A-BRAFV600E cell counts than four weeks of continuous treatments with the same dose in vitro. 9Here, we set out to elucidate the cell-level behaviour related to phenotypic adaptation that gives rise to this macroscopic, population-level result.To investigate the directionality of phenotypic adaptation in the regarded melanoma cells in response to encorafenib exposure and removal, we present four candidate phenotype update strategies for describing cell-level dynamics.These are illustrated in Figure 2a and are mathematically defined in the Methods section.
With the no update strategy, cells never update their phenotype states.With the unbiased strategy, cells propose to update their phenotype states to adjacent states with equal probability, unless at a boundary, and subsequently move to the proposed state.With the semi-biased strategy, cells propose to update their phenotype states to adjacent states with equal probability, but only accept the move if it results in a higher fitness, i.e., net growth rate.This corresponds to the cells having a noisy measurement of the fitness gradient.Finally, with the biased strategy, cells update their phenotype states to states of higher fitness whenever possible.As such, phenotypic adaptation is non-directed in the first two (no update and biased) strategies, but directed in the last two (semi-biased and biased) strategies.This gamut of phenotype progression is highly studied in segments of the mathematical literature 16 but, to the best of our knowledge, has not been applied in the interpretation of experimental data.Moreover, whilst the no update strategy suggests that the presence of drug-resistant cells that exist post treatment arise purely from pre-exisiting phenotypic heterogeneity, the other three strategies allow drug-resistant cells to arise through phenotypic plasticity.
The four update strategies are implemented through a simple computational algorithm (Figure 2b) and produce significantly different phenotype distributions in response to 300 nM encorafenib exposure and removal following continuous and intermittent treatments (Figure 2c).Whilst the non-directed strategies yield phenotype distributions that only moderately change over time in response to intermittent treatments, the directed strategies yield phenotype distributions that markedly move between the drug-sensitive ( = 0) and drug-resistant ( = 1) extrema.Simulation results pertaining to 0 and 300 nM treatment responses have been produced with dose-dependent net growth rates from the fitness matrix (Figure 1d) and the probability factor  = 1 (Figure 2a).Cells have a chance to update their phenotype state  off and  on times per day in drug absence and presence, respectively, where we set  on =  off = 2, ensuring that cells can traverse phenotype space within a week (following observations in Figure 1a).Comparing simulated cell counts in Figure 2d, our model predicts a distinct benefit in using 300 nM encorafenib intermittent over continuous treatments for the directed phenotype update strategies, but not for the non-directed strategies.
Further simulation results with 500 nM demonstrate that only the directed strategies are able to match in vitro data 9 and reproduce cell population dynamics such that 500 nM encorafenib intermittent treatments outperform continuous treatments, in terms of suppressing cell counts (Figure 2).Phenotypedependent net growth rates in response to 500 nM encorafenib are fit to data through a least-squares parameter estimation for each update strategy, where net growth rates are bounded by those for the adjacent drug doses for which data is available i.e., 300 nM and 1000 nM.The number of cellular phenotype updates per day,  on and  off are also optimised.
Taken together, the results in Figure 2c-d suggest that phenotypic adaptation in the regarded WM239A-BRAFV600E cells is directed towards phenotype states of high fitness through semi-biased or biased update strategies (Figure 2e).With the used model and available data, we can not distinguish which of the two directed strategies best describe the model system.However, the conclusion that phenotypic adaptation is directed provides a growth rate data-driven explanation for why intermittent treatments outperform continuous treatments at 500 nM encorafenib.
Moreover, with directed phenotypic adaptation, our model predicts that death events spike once treatments resume after a drug holiday, when the cells have acquired drug-sensitive phenotype states, in the intermittent treatments (Figure 2f).This simulation result is in qualitative agreement with Kavran et al.'s 9 experimental observation that the percentage of propidium iodide (PI) positive cells peak during drug re-challenge, where PI is used as a cell death marker, giving further credence to our model formulation and our simplifying model assumption that positive and negative net growth rates, respectively, can be approximated to yield birth and death events on the cell-level.

Simulations demonstrate the benefits of targeting and leveraging cell plasticity in treatments
Our results beg the question: Given directed phenotypic adaptation, when do intermittent treatments outperform their continuous counterparts?To address this question, we perform a series of simulation experiments in which output dynamics are quantified in response to variations of two model-inputs: the number of times per day that a cell can update its phenotype state (), and the duration of treatment on/off intervals ( ).Next, to enable biologically interpretable and mathematically tractable comparisons of continuous and intermittent treatments, we introduce the effective net growth rate  eff as the average (per cell and day) net growth rate that cells experiences over the course of a simulation.Using the biased update strategy, we derive analytical approximations through ordinary differential equation for  eff in the limit of large total treatment times.For 300 nM and 500 nM encorafenib, simulation results with instantaneous drug on/off switches and analytical expressions both identify pairings of  and  values that outperform continuous treatments, i.e., yield lower effective growth rates than those for fully drug-adapted cells in drug presence.These pairings are highlighted by yellow regions in Figure 3b,e.Throughout Figure 3, we have used phenotype-dependent net growth rates for 500 nM that are obtained through log-linear interpolation between net growth rates for 300 and 1000 nM in the fitness matrix (Figure 1a).
When all cells start in the most drug-sensitive phenotype state, the effective net growth rates increase with the number of cell updates per day for the regarded doses and simulation inputs (Figure 3a-c).In other words intermittent treatment become more effective when  is decreased.As such, our simulations highlight the therapeutic benefit of targeting phenotype adaptation by slowing it down.
Next, motivated by clinical studies in which drug re-challenge show anti-tumour activity in patients who have previously been treated with BRAF-inhibitors, 3,17,18 we investigate the impact of drug-holidays on treatment responses.We thus perform a simulation experiment in which all cells start in the most drug-resistant phenotype state (Figure 3d-f), as would follow from a period of continuous treatment.In such cases, our model produces non-monotonic relationships between effective net growth rates and the number of daily cell updates.This follows from the fact that cells may not have time to traverse phenotype space into states that are susceptible to sudden changes in drug exposure or removal.In an effort to push drug-resistant cells to drug-sensitive states in such cases, we can schedule a drug holiday and delay the treatment restart time.The benefits of delaying intermittent treatments are shown in Figure 3e, where the most effective treatment schedules are those that allow time for cells to traverse phenotype space to reach drug-susceptible states.Our simulation results thereby highlight the benefit of leveraging phenotypic adaptation in treatments.

Discussion
Cell plasticity is not listed as one of the original Hallmarks of Cancer which were famously collated in 2000. 19However, based on the last two decades' research advances, an updated version of the work identifies cancer cells' ability to unlock phenotypic adaptation as an emerging hallmark, intended to stimulate debate amongst researchers and inspire investigations that will improve our understanding of cancer. 8The work presented in this article directly contributes to this understanding by evidencing that phenotypic adaptation in BRAFV600E-mutant melanoma is directed to states of high fitness, and demonstrating that this adaptation enables evasion of drug effects.
Our model emerges from growth rate data analysis of BRAFV600E-mutant melanoma cells exposed to the BRAF-inhibitor encorafenib.Our conclusion that phenotypic adaptation in the regarded cells is directed towards state of high fitness is obtained through a model-selection procedure, in which a gamut of phenotype update strategies on the cell-level are evaluated against dynamic cell count data in response to intermittent and continuous encorafenib treatments.One notable consequence of our results is that dynamic phenotypic adaptation alone suffices to explain why intermittent treatments can outperform continuous treatments.This explanation offers an alternative to variations of the more common inter-clonal competition-models in which two or more distinct subpopulations of cells, such as a drug-resistant and drug-sensitive clones, interact with each other and compete for resources. 20he model developed in this paper is simple by design to enable direct mappings between phenotype states and growth rates.This simplicity, together with the informative results produced by the model, inspires and informs the design of further model-in vitro integrated research towards more detailed quantitative insight into phenotypic adaptation.For instance, data that clearly distinguishes cell division and death events could be used to refine our model to include growth and death rates, as opposed to net growth rates only.Further, our methodology to linearly interpolate between drug-naive and drug-adapted growth rates to assign fitnesses to intermediate phenotype states could be refined by further drug challenge experiments.These would chart out growth rates in phenotype-dose parameter space through cell counts measured at multiple time points after the pre-assay drug incubation, and ideally include pre-assay drug incubation performed at multiple doses.It would also be informative to have data revealing if and when drug-adapted cells can fully recover the growth rates of drug-naive cells.This would allow us to assess dose-dependent velocities of phenotypic adaptation.
Whilst our results clearly indicate that phenotypic adaptation is directed towards state of high fitness, we can not determine if the semi-biased or biased update strategy best describes the regarded cells.This follows from the fact that we only consider mean growth rates in the model and do not have access to data that reveal population-level phenotype distributions.The observation that cell count data alone can not be used to distinguish between semi-biased and biased phenotype adaptation opens up for a broader question: what type of experimental data is needed to identify phenotypic heterogeneity in cell populations and quantify how biased phenotypic adaptation is?To address this question, we suggest moving from our individual-based model to, for instance, a differential equation model in order to allow for rigorous mathematical analysis.
Mathematical models, analysis and simulations, especially those integrated with data, have made significant contributions to our understanding of cancer dynamics and treatment responses. 21,22To highlight a few contributions pertaining to BRAFV600E-mutant melanoma, Gerosa et al. 15 integrated proteomics and modelling to show that exposure to drugs that target the BRAFV600E-MEK-ERK pathway can cause non-genetic drug resistance by signal rewiring in the targeted pathway.Such rewiring was also identified as a mechanism of drug resistance by Fröhlich et al., 13 in an extension of the Gerosa et al. model, together with proteomic, transcriptomic and imaging data, to uncover that drug resistance can be mediated by the co-existence of two functionally distinct channels upstream of ERK, one initiated BRAFV600E monomers and the other by RAS.Whilst these models describe subcellular mechanisms that drive phenotypic adaptation, our model describes how cell-level phenotype traits, in the form of drugdependent growth rates, change in response to BRAF-inhibitor exposure and removal.Comparing all three models which complement each other, ours reveals less detailed subcellular information but is, on the other hand, considerably less data-intensive and therefore more accessible.Other data-integrated mathematical models of melanoma have been used to identify dose combinations of BRAFV600E, MEK, and ERK inhibitors that yield synergistic treatment responses, 23,24 predict initial treatment responses to BRAFV600E-inhibitors in xenografts, 25 and demonstrate phenotypic plasticity and multi-stability. 26Beyond melanoma, cell plasticity, and its impact on treatment responses and the evolution of therapy resistance, has been studied mathematically through structured partial differential equation models 27,28,29,30 which, notably, can be used to suggest model-informed treatment strategies that impede the evolution of resistance.
Treatment scheduling and responses were also investigated in our current study.Our model-generated simulation results suggest that phenotypic adaptation can be both targeted and leveraged in an effort to suppress in vitro cell counts.These results support the development of molecular inhibitors that target mechanisms that drive cell plasticity. 31Furthermore, strategies that leverage phenotypic adaptation, and in other words use cancers' ability to adapt against itself, can implicitly or explicitly be implemented through intermittent treatments and adaptive treatments. 32,33We remark that the treatment-related results presented in this study are based on data-driven models of in vitro systems, and that appropriate model extensions, as well as in vivo or clinical data, are needed to draw any conclusions about the role of phenotypic adaptation in xenografts or clinical tumours.However, one of the major benefits of modelling in vitro systems to study cell plasticity is the access to temporal and easy-to-interpret data that reveal how cancer systems change over time.Such data are paramount to understanding dynamical aspects of cancer.Importantly, our limited understanding of cancer dynamics has been proposed as one of the main factors that hinder the development of efficient targeted therapy protocols 34,35 and whilst most cancer treatments target genomic cancer aberrations without consideration for evolutionary tumour aspects, clinical trials that are informed by mathematical models and account for tumour dynamics have started to emerge. 36he growth rate-to-phenotype modelling pipeline presented in this work can be used to quantify phenotypic adaptation velocities, i.e., directions and rates, in any cell lines for which cell counts can be measured.Our pipeline can thus be used to estimate such velocities in cancerous, non-cancerous and mutated cell lines, and thereby quantitatively assess to what extent phenotypic adaptation stratifies cancer cells, mutated or not, from other cells.This functionality of our pipeline is of significant value to contemporary cancer research, as an emerging body of work has identified non-genetic resistance as a reason for drug resistance in multiple cancer types, including melanoma, 7 neuroblastoma 37 and lung cancer. 38One well-studied example mechanism of non-genetic resistance is the Epithelial-Mesenchymal Transition (EMT), 8,31,39 although a multitude of other mechanisms have been reported. 1,9,40,41,42,43Our growth rate-driven methodology is agnostic to these subcellular mechanisms and can therefore be used to directly examine phenotypic adaptation in its role as an emerging hallmark of cancer in melanoma and beyond.

Using growth rates to parameterise phenotype and drug dose-dependent fitnesses
To assign dose-dependent growth rate values to the intermediate phenotypes  ∈ (0, 1) for all doses included in the fitness matrix (Figure 1d) we linearly interpolate between the extrema  = 0 and  = 1.Since we do not have access to three-day cell count fold changes for 500 nM encorafenib, growth rates at this dose are instead estimated through a least-squares parameter fit evaluated against fold change data for 500 nM continous and intermittent treatments (Figure 2d).In the fit, four free model parameters are estimated: • FM(0, 500 nM): the growth rate for drug-naive cells in 500 nM encorafenib.
•  on : the number of times per day that cells have a chance to phenotype updates their phenotype when the drug is on.
•  off : the number of times per day that cells have a chance to phenotype updates their phenotype when the drug is off.

Modelling cell-level phenotypic adaptation with different update strategies
The four phenotype update (PU) strategies in Figure 2a are formalised below.In order, these are the no, unbiased, semi-biased and biased update strategies.We let  0 denote the current phenotype state of an arbitrary cell.
Probabilities to remain at the current state  0 are implied from the above equations through the relationship ( + ) + ( 0 ) + ( − ) = 1.The value FM(, ) can be read from the fitness matrix in Figure 1d, and denotes the net growth rate of a cell in phenotype state  that is exposed to drug dose .In Eqs.2a-d, (⋅) is the Heaviside function and  denotes the probability at which a cell updates it phenotype state in a random direction at update time (Figure 2a).

Deriving analytical expressions to identify intermittent treatment regimes that outperform their continuous counterparts
To predict when intermittent treatments with regular on/off intervals of equal lengths outperform continuous treatments with twice the drug exposure time, we derive an analytical expression for the solution to the model under the biased phenotype update strategy.We let  denote the interval length (i.e., the amount of time that cells are either on or off drug), and  denote the phenotype adaptation velocity such that /Δ gives the number of phenotype updates per day, where Δ = 0.1 is the step size arising from the choice  = 9 (Figure 1b).Implicit to the biased update rule is that, at any time in a simulation, there is no phenotypic heterogeneity between cells (Figure 2c).As such, all cells have an identical phenotype, denoted by (), the dynamics of which are governed by ,  < 1 and drug is on, −,  > 0 and drug is off, 0 otherwise. ( for drug doses 300 and 500 nM.Over a time period 2 , the cell count  (2 ) for cells on intermittent treatments is given by where the effective growth rate,  eff , is We here assume that cells are initially sensitive so that (0) = 0 and that the adaptation rate is identical in both the sensitisation and desensitization directions so that  =  on =  off from which it follows that (2 ) = 0. Finally, we must consider two cases depending on whether the cells are able to fully traverse phenotype space between the extrema  = 0 and  = 1 within a treatment interval or not.When they do,   ≥ 1. Solving Eq. 5 yields the effective growth rate for cells on intermittent treatment to be Eqs. 7 and 8 respectively correspond to the average growth rate and fitness difference for drug-sensitive and drug-resistant cells.To draw a fair comparison between cells on continuous and intermittent treatments, we compare the long-term effective growth rate of cells under each regime.To avoid effects of the initial drug desensitization for cells undergoing continuous treatments, we consider the effective growth rates to simply be given by FM(1, ).Thus, we expect intermittent treatment to outperform continuous treatment if  eff is lower than FM(1, ).We highlight here that the effective growth rate depends only on the product  , and not each parameter individually.Figure 3 is produced with net growth rates obtained through log-linear interpolation between FM(,300 nM) and FM(,300 nM) so that FM(0,500 nM) = -0.2866day −1 and FM(1,500 nM) = 0.0741 day −1 .Inserting this choice into Eq.6 for   ≥ 1, we find that intermittent treatment will outperform continuous treatment if   < 3.3748 or, equivalently,   < 33.748.

Implementing the discrete adaptive phenotype model and performing simulation experiments
The model is implemented in MATLAB.Information on how to access, run and modify the code files is available on the public GitHub repository https://github.com/SJHamis/phenotype_adaptation.

Figure 1 :Figure 2 :
Figure 1: Empirical growth rates are mapped to phenotype states and give rise to a theoretical model.(a) Growth rates for BRAFV600E-mutant melanoma cells from the WM239A cell-line are shown in response to various doses of the BRAF-inhibitor encorafenib, as measured over 72 hours.The rates are extracted from quadruplicate cell count assays by Kavran et al. 9 Prior to drug challenge, the cells are incubated in 500 nM encorafenib for 0,1,2,3, or 4 weeks (wks), with four replicates per incubation time.(b) Mean growth rates for drug-naive phenotypes (0 wks in a), and drug-adapted phenotypes (1-4 wks in a) are plotted over drug doses.Between the two extreme phenotype states,  intermediate states are introduced (here,  = 9).Marker colours and sizes correspond to phenotype state values and drug doses, respectively.(c) Net growth rates in (b) are mapped onto a discrete phenotype-dose space, as is indicated by marker colours and sizes.(d) Points in phenotype-dose space in (c) are mapped onto a fitness matrix (FM), in which the values of the left-most and right-most columns correspond to the growth rates in for drug-naive and drug-adapted phenotype states, respectively.Other values in the fitness matrix are obtained via linear interpolation between FM( = 0,   ) and FM( = 1,   ) for each of the  doses.Negative net growth rates are highlighted with yellow labels.

Figure 3 :
Figure 3: Simulations demonstrate the benefits of targeting and leveraging cell plasticity in treatments.(a) We consider an initial condition in which all cells are in the most drug-sensitive phenotype state ( = 0) in (a-c).(b) Effective growth rates in response to 8 week intermittent treatments with 300 nM (left) and 500 nM (right) encorafenib are obtained trough simulations with the biased phenotype update strategy.Mean values for 100 simulations are compared to effective growth rates for continuous treatments to show which  −  pairings outperform continuous treatments in terms of suppressing effective growth rates.Analytical results that identify these pairings are shown with dashed lines.(c) The plots show simulated cell count dynamics in response to 1 and 4 cell updates per day and treatment intervals of 1 and 4 days, as marked by stars in (b).These results demonstrate the benefit of targeting phenotypic adaptation in treatments by decreasing the number of cell updates per day.(d) We consider an initial condition in which all cells are in the most drug-resistant phenotype state ( = 1) in (d-e).(e) The experiments in (b) are repeated with the initial condition in (d).(f) The plots show simulated cell count dynamics in response to 1 and 4 cell updates per day and treatment intervals of 1 day, as marked by stars in (d).Treatments commence after a delay period of 1 (left), 4 (middle), or 10 (right) days.The results demonstrate the benefit of leveraging phenotypic adaptation in treatments by delaying treatments and allowing drug-resistant cells to desensitize by traversing phenotype space.
The probabilities of updating to higher ( + ) and lower ( − ) states, if possible with the restriction  ∈ [0, 1], are given below.