Synchronization and Rhythm Transition in a Complex Neuronal Network

Synchronization and rhythm transition in an excitatory-inhibitory balanced cortical neuronal network are investigated in this paper. A small-world neuronal network is performed to be the cortical region of cerebral cortex, which is composed of different types of Izhikevich neurons. The combination of regular spiking (RS) cells, chattering (CH) cells, or mixed RS and CH cells are imitated as excitatory neurons, whereas fast spiking (FS) cells mimic inhibitory neurons. We mainly focus on the effect of different types of neurons on synchronization and rhythm transition of the neuronal network. The simulation results illustrate that it is easier for the neuronal network with CH excitatory neurons to achieve synchronization, and it is also more susceptible to the parameters than the network with RS neurons. Together with the weight of the synaptic connection, the structure of the small-world network is also explored to identify its influence on the synchronization and the rhythm transition. Moreover, we found that the synchronization performance could be improved by increasing the degree of influence among neurons. Importantly, the synchronization states are associated with the rhythm transitions. Especially, the consistency of synchronization of the neuronal network and $\gamma $ band rhythm is illustrated. In addition, our results could have an important significance on further understanding brain rhythms and synchrony in neuroscience.


I. INTRODUCTION
The human brain is the most sophisticated and interactive network with nontrivial topological properties, consisting of 10 11 neurons, with each neuron connecting to more than 10 4 neurons via synapses [1], [2]. Understanding how the human brain, such a complex dynamical system, is capable of processing information accurately in the cortex is a significant issue in the fields of neuroscience, signal transmission, and cognitive informatics [3], [4]. For these reasons, much attention has been paid to brain rhythms observed in scalp electroencephalogram (EEG) and synchronous dynamical processes in complex networks [5], [6].
Electroencephalogram (EEG) signals are recordings of oscillations of brain electric potentials placed on the human scalp and reflect the firing of neurons in the cerebral The associate editor coordinating the review of this manuscript and approving it for publication was Zhanpeng Jin . cortex [7]- [9]. It is shown that the EEG contains different frequency components which show different brain states and cognition. Normally, EEG frequency ranges are categorized as δ = [ 1 to 4 Hz], θ = [4 to 8 Hz], α = [8 to 13 Hz], β = [13 to 30 Hz] oscillation, and high frequencies are referred to as γ = [more than 30 Hz] activity. These different frequency bands may correspond to different behavioral states, cognitive capabilities, or human brain locations. For example, the δ rhythm is dominant during the deep sleep state, whereas the β rhythm correlates with normal wakeful consciousness and emerges when a person is focusing thinking or solving problems. These brain rhythms occur through synchronization between various spiking or bursting characteristics in neural circuits. Neuronal synchronization also plays a key role in transferring and coding mechanisms among different regions of the cerebral cortex. Therefore, it is of great importance to investigate the processing of synchronization and rhythm transition in the brain.
In the past few decades, there have been various techniques used to concern the synchronization of neuronal networks [10]- [13]. On one hand, some scholars have studied how the external time periodic stimulation affected the neuronal synchronization in the brain [14]. On the other hand, considering that time delay is ubiquitous in the real world and is often the source of complex behaviors of the dynamical system, the researchers investigated the synchronization problems for neuronal networks with time-delayed feedback [11], [15]- [17]. In addition, the synaptic plasticity providing the basis for memory and learning of brain, was taken into consideration in many literatures [18]- [20]. However, most of the previous theoretical and numerical works on synchronization mainly based on the single neuron type in the neuronal networks. Due to the complexity of the brain, studying the dynamical behaviors for both cases of excitatory and inhibitory neurons is becoming more and more popular in complex neuronal networks [21], [21]- [24]. The authors have explored the interconnectivity between excitatory-inhibitory neural networks [23], [24], and they also analysed that the strength of inhibitory intra-connectivity plays an essential role in determining the dynamics of Pyramidal Interneuron Network Gamma (PING) mechanism. Here, we will analyse the synchronization and rhythm transition in an excitatory-inhibitory balanced cortical neuronal network.
To better simulate the dynamical system of the human brain, we need to select a suitable model for the local dynamics of a single neuron to offer computational advantages of large-scale brain models. The Hodgkin-Huxley (H-H) model has played a prominent role in computational neuroscience, relying on unprecedented innovation [25]. It consists of four differential equations and up to ten parameters to simulate membrane potential and different ion channels. Despite that it can actually simulate all properties of a neuron, the limitation of this model is that it requires massive computing power and a large number of free parameter estimation. Therefore, many efficient models have been proposed to reduce the complexity of H-H model, such as the FitzHugh-Nagumo model, the integrate-and-fire model, the Morris-Lecar model, the Hindmarsh-Rose model, and so on [26]- [29].
In fact, as a coin has two sides, when reducing the model's complexity of computation, the biophysical accuracy is sacrificed, so fewer firing patterns and properties are produced. To solve these problems, a simple model of a spiking neuron was developed up by Izhikevich [30], [31]. Compared with the H-H model, the Izhikevich model improves computational efficiency, because it only consists of two coupled differential equations, yet is still biologically plausible. Therefore, we chose the Izhikevich model to model the dynamics of a single neuron in this paper [30]- [32].
With the goal of better understanding the dynamics of the neuronal network, special attentions need to be paid to understanding how neurons are connected. In fact, in the human brain, the synaptic connections have been found to have complex network topologies. Many different kinds of network structures are considered to investigate the synchronization and rhythm transition in the neuronal network [12], [13], [18]- [22], [33]- [35]. Qin have studied the vibrational resonance in a randomly network [21]. Kim investigated the burst synchronization in a scale-free neuronal network [20], [35]. In [18], the authors explored the bursting synchronization in a modular neuronal network. Then, the network topologies, including the regular ring network, small-world network, Erdos-Renyi and scale-free network, are compared as the key element affecting the beta-rhythm oscillations and synchronization transition, and it was found that the clustering coefficient played an important role in population synchronization in neuronal networks [6].
Mostly, the small-world network has been considered as the neuronal system since proposed by Watts and Strogatz in 1998 [36], [37]. In addition, studies have proven that the human brain cortex has small-world properties [12], [16], [33], [35], [38]. The small-world model provides a powerful tool to investigate the structures and dynamics of the human brain. A lot of work concerned the neuronal dynamics of small-world networks [16], [38], [39]. Qu and Wang have explored the synchronization and spatiotemporal behaviors of small-world neuronal networks based on different neuron models [13], [22], [33], [34]. Mahsa Khoshkhou's work was based on the RS neurons [6]. Kim's work consisted of excitatory RS pyramidal neurons and inhibitory FS interneurons [22].
However, these previous work have restricted their discussion mostly on the case of one type of excitatory neurons and one type of inhibitory neurons. They have not paid attention to the effects of different types of excitatory and inhibitory neurons on synchronization and rhythm transition in a complex network.
Motivated by the previous works, we explore the rhythm transition and the synchronization of cortical neuronal networks on different types of excitatory and inhibitory neurons. We consider the situation with the regular spiking (RS) pyramidal neurons, chattering (CH) neurons, and hybrid RS and CH neurons being as the excitatory neurons, whereas fast spiking (FS) interneurons being chosen as the inhibitory neurons.
Together with the synaptic coupling, the network property is also considered to be an influence factor of the rhythm types, the synchronous status, and particularly the relationship between synchronization and rhythm transition.
Accordingly, the outline of this paper is organized as follows: In section II, the single Izhikevich neuron model and its individual dynamics are described. A complex small-world neuronal network and the synchrony measurement we used are also introduced in this section. Then, in section III, the effects of many key parameters on synchronization and rhythm transition with excitatory and inhibitory neuronal networks are investigated. There are three different cases to be discussed: (i) the network consists of excitatory regular-spiking (RS) and inhibitory fast-spiking (FS) neurons, (ii) the network is composed of excitatory chattering (CH) and inhibitory fast spiking (FS) neurons, and (iii) the network is constituted of hybrid excitatory regular spiking (RS) and excitatory chattering (CH) neurons, as well as inhibitory fast spiking (FS) neurons. Section III is devoted to presenting and discussing the simulation results about the spatiotemporal order, synchronization, and the rhythm transition. Finally, a brief discussion and conclusion of our work are given in section IV.

