Transient spatiotemporal chaos in the Morris-Lecar neuronal ring network Keegan Keplinger and Renate Wackerbauer Citation: Chaos: An Interdisciplinary Journal of Nonlinear Science 24, 013126 (2014); doi: 10.1063/1.4866974 View online: http://dx.doi.org/10.1063/1.4866974 View Table of Contents: http://scitation.aip.org/content/aip/journal/chaos/24/1?ver=pdfcov Published by the AIP Publishing Articles you may be interested in Synaptic plasticity modulates autonomous transitions between waking and sleep states: Insights from a MorrisLecar model Chaos 21, 043119 (2011); 10.1063/1.3657381 Hyperlabyrinth chaos: From chaotic walks to spatiotemporal chaos Chaos 17, 023110 (2007); 10.1063/1.2721237 Complex dynamics in simple Hopfield neural networks Chaos 16, 033114 (2006); 10.1063/1.2220476 Synchronization and propagation of bursts in networks of coupled map neurons Chaos 16, 013113 (2006); 10.1063/1.2148387 Random dynamics of the Morris–Lecar neural model Chaos 14, 511 (2004); 10.1063/1.1756118

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

CHAOS 24, 013126 (2014)

Transient spatiotemporal chaos in the Morris-Lecar neuronal ring network Keegan Keplingera) and Renate Wackerbauerb) Department of Physics, University of Alaska, Fairbanks, Alaska 99775-5920, USA

(Received 15 September 2013; accepted 14 February 2014; published online 4 March 2014) Transient behavior is thought to play an integral role in brain functionality. Numerical simulations of the firing activity of diffusively coupled, excitable Morris-Lecar neurons reveal transient spatiotemporal chaos in the parameter regime below the saddle-node on invariant circle bifurcation point. The neighborhood of the chaotic saddle is reached through perturbations of the rest state, in which few initially active neurons at an effective spatial distance can initiate spatiotemporal chaos. The system escapes from the neighborhood of the chaotic saddle to either the rest state or to a state of pulse propagation. The lifetime of the chaotic transients is manipulated in a statistical sense through a singular application of a synchronous perturbation to C 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4866974] a group of neurons. V

Transient dynamics may often be more relevant than the asymptotic states of a system in terms of observation and modeling.1 Transient behavior in the brain might play an integral role in perceptual phenomena associated with stimuli and in cognitive function. Brain information processing is observed as sequential switching from one metastable state to another, and the corresponding transients could result from sequential saddles.2 We show that spatiotemporal chaos in a ring network of diffusively coupled excitable Morris-Lecar neurons is transient, and either spontaneously collapses to the rest state, where all neurons are inactive, or to a state of pulse propagation. Such transients due to the existence of a chaotic saddle can be long-lived and demonstrate a dramatic systemintrinsic change in the absence of an external influence, as opposed to transients to an attractor that are typically short and do not exhibit chaotic features.1

I. INTRODUCTION

Spatiotemporal chaos is a generic pattern in extended non-equilibrium systems exhibiting a rapid decay of spatial and temporal correlations.3 Asymptotic spatiotemporal chaos is reported, for example, in fluid experiments,4 chemical reaction-diffusion systems,5 and in cardiac fibrillation.6 Spatiotemporal chaos is transient in systems with a nonattracting chaotic invariant set, such that the spatiotemporal dynamics spontaneously collapses to a coexisting attractor with typically regular dynamics.1 Turbulence in shear flow experiments is transient,7 and “stable chaos” (with a negative Lyapunov exponent) is transient in systems of coupled onedimensional maps.8 Chaotic transients have also been discussed as a mechanism for species extinction in ecology.9,10 In excitable, extended reaction-diffusion systems, transient spatiotemporal chaos was reported in models for semiconductor charge transport,11 for the CO oxidation on a single-crystal Pt surface,12 and for a cubic autocatalytic a)

Electronic mail: [email protected] Electronic mail: [email protected]

b)

1054-1500/2014/24(1)/013126/10/$30.00

reaction.13 In the complex Ginzburg-Landau equation, spatiotemporal chaos is asymptotic for periodic boundary conditions and collapses to an oscillatory medium in the presence of non-periodic (Neumann) boundary conditions.14 The Morris-Lecar (ML) neuron model is an empirically determined ion current model of the barnacle giant muscle fiber.15 It serves as a computationally efficient reduction of the Hodgkin-Huxley neuron model16 and displays class I and class II neural behavior. Various studies considering a single ML-neuron or coupled ML-neurons show dynamical phenomena that resemble (class I) natural neuron behavior. For example, stochastic inputs can provide the high variability in interspike intervals necessary for neural coding;17 few coupled ML-neurons can model various central pattern generators for asymmetric traveling waves in the pyloric network of the lobster stomatogastric ganglion18 and in quadruped locomotion.19 In networks of coupled MLneurons, the constructive role of noise for an optimal neuron response to subthreshold stimuli was reported.20 Spatiotemporal chaos and chaotic itinerancy were discussed by Fujii and Tsuda21 in a network of ML-neurons, diffusively coupled by gap junctions, motivated from physiological observations from the neocortex. This paper explores transient spatiotemporal chaos in a ring network of diffusively coupled (class I) ML neurons. Diffusive coupling of axons is relevant for axo-axonic gap junction networks22 and can be generalized to active cables.23 The model is introduced in Sec. II. The chaotic transients of the neuron activity and their collapse to either a traveling pulse solution or to the rest state are discussed in Sec. III. Immediate pathways from the neuron network in the rest state to the chaotic transient state under singular system perturbation are presented in Sec. IV. In Sec. V, the robustness of the chaotic transients is discussed for various synchronous, initial neuron inputs. II. MODEL

