## Abstract

Spike-timing variability has a large effect on neural information processing. However, for many systems little is known about the noise sources causing the spike-response variability. Here we investigate potential sources of spike-response variability in auditory receptor neurons of locusts, a classic insect model system. At low-spike frequencies, our data show negative interspike-interval (ISI) correlations and ISI distributions that match the inverse Gaussian distribution. These findings can be explained by a white-noise source that interacts with an adaptation current. At higher spike frequencies, more strongly peaked distributions and positive ISI correlations appear, as expected from a canonical model of suprathreshold firing driven by temporally correlated (i.e., colored) noise. Simulations of a minimal conductance-based model of the auditory receptor neuron with stochastic ion channels exclude the delayed rectifier as a possible noise source. Our analysis suggests channel noise from an adaptation current and the receptor or sodium current as main sources for the colored and white noise, respectively. By comparing the ISI statistics with generic models, we find strong evidence for two distinct noise sources. Our approach does not involve any dendritic or somatic recordings that may harm the delicate workings of many sensory systems. It could be applied to various other types of neurons, in which channel noise dominates the fluctuations that shape the neuron's spike statistics.

## Introduction

In many sensory systems, neuronal responses vary significantly when the same stimulus is presented repeatedly (de Ruyter van Steveninck et al., 1997; Kreiman et al., 2000; Ronacher et al., 2004). This response variability has a profound impact on sensory signal processing (Kara et al., 2000; Vogel et al., 2005) and limits the reliability of behavioral responses to external stimuli (Tolhurst et al., 1983; Britten et al., 1996). Noise sources at all levels of neuronal pathways cause this variability (Faisal et al., 2008), e.g., external noise that is amplified during the transduction process (Grewe et al., 2003), synaptic noise resulting from the bombardment of myriads of synapses (Calvin and Stevens, 1967), or noise caused by the stochastic gating of ion channels (“channel noise”) (Neher and Stevens, 1977; White et al., 2000).

The statistics of channel noise depends on the gating dynamics of the ion channels. In particular, slow stochastic adaptation currents generate colored noise with a correlation time scale set by the adaptation time constant (Schwalger et al., 2010). In contrast, fast spike-generating currents lead to approximately white noise. These different time scales have distinct effects on the spiking statistics (Lindner, 2004; Schwalger et al., 2010) as quantified by the distribution and serial correlations of interspike intervals (ISI). Conversely, the ISI statistics points to a certain time scale of the noise and may enable one to identify the dominating noise source for a given cell. This is particularly interesting for sensory neurons where action potentials can be recorded far away from the delicate transduction machinery without damaging it.

Here we investigate spike-timing variability and potential noise sources in a key insect model system—locust auditory receptor neurons. Although these neurons respond quite variable to pure tones with constant amplitude (Schaette et al., 2005), they transmit dynamical stimuli remarkably well (Machens et al., 2001, 2005) based on precise spike timing (Rokem et al., 2006). Locust auditory receptor cells are bipolar neurons that are directly attached to the tympanal membrane (Gray, 1960) (see Fig. 1). The underlying mechanosensory transduction process is well understood in terms of a cascade of linear filters and static nonlinearities (Gollisch et al., 2002; Gollisch and Herz, 2005). Furthermore, these neurons exhibit spike-frequency adaptation (Römer, 1976; Gollisch and Herz, 2004).

To uncover basic properties of the primary noise sources, we first quantify the variability of ISIs by analyzing their statistics during acoustic stimulation with pure tones and compare the results with predictions from canonical models of suprathreshold, i.e., limit-cycle firing. Second, we investigate the distinct effects of channel noise from different ionic currents on the ISI statistics in simulations of conductance-based models of the locust auditory neuron to test the results from the ISI analysis and to study potential biophysical mechanisms. A model variant with both a stochastic adaptation current and a stochastic sodium or receptor current reproduces our experimental findings. Combining experimental data, theory, and simulations enables us to distinguish slow from fast noise solely from the ISI statistics and to put forward strong hypotheses on the noise sources.

## Materials and Methods

#### Electrophysiology

Intracellular recordings were performed from axons of auditory receptor neurons of adult *Locusta migratoria* (seven males and five females). All experiments were performed at room temperature (22−24°C). After decapitation and removal of legs and wings, the abdomen's tip was cut to extract the gut. The ventral side of the animal was then glued to an animal holder using sticking wax. The upper part of the metathorax was removed and the torso was filled with locust saline (Pearson and Robertson, 1981). After exposing the metathoracic ganglion and the auditory nerve, the nerve was fixed with a custom-made forceps. The recordings were performed from the proximal fibers in the auditory nerve. We used standard glass microelectrodes (borosilicate; 1.5 mm outer diameter; World Precision Instruments and GB150F-8P, Science Products) pulled to a resistance of 30–80 Ω (Model P-97; Sutter Instrument Company) and filled with a 1 m KCL solution. The signals from the auditory receptor cells were amplified (BRAMP-01 or SEC-05LX; npi electronic) and recorded by a data acquisition board (PCI-6229; National Instruments) with a sampling rate of 20 kHz. For on-line spike detection, data analysis, and the generation of acoustic stimuli we used the open-source software *RELACS* (www.relacs.net) running on a Debian Linux computer. Acoustic stimuli transmitted by the acquisition board (sampling rate: 120 or 240 kHz) were sent to a custom-made attenuator (ATN-01M; npi electronic) based on the Chrystal CS3310 attenuator chip (Cirrus Logic) for adjusting the sound intensity, amplified (RKB-250; Rotel), and then forwarded to two loudspeakers (DSM 25 FFL-8; Visaton). The loudspeakers were positioned orthogonal to the locust body axis at a distance of 25 cm from the animal. For calibrating the sound intensity we used a 1/2 inch microphone (40AC on a 26AM preamplifier; G.R.A.S.) on a constant voltage amplifier (12AK; G.R.A.S.). To reduce echoes, the experimental setup was lined with sound-attenuating foam.

#### Stimuli

