# Spiral wave dynamics in a neuronal network model

Diogo L. M. Souza1,∗, Fernando S. Borges2,3, Enrique C. Gabrick1,4,5, Lucas E.
Bentivoglio1, Paulo R. Protachevicz6, Vagner dos Santos8,9, Ricardo L. Viana7,
Ibere L. Caldas6, Kelly C. Iarosz1,6,8, Antonio M. Batista1,6,10, Jürgen
Kurths4,5 1 Graduate Program in Sciences, State University of Ponta Grossa,
84030-900, Ponta Grossa, PR, Brazil
2 Department of Physiology and Pharmacology, State University of New York
Downstate Health Sciences University, Brooklyn, New York, USA
3 Center for Mathematics, Computation, and Cognition, Federal University of
ABC, 09606-045 São Bernardo do Campo, SP, Brazil
4 Department of Physics, Humboldt University Berlin, Newtonstraße 15, 12489
Berlin, Germany
5 Potsdam Institute for Climate Impact Research, Telegrafenberg A31, 14473
Potsdam, Germany
6 Institute of Physics, University of São Paulo, 05508-090, São Paulo, SP,
Brazil
7 Department of Physics, Federal University of Paraná, 81530-000, Curitiba,
PR, Brazil
8 Exact and Natural Sciences and Engineering, UNIFATEB University Center,
84266-010, Telêmaco Borba, PR, Brazil.
9 Department of Physics, State University of Ponta Grossa, 84030-900, Ponta
Grossa, PR, Brazil
10 Department of Mathematics and Statistics, State University of Ponta Grossa,
84030-900, Ponta Grossa, PR, Brazil<EMAIL_ADDRESS>
###### Abstract

Spiral waves are spatial-temporal patterns that can emerge in different
systems as heart tissues, chemical oscillators, ecological networks and the
brain. These waves have been identified in the neocortex of turtles, rats, and
humans, particularly during sleep-like states. Although their functions in
cognitive activities remain until now poorly understood, these patterns are
related to cortical activity modulation and contribute to cortical processing.
In this work, we construct a neuronal network layer based on the spatial
distribution of pyramidal neurons. Our main goal is to investigate how local
connectivity and coupling strength are associated with the emergence of spiral
waves. Therefore, we propose a trustworthy method capable of detecting
different wave patterns, based on local and global phase order parameters. As
a result, we find that the range of connection radius ($R$) plays a crucial
role in the appearance of spiral waves. For $R<20$ $\mu$m, only asynchronous
activity is observed due to small number of connections. The coupling strength
($g_{\rm syn}$) greatly influences the pattern transitions for higher $R$,
where spikes and bursts firing patterns can be observed in spiral and non-
spiral waves. Finally, we show that for some values of $R$ and $g_{\rm syn}$
bistable states of wave patterns are obtained.

Keywords: firing patterns, spiral waves, network, neural oscillations,
chimera, pyramidal neurons

##  1 Introduction

It has been identified astonishing spatiotemporal patterns which can be seen
in 2D and 3D networks with non-local interactions, named spiral wave chimera
states. It consists of synchronous oscillators which rotate around a
desynchronized core, which is called phase singularity (PS) [1]. This pattern
is seen in a variety of fields. In heart tissues [2, 3, 4] spiral waves have
been extensively studied due their correlation with some heart issues, such as
ventricular tachycardia and ventricular fibrillation. Davidenko et al. [5]
demonstrated the emergence of spiral waves in sheep and dog epicardial muscle
by means of potentiometric dye with charge-coupled device imaging technology.
The authors showed that sometimes the phase singularity drifted away from its
original position and dissipated in the tissue border, this drift was
associated with a Doppler shift. Supression of spiral waves in the atrial
cardiomyocytes can be done by means of optogenetics [6]. Spiral waves are able
to emerge in ecological network composed of diffusible prey-predator species
locally coupled [7]. Totz et al. [8] observed the emergence of spiral wave
chimeras in large populations of coupled chemical oscillators, they have
studied the motion and splitting of the asynchronous core.

The chimera is defined as the coexistence of coherent and incoherent domains
when similar elements interact [9]. The notion of chimera states has been
discussed in the last decade [10]. Majhi et al. [11] reviewed important
aspects of such dynamics and highlighted the different types, as well as their
relevance in a biophysical context. Chimera states play an important role in
empirical brain networks, specially in the cerebral cortex [12]. It has been
found in networks composed of coupled neurons, such as integrate-and-fire [14,
13], Hodgkin-Huxley [9], Hindmarsh-Rose [15], FitzHugh-Nagumo [17, 16],
Morris-Lecar [18], Nekorkin maps [19] and others. The spatiotemporal patterns,
in which coherence coexists with incoherence, have been reported in local [20,
21] and non-local network configurations [22, 23].

Many studies have improved the understanding of chimera states in neuronal
networks. The parameter regions for chimera states in a two coupling
configuration of identical thermally sensitive Hodgkin-Huxley neurons was
characterized in Ref. [9]. Majhi [24] investigated chimera states in uncoupled
neurons induced by multilayer interactions. They found that the competition
between electrical and chemical synapses can lead to chimeras. Hövel et al.
[25] studied the coexistence of chimera and traveling waves, as well as multi-
stability, in the parameter space of a neuronal system of excitability type I.
Chemical and electrical synapses can affect the features of chimera states and
traveling waves [26]. The coexistence of patterns can emerge due to electrical
and chemical coupling in a multi-weighted of the memristive Fitzhugh-Nagumo
network [16]. It was demonstrated the importance of such interactions to
obtain synchronous and desynchronous behaviours coexisting [27]. Chimera
states were also reported in time-delay networks [28, 29].

