Suppression of absence seizures by using different stimulations in a reduced corticothalamic‑basal ganglion‑pedunculopontine nucleus model

: Coupled neural network models are playing an increasingly important part in the modulation of absence seizures today. However, it is currently unclear how basal ganglia, corticothalamic network and pedunculopontine nucleus can coordinate with each other to develop a whole coupling circuit, theoretically. In addition, it is still difficult to select effective parameters of electrical stimulation on the regulation of absence seizures in clinical trials. Therefore, to develop a coupled model and reduce computation cost, a new model constructed by a simplified basal ganglion, two corticothalamic circuits and a pedunculopontine nucleus was proposed. Further, to seek better inhibition therapy, three electrical stimulations, high frequency stimulation (HFS), 1:0 coordinate reset stimulation (CRS) and 3:2 CRS, were applied to the thalamic reticular nucleus (RE) in the first corticothalamic circuit in the coupled model. The simulation results revealed that increasing the frequency and pulse width of an electrical stimulation within a certain range can also suppress seizures. Under the same parameters of electrical stimulation, the inhibitory effect of HFS on seizures was better than that of 1:0 CRS and 3:2 CRS. The research established a reduced corticothalamic‑basal ganglion‑pedunculopontine nucleus model, which lays a theoretical foundation for future optimal parameters selection of electrical stimulation. We hope that the findings will provide new insights into the role of theoretical models in absence seizures.


Introduction
Absence epilepsy is a common neurological disease of the brain that usually occurs in childhood [1][2][3].It is characterized by a spike wave discharge (SWD) occurring in the electroencephalogram (EEG) with a frequency of 2-4 Hz [4][5][6].Absence epilepsy can lead to cognitive and learning disabilities and temporary loss of consciousness, and these symptoms are related to the abnormality of corticothalamic circuits [7,8].Recent studies have shown that the reticular nucleus of the thalamus acts as a pacemaker in the control of an absence seizure [9][10][11].Similarly, much research evidence shows that the basal ganglia is also involved in the control of absence seizure by the corticothalamic circuit [12][13][14].To more conveniently study how the basal ganglia is involved in controlling absence seizure, Fan et al. [15] simplified the basal ganglia into a theoretical model of 3I:2O modulator to study the regulation of basal ganglia on absence seizure.In addition, some scholars have used the mean field model to analyze the dynamic state of neurons in absence epilepsy.Robinson et al. [16] established a corticothalamic network through the mean field model and successfully simulated the EEG under normal and seizure states.Chen et al. [17] established a basal ganglia-corticothalamic (BGCT) model through the mean model to describe the physiological mechanism of neuronal firing.Jiruska et al. [18] made an in-depth study on the theoretical mechanism of synchronization and desynchronization of epileptic seizures.However, most of the research only focuses on the network composed of a single corticothalamic circuit and basal ganglia or two corticothalamic circuits without considering the influence of the basal ganglia.Therefore, it is necessary to further establish a model composed of two corticothalamic circuits and basal ganglia.
The first motivation of the research was to investigate an expanding model to extend the original single-compartment cortico-thalamic model by coupling other neurons.In addition to the basal ganglia and corticothalamus related to absence seizures, there are also some studies indicating that the pedunculopontine nucleus (PPN) is also involved in the regulation of some neurological diseases [19], and some scholars have established a theoretical model consisting of the corticothalamic circuit, basal ganglia and PPN [20][21][22].Inspired by these models, we introduced a PPN structure into the model to explore the physiological mechanism of seizures.In the human brain, PPN controls movement through gamma amino butyric acid energy (GABAergric) and glutamatergic neurons [23,24].Some experimental studies have shown that PPN is also closely related to basal ganglia.PPN projects signals to the subthalamic nucleus (STN) and globus pallidus internal (GPi) of basal ganglia, and STN also gives feedback signals to the PPN [24,25].Niktarash et al. [26] established a coupling theoretical network between the two based on the projection effect of PPN on GPi.In addition to being connected to the basal ganglia, the experimental study by Wang et al. [27] shows that PPN can project choline to the thalamic reticular nucleus.These findings provide theoretical guidance for us to further improve the pathogenesis of absence seizure.
The second motivation of the research was to identify an optimal suppressive scheme by contrasting different stimulus paradigms and many stimulus parameters.Although there are many clinical trials on neurological diseases, it is difficult to achieve good results through drug treatment for some patients with seizures [28].For a small number of patients with seizures, it can take surgical resection to achieve good treatment results [29].However, some patients with seizures are still difficult to treat.As a new method of treating neurological diseases, electrical stimulation is widely used in clinical treatment and theoretical model research [30][31][32].Wang et al. [9] established a theoretical model to study the effect of electrical stimulation on absence seizures and concluded that electrical stimulation can suppress seizures well.In addition, some scholars have also studied the effect of electrical stimulation on seizures through animal experiments and made in-depth research on the mechanism of electrical stimulation [33][34][35].Recently, coordinated reset stimulation (CRS) was used to study the suppression of seizures as a pulse sequence sent to a specific part of the brain within a specific period of time [36][37][38].Fan et al. [15] used high frequency CRS to stimulate the various structures of the basal ganglia and achieved a good inhibitory effect on seizures.In addition, high frequency pulses are also used in animal experiments to study the relationship between the thalamic reticular nucleus and seizures [39].
However, previous models only considered interactions between a single corticothalamic network and the basal ganglia or two corticothalamic networks.For example, although the 3I:2O model simplifies the basal ganglia, which improves the running time of computing, it cannot reveal more complex neural connection pathways.Therefore, we established a model consisting of two cortical thalamic circuits, PPN and simplified basal ganglia based on the 3I:2O model, which can study more complex neural connection pathways.It is also possible to study the effect of electrical signals passing through the PPN on the electrical signals of the cerebral cortex.Although the model has become more complex, the calculation runs in less time than the non-simplified model.In terms of the effect of electrical stimulation on the brain, previous studies have only explored the inhibitory effect of fixed parameters on seizure.Therefore, solving inhibitory control problems needs much more exploration.
The key contribution of this research has three main creations.First, we simplified the basal ganglia into a modulator, which can efficiently receive signals from other parts.Second, we found that increasing the coupling strength of some pathways can suppress absence seizures.Third, we applied three different forms of electrical stimulation to the thalamic reticular nucleus and explored how the electrical stimulation parameters inhibited SWD under changing conditions.

