## Abstract

We demonstrate the mathematical equivalence of two commonly used forms of firing rate model equations for neural networks. In addition, we show that what is commonly interpreted as the firing rate in one form of model may be better interpreted as a low-pass-filtered firing rate, and we point out a conductance-based firing rate model.

At least since the pioneering work of Wilson & Cowan (1972), it has been common to study neural circuit behavior using rate equations—equations that specify neural activities simply in terms of their rates of firing action potentials, as opposed to spiking models, in which the actual emissions of action potentials, or spikes, are modeled. Rate models can be derived as approximations to spiking models in a variety of ways (Wilson & Cowan, 1972; Mattia & Del Giudice, 2002; Shriki, Hansel, & Sompolinsky, 2003; Ermentrout, 1994; La Camera, Rauch, Luscher, Senn, & Fusi, 2004; Aviel and Gerstner, 2006; Ostojic & Brunel, 2011; reviewed in Ermentrout & Terman, 2010; Gerstner & Kistler, 2002; and Dayan & Abbott, 2001).

**v**-equation and

**r**-equation respectively: Here,

**v**and

**r**are each vectors representing neural activity, with each element representing the activity of one neuron in the modeled circuit.

**v**is commonly thought of as representing voltage, while

**r**is commonly thought of as representing firing rate (probability of spiking per unit time).

**f**(

**x**) is a nonlinear input-output function that acts element-by-element on the elements of

**x**, that is, it has

*i*th element (

**f**(

**x**))

_{i}=

*f*(

*x*) for some nonlinear function of one variable

_{i}*f*.

*f*typically takes such forms as an exponential, a power law, or a sigmoid function, and

*f*(

*v*) is typically regarded as a static nonlinearity converting the voltage of the

_{i}*i*th cell

*v*to the cell's instantaneous firing rate.

_{i}**W**is the matrix of synaptic weights between the neurons in the modeled circuit. and

**I**are the vectors of external inputs to the neurons in the

**v**or

**r**networks, respectively, which may be time dependent. In the appendix, we illustrate a simple heuristic derivation of the

**v**-equation, starting from the biophysical equation for the voltages

**v**. Along the way, we also point to a conductance-based version of the rate equation.

When developing a rate model of a network, it can be unclear which form of equation to use or whether it makes a difference. Here we demonstrate that the choice between equations 1 and 2 makes no difference: the two models are mathematically equivalent, and so will display the same set of behaviors. It has been noted previously (Beer, 2006) that when **I** is constant and **W** is invertible, the two equations are equivalent under the relationship **v** = **Wr** + **I**, . We generalize this result to demonstrate the equivalence of the two equations when **W** is not invertible and inputs may be time dependent.

The **v**-equation is defined when we specify the input across time, , and the initial condition **v**(0); we will call the combination of these and equation 1 a **v**-model. The **r**-equation is defined when we specify **I**(*t*) and **r**(0); we will call the combination of these and equation 2 an **r**-model. We will show that any **v**-model can be mapped to an **r**-model and any **r**-model can be mapped to a **v**-model such that the solutions to equations 1 and 2 satisfy **v** = **Wr** + **I**.

As we will see, the inputs in equivalent models are related by , or . That is, **I** is a low-pass-filtered version of . Note that there is an equivalence class of **I**, parameterized by **I**(0), that all correspond to the same under this equivalence. We assume that the equivalence class has been specified, that is, has been specified (if **I** has been specified, can be found as ). Then a **v**-model is defined by specifying **v**(0), while an **r**-model is defined by specifying the set {**r**(0), **I**(0)}. If **W** is *D* × *D*, then **v**(0) is *D*-dimensional, while {**r**(0), **I**(0)} is 2*D*-dimensional, so we can guess that the map from **r** to **v** takes a D-dimensional space of **r**-models to a single **v**-model, and conversely the map from **v** to **r** takes a single **v**-model back to a D-dimensional space of **r**-models, and we will show that this is true.

**r**evolves according to the

**r**-equation, then

**Wr**+

**I**evolves according to the

**v**-equation. Setting

**v**=

**Wr**+

**I**, we find: Therefore, if

**v**evolves according to the

**v**-equation and

**r**evolves according to the

**r**-equation and

**v**(0) =

**Wr**(0) +

**I**(0), then, since the

**v**-equation propagates

**Wr**+

**I**forward in time,

**v**=

**Wr**+

**I**at all times

*t*> 0. We will thus have established the desired equivalence if we can solve

**v**(0) =

**Wr**(0) +

**I**(0) for any

**v**-model, specified by

**v**(0), or for any

**r**-model, specified by {

**r**(0),

**I**(0)}.