Neuronal networks can exhibit spiral wave spatial patterns when in a sleep-
like state [31, 30]. It is able to modify ongoing activity in the cortex [30,
32]. Local excitatory interactions are of great importance in the generation
of spiral waves [30]. The role of spiral waves in the neocortex is not yet
well understood, thus simulations and experimental methods are crucial to
better understand spiral wave dynamics. The emergence of spiral waves in
computational simulations has been observed in different models throughout the
time. Ma et al. [33] studied the robustness of the spiral wave when noise is
applied to ion channels of the standard Hodking-Huxley model. They showed that
the spiral wave does not sustain its activity, whether the noise increases or
the fraction of active channels is below a threshold. The transmission delay
enhances the coherence of spiral waves in a noisy network [34]. Spiral waves
can be seen in non-locally coupled maps for some coupling parameters [35].
Santos et. al. [36] demonstrated that chimera spiral waves with multi-cores
appear in regular and fractal networks. Their study showed that the initial
condition of the adaptation current has a crucial role in the emergence of
spiral waves.

The formation of spiral waves are not only observed in neuronal simulations,
but also in experimental data [37, 39, 38]. The effects of spiral waves in
cognitive activity are not well defined, although it influences the
frequencies [30]. Prechtl et al. [38] used visual stimuli to register the
appearance of spiral waves in the turtle cortex. Furthermore, Huang et al.
[31] captured the emergence of spiral waves in rat cortex and its influence in
coordinating oscillations of neuronal populations. Spiral waves have been
identified in the neocortex during sleep-like states and pharmacologically
induced oscillations [30]. It is theorized that such spatial pattern plays an
important role in organizing and modulating healthy and pathological brain
activities, for instance, epilepsy [31, 30, 38].

In this work, we propose a neuronal network composed of pyramidal cells, that
are identified as excitatory neurons. Pyramidal neurons are the most common
cells in the brain, and the name is due to its characteristic shape [40]. This
neuron type is also the main one of the excitatory neuron family in the brain
[40]. Our network model is inspired by a two dimensional slice of pyramidal
layer of the hippocampus [41]. The model considered to describe a single
neuron membrane potential is the adaptive exponential integrate-and-fire
(aEIF) [42]. Besides the simplicity of this model and its low computational
cost, it exhibits a great biophysical accuracy [43, 44]. Due to the lack of
reliable tools able to characterize chimera spiral waves, we propose a method
to identify chimeras spiral wave, and analyze their dependence on the coupling
radius ($R$) and the synapse coupling strength ($g_{\rm syn}$). Our results
show that the spiral waves emerge for different combinations of $R$ and
$g_{\rm syn}$. We observe spiral waves when the value of the global phase
order parameter is low ($Z_{\rm g}<0.7$), the local phase order parameter
value is high ($Z_{\rm L}\approx 1$) and the amount of phase singularity is
not greater than a maximum value ($P_{\rm max}\leq 20$). We observe for low
$R$ values the desynchronized patterns. Moreover, there are many pattern
transitions and bistable regions in which spiral and non-spiral waves are
possible for higher $R$ values.

The paper is organized as follows: In Section 2, we present our neuronal
network model. In Section 3, we show the diagnostic tools. Section 4 presents
our results. In the last Section, we draw our conclusions.

##  2 Neuronal Network Model

We consider a network composed of coupled adaptive exponential integrate-and-
fire (aEIF) neurons [42]. The dynamics of the neuronal network is given by

$\displaystyle C\frac{dV_{i}}{dt}$ $\displaystyle=$ $\displaystyle-g_{\rm
L}(V_{i}-E_{\rm L})+g_{\rm L}{\Delta}_{\rm T}\exp\left(\frac{V_{i}-V_{\rm
T}}{{\Delta}_{\rm T}}\right)-w_{i}+I$ $\displaystyle+$ $\displaystyle(V_{\rm
rev}-V_{i})\sum_{j=1}^{N}g_{j}M_{ij},$
$\displaystyle\tau_{w}\frac{dw_{i}}{dt}$ $\displaystyle=$ $\displaystyle
a(V_{i}-E_{\rm L})-w_{i},$ (2) $\displaystyle\tau_{g}\frac{dg_{i}}{dt}$
$\displaystyle=$ $\displaystyle-g_{i}.$ (3)

where $V_{i}$ is the membrane potential, $w_{i}$ is the adaptation variable,
$g_{i}$ is the synaptic conductance of the neuron $i$, $I$ is the injected
current, $C$ is the membrane capacitance, $g_{\rm L}$ is the leak conductance,
$E_{\rm L}$ is the resting potential, $\Delta_{\rm T}$ is the slope factor,
$V_{\rm T}$ is the threshold potential, $V_{\rm rev}$ is the potential
reversal, $\tau_{w}$ is the adaptation time constant, $a$ is the level of
subthreshold adaptation and $\tau_{g}$ is the synaptic time constant. The
parameter $\Delta_{\rm T}$ controls the sharpness of the initial phase of the
spike [45]. This model improves the standard leaky integrate-and-fire by
adding the adaptation current and the exponential term in the membrane
potential. The adaptation current describes the slow activation and
deactivation of some potassium ionic channels [46]. Besides that, a spike
threshold mechanism represented by the exponential term describes the fast
arising when a action potential is generated [43, 42, 47]. For $\Delta_{\rm
T}\rightarrow 0$, the neuron model becomes a standard leaky integrate-and-fire
neuron model [48]. If $V_{i}$ reaches the threshold $V_{\rm{peak}}$ at a
certain time $t$, the following reset conditions are applied:
$V_{i}\rightarrow V_{\rm r}$, $w_{i}\rightarrow w_{\rm r}=w+b$ and
$g_{i}\rightarrow g_{i}+g_{\rm r}$. In our simulations, we consider the
parameter values in Table 1. We constructed a two-dimensional excitatory
neuronal network with 1000 $\mu$m X 1000 $\mu$m dimensions and $N=17,324$
neurons in which the distance between each neuron is 7 $\mu$m in $x$-axis and
8 $\mu$m in $y$-axis [41]. The neuron $i$ is connected with its neighbors that
are within a radius $R$ ($M_{ij}=1$), as displayed in Figure 1(a).

