- 1Department of Theoretical Physics, Kursk State University, Kursk, Russia
- 2Department of Optics and Biophotonics, Saratov State University, Saratov, Russia
- 3Department of Biophysics, Biological Faculty, Lomonosov Moscow State University, Moscow, Russia
- 4Department of Molecular Neurobiology, Institute of Bioorganic Chemistry RAS, Russian Federation, Moscow, Russia
Neuronal firing and neuron-to-neuron synaptic wiring are currently widely described as orchestrated by astrocytes—elaborately ramified glial cells tiling the cortical and hippocampal space into non-overlapping domains, each covering hundreds of individual dendrites and hundreds thousands synapses. A key component to astrocytic signaling is the dynamics of cytosolic Ca2+ which displays multiscale spatiotemporal patterns from short confined elemental Ca2+ events (puffs) to Ca2+ waves expanding through many cells. Here, we synthesize the current understanding of astrocyte morphology, coupling local synaptic activity to astrocytic Ca2+ in perisynaptic astrocytic processes and morphology-defined mechanisms of Ca2+ regulation in a distributed model. To this end, we build simplified realistic data-driven spatial network templates and compile model equations as defined by local cell morphology. The input to the model is spatially uncorrelated stochastic synaptic activity. The proposed modeling approach is validated by statistics of simulated Ca2+ transients at a single cell level. In multicellular templates we observe regular sequences of cell entrainment in Ca2+ waves, as a result of interplay between stochastic input and morphology variability between individual astrocytes. Our approach adds spatial dimension to the existing astrocyte models by employment of realistic morphology while retaining enough flexibility and scalability to be embedded in multiscale heterocellular models of neural tissue. We conclude that the proposed approach provides a useful description of neuron-driven Ca2+-activity in the astrocyte syncytium.
1. Introduction
Astrocytes of the cortical and hippocampal gray matter are important actors in a number of information processing processes, including synaptic plasticity, long-term potentiation, and synchronization of neuronal firing (Haydon, 2001; Lee et al., 2014; De Pitta et al., 2016; Poskanzer and Yuste, 2016) as well as in coupling neuronal activity to blood flow changes (Otsu et al., 2015). Recent evidence converges on a close connection of these functions with whole-brain processes and systemic regulation pathways. Thus, astrocytes respond to and are able to regulate systemic blood pressure (Marina et al., 2020); they significantly (up to 60%) change their volume during sleep or under anesthesia (Xie et al., 2013); astrocytes play an important role in the clearance of beta-amyloids, a process with mechanisms that are now being actively discussed (Iliff et al., 2012; Abbott et al., 2018; Semyachkina-Glushkovskaya et al., 2018; Mestre et al., 2020); both intracellular and network-level activity of astrocytes are significantly different in sleep and during wakefulness, and activates with locomotion (Bojarskaite et al., 2020; Ingiosi et al., 2020; McCauley et al., 2020). It is important to note that many of the mentioned astrocyte functions are not directly related to neural activity, but are governed by their own regulatory pathways (O'Donnell et al., 2015). Some of these functions are tightly linked to dynamic regulation of astrocyte morphology and volume and depend, for example, on the circadian rhythm of aquaporin expression (Hablitz et al., 2020).
In summary, this frames a new mindset for understanding the function of astrocytes and at the same time poses a challenge for modeling studies. Namely, the morphological features should now be considered as a specific control parameter that significantly contribute to the both single-cell dynamics and network activity patterns. This problem breaks down into three specific tasks: (i) to provide tractable, but still biologically reasonable mathematical account for contribution of subcellular morphological features to intracellular calcium dynamics; (ii) to further develop approaches to modeling of Ca2+ dynamics on data-driven irregular structures, both for an individual cell and for a network; (iii) to reveal how realistic morphological features are manifested in the spatiotemporal patterns of the calcium dynamics.
We address these tasks in more detail in the rest of the Introduction.
1.1. Calcium Signaling in Astrocytes
With plasma membranes enriched in a variety of potassium channels and lacking voltage-gated sodium channels, astrocytes are not electrically excitable (Verkhratsky and Nedergaard, 2018). On the other hand, they display a rich repertoire of Ca2+-activity at multiple spatial and temporal scales (Lind et al., 2013; Volterra et al., 2014; Wu et al., 2014; Bindocci et al., 2017). Although astrocytic Ca2+ transients can occur spontaneously, their frequency is modulated by neuronal activity (Stobart et al., 2018), changes in local tissue oxygenation (Mathiesen et al., 2013; Marina et al., 2020), and other factors (Semyanov et al., 2020). As outputs, Ca2+-activity in astrocytes leads to release of signaling molecules: gliotransmitters, such as GABA, D-serine, and glutamate, as well as vasoactive metabolites (Serrano et al., 2006; Henneberger et al., 2010; Bazargani and Attwell, 2016). This has been summarized in a concept of “tripartite synapse,” i.e., sensing of synaptic neurotransmitter release by perisynaptic astrocyte processes, encoding this information in Ca2+ signals and response with secretion of neuroactive molecules (Araque et al., 2014). There is still however an ongoing debate on the mechanisms involved in generating Ca2+ transients in astrocytes and the extent of effect of astrocyte-derived molecules on synaptic plasticity, e.g., on LTP (Fiacco and McCarthy, 2018; Savtchouk and Volterra, 2018).
Recent experimental evidence obtained with genetically encoded or pipette-loaded Ca2+ indicators (Tong et al., 2013; Rungta et al., 2016) heralds functional segregation between the less frequent global internal store-operated Ca2+ transients at the level of cell soma and primary branches, and the more frequent spatially limited microdomain Ca2+ transients in the thin mesh of astrocytic leaflets—ramified nanoscopic processes, also known as perisynaptic processes (PAPs) due to their proximity to synaptic connections between neurons. The transients located in the leaflets primarily rely on influx of Ca2+ through plasma membrane, in part because of the high surface-to-volume ratio in this region and in part because the leaflets are often devoid of organelles including ER (Patrushev et al., 2013) and thus can not support exchange with intracellular stores.
The coupling from synaptic activity to local Ca2+ transients in PAPs and from the latter to global Ca2+ events is an area of active research. As reviewed in Savtchouk and Volterra (2018), early works attributed this to activation of G-protein coupled receptors to glutamate, but later this pathway has been put to question due to apparent lack of mGluR5 receptor expression in adult astrocytes. Alternative sources of microdomain transients have been proposed, such as via TRP channels (Shigetomi et al., 2011), from mitochondria (Agarwal et al., 2017), etc. One plausible alternative causal pathway can be formulated as follows (Rojas et al., 2007; Verkhratsky et al., 2012; Kirischuk et al., 2016; Parpura et al., 2016): neurotransmitters, released from the presynaptic membranes, primarily glutamate, but also GABA, are cleared from the extracellular space by astrocytic transporters utilizing Na+ gradient to drive the neurotransmitters into the cell. This leads to build up of Na+ ions in the cytosol, which can lead to temporary reversal of Na+/Ca2+-exchanger allowing for Ca2+ entry via this transporter. Conceivably, if this local Ca2+ influx happens near the ER and coincides with an increase in inositol trisphosphate (IP3) production by phospholipase C, it can trigger Ca2+-induced release of Ca2+ from intracellular stores via IP3 receptors (IP3Rs) of the ER.
The release of calcium from ER is spatially inhomogeneous due to the non-uniform, clustered, distribution of IP3 receptors (Smith et al., 2009; Taufiq-Ur-Rahman et al., 2009; Ross, 2012), with clusters spaced at about 0.5–5 μm apart. At a detailed level, calcium release from the receptor clusters has a stochastic character. The effect of the stochastic activation of IP3R clusters on the calcium dynamics has been investigated by Shuai and Jung both in point and distributed models (Shuai and Jung, 2002, 2003). In the case of a large enough number of clusters, Ca2+ release events can be averaged to a lumped deterministic description. Particularly, the increase in IP3 level transforms stochastic calcium increases into regular waves.
Recapitulating, calcium signaling mechanisms are inhomogeneous across the cell and depend on local morphological parameters, which has to be taken into account in modeling. It seems practical to introduce a metaparameter to describe the relative inputs of store-related and plasma membrane-related Ca2+ pathways. This metaparameter can reflect local surface-to-volume ratio or the dominant size of processes and can empirically be linked to the astrocyte cytoplasm volume fraction parameter, which can be estimated directly from fluorescent images.
1.2. Cell Morphology and Network Connectivity
Astrocytes have intricate and highly complex morphology, which raises computational issues and demands an elaborate approach to modeling. The contribution of the astrocytic spatial segregation and coupling to brain physiology and functions is still not sufficiently understood, especially taking into account that astrocyte-to-neuron and astrocyte-to-astrocyte interaction mechanisms are diverse and depend on brain region. The existence of intercellular Ca2+ waves traveling across the network of astrocytes suggests a distinct mechanism for long-distance signaling (Cornell-Bell et al., 1990) and plasticity, which operates in parallel to and at much slower time scales than neuronal synaptic transmission (Pirttimaki and Parri, 2013; Sims et al., 2015).
The size of cliques of cortical astrocytes coupled within a local network is estimated around 60–80 cells (Haas et al., 2006; Houades et al., 2006, 2008), but several networks can also connect via a limited number of “hub” astrocytes (Carmignoto, 2000). The implications of inter-astrocyte connectivity have been analyzed in a modeling study by Lallouette et al. (2014) with the main conclusion that sparse short-range connections can promote Ca2+ wave propagation along the network. This allows to conjecture that once initiated, a wave of excitation can propagate over long distances in the brain cortex and affect (activate or inhibit) postsynaptic neurons at distant synaptic terminals, although most Ca2+ events are confined to a single astrocyte spatial domain. Propagating calcium waves can travel distances of more than 100 μm with speed from 7 to 27 μm/s in culture and brain slices (Dani et al., 1992). However, the waves observed in vivo rarely spread more than 80 μm (Hoogland et al., 2009; Brazhe et al., 2013), although this observation can be influenced by imaging protocol, as Kuga and colleagues reported large-scale Ca2+ glissandi in vivo that were only observable under low laser intensity (Kuga et al., 2011).
It follows that for meso-scale problems related to brain tissue physiology, it is computationally cumbersome to build a ground-up model starting from individual processes. We propose a more pragmatic approach based on texture-like volume segmentation to classes such as “soma,” “large branches,” and “gliapil” or a mesh of unresolved thin processes. This rasterization radically simplifies model implementation and scales to large networks. At the same time, by defining morphology-based spatial distribution of a metaparameter, one can study the effects of spatial heterogeneity at different scales. Indeed, the spatial distributions used for simulations are ideally data-driven. Because it is not always possible to infer the astrocytic network structure or even individual domain boundaries from experimental data, and because the networks can be variable anyway, it seems inviting to generate variable astrocytic tilings from images of individual cells.
1.3. Modeling Studies
Models of IP3-mediated Ca2+ oscillations have been extensively reviewed both in general (Dupont et al., 2011) and in application to astrocytes (Riera et al., 2011; Manninen and Havela, 2017; Oschmann et al., 2017a), which included both point- and spatially extended models. In particular, the De Young–Keizer model stemmed several currently popular models of Ca2+ dynamics in literature. This model allows to simulate IP3-sensitive calcium dynamics in cytoplasm and ER occurring at the constant level of IP3 including also a variant of the model with the positive-feedback mechanism of Ca2+ on IP3 production (De Young and Keizer, 1992). Li and Rinzel (1994) reduced De Young–Keizer model to a two-variable system introducing the experimentally observed time scale difference between fast and slow inactivation of IP3 receptor by Ca2+. Adding the dynamics for [IP3] with synthesis dependent on activation of metabotropic glutamate receptors and [Ca2+] degradation leads to a three variable model (Ullah et al., 2006). Also building on Li–Rinzel model and providing a more detailed description of IP3 degradation, De Pitta and co-authors proposed a three-variable model for glutamate-induced intracellular calcium dynamics caused by the synaptic activity in astrocytes (De Pitta et al., 2009).
One of the first models for intercellular propagation of calcium waves has been described in (Sneyd et al., 1994) by adding diffusion of IP3 and cytosolic Ca2+ to the two-pool Ca2+ model. The effect of Ca2+ diffusion rate on spatiotemporal patterns of Ca2+ signaling was studied by Shuai and Jung (2003) in a lattice-based model. Later, Kang and Othmer (2009) regarded networked astroglial Ca2+ signaling in a 2D model using spatial patterns in form of sparsely connected irregularly branching cells with simplified morphology. Both intracellular diffusion of IP3 via gap-junctions and extracellular purinergic signaling was regarded as a mechanism of intercellular communication in Kang and Othmer (2009); Edwards and Gibson (2010); intercellular Ca2+ diffusion was however disregarded in most modeling studies, e.g., Ullah et al., 2006; Kang and Othmer, 2009; Edwards and Gibson, 2010, primarily based on the notion of a much faster diffusion of IP3 than Ca2+, see Allbritton et al., 1992, and small permeability of gap junctions to Ca2+. More recently Savtchenko et al. (2018) suggested an advanced NEURON-based modeling environment for detailed spatially extended models of astrocytes. However, they did not address full calcium dynamics models or morphology-defined variations of mechanisms. Specifically, the relative weights of plasma membrane-dependent mechanisms (IP3 synthesis and Ca2+ influx) and store-dependent mechanisms scale with astrocytic process morphology, as defined by surface to volume ratio, cytoplasm volume fraction and the physical presence of ER in the process. This has been studied in point-models by Oschmann et al. (2017b) and in 1D extended model by Wu et al. (2018). Recently, Brazhe et al. (2018) studied the implications of the spatial segregation between IP3 synthesis and plasma membrane exchange and the IP3-mediated ER exchange in discrete spatial templates of variable complexity.
The tripartite synapse concept and the computational role of astrocytes in neural network activity has early attracted the attention of modeling studies, pioneered by papers by Nadkarni and Jung (2004, 2007). Understanding of the tripartite synapse from the viewpoint of non-linear dynamics and functional models have been developed in works of Postnov et al. (2007, 2008, 2011), and Tewari and Majumdar (2012). Later, tripartite synapses have been adopted in more formal neural network models (Alvarellos-González et al., 2012; Sajedinia and Hélie, 2018; Lenk et al., 2020). Because there is still no consensus based on experimental evidence on mechanisms of Ca2+ transients in PAPs and gliotransmission effects (Savtchouk and Volterra, 2018), it is hard to formulate a comprehensive model that would include all conceivable pathways and still remain tractable. While at this stage refraining from closing the loop from astrocytes to neurons, we believe it is important to understand the spatiotemporal patterning of astrocytic Ca2+ signaling at levels from microdomains to networks.
1.4. The Proposed Modeling Approach
Our main motivation in this study is to learn if uncorrelated background synaptic activity, when sensed by astrocytes, will be shaped into morphology-defined patterns of Ca2+ signaling. We present a model of multi-cellular network of astrocytes based on realistic spatial templates. We start from a single-cell model, which is considerably simpler than in Savtchenko et al. (2018), allowing for smaller computational costs, and move on to connect separate cells together to obtain a network model.
We focus on the implications of the morphology-dependent spatial segregation of the Ca2+ signaling mechanisms between astrocytic leaflets and branches. We follow the lines set out in Brazhe et al. (2018) toward more realistic and larger scale spatial templates, ranging from single astrocytes to networks. In contrast to astrocytes in culture or retina in vivo, cortex astrocytes are non-flat and occupy some volume in 3D space. Nevertheless, we chose to reduce dimensions to 2D and flatten astrocyte images used as spatial templates. One reason for this was to reduce computational cost, especially when addressing network models. Another reason that most existing Ca2+ imaging data are obtained as time series in single focal plane, and the experimentally obtained dynamics is confined to flat 2D anyway. The work of Bindocci et al. (2017) demonstrated the richness of Ca2+ dynamics within the whole astrocytic domain in 3D, but volumetric imaging is not yet widely used in the context of astrocytes. We therefore contemplated that using 2D templates for simulations would not restrict us from observing diverse and physiologically relevant Ca2+ signaling patterns, event if real astrocytes have more degrees of freedom. The rest of the paper is organized as follows: we start from a description of the proposed model in a top-down order: the general concept is followed by proposed algorithm of creating spatial templates for modeling and then continues with description of the differential equations for dynamics of intracellular and ER Ca2+, intracellular IP3, and extracellular glutamate concentrations. Having defined the model, we test its plausibility on single-astrocyte templates and after quantification of Ca2+ event statistics we proceed to behavior of astrocyte networks, where we observe noise-driven regular activation patterns.
2. Model
2.1. Model Design and Overview
In this work we aim to conceptualize our current understanding of spatial organization of the astrocytic Ca2+ dynamics in a form of a spatially detailed model of individual and networked astrocytes excited by stochastic background neuronal activity. In the light of the striking differences between Ca2+ signaling in astrocytic leaflets and thin processes on the one hand and global somatic signaling on the other, we start with segregation of the modeling space into three major classes as shown in Figure 1: astrocyte soma with thick branches (I), a mesh of astrocytic thin processes (II) and extracellular space (III). The continuum between the two extreme classes I and II is defined as local fraction of astrocytic cytoplasm volume (AVF) and a related parameter—local surface-to-volume ratio (SVR) of the astrocytic processes. In extreme class I regions, such as soma, Ca2+ dynamics are dominated by exchange with intracellular stores, and a unit of modeled space (template pixel/voxel) contains only astrocyte, while in extreme class II regions (leaflets), Ca2+ dynamics is dominated by exchange with plasma membrane and each modeled pixel contains a mesh of extremely thin astrocyte processes tangled with neuropil. We thus define a mapping of each pixel in the spatial model template to either class III (no astrocyte) or to a continuous variable between the extreme cases of class I and II with implications in local calcium dynamics and diffusion.
Figure 1. Model structure and molecular mechanisms. (A) Astrocytic network is segmented in three spatial compartments: I—cell bodies and thick branches; II—the mesh of thin branches; III—extracellular space. (B) Model variables (in colored ovals) and main regulatory pathways of intra-astrocyte calcium dynamics.
With regard to the local calcium dynamics, the extreme complexity and sheer number of cellular pathways involved, makes the detailed and comprehensive modeling of every Ca2+-related mechanism extremely challenging. Not surprisingly, there is a substantial body of published models that aim to account for the essential features of calcium dynamics in astrocytes, which do not completely agree with each other (Manninen and Havela, 2017). To choose the best model we build upon a model proposed by Ullah et al. (2006) as a prototype, while other models could fit in the proposed approach as well, for example the “ChI” model by De Pitta et al. (2009), which is similar to that of Ullah et al.
2.1.1. Spatial Structure
To represent astrocyte networks with realistic geometry of the regions I–III, one needs to create such templates algorithmically or, alternatively, obtain them from experimental data. Each of the two variants has its benefits and drawbacks. To provide just two examples, the experiment-based approach was employed in Wallach et al. (2014) and an algorithmic creation of network templates was employed in Postnov et al. (2009). Here, we draw advantages from both approaches by suggesting a simple stochastic data-driven algorithm to create realistic surrogate spatial templates of astrocyte networks. Specifically, we use experimental images of astrocytes obtained from a public database, and arrange network structure using Voronoi partition and simple geometrical transformations (see section 2.2.1).
2.1.2. Neuronal Activity
We assume that astrocytic Ca2+ response to local neuronal activity is primarily driven by the transporter-mediated uptake of neurotransmitters released from presynaptic membranes. One of the possible coupling mechanisms is the reversal of the Na+/Ca2+-exchanger transport due to an increase in [Na+] allowing for a Ca2+ influx. Here, we sacrifice biophysical details in favor of model simplicity and assume that astrocyte calcium dynamics is excited directly by glutamate released from the presynaptic terminals, causing transient fluxes of Ca2+ through the plasma membrane. A typical cortex astrocyte is associated with 300–400 individual dendrites and is in contact with about 104–105 synapses (Bushong et al., 2002; Halassa et al., 2007). Judging by these numbers and taking into account sparsity of neuronal signaling in the cortex it seems safe to treat each pixel in the distributed model template as associated with a single or just a few individual synapses. For as long as we are not focused on information processing in the cortex, we can assume independent stochastic nature of spiking activity in any of the presynaptic units and describe local activity only statistically, neglecting any complex spike timing patterns. Consequently, we describe the synaptic glutamate drive to the model in each pixel as triggered by presynaptic spike trains drawn from independent homogeneous Poisson process ξp(t) with intensity p Hz.
2.2. Astrocyte Network Topology
2.2.1. Data-Driven Network Generation
Astrocytes, like neurons, have complex morphology. Ideally, an algorithm to create spatial templates should provide means to “grow” realistic branching 3D shapes of astrocytes from a set of randomly placed “seed” locations. Indeed, there are many experimental and modeling studies of the branching patterns for various types of neurons (Ascoli et al., 2007; Donohue and Ascoli, 2008; Cuntz et al., 2010; Polavaram et al., 2014), providing means for creation of realistic surrogate shapes of as many neurons as needed. However, unlike neurons, there is less data available on the statistics of astrocyte branching, which makes it harder to create surrogate spatial templates of astrocyte networks. This hindrance can be circumvented by using a public database of microscopic images of cortical and hippocampal astrocytes (Martone et al., 2002, 2008).
To create a library of realistic spatial templates for individual cells, we downloaded a set of 27 fluorescent confocal 3D stacks of hippocampal astrocytes (4-week old rats, microinjection loaded with lucifer yellow in acute slices) (Bushong et al., 2004). The stacks have average lateral resolution of ≈ 0.07 μm/px and vertical (Z-axis) resolution of 0.2 μm; there are 45-60 Z-planes in each stack, thus encompassing the thickness of about 10 μm along the Z-axis. Because our model is set in two-dimensional space, the stacks were flattened along the Z-axis by max-projection, Figure 2. The projections were downsampled 4 × before simulations, resulting in lateral resolution of ≈ 0.28 μm/px. Each of the experimental astrocyte images then serves as a progenitor of randomized offsprings obtained by applying 250 random rotations (from 0° to 360°, shearings and stretchings (within ±20% of original size, uniform distribution), which results in a collection of 6,750 randomized pseudo-experimental astrocyte templates, used to tile the model space. Such data set expansion from a limited number of “real-world” objects is a popular approach in machine learning (Simard et al., 2003; Krizhevsky et al., 2012) helping to prevent overfitting and providing for transformation-invariant feature learning.
Figure 2. Algorithm to create surrogate templates of astrocyte network. First, a set of seeding points on a regular grid (light-gray) is perturbed with random shifts (dark-gray points). Voronoi diagram is then drawn for these points (blue lines). Each patch in the Voronoi partitioning is then filled with the best shape-matching template from an augmented collection of astrocyte images. The lookup collection is created from a set of experimental images taken from CCDB (Martone et al., 2002, 2008) by applying multiple different random rotations and shears to each experimental image.
Inspired by the fact that astrocytes establish distinct non-overlapping territories, we employ an algorithm based on Voronoi partitioning and active contours to tile the model space with astrocytes. First, we create a lattice of “seed points” regularly spaced at some intervals corresponding to average cell density, typically around 50 μm, shown in light gray in Figure 2. The resulting regular grid is then deformed by jittering x and y coordinates of every point by a random displacement value drawn from Gaussian distribution with σ = 10 μm (dark gray points in Figure 2). Different values of spacing and jitter can be used, the ones used here tended to give the most realistic tiling results. Next, a Voronoi diagram, which for each seed point delineates territories closer to it than to any other seed point, is drawn for the jittered points. We then iteratively pick a polygonal area patch from the Voronoi partitioning, look up a template astrocyte from the randomized collection, with a convex hull best matching the shape of the given Voronoi patch, and place this template into the model space. Repeated for all patches in the Voronoi partition, this creates a preliminary tiling with partially overlapping domains of neighboring astrocytes and occasional empty spaces. Next, this draft tiling is optimized with an active deformable model: the perimeter of each cell template is treated as an elastic two-dimensional curve, which optimizes an energy functional designed to promote repulsion between overlapping regions and adhesion between neighboring cells, with a penalization of the major cell shape deformation. After all domain boundaries are settled, the spatial templates are interpolated into the deformed contours. The described process of the network template creation is visualized in Supplementary Video 1.
2.2.2. Computational Design
Our simulations are based on compiling an encoded raster image representation (a template) of the model space to region-specific equations. For the sake of computational simplicity, we use two-dimensional spatial layout—each pixel of the spatial template can be interpreted as a thin slab, occupied either exclusively (e.g., in the soma) or partly by astrocyte cytosol; or as belonging to extracellular space. As follows from this approach, each pixel in the model space has to be assigned to either astrocyte-free space (class III) or astrocyte-occupied space, ranging from class I, astrocyte soma and thick branches, to class II, i.e., elements of volume containing a tangle of thin astrocyte processes and unresolved neuronal structures, e.g., synaptic boutons. We account for a graded transition from thick branches to thin processes to leaflets by introducing a local astrocyte volume fraction (AVF) parameter, which defines the landscape of how much of each pixel volume is occupied by astrocyte in the 3D prototype. AVF here is defined as a ratio between local fluorescence intensity of the template and the intensity at the soma AVF = max(I/Imax, AVFmin). An example of the described mapping from image intensity to AVF is shown in Figures 3A,B, where the colormap in Figure 3B is such that the non-zero blue channel delineates the presence of astrocyte cytoplasm (non-zero AVF), and the intensity of the red channel encodes the AVF value. To describe the relative input of the store-operated calcium flux and plasma membrane flux, we introduce a surface-volume ratio (SVR) parameter, which inversely depends on AVF (Figure 3C). The SVR value is maximal at the edges of the leaflets and minimal in the soma. Accordingly, a simple raster RGB image serves as a spatial template to encode the model space. Specifically, non-zero values in the blue channel define astrocyte-occupied pixels, while intensity in the red channel encodes AVF and ranges from minimal value AVFmin (class II) to 1 (class I). Thus, one can set up computation for a specific spatial template by simply drawing it algorithmically or with an indexed palette using a graphical editor.
Figure 3. Example of AVF color-coding. (A) Maximal projection of a confocal image of a cortical astrocyte. (B) color-coded template ready for simulation; regions with non-zero blue channel delineate astrocyte domain, while intensity of the red channel encodes AVF. (C) Link between color-coded AVF and SVR parameters.
At each integration step the master program module optionally compiles the provided image into a set of equations by mapping each pixel color to equation set following the color-coded dictionary. For each pixel first the point dynamics are applied, i.e., right-hand terms are evaluated. Next, diffusion of ions and molecules and any other short-range interactions is taken into account based on the class of the neighboring pixels. This approach is flexible, but has an overhead of compiling the color-to-equation mapping. To improve the computational performance we employ NVIDIA CUDA, a parallel GPU-based computing technology.
2.3. Intracellular Calcium Dynamics: Principal Quantities and Flows
The model for local Ca2+ dynamics is based on that of Ullah et al. (2006) with a few modifications previously introduced in Brazhe et al. (2018), which sum up to treating ER calcium as a dynamic variable, adding neurotransmitter-dependent calcium influx via plasma membrane, and segregation between thick and thin processes. Below we describe the proposed model, focusing on the differences with the Ullah et al. (2006) model; equations and parameters that are the same here as in the Ullah model are omitted.
The principal variables of the model are (i) the cytosolic calcium concentration , (ii) calcium concentration in the endoplasmic reticulum , (iii) inositol trisphosphate concentration in the cytosol [IP3] and (iv) extracellular glutamate concentration [Glu].
To account for the morphology-based spatial heterogeneity of astrocytes, we introduce a parameter r ∈ (0 < rmin…1]—a scalar quantity, roughly representing local AVF. This parameter also defines a linked parameter s representing local SVR of the astrocyte processes; SVR is inversely related to AVF: s = 1/(1−exp[0.1(r − 0.5)]). SVR scales relative inputs of Ca2+ exchange through plasma membrane and with ER as described below, while AVF scales effective diffusion coefficients for Ca2+ and IP3.
The equation set for the principal variables reads:
where JER is the total flow of calcium ions in exchange between the cytosol and endoplasmic reticulum; Jpm is the total flow of calcium ions through the plasma membrane in exchange between the cytosol and extracellular space; IGlu and ICa stand for glutamate- and calcium-dependent inositol trisphosphate production mechanisms, mediated by phospholipases β and δ; Ieq is a simplified first-order equilibration of inositol trisphosphate concentration to the basal level [IP3]0; [Glu]amb is the ambient concentration of extracellular glutamate, and τGlu is the timescale of its clearance and return to the baseline level; ξp(t) is stochastic source of glutamate from nearby located neuronal synapses triggered by Poisson spike trains in each pixel. Additionally, Jdiff, Idiff, and Gdiff describe the finite-element approximation of diffusion of cytosolic Ca2+, IP3 and extracellular glutamate, respectively and are described below.
The weighting coefficient s accounts for the stratification of intracellular dynamics according to Figure 3: in the leaflets r = rmin ≪ 1 and s ≈ 1, while for deep cytosol locations r = 1 and s ≈ 0. With this we (i) describe that input from all plasma membrane calcium currents is maximal in the leaflets and (ii) assume that endoplasmic reticulum does not invade leaflets much, thus ER exchange is large only in thicker branches and soma. We also assume that all IP3 is produced in the plasma membrane by means of G-protein coupled phospholipase β or Ca2+-dependent phospholipase δ.
2.3.1. Calcium Exchange Between the Cytosol and ER
Total calcium flow across the endoplasmic membrane is composed of IP3R-mediated current JIP3, leak of Ca2+ from endoplasmic reticulum Jleak, and contribution of ER membrane Ca2+ pump Jpump:
Ca2+ current via IP3 receptors JIP3 is modeled in the same way as in Ullah et al. (2006) (Equations 2, 4, 5, 9–12). Two other terms in Equation (5) stand for the leak of calcium from ER, and for the pumping it back, respectively, following Equations (3,6) in Ullah et al. (2006).
2.3.2. Transmembrane Calcium Flows
Transmembrane calcium exchange Jpm consists of three flows:
where Jin describes the sum of background constant Ca2+ influx and agonist-dependent IP3-stimulated Ca2+ influx across plasma membrane from the extracellular space and Jout is an extrusion current (eqs. 7–8 in Ullah et al.); JGlu = γ[Glu] describes the direct effect of extracellular glutamate on additional calcium influx.
2.3.3. Inositol Trisphosphate Turnover
The dynamics for IP3 concentration from Equation (3) has the following terms: first, we use a lumped first-order description of [IP3] equilibration to a resting level [IP3]0:
second, we account for the Ca2+-stimulated IP3 production in the same way as in Ullah et al. (2006), (eq. 14), and third, we account for glutamate-driven IP3 production IGlu following Ullah et al. (2006), (eq. 15).
2.3.4. Synaptic Glutamate Drive
Stochastic glutamate source ξp(t) in each pixel is modeled as quantal release triggered by a spike train drawn from a homogeneous Poisson process with intensity p, which agrees with statistics of neuronal firing (Softky and Koch, 1993). Accordingly, the ξp(t) term in Equation (4) is given by
where A is the instantaneous increase in glutamate release rate associated with each presynaptic event and tk are times of presynaptic spikes in the given pixel following Poisson process with intensity psyn.
2.4. Intracellular Diffusion
Elevated cytoplasmic Ca2+ can remain confined to the spatial domain of a single astrocyte, but can also spread to the neighboring astrocytes (Nedergaard, 1994; Carmignoto, 2000; Falcke, 2004) in a wavelike manner. The involvement of a large number of cells into a wave is still not fully understood though it may be an important aspect of information processing in the brain (Haas et al., 2006). At least two main mechanisms can account for the intercellular wave propagation: (i) secretion and diffusion of extracellular ATP and its action on P2Y receptors on astrocytic membranes and (ii) diffusion of intracellular IP3 and Ca2+ between contacting astrocytic leaflets, via gap junctions. Relative input of the two mechanisms differs across brain regions and for the cortical astrocytes the one mediated by the gap junctions has been reported to prevail (Haas et al., 2006). The current work is therefore focused on the latter mechanism. Accordingly, astrocytes in the model are networked by an analog of gap junctions dispersed over the parts of the cell perimeter and simulated as the connection of areas with low AVF.
Here we employ a rather simplified description of diffusion in the cytoplasm where the region occupied by astrocyte is considered as a continuous space. As corroborated by evidence for autologous gap junctions between the processes of the same astrocyte (Wolff et al., 1998; Nagy and Rash, 2003; Genoud et al., 2015), astrocyte cytosolic volume can be described as a porous sponge-like medium rather than a branched structure or acyclic graph. Thus, possible hindrance to IP3 or Ca2+ diffusion through the intricate mesh of astrocytic processes due to tortuosity and porosity of the astrocytic volume can be accounted for by a simple scaling of the apparent diffusion coefficient. Though an interesting issue, a detailed account for intracellular diffusion and connectivity between neighboring points in an astrocyte is outside the scope of the current study and here we resort to a rather minimalistic description.
Following the study of diffusion coefficients of IP3 and Ca2+ in Xenopus laevis oocytes (Allbritton et al., 1992), most modeling studies assume a much faster diffusion of IP3 (≈ 300 μm2/s) than Ca2+ (≈ 10 μm2/s) due to Ca2+ buffering and often disregard Ca2+ diffusion altogether. However, the effective rate of IP3 diffusion can occur on a much slower timescale (Dickinson et al., 2016), equalizing the signaling range and speed of the two signaling factors. Moreover, gap junctions formed by connexin43, characteristic for astrocytes, are permeable to Ca2+ (De Bock et al., 2012). We thus account for intra- and intercellular diffusion of both IP3 and Ca2+ in the model with the same effective diffusion coefficients. This includes exchange at borders between neighboring astrocytes to imitate the function of gap junctions, which is supported by evidence that IP3 can diffuse through the gap junctions along with Ca2+ (Yule et al., 1996).
Specifically, the diffusive term in Equation (1), e.g., for Ca2+ reads:
where is the diffusion rate defined as the diffusion coefficient for Ca2+ scaled by spatial grid step δx and local AVF value r: , and i enumerates the N nearest neighboring astrocyte-containing units. Here, we regard diffusion in porous media and assume that larger AVF is associated with larger cross-sectional area open for diffusion, as the astrocyte process diameter increases, and simultaneously with less tortuous paths taken up by diffusing molecules as the processes become less entangled. This leads to approximately quadratic scaling of D* with r. The diffusive term for IP3 is defined in a similar way.
Finally, neighboring pixels with different AVF values obviously contain unequal volumes of astrocyte cytoplasm; hence, small concentration changes in areas with high AVF should cause larger diffusion-mediated concentration changes in the neighboring pixels with low AVF. This was accounted for by scaling the concentration rates of change due to diffusive exchange by the ratio of the AVF values of the two neighboring pixels.
2.5. Model Parameters and Numerical Details
The basic set of model parameters is given in Table 1. Only new parameters and values different from that in Ullah et al. (2006) are shown. For convenience, the full set of parameters is provided in Supplementary Table 1. The few values that are different were adjusted in order to provide the reasonable dynamics with the introduced treatment of as a variable in our model and the spatially extended layout. Diffusion coefficient for Ca2+ is taken as a lower-bound estimate in Allbritton et al. (1992). Slow diffusion coefficient of IP3 is based on Dickinson et al. (2016). We also added new parameters, specifically, A, τGlu, and DGlu. The latter was chosen as a small value to describe only minimal spillover from a release site and buffering by binding to transporters. The pair of parameters describing instantaneous glutamate release rate and slower decay could be varied, because it is hard to assess the actual transmitter concentration and decay time as sensed by astrocyte leaflets. Extracellular glutamate transients occurring due do quantal synaptic release as estimated by fluorescent glutamate sensor have decay timescale in close to 100 ms (Jensen et al., 2019), and this value was used for the simulations shown below. This led to local glutamate transients peaking at 1.2 μM and decaying within 200 ms. We note that qualitatively similar Ca2+ signaling dynamics could be obtained with a shorter τGlu value, compensated by higher release rate A.
Numerical integration of the model differential equations is done in an explicit scheme (4th order Runge–Kutta method adopted for stochastic differential equations with a fixed timestep dt = 0.002 s) implemented in AGEOM–CUDA software (Postnov et al., 2012). Spatial grid step was δx = 0.275 μm/pixel for single-cell templates and δx = 0.55 μm/pixel for network templates. For reproducibilty, a reference implementation of spatial template generation and model simulation is available at https://zenodo.org/record/4552726#.YDAz1nUzZQ8 in form of Jupyter notebooks, Python and C code.
To quantify spatiotemporal properties of the simulated Ca2+ dynamics, we examined complementary cumulative distribution functions (CCDF) of areas and lifetimes of individual Ca2+ transients in all single-cell spatial templates to avoid selection bias. Ca2+ transients were thresholded at 25% change from the local baseline level. The resulting contiguous TXY volumes of suprathreshold Ca2+ concentration were treated as discrete events. CCDF of some random variable X is defined as the probability P that X is greater than some x value: . We present the CCDF curves in double logarithmic coordinates to test if they can be approximated by a straight line, implying a power-law behavior. If a given variable is distributed according to a power law with probability density function (PDF) , then the CCDF also has a power-law behavior, but with a smaller exponent .
3. Results
The proposed model, including the modifications to the local calcium dynamics and spatial mapping, was tested in a number of simulation experiments with different parameter settings and different spatial templates (which we call “cells” for shorthand below). To test for agreement between model behavior and the experimentally observed dynamics, first, we looked at the effect of the level of mean neuronal firing rate (local rate of the Poisson point process in terms of our model) on spatio-temporal dynamics of astrocytic calcium, and second, we tested whether the artificial spatial templates could provide for realistic intercellular calcium waves or other collective variants of astrocytic calcium dynamics.
3.1. Wave Patterns in Single-Cell Templates
Figure 4 summarizes simulations of the 27 single-cell templates shown in Figure 2. At low excitation (psyn = 0.005 Hz) most Ca2+ events were spatially confined and tended to start at a small number of sites, as shown with max-span contours and red labeling in Figure 4A (left) for 25 largest events. At higher excitation (psyn = 0.01 Hz) many Ca2+ events spread to occupy the whole cell domain and again tended to initiate at the same sites. Synaptic signaling events were integrated into a spatial glutamate profile as shown in Figure 4A (bottom): local surges of extracellular glutamate are sparse at low excitation, while their instantaneous spatial density increases at high excitation, with a tendency of nearby sparks to blend.
Figure 4. Simulations of single-astrocyte spatial templates. (A) top row: Ca2+ transient initiation sites (red) and maximum span contours at low (left, psyn = 0.005 Hz) and high (right, psyn = 0.01 Hz) synaptic drive, in both cases 25 largest events are shown, at high drive all such events span the whole template; bottom row: snapshots of instantaneous extracellular glutamate concentrations at low and high synaptic drive parameters. Scale bar: 25 μm. (B) Effect of synaptic drive on Ca2+ transient frequency and sizes (n = 27 templates, simulation time 2,500 s after burn-in period of 2,000 s); top row: number of events covering more than 25% of cell area (“large” events) increases with excitation strength (left), number of events covering less than 25% of cell area (“small” events) decreases (right); bottom row: average areas of both “large” (left) and “small” (right) events increases with excitation strength in most templates. (C) Distribution of baseline [Ca2+]i levels with AVF parameter at low and high stimulation drives (left) and an example of spatial distribution (right). Each transparent line corresponds to one template, thick lines: average. (D) Same as in (C), but for [IP3]i . (E) Evolution of a single Ca2+ transient starting in top-right corner of the template; top row: [Ca2+]i , middle row: [IP3]i , bottom row: relative change in [IP3]i .
The tendency, exemplified by a single template in Figure 4A, was supported by the majority of single-cell templates (Figure 4B): the number of events (during 2,500 s simulation time) covering more than 25% of the cell area increased with excitation for nearly all cells except six, which were incapable of generating whole-cell transients at psyn = 0.01 Hz. There were no obvious differences between these cells and all the others in overall morphology or AVF statistics. All cells did generate whole-cell transients at a higher excitation of psyn = 0.02 Hz. Most cells demonstrated a decline in the number of small events, covering less than 25% of the cell area with excitation, as a larger proportion of events was enabled to spread over larger areas, while the rate of event initiation could remain stable. The area of large (> 25%) events increased for all cells which generate large events under low drive conditions except the three, which did not generate large events at all, apart from other cells. The average area of the small (< 25%) events increased with excitation strength for all cells.
Stochastic local glutamate surges initiate two parallel processes: fast localized Ca2+ transients and slower IP3 production. Both integrate over time to steady-state levels of the model variables. Because the relative input of plasma membrane transport is defined by AVF in our model, we expected that the steady state levels and the probability of Ca2+ event initiation should depend on AVF as well. Steady-state values of [Ca2+] and [IP3] decreased with growing AVF (Figures 4C,D), forming an uneven spatial profile. Calcium, as well as IP3 levels, were higher in the periphery and lower near the soma, which is due to Ca2+ entry during the synaptic excitation and due to higher IP3 production in the regions with lower AVF. Different cells varied in steady-state levels of the variables, while intensified stimulation lead on average to a slight elevation of steady-state [Ca2+]i, due to increased Ca2+ entry and did not affect [IP3]i, as the increase in [Ca2+]i was insufficient for activation of PLCδ.
An example of a single large Ca2+-event is shown in Figure 4E. The expanding wave of elevated [Ca2+]i initiates in a small location and spreads over the whole spatial domain (top row). Due to the large range of steady-state [IP3]i concentrations, the event is unclear in absolute [IP3]i values (middle row), but is obvious in the relative scale (bottom row).
Despite the stochastic nature of excitation, Ca2+ activity in most cells is self-organized in a repeated pattern of Ca2+ transient initiation and spreading (Figure 5); event initiation sites were tightly clustered. Interestingly, activation in some clusters lead to spatially confined events, unlinked to activity in the rest of the cell, while transients originating in other sites tended to spread over the whole spatial domain in a repeated fashion. This is illustrated in Figures 5A–D for an example spatial template (see also Supplementary Video 2). This cell is markedly anisotropic, which defines the dominant wave spreading properties. The two active initiation sites labeled as #1 and #2 display different properties: events, starting in the site #1 often spread over large portions of the cell, as shown for a line-scan path in Figure 5B, while line-scan along the path starting in #2 was either activated by a wave coming from #1 or—very locally and with a higher frequency—by confined transients initiating in #2. Averaging small temporal windows around Ca2+ spikes at the origin of path #1 shows that activation along this path is time-locked to activation of the initiation site. On the other hand, averaging similar temporal windows triggered by Ca2+ spike at the origin of path #2 did not reveal any structured activation patterns. Figure 5D shows max-span contours of 25 largest events with their initiation sites mapped in red, as well as 5 peak-delay maps for a repeated pattern of activation. In these maps color indicates delay in seconds between the Ca2+ peak at the initiation site and the Ca2+ peak at each point of the cell. The precise wave initiation site could vary within 3–5μm, but the overall spatiotemporal pattern remained similar, with activation spreading mainly along the long thick processes toward the soma. Additional examples of repeated patterns for other cells are provided in Figure 5E. In some cells preferred wave initiation sites could alternate between two polar positions or a few neighboring regions. There was also some scatter in the maximum delay between the initiation of the wave and its full expansion.
Figure 5. Ca2+ dynamics in single-cell templates self-organizes in repeatable spatiotemporal patterns. (A) A spatial template for simulation and 2 path-scans (#1, purple and #2, yellow) used for rasters in (B); (B) Rasters of [Ca2+]i and relative change in [IP3]i along paths; path starting points are shown as magenta squares. (C) Transient-triggered averages, linked to the [Ca2+]i peaks at the origin of path #1 (left) and to the [Ca2+]i peaks at the origin of path #2 (right); initiation of a [Ca2+]i peak at the origin of path #1 typically leads to a full-cell [Ca2+]i wave, while [Ca2+]i transients emerging at the origin of path #2 remained localized. (D) Left: clustered initiation points (red) and contours of 25 biggest [Ca2+]i waves; right: 5 “large” [Ca2+]i waves with color-coded delay before reaching the peak in [Ca2+]i show a tendency to start in the same area and spread with similar spatiotemporal profiles. (E) Examples of repeated spatiotemporal patterns in 4 other spatial templates: in some cells there was more than one preferred initiation site, the full delay before initiation and waning of the wave also varied.
Though localized Ca2+ events could be initiated in the low-AVF regions due to direct influx through the plasma membrane, these event seeds needed to reach a tip of a thicker branch with higher AVF to be amplified by IP3-mediated ER exchange. Thus, initialization of a global Ca2+ wave critically depends on a coincidence of exactly the right spot in the AVF profile—allowing both for a high enough [IP3]i baseline and sufficient ER exchange—and a wide and long enough cluster of glutamate release due to local increase in synaptic activity. Areas containing only thin processes with low AVF will display only frequent local Ca2+ sparks, unable to invade the neighboring regions and thus will primarily set the baseline levels of [Ca2+]i and [IP3] due to diffusion. Indeed, astrocytes in hippocampal slices display frequent localized Ca2+ events in the cell periphery, often termed “microdomains,” with a characteristic size of high-Ca2+ spots much smaller than the cell size and originating in the thin processes region (Rungta et al., 2016). At the same time, the “lifespan” of the Ca2+ wavefront increases with the already invaded area of the thick branch region due to regenerative Ca2+ release, which in turn relies on background [IP3] level. Thus, the specific topology of the cell template predicts the “hot spots” for the probability of Ca2+ wave seeding.
We next describe statistics of areas and lifetimes of individual Ca2+ transients (thresholded at 25% deviation from baseline). The left column of Figure 6 provides CCDF of areas, covered by individual events, while the right column describes durations of calcium events. Red curves correspond to the case of low neuron activity, the blue ones correspond to high neuron activity. In both cases the CCDF curves are presented in double logarithmic coordinates and can be approximated by a straight line within some ranges of event areas or durations, suggesting a power-law behavior in agreement with experimental data (Wu et al., 2014). After recalculating from CCDF to PDF exponents, the resulting parameters were different from that reported in Wu et al. (2014): for areas, α ≈ 3.3…4.0 in the model vs. α ≈ 2.1…2.4 in cultured astrocytes, and for durations α ≈ 4.5…5.1 in the model vs. α ≈ 1.97…2.16. These discrepancies can be explained by imperfection of the model and the 2D spatial embedding in the model. Defining model parameters that govern the shape of the event size and duration distributions is a potentially interesting outlook for further studies. Increase in synaptic excitation favored larger areas occupied by waves and longer durations of each event. The kinks in the CCDFs for event areas correspond to transition to whole-cell waves, i.e., events covering more than 30–50 μm2 were likely to expand further and cover the whole cell.
Figure 6. Complementary cumulative distribution functions for areas (left) and durations (right) of Ca2+ events in all single-cell templates. Red lines—low excitation (psyn = 0.005 Hz), blue lines—high excitation (psyn = 0.01 Hz). Event area CCDFs: slopes of the fits for the low and high excitations are 3.0 and 2.3, correspondingly; the bend at the large areas corresponds to transition to whole-cell activation. Event duration CCDFs: slopes of the fits are 4.1 and 3.5 for the low and high excitations, correspondingly.
We so far examined Ca2+ transients in 2D spatial templates, which we use to generate astrocytic network in the next section to reduce computational load. We note however that the presented modeling approach is extendable to 3D with minimal modifications. An example of simulation for a 3D single-cell spatial template is presented in Supplementary Video 3. The features reported above for 2D templates remain in 3D: there are Ca2+ transients at different spatial and temporal scales, waves that expand toward soma and vice versa, there are visible repeated patterns of activation.
3.2. Collective Dynamics of Astrocytes
After testing model behavior at microdomain and single-cell level, we turned to ensembles of connected cells. Collective dynamics of astrocytes was simulated using spatial templates containing about 40 cells as shown in Figure 7A; see also Supplementary Video 4. The simulated network is smaller than the size of cliques or networks of connected astrocytes in the neocortex (Houades et al., 2008), but still at the same order of magnitude. Simulations of larger spatial network templates are also possible but are more computationally demanding. Ca2+ activity was not uniform across the spatial template: there were active “hot” brims in the periphery of most of the cell domains and some cells were less active than others, as shown in Figure 7B. We observed Ca2+ transients originating in some astrocytes and expanding to their neighbors in a wavelike manner, an example of such a wave is shown in Figure 7C as a sequence of contours of elevated [Ca2+]i, separated by 1 s.
Figure 7. Ca2+ activity in multicellular template. (A) Layout of individual spatial domains, each cell is color-coded. (B) Activity level is not uniform: color-coded percentage of time that each pixel had [Ca2+]i above 0.1μM. Domain periphery is more active than somatic regions. (C) Example of spatiotemporal evolution of two co-occurring waves. Contours are separated by 1 s, time delay since the first contour is color-coded. (D) Cellwise dynamics: Ca2+ and IP3 values averaged over individual cell domains, line color corresponds to the map in (A), visible are single-cell events as well as packed Ca2+ transients representing multi-cellular waves. (E,F) Repeated patterns of cell activation. (E) Spike-triggered averages of cellwise [Ca2+]i profiles initiated by the cells indicated as #, *, and & in (A,F); activity of the first cell is time-locked to activation of a single other cell, Ca2+ transient the second cell consistently leads to activation of several cells, while the third cell does not participate in repeatable patterns. (F) Local score of activation pattern repeatability based on approach shown in (E) (see text for more details); activation of a subpopulation of cells leads to repeated activation of its neighbors; red regions denote wave initiation sites as in Figure 4A.
To simplify the description of emerging spatiotemporal patterns, each cell can be considered as an element that “fires,” i.e., producing a global whole-cell spike, or remains silent. This cell-wise activity is shown for [Ca2+]i and [IP3]i in Figure 7D, where line colors correspond to domain colors in Figure 7A. We observed individual cell spikes as well as packs of tightly grouped spikes corresponding to multi-cellular waves. A considerable scatter is clear in the baseline levels of IP3, which reflects individual properties of each cell. IP3 concentration peaks are wider and occur with a small delay in comparison to Ca2+ peaks, which reflects the slower kinetics of IP3 production timescale and a lag due to diffusion of IP3 from the periphery to the somatic region.
A large Ca2+ transient in one astrocyte can spread to other cells. In single cell templates we observed self-organized repeated patterns of activation. We were curious, if there will also emerge repeated patterns of intercellular dynamics at a network level. As a simple test for repeatable patterns, we calculated averages of the domain-averaged Ca2+ rasters in a short time window, triggered by Ca2+ spikes in different domains. We collected Ca2+ traces, where each spatial domain served as a large region of interest (ROI), selected one cell as a seed, created time windows around Ca2+ spikes in this cell and averaged Ca2+ dynamics snapshots in all other ROIs across the time windows. In the case of stable activation sequence, the spikes appear with the same time-lag relative to the seeding astrocyte, and are thus visible in the raster plot, while if spiking in other astrocytes is not time-locked to the seeding astrocyte, the average Ca2+ signal will be faint. Examples of such spike-triggered averages are shown in Figure 7E for three cells labeled as “#,” “*,” and “&” in Figure 7A used to center the time windows. Here, Ca2+ spikes in cell “#” was time-locked to activation of a single other cell, activation of cell “*” lead to repeated activation of several other cells with a stable delay, while activation of cell “&” did not repeatedly lead to activation of other cells.
Inspired by this cell-wide “repeatability” measure based on the contrast of spike-triggered averages, we defined a similar score for more a detailed mapping of whether activation in some region repeatedly lead to activation in other areas with stable time lags. To this end, we split the spatial domain into overlapping square windows of size 5 × 5 μm and extracted Ca2+ dynamics from these patches. We then selected the patches where there were more than 10 Ca2+ spikes reaching at least 0.5 μM [Ca2+]i and created spike-triggered 50 s-long averages from the raster of Ca2+ signals in all patches. Percentage of all points with [Ca2+]i > 0.2 μM in such spike-triggered windows was used as the repeatability score for a given patch. The patches were then projected back onto the spatial template, with averaging of the values in overlapping areas between patches. This resulted in an automated mapping of areas leading to repeated downstream activation at the network level, revealing cells, serving as hubs in spreading multicellular events (Figure 7F). As suggested in Brazhe et al. (2018) for a simpler model, some spatial configurations of thick branches and leaflets can trigger persistent pacemaker-like activity, taking over the control of dynamics at the network level. We thus tested if alterations to the spatial template could lead to self-organization in a different spatial pattern of the centers of high repeatability (Figure 8). The cell labeled as “*” in Figure 7 was often activated from its direct neighbor to the right. Unlinking this cell from the network by setting to zero all contacts at domain boundaries (Figure 7A) lead to a change in the repeatability map, decreasing the score for the cell “*” and increasing it for the three cells in the center. Unlinking the cell “*” from the network effectively silenced it, leaving only three cells with relatively high repeatability scores in the template.
Figure 8. Modification of spatial templates changes activation patterns. Color-coding: repeatability score, red: wave initiation sites. (A,B) Unlinking single cells from the neighbors; Blue contours—isolated cells, red regions denote wave initiation sites as in Figure 4A. Isolation of the cell often activated the cell denoted “*” in Figure 7 lowers repeatability of the latter. (B) Isolation effectively silences the cell denoted “*” in Figure 7, which was a center of repeated patterns before. (C,D) Modifications of wave initiation sites. (C) Averaging AVF values within the activation sites increases repeatability of the waves, starting at the cells, isolated in (A,B) and one other cell. (D) Removing cell content from the spatial template within active initiation sites leads to a lower number of wave initiation sites and increased repeatability in several cells.
A subtler alteration of the spatial template can be directed at the sites of frequent Ca2+ transient initiation sites shown in red in Figure 8. Substituting the natural AVF profile at these sites with the average AVF value (Figure 8C) or ablating cell content from these areas (Figure 8D) lead to a dramatic reorganization of the Ca2+ initiation sites and the centers of high repeatability. In both cases the number of initiation sites was reduced, with some of the former sites being silenced and some new sites formed in the neighborhood of the previous sites. Also in both cases the inequality of repeatability score was increased, with a few cells showing very high values. We attribute this to the reduced number of initiation sites, leading to a more repeated activation sequences.
4. Discussion
We proposed a spatially detailed model of astrocytic calcium activity, which reflects current understanding of the two distinct mechanisms of Ca2+ dynamics: excitable IP3-mediated exchange with ER in astrocyte soma and branches and plasma membrane exchange in the fine astrocytic processes and leaflets, sensitive to external conditions. Specifically, we suggest (i) an algorithm for data-based generation of 2D spatial templates matching realistic astrocyte morphology, and (ii) morphology-dependent spatially non-uniform parameter landscape for the calcium dynamics. To this end, we introduce the AVF parameter, which sets locally the relative input of the plasma membrane and ER based pathways and scales effective intracellular diffusion coefficients. The central idea underlying this separation is that astrocytes “sense” synaptic activity with fine processes, and it is where Ca2+ transients are relying on extracellular Ca2+ rather than intracellular stores, and where the bulk of IP3 can be produced, while thicker branches and somata provide the positive feedback gain mechanism for IP3-mediated Ca2+-induced Ca2+ release from ER. This mechanism separation is directly mapped to cell morphology in our approach.
We tested the suggested framework both at individual cell level and for algorithmically created multicellular astrocytic network templates. Our results show that the model is able to reproduce characteristic spatiotemporal patterns of Ca2+ dynamics driven by synaptic activity, represented by spatially uncorrelated point-sources of glutamate release coupled to focal Ca2+ entry, triggered by independent stochastic Poisson spike trains. At the single-cell level, the statistics of Ca2+ event durations and expansion areas turned out to have a power-law distribution resembling experimental data (Wu et al., 2014). Power-law statistics of the Ca2+ transients does not directly follow from the model equations and is an emergent property of the interplay between astrocyte morphology and Ca2+ dynamics.
The presented model is a rather simplified representation of native astrocyte morphology and Ca2+ dynamics. It is less detailed, but also less computationally demanding than the framework proposed by Savtchenko et al. (2018), allowing for simulations on a GPU-equipped laptop rather than on a supercomputer or a cloud. A major simplification of our model, dictating its limitations, is the reduction of real 3D astrocytic morphology to flat 2D spatial templates. The flattening was primarily done for the sake of computational tractability, but also conceptually matches single-plane imaging regime. The emergence of repeated activation patterns should not depend on the embedding space dimensionality, although the repeated propagation patterns can become more complex and elaborate in 3D; one can also expect different expansion rates and initiation probabilities for Ca2+ events in 3D. Notwithstanding, we argue that using 2D patterns can be a useful approximation. First, some astrocyte network systems can be regarded as effectively two-dimensional, e.g., astrocyte cultures or retinal astrocyte networks. Second, “true” astrocyte morphology can be regarded as less than three-dimensional, with the degrees of freedom limited by branching connectivity of the astrocytic processes, creating more or less independent astrocyte lobes and sub-domains. A 2D morphology-based modeling approach has earlier been employed by Kang and Othmer (2009) to study calcium waves in retinal astrocytes. Our interpretation of astrocyte morphology is however different from the one in Kang and Othmer (2009) in several key aspects. First, Kang and Othmer used a GFAP-positive immunostained micrograph of retinal astrocytes, but GFAP stains only a small fraction of astrocyte cytosol volume, while the rest, especially the mesh of thin processes, got excluded from consideration of calcium dynamics, leading for the reconstructed morphology to resemble that of cultured astrocytes (see e.g., Wu et al., 2014 for comparison of morphology). On the other hand, an account for spatial segregation of calcium activity in fine astrocytic processes and thick branches was key to this study. Second, we account for a morphology-based balance between Ca2+ entry via plasma membrane and from ER, which is absent in the work of Kang and Othmer. Third, we use a larger spatial template and describe Ca2+ dynamics on a larger spatial scale. Our model is also more focused on intercellular communication through gap junctions, as we do not account for extracellular purinergic signaling to explain Ca2+ wave propagation, which was done in earlier works focused on retinal astrocytes (Kang and Othmer, 2009; Edwards and Gibson, 2010), keeping in mind that the role of the extracellular pathway can be more pronounced in the retina than in neocortex.
Another limitation comes from our approximation of the 3D mesh of fine astropil sponge by a continuous active medium, parameterized by astrocytic cytosol volume fraction. Effectively, we “glue” together the individual branches and leaflets, assuming that they are at least partly interconnected by autologous gap junctions and branch-to-branch loops. The idea of “loopy” or sponge-like organization of the astropil has experimental support (Wolff et al., 1998; Genoud et al., 2015; Arizono et al., 2020), and we therefore adopt the AVF framework to represent unresolved astrocyte processes, also accounting for the tortuosity of the sponge by AVF-dependent scaling diffusion coefficients.
Indeed, the mentioned simplifications are expected to limit the predictive power of the model with regard to event frequencies, scaling characteristics and propagation speed. On the other hand, the simulated patterns of single-cell calcium transients are qualitatively similar to that observed in a single focal plane, which suggests that this reduction seems to preserve main features of astrocyte dynamics, while it is worth to be investigated further at sub-cellular spatial scales in future work. A proof-of-concept possibility to use the proposed modeling framework to simulate Ca2+ signaling in 3D is presented in Supplementary Video 3.
A notable feature of simulated Ca2+ dynamics in this system is spontaneously emerging stable patterns in both initiation and propagation of calcium transients, in tune to the co-active neuronal and astrocytic cells or repeating sequences of neuronal activation reported in slices (Sasaki et al., 2011, 2014). In agreement with Brazhe et al. (2018) we observed morphology-dependent emergence of hotspots with persistent pacemaker-like activity, taking over control of the dynamics at larger scales. In single-cell templates these preferred initiation sites could lead to activation of either spatially confined microdomains or larger expanding Ca2+ areas, covering up to the whole cell. In multicellular systems we observed self-organized patterns of repeated calcium activity involving multiple cells. All cells in the template sharing the same equations and parameters, local differences in morphology of single astrocytes and geometry of astrocyte-to-astrocyte contacts favored initiation of multicellular Ca2+ waves in some cells, followed by repeated sequences of cell activation as the Ca2+ wave sweeped across the network. Excluding some cells from the original template caused dramatic reshaping of repeating activation patterns; removing or mashing up cell content in event initiation hotspots effectively reduced the number of active initiation sites and led to more stereotyped network activity. Our simulations showed that different cells in the network templates had different levels of activity and could develop patterns of Ca2+ events with preferred directionality. Recently, Wang et al. (2019) have demonstrated heterogeneity of individual neocortical astrocytes with respect to the properties of their Ca2+ activity. They also report a tendency for anatomical directionality during network-wide bursts of Ca2+ activity, accompanying locomotion and presumed to result from adrenergic input from locus coeruleus. In the light of the present study, it is interesting to ask, whether the experimentally observed heterogeneity of astrocytic signaling is defined by different patterning and levels of local synaptic activity or by individual astrocyte morphology, or the interplay of both? At a larger spatial scale it is interesting whether the observed population-wide directionality is a product of the afferentation delays from locus coeruleus or some previously unstudied anatomic directionality of astrocyte morphology? Going even further, it seems an exciting possibility that—if Ca2+ activity indeed affects LTP and Ca2+ activity patterns are shaped by astrocytic morphology template—some part of memory in the neural network can be engraved in the fingerprint of astrocyte morphology. We conclude that with the same parameter set, the specific dynamical regime and role of an individual cell within astrocytic network is to a large extent defined by its morphology and this may have implications for computational performance of the underlying neuronal network.
The concept of spatially segregated oscillator regarded in (Brazhe et al., 2018) attracts growing interest as a new class of bio-inspired dynamic models. For example, Vanag (2020) explores a model of spatially segregated Belousov-Zhabotinsky oscillating reaction, where the catalyst (analogous to calcium induced calcium release in Ca2+ models) is confined to small beads, whereas the rest chemical components of the reaction can diffuse freely, which leads to complex dynamical regimes and interaction between neighboring beads with immobilized catalyst. Another similar concept is that of volume transmission and intercellular communication via diffusing signaling molecules secreted by excitable cells (Sykova and Nicholson, 2008). There also appears to be an interesting connection between the self-organizing Ca2+ signaling in coupled astrocytes and the concept of “reaction-diffusion computing” or “chemical computing” based on oscillatory chemical (e.g., Belousov-Zhabotinsky) systems in coupled reactor volumes, e.g., microemulsions (Adamatzky et al., 2005; Epstein et al., 2012; Showalter and Epstein, 2015; Torbensen et al., 2017; Vanag, 2019; Proskurkin et al., 2020). It is conceivable that astrocytes can provide a substrate for such Ca2+-based reaction-diffusion computing in parallel to, and organizing sparse network-based neuronal connectivity. Interestingly, it is common for self-organizing spatially distributed chemical systems to rely on inhibitory diffusive coupling (Li et al., 2015), while only excitatory coupling of neighboring astrocytes was modeled. Finding a mechanism for negative diffusive coupling, competing with the excitatory one can be an interesting outlook stemming from this analogy. Finally, it is inviting to believe that a combination of “analog” chemical computing approaches and “digital” neuronal ones, as well as a combination of graph-based and neighbor-based connectivity will inspire new algorithms for machine learning.
Data Availability Statement
The model code and data, used for network spatial template generation and simulations are available for download at https://zenodo.org/record/4552726#.YDAz1nUzZQ8 (doi: 10.5281/zenodo.4552726).
Author Contributions
DV: computer simulations, analysis of simulation results, figure preparation. AV: computer simulations, analysis of simulation results, figure preparation, writing the manuscript. DP: program design and architecture for simulations, writing the manuscript. AB: conceptualization of the model, network generation algorithm, figure preparation, writing the manuscript. All authors contributed to the article and approved the submitted version.
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Acknowledgments
DP and AB acknowledge the support from the Russian Foundation for Basic Research, grant 19-515-55016, which covered the development of the modeling approach for spatial segmentation of inter-astrocyte functions. AV and AB acknowledge support from Russian Science Foundation, grant 17-74-20089, which covered development of the numeric algorithm for creation of realistic spatial templates, biophysical interpretation of the model, and the collective astrocyte signaling analysis. AB also acknowledges support from Interdisciplinary Scientific and Educational School of Moscow University Molecular Technologies of the Living Systems and Synthetic Biology, which covered 3D model simulations. DV acknowledges support from the Russian Foundation for Basic Research, grant 16-32-50221 for mobility of young scientists, which covered model interpretation in terms of non-linear dynamics, CUDA-based simulations of the model, and analysis of single-cell dynamics.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fncel.2021.645068/full#supplementary-material
References
Abbott, N. J., Pizzo, M. E., Preston, J. E., Janigro, D., and Thorne, R. G. (2018). The role of brain barriers in fluid movement in the CNS: is there a “glymphatic” system? Acta Neuropathol. 135, 387–407. doi: 10.1007/s00401-018-1812-4
Adamatzky, A., De Lacy Costello, B., and Asai, T. (2005). Reaction-Diffusion Computers. Amsterdam: Elsevier.
Agarwal, A., Wu, P. H., Hughes, E. G., Fukaya, M., Tischfield, M. A., Langseth, A. J., et al. (2017). Transient opening of the mitochondrial permeability. transition pore induces microdomain calcium transients in astrocyte processes. Neuron 93, 587.e7–605.e7. doi: 10.1016/j.neuron.2016.12.034
Allbritton, N. L., Meyer, T., and Stryer, L. (1992). Range of messenger action of calcium ion and inositol 1,4,5-trisphosphate. Science 258, 1812–1815. doi: 10.1126/science.1465619
Alvarellos-González, A., Pazos, A., and Porto-Pazos, A. (2012). Computational models of neuron-astrocyte interactions lead to improved efficacy in the performance of neural networks. Comput. Math. Methods Med. 2012:476324. doi: 10.1155/2012/476324
Araque, A., Carmignoto, G., Haydon, P. G., Oliet, S. H., Robitaille, R., and Volterra, A. (2014). Gliotransmitters travel in time and space. Neuron 81, 728–739. doi: 10.1016/j.neuron.2014.02.007
Arizono, M., Inavalli, V., Panatier, A., Pfeiffer, T., Angibaud, J., Levet, F., et al. (2020). Structural basis of astrocytic Ca. Nat. Commun. 11:1906. doi: 10.1038/s41467-020-15648-4
Ascoli, G. A., Donohue, D. E., and Halavi, M. (2007). Neuromorpho.org: a central resource for neuronal morphologies. J. Neurosci. 27, 9247–9251. doi: 10.1523/JNEUROSCI.2055-07.2007
Bazargani, N., and Attwell, D. (2016). Astrocyte calcium signaling: the third wave. Nat. Neurosci. 19, 182–189. doi: 10.1038/nn.4201
Bindocci, E., Savtchouk, I., Liaudet, N., Becker, D., Carriero, G., and Volterra, A. (2017). Three-dimensional Ca2+ imaging advances understanding of astrocyte biology. Science 356:eaai8185. doi: 10.1126/science.aai8185
Bojarskaite, L., Bjoernstad, D. M., Pettersen, K. H., Cunen, C., Hermansen, G. H., Aabjoersbraaten, K. S., et al. (2020). Astrocytic ca2+ signaling is reduced during sleep and is involved in the regulation of slow wave sleep. Nat. Commun. 11:3240. doi: 10.1038/s41467-020-17062-2
Brazhe, A., Mathiesen, C., and Lauritzen, M. (2013). Multiscale vision model highlights spontaneous glial calcium waves recorded by 2-photon imaging in brain tissue. Neuroimage 68, 192–202. doi: 10.1016/j.neuroimage.2012.11.024
Brazhe, A. R., Postnov, D. E., and Sosnovtseva, O. (2018). Astrocyte calcium signaling: interplay between structural and dynamical patterns. Chaos 28:106320. doi: 10.1063/1.5037153
Bushong, E. A., Martone, M. E., and Ellisman, M. H. (2004). Maturation of astrocyte morphology and the establishment of astrocyte domains during postnatal hippocampal development. Int. J. Dev. Neurosci. 22, 73–86. doi: 10.1016/j.ijdevneu.2003.12.008
Bushong, E. A., Martone, M. E., Jones, Y. Z., and Ellisman, M. H. (2002). Protoplasmic astrocytes in ca1 stratum radiatum occupy separate anatomical domains. J. Neurosci. 22, 183–192. doi: 10.1523/JNEUROSCI.22-01-00183.2002
Carmignoto, G. (2000). Reciprocal communication systems between astrocytes and neurones. Prog. Neurobiol. 62, 561–581. doi: 10.1016/S0301-0082(00)00029-0
Cornell-Bell, A. H., Finkbeiner, S. M., Cooper, M. S., and Smith, S. J. (1990). Glutamate induces calcium waves in cultured astrocytes: long-range glial signaling. Science 247, 470–473. doi: 10.1126/science.1967852
Cuntz, H., Forstner, F., Borst, A., and Hausser, M. (2010). One rule to grow them all: a general theory of neuronal branching and its practical application. PLoS Comput. Biol. 6:e1000877. doi: 10.1371/journal.pcbi.1000877
Dani, J. W., Chernjavsky, A., and Smith, S. J. (1992). Neuronal activity triggers calcium waves in hippocampal astrocyte networks. Neuron 8, 429–440. doi: 10.1016/0896-6273(92)90271-E
De Bock, M., Wang, N., Bol, M., Decrock, E., Ponsaerts, R., Bultynck, G., et al. (2012). Connexin 43 hemichannels contribute to cytoplasmic ca2+ oscillations by providing a bimodal ca2+-dependent ca2+ entry pathway. J. Biol. Chem. 287, 12250–12266. doi: 10.1074/jbc.M111.299610
De Pitta, M., Brunel, N., and Volterra, A. (2016). Astrocytes: Orchestrating synaptic plasticity? Neuroscience. 323, 43–61. doi: 10.1016/j.neuroscience.2015.04.001
De Pitta, M., Goldberg, M., Volman, V., Berry, H., and Ben-Jacob, E. (2009). Glutamate regulation of calcium and ip3 oscillating and pulsating dynamics in astrocytes. J. Biol. Phys. 35, 383–411. doi: 10.1007/s10867-009-9155-y
De Young, G. W., and Keizer, J. (1992). A single-pool inositol 1,4,5-trisphosphate-receptor-based model for agonist-stimulated oscillations in Ca2+ concentration. Proc. Natl. Acad. Sci. U.S.A. 89, 9895–9899. doi: 10.1073/pnas.89.20.9895
Dickinson, G. D., Ellefsen, K. L., Dawson, S. P., Pearson, J. E., and Parker, I. (2016). Hindered cytoplasmic diffusion of inositol trisphosphate restricts its cellular range of action. Sci. Signal. 9:ra108. doi: 10.1126/scisignal.aag1625
Donohue, D. E., and Ascoli, G. A. (2008). A comparative computer simulation of dendritic morphology. PLoS Comput. Biol. 4:e1000089. doi: 10.1371/journal.pcbi.1000089
Dupont, G., Combettes, L., Bird, G. S., and Putney, J. W. (2011). Calcium oscillations. Cold Spring Harb Perspect Biol. 3:a004226. doi: 10.1101/cshperspect.a004226
Edwards, J. R., and Gibson, W. G. (2010). A model for ca2+ waves in networks of glial cells incorporating both intercellular and extracellular communication pathways. J. Theor. Biol. 263, 45–58. doi: 10.1016/j.jtbi.2009.12.002
Epstein, I. R., Vanag, V. K., Balazs, A. C., Kuksenok, O., Dayal, P., and Bhattacharya, A. (2012). Chemical oscillators in structured media. Acc. Chem. Res. 45, 2160–2168. doi: 10.1021/ar200251j
Falcke, M. (2004). Reading the patterns in living cells the physics of Ca2+ signaling. Adv. Phys. 53, 255–440. doi: 10.1080/00018730410001703159
Fiacco, T. A., and McCarthy, K. D. (2018). Multiple lines of evidence indicate that gliotransmission does not occur under physiological conditions. J. Neurosci. 38, 3–13. doi: 10.1523/JNEUROSCI.0016-17.2017
Genoud, C., Houades, V., Kraftsik, R., Welker, E., and Giaume, C. (2015). Proximity of excitatory synapses and astroglial gap junctions in layer IV of the mouse barrel cortex. Neuroscience 291, 241–249. doi: 10.1016/j.neuroscience.2015.01.051
Haas, B., Schipke, C. G., Peters, O., Söhl, G., Willecke, K., and Kettenmann, H. (2006). Activity-dependent ATP-waves in the mouse neocortex are independent from astrocytic calcium waves. Cereb. Cortex 16, 237–246. doi: 10.1093/cercor/bhi101
Hablitz, L. M., Plá, V., Giannetto, M., Vinitsky, H. S., Stæger, F. F., Metcalfe, T., et al. (2020). Circadian control of brain glymphatic and lymphatic fluid flow. Nat. Commun. 11, 1–11. doi: 10.1038/s41467-020-18115-2
Halassa, M. M., Fellin, T., Takano, H., Dong, J. H., and Haydon, P. G. (2007). Synaptic islands defined by the territory of a single astrocyte. J. Neurosci. 27, 6473–6477. doi: 10.1523/JNEUROSCI.1419-07.2007
Haydon, P. G. (2001). Glia: listening and talking to the synapse. Nat. Rev. Neurosci. 2, 185–193. doi: 10.1038/35058528
Henneberger, C., Papouin, T., Oliet, S. H., and Rusakov, D. A. (2010). Long-term potentiation depends on release of D-serine from astrocytes. Nature 463, 232–236. doi: 10.1038/nature08673
Hoogland, T. M., Kuhn, B., Göbel, W., Huang, W., Nakai, J., Helmchen, F., et al. (2009). Radially expanding transglial calcium waves in the intact cerebellum. Proc. Natl. Acad. Sci. U.S.A. 106, 3496–3501. doi: 10.1073/pnas.0809269106
Houades, V., Koulakoff, A., Ezan, P., Seif, I., and Giaume, C. (2008). Gap junction-mediated astrocytic networks in the mouse barrel cortex. J. Neurosci. 28, 5207–5217. doi: 10.1523/JNEUROSCI.5100-07.2008
Houades, V., Rouach, N., Ezan, P., Kirchhoff, F., Koulakoff, A., and Giaume, C. (2006). Shapes of astrocyte networks in the juvenile brain. Neuron Glia Biol. 2, 3–14. doi: 10.1017/S1740925X06000081
Iliff, J. J., Wang, M., Liao, Y., Plogg, B. A., Peng, W., Gundersen, G. A., et al. (2012). A paravascular pathway facilitates csf flow through the brain parenchyma and the clearance of interstitial solutes, including amyloid β. Sci. Transl. Med. 4:147ra111. doi: 10.1126/scitranslmed.3003748
Ingiosi, A. M., Hayworth, C. R., Harvey, D. O., Singletary, K. G., Rempe, M. J., Wisor, J. P., et al. (2020). A role for astroglial calcium in mammalian sleep and sleep regulation. Curr. Biol. 30, 4373–4383.e7. doi: 10.1016/j.cub.2020.08.052
Jensen, T. P., Zheng, K., Cole, N., Marvin, J. S., Looger, L. L., and Rusakov, D. A. (2019). Multiplex imaging relates quantal glutamate release to presynaptic Ca. Nat. Commun. 10:1414. doi: 10.1038/s41467-019-09216-8
Kang, M., and Othmer, H. (2009). Spatiotemporal characteristics of calcium dynamics in astrocytes. Chaos 19:037116. doi: 10.1063/1.3206698
Kirischuk, S., Heja, L., Kardos, J., and Billups, B. (2016). Astrocyte sodium signaling and the regulation of neurotransmission. Glia 64, 1655–1666. doi: 10.1002/glia.22943
Krizhevsky, A., Sutskever, I., and Hinton, G. E. (2012). “ImageNet classification with deep convolutional neural networks,” in Advances in Neural Information Processing Systems 25, eds F. Pereira, C. J. C. Burges, L. Bottou, and K. Q. Weinberger (New York, NY: Curran Associates, Inc.), 1097–1105.
Kuga, N., Sasaki, T., Takahara, Y., Matsuki, N., and Ikegaya, Y. (2011). Large-scale calcium waves traveling through astrocytic networks in vivo. J. Neurosci. 31, 2607–2614. doi: 10.1523/JNEUROSCI.5319-10.2011
Lallouette, J., De Pitta, M., Ben-Jacob, E., and Berry, H. (2014). Sparse short-distance connections enhance calcium wave propagation in a 3D model of astrocyte networks. Front. Comput. Neurosci. 8:45. doi: 10.3389/fncom.2014.00045
Lee, H. S., Ghetti, A., Pinto-Duarte, A., Wang, X., Dziewczapolski, G., Galimi, F., et al. (2014). Astrocytes contribute to gamma oscillations and recognition memory. Proc. Natl. Acad. Sci. U.S.A. 111, E3343–E3352. doi: 10.1073/pnas.1410893111
Lenk, K., Satuvuori, E., Lallouette, J., Ladron-de Guevara, A., Berry, H., and Hyttinen, J. (2020). A computational model of interactions between neuronal and astrocytic networks: the role of astrocytes in the stability of the neuronal firing rate. Front. Comput. Neurosci. 13:92. doi: 10.3389/fncom.2019.00092
Li, N., Tompkins, N., Gonzalez-Ochoa, H., and Fraden, S. (2015). Tunable diffusive lateral inhibition in chemical cells. Eur. Phys. J. E Soft. Matter 38:18. doi: 10.1140/epje/i2015-15018-3
Li, Y., and Rinzel, J. (1994). Equations for insp3 receptor-mediated [ca2+]i oscillations derived from a detailed kinetic model: a hodgkin-huxley like formalism. Proc. Natl. Acad. Sci. U.S.A. 166, 461–473. doi: 10.1006/jtbi.1994.1041
Lind, B. L., Brazhe, A. R., Jessen, S. B., Tan, F. C., and Lauritzen, M. J. (2013). Rapid stimulus-evoked astrocyte ca2+ elevations and hemodynamic responses in mouse somatosensory cortex in vivo. Proc. Natl. Acad. Sci. U.S.A. 110, E4678–E4687. doi: 10.1073/pnas.1310065110
Manninen, T., Havela, R., and andLinne, M.-L. (2017). Reproducibility and comparability of computational models for astrocyte calcium excitability. Front. Neuroinform. 11:11. doi: 10.3389/fninf.2017.00011
Marina, N., Christie, I. N., Korsak, A., Doronin, M., Brazhe, A., Hosford, P. S., et al. (2020). Astrocytes monitor cerebral perfusion and control systemic circulation to maintain brain blood flow. Nat. Commun. 11, 1–9. doi: 10.1038/s41467-019-13956-y
Martone, M. E., Gupta, A., Wong, M., Qian, X., Sosinsky, G., Ludascher, B., et al. (2002). A cell-centered database for electron tomographic data. J. Struct. Biol. 138, 145–155. doi: 10.1016/S1047-8477(02)00006-0
Martone, M. E., Tran, J., Wong, W. W., Sargis, J., Fong, L., Larson, S., et al. (2008). The cell centered database project: an update on building community resources for managing and sharing 3D imaging data. J. Struct. Biol. 161, 220–231. doi: 10.1016/j.jsb.2007.10.003
Mathiesen, C., Brazhe, A., Thomsen, K., and Lauritzen, M. (2013). Spontaneous calcium waves in bergman glia increase with age and hypoxia and may reduce tissue oxygen. J. Cereb. Blood Flow Metab. 33, 161–169. doi: 10.1038/jcbfm.2012.175
McCauley, J. P., Petroccione, M. A., D'Brant, L. Y., Todd, G. C., Affinnih, N., Wisnoski, J. J., et al. (2020). Circadian modulation of neurons and astrocytes controls synaptic plasticity in hippocampal area Ca1. Cell Rep. 33:108255. doi: 10.1016/j.celrep.2020.108255
Mestre, H., Mori, Y., and Nedergaard, M. (2020). The brain's glymphatic system: Current controversies. Trends Neurosci. 43, 458–466. doi: 10.1016/j.tins.2020.04.003
Nadkarni, S., and Jung, P. (2004). Dressed neurons: modeling neural-glial interactions. Phys. Biol. 1, 35–41. doi: 10.1088/1478-3967/1/1/004
Nadkarni, S., and Jung, P. (2007). Modeling synaptic transmission of the tripartite synapse. Phys. Biol. 4:1. doi: 10.1088/1478-3975/4/1/001
Nagy, J. I., and Rash, J. E. (2003). Astrocyte and oligodendrocyte connexins of the glial syncytium in relation to astrocyte anatomical domains and spatial buffering. Cell Commun. Adhes 10, 401–406. doi: 10.1080/cac.10.4-6.401.406
Nedergaard, M. (1994). Direct signaling from astrocytes to neurons in cultures of mammalian brain cells. Science 263, 1768–1771. doi: 10.1126/science.8134839
O'Donnell, J., Ding, F., and Nedergaard, M. (2015). Distinct functional states of astrocytes during sleep and wakefulness: is norepinephrine the master regulator? Curr. Sleep Med. Rep. 1, 1–8. doi: 10.1007/s40675-014-0004-6
Oschmann, F., Berry, H., Obermayer, K., and Lenk, K. (2017a). From in silico astrocyte cell models to neuron-astrocyte network models: a review. Brain Res. Bull. 1–9. doi: 10.1016/j.brainresbull.2017.01.027
Oschmann, F., Mergenthaler, K., Jungnickel, E., and Obermayer, K. (2017b). Spatial separation of two different pathways accounting for the generation of calcium signals in astrocytes. PLoS Comput. Biol. 13:e1005377. doi: 10.1371/journal.pcbi.1005377
Otsu, Y., Couchman, K., Lyons, D. G., Collot, M., Agarwal, A., Mallet, J.-M., et al. (2015). Calcium dynamics in astrocyte processes during neurovascular coupling. Nat. Neurosci. 18, 210–218. doi: 10.1038/nn.3906
Parpura, V., Sekler, I., and Fern, R. (2016). Plasmalemmal and mitochondrial na(+) -ca(2+) exchange in neuroglia. Glia 64, 1646–1654. doi: 10.1002/glia.22975
Patrushev, I., Gavrilov, N., Turlapov, V., and Semyanov, A. (2013). Subcellular location of astrocytic calcium stores favors extrasynaptic neuron-astrocyte communication. Cell Calcium 54, 343–349. doi: 10.1016/j.ceca.2013.08.003
Pirttimaki, T. M., and Parri, H. R. (2013). Astrocyte plasticity: implications for synaptic and neuronal activity. Neuroscientist 19, 604–615. doi: 10.1177/1073858413504999
Polavaram, S., Gillette, T. A., Parekh, R., and Ascoli, G. A. (2014). Statistical analysis and data mining of digital reconstructions of dendritic morphologies. Front. Neuroanat. 8:138. doi: 10.3389/fnana.2014.00138
Poskanzer, K. E., and Yuste, R. (2016). Astrocytes regulate cortical state switching in vivo. Proc. Natl. Acad. Sci. U.S.A. 113, E2675–E2684. doi: 10.1073/pnas.1520759113
Postnov, D., Ryazanova, L., Brazhe, N., Brazhe, A., Maximov, G., Mosekilde, E., et al. (2008). Giant glial cell: new insight through mechanism-based modeling. J. Biol. Phys. 34, 441–457. doi: 10.1007/s10867-008-9070-7
Postnov, D. E., Brazhe, N. A., and Sosnovtseva, O. V. (2011). “Functional modeling of neural-glial interaction,” in Biosimulation in Biomedical Research, Health Care and Drug Development, eds E. Mosekilde, O. Sosnovtseva, and A. Rostami-Hodjegan (Vienna: Springer), 133–151. doi: 10.1007/978-3-7091-0418-7_6
Postnov, D. E., Koreshkov, R. N., Brazhe, N. A., Brazhe, A. R., and Sosnovtseva, O. V. (2009). Dynamical patterns of calcium signaling in a functional model of neuron-astrocyte networks. J. Biol. Phys. 35, 425–445. doi: 10.1007/s10867-009-9156-x
Postnov, D. E., Postnov, D. D., and Zhirin, R. (2012). The “AGEOM_CUDA” Software for Simulation of Oscillatory and Wave Processes in Two-Dimensional Media of Arbitrary Geometry on the Basis of High-Speed Parallel Computing on Graphics Processing Unit Technology CUDA. RF registration certificate #2012610085 from 10.01.2012 (in russian), Moscow.
Postnov, D. E., Ryazanova, L. S., and Sosnovtseva, O. V. (2007). Functional modeling of neural-glial interaction. Biosystems 89, 84–91. doi: 10.1016/j.biosystems.2006.04.012
Proskurkin, I. S., Smelov, P. S., and Vanag, V. K. (2020). Experimental verification of an opto-chemical “neurocomputer”. Phys. Chem. Chem. Phys. 22, 19359–19367. doi: 10.1039/D0CP01858A
Riera, J., Hatanaka, R., Ozaki, T., and Kawashima, R. (2011). Modeling the spontaneous Ca2+ oscillations in astrocytes: Inconsistencies and usefulness. J. Integr. Neurosci. 10, 439–473. doi: 10.1142/S0219635211002877
Rojas, H., Colina, C., Ramos, M., Benaim, G., Jaffe, E. H., Caputo, C., et al. (2007). Na+ entry via glutamate transporter activates the reverse Na+/Ca2+ exchange and triggers Ca(i)2+-induced Ca2+ release in rat cerebellar type-1 astrocytes. J. Neurochem. 100, 1188–1202. doi: 10.1111/j.1471-4159.2006.04303.x
Ross, W. (2012). Understanding calcium waves and sparks in central neurons. Nat. Rev. Neurosci. 13, 157–168. doi: 10.1038/nrn3168
Rungta, R. L., Bernier, L. P., Dissing-Olesen, L., Groten, C. J., LeDue, J. M., Ko, R., et al. (2016). Ca2+ transients in astrocyte fine processes occur via Ca2+ influx in the adult mouse hippocampus. Glia 64, 2093–2103. doi: 10.1002/glia.23042
Sajedinia, Z., and Hélie, S. (2018). A new computational model for astrocytes and their role in biologically realistic neural networks. Comput. Intell. Neurosci. 2018, 1–10. doi: 10.1155/2018/3689487
Sasaki, T., Ishikawa, T., Abe, R., Nakayama, R., Asada, A., Matsuki, N., et al. (2014). Astrocyte calcium signalling orchestrates neuronal synchronization in organotypic hippocampal slices. J. Physiol. 592, 2771–2783. doi: 10.1113/jphysiol.2014.272864
Sasaki, T., Kuga, N., Namiki, S., Matsuki, N., and Ikegaya, Y. (2011). Locally synchronized astrocytes. Cereb. Cortex 21, 1889–1900. doi: 10.1093/cercor/bhq256
Savtchenko, L. P., Bard, L., Jensen, T. P., Reynolds, J. P., Kraev, I., Medvedev, N., et al. (2018). Disentangling astroglial physiology with a realistic cell model in silico. Nat. Commun. 9:3554. doi: 10.1038/s41467-018-05896-w
Savtchouk, I., and Volterra, A. (2018). Gliotransmission: beyond black-and-white. J. Neurosci. 38, 14–25. doi: 10.1523/JNEUROSCI.0017-17.2017
Semyachkina-Glushkovskaya, O., Postnov, D., and Kurths, J. (2018). Blood-brain barrier, lymphatic clearance, and recovery: Ariadne's thread in labyrinths of hypotheses. Int. J. Mol. Sci. 19:3818. doi: 10.3390/ijms19123818
Semyanov, A., Henneberger, C., and Agarwal, A. (2020). Making sense of astrocytic calcium signals - from acquisition to interpretation. Nat. Rev. Neurosci. 21, 551–564. doi: 10.1038/s41583-020-0361-8
Serrano, A., Haddjeri, N., Lacaille, J. C., and Robitaille, R. (2006). Gabaergic network activation of glial cells underlies hippocampal heterosynaptic depression. J. Neurosci. 26, 5370–5382. doi: 10.1523/JNEUROSCI.5255-05.2006
Shigetomi, E., Tong, X., Kwan, K. Y., Corey, D. P., and Khakh, B. S. (2011). Trpa1 channels regulate astrocyte resting calcium and inhibitory synapse efficacy through gat-3. Nat. Neurosci. 15, 70–80. doi: 10.1038/nn.3000
Showalter, K., and Epstein, I. R. (2015). From chemical systems to systems chemistry: patterns in space and time. Chaos 25:097613. doi: 10.1063/1.4918601
Shuai, J. W., and Jung, P. (2002). Stochastic properties of ca(2+) release of inositol 1,4,5-trisphosphate receptor clusters. Biophys. J. 83, 87–97. doi: 10.1016/S0006-3495(02)75151-5
Shuai, J. W., and Jung, P. (2003). Selection of intracellular calcium patterns in a model with clustered ca2+ release channels. Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 67(3 Pt 1):031905. doi: 10.1103/PhysRevE.67.031905
Simard, P. Y., Steinkraus, D., and Platt, J. C. (2003). “Best practices for convolutional neural networks applied to visual document analysis,” in International Conference on Document Analysis and Recognition (Washington, DC: IEEE Computer Society), 958–962.
Sims, R. E., Butcher, J. B., Parri, H. R., and Glazewski, S. (2015). Astrocyte and neuronal plasticity in the somatosensory system. Neural Plast 2015:732014. doi: 10.1155/2015/732014
Smith, I. F., Wiltgen, S. M., Shuai, J., and Parker, I. (2009). Ca2+ puffs originate from preestablished stable clusters of inositol trisphosphate receptors. Sci. Signal. 2:ra77. doi: 10.1126/scisignal.2000466
Sneyd, J., Charles, A. C., and Sanderson, M. J. (1994). A model for the propagation of intercellular calcium waves. Am. J. Physiol. 266(1 Pt 1):C293–C302. doi: 10.1152/ajpcell.1994.266.1.C293
Softky, W. R., and Koch, C. (1993). The highly irregular firing of cortical cells is inconsistent with temporal integration of random EPSPs. J. Neurosci. 13, 334–350. doi: 10.1523/JNEUROSCI.13-01-00334.1993
Stobart, J. L., Ferrari, K. D., Barrett, M., Gluck, C., Stobart, M. J., Zuend, M., et al. (2018). Cortical circuit activity evokes rapid astrocyte calcium signals on a similar timescale to neurons. Neuron 98, 726.e4–735.e4. doi: 10.1016/j.neuron.2018.03.050
Sykova, E., and Nicholson, C. (2008). Diffusion in brain extracellular space. Physiol. Rev. 88, 1277–1340. doi: 10.1152/physrev.00027.2007
Taufiq-Ur-Rahman Skupin, A., Falcke, M., and Taylor, C. W. (2009). Clustering of InsP3 receptors by InsP3 retunes their regulation by InsP3 and Ca2+. Nature 458, 655–659. doi: 10.1038/nature07763
Tewari, S. G., and Majumdar, K. K. (2012). A mathematical model of the tripartite synapse: astrocyte-induced synaptic plasticity. J. Biol. Phys. 38, 465–496. doi: 10.1007/s10867-012-9267-7
Tong, X., Shigetomi, E., Looger, L. L., and Khakh, B. S. (2013). Genetically encoded calcium indicators and astrocyte calcium microdomains. Neuroscientist 19, 274–291. doi: 10.1177/1073858412468794
Torbensen, K., Rossi, F., Ristori, S., and Abou-Hassan, A. (2017). Chemical communication and dynamics of droplet emulsions in networks of belousov-zhabotinsky micro-oscillators produced by microfluidics. Lab Chip 17, 1179–1189. doi: 10.1039/C6LC01583B
Ullah, G., Jung, P., and Cornell-Bell, A. (2006). Anti-phase calcium oscillations in astrocytes via inositol (1, 4, 5)-trisphosphate regeneration. Cell Calcium 39, 197–208. doi: 10.1016/j.ceca.2005.10.009
Vanag, V. K. (2019). Hierarchical network of pulse coupled chemical oscillators with adaptive behavior: chemical neurocomputer. Chaos 29:083104. doi: 10.1063/1.5099979
Vanag, V. K. (2020). Size- and position-dependent bifurcations of chemical microoscillators in confined geometries. Chaos 30:013112. doi: 10.1063/1.5126404
Verkhratsky, A., and Nedergaard, M. (2018). Physiology of astroglia. Physiol. Rev. 98, 239–389. doi: 10.1152/physrev.00042.2016
Verkhratsky, A., Rodriguez, J. J., and Parpura, V. (2012). Calcium signalling in astroglia. Mol. Cell Endocrinol. 353, 45–56. doi: 10.1016/j.mce.2011.08.039
Volterra, A., Liaudet, N., and Savtchouk, I. (2014). Astrocyte Ca2+ signalling: an unexpected complexity. Nat. Rev. Neurosci. 15, 327–335. doi: 10.1038/nrn3725
Wallach, G., Lallouette, J., Herzog, N., De Pitta, M., Ben Jacob, E., Berry, H., et al. (2014). Glutamate mediated astrocytic filtering of neuronal activity. PLoS Comput. Biol. 10:e1003964. doi: 10.1371/journal.pcbi.1003964
Wang, Y., DelRosso, N. V., Vaidyanathan, T. V., Cahill, M. K., Reitman, M. E., Pittolo, S., et al. (2019). Accurate quantification of astrocyte and neurotransmitter fluorescence dynamics for single-cell and population-level physiology. Nat. Neurosci. 22, 1936–1944. doi: 10.1038/s41593-019-0492-2
Wolff, J. R., Stuke, K., Missler, M., Tytko, H., Schwarz, P., Rohlmann, A., et al. (1998). Autocellular coupling by gap junctions in cultured astrocytes: a new view on cellular autoregulation during process formation. Glia 24, 121–140. doi: 10.1002/(SICI)1098-1136(199809)24:1<121::AID-GLIA12>3.0.CO;2-T
Wu, Y.-W., Gordleeva, S., Tang, X., Shih, P.-Y., Dembitskaya, Y., and Semyanov, A. (2018). Morphological profile determines the frequency of spontaneous calcium events in astrocytic processes. bioRxiv. 67, 246–262. doi: 10.1101/410076
Wu, Y. W., Tang, X., Arizono, M., Bannai, H., Shih, P. Y., Dembitskaya, Y., et al. (2014). Spatiotemporal calcium dynamics in single astrocytes and its modulation by neuronal activity. Cell Calcium 55, 119–129. doi: 10.1016/j.ceca.2013.12.006
Xie, L., Kang, H., Xu, Q., Chen, M. J., Liao, Y., Thiyagarajan, M., et al. (2013). Sleep drives metabolite clearance from the adult brain. Science 342, 373–377. doi: 10.1126/science.1241224
Keywords: calcium signaling, cell morphology, noise-driven dynamics, astrocytes, modeling
Citation: Verisokin AY, Verveyko DV, Postnov DE and Brazhe AR (2021) Modeling of Astrocyte Networks: Toward Realistic Topology and Dynamics. Front. Cell. Neurosci. 15:645068. doi: 10.3389/fncel.2021.645068
Received: 22 December 2020; Accepted: 09 February 2021;
Published: 05 March 2021.
Edited by:
Leonid Savtchenko, University College London, United KingdomCopyright © 2021 Verisokin, Verveyko, Postnov and Brazhe. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Alexey R. Brazhe, brazhe@biophys.msu.ru