Generation and annihilation of localized persistent-activity states in a two-population neural-field model

Neural Networks 46 (2013) 75–90

Contents lists available at ScienceDirect

Neural Networks journal homepage: www.elsevier.com/locate/neunet

Generation and annihilation of localized persistent-activity states in a two-population neural-field model M. Yousaf a,b , B. Kriener a,∗ , J. Wyller a , G.T. Einevoll a a

Department of Mathematical Sciences and Technology, Norwegian University of Life Sciences, P.O. Box 5003, N-1432 Ås, Norway

b

Department of Mathematics, COMSATS Institute of Information Technology, Lahore, Pakistan

article

info

Article history: Received 7 February 2013 Received in revised form 25 April 2013 Accepted 26 April 2013 Keywords: Two-population neuronal networks Working memory Integro-differential equations Bumps Persistent states Cortex

abstract We investigate the generation and annihilation of persistent localized activity states, so-called bumps, in response to transient spatiotemporal external input in a two-population neural-field model of the Wilson–Cowan type. Such persistent cortical states have been implicated as a biological substrate for short-term working memory, that is, the ability to store stimulus-related information for a few seconds and discard it once it is no longer relevant. In previous studies of the same model it has been established that the stability of bump states hinges on the relative inhibitory constant τ , i.e., the ratio of the time constants governing the dynamics of the inhibitory and excitatory populations: persistent bump states are typically only stable for values of τ smaller than a critical value τcr . We find here that τ is also a key parameter determining whether a transient input can generate a persistent bump state (in the regime where τ < τcr ) or not. For small values of τ generation of the persistent states is found to depend only on the overall strength of the transient input, i.e., as long as the magnitude and duration of the excitatory transient input are larger and/or long enough, the persistent state will be activated. For higher values of τ we find that only specific combinations of amplitude and duration leads to persistent activation. For the corresponding annihilation process, no such delicate selectivity on the transient input is observed. © 2013 Elsevier Ltd. All rights reserved.

1. Introduction A key ability of the brain is to transiently hold stimulus-related information in so-called working memory, i.e., that information is actively remembered in an online fashion even after the initial stimulus subsides (Goldman-Rakic, 1995; Wang, 2001). In a series of experiments, i.e., delayed-response task studies with awake behaving monkeys (Goldman-Rakic, 1995 and references therein), the persistent activation of groups of neurons in cortical areas such as the prefrontal cortex were identified as a neural correlate underlying this short-term memory. Various candidate models for how cortical networks may generate and sustain the selective activation of subpopulations of neurons have been put forward in the last few decades, amongst them persistent activation by thalamocortical and corticocortical loops, intrinsic cellular bistability, or attractor states of local recurrent networks (see Compte, 2006; Wang, 2001 and references therein). The idea of network attractor states

∗

Corresponding author. Tel.: +47 64965412. E-mail addresses: [email protected] (M. Yousaf), [email protected], [email protected] (B. Kriener), [email protected] (J. Wyller), [email protected] (G.T. Einevoll). 0893-6080/$ – see front matter © 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.neunet.2013.04.012

has in particular inspired many studies in the framework of spiking neuronal networks (e.g., Brunel, 2003; Compte, 2006; Compte, Brunel, Goldman-Rakic, & Wang, 2000; Compte et al., 2003; Renart, Brunel, & Wang, 2003) and neural-field models (e.g., Amari, 1977; Blomquist, Wyller, & Einevoll, 2005; Coombes, 2005; Laing, Troy, Gutkin, & Ermentrout, 2002; Machens & Brody, 2008; Oleynik, Wyller, Tetzlaff, & Einevoll, 2011; Wilson & Cowan, 1973). One pivotal element behind stable persistent activation in network models is sufficient recurrent excitatory feedback to maintain the activation, while inhibition is needed in order to keep the system from entering a state of run-away excitation (see, however, Ermentrout & Drover, 2003; Rubin & Bose, 2004). In neuralfield models the effective synaptic connectivity is therefore often assumed to be of ‘Mexican-hat’ type, i.e., nearby neurons predominantly excite each other, while cells further apart have a net inhibitory interaction, also called lateral inhibition. In combination with the assumption of translation-invariant symmetrical coupling such systems can stably sustain bump attractors, i.e., the persistent activation of a localized subpopulation of the network, e.g., Amari (1977), Blomquist et al. (2005), Laing et al. (2002), Oleynik et al. (2011) and Wilson and Cowan (1973). Neural-field models have provided a powerful tool to investigate the properties of such persistent bumps.

76

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

A large number of modeling studies have focused on investigating generic properties of persistent bumps in the absence of external inputs, in particular their existence and stability (see, e.g., Cressman, Ullah, Ziburkus, Schiff, & Barreto, 2009; Ermentrout & McLeod, 1993; Gray & Robinson, 2009; Guo & Chow, 2005; Kilpatrick & Bressloff, 2010a, 2010b; Pinto & Ermentrout, 2001; Werner & Richter, 2001; Wyller, Blomquist, & Einevoll, 2007). However, a model for working memory must also account for how such persistent-activity states may arise in response to external stimulus presentation. Towards this, several one-population model studies have investigated the effect of external input on bump properties: The seminal paper by Amari (1977) investigated the existence and stability of bumps in a simplified neural-field model of the lateral-inhibition type with stationary and spatially uniform external input. Kishimoto and Amari (1979) subsequently addressed the same problem for a whole set of different firing-rate functions, while Laing et al. (2002) extended the approach to systems allowing for multiple-bump solutions. Folias and Bressloff (2004) also investigated bumps in a onepopulation model with spatially localized external inputs and found that sufficiently large inputs can stabilize bump states, and that reduction of the input amplitude may induce a Hopf instability and the conversion of stable bumps into breather-like oscillatory waves. Rubin and Troy (2004) further showed that bumps can exist and be linearly stable in such systems, even if they have a symmetric off-center synaptic architecture instead of recurrent excitation. To analyze the effect of spatially and temporally localized input on bump formation, Laing and Chow (2001) followed a combined approach of field modeling and spiking neuron simulations. They showed how bump attractors can be persistently activated by targeted transient localized external input, and that such a sustained activation can be switched off again by synchronizing excitatory input. In a very recent paper Marti and Rinzel (2013) moreover study the activation of multiple-bump states in lateral inhibition field models by repeated short subthreshold stimulation with various spatially structured inputs. They find that the resulting bump pattern then depends on which of the possible multistable states is closest to the average stimulus pattern. Several studies also analyzed bump formation in response to external inputs in network simulations of spiking neurons with realistic biophysical parameters. Camperi and Wang (1998), e.g., studied how bumps can be initiated in a ring network with lateral inhibition when the network is stimulated by spatially modulated external input, and can in turn be deactivated by a global inhibitory input. They showed, moreover, how neuron properties, such as cellular bistability, can enhance the stability of such stationary bumps. A similar effect of activation and deactivation of bumps by excitatory stimulation as described by Laing and Chow (2001) was investigated by Gutkin, Laing, Colby, Chow, and Ermentrout (2001) in a spiking neuron network model with lateral inhibition, that used more biologically realistic Hodgkin–Huxley neurons equipped with neurophysiologically motivated synaptic time courses and parameters. In both studies the importance of asynchronous spiking – which corresponds to time-independent rate dynamics in field models – for bump maintenance was stressed, while excess synchrony was shown to work as a deactivation signal. In the line of bump formation studies in spiking neuron networks the study by Rubin and Bose (2004) should not go unmentioned. They demonstrated that by taking into account the intricate interplay between cellular excitability and input structure localized stimulation can lead to bump formation also in networks with purely excitatory distance-dependent coupling. Even if one-population lateral-inhibition type models can reflect some aspects of cortical network behavior, they inherently

cannot capture the individual dynamics of the two basic neuron populations in cortex, i.e., excitatory and inhibitory cells. These populations are present in different fractions (about four to five times more excitatory than inhibitory neurons), and will in general have different dynamical time scales, synaptic footprints, firingrate functions, etc.. Two-population neural-field models are the minimal model for incorporating this aspect of cortical dynamics, and have also been used in the study of existence and stability of persistent bump states in the absence of external inputs (Blomquist et al., 2005; Oleynik et al., 2011). In a recent study we investigated how stationary localized external input modifies the bump structure in such a two-population model (Yousaf, Wyller, Tetzlaff, & Einevoll, 2013). In terms of being a model for working memory, however, the more interesting case is when transient external inputs both generate and annihilate persistent bump activity, and this is the topic of the present study. We generalize a two-population model considered in Blomquist et al. (2005), Oleynik et al. (2011) and Yousaf et al. (2013) by adding terms corresponding to transient external inputs and consider the following Wilson–Cowan-like neural-field model:

