Search

In this article, we introduce a novel moment closure scheme based on concepts from Model Predictive Control (MPC) to accurately describe the time evolution of the statistical moments of the solution of the Chemical Master Equation (CME). The Method of Moments, a set of ordinary differential equations frequently used to consider the firstnmmoments, is generally not closed since lower-order moments depend on higher-order moments. To overcome this limitation, we interpret the moment equations as a nonlinear dynamical system, where the firstnmmoments serve as states and the closing moments serve as control input. We demonstrate the efficacy of our approach using two example systems and show that it outperforms existing closure schemes. For polynomial systems, which encompass all mass-action systems, we provide probability bounds for the error between true and estimated moment trajectories. We achieve this by combining convergence properties of a priori moment estimates from stochastic simulations with guarantees for nonlinear reference tracking MPC. Our proposed method offers an effective solution to accurately predict the time evolution of moments of the CME, which has wide-ranging implications for many fields, including biology, chemistry, and engineering.

MotivationThe Chemical Master Equation is a set of linear differential equations that describes the evolution of the probability distribution on all possible configurations of a (bio-)chemical reaction system. Since the number of configurations and therefore the dimension of the CME rapidly increases with the number of molecules, its applicability is restricted to small systems. A widely applied remedy for this challenge are moment-based approaches which consider the evolution of the first few moments of the distribution as summary statistics for the complete distribution. Here, we investigate the performance of two moment-estimation methods for reaction systems whose equilibrium distributions encounter heavy-tailedness and hence do not possess statistical moments.ResultsWe show that estimation via Stochastic Simulation Algorithm trajectories lose consistency over time and estimated moment values span a wide range of values even for large sample sizes. In comparison, the Method of Moments returns smooth moment estimates but is not able to indicate the nonexistence of the allegedly predicted moments. We furthermore analyze the negative effect of a CME solution’s heavy-tailedness on SSA run times and explain inherent difficulties.While moment estimation techniques are a commonly applied tool in the simulation of (bio-)chemical reaction networks, we conclude that they should be used with care, as neither the system definition nor the moment estimation techniques themselves reliably indicate the potential heavy-tailedness of the CME’s solution.

Many measures have been taken to develop data and modeling standards towards a F.A.I.R. data and model management by the Systems Biology community. Still, there is an ongoing debate about incentives and merits for the individual researcher to make their research results reproducible. Here, we pose the specific question whether reproducible models have higher impact in terms of citations. Therefore, we statistically analyze 328 published models which were recently classified by Tiwari et al. based on their reproducibility. For hypotheses testing, we use a flexible Bayesian approach that provides full distributional information for all quantities of interest and can handle outliers. The results show that in the period from 2013, i.e. 10 years after the introduction of SBML, to 2020, the group of reproducible models is significantly more cited than the non-reproducible group. We show that this effect is not explained by differences in journal impact factors and that it increases with additional standardization of data and error model integration via PEtab. Overall, our statistical analysis demonstrates long-term merits of reproducible modeling for the individual researcher in terms of citations. Moreover, it provides evidences for an increased use of reproducible models in the scientific community.

In recent years, epigenetic memory systems have been developed based on DNA methylation and positive feedback systems. Achieving a robust design for these systems is generally a challenging and multi-factorial task. We developed and validated a novel mathematical model to describe methylation-based epigenetic memory systems that captures switching dynamics of methylation levels and methyltransferase amounts induced by different inputs. A bifurcation analysis shows that the system operates in the bistable range, but in its current setup is not robust to changes in parameters. An expansion of the model captures heterogeneity of cell populations by accounting for distributed cell division rates. Simulations predict that the system is highly sensitive to variations in temperature, which affects cell division and the the efficiency of the zinc finger repressor. A moderate decrease in temperature leads to a highly heterogeneous response to input signals and bistability on a single cell level. The predictions of our model were confirmed by flow-cytometry experiments conducted in this study. Overall, the results of our study give insights into the functional mechanisms of methylation-based memory systems and demonstrate that the switching dynamics can be highly sensitive to experimental conditions.

