Slow negative feedback enhances robustness of square-wave bursting.

Square-wave bursting is an activity pattern common to a variety of neuronal and endocrine cell models that has been linked to central pattern generation for respiration and other physiological functions. Many of the reduced mathematical models that exhibit square-wave bursting yield transitions to an alternative pseudo-plateau bursting pattern with small parameter changes. This susceptibility to activity change could represent a problematic feature in settings where the release events triggered by spike production are necessary for function. In this work, we analyze how model bursting and other activity patterns vary with changes in a timescale associated with the conductance of a fast inward current. Specifically, using numerical simulations and dynamical systems methods, such as fast-slow decomposition and bifurcation and phase-plane analysis, we demonstrate and explain how the presence of a slow negative feedback associated with a gradual reduction of a fast inward current in these models helps to maintain the presence of spikes within the active phases of bursts. Therefore, although such a negative feedback is not necessary for burst production, we find that its presence generates a robustness that may be important for function.

Saddle Slow Manifolds and Canard Orbits in [Formula: see text] and Application to the Full Hodgkin-Huxley Model.

Many physiological phenomena have the property that some variables evolve much faster than others. For example, neuron models typically involve observable differences in time scales. The Hodgkin-Huxley model is well known for explaining the ionic mechanism that generates the action potential in the squid giant axon. Rubin and Wechselberger (Biol. Cybern. 97:5-32, 2007) nondimensionalized this model and obtained a singularly perturbed system with two fast, two slow variables, and an explicit time-scale ratio ε. The dynamics of this system are complex and feature periodic orbits with a series of action potentials separated by small-amplitude oscillations (SAOs); also referred to as mixed-mode oscillations (MMOs). The slow dynamics of this system are organized by two-dimensional locally invariant manifolds called slow manifolds which can be either attracting or of saddle type.In this paper, we introduce a general approach for computing two-dimensional saddle slow manifolds and their stable and unstable fast manifolds. We also develop a technique for detecting and continuing associated canard orbits, which arise from the interaction between attracting and saddle slow manifolds, and provide a mechanism for the organization of SAOs in [Formula: see text]. We first test our approach with an extended four-dimensional normal form of a folded node. Our results demonstrate that our computations give reliable approximations of slow manifolds and canard orbits of this model. Our computational approach is then utilized to investigate the role of saddle slow manifolds and associated canard orbits of the full Hodgkin-Huxley model in organizing MMOs and determining the firing rates of action potentials. For ε sufficiently large, canard orbits are arranged in pairs of twin canard orbits with the same number of SAOs. We illustrate how twin canard orbits partition the attracting slow manifold into a number of ribbons that play the role of sectors of rotations. The upshot is that we are able to unravel the geometry of slow manifolds and associated canard orbits without the need to reduce the model.

Invariant manifolds and global bifurcations.

Invariant manifolds are key objects in describing how trajectories partition the phase spaces of a dynamical system. Examples include stable, unstable, and center manifolds of equilibria and periodic orbits, quasiperiodic invariant tori, and slow manifolds of systems with multiple timescales. Changes in these objects and their intersections with variation of system parameters give rise to global bifurcations. Bifurcation manifolds in the parameter spaces of multi-parameter families of dynamical systems also play a prominent role in dynamical systems theory. Much progress has been made in developing theory and computational methods for invariant manifolds during the past 25 years. This article highlights some of these achievements and remaining open problems.

Adaptive topographies and equilibrium selection in an evolutionary game.

It has long been known in the field of population genetics that adaptive topographies, in which population equilibria maximise mean population fitness for a trait regardless of its genetic bases, do not exist. Whether one chooses to model selection acting on a single locus or multiple loci does matter. In evolutionary game theory, analysis of a simple and general game involving distinct roles for the two players has shown that whether strategies are modelled using a single 'locus' or one 'locus' for each role, the stable population equilibria are unchanged and correspond to the fitness-maximising evolutionary stable strategies of the game. This is curious given the aforementioned population genetical results on the importance of the genetic bases of traits. Here we present a dynamical systems analysis of the game with roles detailing how, while the stable equilibria in this game are unchanged by the number of 'loci' modelled, equilibrium selection may differ under the two modelling approaches.