∂ ue (x, t ) = −ue (x, t ) + ∂t



∞

ωee (x′ − x)  ∞ × Pe (ue (x′ , t ) − θe )dx′ − ωie (x′ − x) −∞

−∞

× Pi (ui (x′ , t ) − θi )dx′ + He (x, t )  ∞ ∂ ωei (x′ − x) τ ui (x, t ) = −ui (x, t ) + ∂t −∞  ∞ × Pe (ue (x′ , t ) − θe )dx′ − ωii (x′ − x)

(1)

−∞

× Pi (ui (x′ , t ) − θi )dx′ + Hi (x, t ). Here ue and ui stand for excitatory and inhibitory activity levels, Pm , with m = e, i, are firing-rate functions, ωmn denotes the distancedependent connectivity strengths, θm and θi are the respective firing-threshold values of the excitatory and inhibitory populations, and τ is the relative inhibition time, i.e., the ratio between the inhibitory and excitatory time constants (Blomquist et al., 2005). Finally, He and Hi represent the spatiotemporal, transient external inputs to the excitatory and inhibitory populations, respectively. Wyller et al. (2007) investigated the Turing instability and pattern formation in this model (1) in the case of no external input, i.e., Hm ≡ 0. They explored the generation and formation of stationary periodic spatial patterns from a spatially homogeneous rest state by Turing-type instabilities, while Blomquist et al. (2005) derived the conditions for existence, uniqueness and stability of stationary symmetric solutions of (1) for Hm ≡ 0. The generic picture consists of two bump-pairs (BPs) in this case, one narrow BP and one broad BP. The narrow BP is generically unstable while the broad BP is stable for small and moderate value of the relative inhibition time τ . The generalization of this work to a wider class of temporal kernels was performed by Oleynik et al. (2011). In Yousaf et al. (2013) the model (1) was investigated for spatially dependent, but temporally constant external input, i.e., Hm (x, t ) = Hm (x). In particular, the existence and stability of bump solutions of (1) were investigated as a function of the amplitude of the external input. At most four BPs were detected, with a maximum of three of them being stable. As a direct extension of Yousaf et al. (2013), we here study the full dynamics of (1) when both spatially and temporally varying external input Hm (m = e, i) are present. In particular, we study the phenomenon of emergence and annihilation of persistent activity bumps in such networks in response to transient external inputs. Such a selective persistent activation and deactivation of

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

neural activity in response to external stimulus cues is observed in the prefrontal cortex during memory tasks (Goldman-Rakic, 1995), and a theoretical understanding of how and under which conditions cortical networks can display such dynamics is thus of great interest (Compte, 2006; Wang, 2001). To our knowledge this is the first study of bump-pair formation and annihilation in response to spatiotemporal external inputs in the context of two-population neural-field models. For conceptual simplicity and mathematical convenience we model the transient inputs as spatiotemporally separable functions, and investigate necessary conditions for generating and annihilating persistent bump states by a combination of analytical and numerical techniques. A first general finding is that a necessary condition to evoke persistent activity in the network is that the excitatory external input must be greater than a particular threshold value. In the ensuing numerical investigations we focus on the case with persistent bump states in the absence of stabilizing stationary external inputs (Blomquist et al., 2005), i.e., the only external input is the transient external input turning the bump state ‘on’ (generation) or ‘off’ (annihilation). In Blomquist et al. (2005) the relative inhibitory time constant τ , cf. (1), was found to be a critical factor in determining bump stability: The generic situation corresponded to two bump states, one with narrow BPs and one with broad BPs. The narrow BPs were always unstable, while the broad BP states were found to be stable for small values of τ , but unstable for τ larger than a critical value τcr . Here, we focus on the generation and annihilation of persistent broad-bump states, and thus τ < τcr , and in particular how these dynamic processes depend on τ . We indeed find the generation and annihilation to strongly depend on the value of τ : For τ much smaller than τcr , generation of the persistent states is found to only depend on the overall strength of the transient input, i.e., as long as the amplitude and duration of the transient input are strong and/or long enough, the persistent state will be activated. For larger values of τ , i.e., τ still smaller than but closer to τcr , the dependence on the amplitude and duration becomes non-monotonic in the sense that only specific combinations of amplitude and duration lead to persistent activation. No such delicate dependence on transient input is found for the annihilation process, however. The paper is organized in the following way: After describing the model in Section 2, we show in Section 3 that the solution of the initial value problem of this system is globally bounded and give a general necessary condition for emergence of persistent activity in the network. In Section 4 we demonstrate numerically the generation and annihilation of persistent-activity states for various types of transient external inputs. Section 5 contains conclusions and discussion of the results. 2. Model The two-population neuronal field model described by (1) can in a more compact form be written as

77

Fig. 1. Diagrammatic view of the two-population model (2) with spatiotemporal external input Hm (x, t ), m = e, i.

and τ is the ratio between inhibitory and excitatory time constants and is therefore called the relative inhibition time (Blomquist et al., 2005). Further, the conversion of activity levels to action-potential firing is done by means of the firing-rate functions Pm (m = e, i). These functions constitute a one-parameter family of smooth and non-decreasing functions mapping the set of real numbers onto the unit interval [0, 1]. The functions Pm (m = e, i) are parameterized by a positive steepness parameter βm (m = e, i). As an example of a firing-rate function Pm , we have Pm (u) =

1 2

(1 + tanh(βm u)).

(4)

For βm → ∞, the firing-rate function Pm approaches the unit step function (Heaviside function) Θ :

Θ ( u) =



0, 1,

u<0 u ≥ 0.

(5)

The functions ωmn (m, n = e, i) in (2) are called connectivity functions. These functions model the synaptic connection strength in the network and are assumed tobe positive, real valued, bounded, ∞ symmetric and normalized, i.e., −∞ ωmn (x)dx = 1. The functions ωmn are parameterized by means of synaptic footprints σmn (m, n = e, i), i.e.,

ωmn (x) =

1

σmn

Φmn (ξmn ),

ξmn =

x

σmn

.

Here Φmn is a non-dimensional scaling function. As an example, with Gaussian connectivity, this scaling function is given as 1

2 Φmn (ξmn ) = √ exp(−ξmn ). π

(7)

The interaction diagram in Fig. 1 illustrates the model (2).

∂ ue = −ue + ωee ⊗ Pe (ue − θe ) − ωie ⊗ Pi (ui − θi ) + He (2a) ∂t ∂ ui τ = −ui + ωei ⊗ Pe (ue − θe ) − ωii ⊗ Pi (ui − θi ) + Hi . (2b) ∂t The operator ⊗ in (2) defines the spatial convolution integral, given

2.1. Spatiotemporal functions describing transient external input

as

Hm (x, t ) = Cm Sm (x) hm (t − tm0 ),

[ωmn ⊗ Pm (um − θm )](x, t )  ∞ = ωmn (x − x′ )Pm (um (x′ , t ) − θm )dx′ .

(3)

−∞

The functions Hm (m = e, i) model the spatiotemporal external inputs, the parameters θe and θi are the firing threshold values,

(6)

We assume the following factorization of the spatiotemporal external input both for generation and annihilation of persistent activity: m = e, i

(8)

where the functions hm (m = e, i) represent the temporal dependence, Sm represent the spatial dependence, and Cm the amplitudes of the external inputs, while tm0 are the onset times for the temporal functions (hm = 0 for t < 0). The spatial functions Sm are continuous and symmetric, i.e., Sm (−x) = Sm (x).

78

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

Fig. 2. Stationary symmetric solutions, i.e., bumps, (solid curves) to the model (2) with no external input (Blomquist et al., 2005). Excitatory and inhibitory bumps are represented by the red and green solid curves, respectively. Top row: Excitatory (a) and inhibitory (b) activity profiles for narrow bump solution. Bottom row: Excitatory (a) and inhibitory (b) activity profiles for broad bump solution. Threshold values: θe = 0.12, θi = 0.08. The connectivity functions are Gaussian, cf. (6), with the following synaptic footprints: σee = 0.35, σei = 0.48, σie = 0.60, σii = 0.69. Pulse-width coordinates a and b are illustrated in panels (c) and (d), respectively. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