For acoustic stimulation, pure tones of different but constant sound intensities (measured in dB sound pressure level; SPL) were presented ipsilaterally to the recorded auditory neuron. The sound frequency was set to the characteristic frequency of the recorded receptor neuron (2–17 kHz, in most cases 4 kHz). All stimuli started and ended with a 2 ms ramp. Short stimuli of 500 ms duration separated by gaps of 1 s were used for most of the analysis. Long stimulus durations of 1 min were used to measure ISI correlations at spike frequencies <50 Hz.

#### Data analysis

For the data analysis, we used *N* = 12 recordings of auditory receptor cells of different animals with stimulus durations of 500 ms (recording duration > 20 min, 30–300 trials per sound intensity) and *N* = 14 recordings of different animals with stimulus durations of 1 min. Due to the firing rate transients induced by spike-frequency adaptation we disregarded the first 200 ms of each trial of 500 ms duration and the first 10 s of each trial evoked by a 1 min pure tone stimulus. Especially at higher sound intensities the data in the remaining analysis windows were not perfectly stationary. We therefore recomputed the ISI statistics after cutting the analysis time window in three. Except for minor quantitative differences this did not change our results so that we report only values based on the full analysis windows. Data analysis was performed with C++ programs, using the numerics library of *RELACS* (www.relacs.net).

ISI variability was quantified by both the coefficient of variation *D* = 〈Δ*T*^{2}〉/(2〈*T*〉^{3}) where 〈*T*〉 is the mean and 〈Δ*T*^{2}〉 = 〈*T*^{2}〉 − 〈*T*〉^{2} the variance of the ISIs. *D* corresponds to the diffusion coefficient of the spike count of the ISI shuffled spike train and will be referred to as the “diffusion coefficient”. For a perfect integrate-and-fire (PIF) model driven by white noise, *D* is proportional to the noise intensity of the driving fluctuations, i.e., to the input noise intensity (Vilela and Lindner, 2009).

We also computed ISI histograms and compared them with probability density functions of simple noise-driven neuron models. The ISI density of a PIF model with white noise driving is given by the following:
with a mean ISI 〈*T*〉 (Gerstein and Mandelbrot, 1964). This probability density also approximates well the ISI density of a PIF neuron with an adaptation current and additive white Gaussian noise of intensity *D* (Schwalger et al., 2010).

In contrast, an approximation of the ISI density of an adapting PIF neuron with a stochastic adaptation current and no white noise driving is given by a probability density resulting from a PIF neuron driven by long-correlated colored noise (Ornstein–Uhlenbeck noise; Lindner, 2004) as follows:
with γ_{1}(*T*) = *T*/τ+*e*^{−T/τ} − 1, γ_{2}(*T*) = 1 − *e*^{−T/τ} and correlation time τ. The parameter ε = σ^{2}/μ^{2} is the variance of the noise rescaled by the square of the mean current μ of the PIF model. The squared coefficient of variation of the density equation (2) is given by the following:
with δ = 〈*T*〉/τ. Given the mean and variance of the ISIs, we can eliminate ε via Equation 3 (the quadratic equation has only one positive solution in ε) and then fit *p _{cn}*(

*T*) to the ISI histograms by a least-square fit of the unknown correlation time constant τ using a simplex algorithm.

To quantify the shape of the ISI densities we used two measures introduced by Schwalger et al. (2010): rescaled versions of the skewness γ* _{s}* = 〈(

*T*− 〈

*T*〉)

^{3}〉/〈Δ

*T*

^{2}〉

^{3/2}and kurtosis γ

*= 〈(*

_{e}*T*− 〈

*T*〉)

^{4}〉/〈Δ

*T*

^{2}〉

^{2}− 3. The rescaled measures are given by the following: and are identical to unity for the inverse Gaussian ISI density (Eq. 1). Values considerably larger than one indicate that the density is more skewed and peaked compared with

*p*(

_{wn}*T*), which is the case for the colored noise distribution

*p*unless τ ≪ 〈

_{cn}*T*〉.

We quantified the correlations among succeeding ISIs by the serial correlation coefficient
For the analysis of the spike train responses evoked by short acoustical stimuli of 500 ms, we calculated ρ_{1} by averaging over all pairs of successive ISIs from all trials, because the number of ISIs within each trial was too small to obtain unbiased estimates for the correlation of each trial separately. For the 1 min stimuli, we divided each trial in sections of 300 ISIs for which we calculated the correlation coefficient according to Equation 6 separately and then averaged over all sections of all trials.

To test the significance of negative serial correlations we determined the distribution of correlations occurring by chance after shuffling the ISIs of each spike train section 2000 times. The measured serial correlation was considered significant if the integral over the shuffled distribution up to this value was smaller than 0.05 (see Fig. 6*B*).

#### Model of the auditory receptors

A schematic drawing of a locust auditory receptor neuron is given in Figure 1*A*. Deflections of the tympanal membrane induced by a sound stimulus cause mechanosensory ion channels to open. The resulting ionic currents activate voltage-dependent ion channels triggering action potentials as well as additional currents mediating spike-frequency adaptation. We separated these processes in two parts: (1) the spike generator with spike-frequency adaptation and (2) the mechanosensory transduction cascade from the sound input to the activation of the mechanosensory currents driving the spike generator.

##### Spike generator with spike-frequency adaptation.

For the simulation of the membrane potential dynamics and the generation of action potentials, we used the single-compartment conductance-based Traub–Miles model (Traub and Miles, 1991) modified by Ermentrout (1998). This model includes an M-type current, i.e., a slow voltage-dependent potassium channel, mediating spike-frequency adaptation. We have chosen this model of mammalian hippocampal pyramidal cells for its well known type-I dynamics. Details of the ion-channel kinetics that might be different in the grasshopper do not matter for our study, since we use this model only for a qualitative but not quantitative exploration of the effects of ion-channel noise.

The membrane potential *V* (in millivolts) is described by the following:
where *C _{m}* = 1 μF/cm

^{2}and

