Traveling pulses in a stochastic neural field model of direction selectivity

Introduction

Continuum neural field models represent the large-scale dynamics of spatially structured networks of neurons in terms of non-linear integro-differential equations, whose associated integral kernels represent the spatial distribution of neuronal synaptic connections (Wilson and Cowan, 1972, 1973; Amari, 1977). As in the case of non-linear partial differential equation (PDE) models of diffusively coupled excitable systems (Keener, 1981; Kuramoto, 1984), non-local neural fields can exhibit a diverse range of spatiotemporal dynamics, including solitary traveling fronts and pulses, stationary pulses, and spatially localized oscillations (breathers), spiral waves, and Turing-like patterns. See, for example, the reviews Ermentrout (1998), Coombes (2005), and Bressloff (2012). In recent years, neural fields have been used to model a wide range of neurobiological phenomena, including wave propagation in cortical slices (Pinto and Ermentrout, 2001; Richardson et al., 2005) and in vivo (Huang et al., 2004), geometric visual hallucinations (Ermentrout and Cowan, 1979; Bressloff et al., 2001), EEG rhythms (Nunez, 1995; Robinson et al., 2001; Liley et al., 2002; Steyn-Ross et al., 2003), orientation tuning in primary visual cortex (V1) (Ben-Yishai et al., 1995; Somers et al., 1995), short term working memory (Camperi and Wang, 1998; Laing and Chow, 2002), control of head direction (Zhang, 1996), direction selectivity (Xie and Giese, 2002), motion perception (Giese, 1999), and binocular rivalry waves (Bressloff and Webber, 2012a). One particularly useful feature of neural fields is that analytical techniques for solving these integro-differential equations can be adapted from previous studies of non-linear PDEs. These include regular and singular perturbation methods, weakly non-linear analysis and pattern formation, symmetric bifurcation theory, Evans functions and wave stability, and homogenization theory (Bressloff, 2012). In particular, we have recently shown how perturbation methods for studying fluctuating fronts in reaction-diffusion PDEs (Schimansky-Geier et al., 1983; de Pasquale et al., 1992; Armero et al., 1998; Sagues et al., 2007) can be extended to the problem of front propagation in stochastic neural fields (Bressloff and Webber, 2012b), and have used this to investigate the effects of noise on binocular rivalry waves (Webber and Bressloff, submitted). Such methods exploit a separation of time-scales in which there is a diffusive-like displacement (wandering) of the front from its uniformly translating position at long time-scales, and fluctuations in the front profile around its instantaneous position at short time-scales.

In this paper, we extend our theory of wave propagation in stochastic neural fields to the case of a neural field that supports traveling pulses rather than fronts. A typical mechanism for generating traveling pulses in an excitatory network is to include some form of slow adaptation, such as spike frequency adaptation (Pinto and Ermentrout, 2001; Coombes and Owen, 2005) or synaptic depression (Kilpatrick and Bressloff, 2010a,b), which suppresses the trailing edge of the wave. One of the motivations for considering excitatory neural fields is that traveling pulses are observed in in vitro cortical slices that have been disinhibited. Here we consider an alternative mechanism for generating pulses, based on asymmetric excitatory/inhibitory synaptic connections. Such a network architecture supports freely propagating pulses without any adaptation, and has been proposed as a simple recurrent mechanism for generating direction selectivity in a network driven by moving stimuli (Mineiro and Zipser, 1998; Xie and Giese, 2002). Most classical models for the direction selectivity of cortical neurons are based on feedforward mechanisms, namely, the linear or non-linear spatiotemporal filtering of afferent thalamo-cortical inputs (Reichardt, 1961; Adelson and Bergen, 1985; Koch and Poggio, 1985; van Santen and Sperling, 1985). Some of these models also involve a combination of lagged (time-delayed) and non-lagged inputs (Saul and Humphrey, 1990; Baker and Bair, 2012). However, there is now considerable experimental data demonstrating that the response of cortical cells is strongly influenced by intracortical circuitry. This has motivated a number of modeling studies that show how direction selectivity can be reproduced by recurrent neural network models with asymmetric lateral excitatory or inhibitory connections and non-direction-selective inputs (Suarez et al., 1995; Maex and Urban, 1996; Mineiro and Zipser, 1998; Xie and Giese, 2002). In this paper, we base our investigation of stochastic traveling pulses on the particular version introduced by Xie and Giese (2002).

The main results of the paper are as follows. We first analyze freely propagating pulses and stimulus-locked pulses in the deterministic case, expanding the analysis of Xie and Giese (2002). In particular, we construct a stability diagram showing the existence and stability of stimulus-locked pulses as a function of stimulus velocity and amplitude. We then turn to a corresponding stochastic version of the model. We show how extrinsic noise in the activity variable leads to a diffusive-like displacement (wandering) of the wave from its uniformly translating position at long time-scales, and fluctuations in the wave profile around its instantaneous position at short time-scales. In the case of freely propagating pulses, the wandering is characterized by pure Brownian motion, whereas in the case of stimulus-locked pulses, it is given by an Ornstein–Uhlenbeck process. This establishes that stimulus-locked pulses are more robust to noise. One major difference between pulses and fronts is that, in principle, noise could significantly affect both the location (center-of-mass) and width of the pulse. We find that fluctuations in the width can be neglected in the case of freely propagating pulses, whereas the saturation of the mean-square displacement of the center-of-mass of the pulse for stimulus-locked pulses means that fluctuations in pulse width can no longer be ignored.

Materials and Methods

Neural Field Model of Direction Selectivity

In this paper we consider a scalar neural field equation of the form

$\begin{array}{cc}\tau \frac{\partial u(x, t)}{\partial t}=-u(x, t)+{\displaystyle {\int }_{-\infty }^{\infty }w}(x-x^{\prime })F(u(x^{\prime }, t))d{x}^{\prime }+h(x, t)& (1)\end{array}$

Here u(x, t) is a measure of activity (current or voltage) within a local population of excitatory and inhibitory neurons at position x ∈ ℝ and time t, τ is a membrane time constant (of order 10 msec), w(x) denotes the spatial distribution of synaptic connections between local populations, F(u) is a non-linear firing rate function and h(x, t) is an external input. (We fix the time-scale by setting τ = 1). F is usually taken to be a sigmoid function

$\begin{array}{cc}F(u)=\frac{1}{1+e^{-\gamma (u-\kappa )}}& (2)\end{array}$

with gain γ and threshold κ. In the high-gain limit γ → ∞, this reduces to the Heaviside function

