Until a few decades ago the quest to better understand high level brain functions, such as learning, memory, and cognition, mainly focused on investigating the rapid, spike-based information processing performed by neurons. Glial cells, and among them astrocytes, were largely regarded as passive support cells, providing neurons with nutrition and structural support without directly participating in information processing functions (Kandel et al., 1991). Over the past 20 years, however, a growing body of evidence has demonstrated that astrocytes do participate in bi-directional signaling with neurons and, therefore, possibly play an important role in shaping communication in the brain (Volterra and Steinhauser, 2004; Perea and Araque, 2005). These findings demand a revision of the traditional neuron-centric model used to explain higher order brain functions to include astrocytes as part of a neuron-glia network model. Within this new framework, signaling includes both, fast spike-based processing and slower modulation mediated by astrocytic elements (Nedergaard and Verkhratsky, 2012).

In this paper, we present a computational modeling framework that spans across several hierarchical layers of the central nervous system (CNS), from the molecular to the synaptic, dendritic, and neuronal levels. The structure of this framework allows us to investigate the influence of glial cells at each of these levels and how they can modulate neuronal communication. Traditionally, high level brain processes are explained using the framework illustrated in red in Figure 1A where molecular dynamics at neuronal synapses are linked to system-level brain functions via synaptic neuronal signaling (Kandel et al., 1991). The new methodology we have chosen to explain and model these complex mechanisms is represented by the combination of the red, neural components, and the green, glial/astrocytic components, which have been added to account for the contribution and influence of astrocytes. Since glutamate is the most important neurotransmitter involved at excitatory synapses and in astrocytic signaling, we focus on modeling the effect of astrocytic glutamate transporters on neuronal spiking within this new framework (Figure 1A).

In the original, neuron-centric framework shown in red, synapses, and neurons are the fundamental building blocks of a network within which information exchange is mediated by molecular mechanisms. The lowest level of hierarchy in Figure 1A comprises elements that interact at the molecular level, such as the three classes of glutamate receptors, ionotropic receptor-channels, AMPAR (α-amino-3-hydroxy-5-methyl-4-isoxazolepropionic acid receptor) and NMDAR (N-Methyl-d-aspartate receptor), and mGluR (metabotropic glutamate receptor). At the synaptic level, ion fluxes through these receptor-channels give rise to synaptic responses of varying time courses and amplitudes depending on channel kinetics. These responses sum in a non-linear spatio-temporal manner and evoke spiking activity in neurons. A network of these neurons constitutes the system that emulates a physiological function.

The box outlined in green shows the inclusion of the previously underappreciated glial components (Nedergaard et al., 2003). One particular type of glial cell, the astrocyte, has thousands of processes and extensions that are often found in close proximity to hundreds of synapses, which it can potentially modulate (Halassa et al., 2007b).

Additionally, it has been discovered that astrocytes express large amounts of neurotransmitter receptors and transporters (Wang and Bordey, 2008). Since a neuron’s pre- and postsynaptic nerve terminals are often ensheathed by astrocytes, which participate in synaptic signal processing, the term tripartite synapse was coined (Araque et al., 1999).

Within a tripartite synapse, biochemical, and morphological studies suggest that excitatory amino-acid transporters (EAATs) expressed on astrocytes are of the type EAAT2. These glutamate transporters maintain low extracellular glutamate concentration, which prevents neurotoxicity in spinal cord, striatum, and hippocampus (Rothstein et al., 1996), and might play a functional role in regulating synaptic currents by clearing glutamate after its synaptic release (Bergles and Jahr, 1998; Diamond, 2005). Neuronal transporters (EAAT3) also take up glutamate from the extracellular space, however, at a significantly lower rate than astrocytes due to their lower expression levels (Rothstein et al., 1996). The role of neuronal transporters could be to limit glutamate spill-over and to slow down glutamate clearance by glial transporters (Diamond, 2001; Scimemi et al., 2009). The importance of astrocytes in the regulation of glutamate uptake, the transformation of glutamate to glutamine (Uwechue et al., 2012) for re-usage in synaptic transmission, and epilepsy pathogenesis has recently been reviewed (Coulter and Eid, 2012).

