There is no design template more important than DNA. Within the sequences of this exquisite substance lie the design plans for each of us and for every living organism. Shaped over billions of years by the creative machinations of evolution, this design template encodes the most complex dynamical systems known to us. Yet, it is only in our lifetimes that we are able to directly edit this template and engineer our own designs. The story that I tell in this article is about our early attempts to design and commission our own control systems in living cells. Guided by what we have learned from controlling man-made systems, we are beginning to develop the theory and methodologies needed to build control systems at the molecular level, an endeavor that is as challenging as it is rewarding. If carried out responsibly, this new ability to reshape the DNA template can have a tremendous benefit for our health and well-being, and will drive major advances in basic science, industrial biotechnology, and medical therapy. In this article, I will take the readers of the proceedings on a journey through the new and promising world of rationally designed genetic control systems. Using a minimum of jargon, I will introduce them to the biological concepts needed to develop an understanding and appreciation of the main design concepts emerging in this nascent area of research. My goal is to convey my own sense of excitement about the possibilities, but it is also to impart a feeling of the opportunities that lay ahead for members of the IEEE to contribute with their own creative ideas to the shaping of this most versatile of design templates, the DNA.

This special issue brings together leading experts in the fields of systems and synthetic biology to provide a quantitative and analytical perspective of the state-of-the-art developments.

Stochastic models of biomolecular reaction networks are commonly employed in systems and synthetic biology to study the effects of stochastic fluctuations emanating from reactions involving species with low copy-numbers. For such models, the Kolmogorov's forward equation is called the chemical master equation (CME), and it is a fundamental system of linear ordinary differential equations (ODEs) that describes the evolution of the probability distribution of the random state-vector representing the copy-numbers of all the reacting species. The size of this system is given by the number of states that are accessible by the chemical system, and for most examples of interest this number is either very large or infinite. Moreover, approximations that reduce the size of the system by retaining only a finite number of important chemical states (e.g. those with non-negligible probability) result in high-dimensional ODE systems, even when the number of reacting species is small. Consequently, accurate numerical solution of the CME is very challenging, despite the linear nature of the underlying ODEs. One often resorts to estimating the solutions via computationally intensive stochastic simulations. The goal of the present paper is to develop a novel deep-learning approach for computing solution statistics of high-dimensional CMEs by reformulating the stochastic dynamics using Kolmogorov's backward equation. The proposed method leverages superior approximation properties of Deep Neural Networks (DNNs) to reliably estimate expectations under the CME solution for several user-defined functions of the state-vector. This method is algorithmically based on reinforcement learning and it only requires a moderate number of stochastic simulations (in comparison to typical simulation-based approaches) to train the "policy function". This allows not just the numerical approximation of various expectations for the CME solution but also of its sensitivities with respect to all the reaction network parameters (e.g. rate constants). We provide four examples to illustrate our methodology and provide several directions for future research.

The development of mechanistic models of biological systems is a central part of Systems Biology. One major challenge in developing these models is the accurate inference of model parameters. In recent years, nested sampling methods have gained increased attention in the Systems Biology community due to the fact that they are parallelizable and provide error estimates with no additional computations. One drawback that severely limits the usability of these methods, however, is that they require the likelihood function to be available, and thus cannot be applied to systems with intractable likelihoods, such as stochastic models. Here we present a likelihood-free nested sampling method for parameter inference which overcomes these drawbacks. This method gives an unbiased estimator of the Bayesian evidence as well as samples from the posterior. We derive a lower bound on the estimators variance which we use to formulate a novel termination criterion for nested sampling. The presented method enables not only the reliable inference of the posterior of parameters for stochastic systems of a size and complexity that is challenging for traditional methods, but it also provides an estimate of the obtained variance. We illustrate our approach by applying it to several realistically sized models with simulated data as well as recently published biological data. We also compare our developed method with the two most popular other likeliood-free approaches: pMCMC and ABC-SMC. The C++ code of the proposed methods, together with test data, is available at the github web page https://github.com/Mijan/LFNS_paper.

Biological systems use a variety of mechanisms to deal with the uncertain nature of their external and internal environments. Two of the most common motifs employed for this purpose are the incoherent feedforward (IFF) and feedback (FB) topologies. Many theoretical and experimental studies suggest that these circuits play very different roles in providing robustness to uncertainty in the cellular environment. Here, we use a control theoretic approach to analyze two common FB and IFF architectures that make use of an intermediary species to achieve regulation. We show the equivalence of both circuits topologies in suppressing static cell-to-cell variations. While both circuits can suppress variations due to input noise, they are ineffective in suppressing inherent chemical reaction stochasticity. Indeed, these circuits realize comparable improvements limited to a modest 25% variance reduction in best case scenarios. Such limitations are attributed to the use of intermediary species in regulation, and as such, they persist even for circuit architectures that combine both IFF and FB features. Intriguingly, while the FB circuits are better suited in dealing with dynamic input variability, the most significant difference between the two topologies lies not in the structural features of the circuits, but in their practical implementation considerations.

