Parameter domains for instability of uniform states in systems with many delays

 We present a computational method for determining regions in parameter space corresponding to linear instability of a spatially uniform steady state solution of any system of two coupled reaction-diffusion equations containing up to four delay terms. At each point in parameter space the required stability properties of the linearised system are found using mainly the Principle of the Argument. The method is first developed for perturbations of a particular wavenumber, and then generalised to allow arbitrary perturbations. Each delay term in the system may be of either a fixed or a distributed type, and spatio-temporal delays are also allowed. Received 19 September 1995; received in revised form 4 September 1996     

On the Population Dynamics of the Malaria Vector

A deterministic differential equation model for the population dynamics of the human malaria vector is derived and studied. Conditions for the existence and stability of a non-zero steady state vector population density are derived. These reveal that a threshold parameter, the vectorial basic reproduction number, exist and the vector can established itself in the community if and only if this parameter exceeds unity. When a non-zero steady state population density exists, it can be stable but it can also be driven to instability via a Hopf Bifurcation to periodic solutions, as a parameter is varied in parameter space. By considering a special case, an asymptotic perturbation analysis is used to derive the amplitude of the oscillating solutions for the full non-linear system. The present modelling exercise and results show that it is possible to study the population dynamics of disease vectors, and hence oscillatory behaviour as it is often observed in most indirectly transmitted infectious diseases of humans, without recourse to external seasonal forcing.     

Analysis of unstable behavior in a mathematical model for erythropoiesis

We consider an age-structured model that describes the regulation of erythropoiesis through the negative feedback loop between erythropoietin and hemoglobin. This model is reduced to a system of two ordinary differential equations with two constant delays for which we show existence of a unique steady state. We determine all instances at which this steady state loses stability via a Hopf bifurcation through a theoretical bifurcation analysis establishing analytical expressions for the scenarios in which they arise. We show examples of supercritical Hopf bifurcations for parameter values estimated according to physiological values for humans found in the literature and present numerical simulations in agreement with the theoretical analysis. We provide a strategy for parameter estimation to match empirical measurements and predict dynamics in experimental settings, and compare existing data on hemoglobin oscillation in rabbits with predictions of our model.     

On the relation of the microtubule length and dynamics: boundary effect and properties of an extended radial array

In interphase cells, microtubules (MT) are long and form extended radial array. The length of individual MTs in living cells exhibits substantial stochastic fluctuations while the average length distribution in a cell remains nearly constant. We present a quantitative model that describes relation of the MT length and dynamics in the steady state in the cell using the minimal set of parameters (cell radius, tubulin concentration, critical concentration for plus end elongation, and the number of nucleation sites). The MT array is approximated as a radial system, where MT minus ends are associated with the nucleation sites on the centrosome, while plus ends grow and shorten. Dynamic instability of MT plus ends is approximated as a random walk process with boundary conditions and the behavior of MT array is quantified using diffusion and drift coefficients (Vorobjev et al., 1997, 1999). We show that establishment of the extended steady-state array could be accomplished solely by the limitation of the MT growth by the cell margin. We determined for the cell radius, tubulin concentration, critical concentration for plus end elongation, and number of nucleation sites the reference point in the parameter space where plus ends of individual MT on average neither elongate nor shorten. In this case average length of MT is equal to the half of cell radius. When any parameter is shifted from its reference value MTs become longer or shorter and consequently acquire positive or negative drift of their ends. In the vicinity of reference point, change in any parameter has major effect on the MT length and rather small effect on the drift. When mean length of the MTs is close to the cell radius the drift of the free plus ends becomes substantial, resulting in processive growth of individual MTs in the internal cytoplasm accompanied by apparent stabilization of the plus ends at the cell margin. Under these conditions small changes in parameters have significant impact on the magnitude of drift. Experimental analysis of the MT plus ends dynamics in different cultured cells shows that in most cases plus ends display positive drift, which, in the framework of the presented model, is in agreement with the simultaneous presence of long MTs.     

A delay reaction-diffusion model of the dynamics of botulinum in fish

A model of interaction between fish and a bacterium (Clostridium botulinum) responsible for avian botulism is introduced, considering diffusion of both fish and bacterium in water. The fish population moves randomly in water. Death fish disintegrate in water, at different locations, causing bacteria to diffuse through water and infect other fish. Existence of uniform steady states is investigated and the linearized stability of the positive uniform steady state is analyzed. A Hopf bifurcation is proved to occur from the uniform steady state when the bifurcation parameter, here the time delay, passes through a critical value and diffusion coefficients satisfy some conditions, that induces time oscillations of the populations. Comments on diffusion-driven instability are provided, and numerical simulations are carried out to illustrate the results.     

Microtubule length and dynamics: Boundary effect and properties of extended radial array

