A Route to Chaotic Behavior of Single Neuron Exposed to External Electromagnetic Radiation

1. Introduction

Dynamic behaviors of single neuron or neural network are essential to understanding complex behaviors in brain or even serious diseases in nervous system. On account of complexity of nervous system, hundreds of equations should be set up to describe neuronal behavior precisely. However, main properties like electrical activities in neurons catch the attention of most researchers, and based on which, several models are established to physically study mode selection of neurons. For example, four-variable Hodgkin-Huxley (HH) equations are usually used to investigate the main properties of neurons via membrane potential (Hodgkin and Huxley, 1990). For further simplification, a three-variable neuron model, called as Hindmarsh-Rose (HR) equation, is reduced from the four-variable HH model (Hindmarsh and Rose, 1982). Besides ordinary differential equations, discrete dynamical systems or called maps, can also be considered as valid phenomenological neuron models to govern evolution of the transmembrane voltage and the dynamics of ionic conductances (Ibarz et al., 2011).

Among all these equations, Fitzhugh-Nagumo (FHN) neuron model shows its validity on oscillatory dynamic behavior in a neuron (Fitzhugh, 1961; Nagumo et al., 1962). The model presents properties of Van der Pol oscillator which can also be analyzed by fast- slow system (Krupa et al., 1997). Many methods are applied for constructing the exact solutions of the FHN equation (Li and Guo, 2006; Dehghan et al., 2010) to study transmission of nerve impulses. FHN equation is usually considered as a non-linear ordinary differential wave equation, and its soliton solutions are detected in different ways (Abbasbandy, 2008; Triki and Wazwaz, 2013). As the most important and complex system in nature, the brain contains a tremendous number of neurons and gliocytes. Collective behaviors of a large set of neurons and dynamics properties of neuron network are of much more importance. Patterns like spiral waves and targeted waves and their breakdown show the complexity of brain and related to horrible disease in neuron system. Sported, stripe, and hexagon patterns are also discovered in a modified FHN model which are very similar to the situation in reaction diffusion system.

Among all the factors affecting dynamic behaviors and pattern formation in neuron model, noise, and magnetic flow etc. are the most important and studied by many researchers. Different dynamical regimes are observed induced by external noise (Garcaojalvo and Schimanskygeier, 1999). Wu measures the pattern transition from subexcitable to excitable media (Ying et al., 2013). Colored noise can enhance the stochastic resonance in FHN neuronal model (Nozaki and Yamamoto, 1998). Noise could also be suppressed by a strong periodic signal (Pankratova et al., 2005). Besides, magnetic flux shows its high affection on collective electrical activities and signals propagation among neurons (Lv et al., 2016). Specifically, mode transition is detected and mismatch of frequency between electromagnetic radiation and the system is found (Ma et al., 2017). Magnetic flux on membrane potential is realized by a mimristor coupling, which leads to non-linear quasi-periodic spatial-temporal patterns (Mvogo et al., 2017).

As a matter of fact, a number of researchers are working on the effects on many facets of biological phenomena of electromagnetic radiation (Ma and Tang, 2017), such as central nervous system (Hossmann and Hermann, 2003), newborn rat cerebellar granule neurons (Lisi et al., 2005), oxidative damage to mitochondrial DNA in primary cultured neurons exposed to 1,800 MHz radio frequency radiation (Xu et al., 2010), cultured hippocampal neurons of rats exposed to microwave radiation (Xu et al., 2006). As to the neuron network, electromagnetic radiation can also play an important role in regulating collective behaviors among a large number of neurons. Usually, chemical and electric synapse is considered as the main type of connection between neurons. But in Ma's view, field coupling also lights the shadow of understanding synchronization problems in neuronal network (Ma and Tang, 2017). Dynamical features of neuron model with electromagnetic radiation considered should be paid more attention.

