Search tips
Search criteria 


Logo of jbiolphysspringer.comThis journalToc AlertsSubmit OnlineOpen Choice
J Biol Phys. 2009 October; 35(4): 425–445.
Published online 2009 June 4. doi:  10.1007/s10867-009-9156-x
PMCID: PMC2750744

Dynamical patterns of calcium signaling in a functional model of neuron–astrocyte networks


We propose a functional mathematical model for neuron-astrocyte networks. The model incorporates elements of the tripartite synapse and the spatial branching structure of coupled astrocytes. We consider glutamate-induced calcium signaling as a specific mode of excitability and transmission in astrocytic–neuronal networks. We reproduce local and global dynamical patterns observed experimentally.

Keywords: Tripartite synapse, Astrocyte, Neuron, Calcium waves, Modeling


In recent years, the importance of neuro-glial interactions for the maintenance of normal function of the central and peripheral nervous systems has been fully recognized [13]. The confinement of the role of glial cells to only the metabolic support and uptake of K + ions and neurotransmitters no longer agrees with data on glial participation in synaptic transmission, long-term potentiation, synaptic plasticity, and the development of neuronal pathologies. As a result, the new concept of the “tripartite synapse” has arisen.

The tripartite synapse includes pre- and postsynaptic neurons and the astrocytic network covering a synapse. A neurotransmitter (usually glutamate, Glu) released from the presynaptic neuron activates not only glutamate ionotropic receptors (i-GluRs) on the postsynaptic membrane but also glutamate metabotropic receptors (m-GluRs) on astrocytes. Interaction with m-GluRs results in IP3 production, subsequent IP3-dependent release of Ca2 + from endoplasmic reticulum (ER) and Ca2 + -dependent Ca2 + release from ER. The initial increase of Ca2 + concentration in the astrocyte cytoplasm depends on the level of synaptic activity. The elevation of Ca2 + above a certain threshold triggers the release of “glial transmitters”: glutamate, adenosine triphosphate (ATP), or D-serine [4]. Glial Glu interacts with i-GluRs and produces additional depolarization of the postsynaptic neuron [5].

Moreover, Ca2 + concentration in the astrocyte can oscillate due to Ca2 + uptake into the ER and the repeated Ca2 + -dependent Ca2 + release. Ca2 + waves can travel not only inside a single astrocyte, but also along a network of astrocytes through gap junctions giving rise to increased Ca2 + concentration in adjacent astrocytes [5, 6]. There are several consequences of this phenomenon. Firstly, Ca2 + waves can return to the place of their origin and can cause the release of glial Glu while the synapse is silent. Secondly, under intensive synaptic activity, the elevation of Ca2 + waves in astrocytes can be significant enough to propagate long distances in the brain and to activate postsynaptic neurons at other synaptic terminals [6]. Indeed, it has been shown that propagating Ca2 + waves can travel distances over 100 μm in both cultures of astrocytes and brain slices [7, 8]. Morphological studies of visual cortex and hippocampus showed the absence of single uncoupled astrocytes, which indicates the functional importance of networks of astrocytes [9, 10]. The number of coupled astrocytes in one network can vary from two to ten cells, and network geometries can be very diverse [11]. However, the important aspect is that one astrocyte network can be connected to another via a common astrocyte [6]. Intracellular Ca2 + waves in networks of astrocytes provide an extraneuronal pathway for long-distance signaling in the brain in parallel with neuronal synaptic transmission [12]. In addition to the above described glutamate-mediated neural–glial interaction, astrocytes can (a) suppress distant synapses via adenosine that accumulates due to the hydrolysis of glial ATP [13, 14] and (b) can increase interneuron excitability and decrease pyramidal neuron excitability due to the interaction of glial ATP with different purinergic receptors [15, 16]. It should be noted that on the neuronal network level the glial ATP-induced depolarization of interneurons leads to GABAa-receptor-mediated synaptic inhibition [15]. Moreover, astrocytes participate in heterosynaptic depression via the release of ATP with its following degradation to adenosine diphosphate (ADP). This prevents synaptic saturation and increases the contrast between activated and nontetanized synapses [17].

At present, most studies are focused on the Glu–Ca2 + –Glu signaling pathway. Its functional implications in astrocytic–neuronal networks, however, can be different. Under normal conditions, they can contribute towards information processing, synapse formation and remodeling, modulation of synaptic transmission and, probably, memory formation. However, hyperactivation of astrocytes can lead to chronic pain, pain sensation [18], the modulation of the strength of seizure-like events under epileptiform neuronal activity [19], or even epilepsy [20]. The spreading of the pain sensation can be caused by the formation of new synapses and is related to self-sustained Ca2 + oscillations in astrocytic networks. Epilepsy also involves the formation of new synapses and is related to the elevated concentration of extracellular Glu [20]. It is proposed that glutamate released from astrocytes may play a crucial role in excessive epileptic synchronous neuronal firing due to propagating Ca2 + waves. (However, there are also other mechanisms, see [19]).

The question is: Which conditions and parameters of astrocytic–neural networks are responsible for the transition from normal to pathological signaling?

The above listed neurophysiological facts and hypotheses provide the basis for the mathematical modeling of astrocytic networks with synaptically bound neurons. It can help us learn more about the functional role of Ca2 + waves in the regulation of synaptic transmission and can help us study the relevant nonlinear mechanisms. The first model of a “dressed neuron” was introduced by Nadkarni and Jung [21, 22]. They suggested a model of a neuron in contact with an astrocyte and studied how the release of glutamate during an action potential triggers the calcium dynamics of astrocytes and, in turn, changes the firing pattern of the same neuron Recently, this model was modified in order to describe the relevant details of synaptic neurotransmitter release [23, 24].

Stamatakis and Mantzaris [25] use the mathematical modeling approach to study the relationship between single cell signal transduction and the wave propagation properties of ATP-coupled astrocytic networks. Two hypotheses on Ca2 + - and IP3-dependent ATP release have been tested. Bennett et al. [26] simulated how properties of purinergic transmission at astrocyte junctions determine the characteristics of Ca2 + propagation in networks of these cells. Gibson et al. [27] constructed a model that connects neural activity and blood flow regulation in arterioles in which astrocytes play a crucial role.

