ROBUSTNESS, STABILITY AND EFFICIENCY OF PHAGE λ GENETIC SWITCH: DYNAMICAL STRUCTURE ANALYSIS

Deciphering the topological landscape of glioma using a network theory framework

Glioma stem cells have been recognized as key players in glioma recurrence and therapeutic resistance, presenting a promising target for novel treatments. However, the limited understanding of the role glioma stem cells play in the glioma hierarchy has drawn controversy and hindered research translation into therapies. Despite significant advances in our understanding of gene regulatory networks, the dynamics of these networks and their implications for glioma remain elusive. This study employs a systemic theoretical perspective to integrate experimental knowledge into a core endogenous network model for glioma, thereby elucidating its energy landscape through network dynamics computation. The model identifies two stable states corresponding to astrocytic-like and oligodendrocytic-like tumor cells, connected by a transition state with the feature of high stemness, which serves as one of the energy barriers between astrocytic-like and oligodendrocytic-like states, indicating the instability of glioma stem cells in vivo. We also obtained various stable states further supporting glioma's multicellular origins and uncovered a group of transition states that could potentially induce tumor heterogeneity and therapeutic resistance. This research proposes that the transition states linking both glioma stable states are central to glioma heterogeneity and therapy resistance. Our approach may contribute to the advancement of glioma therapy by offering a novel perspective on the complex landscape of glioma biology.

Evolution of Telencephalon Anterior-Posterior Patterning through Core Endogenous Network Bifurcation

How did the complex structure of the telencephalon evolve? Existing explanations are based on phenomena and lack a first-principles account. The Darwinian dynamics and endogenous network theory-established decades ago-provides a mathematical and theoretical framework and a general constitutive structure for theory-experiment coupling for answering this question from a first-principles perspective. By revisiting a gene network that explains the anterior-posterior patterning of the vertebrate telencephalon, we found that upon increasing the cooperative effect within this network, fixed points gradually evolve, accompanied by the occurrence of two bifurcations. The dynamic behavior of this network is informed by the knowledge obtained from experiments on telencephalic evolution. Our work provides a quantitative explanation for how telencephalon anterior-posterior patterning evolved from the pre-vertebrate chordate to the vertebrate and provides a series of verifiable predictions from a first-principles perspective.

Uncovering the cell fate decision in lysis-lysogeny transition and stem cell development via Markov state modeling

Understanding the behavior of a complex gene regulatory network is a fundamental but challenging task in systems biology. How to reduce the large number of degrees of freedom of a specific network and identify its main biological pathway is the key issue. In this paper, we utilized the transition path theory (TPT) and Markov state modeling (MSM) framework to numerically study two typical cell fate decision processes: the lysis-lysogeny transition and stem cell development. The application of TPT to the lysis-lysogeny decision-making system reveals that the competitions of CI and Cro dimer binding play the major role in determining the cell fates. We also quantified the transition rates from the lysogeny to lysis state under different conditions. The overall computational results are consistent with biological intuitions but with more detailed information. For the stem cell developmental system, we applied the MSM to reduce the original dynamics to a moderate-size Markov chain. Further spectral analysis showed that the reduced system exhibits nine metastable states, which correspond to the refinement of the five known typical cell types in development. We further investigated the dominant transition pathways corresponding to the cell differentiation, reprogramming, and trans-differentiation. A similar approach can be applied to study other biological systems.

Noise distorts the epigenetic landscape and shapes cell-fate decisions

The Waddington epigenetic landscape has become an iconic representation of the cellular differentiation process. Recent single-cell transcriptomic data provide new opportunities for quantifying this originally conceptual tool, offering insight into the gene regulatory networks underlying cellular development. While many methods for constructing the landscape have been proposed, by far the most commonly employed approach is based on computing the landscape as the negative logarithm of the steady-state probability distribution. Here, we use simple models to highlight the complexities and limitations that arise when reconstructing the potential landscape in the presence of stochastic fluctuations. We consider how the landscape changes in accordance with different stochastic systems and show that it is the subtle interplay between the deterministic and stochastic components of the system that ultimately shapes the landscape. We further discuss how the presence of noise has important implications for the identifiability of the regulatory dynamics from experimental data. A record of this paper's transparent peer review process is included in the supplemental information.

Exact Probability Landscapes of Stochastic Phenotype Switching in Feed-Forward Loops: Phase Diagrams of Multimodality

Feed-forward loops (FFLs) are among the most ubiquitously found motifs of reaction networks in nature. However, little is known about their stochastic behavior and the variety of network phenotypes they can exhibit. In this study, we provide full characterizations of the properties of stochastic multimodality of FFLs, and how switching between different network phenotypes are controlled. We have computed the exact steady-state probability landscapes of all eight types of coherent and incoherent FFLs using the finite-butter Accurate Chemical Master Equation (ACME) algorithm, and quantified the exact topological features of their high-dimensional probability landscapes using persistent homology. Through analysis of the degree of multimodality for each of a set of 10,812 probability landscapes, where each landscape resides over 10⁵–10⁶ microstates, we have constructed comprehensive phase diagrams of all relevant behavior of FFL multimodality over broad ranges of input and regulation intensities, as well as different regimes of promoter binding dynamics. In addition, we have quantified the topological sensitivity of the multimodality of the landscapes to regulation intensities. Our results show that with slow binding and unbinding dynamics of transcription factor to promoter, FFLs exhibit strong stochastic behavior that is very different from what would be inferred from deterministic models. In addition, input intensity play major roles in the phenotypes of FFLs: At weak input intensity, FFL exhibit monomodality, but strong input intensity may result in up to 6 stable phenotypes. Furthermore, we found that gene duplication can enlarge stable regions of specific multimodalities and enrich the phenotypic diversity of FFL networks, providing means for cells toward better adaptation to changing environment. Our results are directly applicable to analysis of behavior of FFLs in biological processes such as stem cell differentiation and for design of synthetic networks when certain phenotypic behavior is desired.

Cell plasticity in cancer cell populations