Reaction networks are systems in which the populations of a finite number of species evolve through predefined interactions. Such networks are found as modeling tools in many biological disciplines such as biochemistry, ecology, epidemiology, immunology, systems biology and synthetic biology. It is now well-established that, for small population sizes, stochastic models for biochemical reaction networks are necessary to capture randomness in the interactions. The tools for analyzing such models, however, still lag far behind their deterministic counterparts. In this paper, we bridge this gap by developing a constructive framework for examining the long-term behavior and stability properties of the reaction dynamics in a stochastic setting. In particular, we address the problems of determining ergodicity of the reaction dynamics, which is analogous to having a globally attracting fixed point for deterministic dynamics. We also examine when the statistical moments of the underlying process remain bounded with time and when they converge to their steady state values. The framework we develop relies on a blend of ideas from probability theory, linear algebra and optimization theory. We demonstrate that the stability properties of a wide class of biological networks can be assessed from our sufficient theoretical conditions that can be recast as efficient and scalable linear programs, well-known for their tractability. It is notably shown that the computational complexity is often linear in the number of species. We illustrate the validity, the efficiency and the wide applicability of our results on several reaction networks arising in biochemistry, systems biology, epidemiology and ecology. The biological implications of the results as well as an example of a non-ergodic biological network are also discussed.

The Chemical Master Equation (CME) is a cornerstone of stochastic analysis and simulation of models of biochemical reaction networks. Yet direct solutions of the CME have remained elusive. Although several approaches overcome the infinite dimensional nature of the CME through projections or other means, a common feature of proposed approaches is their susceptibility to the curse of dimensionality, i.e. the exponential growth in memory and computational requirements in the number of problem dimensions. We present a novel approach that has the potential to "lift" this curse of dimensionality. The approach is based on the use of the recently proposed Quantized Tensor Train (QTT) formatted numerical linear algebra for the low parametric, numerical representation of tensors. The QTT decomposition admits both, algorithms for basic tensor arithmetics with complexity scaling linearly in the dimension (number of species) and sub-linearly in the mode size (maximum copy number), and a numerical tensor rounding procedure which is stable and quasi-optimal. We show how the CME can be represented in QTT format, then use the exponentially-converging hp-discontinuous Galerkin discretization in time to reduce the CME evolution problem to a set of QTT-structured linear equations to be solved at each time step using an algorithm based on Density Matrix Renormalization Group (DMRG) methods from quantum chemistry. Our method automatically adapts the "basis" of the solution at every time step guaranteeing that it is large enough to capture the dynamics of interest but no larger than necessary, as this would increase the computational complexity. Our approach is demonstrated by applying it to three different examples from systems biology: independent birth-death process, an example of enzymatic futile cycle, and a stochastic switch model. The numerical results on these examples demonstrate that the proposed QTT method achieves dramatic speedups and several orders of magnitude storage savings over direct approaches.

Overexpression of the inducible heat shock protein 70, Hsp72, has broadly cytoprotective effects and improves outcome following stroke. A full understanding of how Hsp72 protects cells against injury is elusive, though several distinct mechanisms are implicated. One mechanism is its anti-inflammatory effects. We study the effects of Hsp72 overexpression on activation of the transcription factor NF-κB in microglia combining experimentation and mathematical modeling, using TNFα to stimulate a microglial cell line stably overexpressing Hsp72. We find that Hsp72 overexpression reduces the amount of NF-κB DNA binding activity, activity of the upstream kinase IKK, and amount of IκBα inhibitor phosphorylated following TNFα application. Simulations evaluating several proposed mechanisms suggest that inhibition of IKK activation is an essential component of its regulatory activities. Unexpectedly we find that Hsp72 overexpression reduces the initial amount of the RelA/p65 NF-κB subunit in cells, contributing to the attenuated response. Neither mechanism in isolation, however, is sufficient to attenuate the response, providing evidence that Hsp72 relies upon multiple mechanisms to attenuate NF-κB activation. An additional observation from our study is that the induced expression of IκBα is altered significantly in Hsp72 expressing cells. While the mechanism responsible for this observation is not known, it points to yet another means by which Hsp72 may alter the NF-κB response. This study illustrates the multi-faceted nature of Hsp72 regulation of NF-κB activation in microglia and offers further clues to a novel mechanism by which Hsp72 may protect cells against injury.

