The Epileptor model: a systematic mathematical analysis linked to the dynamics of seizures, refractory status epilepticus and depolarization block.

Abstract One characteristic of epilepsy is the variety of mechanisms leading to the epileptic state, which are still largely unknown. Refractory status epilepticus (RSE) and depolarization block (DB) are other pathological brain activities linked to epilepsy, whose patterns are different and whose mechanisms remain poorly understood. In epileptogenic network modeling, the Epileptor is a generic phenomenological model that has been recently developed to describe the dynamics of seizures. Here, we performed a detailed qualitative analysis of the Epileptor model based on dynamical systems theory and bifurcation analysis, and investigate the dynamic evolution of “normal” activity toward seizures and to the pathological RSE and DB states. The mechanisms of the transition between states are called bifurcations. Our detailed analysis demonstrates that the generic model undergoes different bifurcation types at seizure offset, when varying some selected parameters. We show that the pathological and normal activities can coexist within the same model under some conditions, and demonstrate that there are many pathways leading to and away from these activities. We here archive systematically all behaviors and dynamic regimes of the Epileptor model to serve as a resource in the development of patient-specific brain network models, and more generally in epilepsy research.


Introduction
Epilepsy is a condition of the nervous system in which neuronal populations manifest as repeated epileptic seizures lasting a few minutes. These seizures are spontaneous and commonly accompanied by fast abnormal discharges (10 ms timescale), after which the brain activity slowly returns to normal. Epileptic seizures may present in different forms, as well as transitions from and to these pathological states. Some seizures are controlled by medication, particularly antiepileptic drugs (AEDs). However, there are seizures that last .1 h without returning to baseline, and do not respond to AEDs (31-43%), resulting then in so-called refractory status epilepticus (RSE) (Mayer et al., 2002;Holtkamp et al., 2005;Rossetti et al., 2005). Patients with RSE are at an increased risk of death. Indeed, the short-term mortality rates of RSE are estimated to be between 16% and 39% (Mayer et al., 2002;Holtkamp et al., 2005;Rossetti et al., 2005;Novy et al., 2010).
Other neuroelectric phenomena are linked to epilepsy, among which spreading depression (SD) is the most prominent. SD is characterized by a slowly propagating depolarization wave [or depolarization block (DB)] in neuronal networks, followed by a shutdown of brain activity (Pietrobon and Moskowitz, 2014). DB is a state in which the neuronal membrane is depolarized, but neurons stop firing. Spreading depression may occur during migraine and some seizures (Rogawski, 2008;Pietrobon and Moskowitz, 2014), and was first described by Leo (1944), who observed a depression of electroencephalographic (EEG) activity that moved across the cortex.
The mechanisms underlying the genesis of epilepsy, DB, and RSE are still largely unknown. In this article, we take an integrative approach toward the understanding of these neuroelectric phenomena. Our thinking is inspired by the wish to understand the underlying dynamic underpinnings rather than the biophysiological basis.
In this context, a large number of experimental and computational models have been proposed to clarify the basic mechanisms of seizures and DB. Most computational models rely on biophysically realistic parameters, trying to reproduce experimental data (Kager et al., 2000;Traub et al., 2001;Destexhe, 2008;Cressman et al., 2009;Ullah et al., 2009). Although these models provide important advances in seizures and DB research, they rarely produce general rules. Furthermore, there has been no attempt at evaluating whether or not seizures, RSE, DB, and normal brain activities can coexist, and, if so, under what conditions.
The previously mentioned models are rooted in physiological mechanisms generating a fairly limited range of behaviors. The physiological foundation is critical when the intended therapeutic intervention shall make use of the physiology. Examples include the identification of signaling pathways leading to the control of neurotransmitters linked to excitability such as lamotrigine or topiramate acting on calcium channels. Another type of intervention acts on the brain as a network and harnesses the capacity to modulate networks, such as stimulation, resection, or disconnection. In this case, physiological realism of a mechanism is to be replaced by dynamic realism as the network communication depends more on the type of signal rather than how it is generated. These types of models are phenomenological. A neural mass model of partial seizures called Epileptor was previously developed, which has found many applications in brain network modeling of epilepsy patients and is the network node model used in the European clinical trial EPINOV (www.epinov.com; Jirsa et al., 2014). Given the wide application of this model and its relevance for clinical research, we here archive in detail its dynamic repertoire. The Epileptor comprises one susbsystem (called subsystem 1) with two state variables responsible for generating fast discharges, another subsystem (called subsystem 2) with two state variables generating sharp-wave events (SWEs). The subsystem 1 (fast) and subsystem 2 (slow) are linked to a state variable, z, evolving on a very slow timescale called the permittivity variable (Jirsa et al., 2014). Interestingly, the transition from and to the pathologic states can occur autonomously, under the slow z evolution. Fast discharges and sharp-wave events are pathological features commonly associated with seizures despite their different forms. The goal of this article is, first, to perform a systematic mathematical analysis of the Epileptor; and, second, to determine the range of behaviors present in the Epileptor, with the further reaching goal to ask whether "normal" brain activities, seizures, RSE, and DB can coexist within the same model. To this aim, we present a qualitative analysis of the Epileptor model based on dynamical systems theory and bifurcation analysis.

Materials and Methods
We provide here a detailed bifurcation analysis of the Epileptor model, which is a neural mass model of partial seizures, to analyze seizure dynamics (Jirsa et al., 2014). The Epileptor model consists of a system of coupled nonlinear differential equations with five state variables. It comprises two 2D subsystems and one slow variable z. Subsystems 1 and 2 are responsible for fast discharges and SWEs, respectively (Jirsa et al., 2014). Analysis of the separated subsystems 1 and 2 was performed to provide more in-depth information on the Epileptor dynamics.

Numerical integration
All stochastic simulations are performed with linear additive Gaussian white noise with a zero mean and a variance of 0.0025 using the Euler-Maruyama method. Deterministic simulations are performed using the Runge-Kutta method with a maximal time step of 0.01.

Finding equilibrium points and their stability
In order to understand the Epileptor dynamics, we identify the equilibrium points and investigate their stability. This analysis is performed on the deterministic Epileptor model. We find the equilibrium points Eðx p 1 ; y p 1 ; z p ; x p 2 ; y p 2 Þ by solving _ We obtain a system of algebraic equations that we solve using the Matlab "solve"\ symbolic solver. The variables to solve for are the state variables of the Epileptor model, which are x 1 ; y 1 ; z; x 2 ; y 2 . The inputs of the solve function are the parameter values of the Epileptor model. Here, we analyze the existence of the equilibrium points when varying the m and x 0 parameters. The solutions are the vector equilibrium points. We first determine the equilibrium points of the whole system when the parameter x 0 varies, for two values of the parameter m (m = 0 and m = 0.5). We determine the stability of the equilibrium points E by evaluating the eigenvalues of the Jacobian matrix J at the equilibrium point E, which is stable if all the real parts of the eigenvalues of J are negative (Izhikevich, 2007).

Parameter space of equilibrium points
We explore the stability of the coexisting equilibrium points for a range of m and x 0 values. For each values combination, we determine the equilibrium points using the solve function and their stability using the Jacobian matrix. The complete results are drawn within a (m, x 0 ) parameter space for equilibrium points, which consists of various areas, each of which is characterized by coexisting equilibrium points with different stability. Moreover, when the parameter I ext 2 is varied, the Epileptor behavior changes. We therefore explore two parameter spaces of equilibrium points, for I ext 2 ¼ 0:45 (default Epileptor value) and I ext 2 ¼ 0.

