Phenotype-Structured Modelling of Tumour Plasticity and Drug Resistance under Adaptive Therapy

Full text

Text extracted from the article PDF for reading and search. Tables, figures and equations are in the PDF, which remains the authoritative version.

1. Introduction

Treatment resistance remains a central cause of failure in cancer therapy. A maximum-tolerated-dose strategy can produce rapid tumour reduction, but it also removes drug-sensitive competitors and creates ecological opportunity for resistant cells. Adaptive therapy instead aims to control rather than immediately eliminate the tumour by modulating exposure in response to disease burden. The approach is grounded in evolutionary competition: treatment holidays preserve or restore sensitive cells that can suppress resistant populations when resistance carries a fitness cost (Gatenby et al., 2009; Enriquez-Navas et al., 2016; Zhang et al., 2017) [1, 2, 3]. Clinical translation nevertheless requires models that link measurable tumour states, treatment timing, toxicity and evolutionary change.

Most mathematical descriptions of adaptive therapy divide a tumour into two or a few compartments. These models have clarified competitive release, multidrug scheduling and the stability of treatment cycles (West et al., 2020; Strobl et al., 2021; Vilela et al., 2025) [6, 5, 13]. Their discrete population structure is also convenient for optimal control. However, drug response is often graded rather than binary. Genetically similar cells may occupy intermediate drug-tolerant states, and therapy can induce reversible transcriptional or epigenetic changes. A two-state model can approximate total burden while missing the width, skewness and resistant tail of the underlying response distribution. Phenotype-structured partial differential equations address this limitation by describing a density over a continuous resistance trait. Diffusion represents random phenotypic variation, advection represents directed plasticity, and phenotype-dependent net growth represents selection. Such equations have been applied to treatment optimisation, metastatic resistance and the evolution of heterogeneous cancer populations (Almeida et al., 2019; Padovano and Villa, 2024) [9, 12]. Recent moment-reduction methods make these models more computationally accessible, but a gap remains between the mathematical representation of a full phenotype distribution and a practical feedback rule that uses evolutionary information without sacrificing positivity, mass balance or toxicity safety. A second gap concerns numerical reliability. Small negative densities produced by an unconstrained discretisation are not harmless: they can corrupt tumour burden, mean phenotype and variance, especially when the resistant tail is initially rare. Adaptive switching introduces additional non-smoothness, while non-local logistic competition couples every phenotype cell. A suitable scheme must therefore preserve non-negativity, conserve phenotype-space flux, resolve treatment switches and remain stable across long simulations. Numerical convergence should be demonstrated separately from biological sensitivity.

This study develops a controlled phenotype-structured PDE-ODE model with reversible plasticity, pharmacokinetics and toxicity. It establishes positivity and boundedness, derives exact equations for tumour burden and phenotype moments, identifies the exposure-dependent fitness optimum, formulates a risk-constrained control problem and constructs a positivity-preserving finite-volume method. Reproducible synthetic experiments compare four treatment strategies, assess convergence and quantify parameter uncertainty. The results demonstrate mechanisms and numerical behaviour; they are not presented as patient-specific predictions or evidence of clinical efficacy.

2. Materials and Methods

2.1. Phenotype Space and Modelling Assumptions

The continuous description avoids an artificial sensitive-resistant dichotomy. It permits a tumour to contain intermediate phenotypes and allows treatment to move the centre and width of the distribution. The variable �(�) denotes dimensionless drug exposure and �(�) is a cumulative toxicity index. Treatment intensity �(�)∈[0,1] is the control, with �=0 representing a treatment holiday and �=1 the reference maximum exposure. Four biological assumptions are used. First, sensitive phenotypes have a growth advantage without treatment, represented by a resistance cost. Second, treatment mortality decreases continuously with resistance. Third, random epigenetic or transcriptional changes generate diffusion in phenotype space. Fourth, treatment can bias plastic transitions towards resistance, whereas treatment withdrawal allows partial reversion towards sensitivity. These assumptions distinguish reversible plasticity from irreversible mutation while retaining a tractable model.

