Abstract
In vivo cortical neurons are known to exhibit highly irregular spike patterns. Because the intervals between successive spikes fluctuate greatly, irregular neuronal firing makes it difficult to estimate instantaneous firing rates accurately. If, however, the irregularity of spike timing is decoupled from rate modulations, the estimate of firing rate can be improved. Here, we introduce a novel coding scheme to make the firing irregularity orthogonal to the firing rate in information representation. The scheme is valid if an interspike interval distribution can be well fitted by the gamma distribution and the firing irregularity is constant over time. We investigated in a computational model whether fluctuating external inputs may generate gamma processlike spike outputs, and whether the two quantities are actually decoupled. Wholecell patchclamp recordings of cortical neurons were performed to confirm the predictions of the model. The output spikes were well fitted by the gamma distribution. The firing irregularity remained approximately constant regardless of the firing rate when we injected a balanced input, in which excitatory and inhibitory synapses are activated concurrently while keeping their conductance ratio fixed. The degree of irregular firing depended on the effective reversal potential set by the balance between excitation and inhibition. In contrast, when we modulated conductances out of balance, the irregularity varied with the firing rate. These results indicate that the balanced input may improve the efficiency of neural coding by clamping the firing irregularity of cortical neurons. We demonstrate how this novel coding scheme facilitates stimulus decoding.
 firing irregularity
 neural code
 balanced synaptic input
 brainmachine interface
 gamma distribution
 information geometry