Table 1: Description and values of the parameters in the AEIF system used in the simulations [51]. Parameter | Description | Value
---|---|---
$C$ | Membrane capacitance | 200 pF
$I$ | Constant input current | 500 pA
$N$ | Number of aEIF neurons | 17.324 neurons
$a$ | Subthreshold adaptation | 2.0 nS
$b$ | Triggered adaptation | 70 pA
$\rho$ | Superficial density of neurons | 0.018 neurons/cm2
$R$ | Radius connection | [10,80] $\mu$m
$M_{ij}$ | Adjacency matrix elements | 0 or 1
$w_{i}(0)$ | Initial membrane potential | [0, 70] pA
$V_{i}(0)$ | Initial membrane potential | [-70, -45] mV
$V_{\rm T}$ | Potential threshold | -50 mV
$V_{\rm rev}$ | Reversal potential | 0 mV
$V_{\rm peak}$ | maximum potential | -40 mV
$V_{\rm r}$ | Reset potential | -58 mV
$E_{\rm L}$ | Leak reversal potential | -70 mV
$g_{\rm L}$ | Leak conductance | 12 nS
$g_{\rm syn}$ | Excitatory synaptic conductance | [0.0,0.25] nS
$t_{\rm ini}$ | Initial time in the analyses | 25 s
$t_{\rm fin}$ | Final time in the analyses | 30 s
$\tau_{w}$ | Adaptation time constant | 300 ms
$\tau_{\rm s}$ | Synaptic time constant | 2.728 ms
$\Delta_{\rm T}$ | Slope factor | 2 mV
$d_{\rm x}$ | x-axis distance | 7 $\mu$m
$d_{\rm y}$ | y-axis distance | 8 $\mu$m

##  3 Diagnostic tools

To identify spiking and bursting firing patterns, we utilize the mean
coefficient of variation (CV) of the neuronal inter-spike interval (ISI) [49,
50], that is given by

$\displaystyle{\rm CV}=\frac{{\sigma}_{\rm{ISI}}}{\rm{\overline{ISI}}},$ (4)

where ${\sigma}_{\rm{ISI}}$ is the standard deviation of the ISI normalised by
the mean $\overline{\rm ISI}$ [52]. The ISI is defined by ${\rm
ISI}_{m}=t_{\rm m+1}-t_{\rm m}$, where $\rm t_{m}$ is the time of the $m$-th
neuronal spike [53]. Spike pattern produces ${\rm CV}<0.5$, meanwhile, burst
pattern produces ${\rm CV}\geq 0.5$ [50]. The mean firing rate is calculated
by means of $F=\overline{\rm ISI}^{-1}$ [50].

To determine synchronous behaviour of the whole network, we consider the
global phase order parameter, that is defined as [54]

$Z(t)=\left|\frac{1}{N}\sum_{j=1}^{N}\exp({\rm i}\phi_{j}(t))\right|,$ (5)

where $\phi_{j}$ is the $j$-neuron phase, given by

$\phi_{j}(t)=2\pi m+2\pi\frac{t-t_{j}^{m}}{t_{j}^{m+1}-t_{j}^{m}},$ (6)

$t_{j}^{m}$ corresponds to the time when a $m$-th spike ($m=0,1,2,\dots$) of a
neuron $j$ happens ($t_{j}^{m}<t<t_{j}^{m+1}$). When the network is totally
synchronized $Z(t)=1$, total desynchronized states produce $Z(t)=0$ and when
$0<Z(t)<1$ the network is partially synchronized [54, 53]. To better
understanding of the network dynamics, we calculate the time-average of the
order-parameter, that is given by

$Z_{\rm g}=\frac{1}{t_{\rm fin}-{t_{\rm ini}}}\int_{t_{\rm ini}}^{t_{\rm
fin}}Z(t)dt,$ (7)

where $t_{\rm fin}-t_{\rm ini}$ is the time window. We use a time window equal
to 5 s to compute ${\overline{Z}}$ to avoid transient effects.

We consider the local phase order parameter to study the local synchronization
of each neuron. We separate the 1000 X 1000 $\mu$m network in 25 boxes in each
dimension, totalizing 625 boxes, then each box contains a significant amount
of neurons. The boxes are labeled by the subindex $(i,j)$, where $i$ indicates
the $x$-axis and $j$ the $y$-axis position of the box, $i,j=[1,25]$. The size
of each box is $40\times 40$ $\mu$m and each one of them contains
approximately 28 neurons. The local phase order parameter is given by [55]

$Z_{i,j}(t)=\left|\left[\frac{A^{2}}{N^{2}}\right]\sum_{m=X_{i}}^{X_{i}+40}\sum_{n=Y_{j}}^{Y_{j}+40}\exp({\rm
k}\phi_{m,n})\right|,$ (8)