Solving Winfree's puzzle: the isochrons in the FitzHugh-Nagumo model.

We consider the FitzHugh-Nagumo model, an example of a system with two time scales for which Winfree was unable to determine the overall structure of the isochrons. An isochron is the set of all points in the basin of an attracting periodic orbit that converge to this periodic orbit with the same asymptotic phase. We compute the isochrons as one-dimensional parametrised curves with a method based on the continuation of suitable two-point boundary value problems. This allows us to present in detail the geometry of how the basin of attraction is foliated by isochrons. They exhibit extreme sensitivity and feature sharp turns, which is why Winfree had difficulties finding them. We observe that the sharp turns and sensitivity of the isochrons are associated with the slow-fast nature of the FitzHugh-Nagumo system; more specifically, it occurs near its repelling (unstable) slow manifold.

Continuation-based numerical detection of after-depolarization and spike-adding thresholds.

The changes in neuronal firing pattern are signatures of brain function, and it is of interest to understand how such changes evolve as a function of neuronal biophysical properties. We address this important problem by the analysis and numerical investigation of a class of mechanistic mathematical models. We focus on a hippocampal pyramidal neuron model and study the occurrence of bursting related to the after-depolarization (ADP) that follows a brief current injection. This type of burst is a transient phenomenon that is not amenable to the classical bifurcation analysis done, for example, for periodic bursting oscillators. In this letter, we show how to formulate such transient behavior as a two-point boundary value problem (2PBVP), which can be solved using well-known continuation methods. The 2PBVP is formulated such that the transient response is represented by a finite orbit segment for which onsets of ADP and additional spikes in a burst can be detected as bifurcations during a one-parameter continuation. This in turn provides us with a direct method to approximate the boundaries of regions in a two-parameter plane where certain model behavior of interest occurs. More precisely, we use two-parameter continuation of the detected onset points to identify the boundaries between regions with and without ADP and bursts with different numbers of spikes. Our 2PBVP formulation is a novel approach to parameter sensitivity analysis that can be applied to a wide range of problems.

Geometric analysis of transient bursts.

We consider the effect of a brief stimulation from the rest state of a minimal neuronal model with multiple time scales. Such transient dynamics brings out the intrinsic bursting capabilities of the system. Our main goal is to show that a minimum of three dimensions is enough to generate spike-adding phenomena in transient responses, and that the onset of a new spike can be tracked using existing continuation packages. We take a geometric approach to illustrate how the underlying fast subsystem organises the spike adding in much the same way as for spike adding in periodic bursts, but the bifurcation analysis for spike onset is entirely different. By using a generic model, we further strengthen claims made in our earlier work that our numerical method for spike onset can be used for a broad class of systems.

CROSS-CURRENTS BETWEEN BIOLOGY AND MATHEMATICS: THE CODIMENSION OF PSEUDO-PLATEAU BURSTING.

A great deal of work has gone into classifying bursting oscillations, periodic alternations of spiking and quiescence modeled by fast-slow systems. In such systems, one or more slow variables carry the fast variables through a sequence of bifurcations that mediate transitions between oscillations and steady states. A rigorous classification approach is to characterize the bifurcations found in the neighborhood of a singularity; a measure of the complexity of the bursting oscillation is then given by the smallest codimension of the singularities near which it occurs. Fold/homoclinic bursting, along with most other burst types of interest, has been shown to occur near a singularity of codimension three by examining bifurcations of a cubic Liénard system; hence, these types of bursting have at most codimension three. Modeling and biological considerations suggest that fold/homoclinic bursting should be found near fold/subHopf bursting, a more recently identified burst type whose codimension has not been determined yet. One would expect that fold/subHopf bursting has the same codimension as fold/homoclinic bursting, because models of these two burst types have very similar underlying bifurcation diagrams. However, no codimension-three singularity is known that supports fold/subHopf bursting, which indicates that it may have codimension four. We identify a three-dimensional slice in a partial unfolding of a doubly-degenerate Bodganov-Takens point, and show that this codimension-four singularity gives rise to almost all known types of bursting.