*I*is the driving current passing through the mechanosensory channels of an auditory receptor neuron (in μA/cm

_{R}^{2}). The deterministic ionic currents are given by the following equations (Ermentrout, 1998): In this parametrization the M-type current is mainly activated during action potentials. Subthreshold activated adaptation currents will not change the following results, because we consider exclusively the suprathreshold firing regime far beyond the bifurcation from quiescence to tonic firing (Schwalger et al., 2010).

In the above equations, all numerical values corresponding to voltages and to rate constants are given in millivolts and kilohertz, respectively. To yield a membrane time constant in the submillisecond range (Gollisch and Herz, 2005), the sodium, potassium, and leak conductances were multiplied by ϑ = 20. The membrane time constant for the model using *ḡ _{L}* = 2 mS/cm

^{2}is τ

*=*

_{m}*C*/

_{m}*ḡ*= 0.5 ms. The driving current

_{L}*I*described below was chosen to result in a steady-state firing rate in its saturation level of 180 Hz which is ∼150 Hz smaller than the onset-state firing rate (Benda et al., 2001; Fig. 2

_{R}*B*).

##### Model of the mechanosensory transduction process.

The activating input signal for the mechanosensory transduction process was a pure tone sound pressure wave *s*(*t*) of amplitude *A* and frequency *f _{s}*:

*s*(

*t*) =

*A*sin(2π

*f*). The sound frequency was set to

_{s}t*f*= 4 kHz, a typical characteristic frequency of low-frequency receptors of locusts (Michelsen, 1971). The amplitude

_{s}*A*can be calculated from the sound intensity

*I*(in dB SPL) by

*A*= 20 μPa · 10

^{I}^{/20}. The amplitude of the tympanal deflection depends linearly on the sound pressure for pure tones (Michelsen and Rohrseitz, 1995). Thus, the tympanal vibration is proportional to the sound pressure wave:

*x*(

*t*) =

*c̃*·

*s*(

*t*).

For simulations and estimations of the magnitude of thermal noise arising from the eardrum we also explicitly modeled the tympanum filter. We approximated the tympanum as a damped harmonic oscillator that is driven by sound pressure stimuli *s*(*t*) plus additive thermal noise
where *x*(*t*) is the deflection of the tympanum, ω_{0} = 2π*f _{s}* the Eigenfrequency,

*m*is the mass, and

*B*is the area of the oscillating element.

*B*scales sound pressure into a force driving the oscillator. Despite the periodic driving we assume that dissipation (friction) and fluctuations (noise) are connected via the Einstein relation where

*k*= 4.1 · 10

_{B}T^{−21}kg m

^{2}/s

^{2}at 20°C is the thermal energy. Dividing Equation 8 by

*m*and re-expressing the damping coefficient γ = 2

*m*/τ

*in terms of the damping time constant τ*

_{d}*(≈100 μs; Schiolten et al., 1981; Gollisch and Herz, 2005) results in the following: with α =*

_{d}*B/m*. Stroboscopic measurements estimate the amplitude of the tympanum deflection evoked by a 100 dB SPL tone at 3 kHz to be approximately

*x̄*= 16 μm (Breckow and Sippel, 1985). This can be used to estimate α to be ∼3000 m

^{2}/kg using the gain of the oscillator,

*x*

_{signal}

^{2}〉 =

*x̄*

^{2}/2) divided by the variance 〈

*x*

_{noise}

^{2}〉 caused by the thermal noise is then With the smallest sound intensities used in our experiments of

*I*= 30 dB SPL, corresponding to

*A*= 632 μPa, and a variation of the unknown effective area

*B*between the whole tympanum's size (

*B*

_{tym}≈ 1 mm

^{2}) and a small fraction of this value (

*B*

_{tym}/100), we obtain signal-to-noise ratios between 400 and 4, respectively. Put differently, for plausible values of the effective area, the signal-to-noise ratio is at least half an order of magnitude larger than one. Higher sound intensities add for every 20 dB another factor of 100. For the simulations we varied

*D*

_{tymp}and used for Equation 13 below

*c*= 25 μm

^{−1}and

*x*

_{0}= 0.12 μm.

The cilia of the locust auditory receptor neurons are connected to the tympanum via attachment cells (Gray, 1960). Gollisch and Herz (2005) have shown that the transmitted signal during the transduction process between the mechanical vibration of the eardrum and the accumulation of electrical charge at the receptor cells' membrane undergoes a squaring. Moreover, the spike frequency of spike train responses saturates at stimulus intensities ∼20 dB above threshold intensity (Römer, 1976). This saturation is caused by the saturation of the mechano-electrical transduction process (Wolf and Benda, personal communication).

For a minimal model consistent with this information, we assumed two populations of ion channels in the membrane of the mechanosensory receptor neurons that open either for positive or negative deflections of the eardrum and whose open probabilities follow a Boltzmann distribution. The sum of the two Boltzmann functions of both receptor channel populations comprises, thus, both the quadratic nonlinearity at low and medium sound intensities and saturation of firing rates:
where *p*_{o}^{+}(*x*) and *p*_{o}^{−}(*x*) are the open probabilities for positive and negative tympanal displacements, respectively. Assuming the channel gating acting as a spring following Hooke's law, the open probabilities are given by the following:
with the single channel-gating force *c* indicating the channel's sensitivity to the tympanal deflections and *x*_{0} being the displacement for which the channel open probability is 0.5 (Howard and Hudspeth, 1988; Hudspeth et al., 2000). The parameters *c̃*, *c*, and *x*_{0} influence both the dynamic range of the auditory receptor neuron model and the minimal open probability at the tympanal resting position. We chose *c̃c* = 0.00025 μPa^{−1} and *x*_{0}/*c̃* = 12000 μPa so that the resulting *f-I* curves are comparable to the ones observed for locust auditory receptor cells in particular with a dynamic range that is ∼20 dB wide and with a peak of the diffusion coefficient *D* in the lower dynamic range (Fig. 3).