A. SINGLE NEURONAL MODEL AND ITS INDIVIDUAL DYNAMICS
In this paper, a computationally efficient model, proposed by Izhikevich, is chosen to simulate the individual dynamics of each neuron in the studied network. This model exhibits rich and complex spiking and bursting behaviors of a single neuron with the appropriate choice of parameters, thus allowing us to analyze the different dynamical properties [30], [31]. The model is expressed by a two-dimensional (2-D) system of ordinary differential equations as below: with the auxiliary after-spike resetting: Here, the parameter v denotes the voltage in mV across the neuron membrane. The slow variable u represents the recovery current and accounts for the activation of K + ionic currents and the inactivation of Na + ionic currents. The synaptic current or DC-currents delivered to the neuron are described by I . According to the equation (3), when the membrane voltage v reaches the apex of a spike, V peak = +30 mV, it is reset to c, and the recovery variable u is added by d.
The parameters a, b, c and d are dimensionless, and assigning them various values results in a great variety of firing characteristics and dynamics of classical Izhikevich neurons, including excitatory and inhibitory cortical cells in intracellular recordings. Noting that by adjusting these parameters, more than twenty biological firing dynamics can be achieved, we wonder how different types of classic Izhikevich neurons affect the synchronization and rhythm transition in an excitatory-inhibitory balanced cortical neuronal network.
In our computations, we use these three fundamental classes of firing patterns of neocortical neurons, and the parameter values for them are listed in Table 1.
The synaptic currents or injected DC-currents are passed through the variable I . In terms of its changing, the neuron model produces rich firing dynamics, such as periodic spiking and bursting patterns.
Therefore, for the purpose of getting the detailed rhythm transition process, the bifurcation diagrams of the inter-spike interval (ISI) and the firing rate f with respect to the input currents I for different types of the Izhikevich neurons are shown in Fig. 2 and Fig. 3.
RS neurons, when there is a persistent input current, fire spikes continuously. As I increases, the ISI is nonlinearly reduced and the firing frequency spans the range from 5 Hz to 50 Hz ( Fig. 2(a) and Fig. 3(a)). To better understand CH neurons, Fig. 2(b) and Fig. 3(b) show the firing rate of the inter-burst and the intra-burst. It should not be overlooked that the firing frequency is lower for RS neurons, whereas for the CH neurons, the inter-burst firing frequency can be as high as 50 Hz, and it contributes to the γ oscillations in the brain. Similarly to the RS neuron, the FS neuron can fire in response to an injected current, but the firing frequency is even greater than 300 Hz as shown in Fig. 2(c) and Fig. 3(c).