In spite of the fact that external magnetic flux is introduced into neuron network and it triggers complex patterns, the affection on a single neuron is still worth exploring deeply. Ma (Ma et al., 2017) develops FHN model and chooses a sinusoidal function as the external magnetic flux. There also exists a sinusoidal function as a transmembrane current mapped from external forcing in original model. The fact means the neuron system is driven by a pair of plane waves which leads to abundant non-linear behaviors in system. In order to compare with the results of Ma, we choose the same improved model, even with identical parameter values in reference paper (Ma et al., 2017). Model and numerical methods are explained in section 2. We provide and analysis the main numerical results when varying frequency of external transmembrane current ω in section 3, and when varying the amplitude of external electromagnetic radiation A in section 4, respectively. The conclusions are drawn in section 5.

2. Model Description

The improved FHN model, with magnetic flux considered, is displayed as Equation (1).

$\begin{array}{ll}\{\begin{array}{l}\frac{\text{d}u}{\text{d}t}=-k(u-a)(u-1.0)-uv+I_{st}+k_0\rho (\phi )u\\ \frac{\text{d}v}{\text{d}t}=(\epsilon +\frac{{\mu }_{1}v}{u+{\mu }_2})[-v-ku(u-a-1.0)]\\ \frac{\text{d}\phi }{\text{d}t}=k_{1}u-k_2\phi +{\phi }_{ext}.\end{array}& (1)\end{array}$

The additive magnetic flux is induced by changing the distribution of ionic concentration of lectrolytes. The evolution of magnetic flux across the membrane is described by the third equation in Equation (1), which has great effects on membrane potential u in first equation in Equation (1) by a induced current denoted by the last term k₀ρ(φ)u. The memory conductance ρ(φ) of a memristor controlled by magnetic flux, which is often described by

$\begin{array}{ll}\rho (\phi )=\alpha +3\beta {\phi }^2& (2)\end{array}$

The second v is slow variable for current and magnetic flux across the membrane. The physical meaning of other parameters are present in the reference (Ma et al., 2017). It is important that both the external electromagnetic radiation φ_ext and the transmembrane current I_st are chosen as trigonometric functions like,

$\begin{array}{ll}{\phi }_{ext}=A\text{ cos }2\pi ft& (3)\end{array}$

$\begin{array}{ll}{I}_{st}=I_0\text{ sin }\omega t& (4)\end{array}$

respectively. The system is driven by these two plane waves. As to simulation, we use fourth order Runge-Kutta algorithm and the time step is h = 0.01. We choose the same initial values (u, v, φ) = (0.2, 0.1, 0.8). Values of other parameters are listed here, a = 0.15, μ₁ = 0.2, μ₂ = 0.3, ε = 0.002, α = 0.1, β = 0.2, I₀ = 0.6, k₀ = −1, k₁ = 0.2, k₂ = 1.0.

3. Main Numerical Results As Varying Angular Frequency ω

3.1. Main Types of Motions of Membrane Electrical Behaviors

Varying transmembrane current I₀ could lead to transformation of mode of electrical activities. However, in this section we focus on ω as our control parameter in system. Ma (Ma et al., 2017) already discovered different types of dynamic behavior like periodic and chaotic-like motions and even bursting phenomenon in the sampled time series for membrane potentials. As matter of fact, periodic, quasi-periodic, and chaotic motions are all discovered in this situation. The parameters related to the external electromagnetic radiation are selected by A = 0.1 and f = 0.01 in this section.

In order to specify the different types of behavior, Poincáre section φ = 0 is used. Poincáre section can help us not only estimate the motion types but also decide the value of number N of Period-N when the motion is periodic (Song et al., 2015a,b, 2016). Sampled time series for membrane potential and their Poincáre sections are plotted in Figures 1–3 as ω = 0.1256, ω = 0.4, and ω = 0.6, respectively. According to Poincáre section, three isolated points are displayed in Figure 1B, which means the dynamic behavior when ω = 0.1256 is periodic motion. Closed orbit in Poincáre section is shown in Figure 2B, corresponding to the quasi-periodic motions as ω = 0.4, which is zoomed in Figure 4 more clearly. When ω reaches up to ω = 0.6, the closed orbit is destroyed and certain structure appears in Poincáre section, shown in Figure 3B.

FIGURE 1

Figure 1. Sampled time series (A) and its Poincáre section (B) as ω = 0.1256.

FIGURE 2

Figure 2. Sampled time series (A) and its Poincáre section (B) as ω = 0.4.

FIGURE 3

Figure 3. Sampled time series (A) and its Poincáre section (B) as ω = 0.6.

FIGURE 4

Figure 4. Enlarged view of Figure 2B, A close circle in Poincáre section as ω = 0.4.

Furthermore, if we keep going to increase ω, quasi-periodic, periodic, and chaotic motions can also be captured, which are plotted in Figures 5–7. Discrete points, strange structures, and closed circles are displayed in Poincáre section, respectively.

FIGURE 5

Figure 5. Sampled time series (A) and its Poincáre section (B) as ω = 0.618.

FIGURE 6

Figure 6. Sampled time series (A) and its Poincáre section (B) as ω = 2.0.

FIGURE 7

Figure 7. Sampled time series (A) and its Poincáre section (B) as ω = 4.0.

3.2. Modes Selection to Angular Frequency ω

It seems that motion of system changes irregularly when we increase angular frequency ω. In order to investigate mode selection, mode transition process is studied carefully. Usually, inter-spike interval (ISI) (Duan et al., 2008; Wang et al., 2014) is an effective tool to analysis the bifurcation of system, and it is also shown in reference (Ma et al., 2017). Also, it is easy to use ISI diagram to describe the multi-mode of electrical activities. Hierarchy in ISI diagram means the different discharge process. However, it is difficult to estimate the type of motion in ISI diagram. We propose another diagram to distinguish different types of motion. We plot the second variable v appearing its corresponding Poincáre section under every parameter value. A few discrete points can be observed under certain parameter value if the motion is periodic. At least it is effective to distinguish the periodic motion from quasi-periodic and chaotic motions. We use nomenclature of semi-Poincáre diagram to call this method to describe transition of motion types.

We plot the semi-Poincáre diagram with 0 < ω < 3.0, which is shown in Figure 8. Enlarged view of Figure 8 in different ω intervals is shown in Figure 9. Discrete points and continuous lines distribute alternatively. Even slight change of ω can lead to great change of motion type. Poincáre section as ω is around 0.4 is shown in Figure 10. Discrete points and closed circles appear in Poincáre section alternatively, which means that type of motion switches back and forth between periodic and quasi-periodic motion. But according to Figure 9A the scenario of switching between the two motion types is interrupted by chaotic behaviors. The motion type of system becomes chaotic when ω is >0.778. and it lasts until ω = 2.208. One of ω interval of chaos characterized by continues lines is shown in Figure 9B. In other words, in the interval 0.778 < ω < 2.208, periodic and quasi-periodic motions disappear and only chaotic motion is captured. However,when ω > 2.208 the system goes back to the situation of switching back and forth between periodic and quasi-periodic motions, shown in Figures 9C,D. But the scenario is still very different from that in interval ω < 0.778. specifically, compared to situation of ω < 0.778, quasi-periodic motion dominates the interval ω > 2.208 and the chance of occurrence of periodic motions is small. It is highly likely that quasi-periodic motion is interrupted by periodic motion occasionally.

FIGURE 8

Figure 8. Semi-Poincáre diagram when varying ω during (0, 3.0).

FIGURE 9

Figure 9. Enlarged views of semi-Poincáre diagram in different ω intervals. (A) 0.72 < ω < 0.82, (B) 1.50 < ω < 1.60, (C) 2.20 < ω < 2.30, and (D) 2.9 < ω < 3.0.

FIGURE 10

Figure 10. Poincáre sections at different ω. Transition of motion type when changing ω.

Hopf bifurcation plays the key role in motion type selection. Irreducible frequency introduced into system by Hopf bifurcation leads the system from periodic motion to quasi-periodic motion or from quasi-periodic to chaos. In order to specify the statement, we also provide the frequency spectrum of second variable for the three different motion type when the parameter ω is closed to each other. We select ω = 0.754, ω = 0.755, and ω = 0.78 for periodic, quasi-periodic, and chaotic motion types. One, two, and three irreducible frequencies are discovered in frequency spectrum diagram (see Figure 11). New irreducible frequency generated via Hopf bifurcation brings the discrete discrete spectrum to continuous spectrum of chaotic motion. The route to chaos of this kind by a few Hopf bifurcations is called Ruelle-Takens route to chaos (Ruelle and Takens, 1971). A series of supercritical and subcritical Hopf bifurcations lead to the transition of different motion types alternatively.

FIGURE 11

Figure 11. Frequency spectrum of second variable for the three different motion types. (A) ω = 0.754, (B) ω = 0.755, and (C) ω = 0.78.

4. Main Results as Varying the Amplitude of External Electromagnetic Radiation A

External electromagnetic radiation is chosen as a cosine function φ_ext = A cos 2πft. The amplitude A is also an important control parameter, which decides the state of system when other parameters are fixed. We still use semi-Poincáre diagram to describe the scenario as varying A during 0 < A < 1.0. The situation of A > 1.0 will not shown in this paper because the the system stay in equilibrium state, which means all the variables remain constants.

Three semi-Poincáre diagrams are plotted as ω = 0.618, ω = 0.4, and ω = 0.6, corresponding to periodic, quasi-periodic, and chaotic motion when A is fixed in 0.1, shown in Figures 12–14, respectively.

FIGURE 12

Figure 12. Semi-Poincáre diagram as A through (0, 1.0) when ω = 0.618.

In large span of 0 < A < 1.0, periodic motion is discovered, but is is interrupted by quasi-periodic motion several times, when ω = 0.618 (Figure 12). When A is fixed at 0.4, it seems that quasi-periodic motion is more easy to be captured during the interval 0 < A < 1.0 (Figure 13). Similarly, the chance we find chaotic behavior is much more than periodic and quasi-periodic motions when ω = 0.6, which can be reflected by Figure 14.

FIGURE 13

Figure 13. Semi-Poincáre diagram as A through (0, 1.0) when ω = 0.4.

FIGURE 14

Figure 14. Semi-Poincáre diagram as A through (0, 1.0) when ω = 0.6.

It is worth noting that Ruelle-Takens route to chaos also exists when varying the parameter A. Take the scenario shown in Figure 14 for example, there are only two closed orbits in Poincáre section when A < 0.0019, however, three closed orbits appear when A is >0.0019. The new orbit is growing lager when we increase A. The process of growth of the new orbit is stopped by a Hopf bifurcation and system enters the chaotic state when A is >0.0832. The growth process of orbit in Poincáre section is shown in Figure 15.

FIGURE 15

Figure 15. Process of growth of the orbit in Poincáre section when ω = 0.6.

5. Conclusions and Discussion

In this paper, we study dynamical behavior of a single neuron system driven by two plane waves. One of them is provided by external current denoted by I_st, and another comes from the external electromagnetic radiation φ_ext.

With frequency of electromagnetic radiation φ_ext fixed, system state is chosen from periodic, quasi-periodic, and chaotic motion. Specifically, system behavior switches back and forth between periodic and quasi-periodic when 0 < ω < 0.778. In the interval of 0.778 < ω < 2.208, only chaotic motion is fond. But when ω is beyond 2.208, the system alternates again between periodic and quasi-periodic motions and the chance of finding quasi-periodic motion is much more than 0 < ω < 0.778.

We also study the non-linear behavior of system varying amplitude of electromagnetic radiation A at three ω values. Very different scenarios are displayed. Periodic, quasi-periodic, or chaotic motion type characterize the scenario at different ω value, respectively.

The route to chaos, induced by Hopf bifurcation, is discovered varying both frequency ω of electromagnetic radiation φ_ext and amplitude of electromagnetic radiation A. However, it seems ω plays a much more important role in controlling the system state. In summary, system behavior is very sensitive to frequency of external force, which might be related to the ratio of frequencies of two external plane waves.

Author Contributions

All authors listed have made a substantial, direct and intellectual contribution to the work, and approved it for publication.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

This work is supported by the National Natural Science Foundation of China (Grant Nos. 11472202) and the National Key Basic Research Program of China (973 Program) (No. 2012CB026002).

References

Abbasbandy, S. (2008). Soliton solutions for the Fitzhugh–Nagumo equation with the homotopy analysis method. Appl. Math. Modell. 32, 2706–2714. doi: 10.1016/j.apm.2007.09.019

CrossRef Full Text | Google Scholar

Dehghan, M., Manafian, J., and Saadatmandi, S. (2010). Application of semi analytic methods for the Fitzhugh–Nagumo equation, which models the transmission of nerve impulses. Math. Methods Appl. Sci. 33, 1384–1398. doi: 10.1002/mma.1329

CrossRef Full Text | Google Scholar

Duan, L., Lu, Q., and Wang, Q. (2008). Two-parameter bifurcation analysis of firing activities in the Chay neuronal model. Neurocomputing 72, 341–351. doi: 10.1016/j.neucom.2008.01.019

CrossRef Full Text | Google Scholar

Fitzhugh, R. (1961). Impulses and physiological states in theoretical models of nerve membrane. Biophys. J. 1, 445–466. doi: 10.1016/S0006-3495(61)86902-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Garcaojalvo, J., and Schimanskygeier, L. (1999). Noise induced spiral dynamics in excitable media. Europhys. Lett. 47, 298–303. doi: 10.1209/epl/i1999-00388-9

CrossRef Full Text | Google Scholar

Hindmarsh, J. L., and Rose, R. M. (1982). A model of the nerve impulse using two first-order differential equations. Nature 296, 162–164. doi: 10.1038/296162a0

PubMed Abstract | CrossRef Full Text | Google Scholar

Hodgkin, A. L., and Huxley, A. F. (1990). A quantitative description of membrane current and its application to conduction and excitation in nerve. J. Physiol. 52, 25–71.

PubMed Abstract | Google Scholar

Hossmann, K. A., and Hermann, D. M. (2003). Effects of electromagnetic radiation of mobile phones on the central nervous system. Bioelectromagnetics 24, 49–62. doi: 10.1002/bem.10068

PubMed Abstract | CrossRef Full Text | Google Scholar

Ibarz, B., Casado, J. M., and Sanjuán, M. A. F. (2011). Map-based models in neuronal dynamics. Physics Rep. 501, 1–74. doi: 10.1016/j.physrep.2010.12.003

CrossRef Full Text | Google Scholar

Krupa, M., Sandstede, B., and Szmolyan, P. (1997). Fast and slow waves in the Fitzhugh–Nagumo equation. J. Diff. Equat. 133, 49–97. doi: 10.1006/jdeq.1996.3198

CrossRef Full Text | Google Scholar

Li, H., and Guo, Y. (2006). New exact solutions to the Fitzhugh–Nagumo equation. Appl. Math. Comput. 180, 524–528. doi: 10.1016/j.amc.2005.12.035

CrossRef Full Text | Google Scholar

Lisi, A., Ciotti, M. T., Ledda, M., Pieri, M., Zona, C., Mercanti, D., et al. (2005). Exposure to 50 Hz electromagnetic radiation promote early maturation and differentiation in newborn rat cerebellar granule neurons. J. Cell. Physiol. 204, 532–538. doi: 10.1002/jcp.20322

PubMed Abstract | CrossRef Full Text | Google Scholar

Lv, M., Wang, C., Ren, G., Ma, J., and Song, X. (2016). Model of electrical activity in a neuron under magnetic flow effect. Nonlin. Dyn. 85, 1479–1490. doi: 10.1007/s11071-016-2773-6

CrossRef Full Text | Google Scholar

Ma, J., and Tang, J. (2017). A review for dynamics in neuron and neuronal network. Nonlin. Dyn., 89, 1569–1578. doi: 10.1007/s11071-017-3565-3

CrossRef Full Text | Google Scholar

Ma, J., Wang, Y., Wang, C., Xu, Y., and Ren, G. (2017). Mode selection in electrical activities of myocardial cell exposed to electromagnetic radiation. Chaos Solit. Fract. 99, 219–225. doi: 10.1016/j.chaos.2017.04.016

CrossRef Full Text | Google Scholar

Mvogo, A., Takembo, C. N., Fouda, H. P. E., and Kofane, T. C. (2017). Patterns formation in diffusive excitable systems under magnetic flow effects. Phys. Lett. A 381, 2264–2271. doi: 10.1016/j.physleta.2017.05.020

CrossRef Full Text | Google Scholar

Nagumo, J., Arimoto, S., and Yoshizawa, S. (1962). An active pulse transmission line simulating nerve axon. Proc. IRE 50, 2061–2070. doi: 10.1109/JRPROC.1962.288235

CrossRef Full Text | Google Scholar

Nozaki, D., and Yamamoto, Y. (1998). Enhancement of stochastic resonance in a Fitzhugh–Nagumo neuronal model driven by colored noise. Phys. Lett. A 243, 281–287. doi: 10.1016/S0375-9601(98)00247-3

CrossRef Full Text | Google Scholar

Pankratova, E. V., Polovinkin, A. V., and Spagnolo, B. (2005). Suppression of noise in Fitzhugh–Nagumo model driven by a strong periodic signal. Phys. Lett. A 344, 43–50. doi: 10.1016/j.physleta.2005.05.099

CrossRef Full Text | Google Scholar

Ruelle, D., and Takens, F. (1971). Note concerning our paper: “on the nature of turbulence”. Commun. Math. Phys. 23, 343–344. doi: 10.1007/BF01893621

CrossRef Full Text

Song, Z., Wang, C., and Zhen, B. (2016). Codimension two bifurcation and multistability coexistence in an inertial two neuron system with multiple delays. Nonlin. Dyn. 85, 1–15. doi: 10.1007/s11071-016-2816-z

CrossRef Full Text | Google Scholar

Song, Z., Xu, J., and Zhen, B. (2015a). Multitype activity coexistence in an inertial two neuron system with multiple delays. Int. J. Bifurc. Chaos 25:1530040. doi: 10.1142/S0218127415300402

CrossRef Full Text | Google Scholar

Song, Z., Yang, K., Xu, J., and Wei, Y. (2015b). Multiple pitchfork bifurcation and multi periodicity coexistences in a delay coupled neural oscillator system with inhibitory to inhibitory connection. Commun. Nonlin. Sci. Numer. Simulat. 29, 327–345. doi: 10.1016/j.cnsns.2015.05.013

CrossRef Full Text | Google Scholar

Triki, H., and Wazwaz, A. M. (2013). On soliton solutions for the Fitzhugh–Nagumo equation with time dependent coefficients. Appl. Math. Model. 37, 3821–3828. doi: 10.1016/j.apm.2012.07.031

CrossRef Full Text | Google Scholar

Wang, H. X., Wang, Q. Y., and Zheng, Y. H. (2014). Bifurcation analysis for hindmarsh-rose neuronal model with time-delayed feedback control and application to chaos control. Sci. China Tech. Sci. 57, 872–878. doi: 10.1007/s11431-014-5535-z

CrossRef Full Text | Google Scholar

Xu, S., Ning, W., Xu, Z., Zhou, S., Chiang, H., and Luo, J. (2006). Chronic exposure to GSM 1800 MHz microwaves reduces excitatory synaptic activity in cultured hippocampal neurons. Neurosci. Lett. 398, 253–257. doi: 10.1016/j.neulet.2006.01.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Xu, S., Zhou, Z., Zhang, L., Yu, Z., Zhang, W., Wang, Y., et al. (2010). Exposure to 1800 MHz radiofrequency radiation induces oxidative damage to mitochondrial DNA in primary cultured neurons. Brain Res. 1311, 189–196. doi: 10.1016/j.brainres.2009.10.062

PubMed Abstract | CrossRef Full Text | Google Scholar

Ying, W., Li, J., Liu, S., Pang, J., Du, M., and Pan, L. (2013). Noise induced spatiotemporal patterns in Hodgkin–Huxley neuronal network. Cogn. Neurodyn. 7, 431–440. doi: 10.1007/s11571-013-9245-1

CrossRef Full Text | Google Scholar