2.1.1. Spatial part of external input The following three types of spatial functions are considered: 1. Bump solutions of (2) with Hm (x, t ) = 0 (Blomquist et al., 2005): Se (x) = Wee (x + a) − Wee (x − a) − Wie (x + b) + Wie (x − b)

(9a)

Si (x) = Wei (x + a) − Wei (x − a) − Wii (x + b) + Wii (x − b).

(9b)

Here a and b are the excitatory and inhibitory pulse-width coordinates (see Fig. 2), respectively, and Wmn (x) is the antiderivative given by Wmn (x) =

2t /Tm1 Tm 2(Tm − t )/Tm (Tm − Tm1 ) 0

 hm ( t ) =

0 < t < Tm1 Tm1 ≤ t ≤ Tm otherwise

(13)

where Tm for m = e, i is the time such that hm (t ) > 0, i.e., the input duration, and Tm1 for m = e, i is the time to peak satisfying 0 < Tm1 < Tm . This function has been used previously to model transient input onto neural networks as it in a simple manner allows for probing the effect of varying the rise time of the input and captures measured activity profiles from thalamus to rodent barrel cortex (Pinto, Hartings, Brumberg,  ∞& Simons, 2003). Note also that the function is normalized, i.e., −∞ hm (t )dt = 1.

x



ωmn (y)dy.

(10)

0

Examples of such bump solutions are shown in Fig. 2. 2. Gaussian function:

    2 x Sm (x) = , √ exp − ρm ρm π 1

m = e, i

(11)

with ρm being width parameters (cf. Fig. 3(a)).

    x   , exp −  Sm (x) = 2ρm ρm 

m = e, i

with ρm being width parameters (cf. Fig. 3(b)).

3. General properties of the model In this subsection we elaborate different general properties of the model. First, we prove that the solutions of the initial value problem (2) are uniformly bounded provided both the initial conditions and the external input functions are bounded and continuous. Secondly, we give general conditions for evoking activity in the network. 3.1. Boundedness

3. Exponentially decaying function: 1

2.1.2. Temporal part of external input For the temporal part of the external input hm (t ) we assume the following triangular temporal function, cf. Fig. 3(c):

(12)

The normalization condition imposed on the connectivity functions ωmn for m = e, i, and the fact that the firing-rate functions satisfy 0 ≤ Pm (u) ≤ 1 uniformly for all x and t, lead to 0 ≤ [ωmn ⊗ Pm (um − θm )](x, t ) ≤ 1.

(14)

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

79

where the function γm is given as

γm (t ) =



t

 exp

s



τm

−tm0

hm (s)ds.

(25)

For the triangular temporal function (13), we get

 0,    (0) γm (t ), γm (t ) = γ (0) (Tm1 ) + γm(1) (t ),    m(0) γm (Tm1 ) + γm(1) (Tm ),

t ≤0 0 < t ≤ Tm1 Tm1 < t ≤ Tm t > Tm

(0)

(26)

(1)

where the auxiliary functions γm and γm are given as γm(0) (t ) =

2τm Tm Tm1

 γm(1) (t ) = 2τm 

Fig. 3. Spatial and temporal components of the external input (8) used in the numerical simulations in this paper. (a) Gaussian function given by (11). (b) Exponential function defined by (12). Throughout the present study (and in this figure) the width parameter ρm is set to 0.16 both for m = e and m = i. (c) Triangular temporal function (13) characterized by the peak time Tm1 and total duration time Tm . Note that kernels are normalized to unit area. In the figure Tm1 = 5 and Tm = 10.

Hence we find the bounding inequalities

∂ ue ≤ −ue + 1 + He ∂t ∂ ui −ui − 1 + Hi ≤ τ ≤ −ui + 1 + Hi ∂t −u e − 1 + H e ≤

(15a) (15b)

for the rate of change of the activities ue and ui . Simple computation reveals the bounds me (x, t ) ≤ ue (x, t ) ≤ Me (x, t )

(16a)

mi (x, t ) ≤ ui (x, t ) ≤ Mi (x, t )

(16b)

for ue and ui . Here the lower and upper bounding functions are given as: me (x, t ) = αe (t ) (Ve (x) + 1) − 1 + [αe ∗ He ](x, t )

(17)

Me (x, t ) = αe (t ) (Ve (x) − 1) + 1 + [αe ∗ He ](x, t )

(18)

mi (x, t ) = τ αi (t ) (Vi (x) + 1) − 1 + [αi ∗ Hi ](x, t )

(19)

Mi (x, t ) = τ αi (t ) (Vi (x) − 1) + 1 + [αi ∗ Hi ](x, t )

(20)

where Ve and Vi are the initial conditions of the system (2). The temporal convolution αm ∗ Hm is given as

[αm ∗ Hm ](x, t ) =

t



αm (t − s)Hm (x, s)ds,

m = e, i.

(21)

0

Here the functions αm for m = e, i are given as

αm (t ) =

1

τm

 exp −

t

 (22)

τm

where τm is given as

τm =



1,

τ,

m=e m = i.

(23)

The temporal convolutions that appear in the bounding functions (17) can now be computed for particular choices of Hm . For example, when Hm is given by means of the factorization (8), we get

[αm ∗ Hm ](x, t ) = Cm Sm (x)



t

αm (tm0 − s) hm (s − tm0 ) ds

0

= Cm Sm (x) αm (t − tm0 ) γm (t − tm0 )

(24)





t

t

τm 

exp exp



τm



(t − τm ) + τm



(Tm − t + τm ) + exp

(27) 

Tm1

τm

Tm (Tm − Tm1 )



(Tm1 − Tm − τm )

  . (28)

The corresponding bounding functions me , Me , mi and Mi are plotted in Fig. 4 for an example using the ‘broad’ bump solution (9) as a spatial function together with the triangular temporal function (13). Notice that in case of stationary inputs Hm (x, t ) = Hm (x), the bounds (16) simplify to the bounds deduced in Yousaf et al. (2013). An important result is that the bounding inequalities (16) can be converted to

|ue (x, t ) − Ue (x, t )| ≤ F (t )   t |ui (x, t ) − Ui (x, t )| ≤ F τ

(29a) (29b)

with the function F defined by F (t ) = 1 − exp(−t ).

(30)

Here Ue and Ui are given as Ue (x, t ) = αe (t ) Ve (x) + [αe ∗ He ](x, t )

(31a)

Ui (x, t ) = τ αi (t ) Vi (x) + [αi ∗ Hi ](x, t )

(31b)

which are the solutions of the two-population model (2) without non-local interaction terms. i.e.

∂ ue (x, t ) = −ue (x, t ) + He (x, t ) ∂t ∂ ui ( x , t ) τ = −ui (x, t ) + Hi (x, t ). ∂t

(32a) (32b)

The bounds (29) thus express a comparison result between the temporal evolution of the solutions of the full system (2) and the temporal evolution of the solutions of simplified system (32) (with the same initial conditions). Notice that the derivation holds true for the general case of firing-rate functions, such that Pm (u) ∈ [0, 1], and also in the multidimensional case. The initial value problem of (2) is locally well-posed in the Banach space of bounded continuous functions provided the connectivity functions are absolute integrable and satisfy the Hölder condition, the derivatives of the firing-rate functions are smooth, bounded and the external input functions are continuous in both variables. The proof of this fact is based on the contraction principle for complete metric spaces (Banach’s fixed point theorem) in a way analogous to Potthast and Graben (2009) for the one-population model. Notice that for a unit step-function limit of the firing-rate function, the Banach fixed point theorem is not applicable in the proof of local well-posedness due to the discontinuity of the unit step function, see Potthast and Graben (2009) for details.

80

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

Fig. 4. Example illustration of boundedness of the solutions of the system (2) where the spatial part Sm of Hm is given by the broad bump pair solution (9) with parameters as in Fig. 2, and the temporal part hm by the triangular function (13). First and second rows correspond to t = 3 and t = 20 after stimulation onset t = 0, respectively. Upper and lower bounds for the solutions are shown as black and blue curves, respectively. The red curves represent the excitatory activity (first column) while green curves represent inhibitory activity (second column). Parameters used: τ = 2.5, Tm1 = 2, Tm = 9, and Cm = 4 both for m = e and m = i. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Let us next assume that lim [αm ∗ Hm ](x, t )| = 0.