Diffusively coupled neuron models describe, e.g., axoaxonic gap junction networks22 or the electrophysiology of active cables.23 We consider a ring network of N diffusively 24, 013126-1

C 2014 AIP Publishing LLC V

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

013126-2

K. Keplinger and R. Wackerbauer

Chaos 24, 013126 (2014)

coupled, identical, excitatory Morris-Lecar15 neurons. The state of neuron i (i ¼ 1, …, N) is given by the membrane potential, Vi, and the fraction of open potassium channels, ni, following: I  IL  ICa  IK þ D½Viþ1 þ Vi1  2Vi ; V_ i ¼ C nss  ni ; n_i ¼ sn

(1)

where C is the membrane capacitance, I is the externally applied current, IL is the leak current, ICa is the calcium current, IK is the potassium current, and D is the electrical coupling constant between neurons. nss is the fraction of open potassium channels at steady state and sn is the time constant for opening/closing of potassium channels. The leak (IL), calcium (ICa), and potassium currents (IK) are given by IL ðVÞ ¼ gL ðV  VL Þ; ICa ðVÞ ¼ gCa mss ðV  VCa Þ; IK ðV; nÞ ¼ gK nðV  VK Þ; with leak conductance gL, membrane resting potential VL, maximum calcium conductance gCa and potassium conductance gK, and reversal potential for the calcium current VCa and for the potassium current VK. The fraction of open calcium (potassium) channels at steady state, mss (nss), and the time scale for potassium channel opening, sn, are given by    1 V  V1 ; 1 þ tanh mss ðVÞ ¼ V2 2    1 V  V3 nss ðVÞ ¼ ; 1 þ tanh V4 2 sn ðVÞ ¼

1  ; V  V3 ucosh 2V4

where u is a time constant and V1, V2, V3, V4 are fitting parameters informed by experimental data.15 We adopt widely used parameters for the excitable ML neuron,24 C ¼ 20lF=cm2 ; gK ¼ 8mS=cm2 ; gCa ¼ 4mS=cm2 ; gL ¼ 2mS=cm2 ; VK ¼ 80mV; VCa ¼ 120mV; VL ¼ 60mV, V1 ¼ 1:2mV; 1 V2 ¼ 18mV; V3 ¼ 14:95mV; V4 ¼ 17:4mV, and u ¼ 15 Hz. 1 D ¼ 0:05ðmsÞ , is representative of a range of coupling constants exhibiting spatiotemporally complex behavior. A bifurcation analysis for a single ML neuron (Fig. 1) with bifurcation parameter I reveals a SNIC (saddle-node on invariant circle) bifurcation near I ¼ 38:7 lA=cm2 . The ML neuron exhibits a stable limit cycle and thus oscillatory behavior for applied currents above this bifurcation point and excitable behavior below this bifurcation point, where the neuron dynamics is characterized by three steady states, a stable node (rest state), a saddle point, and an unstable focus. A neuron is called excitable, if its response to a superthreshold perturbation yields an excitation cycle (Fig. 1), which means for the ML neuron that a typical trajectory moves further away from the rest state, passes around the unstable focus, and returns to the stable node. At the

FIG. 1. Bifurcation diagram for the membrane potential (V) as a function of applied current (I). A SNIC bifurcation occurs near I ¼ 38:7 lA=cm2 , and a subcritical Hopf bifurcation occurs near I ¼ 41:4 lA=cm2 . The corresponding fixed points are stable node (full line, filled circle), saddle point (triangle), unstable focus (star), and stable focus (line with squares); the stable limit cycle (line with open circles) and the unstable limit cycle (line with pluses) are represented by their extrema in V. The inset shows a typical excitation cycle (dashed line) in phase space for I ¼ 32 lA=cm2 , together with the cubic V-nullcline and the n-nullcline (full lines), and the three steady states.

subcritical Hopf bifurcation near I ¼ 41:4 lA=cm2 , the unstable focus becomes stable and an unstable limit cycle is born. An in-depth bifurcation analysis of a single ML neuron was performed by Tsumoto et al.24 The ring network of coupled excitable ML neurons is excitable if the excitation can spread between neurons. A numerical analysis shows that the threshold for excitation propagation in the ring network is near I ¼ 28:1 lA=cm2 . This study focuses on the excitatory parameter regime (class I excitability) between the onset of excitation propagation and the SNIC bifurcation, I 2 ½28:1; 38:7 lA=cm2 , where complex spatiotemporal dynamics evolves. The network does not exhibit a chaotic steady state in this parameter regime. III. CHAOTIC TRANSIENTS

Typical spatiotemporal patterns of the excitable ML ring network are shown in Fig. 2 for applied currents in the parameter regime between the onset of excitation propagation and the SNIC bifurcation. Each pattern is characterized with an irregular appearance of neuron groups in the neighborhood of the rest state (white patches) intermingled with active neuron groups, characterized by subtle, persistent ripple-like structures (dark patches), where trajectories spiral around the unstable focus. Pulse trains are frequent for low applied current. As more current is injected into the system, active neuron groups become more dominant, and the (mean) spatial scale of neuron groups close to the rest state decreases. For larger currents (e.g., Figs. 2(c)–2(e)), local oscillatory events (chain of white “beads”) develop, where a group of neighboring neurons undergoes subsequent coupling-induced excitation cycles. The number of subsequent oscillatory events is increasing and the duration of a single oscillatory event (white bead) is decreasing as the SNIC bifurcation is approached. We note that such

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