In this review, we propose a recension of biological observations on plasticity in cancer cell populations and discuss theoretical considerations about their mechanisms.

Host-virus evolutionary dynamics with specialist and generalist infection strategies: Bifurcations, bistability, and chaos

In this work, we have investigated the evolutionary dynamics of a generalist pathogen, e.g., a virus population, that evolves toward specialization in an environment with multiple host types. We have particularly explored under which conditions generalist viral strains may rise in frequency and coexist with specialist strains or even dominate the population. By means of a nonlinear mathematical model and bifurcation analysis, we have determined the theoretical conditions for stability of nine identified equilibria and provided biological interpretation in terms of the infection rates for the viral specialist and generalist strains. By means of a stability diagram, we identified stable fixed points and stable periodic orbits, as well as regions of bistability. For arbitrary biologically feasible initial population sizes, the probability of evolving toward stable solutions is obtained for each point of the analyzed parameter space. This probability map shows combinations of infection rates of the generalist and specialist strains that might lead to equal chances for each type becoming the dominant strategy. Furthermore, we have identified infection rates for which the model predicts the onset of chaotic dynamics. Several degenerate Bogdanov-Takens and zero-Hopf bifurcations are detected along with generalized Hopf and zero-Hopf bifurcations. This manuscript provides additional insights into the dynamical complexity of host-pathogen evolution toward different infection strategies.

Dynamics of gene regulatory networks and their dependence on network topology and quantitative parameters - the case of phage λ

BACKGROUND

Gene regulatory networks can be modelled in various ways depending on the level of detail required and biological questions addressed. One of the earliest formalisms used for modeling is a Boolean network, although these models cannot describe most temporal aspects of a biological system. Differential equation models have also been used to model gene regulatory networks, but these frameworks tend to be too detailed for large models and many quantitative parameters might not be deducible in practice. Hybrid models bridge the gap between these two model classes - these are useful when concentration changes are important while the information about precise concentrations and binding site affinities is partial.

RESULTS

In this paper we study the stable behaviours of phage λ via a hybrid system based model. We identify wild type and mutant behaviours that arise for various orderings of binding site affinities. We propose experiments for detecting these behaviours: we suggest several ways of altering binding affinities with either mutations or genome rearrangements to achieve modified behaviours. The feasibility of these experiments is assessed. The interplay between the qualitative aspects of a network, e.g. network topology, and quantitative parameters, e.g. growth and degradation rates of proteins, is demonstrated. We also provide a software for exploring all feasible states of a hybrid system model and identifying all attractors.

CONCLUSIONS

The behaviours of phage λ are determined mainly by the topology of this network and by the mutual order of binding affinities. Exact affinities and growth and degradation rates of proteins fine tune the system. We show that only two stable behaviours are possible for phage λ if the main constraints of λ switch are preserved - these behaviours correspond to lysis and lysogeny. We identify several variants of both lysis and lysogeny - one wild type and one modified behaviour for each. We elucidate the necessary constraints for binding site affinities to achieve both wild type lysis and lysogeny. Our software is applicable to a wide range of biological models described as a hybrid system.

Quantitative Modelling of the Waddington Epigenetic Landscape

C.H. Waddington introduced the epigenetic landscape as a metaphor to represent cellular decision-making during development. Like a population of balls rolling down a rough hillside, developing cells follow specific trajectories (valleys) and eventually come to rest in one or another low-energy state that represents a mature cell type. Waddington depicted the topography of this landscape as determined by interactions among gene products, thereby connecting genotype to phenotype. In modern terms, each point on the landscape represents a state of the underlying genetic regulatory network, which in turn is described by a gene expression profile. In this chapter we demonstrate how the mathematical formalism of Hopfield networks can be used to model this epigenetic landscape. Hopfield networks are auto-associative artificial neural networks; input patterns are stored as attractors of the network and can be recalled from noisy or incomplete inputs. The resulting models capture the temporal dynamics of a gene regulatory network, yielding quantitative insight into cellular development and phenotype.

Probabilistic control of HIV latency and transactivation by the Tat gene circuit

The reservoir of HIV latently infected cells is the major obstacle for eradication of HIV infection. The "shock-and-kill" strategy proposed earlier aims to reduce the reservoir by activating cells out of latency. While the intracellular HIV Tat gene circuit is known to play important roles in controlling latency and its transactivation in HIV-infected cells, the detailed control mechanisms are not well understood. Here we study the mechanism of probabilistic control of the latent and the transactivated cell phenotypes of HIV-infected cells. We reconstructed the probability landscape, which is the probability distribution of the Tat gene circuit states, by directly computing the exact solution of the underlying chemical master equation. Results show that the Tat circuit exhibits a clear bimodal probability landscape (i.e., there are two distinct probability peaks, one associated with the latent cell phenotype and the other with the transactivated cell phenotype). We explore potential modifications to reactions in the Tat gene circuit for more effective transactivation of latent cells (i.e., the shock-and-kill strategy). Our results suggest that enhancing Tat acetylation can dramatically increase Tat and viral production, while increasing the Tat-transactivation response binding affinity can transactivate latent cells more rapidly than other manipulations. Our results further explored the "block and lock" strategy toward a functional cure for HIV. Overall, our study demonstrates a general approach toward discovery of effective therapeutic strategies and druggable targets by examining control mechanisms of cell phenotype switching via exactly computed probability landscapes of reaction networks.

Endogenous Molecular-Cellular Network Cancer Theory: A Systems Biology Approach