Note that, as expected, a D-dimensional space of **r**-models converges on the same **v**-model. Since {**r**(0), **I**(0)} forms a 2D-dimensional space, which is constrained by the *D*-dimensional equation **v**(0) = **Wr**(0) + **I**(0), the *D*-dimensional subspace of **r**-models {**r**(0), **I**(0)} that satisfy this equation all converge on the same **v**-model.

To go from an **r**-model to a **v**-model is straightforward: we simply set **v**(0) = **Wr**(0) + **I**(0). To go from a **v**-model to an **r**-model, we first define some useful notation:^{1}

- •.
is the null space of

**W**, that is, the subspace of all vectors that**W**maps to 0.**P**_{N}is the projection operator into . - •.
is the subspace perpendicular to . This is the subspace spanned by the rows of

**W**.**P**_{N⊥}is the projection operator into . - •.
is the range of

**W**, that is, the subspace of vectors that can be written**Wx**for some**x**. This is the subspace spanned by the columns of**W**.**P**_{R}is the projection operator into . - •.
is the subspace perpendicular to , also called the left null space.

**P**_{R⊥}is the projection operator into .

For any vector **x**, we define **x**_{N} ≡ **P**_{N}**x**, **x**_{N⊥} ≡ **P**_{N⊥}**x**, **x**_{R} ≡ **P**_{R}**x**, **x**_{R⊥} ≡ **P**_{R⊥}**x**. We rely on the fact that **x** = **x**_{N} + **x**_{N⊥} = **x**_{R} + **x**_{R⊥}.

**v**-model, the equation

**v**(0) =

**Wr**(0) +

**I**(0) has a solution if and only if , which is true if and only if

**v**

_{R⊥}(0) −

**I**

_{R⊥}(0) = 0,

^{2}so we must choose Letting

*D*be the dimension of and

_{R}*D*the dimension of , the fundamental theorem of linear algebra states that

_{N}*D*+

_{R}*D*=

_{N}*D*. So

**I**

_{R⊥}(0) has dimension

*D*. This leaves unspecified

_{N}**I**

_{R}(0), which has dimension

*D*.

_{R}**r**

_{N⊥}(0), we note that the equation

**v**=

**Wr**+

**I**can equivalently be written

**v**=

**Wr**

_{N⊥}+

**I**(because

**Wr**

_{N}= 0, so

**Wr**=

**Wr**

_{N⊥}). That is, knowledge of

**v**specifies only

**r**

_{N⊥}. We define

**W**

^{−1}to be the Moore-Penrose pseudo-inverse of

**W**. This is the matrix that gives the one-to-one mapping of into that inverts the one-to-one mapping of to induced by

**W**, and that maps all vectors in to 0.

^{3}The pseudo-inverse has the property that

**W**

^{−1}

**W**=

**P**

_{N⊥}while

**WW**

^{−1}=

**P**

_{R}. Then we can solve for

**r**

_{N⊥}(0) as This is a

*D*-dimensional equation for the 2

_{R}*D*-dimensional set of unknowns {

_{R}**r**

_{N⊥}(0),

**I**

_{R}(0)}, so it determines

*D*of these parameters and leaves

_{R}*D*free. For example, it could be solved by freely choosing

_{R}**I**

_{R}(0) and then setting

**r**

_{N⊥}(0) =

**W**

^{−1}(

**v**

_{R}(0) −

**I**

_{R}(0)), or by freely choosing

**r**

_{N⊥}(0) and then setting

**I**

_{R}(0) =

**v**

_{R}(0) −

**Wr**

_{N⊥}(0).

Equations 10 and 9 together ensure the equality **v**(0) = **Wr**(0) + **I**(0). Applying **W** to both sides of equation 10 yields **v**_{R}(0) = **Wr**_{N⊥}(0) + **I**_{R}(0) = **Wr**(0) + **I**_{R}(0). This states that the equality holds within the range of **W**; orthogonal to the range of **W**, we have **P**_{R⊥}**Wr** = 0 and **v**_{R⊥}(0) = **I**_{R⊥}(0). Together, these yield **v**(0) = **Wr**(0) + **I**(0).

Finally, we can freely choose **r**_{N}(0), which has no effect on the equation **v**(0) = **Wr**(0) + **I**(0). **r**_{N}(0) has *D _{N}* dimensions, so we have freely chosen

*D*+

_{R}*D*=

_{N}*D*dimensions in finding an

**r**-model that is equivalent to the

**v**-model. That is, we have found a D-dimensional subspace of such

**r**-models—those that satisfy

**v**(0) =

**Wr**(0) +

**I**(0).

