BACKGROUND
the machinery of the eukaryotic cell cycle has been extensively dissected and described, in both simple and complex organisms.
proliferation hinges on the cell’s ability to replicate the genome with high fidelity, segregate the chromosomes equally, and ultimately divide into two genetically identical cells.
a fundamental process in the regulation of dna replication is the step-wise assembly of the pre-replicative complex  at origins of replication.
the pre-rc is a congregation of proteins each performing a specific role.
its formation is facilitated by the six-subunit origin recognition complex , which, in the budding yeast saccharomyces cerevisiae, binds an 11 bp consensus sequence  <cit> .
orc then recruits cdc <dig>  which, like orc, exhibits atpase activity  <cit> .
the co-import of cdt <dig> and the mcm2- <dig> complex  into the nucleus follows  <cit> , and the mcm·cdt <dig> heptamer is then targeted to origins by an interaction between cdt <dig> and orc <dig>  <cit> .
initial loading of an mcm ring at the origin requires cdc <dig> atp-hydrolysis.
reiterative loading of an additional mcm molecule occurs via orc atp-hydrolysis  <cit> , resulting in two rings at each origin  <cit> .
at this point origins are said to be licensed.
in late g <dig> phase, a burst of dbf <dig> synthesis activates the dbf4-dependent kinase cdc <dig> , which then phosphorylates multiple mcm subunits  <cit> , bringing about a conformational change that stimulates mcm helicase activity.
dbf <dig> levels decrease over the course of s-phase and, starting at the metaphase/anaphase transition, dbf <dig> is actively degraded by the anaphase promoting complex  and its activating co-factor, cdc <dig>  <cit> .
in this way, dbf <dig> levels are prevented from rising until the next g1/s transition.
the phosphorylation of mcm by ddk is coincident with the phosphorylation of the protein factors sld <dig> and sld <dig> by clb5-cdc <dig>  a cyclin-dependent kinase  complex, the activity of which rises just prior to s-phase entry.
the sld proteins, once phosphorylated, are stabilized as a complex with the adaptor protein dpb <dig> and the tetrameric gins complex, forming a module that interacts with cdc <dig>  the latter acts as a scaffold for this module, which is then competent to associate with the pre-rc and attract dna polymerase  <cit> .
a recent study shows that the end result is the tight association of cdc <dig>  mcm and gins  with origins, allowing the unwinding of dna and processive replication by dna polymerase  <cit> .
this represents the essential role of cdk in stabilizing polymerase at the moving replication fork and switching the system from a pre-replicative state to a replicative one.
from this point until late in mitosis, cdk levels remain high.
this continued cdk activity prevents re-establishment of pre-rcs at origins that have already fired through a number of mechanisms.
firstly, cdk phosphorylates cdc <dig>  thus causing the scfcdc <dig> complex to target cdc <dig> to the proteasome for degradation  <cit> .
secondly, orc <dig> and orc <dig> are phosphorylated by cdk  <cit> , with the phosphorylation of orc <dig> rendering it refractory to interaction with cdt <dig>  <cit> , thereby preventing further mcm loading.
finally, cdk facilitates the nuclear export of both mcm and cdt <dig>  at different time points.
just prior to initiation, cdt <dig> exits via a cdk-dependent mechanism, while mcm proteins fall off the dna upon fork termination and are then exported in a cdk-dependent manner  <cit> .
thus, while cdk initiates replication, it subsequently prevents pre-rc reassembly.
this illustrates its dual role in triggering initiation through formation of cmg, then preventing re-initiation by inhibiting pre-rc reformation.
mathematical modeling has been successfully used in the past to address various aspects of the cell cycle.
early models  did not incorporate specific biochemical mechanisms; they were hypothetical representations of periodic cellular activity.
as the molecular mechanisms driving the cell cycle were revealed, models appeared that incorporated these findings .
for s. cerevisiae in particular, multiple modeling approaches have been applied, based both on network descriptions  <cit>  and on specific molecular details such as gene expression and biochemical kinetics  <cit>  .
some modeling efforts have been comprehensive, such as the tyson group’s ordinary differential equation -based models  <cit> , while others address specific cell-cycle phenomena, such as the links between cell size and cycle progression  <cit> .
spiesser et al.
<cit>  developed a model of chromosomal replication, which reproduced the spatio-temporal replication profile of yeast chromosomes.
origin firing was also described in  <cit>  wherein the authors used a stochastic model to describe these origin-specific features of replication.
a recent report  <cit>  presented an ode-based model describing the initiation of dna replication, incorporating origin licensing, firing and the network of regulatory phosphorylation events.
the model parameters were partly calibrated against experimental data, but largely selected through an optimization routine designed to attain an idealized function, resulting in a model that is particularly suited to exploring events specifically at the g1/s transition.
here, we present a new model of the initiation of dna replication.
in contrast to the work of brümmer et al.
<cit> , we took a ‘bottom-up’ approach and began by gathering in vivo data for precise protein levels at specific cell cycle time-points, then calibrated our model against these values.
rather than limiting ourselves to the observation of firing near the g1/s transition and fitting to dna-specific replication profiles, we validated our model against the behaviour of the constituent protein complexes throughout the entire cell cycle.
to facilitate the use of our model in a comprehensive description of the cell cycle, we designed it to integrate easily with the model of chen et al.
<cit> .
finally, we validated the model by comparing in silico predictions to experimental observations, using both our own knockdown experiments and results from the literature.
the model presented here consolidates the known interactions between dna replication initiation proteins and the mechanisms that allow them to drive genome duplication.
additionally, regulatory aspects of the system, which ensure that re-replication does not occur, have been modeled.
the model’s behaviour provides a falsifiable hypothesis regarding the dynamics of dna replication initiation.
furthermore, it accurately predicts the phenotypes of known experimental cell cycle mutants as well as those arising through in vivo perturbations to proteins in the network.
because our model is constrained only by the fluctuating levels of replication factors, it provides a unique understanding of the kinetics governing the reactions between them.
successful integration into a whole cell cycle model allows the initiation of dna replication to be explored in a broader quantitative context.
RESULTS
we began construction of our model by identifying the important players in replication initiation and establishing an interaction network, as shown in figure  <dig>  after selecting appropriate descriptions of reaction kinetics, we generated an ode-based model and calibrated the model parameters to in vivo data.
description of model components
the model describes sixteen molecular species  and depends on twenty-four parameters, which characterize the rates of seventeen biochemical processes .
the model describes the following molecular species .
rc <dig> : origin-bound orc
rc2: origin-bound orc associated with cdc6
rc3: origin-bound orc associated with cdc <dig>  with mcm loaded
rc4: origin-bound orc, with mcm loaded
rc5: origin-bound orc, with mcm loaded and dbf <dig> associated
rc6: origin-bound orc, with mcm loaded and dbf <dig> and cdc <dig> associated
rc7: origin-bound phosphorylated orc
fork: the elongation fork, with mcm and cdc <dig> associated
cdc6n: non-chromatin associated nuclear cdc6
dbf4n: non-chromatin associated nuclear dbf4
cdc45n: non-chromatin associated nuclear cdc45
mcmc: cytosolic mcm
cdt1c: cytosolic cdt1
mcm·cdt1n: non-chromatin associated nuclear mcm bound to cdt1
cdt1n: non-chromatin associated nuclear cdt1
mcmn: non-chromatin associated nuclear mcm
the mcm species corresponds to dimers of mcm2- <dig> heterohexamers, as two complexes are loaded at each origin.
similarly, the cdc <dig> species corresponds to a dimer, as described by bowers et al.
<cit> .
we describe concentrations in units of molecules per cell.
reaction events
the seventeen processes that make up the model are shown in table  <dig>  their rates depend on the species concentrations, the model parameters, and on two fixed, time-varying input functions describing the abundance of clb <dig>  and of cdc <dig>
in choosing reaction kinetics, we balanced the complexity of the model against its ability to adequately describe the behaviour of the overall system.
we limited our description of initiation to the interactions between the pre-rc and replisome proteins that we found to be the essential core of the network .
as a result, certain processes were combined into single events, some reactions were presumed irreversible, and only some reaction rates were presumed to have non-linear kinetics.
except for rc <dig>  phosphorylation states are not explicitly described, as we have no data for the individual phosphorylation events.
this is acceptable for our purposes as the lumped function of cdk in each case is consistent with a scenario where the effect of cdk is proportional to its concentration .
additionally, processes that involve multi-protein complexes are represented by a single member – one cmg  complex stabilizes dna polymerase at each replication fork.
of the three protein factors it is comprised of, cdc <dig> is limiting.
although mcm is also included in the gins complex, we model both mcm and cdc <dig> as separate species.
dbf <dig> represents the dbf4-cdc <dig> kinase complex and mcm <dig> represents the mcm2- <dig> helicase.
although the protein factors cdc <dig>  dbp <dig>  sld <dig>  sld <dig> and gins interact to facilitate formation of the pre-initiation complex at origins, we model only cdc <dig>  which is the limiting factor in the cmg complex  <cit> ; ultimately the number of forks fired  is dependent on the cdc <dig> concentration.
we take dbf <dig>  which is the limiting regulatory subunit of cdc <dig>  as representative of active ddk, which is one of the limiting factors in replication initiation  <cit> .
mcm <dig> is used to represent mcm complexes; the mcm <dig> concentration has been reported to approximate the number of total complexes per cell in an asynchronous population  <cit> .
the replication complexes in our models exist only on chromatin and therefore represent the activity of these proteins at the dna as opposed to soluble complexes.
the network shown in figure  <dig> includes both reversible and irreversible reactions as indicated.
association/dissociation reactions are considered reversible, in accordance with a dynamic pre-rc/pre-ic loading mechanism as described above.
in most cases, phosphorylation events are modeled as irreversible, in the absence of identified countervailing enzymes.
we found that it was sufficient to describe most reaction rates by mass action kinetics.
in cases where saturation occurs , we employed michaelis-menten kinetics.
to simplify the description of the phosphorylation of orc by cdk , we do not describe phosporylation and dephosphorylation explicitly, but combine them into a single dephosphorylation event whose rate is inversely proportional to the level of cdk .
we introduce cooperativity in this mechanism to account for multiple phosphorylation events  <cit>  or an additional inhibitory cdk-orc <dig> binding mechanism  <cit> .
the establishment of replication complexes in our model reflects the sequential binding of proteins that constitute the pre-replicative complex.
in some cases the association and dissociation of pre-rc components is reversible.
we treat the loading and maintenance of mcm2- <dig> helicase complexes as a dynamic process, which is dependent on the concentrations of the factors orc, cdc <dig>  cdt <dig> and mcm2- <dig> itself.
a mechanistic model for the dynamic assembly of pre-rcs was first described by the bell lab  <cit> .
the requirement of pre-rc factors for maintenance of helicase-loaded origins in late g <dig> has been further demonstrated by work from these researchers as well as our group  <cit> .
network and differential equations
referring to figure  <dig> and table  <dig>  the dynamics of the system are described as:
drc2dt=v5−v7dforkdt=v14−v15drc3dt=v7−v10dcdc6ndt=v1+v10−v5−v2drc4dt=v10−v12−v11ddbf4ndt=v3+v14−v12−v4drc5dt=v12−v13dcdt1ndt=v7−v8+v9drc6dt=v13−v14dmcmndt=v15+v11+v9−v16drc7dt=v14−v17dmcm·cdt1ndt=v6−v9−v <dig>
the remaining state variables are constrained by the following conservations:
rc1 = rctotal− rc2 − rc3 − rc4 − rc5 − rc6 – rc7cdt1c= cdt1total− cdt1nmcmc= mcmtotal− mcmn− rc3 − rc4 − rc5 − rc6 − fork − mcm·cdt1cdc45n= cdc45total− rc6 − fork
where rctotal, cdt1total, mcmtotal, and cdc45total are the fixed total number of origins, cdt <dig> molecules, mcm2- <dig> complexes, and cdc <dig> dimers, respectively.
these four factors have been shown to be present at constant levels throughout the cell cycle  <cit> .
the value for rctotal used in our model is  <dig>  as described in  <cit> .
system inputs
the biological network responsible for the initiation of dna replication does not oscillate autonomously; it displays periodic behaviour when driven by periodic signals from the cell cycle.
likewise, our model displays oscillations only when driven by periodic forcing input.
in order to facilitate the combination of our model with the cell cycle model of chen et al.
<cit> , we used the simulated profiles of clb <dig> and cdc <dig> from their model as periodic inputs to ours.
cdc <dig> mediates the degradation of dbf <dig> .
clb <dig> is responsible for cdc <dig> degradation , loading of cdc <dig> , nuclear export of free mcm  and cdt <dig> , and phosphorylation of orc <dig> and orc <dig> .
we converted the time-varying profiles of clb <dig> from the chen et al.
<cit>  model to molecules-per-cell units using the genome-wide gfp tagging and localization experiments described in  <cit> .
the profile of cdc <dig> was similarly obtained by scaling to cellular abundance levels reported in another study – while cdc <dig> has been determined to peak at  <dig> copies in a haploid cell, the functional apccdc <dig> level can be estimated by considering the apc cyclosome subunit cdc <dig>  <cit> .
this was reported in different studies to be 593 mol/cell in an asynchronous population  <cit>  and at its maximal value of 750 mol/cell in metaphase  <cit> .
data acquisition
data for cdc <dig> and cdc <dig> levels were obtained from individual isogenic strains in which the open reading frame of the corresponding gene was fused to a sequence encoding a 13myc epitope tag  <cit> .
in figure  <dig>  a representative western blot for cdc45-myc is shown , with the corresponding facs analysis .
the levels of mcm <dig> were determined using an anti-mcm <dig> antibody.
in each of our time course experiments, cells were first arrested in late g <dig> phase with the mating pheromone α-factor and then released synchronously into the cell cycle, as described in methods.
from the literature, we used time course data for chromatin-bound and soluble dbf <dig> and mcm <dig>  from pasero et al.
<cit>  and quantitation of the nuclear fraction of cdt <dig> from tanaka and diffley  <cit> .
in order to convert relative measures of protein abundance to molecule-per-cell numbers, we used scaling factors obtained from the database provided by ghaemmaghami et al.
<cit> .
the data is shown along with a best-fit simulation in figure  <dig>  raw timecourse data can be found in additional files  <dig> and  <dig>  our calculation of molecule-per-cell estimates is demonstrated in additional file  <dig>  the fits in figure  <dig> represent our best solution to a trade-off between quality of fit and model complexity.
we explored the effect of adding additional species and parameters.
these additional features could, in some cases, provide minor improvements to the fit, but our confidence in the parameter estimates  suffered as the complexity of the model grew.
while our model reports proteins and complex concentrations in absolute units of molecules/cell, our accuracy regarding these values is limited to by the literature-reported cellular abundances for the various protein factors.
we have used discretion when inconsistencies arise, choosing the reported values that are best justified by multiple studies.
the network dynamics are based on the relative protein abundances over the time-course .
consequently, changes in the global protein level for a particular factor do not affect the dynamics; instead, they impact our predictions of absolute molecule/cell counts and parameter values.
this is an important consideration regarding the conversion of densitometry readings to absolute values, as the overall levels are ultimately determined by a literature-derived scaling factor.
future experiments will undoubtedly result in improved estimates of protein abundances.
these can be easily incorporated into the model by scaling the parameter values accordingly .
parameter calibration
the model parameters were calibrated using a weighted least-squares comparison with the data described above.
we used a combination of global optimization  and local search  to find the best-fit parameter set shown in table  <dig>  the table also shows the percent error associated with each parameter estimate.
the percent error is the relative size of a 95% confidence interval for the estimate, calculated via the fisher information matrix and the cramer-rao bound  <cit> .
the percentage errors show that some parameters are estimated with high confidence while others are represented with less accuracy.
these values were used to solve the odes in our consensus model.
the percentage error for each parameter is indicated.
parameter values that were well constrained by the data include those specifying the rates of production, degradation and association of cdc <dig>  and dbf <dig>  as well as the rate of origin firing .
this reflects the strong reliability of our data for these two protein factors as well as for the proteins that form the replication complex  that gives rise to active forks.
parameters values in which we have low confidence include those that govern the loading of mcm by cdt <dig> , cdc <dig> dissociation from rc <dig> , and the phosphorylation of orc .
the reversible dissociation of cdc <dig> is needed to accurately fit the data and there is no evidence suggesting that orc-cdc <dig> binding is irreversible.
nevertheless, it is clear that experimental observations specific to this process are required to more precisely estimate this parameter value.
the reaction whereby the mcm·cdt <dig> species loads the mcm complex  is extremely transient  <cit> .
provided parameter k <dig> is sufficiently large, the kinetics of this reaction will be rapid enough to fit the data.
consequently, the data cannot support a precise estimate of the parameter value.
this observation suggests that mcm loading is an extremely rapid biochemical step in pre-rc assembly.
it may point to a role for cdt <dig> in repeatedly targeting mcm complexes to origins throughout g <dig>  such a phenomenon is consistent with the requirement for a dynamic loading mechanism that ensures pre-rc fidelity up until the g1/s transition.
finally, the phosphorylation of orc  contributes to the prevention of repeated origin firing.
however, this mechanism has not been well characterized, and our data is unable to accurately constrain the specifics of this process.
the kinetic rates in this network have not been the subject of prior experiments, but previous reports of protein half-lives are consistent with our predicted parameter values.
drury et al.
<cit>  estimated that cdc <dig> is reduced below the point of detection within  <dig> minutes of s-phase entry, corresponding to a half-life no longer than  <dig>  min.
similarly, cheng et al.
<cit>  reported that dbf <dig> is reduced below visible levels within  <dig> minutes by the apc-dependent pathway, indicating a half-life no longer than 3 min.
our model-based predictions of degradation rates correspond to half-lives of 1 min.
and  <dig>  min.
for cdc <dig> and dbf <dig> respectively, in good agreement with these earlier findings.
figure  <dig> shows the simulated model behaviour for the best-fit parameter set.
some replication complex species – rc <dig>  rc <dig> and rc <dig> – are extremely transient.
their low levels of abundance are shown separately from other rcs, on an appropriate scale.
simulations were carried out in matlab .
perturbations
our initial explorations of the model revealed that the network’s behaviour is particularly sensitive to the abundance of dbf <dig> and cdc <dig> and relatively insensitive to the level of cdt <dig>  we investigated the effects of perturbations by simulating reductions in dbf <dig>  cdt <dig> and cdc <dig>  in our model.
when the cdc <dig> production rate  was reduced to 10% of its nominal  value, persistence of the rc <dig> complex was observed.
similarly, when the dbf <dig> production rate  is reduced by the same relative amount, an accumulation of rc <dig> occurs.
in both cases, the perturbation interferes with pre-initiation complex assembly and blocks the system at the nearest previous persistent rc state; rc <dig> is not persistent since the unloading of mcm causes a rapid transition back to rc <dig>  it is worth noting that because mcm can dissociate from orc , rc <dig> represents a complex containing mcms that will be functionally incorporated into replication forks as opposed to those that loosely associate with origins.
because the timing of our model is fixed, the various state concentrations  indicate the progression from licensing to firing.
a reduction in the fork species compared to the wild-type case suggests a slow-down in s-phase because fewer origins are firing within the prescribed time.
using the peak abundance of the fork species as a measure of replicative efficiency, we saw significant reductions in both simulated knock-downs .
conversely when we simulated the reduction of cdt <dig> abundance to 10% of nominal values, origin firing was only reduced by 23%, suggesting that the network is relatively refractory to depletion of cdt <dig> .
additional file  <dig> shows the levels of the various model components for each perturbation.
to investigate the accuracy of these mathematical predictions, we carried out corresponding wet lab depletion experiments.
reducing dbf <dig> or cdc <dig> concentrations in yeast cells to roughly 90% below normal endogenous levels resulted in a rapid g <dig> phase arrest, evident after  <dig> hours of depletion, as judged by facs analysis indicating the accumulation of cells with 1 c  dna content .
in contrast, a corresponding depletion of cdt <dig> had no appreciable effect, and dna replication defects were still only minimally evident after  <dig> hours of further reduction.
thus, our in silico simulations using our nominal parameter set were predictive of in vivo perturbations.
these experiments were used to validate the model; they were not used for calibration.
the insensitivity to perturbations in cdt <dig> levels is consistent with its apparent excess relative to origins  <cit> , although the number of cdt <dig> molecules that act at each origin has not yet been characterized.
moreover, the mechanism by which cdt <dig> aids in recruiting the helicase molecules to pre-rcs is extremely transient  <cit> .
while many factors are limiting, the system appears to be highly sensitive to the levels of cdc <dig>  due to its low abundance relative to mcm and cdt <dig>  even a moderate depletion of cdc <dig> significantly alters the dynamics of pre-rc loading.
the same is true for dbf <dig>  although in this case its role in activation of the cdc <dig> kinase renders the system highly sensitive to its concentration; firing cannot occur without the dbf4-cdc <dig> complex.
since dbf <dig> is, like cdc <dig>  limiting, flow through the network is blocked when the kinase does not reach a threshold level.
additionally at limiting levels, the number of replication forks produced by our model is significantly reduced, consistent with in vivo reports from the literature showing a lengthening of s-phase  <cit> .
further validation of our model comes from comparison with additional in vivo experiments reported in the literature.
jones et al.
<cit>  showed that the interaction between the mcm complex and dbf <dig> was reduced to half its wild-type level when a dbf <dig> domain that binds mcm <dig> was mutated, impairing s-phase progression.
we mimicked this effect by reducing the rate of association of dbf <dig> with rc <dig>  by 50%, leading to a similar result .
similarly, zou et al.
<cit>  reported that the cdc45- <dig> mutant shows an aberrant growth phenotype at the non-permissive temperature.
this is thought to be due to a disruption of cdc45’s ability to interact with mcm and orc .
as shown in figure 7c, by reducing the rate of cdc <dig> interaction with rc <dig>  by 50%, a marked reduction in the peak abundance of the fork species results, indicative of a slower s-phase, as observed when the mutant was grown at the non-permissive temperature.
the actual reduction in cdc45’s association with the pre-rc due to conformational changes in the mutant might be even more pronounced than a 2-fold reduction.
in any case, our simulation is consistent with cdc45’s origin-initiation role being compromised by impairing its ability to interact with its ligands to form the cmg complex.
linking our dna initiation model to a previously established cell cycle model
our model of the initiation of dna replication only displays oscillatory behaviour when forced with periodic signals from the cell cycle.
by choosing to incorporate signals that correspond to species in the cell cycle model of chen et al.
<cit> , we were able to merge the two models in a straightforward manner.
in the chen model, the initiation of dna replication is represented by a single lumped state variable, called ori.
at the beginning of the cell cycle, ori has value zero.
its rate of growth depends linearly on clb <dig>  when it reaches a threshold value, dna synthesis is presumed to have begun, and triggers an increase in the value of the parameter kmad <dig>  leading to an inactivation of cdc <dig>  which is required for mitotic exit.
this mad2-dependent inhibition of cdc <dig> represents the spindle assembly checkpoint  <cit> , ensuring that cells with replicated dna do not complete mitosis without properly aligning the chromosomes.
when chromosomes have properly aligned on the metaphase plane kmad <dig> drops and cdc <dig> promoted exit from mitosis.
in our model, the level of dna synthesis is represented by the fork species.
to merge the two models, we removed the ori state from the chen model, and instead used the fork species to trigger the change in kmad <dig>  as detailed in methods.
besides “closing the loop” between the two models by incorporating two-way inter-model signalling , we had to deal with a single shared species: both models describe the dynamics of cdc <dig>  we arrived at a merged description of cdc <dig> behaviour by incorporating the dynamics of recognition complex association and dissociation into the chen model’s formulation of cdc <dig> behaviour .
the resulting combined model behaves only marginally differently from either model in isolation, as shown in additional file 6: figure s <dig> and additional file 7: figure s7).
replication metrics of brümmer et al
the recently published model of brümmer et al.
<cit>  also describes the network responsible for the initiation of dna replication.
the  <dig> free parameters of that model were chosen by a combination of fitting and optimization.
the authors used literature-derived data to fix  <dig> of the kinetic parameters.
the remaining  <dig> free parameters were not fit to data, but were selected through a procedure that optimized the coherence of origin firing and minimized re-replication .
while it is impossible to assess the accuracy of the parameter values obtained from this procedure, the resulting idealized model provides a useful starting point for examining how the network structure constrains the system behaviour.
the model of brümmer et al.
focuses on early origin firing and so represents the mechanics of firing at the start of s-phase.
in contrast, our model describes firing dynamics throughout s-phase in order to fit into the broader context of the cell cycle  <cit> .
the parameter set driving our system is not filtered to retain only those that produce replication dynamics consistent with coherent firing just at the g1/s boundary.
rather, the parameters are specified by the actual cellular concentrations of the active protein factors generating replication forks.
while both models incorporate the important role of cdk, brümmer et al.
emphasize the multi-site phosphorylations of several factors involved in mechanisms that minimize potential re-replication.
to this end, brümmer et al.
employed a metric to assess re-replication.
their idealized model exhibits  <dig>  re-replication events per cycle.
applying the same measure to our bottom-up model yields  <dig>  re-replication events per cycle .
the near-zero value obtained by brümmer et al.
is close to their idealized target of zero.
both estimates are consistent with the belief that re-replication occurs in wildtype cells, but at an extremely low rate  <cit> .
because the nature of the dephosphylation of orc  remains uncharacterized, we use a conservative estimate of the number of orc phosphorylation sites.
increasing this number by twofold, consistent with the number of cdk target residues on orc  <cit>  reduces re-replication to a value on the same order of magnitude as brümmer’s value.
thus, both models effectively deal with representing control and prevention of rereplication.
CONCLUSIONS
while our model provides a sound description of the initiation of dna replication, a number of aspects of the network remain unresolved: for example, the kinetics of mcm loading, the mechanism by which cdk phosphorylates orc, and the details of the association of the gins complex, sld <dig>  and sld <dig>  while modeling cdc <dig> captures the events regarding cmg formation at origins, being its limiting factor, a future version of the model could better distinguish the initial cdc <dig> association at origins from subsequent cmg formation.
while this has minimal effect on our network dynamics and no effect on our blocking of re-replication, it would provide a better resolution of events at origins just prior to the g1/s transition.
incorporating timecourse experiment data for levels of a gins complex member would aid in this analysis.
nevertheless, the assumptions made allow us to approximate the aforementioned processes, simplifying the network without losing information about system behaviour at the level we intend to model.
with our nominal parameter set, we observe the system to behave as the ordered accumulation of proteins forming a loading complex at origins throughout the genome.
activation by increasing s-cdk levels and the concentration of dbf <dig>  increase linearly the number of replication forks established as a result of successive forward transition of the various replication complexes .
it should be noted that as is found in vivo, not all origins fire as a consequence of being furnished with mcm-containing pre-rcs.
replication is maximal at the g1/s transition, but continues into s-phase as origin firing is temporally spaced.
this is thought to ensure sufficient time to address any defects in replication and is mediated by the limiting nature of one of the initiation activators, ddk .
our in vivo perturbation of dbf <dig> levels reproduces this consequence and points to other system observations: cdc <dig> levels are intimately controlled by cdk levels to avoid re-replication, however this mechanism is tightly regulated such that clb <dig> levels rising too soon would prevent the assembly of the pre-rc in g <dig>  a feature of the system well documented.
additionally, cdt <dig> appears to act catalytically rather than stoichiometrically given the system is relatively impervious to a reduction of this factor to  <dig> % of its wildtype level.
this might play into its role in chaperoning mcm2- <dig> hexamers to origins, where they are loaded subsequently leading to the release of cdt <dig>  which may then be recycled to mediate the loading of other mcms.
this aspect of the system has not yet been investigated experimentally and would be of future interest.
many human orthologs of the yeast proteins described in our network have been associated with cellular pathologies.
our model is specific to the replication machinery in budding yeast, but the mechanisms driving this process are highly conserved throughout eukarya.
efforts to develop an analogous model in mammalian cells would be useful in understanding and dissecting cell proliferation in humans.
a number of models of the mammalian cell cycle have been proposed  <cit> .
incorrect pre-rc formation has been linked to impaired dna damage repair pathways in humans  <cit> , while both orc <dig>  <cit>  and members of the mcm2- <dig> complex  have been shown to be reliable cancer biomarkers.
recent work by bicknell et al.
<cit>  has shown that point mutations in the human orc <dig>  orc <dig>  orc <dig>  cdt <dig> and cdc <dig> genes are associated with meier-gorlin syndrome, a form of primordial dwarfism, and several of these mutations were determined to interfere with proper pre-rc formation.
these findings highlight the potential utility of in silico mammalian models in further exploring the molecular basis of such disorders.
given that our model shows good predictive capability, it serves not only as an informational tool for yeast biology, but also as a proof of principle for higher order system models.
despite the requirement for a mammalian model to comprehensively verify specific mechanisms, the system of dna replication initiation is conserved well enough that perturbations to proteins such as those described above can, in fact be preliminarily examined.
although previously established replication models  <cit>  consider the ordered timing of origin firing based on genomic replication profiles, our goal was to represent the temporal organization of origin firing as a function of the concentration of active replication species.
a focus on using real protein levels as a determinant of replication dynamics is a novel approach.
when used in concert with models describing genome-level origin characteristics and/or combining the findings with models exploring other cell-cycle modules, future efforts will generate a well-rounded picture of dna replication initiation.

 the machinery of the eukaryotic cell cycle has been extensively dissected and described, in both simple and complex organisms. proliferation hinges on the cell’s ability to replicate the genome with high fidelity, segregate the chromosomes equally, and ultimately divide into two genetically identical cells. a fundamental process in the regulation of dna replication is the step-wise assembly of the pre-replicative complex at origins of replication. the pre-rc is a congregation of proteins each performing a specific role. its formation is facilitated by the six-subunit origin recognition complex, which, in the budding yeast saccharomyces cerevisiae, binds an 11 bp consensus sequence. orc then recruits cdc which, like orc, exhibits atpase activity. the co-import of cdt and the mcm2- complex into the nucleus follows, and the mcm·cdt heptamer is then targeted to origins by an interaction between cdt and orc. initial loading of an mcm ring at the origin requires cdc atp-hydrolysis. reiterative loading of an additional mcm molecule occurs via orc atp-hydrolysis, resulting in two rings at each origin. at this point origins are said to be licensed. in late g phase, a burst of dbf synthesis activates the dbf4-dependent kinase cdc, which then phosphorylates multiple mcm subunits, bringing about a conformational change that stimulates mcm helicase activity. dbf levels decrease over the course of s-phase and, starting at the metaphase/anaphase transition, dbf is actively degraded by the anaphase promoting complex and its activating co-factor, cdc. in this way, dbf levels are prevented from rising until the next g1/s transition. the phosphorylation of mcm by ddk is coincident with the phosphorylation of the protein factors sld and sld by clb5-cdc a cyclin-dependent kinase complex, the activity of which rises just prior to s-phase entry. the sld proteins, once phosphorylated, are stabilized as a complex with the adaptor protein dpb and the tetrameric gins complex, forming a module that interacts with cdc the latter acts as a scaffold for this module, which is then competent to associate with the pre-rc and attract dna polymerase. a recent study shows that the end result is the tight association of cdc mcm and gins with origins, allowing the unwinding of dna and processive replication by dna polymerase. this represents the essential role of cdk in stabilizing polymerase at the moving replication fork and switching the system from a pre-replicative state to a replicative one. from this point until late in mitosis, cdk levels remain high. this continued cdk activity prevents re-establishment of pre-rcs at origins that have already fired through a number of mechanisms. firstly, cdk phosphorylates cdc thus causing the scfcdc complex to target cdc to the proteasome for degradation. secondly, orc and orc are phosphorylated by cdk, with the phosphorylation of orc rendering it refractory to interaction with cdt, thereby preventing further mcm loading. finally, cdk facilitates the nuclear export of both mcm and cdt at different time points. just prior to initiation, cdt exits via a cdk-dependent mechanism, while mcm proteins fall off the dna upon fork termination and are then exported in a cdk-dependent manner. thus, while cdk initiates replication, it subsequently prevents pre-rc reassembly. this illustrates its dual role in triggering initiation through formation of cmg, then preventing re-initiation by inhibiting pre-rc reformation. mathematical modeling has been successfully used in the past to address various aspects of the cell cycle. early models did not incorporate specific biochemical mechanisms; they were hypothetical representations of periodic cellular activity. as the molecular mechanisms driving the cell cycle were revealed, models appeared that incorporated these findings. for s. cerevisiae in particular, multiple modeling approaches have been applied, based both on network descriptions and on specific molecular details such as gene expression and biochemical kinetics. some modeling efforts have been comprehensive, such as the tyson group’s ordinary differential equation -based models, while others address specific cell-cycle phenomena, such as the links between cell size and cycle progression. spiesser et al. developed a model of chromosomal replication, which reproduced the spatio-temporal replication profile of yeast chromosomes. origin firing was also described in wherein the authors used a stochastic model to describe these origin-specific features of replication. a recent report presented an ode-based model describing the initiation of dna replication, incorporating origin licensing, firing and the network of regulatory phosphorylation events. the model parameters were partly calibrated against experimental data, but largely selected through an optimization routine designed to attain an idealized function, resulting in a model that is particularly suited to exploring events specifically at the g1/s transition. here, we present a new model of the initiation of dna replication. in contrast to the work of brümmer et al., we took a ‘bottom-up’ approach and began by gathering in vivo data for precise protein levels at specific cell cycle time-points, then calibrated our model against these values. rather than limiting ourselves to the observation of firing near the g1/s transition and fitting to dna-specific replication profiles, we validated our model against the behaviour of the constituent protein complexes throughout the entire cell cycle. to facilitate the use of our model in a comprehensive description of the cell cycle, we designed it to integrate easily with the model of chen et al.. finally, we validated the model by comparing in silico predictions to experimental observations, using both our own knockdown experiments and results from the literature. the model presented here consolidates the known interactions between dna replication initiation proteins and the mechanisms that allow them to drive genome duplication. additionally, regulatory aspects of the system, which ensure that re-replication does not occur, have been modeled. the model’s behaviour provides a falsifiable hypothesis regarding the dynamics of dna replication initiation. furthermore, it accurately predicts the phenotypes of known experimental cell cycle mutants as well as those arising through in vivo perturbations to proteins in the network. because our model is constrained only by the fluctuating levels of replication factors, it provides a unique understanding of the kinetics governing the reactions between them. successful integration into a whole cell cycle model allows the initiation of dna replication to be explored in a broader quantitative context.  we began construction of our model by identifying the important players in replication initiation and establishing an interaction network, as shown in figure  after selecting appropriate descriptions of reaction kinetics, we generated an ode-based model and calibrated the model parameters to in vivo data. description of model components the model describes sixteen molecular species and depends on twenty-four parameters, which characterize the rates of seventeen biochemical processes. the model describes the following molecular species. rc: origin-bound orc rc2: origin-bound orc associated with cdc6 rc3: origin-bound orc associated with cdc with mcm loaded rc4: origin-bound orc, with mcm loaded rc5: origin-bound orc, with mcm loaded and dbf associated rc6: origin-bound orc, with mcm loaded and dbf and cdc associated rc7: origin-bound phosphorylated orc fork: the elongation fork, with mcm and cdc associated cdc6n: non-chromatin associated nuclear cdc6 dbf4n: non-chromatin associated nuclear dbf4 cdc45n: non-chromatin associated nuclear cdc45 mcmc: cytosolic mcm cdt1c: cytosolic cdt1 mcm·cdt1n: non-chromatin associated nuclear mcm bound to cdt1 cdt1n: non-chromatin associated nuclear cdt1 mcmn: non-chromatin associated nuclear mcm the mcm species corresponds to dimers of mcm2- heterohexamers, as two complexes are loaded at each origin. similarly, the cdc species corresponds to a dimer, as described by bowers et al.. we describe concentrations in units of molecules per cell. reaction events the seventeen processes that make up the model are shown in table  their rates depend on the species concentrations, the model parameters, and on two fixed, time-varying input functions describing the abundance of clb and of cdc in choosing reaction kinetics, we balanced the complexity of the model against its ability to adequately describe the behaviour of the overall system. we limited our description of initiation to the interactions between the pre-rc and replisome proteins that we found to be the essential core of the network. as a result, certain processes were combined into single events, some reactions were presumed irreversible, and only some reaction rates were presumed to have non-linear kinetics. except for rc phosphorylation states are not explicitly described, as we have no data for the individual phosphorylation events. this is acceptable for our purposes as the lumped function of cdk in each case is consistent with a scenario where the effect of cdk is proportional to its concentration. additionally, processes that involve multi-protein complexes are represented by a single member – one cmg complex stabilizes dna polymerase at each replication fork. of the three protein factors it is comprised of, cdc is limiting. although mcm is also included in the gins complex, we model both mcm and cdc as separate species. dbf represents the dbf4-cdc kinase complex and mcm represents the mcm2- helicase. although the protein factors cdc dbp sld sld and gins interact to facilitate formation of the pre-initiation complex at origins, we model only cdc which is the limiting factor in the cmg complex; ultimately the number of forks fired is dependent on the cdc concentration. we take dbf which is the limiting regulatory subunit of cdc as representative of active ddk, which is one of the limiting factors in replication initiation. mcm is used to represent mcm complexes; the mcm concentration has been reported to approximate the number of total complexes per cell in an asynchronous population. the replication complexes in our models exist only on chromatin and therefore represent the activity of these proteins at the dna as opposed to soluble complexes. the network shown in figure  includes both reversible and irreversible reactions as indicated. association/dissociation reactions are considered reversible, in accordance with a dynamic pre-rc/pre-ic loading mechanism as described above. in most cases, phosphorylation events are modeled as irreversible, in the absence of identified countervailing enzymes. we found that it was sufficient to describe most reaction rates by mass action kinetics. in cases where saturation occurs, we employed michaelis-menten kinetics. to simplify the description of the phosphorylation of orc by cdk, we do not describe phosporylation and dephosphorylation explicitly, but combine them into a single dephosphorylation event whose rate is inversely proportional to the level of cdk. we introduce cooperativity in this mechanism to account for multiple phosphorylation events or an additional inhibitory cdk-orc binding mechanism. the establishment of replication complexes in our model reflects the sequential binding of proteins that constitute the pre-replicative complex. in some cases the association and dissociation of pre-rc components is reversible. we treat the loading and maintenance of mcm2- helicase complexes as a dynamic process, which is dependent on the concentrations of the factors orc, cdc cdt and mcm2- itself. a mechanistic model for the dynamic assembly of pre-rcs was first described by the bell lab. the requirement of pre-rc factors for maintenance of helicase-loaded origins in late g has been further demonstrated by work from these researchers as well as our group. network and differential equations referring to figure  and table  the dynamics of the system are described as: drc2dt=v5−v7dforkdt=v14−v15drc3dt=v7−v10dcdc6ndt=v1+v10−v5−v2drc4dt=v10−v12−v11ddbf4ndt=v3+v14−v12−v4drc5dt=v12−v13dcdt1ndt=v7−v8+v9drc6dt=v13−v14dmcmndt=v15+v11+v9−v16drc7dt=v14−v17dmcm·cdt1ndt=v6−v9−v the remaining state variables are constrained by the following conservations: rc1 = rctotal− rc2 − rc3 − rc4 − rc5 − rc6 – rc7cdt1c= cdt1total− cdt1nmcmc= mcmtotal− mcmn− rc3 − rc4 − rc5 − rc6 − fork − mcm·cdt1cdc45n= cdc45total− rc6 − fork where rctotal, cdt1total, mcmtotal, and cdc45total are the fixed total number of origins, cdt molecules, mcm2- complexes, and cdc dimers, respectively. these four factors have been shown to be present at constant levels throughout the cell cycle. the value for rctotal used in our model is as described in. system inputs the biological network responsible for the initiation of dna replication does not oscillate autonomously; it displays periodic behaviour when driven by periodic signals from the cell cycle. likewise, our model displays oscillations only when driven by periodic forcing input. in order to facilitate the combination of our model with the cell cycle model of chen et al., we used the simulated profiles of clb and cdc from their model as periodic inputs to ours. cdc mediates the degradation of dbf. clb is responsible for cdc degradation, loading of cdc, nuclear export of free mcm and cdt, and phosphorylation of orc and orc. we converted the time-varying profiles of clb from the chen et al. model to molecules-per-cell units using the genome-wide gfp tagging and localization experiments described in. the profile of cdc was similarly obtained by scaling to cellular abundance levels reported in another study – while cdc has been determined to peak at copies in a haploid cell, the functional apccdc level can be estimated by considering the apc cyclosome subunit cdc. this was reported in different studies to be 593 mol/cell in an asynchronous population and at its maximal value of 750 mol/cell in metaphase. data acquisition data for cdc and cdc levels were obtained from individual isogenic strains in which the open reading frame of the corresponding gene was fused to a sequence encoding a 13myc epitope tag. in figure  a representative western blot for cdc45-myc is shown, with the corresponding facs analysis. the levels of mcm were determined using an anti-mcm antibody. in each of our time course experiments, cells were first arrested in late g phase with the mating pheromone α-factor and then released synchronously into the cell cycle, as described in methods. from the literature, we used time course data for chromatin-bound and soluble dbf and mcm from pasero et al. and quantitation of the nuclear fraction of cdt from tanaka and diffley. in order to convert relative measures of protein abundance to molecule-per-cell numbers, we used scaling factors obtained from the database provided by ghaemmaghami et al.. the data is shown along with a best-fit simulation in figure  raw timecourse data can be found in additional files and our calculation of molecule-per-cell estimates is demonstrated in additional file the fits in figure  represent our best solution to a trade-off between quality of fit and model complexity. we explored the effect of adding additional species and parameters. these additional features could, in some cases, provide minor improvements to the fit, but our confidence in the parameter estimates suffered as the complexity of the model grew. while our model reports proteins and complex concentrations in absolute units of molecules/cell, our accuracy regarding these values is limited to by the literature-reported cellular abundances for the various protein factors. we have used discretion when inconsistencies arise, choosing the reported values that are best justified by multiple studies. the network dynamics are based on the relative protein abundances over the time-course. consequently, changes in the global protein level for a particular factor do not affect the dynamics; instead, they impact our predictions of absolute molecule/cell counts and parameter values. this is an important consideration regarding the conversion of densitometry readings to absolute values, as the overall levels are ultimately determined by a literature-derived scaling factor. future experiments will undoubtedly result in improved estimates of protein abundances. these can be easily incorporated into the model by scaling the parameter values accordingly. parameter calibration the model parameters were calibrated using a weighted least-squares comparison with the data described above. we used a combination of global optimization and local search to find the best-fit parameter set shown in table  the table also shows the percent error associated with each parameter estimate. the percent error is the relative size of a 95% confidence interval for the estimate, calculated via the fisher information matrix and the cramer-rao bound. the percentage errors show that some parameters are estimated with high confidence while others are represented with less accuracy. these values were used to solve the odes in our consensus model. the percentage error for each parameter is indicated. parameter values that were well constrained by the data include those specifying the rates of production, degradation and association of cdc and dbf as well as the rate of origin firing. this reflects the strong reliability of our data for these two protein factors as well as for the proteins that form the replication complex that gives rise to active forks. parameters values in which we have low confidence include those that govern the loading of mcm by cdt, cdc dissociation from rc, and the phosphorylation of orc. the reversible dissociation of cdc is needed to accurately fit the data and there is no evidence suggesting that orc-cdc binding is irreversible. nevertheless, it is clear that experimental observations specific to this process are required to more precisely estimate this parameter value. the reaction whereby the mcm·cdt species loads the mcm complex is extremely transient. provided parameter k is sufficiently large, the kinetics of this reaction will be rapid enough to fit the data. consequently, the data cannot support a precise estimate of the parameter value. this observation suggests that mcm loading is an extremely rapid biochemical step in pre-rc assembly. it may point to a role for cdt in repeatedly targeting mcm complexes to origins throughout g such a phenomenon is consistent with the requirement for a dynamic loading mechanism that ensures pre-rc fidelity up until the g1/s transition. finally, the phosphorylation of orc contributes to the prevention of repeated origin firing. however, this mechanism has not been well characterized, and our data is unable to accurately constrain the specifics of this process. the kinetic rates in this network have not been the subject of prior experiments, but previous reports of protein half-lives are consistent with our predicted parameter values. drury et al. estimated that cdc is reduced below the point of detection within minutes of s-phase entry, corresponding to a half-life no longer than  min. similarly, cheng et al. reported that dbf is reduced below visible levels within minutes by the apc-dependent pathway, indicating a half-life no longer than 3 min. our model-based predictions of degradation rates correspond to half-lives of 1 min. and  min. for cdc and dbf respectively, in good agreement with these earlier findings. figure  shows the simulated model behaviour for the best-fit parameter set. some replication complex species – rc rc and rc – are extremely transient. their low levels of abundance are shown separately from other rcs, on an appropriate scale. simulations were carried out in matlab. perturbations our initial explorations of the model revealed that the network’s behaviour is particularly sensitive to the abundance of dbf and cdc and relatively insensitive to the level of cdt we investigated the effects of perturbations by simulating reductions in dbf cdt and cdc in our model. when the cdc production rate was reduced to 10% of its nominal value, persistence of the rc complex was observed. similarly, when the dbf production rate is reduced by the same relative amount, an accumulation of rc occurs. in both cases, the perturbation interferes with pre-initiation complex assembly and blocks the system at the nearest previous persistent rc state; rc is not persistent since the unloading of mcm causes a rapid transition back to rc it is worth noting that because mcm can dissociate from orc, rc represents a complex containing mcms that will be functionally incorporated into replication forks as opposed to those that loosely associate with origins. because the timing of our model is fixed, the various state concentrations indicate the progression from licensing to firing. a reduction in the fork species compared to the wild-type case suggests a slow-down in s-phase because fewer origins are firing within the prescribed time. using the peak abundance of the fork species as a measure of replicative efficiency, we saw significant reductions in both simulated knock-downs. conversely when we simulated the reduction of cdt abundance to 10% of nominal values, origin firing was only reduced by 23%, suggesting that the network is relatively refractory to depletion of cdt. additional file shows the levels of the various model components for each perturbation. to investigate the accuracy of these mathematical predictions, we carried out corresponding wet lab depletion experiments. reducing dbf or cdc concentrations in yeast cells to roughly 90% below normal endogenous levels resulted in a rapid g phase arrest, evident after hours of depletion, as judged by facs analysis indicating the accumulation of cells with 1 c dna content. in contrast, a corresponding depletion of cdt had no appreciable effect, and dna replication defects were still only minimally evident after hours of further reduction. thus, our in silico simulations using our nominal parameter set were predictive of in vivo perturbations. these experiments were used to validate the model; they were not used for calibration. the insensitivity to perturbations in cdt levels is consistent with its apparent excess relative to origins, although the number of cdt molecules that act at each origin has not yet been characterized. moreover, the mechanism by which cdt aids in recruiting the helicase molecules to pre-rcs is extremely transient. while many factors are limiting, the system appears to be highly sensitive to the levels of cdc due to its low abundance relative to mcm and cdt even a moderate depletion of cdc significantly alters the dynamics of pre-rc loading. the same is true for dbf although in this case its role in activation of the cdc kinase renders the system highly sensitive to its concentration; firing cannot occur without the dbf4-cdc complex. since dbf is, like cdc limiting, flow through the network is blocked when the kinase does not reach a threshold level. additionally at limiting levels, the number of replication forks produced by our model is significantly reduced, consistent with in vivo reports from the literature showing a lengthening of s-phase. further validation of our model comes from comparison with additional in vivo experiments reported in the literature. jones et al. showed that the interaction between the mcm complex and dbf was reduced to half its wild-type level when a dbf domain that binds mcm was mutated, impairing s-phase progression. we mimicked this effect by reducing the rate of association of dbf with rc by 50%, leading to a similar result. similarly, zou et al. reported that the cdc45- mutant shows an aberrant growth phenotype at the non-permissive temperature. this is thought to be due to a disruption of cdc45’s ability to interact with mcm and orc. as shown in figure 7c, by reducing the rate of cdc interaction with rc by 50%, a marked reduction in the peak abundance of the fork species results, indicative of a slower s-phase, as observed when the mutant was grown at the non-permissive temperature. the actual reduction in cdc45’s association with the pre-rc due to conformational changes in the mutant might be even more pronounced than a 2-fold reduction. in any case, our simulation is consistent with cdc45’s origin-initiation role being compromised by impairing its ability to interact with its ligands to form the cmg complex. linking our dna initiation model to a previously established cell cycle model our model of the initiation of dna replication only displays oscillatory behaviour when forced with periodic signals from the cell cycle. by choosing to incorporate signals that correspond to species in the cell cycle model of chen et al., we were able to merge the two models in a straightforward manner. in the chen model, the initiation of dna replication is represented by a single lumped state variable, called ori. at the beginning of the cell cycle, ori has value zero. its rate of growth depends linearly on clb when it reaches a threshold value, dna synthesis is presumed to have begun, and triggers an increase in the value of the parameter kmad leading to an inactivation of cdc which is required for mitotic exit. this mad2-dependent inhibition of cdc represents the spindle assembly checkpoint, ensuring that cells with replicated dna do not complete mitosis without properly aligning the chromosomes. when chromosomes have properly aligned on the metaphase plane kmad drops and cdc promoted exit from mitosis. in our model, the level of dna synthesis is represented by the fork species. to merge the two models, we removed the ori state from the chen model, and instead used the fork species to trigger the change in kmad as detailed in methods. besides “closing the loop” between the two models by incorporating two-way inter-model signalling, we had to deal with a single shared species: both models describe the dynamics of cdc we arrived at a merged description of cdc behaviour by incorporating the dynamics of recognition complex association and dissociation into the chen model’s formulation of cdc behaviour. the resulting combined model behaves only marginally differently from either model in isolation, as shown in additional file 6: figure s and additional file 7: figure s7). replication metrics of brümmer et al the recently published model of brümmer et al. also describes the network responsible for the initiation of dna replication. the free parameters of that model were chosen by a combination of fitting and optimization. the authors used literature-derived data to fix of the kinetic parameters. the remaining free parameters were not fit to data, but were selected through a procedure that optimized the coherence of origin firing and minimized re-replication. while it is impossible to assess the accuracy of the parameter values obtained from this procedure, the resulting idealized model provides a useful starting point for examining how the network structure constrains the system behaviour. the model of brümmer et al. focuses on early origin firing and so represents the mechanics of firing at the start of s-phase. in contrast, our model describes firing dynamics throughout s-phase in order to fit into the broader context of the cell cycle. the parameter set driving our system is not filtered to retain only those that produce replication dynamics consistent with coherent firing just at the g1/s boundary. rather, the parameters are specified by the actual cellular concentrations of the active protein factors generating replication forks. while both models incorporate the important role of cdk, brümmer et al. emphasize the multi-site phosphorylations of several factors involved in mechanisms that minimize potential re-replication. to this end, brümmer et al. employed a metric to assess re-replication. their idealized model exhibits re-replication events per cycle. applying the same measure to our bottom-up model yields re-replication events per cycle. the near-zero value obtained by brümmer et al. is close to their idealized target of zero. both estimates are consistent with the belief that re-replication occurs in wildtype cells, but at an extremely low rate. because the nature of the dephosphylation of orc remains uncharacterized, we use a conservative estimate of the number of orc phosphorylation sites. increasing this number by twofold, consistent with the number of cdk target residues on orc reduces re-replication to a value on the same order of magnitude as brümmer’s value. thus, both models effectively deal with representing control and prevention of rereplication.  while our model provides a sound description of the initiation of dna replication, a number of aspects of the network remain unresolved: for example, the kinetics of mcm loading, the mechanism by which cdk phosphorylates orc, and the details of the association of the gins complex, sld and sld while modeling cdc captures the events regarding cmg formation at origins, being its limiting factor, a future version of the model could better distinguish the initial cdc association at origins from subsequent cmg formation. while this has minimal effect on our network dynamics and no effect on our blocking of re-replication, it would provide a better resolution of events at origins just prior to the g1/s transition. incorporating timecourse experiment data for levels of a gins complex member would aid in this analysis. nevertheless, the assumptions made allow us to approximate the aforementioned processes, simplifying the network without losing information about system behaviour at the level we intend to model. with our nominal parameter set, we observe the system to behave as the ordered accumulation of proteins forming a loading complex at origins throughout the genome. activation by increasing s-cdk levels and the concentration of dbf increase linearly the number of replication forks established as a result of successive forward transition of the various replication complexes. it should be noted that as is found in vivo, not all origins fire as a consequence of being furnished with mcm-containing pre-rcs. replication is maximal at the g1/s transition, but continues into s-phase as origin firing is temporally spaced. this is thought to ensure sufficient time to address any defects in replication and is mediated by the limiting nature of one of the initiation activators, ddk. our in vivo perturbation of dbf levels reproduces this consequence and points to other system observations: cdc levels are intimately controlled by cdk levels to avoid re-replication, however this mechanism is tightly regulated such that clb levels rising too soon would prevent the assembly of the pre-rc in g a feature of the system well documented. additionally, cdt appears to act catalytically rather than stoichiometrically given the system is relatively impervious to a reduction of this factor to % of its wildtype level. this might play into its role in chaperoning mcm2- hexamers to origins, where they are loaded subsequently leading to the release of cdt which may then be recycled to mediate the loading of other mcms. this aspect of the system has not yet been investigated experimentally and would be of future interest. many human orthologs of the yeast proteins described in our network have been associated with cellular pathologies. our model is specific to the replication machinery in budding yeast, but the mechanisms driving this process are highly conserved throughout eukarya. efforts to develop an analogous model in mammalian cells would be useful in understanding and dissecting cell proliferation in humans. a number of models of the mammalian cell cycle have been proposed. incorrect pre-rc formation has been linked to impaired dna damage repair pathways in humans, while both orc and members of the mcm2- complex have been shown to be reliable cancer biomarkers. recent work by bicknell et al. has shown that point mutations in the human orc orc orc cdt and cdc genes are associated with meier-gorlin syndrome, a form of primordial dwarfism, and several of these mutations were determined to interfere with proper pre-rc formation. these findings highlight the potential utility of in silico mammalian models in further exploring the molecular basis of such disorders. given that our model shows good predictive capability, it serves not only as an informational tool for yeast biology, but also as a proof of principle for higher order system models. despite the requirement for a mammalian model to comprehensively verify specific mechanisms, the system of dna replication initiation is conserved well enough that perturbations to proteins such as those described above can, in fact be preliminarily examined. although previously established replication models consider the ordered timing of origin firing based on genomic replication profiles, our goal was to represent the temporal organization of origin firing as a function of the concentration of active replication species. a focus on using real protein levels as a determinant of replication dynamics is a novel approach. when used in concert with models describing genome-level origin characteristics and/or combining the findings with models exploring other cell-cycle modules, future efforts will generate a well-rounded picture of dna replication initiation. 