where $k=\sqrt{-1}$, $A$ is the total number of boxes, $X_{i}$ is the $x$-axis
position of the box, and $Y_{i}$ is the $y$-axis position of the box. The mean
value of $Z_{i,j}$ is

$\overline{Z}_{i,j}=\frac{1}{t_{\rm fin}-t_{\rm ini}}\int_{t_{\rm
ini}}^{t_{\rm fin}}Z_{i,j}(t)dt,$ (9)

the average of the local phase order parameter of the network is

$Z_{\rm L}=\frac{1}{625}\sum_{i,j}^{25}\overline{Z}_{i,j}.$ (10)

Figure 1(a) displays an illustrative scheme of the connections in our neuronal
network, the neurons connect by means of chemical synapses. There are
connections among neurons within the radius R, otherwise there is no
connection. The radius origin is in the neuron position. Figure 1 (b-d)
exhibit the time evolution of the membrane potential of three different
neurons. In panel (b), for a low coupling strength ${\it g_{\rm syn}}=0.001$
nS, the neurons are not synchronized. Increasing the coupling strength to
${\it g_{\rm syn}}=0.1$ nS, the neurons fire in a synchronized spike pattern,
as shown in panel (c). The panel (d) shows neurons synchronized and firing
according to a burst pattern for ${\it g_{\rm syn}}=0.2$ nS. The values of the
local and global phase order parameters are exhibited in panel (e), where the
circles, squares and triangles correspond to the firing patterns of the panels
(b-d). The neurons in the network synchronize (locally and globally) for small
values of $\it g_{\rm syn}$ and the global order parameter varies much more
than the local order parameter. The panel (f) displays the CV as a function of
$\it g_{\rm syn}$. The CV value below the dashed line indicates the spike fire
patterns. There is a maximum value of coupling strength before the transition
from spike to burst behaviors, $g_{\rm syn}=0.2$ nS.

Figure 1: The panel (a) displays a schematic representation of the two-
dimensional network composed of coupled excitatory neurons. The panels (b-d)
exhibit the time evolution of the membrane potential $V$ for three different
values of $\it g_{\rm syn}$. We observe (b) desynchronyzed spikes for ${\it
g_{\rm syn}}=0.001$ nS (blue circle), (c) synchronized spike for ${\it g_{\rm
syn}}=0.1$ nS (blue square), and (d) synchronized burst for ${\it g_{\rm
syn}}=0.2$ nS (blue triangle). The panel (e) shows the global and local
Kuramoto order parameter as a function of $\it g_{\rm syn}$ and the panel (f)
exhibits CV as a function of $g_{\rm syn}$. In panels (b-f), we consider $\it
R$ = 70 $\mu$m.

##  4 Results

A bidirectional neuronal network is able to exhibit spiral wave activities
[14] with different amounts of phase singularities (PS) and spiral waves. PS
is the spiral core and is composed of asynchronous neurons. Each spiral wave
rotates around its own PS and can share the same core. We observe that the
local order parameter is below 0.7 in the PS. We define a $\overline{Z}_{i,j}$
threshold to characterize the phase singularities. In this work, we focus on
studying the correlation of spiral wave emergence with the network parameter
$R$ and $g_{\rm syn}$. We propose a reliable method to distinguish spiral wave
from other travelling waves in the cortex. We use the diagnostic tools shown
in the previous section.

Figure 2 displays three patterns of waves observable in the network and the
local phase order parameter. The panels (a), (b) and (c) exhibit the phase
(color bar) of each neuron in a time instant and the panels (d), (e) and (f)
show the local order parameter of each wave pattern. The ring waves are shown
in Figure 2(a). The wave pattern produces some dark orange dots in the panel
(d), which are the origin of where the ring waves propagate from. The dark
orange spots are not PS due to the fact that the local order parameter is not
below 0.7. The spiral wave is displayed in the panel (b) and its centre is
shown in the panel (e), by the black region. In the panel (e), we calculate
the local order parameter of the spiral wave, where the black area corresponds
to the phase singularity (PS). The mean local order parameter ($Z_{\rm L}$)
does not decrease significantly with the presence of the phase singularity,
due to the high quantity of local synchronous neurons. For higher $R$ values,
it is possible to observe a burst spiral wave, shown in Figure 2(c). There are
some border abnormalities in this case. In Figure 2(f), the border has low
values of $Z_{\rm L}$, however, we do not count it as PS due the fact that
there is no spiral wave rotating around it. To avoid mistakes in the PS
counting, we discard the borders. High values of $R$ increase the border
effects due to the high connectivity of each neuron, although $R$ influences
the firing pattern.

Figure 2: (a) Spatial pattern of a synchronyzed spikes (${\it R}=55$ $\mu$m
and ${\it Z_{g}}=0.94$), (b) spiral wave pattern (${\it R}=64.5$ $\mu$m and
${\it Z_{\rm g}}=0.15$), (c) burst spiral wave pattern (${\it R}=75$ $\mu$m
and ${\it Z_{\rm g}}=0.36$), (d) local order parameter for (e) local order
parameter of spike spiral wave and (f) local order parameter of burst spiral
wave of each neuron $i$. We consider ${\it g_{\rm syn}}=0.14$ nS.

