Control of coherence resonance by self-induced stochastic resonance in a multiplex neural network
CControl of coherence resonance by self-induced stochastic resonancein a multiplex neural network
Marius E. Yamakou ∗ and J¨urgen Jost
1, 2, † Max-Planck-Institut f¨ur Mathematik in den Naturwissenschaften, Inselstr. 22, 04103 Leipzig, Germany Santa Fe Institute for the Sciences of Complexity, NM 87501, Santa Fe, USA (Dated: May 24, 2019)We consider a two-layer multiplex network of diffusively coupled FitzHugh-Nagumo (FHN) neu-rons in the excitable regime. We show that the phenomenon of coherence resonance (CR) in onelayer cannot only be controlled by the network topology, the intra and inter-layer time-delayedcouplings, but also by another phenomenon, namely, self-induced stochastic resonance (SISR) inthe other layer. Numerical computations show that when the layers are isolated, each of thesenoise-induced phenomena is weakened (strengthened) by a sparser (denser) ring network topology,stronger (weaker) intra-layer coupling forces, and longer (shorter) intra-layer time delays. However,CR shows a much higher sensitivity than SISR to changes in these control parameters. It is alsoshown, in contrast to SISR in a single isolated FHN neuron, that the maximum noise amplitude atwhich SISR occurs in the network of coupled FHN neurons is controllable, especially in the regimeof strong coupling forces and long time delays. In order to use SISR in the first layer of the multiplexnetwork to control CR in the second layer, we first choose the control parameters of the second layerin isolation such that in one case CR is poor and in another case, non-existent. It is then shownthat a pronounced SISR cannot only significantly improve a poor CR, but can also induce a pro-nounced CR, which was non-existent in the isolated second layer. In contrast to strong intra-layercoupling forces, strong inter-layer coupling forces are found to enhance CR. While long inter-layertime delays just as long intra-layer time delays, deteriorates CR. Most importantly, we find that ina strong inter-layer coupling regime, SISR in the first layer performs better than CR in enhancingCR in the second layer. But in a weak inter-layer coupling regime, CR in the first layer performsbetter than SISR in enhancing CR in the second layer. Our results could find novel applications innoisy neural network dynamics and engineering.
I. INTRODUCTION
Noise is ubiquitous and even relatively low noise inten-sities can have significant and counter-intuitive effects onnonlinear dynamical systems. For instance, a randomlyperturbed system in the excitable regime can produceoscillatory responses which possess a high degree of co-herence and yet are completely different from what isobserved in the absence of the random perturbations.In this paper, we focus on two such phenomena: co-herence resonance (CR) [1–4] and self-induced stochasticresonance (SISR) [5–7] in a network of diffusively cou-pled FitzHugh-Nagumo (FHN) neurons; a paradigmaticmodel which describes the excitability and spiking be-havior of neurons [8]. Lee DeVille et al. [9] (see also [10])have analyzed the differences in the mechanisms leadingto CR and SISR in an isolated excitable system. Theyshowed that even though CR and SISR lead to the emer-gence of the same dynamical behavior (i.e., weak-noise-induced coherent oscillations) in excitable systems, theyare fundamentally different in their dynamical and emer-gent nature.CR occurs when the regularity of noise-induced oscilla-tions of an excitable system is a non-monotonic functionof the noise amplitude, and it is optimally correlated at ∗ [email protected] † [email protected] some non-zero value of the noise intensity. Thus, duringCR, there exists a maximal degree of coherence in the os-cillations of the excitable system for some intermediatenoise amplitude. CR is known to be robust in the sensethat the coherence of the noise-induced oscillations is in-sensitive against variations of the timescale separationratio between the fast and slow variables of the excitablesystem and the amplitude of the noise [7, 9]. The crucialcondition necessary for the occurrence of CR is the prox-imity of the system’s parameters to the Hopf bifurcation[2–4, 11, 12] or the saddle-node bifurcation of limit cycles[13–16]. This means that the system must be near butbefore bifurcation thresholds so that a relatively smallnoise amplitude can easily (without overwhelming theentire dynamics) drive the system towards the determin-istic limit cycle which emerges right after the bifurcation.This implies that noise plays a rather passive role in thephenomenon of CR.On the other hand, just like CR, SISR occurs whenan optimal noise amplitude induces coherent oscillationsin a system in the excitable regime (i.e., a regime thatdoes not allow for self-sustained oscillations in absenceof perturbations). The similarity between CR and SISRhowever ends at this point. Unlike CR, SISR is robustto parameter tuning and therefore does not require thesystem’s parameter to be in the immediate neighborhoodof bifurcation thresholds. SISR relies completely on theability of the slow-fast excitable system to match thetimescales of the motion of trajectories along the slow a r X i v : . [ n li n . AO ] M a y anifold and of escape events of these trajectories fromthe slow manifold. The characteristics of the coherentnoise-induced oscillations due to SISR have been shownto depend non-trivially on both the timescale separationratio weakly via its logarithm but more importantly onthe amplitude of the noise [5, 9, 10, 17–20]. This meansthat for SISR, noise plays a more subtle role in the gen-eration of the coherent oscillations than in CR.For the specific case of the FHN model, CR occursonly when the noise term is attached to the slow recov-ery variable equation (one should note that this is nota necessary condition for the occurrence of CR in gen-eral, as CR may still occur in slow-fast systems wherethe noise term is attached to the fast variable equation;see for example [12]). As pointed out earlier, what is cru-cial for CR is the proximity of the system’s parameter tothe bifurcation thresholds. Therefore, for CR, the noiseterm has to be added to the system such that it can di-rectly affect the bifurcation parameter [1–4]. This is whythe noise term is added to the slow variable equation ofthe FHN model in the case of CR.On the other hand, because SISR occurs only in thesingular limit (i.e., the limit as the timescale separationratio tends to zero), it requires the noise term to be at-tached to the fast variable equation of the slow-fast ex-citable system which evolves itself on the fast timescale[5–7, 9, 10, 17–19, 21]. Because SISR is more stable thanCR to parametric perturbations, it might be useful froman engineering point of view to try to control CR usingSISR and not the other way round.One of the most relevant questions today is relatedto the control of noise-induced phenomena in networksof coupled oscillators. Several control schemes of CRbased on time-delayed feedbacks and network topologyand heterogeneity have been studied, in particular withthe FHN model equations. It has been shown that appro-priate selection of the time-delayed feedback parameterscan modulate CR [13, 22–27]. With a particular versionof the FHN model, it has been shown that in a ring oflocally coupled units, time-delay weakens CR while innon-local and global coupling cases, only appropriate de-lay values could strengthen or weaken CR [27]. C. Zhouet at [28] showed that the spatial heterogeneity of the bi-furcation parameter of coupled FHN oscillators can alsocontrol CR. In this control scheme, CR is enhanced inthe network producing oscillations that are more coher-ent than those generated during CR in a single FHN unit.It has also been shown that in a ring of non-locally cou-pled oscillators, time-delays can control the parametersrange in which CR-induced chimera states occur [29].Recently, the use of multiplex network in the controlstrategy of dynamical behaviors of coupled oscillators hasattracted some attention. In a multiplex network, eachtype of interaction between the nodes is described by asingle layer network and the different layers of networksdescribe the different modes of interaction. In such net-works, the layers contain the same number of nodes andthe interaction between the layers are allowed only for replica nodes. In the case of neural networks, the neu-rons can form different layers depending on their connec-tivity through a chemical link or by an ionic channel. Inbrain networks, different regions can be seen connectedby functional and structural neural networks [30–32].Multiplex networks open up new possibilities of con-trol, allowing to regulate nonlinear systems by means ofthe interplay between dynamics and multiplexing. Con-trol mechanisms based on multiplexing have many ad-vantages. In particular, the dynamics of one layer can becontrolled by adjusting the parameters of another layer.This is important from the point of view of engineeringand brain surgery since it is not always possible to di-rectly access the desired layer, while the network withwhich this layer is multiplexed may be accessible andadaptable. The multiplexing of networks has been shownto control many dynamical behaviors including synchro-nization [33–37] and pattern formation [38–44]. However,the control of CR based on the multiplexing of networkshas not been extensively investigated. This alternativecontrol scheme of CR has been only very recently stud-ied in [45]. Here, it is shown that even weak multiplexingcan induce CR in a network of FHN neurons that doesnot demonstrate CR in isolation. Moreover, it has beenshown that the multiplex-induced CR in the layer whichis deterministic in isolation can manifest itself even morestrongly than that in the noisy layer. In [45], each layer ofthe multiplex network considered undergoes CR, as thenoise term in attached only to the slow variable equationsand the bifurcation parameter is set in the non-oscillatoryregime but close to the Hopf bifurcation threshold.The objective of this paper is to propose a new controlstrategy of CR in a two-layer multiplex network of FHNneurons in the excitable regime. Our control strategy isnot only based on the combined effects of time-delayedcouplings, intra-layer network topology, and multiplex-ing, but will also include the phenomenon of SISR. Inshort, we aim at controlling CR with SISR in a time-delayed multiplex network. Our aim and setting are dif-ferent from the one in [45] which considers the control ofCR in one layer of the multiplex with CR in the otherlayer and which ignores time-delays (which are ubiqui-tous in complex systems). The main questions we wantto address in this work are the following:(i) In a two-layer multiplex network, can SISR occurringin one layer be used to control CR in the other layer?(ii) Which phenomenon (i.e., SISR or CR, occurring inone layer) would better enhance CR in the other layer?This paper is organized as follows: In section II, wepresent the model equation and determine the parameterrange of its excitable regime. In section III, we presentour results. In the first part of this section, we considerisolated layers and investigate the effects of time-delayedcouplings and different ring network topologies on CRand SISR. In the second part, we consider the multiplexnetwork of the layers where we present and compare thecontrol schemes of CR based on SISR and on CR. Insection IV, we have the summary and conclusion.2 I. MODEL AND ITS EXCITABLE REGIME
We consider the following two-layer multiplex network,where each layer represents a ring of N diffusively cou-pled FHN neurons in the excitable regime and in thepresence of noise: dv i dt = v i − v i − w i + κ n i + n (cid:80) j = i − n (cid:16) v j ( t − τ ) − v i ( t ) (cid:17) + κ (cid:16) v i ( t − τ ) − v i ( t ) (cid:17) + σ dW i dt ,dw i dt = ε ( v i + α − β w i ) + σ dW i dt ,dv i dt = v i − v i − w i + κ n i + n (cid:80) j = i − n (cid:16) v j ( t − τ ) − v i ( t ) (cid:17) + κ (cid:16) v i ( t − τ ) − v i ( t ) (cid:17) ,dw i dt = ε ( v i + α − β w i ) + σ dW i dt . (1) v i ∈ R and w i ∈ R represent the fast membrane poten-tial variables and the slow recovery current variables inthe first layer, respectively. The index i = 1 , ..., N standsfor the node i in the ring network of N neurons and allindices are modulo N . Similarly, v i ∈ R and w i ∈ R respectively represent the membrane potential and re-covery current variables for the neurons in the secondlayer. 0 < ε (cid:28) < ε (cid:28) β > β > α ∈ (0 ,
1) is a con-stant parameter. For the i th neuron in a ring, n and n represent the number of nearest connected neighbors ineach direction on the ring. In each ring, we can have 3different coupling schemes between the neurons: n , = 1for local coupling, n , = ( N − / N ) for globalcoupling, and 1 < n , < ( N − / n and n acts as independent control parametersfor the ring topology of the underlying networks of themultiplex. The coupling strength within the layers (intra-layer coupling) is measured by κ and κ in the first andsecond layer, respectively. The intra-layer couplings havethe form of a diffusive coupling, that is in each layer, thecoupling term vanishes if v i and v j (resp. v i and v j )are equal. τ and τ are respectively the signal propaga-tion delays within the first and second layer, which wehenceforth refer to as layer 1 and layer 2, respectively.The coupling between the layers (inter-layer coupling) isbidirectional, diffusive, and its strength is measured by κ . τ is the signal propagation delay between layer 1and layer 2 of the multiplex network (see Fig.1). σ , σ ,and σ stand for noise intensities. dW ki dt , k = 1 , , layer 1: κ , τ multiplexing: κ , τ layer 2: κ , τ FIG. 1. (Color online) Schematic diagram showing a multi-plex network of neurons consisting of two layers. Each node ina layer is connected only to the replica node in the other layer.The lower layer (layer 2) in red represents a coupled ring net-work topology, where CR occurs. The upper layer (layer 1) inblue represents a coupled ring network (with possibly a differ-ent ring topology), where either SISR or CR can occur and beused together with the intra and inter time-delayed couplingsparameters { τ , κ , τ , κ , τ , κ } , to control CR in layer 2. uncorrelated Gaussian white noises, the formal derivativeof Brownian motion with mean zero and unit variance.The conditions required for the occurrence of CR andSISR in an isolated FHN neuron are well known [2, 3, 7,9, 19]. An excitable regime is a necessary and commonrequirement for both CR and SISR. In this regime, theisolated FHN neuron has a unique and stable fixed point.Choosing an initial condition in the basin of attractionof this fixed point will result in at most one large non-monotonic excursion into the phase space after which thetrajectory returns to this fixed point and stays there untilthe initial conditions are changed again [19, 46].In any isolated layer ( κ = 0) of the multiplex net-work, we can set an isolated neuron ( κ = 0 or κ = 0)in an excitable regime by fixing the parameters to cer-tain values: We set α = 0 . β and β are chosen such that β > β h ( ε ) and β > β h ( ε ), where β h ( ε ) and β h ( ε )are the Hopf bifurcations values of a neuron in layer 1and layer 2, respectively. For each neuron in layer 1 andin layer 2, we set ε = 0 . ε = 0 .
01 and calcu-late the Hopf bifurcation values as β h ( ε ) = 0 . β h ( ε ) = 0 . β ≤ β h ( ε ) and β ≤ β h ( ε ), an isolated neuron is in theoscillatory regime; a regime that we want to avoid as thecoherent oscillations generated by SISR are due only tothe presence of noise and not because of the occurrence ofa Hopf bifurcation [19]. While the coherent oscillationsgenerated by CR, even though they require proximity tothe Hopf bifurcation threshold, they still need the systemto be in the excitable regime [2].To investigate the effects of time-delayed coupling anddifferent ring network topologies on SISR in layer 1 andCR in layer 2, we have to make sure that for each layer,the entire network is in the excitable regime and not just3he individual neurons. Therefore, each layer should bein the parameter regime where there are no time-delayedcoupling induced oscillations. It is important to iden-tify this regime (if it exists) because certain time-delayedcouplings may induce self-sustained oscillations in a layernetwork, even when all of the individual neurons are inthe excitable regime. A saddle-node bifurcation resultingin a pair of stable and unstable limit cycles may generatethese delayed-coupling induced self-sustained oscillations[47].To ensure that the deterministic isolated layers are ina parameter regime where time-delayed coupling inducedoscillations are absent, we numerically calculate the oscil-lation period, i.e., the interspike interval ( ISI ) of the syn-chronized oscillations for different values of τ and κ inlayer 1 ( ε = 0 .
01) and τ and κ in layer 2 ( ε = 0 . β and β used in the calculations, layer1 and layer 2 showed the same behavior. We thereforeonly show the results of layer 2. Throughout this paper,we will for the purpose of simplicity, set β = β = 0 . ISI is plotted in the ( κ − τ ) and ( N − n )parameter spaces for different values of the Hopf bifur-cation parameter β . In Fig.2 (a) and (b) , β > β h ( ε ).While in Fig.2 (c) and (d) , β ≤ β h ( ε ). We observethat there are no time-delayed coupling oscillations inFig.2 (a) and (b) . In this case, the excitable regimeof each neuron persists in the entire layer network forthe considered delayed-coupling intervals, network size N , and number of neighbors n . The gray region inFig.2 (a) and gray lines in Fig.2 (b) represent the absenceof time-delayed coupling induced oscillations : this is theexcitable regime in which we focus our study and whichconsists of a stable homogeneous steady state given by( v ∗ i , w ∗ i ) = ( − . , − . (c) and (d) , when all the neu-rons in the layer are in the oscillatory regime (i.e., β ≤ β h ( ε )), we instead have the phenomenon of amplitudedeath which is characterized by oscillation quenching in-duced by the presence of delay in the coupling betweenthe neurons in the network [48–51]. Fig.2 (b) and (d) show that the excitable regime of the layer network isindependent of its size N and the number of connectedneighbors n . The oscillatory regime in Fig.2 (c) (redregion) also shows (in Fig.2 (e) ) an independent behav-ior with respect to N and n , except that the period ofthe oscillations increases with decreasing number of con-nected neighbors n . III. METHODS AND RESULTS
The coefficient of variation ( CV ), the correlation time,the power spectral density, and the signal-to-noise ratioare commonly used measures to qualify and quantify theregularity of noise induced oscillations, and therefore ameasure of how pronounced CR and SISR can be at aparticular noise amplitude. From a neurobiological point FIG. 2. (Color online) The
ISI is color coded and correspondsto
ISI + ξ with 0 < ξ (cid:28) (a): the gray region in the ( κ − τ )plane represents the excitable regime (consisting of a stablehomogeneous fixed point) of the isolated layer 2 with theHopf bifurcation parameter β = 0 . > β h ( ε ) = 0 . n = 1. (b): the excitable regime ( κ , τ ) = (0 . , .
0) in (a) is independent of N and n . (c): the red region repre-sents the oscillatory regime β = 0 . < β h ( ε ) and thegray region an excitable regime induced by amplitude death, n = 1. (d): the excitable regime ( κ , τ ) = (0 . , .
0) in (c) is independent of N and n . (e): the oscillatory regime( κ , τ ) = (0 . , .
0) in (c) is independent of N and n . Initialconditions corresponds to a spike for all neurons. Other pa-rameters of the isolated layer 2: N = 25, ε = 0 . α = 0 . σ = 0 . of view, CV is more important because it is related tothe timing precision of the information processing in neu-ral systems [53]. For a Poissonian spike train (rare andincoherent spiking), CV = 1. If CV <
1, the sequencebecomes more coherent, and CV vanishes for a periodicdeterministic spike train. CV values greater than 1 cor-respond to a point process that is more variable than aPoisson process [19, 52]. Because of the importance of CV in neural information processing, we shall use it tocharacterize the regularity of the noise-induced oscilla-tions generated by SISR and CR in our neural network.4he CV of an isolated neuron is defined as [2]: CV = (cid:112) (cid:104) ISI (cid:105) − (cid:104) ISI (cid:105) (cid:104) ISI (cid:105) , (2)where (cid:104) ISI (cid:105) and (cid:104)
ISI (cid:105) represent the mean and themean squared interspike intervals, respectively. Theabove definition of CV is limited to characterizing CRand SISR in an isolated neuron. For a network of cou-pled neurons, CR and SISR can be measured by redefin-ing CV as follows [54]: CV = (cid:113) (cid:104) ISI (cid:105) − (cid:104) ISI (cid:105) (cid:104) ISI (cid:105) , (3)where the extra bar indicates the additional average overneurons in the layer network.In the rest of our numerical calculations, we use thefourth-order Runge-Kutta algorithm for stochastic pro-cesses [55] to integrate over a very long time interval(600000 time units) and then averaging over time, neu-rons ( N = 25), and 5 realizations of each noise amplitude.In the numerical simulations, this long time interval per-mitted us to collect with a small noise amplitude at least70500 ISIs with ε = 0 .
01 used for CR and at least 125ISIs with ε = 0 . (cid:28) ε used for SISR. A. Isolated layers ( κ = 0 ): Network topology andtime-delayed coupling effects on CR and SISR In order to control (strengthen or weaken) CR in layer2 of the multiplex network using SISR (or CR itself) inlayer 1, we should naturally first understand the effectsof the intra-time-delayed couplings and different ring net-work topologies on CR and SISR in the isolated layers.We therefore begin our study with the dynamics of dis-connected ( κ = 0) layers and consider an isolated one-layer network. Coherence Resonance
The conditions required for the occurrence of CR in anisolated ( κ = 0) FHN neuron in layer 2 are: the prox-imity of the parameter to the Hopf bifurcation value anda sufficiently small noise amplitude that does not over-whelm the entire dynamics, and which is also sufficientlylarge so that trajectories do not stick for too long in thebasin of attraction of the stable fixed point to produce aPoisson sequence of spikes [2, 7, 9, 10], that is: β − β h ( ε ) ≤ δ,σ (cid:28) ,σ (cid:38) (1 − w ∗ i − ε β ) / . (4)where δ > σ is attached only to theslow variable equations. It has been shown that CR can be enhanced in a net-work of coupled oscillators by choosing appropriate cou-pling strengths and time delays [54, 56]. But one shouldnote that the response of CR to changes in time-delayedcoupling can vary from one kind of nonlinear oscillator toanother depending on the details of the equation. Thatis, there is no general behavior as stronger coupling forcesmay strengthen CR in a network of one kind of nonlinearoscillators [45], but as we shall see in the current paper,it can instead weaken CR in a network of just slightlymodified equations of the nonlinear oscillators. There-fore, without investigations, one cannot say a priori howCR will behave in response to changes in time-delayedcouplings in a network of a particular nonlinear oscillator,even if this behavior is known for a network of anothernonlinear oscillator, with a slightly modified equation.Fig.3 shows the variation of CV against the noise am-plitude σ in the isolated layer 2 with a locally ( n = 1)coupled ring network topology. In the numerical com-putations, we choose ε = 0 .
01, a standard value used inthe literature of the FHN neuron model. For some valuesof the time delay τ and the coupling force κ , the round-bottom CV -curves show a non-monotonic dependence of CV on σ with a well-pronounced deep minimum. Thenetwork thus exhibits a pronounced CR with coherent os-cillations emerging in a very narrow interval of the noiseamplitudes. For our model equation, the general behav-ior of the CV -curves in response to changes in the time-delayed coupling is the following: CV -curves are shiftedto higher values with their minima occurring at highernoise amplitudes for stronger coupling forces and longertime delays.In Fig.3 (a) , for a fixed and weak coupling force of κ = 0 .
1, all the CV -curves obtained with different timedelays ( τ = 0 . , τ = 2 . , τ = 5 . τ = 10 .
0) show adeep minimum, indicating a high degree of coherence ofthe oscillation due to CR. We notice that for this weakcoupling force, as τ increases, the entire CV -curve isonly slightly shifted up and to the right, thus slightlyweakening CR and increasing the value of the noise am-plitude at which maximal coherence is achieved. For ex-ample, for τ = 0 . CV min = 0 .
13 at σ = 0 . × − and for τ = 10 . CV min = 0 .
20 at σ = 5 . × − .In Fig.3 (b) , we now have a strong coupling force( κ = 1 . (a) , but is now muchstronger, i.e., with strong coupling forces, the time delayshave stronger effects on the CV -curves than with weakcoupling forces. The CV -curves are shifted up to highervalues with the minima occurring at much larger noiseamplitudes than in Fig.3 (a) . For example, for τ = 0 . CV min = 0 .
14 at σ = 0 . × − and for τ = 7 . CV min = 0 .
56 at σ = 4 . × − . We notice that for κ = 1 . τ = 7 .
0, even though the red curve showsa non-monotonic dependence of CV on σ , it still has arelatively high minimum value of CV min = 0 .
56, indicat-ing a poor CR. While for κ = 0 . (a) , τ = 10 . CV min = 0 .
2, indicating high coherence.5
IG. 3. (Color online) Coefficient of variation CV againstnoise amplitude σ of layer 2 in isolation. Increasing thestrength of coupling force and the length of time delay shiftthe CV -curves up and to the right, thus weakening CR. Inparticular, in (b) , for κ = 1 . τ = 7 .
0, the red CV -curve has a minimum of CV min = 0 .
56 at σ = 4 . × − indicating a poor CR. In (d) , for κ = 1 . τ = 7 .
0, thepink CV -curve lies entirely above the line CV = 1 .
0, indi-cating the absence of CR. Parameters of layer 2: N = 25, n = 1, β = 0 . ε = 0 . α = 0 . κ = 0 . In Fig.3 (c) , we now have a fixed and short time delayof τ = 0 . κ . We ob-serve that the coupling force between the neurons doesnot have a significant effect on the coherence of the oscil-lations. The CV -curves at a short time delay of τ = 0 . CV min ≈ . σ ≈ . × − for both weak and strong couplingforces.In Fig.3 (d) , we fix the time delay at a relatively largevalue, τ = 7 .
0, and we notice that the same values of thecoupling force that could not change the values of CV inFig.3 (c) , can now very strongly affect the CV -curves.For long time delays, stronger coupling forces betweenthe neurons tend to decrease the coherence of the oscil-lations. In particular, for τ = 7 . κ = 1 .
5, the pink CV -curve does not show a pronounced non-monotonicbehavior and lies entirely above the line CV = 1 .
0, indi-cating the presence of highly incoherent oscillations andthe absence of CR.Now, we turn our attention to the effects of differentring network topologies on CR. The numerical simula-tions carried out in Fig.3 with a locally ( n = 1) cou-pled ring network topology were correspondingly carriedout with non-locally ( n = 3 and n = 6) and globally ( n = 12) coupled ring network topologies, all with a net-work size of N = 25. These simulations (not shown) de-picted the same qualitative behavior as seen in Fig.3, i.e.,CR deteriorating with stronger coupling forces and longertime delays. For a fixed coupling strength and time de-lay, we determined which of the ring network topologies(local, non-local, or global) will deteriorate CR the most.Fig.4 (a) shows that in a regime of weak coupling force( κ = 0 .
1) and short time delay ( τ = 1 . CV . Inthis regime, a change in the topology does not affect the CV curve and leaves the high coherence of the oscilla-tions unchanged for the network size considered. Here,for all the ring topologies: local ( n = 1), non-local( n = 3 , n = 12), the corresponding CV -curves have the same minimum value of CV min ≈ . σ = 0 . × − . The reason behind this behaviorlies in the fact that we have a weak ( κ →
0) couplinglimit between identical oscillators, in which case there isvery little interaction between oscillators (which behavein this limit as isolated oscillators). And because they areidentical oscillators, they all show a very similar dynamicirrespective of the ring network topology.The effects of the different ring network topologies onthe CV -curve become significant only when we are in aregime of strong coupling forces ( κ > .
5) and long timedelays ( τ > . (b) , we fix the strong couplingforce at κ = 1 . τ = 7 . CV -curves corresponding to dif-ferent ring network topologies are shifted up to higher(compared to the corresponding CV -curves in Fig.4 (a) with a weaker coupling and shorter delay) values as thering network becomes sparser. Firstly, this confirms thefact that irrespective of the ring network topology (lo-cally, non-locally or globally coupled), stronger couplingforces and longer time delays weaken CR which eventu-ally disappears at some sufficiently large values of thesecontrol parameters (see Fig.3 (d) ). Secondly, differentring network topologies can significantly affect the CV -curve only when the coupling forces are strong and thetime delays are long.In Fig.4 (b) , it can be clearly seen that the sparserthe network is, the higher the CV -curve, and hence, theless pronounced is CR. This behavior can be explainedby the fact that in a globally coupled network of iden-tical oscillators with strong coupling forces, the networkrelatively easily synchronizes, with the less coherent os-cillators (having high CV values) synchronizing to adoptthe dynamics of the more coherent ones (having low CV values). This has an overall effect of lowering the aver-aged CV of the entire network thus enhancing CR. Onthe other hand, as the network becomes less dense (i.e.,in a non-locally and especially locally coupled network),all the oscillators in the network can no longer so eas-ily synchronize (in particular, of course, those which arenot connected), hence, the averaged CV of the networkis calculated with high CV values (of less coherent os-6 IG. 4. (Color online) Coefficient of variation CV againstnoise amplitude σ for different (local, nonlocal, global) ringnetwork topologies in layer 2. (a): changing the topology( n ) has no effect on the high coherence of oscillations in aweak coupling and short time delay regime. (b): in a strongcoupling and long time delay regime, the topology can signif-icantly affect the coherence of the oscillations. In this regime,the sparser the network is, the less coherent the oscillationsare. Parameters of layer 2: N = 25, β = 0 . ε = 0 . α = 0 . κ = 0 . cillators which cannot easily synchronize) and low CV values (of more coherent oscillators). This has the over-all effect of shifting the averaged CV to higher values anddeteriorating CR as seen in Fig.4 (b) .A poor or non-existent CR in layer 2 can thereforeeasily be obtained by imposing a strong coupling forcewith a long time delay in a locally coupled ring networktopology. These are the three control parameter regimes( κ = 1 . . τ = 7 . n = 1) in which we will setlayer 2 before controlling (improve or induce) CR withmultiplexing, CR and SISR in layer 1.As we mentioned earlier, it has been shown in [45] withanother version of the FHN neuron model that strongercoupling forces rather shift the CV -curve to lower val-ues, thus enhancing CR. The version of the FHN neuronmodel used here is more general than that in [45] (see also[27] which uses the same version as in [45]) and shows arather opposite effect, i.e., stronger coupling forces dete-riorate CR. The reason behind this difference lies in thedetails of the equations of the two versions used. For ex-ample, in the version used in [27], the network of coupledneurons stays in the excitable regime when the couplingforce is weak and time delay is short. But as the strengthof the coupling force and the length of the time delay in-crease, the network eventually switches from an excitableregime into a delayed-coupling-induced oscillatory regimevia a saddle-node bifurcation. Therefore, increasing thecoupling force tends to bring their version of the FHNmodel closer to the bifurcation threshold, thus enhancingCR. In the current paper, only the slow variable equa-tion of the neurons is different from that in [27] and [45],but one can already see in Fig.1 (a) that increasing thestrength of the coupling force and the length of the timedelay within some interval does not switch the networkout of the excitable regime. As the coupling force be- comes stronger and the time delay longer, the networktends to move rather further away from the oscillatoryregime, i.e., away from the Hopf bifurcation threshold. InFig.2 (c) , we even see that strong coupling and long timedelays can switch the network from an oscillatory regimeinto a delayed-coupling-induced excitable regime; this ef-fect is the exact opposite of that produced by the modelused in [27] and [45]. In general, the further away we arefrom the (Hopf or saddle-node) bifurcation threshold, themore difficult it is for CR to occur. In such situations, weneed stronger noise amplitudes to have CR. This is whywe see in Fig.3 that as κ and τ increase, the CV -curvesachieve their minima at larger values of σ . Self-Induced Stochastic Resonance
All previous works on SISR have considered single iso-lated oscillators. The current work is thus the first to in-vestigate SISR in a network of coupled oscillators. Here,we investigate the effects of time delay ( τ ), couplingstrength ( κ ), and different ring network topologies ( n )on SISR in a layer network of coupled FHN oscillators.We compare the behavior of SISR to that of CR whenthese parameters are changed. Understanding the effectsof τ , κ , and n on SISR in layer 1 will guide us in con-trolling CR in layer 2 with SISR via the multiplexing ofboth layers.In the isolated layer 1 with σ = 0, and in the limitas ε →
0, each neuron equation in this layer reduces tocoupled Langevin equations of the form: dv i = − ∂U i ( v i , w i ) ∂v i dt + σ dW i ,U i ( v i , w i ) = v i − v i + v i w i − κ n i + n (cid:80) j = i − n (cid:16) v i ( t ) v j ( t − τ ) − v i ( t ) (cid:17) , (5)where the interaction potentials U i ( v i , w i ) ( i = 1 , ..., N )viewed as a function of v i with w i nearly constant,are double-well potentials. Fig.5 shows the landscapeof these interaction potentials.When w i < U i ( v i , w i ) is asymmetric with theshallower well at the left. In this case, the neuron is closeto the homogeneous stable fixed point at ( v ∗ i , w ∗ i ) =( − . , − . (cid:52) U li ( w i ) into the right well,see Fig.5 (a) .When w i > U i ( v i , w i ) is also asymmetric. In thiscase, the neuron has spiked and a return to the quiescentstate (the homogeneous stable fixed point) consists ofjumping over the right energy barrier (cid:52) U ri ( w i ) into theleft well, see Fig.5 (b) . When w i = 0, U i ( v i , w i ) issymmetric with (cid:52) U li ( w i ) = (cid:52) U ri ( w i ), in which casethe neuron is half way between the quiescent state andthe spike state, see Fig.5 (c) .7 IG. 5. (Color online) Landscapes of the interaction poten-tial U i ( v i , w i ) in Eq.5 for a locally ( n = 1) coupled ringnetwork topology with the energy barriers indicated in theasymmetric ( w i (cid:54) = 0) cases (a) and (b) and symmetric case (c) . The stronger the intra-layer coupling force κ is, thedeeper the energy barrier functions (cid:52) U li ( w i ) and (cid:52) U ri ( w i )are. In particular, when w i = 0, (cid:52) U li ( w i ) = (cid:52) U ri ( w i ),and both energy barriers achieve the deepest height at thestrongest coupling force κ . The saddle point and the leftand right minima of the interaction potential are located at v i = v ∗ m ( w i ), v i = v ∗ l ( w i ), and v i = v ∗ r ( w i ), respectively,i.e., v ∗ l ( w i ) < v ∗ m ( w i ) < v ∗ r ( w i ) We can also see in Fig.5 that the intra-layer couplingforce κ does not change the symmetry (or asymmetry)of the interaction potential. It only changes the depthof the energy barriers. The stronger κ is, the deeperthe energy barrier functions (cid:52) U li ( w i ) and (cid:52) U ri ( w i ) de-fined in Eq.7 are. Later, this particular behavior shallexplain why SISR is deteriorated by stronger intra-layercoupling forces. We choose the parameters of the coupledneurons in Eq.5 such that they satisfy the conditions nec-essary for the occurrence of SISR. These conditions whichare adapted from those of a single isolated FHN neuron[9, 19] to include the time-delayed coupling between theneurons coupled in a ring network topology are such that: lim ( ε ,σ ) → (0 , σ e ( ε − ) ∈ (cid:16) (cid:52) U li ( w ∗ i ) , F ( κ , τ , n ) (cid:17) , lim ( ε ,σ ) → (0 , σ e ( ε − ) = O (1) ,β − β h ( ε ) > ,σ = 0 , (6) where F ( κ , τ , n ) := (cid:110) ( κ , τ , n ) : (cid:52) U li ( w i ) = (cid:52) U ri ( w i ) (cid:111) , (cid:52) U li ( w i ) := U i (cid:0) v ∗ m ( w i ) , w i (cid:1) − U i (cid:0) v ∗ l ( w i ) , w i (cid:1) , (cid:52) U ri ( w i ) := U i (cid:0) v ∗ m ( w i ) , w i (cid:1) − U i (cid:0) v ∗ r ( w i ) , w i (cid:1) , (7)with v ∗ l,m,r ( w i ) := (cid:110) v i : v i − v i − w i + κ n i + n (cid:88) j = i − n (cid:16) v j ( t − τ ) − v i ( t ) (cid:17) = 0 ,v ∗ l ( w i ) < v ∗ m ( w i ) < v ∗ r ( w i ) (cid:111) . (8)The energy barrier functions (cid:52) U l,ri ( w i ) are obtainedfrom the potential U i ( v i , w i ) by taking the differ-ence between the potential function value at the sad-dle point v ∗ m ( w i ) and at the local minima v ∗ l,r ( w i )of the double potential U i ( v i , w i ), respectively [19]. (cid:52) U li (cid:0) w ∗ i (cid:1) (which has to be crossed to induce a spike)is the value of the left energy barrier function atthe w i -coordinate of the homogeneous stable fixedpoint (cid:0) v ∗ i ( κ , τ , n ) , w ∗ i ( κ , τ , n ) (cid:1) . This is where U li (cid:0) w ∗ i ( κ , τ , n ) (cid:1) gets its κ , τ , and n dependence.We note that n can only change the value of this fixedpoint because v ∗ i and w ∗ i depend on n . But changing n (or even the network size N ) does not take the networkout of the excitable regime as long as β − β h ( ε ) > (b) .Fig.6 shows the variation of CV against the noise am-plitude σ in layer 1 with a locally ( n = 1) coupled ringnetwork topology. In the numerical computations, wechoose ε = 0 . ε = 0 . ε → σ in Eq.(1) which is attached to the slow variable ( w i )equations is set to zero in the conditions necessary forSISR in Eq.(6). We recall that for the FHN model, SISRrequires the noise term to be attached to the fast variable( v i ) equation, while CR requires the noise term to be at-tached to the slow variable ( w i ) equation [9]. Therefore,in our study, the noise amplitude σ will only be used inthe multiplex network setting, when we want to investi-gate the enhancement of CR in layer 2 by the occurrenceof CR in layer 1. This will later allow us to compare theenhancement capabilities of CR and SISR.In Fig.6, the general behavior of SISR is the follow-ing: In the regime of weak coupling forces κ and shorttime delays τ , the CV -curves remain very low, indicat-ing a high degree of coherence of the oscillations. Thishigh degree of coherence is gradually lost as one moves tostronger coupling forces and longer time delay. This gen-eral behavior has also been observed with CR in Fig.3.The reason for this behavior for the case SISR is the sameas for CR, i.e., where larger values of κ and τ bring the8ystem further away from the Hopf bifurcation thresh-old, making the noise-induced oscillations to be less co-herent, as a larger noise amplitude is now required toinduce oscillations. Moreover, the CV -curves of SISR inFig.6 show oscillations that are much more coherent fora wider range of noise amplitude σ than the CV -curvesof CR in Fig.3, which are rather relatively higher withthe minimum occurring at single value of the noise am-plitude σ . This confirms the more robust nature of SISRcompared to CR.The response of the SISR to changes in the couplingforce and time delay in Fig.6 can also be explained interms of the interaction potential U i ( v i , w i ) in Fig.5.We observe in Fig.5 that for a fixed ring network topol-ogy and time delay, as the coupling force κ increases,the energy barriers (cid:52) U li ( w i ) and (cid:52) U ri ( w i ) becomedeeper. In particular, when w i <
0, the trajectory isin the left potential well and as κ becomes stronger(0 . , . , . , . (cid:52) U li ( w i ) be-comes deeper (hence the trajectory at the bottom of thewell get closer to the homogeneous stable fixed point at w ∗ i = − . (cid:52) U li ( w i ) is (in other words, the stronger the cou-pling force κ is), the closer is the trajectory to the stablefixed point and the further away is the system from theoscillatory regime.For the trajectory to jump over this deep left energybarrier, a stronger noise amplitude is thus needed. Thisis why in Fig.6 as κ increases, the left branch of the CV -curve is shifted to the right, meaning that strongernoise amplitudes are required to induced frequent spik-ing (i.e., frequent escape from the deep left energy bar-rier). Furthermore, we can see (from Fig.6 (d) ) that atlonger time delay τ , this effect (the shifting of the leftbranch of the CV -curve to the right) is more pronouncedthan in Fig.6 (c) with a shorter time delay. This is be-cause in Eq.5, the longer the time delay is ( τ (cid:29) (cid:2) v i ( t ) v j ( t − τ ) − v i ( t ) (cid:3) from zero (given that the oscillators are identical), andhence, the stronger the effect of the coupling force κ on the interaction potential, the energy barriers func-tions and consequently, on the CV -curves. Otherwise,if τ →
0, then because the oscillators are identical, (cid:2) v i ( t ) v j ( t − τ ) − v i ( t ) (cid:3) →
0, and κ will have littleeffect on the interaction potential, the energy barriers,and on the CV -curves. This is why the coupling force κ has a stronger effect on SISR only when τ gets longer,and vice versa.In Fig.6 (a) , for a fixed and weak coupling force κ =0 .
1, all the flat-bottom CV -curves obtained with differenttime delays ( τ = 0 . , τ = 5 . , τ = 10 . τ = 20 .
0) showa deep and broad minimum, indicating a high degree ofcoherence for a wide range of the noise amplitude. Wenotice that for this weak coupling force, as τ increases,the CV -curve is not shifted to higher values as in thecase of CR in Fig.3 (a) . Here, only the left branch ofthe CV -curve is significantly shifted to the right. Thismeans that in the weak coupling regime, the coherence of FIG. 6. (Color online) Coefficient of variation CV againstnoise amplitude σ of layer 1 in isolation. Like CR in Fig.3,SISR is weakened by stronger coupling κ and longer delays τ . SISR is, however, much less sensitive to these parametersthan CR. The flat-bottom of the CV -curves indicates thata wider range of noise amplitude can induce coherent oscil-lations during SISR, unlike CR in Fig.3 with round-bottom CV -curves. For the same time-delayed couplings, oscillationdue to SISR are more coherent (lower CV -curves) than thosedue to CR. In (b) and (d) with a strong coupling force andlong time delay regime, the maximum noise amplitude belowwhich SISR occurs is controllable, but fixed in (a) and (c) with a weak coupling force and short time delay regime. Pa-rameters of layer 1: N = 25, n = 1, β = 0 . ε = 0 . α = 0 . κ = 0 . . the oscillations due to SISR is not affected as the delaybecomes longer, but the coherence is achieved only atlarger noise amplitudes σ . In Fig.6 (a) , we have the sameminimum value of CV min ≈ .
012 for: τ = 0 . σ ∈ (cid:0) . × − , . × − (cid:1) ; τ = 5 . σ ∈ (cid:0) . × − , . × − (cid:1) ; τ = 10 . σ ∈ (cid:0) . × − , . × − (cid:1) ; and τ = 20 . σ ∈ (cid:0) . × − , . × − (cid:1) .Notice that the lower bound of the intervals increases asthe delay increases while the upper bounds are almostfixed.In Fig.6 (b) , we switch to a strong coupling ( κ = 1 . CV -curves into asmaller noise interval, while shifting the curves to highervalues, thus decreasing the coherence of the oscillation.Here, we have different noise intervals for different min-ima of CV . We have CV min = 0 .
014 at τ = 0 . σ ∈ (cid:0) . × − , . × − (cid:1) ; CV min = 0 .
071 at τ = 2 . σ ∈ (cid:0) . × − , . × − (cid:1) ; and in the last twocases, the noise intervals in which we have the highestcoherence have contracted to points with CV min = 0 . σ = 2 . × − for τ = 5 .
0; and CV min = 0 .
51 at σ = 8 . × − for τ = 7 . (b) , firstly, we notice that in the strong cou-pling regime, the effect of the time delay on the coher-ence of the oscillations becomes significant unlike whenthe coupling is weak as in Fig.6 (a) . In Fig.6 (b) , at only τ = 7 . CV min is already above 0.5, whereas in theweak coupling regime in Fig.6 (a) , even at τ = 20 .
0, wehave CV min ≈ . (a) and (c) , respectively), theupper bound of the noise interval for which the CV -curves achieve their minima is almost fixed. Here, onlythe lower bound of the noise intervals is shifted to theright. Whereas in the strong coupling and long time de-lay regimes (i.e., in Fig.6 (b) and (d) , respectively), boththe lower and upper bounds of the noise intervals areshifted respectively to the right and to the left as the τ and κ increase. This has the overall effect of shrinkingthe noise interval in which the CV -curves achieve theirminima to a single value of σ .To explain this behavior, we obtain from the first con-dition of Eq.6 the minimum and maximum noise ampli-tudes between which SISR occurs as: σ min = (cid:115) (cid:52) U li (cid:0) w ∗ i ( κ , τ , n ) (cid:1) log e ( ε − ) ,σ max = (cid:115) F (cid:0) κ , τ , n (cid:1) log e ( ε − ) . (9)We observe that σ min depends on the fixed point co-ordinate w ∗ i ( κ , τ , n ) which in turn also depends on κ , τ , and n . Therefore, changing κ , τ , and n willchange the value of w ∗ i ( κ , τ , n ) which will in turnchange the value of σ min via the energy barrier func-tion (cid:52) U li ( w ∗ i ). Numerical computations show that σ min increases as κ and τ increase (see Fig.6) and as n de-creases (see Fig.7 (b) ).On the other boundary, σ max does not depend on thecoordinates of the homogeneous stable fixed point, buton the complicated function F ( κ , τ , n ), completely de-fined by Eq.7 and Eq.8. In Fig.6 (a) and (c) (i.e., inthe weak coupling and short time delay regimes, respec-tively), we notice that σ max = O (10 − ) is almost fixedfor all the values of the time delay and coupling strengthused. This fixation of the upper bound of the noise inter-val in which SISR occurs was already observed in a singleisolated FHN neuron [19]. In the case of a single isolatedFHN neuron, the function in Eq.7 does not depend onany system parameter and takes a simple constant value F = . This implies (for a fixed ε = 0 . σ max = (cid:0) / · log e ( ε − ) (cid:1) / .In the current case of coupled FHN neurons, thefixation of the upper bound of the noise interval in which SISR occurs can therefore only be observed if F ( κ , τ , n ) → C , where C is a constant. In particular,in a weak coupling ( κ →
0) regime and short time delay( τ →
0) regime (or more precisely, v j ( t − τ ) − v i ( t ) → τ →
0, because all the oscillators are identical), F ( κ , τ , n ) → . In these regimes (see Fig.6 (a) and (c) ), we observe that σ max ≈ − , that is the value ob-tained in [19] for the case of a single isolated ( κ = 0)FHN neuron.In the strong ( κ (cid:29)
0) coupling and long ( τ (cid:29) ⇒ v j ( t − τ ) − v i ( t ) (cid:54) = 0) time delay regimes (Fig.6 (b) and (d) ), the function F in Eq.7 is now strongly modified by n and the large values of κ and τ . This is why in theseregimes, the upper bound σ max of the noise interval inwhich SISR occurs is not longer fixed, but shifted to theleft as κ and τ take on larger values. This left-shifting ofthe upper bound of the noise interval cannot occur in anisolated FHN neuron because it owes its occurrence to the(strong) coupling and (long) time delays in a network ofcoupled FHN neurons. Moreover, the squeezing of lowerand upper bounds of the noise interval in which SISRoccurs in Fig.6 (b) and (d) explains the round-bottomshape of the CV -curves (in contrast to the flat-bottom CV -curves in Fig.6 (a) and (c) ).In Fig.6, we have a locally ( n = 1) coupled ring net-work topology. The numerical simulations (not shown)for non-locally ( n = 3 ,
6) and globally ( n = 12) coupledring network topologies displayed the same qualitativebehavior observed in Fig.6 when κ and τ are varied. InFig.7, we fixed κ and τ and we investigate the effectsof different ring network topologies on SISR.Fig.7 (a) shows that in a weak ( κ = 0 .
1) coupling andshort ( τ = 1 .
0) time delay regime, changing the ringnetwork topology has no effects on the CV -curves whichremain low with CV min = 0 .
012 for ( σ min , σ max ) =(4 . × − , . × − ) for the system of N = 25 used.This behavior of SISR is similar to that CR in Fig.4 (a) except that the oscillations are less coherent for CR with CV min = 0 .
14. The reasons for this behavior is the sameas the one given for the case of CR Fig.4 (a) .For the same reasons given for the case of CR inFig.4 (b) , the network topology significantly affects SISRonly when the coupling and the time delay betweenthe neurons take on larger values. This is observed inFig.7 (b) ( κ = 1 . τ = 7 .
0) which shows that thesparser the network is, the higher the CV -curve is, andhence, the less coherent the oscillations due to SISRare. Notice further that in Fig.7 (b) , irrespective of thering network topology, the CV -curves all have a round-bottom shape due to the shrinking (on both ends) of thenoise interval ( σ min , σ max ), caused by the strong cou-pling force and the long time delay.The behavior of the CV -curves in Fig.7 (a) and (b) can also be explained in terms of the interaction poten-tial in Eq.5. When the coupling force becomes weaker κ → τ →
0, the inter-action part of the potential in Eq.5 turns to zero (i.e.,10 n i + n (cid:80) j = i − n (cid:2) v i ( t ) v j ( t − τ ) − v i ( t ) (cid:3) →
0) and differentvalues of n (i.e., different ring network topologies) can-not significantly affect this limit and hence cannot alsosignificantly affect the CV -curves as seen in Fig.7 (a) .Furthermore, in these limits (i.e., κ → τ → CV values inFig.7 (a) .On the other hand, in the strong coupling limit( κ (cid:29)
0) and long time delay limit ( τ (cid:29) κ n i + n (cid:80) j = i − n (cid:2) v i ( t ) v j ( t − τ ) − v i ( t ) (cid:3) , does not longerturn to zero and different values of n can now signif-icantly affect the energy barrier functions. As n de-creases (i.e., as the network becomes sparser), the quan-tify κ n (for a fixed κ (cid:29)
0) increases, which means theenergy barriers become deeper (see Fig.5 with larger val-ues of κ ). And as we saw earlier, to be able to jumpover deep energy barriers, the noise amplitude should belarge enough, and large noise amplitudes deteriorate thecoherence of the oscillations. This is what we observe inFig.7 (b) .Another behavior we can observe is that the CV -curvesof SISR in Fig.6 and Fig.7 are rougher than those ofCR in Fig.3 and Fig.4, even though the coherence of theoscillations due to SISR is higher than those due to CR.The reason for this behavior lies in the fact that CV -curves are computed using the interspike interval of thefast variable ( v i ) equations. Because in the case of SISR,the noise term is attached to fast variable equations, thenoise affect the fast variable v i directly and more thanthe slow variable w i which is indirectly affected. Hencethe roughness of the CV -curves. While in the case ofCR, even though the noise term is instead attached tothe slow variable ( w i ) equations, the CV -curves are stillcomputed using the interspike interval of the fast variablewhich is in this case, not directly affected by the noise.Hence, we have relatively smoother CV -curves with CRthan with SISR. B. Multiplex Network and control of CR
It has been shown in [45] that in a two-layer multiplexnetwork (without time delays), CR in one layer could in-duce and enhance CR in the other layer, even with a weakmultiplexing. Here, we address the question whether acontrol of CR based on multiplexing is still possible whenwe use a different noise-induced phenomenon, namely,SISR instead of CR. Then, we address a second ques-tion: which of the enhancement strategies, CR or SISR,of CR in one layer of the multiplex network is better?To answer these questions, we consider layer 1 andlayer 2, each with a network size of N = 25, coupled FIG. 7. (Color online) Coefficient of variation CV againstnoise amplitude σ for different (local, nonlocal, global) ringnetwork topologies. (a): changing the topology ( n ) has noeffect on the high coherence of oscillations in a weak couplingand short time delay regime. (b): in a strong coupling andlong time delay regime, the ring topology can significantlyaffect the coherence of the oscillations. In this regime, thesparser the network is, the less coherent the oscillations dueto SISR are. Parameters of layer 1: N = 25, β = 0 . ε = 0 . α = 0 . κ = 0 . in a multiplex fashion as in Fig.1. The coupled layersare non-identical, as we set layer 2 in a locally coupledring network topology with a strong intra-layer couplingforce, a long intra-layer time delay, while layer 1 has aglobally coupled ring topology with a weak coupling forceand a short time delay. With this choice of parameters,we set layer 2 in a regime where CR is either poor (i.e.,0 . < CV min < . n = 1, κ = 1 . τ = 7 .
0, see thered curve in Fig.3 (b) ) or non-existent (i.e., CV min > . n = 1, κ = 1 . τ = 7 .
0, see the pink curve inFig.3 (d) ). In layer 1, we set the parameters such that CR(see Fig.4 (a) ) or SISR (see Fig.7 (a) ) is very pronounced,i.e., n = 12, κ = 0 . τ = 1 . τ and the inter-layer coupling force κ that characterizethe coupling between the layers. Numerical investiga-tions have shown that the more pronounced CR or SISRin layer 1 is, the more pronounced CR induced in layer2 would be. Therefore, to maximize the enhancement orinduction of a poor or non-existent CR in layer 2, we usethe intra-layer control parameters of layer 1 which induceone of the most pronounced forms of CR or SISR in thislayer, i.e., κ = 0 . τ = 0 .
5, and n = 12.We now investigate the impacts of the inter-layer cou-pling force κ , the inter-layer time delay τ , and a pro-nounced SISR in layer 1 on a poor or non-existent CR inlayer 2. In Fig.8 (a) , we see that at a weak ( κ = 0 . τ = 0 .
5) inter-layer time delay, the poor CR in layer 2 cannot be en-hanced by a pronounced SISR in layer 1. In isolation( κ =0.0), layer 2 has CV min = 0 .
56 at σ = 4 . × − (see red curve in Fig.8 (a) ), and for a weak inter-layercoupling force ( κ = 0 . CV -curvewith a rather higher minimum value, i.e., CV min = 0 . IG. 8. (Color online) Coefficient of variation CV againstnoise amplitude σ of layer 2 coupled to layer 1 in a multi-plexed network. In (a) and (b) (black curves), when thereis a weak inter-layer coupling ( κ = 0 .
10) and a short inter-layer time delay ( τ = 0 . κ = 0 .
75) (green curvesin (a) and (b) ) and a short inter-layer time delay ( τ = 0 . ε = 0 . κ = 0 . τ = 0 . n = 12, σ = σ , β = 0 . α = 0 . N = 25. Parameters of layer 2: ε = 0 . τ = 7 . n = 1, β = 0 . α = 0 . N = 25, κ = 1 . (a) , κ = 1 . (b) . at σ = 4 . × − (see black curve in Fig.8 (a) ).This inability of SISR to enhance CR in a weak inter-layer coupling regime is also observed in Fig.8 (b) , whereCR is non-existent. Here, layer 1 in isolation ( κ = 0 . CV -curve with CV min = 1 .
14 at σ = 6 . × − (see pink curve) against CV min = 1 .
18 at σ = 3 . × − (see black curve) with a weak coupling force of κ = 0 . τ = 0 .
5. Therefore, at a weakinter-layer coupling force, a pronounced SISR in layer 1,instead of improving CR in layer 2, can actually makeCR even poorer than when layer 2 is isolated.On the other hand, strong inter-layer coupling forcesdrastically reverse the scenario. From the green curvesin Fig.8 (a) and (b) , the inter-layer coupling force takesa larger value, i.e., κ = 0 .
75 and the time delay is keptfixed at τ = 0 .
5. We observe that the CV -curves oflayer 2 drop significantly to low values, indicating anenhancement of a poor CR and an induction of a non-existent CR in layer 2 by a pronounced SISR in layer1. In Fig.8 (a) , the relatively high values of the minimaof the CV -curves of layer 2 (i.e., CV min = 0 .
56 withno ( κ = 0 .
0) multiplexing and CV min = 0 .
67 witha weak ( κ = 0 .
10) multiplexing) are shifted to a lowvalue of CV min = 0 .
21 at σ = 6 . × − , indicating ahigh coherence of the oscillations in layer 2 induced by astrong inter-layer coupling force and a pronounced SISRin layer 1. We can also see in Fig.8 (b) that the very highvalues of the minima of the CV -curves of layer 2 (i.e., CV min = 1 .
14 with κ = 0 . CV min = 1 .
18 with κ = 0 .
10) are shifted to a significantly low value, i.e., CV min = 0 .
23 at σ = 2 . × − by a strong inter-layer coupling force and a pronounced SISR in layer 1.Now, in the last part of our study, we compare theperformance of SISR (in layer 1) in enhancing CR inlayer 2 to the performance of CR (in layer 1) in enhancingCR in layer 2. In the control of CR by SISR (which wehenceforth refer to as SISR-CR control scheme), we setlayer 1 with a pronounced SISR, i.e., we choose a globallycoupled ring network topology ( n = 12), a weak intra-layer coupling force ( κ = 0 . τ = 0 . σ (cid:54) = 0 and σ = 0.In the control of CR by CR (which we henceforth re-fer as to CR-CR control scheme), we set layer 1 suchthat it also has a pronounced CR with of course, thesame ring network topology, the same intra-layer cou-pling force, and the same intra-layer time delay as in theSISR-CR control scheme, i.e., we choose a globally cou-pled ring network topology ( n = 12), a weak intra-layercoupling force ( κ = 0 .
1) and a short intra-layer timedelay ( τ = 0 . σ = 0) on the fast variable ( v i ) equations in layer1 and switch on the noise term ( σ (cid:54) = 0) on the slowvariable ( w i ) equations of this layer. We recall that forthe FHN neuron model, SISR requires the noise term tobe attached to the fast variable equation, while CR re-quires the noise term to be attached to the slow variableequation. Furthermore, because CR (unlike SISR with ε = 0 . ε →
0) to occur, we set the singular parameter ε to itsstandard value ε = 0 . (cid:28) ε . We further note that thecoherence of the noise-induced oscillations due to CR isinsensitive against variations of the timescale separationratio between the fast and slow variables of the excitablesystem [7, 9].In order to control CR in layer 2, we choose an intra-layer coupling force, intra-layer time delay, and a ringnetwork topology such that CR is non-existent in thislayer when it is in isolation, i.e., we choose κ = 1 . τ = 7 . n = 1, see in Fig.3 (d) . In the SISR-CR andCR-CR control schemes, layer 1 and layer 2 are coupledin a multiplex fashion in two different inter-layer couplingregimes. In the first regime, we choose weak inter-layercoupling forces ( κ = 0 . κ = 0 .
4) and a shortinter-layer time delay ( τ = 1 . κ = 0 . κ = 1 .
0) with the same short time delay ( τ = 1 . (a) with κ = 0 . τ = 1 .
0, we see thatthe CR-CR control scheme (blue curve) performs muchbetter than the SISR-CR control scheme (red curve).Here, the CR-CR control scheme provides a minimum CV of CV min = 0 .
21 at σ = 4 . × − against CV min = 1 .
22 at σ = 4 . × − for the SISR-CR controlscheme. This ability of CR (in one layer of a multiplexnetwork) to enhance CR (in the other layer) in a weakinter-layer coupling regime was already shown in [45].Our work further confirms their finding on a different ver-12 IG. 9. (Color online) Coefficient of variation CV againstnoise amplitude σ of layer 2 coupled to layer 1 in a mul-tiplexed network. In (a) and (b) , we have a weak inter-layer coupling regime, where the CR-CR control scheme per-forms better than the SISR-CR control scheme in enhancingCR in layer 2 of the multiplex. In (c) and (d) , we have astrong inter-layer coupling regime, where the SISR-CR con-trol scheme takes over and performs better than the CR-CRcontrol scheme. Parameters of layer 1 with CR: ε = 0 . κ = 0 . τ = 0 . n = 12, σ = σ , β = 0 . α = 0 . N = 25. Parameters of layer 1 with SISR: ε = 0 . κ = 0 . τ = 0 . n = 12, σ = σ , β = 0 . α = 0 . N = 25. Parameters of layer 2: ε = 0 . κ = 1 . τ = 7 . n = 1, β = 0 . α = 0 . N = 25. sion of the FHN model. This behavior is however, absentwith SISR which performs very poorly in the SISR-CRcontrol scheme in a weak inter-layer coupling regime.In Fig.9 (b) , the inter-layer coupling force is increasedto κ = 0 . τ = 1 .
0. The SISR-CR control scheme improveson its performance which is, however, still poorer thanthat of the CR-CR control scheme. Here, the SISR-CRcontrol scheme provides a minimum CV of CV min = 0 . σ = 2 . × − against CV min = 0 .
23 at σ = 4 . × − for the CR-CR control scheme.Surprisingly, in Fig.9 (c) , where we have a strong inter-layer coupling force ( κ = 0 .
5) and the same short inter-layer time delay ( τ = 1 . CV of CV min = 0 .
22 at σ = 1 . × − against CV min = 0 . σ = 2 . × − for the CR-CR control scheme.In Fig.9 (d) , we further increase the strength of the FIG. 10. (Color online) Minimum coefficient of variation( CV min ) against the inter-layer coupling force κ and timedelay τ for layer 2 when coupled to layer 1 in a multiplexednetwork. In (a) and (b) , we have a 3D plot and its 2D pro-jection onto the parameter plane ( κ − τ ) showing the en-hancement performance of the CR-CR control scheme, whichis good even at weak (0 . < κ < .
5) inter-layer couplingforces. In (c) and (d) , we have a 3D plot and its 2D projectiononto the parameter space ( κ − τ ) showing the enhance-ment performance of the SISR-CR control scheme, which ispoor at weak (0 . ≤ κ < .
5) inter-layer coupling forces.But at strong (0 . ≤ κ ≤ .
0) inter-layer coupling forces,the SISR-CR control scheme has a good performance which iseven better than that of the CR-CR control scheme, especiallyat longer inter-layer time delays and weaker noise amplitudes.Parameters of layer 1 with CR: ε = 0 . κ = 0 . τ = 0 . n = 12, σ = σ , β = 0 . α = 0 . N = 25. Parameters oflayer 1 with SISR: ε = 0 . κ = 0 . τ = 0 . n = 12, σ = σ , β = 0 . α = 0 . N = 25. Parameters of layer 2: ε = 0 . κ = 1 . τ = 7 . n = 1, β = 0 . α = 0 . N = 25. inter-layer coupling force to κ = 1 . κ =0 . CV of CV min = 0 .
20 at σ = 9 . × − against CV min =0 .
30 at σ = 5 . × − for the CR-CR control scheme.In Fig.10, we globally compare the performances of theSISR-CR and CR-CR control schemes in the ( κ − τ )parameter space. In Fig.10 (a) and (b) , we have theCR-CR control scheme with a color coded minimum CV ( CV min ). We can see that for very weak inter-layer cou-pling regime (0 ≤ κ < .
2) the CR-CR control schemeperforms poorly as CV min is very high in value (redcolor). But as soon as κ ≥ . CV min drops to dark13lue indicating a good performance of the CR-CR controlscheme at weak inter-layer coupling forces. Furthermore,we observe that as the inter-layer time delay increases, CV min also increases. This effect of the inter-layer timedelay is however more significant in stronger inter-layercoupling regimes i.e., κ ≥ . (c) and (d) we have the SISR-CR controlscheme. We can see that for weak inter-layer couplingregime (0 ≤ κ < . CV min takes red andgreen colors indicating a relatively higher CV min com-pared to those of the CR-CR control scheme. However, inthe strong inter-layer coupling regime (0 . ≤ κ ≤ . CV min takes a dark blue (darker thanthe blue color in the CR-CR control scheme, especiallyat long time delays) color, indicating a high (higher thanthose of the CR-CR control scheme) enhancing perfor-mance of the SISR-CR control scheme. Moreover, wecan also observe (like in the CR-CR control scheme)that, as the inter-layer time delay increases, the perfor-mance of the SISR-CR control scheme reduces with theeffect becoming more significant at stronger ( κ ≥ . (a) - (b) andFig.6 (a) - (b) . IV. SUMMARY AND CONCLUSIONS
In this paper, we have considered a two-layer time-delayed multiplex network of stochastic FHN neurons inthe excitable regime. First, we systematically and in-dependently investigated the impacts of the intra-layertime-delayed couplings and different ring network topolo-gies on the noise-induced phenomena of CR and SISR inthe isolated layers of the multiplex network. We deter-mined the parameter regimes and the network topologywhich enhance or weaken CR and SISR the most. Itwas found that strong intra-layer coupling forces, longintra-layer time delays, and a locally coupled ring net-work topology tend to weaken CR and SISR, with CRshowing a more much sensitive behavior to variations ofthese control parameters than SISR. Furthermore, in theisolated layer network of FHN neurons exhibiting SISR,we found that the maximum noise amplitude at which SISR occurs is not fixed, but controllable, especially in astrong intra-layer coupling force and long intra-layer timedelay regime. This is in contrast to SISR in a single iso-lated FHN neuron, where this maximum noise amplitudeis fixed and uncontrollable.Secondly, we set up two control strategies of CR inlayer 2, based not only on the multiplexing with layer 1,but also on SISR and CR. In one of the control schemes(which we termed as CR-CR control scheme), we have apronounced ( n = 12, κ = 0 . τ = 0 .
5) CR in layer 1which controls CR in layer 2. While in the second controlscheme (SISR-CR control scheme), we have a pronounced( n = 12, κ = 0 . τ = 0 .
5) SISR in layer 1 whichcontrols CR in layer 2. In our control schemes, we haveconsidered layer 2 of the multiplex network with a non-optimal network topology and time-delayed coupling forCR, i.e., the intra-layer time-delayed coupling of layer 2was such that CR is non-existent when the layer is inisolation (i.e., n = 1, κ = 1 . τ = 7 . . ≤ κ < .
5, the CR-CR control scheme couldsignificantly enhance (after inducing a previously non-existent) CR in layer 2. While in the weak inter-layercoupling regime (0 . ≤ κ < . . ≤ κ < .
3) coupling forces, induce CR.But at strong inter-layer coupling force ( κ ≥ . ACKNOWLEDGMENTS
Marius E. Yamakou thanks the Max Planck Institutefor Mathematics in the Sciences, Leipzig, Germany. [1] HuGang et al., Phys. Rev. Lett.
807 (1993).[2] A.S. Pikovsky and J. Kurths, Phys. Rev. Lett. , 775(1997).[3] B. Lindner and L. Schimansky-Geier, Phys. Rev. E ,7270 (1999).[4] B. Lindner et al., Phys. Rep. , 321 (2004).[5] M.I. Freidlin, J. Stat. Phys. , 283 (2001). [6] M.I. Freidlin, Stoch. Dyn. , 261 (2001).[7] C.B. Muratov et al., Physica D , 227 (2005).[8] R. FitzHugh, Biophys. J., , 445 (1961).[9] R.E. Lee DeVille, E. Vanden-Eijnden, C.B. Muratov,Phys. Rev. E , 031105 (2005).[10] C.B. Muratov and E. Vanden-Eijnden, Chaos , 015111(2008).
11] A. Neiman, P.I. Saparin, L. Stone, Phys. Rev. E , 270(1997).[12] V. Beato et al., Philos. Trans. R. Soc. A , 381 (2007).[13] J. Hizanidis and E. Sch¨oll, Phys. Rev. E , 066205(2008).[14] Z.-Q. Liu et al., Physica A , 2642 (2010).[15] B. Jia et al., Chin. Phys. Lett., , 090507 (2011).[16] H.G. Gu et al., Int. J. Mod. Phys. B , 3977 (2011).[17] R.E.L. Deville and E. Vanden-Eijnden, J. Stat. Phys. , 75 (2007).[18] R.E.L. Deville and E. Vanden-Eijnden, Commun. Math.Sci. , 431 (2007).[19] M.E. Yamakou and J. Jost, Nonlinear Dyn. , 2121(2018).[20] M.E. Yamakou and J. Jost, Europhys. Lett. , 18002(2017).[21] J. Shen, L. Chen, and K. Aihara, The 4th InternationalConference on Computational Systems Biology, pp. 251-257 (ISB 2010).[22] N.B. Janson, A.G. Balanov, E. Sch¨oll, Phys. Rev. Lett. , 010601 (2004).[23] J. Hizanidis, A. Balanov,A. Amann, E. Sch¨oll, Phys. Rev.Lett. , 244104 (2006).[24] R. Aust et al., Eur. Phys. J. ST , 77 (2010).[25] P.M. Geffert et al., Eur. Phys. J. B , 291 (2014).[26] V. Semenov et al., Chaos , 033111 (2015).[27] M. Maria et al., Chaos , 101102 (2017).[28] C. Zhou, J. Kurths, B. Hu, Phys. Rev. Lett. , 114320 (2017).[30] M. De Domenico, Gigascience , 1 (2017).[31] A.N. Pisarchik et al., Biol. Cybern. , 397 (2014).[32] A.V. Andreev et al., Chaos Solitons Fractals , 80(2018).[33] L.V. Gambuzza et al., Europhys. Lett. , 20010 (2015).[34] A. Singh et al., Europhys. Lett. , 30010 (2015).[35] I. Leyva et al., Sci. Rep. , 45475 (2017).[36] L. Zhang, A.E. Motter, T. Nishikawa, Phys. Rev. Lett. , 174102 (2017).[37] R.G. Andrzejak et al., Chaos , 053114 (2017).[38] N.E. Kouvaris et al., , Sci. Rep. , 10840 (2015).[39] S. Ghosh and S. Jalan, Int. J. Bifurcation Chaos ,1650120 (2016).[40] S. Ghosh et al., Europhys. Lett. , 60005 (2016).[41] V.A. Maksimenko et al., Phys. Rev. E , 052205 (2016).[42] B.K. Bera et al., Europhys. Lett. , 10001 (2017).[43] A. Bukh et al., Chaos , 111102 (2017).[44] S. Ghosh et al., Chaos Solitons Fractals , 56 (2018).[45] N. Semenova and A. Zakharova, Chaos , 051104(2018).[46] E.M. Izhikevich, Int. J. Bifurc. Chaos , 1171 (2000).[47] E. Sch¨oll et al., Philos. Trans. R. Soc. A , 1079 (2009).[48] D.V.Ramana Reddy, A. Sen, G.L. Johnston, Phys. Rev.Lett. , 5109 (1998).[49] D.V.Ramana Reddy, A. Sen, G.L. Johnston, Phys. Rev.Lett. , 3381 (2000).[50] V. Resmi, G. Ambika, R.E. Amritkar, Phys. Rev. E ,046212 (2011).[51] A. koseska et al., Phys. Rep. , 173 (2013).[52] C. Kurrer and K. Schulten, Phys. Rev. E , 6213 (1995).[53] X. Pei, L. Wilkens, F. Moss, Phys. Rev. Lett. , 4679(1996).[54] M. Masoliver et al., Chaos , 101102 (2017).[55] N. Klasdin, J. Guid. Control Dyn. , 114 (1995). [56] B. Hu and C. Zhou, Phys. Rev. E , R1001(R) (2000)., R1001(R) (2000).