Bifurcation in cell cycle dynamics regulated by p53

We study the regulating mechanism of p53 on the properties of cell cycle dynamics in the light of the proposed model of interacting p53 and cell cycle networks via p53. Irradiation (IR) introduce to p53 compel p53 dynamics to suffer different phases, namely oscillating and oscillation death (stabilized) phases. The IR induced p53 dynamics undergo collapse of oscillation with collapse time \Delta t which depends on IR strength. The stress p53 via IR drive cell cycle molecular species MPF and cyclin dynamics to different states, namely, oscillation death, oscillations of periods, chaotic and sustain oscillation in their bifurcation diagram. We predict that there could be a critical \Delta t_c induced by p53 via IR_c, where, if \Delta t<\Delta t_c the cell cycle may come back to normal state, otherwise it will go to cell cycle arrest (apoptosis).


Introduction
p53 is well known for its abnormally long stability in response to the stress available against genomic integrity [1]. It conglomerated with its negative inhibitor MDM2 in the nucleus due to their strong interaction [2]. When the cell is in stress condition (due to irradiation, stress inducer molecule etc), p53 concentration level rises which leads to cell cycle arrest until repair or doctoring takes place of the impaired DNA. If the repair is not successful the system goes towards the apoptosis [3][4][5][6]. The transcriptional ability of the p53 is kept under controlled level at normal state due to its negative feedback interaction with MDM2 [7]. The hyperbolized concentration of MDM2 helps in degradation of the p53 protein because of its E3-ligase activity, causing adherence of ubiquitin to the lysine rich C-terminal of the p53 molecule [8][9][10]. Introduction of stress in the system is sensed by the activation of ARF protein, initially situated in nucleolar region in the form of nucleophosmin shifts to the nucleoplasm in its independent and active cast, to mark MDM2 for its degradation, thus assisting the p53 stability [11][12][13]. Triggering of p53 in response to stress leads to the expression of several downstream genes apart from the MDM2. p21 protein is one of the most important proteins which is found to be expressed due to p53 accumulation in the cell [14]. p53 acts as a transcription factor for p21. It is also reported that p21 expression is directly proportional to the level of p53 in the system [15]. The role of p21 in model is described briefly as follows. The main component of p53 -MDM2 regulatory network is the feedback loop between p53 and MDM2 [37]. p53 and MDM2 interact to form p53 -MDM2 complex with a rate constant k 17 [37], followed by dissociation of the complex to its respective components with a rate constant k 18 [41,42]. The transcription rate of MDM2 gene to its mRNA (MDM2 -mRNA) is takes place with rate constant k 20 , followed by translation of MDM2 -mRNA to MDM2 with a rate constant k 22 [37,43] and its (MDM2 -mRNA) self-degradation with a rate constant k 21 [44]. The ubiqiutination of MDM2 protein occurs with rate constant k 23 . The p53 synthesis is taken placed with a rate constant k 16 , and gets ubiquitinized at the rate constant k 19 [43]. The DNA damage in system is introduced via irradiation with an estimated rate constant of k 24 [37]. Irradiation is reported to be a major cause of DNA damage. The severity of the DNA damage is depended on the dose of exposure of irradiation [34]. The repair of the DNA damage is then occurred at a rate constant k 25 [45,46]. The activation of ARF due to DNA damage takes place at a rate constant k 26 [37]. Further, this activated form of ARF interacts with MDM2 protein and forms ARF -MDM2 complex with a rate constant k 27 [47]. The degradation of ARF protein is reported to occur at a rate constant k 28 [48]. ARF based degradation of the MDM2 takes place by getting targeted to the complex via proteasome recognition with a much faster rate constant k 29 than individual degradation rates [49]. The p53, being a transcription promoting factor for many of the proteins, also transcribes the gene responsible for the manufacture of p21 protein with a rate constant k 30 as presumed by the approximations made to attain the appropriate oscillations and arrests [14,50,51].
The p21 protein is capable of making complex with the cell division promoting factor MPF with a rate constant k 31 [16,18] with respect to the amount of concerned molecules present in the system [19]. Then the inhibition of MPF, or more appropriately G2 associated Cyclin -Cdk complex, by p21 is approximated with a rate constant k 32 [50,52]. p21 then gets degraded by the virtue of its half-life in the system with a rate constant k 33 [18,23]. The cyclin is assumed to translate at the rate constant k 1 [53]. Further, ubiquitin dependent cyclin degradation or protease independent degradation of the cyclin is reported to happen at a rate constant k Ã [54]. The degradation of the cyclin due to effect of protease activation during cyclin accumulation and interaction between inactive form of MPF with cyclin takes place with a rate constant k 4 [38]. Formation of activated form of MPF (M) occurs due to interaction of cyclin with inactive MPF (M Ã ) with a rate constant k ÃÃ [38,[55][56][57]. Further this activated form of MPF (M) converts to inactivated form (M Ã )with a rate constant k ÃÃÃ [55,57]. The activated form of MPF(M) interact with inactive protease(X Ã ) to generate activated form of protease (X) with a rate constant k ÃÃÃÃ [38,58,59]. The generation of activated form of cyclin protease (X) occurs due to interaction of cyclin protease with inactive X Ã with a rate k ÃÃÃÃÃ [38,55]. The activated form of protease (X) can convert into inactive form (X Ã ) with a rate constant k ÃÃÃÃÃÃ [38,57]. In Fig 1, The blue dots indicates creation and black dots indicates decay of the respective molecular species. The lists of molecular species and biochemical reaction channels involved in this proposed model are listed in Tables 1 and 2 respectively.
The biochemical reaction network shown in Fig 1 are represented by the twenty five reaction channels listed in Table 2, which are participated by thirteen molecular species (Table 1) defined by a vector at any instant of time t, The time evolution of these variables can be translated from the twenty five reaction channels into the following set of nonlinear ordinary differential equations (ODE) based on Mass action law of chemical kinetics, where, the expressions for M Ã and X Ã in the Fig 1 are given by, M Ã = 1-x 10 and X Ã = 1-x 11 . The set of coupled ODEs can be solved using Runge Kutta method of standard numerical integration algorithm [60].

Results and Discussion
We numerically simulate the proposed model and the results demonstrate new phenomena in bifurcation diagram which may be significant to correlate with various experimental situations. The interaction of p53 regulatory network and cell cycle network highlights different form of signal processing between non-identical networks which could be the way of regulating one another. We study the complicated way of this interaction in order to understand some of the basic mechanisms of network interaction.

Dynamics of p53 driven by irradiation
We first present the spatio-temporal behaviour of p53 upon exposure of irradiation in Fig 2. The p53 dynamics maintains minimum concentration level at IR = 0 (normal condition). As IR dose increases p53 start showing damped oscillatory behaviour (Fig 2 second and third panels) indicating stressed behaviour of p53. The increase in IR dose induces increase in time to attain stability of p53 dynamics (amplitude death) indicating increase in unstability of p53 dynamics (Fig 2 third panel). This could be due to the fact that the increase in IR dose may cause high DNA damage leading to more stress in p53. However, if the IR dose is comparatively strong (IR = 5), the damage within the DNA is also high which may cause the collapse of the p53 oscillatory behaviour (Fig 2 fourth panel) and then repaired back the DNA damage to come back to p53 oscillatory condition. We also found that the time of collapse (Δt) increases as IR dose increases (Fig 2 fifth panel) and it becomes difficult to repair back the DNA damage. In general p53 will collapse forever and will not be recovered back if Δt ! 1 (probable case of apoptosis). However, in real situation, one probably can define a critical Δt c such that, if ΔthΔt c , p53 could come back after DNA repair, and otherwise it will go to apoptosis. Nevertheless, it is very difficult to find out this Δt c .
Similarly, we also present the plots of temporal variation of the concentration of MDM2 due to exposure of irradiation in right panels of Fig 2. We observed similar kind of behaviour as  [38][39][40] [16,50] obtained in case of p53 protein dynamics. This is probably due to intercorrelation between p53 and MDM2 in the system via feedback mechanism. It is also noted that corresponding variations in the behaviours of both p53 and MDM2 (as observed by comparing panels in Fig 2) are due to their positive as well as negative feedback regulations prescribed to them.

Phase diagram of p53 compelled by IR
We simulate the maxima of p53 amplitudes after removing the transients as a function of IR (Fig 3) to capture the different phases namely oscillation and oscillation death regimes. The behaviour of Δt as a function of IR follows the functional form Dt ¼ A Bþe ÀIR with the values of A = 6778 and B = 0.00887 (fitting values of the function to the data) (Fig 3 inset). The separation between two phases oscillation death and oscillating regimes are clearly visible after the IR * 3.45 and Δt increases as IR increases.
Generally as Δt ! 1 when IR ! 1, but numerically we approximately found that after IR = R c * 11 Δt become Δt c * 79 hours and becomes constant (Fig 3 inset). This means that for any ΔthΔt c , the p53 can able to recover back to normal stable state by repairing DNA damage, otherwise, the system can't able to come back to normal state, but will go to apoptosis.

Bifurcation in Cyclin regulated by p53
Since cell cycle and p53 regulatory networks are interacted through p21 (Fig 1), the temporal behaviour of cyclin can be regulated by p53 via IR and p21. When IR = 0, the two networks work in normal condition, leaving p53 dynamics at low level (stabilized state) (Fig 2 upper  panel) and sustain oscillation in cyclin dynamics (Fig 4 upper left panel). As IR increases, p53 will get activated through DNA damage giving oscillatory behaviour affecting the dynamics of cyclin. When IR = 0.1, the cyclin dynamics shows chaotic behaviour upto t = 145 hours, and then the dynamics becomes sustain oscillation (Fig 4 second left panel and upper right panel). The chaotic behaviour in cyclin dynamics could due to the sudden activation in p53 dynamics due to IR irradiation. Now as IR increases (IR = 0.5), we get various situations in the cyclin dynamics, namely, the emergence of period two (for t * [10-40] hours), period 3 (for t * [40-85] hours), chaotic regime (for t * [85-175] hours) and sustain oscillation regime (for ti175 hours) (Fig 4 second  right upper panel). Further, as IR increases the emergence of oscillation death regime started to exist in the cyclin dynamics (Fig 4 fourth right panel onwards) and the oscillation death regime become larger. Further increase in IR compels the period 2 and 3 regimes to vanish after some value of IR (IRi9) and the chaotic regime becomes larger.
The perturbation induced by p53 through IR to the cyclin via p21 clearly induces cyclin dynamics to various states shown by the bifurcation diagram (Fig 4 right panels). We also notice that as one decrease or increase to cross over to sustain oscillation, the state just before it is chaotic regime. The emergence of oscillation death regime starts from IRi3 and then switches to sustain oscillation after sometime. This oscillation death regime corresponds to the collapse time due to strong sudden DNA damage. Once the DNA damage is recovered it comes back to sustain oscillation. If the IR is very large then oscillation death regime is large enough that DNA damage can not be repaired back halting the cell cycle permanantly and goes to apoptosis.

Dynamics of MPF regulated by p53
We present the temporal behaviour of MPF regulated by p53 as a function of IR (Fig 5) which induces at different states in MPF shown by bifurcation diagrams. The impact upon the MPF due to p53 via IR is not a direct phenomenon but through p21 molecule in the network. Various studies reported that p21 directly interact with cyclin dependent kinases, which has very important role in the formation of maturation promoting factor (MPF). The interaction of p21 with cdk leads to less availability of cdk due to the formation of MPF. Moreover, various experimental results also reported that p21 directly interacts with MPF [16,18]. It is observed that an IR = 0, the MPF dynamics shows sustain oscillatory behaviour indicating no impact of p53. Further, as IR dose increases the oscillatory behaviour of MPF is abruptly changed inducing different states of MPF as we obtained in the case of cyclin. The increases in IR dose induce different states oscillation death, period 1, 2, 3, chaotic and sustain oscillation regimes indicated by the bifurcation diagram for various IR values. Moreover, as IR increases the width of oscillation death [16] regime also increases and if IR is not strong enough the DNA can able to repair back otherwise the system will go to apoptosis.

Bifurcation in MPF and Cyclin
We study the regulation of cell cycle dynamics by p53 via IR. The maxima values of MPF (MPF M ) and cyclin (Cyc M ) as a function of IR are calculated for a range of time in the range [0, 50] hours (Fig 6). It is observed that for low IR dose, MPF M exhibits chaotic behaviour. Howeover, if IR dose is comparatively high, MPF M becomes almost constant. If the value of IR is  Similarly, one can also observe the IR dependent maxima of cyclin Cyc M in the bifurcation diagram (Fig 6 lower panel)

Conclusion
We study the way how p53, one of the largest hubs in cellular network, regulates and controls cell cycle dynamics. We studied the behaviour of different molecules which are highly involved in the checking of cell cycle at G2 phase driven by p53 via IR. The simulation results of the model provided us to understand the biological phenomenon and mechanism of cell cycle arrest due to DNA damage faced by the cell due to the irradiation. The results we got are closely in agreement with the previous experimental reports [16,17]. Our study suggests that the temporal dynamics of molecular species involved in cell cycle, considered in the model, are controlled by p53. The role of p21 protein in the delay of G2 phase was considered as a cross-talk between p53 regulatory network and cell cycle. The sudden irradiation to the system with high dose induces collapse of the system due to DNA damage, leading to cell cycle arrest. The cell cycle is resumed again to normal situation by repairing back the DNA damage. Moreover, the time of recovery from cell cycle arrest and then resumption of oscillation depends on the amount of dose of IR exposed to the system.
During the process of regulation of cell cyle by p53 via IR we observed the emergence of different periods (1, 2, 3 etc) in the bifurcation diagram of oscillatory dynamics of cell cycle variables (MPF M and Cyc M ) which may have various information of certain biological significance. Further, the dynamics of these variables switched to various states, namely, chaotic, oscillation death (stabilized state), bifurcating to various periods of oscillation and sustain oscillation states during the process of time evolution. These states could be the different phases of the variables to self-recover back to its normal condition from the sudden stress given to the system. However, how these complicated states are used by the system dynamics when the system is perturbed need to be investigated further.
The study also demonstrates the mechanism of cell cycle arrest induced by perturbed p53 via IR indicated by collapse of the oscillation (oscillation death) for certain interval of time (Δt). This collapse time is a function of strength of the perturbation imparted to the system. Our study shows that there is a minimum value of IR = R c , below which the system comes back to its normal state, otherwise the system will go to apoptosis. Our findings will probably be useful for the further study on the impact of p53 on cell cycle checking at G2 phase and related dynamics.