Model and mathematical methods
To analyze the dynamic state of the basal ganglia, corticothalamic circuit and PPN, we built a reduced BGCT-PPN model, as shown in Figure 1.The model was composed of six types of neurons: Excitatory pyramidal neurons (PY1 and PY2); inhibitory interneurons (IN1 and IN2); thalamic relay nuclei (TC1 and TC2); thalamic reticular nucleus (RE1 and RE2); simplified basal ganglia modulator (3I:5O); PPN.The basal ganglia were simplified to a 2I:3O modulator in the previous study [9].In the model of this paper, it's simplified to a modulator with three inputs and five outputs stimulation.The mean field model is described [40][41][42]: where  ∈  ,  ,  ,  ,  ,  ,  ,  , ,  represents the combination of different types of neuron populations, p1 and p2 represent two different excitatory pyramidal neurons in the model, i1 and i2 represent inhibitory interneurons in two different parts, t1 and t2 represent two different thalamic relay nuclei in the model, r1 and r2 represent two different thalamic reticular nuclei in the model, m represents the simplified basal ganglia and ppn is the pedunculopontine nucleus. is the maximum mean discharge rate of neuron population and .y  is the firing threshold of each neuron population. is the standard deviation of the mean discharge threshold.The potential change in the cell can be expressed as follows [43][44][45]: representing the coupling strength of neuron population  to .  is the input transmission rate. is the axon time delay from neuron  to .  represents the axonal and dendritic filtering of the input signal. and  are the attenuation and rise time of the input signal.
where  =  / represents the damping rate of propagation,  represents the characteristic range of axons and  is the rate of signal conduction in the axon.In the paper, we only considered the damping rate of cortical excitatory pyramidal neurons.To facilitate the calculation, we set Mp = Mi and Vp = Vi, where p and i represent excitatory pyramidal neurons and inhibitory interneurons, respectively.
Based on previous research, we combined the basal ganglia, corticothalamic circuit and PPN as a new model, and the mathematical equations were described as follows: The kinetic equations of PY 1 under its own and external stimuli are as follows: The kinetic equation of TC 1 stimulated by the surrounding nervous nuclei is shown below: .