In interphase cells, microtubules (MT) form an extended radial array. The length of individual MTs in living cells exhibits substantial stochastic fluctuations, while the average length distribution in a cell remains nearly constant. We present a quantitative model that describes the relation of the MT length and dynamics in the steady state in the cell using the minimal set of parameters (cell radius, tubulin concentration, critical concentration for plus-end elongation and the number of nucleation sites). The MT array is approximated as a radial system, where minus-ends of MTs are associated with nucleation sites on the centrosome, while plus ends grow and shorten. Dynamic instability of MT plus ends is approximated as a random walk process with boundary conditions; the behavior of an MT array is quantified using diffusion and drift coefficients (Vorobjev et al., 1997; Vorobjev et al., 1999). We show that the establishment of the extended steady-state array could be accomplished solely by the limitation of MT growth by the cell margin. For the cell radius, tubulin concentration, critical concentration for plus-end elongation, and the number of nucleation sites we determined the reference point in the parameter space where plus ends of individual MTs, on average, neither elongate nor shorten. In this case, the average MT length is equal to the half of the cell radius. When any parameter is shifted from its reference value, MTs become longer or shorter and, consequently, acquire a positive or negative drift of their plus ends. In the vicinity of the reference point, a change in any parameter has a major effect on the MT length and a rather small effect on the drift. When the average MT length is close to the cell radius, the drift of free plus ends becomes substantial, resulting in processive growth of individual MTs in the internal cytoplasm, accompanied by the apparent stabilization of plus ends at the cell margin. Under these conditions small changes in parameters have a significant impact on the magnitude of the drift. Experimental analysis of MT plus-end dynamics in different cultured cells shows that, in most cases, plus ends display positive drift, which, in the framework of the presented model, is in agreement with the simultaneous presence of long MTs.     

Indirect Genetic Effects and the Dynamics of Social Interactions

BackgroundIndirect genetic effects (IGEs) occur when genes expressed in one individual alter the expression of traits in social partners. Previous studies focused on the evolutionary consequences and evolutionary dynamics of IGEs, using equilibrium solutions to predict phenotypes in subsequent generations. However, whether or not such steady states may be reached may depend on the dynamics of interactions themselves.ResultsIn our study, we focus on the dynamics of social interactions and indirect genetic effects and investigate how they modify phenotypes over time. Unlike previous IGE studies, we do not analyse evolutionary dynamics; rather we consider within-individual phenotypic changes, also referred to as phenotypic plasticity. We analyse iterative interactions, when individuals interact in a series of discontinuous events, and investigate the stability of steady state solutions and the dependence on model parameters, such as population size, strength, and the nature of interactions. We show that for interactions where a feedback loop occurs, the possible parameter space of interaction strength is fairly limited, affecting the evolutionary consequences of IGEs. We discuss the implications of our results for current IGE model predictions and their limitations.     

Existence of excitation waves for a collection of cardiomyocytes electrically coupled to fibroblasts

We consider mathematical models of a collection of cardiomyocytes (myocardial tissue) coupled to a varying number of fibroblasts. Our aim is to understand how conductivity (δ) and fibroblast density (η) affect the stability of the collection. We provide mathematical and computational arguments indicating that there is a region of instability in the η-δ space. Mathematical arguments, based on a simplified model of the coupled myocyte-fibroblast system, show that for certain parameter choices, a stationary solution cannot exist. Numerical experiments (1D, 2D) are based on a recently developed model of electro-chemical coupling between a human atrial myocyte and a number of associated atrial fibroblasts. The numerical experiments demonstrate that there is a region of instability of the form observed in the simplified model analysis.     

Equilibria and stability of a class of positive feedback loops

Positive feedback loops are common regulatory elements in metabolic and protein signalling pathways. The length of such feedback loops determines stability and sensitivity to network perturbations. Here we provide a mathematical analysis of arbitrary length positive feedback loops with protein production and degradation. These loops serve as an abstraction of typical regulation patterns in protein signalling pathways. We first perform a steady state analysis and, independently of the chain length, identify exactly two steady states that represent either biological activity or inactivity. We thereby provide two formulas for the steady state protein concentrations as a function of feedback length, strength of feedback, as well as protein production and degradation rates. Using a control theory approach, analysing the frequency response of the linearisation of the system and exploiting the Small Gain Theorem, we provide conditions for local stability for both steady states. Our results demonstrate that, under some parameter relationships, once a biological meaningful on steady state arises, it is stable, while the off steady state, where all proteins are inactive, becomes unstable. We apply our results to a three-tier feedback of caspase activation in apoptosis and demonstrate how an intermediary protein in such a loop may be used as a signal amplifier within the cascade. Our results provide a rigorous mathematical analysis of positive feedback chains of arbitrary length, thereby relating pathway structure and stability.     

Dissection of a model for neuronal parabolic bursting