B. NETWORK STRUCTURE
As mentioned in the above section, we consider the Izhikevich model as the local node in the network, and the pulse-coupled neuronal network (PCNN) is used to perform our network and we take the small-world network as the underlying structure [36], [37], [43]. The original small-world model was proposed by Watts and Strogatz [36]. To avoid the isolated vertices, the NW small-world model was brought by Newman and Watts [37].
When the number of network N is large enough and the rewiring probability p is sufficiently small, the properties of the NW model are quite similar to those of the WS model. Recent studies about the neuronal system are also under the NW small-world networks [18], [44]. Hence, in our simulation, we use the Newman-Watts (NW) small-world structure. The NW small-world network starts from a regular ring network with N nodes, and each node is coupled to its 2k nearest neighbors (k on either side). Then, connections between pairs of unconnected nodes are added with the probability p. For p = 0, the network reduces to the originally nearest-neighbor coupled network, and for p = 1, it evolves to the globally connected network. In order to comply with the characteristics of a small-world, we focus on the cases where the shortcut adding probability p is less than 0.2. The dynamics of the neuronal network can be described by the following equations: with the auxiliary after-spike resetting: Here, v i (t) and u i (t) are the state variables of the ith neuron at time t. The neuron index is denoted by i, j = 1, 2, . . . , N , and I ext i describes the external applied current, which is the intensity of the DC-current. The noise current I noise i = Dξ i represents the external or intrinsic fluctuations of the neuron itself, where D is referred to as the intensity of noise. In this paper, we fix the noise intensity for excitatory neurons as D exc = 5 and the inhibitory neurons as D inh = 2, according to the PCNN in [30]. ξ i is a random process without time correction, and it is uniformly distributed on the interval [−1, 1]. I syn i accounts for the total incoming synaptic currents into neuron i from other neurons. Based on the structure of the network, the connectivity matrix (g ij ) N ×N can be defined. Here, g ij represents the matrix element, which takes the following form: if the neuron i is connected with the neuron j (i = j), then g ij = g ji = 1; otherwise, g ij = 0 and g ii = 0 for all i. In addition, (g ij ) N ×N is an asymmetric and irreducible matrix. (w ij ) N ×N describes the coupling strength VOLUME 8, 2020 of the synapse from neuron j to neuron i, and the firing of the jth neuron instantaneously changes the value v i by w ij . The connection weights between inhibitory neurons are always set as w inh = −1, and the connection weights between excitatory neurons are w exc , which is a control parameter in our simulations.
In this paper, we construct a network consisting of N = 1000 neurons in real-time. Motivated by the anatomy of a mammalian cortex, the ratio of excitatory to inhibitory neurons is chosen to be 4 : 1, which means there are 800 excitatory neurons and 200 inhibitory ones. Besides the synaptic inputs, each neuron receives a noisy thalamic input. The RS and CH neurons are used to model the excitatory cells. The FS neurons are used to model the inhibitory cells.

