13 Continuity Equation and the Fokker-Planck Approach 13 Continuity Equation and the Fokker-Planck Approach 13.2 Stochastic spike arrival

13.1 Continuity equation

In a population of neurons, each neuron may be in a different internal state. In this section we derive partial differential equations that describe how the distribution of internal states evolves as a function of time. We start in Section 13.1.1 with a population of integrate-and-fire neurons. Since the state of an integrate-and-fire neuron is characterized by its membrane potential, we describe the dynamics of the population as the evolution of membrane potential densities.

Population activity, $A(t)$ , was introduced in Chapter 7 as the fraction of neurons that fire at time $t$ in a finite population. In this chapter and the next, the population activity is the expected population activity $A(t)\equiv\langle A(t)\rangle$ . Since A(t) is self-averaging, it can also be said that we consider the limit of large populations. The finite-size effects will be discussed in Section 14.6.

The formulation of the dynamics of a population of integrate-and-fire neurons on the level of membrane potential densities has been developed by Abbott and van Vreeswijk (4), Brunel and Hakim (78), Fusi and Mattia (163), Nykamp and Tranchina (367), Brunel (79) and Omurtag et al. (371). The closely related formulation in terms of refractory densities has been studied by Wilson and Cowan (552), Gerstner and van Hemmen (180), Bauer and Pawelzik (41), and Gerstner (183). Generalized density equations have been discussed by Knight (265) and Fourcaud and Brunel (154).

13.1.1 Distribution of membrane potentials

We study a homogeneous population of integrate-and-fire neurons. The internal state of a neuron $i$ is determined by its membrane potential which changes according to

\tau_{m}{{\text{d}}\over{\text{d}}t}u_{i}=f(u_{i})+R\,I_{i}(t)\quad{\rm for}% \quad u_{i}<\theta_{\rm reset}\,.

(13.1)

Here $R$ is the input resistance, $\tau_{m}=R\,C$ the membrane time constant, and $I_{i}(t)$ the total input (external driving current and synaptic input). If $u_{i}\geq\theta_{\rm reset}$ the membrane potential is reset to $u_{i}=u_{r}<\theta_{\rm reset}$ . Here $f(u_{i})$ is an arbitrary function of $u$ ; cf. Eq. (5.2) in Chapter 5. For $f(u_{i})=-(u_{i}-u_{\rm rest})$ the equation reduces to the standard leaky integrate-and-fire model. For the moment we keep the treatment general and restrict it to the leaky integrate-and-fire model only later.

In a population of $N$ integrate-and-fire neurons, we may ask how many of the neurons have at time $t$ a given membrane potential. For $N\to\infty$ the fraction of neurons $i$ with membrane potential $u_{0}<u_{i}(t)\leq u_{0}+\Delta u$ is

\lim_{N\to\infty}\left\{{{\rm neurons~{}with~{}}u_{0}<u_{i}(t)\leq u_{0}+% \Delta u\over N}\right\}=\int_{u_{0}}^{u_{0}+\Delta u}p(u,t)\,{\text{d}}u\,,

(13.2)

where $p(u,t)$ is the membrane potential density; cf. Chapter 8. The integral over this density remains constant over time, i.e.,

\int_{-\infty}^{\theta_{\rm reset}}p(u,t)\,{\text{d}}u=1\,.

(13.3)

The normalization to unity expresses the fact that all neurons have a membrane potential below or equal to threshold.

The aim of this section is to describe the evolution of the density $p(u,t)$ as a function of time. As we will see, the equation that describes the dynamics of $p(u,t)$ is nearly identical to that of a single integrate-and-fire neuron with diffusive noise; cf. Eqs. (8.40) and (8.41) in Chapter 8.

13.1.2 Flux and continuity equation

Let us consider the portion of neurons with a membrane potential between $u_{0}$ and $u_{1}$ ,

{n(u_{0};u_{1})\over N}=\int_{u_{0}}^{u_{1}}p(u^{\prime},t)\,{\text{d}}u^{% \prime}\,.

(13.4)

The fraction of neurons with $u_{0}<u<u_{1}$ increases if neurons enter from below through the boundary $u_{0}$ or from above through the boundary $u_{1}$ ; cf. Fig. 13.1. Since there are many neurons in the population, we expect that in each short time interval $\Delta t$ , many trajectories cross one of the boundaries. The flux $J(u,t)$ is the net fraction of trajectories per unit time that crosses the value $u$ . A positive flux $J(u,t)>0$ is defined as a flux toward increasing values of $u$ . In other words, in a finite population of $N$ neurons, the quantity $N\,J(u_{0},t)\,\Delta t$ describes the number of trajectories that cross in the interval $\Delta t$ the boundary $u_{0}$ from below, minus the number of trajectories crossing $u_{0}$ from above. Note that $u_{0}$ here is an arbitrary value that we chose as the lower bound of the integral in Eq. (13.4) and as such has no physical meaning for the neuron model.

Fig. 13.1: The number of trajectories in the interval $[u_{0},u_{1}]$ changes if one of the trajectories crosses the boundary $u_{0}$ or $u_{1}$ . For a large number of neurons this fact is described by the continuity equation; cf. Eq. (13.6). Schematic figure where only three trajectories are shown.

Since trajectories cannot simply end, a change in the number $n(u_{0};u_{1})$ of trajectories in the interval $u_{0}<u<u_{1}$ can be traced back to the flux of trajectories in and out of that interval. We therefore have the conservation law

{\partial\over\partial t}\int_{u_{0}}^{u_{1}}p(u^{\prime},t)\,{\text{d}}u^{% \prime}=J(u_{0},t)-J(u_{1},t)\qquad\,.

(13.5)

Taking the derivative with respect to the upper boundary $u_{1}$ and changing the name of the variable from $u_{1}$ to $u$ yields the continuity equation,

\qquad{\partial\over\partial t}p(u,t)=-{\partial\over\partial u}J(u,t)\quad{% \rm~{}for~{}}u\neq u_{r}{\rm~{}and~{}}u\neq\theta_{\rm reset}\,,

(13.6)

which expresses the conservation of the number of trajectories. In integrate-and-fire models, however, there are two special voltage values, $u_{r}$ and $\theta_{\rm reset}$ , where the number of trajectories is not conserved, because of the fire-and-reset mechanism.

Since neurons that have fired start a new trajectory at $u_{r}$ , we have a ‘source of new trajectories’ at $u=u_{r}$ , i.e., new trajectories appear in the interval $[u_{r}-\epsilon,u_{r}+\epsilon]$ that have not entered the interval through one of the borders. Adding a term $A(t)\,\delta(u-u_{r})$ on the right-hand side of (13.6) accounts for this source of trajectories. The trajectories that appear at $u_{r}$ disappear at $\theta_{\rm reset}$ , so that we have

\qquad{\partial\over\partial t}p(u,t)=-{\partial\over\partial u}J(u,t)+A(t)\,% \delta(u-u_{r})-A(t)\,\delta(u-\theta_{\rm reset})\,.

(13.7)

The density $p(u,t)$ vanishes for all values $u>\theta_{\rm reset}$ .

The population activity $A(t)$ is, by definition, the fraction of neurons that fire, i.e., those that pass through the threshold. Therefore we find

A(t)=J(\theta_{\rm reset},t)\,.

(13.8)

Equations (13.7) and (13.8) describe the evolution of the the membrane potential densities and the resulting population activity as a function of time. We now specify the neuron model so as to have an explicit expression for the flux.