Several studies demonstrated that astrocytes not only uptake glutamate inside a tripartite synapse, but under certain conditions, can also release glutamate (Araque et al., 1998), through a process termed gliotransmission (Halassa et al., 2007a). This term describes the process of glutamate release from astrocytes due to an increase in intracellular calcium via mGluR-mediated mechanisms in response to neural signaling inside a tripartite synapse (Parpura et al., 1994; Fiacco and McCarthy, 2006; Zur Nieden and Deitmer, 2006).

Because of the astrocyte’ structural characteristics and biochemical signaling mechanisms, these cells may play an important role at the neuron and network levels of the CNS hierarchy (Volterra et al., 2002). However, the exact manner in which astrocytes communicate with neurons in vivo is still unclear. For example, the effect of glutamate uptake in synaptic transmission, and on neuronal spiking, is difficult to verify using state-of-the-art experimental procedures. Furthermore, several experimental and review articles have been published that challenge the tripartite synapse concept (Agulhon et al., 2010; Nedergaard and Verkhratsky, 2012), and raised a number of issue regarding the way astrocytes contribute to synaptic signaling by performing feedforward and/or feedback action through the uptake and release of neuro- and gliotransmitters (Smith, 2010).

One of the primary reasons why these disagreements cannot easily be resolved, is the lack of experimental techniques to directly study astrocytes and their biochemical signaling mechanisms (Nedergaard and Verkhratsky, 2012). Since astrocytes are generally not electrically excitable (do not generate action potentials), techniques for measuring glial cell activity mainly rely on imaging for in vitro studies or genetic manipulations for in vivo experiments. Hence, experiments are often performed under non-physiological conditions (Nedergaard and Verkhratsky, 2012), as in vitro experiments do not represent the natural environment of these cells and their responses might be drastically amplified. Disabling astrocytes through genetic manipulations is equally non-physiological and might result in compensatory effects in vivo (Smith, 2010).

One promising approach that can resolve these uncertainties is to use computational models. Such models provide opportunities to analyze the behavior of the system in response to various stimuli and conditions otherwise difficult to conduct experimentally, either due to the lack of the necessary technology, or difficulties in accessing readouts at the sub-cellular level.

In recent years, many computational modeling groups have demonstrated the glial influence within synaptic, neuronal, or network dynamics. Nadkarni and Jung (2007) have characterized astrocytic effects on spontaneous activity at the postsynaptic level. Silchenko and Tass (2008) have modeled the effects of glutamate release from astrocytes on neuronal depolarization and activity. Somjen et al. (2008) have used computer simulations to demonstrate that neuron-glia interactions are in part mediated by potassium ion fluxes between the two entities. Finally De Pitta et al. (2011) and Wade et al. (2011) have presented a modeling approach with bi-directional communication and synchrony between astrocytes and neuron clusters.

Our modeling approach integrates the dynamics from the molecular to the neuronal level and provides an original framework that allows for a better understanding of the effects of glia in a hierarchical manner. This modeling methodology mimics the structure of the CNS which spans several spatial and temporal scales and forms a hierarchical system from the bottom-up: molecular to synapse level, synapse to neuron level, and neuron to network level, with the inclusion of glial interactions. For example, the activation of receptors at single synapses, which are distributed along the dendritic tree of pyramidal neurons, sum in a non-linear fashion, which changes the membrane potential (Poirazi et al., 2003). These responses influence the network level, along with additional modulatory inputs from other pathways and inhibitory connections.

In this paper, we focus on the role of astrocytic glutamate uptake on synaptic responses and how it can modify neuronal spiking within the context of this modeling framework (Figures 1A,B). We observe that glial glutamate uptake is important to decrease desensitization of ionotropic glutamate receptors resulting in enhanced paired pulse facilitation for very short input intervals. Based on these observations, we propose that the time interval between input pulses and the interactions between receptors and transporters significantly contribute to neuron spiking patterns. Our model provides a unique method to assess critical parameters that can influence neural network behavior in relation to glutamate uptake.

