Model-Based Analysis of Cell Cycle Responses to Dynamically Changing Environments

Cell cycle progression is carefully coordinated with a cell’s intra- and extracellular environment. While some pathways have been identified that communicate information from the environment to the cell cycle, a systematic understanding of how this information is dynamically processed is lacking. We address this by performing dynamic sensitivity analysis of three mathematical models of the cell cycle in Saccharomyces cerevisiae. We demonstrate that these models make broadly consistent qualitative predictions about cell cycle progression under dynamically changing conditions. For example, it is shown that the models predict anticorrelated changes in cell size and cell cycle duration under different environments independently of the growth rate. This prediction is validated by comparison to available literature data. Other consistent patterns emerge, such as widespread nonmonotonic changes in cell size down generations in response to parameter changes. We extend our analysis by investigating glucose signalling to the cell cycle, showing that known regulation of Cln3 translation and Cln1,2 transcription by glucose is sufficient to explain the experimentally observed changes in cell cycle dynamics at different glucose concentrations. Together, these results provide a framework for understanding the complex responses the cell cycle is capable of producing in response to dynamic environments.

T div,1 T div,0 T div,2 T div,3 V div,0 V div,1 V div,2 V div,3 A B Figure S1: Simulating the dynamic response of the cell cycle of a step-change in parameters. The parameter k s,bS in the Barik model undergoes a step-change (A), which changes the volume at division (V div ) and cell cycle duration (T div ) in subsequent generations (B).

Pfeuty
Chen Barik Figure S3: Model sensitivity to changes in the growth rate parameter, µ. In all three models, increasing growth rate leads to larger daughter cells and a reduced duration of G1 (upper panels), consistent with experimental observations. In addition, the models consistently predict that increasing growth rate will monotonically increase the size of daughter cells in subsequent generations until they reach their final size, irrespective of the time at which growth rate is increased (lower panels), with some minor deviations. (mins) Figure S4: Correlation between V dau and T G1 for a range of mutants. Data from [1], after filtering out mutants which displayed changes in growth rate (as described in [1]).  Figure S5: Mode-locking of the Pfeuty model to periodic forcing. (A) shows the predicted change in cell cycle phase in response to a 2.15 minute, 10% increase in the parameter s x,2 (denoted ∆φ(t), red line). This perturbation is repeatedly applied at a period 5 minutes less than the unforced cell cycle period (dashed line). The predicted phase of entrainment is given by the point of intersection of these lines where ∆φ (t) < 0 (red circle). (B) shows the simulated results evaluating the prediction made in (A). The shaded area represents the time at which the perturbation is applied. The dashed vertical line represents the prediction made in (A). The stability of the mode-locking is demonstrated by the consistent phase relationship between the perturbation and the timing of cell division. (min) Figure S6: Linear fit of V dau and T G1 at constant growth rate and different glucose levels. Data from [2]. This quantifies the negative correlation between V dau and T G1 observed as glucose levels change.  Figure S7: Examples of how different parameter combinations can produce the same eventual change in behaviour, but with different dynamic responses. Responses to changes in two pairs of parameters are analysed: k dcmp and k smbM (A, B, C); and k dcm and k sn3 (D, E, F). In each case, the different responses of the individual parameters are combined to give the same eventual change in V div and T div (A, D). However, each case has a distinct dynamic response (compare B,C to E,F).

Model equations
All three models share the basic pattern of growth, budding and division, with cell size (V ) growing exponentially, as given by: The link between growth in cell volume and and the macroscopic observables characterising the cell cycle is then given by (as in the Main Text): This description is coupled to models of the molecular mechanisms of cell cycle progression through thresholds in the concentrations of key components. In each model, a component controls the initiation of S-phase and budding, and this event occurs when the concentration of this component increases through a given threshold. In the case of the Barik model, this threshold is as specified in the original paper [3], while for the Pfeuty and Chen models a threshold was chosen such that the daughter cell volume is an appropriate fraction of the mother cell volume. Similarly, each model includes a component controls the initiation of cytokinesis, and this event occurs when the concentration of this component decreases through a given threshold. All three models already specify such a threshold, and we use these as specified in the original papers. All three models were run with the same basal growth rate (µ = 0.007 min −1 ). Parameters determining how cell size (V ) modulates cell cycle dynamics were rescaled in the Pfeuty and Chen models to make V div the same for all three models (and therefore in common units, of fL), with details described below.

