HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS		QUICK SEARCH:   [advanced] Author: Keyword(s): Year:  Vol:  Page:

This Article Abstract Full Text (PDF) Supplemental Figure All Versions of this Article: 97/6/4079    most recent 01357.2006v1 Alert me when this article is cited Alert me if a correction is posted Citation Map Services Email this article to a friend Similar articles in this journal Similar articles in ISI Web of Science Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by Sussillo, D. Articles by Maass, W. Search for Related Content PubMed PubMed Citation Articles by Sussillo, D. Articles by Maass, W.

Self-Tuning of Neural Circuits Through Short-Term Synaptic Plasticity

¹Center for Theoretical Neuroscience in the Center for Neurobiology and Behavior, Columbia University, New York, New York; and ²Institute for Theoretical Computer Science, Technische Universitaet Graz, Graz, Austria

Submitted 24 December 2006; accepted in final form 28 March 2007

	ABSTRACT

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
Numerous experimental data show that cortical networks of neurons are not silent in the absence of external inputs, but rather maintain a low spontaneous firing activity. This aspect of cortical networks is likely to be important for their computational function, but is hard to reproduce in models of cortical circuits of neurons because the low-activity regime is inherently unstable. Here we show—through theoretical analysis and extensive computer simulations—that short-term synaptic plasticity endows models of cortical circuits with a remarkable stability in the low-activity regime. This short-term plasticity works as a homeostatic mechanism that stabilizes the overall activity level in spite of drastic changes in external inputs and internal circuit properties, while preserving reliable transient responses to signals. The contribution of synaptic dynamics to this stability can be predicted on the basis of general principles from control theory.

	INTRODUCTION

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
Cortical neurons fire not only in response to sensory inputs, but also spontaneously, as part of an ongoing network activity of the neocortex. According to Steriade (2001), the rates of spontaneous firing recorded extracellularly from motor and association areas of awake monkeys and cats are typically in the range of 10–15 Hz. Baddeley et al. (1997) found that cells in area IT (recorded extracellularly in awake macaque) fire at a rate of 14 ± 8.3 Hz during blank-screen viewing. Intracellular recordings from regular-spiking neurons (presumed to be pyramidal cells) in motor, association, primary somatosensory, and visual association areas of awake adult cats yielded a spontaneous firing rate of 9.4 ± 1.7 Hz (Steriade et al. 2001). The functional role of this spontaneous firing of neocortical neurons is unknown. It was previously conjectured that persistent neuronal activity, which is "conspicuously absent in cerebellar and basal ganglia circuits," is an essential feature of the neocortex, which enables it "to incorporate the past into the system's present state" (Buzsáki 2006). It was also previously conjectured to support "flexible cooperation among local and distant cell assembles," which is "believed to underlie the efficacy of cortical performance and is an essential ingredient for cognition" (Buzsáki 2006). On a more technical level, it was pointed out that not only recurrent circuits of neurons in the neocortex, but dynamical systems in general support higher computational performance if they operate not in the ordered (or dissipative) regime, but closer to a critical state—which amounts to a substantial spontaneous firing activity in a neural circuit—so that also small external inputs can cause significant system responses (Buzsáki 2006; Legenstein and Maass 2007).

Spontaneous brain activity is commonly assumed to be generated by a combination of intrinsic electrophysiological properties of single neurons and synaptic interactions in the network (for a recent review see Destexhe and Contreras 2006). It was already pointed out in Abeles (1991) that "it seems to be very difficult to stabilize a network of inhibitory and excitatory neurons at very low levels of activity." Through theoretical analysis and extensive computer simulations we show that the short-term dynamics of synapses—more precisely, the empirically found diversity of mixtures of paired-pulse depression and paired-pulse facilitation that depend on the type (excitatory or inhibitory) of pre- and postsynaptic neurons—endows networks of leaky integrate-and-fire (LIF) neuron models with either conductance-based or current-based synapses, with a remarkable stability of spontaneous network activity, even for rather low average firing rates, such as 10 Hz.

It is well known that synapses in the neocortex do not respond in the same way to each spike in a train of spikes. Rather, the amplitude of the postsynaptic response [excitatory postsynaptic potential (EPSP) or excitatory postsynaptic current (EPSC)] to a presynaptic spike depends on the history of presynaptic firing activity. If the response amplitude for the second spike in a pair of spikes is smaller than that for the first one, one calls this effect paired-pulse depression. If the response amplitude for the second spike is larger, one calls this paired-pulse facilitation. Actually, however, most neocortical synapses exhibit a mixture of depression and facilitation and the amplitude of the postsynaptic response depends in a complex manner on the temporal pattern of several preceding presynaptic spikes. This effect is commonly referred to as synaptic dynamics or short-term plasticity. In contrast to long-term synaptic plasticity this short-term plasticity does not usually depend on the pattern of postsynaptic firing and its effect is not long-lasting (more precisely: it depends only on the temporal pattern of presynaptic spikes during the last few seconds). We refer to Abbott and Nelson (2003), Thomson (2003), and Abbott and Regehr (2004) for recent reviews of biological data and proposed functional roles of short-time synaptic plasticity. The dynamics of cortical synapses can be fitted quite well by the phenomenological model from Markram et al. (1998), where the short-term dynamics of synapses is characterized by three parameters: U (which roughly models the release probability of a synaptic vesicle for the first spike in a train of spikes), D (time constant for recovery from depression), and F (time constant for recovery from facilitation). Hereafter we call this set of parameters UDF. An essentially equivalent model was previously proposed in Varela et al. (1997).

In contrast to the static synapse of an artificial neural network, a dynamic synapse does not just provide a fixed amplification of presynaptic activity, but rather implements a nonlinear filter that could potentially serve a number of different purposes, such as cortical gain control (Abbott et al. 1997), burst detection (Lisman 1997), or temporal integration of information (Haeusler and Maass 2007; Maass et al. 2004). In Tsodyks et al. (1998) it was shown that depressing synapses can create in a mean-field model for a population of excitatory neurons a stable fixed point for its activity and that dynamic synapses can create a variety of rhythmic and irregular activity patterns if one adds a population of inhibitory neurons.

One intriguing aspect of biolocial data on synaptic dynamics is that this short-term dynamics—more precisely, the parameters UDF—differ from synapse to synapse. Moreover, the UDF values form clusters within the three-dimensional (3D) UDF parameter space and the cluster to which the UDF value of a particular synapse belongs depends on the type (excitatory or inhibitory) of its presynaptic neuron and in some cases also on the type of the postsynaptic neuron (Galarreta and Hestrin 1998; Gupta et al. 2000; Markram et al. 1998; Thomson 1997). It has been difficult to identify a possible functional role for the rules that assign a particular region in the 3D UDF parameter space to different types of synaptic connections. Herein we show that the assignment of different ranges of UDF values to different types of synapses makes it possible to implement a self-tuning principle for the firing rate of a cortical circuit. This self-tuning principle enables neural circuits to respond to external perturbations with a characteristic transient response in the firing rate of excitatory neurons and then returns to its previous firing rate within a few hundreds of milliseconds back into a given target range. The firing rate of inhibitory neurons is automatically adjusted by the synaptic dynamics so that it compensates the external perturbation. A similar self-tuning property can be demonstrated for changes within the circuit, such as those caused by long-term synaptic plasticity or changes in the concentration of neurotransmitters.

We show that this self-tuning principle in models of cortical networks of neurons is related to an abstract self-tuning principle recently proposed for much simpler dynamical systems in control theory (Moreau and Sontag 2003a,b). Their abstract self-tuning principle enables some simple two-dimensional (2D) dynamical system to return in a reliable manner to a particular operating regime (e.g., the vicinity of a bifurcation), in spite of external perturbations. They considered the following dynamical system with external input v(t) that has just two variables: x (firing rate) and µ (synaptic strength)

(1)

(2)

