Hypersynchronous (HYP) seizure onset is one of the frequently observed seizure-onset patterns in temporal lobe epileptic animals and patients, often accompanied by hippocampal sclerosis. However, the exact mechanisms and ion dynamics of the transition to HYP seizures remain unclear. Transcranial magneto-acoustic stimulation (TMAS) has recently been proposed as a novel non-invasive brain therapy method to modulate neurological disorders. Therefore, we propose a biophysical computational hippocampal network model to explore the evolution of HYP seizure caused by changes in crucial physiological parameters and design an effective TMAS strategy to modulate HYP seizure onset. We find that the cooperative effects of abnormal glial uptake strength of potassium and excessive bath potassium concentration could produce multiple discharge patterns and result in transitions from the normal state to the HYP seizure state and ultimately to the depolarization block state. Moreover, we find that the pyramidal neuron and the PV+ interneuron in HYP seizure-onset state exhibit saddle-node-on-invariant-circle/saddle homoclinic (SH) and saddle-node/SH at onset/offset bifurcation pairs, respectively. Furthermore, the response of neuronal activities to TMAS of different ultrasonic waveforms revealed that lower sine wave stimulation can increase the latency of HYP seizures and even completely suppress seizures. More importantly, we propose an ultrasonic parameter area that not only effectively regulates epileptic rhythms but also is within the safety limits of ultrasound neuromodulation therapy. Our results may offer a more comprehensive understanding of the mechanisms of HYP seizure and provide a theoretical basis for the application of TMAS in treating specific types of seizures.

HYP seizure onset is characterized as low-frequency, high-amplitude periodic spikes and is thought to result from an imbalance between inhibition and excitation. As a novel neuromodulation, TMAS has attracted widespread attention for its non-invasive regulation of neural activity without tissue damage. However, due to neurons’ intrinsic characteristics and interactions, different stimulation parameters may lead to transitions in neuronal states. Computational models can reproduce phenomena observed from physiological experiments, reveal underlying seizure mechanisms and neuromodulation from the perspective of dynamics, and guide future treatment for epileptic disorders. We complement our network model with variables for the ion concentration to reveal the dynamical mechanisms of the evolution of HYP seizure driven by excess extracellular potassium. Furthermore, we design an effective TMAS modulation strategy to control HYP seizures in the human hippocampus.

Temporal lobe epilepsy in humans is caused by abnormal neural discharges in the focus lesion of the brain, which is often closely associated with the hippocampus. According to electroencephalographic recordings in temporal lobe epileptic patients and animal models, hypersynchronous (HYP) and low-voltage fast (LVF) are identified as the most frequently observed focal seizure-onset patterns, and the LVF onset pattern often follows HYP onset.1–3 High-amplitude periodic spikes that occur at a frequency below 2 Hz are considered electrophysiological characteristics of HYP onset seizure.4 Previous studies and experiments have shown that HYP onset seizure is caused by progressive exhaustion of inhibition and unrestrained enhancement of excitation.5–7 In recent years, it has been documented that the HYP onset pattern can be triggered by optogenetic activation of pyramidal neurons in the vitro 4-aminopyridine model and is mostly accompanied by fast ripples (250–500 Hz).8,9 With the increasing research on neurological diseases, the role of ions in temporal lobe epilepsy is receiving more attention. From the perspective of ion dynamics, there are several possibilities for its contribution to epilepsy: for instance, extracellular potassium accumulation due to blockade of Kv1 potassium channels10 or impaired astrocytic potassium regulation system;11 low magnesium;5 excess intracellular chloride due to abnormal activities of N a + / K + / 2 C l cotransporter-1 (NKCC1) and K + / C l cotransporter-2 (KCC2) cotransporter.12 

There are numerous mathematical models proposed to study the effect of ion dynamics on the transitions from regular neuronal activity to epilepsy. Ullah presented a model containing the following compartments: neurons, glia, and extracellular space to investigate the effects of high potassium, glial cell buffering,13 and cell swell14 on epilepsy. Furthermore, Wei et al.15 used a similar model to reproduce LVF seizure-like events caused primarily by high potassium due to hypoxia and validated the vitro experiments in the Ammon's horn 1 (CA1) region of rat hippocampal slices.16,17 In addition, Gentiletti constructed a small-scale biophysical network to explore neuronal state transitions between the various LVF seizure phases and corresponding changes in Na, K, and Cl ion concentrations.18 Yao used a single-neuron model to study the dynamic mechanisms of seizure and spreading depression based on oxygen and ion concentrations.19 A minimal model called “Epileptor-2” was proposed to indicate that ictal discharges were triggered by potassium accumulation and eliminated by sodium–potassium pump activation.20 However, to the best of our knowledge, in contrast to numerous computational models of LVF onset seizure, only a few were proposed to study mechanisms of HYP seizure onset at the ion concentration level. In particular, Buchin et al.21 presented a hippocampus subiculum network model to study the transition from gamma oscillations to HYP seizure due to the reduced efficacy of the KCC2 cotransporter, but this research mainly focused on transitions between normal and epileptic state rather than on the complete dynamic state transitions.

The existing treatment methods for epilepsy mainly include medication, surgical treatment, and neuromodulation. However, medication lacks brain target specificity, and surgical treatment by implanting electrodes can cause irreversible damage to the brain. Compared with traditional brain stimulation such as deep brain stimulation, optogenetic stimulation, and transcranial magnetic stimulation, transcranial magneto-acoustic stimulation (TMAS) has the advantage of non-invasiveness, high spatial resolution (<3 mm), and high penetration depth.22,23 Therefore, TMAS is gradually playing an important role in treating neurological diseases such as Alzheimer's disease24 and Parkinson's disease.25 As a novel non-invasive neuromodulation, TMAS can generate a stimulating electric current in nerve tissues by Lorentz forces on moving ions produced by the ultrasonic wave in the presence of a magnetic field.26 Previous physiological experiments have shown that low-intensity pulsed ultrasound stimulation can effectively regulate neuronal potassium currents27 and show great potential in the suppression and treatment of epilepsy.28,29 Moreover, Yuan et al.30,31 demonstrated that TMAS could affect neuronal firing rhythms and ion concentrations, which provided a theoretical basis for TMAS in the modulation of HYP seizures. However, the crucial ultrasound stimulation parameters and therapeutic effects of TMAS on epilepsy need further research.

Thus, in this paper, we construct a biophysical network model containing pyramidal neurons, a parvalbumin-positive (PV+) interneuron, and astrocytes to investigate the inherent ionic mechanisms of the complete evolution of HYP seizure and the dynamical properties of bifurcation types at seizure onset/offset. Further, we explore precise control of the HYP epileptic hippocampal network by designing an effective TMAS regulation strategy. Our results may provide a theoretical basis and guidance for diagnosing and treating HYP seizures.

Our biophysical computational hippocampal Ammon's horn 3 (CA3) network model comprises four pyramidal neurons and a PV+ interneuron. Each neuron is embedded with extracellular space and an astrocyte. The extracellular space of each neuron is enclosed in a bath solution that represents blood vessels in vivo. As shown in Fig. 1(a), pyramidal neurons receive excitation from other pyramidal neurons via glutamatergic receptors and inhibition from a PV+ interneuron via gamma-aminobutyric acid type A ( GAB A A ) receptor, whereas the PV+ interneuron only receives glutamate from pyramidal neurons. It is worth mentioning that our network connectivity is all-to-all, so there are synaptic connections between the pyramidal neurons located on the diagonal of Fig. 1(a) (PY1–PY3 and PY2–PY4). For instance, the four lines on the main diagonal that involve PY3 in Fig. 1(a) represent the synaptic connections of PY3 to PY1, PY3 to PV+, PV+ to PY3, and PY1 to PY3.

FIG. 1.

(a) Schematic diagram of the biophysical hippocampal network model. The model consists of four pyramidal neurons (orange) and a PV+ interneuron (green) linked by glutamatergic and GABAergic synapses. Each neuron is embedded with extracellular space and an astrocyte(blue). The extracellular space of each neuron is enclosed in a bath solution (within a blue rounded rectangle). (b) Correlation diagram reflecting the dynamic interaction of the five major ions in our pyramidal neuron model. The ionic dynamics are mainly modulated by ion channels, transporters (KCC2, NKCC1, and Na+/K+ pump) and synaptic receptors (GABAA receptor and glutamate receptor). The different colored arrows represent the direction of the corresponding ion flux.

FIG. 1.