Dynamical systems analysis of spike-adding mechanisms in transient bursts.

Transient bursting behaviour of excitable cells, such as neurons, is a common feature observed experimentally, but theoretically, it is not well understood. We analyse a five-dimensional simplified model of after-depolarisation that exhibits transient bursting behaviour when perturbed with a short current injection. Using one-parameter continuation of the perturbed orbit segment formulated as a well-posed boundary value problem, we show that the spike-adding mechanism is a canard-like transition that has a different character from known mechanisms for periodic burst solutions. The biophysical basis of the model gives a natural time-scale separation, which allows us to explain the spike-adding mechanism using geometric singular perturbation theory, but it does not involve actual bifurcations as for periodic bursts. We show that unstable sheets of the critical manifold, formed by saddle equilibria of the system that only exist in a singular limit, are responsible for the spike-adding transition; the transition is organised by the slow flow on the critical manifold near folds of this manifold. Our analysis shows that the orbit segment during the spike-adding transition includes a fast transition between two unstable sheets of the slow manifold that are of saddle type. We also discuss a different parameter regime where the presence of additional saddle equilibria of the full system alters the spike-adding mechanism.

A unified model of CA1/3 pyramidal cells: an investigation into excitability.

After-depolarisation is a hallmark of excitability in hippocampal pyramidal cells of CA1 and CA3 regions, because it constitutes the subthreshold relation between inward and outward ionic currents. This relationship determines the nominal response to stimuli and provides the necessary conditions for firing a spike or a burst of action potentials. Nevertheless, after-depolarisation is an inherently transient phenomenon that is not very well understood. We study after-depolarisation using a single-compartment pyramidal-cell model based on recent voltage- and current-clamp experimental data. We systematically investigate CA1 and CA3 behaviour and show that changes to maximal conductances of T-type Ca(2+)-current and muscarinic-sensitive and delayed rectifier K(+)-currents are sufficient to switch the behaviour of the model from a CA3 to a CA1 neuron. We use model analysis to define after-depolarisation and bursting threshold. We also explain the influence of particular ionic currents on this phenomenon. This study ends with a sensitivity analysis that demonstrates the influence of specific currents on excitability. Counter-intuitively, we find that a decrease of Na(+)-current could cause an increase in excitability. Our analysis suggests that a change of high-voltage activated Ca(2+)-current can have a similar effect.

Modeling mechanisms of cell secretion.

Secretion is a fundamental cellular process involving the regulated release of intracellular products from cells. Physiological functions such as neurotransmission, or the release of hormones and digestive enzymes, are all governed by cell secretion. Anomalies in the processes involved in secretion contribute to the development and progression of diseases such as diabetes and other hormonal disorders. To unravel the mechanisms that govern such diseases, it is essential to understand how hormones, growth factors and neurotransmitters are synthesized and processed, and how their signals are recognized, amplified and transmitted by intracellular signaling pathways in the target cells. Here, we discuss diverse aspects of the detailed mechanisms involved in secretion based on mathematical models. The models range from stochastic ones describing the trafficking of secretory vesicles to deterministic ones investigating the regulation of cellular processes that underlie hormonal secretion. In all cases, the models are closely related to experimental results and suggest theoretical predictions for the secretion mechanisms.

Full system bifurcation analysis of endocrine bursting models.

Plateau bursting is typical of many electrically excitable cells, such as endocrine cells that secrete hormones and some types of neurons that secrete neurotransmitters. Although in many of these cell types the bursting patterns are regulated by the interplay between voltage-gated calcium channels and calcium-sensitive potassium channels, they can be very different. We investigate so-called square-wave and pseudo-plateau bursting patterns found in endocrine cell models that are characterized by a super- or subcritical Hopf bifurcation in the fast subsystem, respectively. By using the polynomial model of Hindmarsh and Rose (Proceedings of the Royal Society of London B 221 (1222) 87-102), which preserves the main properties of the biophysical class of models that we consider, we perform a detailed bifurcation analysis of the full fast-slow system for both bursting patterns. We find that both cases lead to the same possibility of two routes to bursting, that is, the criticality of the Hopf bifurcation is not relevant for characterizing the route to bursting. The actual route depends on the relative location of the full-system's fixed point with respect to a homoclinic bifurcation of the fast subsystem. Our full-system bifurcation analysis reveals properties of endocrine bursting that are not captured by the standard fast-slow analysis.