The parameter µ₀ is here assumed to be unknown and conditions¹ on the functions f and g are given under which the variable x moves (in the absence of external input) to a steady state x* with = 0 (i.e., µ = µ₀) regardless of the value of this parameter µ₀. The system expressed by Eqs. 1 and 2 is obviously quite far removed from models for recurrent neural networks consisting of excitatory and inhibitory neurons with nonnegative firing rates. However, it is remarkable that its self-tuning property relies on the assumption that the strength µ of the positive feedback in Eq. 2 varies as a function of x, unlike the usual "static" synaptic weights in artificial neural networks. Furthermore, the required assumptions on the functions f and g in Eq. 2 imply that µ decreases when x grows beyond its steady-state value x*. Thus the resulting dynamics of the "synaptic weight" µ has some properties in common with synapses between excitatory neurons in the cortex (which have been found to be usually depressing), although apparently no relation to synaptic dynamics had been intended (or even discussed) in Moreau and Sontag (2003). Herein we show (see the end of RESULTS and the Supplementary material² that this control principle can be generalized and adapted to more complex dynamical systems that consist of two interacting modules (e.g., a population of excitatory neurons and a population of inhibitory neurons) that have profoundly different functional roles. It turns out that an application of the control principle to such a distributed system requires not only an activity-dependent modulation of "weights" µ(x), but different rules for the modulations µ_mn(x) of different weights, that depend on the type n of the source of the feedback µ_mn(x_n)·x_n and the type m of the target of the synaptic connection (where m, n range over excitatory and inhibitory neurons). Thus in applications to models for cortical circuits, this principle requires differential synaptic dynamics for different types of synaptic connections.

More precisely, our theoretical analysis predicts that the UDF parameters that characterize the dynamics of synapses belong to a specific region within the UDF space (which we call the "N-volume") if a strengthening of this type of synaptic connection tends to increase the firing rate of excitatory neurons in the circuit, or into another region ("P-volume") if a strengthening of this type of synaptic connection tends to decrease the firing rate of excitatory neurons. A comparison of these predicted ranges of UDF values with empirically found UDF values gives a quite satisfactory (although not perfect) match. This is remarkable because many aspects that are relevant for the dynamics of cortical circuits of neurons (different input–output behaviors of excitatory and inhibitory neurons), further differentiation of neuronal dynamics according to subclasses of inhibitory (Gupta et al. 2000) and excitatory neurons, subnetworks formed by particular subtypes of neurons (Yoshimura and Callaway 2005; Yoshimura et al. 2005), and modulation of neuronal dynamics through neuromodulators (Gulledge and Jaffe 1998, 2001) are not reflected in the network models consisting of LIF neurons that we examine herein.

Altogether this article provides a novel explanation for the functional role of short-term synaptic dynamics, in particular also for the empirically found diversity of synaptic dynamics in dependency on the type of pre- and postsynaptic neurons. In addition it contributes new ideas and principles for the investigation of a research topic, obviously of primary importance—the remarkable stability of cortical circuits of neurons against numerous perturbations and parameter drifts. This stability has so far not been reproduced with either artificial neural networks or with biologically more realistic models for cortical circuits and thus an understanding of its prerequisites may provide important insight into design principles of cortical networks. The general program to investigate fundamental aspects of the design of neural circuits by focusing on their stability properties has already produced remarkable insight into another type of neural circuit: the pattern-generating pyloric circuit in the crustacean stomatogastric ganglion (for a recent review see Marder and Goaillard 2006). It turns out that this circuit is endowed with intricate homeostatic mechanisms that move a high-dimensional vector of internal parameters (conductances, etc.) to a different setting if one method for producing the desired firing pattern becomes unavailable, such as because of changes in the concentration of neuromodulators. Thus it is essential for the functioning of this pattern-generating circuit that it can produce "multiple solutions to the production of similar behavior" (Marder and Goaillard 2006). We propose that analogous design principles are implemented in cortical circuits of neurons, which have to maintain a stable computational function in spite of drastic changes in the intensity of external inputs and in the concentrations of neuromodulators (that affect the excitability of neurons within the network; Gulledge and Jaffe 1998, 2001). More specifically, we propose that the short-term dynamics of synapses is essentially involved in design principles that produce a stable function of cortical networks.

	METHODS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
Neuron models

We used two leaky integrate-and-fire neuron models (LIF): one with current-based synapses and the other with conductance-based synapses. The current-based LIF neuron is defined by

(3)

where V is the membrane voltage; _mem is the membrane time constant; V_rest is the resting membrane voltage; R_m is the membrane resistance; J_E,j(t) and J_I,j(t) are the K_E and K_I synaptic currents from the excitatory and inhibitory synapses, respectively; and I^ext(t) is an external input composed of mean input current I^inject and a white noise current I^noise(t), with zero mean and SD ^noise. The LIF neuron with conductance-based synapses is defined by

(4)

where C_m is the membrane capacitance; g_leak is the leak conductance of the neuron; and g_E,j(t) and g_I,j(t) are the K_E and K_I synaptic conductances from the excitatory and inhibitory synapses, respectively. The constant E_E is the reversal potential for the excitatory synapses and the constant E_I is the reversal potential for the inhibitory synapses. The dynamics of each synapse over time are defined by

(5)

where j(t) is the instantaneous synaptic current or conductance, depending on the model, and W is the synaptic weight. The currents or conductances decrease exponentially with time constant _syn and increase instantaneously by adding W to the running value of j(t) whenever a spike occurs in the presynaptic neuron at time _sp. Herein we compare the use of static synapses to dynamic synapses and so the synaptic weight W was constant in some simulations and in other simulations allowed to vary according to known biological rules of synaptic plasticity [i.e. W(t)]. In the case of dynamic conductance-based synapses we varied the maximal conductance and in the case of dynamic current-based synapses we varied the maximal current. When static synapses were used (i.e., when W is constant) we use the notation J to denote the synaptic weight and always state whether the simulation used conductance- or current-based synapses. When dynamic synapses were used [i.e., when W(t) is implied] we use the notation µ to denote the synaptic weight.

We chose values for these LIF neurons that qualitatively reflect the high conductance, in vivo UP state measured in in vivo studies of cortical neurons (Destexhe et al. 2003). The membrane resistance was R_m = 10 M, the membrane time constant was _mem = 10 ms, and, in the case of neurons with conductance-based synapses, the capacitance was set to C_m = 10 ms/10 M = 1 nF. The neurons’ resting potential was V_rest = –60 mV and the firing threshold was V_thresh = –50 mV. When the neuron spiked, the membrane voltage was reset to V_reset = –60 mV, where it remained for a refractory period of _refr = 3 ms. The excitatory synaptic time constant _syn was _e = 4 ms; the inhibitory value was _i = 8 ms. In the case of neurons with conductance-based synapses g_leak = 1/R_m = 100 nS and the reversal potentials for the excitatory and inhibitory neurons were E_E = 0 mV and E_I = –80 mV. We modeled background synaptic inputs (from neurons outside of the simulated circuit) through I^ext = I^inject + I^noise with a constant I^inject = 2.455 x 10^–9 A and values of I^noise drawn from a 0-mean Gaussian with ^noise = 6 x 10^–9 A SD. I^ext caused an average membrane potential of –55.4 mV and SD of 4.3 mV in the absence of recurrent synaptic feedback (and when the threshold V_thresh was temporarily raised to 0 mV to prevent spikes). This resulted in a firing rate of approximately 20 Hz for each neuron in the absence of recurrent synaptic input I^rec and was chosen because it is the median value of the firing rates of irregular and tonic discharges reported in in vivo conditions in waking animals (Destexhe et al. 2003; Steriade et al. 2001) and is thus suitable for studying stability involving both up- and downregulation of firing rates.

Spiking network model

In our simulations we created sparse, random-spiking networks of 5,000 neurons, with 4,000 excitatory neurons and 1,000 inhibitory neurons. Neurons were either excitatory or inhibitory. The network models contained either purely LIF neurons with conductance-based synapses or purely LIF neurons with current-based synapses; we refer to the former as "conductance-based networks" and the latter as "current-based networks." We use the term sparse to denote that the number of excitatory connections K_e and inhibitory connections K_i to a given neuron are much less than both the total number of excitatory neurons N_e or inhibitory neurons N_i in the network. We use the term random to denote that the connections between neurons are chosen randomly and follow no other design principle. Following Vogels and Abbott (2005), the neurons were connected randomly to each other with a 2% connection probability. The network connectivity between neuron types is shown in Fig. 1.

View larger version (21K): [in this window] [in a new window]   FIG. 1. Basic network architecture used in mean-field simulations and spiking network simulations. Two populations of neurons connect to each other, one excitatory Exc with firing rate xe and one inhibitory Inh with firing rate xi. Blue denotes excitatory populations or recurrent excitatory connections and red denotes inhibitory populations or inhibitory recurrent connections. If the synaptic connections have a static synaptic weight they are denoted Jmn from presynaptic pool n to postsynaptic pool m; if the synaptic connections have dynamic synaptic weights then they are denoted as µmn. External excitatory input, labeled Iinject + Inoise, is denoted by green. External input in this paper is a constant noisy current with a stationary mean.

We used two equivalent forms of dynamic synapse models for networks with spiking and nonspiking neurons. The first set of expressions constituting the Eq. 6 system is taken from Markram et al. (1998)³ and describes for the case of spiking neurons the magnitude of the synaptic weight µ_k (a current or conductance, depending on the network model) for the kth spike in a spike train with interspike intervals ₁, ₂, ..., _k–1

(6)

where u is the running variable for synaptic utilization and R is the running variable for synaptic availability, with u, R [0, 1]. The constants U, D, and F (for curtness further denoted as UDF or AUDF when the strength scaling factor A is included) represent the release probability for the first spike, the depression time constant, and facilitation time constant, respectively. The synaptic weight for the kth spike is then defined by µ_k (i.e., W := µ_k in Eq. 5) and has an arbitrary strength scale set by A.

An equivalent, continuous mean-field model for dynamic synapses was developed in Tsodyks et al. (1998), which can be used in nonspiking mean-field models with continuous time. The mean-field equations for dynamic synapses are given by

(7)

where x(t) represents the instantaneous firing rate of the presynaptic neuron; U¹ is simply an algebraic term introduced for clarity; and µ, u, R, U, D, and F have the same meaning as in the Eq. 6 system. Note the magnitude µ(t) of the synaptic weight at time t depends not only on the current presynaptic firing rate x(t), but through the auxiliary variables u and R also on the history of this presynaptic firing rate.

The steady-state equations for this continuous dynamic synapse model yield in the case of a constant presynaptic firing rate x(t) = x*

(8)

where the * notation denotes steady-state values of u, U¹, R, and µ. In particular for m, n {e, i} the value µ*_mn (x_n) defined according to Eq. 8 denotes the steady-state synaptic weight from any neuron in population n to any neuron in population m as a function of the presynaptic firing rate x* = x_n.

Linking UDF parameters to network stability

We show that the fundamental property of dynamic synapses for creating a network with a stable firing rate is the relationship between the slopes of the steady-state synaptic weight curves µ*_mn (x_n) for different synapse types. We searched for UDF parameters that generate a negative derivative for the connections E E and I I because for our parameters an increase in the strength of either leads to a higher firing rate x_e of the excitatory neurons and a positive derivative for the steady-state curves for the connections E I and I E because an increase in their strengths leads to a lower firing rate x_e (see Supplemental materials for further explanation). More specifically, we searched for UDF settings that generated either a strictly positive or strictly negative derivative between the frequencies 10 and 100 Hz of the presynaptic firing rate and that additionally were in a physiological range (see Table 1). We uniformly sampled the UDF space in the range U (0, 1], D (0, 1] s, and F (0, 1] s (at increments of 0.014 in each dimension) and found two topologically simple volumes of parameter space that satisfied these requirements. These were designated "P-volume" for the positive derivative and the "N-volume" for the negative derivative. We randomly chose three assignments R1, R2, and R3 to the UDF values of the four synapse types to create three separate dynamic synapse UDF parameters sets. Thus each parameter set R_i (i = 1, 2, 3) consisted of two points from the N-volume, which were randomly assigned as the mean for the UDF-parameters of the synaptic connections E E and I I, and two points from the P-volume, which were also randomly assigned as the mean to the synaptic connections E I and I E. Each individual R_i parameter set was used separately from the others in the network simulations. Additionally, we used experimental data on UDF parameters from Gupta et al. (2000) and Markram et al. (1998) for all four synapse types (see Table 1). All self-tuning simulations were carried out for these parameters based on experimental data, as well as all three parameter assignments R1, R2, and R3 based on a theoretical analysis of desirable network performance. Performance differences between the R1–R3 parameter settings were small and we averaged the responses across the parameters sets. These UDF parameters were fixed for all self-tuning experiments and only the scaling factors A_mn of the dynamic synapses were modified on a per-simulation or per-network basis.

Self-tuning simulation setup

The self-tuning simulations we performed involved measuring and comparing the equilibrated (steady-state) average firing rates of excitatory neurons with static synapses and then with dynamic synapses for three network models: conductance-based networks, current-based networks, and a rate-based mean-field network (see following text for the definition). There were three types of perturbation simulations run, with two parameters perturbed in each case:

1) mean input current strength and SD of current injection

