Oscillations in a model of repression with external control

A mathematical model for control by repression by an extracellular substance is developed, including diffusion and time delays. The model examines how active transport of a nutrient can produce either oscillatory or stable responses depending on a variety of parameters, such as diffusivity, cell size, or nutrient concentration. The system of equations for the mathematical model is reduced to a system of delay differential equations and linear Volterra equations. After linearizing these equations and forming the limiting Volterra equations, the resulting linear system no longer has any spatial dependence. Local stability analysis of the radially symmetric model shows that the system of equations can undergo Hopf bifurcations for certain parameter values, while other ranges of the parameters guarantee asymptotic stability. One numerical study shows that the model can exhibit intracellular biochemical oscillations with increasing extracellular concentrations of the nutrient, which suggests a possible trigger mechanism for morphogenesis.The work of this author was supported in part by NSF grants DMS-8603787 and DMS-8807360The work of this author was supported under the REU program of NSF by grant DMS-8807360The work of this author was supported under the REU program of NSF by grant DMS-8807360     

Steady states and dynamics of urokinase-mediated plasmin activation in silico and in vitro

Plasmin (PLS) and urokinase-type plasminogen activator (UPA) are ubiquitous proteases that regulate the extracellular environment. Although they are secreted in inactive forms, they can activate each other through proteolytic cleavage. This mutual interplay creates the potential for complex dynamics, which we investigated using mathematical modeling and in vitro experiments. We constructed ordinary differential equations to model the conversion of precursor plasminogen into active PLS, and precursor urokinase (scUPA) into active urokinase (tcUPA). Although neither PLS nor UPA exhibits allosteric cooperativity, modeling showed that cooperativity occurred at the system level because of substrate competition. Computational simulations and bifurcation analysis predicted that the system would be bistable over a range of parameters for cooperativity and positive feedback. Cell-free experiments with recombinant proteins tested key predictions of the model. PLS activation in response to scUPA stimulus was found to be cooperative in vitro. Finally, bistability was demonstrated in vitro by the presence of two significantly different steady-state levels of PLS activation for the same levels of stimulus. We conclude that ultrasensitive, bistable activation of UPA-PLS is possible in the presence of substrate competition. An ultrasensitive threshold for activation of PLS and UPA would have ramifications for normal and disease processes, including angiogenesis, metastasis, wound healing, and fibrosis.     

The effects of time delays in a phosphorylation-dephosphorylation pathway

Complex signaling cascades involve many interlocked positive and negative feedback loops which have inherent delays. Modeling these complex cascades often requires a large number of variables and parameters. Delay differential equation models have been helpful in describing inherent time lags and also in reducing the number of governing equations. However the consequences of model reduction via delay differential equations have not been fully explored. In this paper we systematically examine the effect of delays in a complex network of phosphorylation-dephosphorylation cycles (described by Gonze and Goldbeter, J. Theor. Biol., 210, (2001) 167-186), which commonly occur in many biochemical pathways. By introducing delays in the positive and negative regulatory interactions, we show that a delay differential model can indeed reduce the number of cycles actually required to describe the phosphorylation-dephosphorylation pathway. In addition, we find some of the unique properties of the network and a quantitative measure of the minimum number of delay variables required to model the network. These results can be extended for modeling complex signalling cascades.     

Numerical bifurcation analysis of delay differential equations arising from physiological modeling

This paper has a dual purpose. First, we describe numerical methods for continuation and bifurcation analysis of steady state solutions and periodic solutions of systems of delay differential equations with an arbitrary number of fixed, discrete delays. Second, we demonstrate how these methods can be used to obtain insight into complex biological regulatory systems in which interactions occur with time delays: for this, we consider a system of two equations for the plasma glucose and insulin concentrations in a diabetic patient subject to a system of external assistance. The model has two delays: the technological delay of the external system, and the physiological delay of the patient's liver. We compute stability of the steady state solution as a function of two parameters, compare with analytical results and compute several branches of periodic solutions and their stability. These numerical results allow to infer two categories of diabetic patients for which the external system has different efficiency.     

Antagonism and bistability in protein interaction networks