013126-3

K. Keplinger and R. Wackerbauer

FIG. 2. Complex spatiotemporal dynamics of the membrane potential (V) in the ring network of N ¼ 50 Morris-Lecar neurons for an applied current of (a) I ¼ 28:5 lA=cm2 , (b) I ¼ 30 lA=cm2 , (c) I ¼ 32 lA=cm2 , (d) I ¼ 35 lA=cm2 , and (e) I ¼ 38 lA=cm2 . Colors indicate when the membrane potential is close to the resting potential (white) and close to the unstable focus (dark). Initially, 10 randomly selected neurons out of the 50 neurons in the rest state are set to (V, n) ¼ (10 mV, 0). In general, the number of neurons required to likely induce complex spatiotemporal dynamics increases with distance to the SNIC bifurcation point, e.g., at I ¼ 38 lA=cm2 only a single excited neuron would suffice. The parameters used for simulating Eq. (1) are given in Sec. II.

oscillatory events (beads) also exist above the SNIC bifurcation point with a similar duration, but much faster than the period of the stable limit cycle. This finding resembles discussions in the neuroscience literature that axo-axonic gap junctions can induce very fast network oscillations.22,25 Figure 3 shows typical trajectories of the spatiotemporal dynamics from Fig. 2(e) with I slightly below the SNIC bifurcation point (I ¼ 38 lA=cm2 ): fast local oscillations that correspond to a chain of beads (thick solid line), the spiraling around the unstable focus that corresponds to an active neuron in a dark ripple-like structure (thin solid line), and

FIG. 3. Typical trajectories for individual network neurons from Fig. 2(e): fast local oscillations corresponding to a chain of beads (thick solid line), spiraling around the unstable focus corresponding to an active neuron in a dark ripple-like structure (thin solid line), and the dynamics from a neuron in between these two (thin dashed line). The inset enlarges the area near the unstable focus. The cubic V-nullcline and the n-nullcline (thick dashed lines), together with the stable node (dot), the saddle point (triangle), and the unstable focus (star) for an uncoupled ML neuron are added for reference.

Chaos 24, 013126 (2014)

FIG. 4. Convergence behavior of the largest Lyapunov exponent (k) with time (t) for (a) fixed applied current, I ¼ 32 lA=cm2 , and varying network sizes [N ¼ 100 (triangle), N ¼ 500 (solid line), N ¼ 1000 (circle), and N ¼ 2000 (plus)], and (b) fixed network size, N ¼ 100, and varying applied currents [I ¼ 28:5 lA=cm2 (triangles), I ¼ 30 lA=cm2 (solid line), I ¼ 32 lA=cm2 (dashed line), I ¼ 35 lA=cm2 (pluses), and I ¼ 38 lA=cm2 (circles)]. Spatiotemporal chaos is short-lived for I ¼ 28:5 lA=cm2 (triangles, Fig. (b)). The maximum Lyapunov exponent is measured numerically due to Benettin et al.26 with a renormalization time interval of dt ¼ 0.1 ms and an initial trajectory distance of d0 ¼ 0.01 mV.

the dynamics from a neuron in between these two cases (thin dashed line). The complex spatiotemporal dynamics (Fig. 2) is chaotic with a positive maximum Lyapunov exponent. Figure 4 shows the convergence of the finite-time Lyapunov exponent (k) with simulation time, using the numerical method developed by Benettin et al.26 The Lyapunov exponent reaches a unique value for a fixed applied current I and network sizes with a large enough number of neurons (e.g., N ¼ 500, 1000, and 2000; Fig. 4(a)). This indicates that spatiotemporal chaos in the ML ring network is independent of system size above a characteristic length scale, and thus extensive.27 For smaller networks (e.g., N ¼ 100), k is slightly reduced, which reflects that spatiotemporal chaos is still influenced by the periodic boundary conditions and network size. Figure 4(b) shows that the maximum Lyapunov exponent typically decreases as the SNIC bifurcation is approached. After a regime of sustained spatiotemporal chaos, the dynamics collapses suddenly to either the rest state (Fig. 5(a)) or a pulse solution (Fig. 5(b)) throughout the parameter range of interest. Such a sudden collapse points to the existence of a chaotic saddle as a result of the coexisting regular attractors.1 The master stability analysis for transient spatiotemporal chaos in diffusively coupled ordinary differential equations (ODEs)28 can be readily applied to the ML-network in Eq. (1), since the rest state is an attractor on the invariant synchronization manifold (if all neuron states are equal at time t ¼ 0, they are equal for all t > 0). The analysis reveals that the rest state is also attractive in the coupled system, since the eigenvalues of the transverse modes to the synchronization manifold are negative everywhere on the attractor for a Laplacian coupling scheme; this opens up the possibility for the collapse of spatiotemporal chaos to the rest state to occur close to the synchronization manifold, but does not account

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

013126-4

K. Keplinger and R. Wackerbauer

Chaos 24, 013126 (2014)

FIG. 5. The collapse of spatiotemporal chaos to (a) the rest state and (b) a pulse solution. The membrane potential (V) is plotted for a 50-neuron network with I ¼ 30 lA=cm2 and two different initial conditions. Refer to Fig. 2 for technical details.

for the collapse to the pulse state. A long-term numerical analysis reveals that the pulse profile does not change over large time scales, several orders of magnitude larger than the average lifetime of spatiotemporal chaos, which points to the pulse solution being an attractor. The phase coherence of spatiotemporal chaos is measured with Kuramoto’s order parameter.29 The amplitude, R 僆 [0, 1], of the mean field neuron population is defined via Reih ¼