(a) Schematic diagram of the biophysical hippocampal network model. The model consists of four pyramidal neurons (orange) and a PV+ interneuron (green) linked by glutamatergic and GABAergic synapses. Each neuron is embedded with extracellular space and an astrocyte(blue). The extracellular space of each neuron is enclosed in a bath solution (within a blue rounded rectangle). (b) Correlation diagram reflecting the dynamic interaction of the five major ions in our pyramidal neuron model. The ionic dynamics are mainly modulated by ion channels, transporters (KCC2, NKCC1, and Na+/K+ pump) and synaptic receptors (GABAA receptor and glutamate receptor). The different colored arrows represent the direction of the corresponding ion flux.

Close modal
As shown in Eq. (1), the pyramidal neuron is modeled with the Hodgkin–Huxley formulation based on Olufsen.32 PV+ basket neurons are thought to account for 60% of all PV-expressing interneurons, which are the most significant number of neurons in all PV-expressing interneurons.33 Thus, we used PV+ basket interneurons modeled by Wang and Buzsáki,34 as shown in Eq. (2). More details on neuron modeling are available in Refs. 31 and 33,
C d V PY d t = I Na I Nap I K I AHP + I L I syn + I TMAS J pump / γ PY ,
(1)
C d V PV + d t = I Na I K + I L I syn + I TMAS J pump / γ PV + ,
(2)
where C, V, and t denote capacitance density, neuronal membrane potential, and time, respectively. The units of the above variables are μ F / c m 2, mV, and ms. The pyramidal neuron contains sodium voltage-gated current I Na, persistent sodium voltage-gated current I Nap, potassium voltage-gated current I K, calcium-activated potassium or after-hyperpolarization (AHP) current I AHP, leak current I L, synaptic current I syn, and sodium–potassium pump current J pump / γ, whereas the PV+ interneuron does not contain I Nap and I AHP. γ represents the conversion factor from the current density into the ion concentration rate of change. For the pyramidal neuron, we set γ PY = 0.0445 mM / ( C cm ). For the PV+ interneuron, we set γ PV + = 0.0286 mM / ( C cm ).35 The voltage-gated ion channel equations are shown in Eq. (3),36,
I Na = g Na m 3 h ( V E Na ) , I Nap = g Nap m 3 ( V E Na ) , I K = g K n 4 ( V E K ) , I AHP = g AHP [ Ca ] i 1 + [ Ca ] i ( V E K ) , I L = g NaL ( V E Na ) g KL ( V E K ) g ClL ( V E Cl ) ,
(3)
where g Na, g Nap, g K, g AHP, g NaL, g KL, and g ClL represent the maximum conductance corresponding to the voltage-gated ion channel currents. [ Ca ] i denotes intracellular calcium concentration. E Na, E K, and E Cl indicate the reversal potential of sodium, potassium, and chloride, respectively. For the pyramidal neuron, we set g Na = 30 mS / c m 2, g Nap = 1 mS / c m 2, g K = 25 mS / c m 2, and g AHP = 1.5 mS / c m 2. For the PV+ interneuron, we set g Na = 30 mS / c m 2 and g K = 9 mS / c m 2. We assume that the leakage current values of both types of neurons are the same and we set g NaL = 0.175 mS / c m 2, g KL = 0.05 mS / c m 2, and g ClL = 0.05 mS / c m 2. The reversal potential E x (x = Na, K, and Cl) can be calculated from the dynamic intracellular ( [ x ] i ) and extracellular ( [ x ] o ) ionic concentrations by the Nernst equation. In particular, GABA receptors are approximately four times more permeable to chloride than bicarbonate.37 Thus, the effect of chloride and bicarbonate ions on the GABAergic reversal potential E GABA in our model is set to four to one,
E Na = 26.64 ln [ Na ] o [ Na ] i , E K = 26.64 ln [ K ] o [ K ] i , E Cl = 26.64 ln [ Cl ] i [ Cl ] o , E GABA = 26.64 ln 4 [ Cl ] i + [ HC O 3 ] i 4 [ Cl ] o + [ HC O 3 ] o ,
(4)
where [ HC O 3 ] i and [ HC O 3 ] o represent intracellular and extracellular bicarbonate concentrations, respectively. We set [ HC O 3 ] i = 16 mM and [ HC O 3 ] o = 26 mM.
In our model, synaptic currents I syn consisting of inhibitory GABA synaptic current I GABA, excitatory alpha-amino-3-hydroxy-5-methyl-4-isoxazole propionic acid (AMPA) current I AMPA and excitatory N-methyl-D-aspartic acid (NMDA) current I NMDA are described as follows:
I GABA = g GABA , ij s GABA , i ( E GABA V j ) , I AMPA = g AMPA , ij s AMPA , i ( E AMPA V j ) , I NMDA = g NMDA , ij s NMDA , i B ( V j ) ( E NMDA V j ) , B ( V ) = 1 / ( 1 + [ Mg ] o exp ( 0.062 V ) / 3.57 ) , d s x , i d t = 1 + tanh ( V / 4 ) 2 × 1 s x τ R s x τ D ,
(5)
where g x , i j represents the maximum synaptic conductance, where x = GABA, AMPA, or NMDA and the synaptic current is from neuron i to neuron j. s x , i (x = GABA, AMPA, and NMDA) represents the fraction of synaptic gating variable related to the presynaptic neuron i. [ Mg ] o and B ( V j ) denote extracellular magnesium concentration and postsynaptic membrane potential dependence fraction on magnesium, respectively. We set [ Mg ] o = 2 mM. τ R and τ D represent the synaptic rise and decay time constant, respectively. The AMPAergic reversal potential E AMPA and the NMDAergic reversal potential E NMDA are set to zero. For an AMPA-type synapse, we set τ R = 0.1, τ D = 3, g PY PV + = 0.01 mS / c m 2, and g PY PY = 0.04 mS / c m 2. For an NMDA-type synapse, we set τ R = 13, τ D = 150, g PY PV + = 0.3 mS / c m 2, and g PY PY = 3 mS / c m 2. For a GABA-type synapse, we set τ R = 0.3, τ D = 40, and g PV + PY = 0.15 mS / c m 2.
The potassium ( I K , GLU ) and sodium ( I Na , GLU ) currents through the glutamatergic channels depend on voltage with the fractions of the AMPAergic and NMDAergic conductance. The current equations are modified after Chizhov et al. shown in Eq. (6),38 
I K , GLU = 0.2 ( g AMPA s AMPA + g NMDA s NMDA B ( V ) ) ( V E K ) , I Na , GLU = 0.4 ( g AMPA s AMPA + g NMDA s NMDA B ( V ) ) ( V E Na ) .
(6)

We complemented the Hodgkin–Huxley formulation model with the ion concentration using differential equations based on previous studies.39 The specific details are shown in Fig. 1(b). The ionic dynamics of N a +, K +, C l , C a 2 +, and HCO 3 are regulated by ion channels, cotransporters (KCC2 and NKCC1), and synaptic channels. The intracellular and extracellular ion concentrations of each ion type are continuously updated as slow variables.

