Understanding bistability in complex enzyme-driven reaction networks

Dynamical Modeling and Qualitative Analysis of a Delayed Model for CD8 T Cells in Response to Viral Antigens

Although the immune effector CD8 T cells play a crucial role in clearance of viruses, the mechanisms underlying the dynamics of how CD8 T cells respond to viral infection remain largely unexplored. Here, we develop a delayed model that incorporates CD8 T cells and infected cells to investigate the functional role of CD8 T cells in persistent virus infection. Bifurcation analysis reveals that the model has four steady states that can finely divide the progressions of viral infection into four states, and endows the model with bistability that has ability to achieve the switch from one state to another. Furthermore, analytical and numerical methods find that the time delay resulting from incubation period of virus can induce a stable low-infection steady state to be oscillatory, coexisting with a stable high-infection steady state in phase space. In particular, a novel mechanism to achieve the switch between two stable steady states, time-delay-based switch, is proposed, where the initial conditions and other parameters of the model remain unchanged. Moreover, our model predicts that, for a certain range of initial antigen load: 1) under a longer incubation period, the lower the initial antigen load, the easier the virus infection will evolve into severe state; while the higher the initial antigen load, the easier it is for the virus infection to be effectively controlled and 2) only when the incubation period is small, the lower the initial antigen load, the easier it is to effectively control the infection progression. Our results are consistent with multiple experimental observations, which may facilitate the understanding of the dynamical and physiological mechanisms of CD8 T cells in response to viral infections.

Exploiting the structure of biochemical pathways to investigate dynamical properties with neural networks for graphs

MOTIVATION

Dynamical properties of biochemical pathways (BPs) help in understanding the functioning of living cells. Their in silico assessment requires simulating a dynamical system with a large number of parameters such as kinetic constants and species concentrations. Such simulations are based on numerical methods that can be time-expensive for large BPs. Moreover, parameters are often unknown and need to be estimated.

RESULTS

We developed a framework for the prediction of dynamical properties of BPs directly from the structure of their graph representation. We represent BPs as Petri nets, which can be automatically generated, for instance, from standard SBML representations. The core of the framework is a neural network for graphs that extracts relevant information directly from the Petri net structure and exploits them to learn the association with the desired dynamical property. We show experimentally that the proposed approach reliably predicts a range of diverse dynamical properties (robustness, monotonicity, and sensitivity) while being faster than numerical methods at prediction time. In synergy with the neural network models, we propose a methodology based on Petri nets arc knock-out that allows the role of each molecule in the occurrence of a certain dynamical property to be better elucidated. The methodology also provides insights useful for interpreting the predictions made by the model. The results support the conjecture often considered in the context of systems biology that the BP structure plays a primary role in the assessment of its dynamical properties.

AVAILABILITY AND IMPLEMENTATION

https://github.com/marcopodda/petri-bio (code), https://zenodo.org/record/7610382 (data).

Phase plane dynamics of ERK phosphorylation

The extracellular signal-regulated kinase (ERK) controls multiple critical processes in the cell and is deregulated in human cancers, congenital abnormalities, immune diseases, and neurodevelopmental syndromes. Catalytic activity of ERK requires dual phosphorylation by an upstream kinase, in a mechanism that can be described by two sequential Michaelis-Menten steps. The estimation of individual reaction rate constants from kinetic data in the full mechanism has proved challenging. Here, we present an analytically tractable approach to parameter estimation that is based on the phase plane representation of ERK activation and yields two combinations of six reaction rate constants in the detailed mechanism. These combinations correspond to the ratio of the specificities of two consecutive phosphorylations and the probability that monophosphorylated substrate does not dissociate from the enzyme before the second phosphorylation. The presented approach offers a language for comparing the effects of mutations that disrupt ERK activation and function in vivo. As an illustration, we use phase plane representation to analyze dual phosphorylation under heterozygous conditions, when two enzyme variants compete for the same substrate.

An integrated modelling approach for targeted degradation: insights on optimization, data requirements and PKPD predictions from semi- or fully-mechanistic models and exact steady state solutions

The value of an integrated mathematical modelling approach for protein degraders which combines the benefits of traditional turnover models and fully mechanistic models is presented. Firstly, we show how exact solutions of the mechanistic models of monovalent and bivalent degraders can provide insight on the role of each system parameter in driving the pharmacological response. We show how on/off binding rates and degradation rates are related to potency and maximal effect of monovalent degraders, and how such relationship can be used to suggest a compound optimization strategy. Even convoluted exact steady state solutions for bivalent degraders provide insight on the type of observations required to ensure the predictive capacity of a mechanistic approach. Specifically for PROTACs, the structure of the exact steady state solution suggests that the total remaining target at steady state, which is easily accessible experimentally, is insufficient to reconstruct the state of the whole system at equilibrium and observations on different species (such as binary/ternary complexes) are necessary. Secondly, global sensitivity analysis of fully mechanistic models for PROTACs suggests that both target and ligase baselines (actually, their ratio) are the major sources of variability in the response of non-cooperative systems, which speaks to the importance of characterizing their distribution in the target patient population. Finally, we propose a pragmatic modelling approach which incorporates the insights generated with fully mechanistic models into simpler turnover models to improve their predictive ability, hence enabling acceleration of drug discovery programs and increased probability of success in the clinic.

Prevalence of multistability and nonstationarity in driven chemical networks