N 1X ei/j ; N j¼1

FIG. 7. Average transient lifetime (T) for collapse to rest state (plus) and for collapse to pulse state (square) versus (a) network size N for I ¼ 30 lA=cm2 , and versus (b) applied current I for N ¼ 30. The line represents a linear fit of the average transient lifetimes. Each data point was determined from 100 random initial conditions; for clarity, error bars (1 standard deviation) are plotted only for the collapse to the rest state, with similar error bars for the collapse to the pulse state. (c) Percent of finite lifetimes with a collapse to a pulse state, Pp, as a function of applied current I for size N ¼ 20, 30, 50, and 60 (ascending). The dotted line provides visual guidance to identify data points for fixed N. Spatiotemporal chaos collapsed to either the steady-state or the pulse solution for all considered parameters with the exception of I ¼ 38 lA=cm2 , where 8% (2%) of the simulations collapsed to an asymptotic, spatially regular state for N ¼ 20 (N ¼ 30). Missing data points at I ¼ 38 lA=cm2 and N ¼ 50 and N ¼ 60 are due to long transient lifetimes.

(2)

where /j is the phase of neuron j relative to the unstable focus and h is a mean field phase. A homogeneous phase distribution yields R ¼ 0, whereas R ¼ 1 for complete phase coherence. Figure 6 shows that R irregularly fluctuates between these extreme values during the regime of spatiotemporal chaos without a trend towards increased coherence as time progresses; R ¼ 1 is approached several times without a collapse following. Consistent with other cases of transient chaos,1,13 the collapse is spontaneous and not predictable in advance, which supports the existence of a chaotic saddle.

FIG. 6. (a) Regime of spatiotemporal chaos followed by a collapse to the rest state, and (b) corresponding order parameter R [Eq. (2)] for a ML ring network with N ¼ 40 neurons and an applied current of I ¼ 30 lA=cm2 . Refer to Fig. 2 for technical details.

The average lifetime (T) of transient chaotic behavior grows exponentially with network size (Fig. 7(a)) as is typical of transient chaos in reaction-diffusion systems.1,13,30 The lifetimes do not differ significantly from each other in mean or standard deviation whether the collapse is to the rest state or to the pulse solution. This is expected if one assumes a single non-attracting chaotic set, which implies the same escape probability from the neighborhood of the chaotic saddle to either of the two attractors. The exponential growth of the average lifetime with system size is arguably due to the probability of randomly uncorrelated regions generating a global pattern that initiates the collapse of spatiotemporal chaos1,31 and should hold also in the infinite limit. The exponent derived from simulations with small network sizes (Fig. 7(a)) is assumed to differ in the infinite limit, since the range of network sizes that are currently computationally accessible (even with supercomputing power) is limited, and typically below network sizes for which boundary effects are negligible and extensivity is reached.27 Figure 7(b) shows that the average transient lifetime also increases exponentially with the applied current for a fixed network size. In addition, we find that it is increasingly more likely that spatiotemporal chaos collapses to a pulse solution rather than to the rest state when the network size increases (Fig. 7(c)). The probability of collapse to a pulse state, Pp, is slightly enhanced in the neighborhood of I ¼ 28:5 lA=cm2 and I ¼ 35 lA=cm2 for all network sizes under consideration. This is caused by a high chance for fast collapse to the pulse state near I ¼ 28:5 lA=cm2 , superimposed with a monotonously increasing chance for slower collapse to the pulse state with increasing I. As the applied current approaches the SNIC

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

013126-5

K. Keplinger and R. Wackerbauer

FIG. 8. INSs measured immediately34 after the perturbation of the network from the rest state; the perturbation is initiated by a single input neuron in a 100-neuron network and I ¼ 35 lA=cm2 . The state of the input neuron is different for each simulation, in order to map out a rectangular area on the plane of n and V with resolution (DV, Dn) ¼ (2.0 mV, 0.01). The instantaneous state of the network is represented by the color of the pixel: rest state (white), pulse state (gray), and chaotic state (black). Both nullclines (full lines), the stable node (dot), the saddle point (triangle), and the unstable focus (star) for the uncoupled ML neuron are shown.

bifurcation point (I ¼ 38 lA=cm2 in Figs. 7(b) and 7(c)), transient spatiotemporal chaos can additionally collapse to a spatially regular pattern for small network sizes, which happened for 8% (2%) of the simulations in the case of N ¼ 20 (N ¼ 30). This spatiotemporal pattern is characterized by a persistent, spatially localized cluster of rhythmic neuron activity32,33 with a Lyapunov exponent that approaches zero and indicates periodic behavior; it is absent in the presence of noise. IV. PERTURBATIONS ON THE NETWORK AT REST

As demonstrated in Sec. III, the path from transient chaos to an attractor is deterministic and requires no external

Chaos 24, 013126 (2014)