Though astrocytes and their networks are the most popular objects for modeling, glial cells of invertebrates, well studied experimentally, have also received much attention. Postnov et al. [28] incorporated the activation and response pathways of giant glial cells from leeches into a quantitative model and studied common and specific features in comparison with astrocytes. It was thus shown that, in addition to glial activation caused by a slow glutamate-induced production of IP3, the fast increase of extracellular potassium concentration caused by neuronal activity provides depolarization of the glial cell followed by calcium influx and, hence, activation of a glial response.

The considerable diversity of animal and cell-dependent characteristics and the generality of signaling pathways involved in reciprocal neural–glial signaling suggest that a generalized mathematical description might be useful in studying the main operating regimes and underlying nonlinear mechanisms of neural–glial ensembles. Such an attempt was performed by Postnov et al. [29], who suggested the functional (generalized and nondimensional) model of the tripartite synapse. With such an approach, different pathways can be activated or blocked by an appropriate choice of control parameters in the model, thus reproducing specific features of different neural–glial ensembles.

In the present paper, we follow the same approach and develop a functional model of a neuron–astrocyte network that can reproduce typical local and global dynamical patterns observed in biological experiments. We use the concept of the tripartite synapse and apply it to spatially extended neuron-astrocyte networks. The spatial structure of the model mimics the geometry of astrocytic networks in a 2D approximation, while the model remains reasonably simple for theoretical and numerical analysis.

Local neural–glial interaction

Functional model

We follow a functional approach for modeling that usually combines the following steps: (1) to identify physiologically important pathways to be incorporated into the model, (2) to estimate temporal and spatial scales on the basis of experimental data, and (3) to select a minimal number of simple and well-known functional subunits from which to construct the model. Such a functional model cannot provide the quantitative description of specific biological processes because variables and parameters are dimensionless and normalized, and the time scales are fixed as control parameters. However, the model should be able to reproduce the main dynamical patterns and to predict their changes by varying the control parameters.

Although our model is qualitative and dimensionless, to keep connected with the original biological problem, we use neurophysiological terms to describe the meaning of the control parameters and model dynamics. In this way, v1 and v2 for the fast variables of the neuron model represent the transmembrane potentials, Iapp plays the role of the applied current, etc. This will allow us to discuss the model behavior in terms of the processes that we aim to simulate.

Figure Figure11 shows the main elements and signaling pathways of a tripartite synapse (a) and the corresponding feedback loops of the functional model (b). N1, N2, and A denote the presynaptic neuron, the postsynaptic neuron, and the astrocyte, respectively. S represents the synapse itself, IP3 and K represent the two possible modes of glial cell activation, while G1, G2, and ATP indicate different pathways of glial response. Let us briefly describe the mathematical implementation of each unit.

Fig. 1
The main pathways in a tripartite synapse (a) and the corresponding structure of a functional model (b)

Presynaptic neuron N1 It receives external stimuli and activates the synaptic unit S. This unit is described by the well-known FitzHugh–Nagumo model [30]:

equation M1

where v1 and w1 are the fast and the slow variables, ε1 = 0.04 is the time separation parameter, and parameter I1 defines the operating regime. At the selected value I1 = 1.02, the neuron is excitable. To count the various stimuli that may control the operating regime of the presynaptic neuron, an applied current Iapp is introduced.

Synapse S In the framework of a functional approach to modeling, the essential properties of synaptic coupling between two neurons are (a) the delayed response of the postsynaptic neuron activity and (b) the threshold for synapse activation. We describe the synaptic coupling by the first-order differential equation [31]:

equation M2
equation M3

Here, z is the synaptic activation variable and τs describes the time delay. The parameters hs, ss, and ds are responsible for the activation and relaxation of z. When v1 < hs, the synapse is inactive and z  0. Increasing v1 above hs causes a rise in the z variable to z  1 with a rate proportional to 1/τs. When v1 falls below hs, z inactivates back to zero. Once activated, z provides the postsynaptic current Isyn applied to the postsynaptic neuron N2. The factor ks plays the role of conductance, while δGm takes into account the G1 pathway of the glial response (see below). The reference level z0 is calculated on the assumption that, when the presynaptic neuron is silent, then z(t) = z0. This gives:

equation M4

Postsynaptic neuron N2 This neuron is also described by the FitzHugh–Nagumo model:

equation M5

with I2 = 1.02 and ε2 = 0.04. Isyn represents the synaptic current. To implement the glial response pathways G2 and ATP, two additional terms are included. IG2 represents the depolarizing effect of extracellular glutamate (variable Gm) released by the glial cell. IATP qualitatively describes the hyperpolarizing action of adenosine (variable Ga) being the product of ATP hydrolysis in the extracellular space (see details below).

Astrocyte A Calcium dynamics in an astrocyte is modeled using a dimensionless form of the Dupont and Goldbeter model [32, 33] with additional terms describing the contribution from activation pathways K and IP3:

equation M6
equation M7

Here, c denotes the Ca2 + concentration within the astrocyte, ce is the calcium concentration in the internal store (ER), τc, together with the time separation parameter εc, define the characteristic time for Ca2 + oscillations and transients (considerably longer than those for neurons).The term equation M8 represents the total calcium influx into the cytoplasm. Parameter r = 0.31 represents a constant transmembrane current. It controls the initial state of the calcium oscillator with no external influence (at α = 0 and β = 0). The term equation M9 is the implementation of the K activation pathway. Namely, it qualitatively describes the potential-dependent inward calcium current that is activated by the depolarization of the glial cell caused by the elevation of the extracellular potassium concentration. This pathway is approximated by the direct introduction of a recovery variable w2 for the postsynaptic neuron N2. The reference level equation M10 is defined by the condition that the postsynaptic neuron in a resting state will not influence the proceedings.The term βSm qualitatively takes into account the IP3 pathway, which includes the glutamate-induced IP3 production in the glial cell and the IP3-dependent release of Ca2 + from the ER. The addition of such a term into (6), but not into both (6) and (7), is a simplification of functional modeling. The Ca2 + exchange between the cytoplasm and the ER is described by the nonlinear function [32]:

