The contribution of age structure to cell population responses to targeted therapeutics
Abstract
Cells grown in culture act as a model system for analyzing the effects of anticancer compounds, which may affect cell behavior in a cell cycle positiondependent manner. Cell synchronization techniques have been generally employed to minimize the variation in cell cycle position. However, synchronization techniques are cumbersome and imprecise and the agents used to synchronize the cells potentially have other unknown effects on the cells. An alternative approach is to determine the age structure in the population and account for the cell cycle positional effects post hoc. Here we provide a formalism to use quantifiable lifespans from live cell microscopy experiments to parameterize an agestructured model of cell population response.
Keywords: Cell cycle, intermitotic time, renewal equation, exponentially modified gaussian
1 Introduction
When examined individually in time lapse microscopy experiments, cells grown in culture display variability in the length of their cell cycles (between mitotic events), and this variability is represented by intermitotic time (IMT) distributions [21, 34, 35]. These distributions are usually obtained from asynchronously dividing populations of cells, which achieve a steadystate age structure when the population is growing exponentially. In experimental studies that examine the effects of perturbations on cellular proliferation, it is often desirable to determine whether the perturbation is affecting cells in particular stages within their cell cycle, i.e. whether the perturbation has cell cyclespecific effects. Since we do not know where a cell is in the cycle, an alternative is to define the age of a cell as the time elapsed from its last division and deal with this measurable quantity instead of the cell cycle phase. In this paper we adopt this definition and provide a formalism to convert data obtained as IMT distributions to parameterize an agestructured population model and, thus, identifying the contribution of age structure to the response of the cell population to perturbation.
Models of agestructured populations using partial differential equations, such as those originally discribed by A. Lotka, A.G. McKendrick, W.O. Kermack, and F. von Foerster, are well adapted to model the dynamical features of experimental cultures of cells transiting the cell cycle with variable IMT. These models have been widely studied from a mathematical perspective [4, 3, 15, 16, 17, 22, 37, 40, 43, 54], but the application of the model to experimental data has been hampered by an inability to determine the agedependent model parameters. The usual approach for parameter estimation is to solve numerically an inverse problem (see [1, 9, 10, 20, 19, 23, 28, 42, 44, 45, 51, 52, 41] on this question for structured population models), but this requires extensive input data and is specific to a given situation. A much more convenient approach is to assume that the distributed parameters lie in a class of functions with only a few constants (power laws or other forms [24]) and obtain the parameters from the fit of the functions to experimental data.
One function that often provides a better fit to IMT data than others, such as lognormal, inverse normal or gamma functions, is an exponentially modified Gaussian (EMG) [27, 53]. Under conditions in which the IMT distribution can be explained by an EMG model, we submit that the agedependent division rate can be identified as an error function. Starting from this observation, we present a simple method to recover the parameters of this error function from parameters fitted to the experimental IMT data. Once reliably parameterized, the agestructured model can be used to make predictions about cell agedependent effects of perturbations, for example, whether cells arrest during their cell cycle in response to treatment with antiproliferative compounds.
Individual based models (IBM) have also been used to track individual cell behavior in proliferating cell populations [5, 36, 26, 39, 46, 47, 48, 49, 50]. IBM have advantages in their direct connection to discrete events and to cellcell interactions. IBM are usually readily implementable to computer simulations, but may require many repeated runs to access an average outcome. IBM are usually difficult to analyze theoretically with respect to parametric input, particularly with highly sensitive parameters. The advantage of continuum differential equations models, such as the cell age structured models we develop here, are their tractability for theoretical analysis and their reproducibility for specific parametric input. The practical difficulties of determining this parametric input is the problem we focus upon here. IBM and differential equations models should be viewed as companion approaches, which complement and compare their insights and predictions for proliferating cell population behavior, when individual cell variation is a primary consideration.
2 The agestructured model
As with human populations, we can associate an age to each individual cell in a cell population. We define the age of a cell as the time elapsed from its last mitosis. Starting from this definition, we derive an agestructured model giving the dynamics of populations of proliferating and quiescent cells. This model is adapted to investigate the effect of treatment by the drug erlotinib on in vitro PC9 cancer cell lines. It has been recently shown in [53] that the main effect of erlotinib on cancer cells is to induce entry into quiescence. In [53] a system of ordinary differential equations (without age structure) is used to model these experiments. We hypothesize here that erlotinib induced quiescence is linked to the age of the cells involved. Many agestructured models with quiescence can be found in the literature (see [8, 11, 12, 29, 30, 31] for examples). Here we present a simple agestructured model with quiescence which allows us to explain observed delays in response to erlotinib. We start from the observation that there is no effect of the treatment on total population growth during the first twenty hours (see [53] and Figure 7). This time corresponds almost exactly to the minimal age of division observed for PC9 cells. It suggests that erlotinib acts only during a specific phase of the cell cycle, which based on its biological activity would be expected to be in G1. The model we present is based on this idea and considers a fractional rate of cells that become quiescent in an agedependent manner, where the fraction is assumed to reflect the dose of erlotinib used for the treatment.
We start from the the McKendrick–Von Foerster’s model which is widely used to model the cell cycle (see [37, 40, 54] and references). The partial differential equation in this model provides the evolution of the density of cells with age, or “cell cycle phase”, at time To take into account the quiescent cells, we introduce the quantity which representsthe density of quiescent cells at time It evolves according to an ordinary differential equation coupled to the equation on We obtain the system
(1) 
In this model, the proliferaing cells age onetoone with time at speed , and divide with rate When a cell divides at mitosis, it produces two daughters which are either proliferating with age , or quiescent. This is taken into account by the boundary condition at and the first term in the equation for where the parameter represents the fraction of proliferating cells which become quiescent. The coefficient is a death rates. The number of proliferating cells at time with age between and is , and the total number of cells at time is
The division rate has a probabilistic interpretation: the probability that a cell did not divide by age is given by
Since all the proliferating cells must divide at some time by definition, the division rate has to satisfy
(2) 
The model is completed with intial data and We choose the time to be the beginning of the erlotinib treatment, so at this time there are no quiescent cells . The age distribution of the proliferating cells is assumed to be at equilibrium (i.e. see Section 3 for the mathematical definition of ). The experimental values of the total population along time are ploted after normalization by the initial value on a logscale (see Figure 7). Because of this normalization, we consider an initial distribution such that which leads to because of the definition of
We want to compare the solutions of model (1) to the experimental observations presented in [53]. The first step is to estimate the different parameters of the model. To estimate the value of the fraction for different doses of treatment, we use data available in [53]. The fraction of quiescent cells is estimated by examination of whether cells treated with erlotinib divide or not before the end of the experiment. We want to use this experimental fraction to estimate the coefficient of the model. For the sake of simplicity, assume that the death rate is in (1). In this case, at the end of the labeling period , the quantity of labeled quiescent cells corresponds to and the quantity of proliferating cells corresponds to Thus, the fraction of quiescent cells at is
(3) 
Now we compute this quantity from model (1), keeping in mind that we have assumed no mortality. We have
and
Finally we obtain
(4) 
So when there is no death rate, the experimental fraction corresponds exactly to the coefficient In our simulations we consider positive a death rate as suggested by the results in Section 4.2. But because the numerical value we recover for is very small in Figure 6), we can consider that is still wellapproximated by and we use the fractions available in [53] to parameterize
The identification of the coefficients and is much more delicate because is a distributed function. The two following sections are devoted to presenting a method to recover these coefficients from IMT distributions. Then these coefficients are used to compare the solutions of Equation (1) to experimental data in Section 5.
3 Modelling the intermitotic time
Since the age of a cell is defined as the time elapsed from its last mitosis, the IMT of a cell is its age at division. This definition allows us to interpret the IMT distributions in terms of a dynamic agestructured population model.
Experimental IMT distributions can be seen as histograms which represent, for a given population, the density of cells with a certain age of division (see Figure 1 for an example taken from [53]). The age of division is distributed into bins of width For all between and the height of the bar represents the density of cells with an age of division in the window The histogram is normalized to represent a density
(5) 
We briefly explain here the method used in [53] to build IMT histograms (we refer to this paper for more details). The data are obtained using extended temporally resolved automated microscopy (ETRAM) in which cell nuclei are fluorescently labeled, imaged by automated time lapse fluorescence microscopy, and tracked as individual cells from the resultant image stacks. Cells are subjected to various microenvironmental conditions (such as the addition of a drug) at a specific time during image acquisition (set to time and the effect of the perturbation from that point in time is followed across the entire population and individual cells within it. The duration of observation is chosen large enough to observe that almost all cycling cells divide before this final time. The data are organized into the bins to obtain a histogram, which is then normalized to represent a density of cells.
The next step is to describe the IMT distributions in terms of an agestructured model and use this interpretation to parameterize Equation (1). We want to estimate the coefficients and and for this we fix in (1). This corresponds to the situation when there is no erlotinib treatment and thus all the cells are proliferating. We recover in this case the original McKendrick–Von Foerster’s model
(6) 
It is a transport equation which satisfies a maximum principle, namely if the initial distribution is nonnegative (positive) then the distribution remains nonnegative (positive) for all time The solutions of (6) have a remarkable behavior as time evolves, in that the age structure equilibrates no matter what the initial age structure of the population may be. This effect is known as asynchronous exponential growth, and its interpretation is that the population disperses over time to a limiting asymptotic equilibrium age structure where the fraction of the population in any age range satisfies
= a constant independent of the initial age structure [55]. Moreover, the solutions to this equation, as , are known to behave like a separated variables solution, that is, a function of time only a function of age only. More precisely, consider the eigenvalue problem
(7) 
This problem has a unique solution given by
where is the unique value such that
(8) 
and
Then we can prove that, for large times,
(see Appendix A for more details and references). If the population of cells proliferates over a sufficiently long time, we can assume that this asymptotic behavior is reached and use it to investigate the IMT distributions. An experimental observation that the total population (independent of age structure) is growing exponentially is an indicator that the population has effectively reached the equilibrium age distribution, which can be checked using ETRAM for other time series data collection (see Figure 2).
Now we give a continuous expression of the IMT distribution in terms of the agestructured model. The age distribution of the cells relative to time (time of perturbation) is given by a truncation of the equilibrium age distribution
where is a scaling constant. We then follow this age distribution along time and obtain for
According to the agestructured model, and because of the normalization (5), the IMT distribution satisfies
(9) 
by definition of the division rate The fact that no cell can divide in a time less than means
(10) 
Under this condition, for large, the function is close to the function
(11) 
where
(12) 
This convergence is made mathematically precise by the following claim (the proof is given in Appendix B for a more general case in which it is not assumed that is close to the equilibrium age distribution):
Remark 2.
Because in the experimental protocol is chosen large enough to observe no dividing cells at the end, we can approximate by which has a simple expression in terms of the division and death rates (11).
4 From the intermitotic time to the division rate
In this section we explain how we can recover the parameters and from the IMT distribution We first present the method in the case when there is no death This simplification is useful to give an inversion formula which gives the rate in terms of Then we extend the method to include possible death rates.
4.1 The case
In the case the constant is equal to 1 and we have
This expression can be inverted to recover the division rate from the IMT distribution (see [32, 14])
(14) 
We start from the fitting of the experimental IMT distributions, which are observed to be positively skewed. From a fitting procedure for we use Equation (14) to recover the division rate of the McKendrick–Von Foerster’s equation.
First consider as in [32] that the IMT distribution is a shifted gamma function (see also [13] and references therein). Setting
(15) 
we can solve explicitly Equation (14) and we find
(16) 
So by fitting an experimental IMT distribution with a shifted gamma function, we obtain two parameters and which allows reconstruction of the division rate of the renewal equation.
To have a smoother transition at the minimum age of division one can consider a second shifted gamma function
(17) 
Then the corresponding is
(18) 
The different functions (15) to (18) are plotted in Figure 3
It has been observed that an exponentially modified Gaussian (EMG) is often a better model for IMT distributions than the gamma function [27, 53] . An EMG is defined as the convolution of a Gaussian with a decreasing exponential, but after solving it can be written with three parameters as
(19) 
where the (complementary) error function is defined by
Replacing by an EMG in Equation (14), we cannot compute explicitly the expression for But by numerical comparison, we obtain a division rate that is essentially indistinguishable from an error function (see Figure 4).
Instead of fitting the IMT distribution by an EMG and then fitting the corresponding by an error function, we may directly assume that is an error function
(20) 
We can then explicitly derive a new fitting formula for due to Equation (11)
(21) 
where the integral in the exponential can be computed as
(22) 
Using this formula for the IMT instead of an EMG formula, the fitting parameters provide immediately the division rate
When we know the Malthusian growth parameter (population growth rate from experimental data (see Figure 2), it is possible to recover a division rate such that relation (8) is satisfied. This is an important step in order to parameterize and apply an agestructured model to experimental data. We notice the relationship of in terms of is
To take advantage of this relationship between the IMT data, and the malthusian parameter at the equilibrium distribution, we define a new histogram by
(23) 
where is the mean age of the bar. Thus defined, the new histogram incorporates information about and satisfies the relation
(24) 
We then fit this new histogram with the model
(25) 
instead of fitting with Because of the normalization (24), we expect that the fitting provides parameters such that , and this relation can be checked numerically a posteriori.
Method and example.
We divide the process in three steps and illustrate it by an example. Each required fitting step can be performed using the freely available Ezyfit Matlab toolbox
\([<\mbox{\url{http://www.fast.upsud.fr/ezyfit/}}>].\)

Step 1: Determine equilibrium IMT distribution

Step 2: Obtain parameters from model fit to IMT distribution

Verify that the correlation coefficient of the IMT data and the chosen form is close to

Step 3: Parameterize agestructured model
4.2 Introducing a death rate
If we include a death rate in the model, then there does not exist an inversion formula as Equation (14) to recover analytically. So we consider that is a rational fraction given by Equation (16) or (18), or an error function given by Equation (20), and we derive the corresponding model from Equation (11). The death rate will appear in the expression of as an additional fitting parameter and will be determined together with the parameters of by the fitting of experimental IMT distributions. Notice that a fixed rate based on measured data can also be used if available to avoid the addition of a free parameter.
Unlike the case for which we have for a positive death rate that This constant is a function of the fitting parameters (see Equation (12)), but we do not have an analytic expression of this function in general. A solution to this problem is to use the Malthusian parameter which incorporates information about and to define a form as in Section 4.1. Indeed we obtain, because of relation (8),
which does not involve Then, the division rate is obtained from as before using the modified histogram defined from by (23).
Example. We fit the same distribution as in Section 4.1 still considering that is an error function. With a constant death rate we obtain the four parameters model
(26) 
where the integral is given by Equation (22). The fitting provides new parameters for and a positive death rate (see Figure 6). The correlation is slightly better than in Figure 5 and the integral is significantly closer to 1. So we can assume that mortality has to be considered for this cell line.
5 Numerical simulation and discussion
Once the parameters accurately estimated, we use them to compare the behavior of the solution of Equation (1) to the experimental data in [53]. Since the model does not have analytic solutions, we need to perform numerical simulations. Many numerical methods are available for structured population equations (see for instance [2, 7, 6, 18, 25, 33]). Here we use a scheme based on the method of characteristics, as in [7, 6, 33], for its antidissipative properties.
We can see in Figure 7 that the numerical simulations are very similar to the experimental curves. In particular, the delay of twenty hours before the effect of the treatment on the growth of the population is apparent, and this twenty hour period pulses two more times as the population approaches a new equilibrium distribution during the total time of the experiment. We have thus developed an agestructured model that can explain the dynamic effects of erlotinib on PC9 cells, which are intrinsically dependent on the age of proliferating cells. The model is also accurate in a quantitative point of view. The quiescent fractions which are used in the numerical simulations are not chosen arbitrarily to obtain the adequate behavior. They are linked to the treatment dose of erlotinib and estimated from experimental data in a rigorous way so that they are realistic.
The model has been chosen as simple as possible to be able to recover all the parameters from experimental IMT distributions. This simplicity leads to qualitative differences between observed data and simulations in Figure 7. Experimental observation suggests the population reaches equilibrium when treated with drug whereas simulated population size is still increasing at the end of the experiment. An agestructure for quiescent cells together with an agedependent rate of death should be consider to explain this plateau effect. Experimental data also suggest that drug treatment increases population growth rate above that of untreated cells briefly, but model does not capture this. This could be obtained by considering a lower death rate for quiescent cells than for proliferating cells. But for such more complex models, additional experimental data and a new parameter estimation method would be necessary to estimate the death rates.
Conclusion
Linking experimental observation of cell behavior between the singlecell and population scales has recently been described using newly developed mathematical models [53]. However, this approach does not take into account the possible cell agedependent effects of a perturbation on cell behavior, such as would be expected if the effects occur in cells at a specific position in the cell cycle. Since these studies were performed with asynchronously dividing cell populations, it is evident that the mathematical models of these experiments should have age structure as a primary feature. In fact, to fit the data, the authors had to use an artificial time offset to account for the agestructured effects. Here we provide a formalized approach to accurately account for cell agedependent effects on cellular behavior. A major difficulty in the parameterization of agestructured models is the determination of the age dependent division rate. Our study provides a method for the quantitative recovery of this rate by fitting experimental IMT distributions to special forms, such as gamma functions or exponential modified Gaussians. This model, once successfully parameterized, is very useful for simulating and analyzing agedependent phenomena in cell population dynamics. We have presented one such application for the in vitro treatment of cancer cells by erlotinib. This example shows the utility of agestructured population models in explaining the connection of drug therapy to phenomena such as cell cycle phase entry into quiescence. The method we have presented can be implemented readily to many issues in cell population behavior when there is experimental data based on cell age.
Aknowledgments
The research stays of P. Gabriel at Vanderbilt University have been financially supported by a grant of the Fondation Pierre Ledoux Jeunesse Internationale.
Appendix A Convergence to the equilibrium
The longtime behavior can be proven by using semigroups methods [37, 54] or General Relative Entropy techniques [38, 40]. We detail here a result provided by General Relative Entropy.
Consider the unique solution to the adjoint eigenvalue problem
Then the General Relative Entropy method allows us to prove that
where
Appendix B Convergence of the IMT distribution
Theorem 3.
For the sake of simplicity, the proof of Theorem 3 is given in the case But it can be easily adapted to the case with a death rate.
Proof of Theorem 3 in the case .
First we extend the initial distribution and the division rate to negative ages by setting for Since the daughters of the labeled cells are not tracked, the age distribution of labeled cells satisfies a transport equation without boundary condition. Thus it writes, using the characteristic method,
Introducing this expression in the definition of we obtain, with changes of variables,
Since the support of is included in and for due to Assumption (2), we have for all
so writes
Define the primitive
This function is nondecreasing and satifies for and for so we have
Because for we obtain by integration of on
which gives
Since
we conclude that
Then we have
and it ends the proof of Theorem 3.
∎
Footnotes
 footnotemark:
References
 A. S. Ackleh. Parameter identification in sizestructured population models with nonlinear individual rates. Math. Comp. Modelling, 30(910):81 – 92, 1999.
 A. S. Ackleh and K. Ito. An implicit finite difference scheme for the nonlinear sizestructured population model. Numer. Funct. Anal. Optim., 18(910):865–884, 1997.
 B. Ainseba and C. Benosman. Cml dynamics: Optimal control of agestructured stem cell population. Math. Comput. Simul., 81(10):1962–1977, 2011.
 B. Ainseba, D. Picart, and D. Thiéry. An innovative multistage, physiologically structured, population model to understand the european grapevine moth dynamics. J. Math. Anal. Appl., 382(1):34–46, 2011.
 T. Alarcón, H. Byrne, and P. Maini. A cellular automaton model for tumour growth in inhomogeneous environment. Journal of Theoretical Biology, 225(2):257 – 274, 2003.
 O. Angulo and J. LópezMarcos. Numerical integration of fully nonlinear sizestructured population models. Appl. Num. Math., 50(34):291 – 327, 2004.
 O. Angulo and J. C. LópezMarcos. Numerical integration of nonlinear sizestructured population equations. Ecol. Modelling, 133(12):3–14, 2000.
 O. Arino, E. Sánchez, and G. F. Webb. Necessary and sufficient conditions for asynchronous exponential growth in age structured cell populations with quiescence. J. Math. Anal. Appl., 215(2):499–513, 1997.
 H. T. Banks, J. E. Banks, L. K. Dick, and J. D. Stark. Estimation of dynamic rate parameters in insect populations undergoing sublethal exposure to pesticides. Bull. Math. Biol., 69:2139–2180, 2007.
 H. T. Banks, K. L. Sutton, W. C. Thompson, G. Bocharov, M. Doumic, T. Schenkel, J. Arguilaguet, S. Giest, C. Peligero, and A. Meyerhans. A new model for the estimation of cell proliferation dynamics using cfse data. J. of Immunological Methods, To appear, 2011.
 F. Bekkal Brikci, J. Clairambault, and B. Perthame. Analysis of a molecular structured population model with possible polynomial growth for the cell division cycle. Math. Comput. Modelling, 47(78):699–713, 2008.
 F. Bekkal Brikci, J. Clairambault, B. Ribba, and B. Perthame. An ageandcyclinstructured cell population model for healthy and tumoral tissues. J. Math. Biol., 57(1):91–110, 2008.
 S. Bernard, L. PujoMenjouet, and M. C. Mackey. Analysis of cell kinetics using a cell division marker: Mathematical modeling of experimental data. Biophys J., 84(5):3414–3424, May 2003.
 F. Billy, J. Clairambault, O. Fercoq, S. Gaubert, T. Lepoutre, T. Ouillon, and S. Saito. Synchronisation and control of proliferation in cycling cell population models with age structure. Math. Comp. Simul., Accepted, 2011.
 A. Calsina and J. Saldaña. A model of physiologically structured population dynamics with a nonlinear individual growth rate. J. Math. Biol., 33:335–364, 1995.
 P. Clément, H. J. A. M. Heijmans, S. Angenent, C. J. van Duijn, and B. de Pagter. Oneparameter semigroups, volume 5 of CWI Monographs. NorthHolland Publishing Co., Amsterdam, 1987.
 O. Diekmann, H. J. A. M. Heijmans, and H. R. Thieme. On the stability of the cell size distribution. J. Math. Biol., 19:227–248, 1984.
 J. Douglas and F. Milner. Numerical methods for a model of population dynamics. Calcolo, 24:247–254, 1987.
 M. Doumic, P. Maia, and J. Zubelli. On the calibration of a sizestructured population model from experimental data. Acta Biotheoretica, 58:405–413, 2010.
 M. Doumic, B. Perthame, and J. Zubelli. Numerical solution of an inverse problem in sizestructured population dynamics. Inverse Problems, 25(4):045008, 2009.
 M. R. Dowling, D. MilutinoviÄ, and P. D. Hodgkin. Modelling cell lifespan and proliferation: is likelihood to die or to divide independent of age? Journal of The Royal Society Interface, 2(5):517–526, 2005.
 K.J. Engel and R. Nagel. Oneparameter semigroups for linear evolution equations, volume 194 of Graduate Texts in Mathematics. SpringerVerlag, New York, 2000. With contributions by S. Brendle, M. Campiti, T. Hahn, G. Metafune, G. Nickel, D. Pallara, C. Perazzoli, A. Rhandi, S. Romanelli and R. Schnaubelt.
 H. Engl, W. Rundell, and O. Scherzer. A regularization scheme for an inverse problem in agestructured populations. J. Math. Anal. Appl., 182(3):658 – 679, 1994.
 P. Gabriel. The shape of the polymerization rate in the prion equation. Math. Comput. Modelling, 53(78):1451–1456, 2011.
 P. Gabriel and L. M. Tine. Highorder WENO scheme for polymerizationtype equations. ESAIM Proc., 30:54–70, 2010.
 P. Gerlee and A. Anderson. A hybrid cellular automaton model of clonal evolution in cancer: The emergence of the glycolytic phenotype. Journal of Theoretical Biology, 250(4):705 – 722, 2008.
 A. Golubev. Exponentially modified Gaussian (EMG) relevance to distributions related to cell proliferation and differentiation. J. Theor. Biol., 262(2):257–266, 2010.
 M. Gyllenberg, A. Osipov, and L. Päivärinta. The inverse problem of linear agestructured population dynamics. J. Evol. Equ., 2(2):223–239, 2002.
 M. Gyllenberg and G. F. Webb. Agesize structure in populations with quiescence. Math. Biosci., 86(1):67–95, 1987.
 M. Gyllenberg and G. F. Webb. A nonlinear structured population model of tumor growth with quiescence. J. Math. Biol., 28(6):671–694, 1990.
 M. Gyllenberg and G. F. Webb. Quiescence in structured population dynamics: applications to tumor growth. In Math. pop. dyn. (New Brunswick, NJ, 1989), volume 131 of Lecture Notes in Pure and Appl. Math., pages 45–62. Dekker, New York, 1991.
 P. Hinow, S. Wang, C. Arteaga, and G. Webb. A mathematical model separates quantitatively the cytostatic and cytotoxic effects of a HER2 tyrosine kinase inhibitor. Theo. Biol. Med. Modelling, 4(1):14, 2007.
 T. Kostova. An explicit thirdorder numerical method for sizestructured population equations. Numer. Methods Partial Diff. Equ., 19(1):1–21, 2003.
 H. Lee and A. Perelson. Modeling t cell proliferation and death in vitro based on labeling data: Generalizations of the smithâmartin cell cycle model. Bull. Math. Biol., 70:21–44, 2008.
 D. Liu, D. M. Umbach, S. D. Peddada, L. Li, P. W. Crockett, and C. R. Weinberg. A randomperiods model for expression of cellcycle genes. Proceedings of the National Academy of Sciences of the United States of America, 101(19):7240–7245, 2004.
 D. Mallet and L. D. Pillis. A cellular automata model of tumorâimmune system interactions. Journal of Theoretical Biology, 239(3):334 – 350, 2006.
 J. A. J. Metz and O. Diekmann, editors. The dynamics of physiologically structured populations, volume 68 of Lecture Notes in Biomathematics. SpringerVerlag, Berlin, 1986. Papers from the colloquium held in Amsterdam, 1983.
 P. Michel, S. Mischler, and B. Perthame. General relative entropy inequality: an illustration on growth models. J. Math. Pures Appl., 84(9):1235–1260, 2005.
 A. A. Patel, E. T. Gawlinski, S. K. Lemieux, and R. A. Gatenby. A cellular automaton model of early tumor growth and invasion: The effects of native tissue vascularity and increased anaerobic tumor metabolism. Journal of Theoretical Biology, 213(3):315 – 331, 2001.
 B. Perthame. Transport equations in biology. Frontiers in Mathematics. Birkhäuser Verlag, Basel, 2007.
 B. Perthame and J. Zubelli. On the inverse problem for a sizestructured population model. Inverse Problems, 23(3):1037–1052, 2007.
 D. Picart and B. Ainseba. Parameter identification in multistage population dynamics model. Nonlin. Anal.: Real World Appl., 12(6):3315 – 3328, 2011.
 D. Picart, B. Ainseba, and F. Milner. Optimal control problem on insect pest populations. Appl. Math. Lett., 24(7):1160–1164, 2011.
 M. Pilant and W. Rundell. Determining a coefficient in a firstorder hyperbolic equation. SIAM J. Appl. Math., 51(2):pp. 494–506, 1991.
 M. Pilant and W. Rundell. Determining the initial age distribution for an age structured population. Math. Pop. Studies, 3(1):3–20, 1991.
 G. G. Powathil, K. E. Gordon, L. A. Hill, and M. A. Chaplain. Modelling the effects of cellcycle heterogeneity on the response of a solid tumour to chemotherapy: Biological insights from a hybrid multiscale cellular automaton model. Journal of Theoretical Biology, 308(0):1 – 19, 2012.
 V. Quaranta, K. A. Rejniak, P. Gerlee, and A. R. Anderson. Invasion emerges from cancer cell adaptation to competitive microenvironments: Quantitative predictions from multiscale mathematical models. Seminars in Cancer Biology, 18(5):338 – 348, 2008.
 V. Quaranta, A. M. Weaver, P. T. Cummings, and A. R. Anderson. Mathematical modeling of cancer: The future of prognosis and treatment. Clinica Chimica Acta, 357(2):173 – 179, 2005.
 I. RamisConde, D. Drasdo, A. R. Anderson, and M. A. Chaplain. Modeling the influence of the ecadherinÎ²catenin pathway in cancer cell invasion: A multiscale approach. Biophysical Journal, 95(1):155 – 165, 2008.
 B. Ribba, T. Alarcón, K. Marron, P. Maini, and Z. Agur. The use of hybrid cellular automaton models for improving cancer therapy. In P. Sloot, B. Chopard, and A. Hoekstra, editors, Cellular Automata, volume 3305 of Lecture Notes in Computer Science, pages 444–453. Springer Berlin / Heidelberg, 2004.
 W. Rundell. Determining the birth function for an age structured population. Math. Pop. Studies, 1(4):377–395, 1989.
 W. Rundell. Determining the death rate for an agestructured population from census data. SIAM J. Appl. Math., 53(6):pp. 1731–1746, 1993.
 D. R. Tyson, S. P. Garbett, P. L. Frick, and V. Quaranta. A method to quantify populationlevel cell proliferation dynamics from singlecell data. Nature Methods, in press.
 G. F. Webb. Theory of nonlinear agedependent population dynamics, volume 89 of Monographs and Textbooks in Pure and Applied Mathematics. Marcel Dekker Inc., New York, 1985.
 G. F. Webb. An operatortheoretic formulation of asynchronous exponential growth. Trans. Amer. Math. Soc., 303(2):751–763, 1987.