We have obtained new insight into the mechanisms for bursting in a class of theoretical models. We study Plant's model [24] for Aplysia R-15 to illustrate our view of these so-called parabolic bursters, which are characterized by low spike frequency at the beginning and end of a burst. By identifying and analyzing the fast and slow processes we show how they interact mutually to generate spike activity and the slow wave which underlies the burst pattern. Our treatment is essentially the first step of a singular perturbation approach presented from a geometrical viewpoint and carried out numerically with AUTO [12]. We determine the solution sets (steady state and oscillatory) of the fast subsystem with the slow variables treated as parameters. These solutions form the slow manifold over which the slow dynamics then define a burst trajectory. During the silent phase of a burst, the solution trajectory lies approximately on the steady state branch of the slow manifold and during the active phase of spiking, the trajectory sweeps through the oscillation branch. The parabolic nature of bursting arises from the (degenerate) homoclinic transition between the oscillatory branch and the steady state branch. We show that, for some parameter values, the trajectory remains strictly on the steady state branch (to produce a resting steady state or a pure slow wave without spike activity) or strictly in the oscillatory branch (continuous spike activity without silent phases). Plant's model has two slow variables: a calcium conductance and the intracellular free calcium concentration, which activates a potassium conductance. We also show how bursting arises from an alternative mechanism in which calcium inactivates the calcium current and the potassium conductance is insensitive to calcium. These and other biophysical interpretations are discussed.     

Mathematical analysis of an age-structured population model with space-limited recruitment

In this paper, we investigate structured population model of marine invertebrate whose life stage is composed of sessile adults and pelagic larvae, such as barnacles contained in a local habitat. First we formulate the basic model as an Cauchy problem on a Banach space to discuss the existence and uniqueness of non-negative solution. Next we define the basic reproduction number R0 to formulate the invasion condition under which the larvae can successfully settle down in the completely vacant habitat. Subsequently we examine existence and stability of steady states. We show that the trivial steady state is globally asymptotically stable if R0 < or = 1, whereas it is unstable if R0 > 1. Furthermore, we show that a positive (non-trivial) steady state uniquely exists if R0 > 1 and it is locally asymptotically stable as far as absolute value of R0 - 1 is small enough.     

Analysis of Discrete Bioregulatory Networks Using Symbolic Steady States

A discrete model of a biological regulatory network can be represented by a discrete function that contains all available information on interactions between network components and the rules governing the evolution of the network in a finite state space. Since the state space size grows exponentially with the number of network components, analysis of large networks is a complex problem. In this paper, we introduce the notion of symbolic steady state that allows us to identify subnetworks that govern the dynamics of the original network in some region of state space. We state rules to explicitly construct attractors of the system from subnetwork attractors. Using the results, we formulate sufficient conditions for the existence of multiple attractors resp. a cyclic attractor based on the existence of positive resp. negative feedback circuits in the graph representing the structure of the system. In addition, we discuss approaches to finding symbolic steady states. We focus both on dynamics derived via synchronous as well as asynchronous update rules. Lastly, we illustrate the results by analyzing a model of T helper cell differentiation.     

Periodic solutions: a robust numerical method for an S-I-R model of epidemics

 We describe and analyze a numerical method for an S-I-R type epidemic model. We prove that it is unconditionally convergent and that solutions it produces share many qualitative and quantitative properties of the solution of the differential problem being approximated. Finally, we establish explicit sufficient conditions for the unique endemic steady state of the system to be unstable and we use our numerical algorithm to approximate the solution in such cases and discover that it can be periodic, just as suggested by the instability of the endemic steady state. Received: 1 September 1995 / Revised version: 30 April 1997     

Coupled Reversible and Irreversible Bistable Switches Underlying TGFβ-induced Epithelial to Mesenchymal Transition

Epithelial to mesenchymal transition (EMT) plays an important role in embryonic development, tissue regeneration, and cancer metastasis. Whereas several feedback loops have been shown to regulate EMT, it remains elusive how they coordinately modulate EMT response to TGF-β treatment. We construct a mathematical model for the core regulatory network controlling TGF-β-induced EMT. Through deterministic analyses and stochastic simulations, we show that EMT is a sequential two-step program in which an epithelial cell first is converted to partial EMT then to the mesenchymal state, depending on the strength and duration of TGF-β stimulation. Mechanistically the system is governed by coupled reversible and irreversible bistable switches. The SNAIL1/miR-34 double-negative feedback loop is responsible for the reversible switch and regulates the initiation of EMT, whereas the ZEB/miR-200 feedback loop is accountable for the irreversible switch and controls the establishment of the mesenchymal state. Furthermore, an autocrine TGF-β/miR-200 feedback loop makes the second switch irreversible, modulating the maintenance of EMT. Such coupled bistable switches are robust to parameter variation and molecular noise. We provide a mechanistic explanation on multiple experimental observations. The model makes several explicit predictions on hysteretic dynamic behaviors, system response to pulsed stimulation, and various perturbations, which can be straightforwardly tested.     

