13 Continuity Equation and the Fokker-Planck Approach

13.6 Neuronal adaptation and synaptic conductance

In the previous section we analyzed a population of exponential integrate-and-fire neurons driven by diffusive noise with mean μ(t)\mu(t) and variance σ(t)\sigma(t). At a first glance, the results might seem of limited relevance and several concerns may be raised.

(i) What happens if we replace the one-dimensional exponential integrate-and-fire model by a neuron model with adaptation?

(ii) What happens if neurons are embedded into a network with coupling within and between several populations?

(iii) What happens if neurons don’t receive synaptic current pulses that lead to a jump of the membrane potential but rather, more realistically, conductance-based input?

(iv) What happens if the input noise is not ‘white’ but colored?

All of these questions can be answered and will be answered in this section. In fact, all questions relate to the same, bigger picture: suppose we would like to simulate a large network consisting of KK populations. In each population, neurons are described by a multi-dimensional integrate-and-fire model, similar to the adaptive integrate-and-fire models of Chapter 6.

The set of equations (6.1) and (6.2) that controls the dynamics of a single neuron ii is repeated here for convenience

τmduidt\displaystyle\tau_{m}\,\frac{{\text{d}}u_{i}}{{\text{d}}t} =\displaystyle= f(ui)-Rkwk,i+RIi(t)\displaystyle f(u_{i})-R\,\sum_{k}w_{k,i}+R\,I_{i}(t) (13.56)
τkdwk,idt\displaystyle\tau_{k}\,\frac{{\text{d}}w_{k,i}}{{\text{d}}t} =\displaystyle= ak(u-urest)-wk,i+bkτkt(f)δ(t-ti(f)),\displaystyle a_{k}\,(u-u_{\rm rest})-w_{k,i}+{b_{k}\tau_{k}}\sum_{t^{(f)}}% \delta(t-t_{i}^{(f)})\,, (13.57)

where f(u)=-(u-urest)+ΔTexp(u-ϑrhΔT)f(u)=-(u-u_{\rm rest})+{\Delta_{T}}\,\exp\left({u-\vartheta_{rh}\over\Delta_{T% }}\right) is the nonlinear spike generation mechanism of the exponential integrate-and-fire model. The voltage equation (13.56) is complemented by a set of adaptation variable wk,iw_{k,i} which are coupled to the voltage and the spike firing via Eq. (13.57). The first of the above questions concerns the treatment of these slow adaptation variables in the Fokker-Planck framework.

In a network with a more realistic synapse model, the input current IiI_{i} of neuron ii is generated as a conductance change, caused by the spike firings of other neurons jj

Ii(t)=jfgij(t-tj(f))(ui-Eijsyn),I_{i}(t)=\sum_{j}\sum_{f}g_{ij}(t-t_{j}^{(f)})\,(u_{i}-E_{ij}^{\rm syn})\,, (13.58)

where gij(s)g_{ij}(s) for s>0s>0 describes the time course of the conductance change and EijsynE_{ij}^{\rm syn} is the reversal potential of the synapse from jj to ii. Question (ii) concerns the fact that the spike firings are generated by the network dynamics, rather than by a Poisson process. Questions (iii) and (iv) focus on the aspects of conductance (rather than current) input and temporal extension of the conductance pulses.

Let us now discuss each of the four points in turn.

13.6.1 Adaptation currents

In order to keep the treatment simple, we consider the Adaptive Exponential Integrate-and-fire Fire model (AdEx) with a single adaptation current (cf. Chapter 6, Eq. (6.3)). We drop the neuron index ii, and consider, just as in the previous sections, that the stochastic spike arrival can be modeled by a mean μ(t)=RI(t)\mu(t)=R\,\langle I(t)\rangle plus a white-noise term ξ\xi