A protein interaction network (PIN) is a set of proteins that modulate one another's activities by regulated synthesis and degradation, by reversible binding to form complexes, and by catalytic reactions (e.g., phosphorylation and dephosphorylation). Most PINs are so complex that their dynamical characteristics cannot be deduced accurately by intuitive reasoning alone. To predict the properties of such networks, many research groups have turned to mathematical models (differential equations based on standard biochemical rate laws, e.g., mass-action, Michaelis-Menten, Hill). When using Michaelis-Menten rate expressions to model PINs, care must be exercised to avoid making inconsistent assumptions about enzyme-substrate complexes. We show that an appealingly simple model of a PIN that functions as a bistable switch is compromised by neglecting enzyme-substrate intermediates. When the neglected intermediates are put back into the model, bistability of the switch is lost. The theory of chemical reaction networks predicts that bistability can be recovered by adding specific reaction channels to the molecular mechanism. We explore two very different routes to recover bistability. In both cases, we show how to convert the original 'phenomenological' model into a consistent set of mass-action rate laws that retains the desired bistability properties. Once an equivalent model is formulated in terms of elementary chemical reactions, it can be simulated accurately either by deterministic differential equations or by Gillespie's stochastic simulation algorithm.     

A Detailed Modular Analysis of Heat-Shock Protein Dynamics under Acute and Chronic Stress and Its Implication in Anxiety Disorders

Physiological and psychological stresses cause anxiety disorders such as depression and post-traumatic stress disorder (PTSD) and induce drastic changes at a molecular level in the brain. To counteract this stress, the heat-shock protein (HSP) network plays a vital role in restoring the homeostasis of the system. To study the stress-induced dynamics of heat-shock network, we analyzed three modules of the HSP90 network-namely trimerization reactions, phosphorylation-dephosphorylation reactions, and the conversion of HSP90 from an open to a closed conformation-and constructed a corresponding nonlinear differential equation model based on mass action kinetics laws. The kinetic parameters of the model were obtained through global optimization, and sensitivity analyses revealed that the most sensitive parameters are the kinase and phosphatase that drive the phosphorylation-dephosphorylation reactions. Bifurcation analysis carried out with the estimated kinetic parameters of the model with stress as bifurcation parameter revealed the occurrence of "mushroom", a type of complex dynamics in which S-shaped and Z-shaped hysteretic bistable forms are present together. We mapped the molecular events responsible for generating the mushroom dynamics under stress and interpreted the occurrence of the S-shaped hysteresis to a normal level of stress, and the Z-shaped hysteresis to the HSP90 variations under acute and chronic stress in the fear conditioned system, and further, we hypothesized that this can be extended to stress-related disorders such as depression and PTSD in humans. Finally, we studied the effect of parameter variations on the mushroom dynamics to get insight about the role of phosphorylation-dephosphorylation parameters in HSP90 network in bringing about complex dynamics such as isolas, where the stable steady states in a bistable system are isolated and separated from each other and not connected by an unstable steady state.     

BioBayes: a software package for Bayesian inference in systems biology

MOTIVATION: There are several levels of uncertainty involved in the mathematical modelling of biochemical systems. There often may be a degree of uncertainty about the values of kinetic parameters, about the general structure of the model and about the behaviour of biochemical species which cannot be observed directly. The methods of Bayesian inference provide a consistent framework for modelling and predicting in these uncertain conditions. We present a software package for applying the Bayesian inferential methodology to problems in systems biology. RESULTS: Described herein is a software package, BioBayes, which provides a framework for Bayesian parameter estimation and evidential model ranking over models of biochemical systems defined using ordinary differential equations. The package is extensible allowing additional modules to be included by developers. There are no other such packages available which provide this functionality.     

In silico evolution of functional modules in biochemical networks

Understanding the large reaction networks found in biological systems is a daunting task. One approach is to divide a network into more manageable smaller modules, thus simplifying the problem. This is a common strategy used in engineering. However, the process of identifying biological modules is still in its infancy and very little is understood about the range and capabilities of motif structures found in biological modules. In order to delineate these modules, a library of functional motifs has been generated via in silico evolution techniques. On the basis of their functional forms, networks were evolved from four broad areas: oscillators, bistable switches, homeostatic systems and frequency filters. Some of these motifs were constructed from simple mass action kinetics, others were based on Michaelis-Menten kinetics as found in protein/protein networks and the remainder were based on Hill equations as found in gene/protein interaction networks. The purpose of the study is to explore the capabilities of different network architectures and the rich variety of functional forms that can be generated. Ultimately, the library may be used to delineate functional motifs in real biological networks.     

A competition model with dynamically allocated inhibitor production

The chemostat is a basic model for competition in an open system and a model for the laboratory bio-reactor (CSTR). Inhibitors in open systems are studied with a view of detoxification in natural systems and of control in bio-reactors. This study allows the amount of resource devoted to inhibitor production to depend on the state of the system. The feasibility of one dependence is provided by quorum sensing. In contrast to the constant allocation case, a much wider set of outcomes is possible including interior, stable rest points and stable limit cycles. These outcomes are important contrasts to competitive exclusion or bistable attractors that are often the outcomes for competitive systems. The model consists of four non-linear ordinary differential equations and computer software is used for most of the stability calculations.     

