Biological systems are inherently variable, with their dynamics influenced by intrinsic and extrinsic sources. These systems are often only partially characterized, with large uncertainties about specific sources of extrinsic variability and biochemical properties. Moreover, it is not yet well understood how different sources of variability combine and affect biological systems in concert. To successfully design biomedical therapies or synthetic circuits with robust performance, it is crucial to account for uncertainty and effects of variability. Here we introduce an efficient modeling and simulation framework to study systems that are simultaneously subject to multiple sources of variability, and apply it to make design decisions on small genetic networks that play a role of basic design elements of synthetic circuits. Specifically, the framework was used to explore the effect of transcriptional and post-transcriptional autoregulation on fluctuations in protein expression in simple genetic networks. We found that autoregulation could either suppress or increase the output variability, depending on specific noise sources and network parameters. We showed that transcriptional autoregulation was more successful than post-transcriptional in suppressing variability across a wide range of intrinsic and extrinsic magnitudes and sources. We derived the following design principles to guide the design of circuits that best suppress variability: (i) high protein cooperativity and low miRNA cooperativity, (ii) imperfect complementarity between miRNA and mRNA was preferred to perfect complementarity, and (iii) correlated expression of mRNA and miRNA--for example, on the same transcript--was best for suppression of protein variability. Results further showed that correlations in kinetic parameters between cells affected the ability to suppress variability, and that variability in transient states did not necessarily follow the same principles as variability in the steady state. Our model and findings provide a general framework to guide design principles in synthetic biology.