The extraction of novel information from omics data is a challenging task, in particular, since the number of features (e.g. genes) often far exceeds the number of samples. In such a setting, conventional parameter estimation leads to ill-posed optimization problems, and regularization may be required. In addition, outliers can largely impact classification accuracy. Here we introduce ROSIE, an ensemble classification approach, which combines three sparse and robust classification methods for outlier detection and feature selection and further performs a bootstrap-based validity check. Outliers of ROSIE are determined by the rank product test using outlier rankings of all three methods, and important features are selected as features commonly selected by all methods. We apply ROSIE to RNA-Seq data from The Cancer Genome Atlas (TCGA) to classify observations into Triple-Negative Breast Cancer (TNBC) and non-TNBC tissue samples. The pre-processed dataset consists of [Formula: see text] genes and more than [Formula: see text] samples. We demonstrate that ROSIE selects important features and outliers in a robust way. Identified outliers are concordant with the distribution of the commonly selected genes by the three methods, and results are in line with other independent studies. Furthermore, we discuss the association of some of the selected genes with the TNBC subtype in other investigations. In summary, ROSIE constitutes a robust and sparse procedure to identify outliers and important genes through binary classification. Our approach is ad hoc applicable to other datasets, fulfilling the overall goal of simultaneously identifying outliers and candidate disease biomarkers to the targeted in therapy research and personalized medicine frameworks.

Motivation The Chemical Master Equation is a stochastic approach to describe the evolution of a (bio)chemical reaction system. Its solution is a time-dependent probability distribution on all possible configurations of the system. As this number is typically large, the Master Equation is often practically unsolvable. The Method of Moments reduces the system to the evolution of a few moments, which are described by ordinary differential equations. Those equations are not closed, since lower order moments generally depend on higher order moments. Various closure schemes have been suggested to solve this problem. Two major problems with these approaches are first that they are open loop systems, which can diverge from the true solution, and second, some of them are computationally expensive. Results Here we introduce Quasi-Entropy Closure, a moment-closure scheme for the Method of Moments. It estimates higher order moments by reconstructing the distribution that minimizes the distance to a uniform distribution subject to lower order moment constraints. Quasi-Entropy Closure can be regarded as an advancement of Zero-Information Closure, which similarly maximizes the information entropy. Results show that both approaches outperform truncation schemes. Quasi-Entropy Closure is computationally much faster than Zero-Information Closure, although both methods consider solutions on the space of configurations and hence do not completely overcome the curse of dimensionality. In addition, our scheme includes a plausibility check for the existence of a distribution satisfying a given set of moments on the feasible set of configurations. All results are evaluated on different benchmark problems. Supplementary information Supplementary data are available at Bioinformatics online.

The Chemical Master Equation is a standard approach to model biochemical reaction networks. It consists of a system of linear differential equations, in which each state corresponds to a possible configuration of the reaction system, and the solution describes a time-dependent probability distribution over all configurations. The Stochastic Simulation Algorithm (SSA) is a method to simulate sample paths from this stochastic process. Both approaches are only applicable for small systems, characterized by few reactions and small numbers of molecules. For larger systems, the CME is computationally intractable due to a large number of possible configurations, and the SSA suffers from large reaction propensities. In our study, we focus on catalytic reaction systems, in which substrates are converted by catalytic molecules. We present an alternative description of these systems, called SiCaSMA, in which the full system is subdivided into smaller subsystems with one catalyst molecule each. These single catalyst subsystems can be analyzed individually, and their solutions are concatenated to give the solution of the full system. We show the validity of our approach by applying it to two test-bed reaction systems, a reversible switch of a molecule and methyltransferase-mediated DNA methylation.

DNA methylation maintenance by DNMT1 is an essential process in mammals but molecular mechanisms connecting DNA methylation patterns and enzyme activity remain elusive. Here, we systematically analyzed the specificity of DNMT1, revealing a pronounced influence of the DNA sequences flanking the target CpG site on DNMT1 activity. We determined DNMT1 structures in complex with preferred DNA substrates revealing that DNMT1 employs flanking sequence-dependent base flipping mechanisms, with large structural rearrangements of the DNA correlating with low catalytic activity. Moreover, flanking sequences influence the conformational dynamics of the active site and cofactor binding pocket. Importantly, we show that the flanking sequence preferences of DNMT1 highly correlate with genomic methylation in human and mouse cells, and 5-azacytidine triggered DNA demethylation is more pronounced at CpG sites with flanks disfavored by DNMT1. Overall, our findings uncover the intricate interplay between CpG-flanking sequence, DNMT1-mediated base flipping and the dynamic landscape of DNA methylation., DNA methylation is one of the major epigenetic mechanisms that critically influence gene expression, genomic stability and cell differentiation. Here, the authors study DNMT1 in complex with DNA substrates and systematically analyze the mechanism and specificity of DNMT1.