equation M11

To take into account environmental fluctuations of different origin, a noise term with amplitude Da is added to (6). In our model, the noise term provides fluctuations of the excitation threshold for local calcium dynamics. In fact, all relevant environmental changes can contribute to this fluctuation term. For example, it can arise due to a variation in the extracellular neurotransmitter level. The relevant experimental study on cultured astrocyte networks was performed by Jung et al. [34].

Production of secondary mediator IP3 This pathway corresponds to the slow process of IP3 production inside the glial cell. The corresponding equation, in the framework of our functional modeling approach, takes the following form:

equation M12

where the Sm variable qualitatively describes the changing concentration of IP3. It is assumed that this process has a threshold value equation M13 and the steepness of activation is governed by the parameter equation M14. Parameters equation M15 and equation M16 control the overall time scale and deactivation rate, respectively.

Glial responses: G1, G2, and ATP pathways The production of glial glutamate (variable Gm) as well as the release and hydrolysis of ATP (variable Ga) are described by equations of the same form for (2) and (9):

equation M17
equation M18

where the control parameters have a similar functional role as in (2) and (9). The parameter subscripts indicate different underlying processes.As was introduced above in (3) and (5), Gm provides the modulatory effect both on synaptic current Isyn and on glutamate-dependent current IG2 in the postsynaptic neuron.Variations of Ga provide the current IATP (see (5)) which is also applied to the postsynaptic neuron. One should note that in the described local version of a functional model currents IG2 and IATP act just against each other. However, in a spatially extended model (used in Section 3.2) the situation is different due to different time scales and diffusion.A typical set of control parameters is shown in Table 1. The first six parameters (top two rows) describing calcium dynamics are taken from [32] while the remaining parameters are selected to comply qualitatively with the known relations between the processes being modeled. In this way, the reference value for all time scales is equal to 1.0 as given in the left-hand side of the second equation for both neurons. Thus, one can adjust the synaptic delay by varying τs [set membership] [0.5; 10.0]. According to our assumptions, the following relation between the other time scales has to be maintained:

equation M19

Threshold parameters hs, equation M20, equation M21, and equation M22 are selected to distinguish between inactivated and activated states of the variables v1, z, and c, respectively. Parameters ds, equation M23, equation M24, and equation M25 are fixed at typical values estimated from [31]. Below, we use α, β, γ, δ, and η as the main control parameters involved in the regulation of glial activation and response. Iapp is used as a variable external current that triggers spontaneous activity of the presynaptic neuron and, thus, allows us to consider the transient dynamics of the neural–glial ensemble.

Table 1
The typical set of model parameters

Dynamical features

As a first trial, we aim to reproduce specific features of stimulus transmission in a small neural–glial ensemble. For this purpose, we select such a stimulus in the form of a spike train of the presynaptic neuron with a duration of Ttrain = 1,000.0 and with a strength of synaptic transmission of ks = 0.2. Other parameters are fixed according to Table 1 in Section 2.1. For this parameter set, the postsynaptic neuron responds to each presynaptic spike.

To test patterns of calcium dynamics in the astrocyte (variable c), we block the glial response by selecting γ = δ = η = 0. With this, we use α and β as the main control parameters. They characterize the fast (potassium-induced depolarization followed by calcium influx) and the slow (glutamate-dependent IP3 production) activation pathways of the glial cell. Figure 2 shows the stimulus (top two panels) and respective calcium responses at a fixed rate of the fast activation pathway when α = 0.02 and at different strengths β of the slow activation pathway. The simplest response labeled 1α can be observed at vanishingly small β and, thus, is evoked only by the fast activation pathway. It is easy to see that increasing β leads not only to the appearance of additional calcium spikes (from 2 to 6 in the figure), but evokes also a temporal reorganization.

Fig. 2
The representative time courses for model variables illustrate different glial responses to neuronal stimulus at ks = 0.2, Ttrain = 1,000.0, and α = 0.02. Two top panels depict the spike train of ...

Note that the fast (controlled by α) activation pathway is typical for invertebrates, where the glial transmembrane calcium current plays a crucial role. In contrast, the β-controlled slow pathway is assumed to be the main activation pathway for astrocytes. Hence, we focus on β-controlled dynamical patterns which are essential for the functional simulation of astrocyte–neural networks. We visualize a relative contribution of fast and slow activation pathways by means of two-parameter diagrams depicted in Fig. 3. If we block the glial feedback pathways (δ = γ = 0), the increasing α (at small and medium β) adds one calcium spike to the firing patterns. Only at a high enough β the calcium response switches to a five-spike mode. Thus, the fast (potassium-related) activation pathway in our model provides a short-term excitable response from the calcium subsystem but not long-term activity. In contrast, increasing β is able to evoke a calcium response of extended duration. Interestingly, the nonmonotonic behavior of curves in the interval of α [set membership] [0.05; 0.015] reveals the competitive interaction between two mechanisms. For example, in Fig. 3a, at β = 0.0045, the increasing α drives the system initially from six-spike to four-spike behavior, then to five spikes per stimulus, and finally the response pattern contains a single calcium spike at α  0.027. Thus, the slow activation pathway provides a longer and more complex response, but the evolution of calcium patterns can be diverse, depending on the contribution from the fast activation pathway.

Fig. 3
Interplay between fast α and slow β activation pathways. Diagrams show different calcium regimes in response to glial activation. Dynamical patterns are labeled according to the number of calcium spikes. aδ = 0, ...