In light of ever apparent limitation of the current dominant cancer mutation theory, a quantitative hypothesis for cancer genesis and progression, endogenous molecular-cellular network hypothesis has been proposed from the systems biology perspective, now for more than 10 years. It was intended to include both the genetic and epigenetic causes to understand cancer. Its development enters the stage of meaningful interaction with experimental and clinical data and the limitation of the traditional cancer mutation theory becomes more evident. Under this endogenous network hypothesis, we established a core working network of hepatocellular carcinoma (HCC) according to the hypothesis and quantified the working network by a nonlinear dynamical system. We showed that the two stable states of the working network reproduce the main known features of normal liver and HCC at both the modular and molecular levels. Using endogenous network hypothesis and validated working network, we explored genetic mutation pattern in cancer and potential strategies to cure or relieve HCC from a totally new perspective. Patterns of genetic mutations have been traditionally analyzed by posteriori statistical association approaches in light of traditional cancer mutation theory. One may wonder the possibility of a priori determination of any mutation regularity. Here, we found that based on the endogenous network theory the features of genetic mutations in cancers may be predicted without any prior knowledge of mutation propensities. Normal hepatocyte and cancerous hepatocyte stable states, specified by distinct patterns of expressions or activities of proteins in the network, provide means to directly identify a set of most probable genetic mutations and their effects in HCC. As the key proteins and main interactions in the network are conserved through cell types in an organism, similar mutational features may also be found in other cancers. This analysis yielded straightforward and testable predictions on an accumulated and preferred mutation spectrum in normal tissue. The validation of predicted cancer state mutation patterns demonstrates the usefulness and potential of a causal dynamical framework to understand and predict genetic mutations in cancer. We also obtained the following implication related to HCC therapy, (1) specific positive feedback loops are responsible for the maintenance of normal liver and HCC; (2) inhibiting proliferation and inflammation-related positive feedback loops, and simultaneously inducing liver-specific positive feedback loop is predicated as the potential strategy to cure or relieve HCC; (3) the genesis and regression of HCC is asymmetric. In light of the characteristic property of the nonlinear dynamical system, we demonstrate that positive feedback loops must be existed as a simple and general molecular basis for the maintenance of phenotypes such as normal liver and HCC, and regulating the positive feedback loops directly or indirectly provides potential strategies to cure or relieve HCC.

Decoding early myelopoiesis from dynamics of core endogenous network

A decade ago mainstream molecular biologists regarded it impossible or biologically ill-motivated to understand the dynamics of complex biological phenomena, such as cancer genesis and progression, from a network perspective. Indeed, there are numerical difficulties even for those who were determined to explore along this direction. Undeterred, seven years ago a group of Chinese scientists started a program aiming to obtain quantitative connections between tumors and network dynamics. Many interesting results have been obtained. In this paper we wish to test such idea from a different angle: the connection between a normal biological process and the network dynamics. We have taken early myelopoiesis as our biological model. A standard roadmap for the cell-fate diversification during hematopoiesis has already been well established experimentally, yet little was known for its underpinning dynamical mechanisms. Compounding this difficulty there were additional experimental challenges, such as the seemingly conflicting hematopoietic roadmaps and the cell-fate inter-conversion events. With early myeloid cell-fate determination in mind, we constructed a core molecular endogenous network from well-documented gene regulation and signal transduction knowledge. Turning the network into a set of dynamical equations, we found computationally several structurally robust states. Those states nicely correspond to known cell phenotypes. We also found the states connecting those stable states. They reveal the developmental routes-how one stable state would most likely turn into another stable state. Such interconnected network among stable states enabled a natural organization of cell-fates into a multi-stable state landscape. Accordingly, both the myeloid cell phenotypes and the standard roadmap were explained mechanistically in a straightforward manner. Furthermore, recent challenging observations were also explained naturally. Moreover, the landscape visually enables a prediction of a pool of additional cell states and developmental routes, including the non-sequential and cross-branch transitions, which are testable by future experiments. In summary, the endogenous network dynamics provide an integrated quantitative framework to understand the heterogeneity and lineage commitment in myeloid progenitors.

Challenges in structural approaches to cell modeling

Computational modeling is essential for structural characterization of biomolecular mechanisms across the broad spectrum of scales. Adequate understanding of biomolecular mechanisms inherently involves our ability to model them. Structural modeling of individual biomolecules and their interactions has been rapidly progressing. However, in terms of the broader picture, the focus is shifting toward larger systems, up to the level of a cell. Such modeling involves a more dynamic and realistic representation of the interactomes in vivo, in a crowded cellular environment, as well as membranes and membrane proteins, and other cellular components. Structural modeling of a cell complements computational approaches to cellular mechanisms based on differential equations, graph models, and other techniques to model biological networks, imaging data, etc. Structural modeling along with other computational and experimental approaches will provide a fundamental understanding of life at the molecular level and lead to important applications to biology and medicine. A cross section of diverse approaches presented in this review illustrates the developing shift from the structural modeling of individual molecules to that of cell biology. Studies in several related areas are covered: biological networks; automated construction of three-dimensional cell models using experimental data; modeling of protein complexes; prediction of non-specific and transient protein interactions; thermodynamic and kinetic effects of crowding; cellular membrane modeling; and modeling of chromosomes. The review presents an expert opinion on the current state-of-the-art in these various aspects of structural modeling in cellular biology, and the prospects of future developments in this emerging field.

ACCURATE CHEMICAL MASTER EQUATION SOLUTION USING MULTI-FINITE BUFFERS

The discrete chemical master equation (dCME) provides a fundamental framework for studying stochasticity in mesoscopic networks. Because of the multi-scale nature of many networks where reaction rates have large disparity, directly solving dCMEs is intractable due to the exploding size of the state space. It is important to truncate the state space effectively with quantified errors, so accurate solutions can be computed. It is also important to know if all major probabilistic peaks have been computed. Here we introduce the Accurate CME (ACME) algorithm for obtaining direct solutions to dCMEs. With multi-finite buffers for reducing the state space by O(n!), exact steady-state and time-evolving network probability landscapes can be computed. We further describe a theoretical framework of aggregating microstates into a smaller number of macrostates by decomposing a network into independent aggregated birth and death processes, and give an a priori method for rapidly determining steady-state truncation errors. The maximal sizes of the finite buffers for a given error tolerance can also be pre-computed without costly trial solutions of dCMEs. We show exactly computed probability landscapes of three multi-scale networks, namely, a 6-node toggle switch, 11-node phage-lambda epigenetic circuit, and 16-node MAPK cascade network, the latter two with no known solutions. We also show how probabilities of rare events can be computed from first-passage times, another class of unsolved problems challenging for simulation-based techniques due to large separations in time scales. Overall, the ACME method enables accurate and efficient solutions of the dCME for a large class of networks.