MOTIVATION: In the noisy cellular environment, stochastic fluctuations at the molecular level manifest as cell-cell variability at the population level that is quantifiable using high-throughput single-cell measurements. Such variability is rich with information about the cell's underlying gene regulatory networks, their architecture and the parameters of the biochemical reactions at their core. RESULTS: We report a novel method, called Inference for Networks of Stochastic Interactions among Genes using High-Throughput data (INSIGHT), for systematically combining high-throughput time-course flow cytometry measurements with computer-generated stochastic simulations of candidate gene network models to infer the network's stochastic model and all its parameters. By exploiting the mathematical relationships between experimental and simulated population histograms, INSIGHT achieves scalability, efficiency and accuracy while entirely avoiding approximate stochastic methods. We demonstrate our method on a synthetic gene network in bacteria and show that a detailed mechanistic model of this network can be estimated with high accuracy and high efficiency. Our method is completely general and can be used to infer models of signal-activated gene networks in any organism based solely on flow cytometry data and stochastic simulations. AVAILABILITY: A free C source code implementing the INSIGHT algorithm, together with test data is available from the authors. CONTACT: [email protected] SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

SUMMARY: SPSens is a software package for the efficient computation of stochastic parameter sensitivities of biochemical reaction networks. Parameter sensitivity analysis is a valuable tool that can be used to study robustness properties, for drug targeting, and many other purposes. However its application to stochastic models has been limited when Monte Carlo methods are required due to extremely high computational costs. SPSens provides efficient, state of the art sensitivity analysis algorithms in a single software package so that sensitivity analysis can be easily performed on stochastic models of biochemical reaction networks. SPSens implements the algorithms in C and estimates sensitivities with respect to both infinitesimal and finite perturbations to system parameters, in many cases reducing variance by orders of magnitude compared to basic methods. Included among the features of SPSens are serial and parallel command line versions, an interface with Matlab, and several example problems. AVAILABILITY: SPSens is distributed freely under GPL version 3 and can be downloaded from http://sourceforge.net/projects/spsens/. The software can be run on Linux, Mac OS X and Windows platforms.

Although cell polarity is an essential feature of living cells, it is far from being well-understood. Using a combination of computational modeling and biological experiments we closely examine an important prototype of cell polarity: the pheromone-induced formation of the yeast polarisome. Focusing on the role of noise and spatial heterogeneity, we develop and investigate two mechanistic spatial models of polarisome formation, one deterministic and the other stochastic, and compare the contrasting predictions of these two models against experimental phenotypes of wild-type and mutant cells. We find that the stochastic model can more robustly reproduce two fundamental characteristics observed in wild-type cells: a highly polarized phenotype via a mechanism that we refer to as spatial stochastic amplification, and the ability of the polarisome to track a moving pheromone input. Moreover, we find that only the stochastic model can simultaneously reproduce these characteristics of the wild-type phenotype and the multi-polarisome phenotype of a deletion mutant of the scaffolding protein Spa2. Significantly, our analysis also demonstrates that higher levels of stochastic noise results in increased robustness of polarization to parameter variation. Furthermore, our work suggests a novel role for a polarisome protein in the stabilization of actin cables. These findings elucidate the intricate role of spatial stochastic effects in cell polarity, giving support to a cellular model where noise and spatial heterogeneity combine to achieve robust biological function.

BACKGROUND: Activation of the NF-κB transcription factor and its associated gene expression in microglia is a key component in the response to brain injury. Its activation is dynamic and is part of a network of biochemical species with multiple feedback regulatory mechanisms. Mathematical modeling, which has been instrumental for understanding the NF-κB response in other cell types, offers a valuable tool to investigate the regulation of NF-κB activation in microglia at a systems level. RESULTS: We quantify the dynamic response of NF-κB activation and activation of the upstream kinase IKK using ELISA measurements of a microglial cell line following treatment with the pro-inflammatory cytokine TNFα. A new mathematical model is developed based on these data sets using a modular procedure that exploits the feedback structure of the network. We show that the new model requires previously unmodeled dynamics involved in the stimulus-induced degradation of the inhibitor IκBα in order to properly describe microglial NF-κB activation in a statistically consistent manner. This suggests a more prominent role for the ubiquitin-proteasome system in regulating the activation of NF-κB to inflammatory stimuli. We also find that the introduction of nonlinearities in the kinetics of IKK activation and inactivation is essential for proper characterization of transient IKK activity and corresponds to known biological mechanisms. Numerical analyses of the model highlight key regulators of the microglial NF-κB response, as well as those governing IKK activation. Results illustrate the dynamic regulatory mechanisms and the robust yet fragile nature of the negative feedback regulated network. CONCLUSIONS: We have developed a new mathematical model that incorporates previously unmodeled dynamics to characterize the dynamic response of the NF-κB signaling network in microglia. This model is the first of its kind for microglia and provides a tool for the quantitative, systems level study the dynamic cellular response to inflammatory stimuli.

