 Technical advance
 Open Access
 Published:
Different ODE models of tumor growth can deliver similar results
BMC Cancer volume 20, Article number: 226 (2020)
Abstract
Background
Simeoni and colleagues introduced a compartmental model for tumor growth that has proved quite successful in modeling experimental therapeutic regimens in oncology. The model is based on a system of ordinary differential equations (ODEs), and accommodates a lag in therapeutic action through delay compartments. There is some ambiguity in the appropriate number of delay compartments, which we examine in this note.
Methods
We devised an explicit delay differential equation model that reflects the main features of the Simeoni ODE model. We evaluated the original Simeoni model and this adaptation with a sample data set of mammary tumor growth in the FVB/NTg(MMTVneu)202Mul/J mouse model.
Results
The experimental data evinced tumor growth heterogeneity and interindividual diversity in response, which could be accommodated statistically through mixed models. We found little difference in goodness of fit between the original Simeoni model and the delay differential equation model relative to the sample data set.
Conclusions
One should exercise caution if asserting a particular mathematical model uniquely characterizes tumor growth curve data. The Simeoni ODE model of tumor growth is not unique in that alternative models can provide equivalent representations of tumor growth.
Background
Simeoni and colleagues [1,2,3] introduced and developed a pharmacokinetic/pharmacodynamic (PK/PD) model of tumor growth from in vivo animal studies, which could be summarized mathematically by a system of ordinary differential equations (ODEs). The model, and adaptations, have proved quite successful in predicting or modeling tumor growth and efficacy of cancer treatments [4,5,6,7,8,9].
A distinguishing characteristic of the Simeoni tumor growth model is that under chemotherapy, a drug’s action is not instantaneous: rather, tumor cells pass through progressive stages of damage because of the drug’s mechanism of action. This delay in drug action is modeled through a series of delay compartments through which the cells transit before elimination. In this note, we propose an alternative mathematical formulation of the Simeoni model using delay differential equations, without recourse to a series of delay compartments. We describe the two models in the next section, and then compare their performances with experimental data relating to mammary tumor growth in a mouse model. We conclude that one should be cautious if asserting that a particular tumor growth model uniquely characterizes experimental data.
Methods
The Simeoni model
The Simeoni model is depicted in Fig. 1. The central compartment Z_{1} represents the actively growing tumor, which increases according to a growth function TGF. c(t) denotes the concentration of a chemotherapeutic or immunotherapeutic agent in the central compartment Z_{1} at time t. In turn, c(t) induces a fraction of tumor cells to commit to cell death with a killing constant k_{1}. Tumor cells damaged by pharmacological treatment are shunted off successively into the peripheral compartments Z_{2}, Z_{3}, Z_{4} with a rate constant k_{2}, followed by elimination representing cell death. The observed tumor volume is the sum of cells in compartments Z_{1}, Z_{2}, Z_{3}, Z_{4}. The system of differential equations prescribing the Simeoni model is as follows:
with initial conditions Z_{1}(0) = V_{0}, Z_{2}(0) = Z_{3}(0) = Z_{4}(0) = 0. Total tumor volume is
and the tumor growth function TGF(t) is given by
Although we have depicted the Simeoni model with 3 peripheral compartments, we note that the number of peripheral compartments is in general arbitrary.
The model incorporating a delay differential equation is conceptually similar to the Simeoni model, the only difference being that the compartments Z_{2}, Z_{3}, and Z_{4} are replaced by a single compartment that incorporates a delay in elimination (Fig. 2). The system of differential equations describing this model is:
where
TGF(t) and c(t) are as before, V(t) = Z_{1}(t) + Z_{2}(t), and initial conditions are Z_{1}(0) = V_{0}, Z_{2}(0) = 0, and Z_{2}(t) = 0 for t < 0. The variable t_{2} represents the time delay before elimination.
In the absence of chemotherapeutic or immunotherapeutic intervention (i.e., c(t) = 0), tumor growth in the Simeoni formulation reduces to the biphasic process TGF(t), characterized by initial exponential growth followed by linear growth. In particular, in contrast to more classic models of tumor growth, e.g. logistic, Gompertz, Von Bertalanffy, there is no plateau or upper limit of tumor size. One of these alternative models might credibly be more appropriate than TGF in certain experimental settings. The generalized logistic model is defined by
the Gompertz model is given by
and the Von Bertalanffy model is
These three models have analytic solutions [10], which might aid in nonlinear fitting to experimental data.
In vivo tumor growth experiments
A series of experiments was undertaken at PRISM, the goal being to establish a viable and reproducible mouse model for mammary tumors. The FVB/NTg(MMTVneu)202Mul/J model was used for this purpose (Jackson stock #002376). All mice were bred at PRISM, with female mice genotyped by polymerase chain reaction to confirm transgene expression. Focal mammary tumors appear in female mice at about 100 days, reset as day 0 to mark the onset of the experiment. Data from one particular experiment are considered in the present work. In this experiment, 40 mice were randomly assigned to either no treatment (n = 21) or a single dose of cisplatin (5 mg/kg, n = 19) on day 0; all animals were relatively healthy and functioning on day 0, with all tumors smaller than 700 mm^{3} at initiation. The primary experimental outcome consisted of the series of tumor measurements from each animal, which were taken daily (excepting weekends, holidays) by a single individual (TJF) using digital calipers, and tumor volumes were recorded using the approximating formula V = (Length x Width^{2})/2, where L > W. [In this approximating formula, tumor shape is taken as an ellipsoid generated from the rotation of a semiellipsis around its larger axis (length), and a multiplicative constant involving p is ignored.] Mice were euthanized when tumors reached a volume of ~ 2500 mm^{3} (or measured 17 mm in length + width + height), or for complications, e.g., interference, ulceration. The method of euthanization was CO_{2} narcosis followed by cervical dislocation, a method consistent with the 2013 recommendations of the American Veterinary Medical Association Panel on Euthanasia. The protocol for this study was approved by the Institutional Animal Care and Use Committee of Prism, confirmed by the National Institutes of Health Office of Laboratory Animal Welfare Assurance, #D16–00819. Animals were housed in Prism’s animal care facility.
All 40 animals initially randomized to the untreated and the cisplatintreated groups were included in the subsequent analyses. There were a total of 386 observations among the 21 control animals, and 545 observations among the 19 treated animals.
Statistical methods
A nonlinear mixed effects (NLME) model was used to fit growth curves to the longitudinal tumor size data. In our context, a general form for NLME models would be:
where N is the number of individuals, n_{i} is the number of observations for individual i, t is the regression variable time, and y are the observations (tumor volumes). The term f is the structural model, expressed from the systems of ordinary differential equations shown earlier.
The residual error is g(t_{ij}, φ_{i}, ζ)ε_{ij}, where ε_{ij} ∼ N(0, σ^{2}) In our modeling,
we chose the function g to be a linear combination of a constant term and a term proportional to the structural model f with the additional parameters ζ = (a, b), that is, y = f + (a + b ∗ f)e.
The individual parameters φ_{i} are defined as follows:
where m is a pdimensional vector of fixed population parameters, l_{i} is a pvector of random effects, W is the p x p variancecovariance matrix of the random effects, and h is a fixed transformation. We assume that all of the model parameters are lognormally distributed among the individuals in each cohort, so that h(x) = exp(x).
Parameter estimation of the typical population values and the interindividual variability of each model parameter were estimated using Monolix 2019R1 (Lixoft, Antony, France). The Monolix suite of programs implements a stochastic approximation expectation maximization (SAEM) algorithm [11, 12] for maximum likelihood estimation of the model parameters. With the treated animals, following a single bolus injection of cisplatin we took c(t) to represent first order elimination, that is, a constant proportion of the drug is eliminated per unit time.
Model selection and comparison was based on the Akaike information criterion (AIC), the corrected Akaike information criterion (AIC_{c}), and the Bayesian information criterion (BIC). AIC_{c} is a modified version of the AIC, and includes a correction term for small sample sizes [13,14,15]:
where k denotes the number of free parameters, and n is the number of observations. In our setting, we also included variance components from the random effects in the parameter count.
Given a set of candidate models, each with a specific IC (AIC, AIC_{c}, BIC) value, we calculate IC model weights [15–17] for comparative purposes. We first compute for each model j the difference in IC relative to the IC of the best candidate model: D_{j} = IC_{j}  min IC. We then obtain the model weights by transforming to the likelihood scale and normalizing:
These weights are measures of the strength of evidence. They sum to one, and reflect the probability of each model j given the data and the M candidate models [15, 17].
Results
In the experiment outlined above, a total of 21 tumorbearing mice were untreated. Tumor growth over time in this group of animals is shown in Fig. 3.
We use this cohort of controls to investigate various tumor growth functions, namely, generalized logistic, Gompertz, Von Bertalanffny, and Simeoni. We present summary data relative to these fits in Table 1. It turns out that the Simeoni tumor growth function provides the best fit to these data, with minimal AIC, AIC_{c}, and BIC. The weights convey the overwhelming preference for the Simeoni tumor growth function among the candidate functions.
Plots of the individual fits with the Simeoni tumor growth function are shown in Fig. 4. Although the individual growth curves are quite heterogeneous (Fig. 3), the fitted curves quite nicely represent the observed data.
Nineteen tumorbearing mice were treated with cisplatin. Their individual growth patterns are shown in Fig. 5.
We used the Simeoni tumor growth function, and proceeded to fit the Simeoni model with varying numbers of peripheral delay compartments as well as the delay model to these data. Summary statistics are given in Table 2. Among the candidate models, the delay model evinces minimal AIC and AIC_{c} values, whereas the Simeoni model with 2 delay compartments minimizes the BIC. The weights reflect this model uncertainty: one might deem the delay model preferable to the Simeoni  2 and Simeoni  3 models relative to AIC or AIC_{c}, but this order is reversed with BIC. Regardless, the Simeoni  1 model appears inferior to the other candidate models.
Plots of the individual fits with the delay model are shown in Fig. 6. Again, interanimal variability is quite high, especially in terms of tumor response to treatment, with the individual fits modeling the observed data appropriately.
Discussion
The Simeoni tumor growth function switches between exponential and linear growth, but without a plateau phase, in contrast to the classic tumor growth functions. It works quite well in our experimental setting, in which tumorbearing mice were euthanized for ethical reasons when tumor volume reached the maximum allowable tumor mass, and this generally occurred before a plateau phase was reached. We caution, however, that in other experimental or clinical settings, the Simeoni tumor growth function might be supplanted by a more biologically realistic growth law [10, 18]. And, as a reviewer has pointed out, even the classical growth curve models logistic, Gompertz, and von Bertalanffy (our eqs. 5–7) enjoy many variations, specializations, and extensions, which might be judiciously employed to good effect in experimental settings.
The original Simeoni model essentially is a transit compartment model: tumor cells pass through progressive stages of damage because of the drug’s mechanism of action. Dying tumor cells stop proliferating and pass through several stages of dying (depicted by n delay compartments) before eventual death. However, there is some ambiguity over the exact number of delay compartments. With our experimental data, the appropriate number of delay compartments remains somewhat uncertain, relative to goodness of fit. Moreover, we have found that the delay compartments can effectively be replaced by a model incorporating a single peripheral compartment with a delay in elimination; that is, explicit incorporation of a delay term in the Simeoni tumor growth model can operationally replace the multiple transit compartments originally proposed, and effectively model the duration of the death process subsequent to drug action. Delay differential equation models have often been successfully investigated in tumor growth studies [19,20,21,22,23,24], and seem especially applicable in settings where tumor volume decrease appears delayed with respect to observed drug concentrations. With agents such as cisplatin, the delay in cell kill relative to the time course of extracellular exposure likely reflects the kinetics of cellular drug uptake and binding to intracellular targets [25]. Alternatively, more physiologically based models can ameliorate the ambiguity inherent with an arbitrary number of transit compartments [25,26,27,28].
With regard to model selection, we are more interested in the relative performances of the models rather than their absolute AIC, AIC_{c}, or BIC values. In this regard, the weights derived from the information criterion values quantify the probabilities of each model being optimal, given the data and the set of candidate models, hence provide useful criteria for model comparisons. There is little uncertainty about the suitability of the Simeoni growth function in Table 1, whereas the identifiability of the best approximating model in Table 2 is not incontrovertible. Regardless, the delay model is certainly competitive with the Simeoni model, and does provide excellent agreement with the experimental data, even with substantial tumor growth heterogeneity and interanimal diversity in response.
Conclusions
One should exercise caution if asserting a particular mathematical model uniquely characterizes tumor growth curve data. The Simeoni ODE model of tumor growth is not unique in that alternative models can provide equivalent representations of tumor growth.
Availability of data and materials
The datasets analysed in this study are available from the corresponding author on reasonable request.
Abbreviations
 AIC:

Akaike information criterion
 AIC_{c} :

Corrected Akaike information criterion
 BIC:

Bayesian information criterion
 nlme:

Nonlinear mixes effects
 ODEs:

Ordinary differential equations
 pk/pd.:

Pharmacokinetic/pharmacodynamic
 PRISM:

Proteogenomics Research Institute for Systems Medicine
 SAEM:

Stochastic approximation expectation maximization
References
 1.
Simeoni M, Magni P, Cammia C, De Nicolao G, Croci V, Pesenti E, Germani M, Poggesi I, Rocchetti M. Predictive pharmacokineticpharmacodynamic modeling of tumor growth kinetics in xenograft models after administration of anticancer agents. Cancer Res. 2004;64:1094–11.
 2.
Magni P, Simeoni M, Poggesi I, et al. A mathematical model to study the effects of drug administration on tumor growth dynamics. Math Biosci. 2006;200:127–51.
 3.
Simeoni M, De Nicolao G, Magni P, et al. Modeling of human tumor xenographs and dose rationale in oncology. Drug Discov Today Technol. 2013;10:e365–72.
 4.
Rocchetti M, Del Bene F, Germani M, Fiorentini F, Poggesi I, Pesenti E, Magni P, De Nicolao G. Testing additivity of anticancer agents in preclinical studies: a PK/PD modelling approach. Eur J Cancer. 2009;45:3336–46.
 5.
Zhou R, Mazurchuk RV, Tamburlin JH, Harrold JM, Mager DE, Straubinger RM. Differential pharmacodynamic effects of paclitaxel formulations in an intracranial rat brain tumor model. J Pharmacol Exp Ther. 2010;332:479–88.
 6.
Terranova N, Germani M, Del Bene F, Magni P. A predictive pharmacokineticpharmacodynamic model of tumor growth kinetics in xenograft mice after administration of anticancer agents given in combination. Cancer Chemother Pharmacol. 2013;72:471–82.
 7.
Li JY, Ren YP, Yuan Y, Ji SM, Zhou SP, Wang LJ, Mou ZZ, Li L, Lu W, Zhou TY. Preclinical PK/PD model for combined administration of erlotinib and sunitinib in the treatment of A549 human NSCLC xenograft mice. Acta Pharmacol Sin. 2016;37:930–40.
 8.
Pigatto MC, Roman RM, Carrara L, Buffon A, Magni P, Dalla CT. Pharmacokinetic/pharmacodynamic modeling of etoposide tumor growth inhibitory effect in Walker−256 tumorbearing rat model using free intratumoral drug concentrations. Eur J Pharm Sci. 2017;97:70–8.
 9.
Yates JWT, Holt SV, Armelle Logie A, et al. A pharmacokinetic–pharmacodynamic model predicting tumour growth inhibition after intermittent administration with the mTOR kinase inhibitor AZD8055. Br J Pharmacol. 2017;174:2652–61.
 10.
Benzekry S, Lamont C, Beheshti A, et al. Classical mathematical models for description and prediction of experimental tumor growth. PLoS Comput Biol. 2014;10:e1003800.
 11.
Delyon B, Lavielle M, Moulines E. Convergence of a stochastic approximation version of the EM algorithm. Ann Stat. 1999;27:94–128.
 12.
Comets E, Lavenu A, Lavielle M. Parameter estimation in nonlinear mixed effect models using saemix,an R implementation of the SAEM algorithm. Journal of Statistical Software. 2017. https://doi.org/10.18637/jss.v080.i03.
 13.
Sugiura N. Further analysis of the data by Akaike’s information criterion and the finite corrections. Communications in Statistics, Theory and Methods. 1978;7:13–26.
 14.
Hurvich CM, Tsai CL. Model selection for extended quasilikelihood models in small samples. Biometrics. 1995;51:1077–84.
 15.
Burnham KP, Anderson DR. Model selection and multi model inference: A practical informationtheoretic approach. New York: SpringerVerlag; 2002.
 16.
Akaike H. On the likelihood of a time series model. The Statistician. 1974;27:217–35.
 17.
Wagenmakers EJ, Farrell S. AIC model selection using Akaike weights. Psychon Bull Rev. 2004;11:192–6.
 18.
Mould DR, Walz AC, Lave T, Gibbs JP, Frame B. Developing exposure/response models for anticancer drug treatment: special considerations. CPT Pharmacometrics Syst Pharmacol. 2015;4:12–27.
 19.
Villasana M, Radunskaya A. A delay differential equation model for tumor growth. J Math Biol. 2003;47:270–94.
 20.
Cui S, Xu S. Analysis of mathematical models for the growth of tumors with time delays in cell proliferation. J Math Anal Appl. 2007;336:523–41.
 21.
Bodnar M, Piotrowska MJ, Forys U. Gompertz model with delays and treatment: mathematical analysis. Math Biosci Eng. 2013;10:551–63.
 22.
Xu S, Wei X, Zhang F. A timedelayed mathematical model of tumor growth with the effect of a periodic therapy. Comput Math Methods Med. 2016;2016:3643019.
 23.
Mazlan MSA, Rosli N, Azmi NS. Modelling the cancer growth process by stochastic delay differential equations under Verhults and Gompertz’s law. Jurnal Teknologi. 2016;78:77–82.
 24.
Kim PS, Lee PP. Modeling protective antitumor immunity via preventable cancer vaccines using a hybrid agentbased and delay differential equation approach. PLoS Comput Biol. 2012;8:e1002742.
 25.
ElKareh AW, Secomb TW. A mathematical model for cisplatin cellular pharmacodynamics. Neoplasia. 2003;5:161–9.
 26.
Ribba B, Watkin E, Tod M, et al. A model of vascular tumour growth in mice combining longitudinal tumour size data with histological biomarkers. Eur J Cancer. 2011;47:479–90.
 27.
Ribba B, Kaloshi G, Peyre M, et al. A tumor growth inhibition model for lowgrade glioma treated with chemotherapy or radiotherapy. Clin Cancer Res. 2012;18:5071–80.
 28.
Ribba B, Holford NH, Magni P, et al. A review of mixedeffects models of tumor growth and effects of anticancer drug treatment used in population analysis. CPT Pharmacometrics Syst Pharmacol. 2014;3:e113.
Acknowledgments
We than the reviewers for their insightful comments, which led to significant improvements in the manuscript.
Funding
This research is supported in part by grant 1PO1HL119165, R01CA169644 and U01CA193787 from the National Institutes of Health. The funding agency played no role in the design of the study, collection, analysis, and interpretation of data, nor in writing the manuscript.
Author information
Affiliations
Contributions
JES conceived and designed the underlying animal studies, and directed the project. TJF assisted with planning the animal studies and carried out the implementation and data acquisition. JAK performed the mathematical analysis and drafted the manuscript. All authors read, edited, and approved the final manuscript.
Corresponding author
Correspondence to James A. Koziol.
Ethics declarations
Ethics approval and consent to participate
This study was approved by the Institutional Animal Care and Use Committee of PRISM.
Consent for publication
Not applicable.
Competing interests
The authors declare they have no competing interests.
Additional information
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.
About this article
Cite this article
Koziol, J.A., Falls, T.J. & Schnitzer, J.E. Different ODE models of tumor growth can deliver similar results. BMC Cancer 20, 226 (2020). https://doi.org/10.1186/s1288502067030
Received:
Accepted:
Published:
Keywords
 Tumor growth
 Cancer chemotherapy
 Mathematical model
 Ordinary differential equations
 Delay differential equations