Introduction
In vivo cortical neurons exhibit highly irregular firing (Softky and Koch, 1993; Holt et al., 1996; Shadlen and Newsome, 1998), which is often represented as Poissonlike spike trains. Background synaptic input with balanced excitation and inhibition has been suggested to generate the irregular firing, although other mechanisms are also possible (Reyes, 2003; Durstewitz and Gabriel, 2007). Because firing rate is crucial for the information coding by neurons, an observer can correctly interpret the messages of a neuron only with an accurate estimation of the firing rate. Here, the observer can be either external or internal to the brain. However, the estimation is difficult when spike sequence is highly irregular. It is, therefore, unclear whether the irregularity of neuronal firing plays an active role in cortical information processing (Shadlen and Newsome, 1994, 1995; Softky, 1995; Destexhe et al., 2003).
Estimating the firing rate will be easier if irregular neuronal activity can be characterized by a measure that is sensitive only to the firing irregularity, but not to the rate modulation. The coefficient of variation (C_{V}) of interspike intervals (ISIs), a widely used measure for irregular neuronal firing (Softky and Koch, 1993), does not have such a property (Stevens and Zador, 1998; Sakai et al., 1999; Shinomoto et al., 1999). To obtain an alternative measure, Miura et al. (2006a,b) studied a decomposition of fluctuating neuronal firing into firing rate and irregularity when spike generation obeys the gamma process, which is a natural extension of the Poisson process. The gamma process has two parameters: the mean rate ξ and the “shape parameter” κ, which represents how irregular the process is. A formula derived therein to estimate the value of κ from an observed sequence of ISIs does not depend on ξ. In other words, the ξcoordinate is orthogonal to the κ coordinate in the space of information representation.
In this study, we use this κ to characterize the irregular spikes. We attempt to examine whether the brain may actually use the coding scheme based on the gamma process, which facilitates estimation of the parameters by their orthogonalization. The constancy of the firing irregularity can be the key evidence because it is essential for the coding scheme. We use Neuron simulator and the wholecell patchclamp recording technique to inject a fluctuating input current into model and real neurons, respectively. We show that, if excitatory and inhibitory synaptic inputs are covaried in a balanced manner, output spikes are well described by the gamma process and the firing irregularity varies only moderately over time regardless of rate modulations. The degree of the irregularity remains unchanged unless the balance between the excitatory and inhibitory synaptic inputs is changed. We demonstrate that this constancy of the firing irregularity improves the estimation of input firing rate. Accumulating evidence suggests that synaptic inputs are balanced in cerebral cortex (Shu et al., 2003; Haider et al., 2006) and spinal cord (Berg et al., 2007). Our results achieve a novel insight into possible contributions of balanced excitation and inhibition to neural code.
Materials and Methods
Neural code using the gamma process: decoupling of firing rate and irregularity.
We briefly summarize the decomposition scheme proposed previously for spike trains obeying the gamma process (Miura et al., 2006a,b). The scheme enables us to estimate the shape parameter κ that is related to the irregularity of neuronal firing, which generates a series of ISIs according to the following gamma distribution (Cox and Lewis, 1966; Baker and Lemon, 2000; Casella and Berger, 2002; Fellous et al., 2003): where T represents the ISI, ξ denotes the firing rate, and Γ(κ) is the gamma function. The gamma process is a natural extension of the conventional Poisson process, with the shape parameter κ measuring the irregularity of spike trains. When κ = 1, the ISI distribution becomes an exponential distribution that is equivalent to a Poisson process in which the spike train is completely irregular. When κ is large, the gamma distribution can be approximated by a normal distribution, with variance decreasing with increasing κ. In the limit of infinitely large κ, ISIs take only a single value. Thus, κ characterizes the irregularity of spike trains.
Our task is to estimate the sequence of instantaneous firing rates {ξ_{1}, ξ_{2},…, ξ_{N}} that likely produces the observed sequence of ISIs {T_{1}, T_{2},…, T_{N}}. It is in general very difficult to estimate the sequence of firing rates and the shape parameter κ simultaneously. If, however, κ is constant over time, we can use the method of estimating functions (Amari and Kawanabe, 1997) to estimate the value of κ from the observed ISIs without knowing the timedependent firing rates. Let the firing rate ξ be distributed according to k(ξ) that is unknown to the observer. Then, the probability distribution of ISIs is given by the gamma distribution q(T; ξ, κ) weighted by k(ξ): This distribution should fit the experimentally measured distribution of ISIs. The probability that an ISI coincides with T is given as the product of the probability that firing rate is ξ and the probability that T is drawn from the ISI distribution q(T;ξ,κ). Because the same T value can appear at different firing rates, ξ must be integrated over all possible firing rates.
To estimate κ without knowing the firing rates, we introduce an estimating function y(T;κ) that is independent of the firing rates. The function must satisfy the following condition on the expectation value: for arbitrary k(ξ). In practice, we may replace the integration over p(T;k,κ) with the following sample average: which gives a good approximation to Equation 3 if the size of data are sufficiently large. However, for the time being, we use the integral expression for generality. Because Equation 3 holds for arbitrary k(ξ), it should hold for k(ξ) + δη(ξ), where δη(ξ) is an infinitesimal change in the firing rate distribution: Subtracting Equation (3) from Equation (5) yields the following: for arbitrary k(ξ). Equation 6 implies that y(T;κ) must be orthogonal to in the functional space (i.e., the direction in which the ISI distribution changes maximally with changes in the firing rate distribution). Therefore, the estimating function obtained by solving Equation 5 is indeed free from an arbitrary rate change. Therefore, solving Equation 4 using the observed sequence of ISIs gives a rateindependent estimation of κ. In general, it is difficult to find this optimal estimating function because Equation 6 depends on k(ξ) through In fact, in most stochastic models of spike generation, we cannot find y(T;κ) that is independent of k(ξ). However, in the case of the gamma distribution we can derive such an optimal estimating function and find the shape parameter κ̂ that gives the “best fit” for the observed spike train (Equation 15). The derivation is complicated and lengthy, so we do not repeat it here.
Computational models
We considered the following single compartment model described by Destexhe et al. (2001), which successfully recreated the membrane potential of neocortical pyramidal neurons subject to an intense synaptic bombardment: The values of parameters were the same as those in Destexhe et al. (2001) and Destexhe and Pare (1999). The pointconductance model described by Destexhe et al. (2001) was used to generate realistic synaptic background activity. The total synaptic current, I_{syn}, was decomposed into a sum of two independent conductances: where g_{e} and g_{i} are excitatory and inhibitory conductances, respectively. The reversal potentials were set as E_{e} = 0 mV and E_{i} =75 mV. The dynamics of conductances of excitatory and inhibitory synapses, g_{e} and g_{i}, were described by the following Ornstein–Uhlenbeck process (Uhlenbeck and Ornstein, 1930): where η_{e}(t) and η_{i}(t) are normalized Gaussian white noises with zero means. Because η_{e}(t) or η_{i}(t) represents the sum of many independent excitatory or inhibitory synaptic inputs, respectively, each stochastic variable obeys the Poisson distribution in the limit of infinitely many synapses, regardless of the statistics of presynaptic spike trains. Therefore, the variables can be approximated by Gaussian white noises with appropriate means and variances (cf., Cateau and Reyes, 2006). The time constants are τ_{e} = 2.7 ms and τ_{i} = 10.5 ms, g_{e}_{0} and g_{i}_{0} are the average conductances, and σ_{e} and σ_{i} are the SDs. It was suggested previously that the fluctuation amplitude of synaptic conductance is proportional to its average magnitude (Berg et al., 2007). We set the SDs as and so that the fluctuating synaptic conductances may rarely take negative values. Thus, g_{e0} and g_{i0} specify the statistics of the synaptic currents completely.
We may define the effective reversal potential, V_{r}, as the potential at which the net average synaptic current vanishes: By solving the above equation with respect to V_{r}, we obtain In what follows, we use the parameter pair (g_{e0}, V_{r}) instead of (g_{e0}, g_{i}_{0}) to specify the statistics of the external input. The conductance parameter can be reversed as If we increase g_{e0} with V_{r} fixed, g_{i}_{0} also increases in proportion to g_{e0}. In addition, V_{r} is determined by the ratio of g_{e0} and g_{i}_{0}. Thus, g_{e0} gives the scale factor in our model. Note that g_{e} and g_{i} obey stochastic equations, so their values fluctuate over time and samples even if we fix the values of (g_{e0}, V_{r}).
We conducted numerical simulations of the single compartment model by using the Neuron simulator (Carnevale and Hines, 2006), and estimated the firing rate and κ from the obtained data at various values of g_{e0} and V_{r}. In the present simulations, g_{e0} varies from 0.012 to 0.048 μS and V_{r} from −62 to −50 mV. The simulations were performed long enough to ensure the convergence of firing rate and irregularity κ to steady values.
Wholecell patchclamp recordings.
Wistar rats (postnatal days 17–24) were deeply anesthetized with diethyl ether gas and then decapitated. Cortical slices (400 μm thick) were prepared with a microslicer (PRO7; Dosaka, Kyoto, Japan). After 30 min incubation at 31°C and at least 1 h recovery at room temperature, each slice was transferred to a submergedtype recording chamber continuously circulated with normal artificial CSF (ACSF; 32°C) which consisted of (mm) 124 NaCl, 2.5 KCl, 1.2 KH_{2}PO_{4}, 26 NaHCO_{3}, 1.2 MgSO_{4}, 2.5 CaCl_{2}, and 25 dglucose, and was saturated with 95% O_{2} and 5% CO_{2} gas.
Wholecell patchclamp recordings were obtained from pyramidal neurons in rat sensorimotor cortex, using patch pipettes (7–15 MΩ) filled with (in mm) 140 Kgluconate, 2 NaCl, 1 MgCl_{2}, 10 HEPES, 0.2 EGTA, 2 5′ATPNa_{2}, 0.5 GTPNa_{2}, and 10 biocytin, pH 7.4, (Isomura et al., 2003; FujiwaraTsukamoto et al., 2004; Tsubo et al., 2007). The membrane potentials of neurons were recorded with a current clamp amplifier (Axoclamp 2B; Molecular Devices, Union City, CA) in a conventional bridge mode. Recorded signals were digitized at 50 kHz with an analog–digital interface (Digidata 1322A; Molecular Devices). An analog computation amplifier (SM1; Cambridge Conductance, Cambridge, UK) was also used to inject fluctuating conductance input I_{syn} described in the previous section. To block the spontaneous input via AMPA, NMDA, and GABA_{A} receptors, a mixture of 6cyano7nitroquinoxaline2, 3dione (CNQX; 20 μm), dl2amino5phosphonopentanoic acid (dlAP5; 25 μm), and bicuculline methiodide (BMI; 10 μm) was respectively added to the ACSF and bathapplied to the cortical slices. All three antagonists were purchased from Sigma (St. Louis, MO).
In each trial, a fluctuating current I_{syn} in Equation 8 was injected for 10 s by using the dynamicclamp method (Robinson and Kawai, 1993; Sharp et al., 1993). Statistics of the fluctuating input for a trial is specified by the scale factor g_{e0} and the excitation/inhibition balance V_{r} in Equations 10 and 11. Trials were separated by 20 s intertrial intervals. For a given parameter set (g_{e0}, V_{r}), three trials were repeated by using different seeds for random number generation. The median of the firing rates and irregularities κ for three trials were plotted in the figures. All experiments were performed in accordance with animal protocols approved by the Experimental Animal Committee of the RIKEN Institute.
The coefficient of variation.
In some simulations and experiments, we compare our measure κ with a conventional measure for irregular neuronal firing, C_{V} (the coefficient of variation), defined as follows: where T denotes an ISI and the horizontal bars represent averaging over samples (Softky and Koch, 1993). For a completely regular spike train in which all ISIs are the same, C_{V} becomes 0. For highly unpredictable spike trains generated by a Poisson process, it becomes 1.
Instant by instant Bayesian decoding with gamma distributions of ISIs.
As mentioned with Equation 1, our decomposition scheme is valid when ISIs obey the gamma distribution. Based on this probability density, we can calculate the probability that sequence of ISIs {T_{1}, T_{2}, … } are observed for a given temporal profile of firing rate, where T_{i} denotes the ith ISI in the spike train. As we show later, κ can be estimated by using Equation 15 without knowing the firing rate profile. Then, using the estimated value, we can discriminate different stimuli that lead to different temporal profiles of firing rate by comparing the probabilities that these profiles may underlie an observed output spike train.
Consider the simplest case where the firing rate is constant over time, that is, the neuron receives a stationary input. The probability density that a spike train has spikes at {t_{1}, t_{2}, …, t_{n}} until time t is as follows: where t_{i+1} − t_{i} represents the ith ISI. The last integral represents the probability that there is no spike after t_{n} until t. Given the times of spikes {t_{1}, t_{2}, …, t_{n}}, we can compare different firing rates ξ_{1} and ξ_{2} by calculating P({t }, t; ξ_{1}, κ) and P({t }, t; ξ_{2}, κ). The firing rate that gives a larger probability is more probable. From the Bayes' theorem (MacKay, 2003; Gelman et al., 2004), the posterior probability can be written as follows: where we assume that the previous probabilities that stimulus 1 or 2 is presented are even. The posterior probability can be calculated at an arbitrary time t. Usually, the posterior probability increases with time for the correct ξ until it finally reaches unity, because the evidence accumulates as the number of output spikes is increased. A decision can be made, for example, when the posterior probability exceeds some threshold close to unity.
Similarly, we can treat a more complicated case with a timedependent firing rate ξ (t). At a first glance, the calculations of the posterior probability look complicated. However, all we have to do is to rescale the time axis of dynamics from t to τ (Brown et al., 2001; Kass and Ventura, 2001): Then, the previous probabilities are obtained by replacing t with τ in Equation 11.
Results
We proposed a novel scheme to characterize irregular neuronal firing with its decomposition into firing irregularity and firing rate. The scheme is valid if spike sequences are generated by the gamma process, and keeps the firing irregularity constant over time when the firing rate may vary arbitrarily. To test this coding scheme, we injected in vivolike fluctuating synaptic inputs into in vitro and model cortical neurons and induced highly irregular firing in these neurons. We tested whether the ISI histograms can be well fitted by gamma distributions. Then, we investigated the condition in which the firing irregularity remains unchanged or varies only moderately when the output firing rate is changed significantly. We measured the irregularity of neuronal firing using the shape parameter and tested the reliability of this measure in characterizing the firing irregularity.
Estimation of shape parameter κ from a spike train
The shape parameter κ̂ that best fits an observed spike train can be obtained through the method of estimating functions in the informationgeometry (Amari and Kawanabe, 1997; Miura et al., 2006a,b). The optimal parameter value is given as a solution to the following equation (see Materials and Methods): where T_{i} denotes the ith ISI in the spike train and φ(κ) denotes the digamma function defined as Note that the above expression does not depend on the distribution of firing rate k(ξ). Therefore, solving Equation 15 does not require the information about the firing rate. In deriving Equation 15, we constructed pairs of consecutive ISIs, (T_{1}, T_{2}), (T_{3}, T_{4}), etc., and assumed that we can assign the same firing rate to the members of each pair, that is, ξ_{1} = ξ_{2}, ξ_{3} = ξ_{4}, etc. This assumption resulted in the quadratic forms of ISIs inside the logarithm of Equation 15. However, the estimation error of κ is relatively small even if the consecutive firing rates are different (Miura et al., 2006b). If the ISI distribution deviates largely from the gamma distribution, κ does not have a clearcut meaning. However, results of our in vitro and other in vivo (Barbieri et al., 2001; Brown et al., 2001) studies showed that the gamma process always fits the ISIs reasonably well and the deviations from the gamma distributions are small. Therefore, we may still define the firing irregularity with κ̂, or the solution to Equation 15, because such deviations change the value of κ̂ continuously and moderately. Thus, κ is a useful measure for the firing irregularity of in vivo and in vitro data in various practical situations.
The explicit derivation of the optimal estimator is somewhat complicated. However, it is easy to validate it. If we have many, but finite, samples of ISIs, we may replace the averaging over ISIs in Equation 15 with the following integration over the true probability distribution p(T_{i}, T_{i} _{+ 1}): In deriving the last equation, we used the explicit formula of the gamma distribution given in Equation 1 and rescaled the integration variables as ξT_{i} → T_{i} (i = 1,2). The changes in the integration variables eliminate ξ dependence of the integrand other than k(ξ), allowing us to integrate variable ξ using Thus, Equation 15 is justified for the gamma process with an arbitrary k(ξ).
Fitting the ISI distributions of in vitro cortical neurons
To see the relationship between the irregularity of neuronal firing and the value of κ, we constructed artificial spike trains by generating ISIs according to various gamma distributions. Figure 1 displays typical examples of such spike trains at different firing rates and irregularities. Even if the spike counts are the same in different spike trains, their irregularities can be quite different. The value of κ becomes unity for a completely random spike train generated by a Poisson process. As the value is increased, spike trains become more regular (Fig. 1). The value becomes infinitely large for a completely periodic spike train. Thus, κ indicates the regularity of spike trains in a continuous manner.
We injected fluctuating current I_{syn} into cortical neurons by the wholecell patchclamp technique to measure firing rate and κ of output spikes. In each recording trial, the fluctuating input was injected for 10 s. Figure 2, a and c, shows examples of spike trains recorded from two cortical neurons. We examined whether the ISIs obtained in experiments are distributed according to gamma distributions, the most crucial assumption in our scheme. To this end, we estimated the value of κ for each spike train from the ISIs. To investigate the stability of the value of κ, we calculated the value separately in the initial half and last half of a single spike train. The two segments yielded almost the same value. This demonstrates that κ (firing rate as well) remains constant in time as long as the parameters characterizing the fluctuating input (g_{e0} and V_{r}) are fixed. We found that the ISI histograms were well fitted by the gamma distributions, with κ being the average of the two values. Although the spike sequences look highly irregular, both spike trains exhibited κ values larger than unity, demonstrating that they are more regular than the Poisson spike train. Why the gamma distribution fits the ISI data better than the Poisson distribution partly resides in the fact that the former inhibits the occurrence of short ISIs. This property enables us to take the relative refractory period of neuronal firing into account.
Constant firing irregularity under balanced excitation and inhibition
We conducted numerical simulations of a computational model of neocortical pyramidal neurons and calculated the value of κ in the spike output to synaptic inputs obeying various statistics (Materials and Methods). The model we used is the single compartment model described by Destexhe et al. (2001), which successfully recreated the membrane potential fluctuations of in vivo neocortical pyramidal neurons subject to a continuous synaptic bombardment. Figure 3, a and b, plots the firing rate and irregularity κ obtained from the computational model at various values of g_{e0} and V_{r}, respectively, where g_{e0} is the scale factor of excitatory and inhibitory synaptic conductances, and V_{r} the effective synaptic reversal potential in a balanced state. These parameters determine the statistical properties of the fluctuating input current. In particular, Figure 3b reveals that the value of κ depends greatly on that of V_{r}, but not much on that of g_{e0}. Thus, the results of our simulations predict that increasing g_{e0} with V_{r} kept constant increases the firing rate, while keeping the value of κ in a narrow range. In Figure 3b, the value of κ is not strictly constant as g_{e0} is varied. However, the ranges of κ values for different values of V_{r} have no significant mutual overlaps in the range of g_{e0} values tested. Thus, we may say that κ remains nearly constant when the value of V_{r} is fixed and that of g_{e0} is changed. In contrast, the values of κ are significantly different at different values of V_{r}. In fact, as shown later, stimulus decoding is quite successful under the assumption that κ remains constant regardless of cell's firing rate. We note that the explicit range of κ values is not really important for the present coding scheme. Because κ is a dimensionless parameter, only the differences in its value between different input situations are meaningful.
To investigate whether the computational model actually generates ISIs obeying a gamma distribution, we compared the ISI histograms to the gamma distributions with the estimated values of κ. All histograms for V_{r} = −53 mV exhibit single peaks at the mid range of ISIs, and the optimal fitting curves correspond to similar values of κ (Fig. 4a–c). The unimodal shapes of these distributions can be well fitted by the gamma distributions. The histograms for V_{r} = −62 mV are monotonic functions of ISI that peak at the leftmost bin (Fig. 4d–f). Again, the optimal distributions have approximately the same values of κ. In contrast, the ISI distributions for different values of V_{r} show different shapes with different κ values. We may conclude that the gamma distributions well represent the ISI distributions generated by the computational model.
We examined the prediction of the model in in vitro cortical neurons (n = 16). Figure 5 plots the firing rate and κ for two neurons recorded by wholecell patchclamp technique with various g_{e0} and V_{r}. When g_{e0} is increased while keeping V_{r} constant, the firing rate is largely increased, but the manipulation changes the value of κ only moderately. For example, the firing rate of cell 1 at the maximum value of g_{e0} is three times larger than that at the minimum value of g_{e0} for V_{r} = −41 mV. This tendency is consistent with that of the numerical results.
We tested whether the ISI distributions obey gamma distributions in the neurons displayed in Figure 5. Figure 6, a and b, plots the ISI histograms of cell 1 at V_{r} = −41 and −46 mV, respectively. Because the data available for this analysis were limited, we combined the data sets recorded at different values of g_{e0} to construct these histograms. The ISIs were normalized by the corresponding mean ISIs so that the shape of the distributions may be specified only by κ. This manipulation is justified if the ISI distribution represents a gamma distribution (as seen below). The histograms exhibit single peaks and are well represented by gamma distributions. Figure 6, c and d, shows similar plots of the ISIs for cell 2 at V_{r} = −40 and −44 mV, respectively. We find that gamma distributions fit the histograms reasonably well.
Because of neurontoneuron variations in the membrane excitability, the range of V_{r} that produces moderate rate changes differs in different neurons. Nevertheless, the tendency that κ depends on V_{r} rather than on g_{e0} was commonly observed in all recorded neurons. The firing rate and κ plotted for all 16 neurons recorded in this study at various values of g_{e0} and V_{r} confirmed this tendency (Fig. 7). In fact, the mean of the slopes of κ in Figure 7b determined by leastsquare fitting with a straight line is given as −1.4 (μS^{−1}), which is not significantly different from 0 (p > 0.1). However, the mean slope of the firing rate in Figure 5a is given as 540 (Hz/μS), which is significantly different from 0 (p < 0.05), showing an obvious increasing behavior. However, Figure 7, c and d, shows that increasing V_{r} with g_{e0} kept constant increases both firing rate and κ simultaneously. Thus, the spiking irregularity κ can be kept constant only when V_{r} is kept unchanged.
Continuous modulations of firing rate do not affect firing irregularity
In the previous section, we have shown that κ is nearly constant in the condition where the rate of neuronal firing is regarded as constant within each trial, but the rate might vary from trial to trial. However, in a realistic situation, the firing rate can vary over time depending on the stimulus intensity. Results of previous experiments suggested that excitatory and inhibitory synaptic inputs to cortical neurons are always balanced (Shu et al., 2003; Haider et al., 2006; Berg et al., 2007). Although the details of this balanced synaptic input should be further elucidated, it may be the case that the growth of excitation is followed immediately by that of inhibition in the cortical dynamics.
To study such a situation, below we change the firing rate of postsynaptic neuron by modulating the amplitude of g_{e0} (accordingly, that of g_{i}_{0}) periodically, with V_{r} kept at −44 mV. The conductances g_{e}(t) and g_{i}(t) fluctuate around the mean values g_{e0}(t) and g_{i}_{0}(t) obeying Equation 9 (see Materials and Methods), where the mean values were defined as follows: where T, the period of the sinusoidal modulation, was chosen in a range of 100–1600 ms. The values are much longer than the time constants of synaptic conductances in Equation 9 (τ_{e} = 2.7 ms and τ_{i} = 10.5 ms), so the timevarying modulation is sufficiently slow and should not cause the excitation and inhibition to go in and out of balance depending on the phase. The value of g_{e0} was modulated between 0.012 and 0.03 μS. Figure 8 shows the firing rate, κ and C_{V} of cortical neurons recorded with the periodically modulated input. For comparison, we recorded the responses of the same neurons to stationary inputs, in which V_{r} was fixed within a range of −40 to −44 mV and g_{e0} was set to one of the following values: 0.012, 0.018, 0.024, 0.03 μS. We note that the values of firing rate and κ were kept well between the two gray lines that designate the value ranges obtained for stationary inputs with V_{r} = −44 mV. In contrast, the values of C_{V} do not stay within the range delimited by the gray lines. These results demonstrate that κ is robust against the modulations in firing rate, whereas C_{V} is strongly influenced by the modulations.
Why does C_{V} not describe the firing irregularity properly when the firing rate changes with time? In general, C_{V} tends to be large when the firing rate changes with time (Stevens and Zador, 1998; Sakai et al., 1999; Shinomoto et al., 1999, 2005a). In fact, we can mathematically prove that the expectation value of C_{V} is always larger than unity for an inhomogeneous Poisson process with a timevarying firing rate (Tuckwell, 1988; Shinomoto and Tsubo, 2001). To see this clearly, spike trains were generated by Poisson process such that the mean firing rate changes discontinuously at time 0 (Fig. 9). By construction, both κ and C_{V} should be 1 for the (stationary) Poisson processes at negative and positive times. However, if we compute κ and C_{V} for the whole spike train, κ is 1 but C_{V} is This overestimation of C_{V} is attributable to the fact that the SD is calculated around the “global” average of ISIs over the entire spike train. Thus, C_{V} cannot properly capture a local irregularity of spiking when the firing rate is not constant. In contrast, κ is robust against the rate modulation because its effect (i.e., the scaling of neighboring ISIs), is locally cancelled out in the logarithm in Equation 15 (Miura et al., 2006a,b). This property of κ makes it more adequate than C_{V} for measuring the irregularity of neuronal firing.
Decoding signals from neuronal responses
The previous results have shown that a balanced increase of excitatory and inhibitory synaptic inputs keeps the irregularity of postsynaptic spikes approximately constant. In other words, the firing pattern roughly obeys a gamma process in the balanced regime. Below, we demonstrate how this knowledge of the spike generating process may improve our ability to decode signals from neuronal firing. We consider the discrimination of two constant stimuli from output spike trains when their firing rates are different (Fig. 10a). A similar test was considered in Wiener and Richmond (2003). We used the spike data recorded in vitro in Figure 5 to mimic the stimulus discrimination by in vivo neurons. Stimulus 1 or 2 was represented by current injection with g_{e0} = 0.04 μS or 0.025 μS, respectively. In both cases, V_{r} was −41 mV. Note that the output spike trains for the two currents obey the gamma distributions with almost the same κ values (≈2.47), but with different rates. Thus, we can discriminate the inputs according to the probabilities that the observed ISIs are generated from the different distributions. We attempted to discriminate the two stimuli by calculating the posterior probability that given an output spike train stimulus 1 or 2 was presented. We performed these calculations instant by instant based on the gamma distributions (see Materials and Methods). Figure 10b shows the resultant posterior probability of stimulus 1 when it is actually presented. Note that 1 minus this probability gives the posterior probability of stimulus 1 when stimulus 2 is presented. The posterior probability of the correct discrimination on average increased with time or the accumulation of output spikes. The decoding performance was optimal if we knew the optimal value of κ in advance (solid line). If, however, we adjusted the decoder to a Poisson process (κ = 1), the performance was much worse (dotdashline). Thus, we can discriminate the stimuli quite efficiently by estimating κ from the output spikes. The different performances emerge from the fact that the ISI distributions with different κ values have different shapes and assign different posterior probabilities to each value of ISI.
Discussion
We have proposed a novel scheme to characterize irregular neuronal firing based on the fact that the firing irregularity can be decoupled from the firing rate if the spike sequence obeys a gamma distribution. We have tested whether this coding scheme may be valid for cortical neurons by injecting in vivolike fluctuating inputs to them. We have shown that the highly irregular firing of these neurons actually exhibit the ISI distributions obeying the gamma distribution. Moreover, the firing irregularity measured from the ISI distributions remained constant over time. We found that balanced excitatory and inhibitory synaptic input is essential for this coding scheme.
Computational merit of the novel coding scheme
Neuronal information processing significantly relies on the firing rate. For instance, cortical neurons may change their firing rates depending on sensory stimuli from the environment or the motor responses prepared and executed in response to such stimuli. It is, therefore, crucial for decoding neural information to estimate the instantaneous firing rate of temporally localized spikes. If neurons fire rhythmically, we can easily obtain an accurate estimation of the instantaneous firing rate by taking the inverse of an ISI. However, if the spiking pattern is highly irregular, the estimation, or even the definition of the instantaneous firing rate is not trivially easy, and the estimated value may not be fully reliable. The constant firing irregularity shown in this study significantly improves the reliability of the firing rate estimation. If the firing irregularity also varied with time, the estimation of the firing rate would become much more difficult. In fact, it is an illdefined problem to simultaneously estimate two parameters (firing rate and irregularity) from observed data (ISIs). In some cases, the estimation of the firing rate can be very difficult even when the irregularity of neuronal firing remains constant (Koyama and Shinomoto, 2005).
Because the output spike trains of cortical neurons are highly irregular, synaptic input from a population of such neurons to a cortical neuron will exhibit a large amount of noise. The brain may try to keep the irregularity in neuronal firing as constant as possible to decode information from the noisy signals. The instant by instant Bayesian decoding shown in Figure 8 is one example of how this decoding process might take place. The strategy used therein is simple, and is applicable to many practical problems to estimate the firing rate from irregular spike trains: (1) estimate the firing irregularity κ for each neuron of interests in advance, and (2) estimate the firing rate with the help of the estimated value of κ. Note that the estimation of κ can be done whatever the firing rate profile is. At step 1, only knowing the typical value of κ for cortical neurons in advance would suffice for most of the practical situations. For instance, the above procedure may be applied to brainmachine interfaces to improve the performance of signal detection or discrimination. Thus, the use of the firing irregularity κ in conjunction with the mean firing rate provides an efficient method to decode information from neuronal activity.
Balanced synaptic inputs “clamp” the firing irregularity
We have shown that the firing irregularity significantly depends on the effective synaptic reversal potential V_{r}. In this study, the value of V_{r} was kept constant to maintain the firing irregularity at a fixed level. Results of recent experiments suggest that constant V_{r} is indeed the case in cortical networks in vivo. It has been shown that the changes in the average excitatory synaptic conductance are balanced with those of inhibitory ones in cortical and spinal cord neurons (Shu et al., 2003; Haider et al., 2006; Berg et al., 2007). The balanced synaptic inputs keep the membrane potential fluctuating around a constant value of V_{r} and make the irregular firing of cortical neurons describable by the proposed framework. Thus, our results suggest that a functional implication of balanced synaptic input is the decomposition of rate information from the irregularity of neuronal firing.
In the previous experiments, the ratio between the average excitatory and inhibitory conductances varied from neuron to neuron. Because the absolute value of V_{r} can be determined by the ratio in number between the excitatory and inhibitory synapses innervating a neuron, it seems that the individual cortical neurons exhibit different values of V_{r}, or different degrees of the firing irregularity. The balance between excitation and inhibition may be broken when, for example, a neuronal circuit performs certain operations (BorgGraham et al., 1998; Chance et al., 2002; Wehr and Zador, 2003). In that case, the measure κ may detect the shift of the balance as it is sensitive to the ratio between excitation and inhibition. As our measure κ is decoupled from firing rate, the outofbalance state may give κ an opportunity for exploring a role different from that of firing rate in the information representation by neurons.
Why is the irregularity κ kept constant in the balanced regimen? It is well known that choosing every κth spike from a Poisson spike train yields a spike train that obeys a gamma distribution of shape parameter κ. Therefore, the output spike train may obey such a gamma distribution if a neuron fires at every moment when it completes integrating κ input spikes obeying a Poisson process. The number (= κ) of inputs required for spike generation may depend on various parameters such as the amplitude of each input, the time constant of leakage, and, particularly, the distance from the effective reversal potential V_{r} to threshold. Therefore, κ remains unchanged if V_{r} is fixed. Because the observed values of κ are small (∼2), in reality they may represent the required number of spike packets rather than single spikes.
Next, why is κ almost independent of g_{e0}? If we increase the means of g_{e} and g_{i} (i.e., the values of g_{e0} and g_{i0}), their variances are also increased in our model. Increasing the means strengthens the “force” that clamps the membrane voltage at V_{r}, which would lower the firing rate and the irregularity. In contrast, increasing the variances raises the firing rate and the irregularity. The two competing effects as a whole increase the firing rate and maintain κ. Thus, the correlated changes in the mean and SD of synaptic conductance (Shu et al., 2003; Berg et al., 2007) are essential for the constancy of κ.
Does the validity of the firing irregularity strongly depend on the membrane properties of neurons? We examined the constancy of κ in the balanced regime using several different neuron models such as leaky integrateandfire neurons (data not shown). The leaky integrateandfire model is one of the simplest computational neuron models (Koch, 1999; Dayan and Abbott, 2001). The results obtained were essentially the same as those obtained from the realistic neuron model or cortical neurons, implying the generality of our results.
Advantage of measure κ over C_{V}
Previous studies pointed out that the conventional measure C_{V} also remains unchanged in the balanced regime. However, the statistics of input spikes was limited and, more importantly, the firing rate was fixed at constant in these studies (Salinas and Sejnowski, 2000; Chance et al., 2002). Our results have extended these results to a wider class of input statistics, demonstrating a new finding that the firing irregularity is robust against the firing rate changes in the balanced regimen. In contrast, C_{V} is sensitive not only to the firing irregularity but also to the rate change (Fig. 9). In fact, the large C_{V} values observed in vivo, which some attempted to reproduce by in vitro experiments (Stevens and Zador, 1998) or computer simulations (Softky and Koch, 1993; Troyer and Miller, 1997; Shadlen and Newsome, 1998; Sakai et al., 1999; Shinomoto et al., 1999), may be mostly caused by the modulations of firing rate. We propose that κ is more suitable than C_{V} for defining the firing irregularity.
We note that some of the previously proposed measures for the firing irregularity, such as and are not influenced by the rate variation like the present κ (Holt et al., 1996; Shinomoto et al., 2003, 2005b). The derivation of these measures, however, is adhoc and they have no corresponding parameters in statistical models. In contrast, κ has a clearcut meaning based on information mathematics (Bickel et al., 1993; Amari and Kawanabe, 1997; Amari and Nagaoka, 2001; Miura et al., 2006a,b). In fact, κ is an optimal estimator of the shape parameter of the gamma distribution of ISIs. In other words, each value of κ specifies a unique shape of the ISI distribution. Therefore, κ is directly linked to a statistical model and, hence, is useful for decoding.
In a widely accepted view of neuroscience, information is conveyed by the instantaneous firing rates of neurons, and an inhomogeneous Poisson process, where each time bin is independent, is often used for a probabilistic model of the spike generation. However, our results indicate that this process more resembles the gamma process rather than the Poisson process, when cortical neurons receive balanced excitatory and inhibitory synaptic inputs (Shu et al., 2003; Haider et al., 2006; Berg et al., 2007). Other studies performed a model selection using Kolmogolov–Smirnov statistical test on the ISIs recorded from in vivo cortical neurons (Barbieri et al., 2001; Brown et al., 2001). The measurement of L_{V} in in vivo cortical neurons also suggested that the constancy of the irregularity (Shinomoto et al., 2003, 2005b). The results of these studies are consistent with the present ones, confirming that the gamma process fits the ISI distribution reasonably well. However, the previous studies did not address the role of balanced synaptic inputs in keeping the firing irregularity constant. We emphasize that our study has clarified the profound relationship between the gamma process and the constant irregularity (i.e., κ), and the biological conditions to ensure the validity of the gamma process (i.e., balanced synaptic input).
We demonstrated that the shape parameter κ in the gamma process augments the information provided by the firing rate, resulting in more accurate estimation of the posterior probabilities of presented stimuli. This fact may enable us to improve the quality of decoded signals in the brainmachine interface, a recent technological challenge to use brainderived signals for controlling robots and other information processing systems. In summary, we have shown that balanced synaptic inputs clamp the irregularity of neuronal firing against the rate modulation. Our results have revealed a possible functional role of balanced synaptic inputs in the cortical information representation.
Footnotes

This work was supported in part by GrantsinAid from Ministry of Education, Culture, Sports, Science and Technology (1810171, 18020007, 18079003, 17022036). We are grateful to T. Momiyama and A. Momiyama for lending us an SM1. We are grateful to A. Reyes for comments on this manuscript.
 Correspondence should be addressed to Keiji Miura, University of Tokyo, 409 Kibantou, Kashiwanoha, Kashiwa, 2778561 Chiba, Japan. miura{at}mns.k.utokyo.ac.jp