Turing Patterns from Dynamics of Early HIV Infection

We have developed a mathematical model for in-host virus dynamics that includes spatial chemotaxis and diffusion across a two-dimensional surface representing the vaginal or rectal epithelium at primary HIV infection. A linear stability analysis of the steady state solutions identified conditions for Turing instability pattern formation. We have solved the model equations numerically using parameter values obtained from previous experimental results for HIV infections. Simulations of the model for this surface show hot spots of infection. Understanding this localization is an important step in the ability to correctly model early HIV infection. These spatial variations also have implications for the development and effectiveness of microbicides against HIV.     

Coupled Reversible and Irreversible Bistable Switches Underlying TGFβ-induced Epithelial to Mesenchymal Transition

Epithelial to mesenchymal transition (EMT) plays an important role in embryonic development, tissue regeneration, and cancer metastasis. Whereas several feedback loops have been shown to regulate EMT, it remains elusive how they coordinately modulate EMT response to TGF-β treatment. We construct a mathematical model for the core regulatory network controlling TGF-β-induced EMT. Through deterministic analyses and stochastic simulations, we show that EMT is a sequential two-step program in which an epithelial cell first is converted to partial EMT then to the mesenchymal state, depending on the strength and duration of TGF-β stimulation. Mechanistically the system is governed by coupled reversible and irreversible bistable switches. The SNAIL1/miR-34 double-negative feedback loop is responsible for the reversible switch and regulates the initiation of EMT, whereas the ZEB/miR-200 feedback loop is accountable for the irreversible switch and controls the establishment of the mesenchymal state. Furthermore, an autocrine TGF-β/miR-200 feedback loop makes the second switch irreversible, modulating the maintenance of EMT. Such coupled bistable switches are robust to parameter variation and molecular noise. We provide a mechanistic explanation on multiple experimental observations. The model makes several explicit predictions on hysteretic dynamic behaviors, system response to pulsed stimulation, and various perturbations, which can be straightforwardly tested.     

Symmetry-breaking in Chiral Polymerisation

We propose a model for chiral polymerisation and investigate its symmetric and asymmetric solutions. The model has a source species which decays into left- and right-handed types of monomer, each of which can polymerise to form homochiral chains; these chains are susceptible to ‘poisoning’ by the opposite-handed monomer. Homochiral polymers are assumed to influence the proportion of each type of monomer formed from the precursor. We show that for certain parameter values a positive feedback mechanism makes the symmetric steady-state solution unstable. The kinetics of polymer formation are then analysed in the case where the system starts from zero concentrations of monomers and chains. We show that following a long induction time, extremely large concentrations of polymers are formed for a short time, during this time an asymmetry introduced into the system by a random external perturbation may be massively amplified. The system then approaches one of the steady-state solutions described above.     

Rich dynamics of a hepatitis B viral infection model with logistic hepatocyte growth

Chronic hepatitis B virus (HBV) infection is a major cause of human suffering, and a number of mathematical models have examined within-host dynamics of the disease. Most previous HBV infection models have assumed that: (a) hepatocytes regenerate at a constant rate from a source outside the liver; and/or (b) the infection takes place via a mass action process. Assumption (a) contradicts experimental data showing that healthy hepatocytes proliferate at a rate that depends on current liver size relative to some equilibrium mass, while assumption (b) produces a problematic basic reproduction number. Here we replace the constant infusion of healthy hepatocytes with a logistic growth term and the mass action infection term by a standard incidence function; these modifications enrich the dynamics of a well-studied model of HBV pathogenesis. In particular, in addition to disease free and endemic steady states, the system also allows a stable periodic orbit and a steady state at the origin. Since the system is not differentiable at the origin, we use a ratio-dependent transformation to show that there is a region in parameter space where the origin is globally stable. When the basic reproduction number, R ₀, is less than 1, the disease free steady state is stable. When R ₀ > 1 the system can either converge to the chronic steady state, experience sustained oscillations, or approach the origin. We characterize parameter regions for all three situations, identify a Hopf and a homoclinic bifurcation point, and show how they depend on the basic reproduction number and the intrinsic growth rate of hepatocytes.     

On the dynamics of a simple biochemical control circuit

The quantitative dynamics of a biochemical control circuit that regulates enzyme or protein synthesis by end-product feedback is analyzed. We first study a simplified repressible system, which is known to exhibit either a steady state or an oscillatory solution. By showing the analogy of thisn-dimensional system with a time-delay equation for a single variable the mechanism of the self-sustained oscillations becomes transparent. In a more sophisticated system we will find as well either steady state or oscillatory solutions. We determine the role of the parameters with respect to stability and frequency. The most general case will be treated by means of the concept of Lyapunov exponents.