State Space Truncation with Quantified Errors for Accurate Solutions to Discrete Chemical Master Equation

The discrete chemical master equation (dCME) provides a general framework for studying stochasticity in mesoscopic reaction networks. Since its direct solution rapidly becomes intractable due to the increasing size of the state space, truncation of the state space is necessary for solving most dCMEs. It is therefore important to assess the consequences of state space truncations so errors can be quantified and minimized. Here we describe a novel method for state space truncation. By partitioning a reaction network into multiple molecular equivalence groups (MEGs), we truncate the state space by limiting the total molecular copy numbers in each MEG. We further describe a theoretical framework for analysis of the truncation error in the steady-state probability landscape using reflecting boundaries. By aggregating the state space based on the usage of a MEG and constructing an aggregated Markov process, we show that the truncation error of a MEG can be asymptotically bounded by the probability of states on the reflecting boundary of the MEG. Furthermore, truncating states of an arbitrary MEG will not undermine the estimated error of truncating any other MEGs. We then provide an overall error estimate for networks with multiple MEGs. To rapidly determine the appropriate size of an arbitrary MEG, we also introduce an a priori method to estimate the upper bound of its truncation error. This a priori estimate can be rapidly computed from reaction rates of the network, without the need of costly trial solutions of the dCME. As examples, we show results of applying our methods to the four stochastic networks of (1) the birth and death model, (2) the single gene expression model, (3) the genetic toggle switch model, and (4) the phage lambda bistable epigenetic switch model. We demonstrate how truncation errors and steady-state probability landscapes can be computed using different sizes of the MEG(s) and how the results validate our theories. Overall, the novel state space truncation and error analysis methods developed here can be used to ensure accurate direct solutions to the dCME for a large number of stochastic networks.

Not just a colourful metaphor: modelling the landscape of cellular development using Hopfield networks

The epigenetic landscape was introduced by Conrad Waddington as a metaphor of cellular development. Like a ball rolling down a hillside is channelled through a succession of valleys until it reaches the bottom, cells follow specific trajectories from a pluripotent state to a committed state. Transcription factors (TFs) interacting as a network (the gene regulatory network (GRN)) orchestrate this developmental process within each cell. Here, we quantitatively model the epigenetic landscape using a kind of artificial neural network called the Hopfield network (HN). An HN is composed of nodes (genes/TFs) and weighted undirected edges, resulting in a weight matrix (W) that stores interactions among the nodes over the entire network. We used gene co-expression to compute the edge weights. Through W, we then associate an energy score (E) to each input pattern (pattern of co-expression for a specific developmental stage) such that each pattern has a specific E. We propose that, based on the co-expression values stored in W, HN associates lower E values to stable phenotypic states and higher E to transient states. We validate our model using time course gene-expression data sets representing stages of development across 12 biological processes including differentiation of human embryonic stem cells into specialized cells, differentiation of THP1 monocytes to macrophages during immune response and trans-differentiation of epithelial to mesenchymal cells in cancer. We observe that transient states have higher energy than the stable phenotypic states, yielding an arc-shaped trajectory. This relationship was confirmed by perturbation analysis. HNs offer an attractive framework for quantitative modelling of cell differentiation (as a landscape) from empirical data. Using HNs, we identify genes and TFs that drive cell-fate transitions, and gain insight into the global dynamics of GRNs.

Biological Sources of Intrinsic and Extrinsic Noise in cI Expression of Lysogenic Phage Lambda

Genetically identical cells exposed to homogeneous environment can show remarkable phenotypic difference. To predict how phenotype is shaped, understanding of how each factor contributes is required. During gene expression processes, noise could arise either intrinsically in biochemical processes of gene expression or extrinsically from other cellular processes such as cell growth. In this work, important noise sources in gene expression of phage λ lysogen are quantified using models described by stochastic differential equations (SDEs). Results show that DNA looping has sophisticated impacts on gene expression noise: When DNA looping provides autorepression, like in wild type, it reduces noise in the system; When the autorepression is defected as it is in certain mutants, DNA looping increases expression noise. We also study how each gene operator affects the expression noise by changing the binding affinity between the gene and the transcription factor systematically. We find that the system shows extraordinarily large noise when the binding affinity is in certain range, which changes the system from monostable to bistable. In addition, we find that cell growth causes non-negligible noise, which increases with gene expression level. Quantification of noise and identification of new noise sources will provide deeper understanding on how stochasticity impacts phenotype.

Mathematical models of the transitions between endocrine therapy responsive and resistant states in breast cancer

Endocrine therapy, targeting the oestrogen receptor pathway, is the most common treatment for oestrogen receptor-positive breast cancers. Unfortunately, these tumours frequently develop resistance to endocrine therapies. Among the strategies to treat resistant tumours are sequential treatment (in which second-line drugs are used to gain additional responses) and intermittent treatment (in which a 'drug holiday' is imposed between treatments). To gain a more rigorous understanding of the mechanisms underlying these strategies, we present a mathematical model that captures the transitions among three different, experimentally observed, oestrogen-sensitivity phenotypes in breast cancer (sensitive, hypersensitive and independent). To provide a global view of the transitions between these phenotypes, we compute the potential landscape associated with the model. We show how this oestrogen response landscape can be reshaped by population selection, which is a crucial force in promoting acquired resistance. Techniques from statistical physics are used to create a population-level state-transition model from the cellular-level model. We then illustrate how this population-level model can be used to analyse and optimize sequential and intermittent oestrogen-deprivation protocols for breast cancer. The approach used in this study is general and can also be applied to investigate treatment strategies for other types of cancer.

Wright-Fisher dynamics on adaptive landscape