External flows of energy, entropy, and matter can cause sudden transitions in the stability of biological and industrial systems, fundamentally altering their dynamical function. How might we control and design these transitions in chemical reaction networks? Here, we analyze transitions giving rise to complex behavior in random reaction networks subject to external driving forces. In the absence of driving, we characterize the uniqueness of the steady state and identify the percolation of a giant connected component in these networks as the number of reactions increases. When subject to chemical driving (influx and outflux of chemical species), the steady state can undergo bifurcations, leading to multistability or oscillatory dynamics. By quantifying the prevalence of these bifurcations, we show how chemical driving and network sparsity tend to promote the emergence of these complex dynamics and increased rates of entropy production. We show that catalysis also plays an important role in the emergence of complexity, strongly correlating with the prevalence of bifurcations. Our results suggest that coupling a minimal number of chemical signatures with external driving can lead to features present in biochemical processes and abiogenesis.

Trade-offs between number fluctuations and response in nonequilibrium chemical reaction networks

We study the response of chemical reaction networks driven far from equilibrium to logarithmic perturbations of reaction rates. The response of the mean number of a chemical species is observed to be quantitively limited by number fluctuations and the maximum thermodynamic driving force. We prove these trade-offs for linear chemical reaction networks and a class of nonlinear chemical reaction networks with a single chemical species. Numerical results for several model systems support the conclusion that these trade-offs continue to hold for a broad class of chemical reaction networks, though their precise form appears to sensitively depend on the deficiency of the network.

Hematopoietic and Chronic Myeloid Leukemia Stem Cells: Multi-Stability versus Lineage Restriction

There is compelling evidence to support the view that the cell-of-origin for chronic myeloid leukemia is a hematopoietic stem cell. Unlike normal hematopoietic stem cells, the progeny of the leukemia stem cells are predominantly neutrophils during the disease chronic phase and there is a mild anemia. The hallmark oncogene for chronic myeloid leukemia is the BCR-ABLp210 fusion gene. Various studies have excluded a role for BCR-ABLp210 expression in maintaining the population of leukemia stem cells. Studies of BCR-ABLp210 expression in embryonal stem cells that were differentiated into hematopoietic stem cells and of the expression in transgenic mice have revealed that BCR-ABLp210 is able to veer hematopoietic stem and progenitor cells towards a myeloid fate. For the transgenic mice, global changes to the epigenetic landscape were observed. In chronic myeloid leukemia, the ability of the leukemia stem cells to choose from the many fates that are available to normal hematopoietic stem cells appears to be deregulated by BCR-ABLp210 and changes to the epigenome are also important. Even so, we still do not have a precise picture as to why neutrophils are abundantly produced in chronic myeloid leukemia.

A mathematical model captures the role of adenyl cyclase Cyr1 and guanidine exchange factor Ira2 in creating a growth-to-hyphal bistable switch in Candida albicans

Recent biochemical experiments have indicated that in Candida albicans, a commensal fungal pathogen, the Ras signaling pathway plays a significant role in the yeast-to-hyphal transition; specifically, two enzymes in this pathway, Adenyl Cyclase Cyr1 and GTPase activating protein Ira2, facilitate this transition, in the presence of energy sensor ATP. However, the precise mechanism by which protein interactions between Ira2 and Cyr1 and the energy sensor ATP result in the yeast-to-hyphal transition and create a switch-like process are unknown. We propose a new set of biochemical reaction steps that captures all the essential interactions between Ira2, Cyr1, and ATP in the Ras pathway. With the help of chemical reaction network theory, we demonstrate that this set of biochemical reaction steps results in bistability. Further, bifurcation analysis of the differential equations based on this set of reaction steps supports the existence of a bistable switch, and this switch may act as a checkpoint mechanism for the promotion of growth-to-hyphal transition in C. albicans.

Multistationarity in Cyclic Sequestration-Transmutation Networks

We consider a natural class of reaction networks which consist of reactions where either two species can inactivate each other (i.e., sequestration), or some species can be transformed into another (i.e., transmutation), in a way that gives rise to a feedback cycle. We completely characterize the capacity of multistationarity of these networks. This is especially interesting because such networks provide simple examples of "atoms of multistationarity", i.e., minimal networks that can give rise to multiple positive steady states.

A robust method for designing multistable systems by embedding bistable subsystems

Although multistability is an important dynamic property of a wide range of complex systems, it is still a challenge to develop mathematical models for realising high order multistability using realistic regulatory mechanisms. To address this issue, we propose a robust method to develop multistable mathematical models by embedding bistable models together. Using the GATA1-GATA2-PU.1 module in hematopoiesis as the test system, we first develop a tristable model based on two bistable models without any high cooperative coefficients, and then modify the tristable model based on experimentally determined mechanisms. The modified model successfully realises four stable steady states and accurately reflects a recent experimental observation showing four transcriptional states. In addition, we develop a stochastic model, and stochastic simulations successfully realise the experimental observations in single cells. These results suggest that the proposed method is a general approach to develop mathematical models for realising multistability and heterogeneity in complex systems.

A single-cell resolved cell-cell communication model explains lineage commitment in hematopoiesis

Cells do not make fate decisions independently. Arguably, every cell-fate decision occurs in response to environmental signals. In many cases, cell-cell communication alters the dynamics of the internal gene regulatory network of a cell to initiate cell-fate transitions, yet models rarely take this into account. Here, we have developed a multiscale perspective to study the granulocyte-monocyte versus megakaryocyte-erythrocyte fate decisions. This transition is dictated by the GATA1-PU.1 network: a classical example of a bistable cell-fate system. We show that, for a wide range of cell communication topologies, even subtle changes in signaling can have pronounced effects on cell-fate decisions. We go on to show how cell-cell coupling through signaling can spontaneously break the symmetry of a homogenous cell population. Noise, both intrinsic and extrinsic, shapes the decision landscape profoundly, and affects the transcriptional dynamics underlying this important hematopoietic cell-fate decision-making system.