2.2. Phenotype-Structured Tumour Equation

The first term directs the population towards resistance during treatment, while the second represents partial re-sensitisation during a holiday. Because �(0,�)≥0 and �(1,�)≤0, drift points into the phenotype interval at its endpoints.

2.3. Adaptive treatment policies

where the initial 30-day exclusion prevents the baseline state from being classified as immediate progression.

2.4. Positivity and boundedness

The following proposition establishes biological admissibility.

Non-negativity follows from the parabolic comparison principle and the inward-pointing no-flux boundary condition. Bounds (11) follow from comparison with linear scalar equations. These estimates also exclude finite-time blow-up.

2.5. Moment Dynamics and Selection Decomposition

with �?????? collecting boundary terms. Phenotypic diffusion increases variance, whereas directional selection can contract the distribution around a resistant optimum.

2.6. Instantaneous Fitness and Treatment Switching

projected onto [0,1]. When � is small, the maximiser is sensitive because the resistance cost dominates. Sufficiently large exposure moves the optimum towards �=1. This calculation provides the mathematical reason for avoiding unnecessarily high continuous exposure.

in the sense that the selective advantage of low-� phenotypes exceeds the imposed reversion timescale. Condition (23) is sufficient rather than necessary, but it makes explicit why a weak cost of resistance undermines adaptive therapy.

2.7. Risk-Constrained Control Formulation

Equation (9) is a clinically interpretable feedback approximation to (27). It uses measurable summary states rather than the full adjoint trajectory.

2.9. Synthetic Experiment and Reproducibility

The study is a methodological numerical experiment, not a clinical trial. Parameters were selected to produce biologically interpretable competition and are reported in Table 1. The phenotype interval used �=101 nodes and ??????�=0.25 day for the six-year baseline simulations. Time to progression was censored at four years in the uncertainty study. A Latin hypercube design with 160 points varied five influential parameters over prespecified ranges.

3. Results

3.1. Analytical Behaviour

For the baseline functions, the highly resistant boundary retains positive fitness at tolerable exposure, so (33) is not satisfied. Continuous treatment can reduce the initial sensitive peak but cannot eradicate the resistant tail. The model therefore predicts evolutionary replacement rather than durable elimination.

Equation (21) also identifies a moving evolutionary target. At low tumour burden, logistic competition weakens and the cost term ����(1−�/�) becomes larger relative to density-limited proliferation. During sustained exposure, however, the drug term dominates and �⋆ approaches one. A holiday reverses the plastic drift and restores the fitness advantage of sensitive cells. The success of cycling consequently depends on the joint action of competitive suppression and partial phenotypic reversion.

3.2. Comparison of Treatment Protocols

Figure 1 compares the four protocols. Continuous exposure initially reduced burden most rapidly, but mean resistance increased monotonically and progression occurred at 618.5 days. Fixed 28/14 cycling delayed progression to 668.5 days. Burden-adaptive treatment reached 742.75 days, whereas phenotype-aware feedback reached 762.5 days.

Tumour burden, mean phenotype and toxicity under four treatment strategies. All curves are generated from the synthetic baseline parameter set. Fig 1: Baseline dynamics under continuous, fixed and adaptive treatment.

The phenotype-aware policy increased progression time by 23.3% relative to continuous exposure and by 2.7% relative to burden-only feedback. It used 1855.4 treatment-equivalent days, 15.3% fewer than continuous treatment, while its maximum toxicity index was 1.124 compared with 1.222 under continuous exposure.

