Bifurcation and Numerical Simulations of Ca2+ Oscillatory Behavior in Astrocytes

In this paper, the dynamical analysis of Ca2+ oscillations in astrocytes is theoretically investigated by the center manifold theorem and the stability theory of equilibrium point. The global structure of bifurcation and evoked Ca2+ dynamics are presented in a human astrocyte model from a mathematical perspective. Results show that the difference in appearance and disappearance of Ca2+ oscillations is partly due to two subcritical Hopf bifurcation points. In addition, the numerical simulations are performed to further verify the effectiveness of the proposed method.


INTRODUCTION
Ca 2+ as an important second messenger in the cytosol is critical for synaptic neurons and glia cells in the brain [1]. The oscillatory changes in concentration of Ca 2+ are called Ca 2+ oscillations and play an active part in the transmission of chemical and electrical signaling process [2]. Astrocytes comprise approximately 50% of the volume of human brain and exhibit not only neuron-dependent Ca 2+ oscillations but also spontaneous Ca 2+ waves [3]. It was demonstrated that the frequencies and amplitudes of Ca 2+ oscillations play key roles in Ca 2+ signal transduction in the nervous system [4]. Recent results from experiment calcium release-activated calcium channel (CRAC) have shown that it is effective for the control in inhibiting neuronal excitability by enhancing calcium release from astrocytes [5].
It was generally considered that Ca 2+ oscillations in astrocyte take place in response to external stimuli, inducing the release of neuro-active chemicals [6,7]. This view began to change as several lines of evidence indicate that these oscillations can also be formed spontaneously [8]. Nevertheless, the mechanism and functional role involved in these stochastic spontaneous Ca 2+ waves are still not well-understood. Basically, Ca 2+ signal transmission of astrocytes in the brain may vary owing to certain bifurcation principles, and different chemical information is typically characterized by frequency, amplitude, and spatial Ca 2+ propagation [9]. Dynamical mechanisms that underlie the Ca 2+ waves have been investigated from both theoretical and experimental points of view in recent years [10][11][12][13][14][15][16][17][18]. Therefore, the stability and bifurcation analysis are fundamental to investigate the appearance and disappearance of spontaneous Ca 2+ oscillations in astrocytes. In the last decades, existing mathematical models helped explore the possible dynamical mechanism of these oscillatory activities in neuronal excitability [19][20][21][22][23].

STABILITY OF EQUILIBRIUM POINT AND BIFURCATION ANALYSIS
In the present work, we apply an extension of the one-pool model proposed by Lavrentovich and Hemkin as a specific example of the stability of equilibrium point and the bifurcation scenario. This model consists of three main variables: cytosol Ca 2+ concentration (Ca cyt ), Ca 2+ concentration in the endoplasmic reticulum (Ca er ), and IP 3 concentration in cell (IP 3 ). The equations and meanings of each expression in the model are given as follows: where The details of each parameter can be found in Table 1 and [4].

ANALYSIS OF STABILITY AND BIFURCATION OF EQUILIBRIA
In the following, v in is chosen as the bifurcation parameter, corresponding to Ca 2+ inflow into the cytosol through the astrocyte's membrane. For convenience, let x = Ca cyt , y = Ca er , z = IP 3 , and r = v in , we first rewrite model (1) as the following form: The equilibrium of system (2) meets the following equations: Let x 0 , y 0 , and z 0 be the roots of Equation (2) and x 1 = xx 0 , y 1 = yy 0 , and z 1 = zz 0 , we have the following representations: The corresponding equilibrium is (0, 0, 0), and system (4) has the same properties with the equilibrium of system (2). With simple calculation, it is easy to calculate the Jacobian matrix of system (4), a 32 = 0, And one can easily obtain the following characteristic equation: After a simple calculation, we have the following equations: Owing to the meaning of x, y, z and r, special conditions meet the needs whether there exists equilibrium of system (4) when r ∈ [0.02, 0.06]. We consider the Hurwitz matrix using coefficients Q i of the characteristic polynomial: It is easy to verify that the eigenvalues of the linearized system are negative or have a negative real part if the determinants of the three Hurwitz matrices are positive: Consider the stability and bifurcations of system (4) for varying parameter v in in the case of the following Routh-Hurwitz criteria: The corresponding two values can be obtained: After the computation based on the Routh-Hurwitz criteria, when we choose r 1 = 0.02383, As r 2 = 0.05944, It can be seen that all the two values satisfy the Routh-Hurwitz criteria. After using the normal form method, one can easily obtain the following conclusions: (1) r < 0.02383, there is a stable node of system (4); (2) r = 0.02383, and system (4) has a non-hyperbolic equilibrium O 1 = (0.04766, 3.96096098, 0.0153858); (3) 0.02383 < r < 0.05944, system (4) has an equilibrium (saddle); (4) r = 0.05944, and there exists a non-hyperbolic equilibrium O 2 = (0.11886, 0.6665221778, 0.0847979); (5) r > 0.05944, there is a stable node.

