Accurate and fast simulation of channel noise in conductance-based model neurons by diffusion approximation

Abstract

Stochastic channel gating is the major source of intrinsic neuronal noise whose functional consequences at the microcircuit- and network-levels have been only partly explored. A systematic study of this channel noise in large ensembles of biophysically detailed model neurons calls for the availability of fast numerical methods. In fact, exact techniques employ the microscopic simulation of the random opening and closing of individual ion channels, usually based on Markov models, whose computational loads are prohibitive for next generation massive computer models of the brain. In this work, we operatively define a procedure for translating any Markov model describing voltage- or ligand-gated membrane ion-conductances into an effective stochastic version, whose computer simulation is efficient, without compromising accuracy. Our approximation is based on an improved Langevin-like approach, which employs stochastic differential equations and no Montecarlo methods. As opposed to an earlier proposal recently debated in the literature, our approximation reproduces accurately the statistical properties of the exact microscopic simulations, under a variety of conditions, from spontaneous to evoked response features. In addition, our method is not restricted to the Hodgkin-Huxley sodium and potassium currents and is general for a variety of voltage- and ligand-gated ion currents. As a by-product, the analysis of the properties emerging in exact Markov schemes by standard probability calculus enables us for the first time to analytically identify the sources of inaccuracy of the previous proposal, while providing solid ground for its modification and improvement we present here.

Figure 1. Markov kinetic schemes.

In the simplest 2-state kinetics (A), a single channel can be in one of two configurations with only one of them associated to a non-zero conductance (filled grey circle). The kinetic parameters and are rates, as they represent the transition probabilities between states, expressed per time unit. In a more general case, single-channel kinetics is described by an -state scheme. Voltage-gated fast-inactivating sodium (B) and delayed-rectifier potassium channels (C) are two examples, where only one state corresponds to a non-zero channel conductance (filled grey circle). An alternative model for sodium channels (D) (Vandenberg and Bezanilla, 1991) is also shown for comparison. We point out that our method can be applied to any kind of kinetic schemes, where the transition rates are known. For (B–C), each state is identified by an arbitrary name convention (, , , etc.), referring to the underlying mapping of these 8- and 5-state channels into multiple 2-state gated subunits (panel E). Indeed, some -state kinetic schemes may be mapped into, or experimentally identified as, a set of independent 2-state gates: the open state of the full scheme corresponds to all the elementary gates in the open states, simultaneously. For instance, the kinetic scheme (B) could be mapped into a set of four independent 2-state gates (E) (i.e., the familiar activation gates and the inactivation gate of sodium fast-inactivating currents), three of whom are identical.

Figure 2. Steady-state statistical properties of the fraction of open channels , under voltage-clamp.

Panels A–C refer to delayed-rectifier potassium channels (see Fig. 1B and Table 2), whereas panels D–F refer to fast-inactivating sodium channels (see Fig. 1A and Table 2). Black and red dots result from the simulations of the exact kinetic schemes and from our diffusion approximation, respectively. The continuous traces in A,B,D,E are drawn by the analytical expressions derived in the text, and refer to an increasing number of simulated channels (namely, 360, 1800, 3600). The dependence on the membrane-patch voltage is studied for the mean of (A,D) and for its variance (B,E). For an increasing number of channels, the variance decreases, as expected. Panels C,F show the time constant of the best-fit single-exponential, which approximates the covariance of (see Eq. 17). The mismatch between actual best-fit values and the characteristic subunit gating time-constants (, , , shown for comparison), clearly indicates that great care should be taken in deriving accurate Langevin-kind formulations. Panels G–L repeat the very same comparisons presented in panels A–F, for the Langevin-approximation introduced by Fox and coworkers (Fox, 1997; Fox and Lu, 1994): the variance of potassium currents is overestimated (H), whereas the variance of sodium currents is underestimated (K). In addition, the autocorrelation properties are not reproduced correctly (I,L).

Figure 3. Sample time-series of the fraction of open channels , under voltage-clamp ().