Competence is a transiently differentiated state that certain bacterial cells reach when faced with a stressful environment. Entrance into competence can be attributed to the excitability of the dynamics governing the genetic circuit that regulates this cellular behavior. Like many biological behaviors, entrance into competence is a stochastic event. In this case cellular noise is responsible for driving the cell from a vegetative state into competence and back. In this work we present a novel numerical method for the analysis of stochastic biochemical events and use it to study the excitable dynamics responsible for competence in Bacillus subtilis. Starting with a Finite State Projection (FSP) solution of the chemical master equation (CME), we develop efficient numerical tools for accurately computing competence probability. Additionally, we propose a new approach for the sensitivity analysis of stochastic events and utilize it to elucidate the robustness properties of the competence regulatory genetic circuit. We also propose and implement a numerical method to calculate the expected time it takes a cell to return from competence. Although this study is focused on an example of cell-differentiation in Bacillus subtilis, our approach can be applied to a wide range of stochastic phenomena in biological systems.

A central challenge in computational modeling of biological systems is the determination of the model parameters. Typically, only a fraction of the parameters (such as kinetic rate constants) are experimentally measured, while the rest are often fitted. The fitting process is usually based on experimental time course measurements of observables, which are used to assign parameter values that minimize some measure of the error between these measurements and the corresponding model prediction. The measurements, which can come from immunoblotting assays, fluorescent markers, etc., tend to be very noisy and taken at a limited number of time points. In this work we present a new approach to the problem of parameter selection of biological models. We show how one can use a dynamic recursive estimator, known as extended Kalman filter, to arrive at estimates of the model parameters. The proposed method follows. First, we use a variation of the Kalman filter that is particularly well suited to biological applications to obtain a first guess for the unknown parameters. Secondly, we employ an a posteriori identifiability test to check the reliability of the estimates. Finally, we solve an optimization problem to refine the first guess in case it should not be accurate enough. The final estimates are guaranteed to be statistically consistent with the measurements. Furthermore, we show how the same tools can be used to discriminate among alternate models of the same biological process. We demonstrate these ideas by applying our methods to two examples, namely a model of the heat shock response in E. coli, and a model of a synthetic gene regulation system. The methods presented are quite general and may be applied to a wide class of biological systems where noisy measurements are used for parameter estimation or model selection.

Modeling and analysis of complex biological networks presents a number of mathematical challenges. For the models to be useful from a biological standpoint, they must be systematically compared with data. Robustness is a key to biological understanding and proper feedback to guide experiments,including both the deterministic stability and performance properties of models in the presence of parametric uncertainties and their stochastic behavior in the presence of noise. In this paper, we present mathematical and algorithmic tools to address such questions for models that may be nonlinear, hybrid,and stochastic. These tools are rooted in solid mathematical theories, primarily from robust control and dynamical systems, but with important recent developments. They also have the potential for great practical relevance, which we explore through a series of biologically motivated examples.

Biological systems have evolved complex regulatory mechanisms, even in situations where much simpler designs seem to be sufficient for generating nominal functionality. Using module-based analysis coupled with rigorous mathematical comparisons, we propose that in analogy to control engineering architectures, the complexity of cellular systems and the presence of hierarchical modular structures can be attributed to the necessity of achieving robustness. We employ the Escherichia coli heat shock response system, a strongly conserved cellular mechanism, as an example to explore the design principles of such modular architectures. In the heat shock response system, the sigma-factor sigma32 is a central regulator that integrates multiple feedforward and feedback modules. Each of these modules provides a different type of robustness with its inherent tradeoffs in terms of transient response and efficiency. We demonstrate how the overall architecture of the system balances such tradeoffs. An extensive mathematical exploration nevertheless points to the existence of an array of alternative strategies for the existing heat shock response that could exhibit similar behavior. We therefore deduce that the evolutionary constraints facing the system might have steered its architecture toward one of many robustly functional solutions.

In the technological and engineering sciences feedback control methods have been studied and applied extensively to man-made systems. Recently, there has been an increasing interest as to the applicability of feedback and control theory to the life sciences. In a previous paper (H. El-Samad et al., "Calcium homeostasis: a feedback control point of view," Proceedings of the 2000 American Control Conference, Chicago, Il, June 2000), the applicability of control theory to the analysis of calcium homeostasis was investigated. A mathematical model for this homeostatic mechanism was derived. Here, the authors provide an overview of this model, and then use it to study a pathological condition of the calcium regulatory mechanism in dairy cows. The clinical condition studied is referred to as parturient paresis (milk fever) and is the result of severe hypocalcemia at the time of parturition. The authors show that their mathematical feedback model for calcium homeostasis exhibits dynamic behavior that can correspond to milk fever. They then present a rigorous mathematical result that guarantees the presence of milk fever-like dynamics for a range of model parameters.