Here, d [ K ] o d t is the rate of change of extracellular potassium concentration determined by total potassium-dependent channel currents I K , total ( I K , I KL , I AHP , I K , GLU ), N a + / K + pump J pump, K + / C l cotransporter J kcc 2, N a + / K + / 2 C l cotransporter J nkcc 1, extracellular potassium diffusion from a bath solution J diff, glial uptake by astrocyte surrounding the neurons J glia, glia N a + / K + pump J gliapump, and potassium lateral diffusion between neurons J lateral. d [ Na ] i d t is the rate of change of intracellular sodium concentration determined by total sodium-dependent channel currents I Na , total ( I Na , I Nap , I Nal , I Na , GLU ), J pump, and J nkcc 1. d [ Cl ] i d t is the rate of change of intracellular chlorine concentration determined by I ClL, I GABA, J kcc 2, and J nkcc 1. [ K ] i, [ Na ] o, and [ Cl ] o are modeled by linear equations to simplify calculations,39,
d [ K ] o d t = 1 τ ( γ β I K , total 2 β J pump + β J kcc 2 + β J nkcc 1 J diff J glia 2 J gliapump J lateral ) , d [ N a ] i d t = 1 τ ( γ I Na , total 3 J pump J nkcc 1 ) , d [ C l ] i d t = 1 τ ( γ ( I ClL + I GABA ) J kcc 2 2 J nkcc 1 ) , [ K ] i = 140 + ( 18 [ Na ] i ) ( 6 [ Cl ] i ) , [ Na ] o = 144 β ( [ Na ] i 18 ) , [ Cl ] o = 130 β ( [ Cl ] i 6 ) ,
(7)
where τ equals to 1000 and converts seconds to milliseconds. β equals to 7 and represents the ratio of intracellular to extracellular volume. The equations for the change in relevant ion currents are shown as follows:
J pump = ρ pump 1 + exp ( 5.5 [ K ] o ) × 1 1 + exp ( ( 25 [ Na ] i ) / 3 ) , J gliapump = 1 3 × ρ pump 1 + exp ( 5.5 [ K ] o ) × 1 1 + exp ( ( 25 [ Na ] gi ) / 3 ) , J glia = G glia 1 + exp ( ( 18 [ K ] o ) / 2.5 ) , J nkcc 1 = U nkcc 1 26.64 × ( 2 E Cl E Na E K ) 1 + exp ( 16 [ K ] o ) , J kcc 2 = U kcc 2 26.64 × ( E Cl E K ) , J lateral = D E d x 2 ( [ K ] o [ K ] o else ) , J diff = ε ( [ K ] o K bath ) ,
(8)
where ρ pump, G glia, U nkcc 1, and U kcc 2 represent the maximum N a + / K + pump strength, the maximal glial uptake strength of potassium, the maximum NKCC1, and the maximum KCC2 cotransporter strengths, respectively. Potassium lateral diffusion J lateral depends mainly on the difference in potassium concentration between neurons. We assume that the potassium lateral diffusion factor D E = 2.5 × 10 5 c m 2 / s is the same as the potassium diffusion factor in water.40 The mean distance between two neurons d x = 50 μ m is measured in the hippocampus.12  ε and K bath denote the maximal potassium diffusion rate and potassium concentration from a bath solution, respectively. We set parameter ρ pump = 1.25 mM / s, U nkcc 1 = 0.1 mM / s, and U kcc 2 = 0.3 mM / s.
Intracellular calcium concentration [ Ca ] i that primarily affects AHP current is modeled as follows:
d [ Ca ] i d t = 0.002 g Ca ( V E Ca ) 1 + exp ( ( V + 25 ) / 2.5 ) [ Ca ] i 80 ,
(9)
where g Ca and E Ca represent the maximum calcium conductance and reversal potential of calcium, respectively. We set g Ca = 1 mS / c m 2 and E Ca = 120 mV.
Through low-frequency ultrasound, TMAS modulates seizures by producing induced currents in neurons in a static magnetic field.31 Schematic diagram of TMAS is shown in Fig. 2(a). Therefore, based on the previous works,30,41 the TMAS-induced electric current density I TMAS can be added directly to Hodgkin–Huxley models expressed as
I TMAS = σ B 2 U I ρ c 0 ( sin ( 2 π f t ) + b ) ,
(10)
where σ, ρ, c 0, and f represent nerve tissue conductance, tissue density, ultrasound speed, and ultrasound frequency, respectively. Magnetic field intensity B and ultrasound intensity U I will be further analyzed as crucial modulation parameters. The bias parameter b determines the ultrasonic waveforms, which we describe in detail in Fig. 2(b). Here, we set σ = 0.5 Simens / m, ρ = 1120 Kg / m 3, c 0 = 1540 m / s, and f = 0.2 MHz.
FIG. 2.

(a) Schematic diagram of TMAS. TMAS can excite or suppress the neurons through an electric current when different ultrasonic waveforms generated by the transducer stimulate neurons in the presence of a magnetic field. (b) Ultrasonic waveforms with upper sine wave (top), sine wave (center), and lower sine wave (bottom) when the bias ( b ) parameter is 1, 0, and −1, respectively.

FIG. 2.

(a) Schematic diagram of TMAS. TMAS can excite or suppress the neurons through an electric current when different ultrasonic waveforms generated by the transducer stimulate neurons in the presence of a magnetic field. (b) Ultrasonic waveforms with upper sine wave (top), sine wave (center), and lower sine wave (bottom) when the bias ( b ) parameter is 1, 0, and −1, respectively.

Close modal
According to Hendee's theory,42 we can convert ultrasonic intensity to ultrasonic pressure at the lesion using the formula shown in Eq. (11),
P L = 2 U I ρ C 0 ,
(11)
where P L represents the ultrasound pressure at the lesion in MPa.
The mechanical index (MI) indicates the potential biomechanical danger due to ultrasound cavitation. The thermal index (RI) reflects the temperature rise in neural tissue due to the absorption of ultrasound energy. With reference to standard regulation of the United States Food and Drug Administration (FDA),43 combined with our model, the MI and RI formulas are simplified as shown in Eqs. (12) and (13), respectively,
M I = P L C MI × f ,
(12)
T I = U I × S × f C TIS ,
(13)
where S is the ultrasound spot area in cm2, C MI = 1 MPa MH z 1 / 2, and C TIS = 210 mW MHz.