We present a workflow for kinetic modeling of biocatalytic reactions which combines methods from Bayesian learning and uncertainty quantification for model calibration, model selection, evaluation, and model reduction in a consistent statistical frame-work. Our workflow is particularly tailored to sparse data settings in which a considerable variability of the parameters remains after the models have been adapted to available data, a ubiquitous problem in many real-world applications. Our workflow is exemplified on an enzyme-catalyzed two-substrate reaction mechanism describing the symmetric carboligation of 3,5-dimethoxy-benzaldehyde to (R)-3,3',5,5'-tetramethoxybenzoin catalyzed by benzaldehyde lyase from Pseudomonas fluorescens. Results indicate a substrate-dependent inactivation of enzyme, which is in accordance with other recent studies.

Many newly synthesized cellular proteins pass through the Golgi complex from where secretory transport carriers sort them to the plasma membrane and the extracellular environment. The formation of these secretory carriers at the trans-Golgi network is promoted by the protein kinase D (PKD) family of serine/threonine kinases. Here, using mathematical modeling and experimental validation of the PKD activation and substrate phosphorylation kinetics, we reveal that the expression level of the PKD substrate deleted in liver cancer 1 (DLC1), a Rho GTPase-activating protein that is inhibited by PKD-mediated phosphorylation, determines PKD activity at the Golgi membranes. RNAi-mediated depletion of DLC1 reduced PKD activity in a Rho-Rho-associated protein kinase (ROCK)-dependent manner, impaired the exocytosis of the cargo protein horseradish peroxidase, and was associated with the accumulation of the small GTPase RAB6 on Golgi membranes, indicating a protein-trafficking defect. In summary, our findings reveal that DLC1 maintains basal activation of PKD at the Golgi and Golgi secretory activity, in part by down-regulating Rho-ROCK signaling. We propose that PKD senses cytoskeletal changes downstream of DLC1 to coordinate Rho signaling with Golgi secretory function.

Finite mixture models are widely used in the life sciences for data analysis. Yet, the calibration of these models to data is still challenging as the optimization problems are often ill-posed. This holds for censored and uncensored data, and is caused by symmetries and other types of non-identifiabilities. Here, we discuss the problem of parameter estimation and model selection for finite mixture models from a theoretical perspective. We provide a review of the existing literature and illustrate the ill-posedness of the calibration problem for mixtures of uniform distributions and mixtures of normal distributions. Furthermore, we assess the effect of interval censoring on this estimation problem. Interestingly, we find that a proper treatment of censoring can facilitate the estimation of the number of mixture components compared to inference from uncensored data, which is an at first glance surprising result. The aim of the manuscript is to raise awareness of challenges in the calibration of finite mixture models and to provide an overview about available techniques.

Biological systems are intrinsically robust. The system, being an open one, is permanently exposed to perturbations due to external noise from the surrounding environment. In this context, the architecture and information processing capabilities of a signaling network play an important role in attenuating the noise. Multi-layered phosphorylation cascades are architectures prevalent in some of our major signaling pathways. In this work, we investigate the robustness of such signaling cascades with respect to external input variations. We employ local sensitivity and output-variance based sensitivity as measures of robustness for such cascades and observe that the efficiency of high frequency signal attenuation increases with the number of in the cascades. Filtering properties of cascades have been observed previously but not under a very rigorous theoretical framework and not in comparison with other models. This work provides an example of optimality versus robustness tradeoff in the design principles of biological systems.

This study elaborates on the normalization of data from Western blot experiments and its impact on parameter estimation. Western blot data have to be preprocessed appropriately in order to enable comparison across different replicates. This includes a two step normalization procedure, in which the raw signals are normalized to a loading control and additionally to a reference condition. If the signals themselves are normally distributed, the normalized data are described by ratios of normal distributions, which have some peculiarities that can complicate further analysis such as parameter estimation for biochemical network reconstruction. Here we shortly recapitulate some properties of these ratio distributions and conditions for various approximations that facilitate further analysis. We illustrate results on a case study in which Western blot data are used to infer the fold change in a knockdown experiment.