In the models with deterministic receptor current mediated by an infinite population of receptor channels (*N _{R}* → ∞), the dynamics of the open probability of the two kinds of channels (+ or −) is given by the following:
with the receptor channel time constant τ

*= 0.1 ms. For smaller values of τ*

_{r}*the generated action potentials lock on the carrier sound wave, a behavior that has not been observed experimentally, and larger values are not supported by the time scale of ∼0.1 ms of the second linear filter of the locust auditory transduction cascade (Gollisch and Herz, 2005). The mean number of open receptor channels corresponds to*

_{r}*N*=

_{o}*N*

_{R}p^{+}+

*N*

_{R}p^{−}with

*N*being the total number of channels of both receptor channel populations. Note that in this study the total number of simulated receptor channels is always given as the sum of channel numbers of both channel populations and both populations have the same number of channels. The receptor current

_{R}*I*passing a fraction of open channels

_{R}*W*=

*N*/

_{o}*N*is determined by the following: where

_{R}*E*= 0 mV is the reversal potential of the receptor current and

_{R}*ḡ*= 0.6 mS/cm

_{R}^{2}its maximal conductivity. The latter was chosen to ensure an

*f-I*curve sufficiently close to the experimental data.

In the stochastic models, the receptor channel gating was simulated using Gillespie's algorithm (see Simulation of stochastic opening ion channels). The transition rates α* _{r}* and β

*defined in the Markov scheme (17) below are given by the following: For the deterministic models, the current passing the receptor channels is given by the fraction of open channels*

_{r}*W*as defined in Equation 15.

##### Kinetic schemes for the stochastic ion channel models.

To describe spike generation, we used different stochastic models. In the first class of model, one type of ion channel was simulated as a finite population of *N* independent channels mediating a stochastic current while the other channel types corresponded to infinite channel populations mediating deterministic ionic currents. The second class of stochastic model only differed in that we assumed a mixture of two finite channel populations instead of only one.

Stochastic ion channels can be modeled as Markov processes (Neher and Stevens, 1977). Assuming that the channels carrying the receptor current are ion channels consisting of the two states “closed” and “open”, this elementary kinetic scheme can be written as follows:
where α* _{r}* and β

*define the rates (transition probability per unit time) from*

_{r}*r*

_{0}to

*r*

_{1}and vice versa, respectively. According to the Markov assumption the steady-state channel open probability depends only on the state in which the system is and can thus be defined by

*p*

_{∞}= α

*/(α*

_{r}*+ β*

_{r}*). The time required to perform a state transition is τ*

_{r}*= 1/(α*

_{r}*+ β*

_{r}*).*

_{r}Assuming a Markov process for the potassium channel composed of four identical gates *n _{i}* with an opening transition rate α

*and a closing transition rate of β*

_{n}*, the potassium channel kinetic scheme is where*

_{n}*n*

_{4}corresponds to the open state where all four channel gates are open (Neher and Stevens, 1977).

The Markov kinetic scheme for a sodium channel with three activating gates *m _{i}* and one inactivating gate

*h*is expressed by the following: where

*m*

_{3}

*h*

_{1}defines the open state where the three activating gates

*m*and the one inactivating gate

*h*are open (Neher and Stevens, 1977).

For the channels carrying the adaptation current we assumed (as for the receptor channels) a simple two state channel model:
where α* _{w}*(

*V*) =

*w*

_{∞}(

*V*)/τ

*(*

_{w}*V*) and β

*(*

_{w}*V*) = (1 −

*w*

_{∞}(

*V*))/τ

*(*

_{w}*V*) are the transition rates between the closed state w

_{0}and the open state

*w*

_{1}.

##### Simulation of stochastic opening ion channels.

The gating of the stochastic opening channels was simulated using Gillespie's stochastic integration algorithm modified for multiple ion channel populations and channels with multiple states (Gillespie, 1976; Gibson and Bruck, 2000). The algorithm calculates the time interval until the next state transition, determines the reaction type, i.e., which one of the possible state transitions occurs, and updates the number of channels in each possible state accordingly. As the maximum time step we used 1 μs and integrated all deterministic differential equations of the model with a forward Euler step.

## Results

To study the underlying noise source from the ISI statistics we stimulated locust auditory receptor neurons with many repetitions of 500 ms pure tones of constant intensity. Such stimuli evoked a spike response, which adapted over time (Römer, 1976; Benda et al., 2001; Fig. 2*A*,*B*, bottom). The initial high spike frequency directly after the onset of the acoustic stimulus declined over a duration of several tens to hundred milliseconds to its steady-state value. A sufficiently long interval of silence between the successive acoustic stimuli assured that the cell completely recovered from adaptation.

The steady-state spike frequency computed from the final 300 ms of the response increased in a sigmoidal fashion with sound intensity and saturated at frequencies of ∼160–200 Hz (Fig. 3*A*,*B*). Although the sensitivity, i.e., the minimum sound intensity required to evoke a response, varied considerably between the receptor neurons, the overall shape of the *f-I* curves (steady-state spike frequency vs sound intensity) was very similar (Fig. 3*B*).

Spike trains evoked by the tones with constant sound intensity showed notable spike-time variability (Schaette et al., 2005) (Fig. 2*A*,*B*, top) that was most prominent at low spike frequencies (Fig. 2*A*, top). The corresponding ISI histograms were broad at low sound intensities (i.e., low spike frequencies) and got much narrower at higher spike frequencies (Fig. 2*C–E*).

To quantify the dependence of the ISI variability on sound intensity we computed both the coefficient of variation, *CV*, as well as the diffusion coefficient, *D*. The CVs started out with values up to 0.9 for sound intensities at the lower end of the dynamic range (Fig. 3*C*,*D*) where spike frequencies are <50 Hz (compare Fig. 3*B*). With increasing sound intensity the CVs monotonically decreased and saturated at low values of ∼0.1–0.2. The maximum in the CV-versus-relative sound intensity curves visible in a few cells and in the average of 12 cells is not significant (note the large error bars) because of the low number of available ISIs at these low spike frequencies. The diffusion coefficient, *D*, in contrast, showed a pronounced peak at the lower end of the receptor cells' dynamic range (Fig. 3*E*,*F*). The peak height of *D* varied from cell to cell (Mean: 9.36 Hz ± 4.87 Hz SD) while the overall shape of *D* was rather similar (Fig. 3*F*).