Computing parametric sensitivities for oscillators has a now well-understood subtlety associated with the indeterminacy of phase. A less universal, but still vexing, subtlety arises when an oscillator is described by a system of differential equations with "hidden" conservation constraints (HCC's); defined as weighted sums of state variables that are time-invariant. If there are HCC's, as is commonly the case for models of biochemical oscillators but rarely the case for practical circuit oscillators, the now-standard approach to computing parametric sensitivities can yield incorrect results. In addition, the monodromy matrix (the matrix of state sensitivities over one oscillation period), is often defective in a way that interferes with the usual approach to computing oscillator phase noise. In this paper we analyze the HCC case, and show that by augmenting the standard sensitivity approach with explicit HCC's, one can recover the correct parametric sensitivities. In addition, we prove that there is a typically satisfied condition that guarantees that a system with HCCs will have a defective monodromy matrix. A deliberately "flawed" ring oscillator circuit and a cyanobacterial circadian clock biochemical oscillator are used to demonstrate the parametric sensitivity problem and its resolution, and to show the issue of the defective monodromy matrix.

MOTIVATION: The study of complex biological relationships is aided by large and high-dimensional data sets whose analysis often involves dimension reduction to highlight representative or informative directions of variation. In principle, information theory provides a general framework for quantifying complex statistical relationships for dimension reduction. Unfortunately, direct estimation of high-dimensional information theoretic quantities, such as entropy and mutual information (MI), is often unreliable given the relatively small sample sizes available for biological problems. Here, we develop and evaluate a hierarchy of approximations for high-dimensional information theoretic statistics from associated low-order terms, which can be more reliably estimated from limited samples. Due to a relationship between this metric and the minimum spanning tree over a graph representation of the system, we refer to these approximations as MIST (Maximum Information Spanning Trees). RESULTS: The MIST approximations are examined in the context of synthetic networks with analytically computable entropies and using experimental gene expression data as a basis for the classification of multiple cancer types. The approximations result in significantly more accurate estimates of entropy and MI, and also correlate better with biological classification error than direct estimation and another low-order approximation, minimum-redundancy-maximum-relevance (mRMR). AVAILABILITY: Software to compute the entropy approximations described here is available as Supplementary Material. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

The need to determine electrostatic fields in domains bounded by molecular surfaces arises in a number of nanotechnology applications including: biomolecule design, carbon nanotube simulation, and molecular electron transport analysis. Molecular surfaces are typically smooth, without the corners common in electrical interconnect problems, but are often so geometrically complicated that numerical evaluation of the associated electrostatic fields is extremely time-consuming. In this paper we describe and demonstrate a meshless spectrally-accurate integral equation method that only requires a description of the molecular surface in the form of a collection of surface points. Our meshless method is a synthesis of techniques, suitably adapted, including: spherical harmonic surface interpolation, spectral-element-like integral equation discretization, integral desingularization via variable transformation, and matrix-implicit iterative matrix solution. The spectral accuracy of this combined method is verified using analytically solvable sphere and ellipsoid problems, and then its accuracy and efficiency is demonstrated numerically by solving capacitance and coupled Poisson/linearized Poisson-Boltzmann problems associated with a commonly used model of a molecule in solution. The results demonstrate that for a tolerance of 10 −3 this new approach reduces the number of unknowns by as much as two orders of magnitude over the more commonly used flat panel methods.

Mechanism-based chemical kinetic models are increasingly being used to describe biological signaling. Such models serve to encapsulate current understanding of pathways and to enable insight into complex biological processes. One challenge in model development is that, with limited experimental data, multiple models can be consistent with known mechanisms and existing data. Here, we address the problem of model ambiguity by providing a method for designing dynamic stimuli that, in stimulus-response experiments, distinguish among parameterized models with different topologies, i.e., reaction mechanisms, in which only some of the species can be measured. We develop the approach by presenting two formulations of a model-based controller that is used to design the dynamic stimulus. In both formulations, an input signal is designed for each candidate model and parameterization so as to drive the model outputs through a target trajectory. The quality of a model is then assessed by the ability of the corresponding controller, informed by that model, to drive the experimental system. We evaluated our method on models of antibody-ligand binding, mitogen-activated protein kinase (MAPK) phosphorylation and de-phosphorylation, and larger models of the epidermal growth factor receptor (EGFR) pathway. For each of these systems, the controller informed by the correct model is the most successful at designing a stimulus to produce the desired behavior. Using these stimuli we were able to distinguish between models with subtle mechanistic differences or where input and outputs were multiple reactions removed from the model differences. An advantage of this method of model discrimination is that it does not require novel reagents, or altered measurement techniques; the only change to the experiment is the time course of stimulation. Taken together, these results provide a strong basis for using designed input stimuli as a tool for the development of cell signaling models.

Processes that repeat in time, such as the cell cycle, the circadian rhythm, and seasonal variations, are prevalent in biology. Mathematical models can represent our knowledge of the underlying mechanisms, and numerical methods can then facilitate analysis, which forms the foundation for a more integrated understanding as well as for design and intervention. Here, the intracellular molecular network responsible for the mammalian circadian clock system was studied. A new formulation of detailed sensitivity analysis is introduced and applied to elucidate the influence of individual rate processes, represented through their parameters, on network functional characteristics. One of four negative feedback loops in the model, the Per2 loop, was uniquely identified as most responsible for setting the period of oscillation; none of the other feedback loops were found to play as substantial a role. The analysis further suggested that the activity of the kinases CK1delta and CK1varepsilon were well placed within the network such that they could be instrumental in implementing short-term adjustments to the period in the circadian clock system. The numerical results reported here are supported by previously published experimental data.

This paper presents a fast boundary-element method (BEM) algorithm that is well suited for solving electrostatics problems that arise in traditional and bio-microelectromechanical systems (bio-MEMS) design. The algorithm, FFTSVD, is Green's-function-independent for low-frequency kernels and efficient for inhomogeneous problems. FFTSVD is a multiscale algorithm that decomposes the problem domain using an octree and uses sampling to calculate low-rank approximations to dominant source distributions and responses. Long-range interactions at each length scale are computed using the FFT. Computational results illustrate that the FFTSVD algorithm performs better than precorrected-FFT (pFFT)-style algorithms or the multipole-style algorithms in FastCap.

Computer simulation is an important tool for improving our understanding of biomolecule electrostatics, in part to aid in drug design. However, the numerical techniques used in these simulation tools do not exploit fast solver approaches widely used in analyzing integrated circuit interconnects. In this paper we describe one popular formulation used to analyze biomolecule electrostatics, present an integral formulation of the problem, and apply the precorrected-FFT method to accelerate the solution of the integral equations.