Threshold fatigue and information transfer

2.1 Experimental recordings

We used the weakly electric fish Apteronotus leptorhynchus in this study and the experimental protocol has been described in detail previously (Bastian et al. 2002; Chacron 2006). Briefly, animals were immobilized by intramuscular injection of curare and respirated with aerated water at a flow of 10 ml/min. Intracellular recordings from pyramidal cells in the electrosensory lateral line lobe were made with glass micropipettes (resistance 25–40 MΩ) filled with 3 M KCL. Recording sites as determined from surface landmarks and depth were limited to the centrolateral and lateral segments only. Sensory stimulation was achieved by presenting amplitude modulations (AMs) of the animal’s own quasi-sinusoidal electric organ discharge via two electrodes located 20 cm on each side of the animal. The AM consisted of low-pass filtered Gaussian white noise (120 Hz cutoff, 8th order Butterworth filter). All procedures were approved by McGill University's animal care committee.

2.2 Modified Hodgkin–Huxley model

We used a Hodgkin–Huxley model with leak, sodium, and potassium conductances. The model is described by the following equations:

Here V is the membrane potential expressed in mV, g_x is the maximum conductance of channel x and E_x is its reversal potential, ξ(t) is Gaussian white noise with zero mean and autocorrelation function <ξ(t) ξ(t+τ)>=σ₂ δ(τ), I is the applied current, and S(t) is the stimulus which consists of low-passed filtered Gaussian white noise with spectral height α and cutoff frequency f_c. The infinite conductance curves are given by:

Usually, V_1/2 is taken to be constant (Hodgkin and Huxley 1952). We instead make V_1/2 vary with time:

where t_i is the timing of the ith action potential. This implies that V_1/2 is increased by 12 mV immediately after each action potential and is then allowed to relax back to its equilibrium value of −55 mV. The differential equation for V_1/2 is exactly the same as that governing the time varying threshold in more phenomenological models (Geisler and Goldberg 1966; Chacron et al. 2001b; Liu and Wang 2001). This is motivated by the fact that V_1/2 is closely related to the action potential threshold as measured in this study since it is the voltage at which roughly 50% of sodium channels are open.

We used an Euler–Maruyama integration scheme for numerical simulations with time step dt=25 μs. Parameter values used were: g_na=55 μS, g_k=40 μS, g_leak=0.18 μS, E_leak=−70 mV, E_na =40 mV, E_k =−88.5 mV, τ_m =0.02 ms, τ_h =τ_n =0.39 ms, τ_V=40 ms, σ=I=1 nA. The model gave rise to sustained action potential firing in the absence of noise (i.e. σ=0).

2.3 Data analysis

We computed the threshold value for spike time t_i+1, θ_i, as the membrane potential value at the time at which the second derivative of the membrane potential was maximum (see Fig. 1). We note that, while this procedure is different than the one previously proposed (Azouz and Gray 1999, 2000), our results do not qualitatively depend on which of the two algorithms is used to compute the action potential threshold (data not shown).

The interspike interval sequence (ISI), the sequence of times between consecutive action potentials, {I_i}≡{t_i+1−t_i} was obtained from the spike time sequence {t_i}. We computed the ISI serial correlation coefficients ρ_j as:

Where the average <…> is taken over index i. From the ISI sequence, we obtained the instantaneous firing rate sequence as { f_i}≡{1/I_i} and the cross-correlation coefficient between {θ_i} and {f_i} was then computed as:

For signal transmission analysis, the stimulus S(t) was sampled at 2 kHz and we made a binary representation R of the spike train at 2 kHz also. We then computed the stimulus-spike train coherence as (Rieke et al. 1996):

where P_RS ( f ) is the cross-spectrum between the stimulus and the binary sequence and P_RR( f ), P_SS ( f ) are the power spectra of the binary sequence and stimulus, respectively. The coherence ranges between 0 and 1 and measures the correlation between the stimulus and neural response at frequency f. The coherence is furthermore related to the mutual information density between the stimulus and response by: I=−log₂[1−C_RS ( f )] (Rieke et al. 1996). We normalized the mutual information density I by the firing rate to account for its known dependence on firing rate (Borst and Haag 2001).

We used a threshold interspike interval value of 20 ms for classifying spikes as either being part of a burst or not (Gabbiani et al. 1996; Metzner et al. 1998; Chacron et al. 2004c; Oswald et al. 2004). As such, an ISI was considered to be part of a burst only if it was smaller than the threshold and we constructed from the binary sequence R a binary sequence R_bursts with all the spikes that belonged to bursts.

3.1 Sensory neurons display threshold fatigue

Figure 1 illustrates the method to compute the action potential threshold. We took the action potential threshold as the point of greatest curvature of the membrane potential V as determined by the local maximum of the second derivative. We then applied this algorithm to intracellular recordings from pyramidal cells in the electrosensory lateral line lobe of weakly electric fish. An example cell is shown in Fig. 2. Visual inspection of the membrane potential (Fig. 2(a)) shows that the action potential threshold is variable and ranges between −62 and −57 mV for this particular cell (Fig. 2(b)). A more careful visual inspection shows that the action potential threshold actually increases during periods of rapid firing and decreases during periods of silence (Fig. 2(a)). This is most easily seen by comparing the time evolution of the threshold to that of the instantaneous firing rate (Fig. 2(c)). It is seen that periods of elevated firing are associated with higher threshold, leading to a positive correlation coefficient (Fig. 2(d)) as well as negative correlation coefficient in adjacent lags. This neuron also showed a prominent negative ISI correlation coefficient at lag one (ρ₁=−0.2988) (Fig. 2(d), inset).

3.2 Threshold fatigue is most prominent in neurons with high spontaneous firing rates

Previous studies have shown that the pyramidal cell population was heterogeneous both anatomically and physiologically. There is a strong negative correlation between spontaneous firing rate and dendritic morphology: cells with small dendritic trees tend to have large firing rates while cells with large dendritic trees tend to have small firing rates (Bastian and Nguyenkim 2001; Bastian et al. 2004). Figure 3(a) assesses the influence of firing rate on threshold fatigue as quantified by the cross-correlation coefficient between threshold and reciprocal interspike interval (ISI) at lag 0 C₀. A strong positive correlation (R=0.64, p=0.00163, N=21) was observed, meaning that threshold fatigue is more likely to be observed in neurons with large firing rates. Previous modeling studies have shown that threshold fatigue can be associated with negative ISI correlations at lag 1 (i.e. short interspike intervals are on average followed by long ones and vice-versa) (Chacron et al. 2000, 2001b, 2003a; Liu and Wang 2001). We thus plotted the ISI correlation coefficient at lag 1 ρ₁ as function of firing rate (Fig. 3(b)). A strong negative correlation between ρ₁ and firing rate was observed (R = −0.823, p<0.0001, N=21), indicating that cells with high spontaneous firing rates also tended to display negative ISI correlations. Thus, as predicted from theory (Chacron et al. 2004a; Lindner et al. 2005), cells that displayed threshold fatigue also tended to display a negative ISI correlation coefficient at lag 1.

3.3 A biophysical model of threshold fatigue

Although many phenomenological models have proposed variations in the action potential threshold (Geisler and Goldberg 1966; Gestri et al. 1980; Gabbiani and Koch 1996; Chacron et al. 2001b, 2004a), it was only recently shown that more biophysically realistic models such as the Hodgkin–Huxley model displayed threshold variability (Azouz and Gray 1999). We thus tested whether a realistic biophysical model with Hodgkin–Huxley conductances, which has been shown to display threshold variability (Azouz and Gray 1999), could display threshold fatigue. The results are shown in Fig. 4. Although the model does display variability in the action potential threshold (Fig. 4(a,b)), there are no cumulative increases in threshold during rapid firing and there is no significant correlation between the instantaneous firing rate and the threshold (Fig. 4(c,d)). The model furthermore displayed no significant ISI correlations at lag 1 (ρ₁=−0.00152) (Fig. 4(d), inset).

Previous studies have shown that pyramidal cells had a well characterized burst mechanism that relies on an interaction between somatic and dendritic voltage-gated sodium channels (Lemon and Turner 2000; Doiron et al. 2002; Fernandez et al. 2005). Although burst dynamics can sometimes give rise to negative ISI correlations (Chacron et al. 2001a), this particular mechanism does not seem to be sufficient to generate the negative ISI correlations and threshold fatigue seen in our experimental data as revealed by simulations of the model proposed by Doiron et al. (2002) (data not shown).

Thus, we explored additions to the Hodgkin–Huxley model that could reproduce the experimental data. We varied the inflexion point of the activation and inactivation infinite conductance curves of the sodium conductance, V_1/2, in a dynamic manner similar to the threshold in phenomenological models (Chacron et al. 2001b, 2003a; Liu and Wang 2001). The results are shown in Fig. 5. With this minor modification, the model exhibits threshold fatigue: the threshold increases during repetitive firing (Fig. 5(a)) and displays threshold variability that is similar to the experimental data (Fig. 5(b)). There was a strong positive correlation between the instantaneous firing rate and the threshold (Fig. 5(c,d)) as in the experimental data. The model furthermore displayed a significant ISI correlation coefficient at lag 1 (ρ₁=−0.3349) that was similar in magnitude to that seen experimentally (Fig. 5(d), inset). We note that the model does not reproduce the damped oscillation in the cross-correlation observed for the data, such a damped oscillation might be due to underlying network dynamics (Doiron et al. 2003). However, a simple point-process model (i.e. no spatial structure) is sufficient to qualitatively explain the threshold fatigue and negative ISI correlations seen in pyramidal cells.