### ISI distributions

What is the origin of the spike-time variability? To answer this question we first compared the observed ISI statistics with canonical noise-driven neuron models.

Given the high firing rates and the nonexponential but still unimodal ISI distributions (Fig. 2*C–E*), the receptor neurons most likely operate in their suprathreshold periodic spiking regime. This regime is characterized by a limit cycle dynamics that is distorted by some noise. The canonical model for this regime is the PIF model. This is the simplest model capturing the essence of limit cycle dynamics to which all neurons and other models of neuronal spiking converge, provided constant input current and sufficiently weak additive noise. This model is thus in no way specific for the grasshopper auditory receptor neurons, but is rather the most generic possible model for suprathreshold periodic firing. With additive white-noise driving (a signal with equal spectral power at all frequencies) the resulting ISI distribution is described by an “inverse Gaussian” density *p _{wn}* (Eq. 1) (Gerstein and Mandelbrot, 1964), in the following referred to as white-noise density. With additional slow adaptation this still holds true although the parameters of the white-noise density are rescaled (Schwalger et al., 2010). Thus, we expect the ISI density

*p*to describe the measured ISI histograms, if the corresponding noise source is close to a white noise and if the noise is sufficiently weak.

_{wn}The white-noise ISI density Equation 1 is a function solely of the mean and the SD of the ISIs, which can be determined directly from the experimentally measured spike responses. The resulting white-noise ISI density *p _{wn}* matches the ISI histograms of spike responses with low spike frequency (Fig. 4

*A*, dashed line). However, for increasing spike frequencies,

*p*no longer provides a good description of the data (Fig. 4

_{wn}*B*,

*C*, dashed line). The ISI histograms of the auditory receptor cells exhibit a more prominent peak and a heavier tail compared with

*p*.

_{wn}Replacing the additive white noise of the PIF model by a colored noise, in particular an Ornstein–Uhlenbeck noise, results in the colored-noise ISI density *p _{cn}* (Eq. 2) (Lindner, 2004). In contrast to

*p*, the shape of

_{wn}*p*depends on one free parameter, the correlation time constant τ of the colored noise. Fitting

_{cn}*p*to the ISI histograms of the recorded data, the resulting colored-noise ISI distributions describe the data remarkably well at both low and high spike frequencies, although there is only one fit parameter, i.e., τ (Fig. 4, solid line). In the upper range of spike frequencies τ converges to the adaptation time constant, which is in the range of 100 ms in the locust auditory system (Figs. 4

_{cn}*C*, 5

*A*). For lower spike frequencies τ vanishes (Figs. 4

*A*, 5

*A*), so that the colored-noise ISI density in fact approaches the white-noise ISI density (Lindner, 2004). The

*D*values of the Kolmogorov–Smirnov (K-S) test (Fig. 5

*B*) also show that for spike-frequencies below ∼60 Hz

*p*describes the data as well as

_{wn}*p*. At higher spike frequencies, however, the description of the data by

_{cn}*p*becomes worse and

_{wn}*p*is significantly better (Wilcoxon test,

_{cn}*p*< 0.05).

Following Schwalger et al. (2010), we further quantified the transition from a white-noise ISI distribution to a colored-noise ISI distribution with non-zero correlation time by computing rescaled variants of the skewness (α* _{s}*, Eq. 4) and kurtosis (α

*, Eq. 5). These quantities equal unity for an inverse Gaussian distribution. Values greater than one indicate a more skewed distribution with a stronger peak and a heavier tail as is the case for the colored-noise ISI distribution. Indeed, at low spike frequencies, both α*

_{e}*and α*

_{s}*are close to one in agreement with a dominating white-noise source (Fig. 5*

_{e}*C*). With increasing spike frequency both α

*and α*

_{s}*monotonically increase to values larger than one further supporting our observation of a transition from an inverse Gaussian to the more peaked colored-noise distribution arising from a dominating colored-noise source.*

_{e}The shapes of the ISI histograms thus suggest the presence of two different types of noise sources in the auditory receptor neurons of locusts: a white-noise source dominating at low spike frequencies and a colored-noise source, the effect of which is more prominent at higher spike frequencies.

### ISI correlations

Since the receptor neurons exhibit spike-frequency adaptation, we further investigated correlations between successive ISIs that may have resulted from the interaction of the noise with the adaptation dynamics. In particular, we expect *negative* correlations for a deterministic adaptation process interacting with a spike generator driven by a white-noise input (Chacron et al., 2001; Liu and Wang, 2001; Schwalger et al., 2010) and *positive* correlations for a stochastic adaptation process interacting with a deterministic spike generator (Schwalger et al., 2010).

For spike frequencies above ∼50 Hz we found positive correlations between successive ISIs that tend to increase with increasing spike frequency (Fig. 6*A*, black line). The serial correlation coefficients at lag 1, ρ_{1}, were obtained from the responses to the 500 ms stimuli as used before for investigating the ISI distributions. Especially at higher spike frequencies the responses are not stationary (Fig. 2*B*). Over the time of the stimulus the mean spike frequency slowly declines. This nonstationarity introduces positive serial correlations that may obscure correlations arising from a potential colored-noise source. We therefore computed the ISI correlations of a PIF model without spike-frequency adaptation with the time course of the input current adjusted such that the resulting time course of the mean spike frequency matches the experimentally measured one and an additive white-noise that ensured the same CV of the ISIs as observed in the experimental data. The resulting ISI correlations (Fig. 6*A*, gray line) are significantly lower than those obtained from the data. Thus the observed positive ISI correlations are likely to arise in part from a colored-noise source and not solely from the slow decline of the spike frequency.