The other two diagrams in Fig. 3 characterize the changes in the above picture when other control parameters are varied. Panel b represents calcium patterns in response to a stimulus two times shorter (Ttrain = 500.0). One can see that the vertical structure of the diagram is almost unaffected while the β-controlled changes are more pronounced. Namely, the maximal number of calcium spikes is reduced significantly. This indicates that the overall intensity of response of the calcium subsystem is very sensitive to the duration of stimulus. The diagram in panel c shows what happens if one glial feedback pathway is activated by selecting γ = 0.5. With this, the postsynaptic neuron N2 receives the depolarizing input from the glial cell and, thus, a positive feedback involving the fast activation pathway arises. One can see that the zigzag-like segments are shifted to the left. This reflects an easier and more intensive activation of the calcium subsystem at β > 0.002.

Figure 4 represents a variety of firing patterns of postsynaptic neuron N2 when the fast activation pathway is blocked (α = 0). Panels a, b, and c show the time courses of the presynaptic neuron v1, calcium response c, and glutamate production Gm. The other panels show the response of the postsynaptic neuron v2 at different sets of control parameters, namely, panels d, e, and f illustrate how increasing γ supports the firing activity of N2 after the stimulus is terminated. Panels g and h show that the suppression of synaptic strength due to glial feedback can completely block a stimulus transmission. However, interplay between glia-induced activation of the postsynaptic neuron and glia-induced suppression of the synapse results in irregular bursing-like firing patterns (panels i–k). The last panel shows the extreme case when the postsynaptic neuron becomes almost “disconnected” from the presynaptic neuron; however, its activity is maintained due to activated glial feedback at γ = 0.5 and δ = 3.1.

Fig. 4
Top panels: Temporal dynamics of presynaptic neuron activity v1, glial calcium concentration c, and glia-produced mediator Gm. Bottom panels: Evoked activity v2 of the postsynaptic neuron at different γ and δ

To summarize, we have considered the main dynamical features of a functional model of a tripartite synapse (1)–(11). The fast activation pathway can provide rapid, but single-(calcium) spike response, which can be enhanced by increasing glial feedback γ to the postsynaptic neuron. The slow activation pathway is more or less independent of the firing patterns of the postsynaptic neuron and is able to provide a long-term potentiation of the postsynaptic neuron. According to physiological data, the slow activation pathway is expected to be a relevant mechanism for the formation of wave patterns in a spatially extended model. We did not test the ATP response pathway because (a) it has a longer time scale and cannot be activated by a short stimulus and (b) being a heterosynaptic mechanism, it affects other synapses rather than those that have been activated initially.

Global neuron–glia interactions

Spatial astrocytic network: numerical algorithm

To take into account the spatial structure of an astrocyte network, our model should be expanded to represent not only single cell calcium dynamics but also its spatial variations in a network of astrocytes coupled via gap junctions. This problem splits into three tasks: (1) to simulate the spatial geometry that is typical for astrocytic networks, (2) to provide a numerical algorithm capable of working with inhomogeneous and spatially irregular structures, and (3) to modify the model equations in order to describe multi-synapse and multi-neuron structures. Note that in a spatial astrocytic network, the ATP pathway, being the alternative pathway for long-distance heterosynaptic signaling, can be implemented.

Artificial spatial pattern

Calcium wave propagation in a continuous medium is different from propagation along highly branched astrocytic processes. To simulate the spatial structure of an astrocytic ensemble, we propose a two-dimensional algorithm that can also be applied to high-dimensional networks.

For each astrocyte, a number of primary branches (axes) na of initial length L0 are drawn originating from the cell center. The number of axes na (from eight to 13) is chosen randomly for each astrocyte. The axes are distributed uniformly at the interval of [0;2π) with (  π/6; π/6) scatter. On top of each axis, a recursively branching structure is created. For each new tentative branch, it is checked to see whether (a) the branch intersects any previously created branches of the same astrocyte or (b) it appears closer to any other cell center than that of the “ancestor”. If either of these conditions are true, the tentative branch is rejected and a new trial with the same branch length, but with a different angle relative to the “mother” branch on the previous level, has to be made. Effective branching levels rise incrementally by a probability of 10%. When a branch is accepted (i.e., it does not intersect any other branches and is closer to the original cell center than to any other cell centers), some nb (up to 4) new branches (number is randomly chosen) are created on top of it. With a probability of 80% an effective level is added. The length of each new branch is shrunk by an arbitrarily chosen factor and its orientation is adjusted relative to the “mother” branch by a randomly chosen angle within a given scatter range. Recursion stops if the effective level of branching is more than a maximal value nmax randomly chosen a priori for each astrocyte in the range of [5; 13]. Finally, a cell body represented by a circle of a given diameter is added to the cell center. To build a network of astrocytes, a desired number of astrocytes nc is created at random positions on a canvas. The examples of created patterns are shown in Fig. 5.

Fig. 5
Example patterns of simulated astrocytic networks. a A single astrocyte with extensive branching and b a network containing ten cells. The branching parameter and initial branch length are smaller than for a. The arrows indicate some of cell bodies

The described algorithm is a straightforward, rather ad hoc implementation of a generic self-repeating stochastic tree-like branching pattern with parameters, tweaked to have the resulting structures resemble the histology of real astrocytes. Resulting structures were eye compared to astrocyte images known from the literature (see [3541]). The source files for the program used to build the structures are available in [42] and readers interested in the building algorithm are encouraged to refer to it for further information.

Color-coded numerical algorithm

Numerical simulations of a model with arbitrary spatial geometry requires a special numerical algorithm. The spatial pattern created as described above contains two colors denoting whether a point belongs to an astrocytic branch or to the extracellular space. This pattern is updated manually by adding pixels of two other colors representing the presynaptic and postsynaptic neurons, respectively. It is assumed that a synapse is attributed to each postsynaptic neuron. Each color is mapped to one or more model equations. For example, the color of the postsynaptic neuron is attributed to (2) and (5). During the course of the numerical simulation, each integration step includes a loop over all pixels of the spatial color-coded pattern. For each pixel, the colors of its eight nearest neighbors are scanned. According to the color of the neighbor, the appropriate coupling term is added to the equations associated with the current pixel. For example, if the current pixel and its neighbor belong to the astrocyte process, then a diffusion term for calcium is applied, etc.