t →∞

(33)

From the solutions (31), we then conclude that (Ue , Ui ) → (0, 0) as t → ∞ uniformly in x. In this case, we get the asymptotic bounds

|um (x, t → ∞)| ≤ 1.

(34)

This means that the corresponding initial value problem of (2) is globally well-posed under the condition (33) (see Fig. 4). Notice that the condition (33) holds true for the concrete example given by (8)–(13). We immediately get the following theorem from (31): Theorem 1. Assume that the following conditions are fulfilled for x∈R 1. Vm (x) < θm for m = e, i. 2. Hm (x, t ) ≤ θm for m = e, i. Then Um (x, t ) < θm ,

m = e, i

(35)

for t > 0 and x ∈ R Now, let us consider the full system (2) with initial conditions um (x, t )|t =0 = Vm (x),

m = e, i

(36)

and assume that the firing-rate functions are approximated by means of the unit step function, i.e. Pm (u) = Θ (u) for m = e, i. We get the following theorem:

1. Vm (x) < θm , m = e, i. 2. Hm (x, t ) ≤ θm for 0 ≤ t ≤ T ∗ and m = e, i. Then um (x, t ) < θm ,

m = e, i

for x ∈ R and 0 ≤ t ≤ T

For the sake of completeness, the proof of Theorem 2 is given in the Appendix. 3.2. General necessary condition for the emergence of persistent activity From Theorem 2 above we get the following result: If we assume that the initial activity levels are below the firing threshold values θm , and also the spatiotemporal external inputs keep the system uniformly below the threshold values θm for a given time interval, then the dynamical evolution of the activities is described by means of the expressions (31), i.e., the local linear dynamics. The only way to involve the non-linear dependences and excite activity patterns for Vm (x) < θm is thus if one of the external inputs He or Hi becomes sufficiently strong at some point in space and time. We summarize this result in the following theorem: Theorem 3. Assume that the following conditions are fulfilled for x∈R 1. Vm (x) < θm for m = e, i. 2. Hm (x, t ) ≤ θm for m = e, i for 0 ≤ t ≤ tcr . Moreover, we assume that there is a point xcr ∈ R such that

Theorem 2. Assume that the following conditions are fulfilled for x∈R

(37)

∗

• ue (xcr , tcr ) = θe (or ui (xcr , tcr ) = θi ). • He (xcr , tcr ) > θe (or Hi (xcr , tcr ) > θi ).

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

81

with the initial conditions ue |t =tcr < θe ,

ui |t =tcr = θi .

The non-local terms ωem ⊗ Θ (ue − θe ) for m = e, i are zero at t = tcr . By continuity, these terms remain zero for some interval tcr < t < T1∗ . From the Eq. (39a) we hence find that ue is a decreasing function of time t on the interval tcr < t < T1∗ . The non-local terms ωem ⊗ Θ (ue − θe ) = 0 for m = e, i by repeating the same argument for time interval T1∗ < t < T2∗ . Hence, by following the induction principle both terms ωem ⊗ Θ (ue − θe ) for m = e, i are zero for all t ≥ tcr . Thus we get the system

∂ ue = −ue − ωie ⊗ Θ (ui − θi ) ∂t ∂ ui τ = −ui − ωii ⊗ Θ (ui − θi ) + Hi . ∂t

Fig. 5. Visualization of Theorem 3 showing the evolution of the emergence of the activity profile, i.e., t = 0, t < tcr , t = tcr and t > tcr .

Then, the dynamical evolution of the activity levels is resumed by the full model (2) at t = tcr with finite contribution from non-local terms ωmn ⊗ Pm (um − θm ). Proof. By proceeding as in the previous theorem, we conclude that um (x, t ) < θm for 0 ≤ t ≤ tcr . Both um and ∂∂utm for m = e, i are continuous with respect to time at t = tcr . Hence the dynamical evolution of the excitatory and inhibitory activities are described by the linear system (32) at t = tcr . Simple computation now reveals that the speed condition:

∂ ue > 0, ∂t



 ∂ ui or >0 ∂t

(40b)

This means that the excitatory field will die out with time, while the inhibitory drive Hi (x, t ) may cause the formation of an inhibitory activity field. Notice that Theorem 3 is a general result which applies to both strong excitatory and inhibitory external inputs. In the Wilson– Cowan-type two-population model considered here, persistent activity is evoked solely by excitatory external input. This external input must moreover be greater than the corresponding threshold value to allow for the system to assume the persistent state. On the other hand inhibitory activity plays a major role in annihilating the persistent-activity state, for which moreover it is not necessary to have a similar condition as extracted in Theorem 3 (i.e., Hi (x0 , t∗ ) > θi ), i.e., inhibition does not have to be superthreshold to annihilate activity. This is due to the fact that there is already inhibitory activity in the network which requires just a little push to disturb the balance to ultimately annihilate the activity. 4. Results

(38)

is satisfied at (x, t ) = (xcr , tcr ). This means that the activity levels evolve into a regime for which we have to take into account finite contributions from non-local terms.  Hence in order to detect the time at which we get patterns caused by the time-dependent external drive, we solve the equation um (x, t ) = θm subject to the constraint Hm (x, t ) > θm . The solution of this problem is given by a curve in the (x, t )-plane. The time tcr , and corresponding xcr , is given as the minimum value of t in this curve, i.e. tcr < t∗ for which um (x∗ , t∗ ) = θm and Hm (xcr , tcr ) > θm . Graphically, Theorem 3 is demonstrated in Fig. 5, where the red curve represents the initial condition of the system (2). Due to continuity, the solutions (black dotted line) remain lower than threshold value (horizontal dotted green line) with no contribution from non-local terms. The intersection condition um (x, t ) = θm is given by the blue curve and the solution of the full model with the contribution from non-local terms in (2) is given by the black solid curve. Let us now consider a special case where He (x, t ) ≡ 0 and that there is a critical space–time point (xcr , tcr ) for which ui (xcr , tcr ) = θi and Hi (xcr , tcr ) > θi . In this case the evolution equations for t > tcr read

∂ ue = −ue + ωee ⊗ Θ (ue − θe ) − ωie ⊗ Θ (ui − θi ) ∂t ∂ ui τ = −ui + ωei ⊗ Θ (ue − θe ) − ωii ⊗ Θ (ui − θi ) + Hi ∂t

(40a)

(39a) (39b)

In this section we numerically investigate the emergence and annihilation of persistent activity in response to spatiotemporal external input (8). The numerical simulations are based on the numerical code developed and used in Yousaf et al. (2013), now extended to the case of transient spatiotemporal external inputs. The code, based on a fourth-order Runge–Kutta method in time, evaluates the spatial activity profile for each time step, as well as the intersection points (pulse-width coordinates (a, b)) of the activity profiles with the respective threshold values. For a detailed description, see Blomquist et al. (2005) and Yousaf et al. (2013). A first requirement for memory storage is that the persistent bump state is stable in the absence of external input. We thus take advantage of the previous work of Blomquist et al. (2005) which investigated existence and stability of bump states for our model system in this setting: In the generic situation two stationary symmetric solutions, so-called bump-pairs (BPs), were found, cf. (9) and Fig. 2. The narrowest BP, the ‘narrow BP’, was found to be generically unstable, while the ‘broad BP’ was found to be stable for small and moderate values of the relative inhibition time τ (Blomquist et al., 2005; Yousaf et al., 2013). In the stable regime the broad BP was also seen to act as a global attractor for a wide class of spatially localized initial conditions. In the present investigations we thus consider the same model and focus on the situation where the broad BP is stable, i.e., τ -values smaller than the critical value τcr . Further, the model parameters used throughout are the ones used in Blomquist et al. (2005) and Yousaf et al. (2013) (and Fig. 2), i.e., θe = 0.12, θi = 0.08, and Gaussian synaptic footprints where σee = 0.35, σei = 0.48, σie = 0.60, σii = 0.69. This parameter choice implies that τcr = 3.03 (Blomquist et al., 2005).

82

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

b

c