For measuring the ISI correlations at lower spike frequencies we used acoustic stimuli of 1 min duration to obtain sufficiently large ISI statistics (≥300 ISIs per trial). At low sound intensities the mean firing rate remained constant for the whole stimulus duration after an initial transient of a few hundred milliseconds (Fig. 6*B*), which we discarded for the computation of correlations. Thus, the resulting ISI correlations are not obscured by slow changes of the mean spike frequency. Furthermore, the ISI histograms obtained from the spike responses to these long stimuli (Fig. 6*C*) are also well described by the white noise ISI density *p _{wn}* (Eq. 1).

Using the long stimuli we indeed found slightly negative ISI correlations (Fig. 6*A*). We tested the significance of these correlations by first calculating the distribution of serial correlation coefficients occurring by chance after shuffling the ISIs of each spike train of a given sound intensity (Fig. 6*D*, gray bars). The correlation coefficient of the non-shuffled data (Fig. 6*D*, solid line) was considered as significantly negative if the integral over the shuffled distribution (cumulative histogram, Fig. 6*D*, dashed line) up to this value was smaller than 0.05. At spike frequencies in the range of 20–45 Hz most responses indeed showed significant negative ISI correlations (Fig. 6*E*).

The results on ISI correlations exactly match the expectations based on the ISI distributions we measured. At low spike frequencies a white-noise source shapes the ISI distribution and results in interaction with an adaptation process in negative ISI correlations. With increasing spike frequency we observe a transition to a regime that is dominated by colored noise with the respective ISI distributions and the corresponding positive ISI correlations. This suggests that the dominating biophysical mechanisms generating the intrinsic noise change with increasing spike frequency.

### Modeling

To uncover potential biophysical mechanisms causing the observed ISI variability we used a model of the locust auditory receptor neuron and compared the ISI statistics resulting from various possible noise sources with the experimentally measured data. In the model (Fig. 1*A*; see Methods for a detailed description), sound waves cause vibrations of the tympanal membrane. As a result mechanosensory ion channels in the membrane of auditory receptor neurons open and close and thereby transform the mechanical vibrations into electrical currents. The resulting change of the receptor cell's membrane potential activates voltage-dependent ion channels that finally trigger action potentials.

Locust auditory receptor cells differ in their sensitivity from cell to cell but have comparable variability measures, *CV* and *D*, as a function of spike frequency (Fig. 3). This excludes extrinsic noise and tympanal noise (see also below) as possible noise sources. Furthermore, synaptic noise can be ruled out because locust auditory receptors are bipolar neurons that are directly attached to the eardrum (Gray, 1960) and generate action potentials without intermediate synapses. Thermal noise due to the membrane resistance was shown to have only a minor importance in neurons compared with other noise origins (see below). Therefore, we expect that intrinsic noise caused by the random gating of ion channels is the main source of noise (Fig. 1*B*). Especially for small channel populations fluctuations between open and closed states of single ion channels give rise to fluctuations of the ionic currents mediated by the channel population and cause membrane potential fluctuations as well as variability in spike timing.

For a minimal single-compartment conductance-based model of the locust auditory receptor neuron four different types of ionic currents need to be included: (1) the mechanosensory receptor current transducing the sound wave into an electrical signal, the (2) sodium and (3) delayed-rectifier potassium current for action potential generation, and (4) an adaptation current mediating spike-frequency adaptation. In the following simulations we explicitly model the stochastic channel dynamics of each of these currents and investigate whether the resulting stochasticity of the ionic currents might explain the experimentally measured ISI statistics.

#### Single current stochasticity

In a first set of simulations we studied the effects of ion-channel stochasticity for each of the involved ionic currents individually. That is, in each of the simulations either (1) the receptor current, *I _{R}*; (2) the M-type adaptation current,

*I*; or the voltage-dependent (3) sodium,

_{M}*I*, or (4) potassium current,

_{Na}*I*, was modeled as a stochastic current resulting from a finite population of ion channels, while all other ionic currents were modeled as deterministic currents corresponding to an infinite population of ion channels. For each case we varied the number of channels carrying the respective current and at the same time scaled the respective single-channel conductances inversely proportional to the number of channels. This resulted in approximately the same spike frequencies but different ISI variability as quantified by both the coefficient of variation,

_{K}*CV*, and the diffusion coefficient,

*D*.

All stochastic ionic currents resulted in a high *CV* for stimuli close to the firing threshold that declined with increasing sound intensity (data not shown). This finding is consistent with the observed statistics of the spike responses of the auditory receptor cells.

In contrast, *D* showed qualitative differences between the different stochastic ionic currents (Fig. 7). Stochastic receptor channels (Fig. 7*A*) generated a peak in *D* similar to the measured experimental data (Fig. 3*E*,*F*). The peak occurs at approximately the sound intensity inducing a tympanal displacement which causes an opening of the receptor channels with a mean probability of 0.5 (as defined by Eq. 12). For this open probability the variability of the ionic current and thus of the ISIs is maximal (Sigworth, 1980).

Voltage-gated sodium channels are responsible for action potential initiation. At membrane potentials close to the spike threshold, a random opening of channels can easily initiate an action potential (Schneidman et al., 1998). Similar to the model with stochastic receptor channels, sodium channel noise resulted in a peak of *D* in the lower dynamic range of the receptor's response (Fig. 7*B*).

A stochastic delayed-rectifier potassium current responsible for the termination of an action potential, however, resulted in a sigmoidal increase of *D* with increasing sound intensity (Fig. 7*C*). This phenomenon is not in agreement with the experimental data.

For modeling the adaptation current we used a population of voltage-gated M-type potassium channels. These channels are mainly activated during action potentials and slowly deactivate between them. Stochasticity of the adaptation channels again can reproduce the experimentally measured values of *D* with a peak in the lower dynamic range (Fig. 7*D*).