C. MEASUREMENT
Many measures have been used to quantify the dynamics of the network activity in previous studies. Foremost among them is how to characterize the synchronization degree of the network. We adopted the mean field potential (MFP) as a global parameter for visual presentation on the synchronization degree of all neurons. It is defined as the mean membrane potential of all neurons in the network. Since the most significant voltage changes occur during the action potentials, significant MFP deflection can only be expected when most neurons simultaneously produce an action potential at the same time. Only when all neurons fire in a complete coincidence that the maximum amplitude can be achieved, and the MFP reaches the maximum amplitude, which corresponds to population synchronization. However, when the neurons fire with random phase shifts, the MFP will mainly show small deflections [33], [45].
Although the MFP can exhibit the behaviors of the whole network, it can't clearly show the synchronization transition. Therefore, another synchrony measure, an improved method presented in [46], [47], is used to quantify the degree of spiking synchronization in the neuronal network. It involves convolving a Gaussian function with each potential time of every neuron to generate variables v i (t). More details of this method are described as follows: where V (t) is the averaged voltage, defined in equation (9), and N indexes the number of neurons in the network. Additionaly, σ i and σ are the variance of each neuron's voltage and the averaged voltage of the whole network, respectively. In equations (10) and (11), · denotes the time average. Then the synchrony measure S is defined as equation (12). S = 0 means completely asynchronous spiking state, whereas S = 1 is equivalent to fully synchronous activity.
The codes implementing these simulations were written by MATLAB. During the simulations, all the equations were numerically integrated using the Fourth-Order Runge-Kutta Method with a time-step size of 0.2 ms. The simulations of the neuronal network have been improved and after 1000 ms from random initial conditions, when the network reaches a statistically stable state. Throughout this paper, the spatiotemporal patterns of the neuronal network are plotted such that the excitatory cells are given the lowest neuron index toward the bottom of the y-axis whereas the inhibitory cells are given the highest neuron index towards the top of the y-axis.

III. RESULTS
Considering the possibility of a neuron's dynamics, we investigate both the RS and CH cells as the excitatory neurons and the FS cells as the inhibitory neurons. We divide the study into three different cases to discuss: (i) the network consists of excitatory RS neurons and inhibitory FS neurons, (ii) the network is composed of excitatory CH neurons and inhibitory FS neurons, and (iii) the network is constituted of hybrid excitatory RS and CH neurons, as well as inhibitory FS neurons. As for the connection structure, we focus on the effects of the synaptic strength of connections among excitatory neurons, and the critical parameters of the small-world network.

