Modeling Tumor Growth Using Differential Equations

Mathematical oncology applies the tools of applied mathematics, physics, and engineering to understand and predict cancer dynamics. A central problem in this field is characterizing tumor growth kinetics from volumetric measurements over time. Understanding the growth law governing a particular tumor type has direct clinical implications: it determines how quickly untreated disease will progress, informs the timing and intensity of treatment, and constrains the design of dosing protocols for chemotherapy and radiation.

The simplest model is exponential growth: dV/dt = rV, which predicts unlimited growth at a fixed fractional rate r. This is a reasonable approximation during the early avascular phase of tumor development, when the tumor is small enough that nutrient and oxygen delivery by diffusion from surrounding tissue is adequate. However, tumors grow beyond roughly 1-2 mm³ in diameter, they exceed the diffusion limit for oxygen and nutrients. This triggers hypoxia, necrosis in the tumor core, and eventually angiogenesis, the recruitment of new blood vessels that dramatically changes the growth dynamics.

The logistic model introduces a carrying capacity K, representing a maximum sustainable volume: dV/dt = rV(1 − V/K). This is a more realistic approximation for tumors that eventually approach an apparent growth plateau, but it has a characteristic S-shape with an inflection point at V = K/2 that is not always consistent with observations. The Gompertz model, which has been empirically superior for many solid tumor types since Laird's 1964 study, uses: dV/dt = rV·ln(K/V). This produces growth deceleration that begins much earlier than in the logistic model, well below the carrying capacity, consistent with the observation that large tumors grow progressively more slowly not just because they're approaching a fixed limit, but because the growth rate itself continuously declines as volume increases. The von Bertalanffy model, originally developed to describe whole-organism growth, has also been applied to tumors, with a form that balances anabolic (volume-proportional) and catabolic (surface-area-proportional) processes.

We fit four ODE models (exponential, logistic, Gompertz, and von Bertalanffy) to tumor volume time series from the supplementary data of Benzekry et al. (2014, PLOS Computational Biology), a landmark compilation of digitized tumor growth curves from published experimental studies in mouse xenograft models covering breast, lung, colon, and melanoma cancer types.

Model fitting is performed using nonlinear least squares with the Levenberg-Marquardt algorithm, which is robust to poor initial parameter guesses. For each tumor series, we fit all four models and compute: root mean squared error (RMSE) on the fitted data, AIC corrected for small samples (AICc) for model selection that penalizes parameter count, normalized residuals to check model assumptions, and 95% confidence intervals on all parameters estimated by bootstrapping residuals. We define model selection as the fraction of tumor series in which each model achieves the lowest AICc. We also compute hold-out prediction accuracy by fitting to the first 75% of each time series and evaluating on the remaining 25%.