To distinguish between stochastic currents showing similar ISI variability we investigated the shapes of the ISI distributions as well as ISI correlations next. For the analysis of the ISI histogram shapes, we adjusted the size of the stochastic channel populations such that all simulations with the different stochastic ionic currents resulted in the same steady-state spike frequency as well as the same *CV* and *D*. Stochastic simulations of ionic currents with fast gating dynamics, i.e., sodium and receptor current, yielded ISI distributions that are well described by the white-noise density *p _{wn}* (Eq. 1; Fig. 8

*A*). Simulations with stochastic M-type adaptation currents, however, resulted in ISI distributions that are more peaked than

*p*. These distributions can be well described by the colored-noise density

_{wn}*p*(Eq. 2) with a correlation time constant τ which is similar to the adaptation time constant (Fig. 8

_{cn}*B*).

Serial correlations were determined from the ISIs of the simulations used for Figure 7. Both the fast stochastic receptor (Fig. 9*A*) and sodium currents (Fig. 9*B*) generated negative correlations at lag one. The negative correlations are strongest at spike frequencies of ∼20 Hz. In contrast, the slow stochastic adaptation current lead to positive correlation coefficients ρ_{1} that increased with higher spike frequencies (Fig. 9*C*).

In addition to the various channel-noise sources, stochastic responses may result from thermal noise of the eardrum. We simulated such eardrum fluctuation using additive Gaussian white noise on the sound pressure wave. The white noise together with the signal is filtered by the tympanum that we modeled as a damped harmonic oscillator. Sufficiently large noise intensities can lead to ISI statistics that are comparable to those caused by the stochastic receptor or sodium current. However, at these noise intensities the resulting variance of the tympanum deflections is already three orders of magnitude higher than the variance caused by the sound wave. Using the Einstein relation and a minimal effective area of 1% of the total area of the tympanum leads to a rough estimate according to which the variance of the thermal noise is at least half an order of magnitude lower than the signal variance. We have checked that with this value, the effect of thermal noise on the ISI statistics is too small to explain the experimental data.

#### Mixed channel noise

The simulations of single-current noise showed that stochastic M-type adaptation currents can explain both the ISI distributions and positive correlations of the experimentally measured spike responses in the upper dynamic range. In contrast, stochastically gated sodium as well as receptor channels reproduced the ISI distributions in the lower dynamic range. In interaction with the deterministic adaptation current these two channel types also reproduced the observed negative correlations.

This suggests that at least two processes are needed to explain the experimental data: (1) a process generating fast current fluctuations, which dominate at low sound intensities, and (2) a slow adaptation process both for generating negative correlations in interaction with the fast noise and for generating slow current fluctuations, which dominate at high sound intensities. To test this hypothesis we simulated the model with a stochastic M-type adaptation current in combination with stochastic currents mediated by either voltage-dependent sodium or mechanosensory receptor channels.

Both mixed-case models explained the ISI statistics of the auditory receptor cells' spike responses (Fig. 10). The ISI variability quantified by *D* showed a peak at the lower dynamic range of the receptor's response. The coefficient of variation revealed a high variability for low steady-state spike frequencies, which declined with increasing rate (Fig. 11*A*). In addition to the *CV*, both the rescaled skewness α* _{s}* and kurtosis α

*confirmed the transition from an inverse Gaussian to the more strongly peaked colored-noise distribution (Fig. 11*

_{e}*B*). For low spike frequencies, α

*and α*

_{s}*took values close to one indicating a strong white-noise influence. With increasing spike frequency, α*

_{e}*and α*

_{s}*showed a monotonic rise to large values. This implies that the ISI density is more skewed and more peaked compared with an inverse Gaussian, which illustrates the prevalence of colored noise in this firing-rate regime. Furthermore, with increasing spike frequency, ISI correlations exhibited a transition from negative values resulting from predominant white noise driving interacting with the adaptation current to positive values indicating strong impact of colored noise (Fig. 11*

_{e}*C*).

## Discussion

This study aims at characterizing intrinsic neuronal noise sources solely from the resulting ISI statistics. Our data from *in vivo* axonal recordings of locust auditory receptor neurons show that with increasing spike frequency the shape of the ISI histograms changes from an inverse Gaussian to a more peaked probability density with a heavier tail, whereas the corresponding ISI correlations go from negative to positive values.

In response to pure tones the auditory receptor neurons exhibit firing rates of up to a few hundred Hertz (Römer, 1976). Together with the nonexponential ISI distributions (Schaette et al., 2005) this strongly indicates that these neurons operate in the suprathreshold regime, i.e., far beyond a bifurcation from quiescence to tonic (repetitive) firing. The canonical model for such limit-cycle dynamics is the PIF model; to first approximation, any tonically firing, real or model, neuron can be described by this simple model. When driven by additive white Gaussian noise, the PIF's ISI distribution is given by an inverse Gaussian (Eq. 1) (Gerstein and Mandelbrot, 1964). Vice versa, observing an inverse-Gaussian ISI distribution of a periodically spiking neuron in the suprathreshold regime suggests that the dominating noise is Gaussian and effectively white, i.e., has a cutoff frequency much larger than the firing rate.

The locust's receptor neurons also show pronounced spike-frequency adaptation (Römer, 1976; Benda et al., 2001; Gollisch and Herz, 2004). Despite this additional slow dynamics the ISI density of a PIF neuron model with additive white noise can be approximated by an inverse Gaussian (Schwalger et al., 2010), but the ISI correlations become negative (Wang et al., 1998; Chacron et al., 2000; Liu and Wang, 2001; Benda et al., 2010; Schwalger et al., 2010)—a feature that is absent in a white noise-driven PIF model without adaptation. Negative ISI correlations have been reported for many different preparations (for review, see Farkhooi et al., 2009; Avila-Akerberg and Chacron, 2011), and have been suggested to play a role in reducing long-term variability (Ratnam and Nelson, 2000; Chacron et al., 2001, 2007).

At low spike frequencies our recorded ISI distributions agree with inverse Gaussians and ISI correlations were small but significantly negative. Both findings suggest that in this regime the dominating noise source generates white Gaussian noise that interacts with a deterministic adaptation dynamics.