Table 2: Baseline protocol outcomes Protocol Progression; gain; treatment days; maximum toxicity Continuous exposure 618.50 days; 0.0%; 2190.00; 1.222 Fixed 28/14 schedule 668.50 days; 8.1%; 1462.00; 1.025 Burden-adaptive 742.75 days; 20.1%; 1902.00; 1.200 Phenotype-aware adaptive 762.50 days; 23.3%; 1855.37; 1.124 The gain from phenotype information was modest but consistent with the control mechanism. The policy reduced intensity when the mean phenotype was high, thereby avoiding strong selection at times when additional drug mostly benefited resistant cells. When burden increased, the burden signal prevented indefinite undertreatment.

3.3. Evolution of The Phenotype Distribution

Figure 2 shows the complete distribution rather than only its mean. The initially narrow density around �=0.12 broadened under diffusion and then shifted towards resistance. Treatment holidays produced temporary leftward movement, but the distribution did not return to its initial location because resistant phenotypes accumulated across cycles.

Snapshots of the phenotype density under phenotype-aware adaptive therapy. Fig 2: Redistribution of resistance phenotypes during adaptive treatment.

The variance first increased as therapy separated the leading resistant tail from the sensitive core. Later, selection contracted the distribution around high resistance. This non-monotone variance is important: a two-compartment model could reproduce the mean burden but would not represent the changing width or skewness that determines the response to future treatment.

3.4. Numerical Convergence and Conservation

Grid-convergence results for tumour burden and mean resistance. Fig 3: Grid convergence of the positivity-preserving scheme.

The observed spatial rates between the three finest grids were between 1.72 and 2.04, consistent with second-order diffusion and first-order upwind drift becoming small for the smooth baseline solution. No negative cell densities occurred beyond round-off tolerance. The difference between the change in discrete total burden and the integrated reaction contribution remained below 3.1×10−10 per time step, confirming flux cancellation.

Table 3: Grid-convergence results Grid (�,ℎ) Burden error; mean-phenotype error

3.5. Sensitivity to Resistance Cost and Induced Plasticity

The two-parameter surface in Figure 4 demonstrates a competition between the cost of resistance and drug-induced plasticity. Increasing �� prolonged control because resistant cells grew more slowly during holidays. Increasing �on shortened control by moving mass towards resistance even before selection alone would have produced the same shift.

Time to progression across fitness-cost and treatment-induced-drift parameters. Fig 4: Sensitivity of progression time to evolutionary parameters.

The surface was nonlinear. When &#3627408464;&#3627408479;<0.18, increasing the holiday length produced little benefit because resistant phenotypes carried almost no cost. When &#3627408464;&#3627408479;>0.32, the phenotype-aware policy achieved more than three years of control for moderate drift. This interaction indicates that resistance cost cannot be estimated independently of plasticity if adaptive schedules are to be personalised. Local elasticity of progression time ??????prog with respect to parameter &#3627408477; was calculated as &#3627408440;&#3627408477;=&#3627408477; ??????prog ∂??????prog ∂&#3627408477;. (34)

3.6. Parameter Uncertainty

Progression-time distributions across 160 low-discrepancy parameter combinations. Fig 5: Uncertainty in progression time under the four protocols.

Median progression times were 1.67, 1.72, 1.94 and 2.03 years for continuous, fixed, burden-adaptive and phenotype-aware therapy, respectively. The 95% interval for phenotype-aware therapy was 1.40-3.10 years, compared with 1.23-2.37 years under continuous exposure.

Table 4. Uncertainty summary Protocol Median; 2.5%–97.5% quantiles (years) Continuous exposure 1.669; 1.232–2.371 Fixed 28/14 schedule 1.723; 1.255–2.526 Burden-adaptive 1.940; 1.356–2.894 Phenotype-aware adaptive 2.027; 1.403–3.097 The intervals overlap, so the study does not claim universal dominance. Instead, phenotype-aware feedback shifts the distribution towards longer control for the parameter ranges examined. The broad upper tail also shows that accurate estimation of evolutionary parameters is more important than small refinements to the numerical time step.

3.7. Numerical hypothesis comparisons