Modular Response Analysis (MRA) is a method to reconstruct signalling networks from steady-state perturbation data which has frequently been used in different settings. Since these data are usually noisy due to multi-step measurement procedures and biological variability, it is important to investigate the effect of this noise onto network reconstruction. Here we present a systematic study to investigate propagation of noise from concentration measurements to network structures. Therefore, we design an in silico study of the MAPK and the p53 signalling pathways with realistic noise settings. We make use of statistical concepts and measures to evaluate accuracy and precision of individual inferred interactions and resulting network structures. Our results allow to derive clear recommendations to optimize the performance of MRA based network reconstruction: First, large perturbations are favorable in terms of accuracy even for models with non-linear steady-state response curves. Second, a single control measurement for different perturbation experiments seems to be sufficient for network reconstruction, and third, we recommend to execute the MRA workflow with the mean of different replicates for concentration measurements rather than using computationally more involved regression strategies. European Commission Horizon 2020 German Research Foundation via the Cluster of Excellence in Simulation Technology

The relation between design principles of signaling network motifs and their robustness against intrinsic noise still remains illusive. In this work we investigate the role of cascading for coping with intrinsic noise due to stochasticity in molecular reactions. We use stochastic approaches to quantify fluctuations in the terminal kinase of phosphorylation-dephosphorylation cascade motifs and demonstrate that cascading highly affects these fluctuations. We show that this purely stochastic effect can be explained by time-varying sequestration of upstream kinase molecules. In particular, we discuss conditions on time scales and parameter regimes which lead to a reduction of output fluctuations. Our results are put into biological context by adapting rate parameters of our modeling approach to biologically feasible ranges for general binding-unbinding and phosphorylation-dephosphorylation mechanisms. Overall, this study reveals a novel role of stochastic sequestration for dynamic noise filtering in signaling cascade motifs.

New technologies to generate, store and retrieve medical and research data are inducing a rapid change in clinical and translational research and health care. Systems medicine is the interdisciplinary approach wherein physicians and clinical investigators team up with experts from biology, biostatistics, informatics, mathematics and computational modeling to develop methods to use new and stored data to the benefit of the patient. We here provide a critical assessment of the opportunities and challenges arising out of systems approaches in medicine and from this provide a definition of what systems medicine entails. Based on our analysis of current developments in medicine and healthcare and associated research needs, we emphasize the role of systems medicine as a multilevel and multidisciplinary methodological framework for informed data acquisition and interdisciplinary data analysis to extract previously inaccessible knowledge for the benefit of patients. Creating specialised computer systems to integrate big data from multiple sources may help provide tailored future healthcare. A healthcare system that meets the needs of every individual needs the ability to collect and mine data precisely and efficiently. Olaf Wolkenhauer at Rostock University, Germany and co-workers have explored the requirements for developing this ‘systems medicine’ approach. Interdisciplinary teams should be established specifically to develop software, computational models and integrated datasets for use in healthcare. Being able to harness information from huge medical datasets could dramatically improve our understanding of disease, and allow scientists to examine illnesses from the molecular scale to whole body responses. Systems medicine could also allow doctors to investigate an individual's specific needs, trial treatments through computer models, and analyze data captured on patients' smartphones, to provide safer, more personalised heathcare.

In earlier work we have developed the circuit-breaking algorithm. This algorithm uses the topology of a dynamic model for a regulatory network to construct a one-dimensional circuit-characteristic. Zeros of this characteristic correspond to the fixed points of the system. Here we use the circuit-breaking algorithm as a tool to identify feedback circuits that are responsible for fixed point bifurcations, which are a main source of complex dynamic behavior. The methods are applied to and evaluated with two biological network models, a simplified regulatory network for a switch with three interconnected positive circuits, and a more complex network for cell cycle control in Xenopus frog eggs.

Here we present a minimal mathematical model for the sphingomyelin synthase 1 (SMS1) driven conversion of ceramide to sphingomyelin based on chemical reaction kinetics. We demonstrate via mathematical analysis that this model is not able to qualitatively reproduce experimental measurements on lipid compositions after altering SMS1 activity. We prove that a positive feedback mechanism from the products to the reactants of the reaction is one possible model extension to explain these specific experimental data. The proposed mechanism in fact exists in vivo via protein kinase D and the ceramide transfer protein CERT. The model is further evaluated by additional observations from the literature.