Fig. 6. Emergence of persistent activity for our triangular temporal function he (13) and the broad bump (9a) as the spatial part Se of the external input function (8). (a) Depicts the excitatory external input He (x, t ) presented to the system (2) in order to evoke activity. (b) Excitatory population activity ue (x, t ). (c) Inhibitory population activity ui (x, t ). Parameters describing excitatory transient input: te0 = 0, Te1 = 3, Te = 9 and Ce = 4, τ = 2.5, other parameters as in Fig. 2.

4.1. Generation of localized persistent activity Fig. 6 shows an example for how a persistent bump state can be generated from a network without any prior activity (Um (x, t ) ≡ 0, Hm (x, t ) ≡ 0, m = e, i) by a transient excitatory input He (x, t ) > θe (Fig. 6(a)). We observe that a persistent activation remains also after the external input is switched off (Fig. 6(b) and (c)), the hallmark of persistent-activity states in working memory. Next we will explore in detail how the successful generation of bump states depends on the detailed temporal and spatial properties of the transient input. 4.1.1. Dependence on temporal function in external input We first investigate the dependence of the temporal part he (t ) of the external input, keeping the spatial part Se (x) fixed. Specifically we choose this spatial part of the external input to correspond to the excitatory part of the broad bump solution (9a), depicted in Fig. 2. The rationale is that a spatial profile of the external input similar to the spatial profile of the target persistent state should facilitate the generation of this state. In Fig. 7 we show ‘parameter-planes’ of bump-state generation, where throughout the manuscript we use the color code such that red corresponds to a resulting state of persistent activation, while green corresponds to a final state of no network activity. In the figure the amplitude Ce and duration Te of transient input span the axes, and the different panels correspond to different rise times (i.e., time to peak) Te1 and relative inhibitory time constants τ . A strong dependence on the relative inhibition time τ is observed: For τ = 1.5 and τ = 2, i.e., much smaller than the critical time constant τcr = 3.03 (Blomquist et al., 2005), the Ce –Te -planes are very similar, cf. Fig. 7(a)–(f). The separation line between successful generation and failure of generation is then unique and a straight line. We see that for a fixed time to peak Te1 , the input amplitude Ce required for generation increases with the total duration Te . This is due to the normalization of the triangular function he (t ), implying a lower peak of he (t ) for longer duration Te that needs to be compensated for by a larger Ce . Further, as Te1 increases (from left to right), the minimal duration Te , i.e., starting point of the xaxis, increases correspondingly, because the total duration must be equal to or larger than the rise time Te1 . However, the resulting Ce –Te -planes appear the same, with the unphysical part simply cut out (gray-shaded area). When τ increases from 2 to 2.5, there is a dramatic change in the Ce –Te -plane Fig. 7(g)–(l). Now there is a region where for one duration time Te , there are three different Ce -regimes, two of which lead to failure of bump-state generation (marked by A′ and C′ in Fig. 7(g)), and one intermediate Ce -range leading to successful generation (marked by B′ ). The region for successful activation moreover shrinks the closer τ becomes to τcr , see Fig. 7(j)–(l). In order to build a better intuition for this observation, we plot in Fig. 8

the temporal evolution of the system in response to external input in the pulse-width plane at the selected points A, B, C for τ = 2 (Fig. 7(d)), and A′ , B′ , C′ for τ = 2.5 (Fig. 7(g)). The situations A (for τ = 2) and A′ (for τ = 2.5) in Fig. 8 show cases of unsuccessful bump-state generation. In the simulation the system is initially set at a Ce –Te -pair outside the basin of attraction (BOA) of the bump attractor, indicated by the red diamond (I.C.). From there it decays to the trivial fixed-point (a, b) = (0, 0) prior to the arrival of the transient external input at t = 0. For this choice of values for Ce and Te , the external input He (x, t ) is not sufficient to make ue reach the threshold θe , so that there is no further evolution of (a, b) from (0, 0), in line with Theorem 3 in Section 3. The attractor corresponding to the persistent bump state is therefore never reached. For the values indicated by B and B′ in Fig. 8 with slightly larger input amplitude Ce on the other hand, the transient input drives the excitatory population across the threshold. When it reaches a transient pulse-width a close to 0.3, inhibition is activated, as can be seen from the abrupt rise in the inhibitory pulse-width coordinate b. Eventually the system evolves to the bump-state attractor (blue asterisk) and settles there. Persistent activity is thus successfully established. However, for larger values of τ closer to the critical value τcr , i.e., for a slower evolution of the inhibitory dynamics, both the activation of the inhibitory population, and the settling into the bump attractor takes longer with pronounced damped oscillations (‘breather-like activity’) as can be seen from the extended spiral the black curve describes in Fig. 8 (B′ ). For the third set of Ce –Te -values, i.e., even larger input amplitude Ce (indicated by the points C and C′ in Fig. 8), we observe pronounced differences between the two τ -value cases. While for τ = 2 the attractor state is still reached (Fig. 8 (C)), for larger τ that is not the case (Fig. 8 (C′ )). This can be understood as follows: First, for larger input amplitude Ce , more of the excitatory population is activated, i.e., the maximal transient in a becomes larger. This in turn will also recruit more inhibition, leading to a larger maximal transient b. That effect is even amplified for larger τ , because if τ is larger, inhibition will kick in later, allowing a to initially grow even larger than for smaller τ (cf. a-axis in Fig. 8 (C) and (C′ )). This in turn leads to a yet broader recruitment of lateral inhibition, i.e., a larger transient in b, that moreover is longer lasting due to the larger τ (cf. Fig. 8(C′ )). These larger excursions in the (a, b)-plane can eventually lead to an ‘‘overshoot’’, such that the trajectory misses the BOA and decays back to the trivial fixedpoint (a, b) = (0, 0). It is thus the trade-off between the size of the BOA and the increase of longer-lasting, wider transients for larger τ that determines whether the system will settle into the stable bump pair or relax back to the trivial attractor. For completeness we show in Fig. 9 parameter-planes of bumpstate generation with the amplitude Ce and peak time Te1 spanning the two axes, for a set of different duration times Te . Also here,

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

83

Fig. 7. Generation of persistent bump states as a function of amplitude Ce and stimulus duration Te of the excitatory external input for different values of the time to peak Te1 (13). Green color represents failure of activation, while red color indicates successful generation. The gray color represents the region which cannot be investigated since Te must be larger than Te1 . Rows one to four represent the Ce –Te -planes for increasing inhibitory time constant τ for three different peak times Te1 . The spatial part Se (x) of the external input He (x, t ) is chosen to be the broad BP given in (9a). The points marked by A, B, C in (d) and A′ , B′ , C′ in (g) denote the Ce –Te -pairs for which Fig. 8 shows the underlying pulse-width coordinate dynamics. The rest of the parameters are the same as in Fig. 2. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

we see the clear difference between the small-τ (τ = 1.5, 2) and large-τ cases (τ = 2.5, 3), i.e., Fig. 9(a)–(h) versus Fig. 9(i)–(p). For τ = 1.5 and 2 the values of Ce needed to activate the system increase with increasing Te , in agreement with our previous observations. There is, however, hardly any dependence on the time to peak Te1 , i.e., the slope of activation is irrelevant in these cases. For τ = 2.5 and 3, on the other hand, there is a clear dependence on the slope of the external input, in that the region of successful activation is largest for small Te1 , and activation gets harder and eventually fails if Te1 becomes too large, cf. Fig. 9(j)–(l). For these large values of τ , there is, moreover, no activation at all for too small durations Te (Fig. 9(i) and (m)–(o)), as already seen in Fig. 7(g)–(l). So far we have established that for our example recurrent model, the transition between the ‘small-τ ’ regime and the more complex ‘large-τ ’ regime occurs between τ = 2 and τ = 2.5. To investigate this transition in more detail, we show in Fig. 10

the results for bump generation in the Ce –Te plane for three intermediate values of τ , i.e., τ = 2.2, τ = 2.3, and τ = 2.4, with the rise time fixed to Te1 = 3. For τ = 2.2 we already observe a small island for large amplitudes Ce and short durations Te where bump activation fails. This region expands with increasing τ until for τ = 2.4 there is a merging of the non-activation regions, implying a failure of the bump-state activation system when Te is sufficiently small, regardless of the magnitude Ce of the amplitude. The qualitatively different dependence of bump-state generation on the parameters Ce , Te , Te1 for different values of the relative inhibitory time constant τ , is an interesting finding. It implies that the system activation is an ‘all-or-nothing’-effect for small τ , while for larger τ the conditions on the external stimulus to evoke persistent activity are much stricter and thus more specific. This could be advantageous from a coding perspective, in that it means that the system is more selective in terms of relevant external stimuli.