To summarize, we have established the equivalence between **r**-models and **v**-models. For each fixed choice of **W**, τ, and , an **r**-model is specified by {**r**(0), **I**(0)} and equation 2, while a **v**-model is specified by **v**(0) and equation 1. The equivalence is established by setting **v**(0) = **Wr**(0) + **I**(0), which yields a *D*-dimensional subspace of equivalent **r**-models for a given **v**-model. Under this equivalence, **v** obeys equation 1, **r** obeys equation 2, and the two are related at all times by **v** = **Wr** + **I**, with . To go from an **r**-model to its equivalent **v**-model, we simply set **v**(0) = **Wr**(0) + **I**(0). To go from a **v**-model to one of its equivalent **r**-models, we set **I**_{R⊥}(0) = **v**_{R⊥}(0), freely choose **r**_{N}(0), and freely choose {**r**_{N⊥}(0), **I**_{R}(0)} from the *D _{R}*-dimensional subspace of such choices that satisfy

**r**

_{N⊥}(0) =

**W**

^{−1}(

**v**

_{R}(0) −

**I**

_{R}(0)), where

**W**

^{−1}is the pseudoinverse of

**W**.

Finally, note that equation 2 can be written . That is, if we regard **v** as a voltage and *f*(**v**) as a firing rate, as suggested by the “derivation” in the appendix, then **r** is a low-pass-filtered version of the firing rate, just as **I** is a low-pass-filtered version of the input .

## Appendix: Simple “Derivation” of the **v**-Equation

**v**-equation can be “derived” as follows. We start with the equation for the membrane voltage of the

*i*th neuron: where

*C*is the capacitance of the

_{i}*i*th neuron and

*g*is the

_{ij}*j*th conductance onto the neuron, with reversal potential

*E*. We assume that the

_{ij}*g*’s are composed of an intrinsic conductance,

_{ij}*g*, with reversal potential

^{L}_{i}*E*; extrinsic input with reversal potential ; and within-network synaptic conductances, with representing input from neuron

^{L}_{i}*j*with reversal potential . Dividing by ∑

_{k}

*g*and defining τ

_{ik}_{i}(

*t*) =

*C*/∑

_{i}_{k}

*g*gives

_{ik}*r*of neuron

_{j}*j*, with proportionality constant : . This ignores synaptic time courses, among other things. We assume that

*r*is given by the static nonlinearity

_{j}*r*=

_{j}*f*(

*v*) (see Miller & Troyer, 2002; Hansel & van Vreeswijk, 2002; Priebe, Mechler, Carandini, & Ferster, 2004, for such a relationship between firing rate and voltage averaged over a few tens of milliseconds). We assume synapses are either excitatory with reversal potential

_{j}*E*or inhibitory with reversal potential

_{E}*E*, and linearly transform the units of voltage so that

_{I}*E*= 1 and

_{E}*E*= −1. We define . This is now a synaptic weight that is positive for excitatory synapses and negative for inhibitory synapses. We define and define . This yields the conductance-based rate equation, with τ

_{I}_{i}(

*t*) =

*C*/(

_{i}*g*+ ∑

_{i}_{k}|

*W*|

_{ik}*f*(

*v*)).

_{k}*g*is much larger than synaptic and external conductances or if inputs tend to be push-pull, with withdrawal of some inputs compensating for addition of others. We absorb the constant denominator into the definitions of and

^{L}_{i}*W*and note that this also implies that τ

_{ij}_{i}is constant, to arrive finally at the

**v**-equation:

## Acknowledgments

This work was supported by R01-EY11001 from the National Eye Institute and by the Gatsby Charitable Foundation through the Gatsby Initiative in Brain Circuitry at Columbia University.

## Notes

^{1}

If **W** is normal, the eigenvectors are orthogonal, so the null space is precisely the space orthogonal to the range: **P**_{N} = **P**_{R⊥} and **P**_{N⊥} = **P**_{R}. However, if **W** is nonnormal, then vectors orthogonal to the null space can be mapped into the null space; the range always has the dimension of the full space minus the dimension of the null space, but it need not be orthogonal to the null space.

^{2}

**v**=

**Wr**+

**I**can be solved, is true for all time if it is true in the initial condition. We compute: Applying

**P**

_{R⊥}to equation 7 and noting that

**P**

_{R⊥}

**W**= 0, we find If , then

**v**

_{R⊥}(0) −

**I**

_{R⊥}(0) = 0, and hence

**v**

_{R⊥}−

**I**

_{R⊥}= 0 at all subsequent times so at all subsequent times. Note also that for any initial conditions, the condition is true asymptotically as

*t*→ ∞.

^{3}

If the singular value decomposition of a matrix **M** is **M** = **USV**†, where **S** is the diagonal matrix of singular values and **U** and **v** are unitary matrices, then its pseudo-inverse is , where is the pseudoinverse of **S**, obtained by inverting all nonzero singular values in **S**.