Pfeuty model
The Pfeuty model was described in [4], and provides a minimal pseud-biochemical description of cell cycle dynamics. The model equations are: In this model, budding occurs when [X] increases above 0.

Division occurs when [Y ]
decreases below 2.
All parameters were rescaled by dividing by a factor of 2.86 to maintain the dynamics of the model at the common growth rate (µ = 0.007 min −1 ). We note that this does not change the behaviour of the model, but merely amounts to a rescaling of the time coordinate.
The parameter s x,2 was additionally rescaled by dividing by a factor of 20 to put cell size in common units with the Barik model.

Chen model
The Chen model was originally described in [5]. A modified version of this was analysed in [6], and it is this version that we use here. The model equations are: Goldbeter Koshland function: In this model, budding occurs when [Cln2] increases above 0.03. Division occurs when [Clb2] T decreases below 0.05.
The parameters k s,n2 , k s,n2 , k s,b2 , k s,b2 , k s,b5 , k s,b5 , [Bck] 0 , and D n3 were rescaled by dividing by a factor of 33.3 to put cell size in common units with the Barik model.

Barik model
The Barik model was described in [3], and provides a detailed description of cell cycle dynamics in the form of mass-action kinetics. The model equations are: In this model, budding occurs when [ClbS] increases above 37.5nM. Division occurs when [ClbM ] decreases below 12.5nM.

Model parameters of particular interest
While there is not, in general, a straightforward mapping between parameters in different models, there are sets of analogous parameters that represent similar molecular mechanisms.
These include parameters involved in cyclin synthesis and degradation (through the APC), as discussed in the Main Text. Here, we detail which sets of parameters we consider to represent these processes in each model.

Cln3 synthesis
The synthesis of the G1 cyclin Cln3 is represented in both the Chen and Barik models, and is represented by the parameters D n3 and k s,n3 , respectively. The component X in the Pfeuty model plays the role of initiating Start in this model, and its rate of synthesis is controlled by the parameter s x,2 .

Mitotic cyclin synthesis
The synthesis of mitotic cyclin is represented in both the Chen and Barik models, and is represented by the parameters k s,b2 and k s,bM , respectively.

APC synthesis
The APC subunit Cdc20 is represented in the Chen model, and its synthesis is represented by the parameter k s,20 . The APC subunit Cdh1 is represented in the Barik model, and its synthesis is represented by the parameter k s,h1 . Having simulated cell cycle behaviour with the basal and perturbed sets of parameters, it is possible to calculate the sensitivity of the observables to changes in the parameter in question according to: Where ∆Q = Q perturbed − Q basal and ∆k i /k i = 0.001. For a given parameter perturbation, this sensitivity gives a first-order estimate of the change in behaviour: An analogous calculation was performed for calculating the dynamic sensitivities, S Q i k , by performing step changes in parameters at different times during the cell cycle and tracking changes in behaviour down subsequent generations. An illustration of a similar simulation (though with larger changes in parameters) is shown in Supplementary Figure S1.

