Analysis of spatially extended excitable Izhikevich neuron model near instability

The article focuses on the issue of a spatiotemporal excitable biophysical model that describes the propagation of electrical potential called spikes to model the diffusion-induced dynamics based on an analytical development of amplitude equations. Considering the Izhikevich neuron model consisting of coupled systems of ODEs, we demonstrate various results of spatiotemporal architecture (PDEs) using a suitable parameter regime. We analytically perform the saddle-node bifurcation and Hopf bifurcation analysis with bifurcating periodic solutions that show the transition phases in the system dynamics. We study different firing patterns both analytically and numerically by the formation of Riccati differential equation. To examine the characteristics of diffusive instabilities, we use Turing amplitude equations by multiscaling method. The instabilities and Turing bifurcation are established using theoretical analysis and numerical simulations. The spatial dynamics has potential effects on the deterministic system as a result of the diffusive matrices with various couplings, and the coupled oscillators with this nearest-neighbor coupling show synchronization measured by the synchronization factor analysis. Our results qualitatively reproduce different phenomena of the extended excitable system based with an efficient analytical scheme.


Introduction
To understand the various functional mechanisms in different brain areas (especially in cortical region), we definitely need the theoretical modeling and numerical simulations with experimental studies [1,17,24,48]. There are a certain number of models for excitable biophysical systems that are capable to reproduce the observed results by mathematical analysis and that are also computationally efficient [17,18,27,31,34]. The models can generate rich dynamics often found in real neurons for electrophysiological activities. When a neuron receives a certain stimulus, it produces action potential (electrical impulses) [16]. Single neurons produce a spike when they are in the vicinity of a bifurcation point from resting state to firing activity. The following results show how the dynamics can be investigated for different characteristics. We are excited to study an excitable spatiotemporal system that supports the wave propagation of nerve impulses [7,37]. From the mathematical point of view, the model consists of a coupled nonlinear reaction diffusion systems of partial differential equations (PDEs) modeling the propagation of electrical potential. It is derived from a system of coupled ordinary differential equations (ODEs) modeling the cell membrane dynamics. The implementation of the proper analytical and efficient numerical techniques is important to study such spatially extended systems.
Over the last few years, diverse firing activities have been explored in some research articles with the help of bifurcation analysis as the effects of electromagnetic induction on different neuron models. Furthermore, many articles have been carried out on the synchronization behavior of coupled neurons [30,35,44,45]. Besides electromagnetic induction, the researchers have studied the coexistence of firing patterns by considering the asymmetric electrical synaptic connection between two neuronal oscillators [36]. Motivated by the above-mentioned articles and some recent studies, that were investigated on spatiotemporal patterns away from the equilibriums, traveling waves [6,28], spiral waves which frequently occur in biological, chemical and other physical systems [12,22,23,32,33,40,41,47], we introduce the solution of the general amplitude equations and prove its well posedness (for Hopf bifurcation [10,11] and Turing instabilities [39,46]). In our study, we consider zeroflux boundary condition to show that the excitable system rests in an isolated condition [26]. Then, we will examine the corresponding system. We are interested to determine the variations in the firing potential as well as the stability analysis of the spatial system.
In this article, we are concerned with the analysis of spatiotemporal dynamics in the Izhikevich model [17,18,29] that plays a major role in the study of the cortical neurons. The famous mathematician and neuroscientist E. M. Izhikevich proposed an excitable two-dimensional (2D) dynamical model that reproduces spiking-bursting activities [3,18,27,38]. Here, we use this simple model that is biologically plausible and computationally efficient. The stability analysis can be observed from bifurcation scenarios where the predominant parameter is I (injected current stimulus). Then, we focus on the local nonlinear excitations in diffusively coupled Izhikevich model for different parameter regimes. Bifurcation scenarios have been widely used to analyze the dynamical responses that are modulated by uncoupled single model and coupled systems with various injected inputs such as synaptic coupling and current stimulus. We consider three firing regimes: fast spiking, tonic spiking and phasic spiking, and verify the firing activities analytically using the formation of Riccati differential equation [13,25]. Bifurcations are performed to show the transition phases of oscillations. Moreover, the dynamics behind it can generally be described through a phase diagram. We present a methodology to characterize the biophysics of this system, which takes into account the temporal structures of complex firings. This provides us to accurately distinguish the fundamental characteristics of the intrinsic dynamics, and the modeling approaches enrich the spiking network-level activities such as synchronization and its measure with synchronization factor [8,21,43]. For Hopf and Turing instabilities, responsible for close to the homogeneous fluctuations and the emergence of stationary and spatial patterns, respectively, the corresponding amplitude equations are referred as Turing amplitude equation (TAE) [42,49] and complex Ginzburg-Landau equation (CGLE; it is responsible for the onset of Hopf instability in the case of homogeneous fluctuations) [2,19].
We used TAE equation to get the insight into the nonequilibrium phenomenon in spatiotemporal systems. The analysis scheme [9,19,49] for the derivation of amplitude equations is based on the multiple space and time scales. The method provides a Taylorseries expansion of the original nonlinear equations with many power operators, and further, it depends on the expansion of the linear and nonlinear terms of a small perturbation parameter close to the onset of diffusive instability. The small parameter measures the deviation from the diffusive instability region. It is used in our case as the reaction diffusion system has diagonal diffusion matrices. The derivation of the procedure based on Kuramoto's method has been investigated to the systems. The coefficients of the amplitude equations obtained in the present article may have their own values, and it may consider to certain dynamical behavior in the diffusive system of experimental conditions. However, a systematic simulation is required to capture the dynamics of excitable models at a sufficient level of accuracy and it requires the suitable discretization scheme. The effects of diffusion in a biophysical extended system and then collective dynamical formation have been studied. The results show that the reaction diffusion system with different coupling charac-teristics that depend on the control parameters, participates in a collective behavior with diffusive instabilities.
The paper is organized as follows: In Sect. 2, we briefly describe the dynamics of the excitable model equations that shows the electrical potential called as spikes and its bifurcation scenarios. We also verify different types of oscillatory activities analytically. In Sect. 3, we verify our bifurcation results analytically. We demonstrate the Turing instability and spatiotemporal dynamics of the model system with some basic preliminaries. We derive the TAE and CGLE and find the coefficients. We report the main features with numerical simulations and analyzed the results to show the qualitative dynamics. Finally, we conclude the results with a discussion in Sect. 4.