= 𝛼𝛽 −𝑉 (𝑡) + 𝑣 𝜙 (𝑡) − 𝑣 𝑆 𝑉 (𝑡) − 𝑣 𝑆 𝑉 (𝑡 − 𝜏) +𝑃 + 𝑣 𝑆[𝑉 (𝑡)] − (𝛼 + 𝛽)𝑉 . (𝑡). (7)
From Figure 1 and the above equation, the following dynamic equation of state of the thalamic reticular nucleus can be obtained: In the same way, we can obtain the dynamic equation of the discharge of the nucleus of the second corticothalamic circuit in Figure 1: Most of the parameters of the formula are from the papers [17,22].A small part of the parameter was changed within a normal physiological range.The detailed results of the parameters were shown in Tables 1 and 2. To further study the control mechanism of seizures, we applied electrical stimulation to the model TC1.The forms of electrical stimulation are HFS, 1:0 CRS and 3:2 CRS.Their pulse diagram was shown in Figure 2. The mathematical expression of HFS was described as follows [9,49,50]: If r < 0, we set r = 0, and the formula for 1:0 CRS is when L < 0, we set L = 0, and 3:2 CRS can be described as In the HFS formula, A represents the amplitude of the electrical stimulation and T and p0 represent the period and width of the pulse oscillation, respectively.The t and H in the electrical stimulation equation represent the time and heaviside function, respectively.In this paper, we gave stimulation time to the reticular nucleus of thalamus 10 s in the model, according to the electrical stimulation form of Figure 2. In the process of modeling electrical stimulation, we constructed two different square wave functions based on HFS pulse sequence, and then we obtained two forms of CRS by superimposing the two square wave functions with HFS, respectively.Both HFS and CRS are single-phase pulse stimulations, of which CRS is divided into 1:0 CRS and 3:2 CRS in Figure 2 due to the different distribution of pulse sequences.

Dynamic state and frequency analysis of neuron firing
The thalamus plays an important role in controlling seizures through GABAB synapses, especially in the pathway between TC and RE [51][52][53].To verify this effect of the thalamus, the state and frequency of epileptic area was analyzed in Figure 3 by a 2D ( , ) panel.It can be seen from Figure 3(a),(b) that when  was small, only saturation state appeared.When time increased, the discharge frequency of the cerebral cortex increased to the maximum and the saturation state was shown in Figure 3(c).If  continued to increase, SWD appeared, indicating a seizure.It can be seen from Figure 3(b) that the oscillation frequency of the SWD was 2-4 Hz.In addition, the frequencies of simple oscillations were also partially located at 2-4 Hz, but the difference from SWD was shown in Figures 3(d),(e).Simple oscillation had only one peak in one cycle and SWD had two peaks.When  increased to the maximum, there was a low firing state.As shown in Figure 3(f), the firing became a straight line with an oscillation frequency of zero Hz as time increased.From the above analysis and results in Figure 3, it can be known that the model in this paper can successfully simulate four dynamic states of the cerebral cortex.