perturbation. The main purpose of this section is to understand the pathway from the rest state attractor to transient chaos through an initial perturbation. The ring network is perturbed from the rest state by changing the state of a few local neurons, called input neurons, to an arbitrary state (V, n). For more than one input neurons, their states are identical. The corresponding instantaneous network state (INS) is recorded (1) as pulse state, if the extended system immediately approaches the pulse solution; (2) as rest state, if the trajectory immediately approaches the stable node; and (3) as chaotic state, if the trajectory immediately reaches the neighborhood34 of the chaotic saddle, displaying spatiotemporal chaos before collapsing to one of the attractors. The state of the input neuron(s) is different in each simulation, in order to map out a rectangular area on the plane of n and V. The instantaneous network states are recorded as a colored pixel (rest state: white, pulse state: gray, chaotic state: black) at the state of the input neuron(s) to generate a basin-like plot, which relates initial conditions to instantaneous network states and shows the consequences of local interactions of input neurons in a network at rest. The diagram of the INSs is closely related to a two-dimensional slice of the diagram of the basins of attraction, in which the basin boundary (stable manifold of the chaotic saddle) should be a set of vanishing measure.1 For the INSs, a finite neighborhood of the chaotic saddle is considered.34 Figure 8 shows the instantaneous network states for a single input neuron and an intermediate applied current. For most states of the input neuron, the dynamics on the network returned immediately to the rest state. For input neuron states between the nullclines of the uncoupled ML neuron, where dV dn dt > 0 and dt > 0, the dynamics on the ML network typically escapes into the neighborhood of the chaotic saddle or, for larger membrane potential, to the pulse state attractor. In general, such escapes are more likely with increasing applied current as the INS panels in Fig. 9 (d ¼ 0) reveal. For low current I, the region for escapes appears line-like and parallel

FIG. 9. Panel of instantaneous network states as a function of the applied current, I, for a single input neuron (d ¼ 0), and for two input neurons a distance d apart (1  d  6). In each simulation, the two input neurons are in the same state. The INS patterns for d > 6 are very close to the ones for d ¼ 6. Technical details for individual figures are explained in Fig. 8; the three dots (from left to right) represent stable node, saddle point, and unstable focus.

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

013126-6

K. Keplinger and R. Wackerbauer

to the n-nullcline, and the majority of escaping trajectories makes a transition directly to the pulse state. For medium applied currents, the region for escapes is clearly widening, and trajectories predominantly approach the neighborhood of the chaotic saddle. For higher currents, INSs close to the chaotic saddle or at the pulse state are very frequent in this area of phase space, with pulse states becoming more prevalent. This study also confirms the numerical finding in Sec. III that the excitation of a single neuron could not initiate spatiotemporal chaotic patterns for lower currents. Interference between two synchronous initial stimuli is investigated with two input neurons, both of them in the same state and spatially separated in the ring network by a distance of d neurons for 1  d  6. Figure 9 represents a series of INSs for varying input neuron distance (d) and varying applied current (I). The overall trend that increasing applied current is accompanied by more INSs in the neighborhood of the chaotic saddle also holds for two input neurons. For lower currents and d 6¼ 2, the system returns quickly to the rest state for the vast majority of input neuron states. d ¼ 2, however, is an effective distance between input neurons to reach a considerably larger region for escapes, which holds across all currents I; the corresponding boundaries between regions of different INSs can be complex, particularly for low I. For larger currents, INSs in the neighborhood of the chaotic saddle tend to dominate over pulse states for smaller distances d (e.g., d ¼ 1), whereas pulse states tend to dominate over chaotic states for medium distances (e.g., d ¼ 3). The boundaries between different INS regions throughout Fig. 9 typically show more structures on finer scales, appearing as striatal layers parallel to the boundary. For d ¼ 6, the INS pattern for each I is visually identical to the corresponding one for larger distances, d > 6, and mimics the pattern for the single input neuron (d ¼ 0). This indicates that the two input neurons behave as two isolated single input neurons with respect to their immediate pathway from the network at rest. Instantaneous network states resulting from the immediate local interaction of three input neurons are shown in Fig. 10 for varying interneuron distances. The location of the three input neurons is represented by their two nearest distances, ds and dl (dl  ds), between neighboring neurons in the ring network. In many aspects, results from the two-input study can inform the three-input study. At low current (I ¼ 28:5 lA=cm2 , Fig. 10(a)), ds ¼ 2 or dl ¼ 2 is again an effective distance between input neurons in order to reach a considerably larger number of INSs that are not in the rest state. The shapes of the INS patterns do not vary significantly for any given dl as long as ds ¼ 2, and they resemble the corresponding shape for d ¼ 2 in the 2-input study (Fig. 9). The interference of the three input neurons at small distances ds and dl give rise to new shapes of INS patterns for all currents (Figs. 10(a)–10(c)), e.g., (dl, ds) ¼ (2, 1), often exhibiting a higher portion of INSs not in the rest state compared to the study with the two input neurons, e.g., d ¼ 2. In most other cases, the INS patterns are governed by the interference of just two of the three input neurons that have a distance ds, typically ds > 2. Then, the shape of the region in which INSs are not in the rest state (fixed ds, dl  ds) is very

Chaos 24, 013126 (2014)