Formulation of the excitable model and some preliminaries
This article focuses on the complex dynamics of the two variables excitable Izhikevich system over an appropriate range of parameters. The model construction does not account the biophysical behavior; however, the output dynamics is very realistic and biophysically plausible. The time evolution of such a mathematical model is described by the following set of ODEs [17,18] dp dt = (0.04 where p measures the membrane voltage dynamics and q, the recovery variable, that measures the activation of K + and inactivation of Na + ionic currents, with an after-spike resetting constraint, i.e., when the voltage, p reaches the peak value, the following relation is applied: If p ≥ p peak = 30, then p ← c, q ← q + d. p is expressed in mV (millivolt) scale and time t is in ms (millisecond) scale [17]. The values of the parameters a, b, c and d determine spiking-bursting activities. The resting potential depends on the parameter b, indicating the sensitivity of q to the subthreshold fluctuations of the voltage, p. The parameter a measures the timescale of the recovery variable, q. The parameters c and d control the afterspike reset values of p and q, respectively, for the uncoupled dynamics. The function (0.04 p 2 + 5 p + 140) was derived using the spike initiation dynamics of a cortical cell. The different suitable choices of parameters generate various types of oscillations such as regular spiking, various bursting, chattering and fast spiking, often found in neocortical and thalamic neurons [3,20]. The stimulus currents are injected using the constant variable, I . The equilibrium point E = ( p 0 , q 0 ) of system (1) can be found by solving f 1 ( p 0 , q 0 ) = 0 and and q 0 = bp 0 . The Jacobian matrix of system (1) at equilibrium point E = ( p 0 , q 0 ) is given by where J 11 = 0.08 p 0 +5 , J 12 = −1 , J 21 = ab and J 22 = −a. The equilibrium point E = ( p 0 , q 0 ) of system (1) is locally asymptotically stable if, T race(J ) = J 11 + J 22 < 0 and Det (J ) = J 11 J 22 − J 12 J 21 > 0. We have considered the following parameter sets that produce different qualitative firings for the deterministic model [17,18,38] and depend on the dynamics with = 1 (shown in Fig. 1a, b and c for parameter sets I and II ) of system (1). Set I : a = 1, b = 1.5, c = −60, d = 0; Set II : a = −0.02, b = −1, c = −60, d = 8 [17,18,38]. The numerical simulations of the excitable systems are performed using the fourth-order Runge-Kutta method with a time step of 0.001, and the initial conditions are set to p = −63 and q = bp [17,18,38]. The simulation results with a smaller time step do not show any significant differences. Bifurcation diagrams of the dynamical model are computed using the Mat-Cont software package [4].

Oscillatory activities and Bifurcation analysis
We consider three different regimes, i.e., fast spiking, tonic spiking and phasic spiking with the steady state. The bifurcation analysis of system (1) is presented for parameter sets I and II by considering current stimulus, I as a control parameter. System (1) with parameter set I shows a subcritical Hopf bifurcation and a saddle-node bifurcation at I = −65 and −63.4375, respectively (Fig. 1d). For saddle-node bifurcation, as the magnitude of the parameter I changes, a stable equilibrium corresponding to the resting state condition is approached by an unstable fixed point; then, they coalesce and annihilate each other. Since the resting state no longer exists (bistable regime exists, as shown in Fig. 1d) with the increase in I , the neuron starts to fire fast spiking (Fig.  1a, solid blue line). There is a coexistence of resting and spiking states in the case of saddle-node bifurcation. At lower current stimulus (I < 85) with set II, it has two equilibrium points which collide at I = 85 ( Fig. 1e marked by SN1 point) and vanishes together for higher current stimulus I > 85. The system has a saddle-node bifurcation at I = 85. Another interesting behavior is the system changes its stability from stable regime to unstable at I = 78.9975 ( Fig. 1e marked by HB1 point) with the transition from phasic spiking to tonic spiking regime (Fig. 1b, c, solid blue lines). System (1) shows supercritical Hopf bifurcation (HB1) at I = 78.9975. In Fig. 1d, e, dashed and thick blue lines represent the stable and unstable equilibrium regions. Now, we are interested to find the solution of system (1) analytically to verify our numerical simulations, i.e., different types of oscillatory activities. We begin with the recovery variable q and find the solution of the last equation of system (1) by considering p as a constant between t 0 and t, where t 0 indicates the initial time. The second equation of system (1) can be written as: which is a first-order linear ODE. In this context, we can find the solution of Eq. (2) using the integrating factor (I.F.) = e at and the solution is given by which can be written as: For better result, we use small time step, i.e., Δt = t −t 0 is very small and replace t by t 0 + Δt. Now, we find the approximate solution of first equation of system (1) by considering q as a constant over some small time interval as described earlier. The first equation of system (1) can be written as: which is a first-order nonlinear ODE. It is impossible to find the exact solution for most of the nonlinear ODEs, fortunately above equation is a special type of ODE, known as Riccati differential equation [13,25]. We can find the solution of Riccati differential equation using the transformation r = 1 p− p 1 , where p 1 is a particular solution to dp dt . In this context, the particular solution p 1 can be found from dp dt = 0, which leads to Here, we are interested to find only one particular solution, so we choose where B = 5 + 0.08 p 1 . To find the solution of system (3), first we find the solution of above equation and then use the transformation r = 1 p− p 1 . So, the solution of system (3) is given by provided p 1 is real, i.e., B is real. Complex value of B indicates system (1) has no fixed point, and then, let us consider Then, the solution of system (3) is given by . Now, we plot the membrane voltage variable p with same time step as described earlier, i.e., Δt = 0.001, which is shown in Fig. 1 a-c in dotted red lines. It is clear from Fig. 1 a-c that analytical results are in good agreement with the numerical results.

Coupled Izhikevich model with 1D diffusion
The general 1D reaction diffusion system with the Izhikevich model is represented by the following two equations: with the resetting constraint equation mentioned in system (1), where p = p(x, t) and q = q(x, t) are the unknown functions to be evaluated and the subscripts t and x represent the differentiation with respect to these variables. The initial conditions of the coupled PDEs are considered with the known functions p(x, t = 0) and q(x, t = 0) for x ∈ Ω, and the boundary conditions show zero-flux boundary conditions ∂ p ∂n = ∂q ∂n = 0, x ∈ ∂Ω and t > 0, where n measures the outward normal to ∂Ω, the boundary of the interval and domain, and Ω is the bounded interval or square domain for 1D and 2D diffusion, respectively. In the 1D case, the length of the excitable cable is x = 100. D 11 and D 22 are the strengths of the synaptic couplings with positive values. One can explore the similar dynamics with the large finite length of the excitable cable by introducing more integrating space points in the simulation. We use a finite-difference scheme (with forward Euler method and central difference representation) with zero-flux boundary conditions for numerical simulation of a cable of finite length. The time step δt = 0.001, space step δx = 0.5 and the number of nodes N = x/δx = 200 are considered. The second-order partial derivatives ∂ 2 p ∂ x 2 and , respectively, to find the numerical solution of system (4). The initial conditions are considered with suitable periodic perturbations from the initial conditions as in the deterministic case. It is fixed for all the simulations. Biologically, the zero-flux boundary condition shows that the cell membranes are impermeable at the boundaries and it acts as an isolated cable [26].
To find the solution, we deviate system (4) around the uniform steady-state condition and linearize it around the nontrivial fixed point, and we obtain the characteristic equation assuming the particular solution (4), where c.c. stands for the complex conjugate, λ k is the wave length and k is the wave number along the x direction (x is the directional vector). The Jacobian matrix of Eq. (4) at the nontrivial equilibrium point E = ( p 0 , q 0 ) is given by The stability conditions of The eigenvalues of the Jacobian matrix J D are given by Now, we use the analytical approach to find the Hopf bifurcation point with its nature and the bifurcating periodic solution. In Andronov-Hopf bifurcation case [10,11], it produces a limit cycle from a fixed point or equilibrium solution of the dynamical model of ODEs, when the fixed point changes its stability through a pair of purely imaginary eigenvalues. The bifurcation can be supercritical or subcritical, resulting in stable or unstable limit cycle, respectively. For supercritical Andronov-Hopf bifurcation, the stable equilibrium loses its stability and gives birth to a small-amplitude limit cycle attractor. When the magnitude of the parameter increases, the amplitude of the limit cycle increases and it becomes a full-size spiking limit cycle. First, consider model (1)   System (4) in the absence of diffusion shows Hopf bifurcation if Im(λ) = 0 and Re(λ) = 0 with k = 0, i.e., J 11 + J 22 =0, which gives 1 (0.08 p 0 +5)−a = 0, where λ is the eigenvalue of the Jacobian matrix J . Thus, the critical value for the Hopf instability is given by For parameter sets I and II, we obtain the critical value I H = −65 and 78.9975, respectively, which are in good agreement with the numerical simulations. In the following, we study the stability of periodic solutions for Hopf bifurcation. The characteristic equation of the linearized matrix J is For the occurrence of Hopf bifurcation, 1 (0.08 p 0 + 5) − a = 0, Eq. (6) becomes Further simplification of Eq. (7) gives two purely imaginary roots as λ = ±i √ ab − a 2 = ±iβ. Differentiat-ing Eq. (6) with respect to I , we obtain λ = 0.08 (λ + a) dp 0 d I 2λ + a − 1 (0.08 p 0 + 5) .
Suppose η 1 and η 2 be the eigenvectors corresponding to eigenvalues iβ and −iβ, respectively. Now, we define a nonsingular matrix A = col(Re(η 1 ), Im(η 1 )), We have y 1 = p and y 2 = a √ ab−a 2 p − 1 √ ab−a 2 q. Then, system (1) becomes Writing F i , i = 1, 2 for the nonlinear part of Eq. (9), we have F 1 = 0.04y 2 1 and F 2 = 0.04a √ ab−a 2 y 2 1 . Again, we find F 1 11 = 0.08, F 2 11 = 0.08a √ ab−a 2 and all other F i jk are zero. Thus, from the genericity condition, we have n = 1 16 √ ab−a 2 (−F 1 11 F 2 11 ) = −0.04a ab−a 2 . Therefore, n < 0 for a > 0, which indicates that system (1) exhibits subcritical Hopf bifurcation. System (1) shows supercritical Hopf bifurcation if n > 0, i.e., a < 0. In our study, system (1) shows a subcritical and a supercritical Hopf bifurcation for parameter sets I and II, respectively, which also satisfy our numerical result. Now, we study the Turing instability. In the presence of diffusion, Turing instability breaks the spatial symmetry, i.e., the uniform steady state is stable for system (1); however, it is unstable for system (4). Thus, the conditions for Turing instability are given by (i) T race(J ) = J 11 We can find the threshold values of Turing instability and Turing bifurcation using Det (J D ) = 0 and Im(λ)=0, Re(λ)=0 at k = k cr t = 0, respectively. The threshold value of Turing bifurcation parameter is given by If we consider D 22 as a Turing parameter, then the threshold value is given by One can also evaluate the critical value of Turing instability with respect to D 11 . Interestingly, it is clear from the above calculations that synaptic coupling D 22 plays a major role for Turing instability. In the absence of D 22 , only the membrane voltage spatially interacts in the excitable extended systems, and as a result of it, there is no Turing instability in the reaction diffusion system [5,26]. Now, we derive the procedure to obtain the amplitude equations with the generalized TAE equations that explore the characteristics of the diffusive instability. We only consider constant diffusion scheme as it has interesting phenomenon in excitable spatial system. At the threshold of Turing point, the spatial symmetry of the reaction diffusion system is destroyed and arises a stationary form in time and oscillatory dynamics in space. It is noted that in the vicinity of a bifurcation point, the evolution of a spatial system generates critical slowing down, that can be described by an amplitude equation given in the following.

Turing amplitude equation (TAE)
For Turing instability, the dynamics of system (4) can be described by the amplitude equations, which is known as TAE. The general form of TAE is given by [42,49] where η T , g T and D T are constants. The normalized form of this equation with η T = D T = 1 and g T = −1 can be considered as the limiting condition of CGLE equation [49]. The TAE equation is a valid efficient tool in one spatial dimensional (1D) case. Here, we do not consider the case of wave instability [28]. In the following, we find the value of above parameters. At the onset of Turing instability, for critical wave number k cr t , the critical value T is given by

11
, whereJ 11 = J 11 = 0.08 p o + 5 andJ 12 = J 12 = −1. First, we introduce a control parameter μ 1 as normalized deviation from the critical value T , i.e., μ 1 = − T T . The right and left eigenvectors of the matrix L 0 = (J D )| μ 1 =0 are given by P P = P P p P P q = 1 α e ik cr t x and P P * = Using the expressions J 1 i j = d J i j d d dμ 1 | = T and γ i j = −J 0 i j + 4k 2 cr t D i j , the matrices J 1 and γ are given by respectively. Finally, η T , D T and g T are given by [49] and We can observe that effects of diffusion appear in the expressions of all the coefficients of the TAE, since the parameters α and β depend on diffusive couplings. g T = 0 line separates the region between supercritical and subcritical Turing instability. g T < (>)0 indicates supercritical (subcritical) Turing instability. First, we consider parameter set I with I = −68 and D 11 = 1. The deterministic system (1) is bistable around I = −68 (Fig. 1d), i.e., stable and unstable equilibrium branch coexists. The stable equilibrium point E = (−54.43, −81.645) becomes unstable for system (4), i.e., Turing instability occurs for suitable value of diffusion coefficient D 22 . The threshold value for Turing instability is given by (D 22 ) cr t = 11.081 (from Eq. (11)), i.e., stable equilibrium point becomes unstable for D 22 > 11.081. For better understanding, we have shown the real part of the eigenvalue λ for different values of D 22 in Fig. 2a. In Fig. 2a, Re(λ)=0 line separates the Turing and non-Turing domains.
As we already know that, system (4) shows Turing instability at parameter set I, with I = −68 and D 11 = 1 for D 22 > (D 22 ) cr t = 11.081. Now, we are interested to find the nature of the Turing instability, i.e., whether it is supercritical or subcritical. In Fig 2b, we show the value of g T is negative in Turing region. That indicates the instability is supercritical.

Complex Ginzburg-Landau equation (CGLE)
The dynamics of the spatiotemporal system (4) nearly to the Hopf instability can be described by the amplitude equations, which is known as CGLE [14,15,49]. The general form of CGLE is given by where i denotes the imaginary unit and ∇ 2 is the one or two-dimensional Laplacian operator. When the real parameters C 0 , C 1 and C 2 vary, the complex amplitude W shows rich dynamical behavior. First, we have to find the condition for Hopf instability and then the values of C 0 , C 1 and C 2 using Kuramoto's method [19,49]. The manipulations are done using CGLE equation, and only the coefficient C 1 depends on the diffusion. For the case of Hopf instability, T race(J ) = J 11 + J 22 = 0. The critical value H of at the onset of Hopf instability is given by Jacobian matrix J has purely imaginary eigenvalues ±iω 0 , in the case of Hopf instability, which gives ω 2 0 = J 11 J 22 − J 12 J 21 , i.e., To find the values of C 0 , C 1 and C 2 , we introduce a small parameter μ as normalized deviation from the critical value H , i.e., μ = − H H . The matrix J 0 = J | μ=0 is given by The left and right eigenvectors of the matrix J 0 are given by P * = 1 2 −i a ω 0 , 1 b + i a bω 0 and P = Next, we compute P * J 1 P to find the value of constant C 0 in Eq. (16), which gives P * J 1 P = − a+iω 0 2 = σ 1 + iω 1 . Now, Next, we compute P * D P to find the value of constant C 1 in Eq. (16), which gives P * D P = − 1 2 (D 11 + D 22 + i a ω 0 (D 22 − D 11 )) = d 1 + id 2 . Now, To find C 2 = g 2 g 1 , we calculate the vectors H X X and N X X X. H X X and N X X X are quadratic and cubic terms, respectively [14,49]. H X X is described by = 0 as in f 1 ( p, q), no product term of p and q exists. Now, with μ = 0 the first component (H 0 X X) 1 is given by Similarly, the first component (H 0 XX ) 1 is given by 2 .
As f 1 ( p, q) is a quadratic function and f 2 ( p, q) is a linear function, (N X X X) i = 0, for i = 1, 2. Next, we compute the following expressions: and (H 0X Y + ) 2 = 0. Next, we calculate g as g = g 1 + ig 2 = −P * (2H 0 XY 0 + 2H 0X Y + + 3N 0 X XX ), which gives Now, and Using Eqs. (19), (20) and (21) with equating the complex and real parts, we obtain and Finally, The CGLE is valid if g 1 > 0. Hence, the Hopf bifurcation is supercritical if g 1 > 0 and subcritical if g 1 < 0. The phase waves spreading away from the origin of perturbation and toward the origin of perturbation are known as wave (W) and antiwave (AW) respectively [49]. W and AW domains are separated by the curve C 1 + C 2 = 0. W and AW domains are given by C 1 + C 2 < 0 and C 1 + C 2 > 0, respectively.
It is clear from Eq. (22) that the sign of g 1 depends on a, and supercritical Hopf bifurcation exists for a < 0, which is in good agreement with the above analytical and numerical simulations. Now, we consider parameter sets I and II with I = −65 and I = 78.9975, respectively. We derive the following values H = 1, , C 2 = 4.9048, g 1 = 0.0816 > 0 for parameter sets I and II, respectively. In Fig. 2c, we show the value of C 1 + C 2 is positive for different values of D 11 , corresponding to D 22 for set II. Therefore, only antiwave exists.

Synchronization factor
Next, we have considered parameter set II with I = 78 and D 11 = 0.001. At lower diffusion D 22 = 0.00001, system (4) produces a de-synchronized firing pattern (Fig. 3a). For better understanding, we have presented time series of different arbitrarily chosen oscillators in Fig. 3d. At higher diffusion D 22 = 0.01, we observed more synchronized firing pattern (Fig. 3b-e). Finally, at D 22 = 20, all nodes are completely synchronized, i.e., all nodes show same firing behavior. To understand the oscillators in the reaction diffusion coupled systems. At lower diffusion, the system produces a de-synchronized firing pattern and then converges to a synchronized firing pattern for higher diffusion synchronization behavior, we derive a statistical measure of synchronization, i.e., the synchronization factor R, which is given by the following equation The symbol * indicates the mean value of the variable over time. R takes the value between 0 and 1 [8,43], where the value one indicates that all nodes are synchronized. We have shown the effect of the diffusion coefficient D 22 on synchronization factor R in Fig 4 . Further, if we consider a bursting regime, we obtained similar results with no significant changes. All the nodes are not synchronized for weak coupling. However, with the increase in diffusion coefficient D, the value of R tends to 1, which indicates all the nodes are completely synchronized for higher coupling.

Conclusions
The effects of diffusion in a biophysical extended system and then collective dynamical formation have been an ongoing research theme. The diffusion can enhance the spatial characteristics and even induce different diffusive instabilities. We consider diverse excitabil-ities of the deterministic system and verify the analytics of the model using Riccati differential equation. Next, the diffusion-induced characteristics can be evaluated on the basis of amplitude equations that provide a mathematical analysis of reaction diffusion systems nearly to the onset of instabilities [28,49]. We report the effects of the diffusion on the system dynamics in the vicinity of Hopf and Turing instability. The 1D cable, which is in steady state for uncoupled system, can produce irregular oscillatory patterns. We also study the effects of diffusion on synchronization in oscillatory regime. For the TAE equation, it is classified the Turing instability and divided it into two domains. We verify the numerical simulations with proper analytical treatment for different firing patterns, Hopf bifurcation and Turing instability. The diffusive coupling affects the coefficient C 1 , whereas the term C 2 is independent of diffusion for CGLE equation. Note that the value of the diffusion coefficient has strong effects on k T and hence on the diffusive instabilities. The diffusive couplings contribute to all the coefficients of the amplitude equations. The approach demonstrates flexibility of our methodology for various firing regions supported by Hopf bifurcation analysis. Further, it requires the modification to observe the effects on spatially extended fast-slow excitable dynamical system. The present article reports the emergence of two types of bifurcations in the temporal dynamics that control different oscillations exhibited by the system. Neuronal oscillations reflect various rhythmic activities and can be generated either by mechanisms in single neurons or by network of neurons connected in a complex fashion. The results further show the dynamical behavior of the transition phases in single model and the diffusion-induced system by a suitable analytical treatment.