^{1,2}, Hualin Shi

^{1}, Haidong Feng

^{2}and Jin Wang

^{2,3,a)}

### Abstract

The global stability of dynamical systems and networks is still challenging to study. We developed a landscape and flux framework to explore the global stability. The potential landscape is directly linked to the steady state probability distribution of the non-equilibrium dynamical systems which can be used to study the global stability. The steady state probability flux together with the landscape gradient determines the dynamics of the system. The non-zero probability flux implies the breaking down of the detailed balance which is a quantitative signature of the systems being in non-equilibrium states. We investigated the dynamics of several systems from monostability to limit cycle and explored the microscopic origin of the probability flux. We discovered that the origin of the probability flux is due to the non-equilibrium conditions on the concentrations resulting energy input acting like non-equilibrium pump or battery to the system. Another interesting behavior we uncovered is that the probabilistic flux is closely related to the steady state deterministic chemical flux. For the monostable model of the kinetic cycle, the analytical expression of the probabilistic flux is directly related to the deterministic flux, and the later is directly generated by the chemical potential difference from the adenosine triphosphate (ATP) hydrolysis. For the limit cycle of the reversible Schnakenberg model, we also show that the probabilistic flux is correlated to the chemical driving force, as well as the deterministic effective flux. Furthermore, we study the phase coherence of the stochastic oscillation against the energy pump, and argue that larger non-equilibrium pump results faster flux and higher coherence. This leads to higher robustness of the biological oscillations. We also uncovered how fluctuations influence the coherence of the oscillations in two steps: (1) The mild fluctuations influence the coherence of the system mainly through the probability flux while maintaining the regular landscape topography. (2) The larger fluctuations lead to flat landscape and the complete loss of the stability of the whole system.

L.X. thanks Professor Weimou Zheng for helpful discussions. The work was supported by National Science Foundation and National Basic Research Program of China (Grant No. 2007CB814800).

I. INTRODUCTION

II. THEORY ON NON-EQUILIBRIUM POTENTIAL AND FLUX

A. Probabilistic descriptions in non-equilibrium systems

B. The potential landscapes and probabilistic flux for the non-equilibrium dynamical systems and networks

III. RESULTS AND DISCUSSIONS

A. Flux and its origin in a monostable network

1. Deterministic flux of the kinetic cycle

2. Probabilistic flux of the kinetic cycle

B. Flux and the potential landscape in the stable limit cycle

1. Limit cycle in the reversible Schnakenberg model

2. The effective variant form of the flux in the oscillation model

3. Chemical force pushes the flux flowing on the energy landscape

4. Energy pump as the source of global coherence and stability of the oscillations

IV. CONCLUSIONS

### Key Topics

- Fokker Planck equation
- 30.0
- Coherence
- 28.0
- Chemical potential
- 21.0
- Probability theory
- 21.0
- Chemical reactions
- 17.0

##### C12

## Figures

Kinetic cycles of the simplified three species enzyme kinetics. Comparing to the closed system of cycle (a), cycle (b) brings in two more substrates *D* and *E* which can break the detailed balance to generate a non-equilibrium steady state flux, with the reaction: by keeping a non-equilibrium concentration ratio [*D*]/[*E*]. By neglecting the fluctuation of concentrations of *D* and *E*, cycle (b) can be represented by cycle (a) in terms of the pseudo-first-order rate constants: and .

Kinetic cycles of the simplified three species enzyme kinetics. Comparing to the closed system of cycle (a), cycle (b) brings in two more substrates *D* and *E* which can break the detailed balance to generate a non-equilibrium steady state flux, with the reaction: by keeping a non-equilibrium concentration ratio [*D*]/[*E*]. By neglecting the fluctuation of concentrations of *D* and *E*, cycle (b) can be represented by cycle (a) in terms of the pseudo-first-order rate constants: and .

Schematic diagram of the three-species kinetic cycle

Schematic diagram of the three-species kinetic cycle

Scheme for the reversible Schnakenberg model. The dash curve in the box uncovers the self-catalytic mechanism of species *X*, and the dash box shows the key part that generates oscillation while *B* and *A* are energy input and output through the non-equilibrium concentrations.

Scheme for the reversible Schnakenberg model. The dash curve in the box uncovers the self-catalytic mechanism of species *X*, and the dash box shows the key part that generates oscillation while *B* and *A* are energy input and output through the non-equilibrium concentrations.

The stationary solution *y* ^{0} of Eqs. (31) in three different regions of energy pump *B*/*A*. The unstable fixed states are marked as the dash line within 0.1 < *B* < 0.9, that is 2 < *B*/*A* < 18, which refers to a limit cycle around the unstable fixed point. The star on the curve is noted as the Hopf bifurcation point. All parameters are set as in Eqs. (32).

The stationary solution *y* ^{0} of Eqs. (31) in three different regions of energy pump *B*/*A*. The unstable fixed states are marked as the dash line within 0.1 < *B* < 0.9, that is 2 < *B*/*A* < 18, which refers to a limit cycle around the unstable fixed point. The star on the curve is noted as the Hopf bifurcation point. All parameters are set as in Eqs. (32).

