Fractional-Order Temporal Dynamics of Tumour–Immune–Chemotherapy Interactions: A Reproducible Dimensionless Analysis and a Pharmacokinetic Caveat

Abstract Fractional differential equations are increasingly used to represent temporal memory in tumour–immune and pharmacokinetic models, but changing the derivative order without a consistent time scale can make comparisons difficult to interpret. We study a dimensionless, commensurate Caputo system for tumour, immune, normal-cell, and drug states. We establish positivity, boundedness, the tumour-invasion threshold, and local stability conditions. We then use a fractional Adams–Bashforth–Moulton predictor–corrector method with an explicit mesh-convergence study. Two conclusions correct potentially misleading interpretations of the model. First, the equilibrium locations are independent of the fractional order: decreasing the order changes transient time scales but does not by itself alter the model’s asymptotic states. Second, after a unit bolus the drug equation has the exact solution C ( τ ) = E α ( − kτ α ). Its maximum is therefore the initial concentration, not a delayed peak. For 0 1, the algebraic tail makes the infinite-horizon area under the concentration curve diverge. In the illustrative coupled simulations, the apparent late tumour reduction at small α is consequently driven by persistent modelled drug exposure and must not be read as evidence of superior therapy. All parameters are dimensionless and illustrative; no patient data are fitted. Source code and generated tables are included to make every numerical statement reproducible.