ISRN BiomathematicsVolume 2012 (2012), Article ID 546315, 7 pagesdoi:10.5402/2012/546315

Research Article

Transient Periodicity in a Morris-Lecar Neural System

Sreenivasan Rajamoni Nadar¹ and Vikas Rai²

¹National Institute of Mental Health, NIH, Bethesda, MD 20892, USA²Department of Mathematics, Faculty of Sciences, Jazan University, Jazan 45142, Saudi Arabia

Received 6 March 2012; Accepted 17 April 2012

Academic Editors: G. Bocharov, K. Kitano, J. R. C. Piqueira, and M. Santillán

Copyright © 2012 Sreenivasan Rajamoni Nadar and Vikas Rai. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

The dynamical complexity of a system of ordinary differential equations (ODEs) modeling the dynamics of a neuron that interacts with other neurons through on-off excitatory and inhibitory synapses in a neural system was investigated in detail. The model used Morris-Lecar (ML) equations with an additional autonomous variable representing the input from interaction of excitatory neuronal cells with local interneurons. Numerical simulations yielded a rich repertoire of dynamical behavior associated with this three-dimensional system, which included periodic, chaotic oscillation and rare bursts of episodic periodicity called the transient periodicity.

1. Introduction

In this paper, we studied a three-dimensional extension of the Morris-Lecar (ML) system [1, 2]. The ML system is a two-dimensional representation of the four-variable Hodgkin-Huxley system [3]. It assumes that the conductivity of Na⁺ ions is very high. The ML system supports only two types of dynamics, resting and sustained spiking activity and does not exhibit chaotic oscillations. In this paper, we explored the dynamics of a model system that had an extra variable added to the ML system. The model described a condition in which excitatory principal cells interacted with local inhibitory inter neurons. These cells are present in a small volume of physiological space. Hence, the effect of the spatial structure of individual cells on the model dynamics was ignored.

The emergent behavior of a population of neurons was studied by considering population of neurons as the basic unit of an excitatory pathway connecting to inhibitory interneurons. The feedback inhibition present in (1) and (3) is an important part of the system architecture. Thus, by defining a population of neurons as a dynamical system [4], we studied the temporal behavior of neurons interacting with each other as well as with the other population. The - subsystem was governed by a faster time scale compared to that of state variable . All interesting behavior of the model system resulted from an interaction between the faster subsystem and a slower one with the population of inhibitory inter-neurons. The differential equations that represent the dynamical system were derived by considering the behavior of prototypical single neurons (ML neurons).

The mean field equations describing the time evolution of neural system are given as follows: where and are the mean membrane potentials for the excitatory and inhibitory cells. is the fraction of open potassium channels in the population of principal cells. and are the Nernst potentials for the potassium ions and leakage channels. , a voltage dependent constant, is given by where and are the midpoint potentials at which the potassium and calcium currents are half-activated. The and are the total conductances for the population of calcium and potassium channels, respectively. The constant denotes “leakage conductance,” and it measures the effect of the channels that are always open. The parameter represents the relative strength of current into the inhibitory interneurons to the current into principal cells. The parameter , the applied current, is maintained as a constant stimulus. The changing polarity of the mean membrane potential of the principal cells ensured that the mean membrane potential of the inhibitory cells did not grow beyond manageable limit.

The important variables in the model are and are nondimensionalized functions to describe the voltage-regulated calcium and potassium ion channels, respectively. The parameters and are temperature scaling factors. The hyperbolic tangent function describes the collective behavior of a large number of channels. On the right-hand sides of (7) and (8), and are dimensionless synaptic strengths of excitatory principal cells and inhibitory inter-neurons. On the left of these equations are functions that are parameterized by these dimensionless constants. These functions were assumed to be of tangent hyperbolic form. is the steepness parameter for both the functions. and are the thresholds for the synaptic strength functions of the excitatory principal cells and inhibitory interneurons. Both the functions have the same functional properties. The negative sign in (1) signifies an inhibitory effect from the interneurons with mean membrane potential on the principal cells. Parameter values for , , , and are chosen in such a way that nondimensionalized functions and reach their equilibrium values almost instantaneously. The aim of the paper is to report the rich dynamical repertoire displayed by the system of ODEs.

