Emergence of Neuronal Synchronisation in Coupled Areas

One of the most fundamental questions in the field of neuroscience is the emergence of synchronous behaviour in the brain, such as phase, anti-phase, and shift-phase synchronisation. In this work, we investigate how the connectivity between brain areas can influence the phase angle and the neuronal synchronisation. To do this, we consider brain areas connected by means of excitatory and inhibitory synapses, in which the neuron dynamics is given by the adaptive exponential integrate-and-fire model. Our simulations suggest that excitatory and inhibitory connections from one area to another play a crucial role in the emergence of these types of synchronisation. Thus, in the case of unidirectional interaction, we observe that the phase angles of the neurons in the receiver area depend on the excitatory and inhibitory synapses which arrive from the sender area. Moreover, when the neurons in the sender area are synchronised, the phase angle variability of the receiver area can be reduced for some conductance values between the areas. For bidirectional interactions, we find that phase and anti-phase synchronisation can emerge due to excitatory and inhibitory connections. We also verify, for a strong inhibitory-to-excitatory interaction, the existence of silent neuronal activities, namely a large number of excitatory neurons that remain in silence for a long time.

Recently, the synchronisation in neuronal networks in presence of both excitatory and inhibitory synapses has been observed using neuronal models (Bera et al., 2019a;Pal et al., 2021) coupled through hypernetworks (Rakshit et al., 2018a;Bera et al., 2019b) and multiplex configurations (Rakshit et al., 2018b). In cortico-cortical communication, one cortical area can interact with other one by means of excitatory and inhibitory connectivities (Roland et al., 2014;Tamioka et al., 2015;Tovete et al., 2015;D'Souza et al., 2016). In this work, we investigate how the excitatory and inhibitory connectivities from one area to another influence the phase angle and neuronal synchronisation. We consider unidirectional (sender-receiver) and bidirectional interactions between two areas. For the unidirectional interaction and desynchronised neurons in the sender area, we show that the phase angle values and synchronous behaviour of the neurons in the receiver area depend not only on the neuronal dynamics of the sender area, but also on the type of connections between the areas. We find phase, anti-phase, and shift-phase synchronisation in the receiver area when the neurons are synchronised in the sender one. With regard to the bidirectional interaction, we verify phase and anti-phase synchronous behaviour between the areas. The excitatory-toexcitatory (inhibitory-to-inhibitory) and excitatory-to-inhibitory (inhibitory-to-excitatory) connections can induce phase and anti-phase synchronisation between the areas, respectively. For a strong inhibitory-to-excitatory interaction between the areas, a large number of silent excitatory neurons are found in both areas.
The paper is organised as follows. In section 2, we introduce the neuronal network composed of adaptive exponential integrate-and-fire (AEIF) neurons and the diagnostic tools to characterise the synchronous behaviour. Sections 3 and 4 present our results and discussions about the effects of the unidirectional and bidirectional interactions between two areas, respectively, on the neuronal synchronisation and the phase angle. We draw our conclusions in section 5. FIGURE 1 | Connection matrix of two areas, where excitatory and inhibitory neurons are randomly connected. In each area, p ee and p ii are the probabilities of connection between excitatory and inhibitory neurons, respectively. The probability p ei (p ie ) corresponds to the connection from the excitatory (inhibitory) to inhibitory (excitatory) neurons. The superscript "A" is used for the connection probabilities between neurons in different areas.

Network
We build a neuronal network composed of two areas, where each one has a thousand of adaptive exponential integrate-andfire neurons (N = 1, 000) (Brette and Gerstner, 2005). Each area has a fraction of excitatory (P exc = 0.8) and inhibitory (P inh = 0.2) neurons (Noback et al., 2005;di Volo et al., 2019). In each area, the neurons are randomly coupled by means of excitatory and inhibitory connections. The connection is excitatory (inhibitory) when it occurs from an excitatory (inhibitory) neuron. Inside of each area, the probabilities of connections in the same neuronal populations (excitatory or inhibitory) are given by p ee = 0.05 and p ii = 0.2, while between different neuronal populations by p ei = p ie = 0.05 (di Volo et al., 2019). Between the areas, the probabilities are given by p A ee = 0.01 (from excitatory to excitatory neurons), p A ei = 0.05 (from excitatory to inhibitory neurons), p A ii = 0.10 (from inhibitory to inhibitory neurons), and p A ie = 0.05 (from inhibitory to excitatory neurons). Figure 1 shows how the probabilities are distributed in a connection matrix, where k and j correspond to the presynaptic and postsynaptic neurons, respectively. In Figure 1, j, k ∈ [1, 1000] correspond to the neurons in Area 1 and j, k ∈ [1001, 2000] to the neurons in Area 2. For the unidirectional configuration, the connections given by k = [1001,2000] and j = [1, 1000] are not considered. For both unidirectional and bidirectional configuration, we consider only excitatory or inhibitory connections between the areas in each case.
With regard to the coupling intensities, each one is associated with a probability of connection and denoted by g ee , g ii , g ei , g ie , g A ee , g A ei , g A ii , and g A ie . For instance, g ei and g ie are related to the connections with the probabilities p ei and p ie , respectively.

Neuronal Model
The cortex is mainly constituted by excitatory pyramidal neurons and inhibitory interneurons (Atencio and Schreiner, 2008). Excitatory neurons have a relatively lower firing rate than inhibitory ones (Wilson et al., 1994;Inawashiro et al., 1999;Baeg et al., 2001). In the mammalian cortex, excitatory neurons show regular spike (RS), while inhibitory neurons exhibit fast spike (FS) activities (Neske et al., 2015;Wang et al., 2016). In addition, while inhibitory neurons exhibit a negligible adaptation, excitatory neurons show an adaptation mechanism in their firings (Foehring et al., 1991;Mancilla et al., 1998;Hensch and Fagiolini, 2004;Destexhe, 2009;Masia et al., 2018;Borges et al., 2020). The adaptive exponential integrate-and-fire (AEIF) model is able to mimic these different firing patterns, including RS and FS (di Volo et al., 2019). In this work, the dynamics of each neuron j (j = 1, . . . , N) in the network is given by The membrane potential V j and adaptation current w j represent the state of each neuron j. The capacitance membrane is set to C m = 200 pF, the leak conductance to g L = 12 nS, the leak reversal potential to E L = −70 mV, the slope factor to T = 2 mV, and the spike threshold to V T = −50 mV. We consider the injection of current I = 270 pA, which is the intensity above the rheobase current. The application of this constant current allows that the neurons change their potentials from resting potentials to spikes. The level of the subthreshold and triggered adaptation are represented by a j and b j , respectively. We consider inhibitory neurons of fast spiking activities without adaptation (a j = 0 and b j = 0) and excitatory neurons of regular spiking with adaptation mechanisms (a j = [1.9, 2.1] nS and b j = 70 pA). Neuronal adaptation corresponds to the capacity of the neuronal membrane in adapting to its excitability according to the past neuronal activity. A sub-and a triggered-threshold adaptation mechanism can be associate with the parameters a j and b j , respectively. The adaptation current also depends on the adaptation time constant τ w = 300 ms. g j represents the synaptic conductance of each neuron j with an exponential decay associated with the synaptic time constant τ s = 2.728 ms. The connections from excitatory and inhibitory neurons are related to the excitatory and inhibitory matrix, − → M exc and − → M inh , where each matrix element is identified as M exc jk and M inh jk , respectively. A matrix element is equal to 1 when there is a connection from k to j or 0 in the absence of a connection. The excitatory and inhibitory elements of the matrix are associated with the red and blue dots in Figure 1. The chemical current input I chem j arriving on each neuron j is defined by the expression where the excitatory and inhibitory currents are given by and where I xy j and I xy,A j are associated with the excitatory (x=e) or inhibitory connectivity (x=i) arriving at the excitatory (y=e) or inhibitory neurons (y=i). The type of synapse depends on the synaptic reversal potential V REV . We consider the V exc REV = 0 mV for excitatory and V inh REV = −80 mV for inhibitory synapses. N T is the total number of neurons in the network. When the membrane potential of the neuron j is above the threshold V j > V thres (Naud et al., 2008), the state variable is updated by the rule In our simulations, we consider V r = −58 mV. The value of b j depends whether the neuron j is excitatory or inhibitory. Each synaptic current is related to the respective conductance g s . Inside of each area, g s is equal to g ee for synapses between excitatory neurons, g ei for synapses from excitatory to inhibitory neurons, g ii for synapses between inhibitory neurons, and g ie for synapses from inhibitory to excitatory neurons. Between different areas we include the superscript "A." The time delay in the conductance is d exc = 1.5 ms for excitatory connections and d inh = 0.8 ms for inhibitory ones (Borges et al., 2020). Table 1 gives the values of the parameters used in our simulations. We consider g ee = 0.5 nS, g ii = 2 nS, and g ie = 1.5 nS. The areas exhibit synchronous and desynchronous behaviour when uncoupled between them for g ei = 1 nS and g ei = 2 nS, respectively. The area 1 is considered synchronised and desynchronised when uncoupled, while the area 2 is always considered desynchronised when uncoupled.
The initial values of V j are randomly distributed in the interval [−70, −50] mV for all neurons. The initial values of w j are randomly distributed in the interval [0, 300] pA for excitatory neurons and equal to 0 for inhibitory ones. The initial value of g j is equal to 0 for all neurons. To solve the delayed differential equations, we consider that the excitatory and inhibitory neurons  in the network are not spiking before the beginning of the simulation (t = 0). To integrate the set of ordinary differential equations, we use the 4th Runge-Kutta method with the time-step of integration equal to 10 −2 ms.

Synchronisation and Relative Phase Angle
The synchronous behaviour in the network can be identified by means of the complex phase order parameter (Kuramoto, 1984) where R(t) is the amplitude of a centroid phase vector over time.
The phase of each neuron j is obtained through where t j,m corresponds to the time of the m−th spike of the neuron j (t j,m < t < t j,m+1 ) (Rosenblum et al., 1997). We consider that the spike occurs when V j > V thres . The value of R(t) is equal to 0 for completely desynchronised behaviour and equal to 1 for fully synchronised patterns. We calculate the time-average order parameter R (Batista et al., 2017) where t fin − t ini is the time window with t fin = 100 s and t ini = 20 s. The order parameter for each area is given by where A denotes the area number. The mean value of R A (t) (R A ) is computed by Equation (8). The resultant phase angle of each area A is defined as The real R x A and complex R y A components of the order parameter can be described as and respectively. A (t) evolves in the counter-clockwise direction, since each individual neuron evolves in this direction. We define a relative phase angle for each area ′ A (t) (Varela et al., 2001) as The phase of the area 1 changes over time and its relative value ′ 1 (t) is equal to 0. For the area 2, the relative phase angle ′ 2 (t) can change over time. We calculate the mean value of the relative phase angle of the area 2 by means of  To compute the predominant rotation direction of the area 2, we consider the first-order derivative of their relative phase angle, which corresponds to the instantaneous velocity of the relative phase angle (rad/s), that is given bẏ The mean value of the relative velocity of the area 2 is obtained via˙ FIGURE 4 | Unidirectional excitatory interaction between synchronised neurons in the sender area and desynchronised neurons in the receiver area. Synchronised sender area can generate phase and shift-phase synchronisation due to g A ee and g A ei , respectively. (A) The predominant direction of rotation of the relative phase angle of the area 2, (B) the standard deviation of the relative phase angle of the area 2, (C) the mean relative phase angle of the area 2, and (D) the mean order parameter of the area 2. The circle, square, and triangle symbols indicate g A ee = 2.7 nS and g A ei = 0.6 nS, g A ee = 1.5 nS and g A ei = 3 nS, and g A ee = 0.3 nS and g A ei = 5.4 nS, respectively.
The value is close to 0 for non-preponderant direction of rotation and positive (negative) for predominant counter-clockwise (clockwise) direction. The variability of ′ 2 is given by Small and high deviations are given by σ ′ A ≈ 0 and σ ′ A > 0.5, respectively. Figure 2 shows a schematic representation of the mean order parameter of the area 1 and area 2 for Figure 2A desynchronised patterns of both areas out-of-phase, Figure 2B in-phase ( ′ 2 = 0), Figure 2C anti-phase ( ′ 2 = π), and Figure 2D shift-phase synchronisation ( ′ 2 ≈ 4π/3). In Figure 2A, the clockwise (blue) and counter-clockwise (red) arrows indicate that the relative phase angle of area 2 change over time. The amplitude and direction of the relative phase angle of the area A can be described by R A e i ′ A . In our simulations, we observe results in which desynchronised activities can be related to high variability of the relative phase angle of the area 2. When the areas are synchronised, the variability of the relative phase angle can go to 0. The mean value of relative phase angle is efficient for small variability of the relative phase angle. We consider unidirectional interactions with excitatory connections between the areas, where the neurons are synchronised in the sender area and the neurons are initially desynchronised in the receiver area. In the raster plots, the blue and red dots indicate the fires of the inhibitory and excitatory neurons, respectively. Synchronised sender area can actuate on the initially desynchronised receiver area generating (A) phase synchronisation due to g A ee , (B) desynchronisation due to both g A ee and g A ei , and (C) shift-phase synchronisation due to g A ei . In (A), R 1 (black line) and R 2 (green line) exhibit different values due to the fact that there are some neurons with bursting activity in the area 2. The bursts are generated due to the excitatory connections from the area 1 to the excitatory neurons of the area 2. In (C), R 2 is smaller than R 1 due to the excitatory connections from the area 1 to the inhibitory neurons of the area 2.

Excitatory Connections
We analyse a neuronal network separated into two areas (senderreceiver) coupled by means of the excitatory connections. Figure 3 displays a schematic representation of the sender area 1 to receiver area 2 via excitatory connections, that are related to the g A ei and g A ee conductances. Firstly, we consider the case in which the neurons in the sender and receiver area are desynchronised. Through the time evolution of ′ 2 , we verify that, depending on the conductance values, it can occur a relative counterclockwise rotation, clockwise rotation, or neither of them in the area 2. We observe that g A ee contributes to generate positive˙ ′ 2 , while g A ei to negative one. Both areas show order parameters with small values, i.e., the neurons remain desynchronised.
Secondly, we consider that the neurons in the sender area are synchronised while the neurons in the receiver area are initially desynchronised. In Figure 4A, the parameter space g A ee × g A ei exhibits values of˙ ′ 2 in [−0.1, 0.1] rad/s, where the˙ ′ 2 values approximately or less than −0.1 are in a small black diagonal region. For a large set of parameters,˙ ′ 2 is close to zero. Figure 4B displays σ ′ 2 values about 0, except for a small red region where the values are greater than 0.5. We compute the mean relative phase angle of the area 2, as shown in Figure 4C. We verify that g A ee leads the area 2 to a value of the phase angle equal to 0, while g A ei leads to a shifted phase angle, ′ 2 ≈ 4π/3. The stabilisation of the phase angle is associated with the synchronisation of the area 2, as shown in Figure 4D. We see a large region in which R 2 is close to 1, meaning that the neurons in the area 2 are synchronised. Figure 5 displays the raster plot (top), relative phase angle (middle), and order parameter (bottom) for (a) g A ee = 2.7 nS and g A ei = 0.6 nS, (b) g A ee = 1.5 nS and g A ei = 3 nS, and (c) g A ee = 0.3 nS and g A ei = 5.4 nS, that are indicated in Figure 4 through the circle, square, and triangle symbols, respectively. In Figure 5A, we see that ′ 2 is closed to the sender phase angle and the receiver area has synchronised neurons due to a high g A ee conductance. Figure 5B shows that for a combination of g A ee and g A ei , neurons in the area 2 are not synchronised. In Figure 5C, we observe that due to high g A ei conductance, neurons in the area 2 have a shift-phase synchronisation, corresponding to ′ 2 ≈ 4π/3, as indicated in Figure 2D. Figure 6 displays a schematic representation of the sender area to the receiver area via the g A ii and g A ie conductances. For desynchronised neurons in the sender and receiver areas, increasing g A ii , we see that˙ ′ 2 is positive for small values of g A ie , as shown in Figure 7A. On the other hand, g A ie contribute to negative values of˙ ′ 2 . Figure 7B exhibits a high variability (σ ′ 2 > 0.5) of ′ 2 . In Figures 7C,D, we compute the mean order parameters for the areas 1 and 2, respectively. The neurons are desynchronised in the area 1, while in the area 2, we see a small region in the parameter space where there is synchronous FIGURE 7 | Unidirectional inhibitory interaction between two coupled areas with desynchronised neurons. Desynchronised sender area can generate synchronisation and silence on excitatory neurons of the receiver area due to g A ii and g A ie , respectively. (A) The predominant direction of rotation of the relative phase angle of the area 2, (B) the standard deviation of the relative phase angle of the area 2, (C) the mean order parameter of the area 1, and (D) the mean order parameter of the area 2. The circle, square, and triangle symbols indicate g A ii = 3.6 nS and g A ie = 0.3 nS, g A ii = 2 nS and g A ie = 1 nS, and g A ii = 0.4 nS and g A ie = 2.7 nS, respectively.

Inhibitory Connections
behaviour (R 2 > 0.8). The increase of R 2 is due to the g A ii parameter, that induces the inhibition of the inhibitory neurons of the receiver area (circle). For high g A ie values (triangle), a large number of excitatory neurons do not fire. A similar result was reported by Zhou et al. (2010), where silent activities of excitatory neurons were observed due to a strong inhibition. Urban-Ciecko et al. (2015) found which inhibition can silence excitatory synapses in the neocortex. Pals et al. (2020) demonstrated that activity-silence maintenance can be related to a working memory process. The silence of neurons has received great attention in the last years (Mochol et al., 2015;Wiegert et al., 2015;Barbosa et al., 2020;Xu et al., 2020). for the values of g A ii and g A ie indicated by circle, square, and triangle symbols in Figure 7. For g A ii = 3.6 nS and g A ie = 0.3 nS (Figure 8A), the neurons in the area 2 synchronise and ′ 2 denotes a relative counter-clockwise rotation. For g A ii = 2 nS and g A ie = 1 nS (Figure 8B), the neurons in the area 2 are desynchronised and ′ 2 denotes a relative clockwise rotation. In Figure 8C, there is no synchronous behaviour and we see a large number of excitatory neurons that remain in silence for a long time.
We also consider the case in which the sender area has synchronised neurons. Figure 9A shows that the˙ ′ 2 values are positive due to the g A ii conductance with small values of g A ie for small values of g A ie . This result is similar to the situation in which the neurons in the sender area are desynchronised. However, due to the synchronised neurons in the area 1, we verify the existence of regions in the parameter space g A ii × g A ie with small values of the variability, as shown in Figure 9B. σ ′ 2 < 0.5 corresponds to a certain stabilisation of the relative phase angle of the area 2. Figure 9C displays the mean relative phase angle of the area 2. For the lowest σ ′ 2 , we find a region where ′ 2 ≈ π. In Figure 9D, we see that synchronous behaviour in the area 2 for large g A ii and small g A ie values arise, where there is a circle symbol. Partial antiphase synchronisation is observed in the region indicated by the square. The triangle denotes the region in which a high inhibition of the excitatory neurons occurs, and as a consequence a great quantity of excitatory neurons in the receiver area are silenced. conductances. Without an interaction between the areas (g A ee = g A ei = 0), the neurons exhibit desynchronised activities. Figure 11A displays the mean order parameter of the neuronal network. The region, where the circle is located, has a larger value of R, due to the fact that the neurons in the two areas are synchronised, namely phase synchronisation among neurons. In Figure 11B, we verify that the reduction of the variability can be associated with the synchronised activities between the areas. Figures 11C,D shows the mean order parameters of the areas 1 and 2, respectively. We can see that the regions of the small variability of ′ 2 (circle and triangle symbols) correspond to the synchronised activities. For the region with large variability values (square symbol), the neurons of the areas are desynchronised. In the region where the triangle symbol is located, there is an anti-phase synchronisation between the neurons of the areas 1 and 2. Figure 12 displays a schematic representation of the bidirectional configuration interacting through inhibitory connections associated with the g A ii and g A ie conductances. Without interaction between the areas (g A ii = g A ie = 0), the neurons exhibit desynchronised activities. Desynchronised sender area can actuate on the initially desynchronised receiver area generating (A) synchronisation due to g A ii , (B) desynchronisation due to both g A ii and g A ie , and (C) silent activities of excitatory neurons due to g A ie conductances.

Inhibitory Connections
FIGURE 9 | Unidirectional inhibitory interaction between two coupled areas, where the neurons in the sender area are synchronised and the neurons in the receiver area are initially desynchronised. (A) The predominant direction of rotation of the relative phase angle of the area 2, (B) the standard deviation of the relative phase angle of the area 2, (C) the mean relative phase angle of the area 2, (D) the mean order parameter of the area 2. The circle, square, and triangle symbols indicate g A ii = 3.6 nS and g A ie = 0.4 nS, g A ii = 2 nS and g A ie = 0.9 nS, and g A ii = 0.4 nS and g A ie = 3.6 nS, respectively. Synchronised sender area can generate synchronisation and silent activities of excitatory neurons of the receiver area depending on g A ii and g A ie .
In the parameter space g A ii × g A ie , the mean order parameter for the neuronal network ( Figure 13A) shows a region in which the neurons in the areas are synchronised, where a circle symbol is located. For the standard deviation of relative phase angle rotation of the area 2, we identify three regions with values about 0, as shown in Figure 13B. The inhibitory connections are responsible for decreasing the relative phase angle variability (circle, square, and triangle symbols). Figures 13C,D exhibits the mean order parameter of the areas 1 and 2, respectively. The regions of small variability of ′ 2 can correspond to the synchronous behaviour in each area or silence of some excitatory neurons.
In Figure 14, we show the raster plot (top), the relative phase angle (middle), and the order parameter (bottom) for (a) g A ii = 3.6 nS and g A ie = 0.2 nS, (b) g A ii = 1.5 nS and g A ie = 1.0 nS, and (c) g A ii = 0.2 nS and g A ei = 2.5 nS, according to the circle, square, and triangle symbols, respectively, denoted in Figure 13. Figure 14A displays the occurrence of phase synchronisation among the neurons between the areas. Depending on the values of g A ii and g A ie , it is possible to observe partial anti-phase synchronisation and silenced excitatory neurons, as shown in Figures 14B,C. The circle, square, and triangle symbols indicate g A ee = 0.9 nS and g A ei = 0.2 nS, g A ee = 0.5 nS and g A ei = 0.9 nS, and g A ee = 0.1 nS and g A ei = 1.8 nS, respectively. The bidirectional excitatory connectivity can generate phase and anti-phase synchronisation between the two areas due to g A ee and g A ei , respectively.
FIGURE 12 | Schematic representation of bidirectional interaction between the areas by means of inhibitory connections.

CONCLUSIONS
In this work, we investigate the influence of excitatory and inhibitory connections between areas in neuronal synchronous behaviour. We build a network with two areas formed by excitatory and inhibitory neurons. The neuron dynamics is modelled by means of an adaptive exponential integrate-andfire (AEIF) model, that is able to mimic known neuronal activities. We consider unidirectional (sender-receiver) and bidirectional interactions between the areas, as well as different coupling configurations.
In the unidirectional interaction, firstly we analyse the dynamical behaviour of the receiver area with excitatory connections from the sender area. When the neurons in the FIGURE 13 | Bidirectional inhibitory interaction between two areas with desynchronised neurons. (A) The mean order parameter of the neuronal network, (B) the standard deviation of the relative phase angle of the area 2, (C) the mean order parameter of the area 1, and (D) the mean order parameter of the area 2. The circle, square, and triangle symbols indicate g A ii = 3.6 nS and g A ie = 0.2 nS, g A ii = 1.5 nS and g A ie = 1 nS, and g A ii = 0.2 nS and g A ei = 2.5 nS, respectively. The bidirectional inhibitory connectivity can generate phase and anti-phase synchronisation between the two areas depending on g A ii and g A ie . In both areas, silent activities of excitatory neurons are observed for small g A ii and large g A ie values, respectively.
sender area are desynchronised, depending on the conductances values, counter-clockwise and clockwise rotation can arise in the receiver area. For synchronised neurons in the sender area, it is possible to observe phase and shift-phase synchronisation. Secondly, for inhibitory connections from the sender area, we find values of the conductances in which the neurons in the receiver area can be silenced, namely, they do not spike for a long time. The inhibitory connections can also induce synchronous behaviour in the neurons that belong to the receiver area even when the neurons in the sender area are desynchronised. For synchronised or desynchronised neurons in the sender area, the excitatory (inhibitory) connections to the excitatory (inhibitory) neurons in the receiver area generate an increase in the relative phase angle of the receiver area. Otherwise, excitatory (inhibitory) connections from the sender area to inhibitory (excitatory) neurons in the receiver area reduce the relative phase angle of the receiver area. We also verify that the synchronised sender area is more efficient to reduce the variability of the relative phase angle of the receiver area than the desynchronised one.
With regard to bidirectional interactions, the excitatory connections to the excitatory neurons can induce phase synchronisation, while to inhibitory neurons can favour antiphase synchronisation. In our work, the anti-phase mechanism due to the inhibitory connections is similar to the mechanism FIGURE 14 | Raster plot (top), relative phase angle (middle), and order parameter (bottom) for (A) g A ii = 3.6 nS and g A ie = 0.2 nS, (B) g A ii = 1.5 nS and g A ie = 1 nS, and (C) g A ii = 0.2 nS and g A ei = 2.5 nS. We consider bidirectional interaction with inhibitory connections between the areas, where the neurons in the sender and receiver areas are initially desynchronised. Due to the bidirectional inhibitory interaction, (A) phase synchronisation can emerge due to g A ii , (B) anti-phase synchronisation due to a combination of g A ii and g A ie , and (C) silent activities of excitatory neurons due to large values of g A ie .
reported by Kim and Lim (2020). They demonstrated the existence of phase-shift synchronisation among three cluster networks due to inhibitory synaptic coupling. We verify that the inhibitory connections from the areas to the inhibitory neurons of other ones can generate phase synchronisation between them due to the bidirectional interaction. In addition, when the inhibitory connections arriving at the excitatory neurons between the areas are strong, silence activities of the excitatory neurons are observed. Our simulations suggest that the excitatory and inhibitory connections from one area to another play a crucial role in the emergence of phase, anti-phase, shift-phase synchronisation between the neurons in the areas. Our results should be useful to clarify how these types of synchronisation emerge in neuronal areas. For more than two areas, we expect to find phase synchronisation due to g A ee . Nevertheless, we believe that more complex patterns related to the synchronous behaviour will arise. In future works, we plan to study the emergence of neuronal synchronisation in more than 2 coupled brain areas. We will also analyse the influence of different interactions on the neuronal activities as proposed in the model of Potjans and Diesmann (2014).

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
All authors discussed the results and contributed to the final version of the manuscript.

FUNDING
The authors acknowledge the financial support from the Brazilian Federal Agency CNPq, under Grant nos. 137141/2020-3 and 302665/2017-0, and the São Paulo Research Foundation (FAPESP, Brazil), under Grants Nos.