2) excitatory synaptic strength and inhibitory synaptic strength

3) percentage excitatory neuron inactivation (a temporary cessation of function potentially caused by a variety of circumstances such as neuromodulators) and percentage inhibitory neuron inactivation. This resulted in 841 networks for the input perturbation experiments and 900 networks for the synaptic weight perturbation and neuronal inactivation experiments. Each network was generated with random connections and first run in the control case with static synapses. The network was then regenerated with the exact same parameters, but a different random number seed, and run with dynamic synapses. To test the robustness of the ideas presented, each network was tested with four different dynamical synapse parameters sets, three based on theoretical considerations called R1, R2, and R3 and the parameters based on experimental data, thus requiring 841 + (841 x 4 x 3) = 10,933 separate network simulations for the input perturbation experiments and 900 + (900 x 4 x 3) = 11,700 separate network simulations for experiments with synapse strength and neuronal inactivation perturbations. To measure the firing rates we let the control networks with static synapses run for 1,500 ms and the average firing rates of excitatory and inhibitory neurons were measured for the last 1,000 ms. The networks with dynamic synapses ran for 2,000 ms and the steady-state firing rates were measured for the last 1,000 ms, which allowed extra time for the dynamic synapses to equilibrate. The measured firing rates for the theoretically derived UDF values were then averaged across all three dynamic synapse parameters sets R1, R2, and R3. Network simulations were performed using the CSIM neural simulator (Natschläger et al. 2003) and MatLab (The MathWorks, Natick, MA). The simulation time step was 0.1 ms. The mean-field simulations were also performed using MatLab.

Correlation function

Network synchrony plays a large role in the behavior of networks of sparsely connected, spiking neurons (Brunel 2000). We were interested in separating asynchronous irregular activity from other types of network activity. We used the voltage traces of randomly sampled neurons to analyze the network correlation and then applied the standard cross-covariance function. The covariance coefficient for the sequences y and z, with means and , defined at the zeroth lag, is given by

To measure the network correlation, we created two random groups of neurons from the excitatory neuron population (excepting differences in firing rate, the inhibitory neurons behave in same manner as the excitatory neurons). The two groups each consisted of 5% or 250 neurons. The cross-covariance coefficient was computed for every pair of neurons, one from each group (excluding identical neurons, which occurred by chance) and the average cross covariance coefficient was taken as a measure of the network correlation.

Mean-field models for sparsely connected neural circuits

For our analysis of network stability we developed neuronal population models—that is, mean-field firing rate models—for sparse, random networks (van Vreeswijk and Sompolinsky 1996, 1998). They model populations of current-based neurons with static and dynamic synapses and are based on a balance of excitation and inhibition that creates temporally irregular spiking. We assume that there are N_e excitatory neurons and N_i inhibitory neurons. K_mn is the number of synaptic connections from neurons in pool n to neurons in pool m (we assume that K_mn << N_e, N_i). We are interested in networks with randomly chosen sparse synaptic connections because they provide reasonable first approximation models for cortical circuits.

Because the firing rate of a neuron is dependent on both the mean input current and the variance of the input current, we refer to its input–output function as the FMS (firing rate, mean input, SD) surface, which we denote by F. A sparse, random network has a parameter range that allows for asynchronous irregular activity (Brunel 2000; for a review, see Vogels et al. 2005) and the synaptic currents to a given cell are uncorrelated in time for this regime, aside from a potentially fluctuating population firing rate. In the steady-state condition, which is when all the dynamics of the system have been allowed to equilibrate, this allows one to use the time average of the synaptic input for estimating the synaptic input at each point in time. Thus in the asynchronous condition, and assuming Poisson spiking statistics, one can write down the average recurrent synaptic current I^rec and recurrent synaptic input variance ^rec2 for each neuron in their respective pools m (Shriki et al. 2003) as

(9)

where m, n {e, i} represent indices into the post- and presynaptic pools, respectively. So the total average current and total input current variance into a neuron in a pool m is

(10)

where I^inject and ^noise2 represent the mean and variance of external input to the neuron, and I_m^tot and _m^tot2 are the mean and variance of the total current injection to each cell (see METHODS).

To use these formulas to create a working neuron model one can easily sample the FMS surface of a single current-based neuron by injecting current with an arbitrary mean and Gaussian noise component and then use the Eq. 9 system to translate the parameters (J_ee, J_ie, x_e, J_ei, J_ii, x_i) into (I^tot, ^tot). Using this technique we sampled the FMS surface F of a current-based LIF neuron with the same neuronal parameters as the neurons used in our current-based spiking networks. Using the samples generated from this neuron we used smoothing splines to create a nonparametric approximation to F that could be evaluated at any given point to yield a good approximation of the postsynaptic firing rate (Fig. 2). Because the excitatory and inhibitory neurons in our current-based networks have the same neuronal parameters, we used the same neuron model for both types. Finally, the mean-field model for a network of excitatory and inhibitory neurons in an asynchronous network, which also represents a sparse, random network under the same assumptions, is defined by

(11)

where x_e is the average firing rate of the excitatory neurons and x_i is the average firing rate of the inhibitory neurons. For a current-based LIF network with _e, _i < _mem, the correct time constant of the model is the neuronal membrane time constant _mem.

View larger version (36K): [in this window] [in a new window]   FIG. 2. FMS (firing rate, mean current injection, SD of current noise injection) surface surface for leaky integrate-and-fire (LIF) neuron. A: approximated FMS surface of LIF neuron with m = 10 ms, rp = 3 ms, Vrest = –60 mV, and Vreset = –60 mV. Approximated gives a very good fit with mean absolute error 0.2 Hz. Dashed black line shows the firing rate as a function of the external current noise level assuming a constant current injection of 2.455 x 10–9 A. Solid black line shows the firing rate as a function of mean input, with injected noise value of 6 x 10–9 A. Intersection of these 2 lines denotes the external input (mean current and noise level) that every neuron in all simulations of spiking neurons received and resulted in a steady-state firing rate of 20 Hz in the absence of synaptic feedback.

Finally, we introduce dynamic synapses (Eq. 7 system) into the sparse, random mean-field model of the Eq. 11 system. The model is defined by the expressions in Eq. 11 with slightly different definitions for I_m^rec and _m^rec2

(12)

where have replaced the static synaptic weight value J_mn with the dynamic synaptic weight value µ_mn(x_n, t).

	RESULTS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
Conditions on dynamic synapses to create network stability

