Workshop 6: Algebraic Methods in Systems and Evolutionary Biology: Titles & Abstracts
We compare continuous and Boolean dynamics of gene regulatory networks, using simple modules as well as complex networks. We translate the Boolean model for gene activity into continuous dynamics for mRNA and protein concentrations using sigmoidal Hill functions. We establish conditions under which dynamical attractors of the Boolean and the continuous model agree with each other. These conditions depend only on general features such as the ratio of the relevant time scales, the degree of cooperativity of the regulating interactions, the logical structure of the interactions, and the robustness of the dynamical patterns under random time delay.
Over the past five years my group, in collaboration with wet-bench biologists, developed and validated asynchronous Boolean models of several signal transduction networks. Along the way we have encountered obstacles related to the lack of timing knowledge and the large size of the state space. In this talk I will present three methodologies we developed to overcome these obstacles. First, from a comparative analysis of several asynchronous update methods we concluded that updating a single, randomly selected node at each time instant offers the best combination of information and economy.
Second, we developed a two-step network reduction method which was able to reduce the number of variables by 90% in two different systems without affecting their dynamic behaviors. Third, we proposed an integration of Boolean rules into graph theoretical analysis and showed that this semi-structural method can identify critical signal mediators on par with dynamic models.
Gene regulatory networks that operate in the dynamically critical regime (between order and chaos) are optimised with respect to the trade-off between phenotypic robustness and flexibility - a balance that ensures both homeostasis and development. In fact analyses of the gene network architecture and patterns of transcriptome changes in several organisms in the past few years suggest that the gene regulatory networks of living organisms are indeed in the critical regime. But how does a gene regulatory network evolve a structure that confers criticality? While this question has evaded scientists for decades, a related equally fundamental question has over the past years attracted considerable interest: The evolution of evolvability. There is now the consensus that evolvability itself is a selectable trait. Evolvability, similar to criticality, is associated with the trade-off between mutational robustness on the one hand (mutations should not disrupt essential functions) and innovation on the other hand (mutations should alter networks sufficiently to add new functions). In this work I will show that critical dynamics in genetic network models naturally emerge as a robust byproduct of the very same evolutionary processes that select for evolvability -- without fine-tuning of parameters or imposing explicit selection criteria (i.e. arbitrary fitness functions). More specificaly, criticality emerges from the requirement of evolvability in the sense that during evolution, the existing adaptive phenotypes must be preserved while allowing new phenotypes to emerge for the organism to be able to cope with new environmental challenges. Strikingly, the gene networks produced by selecting for evolvability have a structure (topology) that is very similar to the one observed in real organisms, such as Escherichia col, charecterized by the existence of global regulators.
Transcriptional gene regulation is a dynamic process and its proper functioning is essential for all living organisms. By combining the abundant static regulatory data with time series expression data using an Input-Output Hidden Markov model (IOHMM) we were able to reconstruct a dynamic representations for these networks in multiple species. The models lead to testable temporal hypotheses identifying both new regulators and their time of activation. We have recently extended these models, by solving an optimization problem related to graph orientation, to connect signaling and regulatory networks. These reconstructed networks link receptors and proteins that directly interact with the environment to the observed expression outcome. I will discuss the application and experimental validation of predictions made by our methods focusing on stress response in yeast. I would also mention a number of other extensions for integrating microRNAs and discriminative motif search and their applications for studying mice development and response to pathogens in mammalian cells.
Dynamical system models are very commonly used to analyze biological interaction networks, such as the dynamics of concentrations in biochemical reaction networks, the spread of infectious diseases within a population, and the dynamics of species in an ecosystem.
Persistence and permanence are properties of dynamical systems that provide information about the long-term behavior of the system. For example, the persistence property is relevant in deciding if, in the long term, a chemical species will be completely consumed by a reaction network, an infection will die off, or a species in an ecosystem will become extinct.
We prove that two-species mass-action systems derived from weakly reversible networks are both persistent and permanent, for any values of the reaction rate parameters. Moreover, we prove that a larger class of networks, called endotactic networks, also give rise to persistent systems, even if we allow the reaction rate parameters to vary in time. These results also apply to general polynomial dynamical systems and other nonlinear dynamical systems. This is joint work with Fedor Nazarov and Casian Pantea.
Cancer encompasses various diseases associated with loss of cell-cycle control, leading to uncontrolled cell proliferation and/or reduced apoptosis. Cancer is usually caused by malfunction(s) in the cellular signaling pathways. Malfunctions occur in different ways and at different locations in a pathway. Consequently, therapy design should first identify the location and type of malfunction and then arrive at a suitable drug combination. We consider the growth factor (GF) signaling pathways, widely studied in the context of cancer. Interactions between different pathway components are modeled using Boolean logic gates. All possible single malfunctions in the resulting circuit are enumerated and responses of the different malfunctioning circuits to a 'test' input are used to group the malfunctions into classes. Effects of different drugs, targeting different parts of the Boolean circuit, are taken into account in deciding drug efficacy, thereby mapping each malfunction to an appropriate set of drugs.
Standard operator optimization assumes a mathematical model and a family of operators on the model, defines an objective function to measure operator performance, and selects an optimal operator as one minimizing (maximizing) the objective function. With robust optimization the paradigm changes: the model is unknown and belongs to an uncertainty class of operators. If a prior distribution governs the uncertainty class, then an intrinsically Bayesian robust (IBR) operator is one that provides the minimum expectation of the objective function across the uncertainty class. The objective cost of uncertainty for any model in the uncertainty class is measured by the gain in performance from using an optimal operator for the model as opposed to using an IBR operator on the model. The overall cost of uncertainty is measured by the mean objective cost of uncertainty (MOCU) over the uncertainty class. Optimization is constrained by optimizing over the family of operators that are optimal for specific models in the uncertainty class. In this case, we refer to an optimal operator as a model-constrained Bayesian robust (MCBR) operator. Two general kinds of operational intervention have been considered for gene regulatory networks. Structural intervention involves a one-time alteration of the regulatory apparatus, which means a permanent change to the regulatory logic. Stationary intervention involves a time-invariant external control policy whose action at any time point depends on feedback from the system. In the presence of uncertainty, MCBR control has been treated for the stationary control of probabilistic Boolean networks in the framework of an objective cost of control and dynamic programming. This talk discusses IBR structural intervention for PBNs and the corresponding MOCU as it pertains to uncertain gene logic.
The decade of genomic revolution following the human genome's sequencing has produced significant medical advances, and yet again, revealed how complicated human biology is, and how much more remains to be understood. Biology is an extraordinary complicated puzzle; we may know some of its pieces but have no clue how they are assembled to orchestrate the symphony of life, which renders the comprehension and analysis of living systems a major challenge. Recent efforts to create executable models of complex biological phenomena - an approach we call Executable Biology - entail great promise for new scientific discoveries, shading new light on the puzzle of life. At the same time, this new wave of the future forces computer science to stretch far and beyond, and in ways never considered before, in order to deal with the enormous complexity observed in biology. This talk will focus on our recent success stories in using formal methods to model cell fate decisions during development and cancer, and ongoing efforts to develop dedicated tools for biologists to model cellular processes in a visual-friendly way.
A network of biochemical reactions gives rise, under mass-action kinetics, to a polynomial dynamical system, the steady states of which form a real algebraic variety. Methods from computational algebraic geometry have begun to yield biological insights into such networks. Here, we review work on "invariants". An invariant is a polynomial expression in specified state variables that holds in any (positive) steady state of a network. Invariants characterise the behaviour of several model networks and concisely capture their salient behaviour. However, it remains challenging to calculate invariants and we do not yet understand why some networks have biologically-meaningful invariants while others, apparently, do not.
How cells make decisions can be investigated using mathematical models. We describe a procedure to decide whether a model is compatible with steady-state data. This method requires no parameter values-- it is based on techniques from algebraic geometry, linear algebra, and optimization. Cellular decisions also depend on where they occur in the cell (e.g., nucleus or cytoplasm). We also find that cellular information processing can be altered by including spatial organization. Borrowing tools from chemical reaction network theory and dynamical systems, we show that the existence of distinct compartments plays a pivotal role in whether a system is capable of multistationarity.
Balanced minimum evolution (BME) is a statistically consistent distance-based method to reconstruct a phylogenetic tree from an alignment of molecular data. In 2000, Pauplin showed that the BME method is equivalent to optimizing a linear functional over the BME polytope, the convex hull of the BME vectors obtained from Pauplin's formula applied to all binary trees. The BME method is related to the popular Neighbor Joining (NJ) algorithm, now known to be a greedy optimization of the BME principle. In this talk I will elucidate some of the structure of the BME polytope and strengthen the connection between the BME method and NJ Algorithm. I will show that any subtree-prune-regraft move from a binary tree to another binary tree corresponds to an edge of the BME polytope. Moreover, I will describe an entire family of faces parametrized by disjoint clades. Finally, given a phylogenetic tree T, I will show that the BME cone and every NJ cone of T have intersection of positive measure.
This talk reports on our investigations to explore appropriate modelling and analysis techniques for processes evolving simultaneously over time and space, applied to biological systems. Current challenges for modelling in Systems Biology include those associated with issues of complexity and representing systems with multi-scale attributes. A drawback of current modelling approaches, including Petri nets, is their limitation to relatively small networks. We use Stochastic and Continuous Petri Nets to consider continuous time evolution as Markov process or system of Ordinary Differential Equations, and Coloured Petri Nets to statically encode finite discrete space. Combining both concepts yields Coloured Stochastic and Coloured Continuous Petri nets, which allow for directly executable models as well as computational experiments using standard analysis and simulation techniques over very large networks. We illustrate our approach by a couple of case studies, including gradient formation, multistrain bacterial colonies, and planar cell polarity signalling in Drosophila wing.
Systems biology aims to explain how a biological system functions by investigating the interactions of its individual components from a systems perspective. Modeling is a vital tool as it helps to elucidate the underlying mechanisms of the system. Many discrete model types can be translated into the framework of polynomial dynamical systems (PDS), that is, time- and state-discrete dynamical systems over a finite field where the transition function for each variable is given as a polynomial. This allows for using a range of theoretical and computational tools from computer algebra, which results in a powerful computational engine for model construction, parameter estimation, and analysis methods.
This talk will attempt to provide a synthesis of the topics discussed at the workshop and to distill some central themes that point to opportunities and challenges in the field.
The availability of genomes of many closely related bacteria with diverse metabolic capabilities offers the possibility of tracing metabolic evolution on a phylogeny relating the genomes to understand the evolutionary processes and constraints that affect the evolution of metabolic networks. Using simple (independent loss/gain of reactions) or complex (incorporating dependencies among reactions) stochastic models of metabolic evolution, it is possible to study how metabolic networks evolve over time. Here, we describe metabolic network evolution as a discrete space continuous time Markov process and introduce a neighbor-dependent model for the evolution of metabolic networks where the rates with which reactions are added or removed depend on the fraction of neighboring reactions present in the network. The model also allows estimation of the strength of the neighborhood effect during the course of evolution. We present Gibbs samplers for sampling networks at the internal node of a phylogeny and for estimating the parameters of evolution over a phylogeny without exploring the whole search space by iteratively sampling from the conditional distributions of the internal networks and parameters. The samplers are used to estimate the parameters of evolution of metabolic networks of bacteria in the genus Pseudomonas and to infer the metabolic networks of the ancestral pseudomonads. The results suggest that pathway maps that are conserved across the Pseudomonas phylogeny have a stronger neighborhood structure than those which have a variable distribution of reactions across the phylogeny, and that some Pseudomonas lineages are going through genome reduction resulting in the loss of a number of reactions from their metabolic networks.
The logical method proved useful for the modelling of regulatory networks and the analysis of their dynamical properties. It relies on two directed graphs: the regulatory graph, which represents the interactions between regulatory components, each associated with discrete levels of expression (or activity), and the state transition graph, which represents the discrete dynamics of such a model. This discrete modelling framework allows qualitative analyses of the behaviours driven by regulatory networks, based on analytical results or on simulation (i.e. construction of state transition graphs). Although this formalism abstracts and simplifies the biological reality, we still need to cope with challenging issues due to the complexity of ever increasing networks. We present here some results and tools that aim at facilitating the analysis of large networks. In particular, we will show how dynamical properties can be predicted from the presence of particular motifs in the regulatory graph, namely regulatory circuits and combination of such circuits. We will also discuss a method to reduce the model, yet keeping the main features of the dynamics. Finally, we will illustrate these approaches on a generic model of E2F1-dependent apoptosis and cell cycle entries.
Chemical reaction networks taken with mass-action kinetics are dynamical systems governed by polynomial differential equations that arise in systems biology. In general, establishing the existence of (multiple) steady states is challenging, as it requires the solution of a large system of polynomials with unknown coefficients. If, however, the steady state ideal of the system is a binomial ideal, then we show that these questions can be answered easily. This talk focuses on systems with this property, and we say such systems have toric steady states. Our main result gives sufficient conditions for a chemical reaction system to admit toric steady states. Furthermore, we analyze the capacity of such a system to exhibit multiple steady states. An important application concerns the biochemical reaction networks networks that describe the multisite phosphorylation of a protein by a kinase/phosphatase pair in a sequential and distributive mechanism. No prior knowledge of chemical reaction network theory or binomial ideals will be assumed.
This is joint work with Carsten Conradi, Mercedes Pérez Millán, and Alicia Dickenstein.
Regulatory networks of biomolecular interactions in cells govern virtually all cellular behaviors and functions. Modern measurement technologies are being used to generateinformation on many types of interactions, involving transcriptional and microRNA regulatory networks, signaling networks, and cytokine networks. Temporal measurements of gene and protein expression levels and chromatin modifications, coupled with data fusion strategies that incorporate computational predictions of regulatory mechanisms on the basis of other types of information, such as nucleic acid sequence, can be used to constructdynamical system models of these networks. The analysis and simulation of such models, in conjunction with experimental validation, sheds light on biological function and paves the way toward rational and systematic control strategies intended to drive a diseased system toward a desired state by means of targeted interventions. At the same time, such systems approaches permit new biological observables that reflect system-level behavior that cannot be understood by studying individual sets of interactions. Cellular decision making, maintenance of homeostasis and robustness, sensitivity to diverse types of information in the presence of environmental variability, and coordination ofcomplex macroscopic behavior are examples of such emergent systems-level behavior. Information theoretic approaches combined with elements of dynamical systems theory, such as phase transitions and structure dynamics relationships, are promising frameworks for studying fundamental principles governing living systems at all scales of organization.
Boolean models with continuous-time updates can directly represent timing information that determines the dynamics of a gene regulatory network. We study the relation between autonomous Boolean models and differential equation models, emphasizing the advantages of the Boolean approach for certain system types. We first note several crucial features of generic systems that must be implemented within such models. We then consider the experimentally well-studied example of the gene network that controls fly body segmentation. An autonomous Boolean model successfully recapitulates the patterns formed in normal and genetically perturbed fly embryos and permits the derivation of constraints on the time delay parameters for the formation of experimentally observed patterns. The Boolean model captures the essential features of previously presented differential equation models and clarifies the logic associated with different parameter choices.
Discrete models have been used successfully in modeling biological processes such as gene regulatory networks. When certain regulation mechanisms are unknown it is important to be able to identify the best model with the available data. In this context, reverse engineering of finite dynamical systems from partial information is an important problem. In this talk we will present a framework and algorithm to reverse engineer the possible wiring diagrams of a finite dynamical system from data. The algorithm consists on using an ideal of polynomials to encode all possible wiring diagrams, and choose those that are minimal using the primary decomposition. We will also show that these results can be applied to reverse engineer continuous dynamical systems.
The inference of molecular networks is an important problem in systems biology. This includes both the structure of the network in the form of its wiring diagram and its long-term dynamics. While there are many algorithms available that aim to infer network structure from experimental data, less emphasis has been placed on methods that utilize time series data effectively to infer both structure and dynamics. Since the network inference problem is typically underdetermined, it is important to also have the option of incorporating prior biological knowledge about the network into the process along with an effective description of the model search space. Finally, it is important to have an understanding of how a given inference method is affected by experimental and other noise in the data used for this purpose.
In this talk we will introduce a novel inference algorithm within the Boolean polynomial dynamical system (BPDS) framework, to meet all these requirements. The algorithm is able to use time series data, including network perturbations such as knock-out mutants and RNAi experiments. To infer wiring diagrams and dynamical models, it allows for the incorporation of prior biological knowledge while being robust to significant levels of noise in the data used for inference. It uses an evolutionary algorithm for local optimization with an encoding of models as BPDS. We will show how within this BPDS framework it is possible to give an effective representation of the model space to be searched to improve computational performance. We will show a validation of the algorithm with both an in silico network and with microarray expression profiles from a synthetic yeast network.
Consider the case of inferring a gene regulatory/coexpression/metabolic network from biological data in a given species. Recent studies have shown that these networks have many interesting dynamic features such as mutational robustness and evolvability. Every network in a particular species has an associated evolutionary history. Evolution can modify a gene network using mechanisms such as gene duplication, horizontal gene transfer, and neofunctionalization. In addition, sequence evolution can act on individual genes as well as their promoter regions. Thus comparative and evolutionary information may help improve the accuracy of inferred biological networks in individual species. Here we proposed to use Graphical Gaussian model and Graphical Lasso for inference on ancestral network. Here we assume that we have the known species tree. From a systems biology view, one important object is the covariance structure of genes' expression. One popular model for gene expression is the undirected gaussian graphical model. In the gaussian graphical model, the expression measurements are assumed to be guassian $N(\mu, \Sigma)$ with unknown covariance $\Sigma$. Of particular interest are the nonzero entries $K_{i,j}$ in the concentration matrix $K = \Sigma^{-1}$. A nonzero entry $K_{i,j}$ represents a nonzero saturated correlation between the ith, jth genes, after conditioning on the expression of all other genes. Thus, if $K_{i,j} \not = 0$, we may suspect that one gene regulates the other, rather than both being regulated by a confounding factor. Nonzero saturated correlations can be represented as edges in a graph, which we call the coexpression network. Using the Graphical Lasso, we can enforce the sparsity and we can estimate the concentration matrix $\hat{K}$. The glasso is fast in practice, even on thousands of genes, because the objective function is convex of special form and fast coordinate descent methods can be used. We propose an extension of the graphical lasso formulation to infer the ancestral network given the observed species tree and data sets from the current species. We hope this method to apply transcription factor networks on hematopoiesis.
Poster Session
Lake Erie is one of the Great Lakes in North America and has a favorable environment for agriculture. On the other hand, it has witnessed recurrent summertime oxygen depletion and related microbial production of greenhouse gases such as nitrous oxide (N2O). In fact, N2O is an intermediate in denitrification, which is a microbial process of conversion of nitrate (NO3) to nitrogen gas (N2). This poster will introduce the gene regulatory network and its (discrete) mathematical model of Pseudomonas aeruginosa, one of the microbes performing denitrification in Lake Erie. Polynomial Dynamical Systems (PDS) was used to model the network, and the model is analyzed by changing the concentration level of some environmental parameters such as oxygen (O2), nitrate (NO3) and phosphorus (P) to see how these parameters affect the long-run behavior of the network. Analysis is done in Analysis of Dynamic Algebraic Models (ADAM available at http://dvd.vbi.vt.edu/adam.html). This model helps us generate some hypotheses for the reason of accumulation of greenhouse gases in Lake Erie, which is still unknown.
We propose a mathematical model involving discrete dynamical systems to understand the interaction between sea-surface temperature and sea-ice cover over the Arctic Ocean. In particular, we use ideas from dynamical systems such as bifurcations and basins of attraction of equilibria along with basic probability theory to make future projections about the possibility of extinction of Arctic sea-ice cover, one of the top five tipping points in the Earth's climate.
In experiments it has been discovered that the key step of gene expression is a strongly stochastic process. Further studies can then be realized either by conducting even more experiments, or by finding a stochastic modeling framework that captures the underlying features of gene regulation. For such a framework, both, exhaustive simulations and analytic calculations deliver desired results. However, with growing complexity, simulations become increasingly time-consuming and unsuited so that analytic results are preferable.
The Laubenbacher research group at Virginia Tech has developed an easily comprehensible stochastic modeling framework that captures the ubiquitous uncertainty of biochemical networks. Also, the Derrida plot is a well-known indicator for the stability of a Boolean network. The difference between the Hamming distance of any configuration and the Hamming distance of the same configuration after applying given update rules, the so-called Derrida value, is generally small for networks that exhibit stable behavior and big for networks with more chaotic behavior.
In published Boolean networks, one subclass of Boolean functions, the class of nested canalizing functions been found to be chosen particularly often to model update functions. Recently, this subclass has been partitioned further, by characterizing each function by its Hamming weight. Explicit formulas that enable the calculation of Derrida values in the context of the stochastic modeling framework have been found; first, for nested canalizing functions in general, and second, for nested canalizing functions of a particular Hamming weight. An analysis of the derived formulas shows that networks based on nested canalizing functions do not exhibit very stable behavior in general but only those based on functions with a small Hamming weight.
In conclusion, formulas for the Derrida values of any nested canalizing function have been found. Since this makes simulations dispensable, the application of Derrida plots for robustness investigations of complex networks has become much simpler.
Work done in collaboration David Murrugarra1 and Reinhard Laubenbacher.
Iron is essential for the growth and survival of the cells in our body as well as the pathogens attacking them. Lung epithelial cells are a prime target for fungal infection because of constant exposure to airborne pathogens. We present a logical model of iron metabolism in lung epithelial cells exposed to proinflammatory cytokines and the fungus Aspergillus fumigatus. It makes predictions about the way in which lung epithelial cells sequester excess extracellular iron, along with how internal iron is stored and released from the cell. Additionally, it allows for the testing of conditions that are experimentally intractable, a process beneficial to many fields, as novel interactions and relationships can be explored without laboratory experimentation.
Work done in collaboration with Kahmya McAlpina, John Nardinib, Leslie Myintc, and Reinhard Laubenbacher
Modeling stochasticity in gene regulation is an important and complex problem in molecular systems biology. This poster introduces a new modeling framework that incorporates propensity parameters for activation and degradation. Within the discrete paradigm, where genes, proteins, and other components of molecular networks are modeled as discrete variables and interaction among these are given by logical rules representing the biochemical mechanisms governing their regulation; stochasticity is modeled at the biological function level under the assumption that even if the expression levels of the input nodes of an update function guarantee activation or degradation there is a probability that the process will not occur due to stochasticity in the process. This approach allows a finer analysis of discrete models and provides a natural set up for cell population simulations. Applications will be shown using two well-known gene regulatory networks, the outcome of lambda phage infection of bacteria and the p53-mdm2 complex.
Work done in collaboration with Alan Veliz-Cuba, Boris Aguilar, Seda Arat, and Reinhard Laubenbacher.
Many complex diseases such as cancers arise from disruptions in the signaling and regulatory components that maintain homeostasis in living organisms. Because of this, one important focus of research is identifying the effects of dysregulations in the multiple cellular signaling elements involved in these diseases. One recent approach to study dysregulations in these systems is the use of a Boolean network framework. An example of such a system is the T cell large granular lymphocyte (T-LGL) leukemia survival signaling network, in which a Boolean network model was able to successfully predict the key signals responsible for all known deregulations [1].
In the present work we use the T-LGL leukemia network model as a case study and look at the effect of combining different signals and dysregulations in the system. We focus on how the number and accessibility of cell fates (attractors) changes according to the signals and dysregulations present. By using a novel network reduction approach we are able to identify some motifs that stabilize which allow us to predict the cell fates of the system. We find that, depending on the signals and mutations present, some of these motifs stabilize in specific steady states which then determine how many of the original cell fates are accessible by the network. This suggests that the interplay between these motifs and the signals present are the driving forces that make the system choose between the different original cell fates.
Work done in collaboration with Reka Albert.
[1] Zhang et. al. Network Model of Survival Signaling in LGL Leukemia, PNAS 105 (2008).