84

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

Fig. 8. Parametric plots of pulse-widths a (excitatory population) and b (inhibitory population) for the amplitude–duration-time pairs Ce , Te denoted by A, B, and C in Fig. 7(d), and A′ , B′ , and C′ in Fig. 7(g). In the simulation the system is initially set at the Ce –Te -pair indicated by the red diamond (initial condition, ‘I.C.’), outside the basin of attraction of the bump state attractor (blue asterisk), from where it decays (dotted green curve) to the trivial fixed-point (a, b) = (0, 0) (green diamond) prior to the arrival of the transient external input at t = 0. The dynamical evolution after arrival of external input at t = 0 is shown as a solid line. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

4.1.2. Dependence on spatial function in external input We now turn to a more detailed investigation of how the generation of persistent bump states varies with different choices of spatial function Se (x) for the excitatory external input. In Fig. 11 we first show results from comparing use of the excitatory part of the broad BP solution for Se (x) (used in the investigations above) with the situation when the narrow BP solution is used instead. A larger portion of the Ce –Te parameter-plane is shown than in Figs. 7, 9 and 10. While we observe that the parameter-planes are qualitatively similar, a detailed investigation shows that the use of the broadBP solution in the transient input gives a larger range of successful bump generation than the narrow-BP solution: the smallest value of Ce needed for bump generation is 3.10 for the broad BP (occurring for Te = 7.36, Ti = 2; see Fig. 11(b)) while it is 6.20 for the narrow BP (occurring for Te = 10.08, Te1 = 2; see Fig. 11(d)). One reason for the lower threshold value for Ce might be that the broad BP is also the attractor that constitutes the persistent-activity state. However, as seen in Fig. 2 the broad-BP solution also has a higher maximal amplitude than the narrow-BP solution, suggesting another reason why a smaller amplitude Ce is needed in order to generate the persistent bump state. When the spatial part of the input is given by a normalized Gaussian or exponential kernel, cf. (11) and (12), the Ce –Te parameter-planes are seen to be more varied (and in this sense qualitatively similar to the cases using the broad-bump stimulus for slightly smaller inhibitory time constants τ depicted in Fig. 10). Fig. 12(b), (c) and (e), (f) show the corresponding temporal and scaled temporal part of the external input that correspond to the minimal Te (Fig. 12(b) and (e)) and minimal Ce (Fig. 12(c) and (f)), respectively. The minimal Ce –Te values corresponding to these four exemplary spatial profiles are summarized in Table 1. Note, that

Table 1 Minimal values of duration time for lowest amplitude of excitatory input required to generate a persistent bump state for different spatial functions. Peak time fixed at Te1 = 2. Other parameters as in Figs. 11 and 12. Parameters

Broad BP

Narrow BP

Gaussian

Exponential

Ce Te

3.05 7.83

5.42 8.73

0.202 5.04

0.404 7.56

the observed much smaller minimal values of Ce for Gaussian or exponential spatial kernels expectedly simply reflect that the maximal values of Se , due to the choice of normalization, are more than a factor ten larger for these functions compared to the broad-BP function. 4.2. Annihilation of persistent activity Above we investigated the process of generation of persistent bump states and the dependence of this process on the specifics of the spatiotemporal excitatory input functions He . Such activation of persistent states of neuronal populations is one key element of working memory, i.e., the capacity of the system to actively hold a memory item. However, after a while this item might have lost relevance, and the persistent activation needs to be switched off again. Thus, now we investigate what kind of transient external input can annihilate persistent bump states. A first candidate annihilation mechanism is inhibitory input onto the excitatory population, i.e., He (x, t ) < 0. We indeed found that such input easily annihilates persistent bump states in our two-population model (results not shown). However, as discussed in Compte (2006), this mechanism should result in a decreased firing, whereas in experiments an increase in firing rates is actually

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

85

Fig. 9. Generation of persistent bump states as a function of amplitude Ce and peak time Te1 of the excitatory external input for different Te values for the triangular temporal function (13). Green color represents failure of activation, red color successful activation. Rows one to four represent the Ce –Te1 -planes for increasing inhibitory time constant τ for four different durations Te of the external input He (x, t ) each. The spatial part Se (x) of the external input is again chosen as the broad BP given by (9). Remaining parameters as in Fig. 2. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 10. Dependence of persistent bump-states plane (Ce –Te -plane) for three different values of the relative inhibitory time constant τ in the transition range 2 . τ . 2.5 where the system alters its qualitative behavior (see text for details), cf. Fig. 7.

86

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

Fig. 11. Generation of persistent bump states as a function of amplitude Ce and duration time Te of the excitatory external input (8) using narrow and broad bumps (9) as spatial kernels Se (x), respectively. (a) Broad bump, (c) narrow bump. (b) and (d) show the corresponding temporal kernels for the smallest amplitude and shortest duration time sufficient to generate the persistent bump state. The peak time is assumed to be Te1 = 2, τ = 2.5, while the rest of the parameters are the same as used in Fig. 2.

Fig. 12. Generation of persistent bump states as a function of amplitude Ce and duration time Te of the excitatory external input (8) for Gaussian and exponentially decaying (12) functions as spatial kernels. (a) Gaussian function (11) used as spatial function. (b) and (c) show the temporal kernels for the smallest duration time Te and amplitude Ce sufficient to evoke persistent activity. (d) Exponentially decaying function (12) used as spatial function. (e) and (f) show the temporal kernels for the smallest duration time Te and amplitude Ce sufficient to evoke persistent activity. ρm = 0.6, other parameters as in Fig. 11.

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

87

Fig. 13. Illustration of the generation and annihilation of persistent bump states by transient inputs. Triangular temporal function (13) and broad-BP spatial input (9) are used. (a) Excitatory external input onto the excitatory population (He (x, t ) > 0) generating the persistent bump state starting at t = 0. (b) Excitatory external input onto the inhibitory population (Hi (x, t ) > 0) annihilating the persistent bump state starting at t = 30. (c) Excitatory population activity ue (x, t ). (d) Inhibitory population activity ui (x, t ). (e) Corresponding triangular temporal functions of the external input and (f) pulse-width-coordinate evolution. Parameters of transient input are te0 = 0, Te1 = 2, Te = 8, Ce = 3.38, ti0 = 30, Ti1 = 32, Ti = 36, Ci = 2. Further, τ = 2.5, remaining parameters as in Fig. 2.

observed. This rather suggests that the erasure signal is indeed excitatory, and this is what we will focus on here. Again there are two possible mechanisms. The first is excitatory input onto the inhibitory population (which in turn inhibits the excitatory population). Our results above (Fig. 7) suggest a second mechanism, namely that excitatory input onto the excitatory population also can annihilate persistent bump states in cases where there is an upper critical amplitude to evoke activation, such as in Fig. 7(g). It should, however, not be possible for cases such as depicted in Fig. 7(a)–(f). This we indeed find to be the case below. 4.2.1. Excitatory external input onto inhibitory population First we focus on the memory-erasing effect of excitatory external input onto the inhibitory population only, i.e. Hi (x, t ) > 0, He (x, t ) ≡ 0 for t ≥ ti0 . Fig. 13 shows an example of the generation and annihilation of persistent bump states by means of excitatory input on the excitatory (Fig. 13(a)) and inhibitory populations (Fig. 13(b)), respectively. The triangular temporal function he (t ) (13) (Fig. 13(e)) and the broad-BP excitatory activity profile are used in Sm (x). The corresponding triangular pulse-width coordinate evolution are illustrated in Fig. 13(f). Annihilation parameter-planes (Ci –Ti , Ci –Ti1 ) for two different values of the relative inhibitory time constant τ are shown in Fig. 14: τ = 1.5 in the top row (panels a, b), τ = 2.5 in the bottom row (panels c, d). As to be expected, erasure of the persistent state is easier for larger τ , i.e., values of τ closer to the critical values τcr for which the persistent bump state is no longer stable. This can