similar to the shape of that region in the study with the two input neurons (d ¼ ds), besides changes in the ratio of chaotic states to pulse states. This is more obvious in the cases of medium (I ¼ 32 lA=cm2 , Fig. 10(b)) and high currents (I ¼ 35 lA=cm2 , Fig. 10(c)), where such regions are clearly larger than for the small current (I ¼ 28:5 lA=cm2 ) in Fig. 10(a). For example, compare (dl, ds) ¼ (dl, 3) with d ¼ 3 in Fig. 9. These INS patterns (Figs. 9 and 10) from local interactions of few input neurons give insight into the immediate pathways of a randomly perturbed network at rest, typically chosen in the literature and also in Sec. III. Randomly chosen input neurons include various configurations of local input neurons, and the combined effect of the local interferences will result in a higher chance for the system to escape to the pulse state or to the neighborhood of the chaotic saddle. On time scales longer than “immediately after the perturbation,” the spatiotemporal system will mainly reach the state of transient chaos, since a local initiation of the chaotic state somewhere in the network eventually spreads over the entire network wiping out coexisting pulse states. Simulations with N/5 input neurons at randomly chosen locations confirm the prevalent escape to the neighborhood of the chaotic saddle for medium and larger currents, as long as the input neuron can initiate an excitation in the network. For low currents, the INS patterns are similar to the ones for d ¼ 2 in Fig. 9 or ds ¼ 2 in Fig. 10(a), since input neurons with other distances only yield a negligible number of escapes from the rest state. The frequent and fast collapse to the pulse state near I ¼ 28:5 lA=cm2 that was reported in the analysis of the lifetimes (Fig. 7(c)) is consistent with the INS pattern and the input neuron state (V, n) ¼ (10 mV, 0). V. PERTURBATIONS ON THE NETWORK IN THE NEIGHBORHOOD OF THE CHAOTIC SADDLE

The collapse of transient chaos in the Morris-Lecar ring network cannot be predicted, whether it is the time of collapse or the type of state after the collapse. This study is designed to explore the robustness of transient chaos under perturbation, which has consequences for facilitating the collapse in a statistical sense. The ring network is perturbed (at time t ¼ 0) from the initially chaotic state by changing the state of so-called input neurons to an arbitrary state in the plane of n and V, in analogy to Sec. IV. The corresponding INS is recorded at time t ¼ s ¼ 750 ms and can either be still the chaotic state or upon collapse, the rest state or the pulse state. In the absence of a perturbation, the dynamics would still be chaotic at time s. The INSs are recorded in phase space as a colored pixel (rest state: white, pulse state: gray, chaotic state: black) at the state of the input neuron to generate a basin-like plot that represents the escape from the neighborhood of the chaotic saddle under perturbation through synchronous input neurons. Typical patterns of INSs are shown in Fig. 11 for varying applied currents (I) and varying percentages (P) of input neurons in a ring network of size N ¼ 100. The locations of the input neurons are randomly chosen. In the absence of a

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

013126-7

K. Keplinger and R. Wackerbauer

Chaos 24, 013126 (2014)

FIG. 10. Panel of instantaneous network states for three input neurons with varying interneuron distances, ds and dl, for different applied currents, (a) I ¼ 28:5 lA=cm2 , (b) I ¼ 32 lA=cm2 , and (c) I ¼ 35 lA=cm2 . ds and dl represent the shortest and second shortest distance between neighboring neurons. For small distances in a large ring network, ds and dl are interchangeable. In each simulation, the three input neurons are in the same state. Technical details for individual figures are explained in Figs. 8 and 9.

perturbation (P ¼ 0%), the INS is chaotic. A small number of input neurons (P ¼ 1%, P ¼ 10% in Fig. 11) does typically not initiate a collapse, and the INS is chaotic for all or nearly all input neuron states for medium and higher currents. For low currents, the INS patterns become complex with

structure on smaller scales; the INS depends sensitively on the state of the input neuron, and, e.g., the type of collapse cannot be predicted. On the other end, for a high number of input neurons (P ¼ 90%, P ¼ 99% in Fig. 11), the system typically collapses to the rest state for all currents and most

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

013126-8

K. Keplinger and R. Wackerbauer

Chaos 24, 013126 (2014)

FIG. 11. Panel of INSs as a function of the applied current (I) and the percentage (P) of input neurons in a network of size N ¼ 100 in an initially chaotic state. In each simulation, all input neurons are in the same state, given by a point in the plane of n and V with resolution (DV, Dn) ¼ (2.0 mV, 0.01). The instantaneous state of the network is measured at time s ¼ 750 ms and represented by the color of the pixel: rest state (white), pulse state (gray), and chaotic state (black). The initial condition, i.e., the chaotic state and the randomly chosen location of the input neurons, is the same for all simulations within an individual INS pattern. The three dots (from left to right) represent stable node, saddle point, and unstable focus for the uncoupled ML neuron, and both nullclines are given by the full lines.

input neuron states, with the exception of a narrow region parallel to the V-nullcline [e.g., ðI; PÞ ¼ ð30 lA=cm2 ; 90%Þ], where chaotic INSs are prevalent. For medium to higher currents, an additional region around the unstable focus arises, where perturbations do not cause the collapse of spatiotemporal chaos [e.g., ðI; PÞ ¼ ð38 lA=cm2 ; 90%Þ]. The predominant collapse to the rest state is understandable, since a large number of synchronous input neurons can generate sufficiently low gradients in V that cease to excite a neuron. It is also consistent with the observation in a chemical system that continuing spatially homogeneous noise promotes the collapse of spatiotemporal chaos.35 For an intermediate number of input neurons (P ¼ 33%, 50%, and 66% in Fig. 11), the INS patterns become disjected with fractal structures arising. Towards larger currents, chaotic INSs are prevalent, and especially near the unstable focus, collapse becomes increasingly difficult. The details in the complex INS patterns for the lower half of currents can vary considerably for different realizations of chaotic initial state or different random locations of input neurons, but they all show continuous regions for collapse to pulse state or to rest state. The statistical analysis uses 100 INS patterns for each parameter set (I, P), resulting from 10 randomly chosen locations of the input neurons and 10 chaotic initial states that do not collapse before time s in the absence of a perturbation. The details in the 100 INS patterns can vary significantly, particularly for an intermediate number of input neurons, although dominant structures are typically preserved, e.g., the preference for no collapse in the case of input neuron states in the neighborhood of the unstable focus. Figure 12(a) shows the probability of a perturbation-induced collapse, PC, defined as the average fraction of the 100 INS patterns for which trajectories reached the rest state or the pulse state within s ¼ 750 ms of the perturbation. This analysis confirms the trends described in the panel of representative INS patterns (Fig. 11). A perturbation-induced collapse is less likely for increasing applied current but still