Effect of electrical stimulation parameters and the coupling strength of TC1 to RE1 pathway on neuronal dynamics
The electrical stimulation parameters and the coupling strength  were used to explore the suppression of seizures.Figure 4(a),(b) shows the effect of HFS on neuronal dynamics and dominant frequency, respectively.Figure 4(c),(d) shows the effect of stimulation on neuronal dynamics and dominant frequency after applying 1:0 CRS on RE1.The effect of 3:2 CRS on the dynamics and dominant frequency of neurons is shown in Figure 4(e),(f).Figure 4 shows that the SWD areas of HFS, 1:0 CRS and 3:2 CRS account for 3.58, 7.48 and 17.46%, respectively.As the coupling strength and electrical stimulation intensity increased, the SWD area gradually became smaller.It was concluded that the inhibitory effect of HFS on seizures is better than CRS, and the therapeutic effect of 1:0 CRS is better than 3:2 CRS.Therefore, for a better therapeutic effect, we could increase the electrical stimulation intensity and coupling strength in a certain range.
Figure 5 explores the effect of electrical stimulation on the dynamics state of neurons in the ( ,  ) panel. Figure 5(a),(b) shows the dynamic states and frequency analysis after adding HFS.It can be seen from the figure that when the pulse width increased appropriately, the area of SWD area gradually decreases and accounts for 3.36% of the whole dynamic state.Figure 5(c),(d) shows that after applying 1:0 CRS to RE1, we found that its change trend is like that of HFS.In the process of increasing the pulse width, the area of SWD also decreases gradually.The SWD area with 1:0 CRS was obviously larger than the SWD area with HFS. Figure 5(e),(f) shows the effect of 3:2 CRS on RE1.Compared with the above two methods, 3:2 CRS has a relatively weak inhibition effect on SWD, and the SWD area accounts for 42.86%.6(a),(b) shows the suppression effect after adding HFS.When the frequency of electrical stimulation increases, the SWD area is greatly reduced, accounting for only 0.45% of the whole state.The influence of 1:0 CRS on the dynamic state and main frequency of neurons was shown in Figure 6(c),(d).When the stimulation frequency increased, the area of SWD began to decrease significantly, and finally accounted for 12.02% of the whole state area.The impact of 3:2 CRS on neuronal activity was shown in Figure 6(e),(f).Its trend of change was the same as that of HFS and 1:0 CRS, but its SWD area was larger, accounting for 34.01%.

Influence of electrical stimulation parameters on neuron dynamics
Figure 7 shows the inhibitive effect of applying various electrical stimulation to RE1 in the (F, A) panel.Figure 7(a),(b) shows the inhibitive effect of HFS on neuronal activity.When the electrical stimulation intensity continues to increase, the SWD area will gradually decrease and then disappear, and the increase of electrical stimulation frequency will also lead to the decrease of the SWD area and, finally, the percentage of the SWD area is 12.70%.Figure 7(c) shows the relationship between the electrical stimulation intensity and frequency of HFS and the number of SWD.From the figure, we can know that the increase of electrical stimulation intensity and frequency will greatly reduce the number of SWD. Figure 7(d),(e) shows the impact of 1:0 CRS on the state and frequency of neuron firing activity.Its influence on SWD area change is like the dynamic state under the influence of HFS, but the percentage of the SWD area is 25.85%.Figure 7(f) is about the relationship between stimulation intensity and the frequency of 1:0 CRS and SWD quantity.With the increase of electrical stimulation intensity and frequency, the SWD quantity begins to decrease gradually, but the decreasing trend of I: Saturate state II: SWD III: Simple oscillation state IV: Low firing state SWD quantity is not as fast as that under the influence of HFS. Figure 7(g),(h) respectively shows the effect of 3:2 CRS on neuron discharge activity.Its SWD area percentage is 45.58%, which is far weaker than the inhibition effect of the above two kinds of electrical stimulation on SWD. Figure 7(i) is the relationship between the parameter changes of 3:2 CRS and the number of SWD.We can find that the number of SWD will decrease slowly with the increase of electrical stimulation intensity and frequency.Figure 8 represents the effect of electrical stimulation on the dynamic states and firing frequency of neurons in the ( , A) panel.Figure 8(a),(b) shows the effect of HFS on the discharge activities of neurons.With the increase of electrical stimulation intensity and pulse width, the SWD area will decrease greatly, accounting for 9.07% of the whole state area.Figure 8(c) shows the influence of the intensity and frequency of HFS on the number of SWD.With the increase of these two parameters, the number of SWD will decrease greatly, which shows that increasing the intensity and pulse width of electrical stimulation has an obvious effect on inhibiting SWD.The impact of 1:0 CRS on the neuron firing state and frequency was shown in Figure 8(d),(e).The SWD area change was also like the change trend of the SWD area under the HFS.However, the SWD area of 1:0 CRS, a larger one, accounted for 17.01%.Figure 8(f) is the influence of electrical stimulation intensity and pulse width of 1:0 CRS on the number of SWD.The increase of electrical stimulation intensity and pulse width will lead to the decrease of the SWD number, but the inhibition effect of increasing electrical stimulation intensity on SWD is more obvious than that of increasing electrical stimulation pulse width.Figure 8(g),(h) shows the effect of 3:2 CRS on neuronal activity.The change trend of the entire discharge activity was like the previous two, but the SWD area of 3:2 CRS accounts for 29.48%.Figure 8(i) shows the influence of stimulation intensity and pulse width of 3:2 CRS on the SWD number.Although the change trend of the SWD number is like the other electrical stimulation, the decrease trend of the SWD number is slow.Figure 9 shows the effect of electrical stimulation on neuronal dynamics in the ( , F) panel.Figure 9(a),(b) are the effects of HFS on the firing state and the frequency of neurons.Under the action of HFS, with the increase of electrical stimulation frequency and pulse width, the SWD area only accounts for 0.23% of the whole dynamic state.Figure 9(c) is the influence of pulse width and stimulation frequency of HFS on the number of SWD, from which we can know that with the increase of stimulation frequency and pulse width, the number of SWD will obviously decrease.Figure 9(d),(e) shows the impact of 1:0 CRS on neuron firing.The distribution of the SWD area was like Figure 9(a),(b), but the SWD area became larger, which accounts for 32.65%.The effect of 3:2 CRS on neuronal activity was shown in Figure 9(e),(f).Compared with HFS and 1:0 CRS, the SWD area of 3:2 CRS, the largest one accounts for 90.02%.Figure 9(i) is the influence of the pulse width and stimulation frequency of 3:2 CRS on the number of SWD, and its inhibition on the number of SWD is the worst among the three kinds of electrical stimulation.Figure 10 is a summary of the percentage of SWD under the above-mentioned different parameter combinations, from which we can know that HFS is obviously superior to CRS in inhibiting SWD.In addition, among the combination of coupling intensity and electrical stimulation parameters, the most obvious inhibitory effect of HFS on absence seizures is the combination of coupling intensity and electrical stimulation frequency, and the most obvious inhibitory effect of 1:0 CRS on absence seizures is the combination of coupling intensity and pulse width.The best inhibitory effect of 3:2 CRS on absence seizures is the combination of coupling intensity and electrical stimulation intensity.Among the combinations of electrical stimulation parameters, the combination of pulse width and electrical stimulation frequency showed that HFS had the best inhibitory effect on absence seizures.Under the combination of pulse width and electrical stimulation intensity, 1:0 CRS and 3:2 CRS had the most obvious inhibitory effect on absence seizures.