The spindle assembly checkpoint is a conserved signalling pathway that protects genome integrity. Given its central importance, this checkpoint should withstand stochastic fluctuations and environmental perturbations, but the extent of and mechanisms underlying its robustness remain unknown. We probed spindle assembly checkpoint signalling by modulating checkpoint protein abundance and nutrient conditions in fission yeast. For core checkpoint proteins, a mere 20% reduction can suffice to impair signalling, revealing a surprising fragility. Quantification of protein abundance in single cells showed little variability (noise) of critical proteins, explaining why the checkpoint normally functions reliably. Checkpoint-mediated stoichiometric inhibition of the anaphase activator Cdc20 (Slp1 in Schizosaccharomyces pombe) can account for the tolerance towards small fluctuations in protein abundance and explains our observation that some perturbations lead to non-genetic variation in the checkpoint response. Our work highlights low gene expression noise as an important determinant of reliable checkpoint signalling.

In this work we consider convergence of the posterior distribution in a Bayesian learning framework in the large sample size limit. According to the large-sample theory, the posterior distribution converges under mild regularity conditions to a normal distribution with covariance that is proportional to the inverse of the sample size. This result is independent of the choice of prior distribution and holds for structurally identifiable problems, which correspond to well-posed problems in the theory of inverse problems. For structurally non-identifiable problems, which belong to the class of ill-posed optimization problems, concepts from large sample theory cannot directly be applied, since the influence of the prior does not vanish in the large sample limit. A different concept is needed for convergence which includes this influence of the prior, and only few results about convergence of the posterior distribution exist so far in this setting. These are mostly restricted to linear problems and normally distributed noise. We will review these results and present ideas that can be used towards a more general concept of posterior convergence for structurally non-identifiable inverse problems.

Background The description of intracellular processes based on chemical reaction kinetics has become a standard approach in the last decades, and parameter estimation poses several challenges. Sensitivity analysis is a powerful tool in model development that can aid model calibration in various ways. Results can for example be used to simplify the model by elimination or fixation of parameters that have a negligible influence on relevant model outputs. However, models are usually subject to rescaling and normalization to reference experiments, which changes the variance of the output. Thus, the results of the sensitivity analysis may change depending on the choice of these rescaling factors and reference experiments. Although it might intuitively be clear, this fact has not been addressed in the literature so far. Methods In this study we investigate the effect of model rescaling and additional normalization to a reference experiment on the outcome of two different sensitivity analyses. Results are exemplified on a model for the MAPK pathway module in PC-12 cell lines. For this purpose we apply local sensitivity analysis and a global variance-based method based on Sobol sensitivity coefficients, and compare the results for differently scaled and normalized model versions. Results Results indicate that both sensitivity analyses are invariant under simple rescaling of variables and parameters with constant factors, provided that sensitivity coefficients are normalized and that the parameter space is appropriately chosen for Sobol’s method. By contrast, normalization to a reference experiment that also depends on parameters has a large impact on the results of any sensitivity analysis, and in particular complicates the interpretation. Conclusion This work shows that, in order to perform sensitivity analysis, it is necessary to take into account the dependency on parameters of the reference condition when working with normalized model versions.

In earlier work, we have introduced the circuit-breaking algorithm (CBA) for the analysis of intracellular regulation networks. This algorithm uses the network topology to construct a one-dimensional circuit-characteristic whose zeros correspond to the fixed points of the system. In this study, we apply the CBA to monotone systems whose flow preserves a partial order with respect to some orthant cone. We consider relations between stability of fixed points and the derivative of the corresponding zeros of the circuit-characteristic. In particular, we derive sufficient conditions for instability in case of global asymptotic stability of the open-loop system. Furthermore, we fully characterize stability of the fixed points if in addition the system is monotone. Combined with the theory of monotone systems, our results are used to characterize the long-term behavior of two models for different intracellular regulation processes.

We present ideas and concepts towards defining a common framework unifying abstract metrics in order to quantify key features of technical load-bearing structures and biological systems. Our aim is to transfer biological concepts to technical systems at this abstract level rather than on the basis of their outward appearance or actual functionality. This means that the biological concept generators for load-bearing structures do not have to be load-bearing structures themselves but may instead achieve rather different functionalities. We intend to carry out this transference by generalizing graph-based abstractions of both technical and biological worlds to allow comparisons to be made at an abstract level. We focus in particular on the intrinsically competing aims of optimality versus multi-functionality and robustness. In this review, we present initial attempts towards defining suitable quantitative measures for robustness to serve as a common ground for studying technical systems and biological systems simultaneously. We discuss generic properties of a ubiquitous signalling network motif and potential relationships to a minimal model for a robust truss structure. These case studies suggest that topological complexity can serve as a common source for a design that is insensitive to perturbations and thus robust in the measures of both worlds.

