Theoretical analysis of effects of transcranial magneto-acoustical stimulation on neuronal spike-frequency adaptation

Transcranial magneto-acoustical stimulation (TMAS) is a noninvasive technique that has advantages in spatial resolution and penetration depth. It changes the firing properties of neurons through the current generated by focused ultrasound and a static magnetic field. Spike-frequency adaptation is an important dynamic characteristic of neural information processing. To address the effects of TMAS on neural spike-frequency adaptation, this study employs some ultrasound and magnetic field parameters, such as magnetic flux density, ultrasonic intensity, fundamental ultrasonic frequency, modulation frequency, and duty cycle. Using these different ultrasound and magnetic field parameters, membrane potential curves, spike-frequency curves, and adapted onset spike-frequency curves are exhibited and analyzed. The results show that spike-frequency adaptation is strongly dependent on ultrasonic intensity and magnetic flux density and is rarely affected by other parameters. However, modulation frequency and duty cycle influence membrane potentials and spike frequencies to some degree. This study reveals the mechanism of the effects of TMAS on neural spike-frequency adaptation and serves as theoretical guidance for TMAS experiments.


Background
Transcranial magneto-acoustical stimulation (TMAS) is a new brain stimulation technique that has been studied in the neuroscience field. TMAS noninvasively stimulates neurons by generating an electric current with ultrasonic waves and a static magnetic field [1][2][3][4][5]. TMAS has more advantages in spatial resolution and penetrability than transcranial direct current stimulation and repetitive transcranial magnetic stimulation, which have been used in neuropsychiatric disease therapy [6,7]. The high spatial resolution of TMAS benefits from its small ultrasonic spot, which is approximately 2 mm in diameter [8]. Transcranial focused ultrasound stimulation (TFUS) is another noninvasive brain stimulation technique that has been demonstrated to have spatial resolution and penetration equally as good as that of TMAS [9,10]. However, our previous studies have shown that TMAS improves the neuromodulation performance in vivo animal experiments compared with TFUS [11,12].
Spike-frequency adaptation is a significant neurobiological phenomenon in neural information processing, in which the spike frequencies of neurons gradually decline under constant external stimulation [13][14][15]. This effect can be attributed to two common types of adaptation currents, voltage-sensitive potassium currents (M-type currents, I M ) and calcium-sensitive afterhyperpolarization currents (AHP-type currents, I AHP ), which have different properties. Specifically, the M-type ones are active even in the quiescent state, while the AHP-type ones are activated by firing activity. Both types affect the dynamics of neurons, which deserves to discuss [13,[16][17][18][19][20][21]. Previous studies examined spike-frequency adaptation. For example, Benda et al. found that fast spike-frequency adaptation can result in intensity invariance of neural responses [16]. They further studied cases in which spike-frequency adaptation is suitable to be modeled as an adaptation current or a dynamic firing threshold, respectively [17]. Yi et al. concluded that the dynamic characteristics of neurons in electrical fields can be determined by neuronal morphology and neuronal adaptability [18]. Ha et al. investigated the influence of the Ca 2+ -activated chloride channel anoctamin-2 on spike-frequency adaptation in thalamocortical and central nervous system neurons [19,20]. Barraca et al. studied the influence of spike-frequency adaptation on the dynamics of a balanced network of neuronal models [13]. Chang et al. investigated the effects of external electric fields at extremely low frequencies on the adaptability of neurons in the Ermentrout model [21]. However, the influences of TMAS on the spike-frequency adaptation of neurons are still unknown. Thus, we applied TMAS to the Ermentrout model to analyze the effects of TMAS on neuronal spike-frequency adaptation.
The Ermentrout model is a typical neuron model that can be considered as a reduced one-compartment counterpart of the Traub-Miles model [22]. Both M-type current and AHP-current were introduced to the Ermentrout model during the neuronal spike-frequency adaptation simulations [17,23,24]. The membrane potentials of sensory neurons can be changed by external physical stimulation [16,17]. Thus, when the Ermentrout model is stimulated by focused ultrasound waves in a static magnetic field, a stimulus current is added to the neuronal model, and a spike is emitted when the membrane potentials exceed the threshold. A spike-frequency curve is generated by assigning the inverse interspike interval from the membrane potential curve for each spike [25][26][27][28].
In our study, we extend the Ermentrout model by adding the electric current induced by TMAS to the membrane balance equation. We adopt an additional adaptation M-type current to study the firing characteristics and analyze the membrane potential curves as well as the spike-frequency curves as we vary the ultrasound and magnetic field parameters, such as magnetic flux density, ultrasonic intensity, fundamental ultrasonic frequency, modulation frequency, and duty cycles. We also discuss the changes in adapted onset spike-frequency curves under magneto-acoustical parameters based on different initial values of adaptation variables. We investigate the characteristics of spike-frequency adaptation for the Ermentrout neuron under TMAS through computational simulation and provide some discussion.

TMAS principle
A schematic showing the principle of TMAS is presented in Fig. 1a. In TMAS, a low-intensity focused ultrasound beam can drive the charged ions in tissue fluid. Combined with the static magnetic field, Lorentz forces in opposite directions generated on positive and negative charges. The moving charged ions forced by the Lorentz force form an electric current, I ext , which is perpendicular to the direction of the static magnetic field and the ultrasound wave. In this study, we establish the Cartesian coordinate axes, in which ultrasound waves longitudinally propagate along the z-axis, the static magnetic field and the created current density are in the x-axis and the y-axis respectively. The relationship between the electric current density, J y , and the ultrasound and magnetic parameters can be expressed as follows [29,30]: where σ is the tissue conductivity, B x is the magnetic flux density, Ŵ is the ultrasound intensity, ρ is the tissue density, c 0 is the ultrasound speed, and f u is the fundamental ultrasonic frequency. This electric current density, J y , can be represented as the external input current I ext of the Ermentrout neuron under TMAS.
In the simulation, to generate a sodium current and action potentials through stimulation, the fundamental ultrasound wave must be a sine wave with a DC offset [12]. The ultrasonic intensity influence the current density in nerve tissue. In the literature, both continuous and modulated ultrasonic waves are applied to stimulate the motor cortex, and both can safely and repeatedly evoke local field potentials [8]. Thus, we used pulsed ultrasound waves modulated by a square wave and continuous waves in our computational study. These waves obey the following equation: where RF and DC are the modulation frequency and the duty cycle, respectively. Thus, the stimulus current, I ext , generated by TMAS can be expressed as follows: The fixed parameters of this mathematical expression are listed in Table 1 [12].

The Ermentrout model
The Ermentrout model is advantageous for describing voltage and current changes in the neuronal membrane. A circuit diagram of the Ermentrout model with external current, I ext , is exhibited in Fig. 1b. The expression of the membrane potential can be described as follows [22]: where C m is the membrane capacity, I ext is the input stimulus current in tissue fluid, and V is the transmembrane voltage component generated by I ext . The voltagegated currents I Na , I K , and I Ca are given by the following equations: where g Na and E Na are the maximum conductance and the reversal potentials of sodium channels, respectively. Similarly, g K , E K , and g Ca , E Ca represent the same parameters for the potassium and calcium channels, respectively. The kinetics of the gating variables x ∈ m, n, h are given by The leak current can be expressed as follows: where g L and E L represent the maximum conductance and the reversal potential of the chlorine channel, respectively.
The M-type current is a kind of depolarizationactivated and hyperpolarizing potassium current that decays with the activation of muscarinic acetylcholine receptors. It obeys the following expression: where g M is the maximum conductance, and w is the gating variable.
The AHP-type current is activated by the influx of calcium that results from firing and can be expressed as where [Ca 2+ ] is the intracellular calcium concentration.
The kinetics of the M-type current w can be expressed as follows: The Ermentrout model is a reduced one-compartment counterpart of the Traub-Miles model. We introduce both M-type currents and AHP-type currents to establish the modified Ermentrout model for adaptation. However, when we simulated the neuronal spike-frequency adaptation, either an M-type current or an AHP-current was considered. In our simulation, we take the M-type current model as an example. Therefore, we set I AHP of the Ermentrout model constant in this paper. The fixed parameters used for the Ermentrout model in this study are listed in Table 2 [25].
The simulation was performed with MATLAB Simulink software (2018, MathWorks, USA). The main code used for the simulation has been deposited in a publicly available repository. It will be freely available to any scientist wishing to use them for non-commercial purposes.

Adaptation of the Ermentrout model
Adaptation is a significant phenomenon of neurons in neural signal processing. The membrane potentials of sensory neurons can be changed by external physical stimulation. Therefore, we pay close attention to the membrane potential curve of the adaptive model under TMAS. A spike occurs whenever the membrane potential exceeds its upper threshold. We can obtain spikefrequency by assigning each time the inverse inter-spike interval that contains this time. Both the membrane potential curves and spike frequency curves exhibit the adaptation processes. In the spike-frequency curve, two frequencies need special attention. One is the initial spike-frequency, f 0 , which denotes the inverse of the first interspike interval in the membrane potential curve; the other is the steady-state spike-frequency, f ∞ , which represents the stable spike-frequency after an extended period. In addition, the time from the initial to (10) steady-state spike-frequency, which is the settling time, is of interest for spike-frequency adaptation.
Other adaptive characteristics of neurons can be well captured by the onset spike-frequency curve. The onset spike-frequency curve, f 0 (X) , is obtained by mapping the corresponding initial firing frequency f 0 to each input parameter X. Besides, to reveal the onset responses of neurons to stimulus intensity changes, we build the adapted onset spike-frequency curve, f 0 (X, A 0 ) , which consists of different onset spike-frequency curves under different initial stimulus inputs A 0 . The adapted curves indicate how the neuron being in a stable and fixed adaptation state will initially respond to new stimuli. While each neuron has exactly one onset and one steady-state spike-frequency curve, it has for each adaptation varaible a different adapted spike-frequency curve.

Membrane potential curves and spike-frequency curves under different TMAS input parameters
To investigate the effects of different TMAS input signals on the spike-frequency adaptation, we discuss the membrane potential curves and the spike-frequency curves under different TMAS input parameters. As described in the section of methods, five variable parameters ( B x , Ŵ , f u , RF, and DC) of the modulated TMAS influence the spike frequency. To obtain the specific conclusions, we designate one as an adjustable parameter and fix the values of the others as follows: B x = 2 T, Ŵ = 3 W/cm 2 , f u = 500 kHz, RF = 1 Hz and DC = 50% . Then, we discuss the characteristics of the membrane potential curves and the spike-frequency curves for different TMAS parameters.
i. To investigate the influence of TMAS on the neural spike frequency, the membrane potential is simulated based on the Ermentrout model using a range of magnetic flux densities from 0.5 to 8 T. The other TMAS input parameters are set to fixed values. Figure 2a-c show the plots of membrane potentials for magnetic flux densities of 0.5 T, 1 T, and 2 T, respectively. One can see that the spikefrequency ultimately settles into a steady state after a long period in which the intervals tend to be the same value. As shown in Fig. 2d, as B x increases, the spike-frequency curve shifts rightward and vertically toward those with a higher magnetic flux density, which means that the initial spike-frequency increases, and undergoes a longer transient process to attain the higher steady-state spikefrequency. In addition, the onset and steady-state spike frequencies of neurons, as well as the settling time, can be obtained through quantitative calcula-  Fig. 3a-c. Similar to the results exhibited in Fig. 2, the spike-frequency decays with time and ultimately enters to a steady state. We also computed the spike frequencies, and the results are shown in Fig. 3d. We find that the spike-frequency curve shifts rightward and horizontally as the ultrasound intensity increases, which means that as ultrasonic intensity increases, the onset spike-frequency f 0 and the steady-state spike-frequency f ∞ increase, and the settling time becomes longer. The detailed data are shown in Table 4.    Fig. 5a-d, respectively. The spike frequencies settling time of neurons can be obtained through quantitative calculation, as shown in Table 5. Clearly, the interspike interval decreases  as the modulation frequency RF increases. From Fig. 5e, f, we can see that the spike-frequency curves under different modulation frequencies have the same initial spike-frequency, and reach higher steady values more quickly with increasing modulation frequencies. However, when the modulation frequencies are too small, the spike frequency in both the transient and steady states fluctuates with time. In addition, the frequency of the spike-frequency curve increases as the modulation frequency increases, while the amplitude decreases. V. We also investigate the neuronal spike-frequency adaptation under TMAS using various duty cycles, from 10 to 95% . We set the modulation frequency, ultrasonic intensity, magnetic flux density, and fundamental ultrasonic frequency using the fixed values mentioned above. The membrane potential curves with duty cycles of 30% , 50% , 70% , and 90% are shown in Fig. 6a-d. The detailed spike frequencies and settling time with different duty cycles are shown in Table 6. The results indicate that a gap between a cluster decreases as the duty cycle increases. The membrane potentials produce bursting activities with the same initial spike frequency, but the spike frequency within each cluster decreases with both increasing time and increasing duty cycle, as shown in Fig. 6e, f.

Adapted spike-frequency curves with various initial stimulus inputs
Adapted spike-frequency curves f 0 (X, A 0 ) represent the onset response of the neuron to input changes. They can be obtained by taking the different values of A 0 , which represent the initial value of stimulus input. Figure 7a-e shows the adapted spike-frequency curves changing with the magnetic flux density, ultrasonic intensity, fundamental ultrasonic frequency, modulation frequency, and duty cycle, respectively. In these  shown in Fig. 7a, b. When magnetic flux density varies from 0.5 to 4 T, and ultrasonic intensity varies from 0.2 to 2 W/cm 2 , the adapted spike-frequency curves have the same form, and the spike-frequency intervals distribute uniformly, which means the neuron has strong adaptation. However, as these two input parameters further increase, the adapted spike-frequency curves, f 0 (B x , A 0 ) and f 0 (Ŵ, A 0 ) , do not match the onset ones ( A 0 = 0). We can conclude that the neuron has weak adaptation.
Other adapted spike-frequency curves under changing fundamental ultrasonic frequency, modulation frequency, and duty cycle are presented in Fig. 7c-e. We can observe that the initial spike frequencies do not change as the adaptation variables vary, whereas these onset spike frequencies decrease as the initial stimulation input increase. We can conclude that the fundamental ultrasonic frequency, modulation frequency, and duty cycle have no significant influence on the spike-frequency adaptation.

Discussions
Based on the simulation results in the section of results, we analyzed the influence of different TMAS parameters on the adaptation process of the Ermentrout neuron. The conclusions demonstrate that both the initial and steadystate spike frequencies of the neuron increase with the increase of magnetic flux density or ultrasonic intensity of TMAS. According to the physiological mechanism of TMAS, the amplitude of the stimulus input increases as the magnetic flux density and ultrasound intensity increase, thus shortening the time required for membrane potential to reach the threshold. Therefore, both the initial spike frequency and the steady-state spike frequency of the neuron increase. The fundamental ultrasonic frequency of TMAS is in kilohertz level. Such a high-frequency transformation of the input intensity has little influence on the current in the tissue fluid of a neuron. Therefore, changes in the fundamental ultrasonic frequency of TMAS do not affect the neural spike frequency. Besides, we used a square wave to modulate the sinusoidal ultrasonic wave with offset in our simulation. Changes in the modulation frequency and duty cycle alter the current waveform of TMAS and then affect the frequency of the stimulus input acted on neurons. Therefore, changes in the modulation frequency and duty cycle lead to changes in the firing mode and the firing rhythm of the neurons.
We analyze the influence of different TMAS parameters on the spike-frequency adaptation of the Ermentrout neuron based on the adapted spike-frequency curves with various initial stimulus inputs in the section of results. The results show that the fundamental ultrasonic frequency, modulation frequency, and duty cycle do not affect the neural spike-frequency adaptation. The adaptation state usually alters with the slowly changing stimulus input and is subtracted from the stimulus input, so that the average value of the stimulus input changes. As a result, only rapidly fluctuating stimulus signals can be transmitted. In other words, the higher the alternating frequency, the strong adaptation the neuron exhibits. TMAS used high-frequency ultrasound to act on neurons. Therefore, there is no effect on the neural spike frequency adaptation of the fundamental ultrasonic frequency, modulation frequency, or duty cycle. However, the increase of magnetic flux density and ultrasonic intensity directly causes the enhancement of the stimulus strength, and the adaptation current strongly fluctuates between spikes. Therefore, the stronger the TMAS input, the worse the adaptation exhibits. According to Kirchhoff 's current theorem, the dynamic equation of membrane potential is derived from the vector sum of ionic currents. Therefore, with the increase of adaptation current, the onset spike frequency of neurons decreased.
Adaptation, which represents not only the current state of adaptation but also the process of dynamic change, plays a significant role in neural information processing. In particular, a variety of firing currents lead to spikefrequency adaptation and can explain many phenomena, such as weak stimulation and enhanced rapid response to stimulation. Therefore, the adaptation current can be

Conclusions
In this paper, we investigated how TMAS influences neuronal spike-frequency adaptation by exploring different ultrasonic and magnetic field parameters, and furthermore, how it affects neural information transmission. Based on a numerical analysis, the magnetic flux density and ultrasonic intensity are found to be the determining parameters of spike-frequency adaptation. In addition, while the different modulation frequencies and duty cycles may lead to different membrane potentials and spike frequencies, they have no influence on spikefrequency adaptation. Based on the numerical analysis of neuron model, this paper explores the adaptation of  Fig. 7 The adapted onset spike-frequency curves with different initial stimulus inputs neurons to the TMAS input and the effect of the stimulus current on the discharge activity of biological nervous system. These results enrich the theoretical basis of TMAS as a potential clinical treatment of nervous system diseases.