Incorporation of variability into the modeling of viral delays in HIV infection dynamics

We consider classes of functional differential equation models which arise in attempts to describe temporal delays in HIV pathogenesis. In particular, we develop methods for incorporating arbitrary variability (i.e., general probability distributions) for these delays into systems that cannot readily be reduced to a finite number of coupled ordinary differential equations (as is done in the method of stages). We discuss modeling from first principles, introduce several classes of non-linear models (including discrete and distributed delays) and present a discussion of theoretical and computational approaches. We then use the resulting methodology to carry out simulations and perform parameter estimation calculations, fitting the models to a set of experimental data. Results obtained confirm the statistical significance of the presence of delays and the importance of including delays in validating mathematical models with experimental data. We also show that the models are quite sensitive to the mean of the distribution which describes the delay in viral production, whereas the variance of this distribution has relatively little impact.     

Ultrasensitivity in independent multisite systems

Multisite modifications are widely recognized as an essential feature of many switch-like responses in signal transduction. It is usually assumed that the modification of one site directly or indirectly increases the rate of modification of neighboring sites. In this paper we provide a new set of assumptions for a multisite system to become highly ultrasensitive even in the absence of cooperativity or allostery. We assume that the individual sites are modified independently of each other, and that protein activity is an ultrasensitive function of the fraction of modified sites. These assumptions are particularly useful in the context of multisite systems with a large (8+) number of sites. We estimate the apparent Hill coefficient of the dose responses in the sequential and nonsequential cases, highlight their different qualitative properties, and discuss a formula to approximate dose responses in the nonsequential case. As an example we describe a model of bacterial chemotaxis that features robust ultrasensitivity and perfect adaptation over a wide range of ligand concentrations, based on non-allosteric multisite behavior at the level of receptors and flagella. We also include a model of the inactivation of the yeast pheromone protein Ste5 by cell cycle proteins.     

Bistability in the JNK cascade

BACKGROUND: Important signaling properties, like adaptation, oscillations, and bistability, can emerge at the level of relatively simple systems of signaling proteins. Here, we have examined the quantitative properties of one well-studied signaling system, the JNK cascade. We experimentally assessed the response of JNK to a physiological stimulus (progesterone) and a pathological stress (hyperosmolar sorbitol) in Xenopus laevis oocytes, a cell type that is well-suited to the quantitative analysis of cell signaling. Our aim was to determine whether JNK responses are graded (Michaelian) in character; ultrasensitive in character, resembling the responses of cooperative enzymes; or bistable and all-or-none in character. RESULTS: The responses of JNK to both progesterone and sorbitol were found to be essentially all-or-none. Individual oocytes had either very high or very low JNK activities, with few oocytes possessing intermediate levels of JNK activity. Moreover, JNK activation was autocatalytic, indicating that the JNK cascade is either embedded in or downstream of a positive feedback loop. JNK also exhibited hysteresis, a form of biochemical memory, in its response to sorbitol. These findings indicate that the JNK cascade is part of a bistable signaling system in oocytes. CONCLUSIONS: In Xenopus oocytes, JNK responds to physiological and pathological stimuli in an all-or-none manner. The JNK response shows all the hallmarks of a bistable response, including strong positive feedback and hysteresis. Bistability is a recurring theme in the biochemistry of oocyte maturation and early embryogenesis; the Mos/MEK/p42 MAPK cascade also exhibits bistable responses, and the Cdc2/cyclin B system is hypothesized to be bistable as well. However, the mechanisms underpinning the positive feedback and bistability in the three cases are different, suggesting that evolution has repeatedly converged upon bistability as a way of producing digital responses.     

Tuned solutions in dynamic neural fields as building blocks for extended EEG models

The most prominent functional property of cortical neurons in sensory areas are their tuned receptive fields which provide specific responses of the neurons to external stimuli. Tuned neural firing indeed reflects the most basic and best worked out level of cognitive representations. Tuning properties can be dynamic on a short time-scale of fractions of a second. Such dynamic effects have been modeled by localised solutions (also called “bumps” or “peaks”) in dynamic neural fields. In the present work we develop an approximation method to reduce the dynamics of localised activation peaks in systems of n coupled nonlinear d-dimensional neural fields with transmission delays to a small set of delay differential equations for the peak amplitudes and widths only. The method considerably simplifies the analysis of peaked solutions as demonstrated for a two-dimensional example model of neural feature selectivity in the brain. The reduced equations describe the effective interaction between pools of local neurons of several (n) classes that participate in shaping the dynamic receptive field responses. To lowest order they resemble neural mass models as they often form the base of EEG-models. Thereby they provide a link between functional small-scale receptive field models and more coarse-grained EEG-models. More specifically, they connect the dynamics in feature-selective cortical microcircuits to the more abstract local elements used in coarse-grained models. However, beside amplitudes the reduced equations also reflect the sharpness of tuning of the activity in a d-dimensional feature space in response to localised stimuli.     