In many biological processes heterogeneity within clonal cell populations is an important issue. One of the most striking examples is a population of cancer cells in which after a common, identical death signal some cells die whereas others survive. The reason for this heterogeneity is intrinsic and extrinsic noise. In this paper we present a mechanistic multi-scale modeling framework for cell populations, in which the dynamics of every individual cell is captured by a parameter dependent stochastic differential equation (SDE). Heterogeneity among individual cells is accounted for by differences in parameter values, modeling extrinsic influences. Based on the statistical properties of the extrinsic noise and the SDE model for the individual cell, a partial differential equation (PDE) model is derived. This PDE describes the evolution of the population density. To determine the statistics of the extrinsic noise from experimental population data, a density-based statistical data model of the noise-corrupted data is derived. Employing this data model we show that the statistics of the extrinsic can be computed using a convex optimization. This efficient way of assessing the parameters allows for a so far infeasible uncertainty analysis via bootstrapping. To evaluate the proposed method, a model for the caspase activation cascade is considered. It is shown that for known noise properties the unknown parameter densities in this model are well estimated by the proposed method.

Feedback circuits are network motifs that are important for the proper regulation of various cellular processes. Moreover, feedback is related to complex dynamic behavior such as bistability, hysteresis or oscillations, which has been investigated, especially in the context of intracellular networks, in recent years. Consequently, it is convenient to analyze the behavior of biological networks in terms of their feedback structure. While single feedback circuits are well-investigated, a generalization of results to more complex networks is usually not straightforward. On the other hand, various results indicate that network complexity is directly connected to a proper functioning of the cell. Thus, theoretical approaches for the analysis of complex feedback systems are needed in order to understand cellular behavior at a systemic level. This work discusses the role of feedback mechanisms for intracellular regulation processes on different biological examples. We will recapitulate results for single elementary circuits and discuss recent attempts to develop a comprehensive theory that is also applicable for more complex networks of interrelated feedback structures. We will focus mainly on fixed point analysis and discuss concepts for a calculation of those in terms of the underlying graph topology. Copyright © 2011 John Wiley and Sons Asia Pte Ltd and Chinese Automatic Control Society

Motivation: Feedback circuits are important motifs in biological networks and part of virtually all regulation processes that are needed for a reliable functioning of the cell. Mathematically, feedback is connected to complex behavior of the systems, which is often related to bifurcations of fixed points. Therefore, several approaches for the investigation of fixed points in biological networks have been developed in recent years. Many of them assume the fixed point coordinates to be known, and an efficient way to calculate the entire set of fixed points for interrelated feedback structures is highly desirable. Results: In this article, we consider regulatory network models, which are differential equations with an underlying directed graph that illustrates independencies among variables. We introduce the circuit-breaking algorithm (CBA), a method that constructs one-dimensional characteristics for these network models, which inherit important information about the system. In particular, fixed points are related to the zeros of these characteristics. The CBA operates on the graph topology, and results from graph theory are used in order to make calculations efficient. Our framework provides a general scheme for analyzing network models in terms of interrelated feedback circuits. The efficiency of the approach is demonstrated on a model for calcium oscillations based on experiments in hepatocytes, which consists of several interrelated feedback circuits. Contact: radde@ist.uni-stuttgart.de Supplementary information: Supplementary data are available at Bioinformatics online.

Observed phenotypes usually arise from complex networks of interacting cell components. Qualitative information about the structure of these networks is often available, while quantitative information may be partial or absent. It is natural then to ask what, if anything, we can learn about the behaviour of the system solely from its qualitative structure. In this article we review some techniques which can be applied to answer this question, focussing in particular on approaches involving graphical representations of model structure. By applying these techniques to various cellular network examples, we discuss their strengths and limitations, and point to future research directions.

High-throughput techniques allow measurement of hundreds of cell components simultaneously. The inference of interactions between cell components from these experimental data facilitates the understanding of complex regulatory processes. Differential equations have been established to model the dynamic behavior of these regulatory networks quantitatively. Usually traditional regression methods for estimating model parameters fail in this setting, since they overfit the data. This is even the case, if the focus is on modeling subnetworks of, at most, a few tens of components. In a Bayesian learning approach, this problem is avoided by a restriction of the search space with prior probability distributions over model parameters.This paper combines both differential equation models and a Bayesian approach. We model the periodic behavior of proteins involved in the cell cycle of the budding yeast Saccharomyces cerevisiae, with differential equations, which are based on chemical reaction kinetics. One property of these systems is that they usually converge to a steady state, and lots of efforts have been made to explain the observed periodic behavior. We introduce an approach to infer an oscillating network from experimental data. First, an oscillating core network is learned. This is extended by further components by using a Bayesian approach in a second step. A specifically designed hierarchical prior distribution over interaction strengths prevents overfitting, and drives the solutions to sparse networks with only a few significant interactions.We apply our method to a simulated and a real world dataset and reveal main regulatory interactions. Moreover, we are able to reconstruct the dynamic behavior of the network.