Understanding anomalous delays in a model of intracellular calcium dynamics.

In many cell types, oscillations in the concentration of free intracellular calcium ions are used to control a variety of cellular functions. It has been suggested [J. Sneyd et al., "A method for determining the dependence of calcium oscillations on inositol trisphosphate oscillations," Proc. Natl. Acad. Sci. U.S.A. 103, 1675-1680 (2006)] that the mechanisms underlying the generation and control of such oscillations can be determined by means of a simple experiment, whereby a single exogenous pulse of inositol trisphosphate (IP(3)) is applied to the cell. However, more detailed mathematical investigations [M. Domijan et al., "Dynamical probing of the mechanisms underlying calcium oscillations," J. Nonlinear Sci. 16, 483-506 (2006)] have shown that this is not necessarily always true, and that the experimental data are more difficult to interpret than first thought. Here, we use geometric singular perturbation techniques to study the dynamics of models that make different assumptions about the mechanisms underlying the calcium oscillations. In particular, we show how recently developed canard theory for singularly perturbed systems with three or more slow variables [M. Wechselberger, "A propos de canards (Apropos canards)," Preprint, 2010] applies to these calcium models and how the presence of a curve of folded singularities and corresponding canards can result in anomalous delays in the response of these models to a pulse of IP(3).

Mixed-mode oscillations and slow manifolds in the self-coupled FitzHugh-Nagumo system.

We investigate the organization of mixed-mode oscillations in the self-coupled FitzHugh-Nagumo system. These types of oscillations can be explained as a combination of relaxation oscillations and small-amplitude oscillations controlled by canard solutions that are associated with a folded singularity on a critical manifold. The self-coupled FitzHugh-Nagumo system has a cubic critical manifold for a range of parameters, and an associated folded singularity of node-type. Hence, there exist corresponding attracting and repelling slow manifolds that intersect in canard solutions. We present a general technique for the computation of two-dimensional slow manifolds (smooth surfaces). It is based on a boundary value problem approach where the manifolds are computed as one-parameter families of orbit segments. Visualization of the computed surfaces gives unprecedented insight into the geometry of the system. In particular, our techniques allow us to find and visualize canard solutions as the intersection curves of the attracting and repelling slow manifolds.

Resetting behavior in a model of bursting in secretory pituitary cells: distinguishing plateaus from pseudo-plateaus.

We study a recently discovered class of models for plateau bursting, inspired by models for endocrine pituitary cells. In contrast to classical models for fold-homoclinic (square-wave) bursting, the spikes of the active phase are not supported by limit cycles of the frozen fast subsystem, but are transient oscillations generated by unstable limit cycles emanating from a subcritical Hopf bifurcation around a stable steady state. Experimental time courses are suggestive of such fold-subHopf models because the spikes tend to be small and variable in amplitude; we call this pseudo-plateau bursting. We show here that distinct properties of the response to attempted resets from the silent phase to the active phase provide a clearer, qualitative criterion for choosing between the two classes of models. The fold-homoclinic class is characterized by induced active phases that increase towards the duration of the unperturbed active phase as resets are delivered later in the silent phase. For the fold-subHopf class of pseudo-plateau bursting, resetting is difficult and succeeds only in limited windows of the silent phase but, paradoxically, can dramatically exceed the native active phase duration.

Locus of boundary crisis: expect infinitely many gaps.

Boundary crisis is a mechanism for destroying a chaotic attractor when one parameter is varied. In a two-parameter setting the locus of the boundary crisis is associated with curves of homoclinic or heteroclinic bifurcations of periodic saddle points. It is known that this locus has nondifferentiable points. We show here that the locus of boundary crisis is far more complicated than previously reported. It actually contains infinitely many gaps, corresponding to regions (of positive measure) where attractors exist.