The goal of introducing dynamic synapses into the sparse, random model is to set up a stable fixed point at a desired target network firing rate x_e = x_i = x*, so that x_e is robust to various types of perturbations (i.e., the circuit is self-tuning). We view self-tuning as the ability of a system to ultimately return a system value to given range, despite a perturbation in system parameters or input. We define the output of the network to be the average excitatory firing rate x_e, which is the firing rate that is actively controlled, perhaps by increasing or decreasing x_i. Although sparse, random networks also have steady-state firing rates as a function of static synaptic weight values J_mn, we will show that introducing dynamic synapses with relatively general constraints endows the network with much greater robustness under a variety of perturbations. This improved robustness is achieved if the UDF parameters (the parameters of the dynamic synapse model; see METHODS) of the dynamic synapses are chosen to satisfy the following three conditions.⁴

1) Maintain the synaptic weight values if the presynaptic firing rate is equal to the desired firing rate.

2) Increase the strength of the E E and I I synapses and decrease the strength of the E I and I E synapses if the presynaptic firing rate is below the target firing rate.

3) Decrease the strength of the E E and I I synapses and increase the strength of the E I and I E synapses if the presynaptic firing rate is above the target firing rate.

The first constraint involves creating the fixed point at a target firing rate by setting the scaling factor A of the dynamic synapse model according to the formula

(13)

where J_mn represents specific synaptic weight values of a corresponding network with static synapses and R_mn* (x*) and U_mn^1* (x*) are the steady-state values for the internal dynamic synapse variables as a function of presynaptic firing rate, here set to the desired steady-state firing rate (see Eqs. 6 and 7 in METHODS). Equations 6 and 13 ensure that the steady-state dynamic synapse weight is equal to the static synaptic weight [µ_mn* (x*) = J_mn] for a presynaptic firing rate x*. We refer to µ_mn* (x*) as the steady-state value for µ_mn (see METHODS, Eq. 8) and designate the resulting function of x* as the steady-state synaptic weight curve (see Fig. 3). In the asynchronous range, when static synaptic weight values (J_ee, J_ei, J_ie, J_ii) are set for a sparse, random network, the network has a specific steady-state firing rate associated with those synaptic strengths, x* (van Vreeswijk and Sompolinsky 1996, 1998). Because the network already has a fixed point with firing rate x* for these static weights, we simply scale the dynamic synapses to take these static weights at x*. We do this by setting the value of A_mn according to Eq.13, which gives µ_mn* (x*) = A_mnR_mn* (x*)U^1*(x*) = J_mn[R_mn* (x*)U^1*(x*)/R_mn* (x*)U^1*(x*)] = J_mn, thus returning the firing rate x*.

View larger version (11K): [in this window] [in a new window]   FIG. 3. Steady-state synaptic strengths.µmn* are used to control the excitatory firing rate in a network with dynamic synapses. A: synaptic strengths to the excitatory population as a function of presynaptic firing rate. Dashed, straight lines show static synapse values such that xe = xi = 10 Hz (blue, Jee; red, Jei). Thicker, curved lines show the dynamical synapse steady-state strengths tuned such that exactly at 10 Hz (black ""), the static and dynamic synapse values match up [blue, µee(xe); red, µei(xi)] and therefore a network with dynamic synapses also supports the 10-Hz excitatory and inhibitory firing rate of the network with static synapses. B: same as A but for the inhibitory population. Note the inversion of the slopes of the steady-state synaptic strengths. Dashed, straight lines: blue, Jie; red, Jii; thicker, curved lines: blue, µie(xe); red, µii(xi). Arrows in A show how the dynamic synapses will push the excitatory firing rate toward 10 Hz in the event of a system perturbation. If a perturbation causes xe to increase, then µee(xe) decreases (A) and µie(xe) increases (B), thereby causing a decrease in xe and an increase in xi; the latter induces an increase in µei(xi), causing an additional decrease in xe by xi. If a perturbation causes xe to decrease then the same curves explain an increase in xe.

View larger version (38K): [in this window] [in a new window]   FIG. 4. Two volumes in the UDF [designation for collective set of parameters: U (which roughly models the release probability of a synaptic vesicle for the first spike in a train of spikes), D (time constant for recovery from depression), and F (time constant for recovery from facilitation)] parameter space where the steady-state synaptic weights.µmn* (xmn) have either positive or negative derivatives with respect to presynaptic firing rate, which are essential for the theoretical prediction of resulting self-tuning properties of the resulting network. A and B: randomly generated UDF parameters (E E and I I, denoted by blue "*"; and I E and E I denoted by red "*"), which are used to test the theoretically predicted stability properties of UDF values from these 2 volumes. A: "N-Volume" of UDF space (blue) denotes that parameter range in which µ*(xpre) has a negative derivative between the firing rates of 10 and 100 Hz (also in C, E, and G). B: "P-Volume" of UDF space (red) denotes that parameter range in which µ*(xpre) has a positive derivative between the firing rates of 10 and 100 Hz (also in D, F, and H). C and D: sampled UDF values for pyramid to pyramid (blue circles) and pyramid to interneuron (green circles) connections based on means and SDs from Markram et al. (1998). E and F: sampled UDF values for F1 facilitating synapses (green "+"), F2 depressing synapses (cyan "+"), and F3 recovering synapses (black "+") based on means and SDs from Gupta et al.(2000). These are all UDF values of GABAergic synapses whose targets are either excitatory or inhibitory neurons; thus no clear prediction for the ranges of their values can be made through the proposed theoretical analysis. UDF values of F1 and F2 synapses lie close to the border region between the P- and N-volume. G and H: sampled UDF values for medial prefrontal cortex excitatory-to-excitatory synapses (PFCe1, PFCe1a, and PFCe1b samples as yellow""; PFCe2, PFCe2a, and PFCe2b samples as black ""; and PFCe3, PFCe3a, and PFCe3b samples as magenta "") based on means and SDs from Wang et al. (2006). All these UDF values lie in the N-volume, as predicted by the theoretical analysis.

Network parameter ranges

To find ranges of synapse strengths for static synapses that were compatible with low firing rates and low network correlation in current-based or conductance-based networks, we varied the excitatory and inhibitory synaptic weights of conductance-based and current-based networks (see METHODS and Fig. 5). Specifically, we varied two synaptic weight values, one excitatory synaptic weight J_e and one inhibitory synaptic weight J_i, and set the all synaptic weight values such that J_ee = J_ie = J_e and J_ei = J_ii = J_i. So that results between conductance- and current-based networks could be compared we scaled the conductance-based synapses so that the conductance- and current-based synapses had equal strength when the conductance-based synapses were multiplied by their respective driving forces at the neurons’ resting potential. We found a broad parameter range of synaptic weights that could potentially support dynamic synapses as a network feedback mechanism in the asynchronous, irregular regime. However, because we measured the firing rate, network correlation, and coefficient of variation of the interspike interval (ISI CV) for only a 2D manifold of the four-dimensional (4D) space of synaptic weights (J_ee, J_ei, J_ie, J_ii), these parameter maps serve only as a guide with which to start our simulations because the parameter space becomes 4D when the dynamic synapses are allowed to vary the synaptic strengths according to individual dynamic synapse formulas with different AUDF parameters. We always started our networks with synaptic weight initial conditions from these 2D maps, so that the networks started in a setting where the network behavior is asynchronous irregular and the feedback current magnitudes are on par with or smaller than the level of the background current injection. Based on these maps, we used the synaptic weight parameter range shown in Fig. 5, A, D, and G as the starting point for the rest of our conductance- and current-based spiking network simulations, as well as our sparse, random mean-field simulations of the Eq. 11 system, all of which used dynamic synapses. Specifically, J_ee, J_ie [0, 0.1] x 10^–9 A and J_ei, J_ii [0, 1.67] x 10^–9 A for current-based spiking networks and the current-based mean-field model and J_ee, J_ie [0, 1.67] x 10^–9 S and J_ei, J_ii [0, 10] x 10^–9 S for the conductance-based spiking networks and the dynamic synapses were set to these values using Eq. 13.

View larger version (35K): [in this window] [in a new window]   FIG. 5. Parameter search over synaptic weights for conductance-based spiking networks, current-based spiking networks, and mean-field firing-rate networks with static synapses. y-axis is Je = Jee = Jie and the x-axis is Ji = Jei = Jii. A–C: network statistics for conductance-based spiking networks with static synapses. A: network steady-state firing rate. B: coefficient of variation (CV) for conductance-based spiking networks over the range of synaptic weights. C: network correlation of excitatory population as measured as described in METHODS. There was no significant network correlation. D–F: network statistics for current-based spiking networks with static synapses. D: network steady-state firing rate. Five colored circles represent the current-based spiking network settings used in the target firing-rate simulation. E: CV for current-based spiking networks over the range of synaptic weights. F: network correlation of excitatory population. There was no significant correlation over the range of synaptic weights tested. G: network steady-state firing rate for mean-field firing rate model.

Interaction of synaptic dynamics with network activity