Conclusions
Based on the original BGCT model, we constructed a reduced BGCT-PPN model consisting of the PPN, two corticothalamic circuits and basal ganglia.We simplified the basal ganglia into a 3I:5O modulator for reducing the complexity of the model and iteration between neurons.We applied three forms of electrical stimulation to RE1 in the first corticothalamic circuit of the model, and then explored how the SWD region is inhibited by varying the parameters of electrical stimulation.By adjusting the parameters of electrical stimulation and the coupling strength of the model, we found that appropriately increasing the stimulation strength of electrical and the coupling strength of the model can effectively inhibit the SWD area.In addition, the increasing frequency and pulse width of electrical stimulation can also inhibit the SWD region.
The numerical results show that when the electrical stimulation pulse width is 1-3 ms, the electrical stimulation intensity is 0-200 mV, the stimulation frequency is 100-185 Hz and the SWD area will decrease with the increase of electrical stimulation parameters.Second, in different parameter panels, HFS has the most obvious inhibitory effect on SWD.1:0 CRS also has a certain inhibitory effect on the SWD area, and the inhibitory effect is better than 3:2 CRS.These results are of great significance to the application of the mean field model in theoretical experiments and the exploration of how to regulate epilepsy through electrical stimulation parameters.It is hoped that our numerical simulation results will be helpful to clinical experiments.However, the above results are based on theoretical research and may have some limitations, so our next research plan is to test the results of theoretical research through experiments.

Figure 1 .
Figure 1.Reduced BGCT-PPN model.The model consists of four parts, i.e., the first corticothalamic circuit (PY1, IN1, TC1 and RE1), the second corticothalamic circuit (PY2, IN2, TC2 and RE2), 3I:5O modulator and PPN.Red arrow indicates the excitatory stimulation mediated by glutamate as the receptor.Blue solid line and dashed line indicate inhibitory stimulation mediated by GABAA and GABAB as receptors, respectively.

