Abstract
In the intact brain neurons are constantly exposed to intense synaptic activity. This heavy barrage of excitatory and inhibitory inputs was recreated in vitro by injecting a noisy current, generated as an Ornstein–Uhlenbeck process, into the soma of rat neocortical pyramidal cells. The response to such in vivo–like currents was studied systematically by analyzing the time development of the instantaneous spike frequency, and when possible, the stationary mean spike frequency as a function of both the mean and the variance of the input current. All cells responded with an in vivo–like action potential activity with stationary statistics that could be sustained throughout long stimulation intervals (tens of seconds), provided the frequencies were not too high. The temporal evolution of the response revealed the presence of mechanisms of fast and slow spike frequency adaptation, and a medium duration mechanism of facilitation. For strong input currents, the slow adaptation mechanism made the spike frequency response nonstationary. The minimal frequencies that caused strong slow adaptation (a decrease in the spike rate by more than 1 Hz/s), were in the range 30–80 Hz and depended on the pipette solution used. The stationary response function has been fitted by two simple models of integrateandfire neurons endowed with a frequencydependent modification of the input current. This accounts for all the fast and slow mechanisms of adaptation and facilitation that determine the stationary response, and proved necessary to fit the model to the experimental data. The coefficient of variability of the interspike interval was also in part captured by the model neurons, by tuning the parameters of the model to match the mean spike frequencies only. We conclude that the integrateandfire model with spikefrequency–dependent adaptation/facilitation is an adequate model reduction of cortical cells when the mean spikefrequency response to in vivo–like currents with stationary statistics is considered.
INTRODUCTION
Single neuron properties have been thoroughly investigated in the past years showing that neural cells are rich in phenomenology and complex in their structure (see, e.g., Mainen and Sejnowski 1996; McCormick et al. 1985; Rhodes 1999). Even the most detailed stateoftheart model is unable to capture the entire phenomenology observed in the experiments. Such a richness calls for a model reduction that could provide a synthetic description of the response properties of the cells under particular conditions. We studied in vitro those features that are supposedly relevant when the cell is embedded in a large network of interconnected neurons, as it would be in in vivo conditions. The guidelines for selecting the relevant features were dictated by the theoretical framework developed in the last decade to study the dynamic properties of networks of integrateandfire neurons (see Gerstner and Kistler 2002 for a review). This approach was already successful in relating single neuron properties to several in vivo phenomena (Amit and Tsodyks 1991; Brunel 2000b) like the omnipresent spontaneous activity (Amit and Brunel 1997) or the persistent, selective delay activity observed in many areas of the cortex in behaving animals (Amit and Brunel 1997; Wang 2001; Yakovlev et al. 1998). In both cases the recorded spike activity is sustained throughout long intervals: low, spontaneous activity is always present, whereas elevated delay activity can last ≤30s (Fuster 1995). During these intervals the statistics of the total synaptic input are likely to be stationary or quasistationary, that is, to vary on time scales that are much longer than the “reaction” time (Gerstner and Kistler 2002) of the assembly of neurons.
Despite the steadiness of the statistics, the total synaptic current results from a considerably irregular synaptic activity, and hence fluctuates all the time. Neurons on the presynaptic side emit spikes spontaneously at a frequency of a few spikes s^{–}^{1} and the neocortical connectivity is rather high [approximately 10^{4} synapses per neuron (Abeles 1991)]. As a consequence, during the interval between two successive spikes, every neuron integrates hundreds of excitatory and inhibitory postsynaptic potentials that arrive at random times. If the spike activities of the presynaptic neurons are statistically independent (see methods) then the resulting total synaptic conductance can be described as a random walk with a Gaussian distribution and finite time correlation length determined by the time development of unitary synaptic events [an Ornstein–Uhlenbeck process (Tuckwell 1988)]. Interestingly in vitro neocortical neurons produce in vivo–like activity when noisy current waveforms are injected (Destexhe et al. 2001; Mainen and Sejnowski 1995), and there is some preliminary indirect evidence, based on the analysis of the distribution of the membrane potential recorded intracellularly in vivo, that somatic currents could be modeled as an Ornstein–Uhlenbeck process (Destexhe and Paré 1999, 2000). The total synaptic current is also approximately Gaussian (Amit and Tsodyks 1992) (also see discussion and the appendix) and, for fast synaptic currents (AMPA or GABA_{A}), the time correlation length is short compared with other inherent time constants of the neuron (e.g., the membrane time constant). As a consequence the total input current is practically white noise and can be fully characterized by its average m_{I} and SD s_{I}.
In the present work we measured the response function of neocortical pyramidal cells to such a noisy current with stationary statistics. We first studied the time development of the neuronal response and characterized the functional role of the adaptation/facilitation components that determine the statistical properties of stationary responses. These components act on different time scales and, for some input currents, reach a steady regime in which they either disappear or are combined together to determine the stationary statistics of the train of spikes generated by the neuron. We then analyzed quantitatively the stationary responses by exploring systematically the whole parameter space {m_{I}, s_{I}} characterizing the input current and determining the resulting spike frequency f. The measured frequencies have been fitted by the theoretical response functions of two simple models of integrateandfire neurons with a single, effective component of adaptation/facilitation. These response functions have been computed analytically (Fusi and Mattia 1999; Ricciardi 1977) and have a relatively simple form. The value of this data modeling is multiple: besides providing a synthetic and efficient way of describing the whole data set, it allows reliable prediction of the output rate in response to currents not used in the experiment and, hence, to cover a wide class of experiments that study the response of the neuron when moving along some specific trajectory in the parameter space of the input current (see e.g., Chance et al. 2002). Moreover, the choice of the response function of model neurons instead of an arbitrary function gives a direct interpretation of the estimated parameters: they are the effective parameters of a model neuron that can recreate the measured response of pyramidal cells. The knowledge of these parameters also allows one to make quantitative predictions about the global dynamics of a network (in vivo) of this kind of cells. In addition, the response function to inputs with stationary statistics also gives important information about the reaction time of the network (Fourcaud and Brunel 2002; Gerstner and Kistler 2002; Mattia and DelGiudice 2002). If the parameters are tuned to reproduce the mean spike frequencies, the simulated neurons can also recreate the higherorder statistics of the interspike intervals expressing, for instance, the degree of irregularity of the spike train (e.g., the coefficient of variability of the interspike intervals). Finally, the statistics of the effective parameters across different cells provide a quantitative estimate of the heterogeneity of the functional properties of the cells.
METHODS
Experimental preparation and recordings
Parasagittal slices of rat somatosensory cortex (300 μm thick) were prepared from 15 to 40dayold female and male Wistar rats according to the institutional guidelines. The preparation was done in icecold extracellular solution using a Campden vibratome (752M; Campden Instruments, Loughborough, UK). Slices were incubated at 35°C for 25 min and afterward left at room temperature until being transferred to the recording chamber. The cells were visualized by infrared differential interference contrast videomicroscopy using a Newvicon camera (C2400, Hamamatsu City, Japan) and an infrared filter (RG9, Schott Mainz, Germany) mounted on an upright microscope (Axioscope FS, Zeiss, Germany).
We recorded in currentclamp whole cell configuration from the soma of layer 5 regular spiking (McCormick et al. 1985) pyramidal cells. Recordings and stimulations were made with an Axoclamp2A amplifier (Axon Instruments, Burlingame, CA) in combination with Clampex 8 (Axon Instruments). The access resistance and the capacitance were compensated using the bridge balance and the capacitance neutralization after having established the whole cell configuration. The data were lowpass filtered at 2.5 kHz with sampling frequency twice the filter frequency. The temperature of the external solution was 31°C. Neurons were visually identified and some of them were filled with biocytin and then stained according to the avidin–biotin–peroxidase (ABC) procedure (Hsu et al. 1981).
Slices were continuously superfused with an artificial cerebrospinal fluid containing (in mM): 125 NaCl, 25 NaHCO_{3}, 2.5 KCl, 1.25 NaH_{2}PO_{4}, 2 CaCl_{2}, 1 MgCl_{2}, 25 glucose, gassed with 95% O_{2}5% CO_{2}.
The pipette solution for most of the analyzed cells contained (in mM): 110 Kgluconate, 30 KCl, 10 EGTA, 10 HEPES, 4 MgATP, 0.3 Na_{2}GTP, 10 Na_{2}phosphocreatine, pH adjusted to 7.3 with KOH. This solution contains a relatively high concentration of EGTA, which allowed long stable recordings and made the cells produce more consistent responses (see Stimulation protocol and observables). We tried to identify undesired artifacts introduced by EGTA by studying the cell's response when two other pipette solutions were used. In what follows we will refer to the solution with high concentration of EGTA as the EGTA pipette solution. The other two pipette solutions were labeled KMeSO_{4} and KGluc and contained (in mM): KMeSO_{4}: 135 Kmethylsulfate, 20 KCl, 0.08 EGTA, 0.045 CaCl_{2}, 10 HEPES, 4 MgATP, 0.3 Na_{2}GTP, 10 Na_{2}phosphocreatine, pH adjusted to 7.3 with KOH. KGluc: 115 Kgluconate, 20 KCl, 10 HEPES, 4 MgATP, 0.3 Na_{2}GTP, 10 Na_{2}phosphocreatine, pH adjusted to 7.3 with KOH. Pipette solution KGluc was always used with 10 mM biocytin. The measured osmolarity of all three pipette solutions was between 310 and 325 mOsm. Because pyramidal cells of the somatosensory cortex in vitro showed virtually no spontaneous activity, we did not systematically block synaptic input mediated by ligandgated channels. Control experiments with blocked synaptic inputs by adding 50 μM DAPV, 10 μM CNQX, and 10 μM bicuculline to the extracellular solution did not change the spike frequency.
The input resistance of the neurons was calculated from the voltage transients in response to at least three different hyperpolarizing (600ms duration, average of the last 300 ms) current pulses (amplitude, 0.05 nA). The membrane time constant τ_{m} was estimated by injecting brief (0.5 ms) hyperpolarizing current pulses (–2.5 nA) into the soma. From the decaying averaged (n = 50) voltage transient after this current pulse, τ_{m} was obtained from the slope of a straight line fitted through the tail portion of the semilogarithmic plot of the membrane voltage against time (Iansek and Redman 1973).
Model of in vivo–like input current
We assume that a large number of presynaptic neurons emit spikes at random times. On the postsynaptic side this heavy barrage is felt as a total synaptic current I that evolves as a random walk. If the synaptic currents are summed linearly and the different inputs are statistically independent (i.e., the emission of a spike by one presynaptic neuron does not affect or only slightly affects the initiation of an action potential in another presynaptic neuron) then the random walk can be replaced by a smoother version, the Ornstein–Uhlenbeck process (Tuckwell 1988), characterized by a Gaussian distribution (mean m_{I}, SD s_{I}) and by a timecorrelation length τ_{I}. If the synaptic evoked potentials sum linearly and there are N_{e} excitatory AMPA receptors, each activated at a mean rate f_{e} and N_{i} inhibitory GABA_{A} receptors (f_{i}), then (1) (2) For simplicity we assumed that the time courses of AMPA and GABA_{A} receptors are the same (τ_{ampa} = τ_{gaba} = τ_{I}). Ī_{e} (Ī_{i}) is the average peak postsynaptic current evoked by the arrival of a single excitatory (inhibitory) spike. The rise time to this peak current is zero, and then the current decays exponentially: I(t) = Ī_{e,i}e^{–}^{t/}^{τ}^{1}. (See the appendix for more details on how these constants are related to the mean synaptic conductances.) Slow NMDA components can be added to m_{I} only, given that they leave s_{I} almost unaffected (Brunel and Wang 2001). The total synaptic current I(t) can be generated with a single equation (3) which is what was used in the experiment ().
Stimulation protocol and observables
The noisy input current was generated as an Ornstein–Uhlenbeck stochastic process by iterating the following expression (4) where ξ_{t} is a unitary Gauss distributed random variable, updated at every time step. The process was generated and injected at a rate of 5 kHz (Δt = 0.2 ms) and the correlation length τ_{I} was usually 1 ms (for 10 cells we also used τ_{I} = 5 ms). The resulting current I(t) has a stationary Gauss distribution with mean m_{I} = μ_{I}τ_{I} and variance (Cox and Miller 1965).
The space {m_{I}, s_{I}} was systematically explored as follows: data points were collected at fixed s_{I} (ranging from 0 to 500 pA), stepwise increasing m_{I} from a subthreshold value up to nonstationary frequencies. This protocol was used to determine the threshold mean current (the rheobase current). Then, the whole space {m_{I}, s_{I}} was discretized and then explored in random order to prevent correlations between time and one of the two parameters m_{I}, s_{I}. In both protocols, the duration of the stimulation depended on the cell response: it was 10 s long, or shorter if ≥150 spikes were collected. For those cells used to characterize the response over time (see results), a stimulus duration of 10 s was always used. The first transient part of the neuronal response (2 s if the stimulation time was longer than 4 s; 0.5 s otherwise) was discarded when estimating the mean spike frequency and the coefficient of variability (see following text). The intervals between successive stimulations were 50–60 s (Fig. 1).
The mean spike frequency (response) was estimated as the ratio between the total number of action potentials N_{sp} and the stimulus duration T. The confidence intervals (68%) of the experimentally measured frequencies were given by Δ = (Δf^{+} + Δf^{–})/2, with (see, e.g., Meyer 1965) (5) The coefficient of variability (CV) of the interspike interval was estimated as the ratio between the SD and the mean of the interspike intervals.
Particular care was taken to ensure that the response of the cell was consistent throughout the whole recording session. Usually the cells showed a transient phase at the beginning, followed by a long time interval (10–90 min) during which the response was consistent (when the same current was injected, the spike frequency did not change within the statistical error of Eq. 5), and by a final unstable phase. Cells with a stable period shorter than 40 min were not included in the analysis. Response consistency depended on the pipette solution. The cells were classified into three groups depending on the level of consistency: 1) consistent cells: less than a third of the repeated recordings were out of range; 2) poorly consistent cells: more than a third of the repeated recordings were out of range; and 3) clearly inconsistent cells: errors as for poorly consistent cells plus there were clear inconsistencies at low frequencies (e.g., the frequency decreased with s_{I} for some points and increased for others). The number of consistent or poorly consistent cells was (out of total number): 41/44 for EGTA pipette solution, 14/19 for KGluc pipette solution, 6/12 for KMeSO_{4} pipette solution.
Model of the integrateandfire neurons
The subthreshold dynamics of the only dynamic variable, the membrane voltage V, obey where C is the membrane capacitance, L(V) is the leakage, and I is the input current that is integrated until V is driven across the threshold θ and initiates an action potential. After the spike emission, V_{r} is reset to some value V_{r} from which the neuron starts again integrating the input current after a refractory period τ_{r}. For the leakage L(V) we studied two models: 1) the classical leaky integrateandfire neuron (LIF), which has a leakage proportional to the membrane depolarization L(V) = –VC/τ (τ is the membrane time constant); and 2) the constant leakage integrateandfire neuron with a floor (CLIFF), which has a constant leakage L(V) = –λ. For this last neuron, to have the same behavior as for the LIF neuron, it is necessary to impose a rigid lower boundary for the depolarization, which we chose to be the cell resting potential (see Response functions of the model neurons).
Adaptation/facilitation components
The stationary response results from the combination of several dynamic processes occurring on different time scales. Some of them are likely to be transient and to disappear in a few hundreds of milliseconds (see results). Others are longlasting and might reach a steady state on time scales of the order of seconds. When the effects of longlasting processes merge together to produce a neuronal response with stationary statistics, we model them by introducing an additional current I_{α} proportional to the mean output spike frequency f. This current is meant to imitate any combination of adaptation and facilitation processes that determine the neuronal response in stationary condition. The linear dependency on f was suggested by the discrepancies observed between the bestfit theoretical response without adaptation/facilitation and the data. Moreover it is supported by the experiments that focus on specific components of frequencydependent modifications of the input current. For example the calciumdependent potassium current, which is responsible for fast adaptation, is usually modeled as a negative current proportional to the intracellular calcium concentration [Ca^{2}^{+}]_{i}, which in turn is proportional to the spike frequency f. A possible implementation in terms of a detailed spike driven dynamics is as follows (see also Ermentrout 1998; Fuhrmann et al. 2002; Liu and Wang 2001 for related models): calcium (or whatever ion species is responsible for the phenomenon; see Powers et al. 1999; SanchezVives et al. 2000) concentration is increased on every spike emission and then decays exponentially to its resting value (6) where the sum extends over all the spikes emitted by the neuron up to time t (t_{k} is the emission time of the kth spike). If the dynamics of [Ca^{2}^{+}]_{i} are slow compared with the interspike intervals (i.e., τ_{Ca} ≫ 1/f), then [Ca^{2}^{+}]_{i} ≃ α_{Ca}f and the current turns out to be proportional to the spike count in some temporal window, I_{α} = –γα_{Ca}f. Any other model that, in stationary conditions, produces a negative current proportional to the mean spike rate would be equivalent.
Theoretical response functions
If the input current I is Gauss distributed and deltacorrelated, then the equation for V can be written as where ξ is a unitary Gaussdistributed variable that is updated every time step. For simplicity we focus on the CLIFF neuron, the same considerations apply to the LIF neuron. μ and σ^{2} are the instantaneous mean and variance of the change in V(t) per unit time and characterize the statistics of V on short time scales (see, e.g., Cox and Miller 1965). For a deltacorrelated input current the variance of V grows linearly with time on short time scales and is proportional to σ^{2}. For such an input current, the response function can be computed analytically for both models (Fusi and Mattia 1999; Ricciardi 1977) and the expressions for the average firing rate f = Φ(m_{I}, s_{I}) is reported in Table 1. The CV of the interspike intervals can also be computed analytically. The expression can be found in Brunel (2000a) for the LIF neuron and in Fusi and Mattia (1999) and Salinas and Sejnowski (2002) for the CLIFF neuron.
The equations of Table 1 express the mean spike frequency f as a function of the mean μ and variance σ^{2} in unit time of the input current (erf = error function). These parameters are related to m_{I}, s_{I} and the time correlation length τ_{I} as indicated in the bottom row: these expressions are based on the assumption that the input current has only fast synaptic components and hence the time correlation length τ_{I} is much shorter than the typical interspike interval. Indeed when the current of Eq. 3 is injected, the variance of V is (Cox and Miller 1965) which in the limit Δt ≫ τ scales as , as expected from a deltacorrelated process. For a more detailed analysis of the effects of time correlated currents on the response function see Brunel and Sergi (1998) and Fourcaud and Brunel (2002). These effects can be reabsorbed in an effective spikeemission threshold θ, which depends on σ_{I}.
The effects of adaptation/facilitation on the stationary response are introduced as an extra current (m_{I} → m_{I} – αf) proportional to the mean spike frequency. For each combination of m_{I}, s_{I}, the mean output frequency determines the amount of negative current that should be added to m_{I} and that modifies the output frequency. The process is iterated until the equation for the frequency becomes selfconsistent and the frequency that appears in the negative current is the same as the output frequency
The asymptotic stationary response function and the response function with α = 0 for the two model neurons are plotted in Fig. 5 and discussed in Response functions of the model neurons.
Data analysis and fitting procedure
The theoretical f – m_{I} curves have been fitted to the stationary data. The fit was achieved through a Monte Carlo minimization (see e.g., Press et al. 1992) of the distance between the measured mean spike frequency and the theoretical spike frequency predicted by the model (7) with respect to the set of parameters Π = {τ_{r}, Vr, C, α}, plus τ_{m} or λ in case of the LIF or CLIFF neuron, respectively. Δ_{k} represents the confidence intervals defined in Stimulation protocol and observables. Because only 5 of the 6 parameters of the model neurons are independent, we set the threshold arbitrarily at θ = 20 mV. In fact, both response functions (see Table 1) are invariant under the scaling θ → ηθ, V_{r} → ηV_{r}, C → C/η, with η > 0.
The minimum of Eq. 7 with respect to the parameters set Π follows approximately a χ^{2}distribution with N – M degrees of freedom, where N is the number of experimental points, and M = 5 is the number of free parameters. The fit was accepted whenever the probability was greater than 0.1 for consistent cells, and greater than 0.01 for poorly consistent cells.
In some cases we also report the absolute discrepancy d, defined as the average (across all points) absolute difference between the measured and the theoretical frequencies of the bestfit curves. This number is not correlated with the goodnessoffit test and gives a useful indication of the error made in considering a theoretical response function which does not strictly pass the leastsquares test.
RESULTS
Time development of the neuronal response
We measured the cell mean spike frequency in response to noisy currents with stationary statistics (see Fig. 1 and methods for details about the protocol). We identified at least three components of adaptation/facilitation, already known in the literature, which determine the features of the stationary response. What follows in this section is not meant to be an extensive and systematic analysis of all the factors that determine the stationary response. The goal is rather to clearly define what we mean by stationary response, to understand when it is possible to have it, and to summarize euristically all the observable mechanisms that affect the transient on different time scales and that might contribute to determine the stationary response. Such stationary response was usually reached in 1–2 s and then sustained until the end of the stimulation. For strong enough injected currents the cells were unable to sustain the elevated activity imposed by the stimulation and no stationary response was possible. The three components are illustrated in Fig. 2. Their features depend on the pipette solution and are summarized in Table 2. They are as follows.
Initial adaptation (Schwindt et al. 1997): fast and invariably present for all the different preparations, it manifests itself for high enough spike frequencies: when the neuron is injected with a constant current, the second and successive interspike intervals are longer than the first one (see Fig. 2, a3–a5). This component, known in the literature to be attributed to calciumdependent potassium current (see e.g., McCormick et al. 1985; Sah 1996), reaches a steady state in a few spikes and it is usually modeled as negative spike frequency–dependent current, which clearly affects the stationary response. This adaptation component was present for all the cells and all the preparations and in particular it was observed also in the strong presence of EGTA, a calcium chelator. Buffering of calcium concentration transients by EGTA is probably not fast enough to disrupt fast initial adaptation (Smith et al. 1984). The spikefrequency dependency was different depending on the pipette solution. The same effect was produced for lower currents in the case of EGTA pipette solution when compared with the case of KGluc pipette solution: 1.2 ± 0.4 nA (KGluc pipette solution), 0.6 ± 0.2 nA (EGTA pipette solution), corresponding to 36 ± 10 Hz (EGTA pipette solution) and 21 ± 8 Hz (KGluc pipette solution). For KMeSO_{4} pipette solution: 0.8 ± 0.2 nA and 24 ± 7 Hz.
Early facilitation: on a time scale of 1–2 s, the mean spike frequency slowly increases. This form of relatively slow facilitation was invariably present for the different preparations and shown by all of the recorded cells and for all the stimulation strengths (see Fig. 2). It is known in the literature to be attributed to calcium accumulation (Powers et al. 1999) and its effects on the stationary response compete with those of the late adaptation component (described below) when the spike frequency is high enough. For this reason it is difficult to isolate its effects and to characterize its dependency on the preparation. However, it was clear that early facilitation was mostly prominent when KMeSO_{4} pipette solution was used and this would be compatible with the known fact (Zhang et al. 1994) that methylsulfate is the least disruptive to intracellular structures of calcium homeostasis.
Late adaptation: for strong enough input currents the cells were unable to sustain the elevated activity imposed by the stimulation. After 2–3 s, the mean spike frequency estimated on a sliding temporal window was constantly decaying. The threshold current for detecting late adaptation and the rate of frequency decay depended on the preparation. In particular, the same effect (a frequency decay of more than 0.5 spikes/s) required relatively low currents in the case of EGTA pipette solution (0.8 ± 0.3 nA, corresponding on average to approximately 35 Hz for the 12 out of 13 cells that showed the phenomenon), higher for KGluc pipette solution (1.6 ± 0.5 nA, 52 Hz, all 19 cells), and much higher for KMeSO_{4} pipette solution (1.8 nA and 60 Hz) for the few (3/12) cells that showed late adaptation. This kind of slow adaptation has already been studied (Fleidervish et al. 1996; Powers et al. 1999; SanchezVives et al. 2000; Sawczuk et al. 1997; Schwindt et al. 1989) and it is hypothesized to be attributed to slow inactivation of sodium channels (Fleidervish et al. 1996; Powers et al. 1999; Sawczuk et al. 1997) and to an outward Na^{+}dependent K^{+} current (SanchezVives et al. 2000; Schwindt et al. 1989). The spike frequency decay was always accompanied by a progressive drift of the maximal upstroke velocity of the membrane potential, indicating that the spike shape was continuously and slowly broadening. This phenomenon is best displayed by using longer stimuli that elicit high initial spike rates. An example is shown In Fig. 3B), where 60slong stimuli were used.
Many of the recorded cells (19/44) showed an initial doublet of spikes, resembling the beginning of a burst. This doublet sometimes masks the initial adaptation component and makes it difficult to analyze. Notice that all the recorded cells were selected not to be inherently bursting. The doublet was observed for all the preparations for strong enough stimulation (usually stronger than the one that reveals fast spike adaptation) and it is shown in Fig. 2, a4–a5. The threshold current depended on the pipette solution: 1.2 nA for the 4/13 in EGTA pipette solution, about 1.8 nA for 13/19 cells in KGluc pipette solution and for 2/12 cells in KMeSO_{4} pipette solution. Because of its transitory nature, the mechanism responsible for this doublet probably does not affect the stationary response.
Maximal stationary response
The late adaptation component makes the cell unable to sustain a response with stationary statistics above a critical output spike frequency, which depended on the cell and on the pipette solution. Although the observed modifications in the spike mean frequency and in the spike shape were usually quite substantial at the end of long stimulations, all cells could recover during the interstimulus intervals, indicating that the phenomenon is reversible.
For most cells we did not stimulate for very long intervals (≤60 s) as for the cell shown in Fig. 3. Hence we do not have enough data to determine the maximal stationary response for all cells. This analysis would have required the injection of several currents, strong enough to produce nonstationary responses, and it would have limited even further the number of stationary data points. Nonstationary responses are not included in the following analysis and are not modeled.
Usually we restricted the input currents to the range in which the slow frequency decay was below 1 Hz/s. This criterion restricted our analysis of the response function to a limited range of frequencies that depended on the pipette solution. The maximal frequency of this range could be determined for those cells for which the stimulation was strong enough to reveal a frequency decay above 1 Hz/s and it was (see Table 2): 44 ± 12 Hz (n = 11/13) for EGTA pipette solution and 56 ± 9 Hz (n = 5/19) for KGluc pipette solution. In KMeSO_{4} pipette solution the cells did not display maximal stationary response with the exception of only one out of 12, in which case the maximal stationary frequency was 64 Hz.
Experimentally measured response functions
In Fig. 4 we show the response function as measured in the experiment. The measurements are shown in the form of f – m_{I} curves that represent the output frequency f (the response) as a function of the mean current m_{I} at constant SD s_{I.} The CV of the interspike interval is also reported in the same form. The qualitative behavior was the same for all cells: for constant currents (s_{I} = 0), there is obviously no activity for average inputs that are below the rheobase current (i.e., the minimal constant current that makes the neuron fire), whereas the response curve is approximately linear for suprarheobase currents. For noisy inputs (s_{I} > 0) there are essentially two regimes. 1) A subrheobase, fluctuationdominated regime in which the mean current m_{I} is below the rheobase current and hence not sufficient to drive the membrane voltage across the threshold for emitting a spike. In this case the action potentials are sporadically initiated by fluctuations and hence the spike activity is very irregular and the CV is high (Fig. 4C). 2) A suprarheobase, driftdominated regime, in which the membrane potential fluctuates around a rising ramp that leads to the emission threshold at a more regular pace (Fig. 4B). The introduction of noise permits subrheobase activity, which in turn smooths out the response curves at the rheobase. The curves, convex for weak input currents, sometimes bend at high frequencies (see e.g., Fig. 6, B and C) and hence change curvature.
Response functions of the model neurons
The data points have been fitted by the response functions of two simple model neurons sharing a similar qualitative behavior (see methods). The two models differ in the form of the leakage L(V): for the first, classical leaky integrateandfire (LIF) neuron the leakage is proportional to V [L(V) = –VC/τ_{m}, where τ_{m} is the membrane time constant], whereas for the second model—the CLIFF neuron—the leakage is constant [L(V) = –λ] and the input current is integrated linearly. To obtain two qualitatively similar response functions for the two models, it is essential to limit the membrane potential from below in the case of the CLIFF neuron (Fusi and Mattia 1999; Salinas and Sejnowski 2002). Indeed, when both neurons are injected with a subrheobase noisy current, the variance of the fluctuations of the membrane voltage tends to increase linearly with time. However, for the LIF neuron the leakage compensates for this tendency, and, depending on the parameters of the current, can lead to a stationary subrheobase regime in which V fluctuates around an asymptotic value. For the CLIFF neuron this is not possible because the current is integrated linearly and all the fluctuations accumulate: for negative net currents (mean synaptic current minus leakage) the mean voltage decreases and the fluctuations increase boundlessly. In such a subrheobase regime the average spike frequency of the CLIFF neuron is always zero (Fusi and Mattia 1999; Gerstein and Mandelbrot 1964). The rigid barrier that limits the membrane voltage from below permits the neuron to fluctuate steadily around a value near the lower barrier, waiting for a fluctuation that drives V across the threshold. In this case a subrheobase regime with nonzero frequency is possible.
In both neurons the adaptation/facilitation components which determine the asymptotic stationary response were modeled by assuming that the effective mean current driving the cell is reduced by a term proportional to the cell's own spike frequency (m_{I} → m_{I} – αf; see methods). There are several simple and detailed biophysical models that correspond to this phenomenological model: the simple models based on calcium dependent potassium channels (see e.g., Wang 1998) and the more realistic Hodgkin–Huxley type model in Traub and Miles (1991) are but two examples (see also Ermentrout 1998). All these models exploit the fact that the adaptation/facilitation processes are slow compared with the average interspike intervals. Note that these mechanisms do not account for the nonstationarity observed in Fig. 3B, given that our model neuron can always have a stationary response, no matter what the intensity of stimulation.
The response functions corresponding to the two models are plotted in Fig. 5 in the same way as the measured curves are shown in Fig. 4. The main difference between the two model neurons is in the dependency of the spike frequency on s_{I}, the fluctuations amplitude of the noisy current: the CLIFF neuron is more sensitive to large s_{I} and this fact manifests itself in a progressively increasing distance between the f – m_{I} curves. This behavior can be exposed by plotting the spike frequency as a function of s_{I} at constant m_{I} (insets in Fig. 5): the curvature has a different sign for the two neurons. Moreover the response function of the LIF neuron without adaptation/facilitation components is highly nonlinear around the rheobase current for low levels of noise and, for s_{I} = 0, the derivative of the f – m_{I} curve is infinite. These nonlinearities are attenuated by the frequencydependent feedback (Ermentrout 1998).
The CV of the interspike intervals behave in the same way for the two model neurons: in the subrheobase fluctuations dominated regime the CV is close to 1 for a wide range of input currents. This corresponds to the maximal irregularity and the train of spikes has Poisson statistics. As the current drives the neuron toward a driftdominated regime, the CV approaches 0 with a speed that depends on s_{I} (i.e., the fluctuations in the input current).
Fitting the theoretical response functions to the data
Both models could faithfully reproduce most of the cells stationary responses, at least for those neurons that had consistent enough responses to allow for a quantitative analysis (see methods). For each cell we fitted the theoretical response functions to the data points by searching the space of the 5 independent parameters characterizing each model: the capacitance C, the refractory periodτ_{r}, the reset potential V_{r}, the adaptation/facilitation constant α, and the time membrane constant τ_{m} for the LIF neuron and C, τ_{r}, V_{r}, α, and the leakage λ for the CLIFF neuron (see methods). We adopted a Monte Carlo technique that minimizes the mean square distance between the data points and the predicted values, each point being weighted by the inverse of the confidence interval. The χ^{2} test passed either with one model or with the other (or with both) with a probability P > 0.1 for 27 consistent cells out of 37. For all these cells, the EGTA pipette solution was used. For other pipette solutions the consistency was usually poor and the χ^{2} test passed with a probability P > 0.01 for 16 cells out of 24.^{1} The responses could anyway be well approximated by the theoretical functions even when the models did not pass the test: the average discrepancy between the measured and the model frequencies were <2.5 Hz; <1.5 Hz if frequencies below 50 Hz were considered. That is, the average difference between a single data point and its theoretical match was below 3 Hz even when the test did not pass. Examples of fitted response functions are shown in Fig. 6. The LIF model was best suited to reproduce f – m_{I} curves that were almost equally spaced (Fig. 6B), whereas the CLIFF model could capture better those cells that were less sensitive to low levels of noise as in Fig. 6A.
Adaptation/facilitation components of the model dynamics turned out to be essential to fit the response of the models to the data. For the LIF neuron it is necessary to linearize the response function at the rheobase. For the CLIFF neuron, which already has an almost linear f – m_{I} curve at the rheobase, one might wonder if the fitting would have been possible without adaptation/facilitation components. Increasing the membrane capacity C as well as increasing α reduces the slope of the f – m_{I} curves. However, in contrast to adaptation/facilitation, an increased capacity also reduces the effect of fluctuations on the spike frequency and hence diminishes the sensitivity to s_{I}. Thus for the CLIFF neuron, adaptation/facilitation are the only mechanisms that can change the slope of the f – m_{I} curves almost without affecting the distance between f – m_{I} curves that correspond to different values of s_{I} (La Camera et al. 2002).
For high frequencies (>40–50 Hz) there is a preliminary evidence of a phenomenon that is not captured by our model neurons: the f – m_{I} curves corresponding to different s_{I} have a slight tendency to diverge, as if the sensitivity to s_{I} would increase for strong enough stimulation (see, e.g., Fig. 6C). This is not captured by our models of integrateandfire neurons for which asymptotically the f – m_{I} curves always tend to merge. This divergence effect is a source of discrepancy that usually does not compromise the χ^{2} test, but that clearly indicates the activation of an extra process that is not modeled. The increased sensitivity to s_{I} is usually accompanied by an increase in the CV of the interspike intervals (see the previous section) and it might be attributable to the activation of calcium spikes in the distal dendrite (Larkum et al. 1999). The study of this phenomenon will require further investigation.
Some cells (n = 6, all of them in EGTA pipette solution) could be fitted only by searching a region of limit parameter values, letting the reset potential be very close to the threshold. This was the only way to capture sets of f – m_{I} curves that do not tend to converge to a common asymptotic value for elevated currents (Fig. 6D). Although the response function could be described by the theoretical response function, we do not consider the model as appropriate for describing these cells. The theoretical response functions are computed under the assumption that the current is deltacorrelated, and this is not a good approximation when the reset potential is very close to the threshold. Hence we consider these cells belong to a different class that the simple integrateandfire models are probably inappropriate to describe.
We recorded the response of 10 cells to noisy currents with a longer time correlation length (τ_{I} = 5 ms). The f – m_{I} curves were qualitatively the same as in the case of shorter τ_{I} and they could be fitted by the theoretical response functions of Table 1 (not shown), even if they are valid only for deltacorrelated currents. A full account of the analysis of these cells will be reported elsewhere.
Higherorder statistics of the interspike intervals
The main goal of this investigation was to study the neuronal response in terms of mean spike rates. However, it is also interesting to consider the higherorder statistics of the distribution of the interspike intervals. A detailed and systematic analysis would require more data to estimate the higherorder moments of the distribution, although our preliminary analysis of the variance of the interspike intervals already gives interesting indications. In Fig. 7 we plotted, for two cells, the CV of the interspike intervals as a function of the spike frequency. Each curve is obtained by keeping s_{I} fixed and then sweeping along m_{I}. The parameters of the model are tuned to fit only the mean firing frequencies. Then, with the bestfit parameters we computed the CV by means of simulations and compared it to the measured one.
The cell in the top row of Fig. 7 represents the typical case, whereas the one below is our best fit. The agreement between the predictions and the data points is in general reasonably good for most of the points and for most of the cells, even if it would probably not pass a χ^{2} test. The largest discrepancies are observed at high firing frequencies, for elevated values of the variance of the input current. There are at least two reasons for these discrepancies. First, these points might not be completely stationary, given that the firing frequency approaches the maximal allowed value for the cell. As expected, this produces a higher CV than predicted by the theoretical models for which the nonstationarity is not modeled. Second, at high frequencies and high s_{I} we often observed doublets or bursts of spikes (see the bottom left inset in Fig. 7). This is also not captured by the model and increases the measured CV.
Moreover there is often a large discrepancy around the rheobase current, when the amount of noise in the current is small. Near the threshold for spike emission, the response sometimes consisted of sequences of tonic, regular firing interrupted by periods of subthreshold activity (see the bottom right inset in Fig. 7). This behavior, not captured by a simple integrateandfire model, artificially increases the variability of a train of spikes that otherwise would be rather regular.
Effective parameters
The parameters corresponding to the best fit must be considered as effective parameters, that is, as those parameters that provide the best description of the stationary response function. Other observables might not be captured by the same parameters (also see discussion). Fitting the f – m_{I} curves for three different values of s_{I} was already restricting the model parameters to a small region of the parameter space (we had 3 to 5 f – m_{I} curves for each cell). The range in which single parameters can vary and the χ^{2} test still passes depends on the sensitivity of the spike frequency on the different parameters. A rough estimate indicates that the least determined parameters are τ_{r} (up to about ±70% error) and V_{r} (about ±25% error), whereas the other parameters can vary at most in intervals of the order of ±20% (α, CLIFF neuron), ±5% (C, α, τ, LIF neuron; C, λ, CLIFF neuron). The refractory period is clearly the effective parameter to which the response function is least sensitive because it affects mostly the very highfrequency region (f ∼ 1/τ_{r}), usually beyond the observed range of output rates.
The effective parameters for the two model neurons are reported in Table 3 for all the cells that were consistent or poorly consistent and that did pass the χ^{2} test (the number of cells is reported in column N). For the analysis of effective parameters we mainly focused on the cells with EGTA pipette solution because they were stable enough to give consistent responses for a large number of data points. The parameters are reported as an average across all cells with EGTA pipette solution and separately for the other two pipette solutions. Interestingly, the differences between different pipette solutions were usually comparable to the differences across cells (reported as SD).
Two passive parameters characterizing the cell (τ_{m}, C) could be measured directly (see methods), and their values were in general different from the effective parameters. The crosscorrelations were negligible for the CLIFF neuron and weak for the LIF neuron: the correlation between the directly measured capacitance and the effective capacitance of the LIF neuron was 0.59 (0.03 for the CLIFF neuron), and the correlation between the membrane time constants was 0.30. Because real neurons are not integrateandfire neurons, such a mismatch is not very surprising. Indeed, the effective parameters provide a good fit to the mean spike frequencies, while the passive membrane parameters are estimated from very different observables (the subthreshold response to short pulses), in different experimental conditions. No clear patterns between effective or directly measured parameters and the model that could fit the data emerged.
DISCUSSION
The most prominent result of the present work is that simple integrateandfire model neurons can faithfully recreate the response of neocortical pyramidal cells to in vivo–like currents. We also gave a full account of the dependency of the cell response on the fluctuations of the input current (s_{I}), which has usually been overlooked in previous studies of neural response properties. This dependency was recently shown to play an important role when the neuron works in a subrheobase regime (Amit and Brunel 1997; Amit and Tsodyks 1991; Chance et al. 2002; Fusi and Mattia 1999) or when a network of interacting neurons has to respond fast to a time varying input (Rudolph and Destexhe 2001; Silberberg, Bethge, Markram, Tsodyks, and Pawelzik, unpublished observations, 2002).
The agreement between the theoretical and the measured response function is remarkably good for a wide range of statistics of the input currents, on the whole {m_{I}, s_{I}} space and for all the sustainable spike frequencies of the analyzed cells. Hence, our results cover a variety of physiological scenarios that can be studied by measuring the neural response when moving along specific trajectories in the {m_{I}, s_{I}} space. For instance balanced synaptic input, as studied e.g., in Chance et al. 2002, would correspond to increasing the variability of the current s_{I}, while keeping the mean current m_{I} fixed, and the outcome could be predicted by studying the theoretical response function of integrateandfire neurons.
Moreover the knowledge of the response functions to in vivo–like currents allows one to study many dynamic properties of networks of interacting neurons and to relate single neuron properties to the activity observed in in vivo experiments. For instance, imposing the stability of a network state, like spontaneous activity (Amit and Brunel 1997), can constrain the possible architecture of the network, thus providing an indirect measurement of the synaptic efficacies in vivo. This will require a more extensive analysis of different types of cells, in different layers and different areas. The analysis of the dynamic states of a network of neurons will also require an additional study of the dependency of the parameters m_{I}, s_{I} of the current on the mean output rate of the presynaptic neurons. This dependency can be quite complicated, involving for instance the shortterm dynamics of single synaptic responses (Tsodyks et al. 1998), the effects of a particular morphological structure (Rhodes 1999), or calcium dynamics (Larkum et al. 1999). So far the analysis of the in vivo phenomena based on single neuron response functions has been successfully performed with a simple linear relation between both m and and the spike frequency of the presynaptic neurons (Amit and Brunel 1997; Brunel and Wang 2001; Fusi and Mattia 1999; Yakovlev et al. 1998).
Besides these additional studies, many novel issues are raised by this work and remain open for further investigation. Here we focused on a specific observable, the mean spike frequency, to conclude that pyramidal cells respond very much like integrateandfire neurons. This is just a starting point and there are other aspects and observables of the neural activity to be investigated. The analysis of the CV that we present here provides a preliminary indication that integrateandfire neurons can capture more than the mean spike frequency. Because the fluctuations play an important role in determining the exact times of the initiation of the action potential (Mainen and Sejnowski 1995), it would be interesting to see how well the integrateandfire neuron (or related models) can describe the exact spike times under noisy current injection (see Kistler et al. 1997). This issue was beyond the scope of this report, but is currently under investigation (R. Jolivet, A. Rauch, H. R. Lüscher, and W. Gerstner, unpublished observations; M. Larkum, W. Senn, and H. R. Lüscher, unpublished observations).
A more systematic analysis of higherorder statistics of the spike trains and to what extent the exact timing of the spikes can be reproduced by the same simple integrateandfire neurons will be studied elsewhere and will probably require the introductions of new elements in the neuron model. Besides, such a study will certainly further constrain the effective parameters (e.g., the reset potential V_{r} can be determined more precisely by fitting both the mean frequency and the CV). So far, both the LIF and CLIFF neuron models with adaptation/facilitation are good enough to describe the response function of pyramidal cells, at least for the wide, biologically plausible range of input currents that was investigated in our work. Other observables might better expose the differences between the two models.
The response of the cell under strong persistent stimulation has not been modeled here. When neurons cannot sustain the elevated rate to which they would be driven by the input current, it is not possible to define a stationary response function. This behavior has to be incorporated in the model that, so far, can reliably describe the response function of the neuron only when sustained spike frequencies are below the threshold frequencies indicated at the end of Maximal stationary response (30–80 Hz depending on the pipette solution used). Above these frequencies some other mechanism beside fast or slow spike frequency adaptation should be invoked to capture the effective reduction of the spike frequency. The activation of such a mechanism cannot rely only on the number of emitted spikes because this would not account for the behavior illustrated in Fig. 3b, where the cell eventually stops firing. Therefore it is more likely to depend on some dynamic variable related to the average depolarization of the cell. One possibility is the inactivation of sodium channels, which has already been identified as responsible for late adaptation in other works (Fleidervish et al. 1996; Powers et al. 1999; Sawczuk et al. 1997).
An important element of realism that was not considered so far concerns the way synaptic conductances are imitated by the injection of a current into the soma. Generating real synaptic inputs by activating presynaptic neurons is the most realistic (and challenging) way of providing in vivo–like inputs to the neurons. However, it is interesting to study intermediate steps toward a realistic situation. Many investigators (see, e.g., Destexhe et al. 2001 and references therein) imitate synaptic conductances by generating a somatic current that is a product of a Gaussdistributed conductance multiplied by a voltagedependent driving force (the instantaneous membrane voltage minus the reversal potential). We give in the appendix a detailed description of this procedure. This kind of time varying current can be generated in an experiment by using dynamic clamp techniques (see, e.g., Harsch and Robinson 2000). The distribution of the injected current is not distorted much by the voltagedependent driving force, and thus the final total current can still be reasonably approximated by an Ornstein–Uhlenbeck process (see also Amit and Tsodyks 1992). However, to study the effective response function of the neuron, the correct scaling factors for both the mean and the variance of the current should be introduced. These factors depend on the average membrane voltage, which in turn depends on the parameters of the current. As a consequence, the scaling factors in principle might not be the same for all points of the f – m_{I} curves whose shape might be qualitatively different when conductances are considered (Tiesinga et al. 2000). Moreover, conductancebased inputs shorten the effective time correlation length of the membrane depolarization (see, e.g., Paré et al. 1998), which in turn modifies the effective threshold for emitting a spike (Brunel and Sergi 1998). However, the analysis presented in the appendix shows that even when a unique scaling factor for all the input currents is chosen, the distortions of the stationary response functions are negligible and the neuron responds in a similar way to either current injection or conductance drive. There is a preliminary evidence that also transient responses are not substantially modified when the voltagedependent driving force is considered: the differences are quantitative and not qualitative (Fourcaud and Brunel 2002). As a consequence, the effects of the voltagedependent driving force can be compensated by the introduction of some constant factors that multiply the total excitatory and inhibitory currents. These factors might change dramatically the relation between the frequencies of the presynaptic neurons and the statistics of the input current, but they do not affect the results presented in this work, which are more related to the dependency of the output spike frequency on the statistics of the input current.
APPENDIX
Conductance versus current injection
The response function studied in our experiments allows relating the output spike frequency to the statistics of the input current— here described by the pair {m_{I}, s_{I}}—that in turn depends on the mean frequency of the input synaptic events and on their efficacy (see Eqs. 1 and 2). In particular, if the inputs can be grouped into two classes (excitatory and inhibitory), then the pair depends on a quadruplet of parameters: {N_{e}f_{e}, Ī_{e}, N_{i}f_{i}, Ī_{i}} (see methods for the definition of the symbols). Given the response function f = Φ(m_{I}, s_{I}) and the relation between {m_{I}, s_{I}} and the quadruplet mentioned above, one can study the dynamic properties of an arbitrary number of homogeneous populations of neurons. For instance the stability of a single recurrent network of excitatory cells can be analyzed by imposing that the input spike frequency f matches the output spike rate f of any of the cell: f = Φ([m_{I}(f), s_{I}(f)]. In general, what really matters is the relation between the output spike frequency and the input spike frequency of all the presynaptic populations of neurons.
To have a more realistic synaptic drive, one may inject into the soma a current I_{dc}, which is a product of a Gaussdistributed conductance multiplied by a driving force that depends on the instantaneous membrane depolarization (Chance et al. 2002; Destexhe et al. 2001; Harsch and Robinson 2000) (A1) where ḡ_{e}, ḡ_{i} (nS) are the excitatory and inhibitory peak conductances, respectively; V_{e} and V_{i} are the reversal potentials for AMPA and GABA_{A} receptors (70 and –10 mV, respectively); x_{e}_{,}_{i}(t) are Gaussdistributed variables with average μ_{xe,xi}, variance and time correlation length of τ_{x} = 1 ms for both. x_{e,i}(t) can be thought of as the probability of ion channel openings (normalized from zero to ḡ_{e,i}) attributed to a large number of independent excitatory and inhibitory postsynaptic potentials, each characterized by a sharp rise and an exponential decay with a time constant τ_{x}. The current that then actually flows into the cell depends on the depolarization V(t) as expressed in Eq. A1 for I_{dc}. The average and the variance of the AMPA conductance drive are a function of the quadruplet of parameters {N_{e}f_{e},ḡ_{e},N_{i}f_{i},ḡ_{i}} Analogous formulas were used for the inhibitory input (see also methods).
In the case of conductance drive, what is relevant for exploring the collective behavior of a population of neurons is the relation between the output spike frequency and the quadruplet of parameters {N_{e}f_{e},ḡ_{e},N_{i}f_{i},ḡ_{i}}. Each quadruplet defines a possible statistics of the input and the distribution of the current can be determined by running a simulation of a LIF neuron in which the instantaneous current is computed at each time step according to Eq. A1. Given the average m̃ and the SD s̃_{I} of the input current I_{dc}, it is not usually possible to predict correctly the output spike frequency because: 1) the distribution of the current is slightly distorted by the driving force, which depends on the depolarization and the Gaussian approximation is no longer valid; 2) the temporal correlations introduced by the driving force change the temporal statistics of the input current in a way that depends on the neuronal activity (see, e.g., Tiesinga et al. 2000). This means that if one injects a Gaussian current characterized by m_{I} = m̃_{I} and s_{I} = s̃_{I}, in general one would not get the same spike rate as one would in the full simulation with the conductance drive I_{dc}.
However, for each quadruplet {N_{e}f_{e},ḡ_{e},N_{i},f_{i},ḡ_{i}} it is possible to generate a Gaussian current characterized by a pair {m_{I}, s_{I}} that, when injected into a neuron, produces the same spike frequency as in the case of the conductance drive. The mean m_{I} and the variance s^{2}_{I} of that current are given by (A2) where V_{E,I,eff} and are four positive scaling parameters in units of mV and mV^{2}, respectively, which do not depend on the statistics of the current determined by {N_{e}f_{e},ḡ_{e},N_{i}f_{i},ḡ_{i}}, but only on the parameters that characterize the neuronal dynamics. This means that the scaling factors are unique for a wide range of inputs and that it is possible to predict the spike frequency in the case of a conductance drive for each quadruplet {N_{e}f_{e},ḡ_{e},N_{i}f_{i},ḡ_{i}}. The recipe is simple: given the quadruplet one can determine m_{I} and s_{I} according to Eq. A2 and then use the response function measured in the experiment to determine the output firing rate.
To prove that this is possible for a wide range of input currents we explored extensively the space of quadruplets {N_{e}f_{e},ḡ_{e},N_{i}f_{i},ḡ_{i}} and compared the output spike frequency in case of conductance drive and in case of current drive for the same simulated neuron. The neuron's parameters were: τ_{r} = 0, C = 500 pF, θ = 20 mV, V_{r} = 10 mV, τ = 20 ms, α = 0, and resting potential V_{rest} = 0. Our strategy to explore the input space was to keep constant ḡ_{e,i} (nS) on each curve in Fig. A1 [values reported in the plot; the values for the rightmost curve correspond roughly to those used for AMPA and GABA_{A} conductances in Harsch and Robinson (2000)], and then sweeping along the diagonal of the {N_{e}f_{e}, N_{i}f_{i}} plane (i.e., N_{e}f_{e} = N_{i}f_{i}, frequency ranges reported in the plot). For each quadruplet {N_{e}f_{e},ḡ_{e},N_{i}f_{i},ḡ_{i}} we computed the spike frequency in the case of conductance drive [the simulated neuron was injected with I_{dc}(t)] and in the case of the injection of a Gaussian current with m_{I}, s_{I} given by Eq. A2. Note that s_{I} is not constant along each curve contrary to the plots shown in Fig. 6.
We show in Fig. A1 that it is possible to tune the four scaling parameters V_{E}_{,}_{I}_{,eff} and U_{E}_{,}_{I}_{,eff} in such a way that the neuron responds with the same mean spike frequency in the dynamicclamp and in the currentclamp modalities. The unique scaling factors for all the inputs were: V_{E}_{,eff} = 51.1 mV, V_{I}_{,eff} = 27.1 mV, U_{E}_{,eff} = 24.7 mV, and U_{I}_{,eff} = 30.2 mV.
DISCLOSURES
This study was supported by Swiss National Science Foundation Grants 3161335.00 and 3152065234.01 and by Silva Casa Foundation.
Acknowledgments
We are grateful to N. Buchs for having written the software for stimulation and data acquisition. We also thank T. Berger, N. Brunel, P. Del Giudice, A. Destexhe, W. Gerstner, and M. Mattia for stimulating discussions and for many useful remarks on the manuscript and G. Silberberg for helpful discussions about the experimental methods.
Footnotes

↵1 For poorly consistent cells the purely statistical error underestimates the error, and a lower threshold for P was chosen. This threshold corresponds in average to multiplying the statistical error by a factor up to 1.25 at most.
The costs of publication of this article were defrayed in part by the payment of page charges. The article must therefore be hereby marked “advertisement” in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.
↵* A. Rauch and G. La Camera contributed equally to this work.
 Copyright © 2003 by the American Physiological Society