To give an idea of how dynamic synapses interact with the network, we discuss an example in detail. The example in Figs. 6 and 7 shows a burst in a network of current-based spiking neurons followed by a convergence to a steady-state firing rate (a conductance-based example would work nearly identically, except that the effect on the synaptic conductances would be harder to visually interpret as a result of the effect of the driving force on each synapse). We also show a sparse, random mean-field simulation using the same parameters. A network was chosen that normally fires at 20 Hz in the absence of dynamic synapses, J_ee = J_ie = 0.05 x 10^–9 A and J_ei = J_ii = –0.1 x 10^–9 A. The network was initialized with the dynamic synapse parameter initial conditions set to a 5-Hz steady-state firing rate and the desired firing rate for the network was set to 10 Hz, by adjusting the scaling factors A_mn in the Eq. 8 system. The firing rate initial conditions for the network were also set to 5 Hz, as best as possible. Because a firing rate of 5 Hz gives an overall positive network gain when the target firing rate is set to 10 Hz, the firing rate shot up at the beginning of the simulation (Fig. 6, A and B). The dynamic synapses adapted to this new firing rate and brought it back down, ultimately reaching nearly the steady-state of 10 Hz for the excitatory population, whereas the inhibitory population firing rate remained near 20 Hz. Thus we see that an excitatory firing rate very near the target firing rate x_e* was reached from quite different initial conditions and then maintained. During the simulation, the average dynamic synapse weights were modified and caused the changes shown in the network firing rates (Fig. 6, C–J). We can assert that the dynamic synapses caused the change in network firing rate (and not the other way around) because the time constants of the dynamic synapses are much longer than any neuronal time constants. Again there is excellent agreement between the current-based mean-field model and the current-based spiking network simulation.

View larger version (23K): [in this window] [in a new window]   FIG. 6. Example showing the firing rates and synaptic weights over time to demonstrate the effect of dynamic synapses on the network firing rate. Scaling factors Amn were set such that the desired steady-state firing rate is 10 Hz, even though the synaptic weight values [Je, Ji] = [0.5, –1] x 10–10 A gave a steady-state firing rate of 20 Hz without dynamic synapses. Initial conditions of the dynamic synapses were set to values corresponding to 5-Hz steady-state firing rate. A: firing rate curve of mean-field model: xe in blue, xi in red. B: firing-rate curve of network simulation; the firing rate is smoothed with a low-pass filter (box car filter with length equivalent to 20 ms). Two firing rate curves are in good qualitative agreement (the drop in firing rates at 500 ms is only a filtering artifact). C–J: synaptic weights through time as set by the dynamic synapse equations. C and D: µee. E and F: µei. G and H: µie. I and J: µii. C, E, G, and I: mean-field model. D, F, H, and J: spiking network model. Dashed lines represent the values Jee, Jei, Jie, and Jii. Comparison of the panels on the left and right shows a good agreement between the mean-field model and the spiking network.

View larger version (42K): [in this window] [in a new window]   FIG. 7. Detailed analysis of network activity in the network example from Fig. 6. A: spike raster showing the spiking activity of 500 neurons (out of 5,000): 400 excitatory and 100 inhibitory. B–D: data from a single excitatory neuron. B: total synaptic excitatory and inhibitory currents. Blue line shows the total excitatory input over time; red line shows the total inhibitory input over time; black line shows the total amount of current into the neuron excluding the noise injection; green dashed line shows the mean current injection into the neuron, which corresponds to a 20-Hz firing rate in the absence of synaptic feedback. C: individual synaptic currents (each color represents a separate synapse to the neuron). Individual synaptic currents are plotted to give an idea of changes in synaptic weight values. Upward exponential curves represent excitatory currents and downward, inhibitory currents. Peak values of these currents are set by the dynamic synapses and vary roughly over twofold the magnitude of the average current. D: membrane potential of the selected neuron, with action potentials. Thick black line represents the membrane potential of the cell; thin black line represents the average; blue line shows the neuron's firing threshold.

Self-tuning in networks of spiking neurons

To investigate the self-tuning properties of dynamic synapses on spiking networks we asked whether dynamic synapses could help tune a spiking neural network to a steady-state target firing rate even if the external input was perturbed because it occurs, for example, in primary sensory areas of the cortex in response to sensory stimuli. To test this, we used a conductance-based spiking network (J_ee = J_ie = 0.4 x 10^–9 S and J_ei = J_ii = 8.48 x 10^–9 S), a current-based spiking network (J_ee = J_ie = 0.013 x 10^–9 A and J_ei = J_ii = –0.18 x 10^–9 A), and a current-based mean-field model (J_ee = J_ie = 0.013 x 10^–9 A and J_ei = J_ii = –0.18 x 10^–9 A), all with a steady-state firing rate of 10 Hz. We additionally performed all of our self-tuning tests with the current-based mean-field model because it aids us in understanding the self-tuning properties of networks with dynamic synapses and shows, in particular, that the self-tuning properties of a circuit with dynamic synapses are a property of the dynamic synapse model and not something hidden in the many details and parameters of a complicated spiking neural network simulation. Every simulation was performed using both the parameters based on experimental data and the R1–R3 dynamic synapse parameters that were derived from our theoretical considerations (see METHODS). The scaling factors A_mn of the dynamic synapses were set only once such that the steady-state dynamic synapse curves intersected at 10 Hz. Thus each simulation had identical dynamic synapse parameters. We then perturbed in small increments both the mean and SD of the external input from –50 to +50% of their normal values (I^inject = 2.455 x 10^–9 A and noise drawn from a 0-mean Gaussian, ^inject = 6 x 10^–9 A SD). The results, which are the steady-state excitatory firing rates from the last 1 s of each simulation, are shown in Fig. 8. All the networks that had a steady-state firing rate between 9 and 11 Hz (a difference of 1 Hz from the target excitatory firing rate) are shown in blue, networks with firing rates between 8 and 9 and 11 and 12 Hz (a difference between 1 and 2 Hz) are shown in cyan, and networks with firing rates between 7 and 8 and 12 and 13 (a difference between 2 and 3 Hz) are shown in green. In Fig. 8, A, D, and G we show the control case with static synapses for the three different network models for comparison. There were few networks that supported a 10-Hz steady-state firing rate. The maximal firing rates were >40 Hz and many networks were essentially shut down with a firing rate very close to or equal to 0 Hz. Thus the steady-state firing rates of networks with static synapses are quite dependent on mean current injection and injected noise level. When we ran the same simulations with dynamical synapses turned on (Fig. 8, remaining panels), the large firing rates were reduced. We examined two different assignments of UDF values to the four different types of synaptic connections in the network: values that are sufficient for inducing stability according to the previously (see first section of RESULTS) described analytical approach (Fig. 8, B, E, and H) and values taken from experimental data (Fig. 8, C, F, and I). Interestingly, for this perturbation test the dynamic synapse parameters based on experimental data somewhat outperformed the dynamic synapses parameters based on the analytical approach.

View larger version (51K): [in this window] [in a new window]   FIG. 8. Self-tuning with respect to perturbations of the level and variance of external input currents. All networks tested were set to a default steady-state firing rate of 10 Hz. Scaling factors Amn were set only once; thus each set of simulations (a single panel) where dynamic synapses were used had identical synapse parameters. Blue denotes the optimal firing rate of 10 Hz (between 9 and 11 Hz), cyan between 8 and 9 and 11 and 12 Hz, green between 7 and 8 and 12 and 13 Hz, and gray denotes networks between 0 and 1 Hz. A–C: networks with conductance-based synapses. A: response (average excitatory firing rate) of conductance-based networks with static synapses to perturbations of external current mean and noise. B: averaged response of conductance-based networks over dynamic synapse parameters R1–R3. C: response of conductance-based networks with the dynamic synapse parameters based on experimental data. D–F: networks with current-based synapses. D: response of current-based networks with static synapses to perturbations of external current mean and noise strength. E: averaged response of current-based networks over dynamic synapse parameters R1–R3. F: response of current-based networks with the dynamic synapse parameters based on experimental data. G–I: response of current-based mean-field networks. G: response of current-based mean-field model with static synapses to perturbations of external current mean and noise strength. H: averaged response of mean-field networks over dynamic synapse parameters R1–R3. I: response of the mean-field networks with the parameters based on experimental data.

View larger version (39K): [in this window] [in a new window]   FIG. 9. Transient firing rate responses in network models with dynamic synapses. A–C: transient firing rate responses for changes in the mean and variance of external inputs are shown for conductance-based networks with experimentally observed UDF parameters for dynamic synapses (see Table 1). A: firing rate of excitatory neurons during the time interval 1–51 ms (in response to a change in external inputs at time 0). Black circle indicates the specific network that is shown in G. Firing rate of excitatory neurons over periods of 251–300 and 450–500 ms (B and C, respectively). D–F: same for conductance-based networks using R1 dynamic synapse parameters. G: average firing rate of excitatory (blue) and inhibitory (red) neurons as a function of time (averaged with a 20-ms sliding window) for the network represented in A as a black dot. H: same for the network represented in D as a black dot (the same network with static synapses assumed a steady-state firing rate of 19 Hz for both excitatory and inhibitory neuron populations; not shown). One sees that in either case the transient firing rate transmitted substantial information about the new level of external input during a significant amount of time after the change in the external input, suggesting that a network-induced homeostasis of firing rates is compatible with computational processes that process information about the change in external inputs.