τmdudt\displaystyle\tau_{m}\,\frac{{\text{d}}u}{{\text{d}}t} =\displaystyle= -(u-urest)+ΔTexp(u-ϑrhΔT)-Rw+μ(t)+ξ(t)\displaystyle-(u-u_{\rm rest})+{\Delta_{T}}\,\exp\left({u-\vartheta_{rh}\over% \Delta_{T}}\right)-R\,w+\mu(t)+\xi(t) (13.59)
τwdwdt\displaystyle\tau_{w}\,\frac{{\text{d}}w}{{\text{d}}t} =\displaystyle= a(u-urest)-w+bτwt(f)δ(t-t(f)).\displaystyle a\,(u-u_{\rm rest})-w+{b\tau_{w}}\sum_{t^{(f)}}\delta(t-t^{(f)}). (13.60)

The stochastic input drives the neuron into a regime where the voltage fluctuates and the neuron occasionally emits a spike. Let us now assume that the input is stationary, i.e., mean and variance of the input are constant. Suppose furthermore that the time constant τw\tau_{w} of the adaptation variable ww is larger than the membrane time constant and that the increase of ww during spike firing is small: b1b\ll 1. In this case, the fluctuations of the adaptation variable ww around its mean value w\langle w\rangle are small. Therefore, for the solution of the membrane potential density equations p0(u)p_{0}(u) the adaptation variable can be approximated by a constant w0=ww_{0}=\langle w\rangle (431; 186). This separation-of-time-scales approach can also be extended to calculate the steady-state rate and time-dependent response for neurons with biophysically detailed voltage-gated currents (435).

13.6.2 Embedding in a network

Embedding the model neurons into a network consisting of several populations proceeds along the same line of arguments as in Section 13.4 or Chapter 12.

The mean input μ(t)\mu(t) to neuron ii arriving at time tt from population kk is proportional to its activity Ak(t)A_{k}(t). Similarly, the contribution of population kk to the variance σ2\sigma^{2} of the input to neuron ii is also proportional to Ak(t)A_{k}(t); cf. Eq. (13.31). The stationary states and their stability in a network of adaptive model neurons can therefore be analyzed as follows.

(i) Determine the stationary state A0A_{0} self-consistently. To do so, we use the gain function gσ(I0)g_{\sigma}(I_{0}) for our neuron model of choice, where the mean current I0=μ/RI_{0}=\mu/R and the noise level σ\sigma depend on the activity A0A_{0}. The gain function gσ(I0)g_{\sigma}(I_{0}) of adaptive nonlinear integrate-and-fire neurons can be found using the methods discussed above.

(ii) Determine the response to periodic modulation of input current and input variance, G^μ(ω)\hat{G}_{\mu}(\omega) and G^σ(ω)\hat{G}_{\sigma}(\omega), respectively, using the methods discussed above.

(iii) Use the inverse Fourier transform to find the linear response filters Gμ(s)G_{\mu}(s) and Gσ(s)G_{\sigma}(s) that describes the population activity with respect to small perturbations of the input current ΔI(t)\Delta I(t) and to small perturbations in the noise amplitude σ(t)\sigma(t)

A(t)=A0+0Gμ(s)ΔI(t-s)ds+0Gσ(s)Δσ(t-s)ds.A(t)=A_{0}+\int_{0}^{\infty}G_{\mu}(s)\,\Delta I(t-s)\,{\text{d}}s+\int_{0}^{% \infty}G_{\sigma}(s)\,\Delta\sigma(t-s)\,{\text{d}}s\,. (13.61)

(iv) Exploit the fact that the current I(t)I(t) and its fluctuations σ(t)\sigma(t) in a network with self-coupling, Eq. (13.31), are proportional to the current activity A(t)A(t). If we set ΔA(t)=A(t)-A0(t)\Delta A(t)=A(t)-A_{0}(t) we therefore have for the case of a single population feeding its activity back to itself

ΔA(t)=0[JμGμ(s)+JσGσ(s)]ΔA(t-s)ds,\Delta A(t)=\int_{0}^{\infty}[J_{\mu}\,G_{\mu}(s)+J_{\sigma}\,G_{\sigma}(s)]\,% \Delta A(t-s)\,{\text{d}}s, (13.62)

where the constants JμJ_{\mu} and JσJ_{\sigma} depend on the coupling parameters.

(v) Search for solutions ΔA(t)=ϵexp[λ(ω)t]cos(ωt)\Delta A(t)=\epsilon\,\exp[\lambda(\omega)t]\cos(\omega t). The stationary state of asynchronous firing with activity A0A_{0} is unstable, if there exists a frequency for which λ(ω)>0\lambda(\omega)>0.