considerable for a sufficient number of input neurons. For example, at high current, I ¼ 35 lA=cm2 , spatiotemporal chaos collapses in 20% of the simulations for P ¼ 50% input neurons. For a high percentage of input neurons, the collapse is very likely (PC > 0.9) for all currents. For fixed current, the resistance to perturbation-induced collapse reveals a threshold in P, above which the collapse probability strongly increases, in particular for medium and higher currents. The probability (Pp) that a collapse is to a pulse state (Fig. 12(b)) decreases with increasing number of synchronized perturbations (P), independent of the applied current.

FIG. 12. Statistics of the instantaneous network states as a function of the applied current (I) given by (a) the probability PC for a collapse to either the rest state or the pulse state and (b) the probability Pp that a collapse is to a pulse state. The data points for each I are offset for different percentages of input neurons, with P ¼ 10%, 30%, 50%, and 90% (from left to right). The thin dashed lines connect data points with a fixed P. The data point for Pp is undefined for I ¼ 35 lA=cm2 and P ¼ 10%, since PC ¼ 0. The statistics was compiled from 100 simulations for each data point, using 10 different chaotic initial states and 10 different randomly selected locations of input neurons. See Fig. 11 for more technical details.

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

013126-9

K. Keplinger and R. Wackerbauer

VI. CONCLUSIONS

A ring network of diffusively coupled, identical, excitable ML neurons exhibits transient spatiotemporal chaos with a collapse to either the rest state or to a state of pulse propagation. The coexistence of these two escapes from the chaotic saddle differs from other studies on excitable extended systems, where collapse was reported to a single asymptotic state, e.g., to a steady state12,13 or to an oscillatory medium.11 In common with these systems is that the average lifetime of the chaotic transients increases exponentially with the network size. These systems differ, however, by their bifurcation scenarios towards the parameter regime of transient spatiotemporal chaos. The Morris-Lecar network is studied below the SNIC bifurcation point, where the stable limit cycle is absent. The Gray-Scott chemical reactiondiffusion system13 was studied below a subcritical Hopf bifurcation point, and the B€ar-Eiswirth chemical surface reaction-diffusion model was studied in a parameter regime below a saddle-loop bifurcation, where the stable limit cycle originating from a near supercritical Hopf bifurcation disappears.12 The Wacker-Sch€oll model for a layered superconductor has a stable limit cycle that is generated via a supercritical Hopf bifurcation and a diffusion driven Turing instability. These examples show that transient spatiotemporal chaos can exist in systems with various dynamical properties. The pathway from the rest state to transient spatiotemporal chaos is not easily predictable but reveals some trends. The immediate interaction of few, nearby, initially perturbed neurons in a network at rest can yield transient spatiotemporal chaos for a complex set of synchronous initial perturbations. Two neurons represent an effective distance between the perturbations to reach a considerably larger chance for chaos to arise in the network for all applied currents (bifurcation parameter). Another spatial configuration favoring the chaotic state is given by neighboring perturbed neurons in the case of medium to large applied currents. Larger sets of randomly located initial perturbations in a network at rest typically include various configurations of such nearby initial perturbations and result in a higher chance for transient chaos, often revealing the importance of perturbations that are two neurons apart. Manipulating or even controlling the lifetime of transient spatiotemporal chaos is important for various applications.1 We show that the collapse of transient spatiotemporal chaos can be facilitated in a statistical sense, if a singular perturbation is applied to a fraction of neurons when the network dynamics is in the chaotic state. For a high portion of synchronously perturbed neurons, transient chaos typically collapses to the rest state, whereas chaos typically persists until the time of measurement for a small portion of perturbed neurons. For an intermediate portion of perturbed neurons, early collapse is clearly induced, with a higher uncertainty about the type of asymptotic state, i.e., rest state or pulse state. Lifetimes of transient spatiotemporal chaos in excitable systems have been manipulated via persistent noise,35 the addition of nonlocal Laplacian coupling in the network topology,36 or competition between chaotic

Chaos 24, 013126 (2014)

populations.37 The consequences of such type of manipulations for the transient lifetimes in the Morris-Lecar neuron network need to be explored in future studies. Experimental and theoretical studies in recent years have suggested to consider transient states in addition to attractor states in pursuit of a more complete understanding of brain functionality.2 Transient behavior in the brain is thought to play an integral role in perceptual phenomena associated with stimuli and in cognitive function. Brain information processing is observed as a sequential switching of neuron activity among different neurons or neuron groups from one metastable state to another. Such transient and sequential neuron activity can result from a robust heteroclinic cycle, where trajectories can temporarily stay in the neighborhoods of sequential saddles. A recent study of odorevoked neuron activity in the locust antennal lobe reports that short transients play a dominant role in optimal stimulus separation,38 and that these transients might be consistent with the transient and sequential switching of neuron activity due to heteroclinic cycles. Further advances in experimental neuroscience will give insight into the relevance of saddles and related transient dynamics for brain function, and also for possible applications of transient spatiotemporal chaos or transient stable chaos in the neurosciences. ACKNOWLEDGMENTS

This research was based upon work supported by the National Science Foundation under Grant No. PHY-0653086 and by the Arctic Region Supercomputing Center at the University of Alaska Fairbanks. 1