Figure 2 .
Figure 2. Different forms of electrical stimulation.(a) HFS.(b) 1:0 CRS.(c) 3:2 CRS.The pulse width of electrical stimulation can be set as 1-3 ms, the frequency of electrical stimulation can be set between 100-185 Hz and the intensity of electrical stimulation is 0-200 mV.

Figure 4 .
Figure 4. Influence of electrical stimulation on neuron dynamics in ( , A) panel.(a) and (b) represent the influence of HFS on the dynamic state and oscillation frequency of neurons.(c) and (d) represent the influence of 1:0 CRS on the dynamic state and oscillation frequency of neurons.(e) and (f) represent the influence of 3:2 CRS on the dynamic state and oscillation frequency of neurons.During the whole electrical stimulation, we set  = 1 ms,  = 120 s and τ = 80 ms.

Figure 5 .
Figure 5. Influence of electrical stimulation on neuron dynamics in the ( ,  ) panel.(a) and (b) represent the influence of HFS on the dynamic state and oscillation frequency of neurons.(c) and (d) represent the influence of 1:0 CRS on the dynamic state and oscillation frequency of neurons.(e) and (f) represent the influence of 3:2 CRS on the dynamic state and oscillation frequency of neurons.When electrical stimulation was applied, the period and intensity of stimulation were  = 120 s and A = 120 mV.

Figure
Figure6explores the effect of coupling strength and electrical pulse stimulation frequency on the SWD area.Figure6(a),(b) shows the suppression effect after adding HFS.When the frequency of electrical stimulation increases, the SWD area is greatly reduced, accounting for only 0.45% of the whole state.The influence of 1:0 CRS on the dynamic state and main frequency of neurons was shown in Figure6(c),(d).When the stimulation frequency increased, the area of SWD began to decrease significantly, and finally accounted for 12.02% of the whole state area.The impact of 3:2 CRS on neuronal activity was shown in Figure6(e),(f).Its trend of change was the same as that of HFS and 1:0 CRS, but its SWD area was larger, accounting for 34.01%.

Figure 6 .
Figure 6.Influence of electric electrical stimulation on neuron dynamics in the ( , ) panel.(a) and (b) represent the influence of HFS on the dynamic state and oscillation frequency of neurons.(c) and (d) represent the influence of 1:0 CRS on the dynamic state and oscillation frequency of neurons.(e) and (f) represent the influence of 3:2 CRS on the dynamic state and oscillation frequency of neurons.The intensity and pulse width of electrical stimulation are set to A = 120 mV and  = 1 ms, respectively.

Figure 7 .
Figure 7. Influence of electrical stimulation on neuron dynamics in the (F, A) panel.(a) and (b) represent the influence of HFS on the dynamic state and oscillation frequency of neurons.(c) is the relationship between the stimulation intensity and frequency of HFS and the number of SWD.(d) and (e) represent the influence of 1:0 CRS on the dynamic state and oscillation frequency of neurons.(f) is the relationship between the stimulation intensity and frequency of 1:0 CRS and the number of SWD.(g) and (h) represent the influence of 3:2 CRS on the dynamic state and oscillation frequency of neurons.(i) is the relationship between the stimulation intensity and frequency of 3:2 CRS and the number of SWD.The pulse width of the electrical stimulation is  = 1 ms.

I:
Saturate state II: SWD III: Simple oscillation state IV: Low firing state

Figure 8 .Figure 9 .
Figure 8. Influence of electrical stimulation on neuron dynamics in the ( , A) panel.(a) and (b) represent the influence of HFS on the dynamic states and oscillation frequency of neurons.(c) is the relationship between the stimulation intensity and frequency of HFS and the number of SWD.(d) and (e) represent the influence of 1:0 CRS on the dynamic state and oscillation frequency of neurons.(f) is the relationship between the stimulation intensity and frequency of 1:0 CRS and the number of SWD.(g) and (h) represent the influence of 3:2 CRS on the dynamic state and oscillation frequency of neurons.(i) is the relationship between the stimulation intensity and frequency of 3:2 CRS and the number of SWD.The period of electrical stimulation is  = 120 s.

Figure 10 .
Figure 10.Effect of electrical stimulation on SWD percentage under different conditions.

Table 2 .
Strength of coupling.