The results of the synaptic weight perturbations for conductance-based spiking networks, current-based spiking networks, and current-based mean-field models are shown in Fig. 10. In Fig. 10, A, D, and G are the control cases with static synapses, shown for comparison. With static synapses the steady-state firing rates for networks with high excitatory synaptic weights are in excess of 200 Hz. However, as shown in Fig. 10, B, C, E, F, H, and I, in networks with dynamic synapses, regardless of conductance-based (A–C) or current-based (D and E) synapses and regardless of which dynamic synapses parameters are used, the excessive firing rates disappear. Most important, in all but the one case of current-based networks with the dynamic synapse parameters based on experimental data, the number of perturbed networks that nevertheless had firing rates very close to the target firing rate increased dramatically. The continuous mean-field predictions of the dynamic synapses applied to the networks with synaptic weight perturbations are shown in Fig. 10, G–I. The mean-field simulations show excellent performance across nearly the entire range of perturbations and are qualitatively similar to spiking network simulations using dynamic synapses. Network correlation did not play a role in the effect the dynamic synapses had on the networks.

View larger version (41K): [in this window] [in a new window]   FIG. 10. Self-tuning with respect to systematic synaptic strength perturbations. All networks tested were set to a default steady-state firing rate of 10 Hz. Blue denotes the optimal firing rate of 10 Hz (between 9 and 11 Hz), cyan between 8 and 9 and 11 and 12 Hz, green between 7 and 8 and 12 and 13 Hz, and gray denotes networks between 0 and 1 Hz. Synaptic strength scaling factors Amn for all synaptic perturbations used the same rescaling factor. A–C: response (as measured by the steady-state average excitatory firing rate) of the conductance-based spiking networks to systematic synaptic strength perturbations. A: response for conductance-based spiking networks with static synapses. B: response of conductance-based spiking networks with dynamic synapses using the R1–R3 parameters. C: response for conductance-based spiking networks with dynamic synapses using the parameters based on experimental data. D–F: response of current-based spiking networks to systematic synaptic strength perturbations. D: response of current-based networks with static synapses. E: response of current-based networks with dynamic synapses using the R1–R3 parameters. F: response of current-based networks with dynamic synapses using the parameters based on experimental data. G–I: response of current-based mean-field firing rate models to systematic perturbations in synaptic strength. G: response of current-based mean-field models with static synapses. H: response of current-based mean-field models with dynamic synapses using the R1–R3 parameters. I: response of current-based mean-field models with dynamic synapses using the parameters based on experimental data. Mean-field model differs from the current-based network in F; see footnote 5 for discussion.

Why should simply scaling the A_mn for dynamic synapses according to Eq. 13 have a helpful effect on the target firing rate of the network if there is an unknown perturbation of the factor J_mn in that same equation? Assume, for example, that there are static synaptic strengths (J_e*, J_i*) that lead to a target firing rate of x*. Assume further that there is an unknown perturbation of these synaptic strengths that can be written as (J_e* + _e, J_i* – _i). Because we have increased the excitation and decreased the inhibition we should expect an increase in the target firing rate; that is, the perturbed synaptic strengths (J_e* + _e, J_i* – _i) lead to x* + x. If we now examine the same network and synaptic weight perturbation but use dynamic synapses we see a change in behavior. As already outlined we set the A_mn according to Eq. 13, so that µ_me* (x*) = J_e* + _e and µ_mi* (x*) = J_i* – _i because we do not know that there is a perturbation. Now the situation with dynamic synapses is quite different from with static synapses because at x* the dynamic synapses equilibrate to values that support a higher firing rate than x*, specifically x* + x. At a firing rate of x* + x, however, because the derivatives of the dynamic synapses have the correct signs, the network supports a firing rate below x* + x. Thus we can conclude that the firing rate must be between x* and x* + x, which is better than the case of static synapses, where the network firing rate is completely dictated by the perturbed synaptic weights and thus has a firing rate of x* + x. The degree to which the final firing rate will be closer to x* or x* + x depends on the magnitude of these derivatives. We refer to the last part of the Supplementary material for an analytical comparison of network stability with static and dynamic synapses.

As a final self-tuning test we also investigated the effects of systematic changes in the excitability of subpopulations of neurons in the simulated circuit. Diverse changes in the excitability of neurons were previously reported to be caused by neuromodulators (especially dopamine) in various cortical areas (Gulledge and Jaffe 1998, 2001; e.g., see Fig. 2B in first reference). Because such changes in the excitability of neurons cannot be modeled directly in the context of LIF neurons, we simply inactivated (i.e., deleted) randomly selected subsets of neurons in the circuit. For continuity, we chose the same networks as in the input perturbation self-tuning simulation: conductance-based spiking networks (J_ee = J_ie = 0.4 x 10^–9 S and J_ei = J_ii = 8.48 x 10^–9 S), current-based spiking networks (J_ee = J_ie = 0.013 x 10^–9 A and J_ei = J_ii = –0.18 x 10^–9 A), and the current-based mean-field models (J_ee = J_ie = 0.013 x 10^–9 A and J_ei = J_ii = –0.18 x 10^–9 A), all with a steady-state firing rate of 10 Hz. We randomly chose excitatory and inhibitory neurons in these networks to inactivate, with the total number of inactivated neurons between 0 and 70% of their respective pool totals, in small increments. The current-based mean-field model received an analogous perturbation by modifying the K_mn parameter in Eq. 9 according to the neuronal inactivation percentage in the spiking networks. The control simulations with static synapses are shown in Fig. 11 A, D, and G; results with the R1–R3 dynamic synapse parameters are shown in Fig. 11, B, E, and H; and results with the dynamic synapse parameters based on experimental data are shown in Fig. 11, C, F, and I. For the networks tested, the dynamic synapse parameters based on experimental data did not have a significant effect but the R1–R3 dynamic synapse parameters showed improvement over the control, with the vast majority of networks within 2 Hz of the target firing rate of 10 Hz. Again, network correlation effects did not play a role in the results of the simulations.

View larger version (47K): [in this window] [in a new window]   FIG. 11. Sparse, spiking networks demonstrate self-tuning with respect to systematic neuron inactivation. All networks tested were set to a default steady-state firing rate of 10 Hz. Scaling factors Amn were set only once; thus each set of simulations (a single panel) where dynamic synapses were used had identical dynamic synapse parameters. Blue denotes the optimal firing rate of 10 Hz (between 9 and 11 Hz), cyan between 8 and 9 and 11 and 12, green between 7 and 8 and 12 and 13 Hz, and gray denotes networks between 0 and 1 Hz. A–C: response (as measured by the steady-state average excitatory firing rate) of the conductance-based spiking networks to systematic neuronal inactivation. A: response of conductance-based spiking networks with static synapses. B: averaged response of conductance-based spiking networks with dynamic synapses using the R1–R3 parameters. C: response of conductance-based spiking networks with the dynamic synapse parameters based on experimental data. D–F: response (as measured by the steady-state average excitatory firing rate) of the current-based spiking networks to systemic neuronal inactivation. D: response of current-based spiking networks with static synapses. E: averaged response of current-based spiking networks with dynamic synapses using the R1–R3 parameters. F: response of current-based spiking networks with dynamic synapses using the parameters based on experimental data. G–I: response (as measured by the steady-state excitatory firing rate) of the current-based mean-field firing rate models to systemic neuronal inactivation. G: response of current-based mean-field models with static synapses. H: averaged response of current-based mean-field models with dynamic synapses using the R1–R3 dynamic synapse parameters. I: response of current-based mean-field models with dynamic synapses using the parameters based on experimental data.

Implicit in the demonstrated self-tuning properties is the ability to set a specific fixed-point firing rate using dynamic synapses, through use of the scaling in Eq. 13. The network perturbations can then be viewed as a perturbation away from this fixed point. We thus wondered, in principle, how close a spiking network could get to a preprogrammed fixed point, as set by Eq. 13, regardless of what the initial synaptic strengths of the network were. Briefly, we used five networks from the current-based, spiking network synaptic parameter range of Fig. 5D. The networks had J_e and J_i as follows (in nA): blue [0.025, –0.05], red [0.025, –0.15], cyan [0.05, –0.1], magenta [0.075, –0.05], and green [0.075, –0.015], leading to corresponding firing rates of 20, 12, 20, 100, and 20 Hz with static synapses. We then chose as the target firing rates every integral firing rate in the range of 1–120 Hz and used Eq. 13 to appropriately set the dynamic synapses in the five networks. We do not mean to suggest that cortical networks create fixed points at such high firing rates, rather only to test the principle behind Eq. 13. The results, which are the steady-state excitatory firing rates, are shown in Fig. 12A. All networks showed an ability to enhance their firing rates toward the desired target firing rate, although the two networks with the weakest excitatory synapses (denoted by blue and red) were unable to follow the target firing rate closely at >20 Hz. Another complication arose in that the two networks with the strongest excitatory synapses (green and magenta), although able to closely follow the target firing rate, also showed significant network correlation at >20 Hz, shown in Fig. 12B. Despite these caveats, the results clearly demonstrate that the ability to tune the steady-state firing rate of a network is a general property of the dynamic synapse configuration. The results lend further evidence to the idea that self-tuning properties of spiking networks with dynamic synapses are caused by setting a firing rate fixed point, to which the networks attempt to return after a perturbation.