Colored Gaussian noise like Ornstein–Uhlenbeck noise with sufficiently long correlation times results in different ISI statistics. The ISI distribution is more peaked and has longer tails compared with the inverse Gaussian and the noise correlations directly translate into positive correlations between successive ISIs (Lindner, 2004; Schwalger et al., 2010). Positive correlations have been observed in neural systems (Lowen and Teich, 1992; Gabbiani and Krapp, 2006), but the origin of the corresponding colored noise was unknown or hypothesized to be of synaptic origin (Chacron et al., 2001). Recently, Schwalger et al. (2010) demonstrated analytically that stochastic adaptation currents also generate Ornstein–Uhlenbeck noise with a correlation time proportional to the time constant of the adaptation dynamics.

The peaked shape of the ISI distribution as well as the positive ISI correlations observed at higher spike frequencies may thus be caused by a colored-noise source in the adaptation current. The additional match of the correlation time constant of the noise obtained from the ISI distributions (Fig. 5*A*) with the time constant of spike-frequency adaptation (Fig. 2*A*,*B*) supports this hypothesis.

Our simulations of conductance-based models further support our findings that a stochastic adaptation current together with fast channel noise indeed reproduces both the white-noise ISI distributions and negative ISI correlations at low spike frequencies and the transition to colored-noise ISI distributions and positive correlations at higher spike frequencies. The strength of positive ISI correlations and the deviations of the ISI density from an inverse Gaussian grow as the ratio between adaptation time scale and mean ISI becomes larger with increasing firing rate (Schwalger et al., 2010). Additionally, in our model the mean fraction of open adaptation channels increases with spike frequency from values close to zero to at most 10%, as measured in the simulations of the modified Traub–Miles model (Ermentrout, 1998). Because the SD of the fraction of open channels is zero when all channels are closed and monotonically rises until half of the channels are open (a state that is never reached here), the fluctuations of the adaptation current increase with spike frequency and so does the strength of the resulting colored noise. In marked contrast, the level of fluctuations for the fast noise sources saturates at large signal amplitude or high firing rates (simulation data not shown). We note that a stochastic adaptation current could also explain very similar data on CV and ISI correlations found in response to current pulse injections in a looming-sensitive neuron (Gabbiani and Krapp, 2006), in particular, the observed transition from negative ISI correlations at low rate and high CV to positive ISI correlations at high firing rate and low CV.

Which of the many potential adaptation mechanisms, such as M-type or AHP-type currents (Sah and Davies, 2000) or slow inactivation of the sodium current (French, 1987), is causing spike-frequency adaptation and the colored noise in the locust auditory receptor neuron is not known and cannot be deduced by our indirect methods that only use the ISI statistics. For simulations we took M-type currents as a prototypical adaptation mechanism, because in the suprathreshold regime, the exact origin of the mean adaptation effect and the adaptation noise does not matter as long as both are associated with the slow part of the adaptation dynamics (Benda and Herz, 2003; Schwalger et al., 2010).

What, however, is a possible source of the white noise? Excluding thermal noise due to the membrane resistance (Johnston noise; Manwani and Koch, 1999) or fluctuations of the eardrum (see Results) leaves channel noise resulting from random opening and closing of ion channels as the most likely candidate for the white noise (White et al., 2000). In our simulations channel noise originating from the delayed rectifier failed to reproduce the peak in *D*. However, stochastic receptor currents as well as stochastic sodium currents reproduced the experimental data. This agrees with other studies, which demonstrate the dominance of sodium-channel noise over potassium-channel noise in (1) shaping spike-timing variability (Chow and White, 1996; Schneidman et al., 1998), (2) placing limits on axon miniaturization (Faisal et al., 2005), or (3) driving slow membrane potential oscillations (Dorval and White, 2005).

Assuming that the white noise can be mainly attributed to stochastic receptor currents, the simulations allow one to estimate the numbers of ion channels of the receptor current from the ISI statistics. This will be a rough estimate because, for example, decreasing the number of receptor channels can be partially compensated by increasing the number of channels of the adaptation current without significantly changing the ISI statistics. Still, the number of receptor channels we used in the mixed-case model (*N _{R}* = 20; Fig. 11) matches the reported order of magnitude of channel numbers of other mechanosensory transducer currents (

*Drosophila*:

*N*≈ 20 per sensory neuron, Nadrowski et al., 2008; bullfrog:

*N*= 85 ± 45 per hair bundle, Howard and Hudspeth, 1988).

Especially at higher firing rates the ISI variability is quite small in locust auditory receptor cells, both by relative (*CV* ≈ 0.1) and absolute measures (

Direct measurements of the current noise close to the site of receptor current and spike generation potentially damage the delicate sensory structures. In contrast, indirect methods that deduce properties of the transduction cascade from recordings of the auditory nerve that leave the ear intact allow one to precisely investigate functional aspects of sensory transduction and set detailed references for future invasive measurements. For example, using the iso-response method, Gollisch et al. (2002) and Gollisch and Herz (2005) have been able to characterize the auditory transduction cascade of locust auditory receptor neurons and to differentiate between input- and output-driven adaptation processes in these neurons (Gollisch and Herz, 2004).

Using auditory-nerve recordings, we measured ISI statistics far away from the site of spike generation. Using generic models we deduced from these data key properties of the underlying noise sources. This indirect approach for dissecting intrinsic noise sources can be generalized to other neurons as well as to the important case of external synaptic noise.

## Footnotes

This work was supported by BMBF Bernstein Collaboration for Computational Neuroscience 01GQ0722 (K.F.) to A.V.M.H., the Bernstein Center for Computational Neuroscience Munich LMU 01GQ1004A (K.F.) to A.V.M.H., the BMBF Bernstein Center for Computational Neuroscience Berlin HU 01GQ1001A (B.L.), and the BMBF Bernstein Award Computational Neuroscience 01GQ0802 (J.B.).

- Correspondence should be addressed to Jan Benda, Institute for Neurobiology, Eberhard Karls Universität Tübingen, Auf der Morgenstelle 28 E, 72076 Tübingen, Germany. jan.benda{at}uni-tuebingen.de