The synaptic modeling platform EONS/Rhenoms™ (Elementary Objects of the Nervous System) was developed to configure detailed synapses and define the distribution and arrangement of molecular elements within a synaptic environment (Bouteiller et al., 2008). The models and their descriptions are available online at http://synapticmodeling.com/. Each of the individual synaptic elements are described by kinetic schema that characterize the behavior of receptors using parameters fitting the model’s responses to experimental data. The generic synapse model described here comprises presynaptic calcium buffers, voltage-dependent calcium channels, a single vesicular glutamate release site, glutamate diffusion in the synaptic cleft (Savtchenko and Rusakov, 2007), and binding of glutamate to postsynaptic ionotropic and metabotropic glutamate receptors (Ambert et al., 2010; Greget et al., 2011). Glutamate uptake is incorporated in the platform by including a model of a high-affinity glial glutamate transporter (EAAT2/GLT-1) adapted from Bergles et al. (2002), and a neuronal transporter (EAAT3) with parameters adapted from Larsson et al. (2004). Glutamate molecules released from presynaptic vesicles diffuse in the synaptic cleft. Glutamate molecules taken-up by these transporters are subtracted from those available at ionotropic and metabotropic receptors in the postsynaptic density. Numerous studies have quantified the tortuosity of the synaptic cleft and estimated the glutamate diffusion coefficient (Nielsen et al., 2004). The diffusion model used here has been adapted from Savtchenko and Rusakov (2007) to calculate glutamate concentration inside the cleft as a function of the distance of the receptor from the release site using Eq. 1. The glutamate diffusion coefficient is 0.4 μm² ms⁻¹. The total number of transmitter molecules released from the vesicle in the model is 3,000 per release event. The concentration of glutamate inside the cleft is determined by the following equation:

Where Glu, r, and D represent the concentration of glutamate inside the cleft, the radial distance, and the diffusion coefficient, respectively. Q represents the number of glutamate molecules released instantaneously. “δ” stands for the height of the cleft and is maintained constant throughout the simulations at 20 nm. The glutamate profile seen by postsynaptic receptors is shown in Figure 2B.

The concentration levels of glutamate available at AMPAR and NMDAR is determined by their spatial location on the post synaptic density (PSD). The astrocytic glutamate uptake is calculated by the following rate equation for glutamate flux:

where Na2ToH, Na2ToGH, Na2To, Na2ToG are the intermediate states that determine glutamate bound and transported by the transporter (see Figure A3 in Appendix). The values for parameters k₆, k₃ is 6 mM⁻¹ ms⁻¹ and for k₋₆, k₋₃, is 0.5 ms⁻¹. All the other rate constants and ion concentrations are the same as reported in Bergles et al. (2002). The uptake rate in the equation above is for a single transporter channel. The diameter of the postsynaptic disk is 200 nm (Takumi et al., 1999) the height of the synaptic cleft is fixed at 20 nm (Savtchenko and Rusakov, 2007). The astrocyte membrane and its transporters are set at a distance of 400 nm from the release site. The distance of the transporter from the release site does not influence its uptake rate (Diamond, 2005). Transporters are expressed on astrocytes surrounding synapses with densities within a range of 6,500–13,000 μm⁻² (Lehre and Danbolt, 1998; Diamond, 2005). The astrocyte is modeled as a cylindrical surface of height 20 nm and radius 400 nm as depicted by the green ensheathment around the synapse in Figure 2A. Assuming the maximum transporter density (13,000 μm⁻²) and a 50% astrocytic coverage surrounding the synapse, we calculate the number of glutamate transporters per synapse to be 650. For the same conditions but with 50% density, we calculate the number of transporters to be 325. Glutamate uptake from neuronal transporters is calculated in a similar manner to that by glial transporters (Eq. 2) with the kinetic rate constants adapted from Larsson et al. (2004) and a transporter density of 90 μm⁻² (Holmseth et al., 2012). The total glutamate cleared from the receptor vicinity is obtained by multiplying the cleared glutamate rate dGluO/dt, integrating over the time steps of the simulation and multiplying by the number of transporters. This amount of cleared glutamate is then subtracted from the glutamate input concentration available at the receptors.