Conclusion 2:
A subcritical Hopf bifurcation occurs when r passes through r 0 = 0.05944 of system (2). r < r 0 , the equilibrium O 2 is unstable, and system (2) begins to oscillate. r > r 0 , the equilibrium O 2 is stable, and the global oscillations of system (2) vanish.

NUMERICAL SIMULATIONS
In order to investigate the bifurcation phenomenon in different Ca 2+ oscillation patterns, we study the generation process with respect to the parameter v in . The bifurcation diagram of the equilibrium of system (2) in the (Ca cyt , v in )-plane [(Ca cyt , v in )-plane)] is shown in Figures 1A,B. Each point of the curve (solid line) represents a stable equilibrium, and the dashed line represents an unstable equilibrium. The equilibrium undergoes the Hopf bifurcation twice, marked by points HB1 and HB2 with respect to the bifurcation parameter v 1 in = 0.0238 µM/s and v 2 in = 0.0594 µM/s. When v in < v 1 in , there exists stable equilibrium of system (2). As v in increases, the stable equilibrium loses its stability at the point HB1 and returns to being stable at HB2.
In Figure 2, we shall present the time evolutions of cytosol Ca 2+ concentration in this model for different values of the parameter v in by numerical simulation. The left panels represent time series of Ca cyt comparison of parameter v in , and the right panels are the corresponding Ca cyt -Ca er -IP 3 phase portrait. For example, there is a single peak in this type of oscillation for v in = 0.024 µM/s in Figure 2A, and the corresponding 3D phase-space is shown in Figure 2B.
Around v in = 0.033 µM/s, it is seen that the number of peak counts and peak magnitude begin to increase, as shown in Figures 2C,D. Similarly, when v in = 0.052 µM/s, five peaks were obtained (Figures 2E,F). Moreover, it should be mentioned in Figures 2G,E, although the results for peak magnitude look very similar and in agreement with the peak counts, that the oscillatory vibration is significantly different (Figures 2G,H).

CONCLUSION
In this paper, we have theoretically investigated the stability of equilibrium and bifurcation of spontaneous Ca 2+ oscillations with a mathematical model in astrocytes. By choosing the flow of Ca 2+ from the extracellular vesicles through the membrane and into the cytosol as the bifurcation parameter, we conclude that two subcritical Hopf bifurcation points play an important role in the occurrence of Ca 2+ oscillations. By combining the theoretical analysis results in this paper, we numerically gave the Hopf bifurcations, which agree with the theoretical results. Our results may be instructive for better understanding the role of spontaneous Ca 2+ oscillations in astrocytes. Because synchronization of different oscillatory patterns may relate to bifurcation, we will give detailed research in future.

DATA AVAILABILITY STATEMENT
All datasets generated for this study are included in the article/supplementary material.

AUTHOR CONTRIBUTIONS
HZ and MY contributed to the conception and design of the study. HZ organized the literature and wrote the first draft of the manuscript. MY performed the design of figures. All authors contributed to the manuscript revision and read and approved the submitted version.