Bifurcation diagram
One important step in the analysis of a mathematical model is the geometrical analysis of bifurcations, which here was performed on the deterministic Epileptor model. First, we draw a (z, x 1 ) bifurcation diagram of the Epileptor for constant m, here m = 0.5. The (z, x 1 ) bifurcation diagram comprises the equilibrium points of the Epileptor, where x 1 is the first coordinate of the vector of equilibrium points that we determined with respect to z. All equilibrium point solutions were obtained using the Matlab solve symbolic solver, where only the (_ x 1 ; _ y 1 ; _ x 2 ; _ y 2 ) equations were considered for constant z. We discretized the space for 801 points in the z-dimension and performed numerical continuation computing the linear stability using the Jacobian matrix. The (z, x 1 ) curve is divided into different branches, of which each is numerically continued separately representing a different stability type of the equilibrium points solutions. Bifurcations were identified at stability changes for given z-values and are indicated by dots in the corresponding visualizations. As the main parameters m and I ext 2 control the dynamics of the Epileptor, we typically present bifurcation diagrams for varying values of m and I ext 2 . All bifurcation diagrams were verified using XPPAUT and explicit numerical simulations of the model system (typically 2000-4000 units of time, of which the initial transients are removed (Figs. 1, 2), for which representative trajectories were plotted in the corresponding diagrams.

Finding a fast-slow limit cycle
We found a large attractor in the (y 1 , c , z) phase space below the SLE attractor shown in Figure 1 The z equation of the Epileptor model (Jirsa et al., 2014) was introduced as follows: _ z ¼ rðsðx 1 À x 0 Þ À zÞ; 8z: graphically the divergence of the Epileptor by introducing the averaging method, which is used to locate the periodic orbits in the phase space. The periodic orbits are intersections of the z-nullcline and the ,x 1 .-curve.
The z-nullcline is given by _ z ¼ 0, and the ,x 1 .-curve is the average value of x 1 for each z constant written in the following form: where x 1 ¼ f ðt; zÞ: We also determine the periodic orbits by using the Pontryagin's averaging technique (Shilnikov and Kolomiets, 2008). To illustrate this technique, we introduce a slow averaged nullcline written as follows: ,_ z. ¼ rðsð ,x 1 . À x 0 Þ À z À 0:1z 7 Þ if z , 0 rðsð ,x 1 . À x 0 Þ À zÞ if z ! 0 ¼ 0 where periodic orbits are the zeros of the ,_ z.. We determine the stability of periodic orbits using the following derivative: which represents the dynamics of the averaged equation. A periodic orbit is stable when Equation 13 is negative, and graphically if the graph of ,_ z. decreases at the given zero. First, we find the periodic orbits as the parameter x 0 varies in a (x 0 , z) bifurcation diagram of periodic orbits. We determine the bifurcation that results in the appearance or disappearance of the periodic orbits. Second, we explore two parameter spaces of periodic orbits for I ext 2 ¼ 0:45 and I ext 2 ¼ 0.

Stabilizing equilibrium points in the Epileptor model
The parameter x 0 can change the stability of equilibrium points and then can control the Epileptor behavior. Using the fact that a trajectory moves around the equilibrium point when it is stable, we use the bifurcation analysis to find the values of x 0 for which an equilibrium point is stable. There is a qualitative distinction between the stable equilibrium points, which we discuss in more detail later.

Coexisting attractors in the Epileptor model
The Epileptor behavior depends on the stability of its equilibrium points. More, there is a pre-existing limit cycle (LC) attractor, which is a stable periodic orbit. In this case, the trajectory can have a behavior that is controlled by the A, Time series of the Epileptor model (its enlarged view is shown on the right), the first (middle), and second (bottom) subsystem are plotted showing the principal components of a seizure-like event, that is an interictal period with no spikes, emergence of preictal spikes, ictal onset, seizure evolution, and emergence of sharp-wave events toward ictal offset. c , c 1 , and c 2 correspond to Àx 1 1 x 2 , x 1 , and x 2 respectively. B, The trajectory of the whole system is sketched in the (y 1 , c , z) phase space. Seizure offset and ictal onset emerge through the z evolution. Here all stochastic simulations were performed with Gaussian white noise using the Euler-Mayurama method. Main parameters values: m = 0, x 0 = −1.6, and r = 0.00035. Initial conditions are [0 −5 3 0 0 0.01].
Research Article: New Research stability of the equilibrium points or can jump to the limit cycle attractor. Thus, there is a coexistence of two attractors. We identify all coexisting attractors according to the values of the parameters m and x 0 .

Analysis of subsystem 1
We provide a detailed analysis of the subsystem 1 dynamics without coupling (i.e., subsystem 1 is not coupled to subsystem 2).

Subsystem 1 equilibrium points and stability
We analytically find the equilibrium points (x 1 , y 1 ) by solving the following equations: _ x 1 ¼ y 1 À f 1 ðx 1 ; x 2 Þ À z 1 I ext 1 ¼ 0 _ y 1 ¼ c 1 À d 1 x 2 1 À y 1 ¼ 0: We determine the stability of the equilibrium points by evaluating the eigenvalues of the Jacobian matrix J. An equilibrium point is stable if all the real parts of the eigenvalues of J are negative (Izhikevich, 2007).
We graphically find the equilibrium points by intersecting the x 1 -and y 1 -nullclines. We determine the nullclines of subsystem 1 and show how to find the equilibrium points and their stability in a phase plane.

Subsystem 1 bifurcation diagram
When finding the equilibrium points of the subsystem 1, we observe a qualitative change of the phase plane, which is interpreted as a bifurcation. We identify the different types of bifurcations that exist according to z, which is considered as a parameter control, and the parameter m.
Using the bifurcation diagrams, we analyze these bifurcations and the qualitative behavior of the susbsystem 1 with respect to z. The bifurcation analysis is performed on the deterministic subsystem 1. First, we plot the whole bifurcation diagram of the subsystem 1 for m = 0, which consists of two curves. The curves consist of the subsystem 1 equilibrium points for each value of z. When the parameter m changes, the shape of one of the curves changes on an interval of z. As a result, we then plot only this curve in bifurcation diagrams for m = 0 and m = 2, and describe the trajectories behavior for each value of m. Below, we summarize the types of bifurcations that exist with respect to m.

Fast-slow subsystem equations
The equilibrium points of subsystem 1 depend on z. When m = 0 and z = 3.1, three equilibrium points coexist: a saddle, a stable node and a stable focus (see m = 0; see Fig. 39A). When m = 1.5 and z = 3.1, then equilibrium points coexist: a saddle, a stable node, and an unstable focus (see m = 1.5; see Fig. 39A). A stable node and a stable focus are equilibrium points of a resting state and a nonoscillatory state, respectively. When the resting and nonoscillatory states coexist, the trajectories converge to one of them, depending on the initial conditions (see m = 0 and m = 1; see Fig. 39A).
To ensure that trajectories switch between two stable states, we introduce the following equation: r and s are positive constant parameters (r , , 1).z changes the input I ext 1 of the subsystem 1 to Àz 1 I ext 1 . Indeed, when z decreases, then Àz 1 I ext 1 is increased, and only a stable focus exists (see Fig. 39B). If m increases, the stable focus becomes unstable, surrounded by a stable limit cycle. Hence, only a stable limit cycle exists (see Fig. 39B). When z increases, then Àz 1 I ext 1 is reduced, and the stable limit cycle coexists with a stable node (see Fig. 39A). The saddle acts as a separatrix (S) between them. When z further increases, the stable limit cycle disappears through a homoclinic bifurcation, HB (see Fig. 39C), and hence only a stable node exists. Thus, z mimics a slow adaptation of the subsystem 1 to produce resting and oscillatory states, and ways to switch between them. The dynamics of subsystem 1 is fast and z is a slow variable, hence Equations 18-20 represent a fastslow subsystem.
Finding equilibrium points of the fast-slow subsystem We analytically find the equilibrium points by solving Equations 19 and 20, where z is a solution of _ z ¼ 0 (Eq. 18). We determine the stability of the equilibrium points by analyzing the Jacobian matrix J. The equilibrium point is stable if all the real parts of the eigenvalues of J are negative (Izhikevich, 2007). We graphically find the equilibrium points of the fast-slow subsystem by using the (z, x 1 ) bifurcation diagram and the z-nullcline, which is related to the parameter x 0 . We show how the z-nullcline moves in the bifurcation diagram when varying x 0 , and use the bifurcation diagram to find the equilibrium points of the fast-slow subsystem for each value of x 0 .

Finding periodic orbits of the fast-slow subsystem
As follows from the Epileptor analysis, we determine the periodic orbits of the fast-slow subsystem as well as their stability by using the averaging method and the Pontryagin's averaging technique. The slow averaged nullcline is given by the following: x 1 ¼ f ðt; zÞÞ is a solution of the subsystem 1 with z constant.
Here, we plot the graph of ,_ z. to explain how we find the periodic orbits and their stability. The periodic orbits are the zeros of the graph of ,_ z.. It is stable if the graph of ,_ z. decreases at the given zero. More, we find the periodic orbits as x 0 varies and determine the bifurcation leading to their appearance or disappearance.

Stabilizing equilibrium points in the fast-slow subsystem
Our approach to stabilizing the equilibrium points of the fast-slow subsystem is the same as that of the Epileptor model: analyze bifurcation diagrams of the fast-slow subsystem and then find the values of x 0 for which an equilibrium point is stable. We show how stable equilibrium points can correspond to different states.

Coexisting attractors in the fast-slow subsystem
The fast-slow subsystem exhibits a bistability of two attractors as the Epileptor model. The first attractor is the fast-slow limit cycle. The behavior of the second attractor depends on the parameters m and x 0 , which control the stability of the equilibrium point. We discuss the coexisting attractors when m and x 0 have different values.

Analysis of subsystem 2
In this section, we analyze the dynamics of subsystem 2 (without coupling), which generates SWEs. We explore the equilibrium points, evaluating their stability and determining the bifurcation types.

Subsystem 2 equations
The subsystem 2 equations are given by the following: where f 2 ðx 2 Þ ¼ 0 i f x 2 , À0:25 a 2 ðx 2 1 0:25Þ if x 2 ! À0:25 : Subsystem 2 equilibrium points and stability We analytically find the equilibrium points (x 2 , y 2 ) by solving the following equations: We determine the stability of the equilibrium points by evaluating the eigenvalues of the Jacobian matrix J 2 . J 2 is defined at the equilibrium point x 2 , which is stable if all the real parts of the eigenvalues of J 2 are negative (Izhikevich, 2007).
The equilibrium points lie graphically at the intersection of the x 2 -and y 2 -nullclines. We determine the nullclines of the subsystem 2 and show how the equilibrium points and their stability changes as the parameter I ext 2 varies.

Subsystem 2 bifurcation diagram
As the phase plane changes qualitatively when varying the parameter I ext 2 , we discuss this change in a bifurcation diagram where I ext 2 is the parameter control. The analysis is performed on the deterministic subsystem 2.

Results
To get a better understanding of the dynamics of the generation and termination of epileptic seizures, we adopted a computational model that reproduces epileptic activity and perform a detailed analysis using a mathematical approach.

Epileptor dynamics
Epileptor model behavior. The Epileptor equations generate SLEs, which are characterized by an onset and offset (Jirsa et al., 2014). We plot time series of the Epileptor system c , subsystem 1 c 1 , and subsystem 2 c 2 in Figure 1A. We find the major elements of an SLE: onset, timescale, offset of SLEs, and their recurrence. Figure 1A shows that during the ictal phase, fast discharges decrease in frequency with time. Ictal states are separated by a period of normal brain activity (non-ictal state). State variables x 1 and y 1 are responsible for generating fast discharges in the ictal state with a fast timescale ( Fig. 1A; see c 1 ). State variables x 2 and y 2 are responsible for generating the SWEs with an intermediate timescale ( Fig. 1A; see c 2 ). We plot the Epileptor trajectory in a (y 1 , c , z) phase space (Fig. 1B). Ictal and normal states (NSs) coexist and their coexistence necessitates a separation in the state space so-called "separatrix." The ictal onset occurs when the trajectory collides with the separatrix after a transient normal state. The seizure offset occurs when the trajectory collides with the separatrix after a transient ictal state. The separatrix acts as a barrier between ictal and normal states. The slow state variable z is responsible for the alternation of both states under a slow timescale (r ,, 1). Equilibrium points. Using the Matlab solve symbolic solver, we localize the equilibrium points E as x 0 varies in a (x 0 , zp) diagram for m = 0.5 ( Fig. 2A) and m = 0 (Fig.  2B). z p is the third component of the equilibrium points E. For m = 0.5 and m = 0, a unique equilibrium point exists when x 0 = −1.6, while three equilibrium points coexist when x 0 = −0.9 ( Fig. 2A,B). The Jacobian matrix J of the Epileptor at the equilibrium point E is given by the following: The stability of the Epileptor equilibrium points depends on x 0 ( Fig. 2A,B). For x 0 = −1.6 and m = 0 or m = 0.5, the equilibrium point is a saddle. The stable manifold of the saddle equilibrium point corresponds to a separatrix between ictal and normal states (Fig. 1B). The trajectory behavior in Figure 1A corresponds to a recurrent alternation between ictal and normal states, which characterizes SLEs.
• When x 0 = −2.5 and m = 0 or m = 0.5, the equilibrium point is a stable node. The trajectory behavior corresponds to a normal activity (Fig. 2C).
• When x 0 = −0.9 and m = 0, the equilibrium points are: an unstable focus, a saddle, and a stable focus (Fig.  2B). The Epileptor remains in a nonoscillatory state (Fig. 2D).
• When x 0 = −0.9 and m = 0.5, three equilibrium points coexist. The equilibrium points are as follows: one saddle and two unstable foci ( Fig. 2A). The trajectory behavior corresponds to a periodic solution (Fig. 2E). Parameter space of equilibrium points. The stability of equilibrium points changes when we vary m and x 0 ( Fig.  2A,B). Trajectories can exhibit: (1) a normal activity when an equilibrium point is a stable node, (2) a nonoscillatory state when an equilibrium point is a stable focus, (3) SLEs when a unique saddle equilibrium point exists, and (4) a periodic solution when equilibrium points are unstable. To explore the various coexisting equilibrium points, we plot a (m, x 0 ) parameter space of equilibrium points in Figure  3A.
The parameter space comprises nine areas. The equilibrium point is an unstable focus in area 1, a stable node in area 5, and a saddle in area 4. The Epileptor model has three equilibrium points in areas 2, 3, 6, 7, 8, and 9. In area 2, a stable node, a saddle and an unstable focus coexist. In area 3, a stable focus, a saddle, and an unstable focus coexist. In area 6, three saddles coexist. In area 7, two saddles and one unstable focus coexist. In area 9, two unstable foci and one saddle coexist. In area 8, an unstable focus, an unstable node, and a saddle coexist.
We plot a (m, x 0 ) parameter space of equilibrium points for I ext 2 ¼ 0 in Figure 3B. The parameter space comprises seven areas. The equilibrium point is an unstable focus in area 1, a stable node in area 5, a stable focus in area 6, and a saddle in area 7. The Epileptor model has three equilibrium points in areas 2, 3, and 4. In area 2, a stable node, a saddle, and an unstable focus coexist. In area 3, a stable focus, a saddle, and an unstable focus coexist. In area 4, two unstable foci and one saddle coexist.
Main parameters of the Epileptor model. m and x 0 play roles of particular importance on the Epileptor dynamics. We showed above that the equilibrium points of the Epileptor change as x 0 varies, thereby allowing the obtaining of different behaviors of the system, including normal activity, nonoscillatory state (Fig. 2), and the alternation between the normal and ictal states (Fig. 1). m can be considered as a parameter that controls the Epileptor dynamic during the ictal period. The equilibrium points of the Epileptor model are graphically defined by the intersection of the nullclines. The x 1 -nullcline ( _ x 1 ¼ 0) is a straight line for x 1 !0. The sign of the slope of the straight line varies the direction of movement of the state variable x 1 in the phase portrait (8x 1 !0). Then, the stability of the equilibrium points (8x 1 !0) changes when varying the slope, which depends on m, and hence the dynamic of the Epileptor can pass from fast discharges to nonoscillatory state during the ictal period as m varies.
Below, we theoretically demonstrate how the system can switch between normal and epileptic activities by varying x 0 and using bifurcation diagrams. Moreover, we explore the significant role of m in controlling the frequency of discharges and in generating a DB.
Bifurcation diagram of the Epileptor for m = 0.5. The Epileptor can alternate between ictal and normal states (see c ; Fig. 1) when the equilibrium point is a saddle (see x 0 = −1.6; Fig. 2). The alternation is interpreted mathematically as a bifurcation, which is a qualitative change of the system behavior (Izhikevich, 2007). We identify two bifurcation types in an SLE. The first bifurcation corresponds to the transition from normal to ictal activity. The second one corresponds to the transition from ictal to normal activity. We identify bifurcation types and find the equilibrium points in a bifurcation diagram. Using the Matlab solve symbolic solver, we plot a (z, x 1 ) bifurcation diagram of the Epileptor model in Figure 4 for m = 0.5. x 1 is the first component of the equilibrium points E. z acts as a control parameter (r = 0). Figure 4 first shows a Z-shaped curve (Z-curve). Z-lower and Z-upper branches consist of stable nodes and unstable foci, respectively. The Z-middle branch separates Z-lower and Z-upper branches and consists of saddles. The Z-middle and Z-upper branches collide as z increases in a saddle-node bifurcation SN 3 (Fig. 4, inset). The Z-middle and Z-lower branches collide as z decreases in a saddle-node bifurcation SN 1 .
Moreover, there are two branches above the Z-shaped curve (Fig. 4, inset). The lower branch consists of saddles and terminates as z decreases in a saddle-node bifurcation SN 4 . The upper branch comprises two sub-branches: one sub-branch (dashed) consists of unstable foci and another (dash-dotted) consists of saddles. The two subbranches collide as z increases in a saddle-node bifurcation SN 2 . The equilibrium points E of the Epileptor model lie at the intersection of the z-nullcline and the curve of equilibrium points (Fig. 4). The z-nullcline (_ z ¼ 0) is given by the following: and depends on x 0 . The z-nullcline moves downward (x 0 decreases) or upward (x 0 increases) in the bifurcation diagram; hence, it intersects the curve of equilibrium points at different sites. When x 0 = −1.6, then the z-nullcline intersects the Z-middle branch, which consists of saddles (Fig. 4). Consistent with Figure 2A, the equilibrium point is a saddle. The trajectory behavior corresponds to SLEs, which is illustrated in Figure 4 (right). When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation SNs 1 . Then the trajectory switches to the Z-upper branch, which consists of unstable foci surrounded by stable periodic orbits. Hence the trajectory exhibits an oscillatory solution on the Z-upper branch, which approaches the Z-middle branch as z increases. The oscillatory solution is homoclinic to one of the saddles along the Z-middle branch, then it is destroyed as z is increased. The trajectory terminates in a homoclinic bifurcation HB and switches to the Z-lower branch. Z-lower and Z-upper branches correspond to normal and ictal states, respectively. The Z-middle branch acts as a separatrix between ictal and normal states. The transition from normal to ictal states (first bifurcation) occurs through a saddle-node bifurcation, SN 1 . The transition from ictal to normal states (second bifurcation) occurs though a homoclinic bifurcation, HB. The homoclinic bifurcation z(HB) corresponds to an offset threshold and the saddle-node bifurcation zðSN 1 Þ corresponds to an onset threshold. Therefore, the transitions between ictal and normal states occur through a fold/homoclinic bifurcation. The system is bistable on [SN 1 , HB]. Bifurcation diagram when m decreases. To explore the role of m in the dynamics of the Epileptor, we plot a (z, x 1 ) bifurcation diagram, when m decreases. Let m = 0 (Fig. 5), the plot first shows a Z-shaped curve. Z-Lower, Z-middle, and Z-upper branches consist of stable nodes, saddles, and unstable foci, respectively.
The Z-upper branch is surrounded by stable periodic orbits (see maxðx 1 Þ and minðx 1 Þ curves), which terminate as z increases in a saddle-node on invariant circle (SNIC) bifurcation with zðSNICÞ . zðSN 1 Þ. Above the Z-curve, the lower branch consists of saddles and the upper branch comprises two sub-branches (Fig. 5, inset). The first subbranch (solid) consists of stable foci and terminates as z decreases in a SNIC bifurcation. The second sub-branch (dash-dotted) consists of saddles. The two sub-branches collide as z increases in a saddle-node bifurcation SN 2 . We plot a trajectory in Figure 5 (right), which corresponds to SLEs. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddlenode bifurcation SN 1 and the trajectory switches to the Zupper branch exhibiting an oscillatory solution, which terminates as z increases in a SNIC bifurcation (see inset). Then the trajectory is at the upper branch above the Z-curve, which consists of stable foci. The trajectory exhibits a nonoscillatory solution, which terminates as z increases in a saddle-node bifurcation SN 2 and the trajectory switches to the Z-lower branch. The Z-lower branch corresponds to the normal state. The Z-upper branch and Decreasing z, the Z-lower and Z-middle branches collide in an SN 1 bifurcation. Above the Z-curve, lower (dash-dotted) branch consists of saddles, and upper branch is divided into two subbranches: one sub-branch (dashed) consists of unstable foci and another (dash-dotted) of saddles. Increasing z, the two sub-branches collide in an SN 2 bifurcation. The inset is their enlarged view. Decreasing z, upper (dashed) and lower branches above the Z-curve collide in an SN 4 bifurcation. Increasing z, the Z-upper branch and lower branch above collide in an SN 3 bifurcation. The ,x 1 .-curve is the average value of x 1 for each z constant. Let x 0 = −1.6, the z-nullcline (_ z ¼ 0) is at the Z-middle branch. A SLE occurs with a fold/homoclinic bifurcation. A saddle (S) periodic orbit separates the SLE attractor (right) and a stable periodic orbit LC (to the left, final orbit not shown). Deterministic trajectories are plotted on both sides of the separatrix S defining two basins of attraction (indicated by arrows). r = 0.003 for LC and r = 0.0007 for SLEs. the upper branch above the Z-curve correspond to the ictal state (see inset). The transition from normal to ictal states (first bifurcation) occurs through a saddle-node bifurcation, SN 1 . The transition from ictal to normal states (second bifurcation) occurs through a saddle-node bifurcation, SN 2 . Fast discharges characterizing the ictal state correspond to the oscillatory solution, which disappears through a saddle-node on invariant circle bifurcation SNIC. Therefore, the transitions between ictal and normal states are said to be of a fold/circle bifurcation. The system is bistable on [SN 1 , SN 2 ]. We conclude that the bifurcation diagrams of the Epileptor model is a Z-shaped curve. Above, there are the following three branches: two branches consisting of saddles, and one branch consisting of stable foci for Bifurcation diagram for m = 0 and I ext 2 ¼ 0. We now analyze the bifurcation diagram when I ext 2 decreases. Let I ext 2 ¼ 0, we plot a (z, x 1 ) bifurcation diagram for m = 0 in Figure 6, which shows a Z-shaped curve. Z-lower, Z-middle, and Z-upper branches consist of stable nodes, saddles, and unstable foci, respectively. Decreasing z, Z-lower and Z-middle branches collide in a saddle-node bifurcation, SN 1 . Increasing z, Z-upper and Z-middle branches collide in another saddle-node bifurcation, SN 2 . The Z-upper branch terminates as z decreases in a saddle-node bifurcation SN 4 . Moreover, two branches appear as z decreases. The lower branch (dashed) consists of unstable foci, and the upper branch (dash-dotted) consists of saddles. Increasing z, lower and upper branches collide in a saddlenode bifurcation SN 3 . Decreasing z, the upper branch terminates in a saddle-node bifurcation SN 4 . z Figure 5. The Epileptor model bifurcation diagram with respect to the slow variable z (m = 0, I ext2 ¼ 0:45). The Z-lower (solid), Zmiddle (dash-dotted), and Z-upper (dashed) branches consist of stable nodes, saddles, and unstable foci, respectively. Decreasing z, the Z-lower and Z-middle branches collide in an SN 1 bifurcation. Above the Z-curve, lower (dash-dotted) branch consists of saddles, and the upper branch is divided into two sub-branches: one (solid) consists of stable foci and another (dash-dotted) of saddles. Increasing z, the two sub-branches collide in an SN 2 bifurcation. The inset is their enlarged view. Decreasing z, upper (dashed) and lower branches above the Z-curve collide in a SNIC bifurcation. Increasing z, the Z-upper branch and lower branch above collide in an SN 3 bifurcation. The ,x 1 .-curve is the average value of x 1 for each z constant. Let x 0 = −1.6, the z-nullcline (_ z ¼ 0) is at the Z-middle branch. A SLE occurs with a fold/circle bifurcation. A saddle (S) periodic orbit separates the SLE attractor (right) and a stable periodic orbit LC (to the left, final orbit not shown). Deterministic trajectories are plotted on both sides of the separatrix S defining two basins of attraction (indicated by arrows). Let x 0 = −1.6, the z-nullcline is at the Z-middle branch, then the equilibrium point is a saddle. We plot a trajectory in Figure 6, which shows transitions between Z-lower and Z-upper branches. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation, SN 1 , and the trajectory switches to the Z-upper branch, which consists of unstable foci. An unstable focus is surrounded by a stable periodic orbit (see maxðx 1 Þ and minðx 1 Þ curves). Hence, when the trajectory is at the Z-upper branch, it exhibits an oscillatory solution, which terminates as z increases in a homoclinic bifurcation HB. Then the trajectory switches to the Z-lower branch. Z-Lower (solid) and Z-upper (dashed) branches correspond to normal and ictal states, respectively. The transition from normal to ictal states occurs through a saddle-node bifurcation SN 1 . The transition from ictal to normal states occurs through an HB. Thus, the transitions between ictal and normal states occur through a fold/homoclinic bifurcation. The system is bistable on [SN 1 , HB].
Bifurcation diagram when m decreases, for I ext 2 ¼ 0. Let m = −0.5, we plot a (z, x 1 ) bifurcation diagram in Figure 7. The Z-upper branch is divided into two subbranches separated by a Hopf (H) bifurcation (Fig. 7). The first sub-branch (solid) consists of stable foci and terminates as z increases in a saddle-node bifurcation SN 2 . The second sub-branch (dashed) consists of unstable foci and terminates as z decreases in a saddle-node bifurcation SN 4 . An unstable focus is surrounded by a stable periodic orbit (see maxðx 1 Þ and minðx 1 Þ curves), which terminates as z increases in an H bifurcation. When x 0 = −1.6, then the z-nullcline is at the Z-middle branch, hence the equilibrium point is a saddle. We plot a trajectory in Figure 7, which shows transitions between Z-lower and Zupper branches. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation, SN 1 , and the trajectory switches to the Z-upper sub-branch (dashed) exhibiting an oscillatory solution, which terminates as z increases in an H bifurcation. Then, the trajectory is at the Z-upper (solid) sub-branch exhibiting a nonoscillatory solution, which terminates as z increases in a saddle-node bifurcation SN 2 . Z-lower (solid) and Z-upper branches correspond to normal and ictal states, respectively. The oscillatory solution corresponds to fast discharges. The transition from normal to ictal states occurs through a saddle-node bifurcation, SN 1 . The transition from ictal to normal states occurs through a saddle-node bifurcation SN 2 . The oscillatory solution (fast discharges) terminates in an H bifurcation. The transitions between ictal and normal states are said to be of a fold/Hopf bifurcation. The system is bistable on [SN 1 , SN 2 ]. Let m = −1, we plot a (z, x 1 ) bifurcation diagram in Figure 8, which has a Z-shaped curve. The Hopf bifurcation point H disappears, then the Z-upper branch consists of stable foci only, which is the equilibrium point of nonoscillatory state. Increasing z, the Z-upper (solid) branch terminates in a SNIC bifurcation (Fig. 8). Below the Zupper branch, two branches appear as z decreases: the upper (dash-dotted) branch consists of saddles, and the lower (dashed) one consists of unstable foci. An unstable focus is surrounded by a stable periodic orbit (see maxðx 1 Þ and minðx 1 Þ curves), which terminates as z increases in a SNIC bifurcation. Let x 0 = −1.6, the z-nullcline is at the Z-middle branch, then the equilibrium point is a saddle. We plot a trajectory in Figure 8, which shows transitions between Z-lower and Z-upper branches. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation, SN 1 , and the trajectory switches to the Z-upper (solid) branch, which consists of stable foci, thereby exhibiting a nonoscillatory solution. Then, the stable focus disappears as z increases through a saddle-node bifurcation, SN 2 , and the trajectory switches to the Z-lower branch. The transition from Z-lower to Z-upper branches occurs through a saddle-node bifurcation, SN 1 , and the transition from Z-upper to Z-lower branches occurs through a saddle-node bifurcation, SN 2 . The transitions between Zlower and Z-upper branches occur through a fold/fold bifurcation. The system is bistable on [SN 1 , SN 2 ]. Here, the SLE attractor reduces to a periodic switch between a nonoscillatory state and a NS.
We conclude that when I ext 2 ¼ 0, the transitions between Z-lower and Z-upper branches occur through a fold/homoclinic bifurcation, 8m ! 0; a fold/Hopf bifurcation for m = −0.5; and a fold/fold bifurcation, 8m , À 0:5. The system is bistable on [SN 1 , SN 2 ] for m , 0 and on [SN 1 , HB] for m ≥ 0. Moreover, the transitions between ictal and normal states occur through a fold/circle bifurcation when I ext 2 ¼ 0:45, and do not when I ext 2 ¼ 0.

Finding LC: a stable limit cycle
LC in the Epileptor model. We found a large geometrical object in the (y 1 , c , z) phase space below the SLE attractor shown in Figure 1. When using the z equation (Eq. 34), the Epileptor diverges with time ( Fig. 9C). The initial conditions i 2 are on the left of the separatrix (S), which means that i 2 are below the separatrix (S) and the SLE attractor in the phase space. To explain this divergence, we plot the z-nullcline and the ,x 1 .-curve in a (z, x 1 ) bifurcation diagram for m = 0.5 (Fig. 9A). The z-nullcline corresponds to When x 0 = −1.6, then the z-nullcline intersects the ,x 1 .-curve at one periodic orbit labeled as S (Fig. 9A). We consider two different initial conditions i 1 and i 2 such as zði 2 Þ , zðSÞ , zði 1 Þ. Only the trajectory i 1 is plotted in the bifurcation diagram, which exhibits transitions between Z-lower and Z-upper branches through a fold/homoclinic bifurcation. We plot a time series for the initial condition i 2 in Figure 9C, which shows that the trajectory diverges with time.
To stabilize the final state that the Epileptor evolves toward, we modified the z-equation by introducing the following equation: We plot the corresponding z-nullcline (_ z ¼ 0; Eq. 35) in a (z; x 1 ) bifurcation diagram for m = 0.5 (Fig. 9B). When x 0 = −1.6, then the z-nullcline intersects the ,x 1 .-curve at two periodic orbits: S and LC (Fig. 9B). We consider the initial conditions i 1 and i 2 such as zði 2 Þ , zðSÞ , zði 1 Þ. The trajectory i 1 exhibits transitions between Z-lower and Zupper branches, which occur through a fold/homoclinic bifurcation. We plot time series for the initial condition i 2 in Figure 9D, which shows a stable LC. We plot the trajectory as its transients toward LC in Figure 4 and the corresponding time series in Figure 9D. LC is characterized by a large amplitude and a fast-slow invariant cycle (Fig. 9D).  We consider Equation 35 and plot the z-nullcline when x 0 ¼ À0:9 in a (z; x 1 ) bifurcation diagram (Fig. 10). The znullcline intersects three branches of the (z; x 1 ) curve, as follows: lower (dash-dotted) and upper (dashed) branches above the Z-curve, and the Z-upper (dashed) branch. Then three equilibrium points coexist: two unstable foci and one saddle. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits, which are LC and S. Then the trajectory only exhibits an oscillatory solution, which stabilizes at the point LC.
Evolution of periodic orbits. We plot the maxðx 1 Þand minðx 1 Þ-curves to limit the oscillatory solutions (Figs. 5-8). Decreasing z, trajectories converge to LC and the amplitude of periodic orbits increases, leading to LC with a large amplitude. The intersection of the z-nullcline and the ,x 1 .-curve gives rise to S and LC, which are saddle and stable periodic orbits, respectively. Let x 0 = −1.6, two periodic orbits exist for m = 0.5 ( Fig.  4) and m = 0 (Fig. 5). We show only one periodic orbit in Figures 4 and 5, which corresponds to S. The second periodic orbit is shown in Figure 9B, which corresponds to LC. LC and S are stable and saddle periodic orbits, respectively. Trajectories plotted in Figures 4 and 5 exhibit SLEs (right) or converge to LC (left), depending on the initial conditions. SLEs and LC are considered as two attractors. The saddle periodic orbit (S) acts as a separatrix between the basin of attraction of SLEs and the basin of attraction of LC.
The z-nullcline moves downward (x 0 decreases) or upward (x 0 increases) in the bifurcation diagram, hence the existence of periodic orbits depends on x 0 . We find S and LC in a (x 0 , z p ) bifurcation diagram of periodic orbits (Fig.  11A). x 0 acts as a control parameter. Finding S helps us to limit the basin of attraction of SLEs and the basin of attraction of LC. Finding LC helps us to determine the zvalue at which LC stabilizes. Red (1) markers and black dots correspond to LC and S, respectively. Only LC exists for large x 0 . Decreasing x 0 , LC and S coexist, and collide as x 0 is further decreased in a saddle-node of periodic orbits bifurcation (SNPO), then fades. Figure 11B (m = −16) shows that for some x 0 , a third periodic orbit appears, which is labeled as blue top (1) markers. This periodic orbit is stable, and is localized above LC and S (Fig. 11B). Since the amplitude of periodic orbits decreases as z increases, the third periodic orbit is a stable limit cycle with a small amplitude, which we denoted as small limit cycle (SLC). SLC and LC are two stable limit cycles that coexist for some x 0 , and are separated by S. SLC exists after a SNPO bifurcation occurs when m = −16 and disappears when x 0 is further decreased.
Parameter space of periodic orbits. We explore two (m, x 0 ) parameter spaces of periodic orbits: for I ext 2 ¼ 0:45 in Figure 12A and for I ext 2 ¼ 0 in Figure 12B. There are five areas separated by a bold line (SNPO bifurcation). LC exists above, but not below. LC exists in area I and coexists with S in area II. LC and S coexist with SLC in area III. Only SLC exists in area IV, which is localized below a SNPO bifurcation (bold line). Periodic orbits disappear in area V.

SLC dynamics
We used the averaging method to find periodic orbits, which lie at the intersection of the z-nullcline and the ,x 1 .-curve. When I ext 2 ¼ 0, area III shows the coexistence of LC, SLC, and S (Fig. 12B).
The stability of the equilibrium points depends on m and x 0 in area III. We analyze the coexistence of LC and SLC for each stability type of equilibrium points. We plot a (z, x 1 ) bifurcation diagram for m = 1 in Figure 13A. When x 0 ¼ À1:8, then the z-nullcline is at the middle branch and the equilibrium point is a saddle. The z-nullcline intersects the ,x 1 .-curve at three periodic orbits: LC (left), S (middle), and SLC (right). LC and SLC are two stable limit cycles with large and small amplitudes, respectively. S acts as a separatrix between LC and SLC. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation SN 1 and the trajectory switches to the Z-upper branch, which consists of unstable foci. Then the trajectory exhibits an oscillatory solution on the Z-upper branch and continues to SLC at which z stabilizes (Fig. 13A, inset). The transition to the Z-lower branch does not occur even if the equilibrium point is a saddle. In fact, z stabilizes at SLC before the homoclinic bifurcation HB. The trajectory behavior corresponds to a periodic solution.
We plot a (z, x 1 ) bifurcation diagram for m = −2 in Figure  13B. Let x 0 = 0, the equilibrium point is an unstable focus. The z-nullcline intersects the ,x 1 .-curve at three periodic orbits ( Fig. 13B): LC (left), S (middle), and SLC (right). S acts as a separatrix between LC and SLC. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation SN 1 , and the trajectory switches to the Z-upper (solid) branch, which consists of stable foci. Increasing z, the trajectory exhibits a nonoscillatory solution on the Z-upper branch (solid), which terminates in a SNIC bifurcation (Fig.  13B, inset). Then, the trajectory is at the dashed branch exhibiting an oscillatory solution and continues to SLC at which z stabilizes. The transition to the Z-lower branch does not occur and the trajectory behavior corresponds to a periodic solution.
We conclude that when I ext 2 ¼ 0, then LC and SLC coexist in area III and S acts as a separatrix between them. The equilibrium point of the Epileptor model is either a saddle or an unstable focus. Time series of the Epileptor system c , subsystem 1 c 1 , and subsystem 2 c 2 are plotted in Figure 13C when the equilibrium point is a saddle, and in Figure 13D when it is an unstable focus. The SLC patterns depend on the stability of the equilibrium point. The subsystem 2 generates a resting state when the equilibrium point is a saddle (see c 2 ; Fig. 13C) and a spiking state when it is an unstable focus (see c 2 ; Fig. 13D).

Transition from epileptiform fast discharges to a DB
The Epileptor exhibits epileptiform fast discharges in the ictal state, which corresponds to an oscillatory solution at the Z-upper branch of the (z, x 1 ) curve. Let I ext 2 ¼ 0, the Z-upper branch consists of unstable foci when m = 0 ( Fig. 6), and consists of stable foci when m = −1 (Fig. 8).
Recall that the SLE reduces to a periodic switch between a nonoscillatory state and NS when m = −1.
Further decreasing m, the imaginary part of the complex-conjugate eigenvalues goes to zero, and then the Zupper branch only consists of stable nodes. We plot a (z, Figure 14A, which shows a Z-shaped curve. Z-lower (solid), Z-middle (dash-dotted), and Z-upper (solid) branches consist of stable nodes, saddles, and stable nodes, respectively. The Zupper branch terminates as z decreases in a SNIC bifurcation ( Fig. 14A, inset). Below the Z-upper branch (Fig. 14A, inset), lower (dashed) and upper (dash-dotted) branches appear, which approach as m decreases to the Z-upper (solid) branch. The lower branch consists of unstable foci, and the upper branch consists of saddles. Increasing z, the lower and upper branches collide in a saddle-node bifurcation, SN 3 . Note that zðSNICÞ , zðSN 3 Þ , zðSN 1 Þ. Let x 0 = −1.6, the equilibrium point is a saddle. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation, SN 1 , and the trajectory switches to the Z-upper branch, which consists of stable nodes. The trajectory exhibits a silent activity on the Z-upper branch, which terminates as z increases in a saddle-node bifurcation SN 2 . Then, the trajectory switches to the Z-lower branch. The transition between Z-lower and Z-upper branches occur through a fold/fold bifurcation. The system is bistable on [SN 1 , SN 2 ]. Time series of the Epileptor model c , subsystem 1 c 1 , and subsystem 2 c 2 are plotted in Figure 14B. The fast discharges during the ictal period disappear with decreasing m underlying a slow wave. This scenario corresponds to a depolarization block (or excitation block; Izhikevich, 2007), which is indicated by segment number 1 (see c ; Fig. 14B). The NS is indicated by segment number 2. Therefore, the transitions between ictal states and NSs reduce to a periodic transition between DB and NS, which are both slow manifolds.
The DB occurs when I ext 2 ¼ 0 but not when I ext 2 ¼ 0:45. To see this, we plot a (z, x 1 ) bifurcation diagram for m = −8 (I ext 2 ¼ 0:45) in Figure 15A, which is with a Z-shaped curve. Z-lower (solid), Z-middle (dash-dotted), and Zupper (dashed) branches consist of stable nodes, saddles, and unstable foci, respectively. Above the Z-curve, lower branch consists of saddles and upper branch consists of stable nodes (Fig. 15A, inset II). Decreasing z, lower and upper branches collide in a SNIC bifurcation. Let x 0 = −1.6, the equilibrium point is a saddle. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation SN 1 and the trajectory switches to the Z-upper branch exhibiting an oscillatory solution (Fig. 15A, inset I), which terminates as z increases in a SNIC bifurcation (Fig. 15A, inset II). Then, the trajectory exhibits a nonoscillatory solution ( Fig. 15A, inset II), which terminates as z increases in a saddle-node bifurcation, SN 2 (Fig. 15A), and the Epileptor switches to the Z-lower branch. The transitions between Z-lower and Z-upper branches occur through a fold/circle bifurcation. The system is bistable on [SN 1 , SN 2 ]. Time series of the Epileptor model c , subsystem 1 c 1 , and subsystem 2 c 2 are plotted in Figure 15B. The Epileptor exhibits fast discharges in the ictal state (see c ), and hence does not enter into a DB.
To understand why a DB occurs for I ext 2 ¼ 0 and not for I ext 2 ¼ 0:45, we analyze the [SN 1 , SN 2 ] interval, on which the Epileptor model is bistable. We determine the branches of the (z, x 1 ) curve existing on this interval (8x 1 ! 0) for I ext 2 ¼ 0:45 and I ext 2 ¼ 0. Distinct branches coexist on [SN 1 , SN 2 ] depending on m: • When m = 0, then three branches coexist on [SN 1 , SN 2 ], consisting of stable nodes, saddles, and unstable foci for I ext 2 ¼ 0:45 (Fig. 5), and they consist of saddles and unstable foci for I ext 2 ¼ 0 ( Fig. 6).
Therefore, a DB occurs when m = −8 and I ext 2 ¼ 0 because one branch exists on [SN 1 , SN 2 ] (8x ! 0), which consists of stable nodes. The SLE attractor is characterized by transitions between Z-lower and Z-upper branches. The Z-lower branch consists of stable nodes, which are the equilibrium points of the NS. When m = −8 and I ext 2 ¼ 0, the Z-upper branch consists of stable nodes, which are the equilibrium points of the DB. Hence, the SLE attractor reduces as m further decreases to a periodic switch between DB and NS.

Stabilizing equilibrium points in the Epileptor model
The transitions between ictal and normal states of the SLE attractor occur when the equilibrium point is a saddle (x 0 = −1.6). When x 0 increases, the Epileptor stabilizes in the ictal (nonoscillatory) state, and when x 0 decreases, the Epileptor model stabilizes in the normal state.
Stabilizing the equilibrium point of the nonoscillatory state. During the ictal state, epileptiform fast discharges (oscillatory state) disappear through three bifurcation types depending on m and I ext 2 : a homoclinic bifurcation (Figs. 4, 6), a SNIC bifurcation (Fig. 5), or a Hopf bifurcation (Fig. 7). When the equilibrium point is a saddle (x 0 = −1.6), the Epileptor switches differently to the normal  Figure 14. Transition from epileptiform fast discharges to depolarization block. A, Bifurcation diagram (z, x 1 ) of the Epileptor model. The Z-lower (solid), Z-middle (dash-dotted), and Z-upper (solid) branches consist of stable nodes, saddles, and stable nodes, respectively. Decreasing z, the Z-lower and Z-middle branches collide in an SN 1 bifurcation. Increasing z, the Z-upper and Z-middle branches collide in an SN 2 bifurcation. Below the Z-upper branch (see inset), the top (dash-dotted) branch consists of saddles and the bottom (dashed) branch consists of unstable foci. Increasing z, bottom and top branches collide in an SN 3 bifurcation. Decreasing z, the Z-upper (solid) branch and upper (dash-dotted) branch below collide in a SNIC bifurcation. Let x 0 = −1.6, the znullcline (_ z ¼ 0) is at the Z-middle branch. The Z-upper (solid) and Z-lower (solid) branches correspond to DB and NS, respectively. The SLE attractor reduces then to a periodic transition between DB and NS. B, Deterministic time series of the Epileptor system c , subsystem 1 c 1 , and subsystem 2 c 2 . Parameters are: m = − 8, I ext2 = 0, r = 0.0005. state according to the bifurcation type: homoclinic, SNIC, or Hopf. The Epileptor switches to the normal state just after a homoclinic bifurcation, and remains in the nonoscillatory state after SNIC and Hopf bifurcations, thereby exhibiting a nonoscillatory solution. The latter disappears as z increases through a saddle-node bifurcation, SN 2 , and then the Epileptor switches to the normal state. When increasing x 0 , the stability of the equilibrium point changes, and the Epileptor remains in the nonoscillatory state. We illustrate this scenario in Figure 16A-C.
Let m = 0 and x 0 = −0.9 (Fig. 16A), the z-nullcline intersects three branches of the (z, x 1 ) curve: lower and upper (solid) branches above the Z-curve, and the Z-upper (dashed) branch. The z-nullcline intersects the upper branch above Z at C. Then three equilibrium points coexist: an unstable focus, a saddle, and a stable focus (C). When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation, SN 1 , and the trajectory switches to the Z-upper branch exhibiting an oscillatory solution (fast discharges), which terminates as z increases in a SNIC bifurcation. Then, the trajectory exhibits a nonoscillatory solution at the upper (solid) branch above the Z-curve, and continues to C at which z stabilizes. Therefore, the transition to the Z-lower branch does not occur and the Epileptor remains in the nonoscillatory state. We plot this scenario in a (Y, X, z) phase space (see Fig.   28A, top trajectory). Corresponding time series are plotted in Figure 28B. We add numbers to identify trajectory segments. We indicate the normal state by (1), fast discharges by (2), and the final state by (3), which corresponds to the equilibrium point of the nonoscillatory state. Here, the NS exists but not the equilibrium point. Hence, after a transient normal state, the Epileptor exhibits fast discharges which disappear in a SNIC bifurcation, and then remains in the nonoscillatory state.
Let m = −0.5 and x 0 = −0.8 (Fig. 16B), the z-nullcline intersects the Z-upper (solid) sub-branch at C which is a stable focus. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation SN 1 and the trajectory switches to the Z-upper (dashed) sub-branch exhibiting an oscillatory solution (fast discharges), which terminates as z increases in a Hopf bifurcation, H. Then, the trajectory exhibits a nonoscillatory solution at the Z-upper (solid) sub-branch, which consists of stable foci, and continues to C at which z stabilizes. The transition to the Z-lower branch does not occur, and the Epileptor remains in the nonoscillatory state. We plot this scenario in a (Y, X, z) phase space (see Fig. 29A, top trajectory). Corresponding time series are plotted in Figure 29B. We indicate the normal state by (1), fast discharges by (2), and the final state by (3), which corresponds to the equilibrium point of the nonoscillatory state. After a transient normal state, the Epileptor exhibits fast discharges which disappear in a Hopf bifurcation, and then remains in the nonoscillatory state.
Let m = −1, the SLE attractor reduces to transitions between nonoscillatory and normal states when x 0 = −1.6, which occur through a fold/fold bifurcation (Fig. 8). When increasing x 0 , the Epileptor model stabilizes its stable focus, which corresponds to the equilibrium point of the nonoscillatory state. We illustrate this case in Figure 16C. Let x 0 = −0.8, the z-nullcline intersects the Z-upper branch at C, which is a stable focus. When a trajectory is at the Zlower branch, the stable node disappears as z decreases through a saddle-node bifurcation, SN 1 , and the trajectory switches to the Z-upper branch exhibiting a nonoscillatory solution. Then, the trajectory continues to C at which z stabilizes. The transition to the Z-lower branch does not occur and the Epileptor remains in the nonoscillatory state. We plot this scenario in a (Y, X, z) phase space (see Fig. 30A, top trajectory). Corresponding time series are plotted in Figure 30B. We indicate the normal state by (1) and the final state by (3), which corresponds to the equilibrium point of the nonoscillatory state. The imaginary part of the eigenvalues of C is responsible for the oscillations (2) around the stable focus. After a transient normal state, the Epileptor spirals into the equilibrium point of the nonoscillatory state exhibiting damped oscillations, and does not re-enter the normal state.
Stabilizing the equilibrium point of the DB. Figure 14 shows that the SLE attractor reduces to a periodic switch between DB and NS when further decreasing m (I ext 2 ¼ 0). The equilibrium point is a saddle (x 0 = −1.6). When increasing x 0 , the Epileptor model stabilizes its stable node, which corresponds to the equilibrium point of DB. To see this, we plot a (z, x 1 ) bifurcation diagram when x 0 = −0.8 in  Figure 16D. The z-nullcline intersects the Z-upper (solid) branch at C, which is a stable node. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation SN 1 and the trajectory switches to the Z-upper branch and enters into a DB. The trajectory continues as z increases to C at which z stabilizes. Thus, the transition to the Z-lower branch does not occur and the Epileptor remains in DB. We plot this scenario in a (Y, X, z) phase space (see Fig.  27A, top trajectory). Corresponding time series are plotted in Figure 27B. We indicate the normal state by (2) and the final state by (4), which corresponds to the equilibrium point of DB. Then, after a transient NS, the Epileptor enters into DB and then remains in it.
Stabilizing the equilibrium point of the NS. When the equilibrium point is a saddle (x 0 = −1.6), the normal state of the SLE attractor disappears through a saddle-node bifurcation SN 1 . When decreasing x 0 , the Epileptor stabilizes its stable node, which is the equilibrium point of the normal state. To see this, we plot a (z, x 1 ) bifurcation diagram when x 0 = −2.5 in Figure 17. The z-nullcline intersects the Z-lower branch at C, which corresponds to the equilibrium point of the normal state. The ictal state (Z-upper branch) exists but not the equilibrium point, hence the Epileptor remains in the normal state, at C which is a stable node.
When I ext 2 ¼ 0:45, the transition from ictal to normal states occurs according to two scenarios, depending on m.
• m = 0.5 (Fig. 17A): When a trajectory is at the Z-upper branch (dashed), the oscillatory solution terminates as z increases in an HB bifurcation and the trajectory switches to the Z-lower branch.
• m = 0 ( Fig. 17B): When a trajectory is at the Z-upper branch (dashed), the oscillatory solution terminates as z increases in a SNIC bifurcation, then the trajectory exhibits a nonoscillatory solution on the upper branch above the Z-curve, which terminates as z increases in a saddle-node bifurcation SN 2 . Then the trajectory switches to the Z-lower branch.
When I ext 2 ¼ 0, the transition from ictal to normal states occurs according to four scenarios, depending on m.
• For m ≥ 0 (Fig. 17C): The trajectory exhibits an oscillatory solution on the Z-upper (dashed) branch, which terminates as z increases in a homoclinic bifurcation HB, and then it switches to the Z-lower branch.
• For m = −0.5 (Fig. 17D): The trajectory exhibits an oscillatory solution on the Z-upper sub-(dashed) branch, which terminates as z increases in a Hopf bifurcation H. Then the trajectory exhibits a nonoscillatory solution on the Z-upper (solid) sub-branch, which terminates as z increases in a saddle-node bifurcation SN 2 , and the trajectory switches to the Zlower branch.
• For m −1, when a trajectory starts at i with zðSNICÞ , zðiÞ (third scenario), then it exhibits a nonoscillatory solution on the Z-upper branch, which terminates as z increases in a saddle-node bifurcation, SN 2 . Then the trajectory switches to the Z-lower branch (Fig. 17E). When a trajectory starts at i with zðiÞ ,zðSNICÞ (fourth scenario), then it exhibits an oscillatory solution, which terminates as z increases in a SNIC bifurcation. Then, the trajectory is at the Z-upper branch exhibiting a nonoscillatory solution, which terminates as z increases in a saddle-node bifurcation SN 2 . Then, the trajectory switches to the Z-lower branch (Fig. 17F).
For all m, when the trajectory is at the Z-lower branch, which consists of stable nodes, it continues as z decreases to C, at which z stabilizes. Therefore, the transition to the Z-upper branch does not occur and the Epileptor remains in the NS.
Coexistence in the Epileptor model Figure 11, A and B, shows the evolution of periodic orbits as x 0 varies. The periodic orbits LC and S disappear as x 0 decreases through a SNPO bifurcation. The dynamics of SLC when it exists, depends on the equilibrium points stability. Figure 4 shows that the z-nullcline and the ,x 1 .-curve intersect at two points: LC and S. Depending on initial conditions, trajectories either converge to LC or exhibit SLEs. We identify seven types of coexisting attractors in the Epileptor model.
Coexistence of LC and seizures (SLEs). Depending on m and I ext 2 , transitions between ictal and normal states, which characterize the SLE attractor, occur through a  Figure 4. Branches of the (z, x 1 ) curve in B are the same as in Figure 5. Branches of the (z, x 1 ) curve in C are the same as in Figure 6. Branches of the (z, x 1 ) curve in D are the same as in Figure  7. Branches of the (z, x 1 ) curve in E and F are the same as in Figure 8. fold/homoclinic bifurcation, a fold/circle bifurcation, or a fold/Hopf bifurcation. We describe for each bifurcation type the coexistence of LC and the SLE attractor: A fold/homoclinic bifurcation: We plot two trajectories in a (z, x 1 ) bifurcation diagram for m = 0.5 (Fig. 4). The equilibrium point is a saddle (x 0 = −1.6). The behavior of the first trajectory (right) corresponds to SLEs, which occur through a fold/homoclinic bifurcation. The Z-middle branch acts as a separatrix between ictal and normal states. The second trajectory (left) converges to LC. The saddle periodic orbit (S) separates LC and the SLE attractor. We plot LC, SLEs, and S in a (Y, X, z) phase space (Fig. 18A). SLEs and LC coexist and the separatrix between them corresponds to the saddle periodic orbit (S). SLEs lock in LC when z declines to below baseline shift (Fig. 18A). Time series are plotted in Figure 18B for SLEs and in Figure 18C for LC. Figure 6 shows that SLEs occur through a fold/homoclinic bifurcation when I ext 2 ¼ 0. Parameter settings in Figure 6 correspond to area II in Figure 12B. Hence, SLEs and LC coexist and are separated by S. Only the SLE attractor is plotted in Figure 6. We plot LC, SLEs, and S in a (Y, X, z) phase space (Fig. 19A). SLEs and LC coexist and are separated by a separatrix, which corresponds to the saddle periodic orbit (S). Time series are plotted in Figure  19B for SLEs and in Figure 19C for LC.
A fold/circle bifurcation: We plot two trajectories in a (z, x 1 ) bifurcation diagram (Fig. 5). The first trajectory (right) exhibits SLEs, and the second one (left) converges to LC. The equilibrium point is a saddle (x 0 = −1.6). The transitions between ictal and normal states occur through a fold/circle bifurcation. Then SLEs and LC coexist and are separated by a saddle periodic orbit (S). We plot LC, SLEs, and S in a (Y, X, z) phase space (Fig. 20A). SLEs and LC coexist and are separated by a saddle periodic orbit (S), which is the separatrix. Time series are plotted in Figure 20B for SLEs and in Figure 20C for LC.
A fold/Hopf bifurcation: Figure 7 shows that SLEs occur through a fold/Hopf bifurcation when m ¼ À0:5 (I ext 2 ¼ 0). The equilibrium point is a saddle (x 0 = −1.6). Parameter settings in Figure 7 correspond to area II in Figure 12B. Then SLEs and LC coexist and are separated by S. Figure  7 only shows the SLE attractor. We plot LC, SLEs, and S in a (Y, X, z) phase space (Fig. 21A). SLEs and LC coexist and the separatrix between them corresponds to a saddle periodic orbit (S). Time series are plotted in Figure 21B for SLEs and in Figure 21C for LC.
Coexistence of LC and a periodic switch between nonoscillatory state and NS. Figure 8 shows that the SLE attractor reduces to a periodic switch between nonoscillatory state and NS, which occurs through a fold/fold bifurcation when m = −1 (I ext 2 ¼ 0). The equilibrium point is a saddle (x 0 = −1.6). Parameter settings in Figure 8 correspond to area II in Figure 12B. Then LC and the periodic switch between nonoscillatory state and NS coexist, and are separated by S. Only the periodic switch between nonoscillatory state and NS is plotted in Figure 8. We plot LC, this periodic switch, and S in a (Y, X, z) phase space (Fig. 22A). LC and the periodic switch between nonoscillatory state and NS coexist and the separatrix between them corresponds to a saddle periodic orbit (S). Time series are plotted in Figure 22B for this periodic switch and in Figure 22C for LC.
Coexistence of two periodic orbits: LC and SLC. SLC is a stable periodic orbit with a small amplitude. When I ext 2 ¼ 0, then SLC, S, and LC coexist in area III (Fig.  12B). Figure 13 shows that the z-nullcline intersects the ,x 1 .-curve at three points: LC, S, and SLC. Here, we plot the SLC behavior only, which depends on the stability of the equilibrium point. SLC and LC coexist in both cases: when the equilibrium point is a saddle (Fig. 13A) and when it is an unstable focus (Fig. 13B). We plot SLC, LC, and S in a (Y, X, Z) phase space (Figs. 23,24). LC and SLC coexist and are separated by a saddle periodic orbit (S), which corresponds to the separatrix. When the equilibrium point is a saddle, the coexistence of LC and SLC is illustrated in Figure 23A. Time series are plotted in Figure  23B for SLC and in Figure 23C for LC. When the equilibrium point is an unstable focus, the coexistence of LC and SLC is illustrated in Figure 24A. Time series are plotted in Figure 24B for SLC and in Figure 24C for LC.
Coexistence of LC and a periodic switch between DB and NS. A (z, x 1 ) bifurcation diagram is plotted in Figure 14A for m = −8 (I ext 2 ¼ 0). The trajectory characterizes the periodic switch between DB and NS. Parameter settings in Figure 14 correspond to area II in Figure 12B. Thus, LC and the periodic switch between DB and NS coexist, and are separated by a saddle periodic orbit (S). We  Figure 19. A, Coexistence of SLEs and RSE. The simulations are performed without noise. SLEs and a stable LC coexist for m = 0 and x 0 = −1.6 (I ext2 ¼ 0). SLEs occur through a fold/homoclinic bifurcation. The arrows indicate the direction of trajectories. For easier visualization, we plot generalized coordinates (X, Y) corresponding to (À10x 1 1 x 2 ; 5y 1 ) for seizures (top) and to (À0:3x 1 1 x 2 ; 0:06y 1 ) for LC (bottom). LC is characteristic of RSE. B, C, Time series of SLEs (B) and LC (C). Parameter settings correspond to region X in Figure 34 and to region 15 in Figure 35 Figure 22. A, Coexistence of RSE and the periodic switch between nonoscillatory state and NS. The simulations are performed without noise. Here, the SLE attractor reduces to the periodic switch between nonoscillatory state and NS, which coexists with a stable LC for m = −1 and x 0 ¼ À1:6 (I ext2 ¼ 0). This periodic switch occurs through a fold/fold bifurcation. The arrows indicate the direction of trajectories. For easier visualization, we plot generalized coordinates (X, Y) corresponding to (À10x 1 1 x 2 ; 5y 1 ) for the periodic switch between nonoscillatory state and NS (top) and to (À0:3x 1 1 x 2 ; 0:06y 1 ) for LC (bottom). LC is characteristic of RSE. B, C, Time series of this periodic switch (B) and LC (C). Parameter settings correspond to region VIII in Figure 34 and to region 13 in Figure 35. The coexistence of LC and S can be found in area II (Fig. 12B). plot LC, S, and the periodic switch between DB and NS in a (Y, X, z) phase space (Fig. 25A). We observe that DB occurs when the fast manifold of the SLE collapses to a point that traces out a line under the slow z-evolution. The ictal state thus reduces to a silent activity. LC and the periodic switch between DB and NS coexist in the phase space and the separatrix acts as a barrier between them. Time series are plotted in Figure 25B for the periodic switch between DB and NS, and in Figure 25C for LC. DB and NS are indicated by segment numbers 1 and 3, respectively.
Coexistence of LC and a NS. Figure 17 shows that when decreasing x 0 , the z-nullcline is at the Z-lower branch, which consists of equilibrium points of the normal state (stable nodes). Let m = 0, the equilibrium point is a stable node when x 0 is below −2.1 (Fig. 3, area 5). LC and S coexist when x 0 = −2.1 (Fig. 12A, area II). We plot LC, NS, and S in a (Y, X, z) phase space (Fig. 26A). LC and NS coexist in the phase space, and the separatrix acts as a barrier between them. Fast discharges and NS are indicated by segment numbers 1 and 4, respectively. After the transient seizure-like fast discharges, the Epileptor remains in the NS. Time series are plotted in Figure 26B for NS and in Figure 26C for LC.
Coexistence of LC and a DB. Figure 16D shows that when increasing x 0 , the z-nullcline is at the Z-upper branch, which consists of equilibrium points of DB (stable nodes). Parameter settings in Figure 16D correspond to area II in Figure 12B, where LC and S coexist. We plot LC, DB, and S in a (Y, X, z) phase space (Fig. 27A). The transient NS and DB are indicated by segment numbers 2 and 4, respectively. After the transient normal state, the Epileptor model remains in the depolarization block. Time series are plotted in Figure 27B for DB and in Figure 27C for LC.
Coexistence of LC and a nonoscillatory state. Figure  16A-C shows that when increasing x 0 , the equilibrium point C is at the Z-upper branch when I ext 2 ¼ 0 and at the upper branch above the Z-curve when I ext 2 ¼ 0:45. C (stable focus) is the equilibrium point of the nonoscillatory state. Figure 16, A and B, shows that after the transient normal state, and the transient seizure-like fast discharges, the Epileptor remains in the nonoscillatory state. Figure 16C shows that after the transient normal state, the Epileptor spirals into the equilibrium point of the nonoscillatory state exhibiting damped oscillations and then stabilizes at C. Parameter settings in Figure 16A correspond to area II in Figure 12A, and parameter settings in Figure 16, B and C, correspond to area II in Figure 12B, and then LC and S coexist. Only the nonoscillatory state is plotted in Figure 16A-C. We plot LC, the nonoscillatory state, and S in a (Y, X, z) phase space (Figs. 28A, 29A, 30A). LC (below) and the nonoscillatory state (above) coexist. The nonoscillatory state shown in Figures 28A, 29A, and 30A is the same as that shown in Figure 16A-C, respectively. NS and the nonoscillatory state are indicated  Figure 25. A, Coexistence of RSE and a periodic switch between DB and NS. The simulations are performed without noise. Here, the SLE attractor reduces to the periodic switch between DB and NS, which coexists with a stable LC for m = −8 and x 0 = −1.4 (I ext2 ¼ 0). Trajectory segments are numbered in A and B. DB corresponds to the segment 1, and the NS to the segment 3. The arrows indicate the direction of trajectories. For easier visualization, we plot generalized coordinates (X, Y) corresponding to (À10x 1 1 x 2 ; 5y 1 ) for the periodic switch between DB and NS (top) and to (Àx 1 1 x 2 ; 0:3y 1 ) for LC (bottom). LC is characteristic of RSE. B, C, Time series of the periodic switch between DB and NS (B) and LC (C). Parameter settings correspond to region VII in Figure 34 and to region 12 in Figure 35 by segment numbers 1 and 3, respectively. Time series for each case are plotted in Figures 28B, 29B, and 30B for the nonoscillatory states, and in Figures 28C, 29C, and 30C for LC.

Parameter space of equilibrium points and periodic orbits
We identified seven types of coexisting attractors: a coexistence of SLE and LC (Figs. 18-21), a coexistence of LC and periodic switch between nonoscillatory state and NS (Fig. 22), a coexistence of SLC and LC (Figs. 23, 24), a coexistence of LC and periodic switch between DB and NS (Fig. 25), a coexistence of NS and LC (Fig. 26), a coexistence of DB and LC (Fig. 27), and a coexistence of nonoscillatory state and LC . SLEs are characterized by transitions between ictal and normal states. Depending on m, x 0 and I ext 2 , these transitions occur through a fold/homoclinic bifurcation (Figs. 18, 19), a fold/circle bifurcation (Fig.  20), or a fold/Hopf bifurcation (Fig. 21).
To characterize the different types of the coexisting attractors, we explore a (m, x 0 ) parameter space of periodic orbits and equilibrium points in Figures 31 and 32 for I ext 2 ¼ 0:45.
Parameter space of the Epileptor model for I ext 2 ¼ 0:45. The parameter space is divided into two parts which are separated by a bold line (SNPO bifurcation; Fig. 31). Above the bold line, LC exists but does not exist below it. Twelve areas exist. LC exists in area I and coexists with SLC in area II. Only SLC exists in area III. The nonoscillatory state exists in area IV and coexists with LC in area V. Only SLE exists in area VI and coexists with LC in areas VII and VIII. An SLE occurs through a fold/circle bifurcation in areas VI and VII. An SLE occurs through a fold/homoclinic bifurcation in area VIII. NS exists in area IX and coexists with LC in area X. LC coexists with a chaotic state in areas XI and XII.
For each area, the Epileptor model has different equilibrium points, depending on m and x 0 . To explore this further, we determine the equilibrium points and the periodic orbits of each area (I-XII) by using 32, which shows 18 (1-18) areas: (I) LC: Area I in Figure 31 is composed of areas 1, 7, 8, 9, 10, and 14 in Figure 32. The z-nullcline intersects the ,x 1 .-curve at one periodic orbit, which corresponds to LC. The z-nullcline intersects the (z, x 1 ) curve at different equilibrium points. The equilibrium point is an unstable focus in area 1 and a saddle in area 14. The equilibrium points are one saddle and two unstable foci in area 7. The equilibrium points are one unstable focus and two saddles in area 9. The equilibrium points are three saddles in area 10. The equilibrium points are a saddle, an unstable focus, and an unstable node in area 8. Trajectories exhibit only a fast-slow cyclic behavior (LC), which is plotted in Figures 9D and 10.
(II) LC and SLC: Area II in Figure 31 corresponds to area 2 in Figure 32. The z-nullcline intersects the ,x 1 .-curve at three periodic orbits in area 2, which correspond to LC, S, and SLC, respectively. The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point, which is an unstable focus. LC and SLC coexist. The stable manifold of the saddle periodic orbit (S) separates the basin of attraction For easier visualization, we plot generalized coordinates (X, Y) corresponding to (À35x 1 1 x 2 ; 15y 1 ) for DB (top) and to (À0:5x 1 1 x 2 ; 0:1y 1 ) for LC (bottom). LC is characteristic of RSE. B, C, Time series of NS (B) and LC (C). Parameter settings correspond to region X in Figure 31 and to region 18 in Figure 32. The coexistence of LC and S can be found in area II (Fig. 12A).
of LC and the basin of attraction of SLC. Trajectories exhibit two periodic solutions depending on initial conditions. The first solution corresponds to a fast-slow cyclic behavior (LC). The second solution corresponds to a periodic solution with a small amplitude (SLC).
(III) SLC: Area III in Figure 31 corresponds to area 3 in Figure 32. The z-nullcline intersects the ,x 1 .-curve at one periodic orbit in area 3, which corresponds to SLC. Area three is below the SNPO bifurcation (bold line), which means that SLC exists after a SNPO bifurcation. The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point, which is an unstable focus. The trajectories exhibit only a periodic solution with a small amplitude (SLC).
(IV) Nonoscillatory state: Area IV in Figure 31 corresponds to area 4 in Figure 32. Here, the z-nullcline and the ,x 1 .-curve do not intersect, hence periodic orbits LC, S, and SLC do not exist. The z-nullcline intersects the (z, x 1 ) curve at different equilibrium points. The equilibrium points are a stable node, an unstable focus, and a saddle in area 4. Trajectories remain in the nonoscillatory state.
(V) LC and a nonoscillatory state: Area V in Figure 31 is composed of areas 5 and 6 in Figure 32. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits in areas 5 and 6, which correspond to LC and S. The z-nullcline intersects the (z, x 1 ) curve at different equilibrium points. The equilibrium points are a stable node, an unstable focus, and a saddle in area 5. The equilibrium points are a stable focus, an unstable focus, and a saddle in area 6.
Trajectories either exhibit a fast-slow cyclic behavior (LC) or remain in the nonoscillatory state. The coexistence of LC and the nonoscillatory state is plotted in Figure 28A. The equilibrium points belong to area 6 in Figure 32. Time series are plotted in Figures 28B for the nonoscillatory state and in Figure 28C for LC.
(VI) SLE with a fold/circle bifurcation: Area VI in Figure  31 corresponds to area 11 in Figure 32. Here, the z-nullcline and the ,x 1 .-curve do not intersect, hence periodic orbits LC, S, and SLC do not exist. The z-nullcline intersects the Z-middle branch, which consists of saddles. Thus, a unique saddle equilibrium point exists. Since m 0, an SLE occurs through a fold/circle bifurcation.
(VII) LC and SLE with a fold/circle bifurcation: Area VII in Figure 31 corresponds to area 12 in Figure 32. Here, the equilibrium point is a saddle, and an SLE occurs through a fold/circle bifurcation. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits, which correspond to LC and S. Then LC and an SLE with a fold/circle bifurcation coexist and are separated by a saddle periodic orbit (S). This coexistence is plotted in Figure 20A. Time series are plotted in Figure 20B for SLE and in Figure 20C for LC.
(VIII) LC and SLE with a fold/homoclinic bifurcation: Area VIII in Figure 31 corresponds to area 13 in Figure 32. Here, the equilibrium point is a saddle, and an SLE occurs through a fold/homoclinic bifurcation. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits, which correspond to LC and S. Then LC and SLE with a fold/ homoclinic bifurcation coexist and are separated by a saddle periodic orbit (S). This coexistence is plotted in Figure 18A. Time series are plotted in Figure 18B for SLE and in Figure 18C for LC.
(IX) NS: Area IX in Figure 31 corresponds to area 17 in Figure 32. The z-nullcline and the ,x 1 .-curve do not intersect, hence periodic orbits LC, S, and SLC do not exist. The z-nullcline intersects the (z, x 1 ) curve at the Zlower branch, and then the equilibrium point is a stable node in area 17. The Epileptor remains in NS after a transient period, which is shown in Figure 17 (X) LC and NS: Area X in Figure 31 corresponds to area 18 in Figure 32. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits, which correspond to LC and S. The z-nullcline intersects the (z, x 1 ) curve at the Z-lower branch, which is a stable node in area 18. Then LC and NS coexist and are separated by a saddle periodic orbit (S). This coexistence is plotted in Figure 26A. Time series are plotted in Figure 26B for NS and in Figure 26C for LC.
(XI, XII) LC and a chaotic state: Areas XI and XII in Figure  31 correspond to areas 15 and 16 in Figure 32. The z-nullcline intersects the ,x 1 .-curve at three periodic orbits in area 15. Then LC, S, and SLC coexist. S acts as a separatrix between the basin of attraction of LC and the basin of attraction of SLC. Then trajectories exhibit either LC or SLC in area 15. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits in area 16, which correspond to LC and S. The z-nullcline intersects the Z-middle branch which consists of saddles. An SLE occurs through a fold/ homoclinic bifurcation in area 16. Then LC and SLE coexist in area 16. S acts as a separatrix between the basin of attraction of LC and the basin of attraction of SLE in area 16.
We calculate the times between two successive spikes, which are referred to as the interspike intervals (ISIs) with respect to x 0 in Figure 33A for m = 1 and in Figure 33B for m = 1.5. The system is chaotic when the ISI is irregular. Therefore, LC and the chaotic state coexist in areas 15 and 16 and are separated by a saddle periodic orbit (S).
The parameter space of the Epileptor model includes attractors as LC, SLE, and normal state, but neither DB nor SLE with a fold/Hopf bifurcation. In fact, using bifurcation analysis, we demonstrated that DB and Hopf bifurcation appeared when decreasing I ext 2 . To characterize these attractors, we explore a (m, x 0 ) parameter space of periodic orbits and equilibrium points in Figures 34 and  35 for I ext 2 ¼ 0.
Parameter space of the Epileptor model for I ext 2 ¼ 0. The parameter space is divided into two parts, which are separated by a bold line (SNPO bifurcation) (Fig. 34). Above the line, LC exists and does not exist below it. Thirteen areas exist. LC exists in area I and coexists with SLC in area II. Only SLC exists in area III. A DB exists in area IV and coexists with LC in area V. A periodic switch between a DB and a NS exists in area VII and coexists with LC in area VIII. A nonoscillatory state coexists with LC in area VI, and a periodic switch between a nonoscillatory state and a NS coexists with LC in area IX. LC and SLE coexist in areas X and XI. an SLE occurs through a fold/Hopf bifurcation in area X and through a fold/homoclinic bifurcation in area XI. A normal state exists in area XII and coexists with LC in area XIII.
For each area, the Epileptor model has different equilibrium points depending on m and x 0 . To analyze this further, we determine the equilibrium points and the periodic orbits of each area (I-XIII) by using 35, which shows 20 (1-20) areas: (I) LC: Area (I) in Figure 34 is composed of areas 1, 7, and 17 in Figure 35. The z-nullcline intersects the ,x 1 .-curve at one periodic orbit, which is LC. The z-nullcline intersects the (z, x 1 ) curve at different equilibrium points. The equilibrium point is an unstable focus in area 1 and a saddle in area 17. The equilibrium points are one saddle and two unstable foci in area 7. Trajectories exhibit only a fast-slow cyclic behavior (LC).
(II) LC and SLC: Area (II) in Figure 34 is composed of areas 2 and 18 in Figure 35. The z-nullcline intersects the ,x 1 .-curve at three periodic orbits in areas 2 and 18, which correspond to LC, S, and SLC, respectively. Then LC and SLC coexist in areas 2 and 18. The stable manifold of a saddle periodic orbit (S) separates the basin of attraction of LC and the basin of attraction of SLC. The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point, which is an unstable focus in area 2 and a saddle in area 18. The trajectories exhibit two periodic solutions depending on the initial conditions. The first solution corresponds to a fast-slow cyclic behavior (LC). The second solution corresponds to a periodic solution with a small amplitude (SLC). Moreover, the SLC behavior depends on the stability of the equilibrium point. The SLC behavior is plotted in Figure 13A for area 18 and in Figure 13B for area 2. The coexistence of LC and SLC is plotted in a phase space, in Figure 23A for area 18, and in Figure 24A for area 2. Time series are plotted in Figure 23B for SLC and in Figure 23C for LC, for area 18. Time series are plotted in Figure 24B for SLC and in Figure 24C for LC, for area 2.
(III) SLC: Area III in Figure 34 corresponds to area 3 in Figure 35. The z-nullcline intersects the ,x 1 .-curve at one periodic orbit in area 3, which corresponds to SLC. Area 3 is below the SNPO bifurcation (bold line), which means that SLC exists after an SNPO bifurcation occurs. The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point, which is an unstable focus. Trajectories exhibit only a periodic solution with a small amplitude (SLC). The SLC behavior is plotted in Figures 13B and 24B. (IV) DB: Area IV in Figure 34 is composed of areas 4 and 8 in Figure 35. Here, the z-nullcline and the ,x 1 .-curve do not intersect, hence periodic orbits LC, S, and SLC do not exist. The z-nullcline intersects the (z, x 1 ) curve at different equilibrium points. The equilibrium point is a stable node in area 8. The equilibrium points are a stable node, an unstable focus, and a saddle in area 4. The Epileptor remains in DB after a transient period, which is plotted in Figures 16D and 27B.
(V) LC and a DB: Area V in Figure 34 is composed of areas 5 and 9 in Figure 35. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits, which correspond to LC and S. The z-nullcline intersects the (z, x 1 ) curve at different equilibrium points. The equilibrium points are a stable node, an unstable focus, and a saddle in area 5. The equilibrium point is a stable node in area 9. Trajectories either exhibit a fast-slow cyclic behavior (LC) or enter into a DB. The coexistence of LC and DB is plotted in Figure  27A. Time series are plotted in Figure 27B for DB and in Figure 27C for LC.
(VI) LC and a nonoscillatory state: Area VI in Figure 34 is composed of areas 6, 10, and 11 in Figure 35. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits, which correspond to LC and S. The z-nullcline intersects the (z, x 1 ) curve at different equilibrium points. The equilibrium points are a stable focus, an unstable focus, and a saddle in area 6. The equilibrium point is a stable focus in area 10 and a stable node in area 11. Trajectories either exhibit a fast-slow cyclic behavior (LC) or remain in a nonoscillatory state. The coexistence of LC and the nonoscillatory state is plotted in Figures 29A and 30A, which corresponds to area 10 in Figure 35. Time series are plotted in Figures 29B and 30B for the nonoscillatory state, and in Figures 29C and 30C for LC.
(VII) Periodic switch between DB and NS: Area VII in Figure 34 corresponds to area 12 in Figure 35. Here, the z-nullcline and the ,x 1 .-curve do not intersect, hence LC, S, and SLC do not exist. The z-nullcline intersects the Z-middle branch which consists of saddles. Then only a saddle equilibrium point exists in area 12. The ictal state of SLE reduces to DB with decreasing m. Periodic switch between DB and NS occurs through a fold/fold bifurcation, and is plotted in Figures 14 and 25B.
(VIII) LC and a periodic switch between DB and NS: Area VIII in Figure 34 corresponds to area 13 in Figure 35. The z-nullcline intersects the ,x 1 . curve at two periodic orbits, which correspond to LC and S. The equilibrium point is a saddle, and the ictal state of an SLE is reduced to DB, with decreasing m. Periodic switch between DB and NS occurs through a fold/fold bifurcation, and coexists with LC. Both are separated by a saddle periodic orbit (S). The coexistence of LC and periodic switch between  (1) and the nonoscillatory state (final state) by (3). The equilibrium point of nonoscillatory state exists, which is a stable focus. After a transient NS, the Epileptor spirals into the equilibrium point (stable focus) and remains in the nonoscillatory state. The arrows indicate the direction of trajectories. For easier visualization, we plot generalized coordinates (X, Y) corresponding to (À35x 1 1 x 2 ; 15y 1 ) for the nonoscillatory state (top) and to (À0:5x 1 1 x 2 ; 0:1y 1 ) for LC (bottom). LC is characteristic of RSE. B, C, Time series of the nonoscillatory state (B) and LC (C). Parameter settings correspond to region V in Figure 34 and to region 6 in Figure 32 DB and NS is plotted in Figure 25A. Time series are plotted in Figure 25B for periodic switch between DB and NS, and in Figure 25C for LC.
(IX) LC and a periodic switch between nonoscillatory state and NS: Area IX in Figure 34 corresponds to areas 14 in Figure 35. The z-nullcline intersects the ,x 1 . curve at two periodic orbits, which correspond to LC and S. Here, there is a unique saddle equilibrium point. The SLE attractor reduces to a periodic switch between nonoscillatory state and normal state, which occurs through a fold/ fold bifurcation in area 14, and coexists with LC. Both are separated by a saddle periodic orbit (S). This coexistence is plotted in Figure 22A. Time series are plotted in Figure  22B for the periodic switch between nonoscillatory state and NS, and in Figure 22C for LC.
(X) LC and SLE with a fold/Hopf bifurcation: Area X in Figure 34 corresponds to area 15 in Figure 35. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits, which correspond to LC and S. There is a unique saddle equilibrium point, and an SLE occurs through a fold/Hopf bifurcation (Fig. 7). Then LC and SLE with a fold/Hopf bifurcation coexist and are separated by a saddle periodic orbit (S). This coexistence is plotted in Figure 22A. Time series are plotted in Figure 21B for SLE and in Figure 21C for LC.
(XI) LC and SLE with a fold/homoclinic bifurcation: Area XI in Figure 34 corresponds to area 16 in Figure 35. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits, which correspond to LC and S. A unique saddle equilibrium point exists, and an SLE occurs through a fold/homoclinic bifurcation in area 16 (Fig. 6). Then LC and SLE coexist and are separated by a saddle periodic orbit (S). This coexistence is plotted in Figure 19A. Time series are plotted in Figure 19B for SLE and in Figure  19C for LC.
(XII) A NS: Area XII in Figure 34 corresponds to area 19 in Figure 35. Here, the z-nullcline and the ,x 1 .-curve do not intersect, hence periodic orbits LC, S, and SLC do not exist. The z-nullcline intersects the (z, x 1 ) curve at the Zlower branch, and then the equilibrium point is a stable node in area 19. The Epileptor remains in NS after a transient period, which is shown in Figure 17. There are 12 regions separated by a SNPO bifurcation (bold line). LC exists above and it does not exist below. LC exists in area I. Area II shows bistability of LC and SLC. In area III, only SLC exists. Only a nonoscillatory state with damped oscillation exists in area IV, and coexists with LC in area V. Only SLE exists in area VI, and coexists with LC in areas VII and VIII. SLE occurs through a fold/circle bifurcation in areas VI and VII, and through a fold/homoclinic bifurcation in area VIII. NS exists in area IX and coexists with LC in area X. In areas XI and XII, LC coexists with a chaotic state which is periodic in area XI, and it does not in area XII.
(XIII) LC and a NS: Area XIII in Figure 34 corresponds to area 20 in Figure 35. The z-nullcline intersects the ,x 1 .-curve at two periodic orbits, which correspond to LC and S. The z-nullcline intersects the (z, x 1 ) curve at the Z-lower branch, and hence the equilibrium point is a stable node in area 20. Then LC and NS coexist and are separated by a saddle periodic orbit (S). This coexistence is plotted in Figure 26A. Time series are plotted in Figure  26B for NS and in Figure 26C for LC.

Epileptor behaviors as a function of m and x 0
Different areas exist in the parameter space of the Epileptor model as m and x 0 vary (Figs. 31, 34). We show the behavior of the SLE attractor in Figure 36A, a, b, b 1 , and c 1 ; the periodic switch between DB and NS in Figure  36 A, a 1 ; the periodic switch between nonoscillatory state and NS in Figure 36 A, c; the nonoscillatory state in Figure  36B, a; the normal state in Figure 36B, b; s and the chaotic state in Figure 37.
SLE: it occurs through: Figure 36A, a and b, and corresponds to areas VI and VII, respectively, in Figure 31. Figure 36A, b 1 , and corresponds to area XI in Figure 34. Figure 36 A, c 1 , and corresponds to area X in Figure 34.  Figure 32. The Epileptor parameter space of equilibrium points and periodic orbits with respect to the parameters m and x 0 (I ext2 ¼ 0:45). There are 18 areas separated by a SNPO bifurcation (bold line). LC exists above and it does not exist below. The znullcline intersects the (z, x 1 ) curve at different equilibrium points. The equilibrium point is an unstable focus for areas 1 through 3. The equilibrium points are a stable node, an unstable focus and a saddle in areas 4 and 5. The equilibrium points are a stable focus, an unstable focus, and a saddle in area 6. The equilibrium points are one saddle and two unstable foci in area 7. The equilibrium points are a saddle, an unstable focus, and an unstable node in area 8. The equilibrium points are one unstable focus and two saddles in area 9. The equilibrium points are three saddles in area 10. The equilibrium point is a saddle for areas 11 through 16. The equilibrium point is a stable node in areas 17 and 18. In areas 1, 7, 8, 9, 10, and 14, only LC exists. Area 2 presents bistability of LC and SLC. In area 3, only SLC exists. In area 4, only a nonoscillatory state exists. In areas 5 and 6, LC and a nonoscillatory state coexist. Only SLE exists in area 11, and coexists with LC in areas 12 and 13. SLE occurs through a fold/circle bifurcation in areas 11 and 12, and through a fold/homoclinic bifurcation in area 13. NS exists in area 17 and coexists with LC in area 18. In areas 15 and 16, LC coexists with a chaotic state which is periodic in area 15, and it does not in area 16. Periodic switch between DB and NS: The switch occurs through a fold/fold bifurcation in Figure 36A, a 1 , and corresponds to area VII in Figure 34, which shows a periodic switch between DB and NS.

A fold/Hopf bifurcation in
Periodic switch between a nonoscillatory state and NS: The switch occurs through a fold/fold bifurcation in Figure 36A, c, and corresponds to area IX in Figure 34, which shows a periodic switch between nonoscillatory state and NS.
Nonoscillatory state: The Epileptor remains in the nonoscillatory state after a transient seizure-like fast discharges (Fig. 36A). Parameter settings correspond to area V in Figure 31 and to area 6 in Figure 32.
Depolarization block: The Epileptor remains in the DB after a transient NS (Fig. 27A). Parameter settings correspond to area V in Figure 34 and to area 5 in Figure 32.
Normal state: The Epileptor remains in the NS after a transient seizure-like fast discharges (Fig. 36B, b). Parameter settings correspond to area IX in Figure 31 and to area 17 in Figure 32.
Chaotic state: The chaotic state exists in areas XI and XII (I ext 2 ¼ 0:45; Fig. 31). Time series of the Epileptor model c , subsystem 1 c 1 , and subsystem 2 c 2 are plotted in Figure   37. Figure 37A correspond to area XI in Figure 31, and Figure 37B correspond to area XII in Figure 31.

Subsystem 1
Here, we present results on the analysis of the subsystem 1 dynamics without coupling. We analytically determine the equilibrium points and use tables to classify them, according to the trace and the determinant of the Jacobian matrix. We present the different tables as a function of the parameter m.

Analysis of subsystem 1
Subsystem 1 equilibrium points and stability. We analytically find the equilibrium points (x 1 , y 1 ) by solving the following equations:    Figure 35. The Epileptor parameter space of equilibrium points and periodic orbits with respect to the parameters m and x 0 (I ext2 ¼ 0). There are 20 regions separated by a SNPO bifurcation (bold line). LC exists above and it does not exist below. The z-nullcline intersects the (z, x 1 ) curve at different equilibrium points. The equilibrium point is an unstable focus for areas 1 through 3. The equilibrium points are a stable node, an unstable focus, and a saddle in areas 4 and 5. The equilibrium points are a stable focus, an unstable focus, and a saddle in area 6. The equilibrium points are one saddle and two unstable foci in area 7. The equilibrium point is a stable node in areas 8, 9, 11, 19, and 20. The equilibrium point is a stable focus in area 10. The equilibrium point is a saddle for areas 12 through 18. In areas 1, 7, and 17, only LC exists. Areas 2 and 18 present bistability of LC and SLC. In area 3, only SLC exists. In areas 4 and 8, only DB exists. In areas 5 and 9, LC and DB coexist. In areas 6, 10, and 11, LC and nonoscillatory state coexist. Only a periodic switch between DB and NS exists in area 12, and coexists with LC in area 13. A periodic switch between a nonoscillatory state and NS coexists with LC in area 14. LC and SLE coexist in areas 15 and 16. SLE occurs through a fold/Hopf bifurcation in area 15 and through a fold/homoclinic bifurcation in area 16. NS exists in area 19 and coexists with LC in area 20.
c À dx 2 1 À f 1 ðx 1 ; 0Þ À z 1 I ext 1 ¼ 0 where x 1 is a solution of: c À dx 2 1 1 ðm 1 0:6ðz À 4Þ 2 Þx 1 À z 1 I ext 1 ¼ 0 if • If d 1 ! 0, Equation 38 has two solutions, x 1 1 and x 2 1 are written as: Let X = z−4, then z 2 À 1; minðXÞ 1 4 [ ½maxðXÞ 1 4; 1 1½ and X is a solution of We find minðXÞ 1 4 and maxðXÞ 1 4 as m varies in Figure  38 To determine the equilibrium points stability, we analyze the eigenvalues of the following Jacobian matrix J: J is defined at the equilibrium point (x 1 , y 1 ). To find the eigenvalues of J, we solve the characteristic equation: where Tr(J) and Det(J) are the trace and the determinant of the matrix J, respectively. The trace and the determinant of J are given by: DetðJÞ ¼ 2dx 1 À Eðx 1 ; zÞ: An equilibrium point is stable if all the real parts of the eigenvalues of J are negative (Izhikevich, 2007). We can classify the equilibrium points according to the trace and the determinant of J (Izhikevich, 2007). We determine the intervals on which the Tr(J) and Det(J) are negative and positive.
The roots of Tr(J) are as follows: Þ=3a gif x 1 , 0 The roots of Det(J) are as follows: where X is a solution of Equations 40 and 41. We conclude that an equilibrium point x 1 (x 1 ≥ 0) exists (see Eq. 42) if: We find the stability of equilibrium points 8z 2 À 1; zðDet 1 Þ and 8z 2 ½zðDet 2 Þ; 1 1½ in Tables 2-Tables 5,  depending

on m.
We graphically determine the equilibrium points in a phase plane. The equilibrium points lie at the intersection of the x 1 -and y 1 -nullclines. The x 1 -nullcline ( _ x 1 ¼ 0) corresponds to a cubic curve for x 1 , 0 and a straight line for The cubic curve is written in the form: and the straight line is written as follows: y 1 ¼ Àðm 1 0:6ðz À 4Þ 2 Þx 1 1 z À I ext 1 : The y 1 -nullcline ( _ y 1 ¼ 0) corresponds to a parabola given by: We plot the phase plane of the subsystem 1 in Figure  39A for z = 3.1 and in Figure 39B for z = 0. We        . The (x 1 , y 1 ) phase plane of the subsystem 1. Possible intersections of x 1 -(black cubic curve 8x 1 ,0 and dashed straight line 8x 1 !0) and y 1 -(green parabola) nullclines depending on z and m. Trajectories are plotted without noise starting from different initial conditions (black dot). The arrows indicate the direction of trajectories. The equilibrium points where they exist are labeled by red squares for stable nodes, and black stars for saddles. A, z = 3.1, three equilibrium points coexist: a stable node (bottom), a saddle (middle), and a stable focus for m = 0. The stable focus becomes unstable when m = 1 and m = 1.5, surrounded by a stable limit cycle. The limit cycle radius increases as m is increased. B, z = 0, one equilibrium point exists which is an unstable focus surrounded by a stable limit cycle for m = 0, and a stable focus for m = −10. C, m = 1.5, one equilibrium point exists for z = 2.7, and three equilibrium points exist for z = 4: a stable node, a saddle, and an unstable focus. The stable limit cycle does not surround the unstable focus, and then is broken through a homoclinic bifurcation. consider only the negative x 1 -axis for the cubic curve, and plot a straight line of the x 1 -nullcline (8x 1 ! 0) for different m.
When z = 3.1, then z À I ext 1 ¼ 0 (Fig. 39A). We plot a straight line of the x 1 -nullcline for m = 0, m = 1, and m = 1.5. Then the x 1 -and y 1 -nullclines intersect at three points for each m value. We plot four trajectories with different initial conditions (i 1 , i 2 , i 3 , and i 4 ). Trajectories i 1 and i 2 converge to one equilibrium point, which is a stable node. Trajectories i 3 and i 4 converge to one equilibrium point, which is a stable focus when m = 0 (Fig. 39A, m = 0). By increasing (m = 1 and m = 1.5), the stable focus loses its stability, and the equilibrium point is an unstable focus. Trajectories i 3 and i 4 converge to a stable limit cycle, which surrounds the unstable focus. The radius of the limit cycle increases as m increases (Fig. 39A, m = 1 and m = 1.5).
We conclude that when z = 3.1, three equilibrium points coexist. The first one is a stable node and the second one is a stable or an unstable focus, depending on m. The third equilibrium point is a saddle. The stable manifold of the saddle corresponds to a separatrix between the first and second equilibrium point. The subsystem 1 undergoes a supercritical Andronov-Hopf bifurcation, H, as m increases at m(H) = 0.514 [i.e., solution of Tr(J) = 0; Eq. 45].
When z = 0, then z À I ext 1 decreases and the x 1 -nullcline moves downward in the phase plane (m = −10; Fig. 39B). We plot a straight line of the x 1 -nullcline (8x 1 ! 0) for m = 0 and m = −10. Then the x 1 -and y 1 -nullclines intersect at one equilibrium point, which is a stable focus for m = −10 and an unstable focus for m = 0 (Fig. 39B). The stable node and saddle disappear. We plot trajectories in the phase plane, which converge to a stable focus for m = −10 and to a stable limit cycle for m = 0. The limit cycle surrounds the unstable focus. The subsystem 1 undergoes a supercritical Andronov-Hopf bifurcation, H, at m (H) = −8.6 [i.e., solution of Tr(J) = 0; Eq. 45].
We now fix m on two values to identify the equilibrium points as z varies.
• Let m = 0, three equilibrium points coexist when z = 3.1: a stable node, a saddle, and a stable focus (Fig. 39A). When z = 0, then only an unstable focus exists (Fig. 39B). Hence, as z decreases, a stable focus becomes unstable, and the subsystem 1 undergoes a supercritical Andronov-Hopf bifurcation, H. Moreover, the stable node and saddle disappear when z = 0 (Fig. 39B). In fact, the stable node and saddle approach each other as z decreases and coalesce at z % 2.9 (figure not shown), which corresponds to zðx 1 ðDet 1 ÞÞ (see Eq. 48). Decreasing z further, the stable node and saddle disappear through a saddle-node bifurcation.
• Let m = 1.5, three equilibrium points coexist when z = 3.1: a stable node, a saddle, and an unstable focus (Fig. 39A). When z decreases to z % 2.9, then a saddlenode bifurcation occurs (figure not shown). When z = 2.7, the stable node and saddle disappear (Fig. 39C). In contrast, when z = 4 (≥ 3.1), three equilibrium points coexist: a stable node, a saddle, and an unstable focus (Fig. 39C). A stable limit cycle surrounds the unstable focus when z = 3.1 (Fig. 39A), and it does not when z = 4 (Fig. 39C). In fact, the stable limit cycle approaches as z increases to the saddle, touches the saddle, and then terminates. Hence, when z = 4, the stable limit cycle disappears through a homoclinic bifurcation, and the trajectories converge to a stable node.
We conclude that the stability of equilibrium points depends on z and m. In addition, we find three bifurcation types: a saddle-node bifurcation, a homoclinic bifurcation, and a supercritical Hopf bifurcation.
Subsystem 1 bifurcation diagram. We find the equilibrium points as z varies in a (z, x 1 ) bifurcation diagram plotted in Figure 40 for m = 0. The (z, x 1 ) curve comprises two geometrical shapes. The right one comprises two branches: the lower (dash-dotted) branch consists of saddles; and the upper (dashed) branch consists of unstable foci. When decreasing z, the lower and upper branches collide in a saddle-node bifurcation, SN. The left geometrical shape corresponds to a Z-curve known from fast-slow subsystems (Ermentrout and Terman, 2010). Left and right shapes are separated by SN 2 and SN, such that zðSN 2 Þ ¼ I ext 1 1 1 and zðSNÞ ¼ maxðXÞ 1 4 (Eq. 42). The discriminant d 1 of Equation 38 is negative between SN 2 and SN. Lower z SN Figure 40. The (z, x 1 ) subsystem 1 bifurcation diagram with respect to z when m = 0. On the left, the plot shows a Z-shaped curve 8z ,zðSNÞ. Z-lower (solid) and Z-middle (dash-dotted) branches consist of stable nodes and saddles, respectively. Zupper (solid) sub-branch consists of stable foci and Z-upper (dashed) sub-branch consists of unstable foci. The two subbranches are separated by a Hopf bifurcation H. The Z-lower and Z-middle branches collide as z decreases in a saddlenode bifurcation SN 1 . The Z-upper and Z-middle branches collide as z increases in a saddle-node bifurcation SN 2 . The Zshaped curve comprises two branches 8z !zðSNÞ: one (dashed) consists of unstable nodes and another (dash-dotted) consists of saddles. Decreasing z, they collide in a saddlenode bifurcation, SN. and upper branches (right shape) consist of saddles and unstable nodes, respectively 8m. Hence, in what follows, we plot only the Z-curve (left shape) in the next bifurcation diagrams. We plot a (z, x 1 ) bifurcation diagram for m = 0 in Figure 41A. The geometrical shape of the bifurcation diagram corresponds to a Z-curve, which comprises four branches. Z-lower (solid) and Z-middle (dash-dotted) branches consist of stable nodes and saddles, respectively. Decreasing z, Z-lower and Z-middle branches collide in a saddle-node bifurcation point, SN 1 . The Zmiddle branch acts as a separatrix between Z-lower and Z-upper branches.
The Z-upper branch comprises two sub-branches separated by a Hopf bifurcation, H: one sub-branch (solid) consists of stable foci and another (dashed) consists of unstable foci (Fig. 41A). The first sub-branch corresponds to a nonoscillatory state, and the second sub-branch corresponds to an oscillatory state. The subsystem 1 then switches as z decreases from nonoscillatory to oscillatory states.
An unstable focus is surrounded by a stable limit cycle, which is limited by red (maxðx 1 Þ) and green (minðx 1 Þ) curves. The stable limit cycle radius increases as z is decreased. The stable limit cycle terminates in a Hopf bifurcation, H, written in the form:  Figure 41. A, B, The (z; x 1 ) subsystem 1 bifurcation diagram with respect to z (z ,zðSNÞ), when m = 0 (A) and m = 2 (B). Z-lower (solid) and Z-middle (dash-dotted) branches consist of stable nodes and saddles, respectively. Z-upper branch consists of unstable foci (B), and it is divided into two sub-branches separated by a Hopf bifurcation, H (A): one (solid) consists of stable foci and another (dashed) consists of unstable foci. The Z-lower and Z-middle branches collide as z decreases in a saddle-node bifurcation SN 1 . The Z-upper (solid for A and dashed for B) and Z-middle branches collide as z increases in a saddle-node bifurcation SN 2 . A stable limit cycle ends at a Hopf bifurcation (A) and at a homoclinic bifurcation (B). The curves maxðx 1 Þ; minðx 1 Þ and ,x 1 Þ. correspond to the maximum, minimum, and averaged values along periodic orbits.
We plot a (z, x 1 ) bifurcation diagram for m = 2 in Figure  41B. Z-lower (solid) and Z-middle (dash-dotted) branches consist of stable nodes and saddles, respectively. Z-upper (dashed) branch consists only of unstable foci, which are surrounded by stable limit cycles. The stable limit cycle amplitude is limited by red (maxðx 1 Þ) and green (minðx 1 Þ) curves. The stable limit cycle terminates as z increases in a homoclinic bifurcation HB. Decreasing z, Z-lower and Zmiddle branches collide in a saddle-node bifurcation SN 1 .  Increasing z, Z-upper and Z-middle branches collide in a saddle-node bifurcation SN 2 (Fig. 41B, inset). SN 1 and SN 2 correspond to saddle-node bifurcation points, and HB corresponds to a homoclinic bifurcation point. SN 1 and SN 2 points are given by the following: which corresponds to zðx 1 ðDet 1 ÞÞ (see Eq. 48), 8x 1 ,0), where X is a solution of Equations 40 and 41. We find saddle-node bifurcation points SN 1 and SN 2 , and a Hopf bifurcation point, H, in Figure 42 as m varies. Figure 42 shows the following three cases: • For m 0.29, then zðHÞ ,zðSN 1 Þ, and the subsystem Fast-slow subsystem equilibrium points. We analytically find the equilibrium points by solving the Equations 37 and 38, where z is a solution of _ z ¼ 0 (Eq. 18). We determine the stability of the equilibrium points by analyzing the Jacobian matrix J written as: J is defined at the equilibrium point (x 1 , y 1 , z). We find the eigenvalues of J by solving the characteristic equation: detðJðx 1 Þ À l IÞ ¼ l 2 À TrðJÞl 1 DetðJÞ: The trace, Tr(J), and the determinant, Det(J), of J are given by: DetðJÞ ¼ 2dx 1 LðzÞ À Eðx 1 ; zÞLðzÞ 1 rsKðx 1 ; zÞ: We graphically determine the equilibrium points using the nullclines. The x 1 -nullcline corresponds to the Equations 50 and 51, and the y 1 -nullcline corresponds to Equation 52.
The z-nullcline (_ z ¼ 0, Eq. 18) corresponds to a straight line for z ≥ 0 written in the form: and to a curve for z , 0 written as: x 1 ðzÞ ¼ ððz 1 0:1z 7 Þ=sÞ 1 x 0 : The equilibrium points lie at the intersection of the z-, x 1 -, and y 1 -nullclines. We graphically find the equilibrium points using the bifurcation diagram of the subsystem 1 (Fig. 41). We plot the Z-curve using the intersection of the x 1 -and y 1 -nullclines for each z. Therefore, the equilibrium points of the fast-slow subsystem lie at the intersection of the z-nullcline and the Z-curve (Fig. 43). A bifurcation diagram of the subsystem 1 is plotted in Figure 43A for m = 0 and in Figure 43B for m = 2. When x 0 = −1.6, then the znullcline is at the Z-middle branch for m = 0 (Fig. 43A) and m = 2 (Fig. 43B), and then only a saddle exists. When x 0 = −1.9, then the z-nullcline is at the Z-middle branch for m = 2, and the equilibrium point is a saddle (Fig. 43B).
The z-nullcline moves downward (x 0 decreases) or upward (x 0 increases) in the bifurcation diagram, and intersects the (z, x 1 ) curve at different points. To determine these points as x 0 varies, we plot a (x 0 , x 1 ) diagram in Figure 44A for m = 0 and in Figure 44B for m = 2. The curve consists of equilibrium points of the fast-slow subsystem. Blue (bottom solid), black (plus sign markers), red Figure 44. A, B, Finding equilibrium points of the fast-slow subsystem with respect to x 0 , for m = 0 (A) and m = 2 (B). Blue (bottom solid), black (plus sign markers), red (dashed), and green (top solid) branches consist of stable nodes, saddles, unstable foci, and stable foci, respectively. SN 1 and SN 2 correspond to saddle-node bifurcation points and H to a Hopf bifurcation point.
(dashed), and green (top solid) branches consist of stable nodes, saddles, unstable foci, and stable foci, respectively. The branch of stable foci exists for m = 0 but does not for m = 2. SN 1 and SN 2 correspond to saddle-node bifurcation points, and H to a Hopf bifurcation point (Fig. 44).
Fast-slow subsystem behavior. The stability of equilibrium points depends on m and x 0 . Different behavioral patterns are discovered as m and x 0 vary: a resting state, a fold/fold bifurcation, a fold/homoclinic bifurcation, a fold/Hopf bifurcation and a periodic solution (Fig. 45). To observe this dynamical change, we integrated the fastslow subsystem equations for different m and x 0 . We plot time series of the fast-slow subsystem in Figure 45, and trajectories in Figure 46, as m and x 0 vary. NS: We plot a (z, x 1 ) bifurcation diagram for m = 0 in Figure 46D. Let x 0 = −2.6, the z-nullcline is at the Z-lower branch. The equilibrium point is a stable node. When a trajectory is at the Z-upper (dashed) sub-branch, then it exhibits an oscillatory solution, which terminates as z increases in a Hopf bifurcation H. Then, the trajectory exhibits a nonoscillatory solution on the Z-upper (solid) subbranch, which terminates as z increases in a saddle-node bifurcation point, SN 2 , and the trajectory switches to the Z-lower branch. Decreasing z, the trajectory continues to C at which z stabilizes. Time series are plotted in Figure  45, E and E1. If m = 2, then the trajectory exhibits an oscillatory solution on the Z-upper branch, which terminates as z increases in a homoclinic bifurcation and the trajectory switches to the Zlower branch. The trajectory continues as z decreases to C, at which z stabilizes (figure not shown).
Fold/Hopf bifurcation: We plot a (z, x 1 ) bifurcation diagram for m = 0.5 in Figure 46A. When m = 0.5, then zðHÞ.zðSN 1 Þ (Fig. 42). Let x 0 = −1.6, the z-nullcline is at the Z-middle branch. The equilibrium point is a saddle. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation, SN 1 . Then, the trajectory switches to the Zupper (dashed) sub-branch exhibiting an oscillatory solution, which terminates as z increases in a Hopf bifurcation, H. The trajectory exhibits a nonoscillatory solution on the Z-upper (solid) sub-branch, which terminates as z increases in a saddle-node bifurcation, SN 2 , and then switches to the Z-lower branch. Therefore, the transitions between Z-lower and Z-upper branches occur through a fold/Hopf bifurcation. Time series are plotted in Figure 45F.
Fold/fold bifurcation: We plot a (z, x 1 ) bifurcation diagram for m = 0 in Figure 43A. When m = 0, then zðHÞ , zðSN 1 Þ (Fig. 42). Let x 0 = −1.6, the z-nullcline is at the Z-middle branch (Fig. 43A). The equilibrium point is a saddle. If a trajectory is at the Z-upper (dashed) subbranch, then it exhibits an oscillatory solution, which terminates as z increases in a Hopf bifurcation H. The trajectory exhibits a nonoscillatory solution on the Z-upper (solid) sub-branch, which terminates as z increases in a saddle-node bifurcation SN 2 , and then switches to the Zlower branch. Decreasing z, the stable node disappears through a saddle-node bifurcation SN 1 and the trajectory switches to the Z-upper (solid) branch, which consists of stable foci. Therefore, the transitions between Z-lower and Z-upper branches occur through a fold/fold bifurcation. Time series are plotted in Figure 45A. Fold/homoclinic bifurcation: We plot a (z, x 1 ) bifurcation diagram for m = 2 in Figures 43B and 46B. Here, a Hopf bifurcation points does not exist (Fig. 42). Let x 0 = −1.9, the z-nullcline is at the Z-middle branch. The equilibrium point is a saddle. When a trajectory is at the Z-upper (dashed) branch, it exhibits an oscillatory solution, which terminates as z increases in a homoclinic bifurcation HB, and then switches to the Z-lower branch. Decreasing z, the stable node disappears through a saddle-node bifurcation, SN 1 , and the trajectory switches to the Z-upper branch. Therefore, the transitions between Z-lower and Z-upper branches occur through a fold/homoclinic bifurcation. Time series are plotted in Figure 45D.
Periodic solution (SLC 1 ): We plot a (z, x 1 ) bifurcation diagram for m = 2 in Figure 46C. Let x 0 = -1.7, the z-nullcline is at the Z-middle branch. The equilibrium point is a saddle. When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation, SN 1 . Then, the trajectory switches to the Zupper (dashed) branch exhibiting as z increases an oscillatory solution which stabilizes at C. The z-nullcline intersects the ,x 1 .-curve at C, hence the system solution is periodic and corresponds to a stable limit cycle, with a small amplitude denoted by SLC 1 . Time series are plotted in Figure 45C. Figure 45B shows the existence of a periodic solution with a large amplitude when m = 2 and x 0 = −1.6. Here, the equilibrium point is a saddle (Fig. 44B). The solution is a stable limit cycle with a fast-slow cyclic behavior, which corresponds to LC.

Finding LC
LC and SLC 1 are stable limit cycles with large and small amplitudes, respectively. Here we show when periodic orbits exist and how they evolve, determining their stability.
Sensitivity to initial conditions. The Epileptor model is sensitive to initial conditions. In fact, integrating the fast-slow subsystem in Equations 18-20 from an initial condition results in different behaviors as m and x 0 are varied ( Fig. 45A-F, left). We considered another initial condition and used numerical integration techniques to solve the fast-slow subsystem equations. Figure 45A1-F1, (right), shows two different solutions as m and x 0 vary. The first solution corresponds to a stable LC and appears in Figure 45, A1-D1, and F1. The second one corresponds to a normal state, which exists when x 0 = −2.6 ( Fig. 45E 1 ).
We now fix m and x 0 , and compare the solutions of the fast-slow subsystem as the initial conditions vary: • When m = 0 and x 0 = −1.6, two solutions coexist: a fold/fold bifurcation (Fig. 45A) and LC (Fig. 45A 1 ).
Therefore, the fast-slow subsystem is sensitive to initial conditions and exhibits a bistability of LC and SLC 1 , a fold/fold bifurcation, a fold/Hopf bifurcation or a fold/homoclinic bifurcation.
Finding periodic orbits. For x 0 = −2.6, the normal state exists and LC does not (Fig. 45E,E 1 ). However, when m = 2 and x 0 = −1.6, only the LC exists (Fig. 45B,B 1 ). We explain this using the averaging method. We introduce ,x 1 ., which is the average value of x 1 -coordinate associated with the subsystem 1, given by: We plot ,x 1 . in a (z, x 1 ) bifurcation diagram (Fig.  43). Periodic orbits lie at the intersection of the z-nullcline and the ,x 1 .-curve. A periodic orbit branch is referred to the Z-upper branch, which is surrounded by periodic orbits. This branch and the ,x 1 .-curve terminate as z increases in a Hopf bifurcation for m 1 (Fig. 43A) and in a homoclinic bifurcation for m . 1 (Fig. 43B).
Let x 0 = −1.6, Figure 43 shows that the z-nullcline intersects the ,x 1 .-curve at two points: LC and S when m = 0, and only at LC when m = 2. However, LC and S coexist when m = 2 and x 0 = −1.9 (Fig. 43B).
We graphically determine the periodic orbit by using the Pontryagin's averaging technique (Shilnikov and Kolomiets, 2008). We introduce the following slow averaged nullcline: where the periodic orbits correspond to the zeros of the ,_ z.. We determine the slow averaged nullcline by substituting the x 1 -coordinates by their averaged values in the z dynamics Equation 18 and setting the right-hand side of the equation equal to zero (Eq. 65). The slow averaged nullcline depends on x 0 and z. We determine the periodic orbit stability by using the following derivative: that represents the dynamics of the averaged equation. A periodic orbit is stable when Equation 66 is negative and is unstable when Equation 66 is positive.
We plot the ,_ z.-curve in Figure 47A for m = 0 ( Figure  47A,a), and m = 2 ( Figure 47A,b). z is a control parameter. The maximal value of z corresponds to a Hopf bifurcation point H for m 1 and to a homoclinic bifurcation point HB for m . 1. The periodic orbits are the zeros of the ,_ z.. A periodic orbit is stable when the ,_ z.-curve decreases over z at the given zero and is unstable when the ,_ z. curve increases. Stable and saddle periodic orbits are labeled as black circles and squares respectively.
We conclude that LC and S are stable and saddle periodic orbits, respectively (Fig. 43). When m = 2 and x 0 = −1.7, then two zeros of ,_ z. are LC and S. The third zero is a stable periodic orbit, which we denote by SLC 1 .
Evolution of periodic orbits. LC and S are stable and saddle periodic orbits, respectively. SLC 1 is a stable periodic orbit, which exists for m = 2. We now analyze how the periodic orbits evolve when x 0 varies slowly. We plot a (x 0 , z p ) bifurcation diagram of periodic orbits in Figure 47B for m = 0 ( Figure 47B,a), and m = 2 ( Figure 47B, b). We localize the periodic orbits at z p . x 0 is a control parameter.
• Let m = 0, only LC exists [red bottom (1) markers] for large x 0 . Decreasing x 0 , S (black dotted squares) appears. Stable (LC) and saddle (S) periodic orbits approach each other as x 0 decreases, collide in a saddlenode periodic orbit bifurcation, and then fades. LC and S disappear through a SNPO bifurcation (Fig. 47B,a).
• Let m = 2, only LC exists for large x 0 . Decreasing x 0 , S and SLC 1 [blue top (1) markers] appear. Hence, LC and SLC 1 coexist and are separated by S. Decreasing x 0 , SLC 1 disappears. Further decreasing x 0 , LC and S disappear through an SNPO bifurcation (Fig. 47B,b).
We can now explain the solutions found in Figure 45. We plot in the (x 0 , z p ) bifurcation diagram the equilibrium points, labeled as dots (Fig. 47B). Blue (segment 1), black (segment 2), green (segment 3), and purple (segment 4) dots correspond to stable nodes, saddles, stable foci, and unstable foci, respectively. For instance, when m = 0  (1) markers. Equilibrium points are labeled as dots. Blue (segment 1), black (segment 2), green (segment 3), and purple (segment 4) dots correspond to stable nodes, saddles, stable foci and unstable foci, respectively. and x 0 = −2.6, a stable node exists and LC does not exist, whence trajectories exhibit a normal state in Figure 45, E and E1. We discuss all solutions found in Figure 45 in the sections: "Monostability in the fast-slow subsystem" and "Coexistence in the fast-slow subsystem".
Periodic switch between a DB and a NS in the fast-slow subsystem When m 1, a Hopf bifurcation H partitions the Zupper branch into two sub-branches: one consists of an unstable focus (8z ,zðHÞ) and another consists of a stable focus (8z.zðHÞ; Fig. 41A). When decreasing m, the imaginary part of the complex-conjugate eigenvalues corresponding to the stable focus goes to zero, and then the stable focus becomes a stable node. We plot a (z, x 1 ) bifurcation diagram when m = −8 in Figure 48A. The Zlower branch consists of stable nodes, which are the equilibrium points of the NS. The Z-upper branch consists of stable nodes (z.zðSN 1 Þ.zðHÞ) which are responsible for the silent activity, similar to a DB. The z-nullcline is at the Z-middle branch which consists of saddles (x 0 = −1.6). When a trajectory is at the Z-lower branch, the stable node disappears as z decreases through a saddle-node bifurcation SN 1 and then the trajectory enters into a depolarization block, which terminates as z increases in a saddle-node bifurcation, SN 2 . Thus, the transitions between the Z-lower and Z-upper branches are reduced to a periodic switch between a NS and a DB. Time series are plotted in Figure 49C.

Stabilizing equilibrium points in the fast-slow subsystem
The transitions between Z-lower and Z-upper branches occur when the equilibrium point is a saddle (x 0 = −1.6). When x 0 is increased, the fast-slow subsystem stabilizes the equilibrium point of the Z-upper branch, and when x 0 is decreased, the fast-slow subsystem stabilizes the equilibrium point of the Z-lower branch (normal state).
Stabilizing the equilibrium point of the nonoscillatory state. Let x 0 = −0.5. We plot a (z, x 1 ) bifurcation diagram when m = 0 in Figure 48C. The z-nullcline is at the Z-upper branch which consists of a stable focus. Decreasing z, the stable node disappears through a saddle-node bifurcation SN 1 and the fast-slow subsystem switches to the Z-upper branch (dashed) exhibiting an oscillatory solution, which terminates in a Hopf bifurcation, H. The fast-slow subsystem exhibits a nonoscillatory solution at the Z-upper branch (solid) and continues to C, at which z stabilizes. The imaginary part of the eigenvalues of C is responsible for the oscillations around the stable focus, which is the equilibrium point of the nonoscillatory state. Time series are plotted in Figure 49A for the nonoscillatory state.
Stabilizing the equilibrium point of the DB. Let x 0 = −0.5. We plot a (z, x 1 ) bifurcation diagram when m = −8 in Figure 48D The z-nullcline is at the Z-upper branch which consists of stable nodes. Decreasing z, the stable node disappears through a saddle-node bifurcation, SN 1 , then the fast-slow subsystem switches to the Z-upper branch (solid) exhibiting DB and continues to C, at which z stabilizes. The imaginary part of the eigenvalues of C goes to zero, and therefore the oscillations around C disappear. Here, the stable node (C) is the equilibrium point of DB. Time series are plotted in Figure 49D for the DB.
Stabilizing the equilibrium point of the NS. Let x 0 = −2.4. We plot a (z, x 1 ) bifurcation diagram when m = 0 in Figure 48B. The z-nullcline is at the Z-lower branch which consists of stable nodes. When a trajectory is at the Z-upper branch (dashed), it exhibits an oscillatory solution, which terminates as z increases in a Hopf bifurcation, H, and then a nonoscillatory solution on the Z-upper branch (solid), which terminates as z increases in a saddle-node bifurcation SN 2 . The fast-slow subsystem continues to C, at which z stabilizes. The stable node (C) is the equilibrium point of the NS. Time series are plotted in Figure 49B for the NS.
Monostability in the fast-slow subsystem Figure 45 shows that the fast-slow subsystem presents a monostability in two cases: Case 1: Only LC exists. ,_ z. has a simple zero for m = 2 and x 0 = −1.6, which is LC (Fig. 47A,b). Then only LC exists, to which all trajectories converge (Fig. 45,B and B 1 ).
Case 2: Only normal state exists. The z-nullcline and ,x 1 .-curve do not intersect when m = 0 and x 0 = −2.6 ( Fig. 46D). Then LC disappears through an SNPO bifurcation (Fig. 47B,a). As a consequence, only a stable node exists to which all trajectories converge (Fig. 45,E and E 1 ).
Coexistence in the fast-slow subsystem Figure 45 shows that LC coexists with SLE or SLC 1 . Using bifurcation analysis, we identify all coexisting attractors in the fast-slow subsystem.
Coexistence of LC and SLE. A coexistence occurs when ,_ z. has two zeros, which are LC and S. LC coexists with SLE, which occurs through a fold/Hopf bifurcation (Fig. 45,F and F1) or a fold/homoclinic bifurcation (Fig. 45,D and D1). The equilibrium point is a saddle. A saddle periodic orbit (S) acts as a separatrix between LC and SLE (Fig. 47B).
Coexistence of LC and SLC 1 . ,_ z. has three zeros, which are LC, S, and SLC 1 (Fig. 47A,b). LC coexists with a periodic solution (SLC 1 ; Fig. 45C and C 1 ). The equilibrium point is a saddle. A saddle periodic orbit (S) acts as a separatrix between LC and SLC 1 (Fig. 47B,b).
Coexistence of LC and a NS. Figure 47B shows the coexistence of LC and a stable node (blue dots), which is the equilibrium point of the normal state. Therefore, LC and the normal state coexist and are separated by S. Let m = 0 and x 0 = −2.4, time series are plotted in Figure 49D for the normal state and in Figure 49D 1 for LC. The bifurcation diagram is plotted in Figure 48B. Only the normal state is shown.
Coexistence of LC and a nonoscillatory state. (Figure 47B,a), shows the coexistence of LC and a stable focus (green dots), which is the equilibrium point of the nonoscillatory state. Then, LC and the nonoscillatory state coexist and are separated by S. Time series are plotted in Figure 49A for the nonoscillatory state and in Figure  49A 1 for LC. The bifurcation diagram is plotted in Figure  48C. Only the nonoscillatory state is shown.
Coexistence of LC and a periodic switch between nonoscillatory state and NS. When m is decreased and the equilibrium point is a saddle, the fast-slow subsystem switches between the nonoscillatory state and the NS through a fold/fold bifurcation (Fig. 45A). (Figure 47B,a), shows the coexistence of LC and a saddle (black dots), which is the equilibrium point of this periodic switch. Therefore, LC and the periodic switch between nonoscillatory state and NS coexist and are separated by S. Time series are plotted in Figure 45A for the periodic switch between nonoscillatory state and NS, and in Figure 45A 1 for LC.
Coexistence of LC and a periodic switch between DB and NS. When m is further decreased and the equilibrium point is a saddle, the fast-slow subsystem switches between DB and NS (Fig. 48A). Time series are plotted in Figure 49C for the periodic switch between DB and NS, and in Figure 49C1 for LC. Parameters m and x 0 are the same, only initial conditions change. Then, LC and the periodic switch between DB and NS coexist.
Coexistence of LC and a DB. When m is decreased, the fast-slow subsystem enters into DB (Fig. 49C). When x 0 is increased, the fast-slow subsystem remains in DB (Fig. 48D). Time series are plotted in Figure 49D for DB, and in Figure 49D 1 for LC. Hence, LC and DB coexist.

Parameter space of equilibrium points and periodic orbits
To characterize the coexisting attractors, we explore a (m, x 0 ) parameter space of the fast-slow subsystem in Figure 50, using numerical techniques. The parameter space is divided in two parts separated by a bold line (SNPO bifurcation), above which LC exists, while it does exist below it. There are 10 areas. For large values of m and x 0 (area I) only LC exists. The adjacent area II shows bistability of LC and a stable focus. Both attractors are separated by a saddle periodic orbit. In the middle, LC coexists with the SLE attractor separated by a saddle periodic orbit. SLE occurs via a saddle-node/saddle-node bifurcation (area VI) and a saddle-node/homoclinic bifurcation (area VII). In area VIII, the homoclinic bifurcation HB is not completed and gives rise to another (coexisting) stable periodic orbit with a small amplitude, SLC 1 . In area V, DB coexists with LC. In area IX, normal brain activity and LC coexist. In area X, only normal activity exists.
(I) LC: The z-nullcline intersects the ,x 1 .-curve at one point which is LC, and intersects the (z, x 1 ) curve at different equilibrium points. The equilibrium point can be a stable focus, an unstable focus, or a saddle. All trajectories converge to LC. Time series are plotted in Figure 45, b and b 1 , for LC.
(II) Coexistence of LC and a stable focus: The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point, which is a stable focus, and intersects the ,x 1 .-curve at two points: LC and S. The trajectories converge to LC or a stable focus, depending on the initial conditions. A saddle periodic orbit (S) limits the basin of attraction of LC and the basin of attraction of a stable focus. The stable focus is the equilibrium point of the nonoscillatory state. Then, LC and nonoscillatory state coexist and are separated by S.
If a trajectory is in the basin of attraction of the stable focus, then it exhibits a nonoscillatory solution following two scenarios. When the trajectory is at the Z-lower branch, the stable node disappears through a saddlenode bifurcation SN 1 and the trajectory switches to the Zupper branch, which comprises a Hopf bifurcation point, H. If zðHÞ. zðSN 1 Þ (first scenario), then the trajectory exhibits an oscillatory solution, which terminates as z decreases in a Hopf bifurcation, H. Then, the trajectory continues to the stable focus at which z stabilizes. If zðHÞ , zðSN 1 Þ (second scenario), then the trajectory exhibits a nonoscillatory solution, which continues as z decreases to the stable focus at which z stabilizes. The stable focus is the equilibrium point of the nonoscillatory state. The system is bistable on [SN 1 , SN 2 ].
Using analytic techniques, we observe that the stable focus reduces as m decreases to a stable node (i.e., the imaginary part of the complex-conjugate eigenvalues corresponding to the stable focus goes to zero). Here, the stable node is the equilibrium point of the DB, and then DB is present in area II. LC and DB coexist as m decreases. The adjacent area II shows the bistability of LC and a stable focus. Both attractors are separated by a saddle periodic orbit.
In the middle, LC coexists with the SLE attractor separated by a saddle periodic orbit. SLE occurs via a saddle-node/saddlenode bifurcation (area VI) and a saddle-node/homoclinic bifurcation (area VII). In area VIII, the HB is not completed and gives rise to another (coexisting) stable periodic orbit with a small amplitude. In area V, DB coexists with LC. In area IX, normal brain activity and LC coexist. In area X, only normal activity exists.
(V) Coexistence of LC and a periodic switch between DB and NS: The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point, which is a saddle, and intersects the ,x 1 .-curve at two points: LC and S. A Hopf bifurcation point, H, exists for m 1 (Fig. 42). Using analytic techniques (Fig. 42), we can deduce that zðHÞ , zðSN 1 Þ. The transitions between Z-upper and Z-lower branches occur through a fold/fold bifurcation. For small m, the Zupper branch consists of stable nodes, which is the equilibrium point of a DB. The system is bistable on [SN 1 , SN 2 ]. The transitions between Z-upper and Z-lower branches reduce to a periodic switch between DB and NS. The LC and periodic switch between DB and NS coexist and are separated by S.
(VI) Coexistence of LC and SLE with a fold/fold bifurcation: Using numerical techniques, there is a saddle-node bifurcation at both offset and onset seizures. We can deduce that SLE occurs with a fold/fold bifurcation in area VI. The z-nullcline intersects the ,x 1 .-curve at two points: LC and S. As a consequence, LC and SLE with a fold/fold bifurcation coexist and are separated by S. In contrast to area V, the Z-upper branch consists of a stable focus, which is the equilibrium point of the nonoscillatory state. The SLE attractor reduces then to a periodic switch between nonoscillatory state and NS, which occurs through a fold/fold bifurcation, coexisting with LC. Time series are plotted in Figure 45A for a periodic switch between nonoscillatory state and NS, and in Figure 45A 1 for LC.
Using analytic techniques (Fig. 42), we can deduce that as m increases, zðHÞ . zðSN 1 Þ. This means that although there is a saddle-node bifurcation at both offset and onset seizures, there is a Hopf bifurcation occurring during the ictal period before offset seizure. Thus, an SLE occurs as m increases through a fold/Hopf bifurcation in area VI, coexisting with LC. Time series are plotted in Figure 45F for SLE with a fold/Hopf bifurcation, and in Figure 45F1for LC.
(VII) Coexistence of LC and SLE with a fold/homoclinic bifurcation: The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point, which is a saddle, and intersects the ,x 1 .-curve at two points: LC and S. Since m . 1, then the transitions between Z-upper and Z-lower branches occur through a fold/homoclinic bifurcation. For m = 1, a Hopf bifurcation point, H, exists such that zðHÞ . zðSN 1 Þ (Fig. 42). When a trajectory is at the Zupper branch, then it exhibits an oscillatory solution, which touches the Z-middle branch before the Hopf bifurcation, H. Then the transitions between Z-upper and Zlower branches occur through a fold/homoclinic bifurcation for m = 1. LC and a fold/homoclinic bifurcation coexist and are separated by S. Time series are plotted in Figure 45D for fold/homoclinic bifurcation and in Figure  45D 1 for LC.
(VIII) Coexistence of LC and SLC 1 : The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point, which is a saddle, and intersects the ,x 1 .-curve at three points: LC, SLC 1 , and S. Then LC and SLC 1 coexist and are separated by S. If a trajectory is in the basin of attraction of SLC 1 , then it exhibits a periodic solution with a small amplitude. When the trajectory is at the Z-lower branch, the stable node disappears through a saddle-node bifurcation, SN 1 , and the trajectory switches to the Z-upper branch. Then, the trajectory exhibits, as z increases, an oscillatory solution, which stabilizes at SLC 1 . Time series are plotted in Figure 45C for SLC 1 and in Figure 45C 1 for LC.
(IX) Coexistence of LC and a NS: The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point which is a stable node, and intersects the ,x 1 .-curve at two points: LC and S. The stable node is the equilibrium point of the normal state. Then, LC and the normal state coexist and are separated by S. Time series are plotted in Figure  49B) for the normal state and in Figure 49B 1 for LC. The fast-slow subsystem switches to NS according to three scenarios, depending on m. We discuss these scenarios below in area X.
In the following areas, LC and S disappear through a SNPO bifurcation.
(III) DB: The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point, which is a stable node (Fig. 43A). The range of m indicates that a Hopf bifurcation point H exists, with zðHÞ ,zðSN 1 Þ (Fig. 42). When a trajectory is at the Zlower branch, the stable node disappears through a saddle-node bifurcation, SN 1 , and the trajectory switches to the Z-upper branch exhibiting DB and continues to the equilibrium point, at which z stabilizes. Therefore, the fast-slow subsystem remains in the depolarization block.
(IV) Periodic switch between DB and NS: This area exhibits a periodic switch between DB and NS, like area V. In contrast, LC and S disappear through a SNPO bifurcation.
Therefore, the fast-slow subsystem only switches between DB and NS.
(X) NS: The z-nullcline intersects the (z, x 1 ) curve at one equilibrium point, which is a stable node. Time series of this solution are plotted in Figure 45, E and E1. Trajectories switch to the normal state according to three scenarios, depending on m.
• For m , 1 (first scenario), a Hopf bifurcation point, H, exists. If a trajectory is at the Z-upper sub-branch, which consists of unstable foci, it exhibits an oscillatory solution that terminates as z increases in a Hopf bifurcation, H. Then, the trajectory exhibits a nonoscillatory solution on the Z-upper sub-branch, which consists of stable foci. The nonoscillatory solution terminates as z decreases in a saddle-node bifurcation SN 2 and the trajectory switches to the Z-lower branch.
• For m = 1 (second scenario), a Hopf bifurcation point H exists. If a trajectory is at the Z-upper sub-branch, which consists of unstable foci, it exhibits an oscillatory solution, which touches the Z-middle branch before the Hopf bifurcation, H, as z increases. Then the trajectory switches to the Z-lower branch through an HB.
• For m . 1 (third scenario), a Hopf bifurcation point, H, does not exist. If a trajectory is at the Z-upper branch, which consists of unstable foci, then it exhibits an oscillatory solution, which terminates as z increases in a homoclinic bifurcation, and the trajectory switches to the Z-lower branch.
The Z-lower branch consists of stable nodes, which are the equilibrium points of the normal state. At the Z-lower branch, the trajectory continues to a stable node, at which z stabilizes. Therefore, the fast-slow subsystem only remains in the normal state. Time series are plotted in Figure 45, E and E1, for the normal state.
Note that the SNPO curve (bold line) decreases as m increases.

Special cases
The fast-slow subsystem has one equilibrium point, which is a stable focus, an unstable focus, a saddle, or a stable node, depending on m and x 0 . When LC coexists with an attractor and both are separated by S, the trajectories converge to LC or the attractor, depending on initial conditions. When S does not exist (x 0 = −1.6; Fig. 47B,b), then all trajectories converge to LC (Fig. 45B,B1). Figure 47B,a, shows that a stable focus and LC coexist, and S does not for x 0 = −0.1 and x 0 = −0.05. We plot trajectories in a (z, x 1 ) bifurcation diagram for x 0 = −0.1 (Fig.  51A) and x 0 = −0.05 (Fig. 51B). Even if S does not exist, trajectories converge to either LC or a stable focus, depending on the initial conditions.

Subsystem 2
Here we present results of our analysis on the dynamics of subsystem 2 without coupling. We partition the results into two main parts: the equilibrium points and their stability, and bifurcations.

Subsystem 2 equilibrium points
We analytically determine the equilibrium points (x 2 , y 2 ) by solving the following equations: SN. Increasing further I ext 2 , S-middle and S-lower branches collide in a SNIC bifurcation (I ext 2 % 0:385). The limit cycle branch terminates as I ext 2 decreases in a SNIC bifurcation.
We identify three intervals on which the equilibrium points change: • ]-1, SN[: only a stable node exists. • ]SN, SNIC[: a stable node, a saddle, and an unstable focus coexist. Here, the unstable focus is not surrounded by a stable limit cycle.
• ]SNIC,11[: only an unstable focus exists, which is surrounded by a stable limit cycle.
The subsystem 2 is stable in the three intervals. Since the stable limit cycle emerges after a SNIC bifurcation, then there is no bistability in the subsystem 2.
We now link the dynamics of subsystem 2 to its behavior (see c 2 ; Fig. 1, main file) in the Epileptor model, which corresponds to transitions between an oscillatory state and a resting state. To explain these transitions, we plot time series c 2 and equilibrium points of the subsystem 2 in Figure 53. The equilibrium points are labeled as dots. First, we observe that the bifurcation diagram of the subsystem 2 ( Fig. 52B) is repeated over time. Blue, black, and green dots correspond to stable nodes, saddles, and unstable foci, respectively. Red dots, which are surrounded by stable limit cycles, correspond to unstable focus. A SNIC bifurcation occurs at the intersection of blue (stable nodes) and black (saddles) dots. Second, transitions from and to the limit cycle occur through a SNIC bifurcation, and under the evolution of I ext 2 1 0:002g À 0:3ðz À 3:5Þ. In fact, the coupling function g(x 1 ) and z change the input I ext 2 of the subsystem 2 to I ext 2 1 0:002g À 0:3ðz À 3:5Þ (see Epileptor equations). We plot time series of z in Figure 53. Let I ext 2 ¼ 0:45, the subsystem 2 is at the limit cycle branch exhibiting an oscillatory solution (Fig. 52B). When z increases, then I ext 2 1 0:002g À 0:3ðz À 3:5Þ is decreased, and the oscillatory solution terminates in a SNIC bifurcation. Hence, the subsystem 2 exhibits a resting state. When z decreases, then I ext 2 1 0:002g À 0:3ðz À 3:5Þ is increased, and the oscillatory solution emerges through a SNIC bifurcation. We conclude that the transitions between the oscillatory state and the resting state (see c 2 ; Fig. 1, main file) occur under the evolution of g(x 1 ) and z.

Discussion
Epileptic seizures take many forms, as well as transitions from and to these pathological states. The mechanisms involved in epileptic dynamics are assumed to be different according to seizure patterns. We have used a neural mass model of partial seizures the "Epileptor" to analyze seizure dynamics. Using a mathematical approach, we showed that the model can generate different dynamic behaviors including SLEs, as the value of the parameters change. Most computational models tried to reproduce experimental data with biophysical realistic parameters (Kager et al., 2000;Traub et al., 2001;Destexhe, 2008;Cressman et al., 2009;Ullah et al., 2009). In contrast, parameters of the Epileptor model do not directly relate to biological variables. Their meaning can only be inferred based on dynamics. The SLE attractor consists of two states, an oscillatory state (limit cycle) interpreted as ictal activity, and a normal state interpreted as interictal activity. The ictal and normal states coexist and occupy two separated regions in the phase space. The SLE attractor corresponds to a periodic transition between the ictal and normal states. Transitions from ictal to normal states, and vice versa, can autonomously occur under the slow evolution of the variable z and without changing the parameters. The transitions from ictal to normal states occur through a bifurcation of the limit cycle, and the transitions from normal to ictal states occur through a bifurcation of the equilibrium point. The combination of two such bifurcations results in classes of (point-cycle) bursters as proposed by Izhikevich (2000). The Epileptor switches from normal to ictal states through a saddle-node bifurcation (Jirsa et al., 2014). In contrast, to return to the normal state, the Epileptor may undergo different bifurcation types, depending on the value of some selected parameters. There are two main parameters having effects on the Epileptor dynamics. The first, m, controls the Epileptor dynamics during the ictal period, and the second, x 0 , controls the equilibrium points of the Epileptor. Bifurcations that define how the Epileptor qualitatively changes its behavior were identified using the scaling behavior of frequency and amplitude during seizure onset and offset. A taxonomy of seizures (SLEs) assembles the bifurcations into 16 different classes (Jirsa et al., 2014). SLE offsets in the Epileptor show the logarithmic scaling of homoclinic bifurcations that occurred more often (83%) in drug-resistant epileptic patients As a consequence, the fold/homoclinic bifurcation was identified as the predominant class of SLEs (Jirsa et al., 2014). Through a bifurcation diagram, we showed that the Epileptor undergoes a homoclinic bifurcation at seizure offset.
On the other hand, the behavior of SLEs in the 17% of patients with a non fold/homoclinic bifurcation was potentially consistent with two classes: fold/Hopf and fold/fold cycle (Jirsa et al., 2014). We evaluated whether other bifurcations may occur at seizure offset, when changing some parameters of the model. Upon change of the first main parameter, we indeed found that the Epileptor may undergo two different types of bifurcation at seizure offset: Hopf and SNIC bifurcations. Hence, the Epileptor can generate two classes: fold/Hopf and fold/circle. Our analysis demonstrated the existence of another burster, defined by Izhikevich (2000) as a point-point type (i.e., both bifurcations are of equilibrium points), called fold/fold bifurcation. This burster does not pertain to the 16 classes (point-cycle) defined by Izhikevich (2000), and then it does not describe any behavior of SLEs.
We proposed in this study a modification of the slow z variable equation (Eq. 34). This modification resolves the divergences observed in the model behavior for some initial conditions, in particular when the initial z value is negative. More, the modification unraveled the existence of a stable LC whose behavior results from the fast-slow structure of the Epileptor model. To explore this issue theoretically, we explained in this study how the modification of the z variable equation leads the model to generate a fast-slow cyclic behavior, in particular for negative initial z values, using bifurcation analysis and the averaging method. We demonstrated that LC does not destabilize, but disappears through a SNPO bifurcation. Interestingly, LC resembles the RSE and does not exist for certain initial conditions. As a consequence, if we consider LC as a prime candidate for RSE, our analysis of the LC dynamics would contribute to the understanding of the mechanisms underlying RSE and to the prediction of how to treat RSE.
In addition to SLEs and RSE, our analysis unraveled another pathological activity, called DB. Using bifurcation analysis, we explained how the Epileptor enters into DB. Indeed, we showed that reducing the first main parameter, the frequency of oscillations during the ictal period decreases until they disappear and DB occurs. In this case, the SLE attractor reduces to a periodic switch between DB and NS. The Epileptor model can therefore reproduce SLEs, RSE, and a periodic switch between DB and NS. Interestingly, it can remain in DB if we increase the value of the second main parameter m, and in NS if we decrease m.
We predicted that epileptic seizures, refractory status epilepticus, and "normal" brain activity can coexist in the brain, under some conditions. In fact, we demonstrated that the normal and ictal states of the SLE attractor coexist in the phase space, and, for certain conditions, LC also exists below. The transitions between them require a large change of z. The "fold/  Figure 52. Subsystem 2 dynamics. A, The (x 2 , y 2 ) phase plane of the subsystem 2. Possible intersections of x 2 -(cubic curve) and y 2 -nullclines depending on I ext2 . a, I ext2 ¼ 0, one equilibrium point (stable node) exists, to which trajectory converges. b, I ext2 ¼ 0:38, three equilibrium points coexist: a stable node, a saddle, and an unstable focus. Trajectory converges to the stable node. c, I ext2 ¼ 1, one equilibrium point exists that is an unstable focus. Trajectory converges to a stable limit cycle. B, Bifurcation diagram of subsystem 2 with respect to I ext2 . The lower (solid), middle (dash-dotted), and upper (dashed) branches of the S-shaped curve consist of stable nodes, saddles, and unstable foci, respectively. A branch of limit cycles (red) originates at a SNIC bifurcation. The curves above and below correspond to the maximum and minimum values along periodic orbits, respectively.
homoclinic bifurcation" was proposed as the predominant class of SLEs (Jirsa et al., 2014). We demonstrated that the coexistence of LC and SLEs occurs for the predominant class of the SLE attractor, as well as two other classes of the SLE attractor: fold/Hopf and fold/ circle bifurcations. More, we demonstrated that only the SLEs with fold/circle type can exist alone in the phase space, for certain conditions. The Epileptor model presents another type of bistability, which consists in the coexistence of LC and a periodic switch between DB and NS. Under certain conditions, this periodic switch between DB and NS can exist in the phase space without the presence of LC. Therefore, we could deduce that DB, normal brain activity and RSE coexist in the brain. However, epileptic seizures and DB do not coexist in the brain. In fact, the SLEs and periodic switch between DB and NS do not coexist in the phase space. DB was developed during the ictal period when reducing the first parameter value m.
We explored the parameter space of, first, the fastslow subsystem and extend to the whole system Epileptor. For the Epileptor, we explore two parameter spaces, as the parameter control of the subsystem 2 is varied. In one of them, DB appears, when reducing the parameter value. The parameter space can provide pathways to switch among SLEs, RSE, DB, and NS. In fact, there is a region in the parameter space representing the coexistence of LC and DB in the phase space, another where only DB exists, and a third where only LC exists. Based on these parameters, we can then propose how to switch from SLEs to DB to RSE, and how to escape the RSE and DB. Interestingly, there is a region representing the coexistence of LC and NS in the phase space, and another where only the NS exists. Therefore, we could propose how to stay in the NS, and thus how to return to the normal brain activity.
In addition to the stable limit cycle LC, the Epileptor model generates another stable periodic orbit, which we denoted by SLC. This limit cycle has higher frequency and smaller amplitude than LC, and occupies the same region as the SLE attractor. SLC is different from LC; both coexist in the phase space. The transitions between them requires a large change of the z variable. Furthermore, there are two different patterns of SLC depending on the equilibrium point stability of the Epileptor, which can be a saddle as for the SLE attractor, or an unstable focus. SLC mimics the dynamics of epileptic seizures during the ictal period of the SLE attractor. For the SLE attractor, the Epileptor returns to the normal state after a transient epileptic seizures, while for SLC, epileptic seizures do not stop. Hence, we believe that SLC could characterize the status epilepticus, or might have another possible clinical and physiological explanation.
In conclusion, our study provides a better understanding of the Epileptor dynamics, exploring the different behaviors that the model can generate. Three classes can be generated by the Epileptor model, which are: fold/homoclinic, fold/Hopf and fold/circle bifurcations. We predict the coexistence of epileptic seizures (SLEs), RSE, DB, and normal brain activity (NS) in the brain, under a variety of conditions. Furthermore, the model predicts the existence of different paths between these neuroelectric  Figure 53. On the subsystem 2 dynamics. Time series and equilibrium points of subsystem 2 are plotted showing how the subsystem 2 evolves with time when coupled to z and subsystem 1, which is a periodic switch between oscillatory and resting states. On top, time series of Z is plotted. For easier visualization, Z corresponds to z − 2. Blue, black, and green dots correspond to stable nodes, saddles, and unstable foci, respectively. Red dots correspond to unstable foci, which are surrounded by stable limit cycles. Here all simulations were performed without noise.