We used the AMPA receptor model described in detail in Robert and Howe (2003) which represents a 16 states model describing the receptor transitions between resting, desensitized, and conducting open states. Successive binding of two, three, and four glutamate molecules produces conformational changes leading to fast opening and closing of the channel.

The current through the channel is calculated by:

Where the open conducting states evolve as:

where I_AMPA is the current mediated by AMPA receptors, nb_AMPA is the number of AMPA receptors (in this study nb_AMPA is 80), consistent with reported AMPAR numbers between 46 and 147 at CA1 hippocampal synapses (Matsuzaki et al., 2001). g₂, g₃, g₄ are unitary conductances with values 9, 15, and 21 pS associated with the channel in open states when 2, 3, and 4 glutamate molecules are bound respectively. The probabilities for the O₂, O₃, O₄ states are calculated based on ODEs simulated using solvers in SBML™. The derivatives of open states Ȯi (where i = 2, 3, 4) are calculated as a product of matrix M containing other states transition rate constants with input Glu and vector of currents states O_i. V_rev is the reversal potential of the AMPAR (V_rev is 0 mV) and V is the membrane potential that changes dynamically during the simulation. More details on the model can be obtained from Robert and Howe (2003) and for kinetic rate parameters, please, see Section “Appendix.”

Glutamate concentration available at NMDAR after glutamate uptake is provided as input to the NMDAR model represented in a 15-state kinetic scheme, which includes agonist (glutamate) and co-agonist (glycine) binding sites, channel blockers (memantine and magnesium), as well as several antagonist sites. The kinetics of this model are borrowed from Ambert et al. (2010). For validation of the NMDAR model, various protocols were tested. For a single short pulse of glutamate, experimental results reported by Schorge et al. (2005) were used to validate the model. For long or repetitive glutamate inputs to the model, the kinetic parameters were adjusted to properly capture effects of desensitization and to match experimental data from Zhang et al. (2008). The equations to calculate NMDAR-mediated synaptic current are:

where I_NMDA is the current mediated by NMDA receptors. nb_NMDA is set at 20 for this study consistent with observations made in Takumi et al. (1999) and Racca et al. (2000). I_o is the current associated with the open conducting state O(t) calculated using ODEs solved with kinetics described in Ambert et al. (2010). The magnesium concentration in the external solution is set to 1 mM; Ψm is the electrical distance of the magnesium binding site from the outside of the membrane (set at 0.8); R, the molar gas constant (8.31434 J mol^-1 K^-1); F, the Faraday constant (9.64867.104 C mol^-1); T, the absolute temperature (273.15 K); g₁ and g₂ are the conductances associated with the open states when one or two glutamate molecules are bound and are 40 and 247 pS respectively; α = 0.01 is the steepness of the voltage-dependent transition from g₁ to g₂.

The total synaptic current is calculated from the sum of AMPAR and NMDAR currents.

The synaptic currents calculated using Eq. (6) drive the neuron membrane potential and synapses are now approximately acting as current sources (Jaffe and Carnevale, 1999). A CA1 pyramidal cell model with active sodium and potassium channels and a morphology described in Jarsky et al. (2005) was used within the NEURON simulation environment (Hines and Carnevale, 1997). Synapses were placed at 16 random locations in stratum radiatum (middle one-third of the cell, 100–200 μm). The synaptic strength was tuned by a factor of 6 to reach threshold levels for neuronal spiking such that spiking probability was 1.