The arguments in steps (i) to (v) do not rely on the assumption of current-based synapses. In fact, as we will see now, conductance-based synapses can be, in the stationary state, well approximated by an equivalent current-based scheme.

13.6.3 Conductance input vs. current input

Throughout this chapter we assumed that spike firing by a presynaptic neuron jj at time tj(f)t_{j}^{(f)} generates in the postsynaptic neuron ii an excitatory or inhibitory postsynaptic current pulse wijδ(t-tj(f))w_{ij}\,\delta(t-t_{j}^{(f)}). However, synaptic input is more accurately described as a change in conductance g(t-tj(f))g(t-t_{j}^{(f)}), rather than as current injection (126). As mentioned in Eq. (13.58), a time dependent synaptic conductance leads to a total synaptic current into neuron ii

Ii(t)=jfwijgij(t-tj(f))(ui(t)-Esyn),I_{i}(t)=\sum_{j}\sum_{f}w_{ij}g_{ij}(t-t_{j}^{(f)})\,(u_{i}(t)-E_{\rm syn}), (13.63)

which depends on the momentary difference between the reversal potential EsynE_{\rm syn} and the membrane potential ui(t)u_{i}(t) of the postsynaptic neuron. A spike fired by a presynaptic neuron jj at time tj(f)t_{j}^{(f)} can therefore have a bigger or smaller effect, depending on the state of the postsynaptic neuron; cf. Ch. 3.

Nevertheless we will now show that, in the state of stationary asynchronous activity, conductance-based input can be approximated by an effective current input (285; 433; 432; 557). The main effect of conductance-based input is that the membrane time constant of the stochastically driven neuron is shorter than the ‘raw’ passive membrane time constant (50; 126).

To keep the arguments transparent, we consider NEN_{E} excitatory and NIN_{I} inhibitory leaky integrate-and-fire neurons in the subthreshold regime

Cdudt=-gL(u-EL)-gE(t)(u-EE)-gI(t)(u-EI),C{{\text{d}}u\over{\text{d}}t}=-g_{L}\,(u-E_{L})-g_{E}(t)\,(u-E_{E})-g_{I}(t)% \,(u-E_{I}), (13.64)

where CC is the membrane capacity, gLg_{L} the leak conductance and EL,EE,EIE_{L},E_{E},E_{I} are the reversal potentials for leak, excitation, and inhibition, respectively. We assume that input spikes at excitatory synapses lead to an increased conductance

gE(t)=ΔgEjfexp[-(t-tj(f))/τE](t-tj(f)).g_{E}(t)=\Delta g_{E}\,\sum_{j}\sum_{f}\exp[-(t-t_{j}^{(f)})/\tau_{E}]{% \mathcal{H}}(t-t_{j}^{(f)}). (13.65)

with amplitude ΔgE\Delta g_{E} and decay time constant τE\tau_{E}. The Heaviside step function {\mathcal{H}} assures causality in time. The sum over jj runs over all excitatory synapses. Input spikes at inhibitory synapses have a similar effect, but with jump amplitude ΔgI\Delta g_{I} and decay time constant τI\tau_{I}. We assume that excitatory and inhibitory input spikes arrive with a total rate νE\nu_{E} and νI\nu_{I}, respectively. For example, in a population of excitatory and inhibitory neurons the total excitatory rate to neuron ii would be the number of excitatory presynaptic partners CEC_{E} of neuron ii times the typical firing rate of a single excitatory neuron.

Using the methods from Chapter 8 we can calculate the mean excitatory conductance

gE,0=ΔgEνEτE,g_{E,0}=\Delta g_{E}\,\nu_{E}\,\tau_{E}\,, (13.66)

where νE\nu_{E} is the total spike arrival rate at excitatory synapses, and analogously the mean inhibitory conductance. The variance of the conductance is

σE2=0.5(ΔgE)2νEτE.\sigma_{E}^{2}=0.5(\Delta g_{E})^{2}\,\nu_{E}\,\tau_{E}. (13.67)