Modeling the dynamic behaviour of biochemical systems at a molecular level aims at understanding and predicting the interactions of macromolecules inside the cell. Models of small subsystems based on differential equations not only prepare the way for the long-term goal of understanding a whole cell, but are inherently valuable due to their ability to predict the behaviour of the subsystem for varying external conditions or parameters. Nitrogen supply is essential for prokaryotes, thus the nitrogen uptake is an interesting target for model building. The goal is to provide new information about the interactions of the relevant proteins by performing various simulations. A model based on piecewise linear differential equations is formulated for the nitrogen uptake in Corynebacterium glutamicum. We theoretically derive a model for biochemical networks and introduce a general method for the parameter estimation which is also applicable in the case of very short time series. This approach is applied to a special system concerning the nitrogen uptake using Western blot experiments. The equations are set up for the main components of this system, the optimization problem for parameter estimation is formulated and solved, and simulations for the evaluation of the model as well as for predictions are carried out. We show that model building based on differential equations can also, when only a few measurements are performed, lead to a satisfactory model which provides valuable insights into the way it's network components function. For example, we are able to make predictions about the maximal value of the time course as well as the steady-state level of the signal transduction protein GlnK in case of restricted activity of the proteases when considering the transition of nitrogen starvation to nitrogen excess or vice versa.

Biological networks are often modeled by systems of ordinary differential equations. In chemical reaction kinetics, for instance, sigmoid functions represent the regulation of gene expression via transcription factors. Solutions of these models tend to converge to a unique steady state, and feedback control mechanisms are required for a more complex dynamic behavior.This paper focuses on periodic behavior in two-component regulatory networks. Here, a key issue is that oscillations in chemical reaction systems are usually not robust with respect to parameter variations. Small variations lead to bifurcations that change the system's overall qualitative dynamic behavior. This concerns the mechanisms stabilizing periodic behavior in living cells. Using a small sample network, we demonstrate that oscillations can efficiently be stabilized by large time scale differences that correspond to reactions with different velocities. Furthermore, the inclusion of a time delay, reflecting transport and diffusion processes, has a similar effect. This suggests that processes of this kind potentially play a crucial role in biological oscillators.

Microarray chips generate large amounts of data about a cell’s state. In our work we want to analyze these data in order to describe the regulation processes within a cell. Therefore, we build a model which is capable of capturing the most relevant regulating interactions and present an approach how to calculate the parameters for the model from time-series data. This approach uses the discrete approximation method of least squares to solve a data fitting modeling problem. Furthermore, we analyze the features of our proposed system, i.e., which kinds of dynamical behaviors the system is able to show.

Motivation: A quantitative description of interactions between cell components is a major challenge in Computational Biology. As a method of choice, differential equations are used for this purpose, because they provide a detailed insight into the dynamic behavior of the system. In most cases, the number of time points of experimental time series is usually too small to estimate the parameters of a model of a whole gene regulatory network based on differential equations, such that one needs to focus on subnetworks consisting of only a few components. For most approaches, the set of components of the subsystem is given in advance and only the structure has to be estimated. However, the set of components that influence the system significantly are not always known in advance, making a method desirable that determines both, the components that are included into the model and the parameters. Results: We have developed a method that uses gene expression data as well as interaction data between cell components to define a set of genes that we use for our modeling. In a subsequent step, we estimate the parameters of our model of piecewise linear differential equations and evaluate the results simulating the behavior of the system with our model. We have applied our method to the DNA repair system of Mycobacterium tuberculosis. Our analysis predicts that the gene Rv2719c plays an important role in this system. Contact: {radde.gebert}@zpr.uni-koeln.de, chris@lanl.gov