The amplitudes of the concentration of species *X* and *Y* vary with energy pump. Both curves are drawn inside the parameter range of the limit cycle, where *A* = 0.06 and *B* ∈ (0.15, 0.85). All other parameters are the same as Eqs. (32).

The amplitudes of the concentration of species *X* and *Y* vary with energy pump. Both curves are drawn inside the parameter range of the limit cycle, where *A* = 0.06 and *B* ∈ (0.15, 0.85). All other parameters are the same as Eqs. (32).

The non-equilibrium chemical potential difference Δ*G* can generate the non-zero deterministic flux *J* _{ SS }. The parameters we use are all shown in Eqs. (32), while B is chosen within (0, 0.1), and with these parameters this model shows one monostable steady state and leads to a steady state deterministic flux *J* _{ SS }.

The non-equilibrium chemical potential difference Δ*G* can generate the non-zero deterministic flux *J* _{ SS }. The parameters we use are all shown in Eqs. (32), while B is chosen within (0, 0.1), and with these parameters this model shows one monostable steady state and leads to a steady state deterministic flux *J* _{ SS }.

Flux and its direction flowing on the landscape and along the deterministic trajectory. The red arrows in (a), (c), (e), and (g) are vectors of probabilistic flux for different energy pump strengths, and (b), (d), (f), and (h) show their directions. The black solid line refers to the deterministic trajectory. The gradual change color from red to blue represents hill to valley of the potential landscape, which is similar as in Ref. 17.

Flux and its direction flowing on the landscape and along the deterministic trajectory. The red arrows in (a), (c), (e), and (g) are vectors of probabilistic flux for different energy pump strengths, and (b), (d), (f), and (h) show their directions. The black solid line refers to the deterministic trajectory. The gradual change color from red to blue represents hill to valley of the potential landscape, which is similar as in Ref. 17.

Sketch for two different kinds of decomposition of the probabilistic flux vectors. (a) shows how the vector is mapped to the deterministic trajectory and (b) shows that for all the flux vectors that pass through the cut line, their normal component will be summed into current. The start point *O* is set as the coordinate origin in the phase space.

Sketch for two different kinds of decomposition of the probabilistic flux vectors. (a) shows how the vector is mapped to the deterministic trajectory and (b) shows that for all the flux vectors that pass through the cut line, their normal component will be summed into current. The start point *O* is set as the coordinate origin in the phase space.

The relationship of deterministic effective flux and loop integration of probability flux. **Flux** _{ Ave } represents effective deterministic flux and **Flux** _{ IntTra } is loop averaged probability flux which can be calculated by the definition equation (33), and Δ*G* means chemical potential difference and also Gibbs free energy. (a) The two fluxes show the similar behavior with changing of Δ*G*. While the variation range Δ*G* ∈ (10.3, 12.05) is corresponding to the limit cycle range *B* ∈ (0.15, 0.85). (b) We can obtain the direct correlation between probability flux and averaged effective deterministic flux. In inset figure, we perform the linear regression analysis with origin software on two variables, in which *R* − *Square*(*COD*) = 0.53222.

The relationship of deterministic effective flux and loop integration of probability flux. **Flux** _{ Ave } represents effective deterministic flux and **Flux** _{ IntTra } is loop averaged probability flux which can be calculated by the definition equation (33), and Δ*G* means chemical potential difference and also Gibbs free energy. (a) The two fluxes show the similar behavior with changing of Δ*G*. While the variation range Δ*G* ∈ (10.3, 12.05) is corresponding to the limit cycle range *B* ∈ (0.15, 0.85). (b) We can obtain the direct correlation between probability flux and averaged effective deterministic flux. In inset figure, we perform the linear regression analysis with origin software on two variables, in which *R* − *Square*(*COD*) = 0.53222.

The total current passing through the Poincare line *OP* in Fig. 8(b) versus non-equilibrium potential difference Δ*G*. The inse figure comes from linear fitness with *R* − *Square*(*COD*) = 0.93679.

The total current passing through the Poincare line *OP* in Fig. 8(b) versus non-equilibrium potential difference Δ*G*. The inse figure comes from linear fitness with *R* − *Square*(*COD*) = 0.93679.

The relationship of chemical driving force and probability flux under loop integration along with the deterministic trajectory. **Flux** _{ IntTra } is the loop averaged probability flux. **∇μ** _{ IntTra } is loop averaged driving force that calculated in the same way as **Flux** _{ IntTra }. (a) The probabilistic flux and its driving force ∇μ have the similar behavior with respect to the Gibbs free energy Δ*G*, and the variation range of Δ*G* is same as Fig. 9. (b) The probability flux correlates with the chemical potential driving force. The inset figure is the linear regression analysis to two variables, and the coefficient of determination is *R* − *Square* = 0.87872, which means strong linear relation between two variables.