be seen from the smaller area in the parameter-plane for which annihilation fails (red). 4.2.2. Excitatory external input onto excitatory population We next investigate whether successful annihilation of persistent activation also can be attained by excitatory input onto excitatory population. A priori this seems conceivable in the situations observed above where ‘too much’ excitation prevented the generation of a persistent state, e.g., for τ = 2.5 and τ = 3 in Fig. 7. In Fig. 15 we indeed observe that for τ = 2.5 there is a region in the Ce –Te –Te1 space where the system can be successfully annihilated (green), if the amplitude Ce is large enough. For τ = 3, i.e., close to τcr = 3.03, almost no excitatory input is required to annihilate the state. However, for, e.g., τ = 2, such annihilation was found to be impossible, regardless of the (positive) value of Ce . Annihilation by excitatory input to the excitatory population is hence more dependent on the relative inhibitory time constant τ than annihilation due to such input onto the inhibitory population, and can thus be considered a more selective erasure mechanism. 5. Conclusions and discussion The selective persistent activation of populations of neurons in the prefrontal cortex during active memory tasks is one of the best studied neural correlates of a higher cognitive function, the so-called working memory (Goldman-Rakic, 1995). A possible

88

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

Fig. 14. Amplitude Ci versus stimulus duration Ti and the peak time Ti1 for excitatory (Hi (x, t ) > 0) external inputs to the inhibitory population to erase the stable activity for the triangular temporal function (13) for τ = 1.5, 2.5 and τ = 3. Green color indicates successful deactivation (final state of no activity), red failure of deactivation (resulting state of persistent activity). (a) Ci − Ti plane for Ti1 = 3, and (b) Ci − Ti1 plane for fixed duration Ti = 10 with τ = 1.5, while the second row ((c) and (d)) corresponds to τ = 2.5. Note, that ti0 = 30 here, thus the x-axis starts at 30. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

mechanism for such selective activation and deactivation of subpopulations of neural networks is the switching between attractor states in networks, which were found both in spiking neuron (e.g., Brunel, 2003; Camperi & Wang, 1998; Compte, 2006; Compte et al., 2000, 2003; Gutkin et al., 2001; Laing & Chow, 2001; Renart et al., 2003) and neural-field models (e.g., Amari, 1977; Coombes, 2005; Laing & Chow, 2001; Laing et al., 2002; Machens & Brody, 2008; Wilson & Cowan, 1973), by means of external input. Here, we have investigated the effect of spatially dependent and temporally transient external input on bump-pair formation in a two-population Wilson–Cowan type neural-field model, where the firing-rate functions are approximated by means of a unit step function (Heaviside function). This is an extension of our previous work on the same model with temporally constant external input (Yousaf et al., 2013). The results can be summarized as follows: First, we have rigorously shown that the solutions of the initial value problem (2) are bounded for finite external input. These boundedness results also hold true for the multidimensional case and for the general case of firing-rate functions possessing values between zero and one. We moreover derived the general necessary condition for the emergence of persistent activity in response to transient external input for the system (2): A persistent activation of the network can be evoked, only if there is some point (xcr , tcr ) in space and time, such that the excitatory population activity is at threshold, i.e., ue (xcr , tcr ) = θe , and the external input at that point is greater than that threshold value, i.e., He (xcr , tcr ) > θe . We then numerically investigated the emergence and annihilation of persistent activity in more detail for different concrete examples of spatiotemporal external input functions Hm > 0 (m = e, i). In particular, we assumed a factorized form, i.e., Hm (x, t ) = Cm Sm (x) hm (t − tm0 ), that made it possible to study the respective roles of spatial shape and time-course of the external input separately. Initially, if there is no activity in the network and

Hm (x, t ) = 0, the addition of suitable non-zero excitatory external input He (x, t ) > 0 satisfying Theorem 3 can evoke persistent activation that remains stable even if the external input is switched off. For all considered input functions, we numerically derived the amplitude Ce and duration Te (Ce –Te -plane) of the external input required to successfully evoke persistent activity, i.e., push the system to the stable broad bump-pair attractor. As expected, in each case the amplitude Ce and the stimulus duration Te must exceed certain threshold values, that vary for the different choices of spatial Se (x) and temporal functions he (t ) of the external input. We found that the relative inhibitory time constant τ plays a pivotal role in shaping the Ce –Te -plane. If it is much smaller than the critical τcr beyond which the broad bump becomes unstable (Blomquist et al., 2005), activation in terms of Ce and Te is an ‘allor-nothing’ effect. As long as the total resulting amplitude of He is sufficiently large, the attractor state is successfully activated. In this small-τ regime we also found that stimulation of the excitatory population is not a means to deactivate the system, i.e., to erase ongoing working memory. This can then only be achieved by stimulating the inhibitory population such that Hi > 0. The picture becomes however very different if τ becomes larger. In this case it can be that, e.g., too short stimuli cannot activate the system for any magnitude of the amplitude, or that there is only a limited band of amplitudes that lead to successful activation. If the amplitude parameter Ce becomes too large for small duration Te , the system overshoots and, instead of settling into the attractor, leaves the basin of attraction and decays to the trivial ‘no-bump’ attractor. Consequently, also the onset-slope of the stimulus matters in this regime. In this intermediate-τ regime we moreover observed that not only an excitatory stimulation of the inhibitory population, but also stimulation of the excitatory population can serve as an erasure signal. In the model studied

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

89

Fig. 15. Amplitude Ce versus stimulus duration Te and the peak time Te1 for excitatory external inputs to erase the activity for a triangular temporal function (13) and the broad BP as a spatial function Se for τ = 2.5. (a) Ce − Te plane for Te1 = 3, and (b) Ce − Te1 plane for fixed duration Te = 10. Again, te0 = 30, thus the x-axis starts at 30.

here, this is a consequence of the trade-off between recruitment of excitation and inhibition. If too much excitation gets activated, it will in turn activate more lateral inhibition that can eventually lead to a ceasing of activity. It was shown that for spiking network models (Gutkin et al., 2001; Laing & Chow, 2001) excitation of the excitatory population can also lead to annihilation of bump activity. Though we expect that in principle recurrent inhibitory feedback can work as an erasure signal in spiking neuron network realizations of the architecture we analyze here as well, other mechanisms, such as oscillatory instabilities (Laing & Chow, 2001), or transient synchronization leading to too many neurons being in the refractory state (Gutkin et al., 2001; Kumar, Schrader, Aertsen, & Rotter, 2008), as well as shunting inhibition (Gutkin et al., 2001) were shown to be efficient too, irrespective of strong inhibitory feedback. The much stricter constraints external stimuli need to fulfill to persistently activate the system, and the greater flexibility in annihilating this activation if τ is closer to τcr could be useful features from a neural computation perspective, because more selectivity in terms of which stimulus aspects are relevant also makes individual stimuli more informative. Numerical investigation of the transition from the monotonic ‘all-or-nothing’ to the non-monotonic complex Ce –Te –Te1 -dependence of generation and annihilation of the attractor state revealed that it is gradual and smooth. Only a slight change of the relation of the effective inhibitory and excitatory population time constants expressed by τ , for example due to learning or neurosecretion, could thus make the network respond more or less selectively to external stimuli. Future work should be directed to gain a better formal understanding of the general mechanism that determines the activation and deactivation of the attractor state in the two-population Wilson–Cowan-type model. For example, the necessary condition for activation of the attractor state presented here only gives a lower bound for the necessary stimulus size. It is important to identify stimulus classes Hm (x, t ) and firing-rate functions Pm (u(x, t )) such that it is possible also to derive the sufficient condition for the activation of the attractor state. Acknowledgments The authors would like to thank, Dr. Tom Tetzlaff and Dr. Hans Ekkehard Plesser (Norwegian University of Life Sciences), for many fruitful and stimulating discussions during the preparation of this paper. This research has been supported by the Research Council of Norway under the grant 178892 (eScience[eNeuro]) and the BrainScaleS project (269921) in EU FP7. M. Yousaf would like to thank the Higher Education Commission of Pakistan for support. Appendix. Proof of Theorem 2 Proof. The proof proceeds by stepping forward in time. By continuity we have that um (x, t ) < θm for m = e, i on some interval