The mathematical analysis of conductance input proceeds in two steps. First, we write the conductance as the mean gE,0g_{E,0} plus a fluctuating component gE,f(t)=gE(t)-gE,0g_{E,f}(t)=g_{E}(t)-g_{E,0}. This turns Eq. (13.64) into a new equation

Cdudt=-g0(u-μ)-gE,f(t)(u-EE)-gI,f(t),(u-EI),C{{\text{d}}u\over{\text{d}}t}=-g_{0}\,(u-\mu)-g_{E,f}(t)\,(u-E_{E})-g_{I,f}(t% ),(u-E_{I}), (13.68)

with a total conductance g0=gL+gE,0+gI,0g_{0}=g_{L}+g_{E,0}+g_{I,0} and an input-dependent equilibrium potential

μ=gLEL+gE,0EE+gI,0EIg0.\mu={g_{L}\,E_{L}+g_{E,0}\,E_{E}+g_{I,0}\,E_{I}\over g_{0}}\,. (13.69)

We emphasize that Eq. (13.68) looks just like the original equation (13.64). The major differences are, however, that the dynamics in Eq. (13.64) is characterized by a raw membrane time constant C/gLC/g_{L} whereas Eq. (13.68) is controlled by an effective membrane time constant

τeff=Cg0=CgL+gE,0+gI,0,\tau_{\rm eff}={C\over g_{0}}={C\over g_{L}+g_{E,0}+g_{I,0}}, (13.70)

and a mean depolarization μ\mu which acts as an effective equilibrium potential (243).

In the second step we compare the momentary voltage u(t)u(t) with the effective equilibrium potential μ\mu. The fluctuating part of the conductance in Eq. (13.68) can therefore be written as

gE,f(t)(u-EE)=gE,f(t)(μ-EE)+gE,f(t)(u-μ),g_{E,f}(t)\,(u-E_{E})=g_{E,f}(t)\,(\mu-E_{E})+g_{E,f}(t)\,(u-\mu), (13.71)

and similarly for the inhibitory conductance. The second term on the right-hand side of Eq. (13.71) is small compared to the first term and needs to be dropped to arrive at a consistent diffusion approximation (432). The first term on the right-hand side of Eq. (13.71) does not depend on the membrane potential and can therefore be interpreted as the summed effects of postsynaptic current pulses. Thus, in the stationary state, a conductance-based synapse model is well approximated by a current-based model of synaptic input. However, we need to use the effective membrane time constant τeff\tau_{\rm eff} introduced above in Eq. (13.70) in the voltage equation.

Example: Response to conductance-modulating input

Suppose we have a population of uncoupled exponential integrate-and-fire neurons with an adaptation current. Using the methods discussed above in Section 13.5, we can calculate the linear response G(s)G(s) of the population to a short conductance pulse gE(t)g_{E}(t) at an excitatory synapse (Fig. 13.10). The response to some arbitrary time-dependent conductance input can then be predicted by convolving the conductance change with the filter response. The Fourier transform of the filter predicts the response to sinusoidal conductance modulation with period TT (Fig. 13.10).

Fig. 13.10: Conductance input. Frequency response of the population of exponential integrate-and-fire neurons as in Fig. 13.9, but to conductance modulation. For high frequencies, the phase slowly approaches a lag of 270o{}^{o} which can be interpreted as a phase advance of 90o{}^{o}. Solid line: numerical solution. Dashed line: analytical high-frequency response; case (i) and (ii) correspond to low noise and high noise, respectively. adapted from (434).

13.6.4 Colored Noise (*)

In the previous subsection, we replaced synaptic conductance pulses by current input. If the time constants τE\tau_{E} and τI\tau_{I} of excitatory and inhibitory synapses are sufficiently short, we may approximate stochastic spike arrivals by white noise. Some synapse types, such as the NMDA component of excitatory synapses, are, however, rather slow (cf. Chapter 3). As a result of this, a spike that has arrived at an NMDA synapse at time t0t_{0} generates a fluctuation of the input current that persists for tens of milliseconds. Thus the fluctuations in the input exhibit temporal correlations, a situation that is termed colored noise, as opposed to white noise; cf. Chapter 8. Colored noise represents the temporal smoothing caused by slow synapses in a compact form.