Adaptive landscape, proposed by Sewall Wright, has provided a conceptual framework to describe dynamical behaviours. However, it is still a challenge to explicitly construct such a landscape, and apply it to quantify interesting evolutionary processes. This is particularly true for neutral evolution. In this work, the authors study one-dimensional Wright Fisher process, and analytically obtain an adaptive landscape as a potential function. They provide the complete characterisation for dynamical behaviours of all possible mutation rates under the influence of mutation and random drift. This same analysis has been applied to situations with additive selection and random drift for all possible selection rates. The critical state dividing the basins of two stable states is directly obtained by the landscape. In addition, the landscape is able to handle situations with pure random drift, which would be non-normalisable for its stationary distribution. The nature of non-normalisation is from the singularity of adaptive landscape. In addition, they propose a new type of neutral evolution. It has the same probability for all possible states. The new type of neutral evolution describes the non-neutral alleles with 0%. They take the equal effect of mutation and random drift as an example.

Absorbing phenomena and escaping time for Muller's ratchet in adaptive landscape

BACKGROUND

The accumulation of deleterious mutations of a population directly contributes to the fate as to how long the population would exist, a process often described as Muller's ratchet with the absorbing phenomenon. The key to understand this absorbing phenomenon is to characterize the decaying time of the fittest class of the population. Adaptive landscape introduced by Wright, a re-emerging powerful concept in systems biology, is used as a tool to describe biological processes. To our knowledge, the dynamical behaviors for Muller's ratchet over the full parameter regimes are not studied from the point of the adaptive landscape. And the characterization of the absorbing phenomenon is not yet quantitatively obtained without extraneous assumptions as well.

METHODS

We describe how Muller's ratchet can be mapped to the classical Wright-Fisher process in both discrete and continuous manners. Furthermore, we construct the adaptive landscape for the system analytically from the general diffusion equation. The constructed adaptive landscape is independent of the existence and normalization of the stationary distribution. We derive the formula of the single click time in finite and infinite potential barrier for all parameters regimes by mean first passage time.

RESULTS

We describe the dynamical behavior of the population exposed to Muller's ratchet in all parameters regimes by adaptive landscape. The adaptive landscape has rich structures such as finite and infinite potential, real and imaginary fixed points. We give the formula about the single click time with finite and infinite potential. And we find the single click time increases with selection rates and population size increasing, decreases with mutation rates increasing. These results provide a new understanding of infinite potential. We analytically demonstrate the adaptive and unadaptive states for the whole parameters regimes. Interesting issues about the parameters regions with the imaginary fixed points is demonstrated. Most importantly, we find that the absorbing phenomenon is characterized by the adaptive landscape and the single click time without any extraneous assumptions. These results suggest a graphical and quantitative framework to study the absorbing phenomenon.

Perturbation-based analysis and modeling of combinatorial regulation in the yeast sulfur assimilation pathway

Here we establish the utility of a recently described perturbative method to study complex regulatory circuits in vivo. By combining rapid modulation of single TFs under physiological conditions with genome-wide expression analysis, we elucidate several novel regulatory features within the pathways of sulfur assimilation and beyond.

In yeast, the pathways of sulfur assimilation are combinatorially controlled by five transcriptional regulators (three DNA-binding proteins [Met31p, Met32p, and Cbf1p], an activator [Met4p], and a cofactor [Met28p]) and a ubiquitin ligase subunit (Met30p). This regulatory system exerts combinatorial control not only over sulfur assimilation and methionine biosynthesis, but also on many other physiological functions in the cell. Recently we characterized a gene induction system that, upon the addition of an inducer, results in near-immediate transcription of a gene of interest under physiological conditions. We used this to perturb levels of single transcription factors during steady-state growth in chemostats, which facilitated distinction of direct from indirect effects of individual factors dynamically through quantification of the subsequent changes in genome-wide patterns of gene expression. We were able to show directly that Cbf1p acts sometimes as a repressor and sometimes as an activator. We also found circumstances in which Met31p/Met32p function as repressors, as well as those in which they function as activators. We elucidated and numerically modeled feedback relationships among the regulators, notably feedforward regulation of Met32p (but not Met31p) by Met4p that generates dynamic differences in abundance that can account for the differences in function of these two proteins despite their identical binding sites.

Evolutionary mechanisms and diversity in cancer

The recently introduced genome theory of cancer evolution provides a new framework for evolutionary studies on cancer. In particular, the established relationship between the large number of individual molecular mechanisms and the general evolutionary mechanism of cancer calls upon a change in our strategies that have been based on the characterization of common cancer gene mutations and their defined pathways. To further explain the significance of the genome theory of cancer evolution, a brief review will be presented describing the various attempts to illustrate the evolutionary mechanism of cancer, followed by further analysis of some key components of somatic cell evolution, including the diversity of biological systems, the multiple levels of information systems and control systems, the two phases (the punctuated or discontinuous phase and gradual Darwinian stepwise phase) and dynamic patterns of somatic cell evolution where genome replacement is the driving force. By linking various individual molecular mechanisms to the level of genome population diversity and tumorigenicity, the general mechanism of cancer has been identified as the evolutionary mechanism of cancer, which can be summarized by the following three steps including stress-induced genome instability, population diversity or heterogeneity, and genome-mediated macroevolution. Interestingly, the evolutionary mechanism is equal to the collective aggregate of all individual molecular mechanisms. This relationship explains why most of the known molecular mechanisms can contribute to cancer yet there is no single dominant mechanism for the majority of clinical cases. Despite the fact that each molecular mechanism can serve as a system stress and initiate the evolutionary process, to achieve cancer, multiple cycles of genome-mediated macroevolution are required and are a stochastically determined event. Finally, the potential clinical implications of the evolutionary mechanism of cancer are briefly reviewed.

A deterministic map of Waddington's epigenetic landscape for cell fate specification

Background

The image of the "epigenetic landscape", with a series of branching valleys and ridges depicting stable cellular states and the barriers between those states, has been a popular visual metaphor for cell lineage specification - especially in light of the recent discovery that terminally differentiated adult cells can be reprogrammed into pluripotent stem cells or into alternative cell lineages. However the question of whether the epigenetic landscape can be mapped out quantitatively to provide a predictive model of cellular differentiation remains largely unanswered.

Results