Figure 3 exhibits the parameter space of the four diagnostic tools for $R=75$
$\mu$m. The panel (a) shows the parameter space of the global phase order
parameter. The dark regions correspond to the asynchronous activity ($Z_{\rm
g}\in[0,0.7]$) and the brighter regions are the synchronized states of the
whole network. The network shows synchronous behavior for high values of $R$.
The value of $\it g_{\rm syn}$ also plays an important role in the
synchronization. Increasing it, we verify transitions from asynchronous to
synchronous and vice-versa. The local phase order parameter is displayed in
the panel (b). The coupling strength is of great importance for the emergence
of synchronous states. The CV suffers changes with the increase of $R$ and
$g_{\rm syn}$. It is identified transitions from spike-to-burst, shown in the
panel (c). The panel (d) exhibits the average firing rate, there is a quick
increase for high values of $R$ and $g_{\rm syn}$. In some regions of high
$g_{\rm syn}$, it occurs the synchronization of bursting, which is strongly
related to pathological cases, such as epilepsy seizures [50]. The panels (a)
and (b) are quite different in some regions, showing regions in which the
network is globally desynchronized meanwhile local neurons are synchronized,
$Z_{\rm g}<1$ and $Z_{\rm L}\approx 1$. This result characterizes waves in the
network, however, it is not enough to distinguish the wave types.

Figure 3: (a-d) Parameter space of the network. The color code represents the
values of (a) the global order parameter, (b) the local order parameter, (c)
the mean coefficient of variation (CV) of the neuronal inter-spike interval
and (d) the mean firing rate.

It is necessary more conditions to characterize spiral waves. The spiral waves
have a particular characteristic, that is the phase singularities (the centre
which the waves rotate around). However, high amounts of PS are related to
asynchronous activity. We establish a maximum value of PS to avoid
asynchronous activity. We propose the following conditions to detect spiral
waves: (i) the mean global order parameter has to point out asynchronous
activity, $Z_{\rm g}\leq 0.7$, (ii) the mean local order parameter has to
indicate a local-synchronized activity, $Z_{\rm L}\approx 1$, (iii) there must
be at least one phase singularity and PS is characterized by a low local phase
order parameter, $\overline{Z}_{i,j}\leq 0.7$, (iv) the number of PS is
limited to a maximum number ($P_{\rm max}$), $0<{\rm PS}\leq P_{\rm max}$. We
consider $P_{\rm max}=20$.

Conditions (i) and (ii) detect wave activity in the network, however, they do
not describe the type of wave. The main difference among the wave types is
their local phase order parameter pattern. For example, in the ring wave
pattern, the local order parameter does not exhibit PS, as shown in Figure 2
(d). The presence of PS is characteristic of spiral waves, then at least one
focus is necessary, condition (iii). Multi-phase singularities can be observed
and the network becomes asynchronous when there are too many ones. We consider
a maximum number of PS, condition (iv). These conditions all together indicate
spiral wave in the network. The color scheme in Figs. 4(a) and 4(b) show the
transition of patterns for the region of $R=75$ $\mu$m. The yellow points
denote the spiral waves, the gray points correspond to the non-spiral waves,
the blue points indicate the asynchronous activity and the red points show the
synchronized states ($(Z_{\rm g},Z_{\rm L})\approx 1$).

Figure 4: In (a) and (b) we plot, respectively, the global phase order
parameter and the number of phase singularity as function of $g_{\rm syn}$ for
$R=75\mu m$. The colors in (a) and (b) represent the dynamical pattern of the
network.

The pattern parameter space of $R\times g_{\rm syn}$ is shown in Figure 5(a),
where the colors yellow, gray, blue and red correspond to the spiral waves,
non-spiral waves, asynchronous and synchronous patterns, respectively. For
$R\leq 20$ $\mu$m, it only generates asynchronous activity. Hence, the radius
of connection has great influence in the type of pattern. There are some non-
spiral waves regions in the range $40\leq R\leq 35$ $\mu$m and $0.25\leq
g_{\rm syn}\leq 0.05$ nS, indicating a coexistence of synchronous and non-
spiral waves pattern. It is important to highlight that the radius of
connection is about 70 $\mu$m or greater in the human neocortex [41]. The
emergence of burst spiral is intrinsically correlated with long range
connections ($R\geq 60$ $\mu$m) and intense coupling strength ($g_{\rm
syn}\geq 0.12$ nS). The greenish and pinkish hexagons are the starting and
ending point of the hysteresis analysis, respectively. Hysteresis diagrams of
$Z_{\rm g}$ and CV are displayed in Figs. 5(b) and 5(c), respectively.
Asynchronous and synchronous initial conditions are considered in the purple
and green squares. The CV bistable region coincides with one of the regions in
the $Z_{\rm g}$, both are observed in the $g_{\rm syn}$ range [0.11, 0.13] nS.
Two different outcomes are possible, that are burst synchronization for
asynchronous initial conditions and spikes desynchronization for synchronous
initial conditions. There is one more bistable region in the panel (b), in
which the spike desynchronization is easier achieved for synchronous initial
conditions.

Figure 5: (a) Parameter space of the network activity pattern, where each
color represents one type of pattern. The yellow region corresponds to the
spiral wave, the gray area is related to the non-spiral wave (for example ring
pattern), the blue region denotes the desynchronized pattern and the red area
corresponds to the synchronous pattern. The panels (b) and (c) show the
hysteresis diagram for $R=75$ $\mu$m, (b) global order parameter and (c)
coefficient of variation. Asynchronous initial conditions are considered in
the purple squares, meanwhile, synchronous initial conditions are used in the
green squares. The greenish and pinkish hexagons exhibit the points which
start and end the hysteresis cycle, respectively.

##  5 Conclusions