A. SYNCHRONIZATION AND RHYTHM TRANSITION WITH REGULAR SPIKING AND FAST SPIKING NEURONS
In this subsection, the RS neurons are used to model the excitatory neurons, and FS neurons are used to model the inhibitory neurons.
First of all, we studied how the weight of the synaptic connection among the excitatory neurons influences the rhythm transition and collective dynamics. The matrix (w ij ) N ×N indicates the weight of the synaptic connection. As mentioned above, we set the synaptic connection weight of inhibitory neurons to other neurons as w inh = −1. As shown in Fig.4, by changing the synaptic connection weight of excitatory neurons w exc , the neuronal network exhibits various collective behaviors. In this situation, the value of w exc is tuned from 0 to 1. When the excitatory synaptic connection weight is w exc = 0.5, the neurons randomly fire spikes (Fig. 4(a)), and the neurons exhibit collective rhythm behavior in the frequency range corresponding to that of the mammalian cortex in the awake state, such as α band rhythm (10 Hz). As shown in Fig. 4(b), with the synaptic connection weight being increased to w exc = 0.63, the neurons become synchronous. Then we continued to increase the connection weight w exc from 0.63 to 0.8, all neurons fire spikes at the same rhythm (the time interval is about 240 ms) and reach complete synchronization.
Another question we study is how the characteristics of the network influence the collective dynamics and the rhythm transition. Two important parameters of the small-world network which determine the network structure are considered. The one is connection probability p, and the other is the number of nearest neighbors of each node k.   The effect of connection probability p is investigated with other parameters being set as constants. The spatiotemporal patterns of the network in terms of various connection probability p are shown in Fig. 5. The synchronization degree of the system was not distinctly changed as the number of connection probability p was added from 0 to 0.2. Thus, complete synchronization is hard to achieve with middle values of w exc and k.
How the number of nearest neighbors of each node k affects the degree of synchronization and rhythm transition is further studied. The corresponding spatiotemporal patterns with different values of k are shown in Fig. 6. It is demonstrated that the neurons are almost synchronous with k being increased, the neuronal network becomes more and more synchronous. It is shown that the network achieve complete synchronization in the case of k = 350.
In order to further quantify the above results, we calculated the MFP with respect to the above key parameters and the results are shown in Fig. 7. According to the definition of MFP, the maximum amplitude is reached only when all neurons are completely coincident [45]. More intuitively, the larger the amplitude, the more synchronous of the network. Figure 7(a) gives the evolution of the MFP versus the synaptic connection weight w exc . It is obvious that with the connection weight w exc being increased, the neuronal network will show the transition from random firing (w exc < 0.63) to near synchronization (0.63 < w exc < 0.75), then VOLUME 8, 2020  Fig. 4. (b) The evolution of the MFP versus the connection probability p, and other parameters are same as Fig. 5. (c) The evolution of the MFP versus the number of nearest neighbor k of each node, and other parameters are same as Fig. 6.   FIGURE 8. (a) The evolution of the percentage of every rhythm band versus the synaptic connection weight w exc , and other parameters are same as Fig. 4. (b) The evolution of the percentage of every rhythm band versus the connection probability p, and other parameters are same as Fig. 5. (c) The evolution of the percentage of every rhythm band versus the number of nearest neighbor k of each node, and other parameters are same as Fig. 6.
to complete synchronization (w exc > 0.75). As for the effect of the connection probability p, the variation of the MFP is exhibited in Fig. 7(b). The synchronous state changes slowly with p increasing. When p < 0.2, the MFP amplitude has no evident change in its range, which means that the neuronal network is always in the same state. Furthermore, the effect of the number of nearest neighbor k of each node is shown in Fig. 7(c). The neurons randomly fire spikes, and the network is almost synchronous for the value of k being from 0 to 180. Then the network transits to complete synchronization with k > 270.
To profoundly understand how the rhythm transfers, we also exhibit the evolution of the network rhythm with the change of these key parameters, which are presented in Fig. 8. Figure 8(a) provides the evolution of every rhythm band in the neuronal network versus the connection weight w exc . It indicates that the rhythm changes slowly as the connection weight w exc increases from 0 to 0.63. With the synaptic connection weight w exc being further increased, the rhythm evolves considerably. The percentage of γ band rhythm increases dramatically, and other types of rhythm decrease. Around w esc = 0.75, the γ band rhythm becomes dominant and the percentage of θ band rhythm decreases below 0.1 and is stable at this level. However, the other types of rhythm all die out. This means that the main rhythm is γ when the neuronal network arrives complete synchronization. In addition, the percentage of every rhythm band in the neuronal network with various connection probability p is shown in Fig. 8(b). The phenomenon is completely different from the previous one. When the value of k increases, the percentages of α and θ band rhythms fluctuate largely, but the main rhythms are always α and θ band rhythms when p < 0.2. Figure 8(c) gives the evolution of the rhythm with the number of nearest neighbor k of each node in the neuronal network. It is demonstrated that with k being increased, the main rhythm changes from θ band rhythm to α band rhythm, then to γ band rhythm until k > 270.
Summarizing the above results, we find that by increasing the weight of synaptic connections and the number of nearest neighbors of each node, the degree of synchronization can be dramatically improved. Meanwhile, increasing the connection probability p, the synchronization status doesn't change obviously.
Comparing Fig. 7 with Fig. 8, a more exciting discovery is that there is some relationship between the synchronization and rhythm transition. The percentage of every rhythm band changed according to the synchronization degree of the neuronal network. When the neurons randomly fire spikes, the main rhythms are θ and α. Then the neuronal network rapidly arrives at complete synchronization. In addition,   the percentage of γ band rhythm suddenly rises from less than 10% to larger than 80% or more. It is illustrated that the percentage of γ band rhythm increases to the same degree as synchronization. In other words, the positive correlation of synchronization degree and the percentage of γ band rhythm is consistent.