The key parameters of the model are , , and . The parameter ensures that the inherent time scale of the group of inhibitory cells is different from the rest of the system. The parameter measures the strength of feedforward inhibition. When (a fixed threshold potential), the connection between the principal cell and the interneuron is opened. The parameter represents constant external current. It turned out to be useful for studying the effects of disinhibition factors on model dynamics. Yet another key parameter of this model system is , which controls the steepness of the “on-off” switch in and . For a single presynaptic terminal, excitatory or inhibitory, one would expect “all or none” behavior. In that case, one could have a simple representation of the on/off switch via a Heaviside step function. In a population of cells, the connection strengths vary smoothly as distributions of thresholds exist and the values of parameters and represent the end points of this distribution.

Inhibitory inter-neurons are found in many brain regions. The ML system does not consider the effect of the inhibitory inter-neurons on the excitatory principal cells. Since it is based on two primary system variables, the M-L system cannot support deterministic chaos. We explored regular and chaotic dynamics, which resulted from an interaction of principal cells with inhibitory inter-neurons. When the values for parameters of the M-L system of equations (as described above) were suitably chosen, we observed Ca⁺⁺ and K⁺ oscillations. The action of inhibitory inter-neurons on the excitatory system (equations (1) without the last term and (2)) generated chaotic dynamics as the third variable (Z) was added to the system. The base parameter values were chosen in such a way that the M-L system displayed sustained spiking activity, which has a stable limit cycle as the underlying attractor. The remaining parameters were set in such a way that dynamics of the model system unfolded on a chaotic attractor. The base parameter set for which the model system (equations (1), (2) and (3) displayed oscillations is given in Table 1.

Table 1: Parameter values and units for the Morris-Lecar system of equations.

Parameter
Value
Units

0.01
mV

0.15
mV

0.03
mV

0.3
mV

0.0
mV

0.4
mV

0.05
mV

0.7
mV

0.5
mV

0.4
Dimensionless

0.3
μA/cm²

0.15
1/ms

0.238
ohm cm²

1.1
mmho/cm²

2.0
mmho/cm²

1.0
mmho/cm²

1.0
Dimensionless

1.0
Dimensionless

The Morris-Lecar equations are applicable to a spatially isopotential patch of membrane injected with a constant current. Transition from single shot firing to a stable limit cycle (tonic firing) occurred as the excitation current was varied. At very low or very high values of the stimulus, it had a single stable fixed point. The neurons of the inhibitory system provided an additional degree of freedom as it is well known that a minimum of three degrees of freedom are required for a system to support complex dynamics including chaos. We chose the set of values presented in Table 1 for the system parameters as base values. At these values of parameters the system displayed oscillatory behavior. Table 1 gives the chosen parameter values.

2. Simulations

Two-dimensional parameter scans were performed for () parameter space with in the range of 0.1 to 0.15 with a step size of 0.01 and in the range of −0.9 to −0.6 with a step size of 0.05, and the results of these scans are presented in Figure 1. The figure plots points at which the system of differential equations exhibited chaotic solutions and nearby points in the system exhibited other dynamical behavior. We used the characteristic signature of deterministic chaos: sensitive dependence on initial conditions to detect chaos in the system. Chaotic dynamics was detected when two time series (collected after transients were allowed to die out) generated at the same set of parameter values, but, for two different initial conditions which were nearby, did not overlap. Computations of the largest Lyapunov exponents were also performed to confirm the regular or chaotic behavior observed at a chosen set of parameters. Figures 2(a) and 2(b) demonstrate the onset of chaos through period doubling route. We have also examined bifurcation diagrams varying the parameter and confirmed that the system exhibits period doubling route to chaos. This portrays sustained spiking activity displayed by the neuronal system for a given combination of temperature scaling factor and the Nernst potential for potassium in the node. The chaotic oscillations of the system are shown in Figures 3(a) and 3(b), and the largest Lyapunov exponent for the attractor is 0.0621.

Figure 1: Points in the two-dimensional parameter space where deterministic chaos was observed.

(a)

(b)

Figure 2: (a) Bifurcation diagram where is the bifurcation parameter and . (b) Bifurcation diagram where is the bifurcation parameter and .

(a)

(b)

Figure 3: (a) The strange chaotic attractor at and . (b) The time histories for two nearby initial conditions. The nonoverlap between the two signals and positive Lyapunov establishes that the underlying attractor supports chaotic dynamics.

The dynamical behavior of the system was sensitive to changes in both the parameters. It changed from simple periodic oscillations to complex periodic ones and at times transitioned to chaos. Two-dimensional parameter scans showed how sensitive the dynamics of the model system was when it encountered a change in an intrinsic attribute of the system. At certain points (e.g., 0.12, −0.7), the system exhibited bistable behavior; a periodic attractor coexisted with a chaotic attractor. A few points represented bistability, especially those at the outer edges of the parameter regions, and chaotic behavior was observed at discrete points in closed intervals.

The other two-dimensional parameter scan was carried out in () with range from 0.10 to 0.15 with a step size 0.01 and range from 0.3 to 0.8 with a step size 0.1. We found that chaotic behavior exists at (0.1, 0.40), (0.11, 0.3), (0.12, 0.3), (0.13, 0.3), (0.14, 0.4), (0.15, 0.35), and (0.15, 0.4). At a few points, neuronal dynamics displayed resting behavior.

Figures 4(a) and 4(b) show stable limit cycle attractor and corresponding time-series at the point where bi-stability was detected.

(a)

(b)

Figure 4: (a) The period-2 cycle found at and . (b) The corresponding time series obtained for two close initial conditions.

The third parameter scan was performed in () with the following ranges and step sizes. The range for was from 0.165 to 0.25 with step size 0.002, and for the range was from −0.85 to −0.7 with a step size of 0.05. Chaos was observed at the following points: (0.236, −0.85) (0.238, −0.7), (0.238, −0.85), (0.238, −0.80), and (0.24, −0.80). Figures 4 and 5 show period-2 cycle and suggest that the system exhibited period doubling route to chaos.

(a)

(b)

Figure 5: (a) Period-2 cycle. The attractor was detected at and . Other parameters were fixed at their base values. (b) Time series for two nearby initial conditions.

The fourth and the final parameter scan was performed in (), and it showed chaos at the following points: (1, 0.9), (1.0, 1.1), (1.1, 0.8), (1.1, 0.9), (1.1, 1.0), (1.2, 0.8), (1.2, 0.9), (1.2, 1.0), (1.3, 0.8), (1.3, 0.9), and (1.3, 1.0). Figure 6 shows the presence of “edge of chaos” in the neural system. Edge of chaos is a dynamical phenomenon that is characterized by the existence of chaotic behavior in finite but arbitrarily small parameter ranges. This suggests that smooth changes in the strengths of excitatory and inhibitory synapses cause dynamical transitions from chaotic behavior to regular oscillatory dynamics. The model system also show chaotic saddle which is shown in Figure 7.

Figure 6: Points in two-dimensional parameter space () where system dynamics displays chaotic behavior. The system supports regular equilibrium dynamics and sustained periodic oscillations at intervening points.

Figure 7: Two-dimensional view in () plane of the chaotic saddle at .12, .

Results of time series comparisons were verified by computing the largest Lyapunov exponents () using the Wolf et al. algorithm [5]. We present with standard error computed from 20 trials of random initial conditions of the system in Table 2. The choice of the delay for reconstruction of the dynamics from the time series was made by choosing the first minimum of mutual information.

Table 2: Values of the largest Lyapunov exponents.

1
0.9

1
1.1

1.1
0.8

1.1
0.9
0.052 0.005

1.1
1.0
0.069 0.003

1.2
0.8
0.031 0.005

1.2
0.9
0.053 0.006

1.2
1.0

1.3
0.8

1.3
0.9

1.3
1.0

We present an interesting finding of transient periodicity displayed by this neural system, which was detected in the () parameter space which is shown in Figure 8.

(a)

(b)

Figure 8: (a) Two-dimensional view of the chaotic saddle at , in plane. (b) Two signals starting from two nearby initial conditions. The non-overlap between the two (except in the time domain (700–900)) confirms that the dynamical behavior displayed by the system is chaos. The transient periodicity manifests itself in the time domain 700–900.

3. Results

Chaotic dynamics is characterized by sensitive dependence on initial conditions. The error committed in fixing initial conditions grows exponentially with respect to time. The largest Lyapunov exponent is a quantity that characterizes the rate of separation of infinitesimally close trajectories. The time series is termed as chaotic if it returns a positive value for the Lyapunov exponent. Computer simulations of our model system suggested that deterministic chaos manifests itself as short-term recurrent phenomenona.

We discovered rare bursts of “periodic” motion at , (cf. 8b). This phenomenon is known as “transient periodicity” in the dynamical systems literature. This is caused by the presence of semiperiodic saddles (unstable invariant sets) in the chaotic attractor. Motion in the vicinity of these sets is periodic. In the case of EEG time series obtained by recording averaged action potentials of epochs of epileptic seizure [6], burst pattern often contains nearly periodic dynamics. In chaotic systems, such transient periodicity can reflect the existence of semi-periodic saddles contained in the attractor. Motion in the vicinity of such objects has a prominent periodic component as trajectories are temporarily attracted to these neighborhoods before exiting.

4. Discussion

Model simulations and animal experiments suggest that seizure activity occurs when a critical mass of neurons join and get involved in time-linked high frequency discharging [4]. An asynchronous-to-synchronous transition (longer-time state) is caused either by deterministic changes in system parameters or spontaneous dynamical shifts from chaos to order caused by transient periodicity in chaotic saddles. This situation has a parallel in our model system wherein the dynamics switched to sustained spiking activity or to resting behavior from chaotic behavior when small changes in crucial parameters occurred. This dynamical behavior is known as “edge of chaos.” Chaotic saddles are semiperiodic semiattractors that are distinct from the better known attractors in nonlinear dynamical systems which are created during boundary crisis [7]. Presence of these unstable invariant sets generates strong episodic periodic components embedded in chaotic dynamics. This has been found in epidemiological models and also in EEG time series from patients of epilepsy [8, 9].

We envision that chaos is vital for functioning of a healthy brain and synchronization of the neuronal system occurs when all its regions are in transient periodicity represented by chaotic saddles in state space. This is how the intermittent pathology of epileptic seizures is created. It is understood that lack of inhibition [6, 10] caused either by a reduction in the number density of excitatory cells that synapse on the inhibitory cells in the focal region or by the diminished release of neurotransmitter -amino butyric acid in a nearby region leads to seizures (transient periodicity). The ictal period is represented by periodic oscillations and the interictal period by chaotic phase. Stable limit cycles and complex periodic dynamics observed in a narrow range of model parameter, , correspond to a biophysical state of partial dis-inhibition. Biophysical experiments performed previously [11, 12] have attempted to simulate dynamics observed in this range. The transient periodicity, which was detected for higher values of this parameter, has been observed in electrographic recordings from depth and subdural electrodes performed in epilepsy patients [8].

Acknowledgments

The authors thank Dr. Pravitha Ramanand for discussions and Professor Raima Larter for making available a few chapters of B. Speelman's doctoral thesis.

References

C. Morris and H. Lecar, “Voltage oscillations in the barnacle giant muscle fiber,” Biophysical Journal, vol. 35, no. 1, pp. 193–213, 1981.

J. Rinzel and B. Ermentrout, “Analysis of neural excitability and oscillations,” in Methods of Neural Modelling, C. Koch and I. Segev, Eds., pp. 251–292, MIT Press, Cambridge, Mass, USA, 1998.

A. L. Hodgkin and A. F. Huxley, “A quantitative description of membrane current and its application to conduction and excitation in nerve,” The Journal of Physiology, vol. 117, no. 4, pp. 500–544, 1952.

B. Speelman, A dynamical systems approach to the modeling of epileptic seizures [Ph.D. Thesis], Indiana University Perdue University at Indiana Polis, 1997.

A. Wolf, J. B. Swift, H. L. Swinney, and J. A. Vastano, “Determining Lyapunov exponents from a time series,” Physica D, vol. 16, no. 3, pp. 285–317, 1985.

I. Pavlov and D. M. Kullmann, “How much inhibition in an epileptiform burst?” Journal of Physiology, vol. 588, no. 1, pp. 17–18, 2010.

Y. Lai and C. Grebogi, “Converting transient chaos into sustained chaos by feedback control,” Physical Review E, vol. 49, no. 2, pp. 1094–1098, 1994.

W. M. Schaffer, B. E. Kendall, C. W. Tidd, and L. F. Olsen, “Transient periodicity and episodic predictability in biological dynamics,” Mathematical Medicine and Biology, vol. 10, no. 4, pp. 227–247, 1993.

B. E. Kendall, W. M. Schaffer, and C. W. Tidd, “Transient periodicity in chaos,” Physics Letters A, vol. 177, no. 1, pp. 13–20, 1993.

I. Marchionni and G. Maccaferri, “Quantitative dynamics and spatial profile of perisomatic GABAergic input during epileptiform synchronization in the CA1 hippocampus,” Journal of Physiology, vol. 587, no. 23, pp. 5691–5708, 2009.

R. S. Sloviter, “Decreased hippocampal inhibition and a selective loss of interneurons in experimental epilepsy,” Science, vol. 235, no. 4784, pp. 73–76, 1987.

J. W. Swann and J. J. Hablitz, “Cellular abnormalities and synaptic plasticity in seizure disorders of the immature nervous system,” Mental Retardation and Developmental Disabilities Research Reviews, vol. 6, no. 4, pp. 258–267, 2000.