View larger version (15K): [in this window] [in a new window]   FIG. 12. Many spiking networks can follow a target firing rate using dynamic synapses. Five selected current-based spiking networks with dynamic synapses showing the networks’ abilities to follow an excitatory target firing rate using scaling factors Amn as set by the heuristic of Eq. 13. Target firing rate was set to every integral firing rate between 1 and 120 Hz and the A values for each network were calculated based on this target frequency. A: average target firing rate; black dashed line shows the line y = x, which is optimal performance. Mean absolute errors for each network blue, red, cyan, magenta, and green are [20.7, 21.4, 4.5, 2.1, 2.0] Hz, respectively. Networks had Je and Ji as follows (in nA): blue [0.025, –0.05], red [0.025, –0.15], cyan [0.05, –0.1], magenta [0.075, –0.05], and green [0.075, –0.015] and are also denoted in Fig. 5D. Dynamic synapse parameters used were the R1 parameters. B: excitatory network correlation as a function of target firing rate. There was negligible network correlation for the blue, red, and cyan networks and significant network correlation for the magenta and green networks between the firing rates of 20 and 100 Hz.

There is little hope that one can give a complete theoretical justification for the demonstrated self-tuning properties of sparsely connected circuits of spiking neurons with dynamic synapses or for the associated mean-field model defined by Eqs. 7–9. The reason lies in the nonlinearity of these models and the fairly large number of variables that are involved. However, the demonstrated self-tuning properties of these complex systems are clearly related to known self-tuning properties of substantially simpler systems for which a theoretical analysis is feasible. These analytical results will be discussed in the following.

An analytical result on the assignment of differential dynamics to different types of synapses in a distributed circuit

In the first section of RESULTS we proposed a heuristic for choosing the signs of the derivatives of the synaptic weight curves µ_mn* (x_n) at x*. We then demonstrated through numerical simulations of circuits of spiking neurons that this choice was justified. This heuristic can also be justified analytically for a 2D mean-field model, similar to that described in Eq. 11 (assuming dynamic synapses with very small derivatives and instantaneously equilibrating dynamic synapses). One can prove that the derivatives d_mn = (dµ_mn/dx_n)(x*) should have the following signs for the excitatory population to have homeostatic properties with respect to changes in the intensity of external input currents: d_ee < 0, d_ii > 0, d_ie > 0, d_ei < 0 (see Supplementary material for details). These signs agree with the previously proposed heuristic rules for the assignment of depressing and facilitating synapses to the four types of synaptic connections (see first section of RESULTS; note that the negative sign is included in µ_ei and µ_ii). For the excitatory population to also have homeostatic properties with respect to changes in the intensity of external input currents to the inhibitory population the analytically derived condition is somewhat more complicated

where _e is the derivative of the excitatory F–I curve evaluated at the fixed point x*. Taken together these conditions imply that regardless of the size of µ_ee(x*), the E E connection should always be depressing and the E I should always be facilitating to allow the circuit to return to a steady-state firing rate in spite of a changed level of external input. Additionally, if µ_ee(x*) > _e^–1 the I E synapses should be facilitating and the I I synapses should be depressing; otherwise, the I E synapses should be depressing and the I I synapses should be facilitating. The fact that this analytical result does not point to a general rule for the optimal dynamics of synapses from inhibitory neurons (from the perspective of firing rate stability) could be seen as a possible explanation for the diverse distribution of experimentally observed UDF parameters for these synapses (see Fig. 4, E and F). We additionally give justification for our self-tuning heuristic for the parameter range used in our current-based networks in the Supplementary materials.

A theoretical analysis of the impact of changes of synaptic weights on stability properties of networks with static and dynamic synapses

One can prove for the same mean-field model with excitatory and inhibitory neurons that the steady-state firing rate x* is less affected by small changes in synaptic weights (resulting for example from LTP or LTD) if the synaptic connections are dynamic and tuned using an appropriate heuristic. An analytical derivation is provided in the second part of the Supplementary material, with specific reference to the fact that our self-tuning heuristic is appropriate for the parameter regime used in the current-based spiking networks herein. An additional example figure (Fig. S1) is also given.

	DISCUSSION

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
We have shown that the inclusion of more realistic models for synapses, which reflect their experimentally found short-term plasticity, provides models for cortical neural networks with interesting stability properties. In particular, it enables such models to return after a large variety of perturbations to a low but nonzero spontaneous firing rate, even if the perturbation, such as a change in the level of external input, a change in synaptic strength, or even neuron deletion, is of a longer-lasting nature. This surprising stability property of network models with dynamic synapses, which has apparently not been reproduced in models with static synapses, may provide a possible explanation for the surprisingly stable low but nonzero rate of spontaneous firing reported for cortical neurons in a variety of studies. Several other homeostatic processes, such as synaptic scaling, LTP and LTD, and genetic regulation of receptor numbers are likely to support the stability of the rate of spontaneous firing of cortical neurons on a larger timescale (Turrigiano and Nelson 2004). However, none of these mechanisms would apparently be able to support stability on the timescale of seconds, say, in response to a changed level of sensory input or to a changed concentration of neuromodulators in a cortical circuit. On the other hand, we showed in Fig. 9 that the stability property endowed by short-term synaptic plasticity leaves a sufficiently large time window of one to several hundred milliseconds during which the firing rate in the circuit is affected by the external perturbation. Thus information about such external inputs can be transmitted to other cortical circuits and integrated into cortical computation. In fact, the dynamic response to an increased external input shown in Fig. 9H is on a timescale similar to that of typical responses to visual inputs in primary visual cortex (see, e.g., Fig. 2 in Lamme and Roelfsema 2000). Even on a larger timescale, however, the information about the external (or internal) perturbation remains accessible to the neural system, although the firing rate of excitatory neurons is returned through synaptic dynamics to a given target level. The same firing rate of excitatory neurons is then accompanied by a different firing rate of inhibitory neurons (see Fig. 9, G and H) and this rate may therefore contain substantial information about the nature and level of the external perturbation for a much longer time period. This observation also implies that the same firing rate of excitatory neurons can be produced by a virtually infinite combination of concentration levels of neuromodulators and activity levels of inhibitory neurons, analogously to the large sets of different circuit parameters that were reported to support the generation of a target firing pattern in the pyloric circuit in the crustacean stomatogastric ganglion (Marder and Goaillard 2006).

The target value of 10 Hz for the spontaneous firing rate excitatory neurons in cortical circuits on which we have focused in our network simulations was chosen rather arbitrarily, although it lies within a range of spontaneous firing rates reported in a number of experimental studies (see INTRODUCTION). We show in Fig. 12 that a large range of other target values can be achieved by choosing a suitable global scaling parameter for synaptic weights, which might for example be under genetic control. However, to demonstrate stability at really low rates (say, 1 Hz), one needs to add neuron models to the network that fire spontaneously even in the absence of synaptic input (or during a release from inhibition); otherwise, the network cannot recover from a temporary pause in firing and will remain indefinitely silent.

We have shown in this article that the stability of the firing rate of excitatory neurons endowed by short-term synaptic dynamics is rooted in general principles of control theory, recently proposed for much simpler dynamical systems. In fact, activity-dependent rescaling of feedback strength was previously postulated as a powerful mechanism for self-tuning in control theory for purely mathematical reasons (Moreau and Sontag 2003). We demonstrated that this theoretical analysis of stability properties provided by dynamic synapses can be extended to the case of a distributed system that consists of several interacting dynamical modules, such as an interacting population of excitatory and inhibitory neurons in a cortical circuit, each with thousands of neurons. In particular, we have exhibited an analytical method for deriving rules regarding the type of synaptic dynamics that should be assigned to different types of synapses for the purpose of network stability. This theoretical approach may help to provide an explanation for the diversity of the empirically found short-term plasticity of connections between modules with different dynamic roles (such as excitatory and inhibitory neurons). It turns out that this theoretical analysis provides clearer predictions for the dynamics of synapses from excitatory neurons than for GABAergic synapses. This might provide an explanation for the experimentally found diversity in the dynamics of GABAergic synapses.

Finally, we emphasize that the investigation of the relationship between stability properties of networks of spiking neurons and the dynamics of specific types of synapses that we have initiated in this article should be seen only as a first step. In a second step one also needs to investigate through computer simulations more detailed network models consisting of different types of excitatory and especially inhibitory neurons with specific firing properties and more specific connectivity patterns. A quite challenging open problem is whether one can also extend the theoretical analysis of stability properties induced by dynamic synapses to more complex networks.

	GRANTS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
D. Sussillo was supported in part by a Fulbright Research Grant from the U.S. State Department. W. Maass was partially supported by the Austrian Science Fund (Projects P17229 [GenBank] and S9102-N13) and by the European Union Project Fast Analog Computing with Emergent Transient States (FP6-015879).

	ACKNOWLEDGMENTS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
The authors thank Drs. Stephan Haeusler and Eduardo Sontag for fruitful discussions that greatly aided this paper and Dr. Rodney Douglas for very helpful critical comments.

	FOOTNOTES

 
¹ It is shown in Moreau and Sontag (2003) that it suffices to assume that 1) f is decreasing, 2) g is strictly increasing, and 3) there exists some x* such that f(x*) = g(µ₀).

² The online version of this article contains supplemental data.

³ The variable u_k–1 in the third expression in the Eq. 6 system was erroneously replaced by uk in Markram et al. (1998); see the discussion in Maass and Markram (2002).