Here we derive a simple deterministic path-integral quasi-potential, based on the kinetic parameters of a gene network regulating cell fate, and show that this quantity is minimized along a temporal trajectory in the state space of the gene network, thus providing a marker of directionality for cell differentiation processes. We then use the derived quasi-potential as a measure of "elevation" to quantitatively map the epigenetic landscape, on which trajectories flow "downhill" from any location. Stochastic simulations confirm that the elevation of this computed landscape correlates to the likelihood of occurrence of particular cell fates, with well-populated low-lying "valleys" representing stable cellular states and higher "ridges" acting as barriers to transitions between the stable states.

Conclusions

This quantitative map of the epigenetic landscape underlying cell fate choice provides mechanistic insights into the "forces" that direct cellular differentiation in the context of physiological development, as well as during artificially induced cell lineage reprogramming. Our generalized approach to mapping the landscape is applicable to non-gradient gene regulatory systems for which an analytical potential function cannot be derived, and also to high-dimensional gene networks. Rigorous quantification of the gene regulatory circuits that govern cell lineage choice and subsequent mapping of the epigenetic landscape can potentially help identify optimal routes of cell fate reprogramming.

Noise regulation by quorum sensing in low mRNA copy number systems

Background

Cells must face the ubiquitous presence of noise at the level of signaling molecules. The latter constitutes a major challenge for the regulation of cellular functions including communication processes. In the context of prokaryotic communication, the so-called quorum sensing (QS) mechanism relies on small diffusive molecules that are produced and detected by cells. This poses the intriguing question of how bacteria cope with the fluctuations for setting up a reliable information exchange.

Results

We present a stochastic model of gene expression that accounts for the main biochemical processes that describe the QS mechanism close to its activation threshold. Within that framework we study, both numerically and analytically, the role that diffusion plays in the regulation of the dynamics and the fluctuations of signaling molecules. In addition, we unveil the contribution of different sources of noise, intrinsic and transcriptional, in the QS mechanism.

Conclusions

The interplay between noisy sources and the communication process produces a repertoire of dynamics that depends on the diffusion rate. Importantly, the total noise shows a non-monotonic behavior as a function of the diffusion rate. QS systems seems to avoid values of the diffusion that maximize the total noise. These results point towards the direction that bacteria have adapted their communication mechanisms in order to improve the signal-to-noise ratio.

Probability landscape of heritable and robust epigenetic state of lysogeny in phage lambda

Computational studies of biological networks can help to identify components and wirings responsible for observed phenotypes. However, studying stochastic networks controlling many biological processes is challenging. Similar to Schrödinger's equation in quantum mechanics, the chemical master equation (CME) provides a basic framework for understanding stochastic networks. However, except for simple problems, the CME cannot be solved analytically. Here we use a method called discrete chemical master equation (dCME) to compute directly the full steady-state probability landscape of the lysogeny maintenance network in phage lambda from its CME. Results show that wild-type phage lambda can maintain a constant level of repressor over a wide range of repressor degradation rate and is stable against UV irradiation, ensuring heritability of the lysogenic state. Furthermore, it can switch efficiently to the lytic state once repressor degradation increases past a high threshold by a small amount. We find that beyond bistability and nonlinear dimerization, cooperativity between repressors bound to O(R)1 and O(R)2 is required for stable and heritable epigenetic state of lysogeny that can switch efficiently. Mutants of phage lambda lack stability and do not possess a high threshold. Instead, they are leaky and respond to gradual changes in degradation rate. Our computation faithfully reproduces the hair triggers for UV-induced lysis observed in mutants and the limitation in robustness against mutations. The landscape approach computed from dCME is general and can be applied to study broad issues in systems biology.

Computational Cellular Dynamics Based on the Chemical Master Equation: A Challenge for Understanding Complexity

Modern molecular biology has always been a great source of inspiration for computational science. Half a century ago, the challenge from understanding macromolecular dynamics has led the way for computations to be part of the tool set to study molecular biology. Twenty-five years ago, the demand from genome science has inspired an entire generation of computer scientists with an interest in discrete mathematics to join the field that is now called bioinformatics. In this paper, we shall lay out a new mathematical theory for dynamics of biochemical reaction systems in a small volume (i.e., mesoscopic) in terms of a stochastic, discrete-state continuous-time formulation, called the chemical master equation (CME). Similar to the wavefunction in quantum mechanics, the dynamically changing probability landscape associated with the state space provides a fundamental characterization of the biochemical reaction system. The stochastic trajectories of the dynamics are best known through the simulations using the Gillespie algorithm. In contrast to the Metropolis algorithm, this Monte Carlo sampling technique does not follow a process with detailed balance. We shall show several examples how CMEs are used to model cellular biochemical systems. We shall also illustrate the computational challenges involved: multiscale phenomena, the interplay between stochasticity and nonlinearity, and how macroscopic determinism arises from mesoscopic dynamics. We point out recent advances in computing solutions to the CME, including exact solution of the steady state landscape and stochastic differential equations that offer alternatives to the Gilespie algorithm. We argue that the CME is an ideal system from which one can learn to understand "complex behavior" and complexity theory, and from which important biological insight can be gained.

Stochastic probability landscape model for switching efficiency, robustness, and differential threshold for induction of genetic circuit in phage lambda

The genetic switch of phage lambda provides a paradigm for studying developmental biology and cell fate. Although there have been numerous experimental and theoretical studies, the mechanisms for switching efficiency, stability and robustness, maintenance of lysogenic state, and the induction lytic state are not fully understood. In this paper, a new method is adopted to account for the full stochasticity typical in small copy events and the exact steady state probability landscape of the genetic circuit of the switch network is computed. The stability and sensitivity of phage lambda switch and its robustness against perturbations derived from our model for wild type and mutants are consistent with experimental data. Our study has revealed likely mechanism of the phage switch network that was previously unknown.

Global view of bionetwork dynamics: adaptive landscape