Ultrasensitivity in Phosphorylation-Dephosphorylation Cycles with Little Substrate

Cellular decision-making is driven by dynamic behaviours, such as the preparations for sunrise enabled by circadian rhythms and the choice of cell fates enabled by positive feedback. Such behaviours are often built upon ultrasensitive responses where a linear change in input generates a sigmoidal change in output. Phosphorylation-dephosphorylation cycles are one means to generate ultrasensitivity. Using bioinformatics, we show that in vivo levels of kinases and phosphatases frequently exceed the levels of their corresponding substrates in budding yeast. This result is in contrast to the conditions often required by zero-order ultrasensitivity, perhaps the most well known means for how such cycles become ultrasensitive. We therefore introduce a mechanism to generate ultrasensitivity when numbers of enzymes are higher than numbers of substrates. Our model combines distributive and non-distributive actions of the enzymes with two-stage binding and concerted allosteric transitions of the substrate. We use analytical and numerical methods to calculate the Hill number of the response. For a substrate with phosphosites, we find an upper bound of the Hill number of , and so even systems with a single phosphosite can be ultrasensitive. Two-stage binding, where an enzyme must first bind to a binding site on the substrate before it can access the substrate''s phosphosites, allows the enzymes to sequester the substrate. Such sequestration combined with competition for each phosphosite provides an intuitive explanation for the sigmoidal shifts in levels of phosphorylated substrate. Additionally, we find cases for which the response is not monotonic, but shows instead a peak at intermediate levels of input. Given its generality, we expect the mechanism described by our model to often underlay decision-making circuits in eukaryotic cells.

Authors Summary

Dose-response curves are said to be ultrasensitive when they are sigmoidal rather than hyperbolic and often underlay cellular decision-making circuits. Zero-order ultrasensitivity is a well-known mechanism to generate sigmoidal curves in phosphorylation cycles, but one of its assumptions often implies that the substrate is more abundant than the modifying enzymes. We show that this assumption is unlikely to always hold in vivo, and we present a general model that generates ultrasensitivity when the enzymes are in excess of their substrate. The model combines conformational allosteric transitions of the substrate with two-stage binding of the enzymes: the enzymes bind first to a docking site on the substrate and then to the substrate''s phosphosites. Ultrasensitivity is generated because the kinase can bind to the fully phosphorylated form of the substrate (at its docking site) and sequester the substrate away from the phosphatase and, similarly, the phosphatase can bind to the fully dephosphorylated form of the substrate and sequester the substrate away from the kinase. The number of kinase-phosphatase competitions for the substrate determines the degree of ultrasensitivity. Finally, we show that this model can generate non-monotonic responses that peak at intermediate levels of input.     

Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis

We propose a new mathematical model of erythropoiesis that takes a positive feedback of erythrocytes on progenitor apoptosis into account, and incorporates a negative feedback of erythrocytes on progenitor self-renewal. The resulting model is a system of age-structured equations that reduces to a system of delay differential equations where the delays account for progenitor compartment duration and cell cycle length. We compare this model with experimental data on an induced-anemia in mice that exhibit damped oscillations of the hematocrit before it returns to equilibrium. When we assume no self-renewal of progenitors, we obtain an inaccurate fitting of the model with experimental data. Adding self-renewal in the progenitor compartment gives better approximations, with the main features of experimental data correctly fitted. Our results indicate the importance of progenitor self-renewal in the modelling of erythropoiesis. Moreover, the model makes testable predictions on the lifespan of erythrocytes confronted to a severe anemia, and on the progenitors behavior.     

Mathematical and Computational Modeling for Tumor Virotherapy with Mediated Immunity

We propose a new mathematical modeling framework based on partial differential equations to study tumor virotherapy with mediated immunity. The model incorporates both innate and adaptive immune responses and represents the complex interaction among tumor cells, oncolytic viruses, and immune systems on a domain with a moving boundary. Using carefully designed computational methods, we conduct extensive numerical simulation to the model. The results allow us to examine tumor development under a wide range of settings and provide insight into several important aspects of the virotherapy, including the dependence of the efficacy on a few key parameters and the delay in the adaptive immunity. Our findings also suggest possible ways to improve the virotherapy for tumor treatment.