finally, this model can serve as a guide in efforts to generate a comprehensive model of the mammalian cell cycle in order to explore cancer-related phenotypes. the model behaviour conforms to perturbation trials previously reported in the literature, and accurately predicts the results of our own knockdown experiments. eukaryotic cell proliferation involves dna replication, a tightly regulated process mediated by a multitude of protein factors. the model was validated against quantified levels of protein factors over a range of cell cycle timepoints. this study establishes a robust model of the processes driving dna replication initiation. we have developed an ordinary differential equation-based model of the protein-protein interaction network describing replication initiation. furthermore, we successfully incorporated our replication initiation model into an established model of the entire yeast cell cycle, thus providing a comprehensive description of these processes.
 eukaryotic cell proliferation involves dna replication, a tightly regulated process mediated by a multitude of protein factors. in budding yeast, the initiation of replication is facilitated by the heterohexameric origin recognition complex. orc binds to specific origins of replication and then serves as a scaffold for the recruitment of other factors such as cdt cdc the mcm2- complex, cdc and the dbf4-cdc kinase complex. while many of the mechanisms controlling these associations are well documented, mathematical models are needed to explore the network’s dynamic behaviour. we have developed an ordinary differential equation-based model of the protein-protein interaction network describing replication initiation.  the model was validated against quantified levels of protein factors over a range of cell cycle timepoints. using chromatin extracts from synchronized saccharomyces cerevisiae cell cultures, we were able to monitor the in vivo fluctuations of several of the aforementioned proteins, with additional data obtained from the literature. the model behaviour conforms to perturbation trials previously reported in the literature, and accurately predicts the results of our own knockdown experiments. furthermore, we successfully incorporated our replication initiation model into an established model of the entire yeast cell cycle, thus providing a comprehensive description of these processes.  this study establishes a robust model of the processes driving dna replication initiation. the model was validated against observed cell concentrations of the driving factors, and characterizes the interactions between factors implicated in eukaryotic dna replication. finally, this model can serve as a guide in efforts to generate a comprehensive model of the mammalian cell cycle in order to explore cancer-related phenotypes. 