Extension of the model equations

The approach described above implies that model equations should be modified to allow different configurations. Equation (6) should be rewritten as follows:

equation M38
equation M39

where i = 1 ...N enumerates the neurons that are located close to the current point and j = 1 ...M counts the points that belong to the astrocytic network (and thus have the same color). dCa and dIP3 represent calcium and IP3 diffusion rates, respectively. To simplify the implementation of the K pathway, we introduce the Heaviside function equation M40 in the corresponding term. This function can take two values, 0 (when the neuron is silent) or α (when the neuron fires).

To introduce the diffusion of mediators Gm and Ga in the extracellular space, the corresponding equations are added to the initial set:

equation M41
equation M42

Here, N counts the neighboring points that belong to astrocytes, while M counts all points of the extracellular space.

Wave patterns

Local and propagated activation

Depending on the level of neuronal activity, the elevations in intracellular Ca2 + concentration activated by neurotransmitters can either remain restricted to an astrocytic process or can propagate in the form of Ca2 + waves to other astrocytic processes that are in contact with different neurons or astrocytes [43].

We simulate the corresponding environment by setting an active postsynaptic neuron with a synapse at the end of a small network branch. Snapshots of the close vicinity of such a neuron and the representative time courses are given in Fig. Fig.6.6. The right column shows three snapshots, each corresponding to the peak level of astrocyte activity at short (a), medium (b), and long (c) neuronal activity. One can see that there are no activated units in the first case, a few activated network elements (with higher brightness) in the case of medium stimulation, and wave propagation in the case of sufficiently long stimulation.

Fig. 6
(Color online) Right column: The snapshots of a network area in close proximity to a postsynaptic neuron N2 (ac). Left column: The representative time courses for the three network points labeled 1, 2, and 3 as indicated in the snapshot

The three left panels show the time courses for three different locations inside the astrocyte network: subscript 1 denotes the pixel adjacent to the neuron N2, subscript 2 denotes the next pixel to the left, and subscript 3 denotes the first out-of-branch pixel. Each location is characterized by the time course of calcium concentration ci (i = 1,2,3) and of IP3 concentration equation M43 that can diffuse along the branch. The black bar in each panel indicates the time interval of an applied stimulus that evokes synaptic transmission and neuronal firing.

Panel a shows that a short stimulation causes IP3 concentration equation M44 to rise but is unable to activate a calcium response of the adjacent network element. However, a trace of the stimulus is visible in the form of the small but long-term elevation of the c1 variable. Panel b illustrates the case of longer stimulation of the postsynaptic neuron within t [set membership] [200; 520]. Now the IP3 elevation is enough to trigger calcium spiking c1 in the first network element which, in turn, activates the neighboring element c2. However, the activation does not propagate outside the chosen branch. If the stimulation is longer, outer network elements also produce calcium spiking and induce a calcium wave propagated along the astrocytic network (see panel c). In this case, the time course of equation M45 shows a higher elevation, providing a higher level of calcium activity.

The transition from (a) to (b) reproduces the key features of the local model that have been discussed above: there is a threshold level for intensity and/or duration of neuronal activity that is sufficient for glial activation. The transition from the (b) to the (c) state is less obvious and seems to be related to local (adjacent processes) and distant (other network elements) activation of a glial network [43]. One would expect that a wave, once it had appeared in an excitable medium, would propagate over the whole network. However, the situation is different in neuron–astrocyte networks. The open end of an astrocyte process plays the role of a spatial threshold, i.e., a single activated element should provide an excitatory stimulus for up to five connected elements of the network. This threshold is obviously not crossed when local activation occurs, but is exceeded when a calcium wave along the whole network is initiated. The model’s behavior suggests that IP3 diffusion along an astrocytic branch might provide a decrease of the excitation threshold and might thus facilitate wave initiation.

Noise-enhanced wave propagation

The difference in characteristic time scales of neural signaling and glial response suggests that calcium wave propagation in an astrocytic network is much longer than in intraneuron signaling. However, the overall velocity of glial signaling depends on many factors. Note that experimentally observed wave patterns in cultured astrocytes look like a number of moving active spots rather than fully developed wave fronts. There is also experimental evidence that noise itself can initiate waves in cultured astrocytes [34]. Thus, let us consider how wave propagation in our artificial astrocytic network is affected by noise being added to the system.

Figure 7 illustrates the different patterns that arise at different noise intensities. With vanishing noise, the wave front is split according to the network branching structure, but propagates at a constant velocity (not shown). Let us set the noise amplitude at Da = 0.0025. The snapshots in the top row of the figure show the details of the wave propagation. At t = 1,570, the wave front is fully shaped and running toward the astrocyte branches as in the deterministic case. However, a small subthreshold elevation of local calcium concentration considerably increases the local excitability ahead of the wave front. In this area of increased excitability, the noise itself initiates calcium spikes (and waves) with a great probability before the arrival of running wave. At t = 1,880, one can observe many independent wave fragments moving in different directions. The third snapshot shows a few new wave fronts still propagating independently before the whole medium becomes “discharged” and the waves decay. Panels a–d in the figure correspond to different noise amplitudes and show the representative time courses for the four selected points indicated in the first snapshot. The point number indicates the subscript of the ci variables (i = 1...4). At weak noise (panels a and b), calcium spikes come in order: c1 spikes appear first, while c2, c3, c4 are delayed. One can see that, at stronger noise, the time delay between the leading spike c1 and the others decreases, which can be considered as wave acceleration. At Da = 0.01 (panel c) the distant wave propagation is blocked, but all branches of the astrocyte network produce calcium spikes within a short interval of time: calcium spikes c1  c4 appear in a cluster, but in a random order. The further increase of the noise intensity to Da = 0.05 (panel d) destroys the spatially synchronized behavior, i.e., calcium spikes and localized wave fronts appear in an irregular way throughout the network.

