|
|
|
|
Vol. 15, Issue 8, 3841-3862, August 2004
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||





* Department of Biology, Virginia Polytechnic Institute and State University, Blacksburg, Virginia 24061-0406;
Molecular Network Dynamics Research Group of the Hungarian Academy of Sciences and Department of Agricultural and Chemical Technology, Budapest University of Technology and Economics, H-1521 Budapest, Hungary; and
The Rockefeller University, New York, New York 10021
Submitted November 7, 2003;
Revised May 3, 2004;
Accepted May 15, 2004
Monitoring Editor: Thomas Pollard
| ABSTRACT |
|---|
|
|
|---|
| INTRODUCTION |
|---|
|
|
|---|
The cell cycle regulatory system is most fully worked out for budding yeast (Saccharomyces cerevisiae). A hypothetical molecular mechanism for regulating DNA synthesis, bud emergence, mitosis, and cell division in budding yeast is proposed in Figure 1. How does one determine whether such a hypothesis is correct? The classical approach in physical chemistry is to convert the mechanism into a set of kinetic equations (nonlinear ordinary differential equations) and to compare the solutions of these equations to the observed behavior of the chemical reaction system. If a set of rate constants can be found for which the solutions fit the observations, then the mechanism is provisionally confirmed (pending further experimental investigations). If not, inconsistencies identify aspects of the mechanism that require revision and further testing. Although a mechanism can be disproved if it is inconsistent with well-established facts, it can never be proved correct, because new observations may force modifications and additions. Hence, our intention is not to prove that the hypothesis in Figure 1 is "true" but rather only to show that the mechanism is a reasonable approximation to what is going on inside yeast cells.
|
Physical chemists have used this methodology successfully to unravel molecular mechanisms that are not only complicated (many components) but also dynamically complex, involving multiple steady states, oscillations, and chaos (Field et al., 1972
; Gyorgyi and Field, 1992
). We have used this approach to study simpler models of cell cycle regulation in frog eggs (Novak and Tyson, 1993
), fission yeast (Novak et al., 2001
), and budding yeast (Chen et al., 2000
). The model in Chen et al. (2000
) gives an adequate description of the G1-to-S transition, but, since it was published, many more molecular details about the M-to-G1 transition ("exit from mitosis") have come to light. By adding these details to the mechanism in Chen et al. (2000
), we are able to model comprehensively and quantitatively all stages of the chromosome replication-segregation cycle in a eukaryote, based on a mechanism that is almost fully specified at the genetic level.
As explained in MATERIALS AND METHODS, we turn the mechanistic hypothesis (Figure 1) into a mathematical model (Table 1) with a preliminary set of kinetic constants (Table 2). Next, we solve these equations numerically and show, in RESULTS, that the solutions are in accord with the physiological properties of 120 mutant strains of budding yeast out of 131 studied so far (Table 3). There are 11 mutants that the model fails to account for, and these failures identify aspects of the mechanism that need further investigation. We also show how to use the model to interpret existing data, to design new experiments, and to think about the "molecular logic" of cell cycle regulation.
|
|
|
| MATERIALS AND METHODS |
|---|
|
|
|---|
A wiring diagram is a set of boxes (components) interconnected by arrows (reactions). An instantaneous state of the system is a specification of the current concentrations of all its components. Given a state of the system, the chemical reactions (synthesis, degradation, activation, inhibition, binding, and release) indicate how the state will change in the next moment of time. Each reaction proceeds at a rate determined by the state of the system and by kinetic parameters (e.g., rate constants and binding constants). By applying the general principles of biochemical kinetics, we can convert the mechanism in Figure 1 into a set of differential and algebraic equations (Table 1) that determine how the state of the control system evolves in time.
It is important to realize that there is no unique correspondence between a wiring diagram and a set of mathematical equations; the same mechanism can be represented by different forms of equations. For example, the modeler may choose to describe some components by differential equations and others by algebraic equations, and to describe the rates of some reactions by Michaelis-Menten kinetics and others by mass-action kinetics. These choices are somewhat arbitrary, depending on the level of detail desired in the model. Hence, there is a hierarchy of assumptions that go into a model: first, the wiring diagram, which is a hypothesis about how the components of the control system are interconnected; second, the mathematical equations, which are representations of a possible set of kinetic consequences of the mechanism; and third, the assignment of specific values to the rate constants in the equations. Once these choices are made, the equations can be solved numerically and the dynamic behavior of the control system compared with the observed properties of dividing yeast cells. If there are inconsistencies between the behavior of the model and the behavior of the cells, then one usually works backwards through the hierarchy to resolve the discrepancies. First, rate constants are adjusted to try to fit the model to the data. If that proves impossible, then one reconsiders the assumptions made in translating the diagram into equations. If modifications of those assumptions do not help, then one must reconsider the wiring diagram itself. Perhaps there are missing components or interactions that prevent the model from fitting the data; for example, see von Dassow et al. (2000
). By this process of trial, error, and refinement, we have settled temporarily on the wiring diagram (Figure 1), equations (Table 1), and parameter values (Table 2) presented here as representative of the budding yeast cell cycle engine. The model will certainly continue to evolve as it is confronted with new experimental studies of the control system. In our experience, during this evolutionary process, successive models show steady improvements over earlier versions rather than wholesale replacement of one set of equations and parameters with another. We suggest that the modeling process is converging on ever better mathematical approximations of the molecular regulatory system. In the next subsection, we describe in more detail how we choose rate laws and assign values to kinetic constants.
Rate Equations and Kinetic Constants
Although the model as a whole (Figure 1) is complicated, each individual equation (Table 1) is simple. For example, Cln2-dependent kinase activity changes in time due to synthesis and degradation of the Cln2 subunit (assuming Cdc28 is in excess and binds rapidly to cyclins). The rate of Cln2 synthesis, k's,n2 + k''s,n2 · [SBF], includes a constitutive rate, k's,n2 (for simulation of GAL-CLN2 mutants), and a contribution dependent on the relative activity of the transcription factor SBF, controlling expression of the CLN2 gene. The activity of SBF is given by the algebraic equation [SBF] = G(Va,sbf, Vi,sbf, Ja,sbf, Ji,sbf). The "Goldbeter-Koshland function" G(Va, Vi, Ja, Ji) varies sigmoidally between 0 and 1, increasing as a function of Va and decreasing as a function of Vi (Goldbeter and Koshland, 1981
). The sigmoid is steeper (more switch-like) as Ja and Ji are much <1. In the equation for SBF activity, Va,sbf = ka,sbf · (
sbf,n3 · ([Cln3] + [Bck2]) +
sbf,5 · [Clb5]) indicates that SBF is activated by Cln2/Cdc28, Cln3/Cdc28, Bck2, and Clb5/Cdc28, whereas Vi,sbf = k'i,sbf + k''i,sbf · [Clb2] indicates that SBF is inactivated by Clb2/Cdc28 (Amon et al., 1993
). The coefficients
sbf express the relative efficiencies of these components in activating SBF. Our choices for the
values (Table 2) reflect experimental data that Cln3/Cdc28 and Bck2 are about equally efficient in activating SBF because cell sizes are about the same in cln3
and bck2
mutants (Epstein and Cross, 1994
), whereas Cln2/Cdc28 is much less efficient (Dirick et al., 1995
; Stuart and Wittenberg, 1995
). The same experiments indicate that SBF is turned on very abruptly as the cell reaches a critical size, so we choose Ja,sbf and Ji,sbf to be much <1. We have assigned a value of 2 to
sbf,b5, because there is evidence that SBF and MBF are activated by Clb5 and 6 (Schwob and Nasmyth, 1993
). To make cln3
bck2
not rescued by sic1
or GAL-CLB5 (Wijnen and Futcher, 1999
), we set
sbf,b5 < 4.
In our model, SBF turns on in wild-type cells when Cln3 and Bck2 have accumulated in the nucleus to a critical level: [Cln3] + [Bck3] = k'i,sbf/(ka,sbf ·
sbf,n3). Hence, the ratio k'i,sbf/ka,sbf sensitively determines cell size when SBF turns on (Dirick et al., 1995
), and the ratio ka,sbf/k''i,sbf (which sets the Clb2-kinase activity required to inactivate SBF) sensitively determines the timing of Cln2 disappearance in the latter part of the cycle (Amon et al., 1993
). The values given to these parameter ratios in Table 2 were chosen to fit the model to these observations.
Cln2 degradation is assumed to be a simple first-order process, kd,n2 · [Cln2] with half-life = 0.693/kd,n2 = 6 min, in reasonable accord with experiment (Salama et al., 1994
; Barral et al., 1995
). Actually, the Clns are degraded by the SCF-dependent proteasome pathway (Deshaies et al., 1995
; Lanker et al., 1996
), whose first step is phosphorylation of the protein to be degraded. We are assuming (in this version of the model) that the phosphorylation of Cln2, required for recognition by the SCF, happens uniformly throughout the cell cycle.
To fit the model to the total amount of Cln2 observed in an asynchronous culture of budding yeast cells (Cross et al., 2002
), we choose k''s,n2 = 0.15 min-1.
The differential equations for [Clb5] and [Clb2] have terms for synthesis and degradation and for binding to and release from stoichiometric inhibitors Sic1 and Cdc6. The synthesis term has a transcription factor given by a Goldbeter-Koshland function, as for [Cln2]. The degradation term consists of several components describing relative contributions from different pathways: Cdc20/APC and Cdh1/APC for Clb2 (Zachariae and Nasmyth, 1999
; Yeong et al., 2000
), and Cdc20/APC-dependent and -independent pathways for Clb5 (Irniger and Nasmyth, 1997
; Wasch and Cross, 2002
). The rate constants for these pathways are estimated from the stabilities of the Clbs in various phases of the cycle when the different pathways are active. Then, the rate constants for synthesis are estimated from the data of Cross (2002
) on average cyclin contents in an asynchronous culture.
By the same reasoning, one can work through every equation in the model, identifying the molecular steps involved, the assumptions made in defining the rate law, and (where relevant data are available) the numerical values assigned to rate constants. In some cases, directly relevant data are not available; for example, for the binding of Sic1 and Cdc6 to Clb2 and Clb5. Typically, we assume that binding is rapid (the rate constant for association, kas is large compared with rate constants for e.g., synthesis, degradation, and phosphorylation) and the equilibrium binding constant (kas/kdi) is large. Precise values for these parameters are not crucial to the model, with one exception. For reasons indicated in RESULTS, we assume that Cdc6 is an ineffective stoichiometric inhibitor of Clb5 kinase. This assumption is supported by the observation of Archambault et al. (2003
) that Cdc6 does not complex with Clb5 in immunoprecipitation assays.
The rate constants we propose are "effective" values because they capture only the time scales of processes (note that every rate constant has a unit of minute-1). They do not capture the characteristic concentration scales of the rates. When the model was being developed, we did not know the absolute concentrations of most of the proteins in the mechanism, so we expressed the concentration of each protein in an arbitrary unit (au), which may be different for each protein. For example, for Tem1 and Cdc15, the arbitrary units are chosen so that [Tem1]T = [Cdc15]T = 1. Hence, the variables [Tem1] and [Cdc15] represent the fraction of each protein pool in the "active" form. Recently, Ghaemmaghami et al. (2003
) published a comprehensive analysis of protein expression in budding yeast. Their data will allow us to assign nanomolar values to each au in a later version of our model.
Our choices of concentration units are constrained by two facts. 1) For any two proteins involved in a stoichiometric interaction (e.g., Cdc14 and Net1, Clb2 and Sic1, Pds1 and Esp1), we must use the same au for each protein of a pair. 2) Cross et al. (2002
) have measured the concentrations of Clns, Clbs, Sic1, and Cdc28 in an asynchronous culture by tagging them with the same epitope. These data allow us to assign a concentration (in nanomolar) to the au used to express the concentrations of all these proteins. For example, Cross determined that there are 3000 molecules of Cln1 and Cln2 per diploid cell (averaged over the cell cycle). Given the diploid cell volume to be 100 fL, 3000 molecules corresponds to a concentration of 50 nM. (The volume of a diploid cell is calculated from the volume of a haploid cell, which was measured to be 50 fL by Nash et al., 1988
and Yaglom et al., 1995
. A diploid cell is twice as big as a haploid cell (Lorincz and Carter, 1979
) and presumably contains twice the amount of every protein.) From our model, the average value of [Cln2] over a full cycle is 1.25 au; therefore, 1 au of Cln2 (or Clb2 or Sic1) corresponds to a concentration of
40 nM, or 1200 molecules per haploid cell.
We illustrate how to relate rate constants in the model to measurements in a cell culture with a few examples. The rate constant for synthesis of Sic1 (when Swi5 is fully active) is k''s,c1 = 0.12 au/min = 4.8 mM/min. On the other hand, the rate constant for synthesis of Cdc20 (when Mcm1 is fully active), k''s,20 = 0.6 au/min, cannot be expressed in nanomolar per minute until we know the average concentration of Cdc20 in an asynchronous culture, which should correspond to 0.7 au in the model (the average value of [Cdc20]T over a full cycle). The equilibrium binding constant for Clb2 and Sic1 is kas,b2/kdi,b2 = 103 au-1 = 25 nM-1, whereas the binding constant for Cdc14 and Net1, kas,rent/kdi,rent = 200 au-1, cannot be expressed in nanomolar-1 until we know the concentration of Net1 or Cdc14 in a cell. Because Net1 and Cdc14 are expressed in the same au, the ratio [Net1]T/[Cdc14]T = 1.4 is a meaningful prediction of the model: that Net1 is in slight excess (40%) over Cdc14. On the other hand, because Tem1 and Cdc15 are expressed in different au, the ratio [Tem1]T/[Cdc15]T = 1 tells us nothing about the relative amount of these proteins in a cell.
Timing of Cell Cycle Events
Much of the quantitative data available to us refers to the timing of bud emergence, onset of DNA synthesis, and cell separation. To link the output of the model (the temporal evolution of e.g., [Cln2], [Clb5], [Clb2], [Cdc14]) to these events, we use auxiliary variables called [BUD], [ORI], and [SPN].
[BUD] represents proteins that are phosphorylated by Cdc28/cyclin dimers and subsequently initiate a new bud when the phosphorylation state reaches a threshold, [BUD] = 1. The rate of phosphorylation, ks,bud(
bud,n2[Cln2] +
bud,n3[Cln3] +
bud,b5[Clb5]), reflects the fact that bud emergence can be driven by any of the Clns or by Clb5/6 (Cvrckova and Nasmyth, 1993
). The values assigned to the
bring the timing of bud emergence in line with observations in mutants lacking various combinations of these cyclins.
In a similar manner, [ORI] = 1 signals the onset of DNA synthesis. At that time, we invoke the DNA replication/spindle assembly checkpoint by activating Mad2 and Bub2, thereby preventing activation of Cdc20 and Tem1, and keeping cells from mitotic exit.
The checkpoint is lifted when [SPN] = 1, which represents alignment of all chromosomes on the metaphase plate. We assume that exit from mitosis (telophase, cell separation) requires the elimination of Clb2-dependent kinase activity (Zachariae and Nasmyth, 1999
) below a threshold value (Kez = 0.3)
At exit from mitosis, [BUD] and [SPN] are reset to zero. [ORI] is reset to zero only if [Clb2] + [Clb5] drops below another threshold (Kez2 = 0.2), which is our criterion for relicensing origins of replication (Li, 1995
; Su et al., 1995
).
To account for inviability of the triple-cln mutant (cln1
cln2
cln3
), to be described in RESULTS, we are led to assume that when cells grow too large they may not successfully initiate DNA synthesis, and they will be G1 arrested. In the model, we require that a cell must reach [ORI] = 1 before a wild-type cell in the same growth medium has divided twice; if not, the cell is considered G1 arrested, even if [ORI] = 1 sometime later.
Mass at Division
Another commonly observed property of cells is their size at division. An important variable in our model is [mass]. Typically, we record the value of [mass] at cell division in wild-type and mutant cells and compare it with observations such as "these mutant cells are roughly twice the size of wild-type cells." The differential equation for [mass] implies that cells grow exponentially, with mass doubling time (MDT) = 0.693/kg, where kg is the specific growth rate. (MDT = 90 min on glucose medium and 150 min on galactose medium.) When a cell exits mitosis, we divide [mass] asymmetrically between mother and daughter cells, according to a rule described in Chen et al. (2000
).
Notice that [mass] enters the dynamical system as a multiplier of the rates of synthesis of the cyclins. It is based on the assumption, which was made in Chen et al. (2000
) as well, that cyclins are synthesized at a rate proportional to the total number of ribosomes in a cell (which is proportional to cell size) and then accumulate in a constant-volume compartment of the cell (the nucleus). Hence, [Cln3], [Cln2], [Clb5], and [Clb2] represent the concentrations of the cyclins in the nucleus. This assumption (rate of cyclin synthesis proportional to mass) is the simplest mechanism we know to coordinate the growth cycle to the chromosomal replication cycle, so that average cell size is maintained generation after generation (Futcher, 1996
).
Because Sic1 and Cdc14 (Shou et al., 1999
; Edgington and Futcher, 2001
) are found in both nuclear and cytoplasmic compartments, their synthesis rates are not multiplied by [mass]. No additional assumptions are made on the nuclear localization of other components in the regulatory network (e.g., Swi5, Cdc20 are not assumed to be exclusively nuclear). As described in DISCUSSION, as more experimental data become available (Huh et al., 2003
), we will incorporate subcellular localization of proteins in a later version of the model.
Computer Simulation
Using a computer program (http://www.math.pitt.edu/~bard/bardware/binary/winpp.zip) freely available from G. Bard Ermentrout (Department of Mathematics, University of Pittsburgh, Pittsburgh, PA), we solve the equations in Table 1, given the parameter values in Table 2 (which are appropriate for a wild-type cell), for the explicit time dependence of each variable ([Cln2](t), [Clb2](t),...). To compute a solution numerically, we also must assign initial conditions (at t = 0) to all the variables ([Cln2](0), [Clb2](0),...). Initial conditions should reasonably represent an experimental protocol, e.g., the values in Table 2 might represent conditions in a small wild-type cell selected from an elutriating rotor (newborn daughter cell with origin of replication relicensed), but their precise values are not important. The simulation was done for cells growing in culture with MDT = 90 min. Under the asymmetric-division rule of Chen et al. (2000
), we give 0.4587 of the mass-at-division to the daughter, and 0.5413 to the mother; hence, the cycle time for the daughter is 101.2 min and that for the mother is 80.0 min.
To simulate each mutant, we use exactly the same equations (Table 1) and parameter values (Table 2) except for those parameter changes that are dictated by the nature of the mutation. In the simplest case (gene X is deleted), the rate of synthesis of protein X is set to zero. Overexpression is more complicated, because a gene may be overexpressed in different ways: from a plasmid or from genomically integrated copies, with its own or a foreign promoter. In the case of multiple integrated copies under control of the natural promoter, we simply multiply the appropriate k's and k''s parameters by an integer, depending on the number of additional copies. If gene X is constitutively overexpressed (typically using the GAL promoter), then we increase k's,x (the constant rate of synthesis of protein X), and we change the specific growth rate kg to match the new growth medium. The new value assigned to k's,x is somewhat arbitrary, because the rates of expression of proteins from such constructs are unknown and vary for different proteins. For each GAL-construct, we choose a value of k's,x that gives results consistent with phenotypes of all mutants containing this particular GAL-construct (e.g., all mutants containing GAL-CLB2 are assigned the same value of k's,b2 but k's,b5 may be assigned a different value to model GAL-CLB5).
For temperature-sensitive (ts) mutants in Table 3, we assume that the relevant catalytic activity is normal at the permissive temperature and zero at the restrictive temperature. To model "synthetic lethality" of ts alleles (where the double mutant gene1ts gene2ts is inviable at a certain temperature originally permissive to the two single mutants), we assume that the relevant rate constants for gene 1 and gene 2 are reduced to a fraction of the wild-type values at that intermediate temperature.
Assay of Model Robustness
To determine model robustness, the model was transferred to MatLab version 6.5 by using the ode23s integrator. The following rules were then applied to establish whether a given parameter set yielded a "viable" solution. First, it was required that the model executes the following events in order: origin relicensing (due to a drop in [Clb2] + [Clb5] below Kez2); origin activation (due to a subsequent rise in [Clb2] + [Clb5], causing [ORI] to increase above 1); spindle alignment (due to a rise in [Clb2], causing [SPN] to increase above 1); Esp1 activation (the [Esp1] variable going through 0.1, because of Pds1 proteolysis); and finally [Clb2] dropping below a threshold Kez to trigger division. If these "events" did not occur in this order the model was considered inviable. If division occurred in an "unbudded cell" (i.e., when [BUD] < 1), this also was considered inviable. The model also was required to meet a rough criterion for a periodic solution, that the root mean square deviation of all variables was <0.05 at subsequent divisions. Last, if the cell [mass] was >10, this was considered to be an inviable model.
Given these rules, MatLab code was written to systematically vary each of the parameters in the model up to 256-fold up and down, with 1.414-fold increments, and the maximum tolerable variation was recorded for each parameter in each direction. Importantly, only single-parameter variations were tested. Whereas formally this procedure does not establish a "volume" in parameter space, it still gives a sense of the robustness of the model. An accurate and meaningful volume determination runs into difficulties that we have not yet solved; therefore, at present, we are relying on this empirical test.
A similar test was run on several mutants: APC-A (k''a,20 = 0; or equivalently ka,apc = 0); cdh1
(k''d,b2 = 0.001, k''d3,pds = 0.0001; for unknown reasons, the k''d,b2 = 0 model caused an error with the Matlab integrator and so the "trivial" value of 0.001 was used instead), cki
(e.g., all synthesis, association values for Sic1 and Cdc6 set to zero), cln3
, cln2
, clb5
(synthesis values for the relevant cyclins set to zero).
In a similar manner, the robustness of an inviable mutant (e.g., APC-A cdh1
) was assessed by searching for single-parameter increase or decrease that yielded a cycling model from an initially noncycling parameter set. Matlab code and complete data on parameter sensitivity runs are available on request (fcross{at}mail.rockefeller.edu).
| RESULTS |
|---|
|
|
|---|
|
Furthermore, the relative amounts of certain groups of proteins are in rough quantitative agreement with recent measurements by Cross et al. (2002
) and Archambault et al. (2003
). The ratios, for an asynchronous culture with MDT = 90 min, are as follows:
![]() |
The most significant discrepancy between the model and experiment is the ratio of [Sic1] to ([Clb1] + [Clb2] + [Clb5] + [Clb6], 1:11 in experiment and 1:3 in model. Cross's measurements indicate that rapidly growing cells have much less Sic1, on average, than would be expected from the model. By decreasing the rate of synthesis of Sic1, we could get better agreement with the experiments for wild-type cells. But the model must satisfy many other constraints implied by the phenotypes of various mutants. In particular, in the model the triple-cln deletion strain (cln1
cln2
cln3
) would become viable, contrary to observations, if the rate of synthesis of Sic1 were reduced by half. In our simulation of the triple-cln mutant, if Sic1 concentration is too low, then Clb5 activates its own synthesis via SBF/MBF, and the simulated cell is perfectly viable by all our criteria. We do not have any explanation for this discrepancy in Sic1 level between the model and the experiment.
The parameter values in Table 2 have been chosen taking into account the properties of wild-type and mutant cells, and they represent a compromise of many, often competing, observations. To this end we must settle on a "data set" of mutants (Table 3) that will serve as a testing ground for the sufficiency of the model. The parameter values proposed in Table 2 have been selected by a painstaking process of trial-and-error to provide a suitable fit to the full data set.
The Model Conforms to the Phenotypes of >100 Mutant Strains
The wiring diagram in Figure 1 has been composed from evidence provided by the phenotypes of dozens of budding yeast mutants that have been constructed and characterized by deleting or overexpressing each genetic component singly and in multiple combinations. It remains an informal "cartoon" of the molecular regulatory system until it is converted into a precise mathematical model and demonstrated to be consistent with most of the facts about budding yeast mitotic division. In Table 3, we list 131 mutants that have been used to test the model.
For each mutant simulation, we use exactly the same equations (Table 1) and parameter values (Table 2), and we are allowed to change only those parameters that are governed by the nature of the mutation, as described in MATERIALS AND METHODS. We compare the computed behavior of the model with the observed phenotype of the cells. For example, if the mutant is inviable, at what phase of the cell cycle is it blocked? If viable, in what subtle ways does it differ from wild-type: size at onset of DNA synthesis, size at bud emergence, size at division, and duration of G1 phase?
In most cases, 120 of 131 mutants in Table 3, the model agrees well with observations; the 11 exceptions are marked with asterisks in the table.
The Model Is Fully Described on a Web Page
At our Web site, http://mpf.biol.vt.edu, full details about the model and all simulations can be found. On this site, we summarize the basic experimental results on which the wiring diagram (Figure 1) is based. We present simulations of all mutants in Table 3, including the precise parameter values used in each case. We provide facilities for the user to repeat any simulations by using our parameter set or any new choice of parameter values. We also give instructions for how one may modify the wiring diagram, build a revised model, and run new simulations.
How the Model Represents APC Phosphorylation and the FEAR Pathway
In earlier models of cell cycle regulation in frog eggs (Novak and Tyson, 1993
) and fission yeast (Novak et al., 2001
), we proposed the existence of an intermediary enzyme (IE) whose purpose was to introduce a time delay between the activation of cyclin B-dependent kinase at the onset of M phase and the activation of Cdc20 (Fizzy in frog eggs and Slp1 in fission yeast) at the metaphase-to-anaphase transition. There are good reasons to suspect that IE is the APC core complex itself (Cross, 2003
). First, deletion of IE would cause cells to arrest in metaphase, as is the case for most components of the APC. Rudner and Murray (2000
) showed that three components of the APC core, Cdc16, Cdc23, and Cdc27, are phosphorylated by mitotic Cdc28 kinase and that only the phosphorylated forms associate with Cdc20 effectively. They also showed that APC-A mutant cells (where all the possible Cdc28-phosphorylation sites in Cdc16/23/27 are removed by substituting alanine for serine/threonine residues) are viable, with a delay in mitotic exit. This phenotype is consistent with the identification of IE with APC, if we assume that the unphosphorylated form of APC retains some intrinsic ability to activate Cdc20 (simulations not shown). Hence, in the present model, the phosphorylated, active form of IE in our earlier models (IEP) is replaced by APC-P, the phosphorylated state of APC that is necessary for full activity in conjunction with Cdc20. Notice that APC need not be phosphorylated to function in conjunction with Cdh1 (Kramer et al., 2000
; Rudner and Murray, 2000
). The present model is not fully correct in that IEP/APC-P is treated as if it were an enzymatic activator of Cdc20 instead of a binding partner. This problem will be corrected in future versions of the model, by revising the wiring diagram and the differential equations. These changes will undoubtedly require minor revisions of the rate constants, but we do not expect any significant changes in the conclusions reached in this article.
The sequential activation of Cdc20 and Cdh1 suggests that Cdc20/APC-P degrades, directly or indirectly, an inhibitor of Cdh1 (Morgan, 1999
). The inhibitor must act before Cdc14 release from the nucleolus. If the Cdc20-sensitive inhibitor were to act downstream of Cdc14 release, then in the cdc20
strain, Cdc14 would be released by MEN action as soon as the spindle assembly checkpoint is satisfied, which is contrary to observation (Shirayama et al., 1999
). Pds1 is an obvious candidate for this inhibitor: it is degraded by Cdc20 (Yamamoto et al., 1996b
), and overexpression of nondegradable Pds1 delays Clb2 degradation for several hours (Cohen-Fix and Koshland, 1999
). But how might Pds1 degradation relate to Cdc14 release? Cdc14 is released from the nucleolus in two stages (Visintin et al., 2003
): an early, transient stage, via the FEAR pathway (Uhlmann et al., 1999
; Stegmeier et al., 2002
; Yoshida et al., 2002
), involving Esp1 and Cdc5; and a late, sustained stage, involving Cdc15 and other MEN components. By binding and inhibiting Esp1, Pds1 inhibits the FEAR pathway, preventing Cdc14 release (FEAR stands for Cdc-fourteen early anaphase release).
At the time we were formulating this model, the information just cited was not known. In its place, we hypothesized a phosphatase (PPX) that inhibits Cdc14 release by opposing the action of Cdc15 on Net1. The role of Pds1 was to keep PPX level high by inhibiting its degradation by Cdc20/APC-P. Hence, in our model, Pds1 up-regulates PPX, an inhibitor of Cdc14 release. According to the FEAR story, Pds1 down-regulates activators of Cdc14 release. With appropriate choice of kinetic constants, the dynamic consequences of the two formulations should be nearly identical. In a future version of the model, we will replace PPX by a representation of the FEAR pathway.
The Model Is Based on Alternative Stable Steady States Corresponding to G1 and S/G2/M Phases of the Cell Cycle
Although the rigor and precision of the model are essential attributes in its favor, the sheer magnitude of information that comes out of the computer can overwhelm the user. To make sense of this information, we need intuitive ways to understand the model's behavioran intuition disciplined by precise numerical simulations of the equations. We rely on the scheme in Figure 3 (described briefly in Chen et al., 2000
). Fundamental to cell cycle control in budding yeast are the antagonistic relations between B-type cyclins (Clb16, in association with Cdc28), which promote DNA synthesis and mitosis, and G1-stabilizers (Cdh1, Sic1, and the cyclin-dependent kinase inhibitor [CKI]-role of Cdc6), which oppose cell proliferation by degrading Clbs and stoichiometrically inhibiting Clb/Cdc28 complexes (For the CKI-role of Cdc6, see Greenwood et al., 1998
and Calzada et al., 2001
). We shall refer to Sic1 and Cdc6 together as the CKIs.
|
Because Clb-dependent kinases can inactivate Cdh1 (for review, see Zachariae and Nasmyth, 1999
) and destabilize CKIs (Verma et al., 1997
), these two classes of proteins are mutual antagonists (Figure 3). The model is designed to have two coexisting, self-maintaining steady states: a G1 state, when Clbs are scarce and their antagonists (Cdh1 and CKIs) are abundant; and an S/G2/M state, when the reverse is true (Nasmyth, 1996
; Novak and Tyson, 1997
; Novak et al., 1998
; Tyson and Novak, 2001
; Tyson et al., 2001
). When yeast cells are proliferating, the control system is undergoing periodic transitions from the G1 state to the S/G2/M state and back again. These transitions (called START and FINISH) are irreversible and alternating. Once a cell has executed START, it does not normally slip back into G1 phase and does START again. Rather, it must execute a distinctly different transition (FINISH) to return to G1. Likewise, a cell that has executed FINISH does not slip back into mitosis and try to separate its chromosomes a second time. There are exceptions to these rules (endoreplication and meiosis), but they do not nullify the central role played by irreversible, alternating START and FINISH transitions in the cell cycle.
To a first approximation, we view the budding yeast cell cycle as an alternation between these two stable steady states generated by the antagonism between Clb kinases and G1-stabilizers. From the simulation of the wild-type cycle (Figure 2), one can see how the control mechanism shifts from one state to the other, and how the transitions are carried out.
The START transition is facilitated by Cln1,2/Cdc28 complexes, which can phosphorylate and inactivate Cdh1 and CKIs, but they themselves are unopposed by the G1-stabilizers (for reviews, see Schwob et al., 1994
and Peters, 1998
). This transition occurs when the cell has grown large enough to accumulate a critical concentration of Cln3-dependent kinase in the nucleus (Miller and Cross, 2000
; Cross et al., 2002
). Cln3 kinase and a back-up (Bck2) activate SBF and MBF, the transcription factors for Cln1,2 and Clb5,6, so their levels increase. Clb5,6-dependent kinases are inactive due to inhibition by the CKIs, but Cln1,2-dependent kinases are not so inhibited. Cln-dependent kinases depress Cdh1 and CKIs enough to allow the Clb-dependent kinases to assert themselves, switching the control system into the stable S/G2/M state. Once the transition is made, Clb kinases by themselves are able to depress their antagonists without the help of Cln1,2 kinases. Rising activity of Clb1,2/Cdc28 turns off Cln1,2 synthesis, causing Cln-dependent kinase activity to drop. Hence, after doing their job, the START-facilitators disappear.
Cdc20/APC facilitates the FINISH transition. Cdc20 transcription is activated in G2/M phases by the transcription factor complex Mcm1/Fkh2/Ndd1 (Spellman et al., 1998
; Zhu et al., 2000
; Simon et al., 2001
), which is activated in turn by Clb1,2 kinase activity. In addition, APC core proteins (Cdc16, 23, and 27) are phosphorylated by Clb1,2 kinase, which facilitates APC binding with Cdc20 to form an active complex (Rudner and Murray, 2000
). Cdc20/APC-P depresses Clb kinase activity by labeling Clbs for degradation; it also initiates activation of a phosphatase, Cdc14, which reverses the inhibitory effects of Clb/Cdc28 on Cdh1 and CKIs, so the latter two can overpower the Clb kinases. As the G1-stabilizers extinguish Clb kinase activity, the transcription factor Mcm1 turns off and Cdc20 synthesis ceases. Because Cdc20 is an unstable protein, it quickly disappears, preparing the cell for the subsequent START transition.
To consider G1 and S/G2/M as two alternative steady states, however, is an oversimplification. After START, Clb-dependent Cdc28-kinase activity rises in at least two stages. First, a moderate activity, dependent primarily on Clb5,6, is responsible for initiating DNA synthesis and inactivating Cdh1. Then, a high activity, dependent primarily on Clb1,2, is responsible for driving congression of replicated chromosomes to the metaphase plate. During FINISH, Clb-dependent Cdc28-kinase activity drops in at least two stages. The initial drop in activity is dependent primarily on Cdc20-dependent degradation of Clb16, and it coincides with Pds1 degradation and subsequent separation of sister chromatids (anaphase). As mentioned in the previous section, the disappearance of Pds1 also results in the degradation of PPX (or equivalently, activation of the FEAR pathway), which results in activation of Cdh1 and a second, deeper drop in Clb1,2. The drop in Clb-dependent Cdc28-kinase activity, together with the sustained release of Cdc14-phosphatase from RENT complexes, seems to be responsible for the final stages of exit from mitosis: nuclear division (telophase), cell division, and relicensing origins of DNA replication.
S/G2/M is not a unitary state; mutants reveal states of intermediate and high activities of Clb-dependent Cdc28 kinase. For example, clb1
clb2
cells arrest in G2 phase with moderate Clb-dependent kinase activity; CLB2db
cells (Clb2 protein lacks "destruction box" sequences necessary for Cdc20-mediated proteolysis) arrest in telophase with very high Clb kinase activity; and cdc14ts cells (which activate Cdc20 but not Cdh1) arrest in telophase with an intermediate activity of Cdc28 kinase, due mostly to Clb1,2.
How Can Mutant Cells Lacking START or FINISH Facilitators Be Rescued?
The transitions from G1 to S/G2/M and back again, which we refer to as START and FINISH, are induced by helper proteins Cln2 and Cdc20 as described above. The helpers are involved in negative feedback loops, i.e., the processes they induce lead to their own demise. Rising Clb2 activity after START turns off Cln2 synthesis, and falling Clb2 activity after FINISH turns off Cdc20 synthesis. Because Cln2 and Cdc20 are unstable proteins, they disappear rapidly after their job is done. We describe next how the inviability of various mutants lacking the helpers is related to this central scheme and how they can be rescued.
Mutations involving the facilitators of START and FINISH have striking phenotypes. For example, in the absence of Cln-dependent kinase activity (cln1
cln2
cln3
, "triple-cln" for short), START cannot occur and cells arrest in G1 phase (Wittenberg et al., 1990
). In the simulation (Figure 4A), START is delayed many hours, but eventually the cell grows large enough (in the model) for Clb5-dependent kinase activity to drive [ORI] to 1. This is a serious problem for the model. To account for inviability of the triple-cln mutant, in this model we assume that a cell must reach [ORI] = 1 before a wild-type cell in the same growth medium has divided twice; if not, the cell is considered G1 arrested. Although this rule is introduced specifically to account for the phenotype of triple-cln cells, it is applied uniformly to all mutants to decide whether the simulated cells are viable or G1 arrested.
|
In the model, deletion of either SIC1 or CDH1 allows triple-cln mutant cells to undergo START (Figure 4, B and C) before transgressing the just mentioned rule about G1 arrest. The former construct is viable in the model and in reality (Tyers, 1996
). The latter, though able to replicate its DNA (just barely, according to the "rule"), gets stuck in telophase. Schwab et al. (1997
) studied the response of cdh1
cells to pheromone treatment (analogous to triple-cln deletion). Although the cells were sensitive to pheromone (i.e., they did not form colonies), they did not arrest uniformly in G1. A fraction of the cell population proceeded through S phase and arrested with 2C DNA content. Considering that the simulated cells just barely satisfy our condition for entering S phase, we might expect, in a population of real cells that vary around the mean simulated behavior, some cells will arrest with 1C DNA content and some with 2C DNA content, as observed.
In triple-cln sic1
cells, what molecules initiate the START transition? These cells contain modest amounts of Clb5 in G1 phase ([Clb5] = k's,b5/k'd,b5 = 0.08 au). Because Sic1 is missing and Cdc6, we assume, is a poor inhibitor of Clb5/Cdc28, Clb5-dependent kinase activity is not inhibited. Together with Bck2, Clb5 activates SBF and MBF when cells are large enough. At this point, Clb5-kinase activity rises sharply and initiates DNA synthesis. Clb5 also inhibits Cdh1, enabling Clb2 to accumulate and the division cycle to progress normally.
In the case of cln1
cln2
clb5
clb6
, all the possible helpers for START (except for Cln3) are eliminated, and the cell arrests in G1 (Schwob and Nasmyth, 1993
). In the model, this mutant can be rescued only be restoring CLN2 or CLB5. Deletion of Sic1 (or CKI altogether) is not enough for rescue, because Cdh1 is active, and it keeps the cell arrested in G1. Deletion of Cdh1 (or addition of GAL-CLB2) allows DNA synthesis to initiate and exit of mitosis to occur, but these quintuple mutant cells are predicted to be inviable because they do not form buds (See Figure S1 in Online Supplement A).
Exit from mitosis is initiated by Cdc20, and, as expected, cdc20
is lethal (Lim et al., 1998
; Shirayama et al., 1998
), with cells arrested in metaphase. Cdc20 plays a role in the degradation of Clb2, Clb5, Pds1 (Yamamoto et al., 1996a
; Shirayama et al., 1999
; Yeong et al., 2000
) and our hypothetical PPX. The double mutants cdc20
pds1
and cdc20
clb5
arrest in telophase and metaphase, respectively (Figure 4, D and E), but the triple mutant cdc20
pds1
clb5
is viable (Shirayama et al., 1999
; Figure 4F).
What is the helper for the FINISH transition of the triple mutant cdc20
pds1
clb5
? Shirayama et al. (1999
) have shown that APC/Cdh1 activity is essential for the viability of this mutant. In the model, the negative feedback loop responsible for exiting mitosis is the following. Clb2 activates chromosome alignment on the mitotic spindle (SPN). When SPN = 1, Bub2 is inhibited and Cdc15 is activated, which in turn activates Cdc14 and Cdh1, causing Clb2 level to decrease. The CKIs are restored, and the cell returns to G1. This interpretation predicts that bub2
would render cdc20
pds1
clb5
inviable: because the negative feedback is broken, the quadruple mutant would arrest in the G1 phase.
How Can Mutant Cells Lacking G1-Stabilizers Be Rescued?
Removing G1-stabilizers (Cdh1 and/or CKIs) causes problems at START and FINISH.
Although cdh1
and sic1
are both viable, they show synthetic lethality with viable mutations that have higher than normal Clb2-kinase activity. For example, both cdh1
and sic1
show synthetic lethality with cdc14ts at 29°C (Yuste-Rojas and Cross, 2000
).
Deletion of all three G1-stabilizers (sic1
cdc6
2-49 cdh1
, "triple-antagonist" for short) abolishes the stable G1 steady state. cdc6
2-49 encodes a truncated Cdc6 protein that has lost its role as a CKI but retains its role as a DNA licensing factor. In simulations, the triple-antagonist arrests in telophase, with intermediate activity of Clb2-kinase (Figure 4G). However, in reality, although these cells are inviable (Cross, 2003
), they are not arrested in telophase (Archambault et al., 2003
). In the latter report, the authors showed that these cells have an unusual phenotype, they are able to undergo DNA synthesis and nuclear division but not cell division, and arrest as a "four-cell body object" with 4C DNA content. (We will describe in more detail below the problems of our model with this mutant.)
Although the model is not able to describe the phenotype of the triple-antagonist correctly, it is able to simulate the rescue of this mutant by overproduction of Cdc20 (with GALL-CDC20) as reported by Cross (2003
) (their Supplementary Figure 4). In this case, Cdc20 serves as the G1-stabilizer. Because it is synthesized constitutively, Cdc20 does not disappear from these cells as Clb2 is degraded during FINISH (Figure 4H). SBF and MBF turn on as usual, but Clb5 does not accumulate because it continues to be degraded by Cdc20. The cell is delayed in G1 until Clb5 kinase eventually triggers DNA synthesis. At that point, Cdc20 is inactivated by Mad2, allowing Clb2 to rise and drive the cell into mitosis. When chromosomes are aligned, [SPN] = 1, the inhibition on Cdc20 is removed, and Clb2 is degraded, returning the cell to G1. This interpretation predicts that mad2
would render sic1
cdc6
2-49 cdh1
GALL-CDC20 inviable.
The Model Predicts Phenotypes of Novel Mutants
As we have already pointed out with regard to (cln1
cln2
clb5
clb6
cdh1
), (cdc20
pds1
clb5
bub2
), and (sic1
cdc6
2-49 cdh1
GALL-CDC20 mad2
), we can use the model to predict phenotypes of mutants that have not yet been investigated experimentally, to test crucial features of the model.
For example, we assume that the Cdc6 is a much weaker inhibitor of Clb5,6 kinases than is Sic1 and that Clb5,6 kinases are not able to phosphorylate and destabilize Cdc6 as efficiently as Clb1,2 kinases. The first assumption that Cdc6 is a weak inhibitor of Clb5 kinase is based on the following observations. 1) sic1
cells have very short G1 period (Schneider et al., 1996
), initiating DNA synthesis before SBF/MBF activation and budding. For DNA replication to occur early in those mutant cells, the high concentration of Cdc6 in G1 must not be able to inhibit even the low level of Clb5,6 kinases generated by MBF-independent transcription. Hence, Cdc6 cannot be a strong inhibitor of Clb5,6. 2) GAL-CLB5 cells do not show advancement in DNA synthesis (Schwob et al., 1994
), indicating that Clb5/Cdc28 must be effectively inhibited by Sic1. This assumption was made independently of experiments published recently by Archambault et al. (2003
), who showed that Sic1 coimmunoprecipitates with Clb5 but Cdc6 does not.
The second assumption, that Clb5 kinase is less able to phosphorylate and inactivate Cdc6, is inferred from the large size of cln1
cln2
cln3
sic1
cells. As described in 1), Clb5,6 kinases are able to initiate DNA synthesis early and to inactivate Cdh1 in the quadruple mutant, but Clb2 is still inhibited by Cdc6. The cell has to grow to very large size to accumulate enough Clb5 to inactivate Cdc6. Only then can Clb2 kinase activity rise, driving the cell into mitosis.
If it is true that Cdc6 is a weaker inhibitor to Clb5 kinase than is Sic1, then whenever the activity of Clb5 kinase is important, sic1
will have very different effects than cdc6
2-49. For example, triple-cln sic1
is viable (Tyers, 1996
), but the model predicts that triple-cln cdc6
2-49 will remain arrested in G1. Similarly, the model predicts (Table 4, row i) that the inviability of GAL-CLB5-db
(Jacobson et al., 2000
) can be rescued by overproducing Sic1 but not Cdc6. On the other hand, because Sic1 and Cdc6 are both strong inhibitors of Clb2, the inviability of Clb2-db
(Wasch and Cross, 2002
) should be rescued by overproduction of either Sic1 or Cdc6 (Table 4, row ii), as confirmed recently by Archambault et al. (2003
).
|
Viability of the triple mutant cdc20
pds1
clb5
(Shirayama et al., 1999
) depends on a feedback loop that activates Cdh1 at the end of the cycle. If we perturb elements in this loop (Table 4, row iii), we are likely to render cdc20
pds1
clb5
inviable. On the other hand, the viability of cdc20
pds1
clb5
does not depend on CKIs.
As shown in Figure 3, Cdc20 helps the FINISH transition by degrading Clb kinases and by activating the CKIs and Cdh1 (through Cdc14). Hence, in row iv (Table 4), the model predicts that inviability of the double mutant cdc20
pds1
(telophase arrest, with Cdc14 released from the nucleolus; Shirayama et al., 1999
) could be rescued by an extra copy of genomic CDC15. The double mutant also can be rescued by TAB6-1 (a dominant Cdc14 mutation where Cdc14 binding to NET1 is reduced; Shou et al., 2001
). In TAB6-1 cells, less Cdc14 is sequestered in the nucleolus, making it easier for the triple mutant (cdc20
pds1
TAB6-1) to exit from mitosis even in the presence of Clb5.
Row v (Table 4) shows when Cdc20 is crippled, as in APC-A mutants, Cdh1 is more important than the CKIs in getting out of M phase. Row vi (Table 4) suggests that ability of Cdc20 overexpression to rescue the lethality of the triple-antagonist critically depends on the checkpoint mechanism. Row vii (Table 4) indicates that net1ts can retain viability in nocodazole if Clbs are overexpressed but not if CKIs are deleted.
The Model Predicts Effective Values of Rate Constants
Table 2 makes
100 predictions about the rates of component biochemical processes involved in cell cycle control, e.g., synthesis and degradation rates, phosphorylation and dephosphorylation rates, relative activities of cyclins with overlapping substrate specificities. Some of these rates are well established experimentally (e.g., protein half-lives are easily measured), whereas most others were completely unknown. As biochemists seek to measure these rate constants and binding constants directly, our estimates will serve as landmarks for relating in vitro measurements to in vivo activities. Some of these predictions are described in detail in Online Supplement A, along with a few supporting observations.
Inconsistencies between Simulations and Observations Identify Problems in the Model
Although the model gives a quantitatively accurate representation of many aspects of the budding yeast cell cycle, it fails to account for certain properties of 11 mutants in our test series (Table 3). Although these failures (Table 5) may reflect faulty parameter settings (Table 2) or mistranslations of the mechanism into equations (Table 1), careful consideration of the inconsistencies (see Online Supplement B) indicates that there are problems in the wiring diagram itself (Figure 1). Identifying these problems suggests ways to improve the model in the future.
|
The most troublesome mutants for the model are the triple mutant cdh1
sic1
cdc6
2-49 and the double mutant swi5
cdh1
strains, which have very similar phenotypes (Table 5).
The triple-antagonist strain (cdh1
sic1
cdc6
2-49, where all three G1-stabilizers are deleted) is inviable (Archambault et al., 2003
), but the cells are not telophase arrested. When a triple-antagonist GAL-SIC1 (integrated) strain is transferred from galactose to glucose medium, the mutant cells reproducibly accumulate as four-cell body objects that are resistant to sonication. The cells are able to undergo DNA synthesis and nuclear division but not cell division, and they arrest as binucleate cells with 4C DNA content. In our simulation, this mutant would be unbudded, arrested in telophase with 2C DNA content.
We do not know how to interpret the phenotype of the triple-antagonist and how to model it properly. There are three problems. First, how can the mutant exit from mitosis? What drives Clb2 activity below the threshold for nuclear division? Degradation by Cdc20 is not enough, because cdc14
(which would be equivalent to the triple-antagonist) arrests in telophase. The observed phenotype of the triple-antagonist suggests that our requirement for nuclear division is incorrect and that Cdc14 must have additional roles besides activating CKI and Cdh1.
Perhaps it is the rise in Cdc14 phosphatase activity, rather than (or in addition to) the fall of Clb2 kinase activity, that triggers nuclear division (see Online Supplement C). Suppose that nuclear division is triggered when the phosphorylation state of a target protein drops below a critical value. If the target protein is phosphorylated by Clb2 kinase and dephosphorylated by Cdc14 phosphatase, then the triple-antagonist cells, which have high Cdc14 activity, could exit from mitosis. However, there are other problems.
Even if they could execute nuclear division, they would have trouble forming a bud in the next cell cycle. The persistent high Clb2 kinase activity after nuclear division would prevent SBF activation and Cln2 accumulation, hence in simulation [BUD]max = 0.25 (it never reaches 1). However, triple-antagonist cells do form buds. It may be that Clb2 kinases are able to trigger bud formation, albeit at very low efficiencies compared with Cln2 kinases.
A third problem for the triple-antagonist cells is that high Clb2 kinase ac