Y. Lai and T. Tel, Transient Chaos: Complex Dynamics on Finite Time Scales, Applied Mathematical Sciences (Springer, New York, 2011). 2 M. I. Rabinovich and P. Varona, Front. Comp. Neurosci. 5, Art. No. 24, 1–10 (2011). 3 M. C. Cross and P. C. Hohenberg, Rev. Mod. Phys. 65, 851 (1993). 4 K. E. Daniels and E. Bodenschatz, Phys. Rev. Lett. 88, 034501 (2002). 5 Q. Ouyang and J.-M. Flesselles, Nature 379, 143 (1996). 6 Z. Qu, J. N. Weiss, and A. Garfinkel, Phys. Rev. Lett. 78, 1387 (1997). 7 B. Hof, J. Westerweel, T. M. Schneider, and B. Eckhardt, Nature 443, 59 (2006). 8 A. Politi, R. Livi, G.-L. Oppo, and R. Kapral, Europhys. Lett. 22, 571 (1993). 9 K. McCann and P. Yodzis, Am. Nat. 144, 873 (1994). 10 L. Shulenburger, Y.-C. Lai, T. Yalc inkaya, and R. D. Holt, Phys. Lett. A 260, 156 (1999). 11 A. Wacker, S. Bose, and E. Sch€ oll, Europhys. Lett. 31, 257 (1995). 12 M. C. Strain and H. S. Greenside, Phys. Rev. Lett. 80, 2306 (1998). 13 R. Wackerbauer and K. Showalter, Phys. Rev. Lett. 91, 174103 (2003). 14 S. Houghton, E. Knobloch, S. Tobias, and M. Proctor, Phys. Lett. A 374, 2030 (2010). 15 C. Morris and H. Lecar, Biophys. J. 35, 193 (1981). 16 A. L. Hodgkin and A. F. Huxley, J. Physiol. 117, 500 (1952). 17 B. S. Gutkin and G. B. Ermentrout, Neural Comput. 10, 1047 (1998). 18 S. Still and G. Le Masson, Neurocomputing 26, 533 (1999). 19 P.-L. Buono and M. Golubitsky, J. Math. Biol. 42, 291 (2001). 20 N. Montejo, M. Lorenzo, V. Perez-Villar, and V. Perez-Mu~ nuzuri, Phys. Rev. E 72, 011902 (2005). 21 H. Fujii and I. Tsuda, Neurocomputing 58–60, 151 (2004). 22 D. Feldmeyer and J. L€ ubke, New Aspects of Axonal Structure and Function (Springer e-Books, Springer, New York, 2010). 23 G. Ermentrout and D. Terman, Mathematical Foundations of Neuroscience, Interdisciplinary Applied Mathematics (Springer, New York, 2010).

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

013126-10 24

K. Keplinger and R. Wackerbauer

K. Tsumoto, H. Kitajima, T. Yoshinaga, K. Aihara, and H. Kawakami, Neurocomputing 69, 293 (2006). 25 A. Draguhn, R. D. Traub, D. Schmitz, and J. G. R. Jefferys, Nature 394, 189 (1998). 26 G. Benettin, L. Galgani, and J.-M. Strelcyn, Phys. Rev. A 14, 2338 (1976). 27 D. Stahlke and R. Wackerbauer, Phys. Rev. E 83, 046204 (2011). 28 R. Wackerbauer, Phys. Rev. E 76, 056207 (2007). 29 Y. Kuramoto and I. Nishikawa, J. Stat. Phys. 49, 569 (1987). 30 D. Stahlke and R. Wackerbauer, Phys. Rev. E 80, 056211 (2009). 31 J. P. Crutchfield and K. Kaneko, Phys. Rev. Lett. 60, 2715 (1988). 32 M. B€ar, M. Hildebrand, M. Eiswirth, M. Falke, H. Engel, and M. Neufeld, Chaos 4, 499 (1994). 33 Such regular patterns for small system sizes have been noticed in other studies,13,27,32 typically close to a bifurcation point, and influenced by boundary condition, network size and spatial scale of the dynamics. They are not considered further in this study.

Chaos 24, 013126 (2014) 34

Numerical analysis shows that the system either immediately goes extinct (no propagation), develops immediately into a pulse, or has a non-pulse activity that spreads throughout the network and always ends up to generate the chaotic pattern. The chaotic pattern however may not live long, so the trajectory may not get close to the chaotic saddle. We distinguish these three cases by evaluating the activity at a test neuron that is 25 neurons away from the input neuron. If the excitation does not reach the test neuron, the system went extinct. For a pulse solution, the integrated activity (over 1000 ms) at the test neuron is small in comparison to the chaotic state, since the pulse just travels through the test neuron and causes it to be active for a short period of time only. 35 R. Wackerbauer and S. Kobayashi, Phys. Rev. E 75, 066209 (2007). 36 S. Yonker and R. Wackerbauer, Phys. Rev. E 73, 026218 (2006). 37 R. Wackerbauer, H. Sun, and K. Showalter, Phys. Rev. Lett. 84, 5018 (2000). 38 O. Mazor and G. Laurent, Neuron 48, 661 (2005).

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 158.125.231.154 On: Mon, 24 Nov 2014 17:04:33

Transient spatiotemporal chaos in the Morris-Lecar neuronal ring network.

Transient behavior is thought to play an integral role in brain functionality. Numerical simulations of the firing activity of diffusively coupled, ex...
3MB Sizes 0 Downloads 3 Views