Characterising monotonic and nonmonotonic sensitivity responses
In the Main Text, the pervasiveness of nonmonotonic changes in T G1 down generations in response to a step-change in parameters is discussed (Figure 4). For a perturbation of a parameter k at time t, the sensitivity of T G1 down generations is characterised by the sequence: In order to quantify montonicity, the sequence of S changes in behaviour are given by (following [7,8]): Therefore, for two parameters with linearly independent sensitivity vectors ( the parameter perturbation required to obtain the eventual change in behaviour (δV div , δT div ) T is given by: The significance of this is illustrated in Supplementary Figure S7.
From the preceding analysis, it is clear that the cell cycle models considered are capable of a multitude of different dynamic responses. One straightforward way of simplifying and comparing these results is to consider what the average response is to a particular perturbation, and how variable that response is. This is sensible in the context of perturbations of populations of cells, where the timing of an external perturbation relative to the cell cycle phase of a particular cell is essentially random. The average change in daughter cell size in generation i is denoted byδ V dau,i and given by: The distribution of daughter cell sizes around this average is characterised by the standard deviation, given by σ V dau,i : Under constant conditions, the volumes of a cell at budding and division are given by: These are related to the volume of the daughter cell at birth: Substituting, we get: Differentiating with respect to a generic parameter, k, fixing dµ/dk = 0, and noting Noting that the fraction of mother cell volume donated to the daughter cell at division, f , is given by: we obtain:

Phase responses and fraction of mass donated to daughter cells
In this section, we use the phenomenological model to relate changes in cell cycle phase (i.e. changes in the timing of cell cycle events relative to a reference case) to changes in cell volume at budding and division at constant growth rate.
We consider a cell cycle which initially has a period T 0 , a size at division V div,0 , and a daughter size V dau,0 , meaning the resulting fraction of volume given to the daughter cell is f 0 = V dau,0 /V div,0 . Two perturbations are applied which result in changes in these characteristics down generations. The first set of characteristics are labelled T 1 , T 2 , T 3 , ..., V div,1 , V div,2 , V div,3 , ..., V dau,1 , V dau,2 , V dau,3 , ..., giving daughter fractions f 1 , f 2 , f 3 , ..., while the second are labelled τ 1 , τ 2 , τ 3 , ..., ν div,1 , ν div,2 , ν div,3 , ..., and ν dau,1 , ν dau,2 , ν dau,3 , ..., giving daughter fractions ρ 1 , ρ 2 , ρ 3 , .... Initially the cells are equally sized: For subsequent generations, we have: Thus after n generations, using the equality from Equation 19: Taking logarithms and rearranging: Since the perturbations applied are, by assumption, temporary, we can take the limit as n → ∞, so that V dau,n → ν dau,n , and as defined in Equation 12 (Main Text). This gives: This converges, since f n → ρ n as n → ∞. In practice, for parameter changes that last less than a cell cycle period, this limit converges very rapidly (within a few generations). When the comparison case is simply the eventual behaviour of the cell cycle under the basal parameter set, we have ρ i = f 0 , and f i = f 0 + ∆f i , giving Equation 14 (Main Text), i.e.: This demonstrates that, within the framework of this phenomenological model, there is a strict relationship between changes in the timing of cell cycle events and the fraction of mass donated to the daughter cell upon division.

Inferring the correspondence between cell cycle characteristics of populations and single cells
In [9], cells were grown under 6 different nutrient limitations, each at 6 different growth rates, and cell cycle characteristics were measured, and presented in the form of population averages. These include the fraction of cells in G1 (denoted F G1 ), and the average cell size in the population (denotedV ). As discussed in the Main Text, the models consistently predicted that, at constant growth rate, changes in V dau and T G1 should be negatively correlated with one another. This prediction was validated under a range of nutrient and genetic perturbations using data from [2], as presented in Figure 2. However, in order to use averaged population data from [9] to provide further validation, a correspondence between the population properties (F G1 andV ) and single-cell properties (T G1 and V dau ) must be found. We begin by relating F G1 to T G1 . F G1 is given by: Here, T G1,dau , T div,dau , T G1,moth , and T div,moth denote the durations of the G1 phases, and the durations of the entire cell cycle, in daughter and mother cells, respectively. In addition, F dau and F moth denote the fractions of daughter and mother cells in the population, respectively (F dau + F moth = 1). Noting that mother cells spend little time in G1, we can approximate F G1 as: This demonstrates the expected proportional relationship between T G1,dau and F G1 (as has been noted previously [10]).
In the case of theV , we have:V The average daughter cell size is given by: For a mother cell, we assume that they are born at size V bud and produce daughter cells at intervals of T S/G2/M (i.e. their G1 phases are of negligible duration). Furthermore, we approximate these values by their values in daughter cells. Thus, the averaged mother cell size is given by: