Dynamic Gain Decomposition Reveals Functional Effects of Dendrites, Ion Channels, and Input Statistics in Population Coding

Decomposing dynamic gain provides subcellular resolution for the analysis of dynamic population coding

Neurons that receive a common feed-forward input constitute a population that encodes this input into changes of their population firing rate Figure 1A₁. Their encoding capability can be quantified with a linear filter, i.e., the dynamic gain function G(f). The dynamic gain’s magnitude and phase capture, in a spectrally resolved manner, how modulations of the input current cause scaled (magnitude) and delayed (phase) modulations of the population’s firing rate. Because dynamic gain reflects population-level function but is also sensitive to cellular properties, it can serve as a central link between cell-level and network-level function (Brunel and Wang, 2003). With the following decomposition of the dynamic gain, we make this link more accessible.

• Download figure
• Open in new tab
• Download powerpoint

Figure 1.

Physical signal processing steps guide dynamic gain decomposition. A₁, A population encoding scheme. Feed-forward input (green) diverges onto a recurrently connected neuronal population. The dominating, asynchronous, recurrent input causes weakly correlated AP patterns. Nevertheless, the population rate can faithfully reflect changes in the common input. A downstream neuron can detect this signal, filtered through the population’s transfer function, the dynamic gain G(f). A₂, Experimentally determined dynamic gain function of a layer 5 pyramidal cell, obtained from 2,330 APs (data from Revah et al., 2019). Note the wide bandwidth, and the maxima and minima. Confidence interval (gray band) and noise floor (dashed) are obtained by bootstrapping (see Materials and Methods). B₁, Evolution of neuron models’ phase plots from point neurons without initiation dynamics to Eyal’s multi-compartment model. Locations of local minima during AP initiation are marked with dots and denoted as V_loc. Normalized dynamic gain functions of the first three model variants are shown next to the phase plots. B₂, Schematic representation of neuronal signal transformation. Distributed synaptic inputs drive current into the soma and axon; the axonal voltage (V_a) changes. For real neurons and biophysically plausible multi-compartment models, various parameters shape these transformations, attribution of dynamic gain changes to specific parameters is a challenge. C₁, Globally, Eyal’s model displays irregular firing patterns in response to stochastic stimuli. Locally, at AP onset, the initiation dynamics varies, signified by the large spread in the delay between the last positive crossing of V_loc and the AP detection voltage (0 mV). C₂, Decomposition of AP initiation process. Somatic input is first transformed into voltage fluctuations at the AP initiation site (V_a). V_a(t) contains the subthreshold dynamics related to the firing pattern (time intervals of neighboring spikes when their voltages reach V_loc), and the suprathreshold dynamics related to the precise, sub-millisecond AP timing (determined by the spike initiation dynamics after V_loc). This decomposition can disentangle the functional effects of neuronal properties on population encoding.

The decomposition is motivated by the simplifications that underlie the lineage of widely used simple neuron models (Fig. 1B₁). The simplest of these models, the LIF neuron 1, abstracts all suprathreshold dynamics as infinitely fast, and all subthreshold dynamics as a low-pass filter, i.e., an impedance Z(f) = R_mem/(1 + iωτ_mem). We can reformulate the LIF model’s dynamic gain function G(f) as a multiplication of a current-to-voltage transformation, captured by the impedance Z(f), and a voltage-to-firing rate transformation, the spike gain G_sp(f). To this end, we expand 6, the definition of the dynamic gain, with the Fourier transformation of the voltage F(V) :G(f)=|F(ν)F(I)|=|F(V)F(I)|⋅|F(ν)F(V)|=Z(f)⋅Gsp(f). (11)This first decomposition step is exact in LIF neurons, but how does it generalize to models with instrinsic AP initiation dynamics?

The EIF model adds an exponentially increasing current term (Eq. 2). In the subthreshold region, the voltage to current transformation becomes voltage dependent. The negative slope of the phase plot flattens, until the local minimum at V_loc (Fig. 2A). Following the analysis as described in the methods, we derived the average current-to-voltage transformation of the EIF model. We found that, in the limit of low firing rates, this effective impedance converges to the impedance of a LIF model. Specifically, under stimulation conditions for which the EIF model’s voltage attains an average value of 〈V〉, the EIF model’s effective impedance converges to the impedance of a LIF model with an effective resistance R′=R/(exp((ϑ−⟨V⟩)/ΔT)−1) in the limit of low firing rates (Fig. 2C, dark blue).

• Download figure
• Open in new tab
• Download powerpoint

Figure 2.

Dynamic gain decomposition of an EIF model retrieves basic building blocks corresponding to simpler LIF model variants. A, Five phaseplots that describe the dynamics of Equation 2 of the EIF model for four different current levels (I in pA: 139.3 (dark blue), 143.4, 146.3, 147.7, and 151.5 (red)). The inset magnifies the region around the local minimum. The arrows indicate the AP detection voltage used for the dynamic gain (0 mV) and the zero-delay dynamic gain (−45 mV). B, When driven with fluctuating currents with averages as in A, a standard deviation of σ_I = 15 pA, and a correlation time τ_corr = 25 ms, the EIF model yields the dynamic gain displayed here in continuous lines. The corresponding zero-delay dynamic gain curves (see Materials and Methods) are shown with dashed lines. The number of APs varies between the simulations (720,000 APs for 5 Hz (red), 120,000 APs for 1 Hz (dark blue), and 50,000 APs for the others, shown in fainter shades). Note how the steep decay of the dynamic gain curves turns to a flat high-frequency behavior for the zero-delay dynamic gain. C, The effective impedances for the simulations from B are displayed with continuous lines (calculated as described in Materials and Methods). The average voltage during these simulations were, in mV, −49.45, −48.86, −48.54, −48.41, and −48.17. The linear approximation of Equation 2, shown by the phaseplots in A, at these voltages are the first-order approximation of the current I to voltage V transformation. As described in the first section of the results, these linear approximations correspond to LIF models with modified values for membrane time constant and resistance. The impedances of these LIF models are shown with dashed lines in the color corresponding to the simulation that is approximated. Toward smaller and smaller firing rates, the numerically determined effective impedances tend toward the corresponding analytical impedances of the linear approximations around the average voltage.

• Download figure
• Open in new tab
• Download powerpoint

Figure 3.

Eyal’s model can realize high-bandwidth encoding only with a type II, but not a type I excitability. B, F–I curves with the rheobase thresholds aligned to 0 pA. C, To switch the excitability of Eyal’s model from type II to type I, the voltage dependence of sodium channel activation and inactivation is shifted 10 mV toward more hyperpolarized potentials. V_a denotes the voltage at the AP initiation site. D, 2D firing rate surface of Eyal’s type II model as a function of mean and std of stochastic stimulus (τ = 5 ms). The 5 Hz iso-firing rate line is labeled black. E, σ_I–μ_I and F CV_ISI–μ_I relation of the type II model at 5 Hz firing rate. Colored dots indicate three example working points, ranging from nearly mean-driven (orange) to fluctuation-driven (blue), and to very fluctuation-driven (green). Corresponding ISI distributions are shown in the inset panel. G, Dynamic gain functions of the three example working points. High-bandwidth encoding is realized only when the neuron model is fluctuation-driven (blue). H, through J as E through G, but for the type I model. J, The type I model does not realize high-bandwidth encoding. Note that an increased σ_I enhances the dynamic gain in the high-frequency region. Nevertheless, the cutoff frequency always remains below 16 Hz. The zero-delay gains of the type I model are shown in the Extended Data Figure 3-1 in Zhang et al. (2023).

Unlike the instantaneous AP “firing” of the LIF model, the more realistic AP initiation of the EIF model adds a suprathreshold component to the dynamic gain, which we can also extract. Beyond the local minimum, the phase plot’s slope dV/dt increases approximately exponentially with voltage (Fig. 2A). APs are now registered when the voltage reaches an AP detection threshold, with is much more depolarized than V_loc. Close to detection at V_detect, the dynamics is dominated by the intrinsic currents, whose magnitude far exceeds that of the extrinsic stimulus. As a consequence, the dynamics varies very little, leading to a very stereotypical AP shape in this voltage region. In contrast, early after passing V_loc, the dynamics is co-determined by intrinsic AP initiation current and the extrinsic stochastic stimulus of similar magnitude. This causes variable initiation delays preceding the rapid, stereotypical AP rise (see example in Fig. 1C₁). The addition of the intrinsic dynamics thus causes a variable AP onset, which drastically changes the high-frequency limit of G(f) from a constant value for LIF (Fig. 5F) to a power law decay for EIF (Fig. 2B; Brunel et al., 2001; Fourcaud-Trocmé et al., 2003). The addition of the intrinsic dynamics causes this variable AP onset, which drastically changes the high-frequency limit of G(f) from a constant value for LIF (Fig. 5F) to a power law decay for EIF (Fig. 2B; Brunel et al., 2001; Fourcaud-Trocmé et al., 2003).

The simple structure of the EIF model allowed Fourcaud-Trocmé and colleagues (Fourcaud-Trocmé et al., 2003) to relate the voltage dependence of the initiation current to the bandwidth of G(f). Here we propose to study this relation with a simple, phenomenological approach that is readily generalized to models of arbitrary complexity. Specifically, we will lower V_detect toward V_loc. At first, as long as we stay in a voltage region, where intrinsic currents dominate, the change of detection voltage will lead to small shifts in AP detection times that affect all APs to a very similar degree. However, when the lowered V_detect reaches the peri-threshold region, where the stochastic extrinsic stimulus dominates the dynamics, the different APs are shifted differently and we approach the behavior of a LIF-like model with a hard threshold. Figure 2B shows the “zero-delay” dynamic curves, resulting for V_detect = V_loc. Just as for LIF models, they are constant in the high-frequency limit. Therefore, the gain decay of the EIF model is closely associated with the suprathreshold dynamics after V_loc. In the general case, we obtain this gain decay as the ratio of two dynamic gain curves, G⁰(f)/G(f), where the “zero-delay” dynamic gain G⁰(f) is obtained with “zero-delay” APs detected already at V_loc, while G(f) is obtained with the conventional, much more depolarized AP detection threshold (Fig. 2B). In other words, we relate the encoding capability of a LIF-like model version to that of the full model and capture the result in the spectrally resolved gain decay.

The introduction of other ion channel types marks the next level of model complexity relating to impedance and spike gain, represented by conductance-based models, such as the Wang–Buszáki model (Xiao-Jing and György, 1996). These additional dynamical variables enable AP repolarization and richer neuronal dynamics, introducing different dynamical bifurcations at the AP threshold. This allows for a variety of firing patterns and consequently for different spike gain shapes. In spatially extended multi-compartment models (e.g., Eyal’s model Eyal et al., 2014), currents and voltage gradients between compartments increase the complexity even further. The presence of various ion channels together with the extended morphology results in a more complex current-to-voltage transformation that we describe with an effective impedance Z_eff(f), which also includes the stimulus filtering along the path from the somatic stimulus source to the axonal AP initiation site.

In summary, our decomposition approach establishes three G(f) components: (1) effective impedance Z_eff(f), (2) zero-delay spike gain Gsp0(f) and (3) dynamic gain decay G⁰(f)/G(f). The first component captures the subthreshold transformation of current fluctuations into voltage fluctuations. The third component captures the impact of variable AP initiation on the precise AP timing, reflecting the variability in the suprathreshold AP shape (Fig. 1C₁). The remaining zero-delay spike gain describes how voltage fluctuations are transformed into fluctuations of the firing rate, reflecting intrinsic frequency preferences of the neuronal dynamics, detailed in the following sections. Dynamic gain is the product of these transformations:G(f)=Zeff(f)⋅Gsp0(f)⋅G(f)G0(f). (12)For such a decomposition, we only require the waveforms of stimulus and membrane voltage (Materials and Methods). The decomposition can therefore be performed on neuron models of arbitrary complexity and even on recordings from real neurons (Fig. 1B₂). In the next sections, we will analyze the functional effects of biophysical and morphological parameters using experimental data and different variants of Eyal’s model, which feature a relatively wide encoding bandwidth (Eyal et al., 2014). But first, Eyal’s model will be characterized and studied at different working points.

Population encoding performance depends on the neuron’s excitability and working point

We first examined the dynamical properties of the multi-compartment model characterized by Eyal and colleagues (Eyal et al., 2014), specifically the neuron model consisting of only a soma and an axon but no dendrite. In response to various constant currents injected into the soma, the firing rate ν of the original model displays a large discontinuity of about 32 Hz upon reaching the rheobase, indicating that it is a type II model. For comparison, we devised a model variant that shares its morphology and is equipped with the same ion channel types at the same densities and with the same voltage sensitivities. However, in the new model variant, sodium channels activate and inactivate with a 10 mV shift toward more hyperpolarized potentials (Fig. 3C). This model variant has a continuous F-I curve, indicating it is a type I model (Fig. 3B, rheobase currents aligned at 0 pA).

We next examined the models’ response to stochastic input with a correlation time τ = 5 ms. Given that the firing rate has a strong impact on the population encoding bandwidth (Fourcaud-Trocmé et al., 2003), the firing rate is often fixed when different models or stimulus conditions are compared (Lindner and Schimansky-Geier, 2001). Here, probing the two-dimensional firing rate surface spanned by stimulus mean μ_I and standard deviation σ_I, we identified the iso-firing rate line at 5 Hz (Fig. 3D). Along this line, the neuron’s firing regime changes from nearly mean-driven to strongly fluctuation-driven. The curvature of the iso-5 Hz curve in the σ_I-μ_I plane differs between the two excitability types. It is slightly concave for type I, and slightly convex for type II excitability (Fig. 3E and H). The firing irregularity, quantified by the coefficient of variation of inter-spike intervals CV_ISI, differs more noticeably between the model variants. While the type I model’s firing irregularity increases monotonously with σ_I (Fig. 3I), the type II model displays a minimum CV_ISI for intermediate σ_I values. Toward more mean-driven conditions, CV_ISI increases strongly (Fig. 3F), because the intrinsic firing rate of 32 Hz is much higher than the 5 Hz target rate. In the presence of small input fluctuations around a near-rheobase average input, this low target rate is realized by bursts of near-regular rapid firing separated by long intervals. This is qualitatively very different from the irregular firing produced at working points further away from the mean-driven regime. For both model variants, we chose three working points, i.e., μ_I-σ_I combinations, covering the range from mean-driven to strongly fluctuation-driven (colored dots in Fig. 3E and H). The inter-spike interval (ISI) distributions at those working points (insets in Fig. 3F and I) further illustrate the very regular firing pattern of the type II model close to the mean-driven condition. As the input fluctuations increase, the ISI distributions become more similar.

When we calculated the dynamic gain functions for the two models at the respective working points, the results were again comparable for the strongly fluctuation-driven cases, but strikingly different for other working points. When the original type II model is mean-driven (σ_I = 5.13 pA), the dynamic gain function appears to have a low bandwidth, dropping to 70% magnitude already at 10 Hz (orange in Fig. 3G). At intermediate and high frequencies, the shape of the dynamic gain function is dominated by a strong resonance around 40 Hz, which mirrors the peak in the ISI distribution around 25 ms. Increasing σ_I to 21.13 pA significantly broadens the resonance peak and increases the cutoff frequency of the dynamic gain function, which now resembles the gain curve shown in Eyal et al. (2014) (blue line). The corresponding ISI distribution is also substantially flattened, indicating the transition from mean-driven to fluctuation-driven. Interestingly, a further increase to σ_I = 80.26 pA, does not enhance the population encoding ability further (green line). Instead, at this strongly fluctuation-driven working point, the dynamic gain function is dominated by a low-pass filter in the low-frequency region with an apparent cutoff frequency of around 12 Hz. Toward higher frequencies, a shoulder appears and results in a local maximum around 60 Hz. The dynamic gain in the high-frequency region exceeds the other two gain curves. In parallel, the ISI distribution narrows again and moves closer to zero. The most probable ISI now is approximately 70 Hz, corresponding to the shoulder in the dynamic gain function. These results demonstrate that the encoding ability critically depends on the working point; the high encoding bandwidth reported for this type II model, is realized only in a fluctuation-driven regime with intermediate σ_I.

The dynamic gain curves of the type I model behave strikingly different. None of the three working points lead to high-bandwidth encoding (Fig. 3J). When mean-driven (σ_I = 4.54 pA), the dynamic gain function has a cutoff frequency below 30 Hz, and decays with a slope of −1 in the log–log scale (orange dashed line), similar to the EIF model with standard AP initiation dynamics (Fourcaud-Trocmé et al., 2003). Increasing σ_I to 22.90 pA to reach the fluctuation-driven regime, the dynamic gain is enhanced in the high-frequency region, while the cutoff frequency becomes smaller (blue dashed line). Further increasing σ_I to 109.31 pA leads to a plateau in the low-frequency region with an even lower cutoff frequency (green dash line). The dynamic gain function is larger in the high frequency with a shoulder around 100 Hz, similar to the type II model’s gain curve at the very fluctuation-driven working point (compare two green curves in F and I). Again, the ISI distribution features a peak around the time interval that corresponds to the shoulder’s frequency range. These results demonstrate that the type I model cannot reproduce high-frequency encoding throughout the biophysically plausible range of working points. In summary, we found that the encoding bandwidth can be increased 10-fold, not by manipulating morphology or ion channel voltage sensitivity, but simply by increasing the potassium current around threshold to achieve a type II excitability. We next use the dynamic gain decomposition to investigate how this occurs.

AP initiation dynamics is not the main bandwidth limitation in Eyal’s model

To decompose the models’ dynamic gains curves according to Equation 12, we begin with the dynamic gain decay. This component elucidates the contribution of AP initiation dynamics on dynamic gain. To prepare this analysis, we first take a closer look at the AP initiation dynamics. We previously demonstrated that high-bandwidth encoding is closely associated with the voltage sensitivity of intrinsic AP initiation dynamics, especially in the voltage region close to the local minimum of the phase plot (Zhang et al., 2022). Using the voltage V_a, recorded at the AP initiation site, we compared the phase plots of type I and type II models with their local minima (V_loc) aligned at 0 mV (Fig. 4A). Although the voltage sensitivity of the sodium channel dynamics is identical for the two model variants, the voltage derivative of the type II model rises substantially more quickly out of V_loc. The slower initiation of the type I model is caused by its 16 mV more hyperpolarized V_loc. The sodium channels’ voltage dependence was shifted by only 10 mV, meaning that AP initiation proceeds at a voltage range with fewer activated sodium channels and of course even less potassium channels. As a consequence, the type I model’s AP initiation current displays a weaker voltage sensitivity.

When we gradually lowered the AP detection voltage V_detect down to V_loc, we obtained new, earlier AP times, but the number of APs stayed constant because we considered only those last positive crossings of V_detect that lead to fully developed APs. The resulting neuron model is similar to a LIF model with a hard threshold at V_loc, and indeed the redefined AP times lead to a zero-delay dynamic gain function that is almost flat in the high-frequency region as expected for a LIF model (continuous light blue curve in Fig. 4B). As we detect the APs earlier, we limit the voltage range, across which the initiation dynamics and extrinsic stimulus can impact the precise AP time. Comparing the dynamic gain functions for two detection thresholds, we reveal how the interplay between intrinsic and extrinsic currents in that interval (Vdetect1,Vdetect2 ) shape the dynamic gain. We find that lowering V_detect all the way from −3 mV down to −48 mV has very limited impact on the dynamic gain function. It only slightly enhances the encoding at large frequencies (above 500 Hz). However, a further decrease of V_detect by 5–10 mV drastically improves high-frequency encoding and has sizeable effects at frequencies down to 50 Hz, but not below. These results demonstrate that the AP initiation dynamics impacts the decay of the dynamic gain function mostly in the limit of high frequencies. It does not, however, dominate itsdecay between 40 and 200 Hz, the region where it falls from its peak to its cutoff frequency. In summary, for Eyal’s model the encoding bandwidth is governed by subthreshold factors, while dynamic gain at very high frequencies is determined by the AP initiation dynamics. Consequently, the AP initiation dynamics is not the main determinant of the encoding bandwidth.

• Download figure
• Open in new tab
• Download powerpoint

Figure 4.

Dynamic gain decomposition based on the decomposition of the AP initiation process. A, Insets show the AP phase plots in response to constant inputs just above rheobase. Arrows indicate the positions of the local minima V_loc of type II (−58.05 mV) and type I (−74.23 mV). The region at the beginning of AP initiation is magnified with voltages shifted to align the V_loc values to 0 mV. B, Reducing the AP detection voltage V_detect toward V_loc enhances the dynamic gain in the high-frequency region, exemplified here for the type II model. The zero-delay dynamic gain function (light blue continuous line) behaves similar to that of a LIF-like model. Further reducing the detection threshold reduces the dynamic gain in the low-frequency region (light blue dash line). Note that we only include threshold crossings that continue to a full AP. C, D, and E, The dynamic gain functions from Figure 3F and I are decomposed into effective impedance (C), zero-delay spike gain (D), and dynamic gain decay (E). Inset panels in E are the AP initiation delay distributions at corresponding working points, color code as in Figure 3. Effective impedance captures the transformation from somatic current to axonal voltage V_a. Zero-delay spike gain is the ratio between zero-delay dynamic gain and effective impedance. Dynamic gain decay is the ratio between original dynamic gain and zero-delay dynamic gain (see Materials and Methods). AP initiation delay is the time interval between voltage crossing V_detect (original AP time) and the last previous positive crossing of V_loc.

Further decreasing the detection threshold cuts into the subthreshold dynamics. The resulting changes in AP times can be substantial and affect the firing pattern and thereby also the dynamic gain in lower frequency regions (light blue dashed line in B). We can conclude that the type II model’s near constant encoding capability between 5 and 50 Hz, does not originate from its fast AP initiation after V_loc. Instead, it is determined by the subthreshold dynamics before V_loc.

Even the type I model’s slower AP initiation dynamics does not limit encoding between 5 and 50 Hz frequency for the fluctuation-driven regimes (blue and green). This is evident from the results of the gain decay analysis applied to the type I model (Extended Data Fig. 3-1 in Zhang et al., 2023). It is also reflected in Figure 4E, indicating that only in the nearly mean-driven regime (orange), does the bandwidth of the gain decay limit the overall dynamic gain bandwidth.

To quantify the voltage trajectories’ variability during AP initiation, we determined the time required for each AP, to progress from the last positive crossing of V_loc to V_detect = −3 mV (see Materials and Methods). This AP initiation delay is a random variable. Its statistics are determined by the interplay between intrinsic AP initiation dynamics and extrinsic stimulus fluctuations. The initiation delay’s variance originates mainly from very variable dynamics very close to V_loc, where the intrinsic currents are lowest and initiation proceeds slowest. For the type II model, 90 % of the total initiation delay are spent on crossing the first 10 mV after V_loc, even though they represent only 18% of the voltage interval. These 10 mV also contribute 99 % of the delay’s standard deviation. Figure 4B illustrates that these first 10 mV, from −58 to −48 mV, also account for the bulk of the dynamic gain decay associated with the AP initiation delay. We can use the dynamic gain decay G(f)/G⁰(f) (Fig. 4E), together with the distributions of AP initiation delay, to quantify the impact of AP initiation dynamics on high-frequency encoding at the three working points in Figure 3G and J. For the type II model, when the activity is mean-driven (orange), the AP initiation delay distribution is relatively flat with a most probable delay as large as 4 ms. Increasing σ_I reduces the mean of the distribution toward 0 ms (blue and green), indicating that the AP initiation delay is no longer limited by the intrinsic AP initiation dynamics. Instead, larger depolarizing inputs increase dV/dt, and accelerate the escape from the slow AP initiation region near V_loc. This faster initiation is paralleled by a rightward shift in the dynamic gain decay. The bandwidth related to the initiation dynamics increases more than three-fold. For the type I model, similar transitions of the dynamic gain decays and AP initiation delay distributions can be observed when increasing σ_I. While the two model types behave similarly when driven with very large σ_I, toward the mean-driven working points the type I model’s AP initiation delay distribution is substantially more flattened, and consequently, the dynamic gain decay has a much lower bandwidth. We attribute this to the different firing patterns close to the mean-driven condition. The type I model produces near-regular firing without the high-frequency bursts of its type II relative. Therefore, the average AP initiation delay would approach 200 ms when σ_I decreases to zero and the firing rate is kept at 5 Hz. For such initiation delays exceeding the input correlation time, the input fluctuations can cause particularly large variability and consequentially deteriorate encoding precision.

Low-frequency encoding is controlled by effective impedance and firing pattern preferences

The subthreshold part of the encoding process, which is described by the zero-delay dynamic gain G⁰(f), can be further decomposed into effective impedance Z_eff(f) and zero-delay spike gain Gsp0(f) (Fig. 4C and D). The effective impedance is calculated as the Fourier transform of output voltage divided by the Fourier transform of input current (see Materials and Methods). The zero-delay spike gain is the ratio of zero-delay dynamic gain and effective impedance: Gsp0(f)=G0(f)/Zeff(f) .

For the type II model, the effective impedance at low frequencies varies substantially with the working point, increasing more than 10-fold at 1 Hz (Fig. 4C, top panel). In a completely passive model without voltage-dependent conductance, the impedance does not depend on the input statistics. The observed changes in effective impedance, therefore, result from varying degrees of ion channel activity. At the mean-driven working point, the subthreshold voltage fluctuates just below V_loc, where some potassium channels are already activated. Depolarization-activated potassium current opposes further depolarization, lowering the impedance. This effect is strong for low frequencies, where the potassium channel activation can follow the input changes, which explains the large drop in the orange impedance curve. At higher frequencies, the potassium channels do not actively oppose depolarization. In the limit of frequencies much higher than their activation time constant, the channels present merely a passive leak and their influence on impedance is much smaller, such that the effective impedance decays with the slope of the passive model. The impedance peak around the frequency of the intrinsic firing rate probably results from subthreshold resonance.

Increasing σ_I and lowering μ_I shifts the average voltage to regions where fewer potassium channels are activated. This increases the effective impedance at low frequency (Fig. 4C, blue continuous line) and reduces the resonance around 40 Hz. Together, these changes contribute to the 100 Hz wide bandwidth of Z_eff at the intermediate, fluctuation-driven working point. At the very fluctuation-driven working point, the subthreshold fluctuations are even less shaped by potassium current. The effective impedance in this condition becomes a low-pass filter, determined by the passive neuronal properties (Fig. 4C, green continuous line). The type I model, in comparison, has less active ion conductances in the subthreshold range. Therefore, Z_eff(f) is far less sensitive to input conditions. At both fluctuation-driven working points, Z_eff(f) behaves similarly to a passive filter with a 10 Hz cutoff (Fig. 4C, blue and green dashed lines).

At first glance, the second subthreshold component, the zero-delay spike gain varies similarly for both models as the working points are changed. It attains larger values close to mean-driven conditions and drops with increasing σ_I (Fig. 4D). This overall behavior is dictated by the condition of constant firing rates. As μ_I decreases and σ_I increases, larger and larger voltage fluctuations generate the same firing rate, leading to a decrease in zero-delay spike gain. Another general feature of all zero-delay spike gain curves is their monotonic, power law increase in the high-frequency limit. It compensates the decay in the effective impedance to reproduce the rather flat zero-delay dynamic gain (Fig. 4D).

Besides these general trends, the shape of the zero-delay spike gain changes between working points, in particular in the low-frequency region. For the type II model, a high plateau forms toward the mean-driven condition (orange). Conversely, a sag forms at the extremely fluctuation-driven working point (green). We will see below, that the zero-delay spike gain is the one component, that is most sensitive to the firing pattern. It represents the neuron’s propensity to turn voltage fluctuations into firing rate fluctuations based on the intrinsic dynamics. It could therefore be considered surprising, that the zero-delay spike gains of these two working points show opposite shapes. After all, both working points are characterized by bursty firing patterns, evident from the high CV_ISI values and the ISI distributions (Fig. 3F). The difference lies in the source of the bursts. At the mean-driven working point, relatively small amplitude low-frequency components drive bursts fired at the intrinsic frequency. At the very fluctuation-driven working point, large fluctuations of the extrinsic stimulus cause bursts with an intra-burst firing rate that is almost twice as high. The increased stimulus fluctuation size, together with a lower number of APs within bursts causes the zero-delay spike gain to drop at low input frequencies. This will be explored further when we analyze how repetitively fired, intra-burst APs contribute to the dynamic gain.

The decomposition of G(f) curves across a large range of firing regimes has shown that under most conditions the initiation dynamics does not limit the overall encoding bandwidth, not even in the type I model with slower suprathreshold dynamics. The variable G(f) curves of the type II model and its ability for high-bandwidth encoding are explained by a strong stimulus dependence of Z_eff(f) and Gsp0(f) . Toward large fluctuations, both model types behave very similarly, because the extrinsic stimulus fluctuations dominate the entire signal transformation process, leading to similar Z_eff(f), Gsp0(f) and G(f)/G⁰(f). Decomposing the subthreshold and suprathreshold impact on dynamic gain, we can also aid the attribution of G(f) changes to a particular parameter change, as we next study the influence of input correlations and neuron morphology.

Fixing the suprathreshold impact reveals the influence of input correlations on dynamic gain

In the previous section, we studied the encoding abilities of Eyal’s type I and II models at working points ranging from mean-driven to strongly fluctuation-driven and found that only the type II model could support high-bandwidth encoding. Here we continue to examine this type II model, specifically how its encoding depends on stimulus correlation times (τ). When τ is increased, cortical neurons, but also LIF models, show enhanced G(f) at high frequencies, as compared to low frequencies (Brunel et al., 2001; Tchumatchenko et al., 2011; Lazarov et al., 2018; Merino et al., 2021). We have called this type of input dependence the “Brunel effect” and considered it an important feature of real neurons. For our study, we use the dynamic gain function in Figure 3G (blue curve, τ = 5 ms) as a reference, and compare it to results obtained with more slowly fluctuating input with τ = 50 ms.

When we set out to define the exact working points for this juxtaposition, we concluded that CV_ISI is of limited practicality, because its non-monotonic CV_ISI-μ_I relation in the type II model (Fig. 3F) does not allow a unique definition of working points. The dynamic gain decomposition motivates two other criteria for working point selection. While a constant CV_ISI aims to obtain a comparable output statistics, i.e., a comparable firing pattern, we propose to also study iso-fluctuation and iso-delay working points. They aim at a comparable subthreshold and suprathreshold dynamics, respectively. The iso-fluctuation working point is found by choosing σ_I and μ_I to fix σ_V, the standard deviation of subthreshold voltage fluctuations. This is closely related to criteria used in several experimental studies, where σ_V and the target firing rate have been used to define the input parameters (Köndgen et al., 2008; Tchumatchenko and Wolf, 2011; Ilin et al., 2013). The iso-delay criterion controls the suprathreshold impact by fixing the average AP initiation delay. Because the intrinsic AP initiation dynamics changes not much between different working points, fixing the average initiation delay for inputs with different τ values, effectively fixes the impact of the extrinsic currents on AP initiation delay. Figure 5, A and B, shows the AP initiation delay distributions, ISI distributions and, as an inset, the σ_V, of the working point with τ = 5 ms, and its iso-fluctuation and iso-delay counterparts with τ = 50 ms.

• Download figure
• Open in new tab
• Download powerpoint

Figure 5.

Brunel effect evaluated with two criteria fixing the sub- and suprathreshold impact on dynamic gain. A, AP initiation delay distributions of three working points of Eyal’s type II model. Two share the same mean delay, denoted as the iso-delay working points (the thin blue curve with τ = 5 ms and the wide black curve with τ = 50 ms). In contrast, the thin and wide blue curves share the same σ_V, denoted as the iso-fluctuation working points. This relation is indicated in the inset in B, showing the σ_V–σ_I relation for τ = 5 and 50 ms. The filled circles represent the three working points used from A to E. B–E, ISI distributions, dynamic gain and in lighter colors the corresponding zero-delay dynamic gain functions (C), effective impedances (D), and zero-delay spike gain functions (E) of the three working points above. Brunel effect, and the constituting sub- and suprathreshold contribution are different when evaluated at the iso-delay and iso-fluctuation working points. F, For a LIF model (Materials and Methods) the iso-delay condition cannot be defined by a delay between threshold crossing and reaching a AP detection threshold, because the two coincide. However, this panel demonstrates, that an equivalent condition can be defined by matching the rate of change of the membrane voltage at the threshold crossing. If μ_I and σ_I are chosen to keep the firing rate and the average dV/dt at threshold (see inset), then the high-frequency limit of the dynamic gain does not change. This is analogous to the behavior of the zero-delay dynamic gains of the iso-delay conditions in C.

Comparing the dynamic gain functions at the iso-fluctuation working points (Fig. 5C, thin and wide blue lines), we find that increasing τ merely introduces a resonance around 50 Hz, consistent with the ISI distribution peak around 20 milliseconds (Fig. 5B). In the low-frequency region, the two gain curves basically coincide, suggesting that the iso-fluctuation criterion indeed equalizes the subthreshold dynamics across the two input conditions. Because the dynamic gain curves in the high-frequency region are also nearly identical, one could conclude that the Brunel effect is absent under iso-fluctuation conditions. However, a comparison of the zero-delay dynamic gain functions (Fig. 5C, wide and thin light blue curves), does show improved high-frequency encoding, for τ = 50 ms as compared to τ = 5 ms. This means that the dynamic gain decay for τ = 50 ms is larger, which is consistent with the more variable AP initiation delay (Fig. 5A).

Next, we examined the Brunel effect at the iso-delay working points. With the same average AP initiation delay, the overall shapes of the AP initiation distributions are also similar (Fig. 5A thin blue line and gray line). Achieving similar initiation delays with much longer input correlations requires a larger σ_I. σ_V is even more than quadrupled from 6 to 26 mV (inset in Fig. 5B). But even the drastic increase in voltage fluctuations shifts the ISI distribution only slightly to smaller values (Fig. 5B wide lines). The large increase in σ_I for τ = 50 ms results in a strong reduction of the dynamic gain function in the low-frequency region (Fig. 5C). It now lies below the dynamic gain function for τ = 5 ms. The high-frequency limits of the two G(f) curves are close to each other, and even the zero-delay gain functions G⁰(f), now attain similar values in the high-frequency limit. We conclude that the type II model reproduces the Brunel effect at the iso-delay working points, which fix the impact of suprathreshold dynamics on population encoding. In this condition, the relative improvement of high-frequency encoding observed in the normalized dynamic gain functions is realized by reducing the unnormalized dynamic gain in the low-frequency region when τ is larger. At iso-fluctuation working points, on the other hand, the low-frequency encoding is nearly constant and no Brunel effect manifests in the dynamic gain curves because the stronger dynamic gain decay cancels increased zero-delay dynamic gains at τ = 50 ms.

We can now use the dynamic gain decomposition to check whether iso-fluctuation and iso-delay criteria indeed standardize sub- and suprathreshold contributions, and to further understand the source of the Brunel effect (Fig. 5D and E). We find that the effective impedances and zero-delay spike gains are very similar at the iso-fluctuation working points, in accordance with its purpose of fixing subthreshold dynamics. Within the identical subthreshold voltage fluctuation range, the active ion conductances recruited for stimulus filtering are close to each other. The remaining difference probably stems from the slightly more depolarized average voltage for τ = 50 ms, leading to slightly higher potassium current activation. Two zero-delay spike gain curves are almost paralleled with each other, indicating that the voltage frequency components are transformed into firing frequency components in similar ways for both correlation times. However, at the iso-delay working points, the effective impedance and spike gain in the low-frequency region show large but opposite changes in response to increasing τ. While the effective impedance at 1 Hz is more than 3 times larger for τ = 50 ms, the spike gain at 1 Hz is more than 5 times larger for τ = 5 ms. The difference in effective impedance is readily explained by the difference in average potassium channel conductance. σ_V is about 6 mV for τ = 5 ms, and increases to about 26 mV for τ = 50 ms. The substantially larger fluctuations into the hyperpolarized range cause ion channel deactivation, a lower membrane conductance and thereby a larger effective impedance. This is reflected in the much lower cutoff frequency of 8.3 Hz for τ = 50 ms as compared to the 56 Hz for τ = 5 ms. Low-pass effective impedance keeps more low-frequency components of stimulus in the voltage fluctuations. As a result, the zero-delay spike gain at the iso-delay working point is drastically reduced in the low-frequency region (Fig. 5E, black curve). The decomposition of zero-delay dynamic gain functions implies two different ways to realize high-bandwidth encoding for two correlation times of input, by influencing the firing patterns emitted by the neuron. How that impacts the dynamic gain will become clear, when we next study the differential encoding capacity of fluctuation-driven APs and repetitively firing APs.

Decomposing dynamic gain functions into repetitive firing and individual AP components

Our previous arguments suggest that the firing pattern, in itself, is an important contributor to the shape of the dynamic gain function. It could even be thought that the type II model is capable of high-bandwidth encoding in part due to its ability to fire high-frequency repetitive APs. We next attempt to analyze this more stringently, by decomposing the dynamic gain function into two parts, contributed by two AP populations (Fig. 6A).

• Download figure
• Open in new tab
• Download powerpoint

Figure 6.

Intrinsic high-frequency firing undermines high-bandwidth encoding. A, A short sequence of APs from the simulations in Figure 5 (τ = 50 ms, iso-delay), shown in a raster plot. The APs are classified into two groups according to the preceding ISI (ISI ≥50 ms long, black dashes, ISI <50 ms shorter, gray dashes). Separate dynamic gain functions of the iso-delay and iso-fluctuation dynamic gain functions in Figure 5C are calculated for both spike classes and shown here. B, The dynamic gain curves derived from fluctuation-driven APs. For the slower fluctuations, high-frequency input components are encoded better than low-frequency components as a consequence of resonant firing. C, The dynamic gain curves derived from APs fired in close succession, e.g., repetitively fired APs. Resonance effects cause the dips and peaks (see main text). D and E, The AP initiation delay distributions of the two spike classes are very similar in the case of iso-delay, and hence very similar to the joint distribution in Figure 5A. The delay distributions in the iso-fluctuation case show a much longer tail. Because the type I model produces much less bursting, the detrimental effects on encoding are weaker. The corresponding data are shown in Extended Data Figure 6-1 in Zhang et al. (2023).

One population is fluctuation-driven APs, defined as those fired after at least 50 ms of silence. These APs are likely to be the first one, or the only one, fired in a voltage upstroke, while the rest of APs are likely part of a burst later in an upstroke, denoted as repetitive firing APs. Figure 6A is an illustration of the AP classification, the relative fraction of the two components is stated in Figure 6B and C between correlation times. When τ = 50 ms, the ISI distribution is significantly more peaked around 25 ms (Fig. 5B), and more than half of the APs are classified as repetitively firing APs. This fraction drops to 20% for τ = 5 ms, because bursts are terminated by the rapid V_a fluctuations.

We now calculated dynamic gain functions based on APs from only fluctuation-driven or only repetitive firing APs. When their phases are taken into account, these two dynamic gains can be summed to yield the curves in Figure 5C (black and thin blue). The fluctuation-driven APs alone encode the input with a higher bandwidth, as compared to all APs together (thin blue line in Fig. 6B vs iso-delay working points in Fig. 5C). The corresponding cutoff frequencies are 133 Hz vs 105 Hz for τ = 5 ms, 234 Hz vs 166 Hz for τ = 50 ms. The shapes of these gain curves reflect the broad 20–50 Hz frequency encoding preference of the fluctuation-driven APs. When the intrinsic preference is overridden by very strong external fluctuations, this preference is suppressed (Extended Data Fig. 6-1A and 1C in Zhang et al., 2023).

The repetitive firing APs alone lead to dynamic gain curves with very pronounced resonance around the preferred firing rate and a pronounced dip around half of that frequency (Fig. 6C). The dip is a direct consequence of the resonance because APs locked to an input component at the resonance frequency and fired in a short burst will be locked worse than random to half of that frequency. The dynamic gain functions derived from repetitive APs alone have a comparatively narrow bandwidth, especially for τ = 5 ms, where the cutoff frequency is below 8 Hz. In conclusion, the two dynamic gain functions for τ = 5 and 50 ms at the iso-delay working points realize high-bandwidth encoding in different ways. One reduces the fraction of repetitive firing APs, while the other forms longer bursts of repetitive firing APs to weaken its low-pass filter effect. In either case, our results clearly show that the high encoding bandwidth of Eyal’s model with type II excitability is not a direct result of the burst firing. APs fired within bursts only contribute the pronounced peaks at the preferred firing rate. We can also use this AP classification strategy to further investigate the variation of spike gain shape across working points. Low frequencies of only a few Hertz were particularly well represented in the zero-delay dynamic gain curve, when the neuron was nearly mean-driven (Fig. 4D, orange). This is caused by a relatively large fraction of repetitive firing APs (Extended Data Fig. 6-1B in Zhang et al., 2023). They occur in long bursts and thereby limit the encoding of intermediate frequencies. This can be understood when considering the example of a 100 ms long burst, for which first and last APs fall on diametrical phases of a 5 Hz stimulus component and thereby limit encoding of this frequency. In contrast, the fluctuation-driven APs, do not show such a preference for very low frequencies, but they make up only a small fraction of all APs at the nearly mean-driven working point.

Dendrites enhance high-bandwidth encoding by suppressing low-frequency effective impedance

As an application of the dynamic gain decomposition, we disentangled the complicated effects of dendritic morphology on each stage of information processing. This analysis elucidates the impact of dendrites on population encoding. We took the two model variants with a median and a large dendrite from Eyal et al. (2014) for comparison (see Materials and Methods), both display type II excitability (see inset Fig. 7A). We first examined the impact of dendrite size on intrinsic AP initiation dynamics. Increasing the dendrite size increases the lateral current from the initiation site toward the soma, due to the larger somato-dendritic current sink. V_loc is shifted to slightly more depolarized values from −58.05 to −54.72 and −53.61 mV. As a result, the AP onset is shifted toward voltages at which the sodium channels have steeper voltage dependence. Consequently, at later stages, e.g., 10 mV/ms, the slope of the phase plot is actually larger for models with larger dendrites, as has been previously reported in Eyal et al. (2014) (Fig. 7A). However, when aligning V_loc to (0 mV, 0 mV/ms), we found that a larger dendrite does not accelerate AP initiation close to V_loc, instead, the local slope of the phase plot is reduced. A first evaluation of the subthreshold effects of the dendrites can be gained from the σ_V–σ_I relation (Fig. 7B, inset panel). The large current sink of the dendrite drastically decreases the impedance of the neuron, reflected by the reduced slopes. Compared to this drastic impedance effect, the dendrite’s impact on the suprathreshold dynamics appears to be subtle. Therefore, to compare the population encoding capabilities of the three models, we chose iso-delay working points to harmonize the suprathreshold contribution and focus on different subthreshold contributions of the dendrite.

• Download figure
• Open in new tab
• Download powerpoint

Figure 7.

A larger dendrite shunts low-frequency inputs, but improves low-frequency spike gain. A, Phase plots of three model variants with different dendrite sizes, V_loc aligned at (0 mV, 0 mV/ms). The inset shows the F–I curves of Eyal’s type II models with different dendrite sizes. Larger dendrites increase the minimal firing rate. B, AP initiation delay distributions at the iso-delay working points (ν = 5 Hz). Inset panel shows corresponding σ_V on σ_V–σ_I relations. C, ISI distributions at the iso-delay working points. D, E, and F, Dynamic gain functions and in lighter colors the corresponding zero-delay dynamic gain functions (D), effective impedances (E) and zero-delay spike gain (F) of the three model variants calculated at iso-delay working points (upper panels), together with their normalized gain curves (lower panels). See the Extended Data Figure 7-1 in Zhang et al. (2023) for the results with τ = 50 ms.

Taking the dynamic gain function of τ = 5 ms in Figure 5 for reference, we fixed the average AP initiation delay for the three model variants, leading to very similar distribution shapes (Fig. 7B). Under this criterion, σ_V is smaller for the neuron model with a dendrite (scattered dots in inset panel), and the three ISI distributions are similar to each other (Fig. 7C). Compared to the dynamic gain function without the dendrite (blue), the other two gain curves (red and green) have smaller dynamic gain below 200 Hz. In the high-frequency region, gain curves decay with similar trends, and corresponding zero-delay dynamic gain functions are also close to each other (Fig. 7D, top panel). Normalizing the dynamic gain functions, we observed an enhancement of high-frequency encoding with a larger dendrite (lower panel in D).

Decomposing the zero-delay dynamic gain functions into effective impedances and zero-delay spike gains, we found that adding a dendrite reduces the effective impedance mostly for low frequencies, but not in the high-frequency limit. (Fig. 7E, top panel). Low-frequency components are suppressed more strongly, because they charge a larger portion of the dendrite, while high frequencies experience much stronger spatial filtering and charge only the proximal dendrite. When the effective impedances are normalized (lower panel in E), the dendrites’ suppression of low frequencies appears as a boost of high-frequency representation in the voltage. Interestingly, a further substantial increase in the dendrite size hardly affects the shape of the effective impedance. The zero-delay spike gain curves are relatively flat in the low-frequency region, and they increase in similar trends in the high-frequency region (Fig. 7F, top panel). The curves for the two dendrite-bearing models are almost identical (red and green), while the zero-delay gain curve of the original model has lower values (blue) because it produces very similar firing patterns (Fig. 7C) with approximately three times larger subthreshold voltage fluctuations (inset in Fig. 7B). The normalized curves in the lower panel show that, the dynamic gain enhancement caused by the effective impedance is slightly undermined at the stage of zero-delay spike gain (from blue to red, and green). Taken together, the improved high-frequency encoding in the presence of a dendrite is primarily due to a suppression of low frequencies in the effective impedance, with a minor effect of the accelerated AP initiation, leading to weaker dynamic gain decay. Studying the effect of the dendrite under slower input correlations, led to the same conclusion. In particular, the dendrite did not change the Brunel effect (Extended Data Fig. 7-1F in Zhang et al., 2023).

Dynamic gain decomposition provides subcellular dissection of a pathophysiological insult on population coding

The dynamic gain decomposition is designed to be applicable to any neuron, independent of its complexity. It should be applicable even to recordings from real neurons, provided that the recorded somatic voltage contains sufficient information. For the subthreshold analysis, this is very likely, but the utility of somatic recordings for initiation-related analysis is not clear a priori. We therefore set out to apply the dynamic gain decomposition to experimental data and chose a specific data set from a previous study for two reasons. First, the two groups in the data set had been reported to have different input resistances, prompting us to expect different effective impedances. Second, the treatment group, but not the control group, contained neurons with very different dynamic gain curves, which we would like to understand better by decomposition.

The data in question were recorded from layer 5 pyramidal neurons in coronal slices of mouse somatosensory cortex and originally published in Revah et al. (2019). Only slices in the treatment group underwent two brief hypoxic episodes that induced spreading depolarization (SD). Unlike simulations, experiments record voltage at the soma and not the initiation site. We showed earlier, that the transfer impedance from soma to initiation site does pose an additional filter for higher input frequencies. However, for distances up to 50 μm, this effect was very small (Zhang et al., 2022), and hence we speculate that using the somatic voltage will still allow decomposition.

We decomposed the dynamic gain functions of each pyramidal neuron into effective impedance and spike gain. For two, very differently affected neurons from the treatment group of Revah et al. (2019), the results are shown in Figure 8. Neuron 1 appears very similar to control neurons, it displays encoding with a bandwidth of approximately 350 Hz (Fig. 8A). Its dynamic gain’s confidence interval widens substantially, as the gain curve drops below the noise floor. In comparison, the dynamic gain of neuron 2 (Fig. 8B) is substantially larger in the low-frequency region, decays more steeply at intermediate frequencies, and displays a lower cutoff frequency of around 200 Hz. But surprisingly, above 400 Hz, its dynamic gain reemerges above the noise floor and increases with a smaller bootstrapping confidence interval as compared to neuron 1. Also note that the absolute magnitude of the dynamic gain, for instance at 200 Hz, is higher in neuron 2, although its cutoff frequency is lower. The underlying properties that cause these differences became clear, when we decomposed the dynamic gain into the sub- and suprathreshold components.

• Download figure
• Open in new tab
• Download powerpoint

Figure 8.

Ischemic insult slows AP initiation. A, Dynamic gain function G(f) (upper panel) of one pyramidal neuron from the treatment group in Revah et al. (2019), and the associated gain curves recalculated with lower AP detection thresholds ranging from 0 mV to −50 mV (lower panel). Inset panels are the phase plots estimated by averaging voltage derivatives at each voltage value. Note the wide bandwidth of the original gain function. Gray areas are the 95th percentile bootstrapping confidence intervals. Thin line represents the 95th percentile noise floor (see Materials and Methods). B, as A, but for another neuron from the treatment group with drastically different gain. Firing rates are 4.10 and 4.68 Hz in A and B respectively. C and D, Effective impedances and zero-delay spike gains of the dynamic gain functions in A and B. Gsp0 is the ratio between zero-delay dynamic gain (V_detect = −50 mV in C and D) and effective impedance. E, Grand average of G(f), Z_eff, and G_sp(f) for control neurons (CT, n = 9) and neurons undergoing brief hypoxia and SD (n = 15). Lines and shades represent means and their standard error. G(f) and G_sp(f) are displayed for the frequency region in which all individual traces are significant, i.e., above the respective noise floor. F, G(f) at 200 Hz does not differ between groups. Due to the large Z_eff differences between the groups, the spike gains at 200 Hz are significantly different (t-test, *p = 0.048).

As a first step, we estimated the phase plot from the recordings (see Materials and Methods) by averaging voltage derivatives at different voltages (inset panels in A and B). This reveals a difference in threshold voltage, but most importantly, a large difference in the amount of voltage-sensitive depolarizing currents between the two neurons. For neuron 1, the voltage derivative changed little for voltages below −40 mV, then suddenly increased to 30 mV/ms within 8 mV, reflecting large AP initiation currents. Neuron 2, in contrast, had smaller depolarizing currents and a higher threshold. The position of the local minimum V_loc was around −30 mV for neuron 2 and, although it could not be discerned precisely, it was likely just below −40 mV for neuron 1.

In the second step, we reduced the AP initiation delay by lowering the AP detection threshold V_detect analogous to the analysis of Eyal’s model in Figure 4B. The recalculated dynamic gain functions for the various V_detect values are given in the lower panels of Figure 8A and B. For neuron 1, the gain curves changed abruptly, once V_detect reached −40 mV and already at −45 mV a saturation was reached. Apparently, only this narrow suprathreshold voltage range contributed to the uncertainty in AP timing. This corresponds to the steep change in the phase plot (inset). The initiation-related decay of the dynamic gain affected only frequency components above 300–350 Hz. When the AP was detected around V_loc, the dynamic gain curve resembled that of a LIF neuron model without voltage-dependent initiation and with a hard threshold Brunel et al. (2001). This result replicates the simulations shown in Figure 4B, albeit with an almost four times higher bandwidth of the dynamic gain decay, as compared to Eyal’s type II model.

In neuron 2, the recalculated dynamic gain curves behaved differently. The transition to LIF-like behavior was much more gradual (Fig. 8B, colored traces). The uncertainty in AP timing accumulated during the transition through a much larger voltage range from above −25 mV down to −40 mV, mirroring the more gradual change in phase plot slope (inset). Another important difference is the affected frequency range. In neuron 2, suprathreshold delay variability already reduced the encoding for frequencies above 90–100 Hz, similar to what we observed for Eyal’s model.

Our insights into the different suprathreshold dynamics are complemented by an analysis of the subthreshold contribution to dynamic gain. We calculated the effective impedance as the ratio of Fourier transform of voltage and Fourier transform of current (see Materials and Methods). The two neurons’ effective impedances and zero-delay spike gains are juxtaposed in Figure 8C and D. This reveals the considerably larger effective impedance of neuron 2 as the origin of this neuron’s large absolute dynamic gain value. The fact that this effective impedance rises for frequencies above 30 Hz is puzzling and we can only speculate about a contribution of sodium channels, possibly due to the increased threshold in neuron 2. For neuron 1, the effective impedance approximates the shape one would expect from a passive soma with dendrites. It is dominated by a cutoff of 12.8 Hz, i.e., a membrane time constant of 12.5 ms. Beyond the cutoff, the slope is less negative than −1, most likely due to the presence of dendrites. The zero-delay spike gains of the two neurons are rather similar at low frequencies, but above 20 Hz, they deviate increasingly. Because these curves are obtained by dividing zero-delay dynamic gain curves by the effective impedance, this deviation originates in the unusual impedance increase in neuron 2.

In summary, the differences between the two neurons’ dynamic gain curves can be attributed to subthreshold and suprathreshold contributions as follows: the difference in absolute values is caused by the higher effective impedance of neuron 2. The different locations of the dynamic gain drop, and consequently the different bandwidths are explained by the substantially slower AP initiation in neuron 2. It causes an earlier drop in dynamic gain due to larger initiation delay variability. The remaining difference is the striking rise in the high-frequency region of the second neuron’s dynamic gain. We suspect that at the high frequencies the extrinsic, stimulation-induced currents dominate the second neuron’s weaker intrinsic currents. Consequently, the exact time of threshold crossing is influenced by high-frequency input fluctuations, similar to the situation in a LIF neuron model. This is the reason for the rising dynamic gain and the relatively narrow confidence interval. In the context of the original study (Revah et al., 2019), it is interesting to note, that cells with such a weak intrinsic AP initiation appeared only after hypoxia and SD, which also compromised the molecular integrity of the AIS.

Applying dynamic gain decomposition to all the cells in the original data set, we could further investigate the previously described effect of hypoxia-induced SD on the dynamic gain. The grand averages for dynamic gain 〈G(f)〉, effective impedance 〈Z_eff(f)〉 and spike gain 〈G_sp(f)〉 were calculated over the frequency range, where data were available from all cells. For 〈G(f)〉 and 〈G_sp(f)〉 that only includes the range up to 300 Hz, at which point at least one cell’s G(f) touched the noise floor. The grand averages, together with their standard errors are shown in Figure 8E. They show clearly that the average effective impedance is increased in the treated cells, and that this difference dominates the unnormalized dynamic gain curves. While the dynamic gain of the treated cells appears to drop earlier, between 100 and 200 Hz as compared to the control cells >300 Hz, this difference does not manifest in a statistically significant difference. The dynamic gain values at 200 Hz are spread out over similar ranges for both groups (Fig. 8F). In contrast, the impedance-corrected spike gain values are significantly different (Student’s t-test, p = 0.048). These examples show that all aspects of dynamic gain decomposition can be applied to experimental data and serve to suggest biophysical parameters as the basis of individual dynamic gain features.