In this paper we focused to study the role of glutamate uptake mediated by glutamate transporters present on astrocytic membranes surrounding CA1 hippocampal synapses on synaptic transmission and neuronal spiking. First we showed that the model was able to reproduce previously observed phenomena regarding the role of glutamate uptake on synaptic transmission by modulating responses mediated by AMPA and NMDA receptors. This was achieved by removing glutamate molecules that are bound to or transported by glutamate transporters from those available at the levels of the receptors. Previous studies have emphasized the importance of taking into account both glutamate diffusion and binding to transporters to determine changes in the decay kinetics of synaptic glutamate concentration (Wadiche et al., 1995). In agreement with this study, glutamate is removed from synaptic cleft by diffusion, and binding to and transport by glutamate transporters.

Some of the phenomena demonstrated here under the conditions we used are summarized:

We also underscore the diversity of synapses that arises from the variability in spatial location of the receptors, which we are investigating in a separate study. It is interesting to note that most of the receptors, such as NMDARs and mGluRs, have many modulatory sites, which consequently increase the number of parameters needed in the simulation. Our highly configurable geometric synapse model allows for the exploration of various parameters that influence sub-cellular/molecular level interactions and their direct or indirect influence on the synapse and neuron levels. By linking such complex unified model of a synapse to morphologically realistic models of neuron within the NEURON simulation environment, we can investigate this complexity in an orderly and hierarchical fashion. This modeling effort allows for the investigation of key phenomena that are otherwise difficult to explore through mainstream reductionist modeling approaches. The main technical drawback of this kind of approach is the computational overhead involved. The modeling paradigm itself is complex in its nature due to (i) the level of parametric details and (ii) the time scale of processing of some elementary models, which in some cases takes place within tenths of microseconds, thereby slowing down the entire system.

The astrocyte model presented here is not a complete model, and astrocytes are known for their role in influencing synaptic transmission beyond glutamate uptake and clearance. This model needs to be expanded to incorporate other important features, such as direct neurotransmitter release from astrocytes and signaling to neurons. All synapses are wrapped differently by astrocyte processes, covering smaller or larger areas with different levels of transporter expression.

We have assumed for simplicity a constant wrapping for each individual synapses, which may not be necessarily true in vivo. Parametric models, where the model behavior is explained by a set of parameters are in general limited by the scope of available experimental evidence. However, with a parametric modeling paradigm, we can test the reliability and sensitivity of these parameters. The observations described through our simulation studies are still preliminary and the modeling architecture established here will enable us to further investigate the effects caused by changing the amount of ensheathment around the synapses, as the density of EAATs in both glia and neurons appears to play an influential role in shaping synaptic glutamate concentration profile and its functional consequences. Including other details of a tortuous path for glutamate and including extrasynaptic NMDARs, may also affect synaptic responses. Future modeling efforts will be directed toward investigating the hierarchical effects of astrocytes on sub-populations of neurons and synapses contacted by astrocyte processes, by incorporating geometry-related considerations. Here, we demonstrated the hierarchical link between synaptic currents to spike generation while also taking into account astrocytic glutamate uptake effects on molecular elements. This is a novel approach and one of the few times such a link has been shown. These results may have significant implications for understanding glial cell effects on nerve cell membrane potential and thus, nerve cell spiking, i.e., neuronal information flow. These results also are important because they strongly suggest that glial cell uptake of synaptic glutamate during neuron-to-neuron synaptic transmission should influence spike-dependent processes that are relevant to secondary messenger pathways and other long-term effects. This comprehensive framework will allow investigating complex mechanisms within large neuron/glia networks, including neurodegenerative diseases and their underlying processes. We can identify the most sensitive parameters in neuron glial interactions and develop different testing paradigms to understand the molecular basis of diseases associated with astrocytic dysfunction.

Jean-Marie C. Bouteiller, Michel Baudry and Theodore W. Berger have a conflict of interest. The University of Southern California holds an equity interest in Rhenovia Pharma and also has received licensing income from Rhenovia Pharma.