In this work, we investigate the emergence of spiral wave patterns within a
bidirectional network composed of pyramidal neuron cells. The individual
dynamics of each neuron is described by the adaptive exponential integrate-
and-fire model. The roles of spiral waves in cognitive activities is not fully
understood. In this way, comprehending how network parameters influence the
emergence of spirals and developing the capability to detect them is crucial
for advancing our understanding of their role in cognitive activities. Hence,
we investigate how network parameters influence the emergence of spirals.
Additionally, we introduce a novel method capable of identifying chimera
spiral waves.

We show that the network composed of coupled pyramidal neurons generates
different spatial-temporal patterns. The local and global phase order
parameters together have shown to be great tools of wave detection. We observe
that wave activity generate a specific range of values for both mean order
parameters. Locally the neurons are synchronized for almost all values of $R$
and $g_{\rm syn}$, indicating that waves do not disrupt the local dynamics.
The waves greatly affect the global dynamics and makes the global order
parameter goes to values close to 0, characterizing desynchronization in the
network. The low impact that waves have in local synchronization is due to the
firing of the closest neighbors when the travelling wave arrives at the
neuron. In the entire network, the passage of wave helps the local neurons to
fire meanwhile the neurons that are not at the wave location do not receive a
firing support. These characteristics hold great power to identify waves in
the network, thus we suggest a diagnostic capable of identifying waves. They
are identified by a low global and high local order parameters. Applying this
criteria, we identify some waves. Mostly spike waves are found for radius
lower than 55 $\mu$m. The type of wave is not pointed out by this diagnostic,
then more criteria to identify spiral waves are necessary.