All simulations are carried out in MATLAB, version R2021b (https://ww2.mathworks.cn). Differential equations are solved with the fourth-order Runge–Kutta method. Taking the high-frequency characteristics of ultrasound into consideration, we set a fixed time step t = 0.0005 ms. Taking the time step and ultrasound frequency into Eq. (10), we obtain the sine function with a time step 2 π f t = 0.2 π. Ten time steps are used to simulate a complete sine cycle, which ensures the precision of TMAS. The Parallel Computing Toolbox in MATLAB is used to speed up computation when we solve equations in parameter space. In addition, the initial conditions for membrane potentials and ion concentrations are determined based on the normal distribution with the following means: V PY = 67 mV, V PV + = 65 mV, [ K ] o = 4 mM, [ Na ] i = 18 mM, [ Cl ] i = 6 mM, and variance at 5% of the corresponding mean.

In this section, we focus on the complete state transitions and the corresponding ionic concentration dynamics due to changes in bath potassium concentration ( K bath ) and the glial uptake strength of potassium ( G glia ) parameters which may be associated with HYP seizures observed in hippocampal slices in high K + medium.44 Alterations in physiologic parameters K bath and G glia cause rapid changes in [ K ] o within pyramidal neurons and PV+ interneuron, accompanied by subsequent changes in the other ion concentrations ( [ Na ] i , [ Cl ] i , [ Ca ] i ), and ultimately induction of state transitions. As shown in Fig. 3, two-parameter space diagrams for parameters K bath and G glia are elaborated to study the quantitative effects on the network's state transitions. There are six states distinguished by discharge amplitude and interspike interval (ISI) including the resting state (1), fast-spiking state (2), synchronous bursting state (3), HYP seizure-onset state (4), sustained ictal activity state (5), and sustained depolarization block state (6). We can see that elevating K bath can gradually cause neurons to transition from normal discharges (resting state or fast-spiking state) to HYP seizures (synchronous bursting state, HYP seizure-onset state or sustained ictal activity) and eventually to the depolarization block state. However, increasing G glia improves the glial uptake strength of potassium and effectively reduces [ K ] o, converting HYP seizure onset into synchronous bursting discharges or fast-spiking state.

FIG. 3.

Two-parameter space diagrams for potassium from bath solution ( K bath ) and glial uptake strength of potassium ( G glia ). The six states are distinguished by different colors and represent the resting state (1), fast-spiking state (2), synchronous bursting state (3), HYP seizure-onset state (4), sustained ictal activity state (5), and depolarization block state (6), respectively.

FIG. 3.

Two-parameter space diagrams for potassium from bath solution ( K bath ) and glial uptake strength of potassium ( G glia ). The six states are distinguished by different colors and represent the resting state (1), fast-spiking state (2), synchronous bursting state (3), HYP seizure-onset state (4), sustained ictal activity state (5), and depolarization block state (6), respectively.

Close modal

To further understand the role of K bath in the neural discharge patterns and ionic concentration dynamics in the network model, we gradually increase K bath and fix Gglia = 10 mM/s. As a result of GABA inhibitory signaling from the interneuron, all four pyramidal neurons have essentially the same discharge pattern after entering the steady state. Therefore, we only present one of the pyramidal neurons to illustrate the differences between the two types of neurons in different states. Figure 4 shows the time series of the membrane potentials and ion concentrations in two classes of neurons vs different parameters K bath. As shown in Fig. 4(a1), all neurons remain at a resting state and ion concentrations remain stable when Kbath = 4 mM/s. As K bath increases, slight increases in [ K ] o make the neurons more excitable, resulting in the network transition to a fast-spiking state [Fig. 4(a2)]. According to the first and third panels in Fig. 4(a2), it can be clearly seen that the PV+ interneuron discharges faster than the pyramidal neuron. When K bath further increases, the pyramidal neuron and PV+ interneuron appear in a synchronized bursting state [Fig. 4(a3)], and their ion concentrations begin to fluctuate dramatically within a certain range. When K bath approximately approaches 16 mM, the pyramidal neuron remains bursting, while the PV+ interneuron partially enters a phase of depolarization block [Fig. 4(a4)]. This phenomenon is defined as the HYP seizure-onset state, which is similar to physiological experiments reported by Zhang et al.6 and Cammarota et al.45 The phase of depolarization block in PV+ interneuron is mediated by the inactivation of sodium channels due to extracellular potassium accumulation, which leads to a decrease in inhibitory GABA currents and more excitation of pyramidal neurons. Further increasing K bath, the neural network is considered to enter a sustained ictal activity [Fig. 4(a5)]. In this state, the membrane potentials and ion concentrations of the pyramidal neuron enter a state of high-frequency slight oscillation. However, the interneuron remains in the sustained depolarization block and its extracellular potassium ion concentration remains constant at a high value. Eventually, all neurons enter a complete depolarization block state when K bath reaches a higher level [Fig. 4(a6)].

FIG. 4.

Time series of pyramidal neuron membrane potentials ( V PY ) and ion concentrations ( Io n PY ) AND PV+ interneuron membrane potentials ( V PV ) and ion concentrations ( Io n PV ) corresponding to the six firing states in Fig. 3: (a1) rest state, K bath = 4 mM; (a2) fast-spiking state, K bath = 10 mM; (a3) synchronous bursting state, K bath = 15 mM; (a4) HYP seizure-onset state, K bath = 16 mM; (a5) sustained ictal activity, K bath = 34 mM; and (a6) sustained depolarization block state, K bath = 40 mM. Dynamic ion concentrations in both types of neurons contain intracellular sodium ( N a i ), extracellular potassium ( K o ), and intracellular chloride ( C l i ).

FIG. 4.

Time series of pyramidal neuron membrane potentials ( V PY ) and ion concentrations ( Io n PY ) AND PV+ interneuron membrane potentials ( V PV ) and ion concentrations ( Io n PV ) corresponding to the six firing states in Fig. 3: (a1) rest state, K bath = 4 mM; (a2) fast-spiking state, K bath = 10 mM; (a3) synchronous bursting state, K bath = 15 mM; (a4) HYP seizure-onset state, K bath = 16 mM; (a5) sustained ictal activity, K bath = 34 mM; and (a6) sustained depolarization block state, K bath = 40 mM. Dynamic ion concentrations in both types of neurons contain intracellular sodium ( N a i ), extracellular potassium ( K o ), and intracellular chloride ( C l i ).

Close modal

We further analyze the HYP seizure mechanism from the perspective of nonlinear dynamics combined with ion concentration dynamics. Complex oscillations such as bursting are triggered by fast/slow systems activity, where the slow subsystem drives the fast subsystem to transition between different states.46 In our model, the neuronal membrane potential timescale is milliseconds as a variable in the fast subsystem, while the ion concentration timescale is seconds as a variable in the slow subsystem. Figure 5 illustrates the evolution of the three-dimensional phase space of the two types of neurons in different firing patterns. From a dynamical perspective, neuronal resting and fast-spiking states are monostable. They correspond to a stable fixed point [Fig. 5(a1)] and a limit cycle [Fig. 5(a2)] in the phase space, respectively. The transition from the resting state to the fast-spiking state involves a saddle-node-on-invariant-circle (SNIC) bifurcation. After entering HYP seizures, the system switches between the quiescent and oscillatory behaviors. Jirsa et al. proposed a dynamical taxonomy of seizure-like events that systematically classified epileptic seizures into 16 different onset/offset bifurcation pairs.47 Saggio et al. provided a dynamical framework that can easily identify the features of oscillation onset/offset bifurcation pair by neuronal discharge amplitude, ISI, baseline jump, and flow changes in the phase space.48,49 As shown in Fig. 5(a3), the phase space of the pyramidal neuron starts at a stable fixed point corresponding to a quiescent state. The reduction of slow-variable [ K ] o drive trajectory to slide along the curve with a black arrow, where the black arrow indicates the direction of movement. When [ K ] o approximately approaches 13 mM, the potassium gate variable is activated and the trajectory transits to the stable limit cycle corresponding to the oscillatory state by crossing an SNIC onset bifurcation. Then, [ K ] o begins to gradually increase and the trajectory displays the large amplitude oscillations forming a tube-like structure. As [ K ] o approximately approaches 15 mM, the system reaches an offset bifurcation and the oscillation ends. At this time, the trajectory jumps back to a stable fixed point by crossing an SH offset bifurcation. According to the dynamical taxonomy of Izhikevich, this would represent a “circle/homoclinic” bursting.46 However, the phase space of the PV+ interneuron seems to be different. The trajectory starts from a stable point and undergoes small oscillations triggered by the excitation of pyramidal neurons. Then, the trajectory jumps up to the stable limit cycle by crossing an SN onset bifurcation. This leads to a baseline jump in the amplitude of PV+ interneuron from a quiescent state to an oscillatory state [third panel in Fig. 4(a3)]. After crossing the onset bifurcation, the trajectory exists in the form of limit cycles. Eventually, the trajectory jumps back to a stable fixed point by crossing an SH offset bifurcation. This leads to a sudden jump back to the baseline from an oscillatory state to a quiescent state [third panel in Fig. 4(a3)]. Similarly, the pyramidal neuron and PV+ interneuron in HYP seizure-onset state show SNIC/SH and SN/SH bifurcations pair at seizure onset/offset, respectively [Fig. 5(a4)]. It is worth noting that the trajectory of PV+ interneuron in HYP seizure-onset state shrinks from the limit cycle to a point corresponding to depolarization block via a supercritical Hopf (SupH) offset when [ K ] o approximately approaches 15.2 mM. Then, the neuronal amplitude increases again and the trajectory transit from the point to the limit cycle via a supercritical Hopf (SupH) offset. However, the final bifurcation offset is SH. Hence, we still classify this bifurcation pair type as an SN/SH. When [ K ] o increases beyond a certain threshold, the neurons transition to a depolarizing block state, causing the system to re-enter monostability via SupH bifurcation. It is just like the phase space of the PV+ interneuron shown in Fig. 5(a5) and two types of neurons in Fig. 5(a6).

FIG. 5.

The evolution of the three-dimensional phase space with V, K o, and potassium gate variable ( n ) corresponding to four firing states in Fig. 3: (a1) rest state, K bath = 4 mM; (a2) fast-spiking state, K bath = 10 mM; (a3) synchronous bursting state, K bath = 15 mM; (a4) HYP seizure-onset state, K bath = 16 mM; (a5) sustained ictal activity, K bath = 34 mM; and (a6) sustained depolarization block state, K bath = 40 mM. Red and blue lines represent the trajectories of pyramidal neuron and PV+ interneuron, respectively. The black arrow indicates the direction of movement.

FIG. 5.

The evolution of the three-dimensional phase space with V, K o, and potassium gate variable ( n ) corresponding to four firing states in Fig. 3: (a1) rest state, K bath = 4 mM; (a2) fast-spiking state, K bath = 10 mM; (a3) synchronous bursting state, K bath = 15 mM; (a4) HYP seizure-onset state, K bath = 16 mM; (a5) sustained ictal activity, K bath = 34 mM; and (a6) sustained depolarization block state, K bath = 40 mM. Red and blue lines represent the trajectories of pyramidal neuron and PV+ interneuron, respectively. The black arrow indicates the direction of movement.

Close modal

In this section, we investigate the therapeutic effect of TMAS on HYP seizure onset where different ultrasonic waveforms are chosen to target the pyramid neurons and PV+ interneuron. After the system is in a homeostatic HYP seizure onset for 7 s, continuous TMAS currents are simultaneously applied to both types of neurons. Figures 6(a1) and 6(a2) show the curve change of the pyramidal neurons and PV+ interneuron firing rates (FR) under different ultrasound waveforms with ultrasound intensity on HYP seizure onset, respectively. We find that the pyramidal neuron and PV+ interneuron firing rates gradually decrease with the increase of U I under lower sine wave stimulation, even reaching zero when ultrasound intensity exceeds 0.5 W/cm2. However, pyramidal neuron and PV+ interneuron firing rates gradually increase under upper sine wave stimulation. Interestingly, two types of neuronal firing rates are basically unchanged under sine wave stimulation and remain similar to those in the absence of stimulation. In Figs. 6(b1)6(b3), we present the time series of the membrane potentials and ion concentrations in two classes of neurons under the upper lower wave, sine wave, and upper sine wave stimulation when the parameters B = 0.9 T and U I = 0.7 W / c m 2. The lower sine wave can inhibit neuron discharges and lead to the network transition from HYP seizure-onset state to rest state. In contrast, the upper sine wave can promote neuron discharges and increase the duration of HYP seizure. The sine wave seems to have no effect on neuronal activity, possibly because the promoting and inhibiting effects of the positive and negative phases of the sine wave cancel each other at the microsecond level, and it is ineffective in changing membrane potentials (at the millisecond level) and ion concentrations (at the second level) on a much larger time scale.

FIG. 6.

The therapeutic effect of TMAS on the HYP seizure onset [Fig. 4(a4), K bath = 16 mM] under different ultrasonic waveforms. (a1) and (a2) are the curves of a pyramidal neuron and PV+ interneuron firing rate (FR) vs TMAS ultrasound intensity, respectively. (b1)–(b3) are the time series of membrane potentials and ion concentrations in the pyramidal neuron and the PV+ interneuron under upper sine wave, sine wave, and lower sine wave stimulation with U I = 0.7 W / c m 2, B = 0.9 T, respectively.

FIG. 6.

The therapeutic effect of TMAS on the HYP seizure onset [Fig. 4(a4), K bath = 16 mM] under different ultrasonic waveforms. (a1) and (a2) are the curves of a pyramidal neuron and PV+ interneuron firing rate (FR) vs TMAS ultrasound intensity, respectively. (b1)–(b3) are the time series of membrane potentials and ion concentrations in the pyramidal neuron and the PV+ interneuron under upper sine wave, sine wave, and lower sine wave stimulation with U I = 0.7 W / c m 2, B = 0.9 T, respectively.

Close modal

After exploring the previous results, we choose the continuous lower sine wave ultrasound as the desired TMAS waveform. Below, we examine the threshold of different stimulation parameters for the HYP seizure onset. Compared with the neuronal firing rate which mainly indicates the overall activity level and responsiveness of neurons, ISI can accurately provide information about changes in individual neurons’ discharge frequency and bifurcations of discharge patterns from a dynamical perspective. Therefore, we calculate the ISI of neurons as a metric to study the effect of TMAS on neuronal spike frequency at different magnetic field intensities. It is observed from Fig. 7 that the larger value of the ISI (>1 s) represents the interval between two seizures, while the smaller value represents the time interval between two neuronal spikes in a seizure. Therefore, we hypothesize that the maximal ISI of pyramidal neurons after TMAS is an indicator of epileptic latency. As shown in Fig. 7, we find the maximum ISI of the pyramidal neuron and interneuron gradually increases before reaching the bifurcation point where B = 0.85 T. It is worth mentioning that the maximum ISI of pyramidal neurons exceeds five seconds when 0.6 T B 0.8 T. After crossing the bifurcation point, the ISI of the pyramidal neuron and interneuron suddenly drops to zero. We quantitatively identified a bifurcation point for state transitions by varying the ultrasound parameter B. These results indicate that lower sine wave TMAS can increase the latency of HYP seizures and eventually result in the transition to the resting state at sufficient intensity.

FIG. 7.

ISI bifurcation diagrams of pyramidal neurons (a) and PV+ interneuron (b) vs magnetic field intensity under lower sine wave stimulation when U I = 0.5 W / c m 2, respectively. The value of ISI equals to zero when 0.85 T B 1 T.

FIG. 7.

ISI bifurcation diagrams of pyramidal neurons (a) and PV+ interneuron (b) vs magnetic field intensity under lower sine wave stimulation when U I = 0.5 W / c m 2, respectively. The value of ISI equals to zero when 0.85 T B 1 T.

Close modal

Next, we investigate the state transitions and the corresponding ionic concentration dynamics due to changes in crucial ultrasound parameters including B and U I under lower sine wave stimulation. The dependence of the control effect on B and U I is shown in Fig. 8, distinguished by neuronal discharge amplitude and ISI. Here, the red parameter area indicates epileptiform discharges and is considered as an invalid control. It can be observed in Fig. 8(b1) that the neurons still undergo HYP seizures due to insufficient stimulation when U I = 0.1 W / c m 2 and B = 0.9 T. As U I and B increase together, the latency of HYP seizures gradually exceeds 5 s, and, therefore, the yellow parameter area is classified as a partially invalid control. Compared with those in Fig. 8(b1), the latency of HYP seizure in Fig. 8(b2) significantly increases and HYP seizure duration in Fig. 8(b2) significantly decreases at the same regulatory time. When B and U I are further increased, the blue parameter area shows that neurons after TMAS enter the resting state, which is considered as a valid control. As shown in Fig. 8(b3), the membrane potentials of pyramidal neuron and interneuron were essentially stabilized at −67 and −70 mV, respectively, similar to the membrane potentials of neurons at rest rather than depolarization block. Meanwhile, we find decay of [ Na ] i and [ K ] o in two classes of neurons, which may result from activation of the N a + / K + pump.

FIG. 8.

The therapeutic effect of lower sine wave TMAS on HYP seizure onset [Fig. 4(a4), K bath = 16 mM] as ultrasound intensity ( U I ) and magnetic field intensity ( B ) are varied. (a) Two-parameter space diagrams for U I and B. The three states are distinguished by different colors and represent the invalid control state, partial invalid control state, and valid control state, respectively. The time series of membrane potentials and ion concentrations in the pyramidal neuron and the PV+ interneuron under lower sine wave stimulation vs ultrasound intensity corresponding to the three control states: (b1) U I = 0.1 W / c m 2, B = 0.9 T; (b2) U I = 0.4 W / c m 2, B = 0.9 T; and (b3) U I = 0.7 W / c m 2, B = 0.9 T.

FIG. 8.

The therapeutic effect of lower sine wave TMAS on HYP seizure onset [Fig. 4(a4), K bath = 16 mM] as ultrasound intensity ( U I ) and magnetic field intensity ( B ) are varied. (a) Two-parameter space diagrams for U I and B. The three states are distinguished by different colors and represent the invalid control state, partial invalid control state, and valid control state, respectively. The time series of membrane potentials and ion concentrations in the pyramidal neuron and the PV+ interneuron under lower sine wave stimulation vs ultrasound intensity corresponding to the three control states: (b1) U I = 0.1 W / c m 2, B = 0.9 T; (b2) U I = 0.4 W / c m 2, B = 0.9 T; and (b3) U I = 0.7 W / c m 2, B = 0.9 T.

Close modal

In this paper, we propose a biophysical network to study state transitions from the normal state to the HYP seizure state (and ultimately to the depolarization block state) and explore the oscillatory dynamics of ion concentrations in different states, which is the first time in the model network to reproduce the complete evolution of HYP seizure driven by excess extracellular potassium. The results imply that HYP seizures are triggered by a dynamic decrease in inhibitory GAB A A signaling due to excess bath potassium concentration and abnormal glial uptake strength of potassium. Moreover, we find that the pyramidal neuron and PV+ interneuron in HYP seizure-onset state exhibit SNIC/SH and SN/SH bifurcation pairs at seizure onset/offset, respectively. On this basis, we study the therapeutic effects of TMAS in modulating HYP seizure onset under different ultrasonic waveforms and eventually find that negative waves can effectively increase the latency of HYP seizures and even completely suppress epilepsy by the metric of firing rate. Further, we focus on the modulation of HYP seizure-onset thresholds by lower sine wave TMAS and find that the greater the magnetic field strength and ultrasonic intensity, the better the control effect. More importantly, we show the parameter area that efficiently modulates HYP seizure onset.

Next, we discuss two ionic factors that may contribute to HYP seizure. First, the PV+ interneuron gradually enters a phase of depolarization block as [ K ] o exceeds 16 mM, which leads to more excitation of pyramidal neurons due to the lack of inhibitory GABA currents. When the interneuron enters the sustained depolarization block, the system converts to the sustained ictal activity state because the pyramidal neurons only receive excitatory glutamate signals from other pyramidal neurons rather than inhibitory GABA signals from the PV+ interneuron. In vitro studies in the isolated rat hippocampus in 4-AP or low-Mg medium have also shown that pyramidal neurons are recruited to ictal discharges when GABAergic interneurons undergo a depolarizing block phase.16,45 Interestingly, in our model, the interneuron is more likely to enter a depolarizing block than pyramidal neurons, which is also consistent with physiological experiments.50 For this phenomenon, we are inclined to hypothesize that the AHP current in pyramidal neurons prevents the depolarization block. Second, HYP seizures may be related to intracellular chlorine accumulation. Work in temporal lobe slices indicated that cotransporter NKCC1 activity is upregulated and cotransporter KCC2 activity is downregulated in epileptic patients with hippocampal sclerosis.12,51 In our study, it can be seen in Fig. 4 that chloride ions gradually influx during seizures in response to cotransporter NKCC1 and KCC2 activities. Then, according to Eq. (4), we can find a gradual positive shift in the GABAA reversal potential, which is closely linked to postsynaptic GABAA receptor activation in pyramidal neurons, resulting in dynamically weakening inhibitory GAB A A signaling. In addition, Weiss suggested that HYP onset may be caused by depolarizing inhibitory postsynaptic currents due to chloride ion dysregulation.52 Furthermore, intracellular recordings in rat perirhinal cortices during continuous 4-AP application indicated that HYP ictal onset may be triggered by diminished GABAA signaling due to extracellular potassium accumulation.10 The above two studies validate our proposed mechanism of HYP seizure onset.

Numerous biophysical models using potassium ions as a slow variable in dynamic analyses have shown that the single pyramidal neuron transitions to seizure through a SNIC bifurcation and transitions to depolarization block through a Hopf bifurcation.53,54 However, many phenomenological models with few variables and parameters are characterized by SN/SH bifurcation pairs at the seizure onset/offset.47,55 In our model, the pyramidal neuron and the PV+ interneuron in the HYP seizure-onset state exhibit SNIC/SH and SN/SH bifurcation pairs at seizure onset/offset, respectively. Meanwhile, the behaviors of the pyramidal neuron and PV+ interneuron are similar to seizure-like events dynamics observed in a simplified single-neuron model,56 corresponding to paths “c6” and “c2” in existing dynamic frameworks, respectively.48 We emphasized the bifurcation of onset and offset pairs at specific stages of HYP seizure based on invariant dynamics characterizations, which provide new insight into epileptic mechanisms and treatment. However, due to random noise and ambiguous data, the classifications of epileptic bifurcations in clinical conditions have been a challenge. Therefore, effective identification and quantification of dynamic features may require more work in the future.

As a novel neuromodulation, the safety of TMAS must be taken into consideration. Previous research has shown that high-intensity focused ultrasound causes irreversible nerve tissue necrosis due to thermal injury caused by prolonged exposure, whereas lower-intensity focused ultrasound can activate or inhibit bioelectric activity without any concomitant brain damage.57,58 Therefore, three medical ultrasound metrics including spatial-peak temporal-average intensity ( I spta ), MI, and RI are calculated to estimate the possible biological effects of ultrasound. It is worth noting that the US FDA allows diagnostic systems with an upper limit of 0.72 W/cm2 in I spta, 1.9 in MI, and 6 in RI, which is also used as a reference for ultrasound neuromodulation parameters.59 Hence, in our model, we set the maximum U I to 0.7 W/cm2 to minimize possible thermal effects and comply with diagnostic regulations of I spta. We can calculate p L = 0.155 MPa, S = 7.065 × 10 2 c m 2 with the ultrasonic parameters U I = 0.7 W / c m 2, f = 0.2 MHz, and the diameter of the ultrasonic spot equals to 3 mm. Taking the above results into Eqs. (11) and (12), we obtain M I = 0.347 and R I = 0.047, which are well below the FDA recommended limit that can cause thermal or mechanical damage to nerve tissue. Meanwhile, for the sake of safety, our magnetic field intensity ( B 1 T ) is also within the parameter range required for diagnosis. In addition, ultrasound frequency is also a safety indicator of great concern. The frequency range of diagnostic medical ultrasound is between 1 and 15 MHz, while the frequency range of therapeutic ultrasound is in kHz.60,61 Low-frequency ultrasound (f < 0.65 MHz) can effectively reduce ultrasonic attenuation in the scalp and skull and is widely used to stimulate brain circuits and modulate neural activity.23 In our model, ultrasound frequency (f = 0.2 MHz) also adheres to low-frequency ultrasound standards. Overall, our TMAS strategy is capable of safely and reliably modulating epilepsy.

There are two potential limitations in our model. First, as observed in Fig. 5, ion concentrations such as [ K ] o and [ Na ] i do not return to the baseline level when neurons enter the quiescence state due to the termination of neuronal discharges or modulation by TMAS. We find out that this may be the consequence of the delayed action of the glial cells and excessive bath potassium concentration in our model. Another limitation is that we focus mainly on ultrasound’s temporal rather than spatial properties. Therefore, we did not consider the distribution of induced TMAS current in the brain. Our subsequent work may need to consider the effect of ultrasound attenuation on the skull and brain tissue.

Although our model focuses on a small-scale network of five neurons in the hippocampal CA3 region, we chose two classes of neurons including four pyramidal neurons and a PV+ interneuron, and complemented our neuron model with expressions for the ion concentration. First, the physiological ratio of pyramidal neurons to interneurons in the hippocampus is four to one. Second, some studies indicate that HYP seizure onset may depend on the preponderant involvement of principal (glutamatergic) neuron networks.7,8 Based on the two reasons mentioned above, we propose a quantitative configuration of these two types of neurons as a minimum model. In our paper, we concentrated on the transition of neuronal states induced by micro-level ionic concentration changes rather than the neuronal clustering activities in complex networks. Our results show HYP epilepsy evolution due to the depolarization block of PV+ interneurons and hyperexcitability of pyramidal neurons, which reproduce phenomena observed from physiological experiments. Inspired by in vitro experiments,50 our future work may construct a larger hippocampal CA3 region network where more GABAergic interneuron subtypes such as the cholecystokin interneuron and the somatostatin interneuron are added to study the effects of various interneurons on HYP seizures dynamics. In addition, another step may be to expand the model to include other brain regions to clarify the flow of information in pathological networks by nonlinear methods. Interestingly, according to Fig. 6, we find that the upper sine wave increases neuron firing rate while the lower sine wave decreases neuron firing rate. The different effects of ultrasound waveforms on neurons provide us with some inspiration. For neurological disorders such as epilepsy caused by neuronal over-excitation, we can choose the lower sine wave to inhibit neuron discharges. On the other hand, for neurological disorders such as major depressive disorder caused by neuronal over-inhibition, we can choose the upper sine wave to promote neuron discharges. Meanwhile, based on the high spatial resolution of TMAS, we can further explore potential target brain regions where TMAS can effectively modulate HYP seizure.

This work was supported by the National Natural Science Foundation of China (NNSFC) (Grant Nos. 12102014, 11932003, 32271361, and 12202022) and the National Key Research and Development Program of China (Grant Nos. 2021YFA1000200 and 2021YFC1000202). We thank the Figdraw platform (https://www.figdraw.com) for the visualization of Fig. 1.

The authors have no conflicts to disclose.

Zhiyuan Ma: Conceptualization (lead); Data curation (lead); Writing – original draft (lead); Writing – review & editing (lead). Yuejuan Xu: Data curation (equal); Writing – review & editing (equal). Gerold Baier: Conceptualization (equal); Writing – review & editing (equal). Youjun Liu: Writing – review & editing (equal). Bao Li: Conceptualization (equal). Liyuan Zhang: Conceptualization (equal); Data curation (equal); Writing – original draft (supporting).

The data that support the findings of this study are available from the corresponding author upon reasonable request.

1
P.
Salami
,
M.
Lévesque
,
J.
Gotman
, and
M.
Avoli
, “
Distinct EEG seizure patterns reflect different seizure generation mechanisms
,”
J. Neurophysiol.
113
(
7
),
2840
2844
(
2015
).
2
P.
Perucca
,
F.
Dubeau
, and
J.
Gotman
, “
Intracranial electroencephalographic seizure-onset patterns: Effect of underlying pathology
,”
Brain
137
,
183
196
(
2014
).
3
A.
Bragin
,
J. J.
Engel
,
C. L.
Wilson
,
E.
Vizentin
, and
G. W.
Mathern
, “
Electrophysiologic analysis of a chronic seizure model after unilateral hippocampal KA injection
,”
Epilepsia
40
(
9
),
1210
1221
(
1999
).
4
A. L.
Velasco
,
C. L.
Wilson
,
T. L.
Babb
, and
J. J.
Engel
, “
Functional and anatomic correlates of two frequently observed temporal lobe seizure-onset patterns
,”
Neural Plast.
7
(
1–2
),
49
63
(
2000
).
5
M.
Derchansky
,
S. S.
Jahromi
,
M.
Mamani
,
D. S.
Shin
,
A.
Sik
, and
P. L.
Carlen
, “
Transition to seizures in the isolated immature mouse hippocampus: A switch from dominant phasic inhibition to dominant phasic excitation
,”
J. Physiol.
586
(
2
),
477
494
(
2008
).
6
Z. J.
Zhang
,
J.
Koifman
,
D. S.
Shin
,
H.
Ye
,
C. M.
Florez
,
L.
Zhang
,
T. A.
Valiante
, and
P. L.
Carlen
, “
Transition to seizure: Ictal discharge is preceded by exhausted presynaptic GABA release in the hippocampal CA3 region
,”
J. Neurosci.
32
(
7
),
2499
2512
(
2012
).
7
M.
Avoli
,
M.
de Curtis
,
V.
Gnatkovsky
,
J.
Gotman
,
R.
Köhling
,
M.
Lévesque
,
F.
Manseau
,
Z.
Shiri
, and
S.
Williams
, “
Specific imbalance of excitatory/inhibitory signaling establishes seizure onset pattern in temporal lobe epilepsy
,”
J. Neurophysiol.
115
(
6
),
3229
3237
(
2016
).
8
Z.
Shiri
,
F.
Manseau
,
M.
Lévesque
,
S.
Williams
, and
M.
Avoli
, “
Activation of specific neuronal networks leads to different seizure onset types
,”
Ann. Neurol.
79
(
3
),
354
365
(
2016
).
9
M.
Lévesque
,
P.
Salami
,
J.
Gotman
, and
M.
Avoli
, “
Two seizure-onset types reveal specific patterns of high-frequency oscillations in a model of temporal lobe epilepsy
,”
J. Neurosci.
32
(
38
),
13264
13272
(
2012
).
10
R.
Köhling
,
M.
D’Antuono
,
R.
Benini
,
P.
de Guzman
, and
M.
Avoli
, “
Hypersynchronous ictal onset in the perirhinal cortex results from dynamic weakening in inhibition
,”
Neurobiol. Dis.
87
,
1
10
(
2016
).
11
F.
Wang
,
X.
Qi
,
J.
Zhang
, and
J.
Huang
, “
Astrocytic modulation of potassium under seizures
,”
Neural Regen. Res.
15
(
6
),
980
987
(
2020
).
12
G.
Huberfeld
,
L.
Wittner
,
S.
Clemenceau
,
M.
Baulac
,
K.
Kaila
,
R.
Miles
, and
C.
Rivera
, “
Perturbed chloride homeostasis and GABAergic signaling in human temporal lobe epilepsy
,”
J. Neurosci.
27
(
37
),
9866
9873
(
2007
).
13
G.
Ullah
,
J. R.
Cressman
, Jr.
,
E.
Barreto
, and
S. J.
Schiff
, “
The influence of sodium and potassium dynamics on excitability, seizures, and the stability of persistent states: II.: Network and glial dynamics
,”
J. Comput. Neurosci.
26
(
2
),
171
183
(
2009
).
14
G.
Ullah
,
Y.
Wei
,
M. A.
Dahlem
,
M.
Wechselberger
, and
S. J.
Schiff
, “
The role of cell volume in the dynamics of seizure, spreading depression, and anoxic depolarization
,”
PLoS Comput. Biol.
11
(
8
),
e1004414
(
2015
).
15
Y.
Wei
,
G.
Ullah
,
J.
Ingram
, and
S. J.
Schiff
, “
Oxygen and seizure dynamics: II.: Computational modeling
,”
J. Neurophysiol.
112
(
2
),
213
223
(
2014
).
16
J.
Ziburkus
,
J. R.
Cressman
,
E.
Barreto
, and
S. J.
Schiff
, “
Interneuron and pyramidal cell interplay during in vitro seizure-like events
,”
J. Neurophysiol.
95
(
2
),
3948
(
2006
).
17
J.
Žiburkus
,
J. R.
Cressman
, and
S. J.
Schiff
, “
Seizures as imbalanced up states: Excitatory and inhibitory conductances during seizure-like events
,”
J. Neurophysiol.
109
(
5
),
1296
1306
(
2013
).
18
D.
Gentiletti
,
M.
de Curtis
,
V.
Gnatkovsky
, and
P.
Suffczynski
, “
Focal seizures are organized by feedback between neural activity and ion concentration changes
,”
Elife
11
,
e68541
(
2022
).
19
C.
Yao
,
Z.
He
,
T.
Nakano
, and
J.
Shuai
, “
Spiking patterns of a neuron model to stimulus: Rich dynamics and oxygen's role
,”
Chaos
28
(
8
),
83112
(
2018
).
20
A. V.
Chizhov
,
A. V.
Zefirov
,
D. V.
Amakhin
,
E. Y.
Smirnova
, and
A. V.
Zaitsev
, “
Minimal model of interictal and ictal discharges ‘Epileptor-2,’
PLoS Comput. Biol.
14
(
5
),
e1006186
(
2018
).
21
A.
Buchin
,
A.
Chizhov
,
G.
Huberfeld
,
R.
Miles
, and
B. S.
Gutkin
, “
Reduced efficacy of the KCC2 cotransporter promotes epileptic oscillations in a subiculum network model
,”
J. Neurosci.
36
(
46
),
11619
11633
(
2016
).
22
Y.
Yuan
,
Y.
Chen
, and
X.
Li
, “
Noninvasive transcranial focused ultrasonic-magnetic stimulation for modulating brain oscillatory activity
,”
Eur. Phys. J. Appl. Phys.
73
(
2
),
21201
(
2016
).
23
Y.
Tufail
,
A.
Matyushov
,
N.
Baldwin
,
M. L.
Tauchmann
,
J.
Georges
,
A.
Yoshihiro
,
S. I.
Tillery
, and
W. J.
Tyler
, “
Transcranial pulsed ultrasound stimulates intact brain circuits
,”
Neuron
66
(
5
),
681
694
(
2010
).
24
F.
Chu
,
R.
Tan
,
X.
Wang
,
X.
Zhou
,
R.
Ma
,
X.
Ma
,
Y.
Li
,
R.
Liu
,
C.
Zhang
,
X.
Liu
,
T.
Yin
, and
Z.
Liu
, “
Transcranial magneto-acoustic stimulation attenuates synaptic plasticity impairment through the activation of piezo1 in Alzheimer’s disease mouse model
,”
Research
6
,
130
(
2023
).
25
Y.
Wang
,
L.
Feng
,
S.
Liu
,
X.
Zhou
,
T.
Yin
,
Z.
Liu
, and
Z.
Yang
, “
Transcranial magneto-acoustic stimulation improves neuroplasticity in hippocampus of Parkinson’s disease model mice
,”
Neurotherapeutics
16
(
4
),
1210
1224
(
2019
).
26
S. J.
Norton
, “
Can ultrasound be used to stimulate nerve tissue?
,”
Biomed. Eng. Online
2
(
6
),
6
(
2003
).
27
Z.
Lin
,
X.
Huang
,
W.
Zhou
,
W.
Zhang
,
Y.
Liu
,
T.
Bian
,
L.
Niu
,
L.
Meng
, and
Y.
Guo
, “
Ultrasound stimulation modulates voltage-gated potassium currents associated with action potential shape in hippocampal CA1 pyramidal neurons
,”
Front. Pharmacol.
10
,
544
(
2019
).
28
Z.
Zhang
,
W.
Qiu
,
H.
Gong
,
G.
Li
,
Q.
Jiang
,
P.
Liang
,
H.
Zheng
, and
P.
Zhang
, “
Low-intensity ultrasound suppresses low-Mg(2+)-induced epileptiform discharges in juvenile mouse hippocampal slices
,”
J. Neural Eng.
16
(
3
),
36006
(
2019
).
29
Z.
Lin
,
L.
Meng
,
J.
Zou
,
W.
Zhou
,
X.
Huang
,
S.
Xue
,
T.
Bian
,
T.
Yuan
,
L.
Niu
,
Y.
Guo
, and
H.
Zheng
, “
Non-invasive ultrasonic neuromodulation of neuronal excitability for treatment of epilepsy
,”
Theranostics
10
(
12
),
5514
5526
(
2020
).
30
Y.
Yuan
,
Y.
Chen
, and
X.
Li
, “
Theoretical analysis of transcranial magneto-acoustical stimulation with Hodgkin-Huxley neuron model
,”
Front. Comput. Neurosci.
10
,
35
(
2016
).
31
Y.
Yuan
,
N.
Pang
,
Y. D.
Chen
,
Y.
Wang
, and
X. L.
Li
, “
Theoretical analysis of the effects of transcranial magneto-acoustical stimulation on neuronal firing rhythm and Ca2+ concentration with Chay neuron model
,”
Biomed. Phys. Eng. Express
3
,
055006
(
2017
).
32
M. S.
Olufsen
,
M. A.
Whittington
,
M.
Camperi
, and
N.
Kopell
, “
New roles for the gamma rhythm: Population tuning and preprocessing for the beta rhythm
,”
J. Comput. Neurosci.
14
(
1
),
33
54
(
2003
).
33
M. J.
Bezaire
and
I.
Soltesz
, “
Quantitative assessment of CA1 local circuits: Knowledge base for interneuron-pyramidal cell connectivity
,”
Hippocampus
23
(
9
),
751
785
(
2013
).
34
X. J.
Wang
and
G.
Buzsáki
, “
Gamma oscillation by synaptic inhibition in a hippocampal interneuronal network model
,”
J. Neurosci.
16
(
20
),
6402
6413
(
1996
).
35
M.
Desroches
,
O.
Faugeras
,
M.
Krupa
, and
M.
Mantegazza
, “
Modeling cortical spreading depression induced by the hyperactivity of interneurons
,”
J. Comput. Neurosci.
47
(
2–3
),
125
140
(
2019
).
36
J. R.
Cressman
,
G.
Ullah
,
J.
Ziburkus
,
S. J.
Schiff
, and
E.
Barreto
, “
The influence of sodium and potassium dynamics on excitability, seizures, and the stability of persistent states: I.: Single neuron dynamics
,”
J. Comput. Neurosci.
26
(
2
),
159
170
(
2009
).
37
N.
Doyon
,
S. A.
Prescott
,
A.
Castonguay
,
A. G.
Godin
,
H.
Kröger
, and
Y.
De Koninck
, “
Efficacy of synaptic inhibition depends on multiple, dynamically interacting mechanisms implicated in chloride homeostasis
,”
PLoS Comput. Biol.
7
(
9
),
e1002149
(
2011
).
38
A. V.
Chizhov
,
D. V.
Amakhin
, and
A. V.
Zaitsev
, “
Mathematical model of Na-K-Cl homeostasis in ictal and interictal discharges
,”
PLoS One
14
(
3
),
e213904
(
2019
).
39
Y.
Wei
,
G.
Ullah
, and
S. J.
Schiff
, “
Unification of neuronal spikes, seizures, and spreading depression
,”
J. Neurosci.
34
(
35
),
11733
11743
(
2014
).
40
R. S.
Fisher
,
T. A.
Pedley
, and
D. A.
Prince
, “
Kinetics of potassium movement in normal cortex
,”
Brain Res.
101
(
2
),
223
237
(
1976
).
41
A.
Montalibet
,
J.
Jossinet
,
A.
Matias
, and
D.
Cathignol
, “
Electric current generated by ultrasonically induced Lorentz force in biological media
,”
Med. Biol. Eng. Comput.
39
(
1
),
15
20
(
2001
).
42
W. R.
Hendee
,
E. R.
Ritenour
, and
K. R.
Hoffmann
, “
Medical imaging physics, fourth edition
,”
Med. Phys.
30
(
4
),
730
(
2003
).
43
J. G.
Abbott
, “
Rationale and derivation of MI and TI–a review
,”
Ultrasound Med. Biol.
25
(
3
),
431
441
(
1999
).
44
S. F.
Traynelis
and
R.
Dingledine
, “
Potassium-induced spontaneous electrographic seizures in the rat hippocampal slice
,”
J. Neurophysiol.
59
(
1
),
259
276
(
1988
).
45
M.
Cammarota
,
G.
Losi
,
A.
Chiavegato
,
M.
Zonta
, and
G.
Carmignoto
, “
Fast spiking interneuron control of seizure propagation in a cortical slice model of focal epilepsy
,”
J. Physiol.
591
(
4
),
807
822
(
2013
).
46
E. M.
Izhikevich
,
Dynamical Systems in Neuroscience: The Geometry of Excitability and Bursting
(
MIT Press
,
2007
).
47
V. K.
Jirsa
,
W. C.
Stacey
,
P. P.
Quilichini
,
A. I.
Ivanov
, and
C.
Bernard
, “
On the nature of seizure dynamics
,”
Brain
137
(
Pt 8
),
2210
2230
(
2014
).
48
M. L.
Saggio
,
A.
Spiegler
,
C.
Bernard
, and
V. K.
Jirsa
, “
Fast-slow bursters in the unfolding of a high codimension singularity and the ultra-slow transitions of classes
,”
J. Math. Neurosci.
7
(
1
),
7
(
2017
).
49
M. L.
Saggio
,
D.
Crisp
,
J. M.
Scott
,
P.
Karoly
,
L.
Kuhlmann
,
M.
Nakatani
,
T.
Murai
,
M.
Dümpelmann
,
A.
Schulze-Bonhage
,
A.
Ikeda
,
M.
Cook
,
S. V.
Gliske
,
J.
Lin
,
C.
Bernard
,
V.
Jirsa
, and
W. C.
Stacey
, “
A taxonomy of seizure dynamotypes
,”
Elife
9
,
e55632
(
2020
).
50
P.
Scalmani
,
R.
Paterra
,
M.
Mantegazza
,
M.
Avoli
, and
M.
de Curtis
, “
Involvement of GABAergic interneuron subtypes in 4-aminopyridine-induced seizure-like events in mouse entorhinal cortex in vitro
,”
J. Neurosci.
43
(
11
),
1987
2001
(
2023
).
51
A.
Muñoz
,
P.
Mendez
,
J.
Defelipe
, and
F. J.
Alvarez-Leefmans
, “
Cation-chloride cotransporters and GABA-ergic innervation in the human epileptic hippocampus
,”
Epilepsia
48
(
4
),
663
673
(
2007
).
52
S. A.
Weiss
, “
Chloride ion dysregulation in epileptogenic neuronal networks
,”
Neurobiol. Dis.
177
,
106000
(
2023
).
53
E.
Barreto
and
J. R.
Cressman
, “
Ion concentration dynamics as a mechanism for neuronal bursting
,”
J. Biol. Phys.
37
(
3
),
361
373
(
2011
).
54
K.
Ma
,
H.
Gu
, and
Z.
Zhao
, “
Fast-slow variable dissection with two slow variables: A case study on bifurcations underlying bursting for seizure and spreading depression
,”
Int. J. Bifurcation Chaos
31
(
06
), 2150096 (
2021
).
55
K. E.
Houssaini
,
A. I.
Ivanov
,
C.
Bernard
, and
V. K.
Jirsa
, “
Seizures, refractory status epilepticus, and depolarization block as endogenous brain activities
,”
Phys. Rev. E
91
(
1
),
010701
(
2015
).
56
D.
Depannemaecker
,
A.
Ivanov
,
D.
Lillo
,
L.
Spek
,
C.
Bernard
, and
V.
Jirsa
, “
A unified physiological framework of transitions between seizures, sustained ictal activity and depolarization block at the single neuron level
,”
J. Comput. Neurosci.
50
(
1
),
33
49
(
2022
).
57
G.
ter Haar
, “
Therapeutic applications of ultrasound
,”
Prog. Biophys. Mol. Biol.
93
(
1–3
),
111
129
(
2007
).
58
H.
Baek
,
K. J.
Pahk
, and
H.
Kim
, “
A review of low-intensity focused ultrasound for neuromodulation
,”
Biomed. Eng. Lett.
7
(
2
),
135
142
(
2017
).
59
T. R.
Nelson
,
J. B.
Fowlkes
,
J. S.
Abramowicz
, and
C. C.
Church
, “
Ultrasound biosafety considerations for the practicing sonographer and sonologist
,”
J. Ultrasound Med.
28
(
2
),
139
150
(
2009
).
60
W. J.
O’Brien
, “
Ultrasound-biophysics mechanisms
,”
Prog. Biophys. Mol. Biol.
93
(
1–3
),
212
255
(
2007
).
61
P.
Bowary
and
B. D.
Greenberg
, “
Noninvasive focused ultrasound for neuromodulation: A review
,”
Psychiatr. Clin. North Am.
41
(
3
),
505
514
(
2018
).