There are two different approaches to colored noise in the membrane potential density equations.

The first approach is to approximate colored noise by white noise and replace the temporal smoothing by a broad distribution of delays. In order to keep the treatment transparent, let us assume a current-based description of synaptic input. The mean input to a neuron ii in population nn arising from other populations kk is

Ii(t)=kCnkwnk0αnk(s)A(t-s)ds,I_{i}(t)=\sum_{k}C_{nk}w_{nk}\int_{0}^{\infty}\alpha_{nk}(s)\,A(t-s)\,{\text{d% }}s\,, (13.72)

where CnkC_{nk} is the number of presynaptic partner neurons in population kk, wnkw_{nk} is the typical weight of a connection from kk to a neuron in population nn and αnk(s)\alpha_{nk}(s) is the time course of a synaptic current pulse caused by a spike fired in population kk at s=0s=0. Suppose CnkC_{nk} is a large number, e.g., Cnk=C_{nk}=1 000, but the population kk itself is at least 10 times larger (e.g., Nk=N_{k}=10 000) so that the connectivity Cnk/Nk1C_{nk}/N_{k}\ll 1. We now replace the broad current pulses α(s)\alpha(s) by short pulses qδ(s-Δ)q\,\delta(s-\Delta) where δ\delta denotes the Dirac δ\delta-function and q=0α(s)dsq=\int_{0}^{\infty}\alpha(s)\,{\text{d}}s and Δ>0\Delta>0 is a transmission delay. For each of the CnkC_{nk} connections we randomly draw the delay Δ\Delta from a distribution p(Δ)=α(Δ)/qp(\Delta)=\alpha(\Delta)/q. Because of the low connectivity, we may assume that the firing of different neurons is uncorrelated. The mean input current Ii(t)I_{i}(t) to neuron ii is then given again by Eq. (13.72), with fluctuations around the mean that are approximately white, because each spike arrival causes only a momentary current pulse. The broad distribution of delays stabilizes the stationary state of asynchronous firing (79).

The second approach consists in an explicit model of the synaptic current variables. To keep the treatment transparent and minimize the number of indices, we focus on a single population coupled to itself and suppose that the synaptic current pulses are exponential α(s)=(q/τq)exp(-s/τq)\alpha(s)=(q/\tau_{q})\,\exp(-s/\tau_{q}). The driving current of a neuron ii in a population nn arising from spikes of the same population is then described by the differential equation

dIidt=-Iiτq+CnnwnnqτqAn(t){{\text{d}}I_{i}\over{\text{d}}t}=-{I_{i}\over\tau_{q}}+C_{nn}\,w_{nn}{q\over% \tau_{q}}\,A_{n}(t) (13.73)

which we can verify by taking the temporal derivative of Eq. (13.72). As before the population activity An(t)A_{n}(t) can be decomposed into a mean μ(t)\mu(t) (which is the same for all neurons) and a fluctuating part ξi(t)\xi_{i}(t) with white-noise characteristics:

dIidt=-Iiτq+μ(t)+ξi(t).{{\text{d}}I_{i}\over{\text{d}}t}=-{I_{i}\over\tau_{q}}+\mu(t)+\xi_{i}(t)\,. (13.74)

However, Eq. (13.74) does not lead directly to spike firing but needs to be combined with the differential equation for the voltage

τmduidt=f(ui)+RIi(t),\tau_{m}\,\frac{{\text{d}}u_{i}}{{\text{d}}t}=f(u_{i})+R\,I_{i}(t)\,, (13.75)

and the reset condition: if u=θresetu=\theta_{\rm reset} then u=uru=u_{r}. Since we now have two coupled differential equations, the momentary state of a population of NN neurons is described by a two-dimensional density p(u,I)p(u,I).

We recall that, in the case of white noise, the membrane potential density at threshold vanishes. The main insight for the mathematical treatment of the membrane potential density equations in two dimensions is that the density at threshold p(θreset,I(t))p(\theta_{\rm reset},I(t)) is finite, whenever the momentary slope of the voltage du/dtRI(t)+f(θreset){\text{d}}u/{\text{d}}t\propto R\,I(t)+f(\theta_{\rm reset}) is positive (154).