The identification of spiral waves can not be done by using only the condition
proposed earlier, due to the fact that it does not specify the wave type. The
spiral waves generate a specific abnormality in the phase producing very low
local order parameter in it. The neurons localized in this abnormality are
desynchronized and are the core of the spiral, called phase singularity (PS).
The neurons in the PS exhibit a low value in the local order parameter
inferior to 0.7, then we propose that PS are characterized by
$\rm\overline{Z}_{i,j}<0.7$. It is not necessary only one PS to exist spiral
waves, multi-phase singularies are possible. However, a great amount of PSs
promote desynchronized activity. We define a maximum number of PSs in the
network, $P_{\rm max}=20$. We propose four conditions which must be met to
identify spiral waves, these conditions regard the range of both phase order
parameters and the amount of PS. The synchronized states are defined by high
global and local order parameters, non-spiral waves are characterized by
$(Z_{\rm L}\approx 1$, $Z_{\rm g}<0.7$ and PS = 0, spiral waves exhibit the
same conditions of non-spiral waves with the addition of PS $>0$ and
desynchronized states are characterized by $Z_{\rm L}\approx 0$ and $Z_{\rm
g}\approx 0$.

The spiral waves are generated by different combinations of $R$ and $g_{\rm
syn}$. This wave pattern is not only restricted for spike firing, but also for
burst spiral. With regard to the burst spiral pattern, it is not as easily
observed as spike spiral. The burst spiral appearance is restricted for
specific values of $R$ and coupling strength, for example $R=75$ $\mu$m and
$g_{\rm syn}=0.14$ nS. The emergence of non-spiral wave is related to the
initial conditions. To understand the initial conditions influence in the
system, we analyze the hysteresis diagram for a constant radius and vary the
$g_{\rm syn}$. Desynchronous behavior is easier to achieve if the initial
conditions are desynchronized. We find a region of bistability in the range
$g_{\rm syn}=$[0.11, 0.13] nS. Depending on the initial conditions, it is
possible to verify the existence of spike or burst waves. The identification
of bistable state is crucial due to its strong correlation with epileptic
seizures [56].

## Acknowledgments

We acknowledge the support from CAPES (Finance Code 88887.849164/ 2023-00;
88881.846051/2023-01; 88881.143103/2017-01), CNPq (403120/2021-7;
301019/2019-3), and FAPESP (Grant Nos. 2018/03211-6; 2020/04624-2;
2022/04251-7).

## References

## References

  * [1] Li Teng-Chao, Pan De-Bei, Zhou Kuangshi, Jiang Ruhong, Jiang Chenyang, Zheng Bo and Zhang Hong 2018 Jacobian-determinant method of identifying phase singularity during reentry. Phys. Rev. E, v. 98, 062405.
  * [2] Ding Q, Wu Y, Hu Y, Liu C, Hu X and Jia Y 2023 Tracing the elimination of reentry spiral waves in defibrillation: Temperature effects. Chaos, Solitons and Fractals,v. 174, 113760.
  * [3] Holden A 1997 The restless heart of a spiral. Nature,v. 387, 655–657.
  * [4] Galappaththige S K, Pathmanathan P, Bishop M J and Gray R A 2019 Effect of Heart Structure on Ventricular Fibrillation in the Rabbit: A Simulation Study. Front. Physiol., 10:564.
  * [5] Davidenko J, Pertsov A, Salomonsz R, Baxter W and Jalife J 1992 Stationary and drifting spiral waves of excitation in isolated cardiac muscle. Nature,355, 349–351.
  * [6] Bingen B O, Engels M C, Schalij M J, Jangsangthong W, Neshati Z, Feola I, Ypey D L, Askar S F A, Panfilov A V, Pinappels D A and Vries A A F 2014 Light-induced termination of spiral wave arrhythmias by optogenetic engineering of atrial cardiomyocytes. Cardiovascular Research,v. 104, 194-205.
  * [7] Kundu S, Majhi S, Muruganandam P and Ghosh D 2018 Diffusion induced spiral wave chimeras in ecological system. Eur. Phys. J. Spec. Top.,227, 983–993.
  * [8] Totz J F, Rode J, Tinsley M R, Showalter K and Engel H 2018 Spiral wave chimera states in large populations of coupled chemical oscillators. Nature Phys, 282–285.
  * [9] Glaze T A, Lewis S and Bahar S. 2016 Chimera states in a Hodgkin-Huxley model of thermally sensitive neurons. Chaos,26, 083119.
  * [10] Haugland S W 2021 The changing notion of chimera states, a critical review. Journal of Physics: Complexity,2, 032001.
  * [11] Majhi S, Bera B K, Ghosh D and Perc M 2019 Chimera states in neuronal networks: a review. Physics of Life Reviews,28, 100-121.
  * [12] Wang Z and Li1u Z 2020 A brief review of chimera state in empirical brain networks. Frontiers in Physiology,11, 724.
  * [13] Provata A and Venetis I 2020 Chimera states in leaky integrate-and-fire dynamics with power law coupling. The European Physical Journal B, 93, 160.
  * [14] Santos M S, Protachevicz P R, Caldas I L, Iarosz K C, Viana R L, Szezech Jr J D, de Souza S L T and Batista A M 2021 Spiral wave chimera states in regular and fractal neuronal networks. Journal of Physics: Complexity,2, 015006.
  * [15] Santos M S, Szezech Jr. J D, Borges F S, Iarosz K C, Caldas I L, Batista A M, Viana R L and Kurths J 2017 Chimera-like states in a neuronal network model of the cat brain. Chaos, Solitons and Fractals,101, 86-91.
  * [16] Ramados J, Aghababaei S, Parastesh F, Rajagopal K, Jafari S and Hussain I 2021 Chimera state in the network of fractional-order FitzHugh Nagumo neurons. Complexity, 2437737.
  * [17] Zhang Z and Dai L 2022 Effects of synaptic pruning on phase synchronization in chimera states of neural network. Applied Sciences,12, 1924.
  * [18] Wu X and Ma J 2013 The formation mechanism of defects, spiral wave in the network of neurons. Plos One, 8(1), e55403.
  * [19] Burk A, Strelkova G and Anishchenko V 2018 Spiral wave patterns in a two-dimensional lattice of nonlocally coupled maps modeling neural activity. Chaos, Solitions and Fractals, 120, 75-82.
  * [20] Santos M S, Protachevicz P R, Iarosz K C, Caldas I L, Viana R L, Borges F S, Ren H-P, Szezech Jr. J D, Batista A M and Grebogi C 2019 Spike-burst chimera states in an adaptive exponential integrate-and-fire neuronal network. Chaos, 29, 043106.
  * [21] Santos V, Borges F S, Iarosz K C, Caldas I L, Szezech J D, Viana R L, Baptista M S and Batista A M 2020 Basin of attraction for chimera states in a network of Rossler oscillators. Chaos, 30, 083115.
  * [22] Hizanidis J, Kanas V, Bezerianos A and Bountis T 2014 Chimera states in network of nonlocally coupled Hindmarsh-Rose neuron models. International Journal of Bifurcation and Chaos, 24(3), 1450030.
  * [23] Tsigkri-DeSmedt N D, Hizanidis J, Hövel P and Provata A 2015 Multi-chimera states in the leaky integrate-and-fire model. Procedia Computer Science, 66, 13-22.
  * [24] Majhi S, Perc M and Ghosh D 2016 Chimera states in uncoupled neurons induced bu a multilayer structure. Scientific Reports, 6, 39033.
  * [25] Hövel P, Vülling A, Omelchenko I and Hizanidis J 2016 Chimera states in neuronal systems of excitability tyoe-I. In: Battiston, S., De Pellegrini, F., Caldarelli, G., Merelli, E. (eds) Proceedings of ECCS 2014. Springer Proceedings in Complexity. Springer, Cham.
  * [26] Calim A, Torres J J, Ozer Mand Muhammet U 2020 Chimera states in hybrid coupled neuron populations, Neural Network, 126, 108-117.
  * [27] Hussian I, Jafari S, Perc M and Ghosh D 2021 Chimera states in a multi-weighted neuronal network. Physics Letters A, 424, 12747.
  * [28] Lucchetti A, Jensen M H and Heltberg M L 2021 Emergence of chimera states in a neuronal model of delayed oscillators. Physical Review Research, 3, 033041.
  * [29] Nikitin D, Omelchenko I, Zakharova A, Avetyan M, Fradkov A L and Schöll E 2022 Complex partial synchronization patterns in network of delay-coupled neurons. Philosophical Transactions of the Royal Society A, 377, 20180128.
  * [30] Huang X, Xu W, Liang J, Takagaki K, Gao X and Wul J 2010 Spiral Wave Dynamics in Neocortex. Neuron, 68, 978-990.
  * [31] Huang X, Troy W C, Yang Q, Ma H, Laing C R, Schiff S J and Wu J-Y 2004 Spiral waves in disinhibited mammalian neocortex. The Journal of Neuroscience, 24(44), 9897-9902.
  * [32] Xu Y, Long Y, Feng J and Gong P 2023 Interacting spiral wave patterns underlie complex brain dynamics and are related to cognitive processing. Nature human behaviour, 7, 1196–1215.
  * [33] Ma J, Huang L, Tang J, Ying H and Jin W 2012 Spiral wave death, breakup induced by ion channel poisoning on regular Hodgkin–Huxley neuronal networks. Communications in Nonlinear Science and Numerical Simulation, v. 17, 4281-4293.
  * [34] Wang Q, Perc M, Duan Z and Chen G 2008 Delay-enhanced coherence of spiral waves in noisy Hodgkin–Huxley neuronal networks. Physics Letters A, v. 372, 5681-5687.
  * [35] Bukh A, Strelkova G and Anishchenko V 2019 Spiral wave patterns in a two-dimensional lattice of nonlocally coupled maps modeling neural activity. Chaos, Solitons and Fractals, v. 120, 75-82.
  * [36] Santos M S, Protachevicz P R, Caldas I L, Iarosz K C, Viana R L, Szezech J D, Souza S L T and Batista A M 2021 Spiral wave chimera states in regular and fractal neuronal networks. Journal of Physics: Complexity, v. 2, 015006.
  * [37] Uzelac I, Iravanian S, Bhatia N and Fenton F H 2022 Spiral wave breakup: Optical mapping in an explanted human heart shows the transition from ventricular achycardia to ventricular fibrillation and self-termination. Heart Rhythm, 19, 1914-1915.
  * [38] Prechtl J C P, Cohens L B, Pesaran B, Mitra P P and Leinfeld D K 1997 Visual stimuli induce waves of electrical activity in turtle cortex. Proceedings of the National Academy of Sciences, v. 94, 7621-7626.
  * [39] Townsend R G, Solomon S S, Chen S C, Pietersen A N J, Martin P R, Solomon S G and Gong P 2015 Emergence of Complex Wave Patterns in Primate Cerebral Cortex. The Journal of Neuroscience, v. 35(11), 4657– 4662.
  * [40] Bekkers J M 2011 Pyramidal neurons.Current Biology Quick Guide, 21, PR975.
  * [41] Bezaire M J, Raikov I, Burk K, Vyas D and Soltesz I 2016 Interneuronal mechnisms of hippocampal theta oscillations in a full-scale model of the rodent CA1 circuit. eLife, 5e18566.
  * [42] Brette R and Gerstner W 2005 Adaptive exponential integrate-and-fire model as an effective description of neuronal activity. Journal of Neurophysiology, 94, 3637-3642.
  * [43] Gerstner W, Kistler W M, Naud R and Paninski L 2014 Neuronal Dynamics: From Single Neurons to Networks and models of Cognition. 1st ed. New York: Cambridge.
  * [44] Hertäg L, Hass J, Golovko T and Durstewitz D 2012 An approximation to the adaptive exponential integrate-and-fire neuron model allows fast and predictive fitting to physiological data. Frontiers in Computational Neuroscience, 6, 62, 1.
  * [45] Badel L, Lefort S, Brette R, Petersen C C H, Gerstner W and Richardson M J E 2008 Dynamic IV curves are reliable predictors of naturalistic pyramidal-neuron voltage traces. Am. Physiological Soc., 99, 656-666.
  * [46] Brette R 2015 What is the most realistic single- compartment model of spike initiation? Plos computational Biology, v. 11.
  * [47] Naud R, Marcille N, Clopath C and Gerstner W 2008 Firing patterns in the adaptive exponential integrate-and-fire model. Biological Cybernetics, 99, 335–347.
  * [48] Clopath C, Jolivet R, Rauch A, Lüscher H R and Gerstner W 2007 Predicting neuronal activity with simple models of the threshold type: Adaptive Exponential Integrate-and-Fire model with two compartments. Neurocomputing, 70 (10-12), 1668-1673.
  * [49] Borges F S, Protachevicz P R, Lameu E L, Bonetti R C, Iarosz K C, Caldas I L, Baptista M S, Batista A M 2017 Synchronised firing patterns in a random network of adaptive exponential integrate-and-fire neuron model. Neural Networks, 90, 1-7.
  * [50] Protachevicz P R, Borges F S, Lameu E L, Ji P, Iarosz K C, Kihara A H, Caldas I L, Szezech Jr. J D, Baptista M S, Macau E E N, Antonopoulos C G, Batista A M and Kurths J 2019 Bistable firing pattern in a neural network model. Frontiers in Computational Neuroscience, 13, 1-8.
  * [51] Naud R, Marcille N, Clopath C and Gerstner W 2008 Firing patterns in the adaptive exponential integrate-and-fire model. Biological Cybernetics, 99, 335-347.
  * [52] Gabbiani F and Koch C 1998 Principles of spike train analysis, in Methods in Neuronal Modeling: From Ions to Networks, 2nd Edn, eds C. Koch and I. Segev (Cambridge, MA: MIT Press), 313-360.
  * [53] Protachevicz P R, Borges R R, Reis A S, Borges F S, Iarosz K C, Caldas I L, Lameu E L, Macau E E N, Viana R L, Sokolov I M, Ferrari F A S, Kurths J, Batista A M, Lo C-Y, He Y, Lin C-P 2018 Synchronous behaviour in network model based on human cortico-cortical connections. Physiological Measurement, 39, 074006.
  * [54] Kuramoto Y 1984 Chemical Oscillations, Waves, and Turbulence. Berlin: Springer-Verlag.
  * [55] Santos M S, Protachevicz P R, Caldas I L, Iarosz K C, Viana R L, Szezech Jr. J D, Souza S L T and Batista A M 2020 Spiral wave chimera states in regular and fractal neuronal networks. Journal of Physics: Complexity, 2, 015006.
  * [56] Borges F S, Protachevicz P R, Souza D L M, Bittencourt C F, Gabrick E C, Bentivoglio L E, Szezech J D Jr., Batista A M, Caldas I L, Dura-Bernal S and Pena R F O 2023 The Roles of Potassium and Calcium Currents in the Bistable Firing Transition. Brain Sci., 13, 1347.