0 ≤ t ≤ t1 (t1 < T ∗ ). The non-local terms ωmn ⊗ Pm (um − θm ) in system (2) do not give any contribution on this time interval and hence the solution in this case is given by the expressions (31). These expressions can be extended to t = t1 . Thus from Theorem 1, we have that ue (x, t1 ) < θe and ui (x, t1 ) < θi . We then repeat the argument for ue (x, t1 ) and ui (x, t1 ) as initial conditions, continuing in time implies that ue (x, t ) < θe and ui (x, t ) < θi for some time interval t1 ≤ t < t2 (with t2 < T ∗ ). Hence, there will also be no contribution from non-local terms on that time interval from which it follows that the solution of (2) is given as Ue (x, t ) = αe (t − t1 )Ve (x) +

t



αe (t − s)He (x, s)ds t1

Ui (x, t ) = αi (t − t1 )Vi (x) +



t

αi (t − s)Hi (x, s)ds t1

on the interval t1 ≤ t < t2 . Again by Theorem 1, we can extend the solution to t = t2 and we get um (x, t2 ) < θm for m = e, i. By proceeding in a similar way, um (x, t2 ) can be assumed as initial conditions for the system (2) and we get the solution on the time interval t2 ≤ t < t3 with the property um (x, t ) < θm for m = e, i. Hence, the time interval 0 ≤ t < T ∗ can be divided into closed bounded subintervals [tn , tn+1 ], with n = 0, 1, 2, 3, . . . , N − 1, t0 ≡ 0 and tN = T ∗ . The solution corresponding to these time intervals is given by (31). Moreover the solutions are continuously glued to together at each division point tn . Thus, the solutions satisfy um (x, t ) < θm for 0 ≤ t ≤ T ∗ . 

References Amari, S. (1977). Dynamics of pattern formation in lateral-inhibition type neural fields. Biological Cybernetics, 27, 77–87. Blomquist, P., Wyller, J., & Einevoll, G. T. (2005). Localized activity patterns in two population neuronal networks. Physica D, 206, 180–212. Brunel, N. (2003). Dynamics and plasticity of stimulus-selective persistent activity in cortical network models. Cerebral Cortex, 13, 1151–1161. Camperi, M., & Wang, X.-J. (1998). A model of visuospatial working memory in prefrontal cortex: recurrent network and cellular bistability. Journal of Computational Neuroscience, 5, 383–405. Compte, A. (2006). Computational and in vitro studies of persistent activity: edging towards cellular and synaptic mechanisms of working memory. Neuroscience, 139, 135–151. Compte, A., Brunel, N., Goldman-Rakic, P. S., & Wang, X. J. (2000). Synaptic mechanisms and network dynamics underlying spatial working memory in a cortical network model. Cerebral Cortex, 10, 910–923. Compte, A., Constantinidis, C. T. J., Raghavachari, S., Chafee, M. V., Goldman-Rakic, P. S., & Wang, X, J. (2003). Temporally irregular mnemonic persistent activity in prefrontal neurons of monkeys during a delayed response task. Journal of Neurophysiology, 90, 3441–3454. Coombes, S. (2005). Waves, bumps, and patterns in neural field theories. Biological Cybernetics, 93, 91–108. Cressman, J. R., Ullah, G., Ziburkus, J., Schiff, S. J., & Barreto, E. (2009). The influence of sodium and potassium dynamics on excitability, seizures and the stability of persistent states: I. Single neuron dynamics. Journal of Computational Neuroscience, 26, 159–170.

90

M. Yousaf et al. / Neural Networks 46 (2013) 75–90

Ermentrout, B., & Drover, J. D. (2003). Nonlinear coupling near a degenerate Hopf (Bautin) bifurcation. SIAM Journal on Applied Mathematics, 63(5), 1627–1647. Ermentrout, G. B., & McLeod, J. B. (1993). Existence and uniqueness of travelling waves for a neural network. Proceedings of the Royal Society of Edinburgh Section A, 123, 461–478. Folias, S. E., & Bressloff, P. C. (2004). Breathing pulses in an excitatory neural network. SIAM Journal on Applied Dynamical Systems, 3, 378–407. Goldman-Rakic, P. S. (1995). Cellular basis of working memory. Neuron, 14, 477–485. Gray, R. T., & Robinson, P. A. (2009). Stability and structural constraints of random brain networks with excitatory and inhibitory neural populations. Journal of Computational Neuroscience, 27, 81–101. Guo, Y., & Chow, C. (2005). Existence and stability of standing pulses in neural networks: II. Stability. SIAM Journal on Applied Dynamical Systems, 4, 249–281. Gutkin, B. S., Laing, C. R., Colby, C. L., Chow, C. C., & Ermentrout, G. B. (2001). Turning on and off with excitation: the role of spike-timing asynchrony and synchrony in sustained neural activity. Journal of Computational Neuroscience, 11, 121–134. Kilpatrick, Z. P., & Bressloff, P. C. (2010a). Effects of synaptic depression and adaptation on spatio-temporal dynamics of an excitatory neuronal network. Physica D, 239, 547–560. Kilpatrick, Z. P., & Bressloff, P. C. (2010b). Stability of bumps in piecewise smooth neural fields with nonlinear adaptation. Physica D, 239, 1048–1060. Kishimoto, K., & Amari, S. (1979). Existence and stability of local excitations in homogenious neural fields. Journal of Mathematical Biology, 7, 303–318. Kumar, A., Schrader, S., Aertsen, A., & Rotter, S. (2008). The high-conductance state of cortical networks. Neural Computation, 20, 1–43. Laing, C. R., & Chow, C. C. (2001). Stationary bumps in networks of spiking neurons. Neural Computation, 13, 1473–1494. Laing, C. R., Troy, W. C., Gutkin, B., & Ermentrout, G. B. (2002). Multiple bumps in a neuronal model of working memory. SIAM Journal on Applied Mathematics, 63, 62–97. Machens, C., & Brody, C. (2008). Design of continuous attractor networks with monotonic tuning using a symmetry principle. Neural Computation, 20, 452–485.

Marti, D., & Rinzel, J. (2013). Dynamics of feature categorization. Neural Computation, 25, 1–45. Oleynik, A., Wyller, J., Tetzlaff, T., & Einevoll, G. T. (2011). Stability of bumps in a twopopulation neural-field model with quasi-power temporal kernels. Nonlinear Analysis: Real World Applications, 12, 3073–3094. Pinto, D. J., & Ermentrout, G. B. (2001). Spatially structured activity in synaptically coupled neuronal networks: II. Lateral inhibition and standing pulses. SIAM Journal on Applied Mathematics, 62, 226–243. Pinto, D. J., Hartings, J. A., Brumberg, J. C., & Simons, D. J. (2003). Cortical damping: analysis of thalamocortical response transformations in rodent barrel cortex. Cerebral Cortex, 33–44. Potthast, R., & Beim Graben, P. (2009). Existence and properties of solutions for neural field equations. Mathematical Methods in the Applied Sciences, 33(8), 935–949. Renart, A., Brunel, N., & Wang, X. J. (2003). Mean-field theory of irregularly spiking neuronal populations and working memory in recurrent cortical networks. In Computational neuroscience: a comprehensive approach (pp. 431–490). Boca Raton FD, USA: CRC Press (Chapter 15). Rubin, J., & Bose, A. (2004). Localized activity patterns in excitatory neuronal networks. Network: Computation in Neural Systems, 15, 133–158. Rubin, J. E., & Troy, W. C. (2004). Sustained spatial patterns of activity in neuronal populations without recurrent excitation. SIAM Journal on Applied Mathematics, 64, 1609–1635. Wang, X. J. (2001). Synaptic reverberation underlying mnemonic persistent activity. Trends in Neurosciences, 24, 455–463. Werner, H., & Richter, T. (2001). Circular stationary solutions in two-dimensional neural fields. Biological Cybernetics, 85, 211–217. Wilson, H. R., & Cowan, J. D. (1973). A mathematical theory of the functional dynamics of cortical and thalamic nervous tissue. Kybernetik, 13, 55–80. Wyller, J., Blomquist, P., & Einevoll, G. T. (2007). During instability and pattern formation in a two-population neuronal network model. Physica D, 22, 75–93. Yousaf, M., Wyller, J., Tetzlaff, T., & Einevoll, G. T. (2013). Effect of localized input on bump solutions in a two-population neural-field model. Nonlinear Analysis: Real World Applications, 14, 997–1025.