Summary : We present a new C implementation of an advanced Markov chain Monte Carlo (MCMC) method for the sampling of ordinary differential equation (ode) model parameters. The software mcmc_clib uses the simplified manifold Metropolis-adjusted Langevin algorithm (SMMALA), which is locally adaptive; it uses the parameter manifold’s geometry (the Fisher information) to make efficient moves. This adaptation does not diminish with MC length, which is highly advantageous compared with adaptive Metropolis techniques when the parameters have large correlations and/or posteriors substantially differ from multivariate Gaussians. The software is standalone (not a toolbox), though dependencies include the GNU scientific library and sundials libraries for ode integration and sensitivity analysis. Availability and implementation : The source code and binary files are freely available for download at http://a-kramer.github.io/mcmc_clib/ . This also includes example files and data. A detailed documentation, an example model and user manual are provided with the software. Contact : andrei.kramer@ist.uni-stuttgart.de

Background Parameter estimation for differential equation models of intracellular processes is a highly relevant bu challenging task. The available experimental data do not usually contain enough information to identify all parameters uniquely, resulting in ill-posed estimation problems with often highly correlated parameters. Sampling-based Bayesian statistical approaches are appropriate for tackling this problem. The samples are typically generated via Markov chain Monte Carlo, however such methods are computationally expensive and their convergence may be slow, especially if there are strong correlations between parameters. Monte Carlo methods based on Euclidean or Riemannian Hamiltonian dynamics have been shown to outperform other samplers by making proposal moves that take the local sensitivities of the system’s states into account and accepting these moves with high probability. However, the high computational cost involved with calculating the Hamiltonian trajectories prevents their widespread use for all but the smallest differential equation models. The further development of efficient sampling algorithms is therefore an important step towards improving the statistical analysis of predictive models of intracellular processes. Results We show how state of the art Hamiltonian Monte Carlo methods may be significantly improved for steady state dynamical models. We present a novel approach for efficiently calculating the required geometric quantities by tracking steady states across the Hamiltonian trajectories using a Newton-Raphson method and employing local sensitivity information. Using our approach, we compare both Euclidean and Riemannian versions of Hamiltonian Monte Carlo on three models for intracellular processes with real data and demonstrate at least an order of magnitude improvement in the effective sampling speed. We further demonstrate the wider applicability of our approach to other gradient based MCMC methods, such as those based on Langevin diffusions. Conclusion Our approach is strictly benefitial in all test cases. The Matlab sources implementing our MCMC methodology is available from https://github.com/a-kramer/ode_rmhmc. Electronic supplementary material The online version of this article (doi:10.1186/1471-2105-15-253) contains supplementary material, which is available to authorized users.

We present a visual analytics system that supports an uncertainty-aware analysis of static and dynamic attributes of biochemical reaction networks (BRNs). These are often described by mathematical models, such as ordinary differential equations (ODEs), which enable the integration of a multitude of different data and data types using parameter estimation. Due to the limited amount of data, parameter estimation does not necessarily yield a single point in parameter space and many attributes of the model remain uncertain. Our system visualizes the model as a graph, where the statistics of the attributes are mapped to the color of edges and vertices. The graph view is combined with several linked views such as lineplots, scatterplots, and correlation matrices, to support the identification of uncertainties and the analysis of their mutual dependencies as well as their time dependencies. To assess the utility of the individual visualization approaches and multiple linked views, a qualitative user study with domain experts was performed. We found that all users were able to process analysis tasks using our system.

Motivation: Experiment design strategies for biomedical models with the purpose of parameter estimation or model discrimination are in the focus of intense research. Experimental limitations such as sparse and noisy data result in unidentifiable parameters and render-related design tasks challenging problems. Often, the temporal resolution of data is a limiting factor and the amount of possible experimental interventions is finite. To address this issue, we propose a Bayesian experiment design algorithm to minimize the prediction uncertainty for a given set of experiments and compare it to traditional A-optimal design. Results: In an in depth numerical study involving an ordinary differential equation model of the trans-Golgi network with 12 partly non-identifiable parameters, we minimized the prediction uncertainty efficiently for predefined scenarios. The introduced method results in twice the prediction precision as the same amount of A-optimal designed experiments while introducing a useful stopping criterion. The simulation intensity of the algorithm's major design step is thereby reasonably affordable. Besides smaller variances in the predicted trajectories compared with Fisher design, we could also achieve smaller parameter posterior distribution entropies, rendering this method superior to A-optimal Fisher design also in the parameter space. Availability: Necessary software/toolbox information are available in the supplementary material. The project script including example data can be downloaded from http://www.ist.uni-stuttgart.de/%7eweber/BayesFisher2012. Contact: patrick.weber@ist.uni-stuttgart.de Supplementary Information: Supplementary data are available at Bioinformatics online.