Panels A–C refer to delayed-rectifier potassium channels (see Fig. 1B and Table 2), and panels D–F to fast-inactivating sodium channels (see Fig. 1A and Table 2). Black and red traces and dots result from the simulations of the exact kinetic schemes and from our diffusion approximation, respectively. The continuous traces in A,D are steady-state realisations of the fraction of open potassium and open sodium channels, respectively, while panels B,E display the amplitude histogram. Under the conditions considered here (360 potassium and 1200 sodium channels), the Gauss-distributed effective stochastic process approximates well the microscopic model. Panels C,F report the autocorrelation function of channel noise fluctuations, demonstrating an excellent agreement of the effective and microscopic simulations (see also Fig. 2C,F). Panels G–L repeat the same comparisons presented in panels A–F, for the Langevin-approximation introduced by Fox and coworkers (Fox, 1997; Fox and Lu, 1994). As in Fig. 2H,K the variance of potassium currents is overestimated (G–H) while the variance of sodium currents is underestimated (J–K). In addition, the autocorrelation properties are not reproduced correctly (I,L). Additional simulations, for distinct values of the holding membrane potential, are provided as Supporting Information (Figures 5–10 in Text S1).

Figure 4. Spontaneous firing in the microscopic and effective models.

When weakly depolarising DC currents (A, ) are applied to both the microscopic (black sample trace) and the effective models (red sample trace), the increase in channel noise variances (see Fig. 2C,F) induces a highly irregular spontaneous emission of action potentials, with qualitatively very similar properties. In these simulations, both length and diameter of the neuron are set to , and the single channel conductance for both sodium and potassium channels is . Panels B,C show respectively the CV of the ISI distribution and the mean firing rate as a function of cell diameter: results are reported for the microscopic, effective and Fox's models (black, red and blue traces, respectively). The results of panels B,C refer to spontaneous activity (i.e., no injected current) with neuron length held fixed at the value .

Figure 5. Comparison of firing efficacy, latency and jitter of a sharp current pulse.

Panels A, B and C display the firing efficacy, the average latency and the jitter of the evoked responses, respectively, after the application of a monophasic stimulus of duration repeated for 10000 trials. Black and red traces and dots result from the simulations of the exact kinetic schemes and from our diffusion approximation, and in blue we indicate the results from the simulation of the Langevin-approximation introduced by Fox. Panel D shows the distribution of spike occurrence times, evoked by a biphasic stimulus over 10000 trials: the duration and amplitude of the preconditioning part are and , respectively, the duration and amplitude of the second part are and . In all panels, the neuron is simulated as a single cylindrical compartment of length and diameter equal to and single channel conductances equal to , for both sodium and potassium channels. The integration time step was set to .

Figure 6. Raster plots and peristimulus time histograms (PSTH) for the timing of spiking responses to repeated identical DC pulses (A) and fluctuating currents (B).

Red traces and markers refer to Montecarlo microscopic simulations of the full model, while black traces and markers refer to the effective model. The values of reliability () and precision () are in accordance with those measured in in vitro experiments. In particular, in panel A: , for the microscopic model, , for the effective model. Panel B: , for the microscopic model, , for the effective model. The DC pulse has an amplitude of , whereas the noisy stimulus is the same realisation of an Ornstein-Uhlenbeck's process, with mean and standard deviation set to , and with autocorrelation time length set to .

Figure 7. Frequency-current () response curves.

Mean firing rate, in response to a DC current injection, studied for increasing stimulus intensities in both Montecarlo microscopic (black trace) and effective model (red trace) simulations. Single-channel conductance for both sodium and potassium channels set to .

Figure 8. Voltage power spectral densities of subthreshold membrane potential trajectories.

Comparison between the microscopic (thick shaded lines) and the effective (thin solid lines) models. of simulated recordings of the membrane potential were obtained under weak holding currents ({), resulting in membrane potential traces fluctuating around an offset ({). Rare spontaneous spikes were removed from the analysis, excluding the preceding and the following each spike. The spectra have been obtained by applying the Welch method, on moving windows of duration and overlapping by , and subsequently averaging the results.