For each Latin hypercube point, paired differences in progression time were calculated because every protocol used the same parameter vector. The phenotype-aware policy exceeded continuous exposure in 153 of 160 simulations and exceeded burden-only feedback in 119 of 160 simulations. A paired Wilcoxon signed-rank calculation gave &#3627408477;<0.001 for phenotype-aware versus continuous treatment. This numerical significance describes robustness across simulated parameter sets; it is not a clinical efficacy test. The effect size against continuous therapy was summarised by ????????????=median(??????PA−??????CT??????CT), (35)

which equalled approximately 0.19. The corresponding treatment-use contrast was negative, showing that improved control did not arise from systematically administering more drug.

4. Discussion

The central contribution is a phenotype-continuous adaptive-therapy model in which tumour burden, drug exposure, toxicity and reversible plasticity interact. Existing adaptive-therapy models often divide the tumour into a small number of sensitive and resistant populations (Cunningham et al., 2018; West et al., 2018; Greene and Gevertz, 2021) [4, 7, 22]. Such models are valuable for control design, but the binary representation forces all intermediate states into predefined compartments. Recent phenotype-structured formulations provide richer evolutionary descriptions (Lorenzi et al., 2016; Almeida et al., 2019; Padovano and Villa, 2024) [10, 9, 12], while recent moment methods reduce their computational burden (Villa et al., 2025). The present framework connects these directions to an implementable feedback strategy and a positivity-preserving numerical method.

The analysis clarifies why burden-only feedback can be incomplete. Two tumours with equal burden may have different phenotype distributions and therefore different future responses. Equation (16) shows that burden depends on the averaged proliferation and drug-kill functions. Equation (17) then shows that the same treatment simultaneously changes the mean phenotype through directed plasticity and differential selection. A policy based only on &#3627408449;(&#3627408481;) observes the consequence but not the evolving cause of resistance.

The numerical comparison supports the ecological rationale of adaptive therapy and competitive suppression developed in earlier experimental and mathematical work (Enriquez-Navas et al., 2016; Hansen et al., 2017; Gallaher et al., 2018; Hansen and Read, 2020) [20, 2, 18, 17]. Continuous exposure produced the quickest initial decline but the strongest long-term selection for resistance. Fixed cycling reduced drug use and modestly delayed progression. Burden-adaptive treatment performed better because holidays restored competition from sensitive phenotypes. Phenotype-aware feedback produced the longest baseline and median uncertainty-study control, while reducing cumulative exposure relative to continuous treatment. Its state-dependent logic complements recent optimisation and cycle-stability studies (Vilela et al., 2025; Li et al., 2025; Gallagher et al., 2025) [13, 15].

The result should be interpreted cautiously. The phenotype-aware advantage over burden feedback was small in the baseline case, and uncertainty intervals overlapped. This is scientifically useful: measuring phenotype is justified only when its information changes a treatment decision enough to offset measurement cost and uncertainty. The model suggests that phenotype information is most valuable when treatment-induced drift is appreciable, resistance has a measurable fitness cost, and tumour burden alone responds slowly. The mathematical decomposition also identifies candidate biomarkers. The mean phenotype &#3627409159; represents average resistance, while ?????? measures heterogeneity. A broad distribution with moderate mean may contain a clinically relevant resistant tail even when the average response remains favourable. This emphasis is consistent with evidence that switching, treatment-induced plasticity and non-genetic mechanisms can create reversible resistance (Pisco and Huang, 2015; Gunnarsson et al., 2020; Marine et al., 2020; Shen and Clairambault, 2020; Pressley et al., 2021) [25, 23,26, 27, 24]. Longitudinal single-cell, circulating tumour DNA, ex vivo drug-response or imaging-derived measurements could be mapped to moments or basis coefficients of ??????(&#3627408481;,&#3627408485;). The full density need not be observed at every clinical visit if a filtering method reconstructs it from sparse observations.