The relationship of chemical driving force and probability flux under loop integration along with the deterministic trajectory. **Flux** _{ IntTra } is the loop averaged probability flux. **∇μ** _{ IntTra } is loop averaged driving force that calculated in the same way as **Flux** _{ IntTra }. (a) The probabilistic flux and its driving force ∇μ have the similar behavior with respect to the Gibbs free energy Δ*G*, and the variation range of Δ*G* is same as Fig. 9. (b) The probability flux correlates with the chemical potential driving force. The inset figure is the linear regression analysis to two variables, and the coefficient of determination is *R* − *Square* = 0.87872, which means strong linear relation between two variables.

(a) The period of oscillation *X* varies with energy pump. (b) Comparison between the period and the deterministic effective flux **Flux** _{ Ave } (same as figures above).

(a) The period of oscillation *X* varies with energy pump. (b) Comparison between the period and the deterministic effective flux **Flux** _{ Ave } (same as figures above).

Phase coherence varies with fluctuation and Gibbs free energy. (a) Phase coherence versus system volume. The axis of 1/*V* is set as logarithmic coordinate. We choose the variation range of 1/*V* as (0.0008, 0.2) which includes the value 0.005 that we used for all above figures. Here *B* = 0.5, corresponding to Δ*G* = 11.5, can ensure that the Schnakenberg model behaves as stochastic oscillation. (b) Phase coherence versus the Gibbs free energy Δ*G*. Here 1/*V* = 0.005 is fixed, and the variation region of Δ*G* is (10.2, 12.2), same as above, which keeps the model behaving as oscillation. Both vertical coordinates of the two curves use the same scale.

Phase coherence varies with fluctuation and Gibbs free energy. (a) Phase coherence versus system volume. The axis of 1/*V* is set as logarithmic coordinate. We choose the variation range of 1/*V* as (0.0008, 0.2) which includes the value 0.005 that we used for all above figures. Here *B* = 0.5, corresponding to Δ*G* = 11.5, can ensure that the Schnakenberg model behaves as stochastic oscillation. (b) Phase coherence versus the Gibbs free energy Δ*G*. Here 1/*V* = 0.005 is fixed, and the variation region of Δ*G* is (10.2, 12.2), same as above, which keeps the model behaving as oscillation. Both vertical coordinates of the two curves use the same scale.

Landscape barrier height versus 1/*V* and Gibbs free energy. Barrier height of the potential landscape is calculated from the landscape top (red center or edge) to the bottom of the oscillation ring (blue valley) shown in Fig. 7. (a) is drawn when 1/*V* ∈ (0.0008, 0.01) with fixed energy pump *B* = 0.5, corresponding to Δ*G* = 11.5. Based on the logarithm data of 1/*V*, we obtain the linear fitness curve with *R* − *Square*(*COD*) = 0.86391. (b) is drawn when Gibbs free energy Δ*G* ∈ (10.2, 12.2) with 1/*V* = 0.05, with linear fitness *R* − *Square*(*COD*) = 0.65715. Both vertical coordinates of the two curves use the same scale.

Landscape barrier height versus 1/*V* and Gibbs free energy. Barrier height of the potential landscape is calculated from the landscape top (red center or edge) to the bottom of the oscillation ring (blue valley) shown in Fig. 7. (a) is drawn when 1/*V* ∈ (0.0008, 0.01) with fixed energy pump *B* = 0.5, corresponding to Δ*G* = 11.5. Based on the logarithm data of 1/*V*, we obtain the linear fitness curve with *R* − *Square*(*COD*) = 0.86391. (b) is drawn when Gibbs free energy Δ*G* ∈ (10.2, 12.2) with 1/*V* = 0.05, with linear fitness *R* − *Square*(*COD*) = 0.65715. Both vertical coordinates of the two curves use the same scale.

Ratio of diffusion force over probability flux force versus 1/*V* and Gibbs free energy Δ*G*. The ratio is defined by the ratio of the most probable magnitude of the diffusion force and the most probable probability magnitude of the flux force . (a) The ratio is drawn against the fluctuation magnitude in the loglog coordinate with 1/*V* ∈ (0.0008, 0.01) and *B* = 0.5. Based on the double logarithm data of 1/*V*, we obtain the linear fitness curve with high score *R* − *Square*(*COD*) = 0.97935. (b) The curve is drawn versus Gibbs free energy with the same variation as above and the fluctuation magnitude is fixed at 1/*V* = 0.005, which also shows good linear trend with the linear fitness score *R* − *Square*(*COD*) = 0.96207.

Ratio of diffusion force over probability flux force versus 1/*V* and Gibbs free energy Δ*G*. The ratio is defined by the ratio of the most probable magnitude of the diffusion force and the most probable probability magnitude of the flux force . (a) The ratio is drawn against the fluctuation magnitude in the loglog coordinate with 1/*V* ∈ (0.0008, 0.01) and *B* = 0.5. Based on the double logarithm data of 1/*V*, we obtain the linear fitness curve with high score *R* − *Square*(*COD*) = 0.97935. (b) The curve is drawn versus Gibbs free energy with the same variation as above and the fluctuation magnitude is fixed at 1/*V* = 0.005, which also shows good linear trend with the linear fitness score *R* − *Square*(*COD*) = 0.96207.

Sketch map for the definition of phase coherence.

Sketch map for the definition of phase coherence.

Article metrics loading...

Full text loading...

Commenting has been disabled for this content