⁴ See the Supplemental materials for a theoretical justification of these choices.

⁵ In a spiking network without dynamic synapses there is a distribution of firing rates, whose mean we have used as the steady-state firing rate of the network. In these networks even the tails of the neuronal firing rate distribution have the same synaptic weight values J_mn. However, in the spiking network with dynamic synapses turned on, the lower part of the firing rate distribution will typically have synapses that are stronger than the set value µ_mn(x*) = J_mn for the target firing rate and the higher part of the distribution will typically have synapses that are weaker than J_mn. Because the steady-state dynamic synapse efficacy curves were constructed only to make sure the derivatives had the right signs there is no reason that the expected value of the steady-state synaptic weight across all synapses from neuron type n to m should exactly equal the mean-field average [i.e., that P(x_n) µ_mn (x_n) = J_mn, where P(x_n) is the distribution of firing rates for neuron type n]. Furthermore, changing the target firing rate may change the shape of the steady-state firing rate distribution of the network, thus making it effectively impossible to cancel this effect if one desires to tune the same network to different firing rates. Overall, this effect results in no more than a 1- to 2-Hz difference between the mean-field prediction of the effect of dynamic synapses and the spiking network results with dynamic synapses.

Address for reprint requests and other correspondence: W. Maass, Institute for Theoretical Computer Science, Technische Universitaet Graz, A-8010 Graz, Austria (E-mail: maass{at}igi.tugraz.at)

	REFERENCES

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
Abbott LF, Nelson SB. Temporal dynamics of biological synapses. In: The Handbook of Brain Theory and Neural Networks, edited by Arbib MA. Cambridge, MA: The MIT Press, 2003, p. 1156–1158.

Abbott LF, Regehr WG. Synaptic computation. Nature 431: 796–803, 2004.[CrossRef][Medline]

Abbott LF, Varela JA, Sen K, Nelson SB. Synaptic depression and cortical gain control. Science 275: 220–224, 1997.[CrossRef][Web of Science][Medline]

Abeles M. Corticonics: Neural Circuits of the Cerebral Cortex. Cambridge, UK: Cambridge Univ. Press, 1991.

Baddeley R, Abbott LF, Booth MC, Sengpiel F, Freeman T, Wakeman EA, Rolls ET. Responses of neurons in primary and inferior temporal visual cortices to natural scenes. Proc Biol Sci 264: 1775–1783, 1997.[Abstract/Free Full Text]

Brunel N. Dynamics of sparsely connected networks of excitatory and inhibitory spiking neurons. J Comput Neurosci 8: 183–208, 2000.[CrossRef][Web of Science][Medline]

Buzsáki G. Rhythms of the Brain. New York: Oxford Univ. Press, 2006.

Destexhe A, Contreras D. Neuronal computations with stochastic network states. Science 314: 85–90, 2006.[Abstract/Free Full Text]

Destexhe A, Rudolph M, Paré D. The high-conductance state of neocortical neurons in vivo. Nat Rev 4: 739–751, 2003.[CrossRef]

Galarreta M, Hestrin S. Frequency-dependent synaptic depression and the balance of excitation and inhibition in the neocortex. Nat Neurosci 1: 587–594, 1998.[CrossRef][Web of Science][Medline]

Gulledge AT, Jaffe DB. Dopamine decreases the excitability of layer V pyramidal cells in the rat prefrontal cortex. J Neurosci 18: 9139–9151, 1998.[Abstract/Free Full Text]

Gulledge AT, Jaffe DB. Multiple effects of dopamine on layer V pyramidal cell excitability in rat prefrontal cortex. J Neurophysiol 86: 586–595, 2001.[Abstract/Free Full Text]

Gupta A, Wang Y, Markram H. Organizing principles for a diversity of GABAergic interneurons and synapses in the neocortex. Science 287: 273–278, 2000.[Abstract/Free Full Text]

Haeusler S, Maass W. A statistical analysis of information-processing properties of lamina-specific cortical microcircuit models. Cereb Cortex 17: 149–162, 2007.[Abstract/Free Full Text]

Lamme VA, Roelfsema PR. The distinct modes of vision offered by feedforward and recurrent processing. Trends Neurosci 23: 571–579, 2000.[CrossRef][Web of Science][Medline]

Legenstein R, Maass W. Edge of chaos and prediction of computational performance for neural microcircuit models. Neural Networks In press.

Lisman JE. Bursts as a unit of neural information: making unreliable synapses reliable. Trends Neurosci 20: 38–43, 1997.[CrossRef][Web of Science][Medline]

Maass W, Markram H. Synapses as dynamic memory buffers. Neural Networks 15: 155–161, 2002.[CrossRef][Web of Science][Medline]

Maass W, Natschlager T, Markram H. Fading memory and kernel properties of generic cortical microcircuit models. J Physiol (Paris) 98: 315–330, 2004.[CrossRef][Medline]

Marder E, Goaillard JM. Variability, compensation and homeostasis in neuron and network function. Nat Rev 7: 563–574, 2006.[CrossRef]

Markram H, Wang Y, Tsodyks M. Differential signaling via the same axon of neocortical pyramidal neurons. Proc Natl Acad Sci USA 95: 5323–5328, 1998.[Abstract/Free Full Text]

Moreau L, Sontag E. Balancing at the border of instability. Phys Rev 68: 020901, 2003a.

Moreau L, Sontag E. Feedback tuning of bifurcations. Syst Control Lett 50: 229–239, 2003b.[CrossRef]

Rainer G, Lee H, Logothetis NK. The effect of learning on the function of monkey extrastriate visual cortex. PLoS Biol 2: e44, 2004.[CrossRef][Medline]

Shriki O, Hansel D, Sompolinsky H. Rate models for conductance-based cortical neuronal networks. Neural Comput 15: 1809–1841, 2003.[CrossRef][Web of Science][Medline]

Steriade M. The Intact and Sliced Brain. Cambridge, MA: The MIT Press, 2001.

Steriade M, Timofeev I, Grenier F. Natural waking and sleep states: a view from inside neocortical neurons. J Neurophysiol 85: 1969–1985, 2001.[Abstract/Free Full Text]

Thomson AM. Activity-dependent properties of synaptic transmission at two classes of connections made by rat neocortical pyramidal axons in vitro. J Physiol 502: 131–147, 1997.[Abstract/Free Full Text]

Thomson AM. Presynaptic frequency- and pattern-dependent filtering. J Comput Neurosci 15: 159–202, 2003.[CrossRef][Web of Science][Medline]

Tsodyks M, Pawelzik K, Markram H. Neural networks with dynamic synapses. Neural Comput 10: 821–835, 1998.[CrossRef][Web of Science][Medline]

Turrigiano GG, Nelson SB. Homeostatic plasticity in the developing nervous system. Nat Rev 5: 97–107, 2004.[CrossRef]

van Vreeswijk C, Sompolinsky H. Chaos in neuronal networks with balanced excitatory and inhibitory activity. Science 274: 1724–1726, 1996.[Abstract/Free Full Text]

van Vreeswijk C, Sompolinsky H. Chaotic balanced state in a model of cortical circuits. Neural Comput 10: 1321–1371, 1998.[CrossRef][Web of Science][Medline]

Varela JA, Sen K, Gibson J, Fost J, Abbott LF, Nelson SB. A quantitative description of short-term plasticity at excitatory synapses in layer 2/3 of rat primary visual cortex. J Neurosci 17: 7926–7940, 1997.[Abstract/Free Full Text]

Vogels TP, Abbott LF. Signal propagation and logic gating in networks of integrate-and-fire neurons. J Neurosci 25: 10786–10795, 2005.[Abstract/Free Full Text]

Vogels TP, Rajan K, Abbott LF. Neural network dynamics. Annu Rev Neurosci 28: 357–376, 2005.[CrossRef][Web of Science][Medline]

Wang Y, Markram H, Goodman PH, Berger TK, Ma J, Goldman-Rakic PS. Heterogeneity in the pyramidal network of the medial prefrontal cortex. Nat Neurosci 9: 534–542, 2006.[CrossRef][Web of Science][Medline]

Yoshimura Y, Callaway EM. Fine-scale specificity of cortical networks depends on inhibitory cell type and connectivity. Nat Neurosci 8: 1552–1559, 2005.[CrossRef][Web of Science][Medline]

Yoshimura Y, Dantzker JL, Callaway EM. Excitatory cortical neurons form fine-scale functional networks. Nature 433: 868–873, 2005.[CrossRef][Medline]

This article has been cited by other articles:

				P. J. Sjostrom, E. A. Rancz, A. Roth, and M. Hausser Dendritic Excitability and Synaptic Plasticity Physiol Rev, April 1, 2008; 88(2): 769 - 840. [Abstract] [Full Text] [PDF]

This Article Abstract Full Text (PDF) Supplemental Figure All Versions of this Article: 97/6/4079    most recent 01357.2006v1 Alert me when this article is cited Alert me if a correction is posted Citation Map Services Email this article to a friend Similar articles in this journal Similar articles in ISI Web of Science Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by Sussillo, D. Articles by Maass, W. Search for Related Content PubMed PubMed Citation Articles by Sussillo, D. Articles by Maass, W.