The toxicity variable is deliberately generic. It could represent a composite adverse-event burden or a specific laboratory marker after recalibration. The state constraint in (25) makes toxicity mathematically different from a soft penalty: no improvement in predicted tumour control is permitted to violate the safety boundary. Recent adaptive-therapy studies have begun to incorporate toxicity and discrete monitoring (Gevertz et al., 2026) [16]. The present PDE framework adds continuous heterogeneity and plasticity to that discussion.

The numerical method was designed around biological invariants and standard finite-volume principles (Eymard et al., 2000; LeVeque, 2002; Hundsdorfer and Verwer, 2003) [32, 33, 31]. Negative tumour densities are meaningless and can corrupt mean-trait calculations. Conservative fluxes ensure that phenotype redistribution does not artificially create or destroy tumour mass. Implicit diffusion avoids a restrictive parabolic time-step condition, while upwind drift maintains monotonicity. More accurate high-order reconstructions could be introduced, but only with limiters that preserve positivity and total burden. The continuous-trait formulation also sits within the broader theory of structured transport and selection equations (Perthame, 2007; Lorenzi et al., 2017) [29, 28], while its reaction-diffusion interpretation is related to classical mathematical-biology frameworks (Murray, 2003) [30].

Several limitations remain. First, the parameters are synthetic and the numerical results demonstrate mechanism rather than clinical effectiveness. Second, phenotype is one-dimensional. Real resistance involves multiple, sometimes non-aligned mechanisms, including multidrug interactions (West et al., 2020) [6]. A multidimensional trait model would be more expressive but raises severe identifiability and computational challenges. Third, the logistic competition term assumes a common carrying capacity and does not include space, immune response or metastatic migration. Fourth, drug pharmacokinetics are represented by one compartment. Finally, the feedback policy uses the model state without observation error or delay. Patient-specific biomarker dynamics (Brady-Nicholls et al., 2020) [19], therapy optimisation experiments (Carrère, 2017) [21], and the open translational questions catalogued by West et al. (2023) [8] indicate where future calibration and validation should concentrate.

Other extensions include spatial tumour structure, immune interactions, multiple drugs, metastasis and evolutionary rescue. The model could also be coupled to a moment closure, using the full PDE only when an a posteriori indicator detects excessive skewness or multimodality. Such an adaptive model hierarchy would preserve interpretability while reducing the cost of patient-specific optimisation.

5. Conclusion

This study developed a non-local phenotype-structured model of tumour plasticity and resistance under adaptive therapy. The model combines phenotype diffusion, treatment-dependent drift, ecological competition, pharmacokinetics and toxicity. Positivity and boundedness were established, and moment equations separated the effects of plasticity, selection and treatment on mean resistance and heterogeneity.

A conservative positivity-preserving finite-volume method was constructed for the coupled PDE-ODE system. Synthetic experiments showed that continuous exposure accelerated resistant replacement, while adaptive schedules delayed progression. Phenotype-aware feedback achieved the longest baseline control, reduced exposure relative to continuous therapy and retained its advantage across most uncertainty samples.

The novelty lies in using the evolving phenotype distribution directly in a clinically interpretable feedback rule while retaining mathematical and numerical guarantees. The framework does not provide a clinical recommendation without calibration. It provides a rigorous basis for studying when phenotype measurements add decision value beyond tumour burden and how evolutionary information can be incorporated into safe adaptive treatment.

Acknowledgement

The author acknowledges colleagues in applied mathematics, numerical analysis and mathematical oncology whose scholarly discussions motivated the development of this framework.

Data and Code Availability No patient-level, personal or confidential data were used. All findings are generated from synthetic parameter sets reported in the manuscript. The numerical algorithm, parameter values and random seed are described sufficiently to reproduce the computational experiment.

Ethical Approval Ethical approval was not required because the study used no human participants, animal subjects or identifiable clinical data.

Conflict of Interest

The author declares no competing financial or non-financial interests.

Funding

No external funding was received for this research.