Based on recent work, I will give a nontechnical brief review of a powerful quantitative concept in biology, adaptive landscape, initially proposed by S. Wright over 70 years ago, reintroduced by one of the founders of molecular biology and by others in different biological contexts, but apparently forgotten by modern biologists for many years. Nevertheless, this concept finds an increasingly important role in the development of systems biology and bionetwork dynamics modeling, from phage lambda genetic switch to endogenous network for cancer genesis and progression. It is an ideal quantification to describe the robustness and stability of bionetworks. Here, I will first introduce five landmark proposals in biology on this concept, to demonstrate an important common thread in theoretical biology. Then I will discuss a few recent results, focusing on the studies showing theoretical consistency of adaptive landscape. From the perspective of a working scientist and of what is needed logically for a dynamical theory when confronting empirical data, the adaptive landscape is useful both metaphorically and quantitatively, and has captured an essential aspect of biological dynamical processes. Though at the theoretical level the adaptive landscape must exist and it can be used across hierarchical boundaries in biology, many associated issues are indeed vague in their initial formulations and their quantitative realizations are not easy, and are good research topics for quantitative biologists. I will discuss three types of open problems associated with the adaptive landscape in a broader perspective.

Optimal enumeration of state space of finitely buffered stochastic molecular networks and exact computation of steady state landscape probability

BACKGROUND

Stochasticity plays important roles in many molecular networks when molecular concentrations are in the range of 0.1 muM to 10nM (about 100 to 10 copies in a cell). The chemical master equation provides a fundamental framework for studying these networks, and the time-varying landscape probability distribution over the full microstates, i.e., the combination of copy numbers of molecular species, provide a full characterization of the network dynamics. A complete characterization of the space of the microstates is a prerequisite for obtaining the full landscape probability distribution of a network. However, there are neither closed-form solutions nor algorithms fully describing all microstates for a given molecular network.

RESULTS

We have developed an algorithm that can exhaustively enumerate the microstates of a molecular network of small copy numbers under the condition that the net gain in newly synthesized molecules is smaller than a predefined limit. We also describe a simple method for computing the exact mean or steady state landscape probability distribution over microstates. We show how the full landscape probability for the gene networks of the self-regulating gene and the toggle-switch in the steady state can be fully characterized. We also give an example using the MAPK cascade network. Data and server will be available at URL: http://scsb.sjtu.edu.cn/statespace.

CONCLUSION

Our algorithm works for networks of small copy numbers buffered with a finite copy number of net molecules that can be synthesized, regardless of the reaction stoichiometry, and is optimal in both storage and time complexity. The algorithm can also be used to calculate the rates of all transitions between microstates from given reactions and reaction rates. The buffer size is limited by the available memory or disk storage. Our algorithm is applicable to a class of biological networks when the copy numbers of molecules are small and the network is closed, or the network is open but the net gain in newly synthesized molecules does not exceed a predefined buffer capacity. For these networks, our method allows full stochastic characterization of the mean landscape probability distribution, and the steady state when it exists.

Reaction coordinates for the flipping of genetic switches

We present a detailed analysis, based on the forward flux sampling simulation method, of the switching dynamics and stability of two models of genetic toggle switches, consisting of two mutually repressing genes encoding transcription factors (TFs); in one model (the exclusive switch), the two transcription factors mutually exclude each other's binding, while in the other model (general switch), the two TFs can bind simultaneously to the shared operator region. We assess the role of two pairs of reactions that influence the stability of these switches: TF-TF homodimerization and TF-DNA association/dissociation. In both cases, the switch flipping rate increases with the rate of TF dimerization, while it decreases with the rate of TF-operator binding. We factorize the flipping rate k into the product of the probability rho(q*) of finding the system at the dividing surface (separatrix) between the two stable states, and a kinetic prefactor R. In the case of the exclusive switch, the rate of TF-operator binding affects both rho(q*) and R, while the rate of TF dimerization affects only R. The general switch displays a higher flipping rate than the exclusive switch, and both TF-operator binding and TF dimerization affect k, R, and rho(q*). To elucidate this, we analyze the transition state ensemble. For the exclusive switch, the transition state ensemble is strongly affected by the rate of TF-operator binding, but unaffected by varying the rate of TF-TF binding. Thus, varying the rate of TF-operator binding can drastically change the pathway of switching, while changing the rate of dimerization changes the switching rate without altering the mechanism. The switching pathways of the general switch are highly robust to changes in the rate constants of both TF-operator and TF-TF binding, even though these rate constants do affect the flipping rate; this feature is unique for nonequilibrium systems.

Dynamics and Intrinsic Statistical Fluctuations of a Gene Switch

We studied a simple gene regulatory network, the toggle switch. Specifically, we examined the means and statistical fluctuations in numbers of proteins. We found that when omega, the ratio of rates of protein-gene unbinding to protein degradation, was between approximately 10(-3) and approximately 10, the fluctuations were much larger than those we would have expected from Poisson statistics. In addition, we examined characteristic time values for system relaxation and found both that they increased with omega and that they have significant phase transition effects, with a secondary time scale appearing near the boundary between bistable and other phases. Last, we discuss the bistability of the toggle switch.

Nonlinear regulation enhances the phenotypic expression of trans-acting genetic polymorphisms

BACKGROUND

Genetic variation explains a considerable part of observed phenotypic variation in gene expression networks. This variation has been shown to be located both locally (cis) and distally (trans) to the genes being measured. Here we explore to which degree the phenotypic manifestation of local and distant polymorphisms is a dynamic feature of regulatory design.

RESULTS

By combining mathematical models of gene expression networks with genetic maps and linkage analysis we find that very different network structures and regulatory motifs give similar cis/trans linkage patterns. However, when the shape of the cis-regulatory input functions is more nonlinear or threshold-like, we observe for all networks a dramatic increase in the phenotypic expression of distant compared to local polymorphisms under otherwise equal conditions.

CONCLUSION

Our findings indicate that genetic variation affecting the form of cis-regulatory input functions may reshape the genotype-phenotype map by changing the relative importance of cis and trans variation. Our approach combining nonlinear dynamic models with statistical genetics opens up for a systematic investigation of how functional genetic variation is translated into phenotypic variation under various systemic conditions.