3.4 Effects of firing rate on threshold fatigue in the model

In order to vary the firing rate in the model with threshold fatigue, we varied the bias current I. Figure 6 shows the effects of varying the bias current on the ISI correlation coefficient at lag 1, and the cross-correlation coefficient between the instantaneous firing rate and the action potential threshold. Both the ISI correlation coefficient (Fig. 6(a)) and the cross-correlation coefficient (Fig. 6(b)) increased in magnitude as a function of firing rate thereby reproducing experimental results (Fig. 3). As such, our model qualitatively explains the variations seen in experimental data and predicts that these are predominantly due to differences in firing rate amongst pyramidal cells. The effects seen can be explained as follows: reducing the bias current reduces the firing rate, thereby increasing the mean ISI. There is therefore less accumulation in threshold during rapid firing, which leads to a decrease in the ISI correlation coefficient (ρ₁) and the cross-correlation coefficient (C₀). The decrease in ρ₁ as a function of firing rate is similar to previous results obtained with phenomenological models (Liu and Wang 2001; Chacron et al. 2003a).

3.5 Threshold fatigue regularizes firing and confers robustness to noise

We quantified the effects of threshold fatigue by comparing the performance of the model with dynamic V_1/2 to that of the model with constant V_1/2, henceforth referred to as the model without threshold fatigue. Overall, the model with threshold fatigue produced spike trains that were more regular than those produced by the model without threshold fatigue as seen by looking at the ISI distribution (Fig. 7(a)). The coefficient of variation (CV, standard deviation to mean ratio of the ISI distribution) was 0.21 with threshold fatigue and 0.44 without. The power spectra of spike trains (Fig. 7(b)) also indicate less output noise in the case of threshold fatigue: power at low frequencies is reduced and the peak at the neuron’s firing rate (about 30 Hz) becomes much sharper. These effects are a consequence of diminished ISI variability (lower CV) and negative correlations among ISIs ${\sum }_{i=1}^{\infty }{\rho }_j<0$ both caused by threshold fatigue. In particular, in the low frequency limit it is known that $S(f=0)={\text{CV}}^{2}r\left[1+2{\sum }_{i=1}^{\infty }{\rho }_j\right]$ for a stationary spike train (Cox and Lewis 1966; Holden 1976), hence both a decrease in CV and the presence of negative correlations lead to a drop in low-frequency power. We note that these results are similar to previous ones using more phenomenological models (Chacron et al. 2004a,b; Lindner et al. 2005) as well as experimental results (Chacron et al. 2005b).

We also computed the area A( f ) between 0 and frequency f under each power spectrum. The areas computed with and without threshold fatigue are equal for f >60 Hz, confirming that this is indeed a reshaping of the power spectrum (Fig. 7(b), inset) (i.e. the total power is conserved). The spike count distribution (Barlow and Levick 1969a,b; Teich and Khanna 1985) also displayed less variance with threshold fatigue (Fig. 7(c)). We also varied the noise standard deviation σ and found that spike train variability, as quantified by CV, was generally lower in the presence of threshold fatigue (Fig. 7(d)). As such, the model with threshold fatigue is more robust to noise than the model without threshold fatigue.

3.6 Threshold fatigue and information transmission

We finally quantified the effects of threshold fatigue on information transmission. However, the models with and without threshold fatigue gave rise to very different ISI probability densities (Fig. 7(a)). In particular, the probability density of high frequency ISIs (<20 ms) was lower for the model with threshold fatigue than for the model without threshold fatigue. Since previous studies have shown that high frequency ISIs coded for low frequency stimuli (Oswald et al. 2004), we compared information transmission by high frequency ISIs with and without threshold fatigue. The results are shown in Fig. 8 were the information densities normalized by the event rates are compared. It is seen that, for low (<20 Hz) and higher frequencies (38–60 Hz), the model with threshold fatigue had a larger information density than the model without threshold fatigue. The situation was opposite however in the intermediate frequency range (20–38 Hz). A stimulus with power only in this range would thus be best transmitted with a neuron without threshold fatigue as predicted from theory (Lindner et al. 2005). This confirms that the noise shaping seen in the power spectrum does translate to increased information transmission. Though threshold fatigue will give rise to less high frequency ISIs, each one will transmit a greater amount of information as there is less noise.