This article has an associated ‘The people behind the papers’ interview.

Summary: Through theory and computational modeling, cell-cell communication is revealed to be a crucial and under-appreciated determinant of cell-fate decision-making during hematopoiesis.

Bistability in fatty-acid oxidation resulting from substrate inhibition

In this study we demonstrated through analytic considerations and numerical studies that the mitochondrial fatty-acid β-oxidation can exhibit bistable-hysteresis behavior. In an experimentally validated computational model we identified a specific region in the parameter space in which two distinct stable and one unstable steady state could be attained with different fluxes. The two stable states were referred to as low-flux (disease) and high-flux (healthy) state. By a modular kinetic approach we traced the origin and causes of the bistability back to the distributive kinetics and the conservation of CoA, in particular in the last rounds of the β-oxidation. We then extended the model to investigate various interventions that may confer health benefits by activating the pathway, including (i) activation of the last enzyme MCKAT via its endogenous regulator p46-SHC protein, (ii) addition of a thioesterase (an acyl-CoA hydrolysing enzyme) as a safety valve, and (iii) concomitant activation of a number of upstream and downstream enzymes by short-chain fatty-acids (SCFA), metabolites that are produced from nutritional fibers in the gut. A high concentration of SCFAs, thioesterase activity, and inhibition of the p46Shc protein led to a disappearance of the bistability, leaving only the high-flux state. A better understanding of the switch behavior of the mitochondrial fatty-acid oxidation process between a low- and a high-flux state may lead to dietary and pharmacological intervention in the treatment or prevention of obesity and or non-alcoholic fatty-liver disease.

SWITCHES: Searchable web interface for topologies of CHEmical switches

 