B. SYNCHRONIZATION AND RHYTHM TRANSITION WITH CHATTERING AND FAST SPIKING NEURONS
In this subsection, the chattering (CH) neurons are used to model the excitatory neurons, and fast spiking (FS) neurons are used to model the inhibitory neurons. Similar to the former subsection, how the weight of the synaptic connection w exc , the connection probability p, and the number of nearest neighbor of each node k in the small-world network influence the synchronization and rhythm transition will be discussed.
At first, the connection probability p and the number of nearest neighbors of each node k are fixed as p = 0.1 and k = 200. The weight of the synaptic connection w exc is taken as the control parameter. The spatiotemporal patterns of the neuronal network with respect to the synaptic connection weight w exc are depicted in Fig. 9. When the excitatory synaptic weight w exc = 0.3, the neurons casually fire spikes as shown in Fig. 9(a). Further increasing w exc causes the neurons in the network to spike with the same rhythm and  Fig. 9. (b) The evolution of the MFP versus the connection probability p, and other parameters are same as Fig. 10. c The evolution of the MFP versus the number of nearest neighbor k of each node, and other parameters are same as Fig. 11.   FIGURE 13. (a) The evolution of the percentage of every rhythm band versus the connection weight w exc , and other parameters are same as Fig. 9. (b) The evolution of the percentage of every rhythm band versus the connection probability p, and other parameters are same as Fig. 10. (c) The evolution of the percentage of every rhythm band versus the number of nearest neighbor k of each node, and other parameters are same as Fig. 11. the complete synchronization is eventually achieved when w exc = 0.4, as shown in Fig. 9(c).
Secondly, the effect of the connection probability p is presented in Fig. 10. Other parameters are set as constants when p is tuned. The weight of synaptic connections among excitatory neurons w exc is fixed as 0.35, and the number of nearest neighbors of each node is fixed as k = 200. Like the above results, it is illustrated that the network is in completely synchronization status if the connection probability is set higher, such as p = 0.2 (Fig. 10(c)), and the neurons fire with the same rhythm (time interval is about 240 ms).
Finally, how the number of nearest neighbor of each node k affects the synchronization and rhythm transition is also considered. The other parameters are set as w exc = 0.35 and p = 0.1. The spatiotemporal patterns of the neuronal network with different values of k are shown in Fig. 11. It is demonstrated that as the value of k increases, the synchronization degree of the network increases. At the same time, all the neurons are self-organized into assemblies and exhibit collective rhythmic behavior.
The evolution of the MFP with the change of these parameters in Figure. 12 helps us to quantify the previous results. Compared to the previous study, it is proven that as the weight of synaptic connection w exc increases from 0.2 to 0.4, the connection probability p increases from 0.05 to 0.15 and the number of nearest neighbors of each node k increases from 100 to 200, the MFP increases dramatically, which means the system becomes more synchronous. Further increasing w exc , p and k causes the MFP to reach the maximum and the complete synchronization is achieved.
The phenomenon that interests us is the connection probability p. In contrast to the situation of RS cells being the excitatory neurons, the connection probability p also significantly influences the system synchronization, as Fig. 12(b) shows. As the value of p continues to increase, the system becomes more and more synchronous. This means that the neuronal network with chattering excitatory neurons can arrive synchronization more easily and the system is more sensitive to all the parameters.
To better understand the rhythm transition, the percentages of all types of rhythm are presented in Fig. 13. It is indicated that the γ band rhythm is always in the dominant position. Furthermore, when the degree of synchronization is enhanced with parameters of great value, the percentage of the γ band rhythm also increases significantly to 80% or even higher whereas the ratios of the other types of rhythm become lower and maintain at a relatively low level. Especially for the connection probability p, the results differ entirely from the above subsection. When the value of p increases, the system becomes more synchronous. As a result, the percentage of γ band rhythm also becomes larger and more stable at this level.
As proved by the above results, synchronization and rhythm transition are closely related. The consistency of the synchronization degree of the neuronal network and the γ band rhythm are illustrated by comparing Fig. 12 with Fig. 13. We also found that a model utilizing CH cells as excitatory neurons was easier to synchronize than models employing RS neurons and the system was more sensitive to all parameters.