$\begin{array}{cc}F(u)\to H(u-\kappa )=\{\begin{array}{cc}1& ifu>\kappa \\ 0& ifu\le \kappa .\end{array}& (3)\end{array}$

The function w(x − x') represents the distribution of synaptic weights from the local population at x' to the population at x. Usually, w is taken to be a symmetric or even function such that w(x) = w(−x). A common choice for the weight distribution is a “Mexican hat” function, with a center excitatory region surrounded by flanking inhibitory regions. As originally shown by Amari (1977), symmetric Mexican hat functions tend to support stationary activity “bumps.” Following Xie and Giese (2002), however, we will use an asymmetric Mexican hat function whose maximum is offset by an amount x₀, that is w(x−x₀) = w(−[x−x₀]); the resulting neural field then supports freely propagating pulses that depend on the degree of offset. Note that such a choice should be contrasted with a symmetric function w with peaks offset from zero see e.g., (Hutt and Atay, 2005). In the case of exponential functions, w takes the form (see Figure 1)

$\begin{array}{cc}w(x)=a_e{\text{e}}^{-{\sigma }_e|x-x_0|}-a_i{\text{e}}^{-{\sigma }_i|x-x_0|}, & (4)\end{array}$

FIGURE 1

Figure 1. Plot of weight distribution for various values of the shift x₀. Here a_e = 5, a_i = 1, σ_e = 0.42, and σ_i = 0.1.

where a_e > a_i and σ_e > σ_i. Setting x₀ = 0 recovers the standard Mexican hat function. Note that one could equally use other functions such as a difference-of-Gaussians without changing the main results of the paper; the advantage of exponentials is that one can carry out explicit calculations.

Finally, the external input h(x, t) consists of two components:

$\begin{array}{cc}h(x, t)=I(x-vt)+\sqrt{ϵ}g(u(x, t))\xi (x, t).& (5)\end{array}$

Here I(x − vt) represents an external pulse-like stimulus moving with constant speed v and amplitude I₀, whereas the second term represents an extrinsic, multiplicative noise source. In particular, ξ(x, t) is a Gaussian process with zero mean and two-point correlations

$\begin{array}{cc}\langle \eta (x, t)\eta (x^{\prime }, t^{\prime })\rangle =2C(|x-x^{\prime }|/\lambda )\delta (t-t^{\prime }).& (6)\end{array}$

Thus the noise is white in time and colored in space with correlation length λ. Formally speaking, η(x, t)dt = dW(x, t) where dW(x, t) is a corresponding Wiener process. The amplitude of the noise is determined by the parameter ∈, and the function g(u) incorporates any activity-dependence. Note that (Xie and Giese, 2002) only considered the deterministic case (∈ = 0). They showed how the deterministic neural field supports freely propagating pulses of fixed speed c when I₀ = 0. This then provides a mechanism for direction selectivity, since these pulses can lock to a moving stimulus of speed v provided that |c − v| is sufficiently small; the range of locking depends on the amplitude I₀. In this paper, we develop a more systematic analysis of stimulus-locking in the absence of noise, and then investigate the effects of noise on both freely propagating and stimulus-locked pulses.

Results

Determinstic Neural Field

We begin by analyzing traveling pulse solutions of the neural field Equation (1) in the absence of noise (∈ = 0). Following the original formulation of Amari (1977), we investigate the existence and stability of traveling pulses by setting the firing rate function to be the Heaviside (Equation 3).

Freely propagating pulses

For the moment, suppose that there are no external inputs so that h(x, t) = 0 in Equation (1). A traveling pulse of velocity c is then defined according to u(x, t) = U(ξ), with ξ = x − ct a traveling wave coordinate such that lim_ξ→±∞ U(ξ) = 0. Moreover, the wave profile is restricted to be super threshold in a connected interval of width d. Since the neural field is equivariant with respect to uniform translations (in the absence of external stimuli), we choose the two threshold crossing points to be

$\begin{array}{cc}U(0)=\kappa ,U(d)=\kappa .& (7)\end{array}$

Thus, U(ξ) > κ for 0 < ξ < d, U(ξ) < κ for ξ < 0, and ξ > d. It turns out the wave travels in the same direction as the offset so we restrict ourselves to the case x₀ > 0 and c > 0. Substituting the traveling pulse solution into Equation (1) gives

$\begin{array}{cc}-c\frac{\partial U(\xi )}{\partial \xi }=-U(\xi )+{\displaystyle {\int }_0^{d}w}(\xi -{\xi }^{\prime })d{\xi }^{\prime }& (8)\end{array}$

Multiplying both sides by e^−ξ/c and integrating gives the following equation for the wave solution:

$\begin{array}{cc}U(\xi )=\frac{e^{\xi /c}}{c}{\displaystyle {\int }_{\xi }^{\infty }\text{}}W({\xi }^{\prime }){\text{e}}^{-{\xi }^{\prime }/c}\text{d}{\xi }^{\prime },& (9)\end{array}$

where

$W(\xi )\equiv {\displaystyle {\int }_{\xi -d}^{\xi }\text{}}w(x)\text{d}x.$

It is convenient to express the weight function in piecewise form as follows:

$\begin{array}{cc}w(x)=\{\begin{array}{ll}{a}_e{\text{e}}^{-{\sigma }_e(x-x_0)}-a_i{\text{e}}^{-{\sigma }_i(x-x_0)}\equiv w_1(x),\hfill & \text{if }x\ge x_0\hfill \\ a_e{\text{e}}^{{\sigma }_e(x-x_0)}-a_i{\text{e}}^{{\sigma }_i(x-x_0)}\equiv w_2(x),\hfill & \text{if }x\le x_0\hfill \end{array}& (10)\end{array}$

We then obtain a piecewise expression for W(ξ) of the form

$\begin{array}{cc}W(\xi )=\{\begin{array}{ll}{W}_3(\xi )\equiv {\displaystyle {\int }_{\xi -d}^{\xi }\text{}}{w}_2(x)\text{d}x,\hfill & \text{if }\xi \le x_0\hfill \\ W_2(\xi )\equiv {\displaystyle {\int }_{\xi -d}^{x_0}\text{}}{w}_2(x)\text{d}x\hfill & \hfill \\ +{\displaystyle {\int }_{x_0}^{\xi }\text{}}{w}_1(x)\text{d}x,\hfill & \text{if }{x}_0\le \xi \le x_0+d\hfill \\ W_1(\xi )\equiv {\displaystyle {\int }_{\xi -d}^{\xi }\text{}}{w}_1(x)\text{d}x,\hfill & \text{if }\xi \ge x_0+d\hfill \end{array}& (11)\end{array}$

We then have

$U(\xi )=\{\begin{array}{ll}\frac{1}{c}{\text{e}}^{\xi /c}(M_3(\xi )+M_1(x_0+d)\hfill & \hfill \\ +M_2(x_0)),\hfill & \text{if }\xi \le x_0\hfill \\ \frac{1}{c}{\text{e}}^{\xi /c}\left(M_2(\xi )+M_1(x_0+d)\right),\hfill & \text{if }{x}_0\le \xi \le x_0+d\hfill \\ \frac{1}{c}{\text{e}}^{\xi /c}{M}_1(\xi ),\hfill & \text{if }\xi \ge x_0+d\hfill \end{array}.$

where

$\begin{array}{cc}{M}_n(\xi )={\displaystyle {\int }_{\xi }^{{\xi }_n}{W}_n}({\xi }^{\prime }){\text{e}}^{-{\xi }^{\prime }/c}d{\xi }^{\prime }& (12)\end{array}$

with ξ₁ = ∞, ξ₂ = x₀ + d and ξ₃ = x₀.

Having obtained the piecewise wave profile U(ξ), the threshold conditions (Equation 7) can now be used to determine the pulse speed c and width d; the resulting transcendental equations have to be solved numerically. Figure 2 shows solutions for the pulse speed and width as functions of the threshold. It turns out that the solution with slower speed (and larger width) is stable (see below). This differs from traveling pulse solutions found in adaptive neural fields, where the faster wave (with larger width) tends to be stable (Pinto and Ermentrout, 2001; Kilpatrick and Bressloff, 2010a). Figure 3A shows a typical pulse waveform and Figure 3B shows a numerical simulation of the neural field Equation (1) using the wave solution as the initial condition. The pulse propagates at the predicted speed without changing shape significantly. This occurs because the parameters were chosen to make the pulse solution linearly stable.

FIGURE 2

Figure 2. (A) Plots of pulse speed c and (B) pulse width d as a function of the threshold κ. Weight parameters are as in Figure 1 with offset x₀ = 3. Stable (unstable) branches are indicated by black (gray) curves.

FIGURE 3

Figure 3. (A) Plot of traveling wave profile U(ξ) obtained analytically. Same parameters as Figure 2 for threshold κ = 4. (B) Spacetime plot of a traveling pulse using the profile of (A) as the initial condition. High (low) activity indicated by light (gray).

Stability

In order to determine the linear stability of a traveling pulse solution U(ξ) in the moving frame, we linearize Equation (1) with h(x, t) = 0 by setting

$U(\xi , t)=U(\xi )+\phi (\xi , t),$

and Taylor expanding to first order in ϕ. This gives

$\begin{array}{cc}\begin{array}{c}\frac{\partial \phi (\xi , t)}{\partial t}=c\frac{\partial \phi (\xi , t)}{\partial \xi }-\phi (\xi , t)+{\displaystyle {\int }_{-\infty }^{\infty }\text{}}w(\xi -y)\\ \times F^{\prime }(U(y))\phi (y, t)\text{d}y.\end{array}& (13)\end{array}$

In the case of the Heaviside rate function (Equation 3), we have

$\begin{array}{cc}{F}^{\prime }(U(\xi ))=\frac{\delta (\xi )}{|U^{\prime }(0)|}+\frac{\delta (\xi -d)}{|U^{\prime }(d)|}.& (14)\end{array}$

Moreover, differentiating Equation (9) with respect to ξ shows that

$U^{\prime }(\xi )=\frac{1}{c}(U(\xi )-W(\xi )).$

Substituting the previous two results into Equation (13) gives

$\begin{array}{cc}\begin{array}{l}\frac{\partial \phi (\xi ,t)}{\partial t}=ℒ\phi (\xi ,t)\\ \equiv c\frac{\partial \phi (\xi ,t)}{\partial \xi }-\phi (\xi ,t)+\frac{c\phi (0,t)}{|\kappa -W(0)|}w(\xi )\\ +\frac{c\phi (d,t)}{|\kappa -W(d)|}w(\xi -d).\end{array}& (15)\end{array}$

where κ is the threshold. Looking for solutions of the form

$\begin{array}{cc}\phi (\xi ,t)={\text{e}}^{\lambda t}\phi (\xi ).& (16)\end{array}$

then leads to the spectral problem

$\begin{array}{cc}ℒ\phi (\xi ,t)=\lambda \phi (\xi ,t).& (17)\end{array}$

We take the linear operator ℒ to act on a Banach space ℬ of continuous, bounded functions ψ(ξ) that are defined for ξ ∈ ℝ, and that decay exponentially as ξ→±∞. Let σ(ℒ) denote the spectrum of the linear operator ℒ, and define the associated resolvent operator according to ℛ_λ ≡(ℒ − λ I)⁻¹, where I is the identity operator. The spectrum can be defined as those values of λ for which ${\mathcal{T}}_{\lambda }\equiv ℒ-\lambda I$ is not bijective. The spectrum is composed of three disjoint sets, the point or discrete spectrum, the residual spectrum, and the continuous spectrum. The point spectrum is defined as the values of λ (eigenvalues) for which the resolvent does not exist. The residual spectrum are the spectral values for which the resolvent exists but is not defined on a dense subset of ℬ. The continuous spectrum are the spectral values for which the resolvent exists and is densely defined but is unbounded (Kreyszig, 1978). Given these definitions, the traveling pulse is said to be linearly stable if (1) Re(λ) < 0 for all λ ın σ(ℒ), λ ≠ 0 and (2) the zero eigenvalue is simple. The existence of a zero eigenvalue with corresponding eigenfunction ϕ(ξ, t) = U′(ξ) reflects translation invariance, and immediately follows from differentiating (Equation 9) with respect to ξ.

We first consider the discrete spectrum by solving the eigenvalue equation

$\frac{d\psi (\xi )}{d\xi }-\frac{\lambda +1}{c}\psi (\xi )+K_0\psi (0)w(\xi )+K_d\psi (d)w(\xi -d)=0,$

where we have introduced the constants

$K_0=\frac{1}{|\kappa -W(0)|}\text{ and }{K}_d=\frac{1}{|\kappa -W(d)|}.$

Multiplying both sides by e^{−(λ + 1)ξ/c} and integrating gives

$\begin{array}{c}\psi (\xi )=K_0\psi (0){\displaystyle {\int }_{\xi }^{\infty }\text{}}w(y){\text{e}}^{(\lambda +1)(\xi -y)/c}\text{d}y\\ +K_d\psi (d){\displaystyle {\int }_{\xi -d}^{\infty }\text{}}w(y){\text{e}}^{(\lambda +1)(\xi -d-y)/c}\text{d}y,\end{array}$

which can be rewritten in the more compact form

$\begin{array}{cc}\psi (\xi )=K_0\psi (0)(w\times P_{\lambda })(\xi )+K_d\psi (d)(w\times P_{\lambda })(\xi -d)& (18)\end{array}$

with

$\begin{array}{cc}(w\times P_{\lambda })(\xi )={\displaystyle {\int }_{-\infty }^{\infty }\text{}}w(y)P_{\lambda }(\xi -y)\text{d}y,P_{\lambda }(\xi )=H(-\xi ){\text{e}}^{(\lambda +1)\xi /c}& (19)\end{array}$

The eigenvalues are now determined by imposing self-consistency at ξ = 0 and ξ = d. Setting ξ = 0 and ξ = d in Equation (18) leads to the vector equation

$\begin{array}{cc}\left[\begin{array}{ll}{K}_0(w\times P_{\lambda })(0)-1\hfill & K_d(w\times P_{\lambda })(-d)\hfill \\ \hfill & \hfill \\ K_0(w\times P_{\lambda })(d)\hfill & K_d(w\times P_{\lambda })(0)-1\hfill \end{array}\right]\left[\begin{array}{l}\psi (0)\hfill \\ \psi (d)\hfill \end{array}\right]=0& (20)\end{array}$

This has a non-trivial solution if and only if the determinant of the matrix is zero. The determinant expressed as a function of λ, ℰ(λ), is a complex analytic function known as the Evans function:

$\begin{array}{cc}\begin{array}{l}ℰ(\lambda )=\left[K_0(w\times P_{\lambda })(0)-1\right]\left[K_d(w\times P_{\lambda })(0)-1\right]\\ -K_0{K}_d(w\times P_{\lambda })(d)(w\times P_{\lambda })(-d).\end{array}& (21)\end{array}$

Thus, the zeros of the Evans function determine the discrete spectrum of the linear operator formed by linearizing the neural field equation about the pulse solution. Evans functions were originally introduced within the context of the stability of solitary pulses in diffusive Hodgkin–Huxley type equations describing action potential propagation in nerve axons (Evans, 1975). Since then the Evans function construction has been extended to a wide range of PDEs, see the review (Sandstede, 2002). It has also recently been applied to neural field equations (Zhang, 2003; Coombes and Owen, 2004; Rubin, 2004; Folias and Bressloff, 2005; Pinto et al., 2005; Sandstede, 2007) and more general non-local problems (Kapitula et al., 2004). An example plot of the real and imaginary parts of ℰ(λ) = 0 on the complex plane is shown in Figure 4. It can be seen that there is a zero eigenvalue and one negative real eigenvalue, indicating that the corresponding traveling pulse is linearly stable.

FIGURE 4

Figure 4. Graphs of the zero sets of the real (dark curves) and imaginary (light curves) parts of the Evans function determining the stability of a freely propagating pulse; intersection points (filled circles) indicate eigenvalues. The line Im λ = −1 indicates the essential spectrum. Same parameter values as Figure 3.

To find the essential spectrum, which is the union of the residual and continuous spectra, we will derive an explicit expression for the resolvent ℛ_λ. We start by writing an inhomogeneous equation of the form ${\mathcal{T}}_{\lambda }\psi (\xi )=h(\xi ),$, where h(ξ) represents a general function from the Banach space ℬ. This can be manipulated as before to give

$\begin{array}{cc}\begin{array}{l}\psi (\xi )=-(h\times P_{\lambda })(\xi )+K_0\psi (0)(w\times P_{\lambda })(\xi )\\ \text{ +}{K}_d\psi (d)(w\times P_{\lambda })(\xi -d)\end{array}& (22)\end{array}$

By evaluating at ξ = 0 and ξ = d as before, we arrive at the following vector equation, which differs from (Equation 20) only on the right-hand side:

$\begin{array}{l}\left[\begin{array}{ll}{K}_0(w\times P_{\lambda })(0)-1\hfill & K_d(w\times P_{\lambda })(-d)\hfill \\ \hfill & \hfill \\ K_0(w\times P_{\lambda })(d)\hfill & K_d(w\times P_{\lambda })(0)-1\hfill \end{array}\right]\left[\begin{array}{l}\psi (0)\hfill \\ \psi (d)\hfill \end{array}\right]\\ =\left[\begin{array}{l}(h\times P_{\lambda })(0)\hfill \\ (h\times P_{\lambda })(d)\hfill \end{array}\right].\end{array}$

Since we are looking for spectral values outside the discrete spectrum, the determinant of the matrix satisfies ℰ(λ)≠ 0. Therefore, multiplying both sides by the inverse matrix yields expressions for ψ(0) and ψ(d) in terms of h:

$\psi (0)=\frac{S_{0}h}{ℰ(\lambda )},\psi (d)=\frac{S_{d}h}{ℰ(\lambda )}$

where

$\begin{array}{l}{S}_{0}h=\left(K_d(w\times P_{\lambda })(0)-1\right)(h\times P_{\lambda })(0)\\ -K_d(w\times P_{\lambda })(-d)(h\times P_{\lambda })(d)\end{array}$

and

$\begin{array}{l}{S}_{d}h=-K_0(w\times P_{\lambda })(d)(h\times P_{\lambda })(0)\\ +\left(K_0(w\times P_{\lambda })(0)-1\right)(h\times P_{\lambda })(d).\end{array}$

Substituting into Equation (22) gives the following expression for the resolvent operator, of the form ℛ_λh = ϕ:

$\begin{array}{cc}\begin{array}{l}-(h\times P_{\lambda })(\xi )+\frac{K_0{S}_{0}h}{ℰ(\lambda )}(w\times P_{\lambda })(\xi )\\ +\frac{K_d{S}_{d}h}{ℰ(\lambda )}(w\times P_{\lambda })(\xi -d)=\psi (\xi )\end{array}& (23)\end{array}$

The resolvent is well-defined for all h in ℬ, so the residual spectrum of ℒ is empty. To find the continuous spectrum, we Fourier transform Equation (23):

$\begin{array}{cc}\begin{array}{l}-\widehat{h}(k){\widehat{P}}_{\lambda }(k)+\frac{K_0{S}_{0}h}{\epsilon (\lambda )}\widehat{w}(k){\widehat{P}}_{\lambda }(k)\\ +\frac{K_d{S}_{d}h}{\epsilon (\lambda )}\widehat{w}(k){\widehat{P}}_{\lambda }(k){\text{e}}^{-2\pi idk}=\widehat{\psi }(k)\end{array}& (24)\end{array}$

It follows that the resolvent operator is unbounded when ${\widehat{P}}_{\lambda }$ is unbounded. Equation (19) implies that

$\begin{array}{cc}{\widehat{P}}_{\lambda }(k)={\displaystyle {\int }_{-\infty }^{\infty }\text{}}{P}_{\lambda }(\xi ){\text{e}}^{2\pi ik\xi }\text{\hspace{0.17em}d}\xi =\frac{1}{\frac{\lambda +1}{c}+2\pi ik}.& (25)\end{array}$

Hence, ${\widehat{P}}_{\lambda }$ is unbounded for λ = −1−2π ikc so that the continuous spectrum of ℒ is a vertical line in the complex plane at Re(λ) = −1. Since Re(λ) < 0, the continuous spectrum will not make any pulse solution of our model unstable.

Stimulus-locked pulses

Now suppose that the neural field is driven by a moving external pulse stimulus of speed v so that Equation (1) becomes

$\begin{array}{cc}\begin{array}{c}\frac{\partial u(x,t)}{\partial t}=-u(x,t)+{\displaystyle {\int }_{-\infty }^{\infty }\text{}}w\left(x-y\right)H\left(u\left(y,t\right)-\kappa \right)\text{d}y\\ +I(x-vt).\end{array}& (26)\end{array}$

In order to study the existence of stimulus-locked pulses, we will define a “stimulus coordinate” ξ = x − vt and look for pulse solutions that move at the same speed as the stimulus, that is, u(x, t) = U(ξ) with

$\begin{array}{cc}-v\frac{\partial U(\xi )}{\partial \xi }=-U(\xi )+{\displaystyle {\int }_{-\infty }^{\infty }\text{}}w(\xi -y)H(U(y)-\kappa )\text{d}y+I(\xi ).& (27)\end{array}$

For concreteness, the stimulus will be represented by a rectangular wave of amplitude I₀ and width d, defined formally as

$I(\xi )=\{\begin{array}{ll}{I}_0,\hfill & \text{if }0\le \xi \le d\hfill \\ 0,\hfill & \text{if }\xi <0\text{ or }\xi >d.\hfill \end{array}$

Since translation invariance no longer holds, it is necessary to determine both threshold crossing points, which we denote by ξ = d₁ and ξ = d₂. Proceeding in a similar fashion to the case of freely propagating pulses, we find that

$U(\xi )=\frac{{\text{e}}^{\xi /v}}{v}{\displaystyle {\int }_{\xi }^{z_0}\text{}}{\text{e}}^{-y/v}W(y)\text{\hspace{0.17em}d}y+\frac{{\text{e}}^{\xi /v}}{v}{\displaystyle {\int }_{\xi }^{z_0}\text{}}{\text{e}}^{-y/v}I(y)\text{\hspace{0.17em}d}y,$

where z₀ = ∞ if v > 0, z₀ = −∞ if v < 0, and

$W(\xi )\equiv {\displaystyle {\int }_{\xi -d_2}^{\xi -d_1}\text{}}w(x)\text{\hspace{0.17em}d}x.$

The latter can be expressed in the piecewise form

$\begin{array}{cc}W(\xi )=\{\begin{array}{ll}{W}_3(\xi )\equiv {\displaystyle {\int }_{\xi -d_2}^{\xi -d_1}\text{}}{w}_2(x)\text{\hspace{0.17em}d}x,\hfill & \text{if }\xi \le x_0+d_1\hfill \\ W_2(\xi )\equiv {\displaystyle {\int }_{\xi -d_2}^{x_0}\text{}}{w}_2(x)\text{\hspace{0.17em}d}x\hfill & \hfill \\ +{\displaystyle {\int }_{x_0}^{\xi -d_1}\text{}}{w}_1(x)\text{\hspace{0.17em}d}x,\hfill & \text{if }{x}_0+d_1\le \xi \le x_0+d_2\hfill \\ W_1(\xi )={\displaystyle {\int }_{\xi -d_2}^{\xi -d_1}\text{}}{w}_1(x)\text{\hspace{0.17em}d}x,\hfill & \text{if }\xi \ge x_0+d_2\hfill \end{array}& (28)\end{array}$

where w₁ and w₂ are defined as in Equation (10). After evaluating the integrals along similar lines to section “Neural Field Model of Direction Selectivity,” we obtain the following expressions for the pulse solution, defined independently for positive and negative stimulus directions:

v > 0:

$U(\xi )=\{\begin{array}{ll}{U}_3(\xi ),\hfill & \text{if }\xi \le x_0+d_1\hfill \\ U_2(\xi ),\hfill & \text{if }{x}_0+d_1\le \xi \le x_0+d_2,\hfill \\ U_1(\xi ),\hfill & \text{if }\xi \ge x_0+d_2,\hfill \end{array}$

with

$\begin{array}{l}{U}_3(\xi )=\frac{1}{v}{\text{e}}^{\xi /v}\left(M_3(\xi )+M_1(x_0+d_2)+M_2(x_0+d_1)\right)+Z(\xi )\\ U_2(\xi )=\frac{1}{v}{\text{e}}^{\xi /v}\left(M_2(\xi )+M_1(x_0+d_2)\right)+Z(\xi )\\ U_1(\xi )=\frac{1}{v}{\text{e}}^{\xi /v}{M}_1(\xi )+Z(\xi ),\\ M_n(\xi )={\displaystyle {\int }_{\xi }^{{\xi }_n}{W}_n}({\xi }^{\prime }){\text{e}}^{-{\xi }^{\prime }/v}d{\xi }^{\prime }\end{array}$

for ξ₁ = ∞, ξ₂ = x₀ + d₂, ξ₃ = x₀ + d₁, and

$Z(\xi )=\{\begin{array}{ll}\left({\text{e}}^{\xi /v}-{\text{e}}^{(\xi -d)/v}\right)I_0,\hfill & \text{if }\xi <0\hfill \\ \left(1-{\text{e}}^{(\xi -d)/v}\right)I_0,\hfill & \text{if }0\le \xi \le d.\hfill \\ 0,\hfill & \text{if }\xi >d\hfill \end{array}$

v < 0:

$U(\xi )=\{\begin{array}{ll}{U}_3(\xi ),\hfill & \text{if }\xi \le x_0+d_1\hfill \\ U_2(\xi ),\hfill & \text{if }{x}_0+d_1\le \xi \le x_0+d_2,\hfill \\ U_1(\xi ),\hfill & \text{if }\xi \ge x_0+d_2\hfill \end{array}$

with

$\begin{array}{l}{U}_3(\xi )=-\frac{1}{v}{\text{e}}^{\xi /v}{N}_3(\xi )+Z(\xi )\\ U_2(\xi )=-\frac{1}{v}{\text{e}}^{\xi /v}\left(N_2(\xi )+N_3(x_0+d_1)\right)+Z(\xi )\\ U_1(\xi )=-\frac{1}{v}{\text{e}}^{\xi /v}\left(N_1(\xi )+N_2(x_0+d_1)+N_3(x_0+d_1)\right)+Z(\xi )\\ N_n(\xi )={\displaystyle {\int }_{{\xi }_n}^{\xi }{W}_n}({\xi }^{\prime }){\text{e}}^{-{\xi }^{\prime }/v}d{\xi }^{\prime }\end{array}$

for ξ₃ = −∞, ξ₂ = x₀ + d₁, ξ₁ = x₀ + d₂, and

$Z(\xi )=\{\begin{array}{ll}0,\hfill & \text{if }\xi <0\hfill \\ \left(1-{\text{e}}^{\xi /v}\right)I_0,\hfill & \text{if }0\le \xi \le d.\hfill \\ \left({\text{e}}^{(\xi -d)/v}-{\text{e}}^{\xi /v}\right)I_0,\hfill & \text{if }\xi >d\hfill \end{array}$

The threshold crossing points (d₁ and d₂) are determined in the same way the pulse speed and width were determined in the no-stimulus case, which is by numerically solving a system of two transcendental equations. The first equation is given by U₃(d₁) = κ. The second equation is U₃(d₂) = κ if d₂ < x₀ + d₁, else it is given by U₂(d₂) = κ. Figure 5 shows a plot of d₁ (black curves) and d₂ (gray curves) vs. the threshold κ. It can be seen that for a certain range of thresholds there exists more than one stable/unstable pair of pulses. Figure 6 shows the linear stability and the number of solutions for different combinations of stimulus speed (v) and strength (I₀). The offset x₀ = 3 and the corresponding spontaneous wave speed is c = 4. (Note that for smaller offsets x₀ and thus smaller wave speeds c, one finds stimulus-locked waves for negative values of v). The stability of solutions in the presence of a stimulus is determined in much the same way as without a stimulus. We again define u(x, t) = U(ξ) + ϕ(ξ,t) and look at the behavior of the perturbations described by ϕ(ξ, t). Substituting into Equation (26), the stimulus term drops out when we perform the linearization, so that

$\begin{array}{c}\frac{\partial \phi (\xi ,t)}{\partial t}=v\frac{\partial \phi (\xi ,t)}{\partial \xi }-\phi (\xi ,t)+{\displaystyle {\int }_{-\infty }^{\infty }\text{}}w(\xi -y)\\ \times F^{\prime }(U(y))\phi (y,t)\text{\hspace{0.17em}d}y.\end{array}$

FIGURE 5

Figure 5. Plot of leading threshold crossing position d₂ (gray curves) and trailing threshold crossing position d₁ (black curves) of stimulus-locked pulses as a function of threshold κ. There exists at most one stable pulse (indicated by arrows) and up to three unstable pulses. Weight parameters as in Figure 2 and κ = 4. Stimulus parameters are d = 5, I₀ = 5, and v = 5.

FIGURE 6

Figure 6. Stability diagram for stimulus-locked pulses in (v, I₀)-parameter space. Weight parameters as in Figure 2, κ = 4, and d = 5. Emerging from the stable pulse solution when I₀ = 0 is a tongue consisting of a stable/unstable pair of pulses (light gray). Similarly, emerging from the unstable solution when I₀ = 0 is a tongue consisting of two unstable pulses (medium gray). As I₀ increases within a tongue an unstable pulse can disappear due to the development of multiple super-threshold regions (indicated by solid curves). All solutions coexist when tongues overlap (dark gray).

Setting ϕ(ξ, t) = e^λtϕ ultimately yields the spectral problem

$\begin{array}{cc}\begin{array}{l}\lambda \phi (\xi )\equiv ℒ\phi (\xi )\\ =v\frac{\partial \phi (\xi ,t)}{\partial \xi }-\phi (\xi ,t)\\ +\frac{|v|\phi (d_1)}{|\kappa -W(d_1)-I(d_1)|}w(\xi -d_1)\\ +\frac{|v|\phi (d_2)}{|\kappa -W(d_2)-I(d_2)|}w(\xi -d_2).\end{array}& (29)\end{array}$

The corresponding Evans function is now

$\begin{array}{cc}\begin{array}{l}ℰ(\lambda )=\left[K_1(w\times P_{\lambda })(0)-1\right]\left[K_2(w\times P_{\lambda })(0)-1\right]\\ -K_1{K}_2(w\times P_{\lambda })(d_2-d_1)(w\times P_{\lambda })(d_1-d_2),\end{array}& (30)\end{array}$

where

$K_n=\frac{\text{sgn}(v)}{|\kappa -W(d_n)-I(d_n)|},n=1,2,$

and

$P_{\lambda }(\xi )=H(-\text{sgn}(v)\xi ){\text{e}}^{(\lambda +1)\xi /v}.$

It is easy to establish as before that the residual spectrum is empty and that the continuous spectrum consists of a vertical line in the complex plane at Re(λ) = −1. So the stability is again determined only by the discrete spectrum, which consists of the zeros of the Evans function. Figure 7A shows an example of a numerical simulation of a stable stimulus-locked pulse solution of Equation (26) with the analytical pulse solution U(x) as an initial condition. Figure 7B shows the same simulation except with the zero initial condition u(x,0) = 0. It can be seen from Figure 7C that both initial conditions converge to the same pulse profile.

FIGURE 7

Figure 7. Space-time plots of a stimulus-locked pulse for (A) u(x, 0) = U(x) and (B) u(x, 0) = 0. Both initial conditions converge to the same wave profile in the large t limit as indicated in (C) for t = 50. Weight parameters as in Figure 2 and κ = 4. Stimulus parameters are d = 5, I₀ = 8, and v = 3.

Stochastic Neural Field

Several recent studies have considered stochastic versions of neural field equations that are based on a corresponding Langevin equation formulation (Brackley and Turner, 2007; Hutt et al., 2008; Faugeras et al., 2009; Bressloff and Webber, 2012b). Motivated by these examples, we consider the following Langevin equation (or stochastic PDE) for the stochastic activity variable U(x, t), which is a rewriting of Equation (1) with h(x, t) given by Equation (5) for I₀ = 0 and ϵ > 0:

$\begin{array}{cc}\begin{array}{c}\text{d}U(x,t)=\left[-U(x,t)+{\displaystyle {\int }_{-\infty }^{\infty }\text{}}w(x-y)F(U(y,t))\text{d}y\right]\text{d}t\\ +{ϵ}^{1/2}g(U(x,t))\text{d}W(x,t),\end{array}& (31)\end{array}$

where dW(x, t) is an independent Wiener process with zero mean and correlation given by

$\begin{array}{cc}\langle \text{d}W(x,t)\text{d}W(x^{\prime },t^{\prime })\rangle =2C(|x-x^{\prime }|/\lambda )\delta (t-t^{\prime })dtd{t}^{\prime }.& (32)\end{array}$

Here λ is the spatial correlation length of the noise such that C(x/λ)→ δ(x) in the limit λ → 0, and ϵ determines the strength of the noise, which is assumed to be weak. For the sake of generality, we take the noise to be multiplicative rather than additive; however, the main results of the paper hold for both. Following standard formulations of Langevin equations (Gardiner, 2009), the multiplicative noise term is taken to be of Stratonovich form in the case of extrinsic noise. Note, however, that an alternative formulation of stochastic neural field theory has been developed in terms of a neural master equation (Buice and Cowan, 2007; Bressloff, 2009, 2010; Buice et al., 2010), in which the underlying deterministic equations are recovered in the thermodynamic limit N→ ∞, where N is a measure of the system size of each local population. In the case of large but finite N, a Kramers-Moyal expansion of the master equation yields a Langevin neural field equation with multiplicative noise of the Ito form Bressloff (2009, 2010). Multiplicative noise in the Stratonovich sense causes a shift in the speed and width of the pulse. This happens because 〈 g(U)dW〉≠ 0, even though 〈 dW〉 = 0. We can use Novikov's theorem (Novikov, 1965) to calculate the former average:

${ϵ}^{1/2}\langle g(U)\text{d}W\rangle =ϵC(0)\langle g^{\prime }(U)g(U)\rangle \text{d}t.$

The average can also be calculated by Fourier transforming Equation (31) and taking averages using the corresponding Fokker–Planck equation (Armero et al., 1998; Bressloff and Webber, 2012b). In the limit that λ approaches 0, we set C(0) → 1/Δx, where Δx is a lattice cut-off that can be identified with the spatial discretization step size in numerical simulations (Bressloff and Webber, 2012b). Following Ref. Armero et al. (1998), we rewrite Equation (31) so the fluctuating term has zero mean:

$\begin{array}{cc}\begin{array}{c}\text{d}U(x,t)=\left[h(U(x,t))+{\displaystyle {\int }_{-\infty }^{\infty }\text{}}w(x-y)F(U(y,t))\text{d}y\right]\text{d}t\\ +{ϵ}^{1/2}\text{d}R(U,x,t),\end{array}& (33)\end{array}$

where

$\begin{array}{cc}h(U)=-U+ϵC(0)g^{\prime }(U)g(U)& (34)\end{array}$

and

$\begin{array}{cc}\text{d}R(U,x,t)=g(U)\text{d}W(x,t)-{ϵ}^{1/2}C(0)g^{\prime }(U)g(U)\text{d}t.& (35)\end{array}$

The stochastic process R has zero mean and correlation

$\begin{array}{cc}\begin{array}{l}\langle dR(U,x,t)dR(U,x^{\prime },t^{\prime })\rangle \\ =\langle g(U(x,t))\text{d}W(x,t)g(U(x^{\prime },t^{\prime }))\text{d}W(x^{\prime },t^{\prime })\rangle +\mathcal{O}({ϵ}^{1/2}).\end{array}& (36)\end{array}$

Separation of time-scales

The effects of additive or multiplicative extrinsic noise on traveling waves can be analyzed using multiple time-scale methods originally developed for reaction-diffusion equations (Schimansky-Geier et al., 1983; de Pasquale et al., 1992; Armero et al., 1998; Sagues et al., 2007), which were recently extended to neural field equations in Ref. Bressloff and Webber (2012b). The main idea is to assume that the fluctuating term generates two distinct phenomena that occur on different time-scales: a diffusive-like displacement of the traveling wave from its uniformly translating position at long time-scales, and fluctuations in the wave profile around its instantaneous position at short time-scales. It is important to point out that, in contrast to traveling front solutions of scalar neural field equations (Bressloff and Webber, 2012b), we are now considering traveling pulse solutions. Thus in addition to the center-of-mass of the traveling pulse wave, which moves with speed c in the absence of noise, there is an additional degree of freedom corresponding to the “width” of the pulse. (In the case of a Heaviside rate function, the width Δ is determined by the threshold crossing points). For simplicity, we assume that the width of the wave is only weakly affected by the noise; this is consistent with what is found numerically. We now express the solution U of Equation (33) as a combination of a fixed wave profile U₀ that is displaced by an amount Δ(t) from its uniformly translating position ξ = x−c_ϵt, where c_ϵ is a noise-dependent speed, and a time-dependent fluctuation Φ in the wave shape about its instantaneous position:

$\begin{array}{cc}U(x,t)=U_0(\xi -\Delta (t))+{ϵ}^{1/2}\Phi (\xi -\Delta (t),t).& (37)\end{array}$

The wave profile U₀ and associated wave speed/width c_ϵ, Δ_ϵ are obtained by solving the modified deterministic equation

$\begin{array}{cc}-c_{ϵ}\frac{\text{d}{U}_0(\xi )}{\text{d}\xi }-h(U_0(\xi ))={\displaystyle {\int }_{\infty }^{\infty }\text{}}w(\xi -{\xi }^{\prime })F(U_0({\xi }^{\prime }))\text{d}{\xi }^{\prime }.& (38)\end{array}$

The results depend on ϵ due to the ϵ-dependence of h. Equation (38) is chosen so that that to leading order, the stochastic variable Δ(t) undergoes unbiased Brownian motion with a diffusion coefficient D(ϵ) = $\mathcal{O}$(ϵ) (see below). The next step is to substitute the decomposition Equation (37) into (33) and expand to first order in $\mathcal{O}$(ε^1/2):

$\begin{array}{l}-[c_{\epsilon }+\dot{\Delta }]U_{0^{\prime }}({\xi }_{\Delta })dt+{ϵ}^{1/2}\left[d\Phi ({\xi }_{\Delta },t)-[c_{ϵ}+\dot{\Delta }]{\Phi }^{\prime }({\xi }_{\Delta },t)dt\right]\\ =h(U_0({\xi }_{\Delta }))dt+{ϵ}^{1/2}{h}^{\prime }(U_0({\xi }_{\Delta }))\Phi ({\xi }_{\Delta },t)dt\\ +{\displaystyle {\int }_{-\infty }^{\infty }w}(\xi -{\xi }^{\prime })F(U_0({\xi }_{{\Delta }^{\prime }}))d{\xi }^{\prime }dt\\ +{ϵ}^{1/2}{\displaystyle {\int }_{-\infty }^{\infty }w}(\xi -{\xi }^{\prime })F^{\prime }(U_0({\xi }_{{\Delta }^{\prime }}))\Phi ({\xi }_{{\Delta }^{\prime }},t)d{\xi }^{\prime }dt\\ +{ϵ}^{1/2}dR(U_0({\xi }_{\Delta }),\xi ,t)+\mathcal{O}(ϵ).\end{array}$

where we have set ξ_Δ = ξ − Δ(t) and ξ_Δ' = ξ' − Δ(t). We now use Equation (38) for U₀, after shifting ξ→ ξ − Δ(t), to eliminate terms and then divide through by $\sqrt{\in }$. This gives the inhomogeneous equation to $\mathcal{O}$(ϵ^1/2)

$\begin{array}{cc}\begin{array}{c}d\Phi ({\xi }_{\Delta },t)-\widehat{L}\Phi ({\xi }_{\Delta },t)\text{d}t={ϵ}^{-\frac{1}{2}}{U^{\prime }}_0({\xi }_{\Delta })d\Delta (t)\\ +dR(U_0({\xi }_{\Delta }),\xi ,t)\end{array}& (39)\end{array}$

where the non-self-adjoint linear operator

$\begin{array}{cc}\begin{array}{l}\widehat{L}A(\xi )\equiv c_{ϵ}\frac{\text{d}A(\xi )}{\text{d}\xi }+h^{\prime }(U_0(\xi ))A(\xi )\\ +{\displaystyle {\int }_{-\infty }^{\infty }\text{}}w(\xi -{\xi }^{\prime })F^{\prime }(U_0({\xi }^{\prime }))A({\xi }^{\prime })\text{d}{\xi }^{\prime }\end{array}& (40)\end{array}$

is defined for all functions A(ξ) in ℒ₂(ℝ). Note that for all terms in Equation (40) to be of the same order we have taken Δ(t) = $\mathcal{O}$(ε^1/2). It then follows that U₀(ξ − Δ(t)) = U₀(ξ) + $\mathcal{O}$(ε^1/2) etc., and Equation (39) reduces to

$\begin{array}{cc}d\Phi (\xi ,t)-\widehat{L}\Phi (\xi ,t)\text{d}t={ϵ}^{-\frac{1}{2}}{U^{\prime }}_0(\xi )d\Delta (t)+d{R}_u(U_0(\xi ),\xi ,t)& (41)\end{array}$

If U₀(ξ) were a traveling front solution of a neural field model with a symmetric, excitatory weight distribution w, then it could be proven that the operator $\widehat{L}$ has a 1D null space spanned by U₀'(ξ) (Ermentrout and McLeod, 1993)). We will assume that such a result carries over to traveling pulse solutions of a neural field with w given by an asymmetric Mexican hat function; the fact that U₀'(ξ) belongs to the null space follows immediately from differentiating Equation (38) with respect to ξ. We then have the solvability condition for the existence of a non-trivial bounded solution of Equation (41), namely, that the inhomogeneous part is orthogonal to all elements of the null space of the adjoint operator ${\widehat{L}}^{\ast }$. The latter is defined with respect to the inner product

${\displaystyle {\int }_{-\infty }^{\infty }\text{}}B(\xi )\widehat{L}A(\xi )\text{\hspace{0.17em}d}\xi ={\displaystyle {\int }_{-\infty }^{\infty }\text{}}\left[{\widehat{L}}^{\ast }B(\xi )\right]A(\xi )\text{\hspace{0.17em}d}\xi .$

Integrating by parts and using (Equation 14) leads to

$\begin{array}{cc}\begin{array}{l}{\widehat{L}}^{\ast }B(\xi )=-c_{ϵ}\frac{\text{d}B(\xi )}{\text{d}\xi }+h^{\prime }(U_0(\xi ))B(\xi )\\ +F^{\prime }(U_0(\xi )){\displaystyle {\int }_{-\infty }^{\infty }w}({\xi }^{\prime }-\xi )B({\xi }^{\prime })\text{d}{\xi }^{\prime }.\end{array}& (42)\end{array}$

We will assume that the null space of the adjoint operator ${\widehat{L}}^{\ast }$ is also one-dimensional and is spanned by some yet to be determined function $\mathcal{V}(\xi )$. (In the case of a Heaviside firing function, we will determine the null space explicitly). Hence, we can write the solvability condition as

${\displaystyle {\int }_{-\infty }^{\infty }\mathcal{V}}(\xi )\left[U_{0^{\prime }}(\xi )\text{d}\Delta (t)+{ϵ}^{1/2}\text{d}R(U_0(\xi ),\xi ,t)\right]\text{\hspace{0.17em}d}\xi =0.$

which leads directly to the stochastic differential equation

$\text{d}\Delta (t)=-{ϵ}^{1/2}\frac{{\displaystyle {\int }_{-\infty }^{\infty }\text{}}\mathcal{V}(\xi )\text{d}R(U_0,\xi ,t)\text{\hspace{0.17em}d}\xi }{{\displaystyle {\int }_{-\infty }^{\infty }\text{}}\mathcal{V}(\xi )U_{0^{\prime }}(\xi )\text{d}\xi }.$

Using the lowest order approximations dR(U₀,ξ, t) = g(U₀(ξ))dW(ξ, t), we deduce that [for Δ(0) = 0]

$\begin{array}{cc}\langle \Delta (t)\rangle =0,\langle \Delta {(t)}^2\rangle =2D(\epsilon )t,& (43)\end{array}$

where D(ϵ) is the effective diffusivity

$\begin{array}{cc}D(ϵ)=ϵ\frac{{\displaystyle {\int }_{-\infty }^{\infty }\text{}}{\mathcal{V}}^2(\xi )g{(U_0(\xi ))}^2\text{\hspace{0.17em}d}\xi }{{\left[{\displaystyle {\int }_{-\infty }^{\infty }\text{}}\mathcal{V}(\xi )U_{0^{\prime }}(\xi )\text{\hspace{0.17em}d}\xi \right]}^2}.& (44)\end{array}$

Explicit results for heaviside rate function

In order to illustrate the above analysis, we consider a particular example where the mean speed c_ϵ and diffusion coefficient D(ϵ) can be calculated explicitly. That is, set g(U) = g₀U for the multiplicative noise term and take F(U) = H(u − κ). (The constant g₀ has units of $\sqrt{\text{length}/\text{time}}$). Note that the choice for g(U) can be interpreted physiologically in terms of an effective modification in the membrane time constant of neurons due to stochastic background synaptic activity (Bernander et al., 1991; Rapp et al., 1992; Bressloff, 1994). The deterministic Equation (38) for U₀ then reduces to

$\begin{array}{cc}-\frac{\text{d}{U}_0(\xi )}{\text{d}\xi }+\Gamma (ϵ)U_0(\xi )=\frac{1}{c_{ϵ}}{\displaystyle {\int }_{\infty }^{\infty }\text{}}w(\xi -{\xi }^{\prime })H(U_0({\xi }^{\prime })-\kappa )\text{d}{\xi }^{\prime },& (45)\end{array}$

where Γ(ϵ) = (1 − ϵ C(0)g²₀)/c_ϵ. Hence,

$\begin{array}{cc}U(\xi )=\frac{e^{\Gamma \xi }}{c_{ϵ}}{\displaystyle {\int }_{\xi }^{\infty }W}({\xi }^{\prime }){\text{e}}^{-\Gamma {\xi }^{\prime }}\text{\hspace{0.17em}d}{\xi }^{\prime }.& (46)\end{array}$

The deterministic pulse profile can be evaluated along identical lines to section “Neural Field Model of Direction Selectivity.” In order to calculate the diffusion coefficient, it is first necessary to determine the null vector $\mathcal{V}(\xi )$ of the adjoint linear operator ${\widehat{L}}^{\ast }$. Substituting F(U) = H(U − κ) and g(U) = g₀U into Equation (42) shows that

$\begin{array}{cc}\begin{array}{l}\frac{d\mathcal{V}(\xi )}{d\xi }+\Gamma (ϵ)\mathcal{V}(\xi )=\frac{\delta (\xi )}{c|{U^{\prime }}_0(0)|}{\displaystyle {\int }_{-\infty }^{\infty }w}(z)\mathcal{V}(z)dz\\ +\frac{\delta (\xi -\Delta )}{c|{U^{\prime }}_0(\Delta )|}{\displaystyle {\int }_{-\infty }^{\infty }w}(z-\Delta )\mathcal{V}(z)dz.\end{array}& (47)\end{array}$

Proceeding along similar lines to Bressloff (2001) and Kilpatrick et al. (2008), we make the ansatz that

$\begin{array}{cc}\mathcal{V}(\xi )=AH(\xi ){\text{e}}^{-\Gamma \xi }+BH(\xi -\Delta ){\text{e}}^{-\Gamma (\xi -\Delta )}.& (48)\end{array}$

Substituting into Equation (47) shows that

$A=\frac{1}{|U_{0^{\prime }}(0)|}[Ab(0)+Bb(\Delta )],B=\frac{1}{|U_{0^{\prime }}(\Delta )|}[Ab(-\Delta )+Bb(0)]$

where

$\begin{array}{cc}b(z)\equiv \frac{1}{c}{\displaystyle {\int }_z^{\infty }{\text{e}}^{-\Gamma ({\xi }^{\prime }-z)}}w({\xi }^{\prime })\text{\hspace{0.17em}d}{\xi }^{\prime }.& (49)\end{array}$

Differentiating Equation (46) shows that U'(ξ) = b(ξ) − b(ξ − Δ), so that we obtain the vector equation

$\left[\begin{array}{cc}\frac{b(0)}{b(0)-b(-\Delta )}-1& \frac{b(\Delta )}{b(0)-b(-\Delta )}\\ \frac{b(-\Delta )}{b(0)-b(\Delta )}& \frac{b(0)}{b(0)-b(\Delta )}-1\end{array}\right]\left[\begin{array}{l}A\hfill \\ B\hfill \end{array}\right]=0$

The matrix has rank 1, confirming that the linear operator ${\widehat{L}}^{\ast }$ has a 1D null-space. The latter is spanned by the function

$\begin{array}{cc}\mathcal{V}(\xi )=b(\Delta )H(\xi ){\text{e}}^{-\Gamma \xi }-b(-\Delta )H(\xi -\Delta ){\text{e}}^{-\Gamma (\xi -\Delta )}.& (50)\end{array}$

In Figure 8 we show the temporal evolution of a freely propagating stochastic traveling pulse, which is obtained by numerically solving the Langevin Equation (31) for F(U) = H(U − κ), g(U) = U and the asymmetric difference-of-exponentials (Equation 4). Note that the location of the stochastic wave appears to coincide with the underlying mean solution. However, over longer time-scales the wandering of the pulse about its mean position would be seen. In Figure 9 we plot the mean position $\overline{X}(t)$ and variance σ²_X(t) of the leading and trailing edges of the pulse as a function of t. It can be seen that they all vary linearly with t, consistent with the assumption that there is a diffusive-like displacement of the center-of-mass of the pulse from its uniformly translating position at long time-scales. The slopes of these curves then determine the effective wave speed and diffusion coefficient according to $\overline{X}(t)$ ~ c_εt and σ²_X(t)~ 2 D(ε) t. Both the leading and trailing edges exhibit the same speeds and diffusivities (after a transient phase). The transients are caused by fluctuations in the mean width of the pulse which can be neglected for large t, where the difference in the size of fluctuations of the leading and trailing edges can be neglected.

FIGURE 8

Figure 8. Numerical simulation of freely propagating pulse solution of the stochastic neural field Equation (31) for a Heaviside rate function F(U) = H(U − κ) with κ = 4, and weight function (Equation 4) with a_e = 5, a_i = 1, σ_e = 0.42, σ_i = 0.1, and x₀ = 3. The multiplicative noise is taken to be g(U) = U, the noise strength is ϵ = 0.005, and C(0) = 10. The wave profile is shown at successive times (A) t = 0 (B) t = 12 (C) t = 18, and (D) t = 24, with the initial profile at t = 0 given byequation U₀. In numerical simulations we take the discrete space and time steps Δx = 0.1, Δt = 0.01. The deterministic part U₀ of the stochastic wave is shown by the black curves.

FIGURE 9

Figure 9. (A) Plot of mean position $\overline{X}(t)$ of leading (blue) and trailing (black) edges of pulse as a function of time t averaged over N = 4096 trials. (B) Corresponding plots of the variance σ²_X(t). Same parameter values as Figure 8.

In order to find the mean location of the leading or trailing edge of the pulse as a function of time, we numerically carry out a large number of level set position measurements. That is, we determine the positions X_a(t) such that U(X_a(t), t) = a, for various level set values a and then define the mean location to be $\overline{X}(t)=\mathbb{E}[X_a(t)]$, where the expectation is first taken with respect to the sampled values a and then averaged over N trials. The corresponding variance is given by ${\sigma }_X^2(t)=\mathbb{E}[{(X_a(t)-\overline{X}(t))}^2]$. In order to compare the numerical results with our theoretical analysis, we assume that X_a(t) = Δ(t) + c_ϵt +x_a(0) for each a on either the leading or trailing edge. It then follows that $\overline{X}(t)$ = c_ϵt + x_a(0) and σ²_X(t) = 〈 Δ(t)²〉. In Figure 10 we plot the numerically estimated diffusion coefficient for various values of the threshold κ and compare these to the corresponding theoretical curves obtained using the above analysis. It can be seen that there is excellent agreement with our theoretical predictions. Finally, note that we can also use the level set data to estimate fluctuations in the width of the pulse. Suppose that X_d(t) and Y_d(t) denote the threshold crossing points of the leading and trailing edges of the pulse at time t. Then the stochastic width of the pulse can be defined according to D(t) = X_d(t) − Y_d(t). We find that after a transient phase, 〈 D(t)²〉 − 〈 D(t)〉 ² « σ²_X(t).

FIGURE 10

Figure 10. Plot of diffusion coefficient D(ε) as a function of threshold κ. Numerical results (“+” for leading edge,“X” for trailing edge) are obtained by averaging over N = 4096 trials starting from the initial condition given by U₀. Corresponding theoretical predictions (solid curves) for D(ε) are based on Equation (44). Other parameters as in Figure 8.

Stimulus-locked pulses

We now add a stimulus term I to the stochastic neural field Equation (33), that is

$\begin{array}{cc}\begin{array}{c}\text{d}U(x,t)=\left[h(U(x,t))+{\displaystyle {\int }_{-\infty }^{\infty }\text{}}w(x-y)F(U(y,t))\text{\hspace{0.17em}d}y\right]\text{d}t\\ +I(x-vt)dt+{ϵ}^{1/2}\text{d}R(U,x,t),\end{array}& (51)\end{array}$

where the stimulus is again a rectangular wave of amplitude I₀ and width d, moving with speed v. Here h and dR are defined by Equations (34) and (35). The stochastic activity variable is now decomposed according to Equation (37), with ξ = x − vt, and the modified deterministic equation

$\begin{array}{cc}-v\frac{\text{d}{U}_0}{\text{d}\xi }-h(U_0(\xi ))-I(\xi )={\displaystyle {\int }_{\infty }^{\infty }\text{}}w(\xi -{\xi }^{\prime })F(U_0({\xi }^{\prime }))\text{\hspace{0.17em}d}{\xi }^{\prime }.& (52)\end{array}$

Through a similar process as in the previous section, we expand to $\mathcal{O}$(ϵ^1/2) to obtain the inhomogeneous equation

$\begin{array}{cc}\begin{array}{c}d\Phi (\xi ,t)-\widehat{L}\Phi (\xi ,t)\text{d}t=-{ϵ}^{-1/2}\left[U_{0^{\prime }}(\xi )-I^{\prime }(\xi )\right]\text{d}\Delta (t)\\ -\text{d}R(U_0(\xi ),\xi ,t)+\mathcal{O}({ϵ}^{1/2}),\end{array}& (53)\end{array}$

where $\widehat{L}$ is defined as in Equation (40) but with c_ϵ → v. The solvability condition is now

$\begin{array}{cc}\begin{array}{l}{\displaystyle {\int }_{-\infty }^{\infty }\mathcal{V}}(\xi )[U_{0^{\prime }}(\xi )\text{d}\Delta (t)+I^{\prime }(\xi )d\Delta (t)\\ +{ϵ}^{1/2}\text{d}R(U_0(\xi ),\xi ,t)]\text{d}\xi =0.\end{array}& (54)\end{array}$

This can be manipulated to give, to leading order, the Ornstein–Uhlenbeck equation (Gardiner, 2009):

$\begin{array}{cc}\text{d}\Delta (t)+A\Delta (t)\text{d}t=\text{d}\widehat{W}(t),& (55)\end{array}$

where

$A=\frac{{\displaystyle {\int }_{-\infty }^{\infty }\text{}}\mathcal{V}(\xi )I^{\prime }(\xi )\text{d}\xi }{{\displaystyle {\int }_{-\infty }^{\infty }\text{}}\mathcal{V}(\xi )U_{0^{\prime }}(\xi )\text{d}\xi },$

and

$\text{d}\widehat{W}(t)=-{ϵ}^{1/2}{g}_0\frac{{\displaystyle {\int }_{-\infty }^{\infty }\text{}}\mathcal{V}(\xi )U_0(\xi )\text{d}W(\xi ,t)\text{\hspace{0.17em}d}\xi }{{\displaystyle {\int }_{-\infty }^{\infty }\text{}}\mathcal{V}(\xi )U_{0^{\prime }}(\xi )\text{d}\xi }.$

Solving the stochastic differential equation in Equation (55) and taking averages shows that 〈Δ(t)〉 = Δ(0)e ^{− At} and

$\begin{array}{cc}\langle \Delta {(t)}^2\rangle -{\langle \Delta (t)\rangle }^2\approx \frac{D(ϵ)}{A}\left[1-{\text{e}}^{-2At}\right],& (56)\end{array}$

where D(ϵ) is given by Equation (44) except for a modified null vector $\mathcal{V}(\xi )$. Thus the variance of Δ(t) approaches D(ϵ)/A in the large t limit.

As in the case of freely propagating pulses, we can explicitly solve for $\mathcal{V}(\in )$ and thus calculate the diffusion coefficient D(ϵ) when F(U) = H(U − κ) and g(U) = U. Since the steps are similar to the previous case, we simply present our results here. In Figure 11 we show the temporal evolution of a single stimulus-locked front, which is obtained by numerically solving the Langevin Equation (51) for F(U) = H(U − κ), g(U) = U and the weight distribution (Equation 4). The external input is taken to be a square pulse of amplitude I₀ = 5, width d = 5, and speed v = 5. Next we determine the mean $\overline{X}(t)$ and variance σ²_X(t) of the position of the leading and trailing edges by averaging over level sets along identical lines to the freely-propagating case. The results are shown in Figure 12. It can be seen that, as predicted by the analysis, $\overline{X}(t)$ varies linearly with t with a slope equal to the stimulus speed v = 5. Moreover, the variance σ²_X(t) approaches a constant value as t → ∞ for both the trailing and leading edges. Thus, we find that stimulus-locked pulses are much more robust to noise than freely propagating pulses, since the variance of the mean position of the leading and trailing edges saturate as t→ ∞. Consequently, stimulus locking persists in the presence of noise over most of the parameter range for which stimulus locking is predicted to occur. However, the trailing edge has an asymptotic variance that is at least an order of magnitude larger than the leading edge, which implies that fluctuations in the width of the pulse can no longer be neglected.

FIGURE 11

Figure 11. Numerical simulation showing the propagation of a stimulus-locked pulse solution of the stochastic neural field Equation (51). The external input is taken to be a square pulse with amplitude I₀ = 5, width d = 5, and speed v = 5. All other parameters are as in Figure 8. The wave profile is shown at successive times (A) t = 0 (B) t = 6 (C) t = 12, and (D) t = 24, with the initial profile at t = 0 given by the solution U₀. In numerical simulations we take the discrete space and time steps Δ x = 0.1,Δ t = 0.01.

FIGURE 12

Figure 12. (A) Plot of mean position $\overline{X}(t)$ of leading (blue) and trailing (black) edges of stimulus-locked pulse as a function of time t averaged over N = 1000 trials. (B,C) Corresponding plots of the variance σ²_X(t) of the leading and trailing edges. Same parameter values as Figure 11.

Discussion

In this paper we have explored the effects of extrinsic noise on propagating pulses in a one-dimensional scalar neural field with asymmetric weights. Such a network has previously been proposed as a continuum model of direction selectivity. We have shown that the effects of noise on the wandering of the mean front position depends on properties of the underlying deterministic pulse. In the case of a freely propagating pulse, we find diffusive wandering with the mean square displacement growing linearly with time t. Moreover, in the large time limit, fluctuations in the width of the pulse can be neglected. On the other hand, if the pulse is locked to a moving pulse-like stimulus, then the wandering is described by an Ornstein–Uhlenbeck process and the mean square displacement saturates in the long time limit. However, we find that fluctuations in the pulse width can no longer be ignored.

In summary, this paper further illustrates how methods developed for studying wave propagation in stochastic PDEs can be adapted to study wave propagation in stochastic neural fields. As we have previously found for fronts, stimulus-locked waves are more robust to noise, which is a desirable property of a network performing some form of stimulus-processing such as direction selectivity.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

References

Adelson, E. H., and Bergen, J. R. (1985). Spatiotemporal energy models for the perception of motion. J. Opt. Soc. Am. 2, 284–299.

Pubmed Abstract | Pubmed Full Text

Amari, S. (1977). Dynamics of pattern formation in lateral inhibition type neural fields. Biol. Cybern. 27, 77–87.

Pubmed Abstract | Pubmed Full Text

Armero, J., Casademunt, J., Ramirez-Piscina, L., and Sancho, J. M. (1998). Ballistic and diffusive corrections to front propagation in the presence of multiplicative noise. Phys. Rev. E 58, 5494–5500.

Baker, P. M., and Bair, W. (2012). Inter-neuronal correlation distinguishes mechanisms of direction selectivity in cortical circuit models. J. Neurosci. 32, 8800–8816.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Ben-Yishai, R., Bar-Or, R. L., and Sompolinsky, H. (1995). Theory of orientation tuning in visual cortex. Proc. Natl. Acad. Sci. U.S.A. 92, 3844–3848.

Pubmed Abstract | Pubmed Full Text

Bernander, O., Douglas, R. J., Martin, K. A. C., and Koch, C. (1991). Synaptic background activity influences spatiotemporal integration in single pyramidal cells. Proc. Natl. Acad. Sci. U.S.A. 88, 11569–11573.

Pubmed Abstract | Pubmed Full Text

Brackley, C. A., and Turner, M. S. (2007). Random fluctuations of the firing rate function in a continuum neural field model. Phys. Rev. E 75, 041913.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Bressloff, P. C. (1994). A green's function approach to analysing the effects of random synaptic background activity in a model neural network. J. Phys. A 27, 4097–4113.

Bressloff, P. C. (2001). Traveling fronts and wave propagation failure in an inhomogeneous neural network. Physica D 155, 83–100.

Bressloff, P. C. (2009). Stochastic neural field theory and the system-size expansion. SIAM J. Appl. Math 70, 1488–1521.

Bressloff, P. C. (2010). Metastable states and quasicycles in a stochastic Wilson-Cowan model of neuronal population dynamics. Phys. Rev. E 85, 051903.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Bressloff, P. C. (2012). Spatiotemporal dynamics of continuum neural fields. J. Phys. A 45, 033001.

Bressloff, P. C., Cowan, J. D., Golubitsky, M., Thomas, P. J., and Wiener, M. (2001). Geometric visual hallucinations, Euclidean symmetry and the functional architecture of striate cortex. Phil. Trans. R. Soc. Lond. B 356, 299–330.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Bressloff, P. C., and Webber, M. A. (2012a). Neural field model of binocular rivalry waves. J. Comput. Neurosci. 32, 233–252.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Bressloff, P. C., and Webber, M. A. (2012b). Front propagation in stochastic neural fields. SIAM J. Appl. Dyn. Syst. 11, 708–740.

Buice, M., and Cowan, J. D. (2007). Field-theoretic approach to fluctuation effects in neural networks. Phys. Rev. E 75, 051919.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Buice, M., Cowan, J. D., and Chow, C. C. (2010). Systematic fluctuation expansion for neural network activity equations. Neural Comput. 22, 377–426.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Camperi, M., and Wang, X.-J. (1998). A model of visuospatial short-term memory in prefrontal cortex: recurrent network and cellular bistability. J. Comput. Neurosci. 5, 383–405.

Pubmed Abstract | Pubmed Full Text

Coombes, S. (2005). Waves, bumps and patterns in neural field theories. Biol. Cybern. 93, 91–108.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Coombes, S., and Owen, M. R. (2004). Evans functions for integral neural field equations with Heaviside firing rate function. SIAM J. Appl. Dyn. Syst. 4, 574–600.

Coombes, S., and Owen, M. R. (2005). Bumps, breathers, and waves in a neural network with spike frequency adaptation. Phys. Rev. Lett. 94, 148102.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

de Pasquale, F., Gorecki, J., and Poielawski, J. (1992). On the stochastic correlations in a randomly perturbed chemical front. J. Phys. A 25, 433.

Ermentrout, G. B. (1998). Neural networks as spatio-temporal pattern-forming systems. Rep. Prog. Phys. 61, 353–430.

Ermentrout, G. B., and Cowan, J. (1979). A mathematical theory of visual hallucination patterns. Bio. Cybern. 34, 137–150.

Pubmed Abstract | Pubmed Full Text

Ermentrout, G. B., and McLeod, J. B. (1993). Existence and uniqueness of travelling waves for a neural network. Proc. R. Soc. Edinburgh A 123, 461–478.

Evans, J. (1975). Nerve axon equations iv: the stable and unstable impulse. Indiana Univ. Math. J. 24, 1169–1190.

Faugeras, O., Touboul, J., and Cessac, B. (2009). A constructive mean–field analysis of multi–population neural networks with random synaptic weights and stochastic inputs. Front. Compt. Neurosci. 3:1. doi: 10.3389/neuro.10.001.2009

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Folias, S. E., and Bressloff, P. C. (2005). Stimulus–locked traveling pulses and breathers in an excitatory neural network. SIAM J. Appl. Math. 65, 2067–2092.

Gardiner, C. W. (2009). Handbook of Stochastic Methods, 4th Edn. Berlin: Springer.

Giese, M. A. (1999). Neural Field Theory for Motion Perception. Dordrecht: Kluwer Academic Publishers.

Huang, X., Troy, W. C., Yang, Q., Ma, H., Laing, C. R., Schiff, S. J., et al. (2004). Spiral waves in disinhibited mammalian neocortex. J. Neurosci. 24, 9897–9902.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Hutt, A., and Atay, F. (2005). Analysis of nonlocal neural fields for both general and gamma-distributed connectivities. Physica D 203, 30–54.

Hutt, A., Longtin, A., and Schimansky-Geier, L. (2008). Additive noise-induces turing transitions in spatial systems with application to neural fields and the swift-hohenberg equation. Physica D 237, 755–773.

Kapitula, T., Kutz, N., and Sandstede, B. (2004). The Evans function for nonlocal equations. Indiana Univ. Math. J. 53, 1095–1126.

Keener, J. P. (1981). Waves in excitable media. SIAM J. Appl. Math. 39, 528–548.

Kilpatrick, Z. P., and Bressloff, P. C. (2010a). Effects of synaptic depression and adaptation on spatiotemporal dynamics of an excitatory neuronal network. Physica D 239, 547–560.

Kilpatrick, Z. P., and Bressloff, P. C. (2010b). Spatially structured oscillations in a two-dimensional neuronal network with synaptic depression. J. Comput. Neurosci. 28, 193–209.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Kilpatrick, Z. P., Folias, S. E., and Bressloff, P. C. (2008). Traveling pulses and wave propagation failure in inhomogeneous neural media. SIAM J. Appl. Dyn. Syst. 7, 161–185.

Koch, C., and Poggio, T. (1985). “Orientation selectivity in the visual cortex?” in Models of the Visual Cortex, eds D. Rose and V. G. Dobson (New York, NY: Wiley), 408–419.

Kreyszig, E. (1978). Introductory Functional Analysis with Applications. New York, NY: Wiley.

Kuramoto, Y. (1984). Chemical Oscillations, Waves and Turbulence. New York, NY: Springer-Verlag.

Laing, C. R., and Chow, C. C. (2002). A spiking neuron model for binocular rivalry. J. Comput. Neurosci. 12, 39–53.

Pubmed Abstract | Pubmed Full Text

Liley, D. J. T., Cadusch, P. J., and Dafilis, M. P. (2002). A spatially continuous mean field theory of electrocortical activity. Network 13, 67–113.

Pubmed Abstract | Pubmed Full Text

Maex, R., and Urban, G. A. (1996). Model circuit of spiking neurons generating directional selectivity in simple cells. J. Neurophysiol. 75, 1515–1545.

Pubmed Abstract | Pubmed Full Text

Mineiro, P., and Zipser, D. (1998). Analysis of direction selectivity arising from recurrent cortical interactions. Neural Comput. 10, 353–371.

Pubmed Abstract | Pubmed Full Text

Novikov, E. A. (1965). Functionals and the random-force method in turbulence theory. Sov. Phys. JETP 20, 1290.

Nunez, P. I. (1995). Neocortical Dynamics and Human EEG Rhythms. New York, NY: Oxford University Press.

Pinto, D., and Ermentrout, G. B. (2001). Spatially structured activity in synaptically coupled neuronal networks: I. Traveling fronts and pulses. SIAM J. Appl. Math 62, 206–225.

Pinto, D., Jackson, R. K., and Wayne, C. E. (2005). Existence and stability of traveling pulses in a continuous neuronal network. SIAM J. Appl. Dyn. Syst. 4, 954–984.

Rapp, M., Yarom, Y., and Segev, I. (1992). The impact of parallel fiber background activity on the cable properties of cerebellar purkinje cells. Neural Comput. 4, 518–533.

Reichardt, W. (1961). “Autocorrelation, a principle for the evaluation of sensory information by the central nervous system,” in Sensory Communication, ed A. Rosenblith (Cambridge, MA: MIT Press), 303–317.

Richardson, K. A., Schiff, S. J., and Gluckman, B. J. (2005). Control of traveling waves in the mammalian cortex. Phys. Rev. Lett. 94, 028103.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Robinson, P. A., Rennie, C. J., Wright, J. J., Bahramali, H., Gordon, E., and Rowe, D. I. (2001). Prediction of electroencephalographic spectra from neurophysiology. Phys. Rev. E 63, 021903.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Rubin, J. E. (2004). A nonlocal eigenvalue problem for the stability of a traveling wave in a neuronal medium. Discrete Contin. Dyn. Syst. 10, 925–940.

Sagues, F., Sancho, J. M., and Garcia-Ojalvo, J. (2007). Spatiotemporal order out of noise. Rev. Mod. Phys. 79, 829–882.

Sandstede, B. (2002). “Stability of travelling waves,” in Handbook of Dynamical Systems, Vol. 2, ed B. Fiedler (Amsterdam: North-Holland), 983–1055.

Sandstede, B. (2007). Evans functions and nonlinear stability of traveling waves in neuronal network models. Int. J. Bifurcation Chaos Appl. Sci. Eng. 17, 2693–2704.

Saul, A. B., and Humphrey, A. L. (1990). Spatial and temporal response properties of lagged and nonlagged cellscells in cat lateral geniculate nucleus. J. Neurophysiol. 64, 206–224.

Pubmed Abstract | Pubmed Full Text

Schimansky-Geier, L., Mikhailov, A. S., and Ebeling, W. (1983). Effects of fluctuations on plane front propagation in bistable nonequilibrium systems. Ann. Phys. 40, 277.

Somers, D. C., Nelson, S., and Sur, M. (1995). An emergent model of orientation selectivity in cat visual cortical simple cells. J. Neurosci. 15, 5448–5465.

Pubmed Abstract | Pubmed Full Text

Steyn-Ross, M. L., Steyn-Ross, D. A., Sleigh, J. W., and Whiting, D. R. (2003). Theoretical predictions for spatial covariance of the electroencephalographic signal during the anesthetic-induced phase transition: increased correlation length and emergence of spatial self-organization. Phys. Rev. E 68, 021902.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Suarez, H., Koch, C., and Douglas, R. J. (1995). Modeling direction selectivity of simple cells in striate visual cortex using the canonical microcircuit. J. Neurosci. 15, 6700–6719.

Pubmed Abstract | Pubmed Full Text

van Santen, J. P., and Sperling, G. (1985). Elaborated reichardt detectors. J. Opt. Soc. Am. A 256, 300–320.

Pubmed Abstract | Pubmed Full Text

Wilson, H. R., and Cowan, J. D. (1972). Excitatory and inhibitory interactions in localized populations of model neurons. Biophys. J. 12, 1–23.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Wilson, H. R., and Cowan, J. D. (1973). A mathematical theory of the functional dynamics of cortical and thalamic nervous tissue. Kybernetik 13, 55–80.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Xie, X., and Giese, M. A. (2002). Nonlinear dynamics of direction-selective recurrent neural media. Phys. Rev. E 65, 051904.

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Zhang, K. (1996). Representation of spatial orientation by the intrinsic dynamics of the head-direction cell ensemble: a theory. J. Neurosci. 16, 2112–2126.

Pubmed Abstract | Pubmed Full Text

Zhang, L. (2003). On the stability of traveling wave solutions in synaptically coupled neuronal networks. Diff. Integral Equat. 16, 513–536.