Noise-induced canard and mixed-mode oscillations in large stochastic networks with multiple timescales
NNOISE-INDUCED CANARD AND MIXED-MODE OSCILLATIONSIN LARGE STOCHASTIC NETWORKS WITH MULTIPLETIMESCALES
JONATHAN TOUBOUL ‡§∗ , MARTIN KRUPA † , AND
MATHIEU DESROCHES † Abstract.
We investigate the dynamics of large stochastic networks with different timescalesand nonlinear mean-field interactions. After deriving the limit equations for a general class of networkmodels, we apply our results to the celebrated Wilson-Cowan system with two populations with orwithout slow adaptation, paradigmatic example of nonlinear mean-field network. This system hasthe property that the dynamics of the mean of the solution exactly satisfies an ODE. This reductionallows to show that in the mean-field limit and in multiple populations with multiple timescales, noiseinduces canard explosions and Mixed-Mode Oscillations on the mean of the solution. This sheds newlight on the qualitative effects of noise and sensitivity to precise noise values in large stochasticnetworks. We further investigate finite-sized networks and show that systematic differences with themean-field limits arise in bistable regimes (where random switches between different attractors occur)or in mixed-mode oscillations, were the finite-size effects induce early jumps due to the sensitivity ofthe attractor.
Key words. mean field equations, neural mass models, bifurcations, noise, dynamical systems
AMS subject classifications.
Introduction.
At functional scales, cortical networks are composed of a largenumber of cells subject to an intense noise. Such networks process and encode in-formation collectively, and in the variety of information transmitted, synchronousactivity and periodic behaviors are ubiquitous [18, 10], both in physiological (as-sociated to memory and attention) and pathological (epilepsy) states. Oscillatorypatterns observed can be extremely complex. Theoretical descriptions of such com-plex oscillations have been satisfactorily performed in deterministic models of singleneurons, and are generally attributed to the presence of several very different timescales in the dynamics. The analysis of these systems shed light on important andsensitive dynamical scenarios such as canard explosions (see for example [20]). Thequestion we address in this work is whether such oscillatory patterns can arise in thecontext of large-scale neuronal aggregates, including a very large number of neurons,strongly interconnected and subject to noise. In other words, we will be interested incomplex oscillatory patterns that can appear in large scale neuronal assemblies in thepresence of noise. As was the case of deterministic small dimension model, we willspecifically focus on the canard phenomenon and the related phenomenon of
MixedMode Oscillations (MMOs).Canard phenomenon, known also as canard explosion is a transition from a small,Hopf type oscillation to a relaxation oscillation, occurring upon variation of a parame-ter. This transition has been first found in the context of the van der Pol equation [4]and soon after in numerous models of phenomena occurring in engineering and inchemical reactions [9]. A common feature of all these models is the presence of timescale separation (slow and fast variables). A particular feature of canard explosion is ‡ The Mathematical Neurosciences Laboratory, Center for Interdisciplinary Research in Biology(CNRS UMR 7241, INSERM U1050, UPMC ED 158, MEMOLIFE PSL*) Coll`ege de France, 11place Marcelin Berthelot, 75005 Paris, § BANG Laboratory, Inria Paris-Rocquencourt. ∗ [email protected] † SISYPHE Laboratory, Inria Paris-Rocquencourt.1 a r X i v : . [ m a t h . D S ] F e b M. DESROCHES, M. KRUPA & J. TOUBOUL that it takes place in a very small parameter interval. For the van der Pol system,where the ratio of the timescales is given by a small parameter ε , the width of thisparameter interval can be estimated by O (exp( − c/ε ), where c > canard cycles whose strikingfeature is that they contain segments following closely a stable slow manifold andsubsequently an unstable slow manifold.The work on canard explosion led to investigations of slow/fast systems in threedimensions, with two slow and one fast variables. In the context of these systemsthe concept of a canard solution or simply canard has been introduced, as a solutionpassing from a stable to an unstable slow manifold [3, 23, 24]. Canards arise near socalled folded singularities of which the best known is folded node [3, 24, 30]. Unlike insystems with one slow variable canards occur robustly (in O (1) parameter regions).Related to canards are Mixed Mode Oscillations (MMOs), which are trajectories thatcombine small oscillations and large oscillations of relaxation type, both recurring inan alternating manner. Recently there has been a lot of interest in MMOs that arisedue to the presence of canards, starting with the work of Milik, Szmolyan, Loeffelmannand Groeller [23]. The small oscillations arise during the passage of the trajectoriesnear a folded node or a related folded singularity. The dynamics near the foldedsingularity is transient, yet recurrent: the trajectories return to the neighborhood ofthe folded singularity by way of a global return mechanism [8].Due to the sensitivity of this phenomena, one may expect that canard explosionsdisappear in the presence of noise. This is implied by the results of [5], who studya system with small noise. Using a pathwise approach, they show that small noiseperturbations dramatically disturb the behavior of the slow-fast system they consider,hiding small oscillations in the random fluctuation and having a positive probabilityto escape from the canard solution. In this picture, it seems relatively hopeless toobtain at the levels of noise typical of the cortical activity, and at the dimension ofthe system, collective canard phenomena and individual canards with a high degree ofreproducibility and regularity. However, in large scale stochastic systems, averagingeffects can take place, regularizing the probability distribution of the macroscopicactivity, as well as that of individual cells. These macroscopic behaviors and theaveraging effects arising at the level of large networks are mathematically characterizedby limit systems, the mean-field limits . In the present manuscript we show that canardexplosion and complex oscillatory pattern appear in the presence of (not necessarilysmall) noise in these mean-field limits. Moreover, we show that such phenomena mayarise through noise induced bifurcations.Usual computational neurosciences descriptions of large scale networks rely onheuristic models which were derived in the seminal works of Wilson, Cowan andAmari [1, 2, 31, 32]. These models implicitly make the assumption that noise cancelsout through averaging effects in the asymptotic regime where the number of neuronsgoes to infinity. The analysis of these heuristic models successfully accounted for anumber of cortical phenomena such as spatio-temporal pattern formation in spatiallyextended models [13, 14, 7].Indeed, in order to uncover, in large scale noisy system, the presence of fine effectssuch as complex oscillatory patterns or canard explosions, one needs to very preciselydescribe the limits of the network activity as the number of neurons increases. Thisis why we will apply a recent developed approach based on the theory of interactingparticle systems to our networks. This approach was used in a very general context in-cluding nonlinear models of spiking neurons, spatial extension and propagation delays tochastic Networks with multiple timescales
1. Network models and their mean field limits.
This section introduces ageneral model of multi-population neuronal network and derives the related macro-scopic mean-field equations arising as a limit when the network size goes to infinity.
Models.
We are interested in the large scale behavior arising from the nonlinearcoupling of a large number N of stochastic diffusion processes representing the mem-brane potential of neurons. Each neuron belongs to a specific population α ∈ { · · · P } ,characterizing the dynamics of neurons and their time constants, and the interactionwith other neurons. Each population is composed of a large number of neurons de-noted ( N α , α = 1 · · · P ). The state of neuron i in population α is described by avariable X it ∈ E = R d typically describing the membrane potential of the cell anddifferent ionic concentrations. The intrinsic dynamics of this variable (disregarding M. DESROCHES, M. KRUPA & J. TOUBOUL input and interactions with other neurons) is governed by a nonlinear function f α anda time constant τ α .The action of neurons of population β onto neurons of population α will be mod-eled through a generalized version of Wilson-Cowan system, as a nonlinear sigmoidaltransformation of a the average of the activity of all cells, weighted by the typicalsynaptic connectivity J αβ . Beyond network interactions, neurons receive inputs cor-responding to the sum of a deterministic current I αt only depending on the populationthey belong to, and stochastic current driven by independent stochastic process ξ it withcommon law (denoted with a slight abuse of notations ξ αt ), that will be considered tobe Ornstein-Uhlenbeck processes, taking into account the synaptic filtering. In theirmost general form, the class of models analyzed here satisfy the network equations: τ α d X it d t = f α ( X it ) + S α X it , P (cid:88) β =1 J αβ N β (cid:88) p ( j )= β X jt + I αt + ξ it . (1.1)where S α ( x, y ) represents the effect of an input current y onto a cell of population α with potential x , (cid:80) Pβ =1 J αβ N β (cid:80) p ( j )= β X jt the total current received through the net-work, and p : N (cid:55)→ { · · · P } is the population function associating to a neuron j thepopulation β it belongs to.This class of models is particularly interesting since it allows dealing with a num-ber of situations arising in neurosciences. A prominent example of this class that willbe useful in our developments is the Wilson-Cowan model, or activity-based firingrate model [31, 32]. This model describes the mean firing-rate variable of a neuron (areal variable) of neuron i in the population α through the differential equation: τ α dX it dt = − X it + S α (cid:16) P (cid:88) β =1 J αβ N β (cid:88) p ( j )= β X jt + I αt + ξ it (cid:17) (1.2)In this model, the variable X it represents the activity of the cell, and the right handside is the incoming firing-rate to the cell i , filtered through the voltage to rate function S α characterizing the integration properties of the synapse when subject to a inputcurrent (see [13]).Another model of this class allows taking into account finer phenomena arisingin the cell’s activity motivating the description of additional variables, such as slowadaptation currents [11] or slow modulation currents [17]. These models correspond tominimal descriptions of systems with slow adaptation proposed by Izhikevich in [17]:in two-dimensional systems, a variety of behavior can be accounted for using a singleslow adaptation. We slightly generalize this model and define a multi-populationnetworks with mean-field slow modulation current. These systems are described bythe following equations, for neuron i in population α : τ α d X it d t = − X it + S α ( (cid:80) Pβ =1 J αβ N β (cid:80) p ( j )= β X jt + λ α U it + I αt + ξ i, t ) d U it d t = ε α ( k + γU it − (cid:80) Pβ =1 1 N β (cid:80) p ( j )= β X jt ) (1.3)Beyond these simple models that will be exactly reducible to ODEs, the classof models includes classical models in neurosciences such as the Fitzhugh-Nagumomodel ( X ∈ R , f α has a cubic and an affine component, and S ( x, y ) is either a linearfunction when considering electrical synapses or a sigmoid transform when considering tochastic Networks with multiple timescales S . We now derive the limitsof these stochastic networks equations as the number of neurons go to infinity. Macroscopic Limits.
Based on the general network given by equations (1.1), wederive the macroscopic mean-field limit under a few assumptions on the dynamics. Forthe sake of simplicity, we denote by f ( x ) : E P (cid:55)→ E P (resp. S ( x, y ) : E P × E P (cid:55)→ E P )the function with component α equal to f α ( x α ) (resp. S α ( x α , y α )). We make thefollowing assumptions:(H1). f is K -Lipschitz continuous(H2). f satisfies the linear growth condition: f ( x ) ≤ C (1 + | x | ).(H3). S is L -Lipschitz continuous in both variables(H4). S is bounded with supremum denoted (cid:107) S (cid:107) ∞ . Remarks. • Rigorously, in the Fitzhugh-Nagumo or Hodgkin-Huxley system, f is not globallyLipschitz continuous and does not satisfies a linear growth condition. However, thesefunctions are locally Lipschitz continuous, and satisfy the estimate x T f ( x ) ≤ K (1 + | x | ) . These two assumptions allow to demonstrate the same kind of properties,by using Itˆo formula and truncation method. This was done in [28] in a moregeneral context, and the interested reader may readily apply this generalization toour system. • In these assumptions, we consider without loss of generality the same constants K for the sake of simplicity. • All these assumptions are obviously satisfied for the Wilson-Cowan system.
We will show that in the mean-field limit, the propagation of chaos applies andthat the state of all neurons in that limit satisfy the implicit equation:d X αt d t = f α ( X αt ) + S α X αt , P (cid:88) β =1 J αβ E (cid:104) X βt (cid:105) + I αt + ξ αt . (1.4)which can be written in vector form:d X t d t = f ( X t ) + S ( X t , J · E [ X t ] + I t + ξ t ) (1.5)where J = ( J αβ ) α,β ∈{ ··· P } ∈ R P × P and I t = ( I αt ) α =1 ··· P , ξ t = ( ξ αt , α ∈ { · · · P } ).This equation is not a usual stochastic differential equation in that in involvesthe expectation of the solution. It is rather an implicit equation on the space ofprobability distributions. The first question that may arise in this type of unusualstochastic equation is its the well-posedness. This is subject to the following: Proposition 1.1.
There exists a unique, square integrable, strong solution tothe equations (1.4) starting from a square integrable initial condition.Proof . The existence and uniqueness of solutions is demonstrated using a routinecontraction argument. Let us fix a square integrable initial condition random variable X ∈ E P and define the map Φ acting on the space of square integrable stochasticprocesses M ([0 , T ] , E P ):Φ : X (cid:55)→ (cid:18) X + (cid:90) t (cid:16) f ( X s ) + S (cid:0) X s , J · E [ X s ] + I s + ξ s (cid:1)(cid:17) ds (cid:19) t ∈ [0 ,T ] . M. DESROCHES, M. KRUPA & J. TOUBOUL
Fixed points of Φ are exactly the solutions of the mean-field equations. Moreover, for X a square integrable, we have Φ( X ) ∈ M ([0 , T ] , E P ). Indeed, thanks to Cauchy-Schwarz inequality, we have: E (cid:34) sup t ∈ [0 ,T ] | Φ( X ) t | (cid:35) ≤ (cid:26) E (cid:2) | X | (cid:3) + T (cid:90) T C (1 + E (cid:2) | X s | (cid:3) ) + (cid:107) S (cid:107) ∞ ds (cid:27) < ∞ . Fixing a square integrable process X ∈ M ([0 , T ] , E P ), we define the sequence ofsquare integrable processes X n +1 = Φ( X n ). We note D nt = E (cid:104) sup t ∈ [0 ,T ] | X n +1 t − X nt | (cid:105) , (cid:107) J (cid:107) the operator norm of J and Z kt = J · E [ X s ] + I s + ξ s . We have: D nt ≤ T E (cid:40) sup s ∈ [0 ,t ] (cid:90) s | f ( X nu ) − f ( X n − u ) | + | S ( X nu , Z nu ) − S ( X n − u , Z nu ) | + | S ( X n − u , Z nu ) − S ( X n − u , Z n − u ) | du (cid:41) ≤ T E (cid:40) sup s ∈ [0 ,t ] (cid:90) s ( K + L ) | X nu − X n − u | + L | J · E (cid:2) X nu − X n − u (cid:3) | du (cid:41) ≤ T E (cid:40) sup s ∈ [0 ,t ] (cid:90) s ( K + L ) | X nu − X n − u | + L (cid:107) J (cid:107) E (cid:2) | X nu − X n − u | (cid:3) du (cid:41) ≤ T ( K + L (1 + (cid:107) J (cid:107) )) (cid:90) t D n − s ds We hence obtain by immediate recursion: D nt ≤ ( K (cid:48) ) n T n n ! D T with K (cid:48) = 3 T ( K + L (1+ (cid:107) J (cid:107) )). D T is a finite quantity because of the boundedness of E (cid:2) | Φ( X ) | (cid:3) . From this inequality, routine methods using the Bienaym´e-Chebyshevinequality allow concluding on the existence and uniqueness of solutions of the mean-field equation (see e.g. [19]).The next theorem proves that the unique solution to the mean-field equations willdescribe exactly the behavior of the network when the network size goes to infinity(i.e. when all N α → ∞ ), and that all neurons will be independent provided thattheir initial condition is. This latter property is referred as the propagation of chaosproperty . Theorem 1.2.
If the initial conditions to the network equations (1.1) are inde-pendent, and identically distributed population per population, then in the mean-fieldlimit, the propagation of chaos applies and the laws of the network solutions convergetowards that of the mean-field equations (1.4) , in the sense that for any i ∈ N thereexists ¯ X it with law given by the solution of the mean-field equations and such that: E (cid:34) sup t ∈ [0 ,T ] | X it − ¯ X it | (cid:35) ≤ C ( T )min α √ N α Proof . The proof of the theorem uses the usual coupling argument popularized bySznitman [26] consisting in identifying an almost sure limit for the network equations. tochastic Networks with multiple timescales i of population α , we define the coupled process ¯ X it solution of equation (1.4) with initial condition X i identical to the initial conditionof neuron i in the finite network and driven by the same noise ξ it : (cid:40) d ¯ X it d t = f α ( ¯ X it ) + S α (cid:16) ¯ X it , (cid:80) Pβ =1 J αβ E (cid:104) ¯ X βt (cid:105) + I αt + ξ it (cid:17) ¯ X i = X i We denote by ˜ J the maximal value of the coupling strengths J αβ . We aim at show-ing that X i almost surely converges towards ¯ X i , and to this purpose analyze thedistance between the two processes, D αt = E (cid:104) sup t ∈ [0 ,T ] | X i,Nt − ¯ X it | (cid:105) . Because of in-terchangeability in the network, that quantity does not depend on the specific neuron i we consider but only on the neuronal population α it belongs to. For the sake ofsimplicity, we note: ¯ Z αt = (cid:80) Pβ =1 J αβ E (cid:104) ¯ X βt (cid:105) + I αt + ξ it ¯ Z N,αt = (cid:80) Pβ =1 J αβ N β (cid:80) p ( j )= β ¯ X jt + I αt + ξ it Z N,αt = (cid:80) Pβ =1 J αβ N β (cid:80) p ( j )= β E (cid:104) X jt (cid:105) + I αt + ξ it . In order to analyze this distance, we decompose it into the sum of four simpler termsthat are more easily controllable: | X it − ¯ X it | ≤ (cid:90) t | f α ( X is ) − f α ( ¯ X is ) | + | S ( X is , Z α,Ns ) − S ( ¯ X is , Z α,Ns ) | + | S ( ¯ X is , Z α,Ns ) − S ( ¯ X is , ¯ Z α,Ns ) | + | S ( X is , ¯ Z α,Ns ) − S ( ¯ X is , ¯ Z αs ) | ds ≤ ( K + L ) (cid:90) t | X is − ¯ X is | + L ˜ J P (cid:88) β =1 N β (cid:88) p ( j )= β | X jt − ¯ X jt | + | (cid:88) p ( j )= β ¯ X jt − E (cid:104) ¯ X jt (cid:105) | ds readily yielding: D αt ≤ (cid:90) t ( K + L ) D αs + L ˜ J P (cid:88) β =1 D βs + L ˜ J P (cid:88) β =1 N β E (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (cid:88) p ( j )= β ¯ X js − E (cid:2) ¯ X js (cid:3)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ds and hence for D t = (cid:80) Pα =1 D αt we obtain: D t ≤ (cid:90) t ( K + L + P L ˜ J ) D s + P L ˜ J P (cid:88) β =1 N β E (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (cid:88) p ( j )= β ¯ X ju − E (cid:2) ¯ X ju (cid:3)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ds (1.6)Let us now evaluate the second term of the righthand side. We have, using Cauchy-Schwarz inequality: E (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:88) j ∈ β ¯ X js − E (cid:2) ¯ X js (cid:3)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ≤ E (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:88) j ∈ β ¯ X ju − E (cid:2) ¯ X ju (cid:3)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) / The processes ( ¯ X js − E (cid:2) ¯ X js (cid:3) ) are centered and independent since they are driven byindependent processes W jt and ξ jt . Developing the square of the sum, we are hence M. DESROCHES, M. KRUPA & J. TOUBOUL left with less that N β non-zero terms which are equal the variance of ¯ X βt , namely E (cid:2) ( ¯ X js − E (cid:2) ¯ X js (cid:3) ) (cid:3) which is a bounded quantity as a result of proposition 1.1. Letus denote by C this bound. Returning to equation (1.6), we obtain: D t ≤ (cid:90) t ( K + L + P L ˜ J ) D s ds + T P L ˜ J (cid:115) C min β N β and Gronwall’s lemma implies that: D t ≤ K min β (cid:112) N β e K t with K = ( K + L + P L ˜ J ) and K = P L ˜ J √ CT . This quantity tends to zero whenall N α tend to zero (i.e. when min β N β → ∞ ), which implies the almost sure conver-gence ( X it ) towards ( ¯ X it ), and hence convergence in law of ( X it ) towards the solutionof the mean-field equations. Moreover, considering a finite number l of neurons in thenetwork ( i , · · · , i l ), it is easy to see using the above result that ( X i t , · · · , X i l t ) con-verge almost surely towards ( ¯ X i t , · · · , ¯ X i l t ) which is a set of independent processes.We therefore have propagation of chaos, i.e. that in the limit where all N β tend toinfinity, the state of these neurons are independent for all times. Remark.
Note that this proof has been performed in the absence of Brownian additivenoise in contrast with usual approaches. Yet, the noisy input received by each cell allows torecover the same kind of properties. We also emphasize the fact that the proof provided herereadily extends to networks with additive noise.
As announced, the results of proposition 1.1 and theorem 1.2 readily apply to theWilson-Cowan model, as we now show:
Corollary 1.3.
In the mean-field limit where all population sizes N α tend toinfinity, the Wilson-Cowan system (1.2) in the mean-field limit satisfies the propaga-tion of chaos property and the law of neurons of population α asymptotically satisfythe mean-field equations: τ α dX αt dt = − X αt + S (cid:16) P (cid:88) β =1 J αβ µ β ( t ) + I αt + ξ αt (cid:17) (1.7) where we denoted µ β ( t ) = E (cid:104) X βt (cid:105) . These quantities ( µ β ( t )) β =1 ··· P satisfy a closedequation: τ α ˙ µ α = − µ α ( t ) + G α (cid:16) P (cid:88) β =1 J αβ µ β ( t ) + I αt (cid:17) (1.8) where G α := x (cid:55)→ E [ S α ( x + ξ αt )] . (1.9) In the case of slow modulation currents (1.3) , the mean-field equations read: (cid:40) τ α d X αt d t = − X αt + S ( (cid:80) Pβ =1 J αβ µ α ( t ) + λ α U t + I αt + ξ αt ) d U t d t = ε α ( k + γU t − (cid:80) Pβ =1 µ β ( t )) (1.10) tochastic Networks with multiple timescales and in that case, U t is a deterministic variable satisfying an ODE and again, theexpectation µ α ( t ) of the process X αt satisfies a closed set of ordinary differential equa-tions: (cid:40) τ α d µ α ( t ) d t = − µ α ( t ) + G α ( (cid:80) Pβ =1 µ β ( t ) + λ α U t + I αt ) d U t d t = ε α ( k + γU t − (cid:80) Pβ =1 µ β ( t )) . (1.11) Proof . Both cases are can be readily set in the general framework of the currentsection, hence proposition 1.1 and theorem 1.2 apply, ensuring propagation of chaosand yielding the announced mean field limit. The averaged equation is trivial to derivetaking the expectation of the process.We therefore proved that the systems under consideration here enjoy a consider-able dimension reduction: the original system is a very large network of nonlinearlyinteracting stochastic processes with random currents, in which the mean is not so-lution of an ordinary differential equation, in contrast to the mean-field equationsdescribing the infinite-size limit which gives rise to a closed-form ODE for the mean.The global dynamics of the system can hence be very precisely understood byanalyzing the dynamics of the obtained system of differential equations. These involvean unspecified function, G α , which is the expectation of our sigmoidal transform S α ( x + ξ αt ) with respect to ξ α . In detail, if ξ αt has a density p αt , this function writes: G α ( x ) = (cid:90) R S ( x + y ) dp αt ( y ) , coupling the dynamics of the mean equation to the distribution of the noisy currentreceived by each cell. In what follows, we will consider that the noisy current is astationary Ornstein-Uhlenbeck process, i.e. are centered, Gaussian processes withconstant standard deviation denoted σ α . In the case where S α ( x ) = erf( g α x ), asimple change of variables and integration by parts (see e.g. [29]) gives the function G α in closed form: G α ( x ) = erf( g α x (cid:112) g α σ α ) . (1.12)The noise level σ α is hence a parameter of the vector field, and may lead to noise-driven bifurcations .From these fully determined systems, we are now in a position to use classicaltools from the analysis of slow-fast systems to uncover the presence of canards andmixed-mode oscillations in macroscopic activity in relationship to the noise level σ α .
2. Canards and MMO in the mean-field equations.
We now analyze thelimit systems (1.8) and (1.11) with a particular focus on the effects of noise on thedynamics. We start by showing the role of noise in the emergence of canard explosion,prior to extend the model to slow modulation currents systems where the additionalvariable gives rise to noise-induced MMOs.
In this section, we study the mean-field equations (1.8) arising fromWilson-Cowan system with two populations. The network equations are given by (1.2),with P = 2 populations evolving on vastly different timescales τ = ε and τ = 1, I = ze , I = 0. Moreover, in order to reduce the system to closed-form mean equa-tions, we choose S α to be erf functions and ξ α to be a stationary Ornstein-Uhlenbeck0 M. DESROCHES, M. KRUPA & J. TOUBOUL process with standard deviation σ α . As a result of the previous section, we know thatthe associated mean-field equations (1.8) have solutions whose mean satisfy a closedset of ODEs given by the equations: (cid:40) ε ˙ u = − u + G ( J u + J u + z e )˙ u = − u + G ( J u + J u ) (2.1)In these equations, the two different effective nonlinearities G and G are of theform (1.12) and therefore depend on the noise intensities σ α .Let us now analyze the dynamics of these equations. We start by analyzing thebehavior of the solutions upon variation of the parameter z e , corresponding to fixedapplied current to the first population. The behavior of the system as a function ofcurrent levels z e is particularly important in our context for its relationship with theslow modulation currents as will be made explicit in the following section. Moreover,this parameter plays the important role of moving the system from the excitableto the spiking regime. It is well known that such a parameter transition, undersome degeneracy assumptions, can lead to a canard explosion, i.e. a transition fromsmall amplitude oscillations to relaxation oscillations, see for example [21]. And inthe present system, canard explosions can be found by varying the parameter z e are observed (see bifurcation diagram generated with AUTO in Fig. 2.1(a)). Thediagrams suggest that the shown canard explosions are non-degenerate (see [21] for alist of non-degeneracy conditions). It would be straightforward but tedious to verifythese conditions and we omit such computations in this paper. This explosion informsus of the presence of canard explosions, upon variation of the input, in the Wilson-Cowan system (2.1), in the presence of noise.The question that arises is how this explosion depends on the amplitude of thenoise neurons are submitted to. This can be addressed here using usual bifurcationand singular perturbation methods since the noise appears as a parameter of theODEs, acting on the maximal slope of the sigmoids. Numerical simulation of thebifurcation diagram shows indeed that the presence of the canard explosion dependson the level of noise, and in particular that noise can induce canard explosions. Asan example of such noise induced canard explosion , we provide in Fig. 2.1(b) theAUTO bifurcation diagram of equations (2.1) upon variation of σ the noise receivedby neurons in the fast population. Therefore, noise can also bring the system fromthe region of excitability to the region of spiking, through a non-degenerate canardexplosion.As a side result, this analysis confirms the direct relationship between levels ofnoise and synchronized oscillations in the system [29, 27]. This is an extremely sur-prising phenomena for which no microscopic interpretation has been provided thusfar, and which in the context of slow-fast Wilson-Cowan system under consideration,appears even more surprising: indeed, the presence of this explosion shows that verysmall changes in the value of the variance of the noise yield dramatic qualitativeswitches on the macroscopic dynamics. Moreover, since we have proved that all neu-rons have the same probability distribution population-wise, the finding shows thatthe mean presents a sudden switch to a large-amplitude periodic orbit, correspond-ing to all neurons presenting an extreme change upon very small variation of thenoise. More than noise-induced oscillations, the phenomenon shows how sensitive thedependence upon noise levels is in large-scale systems.Now that the two-population system has been analyzed as a function of parame-ters z e and σ , we proceed by adding to the same system a slow adaptation current. tochastic Networks with multiple timescales u u C (b) (cid:1) (cid:6)(cid:1)(cid:5) (cid:1) (cid:6) (cid:1) (cid:5)(cid:1)(cid:5) (cid:1) (cid:5)(cid:2)(cid:2)(cid:1)(cid:5)(cid:3)(cid:3)(cid:1)(cid:5)(cid:4) z e u HBHB • • (a) (a) Bifurcation Diagram in z e v u HBHB • • (a) σ u u C (b) (b) Bifurcation Diagram in σ Figure 2.1 . Canard explosion in the 2D WC model, for a variation of parameter z e (toppanels) and v (bottom panels). Shown are the bifurcation branch with the two explosions (leftpanels), cycles corresponding to the canard explosion highlighted by a dashed box on the bifurcationdiagram (right panels). We now en-rich the network model analyzed in the previous section by considering slow adapta-tion currents. The underlying microscopic equation correspond to a two-populationsWilson-Cowan system with slow adaptation as introduced in equation (1.3). In orderto simplify the model, we fix λ = 1, λ = 0 and I = I = 0. The results of section 1ensure that the propagation of chaos occur and that the network converges towardsmean-field equations, given by equation (1.10), and that the mean of the solution tothe mean-field equation satisfies an equation of type (1.11), which in our two pop-ulations setting with one fast and one slow population and parameters provided insection 2.1, read: ε ˙ u = − u + G ( J u + J u + z e )˙ u = − u + G ( J u + J u )˙ z e = k + γz e − u − u , (2.2)where the slow adaptation current (noted U in the general equation (1.11)) is denotedby z e . This notation clearly shows that equation (2.2) can also be seen as general-ization of the two population Wilson-Cowan system (2.1) obtained by turning theparameter z e into a dynamic variable. In these equations, the new parameter k con-trols the position of the equilibria of (2.2). One of the most common routes to MMOsis the passage through a so-called Folded Saddle Node of type II (FSN II, see [22]).This transition is similar to canard explosion in the sense that it corresponds to apassage from a regime of a stable excitable equilibrium to MMOs or relaxation oscil-lators. We now investigate the presence of FSN II in our system (2.2) by considering2 M. DESROCHES, M. KRUPA & J. TOUBOUL
Figure 2.2 . The top left panel shows a two-parameter bifurcation diagram of the 3D WilsonCowan system (2.2) in the ( k, σ ) -plane, with Hopf (black) and period-doubling (blue) bifurcationbranches.The red curve is the locus of FSN II bifurcation. All other panels show the long-termbehavior of the system for particular values of k with fixed σ = 0 (panels (a) to (f)) and forparticular values of σ with fixed k = − (panels A to F), respectively. the reduced system, obtained by setting ε = 0:0 = − u + G ( J u + J u + z e ) (2.3)˙ u = − tu + G ( J u + J u ) (2.4)˙ z e = k + γz e − u − u (2.5)The constraint (2.3) is eliminated by solving this algebraic equation in u . Thisprovides an expression of u as a function of ( u , z e ): u = 1 J (cid:0) G − ( u ) − J u − z e (cid:1) =: F ( u , z e , σ )and G − is expressed in terms of the inverse of the erf function and of σ , dependencethat is made explicit in F . We can now rewrite (2.3)-(2.5) in the form: ∂F∂u ˙ u = − F ( u , z e , σ ) + G ( J u + J F ( u , z e , σ )) − ∂F∂z e ( k + γz e − u − F ( u , z e , σ )) (2.6)˙ z e = k + γz e − u − F ( u , z e , σ ) (2.7)Now the fold line for system (2.6)-(2.7) corresponds to the equation ∂F∂u = 0 andfolded singularities correspond to points on the fold line for which the RHS of (2.6)is 0. An FSN II point is a point on the fold line for which the RHS of both (2.6)and (2.7) are zero. Such points are of codimension 1 in the parameter space, and,under some non-degeneracy conditions, correspond to a passage of a true equilibriumof (2.6)-(2.7) through the fold line [22]. As a parameter goes through an FSN II pointthe system passes from a regime of excitability (stable equilibrium near the fold) to a tochastic Networks with multiple timescales tu (b) u u z e C F + F − • p fn MMO (a)
Figure 2.3 . Panel (a) shows the critical manifold C of system (2.2) (green surface) viewed inthe three-dimensional phase space ( u , y, z e ) . Also shown are both the fold curves F ± (black lines)and an MMO of type , solution to the problem, whose existence can be explained by the presenceof a folded node p fn (black dot) in the system. Panel (b) shows the time profile of u for this MMO. regime of MMOs or relaxation , depending on the global part of the flow [8]. On thelevel of folded singularities, FSN II corresponds to a passage from folded saddle tofolded node, and folded nodes are associated to the occurrence of MMOs [30, 8, 12].Note eventually, as demonstrated in [16], the presence of a curve of Hopf bifurcationsin the vicinity of an FSN II point.This is precisely the scenario of transition occurring in out system. Indeed, inFigure 2.2, we display in the upper left corner the FSN II curve (shown in red),computed using MATLAB based on the above formulae. Nearby (in black) is the curveof Hopf bifurcations (HB) for ε = 0 . k , is highlighted in figure 2.4. In particular, we find a cascade of period-doubling bifurcations (PD) (best seen on panel (b), which is an enlargement of panel(a)) and a family of complicated closed branches or isolas (we computed four of them).This structure is compatible with the analysis done in [16] on a minimal system of asingular Hopf bifurcation.Changes in the parameters ( k, σ ) correspond to paths in this diagrams governingthe bifurcation scenarios. Two particular paths are displayed in Figure 2.2: (i) varying k while keeping σ = 2 fixed (small letters) and (ii) fixed k = − σ varying(capital letters), and a phase space view of an MMO trajectory is given in Figure 2.3.Multiple MMO traces corresponding to these two variations of parameters highly thetypical traces found inside the closed curve corresponding to FSN II. The scenario (i)will be explored numerically in section 3 in the network equations, and scenario (ii)highlights the fact that for k fixed and σ increasing from 0 to a positive value, weobserve a transition from excitability to MMOs as a noise induced bifurcation.Therefore, the analysis of the system of ODEs arising from taking the average ofthe solution to the mean-field equations in Wilson-Cowan systems with slow adapta-tion currents yields to the conclusion that (i) qualitative behaviors sensitively dependupon noise levels, as we already learnt from the analysis of section 2.1 and (ii) thatvery delicate phenomena such as canard-induced MMOs persist in the presence of In fact other types of dynamics could occur, but a passage to MMOs or relaxation is what isoften observed. M. DESROCHES, M. KRUPA & J. TOUBOUL (cid:1) (cid:3) (cid:1) (cid:2)(cid:1)(cid:6)(cid:4) (cid:1) (cid:2)(cid:1)(cid:4) (cid:1) (cid:2)(cid:1)(cid:3)(cid:4) (cid:1) (cid:2)(cid:2)(cid:1)(cid:4)(cid:2)(cid:1)(cid:5)(cid:2)(cid:1)(cid:6)(cid:2)(cid:1)(cid:7)(cid:2)(cid:1)(cid:8)(cid:3) ku (b) • HB branch from HBbranch from PDisola of MMOs (cid:1) (cid:8) (cid:1) (cid:7) (cid:1) (cid:6) (cid:1) (cid:5) (cid:1) (cid:4) (cid:1) (cid:3)(cid:2)(cid:2)(cid:1)(cid:7)(cid:3)(cid:3)(cid:1)(cid:7)(cid:4) ku HB • (a) • HB Figure 2.4 . Bifurcation diagram of the original system (2.2) as a function of parameter k .Panel (a) shows a large view of the branch of equilibria (black) and the branches of periodic solutions(blue) born at two Hopf bifurcations (black dot) that are symmetric with respect to k = 0 . Panel (b)is an enlarged view in the parameter region where additional branches of periodic solutions exist,both branches born at multiple period-doubling bifurcations ( PD ) as well as closed branches (isolas).The vertical axis correspond to the computed L -norm of the solutions. possibly large noise. This phenomenon is, again, very surprising, and may only beunderstood in the context of collective effects in very large networks. Indeed, largenoise of course will destroy the fine structure of small oscillations in regimes wherethe noiseless system presents mixed-mode oscillations, but we observe here that thesepersist for large noise as a collective effect. Even more interesting is that noise ac-tually induces the apparition of such solutions, and that the shape of the solutionsensitively depends on the precise noise levels neurons are submitted to.
3. Back to the network dynamics.
Thus far, we showed that in the limit N → ∞ , the network equations (1.2) can be reduced to mean-field equations whosemean satisfy an ordinary differential equation where the noise parameter, appearingas a parameter of the limit equation, is related to the presence of canards and MMOs.We now return to the original network equations (1.2) with finite N . We show thatthe predicted presence of complex oscillations as well as canard explosion are generallyrecovered in the simulations of large but finite networks. In certain parameter regions,systematic errors related to the finiteness of the number of neurons considered howeverappear. These finite-size effects, corresponding to the presence of random MMOs intwo-populations networks (in which cases no MMO is present in the mean-field limit)are analyzed in view of the competition between multiple attractors of the mean-fieldlimit and finite-size fluctuations around that limit. In three populations networks, theMMO solutions are generally very fragile close from the transitions, and the finitenessof the number of neurons generally induces early jumps shortening the duration ofthe small oscillations phase.The simulation of large networks presents a number of theoretical and technicalobstructions. Indeed, at least three difficulties arise: the stiffness of the solutionsand the large periods of relaxation oscillations related to the slow-fast nature of theequations, the dimensionality issue related to the number of neurons considered ( N ),and the stochastic nature of the trajectories. The stiffness issue is generally takencare by using refined simulation algorithms, that generally do not have a stochasticcounterpart. Here, accurate simulations will necessitate to use small timesteps dt and large simulation times T (due to the typical periods of relaxation oscillations), tochastic Networks with multiple timescales T /dt , and the number of independent random variables we need to draw during thesimulation will increase proportionally to
N T /dt , which can raise concerns about thefact that the correlations of the pseudo-random sequence of random numbers maysurface. Moreover, with large simulation times, error accumulation will arise, pro-portionally to the computation time. With these limitations in mind, we developeda vectorized simulation algorithm based on Euler-Maruyama scheme implemented inMatlab R (cid:13) , which remains relatively efficient (computation times are almost indepen-dent of the network size) at the levels of precision needed. Typically, the simulation ofa N = 2 000 neurons network over a time interval T = 1 500 with timestep dt = 0 . s on a Macbook Pro with processor 2.3 GHz intel core i7 with 16GB 1600MHz DDR3. The specific development of computational algorithm, precise analysisof the error ad limitations raised by the number of random numbers drawn are notin the scope of the present paper, and we are confident that our implementation,with proper choices of timestep, total times and number of neurons, provide accuratedescriptions of the network behavior in the examples provided below. Westart by analyzing the finite-size effects in our two population network. For our choiceof parameters, we have seen that a subcritical Hopf bifurcation and canard explosionoccurs when varying the intensity of the noise (see Fig. 3.1(a)). The full bifurcationdiagram, plotted in Fig. 2.1 displayed a canard explosion with a branch of limit cyclespersisting until very large values of the noise standard deviations (above 45). Suchstandard deviations correspond to extremely large noise that are generally very hardto simulate accurately and seem unreasonably large in the Wilson-Cowan model. Arestricted bifurcation diagram is displayed in Fig. 3.1(a) and shows that, as a functionof noise levels, the system either presents a unique stable stationary distribution (blueregion), an unstable fixed point and a stable periodic orbit corresponding to relaxationoscillations for large noise (pink region), and an intermediate noise region (yellow)where we have co-existence of a stable stationary distribution and a stable periodicorbit, separated by an unstable periodic orbit. This parameter region is referred toas the bistability regime.In the case where a unique stable solutions exists for the limit averaged equa-tion (pink or blue regions), network trajectories will display random fluctuations wellapproximating the attractor. This is what is plotted in Fig. 3.1, both in the caseof stationary and periodic solutions. The canard explosion, displayed in section 2,corresponds to the sharp transition from fixed point to relaxation oscillations, anddue to the fact that the canard explosion occurs through subcritical Hopf bifurcationprevents from observing the full canard explosion and canard cycles, that correspondto unstable trajectories of the limit ODE and of the network equations.In the bistability regime, trajectories are more complex and finite-size effects willresult in random switches between the two attractors. Typical shape of the phaseplane, displayed in Fig. 3.2(a), show the coexistence of relaxation oscillations andstable fixed points. The separatrix between the two attractors is the unstable limitcycle arising from the Hopf bifurcation, winding around the fixed point. Incidentally,we observe that the unstable periodic orbit and the relaxation oscillations cycles areextremely close (they are arbitrarily close when the intensity of the noise approachesthe value corresponding to the fold of limit cycles). Therefore, in this parameter re-gion, random switches from relaxation oscillations to random fluctuations around thestable fixed point will occur, phenomenon that we term random attractor switching.6
M. DESROCHES, M. KRUPA & J. TOUBOUL
HFLC (a) Bifurcation Diagram (cid:9) (cid:9) (cid:9) (cid:9) (b) Stationary regime, v = 0 . (cid:9) (cid:9) (cid:9) (cid:9) (c) Periodic Regime, v = 2 Figure 3.1 . Network Simulations for the Wilson-Cowan two populations system, whose asymp-totic dynamics was analyzed in 2.1. (a) zoomed bifurcation diagram for relevant noise intensities,(b) and (c) are examples of trajectories of the network equations in the cyan and pink regions re-spectively. (b): stationary solutions, v = 0 . : (left) 10 randomly chosen trajectories of the network(blue: population 1, red: population 2) display random fluctuations of the around the predictedstationary average (cyan: population 1 and black: population 2) computed through the limit ODE.(right) empirical average vs theoretical limit average, show a very precise fit. (c) Same simula-tion and color code as in (b) but for v = 2 corresponding to the oscillatory regime (pink region).Simulation parameters: neurons per population, Euler-Maruyama scheme with dt = 0 . andrandomly chosen initial conditions. Similarly to what was shown in [6], we may interpret the switching trajectoriesas random MMOs. However, in contrast to their setting, the noise perturbationinvolved in our case are not driven by white noise with constant standard deviation.As we analyze in a forthcoming study, the finite-size perturbations correspond toGaussian perturbations with standard deviation function of the solution of the mean-field equations. The nature of this process implies that the random switches donot occur completely randomly, but depend on the size of the fluctuations, which isgoverned by a function of the mean-field equation. A consequence will be, as presentedin the example plotted in Fig. 3.2(b), that the system has a bias to switch in phasewith the virtual periodic attractor. This surprising phenomenon will be analyzed ingreater detail in a forthcoming study.Let us eventually discuss the occurrence of switches as a function of the size of thenetwork and of the parameters. We observe that the attraction bassin of the stationarysolution becomes increasingly large as noise approaches the value corresponding to thefold of limit cycles, and increasingly small when approaching the Hopf bifurcation.Close from the value of the fold of limit cycles, random switches from the periodic tothe stationary solutions will become very probable and switches from the stationary tothe periodic solutions are less probable, implying that the averaged time spent in theneighborhood of the stationary solution will become larger as we approach the fold oflimit cycle. A symmetrical phenomenon will arise close from the Hopf bifurcation. Inorder to illustrate this phenomenon, we display in Fig. 3.3 the average time spent in a tochastic Networks with multiple timescales -0.500.511.522.5U2 0 0.5 1 1.5 2U1 (a) Phase plane, v = 1 . v = 1 . Figure 3.2 . Bistable regime in the two-populations Wilson-Cowan system. (a) Geometry of thephase plane of the limit ODE, shows the presence of a relaxation oscillations cycle (black) and anunstable cycle (blue) delineating the attraction basin of the stable fixed point. Red and Green curvesare the nullclines of the system, the triangle (square) depicts a saddle (unstable focus) fixed point.(b) depicts a typical trajectory of random attractor switching, with neurons per population and v = 1 . , same color code as in Fig. 3.1. Dotted lines depict the location of the stationary mean forboth populations. neighborhood of the stationary solution (actually inside the deterministic attractionbassin of the stable fixed point depicted in blue in Fig. 3.2(a)) as a function of thenoise level, averaged across 20 trajectories. These switching probabilities also dependon the size of the network. Indeed, as the number of neurons tends to infinity, thatprobability tends to zero, and the trajectory will be deterministically locked to thatpredicted by the mean-field theory and the initial condition, and the smaller thenetwork size, the more likely the switches will occur. Figure 3.3 . Proportion of time along trajectories spent in the vicinity of the stationary solutionas a function of noise levels. Colors correspond to the regions in the bifurcation diagram 3.1(a).
We therefore conclude from the present analysis that finite-size effects barelyaffect the solution, and that the mean-field model accurately describes the behavior8
M. DESROCHES, M. KRUPA & J. TOUBOUL of the network in most of the parameter space, except during the canard explosion,in an exponentially small region of parameters, and in that region, random switchesbetween multiple attractors may occur, as a function of the network size and thegeometry of the limit dynamical system.
We now investigate the behavior of finite-sized networks in the 2 populationsWilson-Cowan system with slow adaptation currents analyzed in the mean-field limitin section 2.2. The analysis the resulting mean-field system showed that the systempresented a FSN II bifurcation related to the presence of mixed-mode oscillations, theshape of which sensitively depended of the noise levels σ and k , as depicted in Fig. 2.2.The reduction obtained of these equations allowed analyze the behavior for extremelyhigh values of the noise which are not accessible through numerical simulations of thenoisy network, chiefly due to precision constraints as we mentioned.We now confront the different behaviors observed in section 2.2 to simulationsof the noisy network equations with the same parameters. We expect in this casethat the trajectories obtained in the MMO regime to be less resilient to perturbationsthan the stationary and periodic orbits analyzed in the two-populations network,due to the sensitivity of the attractor in particular in the small oscillations phase.In contrast with the analysis of the two-populations Wilson-Cowan model, we willtherefore observe more drastic finite-size effects even in regimes where there exists aunique attractor, since it is the very structure of the attractor that presents sensitivityto perturbations. These will lead to escape the attractor during the small oscillationsphase, or early jumps .This is what we observe in Fig. 3.4 where we reproduce the bifurcation scenariodisplayed in Fig. 2.2 for fixed noise σ = 2 and different values of the parameter k (the only scenario we did not display was (B), which requires very simulations on verylong time intervals due to the very large period of the mixed-mode oscillation). In allcases, we recover the qualitative nature of the solutions predicted in the infinite-sizelimit, and early jumps, corresponding to systematic random escape of the attractor,occur due to finite-size effects. Moreover, we observe that the small oscillations tendto be hidden by the noise related to the finite-size effects. These two finite-sizeeffects of course disappear when considering larger networks and one recovers withprecision the attractor predicted by the theory, as we show in Fig. 3.5 for a network of N = 10 000 neurons. But as soon as the number of neurons is smaller, the network willsystematically fail to follow precisely the predicted MMO and will generally escapethe attractor prior than the infinite-sized system.
4. Conclusion.
In a number of applied domains such as physics, chemistry andneurosciences, the emergence of complex oscillations such as mixed-mode oscillationshave been theoretically related to canard explosions and multiple timescales dynamicsof deterministic dynamical systems (see [12] for a review). But the mechanisms un-derlying the existence of such orbits generally involve a number of interacting agents,for instance neurons, that are characterized by a stochastic activity. It was henceimportant in that context to show that such orbits persist in these contexts, andto identify the mechanisms accounting for the emergence of such collective behav-iors. At the level of one single cell with stochastic perturbations, studies showed thatsystems presenting canard explosions displayed random MMOs, with the effect of ran-domly modifying the period of the oscillations and changing the qualitative natureof the solution even in small noise [6]. In this paper, the authors analyze stochas-tic perturbations of Fitzhugh-Nagumo model in the relaxation oscillations regime tochastic Networks with multiple timescales (a) k = − . (b) k = − . (c) k = − . (d) k = − (e) k = − Figure 3.4 . Solution to the mean-field equations (1.7) (red and blue lines) and empiricalaverages of the network equation for N = 2 000 neurons in the cases (a), (c), (d), (e) and (f) ofFig. 2.2. In each case, we observe that qualitative behaviors are recovered in the network equations,and that early jumps occur systematically. (a) k = − . v = 2 (b) k = − . Figure 3.5 . Same quantities as in Fig. 3.4 for N = 10 000 neurons, (a) same parameters as inFig. 2.2(d) and (b): same parameters as in Fig. 2.2 with k = − . , trajectory presenting a numberof visible small oscillations. concluded indeed noise resulted in sensible modifications of the trajectories with ran-dom interspike intervals and trajectories similar to MMOs with a random numberof small oscillations that asymptotically distribute as a geometric random variable.As of mixed-mode oscillations, Berglund and collaborators [5] analyzed the effect ofGaussian noise on a 3 dimensional system with different timescales known to featuremixed-mode oscillations and show that noise changes the global dynamics and mixedmode patterns, and that there exists a maximal amount of noise above which smalloscillations are hidden by noise. These studies pointed towards the fact that otherphenomena, more resilient to noise, took place in large noisy systems giving rise to0 M. DESROCHES, M. KRUPA & J. TOUBOUL such complex oscillations in the presence of multiple timescales.Here, we revisited these questions and proposed that self-averaging of the trajec-tories in large networks could provide a simple generic explanation to the emergence ofcomplex oscillations in large-scale noisy networks. We based our analysis on a classicalmodel of neuronal network, the Wilson-Cowan system, which describes heuristicallythe dynamics of firing-rates of neurons. This system presented two advantages inour view. First an unsolved mathematical problem was to rigorously derive asymp-totic limits of these networks with nonlinear mean-field interactions as the numberof neurons tends to infinity. The nature of is indeed very different to usual networkequations analyzed in the domain of mean-field limits and propagation of chaos as de-veloped in gas dynamics and extended to neurosciences [26, 28], since the interactingterm is a nonlinear function of a mean-field term, and moreover that the neurons arenot governed by a stochastic differential equation but rather by a differential equationwith stochastic coefficients ξ it . This lead us to extend the theory of mean-field limits towards this new kind of equations, which we developed in a relatively general fashionin section 1. Applied to Wilson-Cowan equations, this result allowed us to exactlyreduce the network equations to a multiple timescales ODE acting on the mean of thesolution. Interestingly, these equations reduce to a modified Wilson-Cowan systemwhere the noise level is embedded in the nonlinearities.The analysis of the mean equation developed in section 2 had at least two interests.First, the fact that the asymptotic equation on the mean of the solution is itselfsolution of a Wilson-Cowan type of equation allowed to go deeper in the understandingof these models. And indeed, the behavior of Wilson-Cowan systems in the presenceof multiple timescales was still an open problem important to solve since the originalWilson-Cowan model introduced in 1972 in [31] had an intrinsic slow-fast structurethat was never explored as such in the literature. Second, this analysis allowed tounderstand the qualitative role of noise in the dynamics. As was done in [29, 27],we demonstrated here again that the noise did not have a pure disturbing effect onthe trajectories but rather a dramatic qualitative effect on the form of the solutions,governing transitions from stationary to periodic behaviors and complex oscillatorybehaviors.Of course, the same study would hold in the case of the simpler Amari-Hopfieldmodel analyzed in [29], where the interaction between neurons is the sum of the firing-rate seen as a sigmoidal transform of the voltage. These models do not take intoaccount the nonlinearities in the synaptic integration taken into account in Wilson-Cowan systems. For such networks, the mean-field equations has Gaussian asymptoticsolutions, whose mean and standard deviation satisfy a closed set of ordinary differ-ential equations [29]. Our analysis can easily be extended to such systems, where thelimit equations obtained can be shown to have Gaussian solution whose mean satisfythe set of ODEs: τ α d µ α d t = − µ α + P (cid:88) β =1 J αβ G αβ ( µ β , σ α )where G αβ is similar to our effective sigmoid function 1.12 and embeds dependenceon the noise level σ α . In the presence of multiple timescales, we expect to observethe same kind of phenomena, and a similar dependence on noise levels, and these The proof uses similar tools as usual analysis.tochastic Networks with multiple timescales
REFERENCES[1]
S. Amari , Characteristics of random nets of analog neuron-like elements , Syst. Man Cybernet.SMC-2, (1972).[2]
S.-I. Amari , Dynamics of pattern formation in lateral-inhibition type neural fields , BiologicalCybernetics, 27 (1977), pp. 77–87.[3]
E. Benoˆıt , Canards et enlacements , Publ. Math. IHES, 72 (1990), pp. 63–91.[4]
E. Benoˆıt, J.-L. Callot, F. Diener, and M. Diener , Chasse au canard , Collect. Math., 31(1981), pp. 37–119.[5]
N. Berglund, B. Gentz, and C. Kuehn , Hunting french ducks in a noisy environment ,Journal of Differential Equations, (2012).[6]
N. Berglund and D. Landon , Mixed-mode oscillations and interspike interval statistics inthe stochastic fitzhugh–nagumo model , Nonlinearity, 25 (2012), p. 2303.[7]
P.C. Bressloff, J.D. Cowan, M. Golubitsky, P.J. Thomas, and M.C. Wiener , WhatGeometric Visual Hallucinations Tell Us about the Visual Cortex , Neural Computation,14 (2002), pp. 473–491.[8]
M. Brøns, M. Krupa, and M. Wechselberger , Mixed-Mode Oscillations Due to GeneralizedCanard Phenomenon , in Bifurcation Theory and spatio-temporal pattern formation, FieldsInstitute Communications, 2006, pp. 39– 64.[9]
M. Brøns , Bifurcations and instabilities in the Greitzer model for compressor system surge ,Math. Eng. Ind., 2 (1988), pp. 51–63.[10]
Gyorgy Buzsaki , Rhythms of the Brain , Oxford University Press, USA, 2009.[11]
R. Curtu and J. Rubin , Interaction of canard and singular hopf mechanisms in a neuralmodel , SIAM Journal on Applied Dynamical Systems, 10 (2011), p. 1443.[12]
M. Desroches, J. Guckenheimer, B. Krauskopf, C. Kuehn, H. M Osinga, and M. Wech-selberger , Mixed-mode oscillations with multiple time scales , SIAM Review, 54 (2012),pp. 211–288. M. DESROCHES, M. KRUPA & J. TOUBOUL[13]
B. Ermentrout , Neural networks as spatio-temporal pattern-forming systems , Reports onprogress in physics, 61 (1998), p. 353.[14]
GB Ermentrout and J.D. Cowan , Temporal oscillations in neuronal nets , Journal of math-ematical biology, 7 (1979), pp. 265–280.[15]
B. Fernandez and S. M´el´eard , A hilbertian approach for fluctuations on the mckean-vlasovmodel , Stochastic processes and their applications, 71 (1997), pp. 33–53.[16]
J. Guckenheimer , Singular Hopf bifurcation in systems with two slow variables , SIADS, 7(2008), pp. 1355–1377.[17]
E.M. Izhikevich , Neural excitability, spiking and bursting , International Journal of Bifurcationand Chaos, 10 (2000), pp. 1171–1266.[18]
Eric R Kandel, James H Schwartz, Thomas M Jessell, et al. , Principles of neural science ,vol. 4, McGraw-Hill New York, 2000.[19]
I. Karatzas and S. Shreve , Brownian motion and stochatic calculus , Springer, 1987.[20]
M. Krupa, N. Popovi´c, N. Kopell, and H.G. Rotstein , Mixed-mode oscillations in a threetime-scale model for the dopaminergic neuron , Chaos: An Interdisciplinary Journal ofNonlinear Science, 18 (2008), pp. 015106–015106.[21]
M. Krupa and P. Szmolyan , Relaxation oscillations and canard explosion , JDE, 174 (2001),pp. 312–368.[22]
Martin Krupa and Martin Wechselberger , Local analysis near a folded saddle-node sin-gularity , Journal of Differential Equations, 248 (2010), pp. 2841–2888.[23]
A. Milik, P. Szmolyan, H. Loeffelmann, and E. Groeller , Geometry of mixed-mode os-cillations in the 3d autocatalator , Int. J. Bifur. Chaos, 8 (1998), pp. 505–519.[24]
P. Szmolyan and M. Wechselberger , Canards in R , JDE, 177 (2001), pp. 419–453.[25] AS Sznitman , Nonlinear reflecting diffusion process, and the propagation of chaos and fluctu-ations associated , Journal of Functional Analysis, 56 (1984), pp. 311–336.[26] ,
Topics in propagation of chaos , Ecole d’Et´e de Probabilit´es de Saint-Flour XIX, (1989),pp. 165–251.[27]
Jonathan Touboul , On the dynamics of mean-field equations for stochastic neural fields withdelays .[28] ,
Propagation of chaos in neural fields , (submitted), (2011).[29]
Jonathan Touboul, Geoffroy Hermann, and Olivier Faugeras , Noise-induced behaviorsin neural mean field dynamics , SIAM J. on Dynamical Systems, 11 (2011).[30]
M. Wechselberger , Existence and Bifurcations of Canards in R in the Case of the FoldedNode , SIADS, 4 (2005), pp. 101–139.[31] H.R. Wilson and J.D. Cowan , Excitatory and inhibitory interactions in localized populationsof model neurons , Biophys. J., 12 (1972), pp. 1–24.[32]