Accelerated search for biomolecular network models to interpret high-throughput experimental data

BACKGROUND

The functions of human cells are carried out by biomolecular networks, which include proteins, genes, and regulatory sites within DNA that encode and control protein expression. Models of biomolecular network structure and dynamics can be inferred from high-throughput measurements of gene and protein expression. We build on our previously developed fuzzy logic method for bridging quantitative and qualitative biological data to address the challenges of noisy, low resolution high-throughput measurements, i.e., from gene expression microarrays. We employ an evolutionary search algorithm to accelerate the search for hypothetical fuzzy biomolecular network models consistent with a biological data set. We also develop a method to estimate the probability of a potential network model fitting a set of data by chance. The resulting metric provides an estimate of both model quality and dataset quality, identifying data that are too noisy to identify meaningful correlations between the measured variables.

RESULTS

Optimal parameters for the evolutionary search were identified based on artificial data, and the algorithm showed scalable and consistent performance for as many as 150 variables. The method was tested on previously published human cell cycle gene expression microarray data sets. The evolutionary search method was found to converge to the results of exhaustive search. The randomized evolutionary search was able to converge on a set of similar best-fitting network models on different training data sets after 30 generations running 30 models per generation. Consistent results were found regardless of which of the published data sets were used to train or verify the quantitative predictions of the best-fitting models for cell cycle gene dynamics.

CONCLUSION

Our results demonstrate the capability of scalable evolutionary search for fuzzy network models to address the problem of inferring models based on complex, noisy biomolecular data sets. This approach yields multiple alternative models that are consistent with the data, yielding a constrained set of hypotheses that can be used to optimally design subsequent experiments.

Single molecule dynamics and statistical fluctuations of gene regulatory networks: A repressilator

The authors developed a time dependent method to study the single molecule dynamics of a simple gene regulatory network: a repressilator with three genes mutually repressing each other. They quantitatively characterize the time evolution dynamics of the repressilator. Furthermore, they study purely dynamical issues such as statistical fluctuations and noise evolution. They illustrated some important features of the biological network such as monostability, spirals, and limit cycle oscillation. Explicit time dependent Fano factors which describe noise evolution and show statistical fluctuations out of equilibrium can be significant and far from the Poisson distribution. They explore the phase space and the interrelationships among fluctuations, order, amplitude, and period of oscillations of the repressilators. The authors found that repressilators follow ordered limit cycle orbits and are more likely to appear in the lower fluctuating regions. The amplitude of the repressilators increases as the suppressing of the genes decreases and production of proteins increases. The oscillation period of the repressilators decreases as the suppressing of the genes decreases and production of proteins increases.

Highly designable phenotypes and mutational buffers emerge from a systematic mapping between network topology and dynamic output

Deciphering the design principles for regulatory networks is fundamental to an understanding of biological systems. We have explored the mapping from the space of network topologies to the space of dynamical phenotypes for small networks. Using exhaustive enumeration of a simple model of three- and four-node networks, we demonstrate that certain dynamical phenotypes can be generated by an atypically broad spectrum of network topologies. Such dynamical outputs are highly designable, much like certain protein structures can be designed by an unusually broad spectrum of sequences. The network topologies that encode a highly designable dynamical phenotype possess two classes of connections: a fully conserved core of dedicated connections that encodes the stable dynamical phenotype and a partially conserved set of variable connections that controls the transient dynamical flow. By comparing the topologies and dynamics of the three- and four-node network ensembles, we observe a large number of instances of the phenomenon of "mutational buffering," whereby addition of a fourth node suppresses phenotypic variation amongst a set of three-node networks.

Cycle kinetics, steady state thermodynamics and motors-a paradigm for living matter physics

An integration of the stochastic mathematical models for motor proteins with Hill's steady state thermodynamics yields a rather comprehensive theory for molecular motors as open systems in the nonequilibrium steady state. This theory, a natural extension of Gibbs' approach to isothermal molecular systems in equilibrium, is compared with other existing theories with dissipative structures and dynamics. The theory of molecular motors might be considered as an archetype for studying more complex open biological systems such as biochemical reaction networks inside living cells.

Nonlinear protein degradation and the function of genetic circuits

The functions of most genetic circuits require a sufficient degree of cooperativity in the circuit components. Although mechanisms of cooperativity have been studied most extensively in the context of transcriptional initiation control, cooperativity from other processes involved in the operation of the circuits can also play important roles. In this work, we examine a simple kinetic source of cooperativity stemming from the nonlinear degradation of multimeric proteins. Ample experimental evidence suggests that protein subunits can degrade less rapidly when associated in multimeric complexes, an effect we refer to as "cooperative stability." For dimeric transcription factors, this effect leads to a concentration-dependence in the degradation rate because monomers, which are predominant at low concentrations, will be more rapidly degraded. Thus, cooperative stability can effectively widen the accessible range of protein levels in vivo. Through theoretical analysis of two exemplary genetic circuits in bacteria, we show that such an increased range is important for the robust operation of genetic circuits as well as their evolvability. Our calculations demonstrate that a few-fold difference between the degradation rate of monomers and dimers can already enhance the function of these circuits substantially. We discuss molecular mechanisms of cooperative stability and their occurrence in natural or engineered systems. Our results suggest that cooperative stability needs to be considered explicitly and characterized quantitatively in any systematic experimental or theoretical study of gene circuits.

Noise in a Small Genetic Circuit that Undergoes Bifurcation

Based on the consideration of Boolean dynamics, it has been hypothesized that cell types may correspond to alternative attractors of a gene regulatory network. Recent stochastic Boolean network analysis, however, raised the important question concerning the stability of such attractors. In this paper a detailed numerical analysis is performed within the framework of Langevin dynamics. While the present results confirm that the noise is indeed an important dynamical element, the cell type as represented by attractors can still be a viable hypothesis. It is found that the stability of an attractor depends on the strength of noise related to the distance of the system to the bifurcation point and it can be exponentially stable depending on biological parameters.