C. SYNCHRONIZATION AND RHYTHM TRANSITION WITH HYBRID EXCITATORY NEURONS AND FAST SPIKING INHIBITORY NEURONS
In fact, in the real cerebral cortex of mammals, it is almost impossible to have only one type of excitatory neurons. Therefore, the rhythm transition and dynamics based on mixed excitatory neurons are discussed in this subsection.
The mixed RS and CH cells are used to model the excitatory neurons and the FS cells are still chosen as the inhibitory neurons. In this part, we mainly wonder how the number of RS and CH neurons affects the dynamics of the entire neuronal network.
For the convenience of intuitive display, the widely applied synchrony measure S is used to get the synchronization status of the network. Different from the global parameter MFP, S is simple and it is the application of variance analysis with the time of each potential for every neuron, which lies in [0, 1]. Due to the difference of the excitatory neurons, the entire system should synchronize if S > 0.6. Similar to the above subsections, we also identify the effects of the synaptic connection weight w exc , the connection probability p, and the number of nearest neighbors k of each node. We also calculate the synchrony measure MFP versus the considered parameters and present them in Fig. 12. In order to identify the effect of the number of different excitatory neurons, the number of RS neurons and one of the parameters considered above are considered as two control parameters.
Firstly, The number of RS neurons and the weight of synaptic connection w exc are controlled and the result is shown in Fig. 14(a). The small-world parameters p and k are fixed as p = 0.1 and k = 200. We set the values of p and k to be middle such that the neuronal network is in the asynchronous state initially according to Fig. 5. As the number of RS neurons increases, the synchrony measure S gets lower and lower. In other words, the larger the number of CH neurons, the better the synchronization degree of the network. It is also shown that combined with an increase of w exc , the degree of synchronization rises in a curve. This result confirms our previous conclusions. Both the synaptic connection weight and the CH excitatory neurons help to improve the synchronization degree of the network.
So what happens when the structure of the small-word neuronal network changes?
Then the synchrony measure S as a function of the number of RS neurons and the connection probability p is revealed; see Fig. 14(b). we realize that p has a relatively slight impact on the synchronization degree. Similar to the results presented in Fig. 7(b) and Fig. 12(b), when the number of RS neurons is greater than the number of CH neurons, there is no significant change in the degree of synchronization, regardless of the value of p. However, when the number of CH neurons is dominant, no matter how much p is, the degree of synchronization is stable at a relatively high level.
Furthermore, we investigate the effect of the number of RS neurons and the number of nearest neighbors k. As shown in Fig. 14(c), it also proves that the number of CH neurons and the number of nearest neighbors k are advantageous to the network synchronization. As the number of nearest neighbors increases, the synchronization degree increases nonlinearly, which is similar to those of the weight of synaptic connection w exc . This indicates that the synchronization degree of the network can be improved by increasing the control parameters to finally stable value.
Aside from this, we also discover that the ratio of RS neurons to CH neurons around 1 is a special case. In this case, the synchronous degree of the network transits first fall from 0.6 to 0.4 and then rise from bad to good, and finally stabilize at a high level, see Fig. 14(a) and Fig. 14(c). This phenomenon may have important implications for the study of the γ mechanism, called ''PING'' in the pyramidal neuronal network [48], [49].