Fig. 7
(Color online) Top panel: Snapshots of wave patterns in the astrocytic network at Da = 0.0025. Bottom panels: Representative time courses for the selected points (indicated in the left snapshot) at different noise amplitudes: aDa = 0.0025, ...

The results described above suggest that cases (c) and (d) are not physiologically relevant, while cases (a) and (b) are of most interest: On the one hand, noise brings irregularity to the spatial structure of calcium waves, but, on the other hand, these spatially irregular waves still occur in response to neuronal activity, i.e. the causality of the signaling process is still preserved. Note that the propagation of such noise-induced irregular waves is much faster than in the deterministic case.

Heterosynaptic suppression

More than a decade ago, it was found that activity in a synaptic pathway can regulate the relative strength of neighboring synapses [44]. More recently, the relative contributions of different glial-produced mediators have been discussed [14]. By releasing glutamate and ATP, astrocytes provide spatially and temporally balanced excitation and inhibition to coordinate neuronal and synaptic networks.

The ATP pathway in our model describes qualitatively a slow process of ATP production and its hydrolysis to adenosine which can spread over long distances by means of diffusion.

To illustrate the effect of heterosynaptic inhibition, we implement a small neural–glial ensemble (Fig. 8). Two independent presynaptic neurons N1 and N3 have excitatory synaptic connections S1 and S3 with one postsynaptic neuron N2. It is assumed that only S3 can activate the adjacent astrocyte, while the synaptic strength of S1 can be reduced by elevated adenosine concentration Ga introduced as (ks  ξaGa) instead of (ks  δGm) in (2). The spatial spreading of Ga is implemented via diffusion over the spatial network. Patterns of Ga elevation similar to the snapshots in Fig. 7, but outside the network branches, can be detected.

Fig. 8
A schematic representation of ATP-mediated intersynaptic interaction

Figure 9 represents the time courses v1, v3, and v2 for all three neurons. To distinguish the two presynaptic signals, we introduce different bursting patterns for N1 and N3. Additionally, N3 starts to fire, after some delay, at t = 500. As is clearly seen, the postsynaptic neuron N2 initially transmits the firing pattern of the N1 neuron. When N3 becomes active, N2’s firing pattern clearly shows a mixture of both signals. However, after some time delay which is necessary for Ga propagation from S3 to S1, the synaptic transmission via S1 becomes blocked and the postsynaptic neuron N2 gets only the signal from N3 (the second half of the time course for v2).

Fig. 9
The time courses for variables v1, v3, and v2 represent the temporal firing patterns for N1, N3, and N2 neurons according to Fig. 8


To summarize, we have proposed a functional mathematical model for neuron–astrocyte networks. Despite being qualitative and simplified, it nevertheless reproduces the most typical glial responses and patterns of signal transmission.

At a local level, our model takes into account the different signaling pathways of activation and response. The fast activation pathway can provide a rapid, but single-spike, glial response that can be enhanced by increasing glial feedback to the postsynaptic neuron. The slow activation pathway is less dependent on the firing patterns of the postsynaptic neuron and is able to provide a long-term potentiation of the postsynaptic neuron. According to physiological data, the slow activation pathway is expected to be a relevant mechanism for wave formation in the spatially extended model.

We have suggested a relatively simple numerical algorithm for the creation of surrogate planar patterns that resemble the geometry of astrocyte networks. Moreover, we have developed a numerical implementation for the modeling of spatially extended neural–glial ensembles. With this, we have simulated a number of dynamical patterns of calcium signaling which are assumed to be essential features of neural–glial ensembles, namely (a) the local and propagated calcium activation in astrocyte networks, (b) a possible acceleration of excitation transmission by noise, and (c) the ATP-mediated heterosynaptic suppression effect.

We also address here the problem of model validation and generality of the results. When a model is quantitative, the best way to verify it is by direct comparison with experimentally recorded time courses for the relevant variables. This is not applicable to any functional model which is qualitative and nondimensional. Considering the results of such qualitative numeric simulations, one should keep in mind the problems of structural stability of the model’s regimes, namely (a) Can the obtained results be observed over a reasonably wide range of control parameters rather than only at specifically tuned combinations? (b) How will the temporal irregularity of a stimulus affect the response of the model? (c) How does the specific spatial structure affect observed wave patterns?

We have tested the sensitivity of the local model to a specific choice of control parameters. From the diagrams in Fig. 3, one can see that the same patterns of response can be obtained at different combinations of α and β. The same is true when varying other control parameters; thus, from a mathematical viewpoint, the model is structurally stable.

It is known that applied noise can affect the stability of the system. For example, the region at α = 0.01 and β = 0.03 in Fig. 3c, where a few bifurcational transitions between different dynamical regimes occur in a narrow range of control parameters, is very sensitive to fluctuations of different origin. Here, however, noise does not actually affect the pathways themselves, so the main features of system response are preserved. We have also recently studied [45] what happens if the presynaptic neuron receives a noisy stimulus rather than constant excitatory input. In this case, the neurons produce an irregular sequence of spikes, and the complex dynamical responses of the model are washed out, but the activation and response pathways remain in charge of the produced patterns. This is not surprising, bearing in mind that both the synapse and the slow activation pathway actually integrate presynaptic spike trains and respond to a mean firing rate rather than to individual firing events.

In Section 3.2, we have illustrated some examples of neural–glial communications. We have shown that our model is capable of reproducing the main responses and of explaining the underlying nonlinear mechanisms. In this sense, the discussed results are quite generic. Moreover, the model can be used to predict dynamical patterns and structures which have not yet been observed experimentally. For example, the mechanism discussed in Section 3.2.1 requires a spatial configuration where narrow astrocyte processes sharply open to broader branches. Specially designed experiments could target this question. Moreover, we believe that the proposed technique can be useful in the investigation of nonlinear and biological mechanisms underlying the formation of dynamical and noise-induced patterns in neuron–astrocyte networks.