Bistable biochemical switches are key motifs in cellular state decisions and long-term storage of cellular 'memory'. There are a few known biological switches that have been well characterized, however these examples are insufficient for systematic surveys of properties of these important systems. Here we present a resource of all possible bistable biochemical reaction networks with up to 6 reactions between 3 molecules, and 3 reactions between 4 molecules. Over 35,000 reaction topologies were constructed by identifying unique combinations of reactions between a fixed number of molecules. Then, these topologies were populated with rates within a biologically realistic range. The Searchable Web Interface for Topologies of CHEmical Switches (SWITCHES, https://switches.ncbs.res.in) provides a bistability and parameter analysis of over 7 million models from this systematic survey of chemical reaction space. This database will be useful for theoreticians interested in analyzing stability in chemical systems and also experimentalists for creating robust synthetic biological switches.

AVAILABILITY AND IMPLEMENTATION

Freely available on the web at https://switches.ncbs.res.in. Website implemented in PHP, MariaDB, Graphviz, and Apache, with all major browsers supported.

Dynamic bistable switches enhance robustness and accuracy of cell cycle transitions

Bistability is a common mechanism to ensure robust and irreversible cell cycle transitions. Whenever biological parameters or external conditions change such that a threshold is crossed, the system abruptly switches between different cell cycle states. Experimental studies have uncovered mechanisms that can make the shape of the bistable response curve change dynamically in time. Here, we show how such a dynamically changing bistable switch can provide a cell with better control over the timing of cell cycle transitions. Moreover, cell cycle oscillations built on bistable switches are more robust when the bistability is modulated in time. Our results are not specific to cell cycle models and may apply to other bistable systems in which the bistable response curve is time-dependent.

Many systems in nature show bistability, which means they can evolve to one of two stable steady states under exactly the same conditions. Which state they evolve to depends on where the system comes from. Such bistability underlies the switching behavior that is essential for cells to progress in the cell division cycle. A quick switch happens when the cell jumps from one steady state to another steady state. Typical of this switching behavior is its robustness and irreversibility. In this paper, we expand this viewpoint of the dynamics of the cell cycle by considering bistable switches which themselves are changing in time. This gives the cell an extra layer of control over transitions both in time and in space, and can make those transitions more robust. Such dynamically changing bistability can appear very naturally. We show this in a model of mitotic entry, in which we include a nuclear and cytoplasmic compartment. The activity of a crucial cell cycle protein follows a bistable switch in each compartment, but the shape of its response is changing in time as proteins are imported into and exported from the nucleus.

Entropic analysis of bistability and the general evolution criterion

An open flow bistable Schlögl model in a well-stirred isothermal reaction volume V showing the balance of the entropy production and the exchange entropy on the non-equilibrium stationary states.

On the definition of a self-sustaining chemical reaction system and its role in heredity

Background

The ability to self-sustain is one of the essential properties of life. However, a consistent and satisfying definition of self-sustainability is still missing. Currently, self-sustainability refers to either “no-intervention by a higher entity” or “regeneration of all the system’s components”. How to connect self-sustainability with heredity, another essential of life, is another problem, as they are often considered to be independent of each other. Last but not least, current definitions of self-sustainability failed to provide a practical method to empirically discern whether a chemical system is self-sustaining or not.

Results

Here I propose a definition of self-sustainability. It takes into account the chemical reaction network itself and the external environment which is simplified as a continuous-flow stirred tank reactor. One distinct property of self-sustaining systems is that the system can only proceed if molecular triggers (or called, seeds) are present initially. The molecular triggers are able to establish the whole system, indicating that they carry the preliminary heredity of the system. Consequently, life and a large group of fires (and other dissipative systems) can be distinguished. Besides, the general properties and various real-life examples of self-sustaining systems discussed here together indicate that self-sustaining systems are not uncommon.

Conclusions

The definition I proposed here naturally connects self-sustainability with heredity. As this definition involves the continuous-flow stirred tank reactor, it gives a simple way to empirically test whether a system is self-sustaining or not. Moreover, the general properties and various real-life examples of self-sustaining systems discussed here provide practical guidance on how to construct and detect such systems in real biology and chemistry.

Reviewers

This article was reviewed by Wentao Ma and David Baum.

Stochastically modeled weakly reversible reaction networks with a single linkage class

It has been known for nearly a decade that deterministically modeled reaction networks that are weakly reversible and consist of a single linkage class have trajectories that are bounded from both above and below by positive constants (so long as the initial condition has strictly positive components). It is conjectured that the stochastically modeled analogs of these systems are positive recurrent. We prove this conjecture in the affirmative under the following additional assumptions: (i) the system is binary, and (ii) for each species, there is a complex (vertex in the associated reaction diagram) that is a multiple of that species. To show this result, a new proof technique is developed in which we study the recurrence properties of the n-step embedded discrete-time Markov chain.

Periodic propagating waves coordinate RhoGTPase network dynamics at the leading and trailing edges during cell migration

Migrating cells need to coordinate distinct leading and trailing edge dynamics but the underlying mechanisms are unclear. Here, we combine experiments and mathematical modeling to elaborate the minimal autonomous biochemical machinery necessary and sufficient for this dynamic coordination and cell movement. RhoA activates Rac1 via DIA and inhibits Rac1 via ROCK, while Rac1 inhibits RhoA through PAK. Our data suggest that in motile, polarized cells, RhoA–ROCK interactions prevail at the rear, whereas RhoA-DIA interactions dominate at the front where Rac1/Rho oscillations drive protrusions and retractions. At the rear, high RhoA and low Rac1 activities are maintained until a wave of oscillatory GTPase activities from the cell front reaches the rear, inducing transient GTPase oscillations and RhoA activity spikes. After the rear retracts, the initial GTPase pattern resumes. Our findings show how periodic, propagating GTPase waves coordinate distinct GTPase patterns at the leading and trailing edge dynamics in moving cells.

Detection of biological switches using the method of Gröebner bases

BACKGROUND

Bistability and ability to switch between two stable states is the hallmark of cellular responses. Cellular signaling pathways often contain bistable switches that regulate the transmission of the extracellular information to the nucleus where important biological functions are executed.

RESULTS

In this work we show how the method of Gröebner bases can be used to detect bistability and output switchability. The method of Gröebner bases can be seen as a multivariate, non-linear generalization of the Gaussian elimination for linear systems which conveniently seperates the variables and drastically simplifies the simultaneous solution of polynomial equations. A necessary condition for fixed-point state bistability is for the Gröbner basis to have three distinct solutions for the state. A sufficient condition is provided by the eigenvalues of the local Jacobians. We also introduce the concept of output switchability which is defined as the ability of an output of a bistable system to switch between two different stable steady-state values. It is shown that bistability does not necessarily guarantee switchability of every state variable of the system. We further show that, for a bistable system, the necessary conditions for output switchability can be derived using the Gröebner basis. The theoretical results are incorporated into an analysis procedure and applied to several systems including the AKT (Protein kinase B), RAS (Rat Sarcoma) and MAPK (Mitogen-activated protein kinase) signal transduction pathways. Results demonstrate that the Gröebner bases can be conveniently used to analyze biological switches by simultaneously detecting bistability and output switchability.

CONCLUSION

The Gröebner bases provides a novel methodology to analyze bistability. Results clarify the distinction between bistability and output switchability which is lacking in the literature. We have shown that theoretically, it is possible to have an output subspace of an n-dimensional bistable system where certain variables cannot switch. It is possible to construct such systems as we have done with two reaction networks.

Hysteretic Genetic Circuit for Detection of Proteasomal Degradation in Mammalian Cells

Synthetic hysteretic mammalian gene circuits generating sustained cellular responses to transient perturbations provide important tools to investigate complex cellular behaviors and reprogram cells for a variety of applications, ranging from protein production to cell fate decisions. The design rules of synthetic gene circuits with controlled hysteretic behaviors, however, remain uncharacterized. To identify the criteria for achieving predictable control of hysteresis, we built a genetic circuit for detection of proteasomal degradation (Hys-Deg). The Hys-Deg circuit is based on a tetracycline-controlled transactivator (tTA) variant engineered to interface with the ubiquitin proteasome system (UPS). The tTA variant activates its own expression, generating a positive feedback loop that is triggered by expression of another tTA gene that is constitutively regulated. Guided by predictive modeling, we characterized the hysteretic response of the Hys-Deg circuit. We demonstrated that control of the hysteretic response is achieved by modulating the ratio of expression of constitutive to inducible tTA. We also showed that the system can be finely tuned through dosage of the inducer tetracycline to calibrate the circuit for detection of the desired levels of UPS activation. This study establishes the design rules for building a hysteretic genetic circuit with an autoregulatory feedback loop and provides a synthetic memory module that could be easily integrated into regulatory gene networks to study and engineer complex cellular behaviors.

Dichotomous Nature of Bistability Generated by Negative Cooperativity in Receptor-Ligand Binding

Positive cooperativity in receptor-ligand binding plays an important role in cell signaling as it generates an ultrasensitive response, a requirement for nonlinear phenomena such as bistability and oscillations in feedback regulated reaction networks. On the other hand, negative cooperativity typically produces a hyperbolic response and is thus less explored. However, recently negative cooperativity was shown to generate an ultrasensitive response under the condition of strong ligand affinity. In this work, we have used mathematical modeling to investigate the effect of negative cooperativity in receptor-ligand interaction on the bistability in a positive feedback regulatory motif. We systematically investigated the effect of negative cooperativity, modifying the two equilibrium constants of the receptor-ligand binding, on the robustness and tunability of bistability. We show that in the regime where negative cooperativity exhibits robust bistability, positive cooperativity results in poor bistability and vice versa. Further we find that the robustness and tunability of bistability depend crucially on the stability of singly and doubly engaged receptors. Our modeling highlights the ability of negative cooperativity to produce complex phenomena with potential applications in designing synthetic devices or in explaining experimental observations in cell biology.

Genetically Encoded Fluorescent Biosensors Illuminate the Spatiotemporal Regulation of Signaling Networks

Cellular signaling networks are the foundation which determines the fate and function of cells as they respond to various cues and stimuli. The discovery of fluorescent proteins over 25 years ago enabled the development of a diverse array of genetically encodable fluorescent biosensors that are capable of measuring the spatiotemporal dynamics of signal transduction pathways in live cells. In an effort to encapsulate the breadth over which fluorescent biosensors have expanded, we endeavored to assemble a comprehensive list of published engineered biosensors, and we discuss many of the molecular designs utilized in their development. Then, we review how the high temporal and spatial resolution afforded by fluorescent biosensors has aided our understanding of the spatiotemporal regulation of signaling networks at the cellular and subcellular level. Finally, we highlight some emerging areas of research in both biosensor design and applications that are on the forefront of biosensor development.

A mechanism for bistability in glycosylation

Glycosyltransferases are a class of enzymes that catalyse the posttranslational modification of proteins to produce a large number of glycoconjugate acceptors from a limited number of nucleotide-sugar donors. The products of one glycosyltransferase can be the substrates of several other enzymes, causing a combinatorial explosion in the number of possible glycan products. The kinetic behaviour of systems where multiple acceptor substrates compete for a single enzyme is presented, and the case in which high concentrations of an acceptor substrate are inhibitory as a result of abortive complex formation, is shown to result in non-Michaelian kinetics that can lead to bistability in an open system. A kinetic mechanism is proposed that is consistent with the available experimental evidence and provides a possible explanation for conflicting observations on the β-1,4-galactosyltransferases. Abrupt switching between steady states in networks of glycosyltransferase-catalysed reactions may account for the observed changes in glycosyl-epitopes in cancer cells.

Subharmonics and Chaos in Simple Periodically Forced Biomolecular Models

This article uncovers a remarkable behavior in two biochemical systems that commonly appear as components of signal transduction pathways in systems biology. These systems have globally attracting steady states when unforced, so they might have been considered uninteresting from a dynamical standpoint. However, when subject to a periodic excitation, strange attractors arise via a period-doubling cascade. Quantitative analyses of the corresponding discrete chaotic trajectories are conducted numerically by computing largest Lyapunov exponents, power spectra, and autocorrelation functions. To gain insight into the geometry of the strange attractors, the phase portraits of the corresponding iterated maps are interpreted as scatter plots for which marginal distributions are additionally evaluated. The lack of entrainment to external oscillations, in even the simplest biochemical networks, represents a level of additional complexity in molecular biology, which has previously been insufficiently recognized but is plausibly biologically important.

Memory functions reveal structural properties of gene regulatory networks

Gene regulatory networks (GRNs) control cellular function and decision making during tissue development and homeostasis. Mathematical tools based on dynamical systems theory are often used to model these networks, but the size and complexity of these models mean that their behaviour is not always intuitive and the underlying mechanisms can be difficult to decipher. For this reason, methods that simplify and aid exploration of complex networks are necessary. To this end we develop a broadly applicable form of the Zwanzig-Mori projection. By first converting a thermodynamic state ensemble model of gene regulation into mass action reactions we derive a general method that produces a set of time evolution equations for a subset of components of a network. The influence of the rest of the network, the bulk, is captured by memory functions that describe how the subnetwork reacts to its own past state via components in the bulk. These memory functions provide probes of near-steady state dynamics, revealing information not easily accessible otherwise. We illustrate the method on a simple cross-repressive transcriptional motif to show that memory functions not only simplify the analysis of the subnetwork but also have a natural interpretation. We then apply the approach to a GRN from the vertebrate neural tube, a well characterised developmental transcriptional network composed of four interacting transcription factors. The memory functions reveal the function of specific links within the neural tube network and identify features of the regulatory structure that specifically increase the robustness of the network to initial conditions. Taken together, the study provides evidence that Zwanzig-Mori projections offer powerful and effective tools for simplifying and exploring the behaviour of GRNs.

Gene regulatory networks are essential for cell fate specification and function. But the recursive links that comprise these networks often make determining their properties and behaviour complicated. Computational models of these networks can also be difficult to decipher. To reduce the complexity of such models we employ a Zwanzig-Mori projection approach. This allows a system of ordinary differential equations, representing a network, to be reduced to an arbitrary subnetwork consisting of part of the initial network, with the rest of the network (bulk) captured by memory functions. These memory functions account for the bulk by describing signals that return to the subnetwork after some time, having passed through the bulk. We show how this approach can be used to simplify analysis and to probe the behaviour of a gene regulatory network. Applying the method to a transcriptional network in the vertebrate neural tube reveals previously unappreciated properties of the network. By taking advantage of the structure of the memory functions we identify interactions within the network that are unnecessary for sustaining correct patterning. Upon further investigation we find that these interactions are important for conferring robustness to variation in initial conditions. Taken together we demonstrate the validity and applicability of the Zwanzig-Mori projection approach to gene regulatory networks.

Open flow non-enzymatic template catalysis and replication

Fully reversible template assisted ligation and replication leads to high template yields in open flow reactors.

Flows, scaling, and the control of moment hierarchies for stochastic chemical reaction networks

Stochastic chemical reaction networks (CRNs) are complex systems that combine the features of concurrent transformation of multiple variables in each elementary reaction event and nonlinear relations between states and their rates of change. Most general results concerning CRNs are limited to restricted cases where a topological characteristic known as deficiency takes a value 0 or 1, implying uniqueness and positivity of steady states and surprising, low-information forms for their associated probability distributions. Here we derive equations of motion for fluctuation moments at all orders for stochastic CRNs at general deficiency. We show, for the standard base case of proportional sampling without replacement (which underlies the mass-action rate law), that the generator of the stochastic process acts on the hierarchy of factorial moments with a finite representation. Whereas simulation of high-order moments for many-particle systems is costly, this representation reduces the solution of moment hierarchies to a complexity comparable to solving a heat equation. At steady states, moment hierarchies for finite CRNs interpolate between low-order and high-order scaling regimes, which may be approximated separately by distributions similar to those for deficiency-zero networks and connected through matched asymptotic expansions. In CRNs with multiple stable or metastable steady states, boundedness of high-order moments provides the starting condition for recursive solution downward to low-order moments, reversing the order usually used to solve moment hierarchies. A basis for a subset of network flows defined by having the same mean-regressing property as the flows in deficiency-zero networks gives the leading contribution to low-order moments in CRNs at general deficiency, in a 1/n expansion in large particle numbers. Our results give a physical picture of the different informational roles of mean-regressing and non-mean-regressing flows and clarify the dynamical meaning of deficiency not only for first-moment conditions but for all orders in fluctuations.

Analysis of network motifs in cellular regulation: Structural similarities, input-output relations and signal integration

Much of the complexity of regulatory networks derives from the necessity to integrate multiple signals and to avoid malfunction due to cross-talk or harmful perturbations. Hence, one may expect that the input-output behavior of larger networks is not necessarily more complex than that of smaller network motifs which suggests that both can, under certain conditions, be described by similar equations. In this review, we illustrate this approach by discussing the similarities that exist in the steady state descriptions of a simple bimolecular reaction, covalent modification cycles and bacterial two-component systems. Interestingly, in all three systems fundamental input-output characteristics such as thresholds, ultrasensitivity or concentration robustness are described by structurally similar equations. Depending on the system the meaning of the parameters can differ ranging from protein concentrations and affinity constants to complex parameter combinations which allows for a quantitative understanding of signal integration in these systems. We argue that this approach may also be extended to larger regulatory networks.

Core signalling motif displaying multistability through multi-state enzymes

Bistability, and more generally multistability, is a key system dynamics feature enabling decision-making and memory in cells. Deciphering the molecular determinants of multistability is thus crucial for a better understanding of cellular pathways and their (re)engineering in synthetic biology. Here, we show that a key motif found predominantly in eukaryotic signalling systems, namely a futile signalling cycle, can display bistability when featuring a two-state kinase. We provide necessary and sufficient mathematical conditions on the kinetic parameters of this motif that guarantee the existence of multiple steady states. These conditions foster the intuition that bistability arises as a consequence of competition between the two states of the kinase. Extending from this result, we find that increasing the number of kinase states linearly translates into an increase in the number of steady states in the system. These findings reveal, to our knowledge, a new mechanism for the generation of bistability and multistability in cellular signalling systems. Further the futile cycle featuring a two-state kinase is among the smallest bistable signalling motifs. We show that multi-state kinases and the described competition-based motif are part of several natural signalling systems and thereby could enable them to implement complex information processing through multistability. These results indicate that multi-state kinases in signalling systems are readily exploited by natural evolution and could equally be used by synthetic approaches for the generation of multistable information processing systems at the cellular level.

Cell-cycle transitions: a common role for stoichiometric inhibitors

The abrupt and irreversible transitions that drive cells through the DNA replication-division cycle are governed by molecular mechanisms that function as bistable “toggle” switches. A common theme of these switches is a network motif consisting of a “beleaguered” enzyme and its “domineering” substrate, locked in a feedback amplification loop.

The cell division cycle is the process by which eukaryotic cells replicate their chromosomes and partition them to two daughter cells. To maintain the integrity of the genome, proliferating cells must be able to block progression through the division cycle at key transition points (called “checkpoints”) if there have been problems in the replication of the chromosomes or their biorientation on the mitotic spindle. These checkpoints are governed by protein-interaction networks, composed of phase-specific cell-cycle activators and inhibitors. Examples include Cdk1:Clb5 and its inhibitor Sic1 at the G1/S checkpoint in budding yeast, APC:Cdc20 and its inhibitor MCC at the mitotic checkpoint, and PP2A:B55 and its inhibitor, alpha-endosulfine, at the mitotic-exit checkpoint. Each of these inhibitors is a substrate as well as a stoichiometric inhibitor of the cell-cycle activator. Because the production of each inhibitor is promoted by a regulatory protein that is itself inhibited by the cell-cycle activator, their interaction network presents a regulatory motif characteristic of a “feedback-amplified domineering substrate” (FADS). We describe how the FADS motif responds to signals in the manner of a bistable toggle switch, and then we discuss how this toggle switch accounts for the abrupt and irreversible nature of three specific cell-cycle checkpoints.

Bistability and Bifurcation in Minimal Self-Replication and Nonenzymatic Catalytic Networks

Bistability and bifurcation, found in a wide range of biochemical networks, are central to the proper function of living systems. We investigate herein recent model systems that show bistable behavior based on nonenzymatic self-replication reactions. Such models were used before to investigate catalytic growth, chemical logic operations, and additional processes of self-organization leading to complexification. By solving for their steady-state solutions by using various analytical and numerical methods, we analyze how and when these systems yield bistability and bifurcation and discover specific cases and conditions producing bistability. We demonstrate that the onset of bistability requires at least second-order catalysis and results from a mismatch between the various forward and reverse processes. Our findings may have far-reaching implications in understanding early evolutionary processes of complexification, emergence, and potentially the origin of life.

A new efficient approach to fit stochastic models on the basis of high-throughput experimental data using a model of IRF7 gene expression as case study

Background

Mathematical models are used to gain an integrative understanding of biochemical processes and networks. Commonly the models are based on deterministic ordinary differential equations. When molecular counts are low, stochastic formalisms like Monte Carlo simulations are more appropriate and well established. However, compared to the wealth of computational methods used to fit and analyze deterministic models, there is only little available to quantify the exactness of the fit of stochastic models compared to experimental data or to analyze different aspects of the modeling results.

Results

Here, we developed a method to fit stochastic simulations to experimental high-throughput data, meaning data that exhibits distributions. The method uses a comparison of the probability density functions that are computed based on Monte Carlo simulations and the experimental data. Multiple parameter values are iteratively evaluated using optimization routines. The method improves its performance by selecting parameters values after comparing the similitude between the deterministic stability of the system and the modes in the experimental data distribution. As a case study we fitted a model of the IRF7 gene expression circuit to time-course experimental data obtained by flow cytometry. IRF7 shows bimodal dynamics upon IFN stimulation. This dynamics occurs due to the switching between active and basal states of the IRF7 promoter. However, the exact molecular mechanisms responsible for the bimodality of IRF7 is not fully understood.

Conclusions

Our results allow us to conclude that the activation of the IRF7 promoter by the combination of IRF7 and ISGF3 is sufficient to explain the observed bimodal dynamics.

Electronic supplementary material

The online version of this article (doi:10.1186/s12918-017-0406-4) contains supplementary material, which is available to authorized users.

Stoichiometric network analysis of spontaneous mirror symmetry breaking in chemical reactions

Stoichiometric network analysis (SNA) is used to study spontaneous mirror symmetry breaking in chemical reaction schemes.

Splitting nodes and linking channels: A method for assembling biocircuits from stochastic elementary units

Akin to electric circuits, we construct biocircuits that are manipulated by cutting and assembling channels through which stochastic information flows. This diagrammatic manipulation allows us to create a method which constructs networks by joining building blocks selected so that (a) they cover only basic processes; (b) it is scalable to large networks; (c) the mean and variance-covariance from the Pauli master equation form a closed system; and (d) given the initial probability distribution, no special boundary conditions are necessary to solve the master equation. The method aims to help with both designing new synthetic signaling pathways and quantifying naturally existing regulatory networks.

Mathematical and Statistical Techniques for Systems Medicine: The Wnt Signaling Pathway as a Case Study

The last decade has seen an explosion in models that describe phenomena in systems medicine. Such models are especially useful for studying signaling pathways, such as the Wnt pathway. In this chapter we use the Wnt pathway to showcase current mathematical and statistical techniques that enable modelers to gain insight into (models of) gene regulation and generate testable predictions. We introduce a range of modeling frameworks, but focus on ordinary differential equation (ODE) models since they remain the most widely used approach in systems biology and medicine and continue to offer great potential. We present methods for the analysis of a single model, comprising applications of standard dynamical systems approaches such as nondimensionalization, steady state, asymptotic and sensitivity analysis, and more recent statistical and algebraic approaches to compare models with data. We present parameter estimation and model comparison techniques, focusing on Bayesian analysis and coplanarity via algebraic geometry. Our intention is that this (non-exhaustive) review may serve as a useful starting point for the analysis of models in systems medicine.

Implications of network topology on stability

In analogy to chemical reaction networks, I demonstrate the utility of expressing the governing equations of an arbitrary dynamical system (interaction network) as sums of real functions (generalized reactions) multiplied by real scalars (generalized stoichiometries) for analysis of its stability. The reaction stoichiometries and first derivatives define the network's "influence topology", a signed directed bipartite graph. Parameter reduction of the influence topology permits simplified expression of the principal minors (sums of products of non-overlapping bipartite cycles) and Hurwitz determinants (sums of products of the principal minors or the bipartite cycles directly) for assessing the network's steady state stability. Visualization of the Hurwitz determinants over the reduced parameters defines the network's stability phase space, delimiting the range of its dynamics (specifically, the possible numbers of unstable roots at each steady state solution). Any further explicit algebraic specification of the network will project onto this stability phase space. Stability analysis via this hierarchical approach is demonstrated on classical networks from multiple fields.

Laplacian Dynamics with Synthesis and Degradation

Analyzing qualitative behaviors of biochemical reactions using its associated network structure has proven useful in diverse branches of biology. As an extension of our previous work, we introduce a graph-based framework to calculate steady state solutions of biochemical reaction networks with synthesis and degradation. Our approach is based on a labeled directed graph G and the associated system of linear non-homogeneous differential equations with first-order degradation and zeroth-order synthesis. We also present a theorem which provides necessary and sufficient conditions for the dynamics to engender a unique stable steady state. Although the dynamics are linear, one can apply this framework to nonlinear systems by encoding nonlinearity into the edge labels. We answer an open question from our previous work concerning the non-positiveness of the elements in the inverse of a perturbed Laplacian matrix. Moreover, we provide a graph theoretical framework for the computation of the inverse of such a matrix. This also completes our previous framework and makes it purely graph theoretical. Lastly, we demonstrate the utility of this framework by applying it to a mathematical model of insulin secretion through ion channels in pancreatic β-cells.

Parameter-free methods distinguish Wnt pathway models and guide design of experiments

The canonical Wnt signaling pathway, mediated by β-catenin, is crucially involved in development, adult stem cell tissue maintenance, and a host of diseases including cancer. We analyze existing mathematical models of Wnt and compare them to a new Wnt signaling model that targets spatial localization; our aim is to distinguish between the models and distill biological insight from them. Using Bayesian methods we infer parameters for each model from mammalian Wnt signaling data and find that all models can fit this time course. We appeal to algebraic methods (concepts from chemical reaction network theory and matroid theory) to analyze the models without recourse to specific parameter values. These approaches provide insight into aspects of Wnt regulation: the new model, via control of shuttling and degradation parameters, permits multiple stable steady states corresponding to stem-like vs. committed cell states in the differentiation hierarchy. Our analysis also identifies groups of variables that should be measured to fully characterize and discriminate between competing models, and thus serves as a guide for performing minimal experiments for model comparison.

Sharper graph-theoretical conditions for the stabilization of complex reaction networks

Across the landscape of all possible chemical reaction networks there is a surprising degree of stable behavior, despite what might be substantial complexity and nonlinearity in the governing differential equations. At the same time there are reaction networks, in particular those that arise in biology, for which richer behavior is exhibited. Thus, it is of interest to understand network-structural features whose presence enforces dull, stable behavior and whose absence permits the dynamical richness that might be necessary for life. We present conditions on a network's Species-Reaction Graph that ensure a high degree of stable behavior, so long as the kinetic rate functions satisfy certain weak and natural constraints. These graph-theoretical conditions are considerably more incisive than those reported earlier.

A technique for determining the signs of sensitivities of steady states in chemical reaction networks

This paper studies the direction of change of steady states to parameter perturbations in chemical reaction networks, and, in particular, to changes in conserved quantities. Theoretical considerations lead to the formulation of a computational procedure that provides a set of possible signs of such sensitivities. The procedure is purely algebraic and combinatorial, only using information on stoichiometry, and is independent of the values of kinetic constants. Three examples of important intracellular signal transduction models are worked out as an illustration. In these examples, the set of signs found is minimal, but there is no general guarantee that the set found will always be minimal in other examples. The paper also briefly discusses the relationship of the sign problem to the question of uniqueness of steady states in stoichiometry classes.

A structural classification of candidate oscillatory and multistationary biochemical systems

Molecular systems are uncertain: The variability of reaction parameters and the presence of unknown interactions can weaken the predictive capacity of solid mathematical models. However, strong conclusions on the admissible dynamic behaviors of a model can often be achieved without detailed knowledge of its specific parameters. In systems with a sign-definite Jacobian, for instance, cycle-based criteria related to the famous Thomas' conjectures have been largely used to characterize oscillatory and multistationary dynamic outcomes. We build on the rich literature focused on the identification of potential oscillatory and multistationary behaviors using parameter-free criteria. We propose a classification for sign-definite non-autocatalytic biochemical networks, which summarizes several existing results in the literature. We call weak (strong) candidate oscillators systems which can possibly (exclusively) transition to instability due to the presence of a complex pair of eigenvalues, while we call weak (strong) candidate multistationary systems those which can possibly (exclusively) transition to instability due to the presence of a real eigenvalue. For each category, we provide a characterization based on the exclusive or simultaneous presence of positive and negative cycles in the associated sign graph. Most realistic examples of biochemical networks fall in the gray area of systems in which both positive and negative cycles are present: Therefore, both oscillatory and bistable behaviors are in principle possible. However, many canonical example circuits exhibiting oscillations or bistability fall in the categories of strong candidate oscillators/multistationary systems, in agreement with our results.

Morphisms of reaction networks that couple structure to function

BACKGROUND

The mechanisms underlying complex biological systems are routinely represented as networks. Network kinetics is widely studied, and so is the connection between network structure and behavior. However, similarity of mechanism is better revealed by relationships between network structures.

RESULTS

We define morphisms (mappings) between reaction networks that establish structural connections between them. Some morphisms imply kinetic similarity, and yet their properties can be checked statically on the structure of the networks. In particular we can determine statically that a complex network will emulate a simpler network: it will reproduce its kinetics for all corresponding choices of reaction rates and initial conditions. We use this property to relate the kinetics of many common biological networks of different sizes, also relating them to a fundamental population algorithm.

CONCLUSIONS

Structural similarity between reaction networks can be revealed by network morphisms, elucidating mechanistic and functional aspects of complex networks in terms of simpler networks.