IV. CONCLUSION
In this paper, we proposed an excitatory-inhibitory balanced cortical neural network based on different types of Izhikevich neurons to investigate the synchronization and rhythm transitions. The neuronal network exhibited different synchronization status and rhythm transitions by tuning the weight of the synaptic connections, the synaptic connection probability, and the number of nearest neighbors of each node.
We divided the study into three different situations to discuss.
In the first place, we used RS neurons as the excitatory neurons and FS neurons as the inhibitory neurons in the network. It was illustrated that increasing the weight of the synaptic connection and the number of nearest neighbors for each node, the synchronization degree dramatically increased. However, it was difficult to improve the synchronization degree of the neuronal network using the synaptic connection probability. It was proved that the percentage of the γ band rhythm was increasing with the increase of the synchronization degree.
In the second place, CH neurons were used to emulate the excitatory neurons and FS neurons were chosen as the inhibitory neurons in the network. By increasing the weight of the synaptic connection, the synaptic connection probability, and the number of nearest neighbors of each node, the neuronal network became more and more synchronous. That is different from the previous situation, which means that the neuronal network utilizing CH neurons can achieve synchronization more efficiently, and the network was more sensitive to all the parameters.
Mixing the above two cases, hybrid RS and CH neurons were used to simulate the excitatory neurons, whereas FS neurons were always acting as the inhibitory neurons in the network. For increasing the coupling strength among synapses, or increasing the number of synaptic edges in the network, it facilitates the network to synchronize. Whereas the difference is that the change of p has less influence on the degree of synchronization compared with k and the number of CH neurons. When N CH is more than N RS , no matter how large p is, the entire system will maintain a high synchronization state. It means that the neural network with chattering neurons can arrive synchronization easier than the case of regular spiking neurons, and it is more susceptible to the parameters.
So what are the causes of these phenomena? We attribute the leading cause of these phenomena to the firing mechanisms of different neurons. RS neuron exhibits spiking, and it has a long interspike interval, so the firing frequency is lower, as shown in Fig. 2(a) and Fig. 3(a). However, CH neuron shows bursting. The firing frequency is high inside the burst even though the firing time is averaged when we calculate S, as shown in Fig. 2(b) and Fig. 3(b). Similarly, for FS neuron, it has a high firing frequency. Comparing these three different types of neurons, we believe that it is more easier for CH and FS neurons to match the firing frequency instead of RS and FS neurons, which is just in line with our experiments. Thus, the neuronal network with CH neurons can arrive synchronization more efficiently than that with RS neurons.
For increasing the number of the nearest neighbors k in the network, it facilitates the network to synchronize. Whereas the difference is that the change of p has less influence on the degree of synchronization compared with k. This is because our value of k is large enough and the value of p is less than 0.2. The change of p has a small effect on the number of synaptic connections.
Currently, we have evidence that γ oscillations are generated by synchronous behaviors of fast spiking (FS) interneurons [50]. Furthermore, the relationship between γ oscillations and synchronization may have an important implication for the study of the Pyramidal Interneurons Network Gamma (PING) mechanism in E-I networks [51]. Although there are important discoveries revealed by these studies, there are also limitations in our model. Since our work mainly aims to study how dynamical responses depend on excitatory and inhibitory neurons, the electrical synapses and chemical synapses are not specifically discussed. Therefore, we will further explore the case of electrical synapses and chemical synapses in the future. In summary, our results might have significant implications for a better understanding of the brain rhythms, cognitive informatics, and synchrony dynamics in neuroscience. YUAN