This work was partly supported by the European Commission (NoE BioSim, Contract No. 4SHB-CT-2004-005137). N.B and A.B. acknowledge support from the Lundbeck Foundation and O.S. acknowledges the Skou grant from Forskningsraadet for Natur og Univers. D.P. and R.K. acknowledge support from the RFBR grant 090201049.

Contributor Information

D. E. Postnov, ur.tennur.uss.soahc@vontsop.

R. N. Koreshkov, ur.liam@9891hserok.

N. A. Brazhe, ur.usm.syhpoib@anu.

A. R. Brazhe, ur.usm.syhpoib@ehzarb.

O. V. Sosnovtseva, kd.utd.kisyf@aglo.


1. Haydon PG. Glia: listening and talking to the synapse. Nat. Rev. 2001;2:185–193. doi: 10.1038/35058528. [PubMed] [Cross Ref]
2. Fields RD, Stevens-Graham B. New insights into neuron–glia communication Science 2002. 298556–562.5622002Sci...298..556F10.1126/science.298.5593.556 [PMC free article] [PubMed] [Cross Ref]
3. Bonvento G, Giaume C, Lorenceau J. Neuron–glia interactions: from physiology to behavior. J. Physiol. 2002;96:167–168.
4. Fiacco TA. Advances in understanding new roles for astrocytes in the modulation of neuronal activity. Physiol. News. 2008;72:18–20.
5. Nedergaard M. Direct signaling from astrocytes to neurons in cultures of mammalian brain cells Science 1994. 2631768–1771.17711994Sci...263.1768N10.1126/science.8134839 [PubMed] [Cross Ref]
6. Carmignoto G. Reciprocal communication systems between astrocytes and neurons. Prog. Neurobiol. 2000;62:561–581. doi: 10.1016/S0301-0082(00)00029-0. [PubMed] [Cross Ref]
7. Finkbeiner S. Calcium waves in astrocytes—filling in the gaps. Neuron. 1992;8:1101–1108. doi: 10.1016/0896-6273(92)90131-V. [PubMed] [Cross Ref]
8. Dani JW, Chernjavsky A, Smith SJ. Neuronal activity triggers calcium waves in hippocampal astrocyte network. Neuron. 1992;8:429–440. doi: 10.1016/0896-6273(92)90271-E. [PubMed] [Cross Ref]
9. Butt AM, Ransom BR. Morphology of astrocytes and oligodendrocytes during development in the intact optic rat nerve. J. Comp. Neurol. 1993;338:141–158. doi: 10.1002/cne.903380110. [PubMed] [Cross Ref]
10. Pasti L, Pozzan T, Carmignotto G. Intracellular calcium oscillations in astrocytes: a highly plastic, bidirectional form of communication between neurons and astrocytes in situ. J. Neurosci. 1997;17:7817–7830. [PubMed]
11. Lee SH, Kim WT, Cornell-Bell AH, Sontheimer H. Astrocytes exhibit regional specificity in gap-junctional coupling. Glia. 1994;11:315–325. doi: 10.1002/glia.440110404. [PubMed] [Cross Ref]
12. Cornell-Bell AH, Sontheimer H, Cooper SM, Smith SJ. Glutamate induces calcium waves in cultured astrocytes: long-range glial signaling Science 1990. 247470–473.4731990Sci...247..470C10.1126/science.1967852 [PubMed] [Cross Ref]
13. Newman EA. Glial cell inhibition of neurons by release of ATP. J. Neurosci. 2003;23:1659–1666. [PMC free article] [PubMed]
14. Fellin T, Pascual O, Haydon PG. Astrocytes coordinate synaptic networks: balanced excitation and inhibition. Physiology. 2006;21:208–215. doi: 10.1152/physiol.00161.2005. [PubMed] [Cross Ref]
15. Bowser DN, Khakh BS. ATP excites interneurons and astrocytes to increase synaptic inhibition in neuronal networks. J. Neurosci. 2004;24:8606–8620. doi: 10.1523/JNEUROSCI.2660-04.2004. [PubMed] [Cross Ref]
16. Koizumi S, Fujishita K, Tsuda M, Shigemoto-Mogami Y, Inoue K. Dynamic inhibition of excitatory synaptic transmission by astrocyte-derived ATP in hyppocampal cultures Proc. Natl. Acad. Sci. U. S. A. 2003. 10011023–11028.110282003PNAS..10011023K10.1073/pnas.1834448100 [PubMed] [Cross Ref]
17. Serrano A, Naddjeri N, Lacaille J-C, Robitaille R. GABAergic network activation of glial cells underlies hippocampal heterosynaptic depression. J. Neurosci. 2006;26:5370–5382. doi: 10.1523/JNEUROSCI.5255-05.2006. [PubMed] [Cross Ref]
18. Hansson E. Could chronic pain and spread of pain sensation be induced and maintained by glial activation? Acta Physiol. (Oxf.) 2006;187:321–327. doi: 10.1111/j.1748-1716.2006.01568.x. [PubMed] [Cross Ref]
19. Fellin T, Gomez-Gonzalo M, Gobbo S, Carmignoto G, Haydon PG. Astrocytic glutamate is not necessary for the generation of epileptiform neuronal activity in hippocampal slices. J. Neurosci. 2006;26:9312–9322. doi: 10.1523/JNEUROSCI.2836-06.2006. [PubMed] [Cross Ref]
20. Tian G-F, Azmi H, Takano T, Xu Q, Peng W, Lin J, Oberheim NA, Lou N, Zeilke R, Kang R, Nedergaard M. An astrocyte basis of epilepsy. Nat. Med. 2005;11:973–981. [PMC free article] [PubMed]
21. Nadkarni S, Jung P. Spontaneous oscillations of dressed neurons: a new mechanism for epilepsy? Phys. Rev. Lett. 2003. 91268101.2003PhRvL..91z8101N10.1103/PhysRevLett.91.268101 [PubMed] [Cross Ref]
22. Nadkarni S, Jung P. Dressed neurons: modeling neural–glia interactions Phys. Biol. 2004. 135–41.412004PhBio...1...35N10.1088/1478-3967/1/1/004 [PubMed] [Cross Ref]
23. Nadkarni S, Jung P. Modeling synaptic transmission of the tripartite synapse Phys. Biol. 2007. 41–9.92007PhBio...4....1N10.1088/1478-3975/4/1/001 [PubMed] [Cross Ref]
24. Nadkarni S, Jung P, Levine H. Astrocytes optimize the synaptic transmission of information PLoS Comput. Biol. 2008. 451000088.240749310.1371/journal.pcbi.1000088 [PMC free article] [PubMed] [Cross Ref]
25. Stamatakis M, Mantzaris NV. Modeling of ATP-mediated signal transduction and wave propagation in astrocytic cellular networks J. Theor. Biol. 2006. 241649–668.668225491410.1016/j.jtbi.2006.01.002 [PubMed] [Cross Ref]
26. Bennett MR, Farnell L, Gibson WG. A quantitative model of purinergic junctional transmission of calcium waves in astrocyte networks. Biophys. J. 2005;89:2235–2250. doi: 10.1529/biophysj.105.062968. [PubMed] [Cross Ref]
27. Gibson WG, Farnell L, Bennett MR. A computational model relating changes in cerebral blood volume to synaptic activity in neurons. Neurocomputing. 2007;70:1674–1679. doi: 10.1016/j.neucom.2006.10.071. [Cross Ref]
28. Postnov DE, Ryazanova LS, Brazhe NA, Brazhe AR, Maximov GV, Mosekilde E. Giant glial cell: new insight throuth mechanism-based modeling. J. Biol. Phys. 2008;34:441–457. doi: 10.1007/s10867-008-9070-7. [PMC free article] [PubMed] [Cross Ref]
29. Postnov DE, Ryazanova LS, Sosnovtseva OV. Functional modeling of neural–glial interaction. Biosystems. 2007;89(1–3):84–91. doi: 10.1016/j.biosystems.2006.04.012. [PubMed] [Cross Ref]
30. FitzHugh RA. Impulses and physiological states in theoretical models of nerve membrane Biophys. J. 1961. 1445–446.4461961BpJ.....1..445F10.1016/S0006-3495(61)86902-6 [PubMed] [Cross Ref]
31. Kopell N, Ermentrout GB, Whittington MA, Traub RD. Gamma rhythms and beta rhythms have different synchronization properties Proc. Natl. Acad. Sci. U. S. A. 2000. 971867–1872.18722000PNAS...97.1867K10.1073/pnas.97.4.1867 [PubMed] [Cross Ref]
32. Keener J, Sneyd J. Mathematical Physiology. New York: Springer; 1998.
33. Dupont G, Goldbeter A. One-pool model for Ca2+ oscillations involving Ca2+ and inositol 1,4,5-triphosphate as co-agonists for Ca2+ release. Cell Calcium. 1993;14:311–322. doi: 10.1016/0143-4160(93)90052-8. [PubMed] [Cross Ref]
34. Jung P, Cornell-Bell A, Madden KS, Moss F. Noise-induced spiral waves in astrocyte syncytia show evidence of self-organized criticality. J. Neurophysiol. 1998;79:1098–1101. [PubMed]
35. Zhu X, Bergles DE, Nishiyama A. NG2 cells generate both oligodendrocytes and gray matter astrocytes. Development. 2008;135:145–157. doi: 10.1242/dev.004895. [PubMed] [Cross Ref]
36. Katyal S, El-Khamisy SF, Russell HR, Li Y, Ju L, Caldecott KW, McKinnon PJ. TDP1 facilitates chromosomal single-strand break repair in neurons and is neuroprotective in vivo. EMBO J. 2007;26:4720–4731. doi: 10.1038/sj.emboj.7601869. [PubMed] [Cross Ref]
37. Martin A, Hofmann HD, Kirsch M. Glial reactivity in ciliary neurotrophic factor-deficient mice after optic nerve lesion. J. Neurosci. 2003;23(13):5416–5424. [PubMed]
38. Bushong EA, Martone ME, Ellisman MH. Maturation of astrocyte morphology and the establishment of astrocyte domains during postnatal hippocampal development. Int. J. Dev. Neurosci. 2004;22:73–86. doi: 10.1016/j.ijdevneu.2003.12.008. [PubMed] [Cross Ref]
39. Bushong EA, Martone ME, Jones YZ, Ellisman MH. Protoplasmic astrocytes in CA1 stratum radiatum occupy separate anatomical domains. J. Neurosci. 2002;22(1):183–192. [PubMed]
40. Mi H, Barres BA. Purification and characterization of astrocyte precursor cells in the developing rat optic nerve. J. Neurosci. 1999;19(3):1049–1061. [PubMed]
41. Lee MY, Deller T, Kirsch M, Frotscher M, Hofmann HD. Differential regulation of ciliary neurotrophic factor (CNTF) and CNTF receptor alpha expression in astrocytes and neurons of the fascia dentata after entorhinal cortex lesion. J. Neurosci. 1997;17(3):1137–1146. [PubMed]
43. Fellin T, Carmignoto G. Neurone-to-astrocyte signaling in the brain represents a distinct multifunctional unit. J. Physiol. 2004;559(1):3–15. doi: 10.1113/jphysiol.2004.063214. [PubMed] [Cross Ref]
44. Manzoni OJ, Manabe T, Nikoll RA. Release of adenosine by activation of NMDA receptor in the hippocampus Science 1994. 2652098–2101.21011994Sci...265.2098M10.1126/science.7916485 [PubMed] [Cross Ref]
45. Ryazanova, L., Trenikhina, Y., Zhirin, R., Postnov, D.: Noise-induced firing patterns in generalized neuron model with subthreshold oscillations. Proc. SPIE 6436, 64360W (2007)

Articles from Journal of Biological Physics are provided here courtesy of Springer Science+Business Media B.V.