Iterative approach to model identification of biological networks

Integrating a tailored recurrent neural network with Bayesian experimental design to optimize microbial community functions

Microbiomes interact dynamically with their environment to perform exploitable functions such as production of valuable metabolites and degradation of toxic metabolites for a wide range of applications in human health, agriculture, and environmental cleanup. Developing computational models to predict the key bacterial species and environmental factors to build and optimize such functions are crucial to accelerate microbial community engineering. However, there is an unknown web of interactions that determine the highly complex and dynamic behavior of these systems, which precludes the development of models based on known mechanisms. By contrast, entirely data-driven machine learning models can produce physically unrealistic predictions and often require significant amounts of experimental data to learn system behavior. We develop a physically-constrained recurrent neural network that preserves model flexibility but is constrained to produce physically consistent predictions and show that it can outperform existing machine learning methods in the prediction of certain experimentally measured species abundance and metabolite concentrations. Further, we present a closed-loop, Bayesian experimental design algorithm to guide data collection by selecting experimental conditions that simultaneously maximize information gain and target microbial community functions. Using a bioreactor case study, we demonstrate how the proposed framework can be used to efficiently navigate a large design space to identify optimal operating conditions. The proposed methodology offers a flexible machine learning approach specifically tailored to optimize microbiome target functions through the sequential design of informative experiments that seek to explore and exploit community functions.

Analysis and design of single-cell experiments to harvest fluctuation information while rejecting measurement noise

Introduction: Despite continued technological improvements, measurement errors always reduce or distort the information that any real experiment can provide to quantify cellular dynamics. This problem is particularly serious for cell signaling studies to quantify heterogeneity in single-cell gene regulation, where important RNA and protein copy numbers are themselves subject to the inherently random fluctuations of biochemical reactions. Until now, it has not been clear how measurement noise should be managed in addition to other experiment design variables (e.g., sampling size, measurement times, or perturbation levels) to ensure that collected data will provide useful insights on signaling or gene expression mechanisms of interest. Methods: We propose a computational framework that takes explicit consideration of measurement errors to analyze single-cell observations, and we derive Fisher Information Matrix (FIM)-based criteria to quantify the information value of distorted experiments. Results and Discussion: We apply this framework to analyze multiple models in the context of simulated and experimental single-cell data for a reporter gene controlled by an HIV promoter. We show that the proposed approach quantitatively predicts how different types of measurement distortions affect the accuracy and precision of model identification, and we demonstrate that the effects of these distortions can be mitigated through explicit consideration during model inference. We conclude that this reformulation of the FIM could be used effectively to design single-cell experiments to optimally harvest fluctuation information while mitigating the effects of image distortion.

Differential methods for assessing sensitivity in biological models

Differential sensitivity analysis is indispensable in fitting parameters, understanding uncertainty, and forecasting the results of both thought and lab experiments. Although there are many methods currently available for performing differential sensitivity analysis of biological models, it can be difficult to determine which method is best suited for a particular model. In this paper, we explain a variety of differential sensitivity methods and assess their value in some typical biological models. First, we explain the mathematical basis for three numerical methods: adjoint sensitivity analysis, complex perturbation sensitivity analysis, and forward mode sensitivity analysis. We then carry out four instructive case studies. (a) The CARRGO model for tumor-immune interaction highlights the additional information that differential sensitivity analysis provides beyond traditional naive sensitivity methods, (b) the deterministic SIR model demonstrates the value of using second-order sensitivity in refining model predictions, (c) the stochastic SIR model shows how differential sensitivity can be attacked in stochastic modeling, and (d) a discrete birth-death-migration model illustrates how the complex perturbation method of differential sensitivity can be generalized to a broader range of biological models. Finally, we compare the speed, accuracy, and ease of use of these methods. We find that forward mode automatic differentiation has the quickest computational time, while the complex perturbation method is the simplest to implement and the most generalizable.

Clustering Based Parameter Estimation of Thyroid Hormone Pathway

An ordinary differential equation (ODE) model of the working of the thyroid system for euthyroidism has been presented. As clinical data for thyroid hormones is relatively scarce, such modelling offers potential benefits over wet lab procedures. Genetic algorithms developed for determining of parameters of the ODE system using the available data have been presented and evaluated. This approach enables subject specific parameter estimation towards characterisation of individual thyroid operation. Initially, a simple steady state model was used. Later a cosinor model for the circadian variation of thyroid hormones was used to obtain more reliable results, as indicated through sensitivity analysis in conjunction with other statistical methods. Our parameter determination method has been tested on groups of patients with similar observed values of thyroid stimulating hormone (TSH), free T ₃ and free T ₄ (identified through clustering) to determine their parameter values jointly. This approach appears to produce parameter sets with lower variation than parameters determined independently, thus leading to better parameter determination.

Fisher information for optimal planning of X-ray diffraction experiments

Fisher information is a powerful mathematical tool suitable for quantification of data `informativity' and optimization of the experimental setup and measurement conditions. Here, it is applied to X-ray diffraction and an informational approach to choosing the optimal measurement configuration is proposed. The core idea is maximization of the information which can be extracted from the measured data set by the selected analysis technique, over the sets of accessible reflections and measurement geometries. The developed approach is applied to high-resolution X-ray diffraction measurements and microstructure analysis of multilayer samples, and its efficiency and consistency are demonstrated with the results of more straightforward Monte Carlo simulations.

Development of a hybrid model for a partially known intracellular signaling pathway through correction term estimation and neural network modeling

Developing an accurate first-principle model is an important step in employing systems biology approaches to analyze an intracellular signaling pathway. However, an accurate first-principle model is difficult to be developed since it requires in-depth mechanistic understandings of the signaling pathway. Since underlying mechanisms such as the reaction network structure are not fully understood, significant discrepancy exists between predicted and actual signaling dynamics. Motivated by these considerations, this work proposes a hybrid modeling approach that combines a first-principle model and an artificial neural network (ANN) model so that predictions of the hybrid model surpass those of the original model. First, the proposed approach determines an optimal subset of model states whose dynamics should be corrected by the ANN by examining the correlation between each state and outputs through relative order. Second, an L2-regularized least-squares problem is solved to infer values of the correction terms that are necessary to minimize the discrepancy between the model predictions and available measurements. Third, an ANN is developed to generalize relationships between the values of the correction terms and the system dynamics. Lastly, the original first-principle model is coupled with the developed ANN to finalize the hybrid model development so that the model will possess generalized prediction capabilities while retaining the model interpretability. We have successfully validated the proposed methodology with two case studies, simplified apoptosis and lipopolysaccharide-induced NFκB signaling pathways, to develop hybrid models with in silico and in vitro measurements, respectively.

An intracellular signaling pathway is often represented by a set of nonlinear ordinary differential equations, which translate our current knowledge about the signaling pathway into a testable mathematical model. However, predictions from such models are often subject to high uncertainty since many signaling pathways are only partially known beforehand. In this study, we propose a systematic approach to develop a hybrid model to improve model accuracy by combining machine learning and the first-principle modeling. Specifically, model correction terms are learned from discrepancy between model predictions and measurements, and these terms are added to the first-principle model to enhance the prediction accuracy. Once these correction terms are learned from the data, an artificial neural network (ANN) model is developed to find an empirical relation between the model and the correction terms so that the developed ANN can be used to posses improved predictive capabilities even in new operating conditions (i.e., generalizability). The final hybrid model is then constructed by coupling the first-principle model with the developed ANN.

Application of Systems Engineering Principles and Techniques in Biological Big Data Analytics: A Review

In the past few decades, we have witnessed tremendous advancements in biology, life sciences and healthcare. These advancements are due in no small part to the big data made available by various high-throughput technologies, the ever-advancing computing power, and the algorithmic advancements in machine learning. Specifically, big data analytics such as statistical and machine learning has become an essential tool in these rapidly developing fields. As a result, the subject has drawn increased attention and many review papers have been published in just the past few years on the subject. Different from all existing reviews, this work focuses on the application of systems, engineering principles and techniques in addressing some of the common challenges in big data analytics for biological, biomedical and healthcare applications. Specifically, this review focuses on the following three key areas in biological big data analytics where systems engineering principles and techniques have been playing important roles: the principle of parsimony in addressing overfitting, the dynamic analysis of biological data, and the role of domain knowledge in biological data analytics.

A scalable method for parameter identification in kinetic models of metabolism using steady-state data

MOTIVATION

In kinetic models of metabolism, the parameter values determine the dynamic behaviour predicted by these models. Estimating parameters from in vivo experimental data require the parameters to be structurally identifiable, and the data to be informative enough to estimate these parameters. Existing methods to determine the structural identifiability of parameters in kinetic models of metabolism can only be applied to models of small metabolic networks due to their computational complexity. Additionally, a priori experimental design, a necessity to obtain informative data for parameter estimation, also does not account for using steady-state data to estimate parameters in kinetic models.

RESULTS

Here, we present a scalable methodology to structurally identify parameters for each flux in a kinetic model of metabolism based on the availability of steady-state data. In doing so, we also address the issue of determining the number and nature of experiments for generating steady-state data to estimate these parameters. By using a small metabolic network as an example, we show that most parameters in fluxes expressed by mechanistic enzyme kinetic rate laws can be identified using steady-state data, and the steady-state data required for their estimation can be obtained from selective experiments involving both substrate and enzyme level perturbations. The methodology can be used in combination with other identifiability and experimental design algorithms that use dynamic data to determine the most informative experiments requiring the least resources to perform.

AVAILABILITY AND IMPLEMENTATION

https://github.com/LMSE/ident.

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

Identification of a time-varying intracellular signalling model through data clustering and parameter selection: application to NF-[inline-formula removed]B signalling pathway induced by LPS in the presence of BFA

Developing a model for a signalling pathway requires several iterations of experimentation and model refinement to obtain an accurate model. However, the implementation of such an approach to model a signalling pathway induced by a poorly-known stimulus can become labour intensive because only limited information on the pathway is available beforehand to formulate an initial model. Therefore, a large number of iterations are required since the initial model is likely to be erroneous. In this work, a numerical scheme is proposed to construct a time-varying model for a signalling pathway induced by a poorly-known stimulus when its nominal model is available in the literature. Here, the nominal model refers to one that describes the signalling dynamics under a well-characterised stimulus. First, global sensitivity analysis is implemented on the nominal model to identify the most important parameters, which are assumed to be piecewise constants. Second, measurement data are clustered to determine temporal subdomains where the parameters take different values. Finally, a least-squares problem is solved to estimate the parameter values in each temporal subdomain. The effectiveness of this approach is illustrated by developing a time-varying model for NF-[inline-formula removed]B signalling dynamics induced by lipopolysaccharide in the presence of brefeldin A.

Parameter estimation in models of biological oscillators: an automated regularised estimation approach

BACKGROUND

Dynamic modelling is a core element in the systems biology approach to understanding complex biosystems. Here, we consider the problem of parameter estimation in models of biological oscillators described by deterministic nonlinear differential equations. These problems can be extremely challenging due to several common pitfalls: (i) a lack of prior knowledge about parameters (i.e. massive search spaces), (ii) convergence to local optima (due to multimodality of the cost function), (iii) overfitting (fitting the noise instead of the signal) and (iv) a lack of identifiability. As a consequence, the use of standard estimation methods (such as gradient-based local ones) will often result in wrong solutions. Overfitting can be particularly problematic, since it produces very good calibrations, giving the impression of an excellent result. However, overfitted models exhibit poor predictive power. Here, we present a novel automated approach to overcome these pitfalls. Its workflow makes use of two sequential optimisation steps incorporating three key algorithms: (1) sampling strategies to systematically tighten the parameter bounds reducing the search space, (2) efficient global optimisation to avoid convergence to local solutions, (3) an advanced regularisation technique to fight overfitting. In addition, this workflow incorporates tests for structural and practical identifiability.

RESULTS

We successfully evaluate this novel approach considering four difficult case studies regarding the calibration of well-known biological oscillators (Goodwin, FitzHugh-Nagumo, Repressilator and a metabolic oscillator). In contrast, we show how local gradient-based approaches, even if used in multi-start fashion, are unable to avoid the above-mentioned pitfalls.

CONCLUSIONS

Our approach results in more efficient estimations (thanks to the bounding strategy) which are able to escape convergence to local optima (thanks to the global optimisation approach). Further, the use of regularisation allows us to avoid overfitting, resulting in more generalisable calibrated models (i.e. models with greater predictive power).

A parallel metaheuristic for large mixed-integer dynamic optimization problems, with applications in computational biology

Background

We consider a general class of global optimization problems dealing with nonlinear dynamic models. Although this class is relevant to many areas of science and engineering, here we are interested in applying this framework to the reverse engineering problem in computational systems biology, which yields very large mixed-integer dynamic optimization (MIDO) problems. In particular, we consider the framework of logic-based ordinary differential equations (ODEs).

Methods

We present saCeSS2, a parallel method for the solution of this class of problems. This method is based on an parallel cooperative scatter search metaheuristic, with new mechanisms of self-adaptation and specific extensions to handle large mixed-integer problems. We have paid special attention to the avoidance of convergence stagnation using adaptive cooperation strategies tailored to this class of problems.

Results

We illustrate its performance with a set of three very challenging case studies from the domain of dynamic modelling of cell signaling. The simpler case study considers a synthetic signaling pathway and has 84 continuous and 34 binary decision variables. A second case study considers the dynamic modeling of signaling in liver cancer using high-throughput data, and has 135 continuous and 109 binaries decision variables. The third case study is an extremely difficult problem related with breast cancer, involving 690 continuous and 138 binary decision variables. We report computational results obtained in different infrastructures, including a local cluster, a large supercomputer and a public cloud platform. Interestingly, the results show how the cooperation of individual parallel searches modifies the systemic properties of the sequential algorithm, achieving superlinear speedups compared to an individual search (e.g. speedups of 15 with 10 cores), and significantly improving (above a 60%) the performance with respect to a non-cooperative parallel scheme. The scalability of the method is also good (tests were performed using up to 300 cores).

Conclusions

These results demonstrate that saCeSS2 can be used to successfully reverse engineer large dynamic models of complex biological pathways. Further, these results open up new possibilities for other MIDO-based large-scale applications in the life sciences such as metabolic engineering, synthetic biology, drug scheduling.

Predictive Modeling of a Batch Filter Mating Process

Quantitative characterizations of horizontal gene transfer are needed to accurately describe gene transfer processes in natural and engineered systems. A number of approaches to the quantitative description of plasmid conjugation have appeared in the literature. In this study, we seek to extend that work, motivated by the question of whether a mathematical model can accurately predict growth and conjugation dynamics in a batch process. We used flow cytometry to make time-point observations of a filter-associated mating between two E. coli strains, and fit ordinary differential equation models to the data. A model comparison analysis identified the model formulation that is best supported by the data. Identifiability analysis revealed that the parameters were estimated with acceptable accuracy. The predictive power of the model was assessed by comparison with test data that demanded extrapolation from the training experiments. This study represents the first attempt to assess the quality of model predictions for plasmid conjugation. Our successful application of this approach lays a foundation for predictive modeling that can be used both in the study of natural plasmid transmission and in model-based design of engineering approaches that employ conjugation, such as plasmid-mediated bioaugmentation.

Parameter estimation in large-scale systems biology models: a parallel and self-adaptive cooperative strategy

Background

The development of large-scale kinetic models is one of the current key issues in computational systems biology and bioinformatics. Here we consider the problem of parameter estimation in nonlinear dynamic models. Global optimization methods can be used to solve this type of problems but the associated computational cost is very large. Moreover, many of these methods need the tuning of a number of adjustable search parameters, requiring a number of initial exploratory runs and therefore further increasing the computation times.

Here we present a novel parallel method, self-adaptive cooperative enhanced scatter search (saCeSS), to accelerate the solution of this class of problems. The method is based on the scatter search optimization metaheuristic and incorporates several key new mechanisms: (i) asynchronous cooperation between parallel processes, (ii) coarse and fine-grained parallelism, and (iii) self-tuning strategies.

Results

The performance and robustness of saCeSS is illustrated by solving a set of challenging parameter estimation problems, including medium and large-scale kinetic models of the bacterium E. coli, bakerés yeast S. cerevisiae, the vinegar fly D. melanogaster, Chinese Hamster Ovary cells, and a generic signal transduction network.

The results consistently show that saCeSS is a robust and efficient method, allowing very significant reduction of computation times with respect to several previous state of the art methods (from days to minutes, in several cases) even when only a small number of processors is used.

Conclusions

The new parallel cooperative method presented here allows the solution of medium and large scale parameter estimation problems in reasonable computation times and with small hardware requirements. Further, the method includes self-tuning mechanisms which facilitate its use by non-experts. We believe that this new method can play a key role in the development of large-scale and even whole-cell dynamic models.

Electronic supplementary material

The online version of this article (doi:10.1186/s12859-016-1452-4) contains supplementary material, which is available to authorized users.

Structured model and parameter estimation in plant cell cultures of Thevetia peruviana

In this work, a mechanistic model for predicting the dynamic behavior of extracellular and intracellular nutrients, biomass production, and the main metabolites involved in the central carbon metabolism in plant cell cultures of Thevetia peruviana is presented. The proposed model is the first mechanistic model implemented for plant cell cultures of this species, and includes 28 metabolites, 33 metabolic reactions, and 61 parameters. Given the over-parametrization of the model, its nonlinear nature and the strong correlation among the effects of the parameters, a parameter estimation routine based on identifiability analysis was implemented. This routine reduces the parameter's search space by selecting the most sensitive and linearly independent parameters. Results have shown that only 19 parameters are identifiable. Finally, the model was used for analyzing the fluxes distribution in plant cell cultures of T. peruviana. This analysis shows high uptake of phosphates and parallel uptake of glucose and fructose. Furthermore, it has pointed out the main central carbon metabolism routes for promoting biomass production in this cell culture.

Controlling the time evolution of mAb N-linked glycosylation - Part II: Model-based predictions

N-linked glycosylation is known to be a crucial factor for the therapeutic efficacy and safety of monoclonal antibodies (mAbs) and many other glycoproteins. The nontemplate process of glycosylation is influenced by external factors which have to be tightly controlled during the manufacturing process. In order to describe and predict mAb N-linked glycosylation patterns in a CHO-S cell fed-batch process, an existing dynamic mathematical model has been refined and coupled to an unstructured metabolic model. High-throughput cell culture experiments carried out in miniaturized bioreactors in combination with intracellular measurements of nucleotide sugars were used to tune the parameter configuration of the coupled models as a function of extracellular pH, manganese and galactose addition. The proposed modeling framework is able to predict the time evolution of N-linked glycosylation patterns during a fed-batch process as a function of time as well as the manipulated variables. A constant and varying mAb N-linked glycosylation pattern throughout the culture were chosen to demonstrate the predictive capability of the modeling framework, which is able to quantify the interconnected influence of media components and cell culture conditions. Such a model-based evaluation of feeding regimes using high-throughput tools and mathematical models gives rise to a more rational way to control and design cell culture processes with defined glycosylation patterns. © 2016 American Institute of Chemical Engineers Biotechnol. Prog., 32:1135-1148, 2016.

Efficient Optimization of Stimuli for Model-Based Design of Experiments to Resolve Dynamical Uncertainty

This model-based design of experiments (MBDOE) method determines the input magnitudes of an experimental stimuli to apply and the associated measurements that should be taken to optimally constrain the uncertain dynamics of a biological system under study. The ideal global solution for this experiment design problem is generally computationally intractable because of parametric uncertainties in the mathematical model of the biological system. Others have addressed this issue by limiting the solution to a local estimate of the model parameters. Here we present an approach that is independent of the local parameter constraint. This approach is made computationally efficient and tractable by the use of: (1) sparse grid interpolation that approximates the biological system dynamics, (2) representative parameters that uniformly represent the data-consistent dynamical space, and (3) probability weights of the represented experimentally distinguishable dynamics. Our approach identifies data-consistent representative parameters using sparse grid interpolants, constructs the optimal input sequence from a greedy search, and defines the associated optimal measurements using a scenario tree. We explore the optimality of this MBDOE algorithm using a 3-dimensional Hes1 model and a 19-dimensional T-cell receptor model. The 19-dimensional T-cell model also demonstrates the MBDOE algorithm’s scalability to higher dimensions. In both cases, the dynamical uncertainty region that bounds the trajectories of the target system states were reduced by as much as 86% and 99% respectively after completing the designed experiments in silico. Our results suggest that for resolving dynamical uncertainty, the ability to design an input sequence paired with its associated measurements is particularly important when limited by the number of measurements.

Many mathematical models that have been developed for biological systems are limited because the complex systems are not well understood, the parameters are not known, and available data is limited and noisy. On the other hand, experiments to support model development are limited in terms of costs and time, feasible inputs and feasible measurements. MBDOE combines the mathematical models with experiment design to strategically design optimal experiments to obtain data that will contribute to the understanding of the systems. Our approach extends current capabilities of existing MBDOE techniques to make them more useful for scientists to resolve the trajectories of the system under study. It identifies the optimal conditions for stimuli and measurements that yield the most information about the system given the practical limitations. Exploration of the input space is not a trivial extension to MBDOE methods used for determining optimal measurements due to the nonlinear nature of many biological system models. The exploration of the system dynamics elicited by different inputs requires a computationally efficient and tractable approach. Our approach plans optimal experiments to reduce dynamical uncertainty in the output of selected target states of the biological system.

BioPreDyn-bench: a suite of benchmark problems for dynamic modelling in systems biology

Background

Dynamic modelling is one of the cornerstones of systems biology. Many research efforts are currently being invested in the development and exploitation of large-scale kinetic models. The associated problems of parameter estimation (model calibration) and optimal experimental design are particularly challenging. The community has already developed many methods and software packages which aim to facilitate these tasks. However, there is a lack of suitable benchmark problems which allow a fair and systematic evaluation and comparison of these contributions.

Results

Here we present BioPreDyn-bench, a set of challenging parameter estimation problems which aspire to serve as reference test cases in this area. This set comprises six problems including medium and large-scale kinetic models of the bacterium E. coli, baker’s yeast S. cerevisiae, the vinegar fly D. melanogaster, Chinese Hamster Ovary cells, and a generic signal transduction network. The level of description includes metabolism, transcription, signal transduction, and development. For each problem we provide (i) a basic description and formulation, (ii) implementations ready-to-run in several formats, (iii) computational results obtained with specific solvers, (iv) a basic analysis and interpretation.

Conclusions

This suite of benchmark problems can be readily used to evaluate and compare parameter estimation methods. Further, it can also be used to build test problems for sensitivity and identifiability analysis, model reduction and optimal experimental design methods. The suite, including codes and documentation, can be freely downloaded from the BioPreDyn-bench website, https://sites.google.com/site/biopredynbenchmarks/.

Electronic supplementary material

The online version of this article (doi:10.1186/s12918-015-0144-4) contains supplementary material, which is available to authorized users.

Individualized model discovery: the case of anemia patients

The universal sequel to chronic kidney condition (CKD) is anemia. Patients of anemia have kidneys that are incapable of performing certain basic functions such as sensing of oxygen levels to secrete erythropoietin when red blood cell counts are low. Under such conditions, external administration of human recombinant erythropoietin (EPO) is administered as alternative to improve conditions of CKD patients by increasing their hemoglobin (Hb) levels to a given therapeutic range. Presently, EPO dosing strategies extensively depend on packet inserts and on "average" responses to the medication from previous patients. Clearly dosage strategies based on these approaches are, at best, nonoptimal to EPO medication and potentially dangerous to patients that do not adhere to the notion of expected "average" response. In this work, a technique called semi-blind robust identification is provided to uniquely identify models of the individual patients of anemia based on their actual Hb responses and EPO administration. Using the a priori information and the measured input-output data of the individual patients, the procedure identifies a unique model consisting of a nominal model and the associated model uncertainty for the patients. By incorporating the effects of unknown system initial conditions, considerably small measurement samples can be used in the modeling process.

Approximate Bayesian computation schemes for parameter inference of discrete stochastic models using simulated likelihood density

BACKGROUND

Mathematical modeling is an important tool in systems biology to study the dynamic property of complex biological systems. However, one of the major challenges in systems biology is how to infer unknown parameters in mathematical models based on the experimental data sets, in particular, when the data are sparse and the regulatory network is stochastic.

RESULTS

To address this issue, this work proposed a new algorithm to estimate parameters in stochastic models using simulated likelihood density in the framework of approximate Bayesian computation. Two stochastic models were used to demonstrate the efficiency and effectiveness of the proposed method. In addition, we designed another algorithm based on a novel objective function to measure the accuracy of stochastic simulations.

CONCLUSIONS

Simulation results suggest that the usage of simulated likelihood density improves the accuracy of estimates substantially. When the error is measured at each observation time point individually, the estimated parameters have better accuracy than those obtained by a published method in which the error is measured using simulations over the entire observation time period.

A Bayesian active learning strategy for sequential experimental design in systems biology

BackgroundDynamical models used in systems biology involve unknown kinetic parameters. Setting these parameters is a bottleneck in many modeling projects. This motivates the estimation of these parameters from empirical data. However, this estimation problem has its own difficulties, the most important one being strong ill-conditionedness. In this context, optimizing experiments to be conducted in order to better estimate a system¿s parameters provides a promising direction to alleviate the difficulty of the task.ResultsBorrowing ideas from Bayesian experimental design and active learning, we propose a new strategy for optimal experimental design in the context of kinetic parameter estimation in systems biology. We describe algorithmic choices that allow to implement this method in a computationally tractable way and make it fully automatic. Based on simulation, we show that it outperforms alternative baseline strategies, and demonstrate the benefit to consider multiple posterior modes of the likelihood landscape, as opposed to traditional schemes based on local and Gaussian approximations.ConclusionThis analysis demonstrates that our new, fully automatic Bayesian optimal experimental design strategy has the potential to support the design of experiments for kinetic parameter estimation in systems biology.

Ensemble inference and inferability of gene regulatory networks

The inference of gene regulatory network (GRN) from gene expression data is an unsolved problem of great importance. This inference has been stated, though not proven, to be underdetermined implying that there could be many equivalent (indistinguishable) solutions. Motivated by this fundamental limitation, we have developed new framework and algorithm, called TRaCE, for the ensemble inference of GRNs. The ensemble corresponds to the inherent uncertainty associated with discriminating direct and indirect gene regulations from steady-state data of gene knock-out (KO) experiments. We applied TRaCE to analyze the inferability of random GRNs and the GRNs of E. coli and yeast from single- and double-gene KO experiments. The results showed that, with the exception of networks with very few edges, GRNs are typically not inferable even when the data are ideal (unbiased and noise-free). Finally, we compared the performance of TRaCE with top performing methods of DREAM4 in silico network inference challenge.

Reverse-engineering post-transcriptional regulation of gap genes in Drosophila melanogaster

Systems biology proceeds through repeated cycles of experiment and modeling. One way to implement this is reverse engineering, where models are fit to data to infer and analyse regulatory mechanisms. This requires rigorous methods to determine whether model parameters can be properly identified. Applying such methods in a complex biological context remains challenging. We use reverse engineering to study post-transcriptional regulation in pattern formation. As a case study, we analyse expression of the gap genes Krüppel, knirps, and giant in Drosophila melanogaster. We use detailed, quantitative datasets of gap gene mRNA and protein expression to solve and fit a model of post-transcriptional regulation, and establish its structural and practical identifiability. Our results demonstrate that post-transcriptional regulation is not required for patterning in this system, but is necessary for proper control of protein levels. Our work demonstrates that the uniqueness and specificity of a fitted model can be rigorously determined in the context of spatio-temporal pattern formation. This greatly increases the potential of reverse engineering for the study of development and other, similarly complex, biological processes.

The analysis of pattern-forming gene networks is largely focussed on transcriptional regulation. However, post-transcriptional events, such as translation and regulation of protein stability also play important roles in the establishment of protein expression patterns and levels. In this study, we use a reverse-engineering approach—fitting mathematical models to quantitative expression data—to analyse post-transcriptional regulation of the Drosophila gap genes Krüppel, knirps and giant, involved in segment determination during early embryogenesis. Rigorous fitting requires us to establish whether our models provide a robust and unique solution. We demonstrate, for the first time, that this can be done in the context of a complex spatio-temporal regulatory system. This is an important methodological advance for reverse-engineering developmental processes. Our results indicate that post-transcriptional regulation is not required for pattern formation, but is necessary for proper regulation of gap protein levels. Specifically, we predict that translation rates must be tuned for rapid early accumulation, and protein stability must be increased for persistence of high protein levels at late stages of gap gene expression.

Simultaneous model discrimination and parameter estimation in dynamic models of cellular systems

Background

Model development is a key task in systems biology, which typically starts from an initial model candidate and, involving an iterative cycle of hypotheses-driven model modifications, leads to new experimentation and subsequent model identification steps. The final product of this cycle is a satisfactory refined model of the biological phenomena under study. During such iterative model development, researchers frequently propose a set of model candidates from which the best alternative must be selected. Here we consider this problem of model selection and formulate it as a simultaneous model selection and parameter identification problem. More precisely, we consider a general mixed-integer nonlinear programming (MINLP) formulation for model selection and identification, with emphasis on dynamic models consisting of sets of either ODEs (ordinary differential equations) or DAEs (differential algebraic equations).

Results

We solved the MINLP formulation for model selection and identification using an algorithm based on Scatter Search (SS). We illustrate the capabilities and efficiency of the proposed strategy with a case study considering the KdpD/KdpE system regulating potassium homeostasis in Escherichia coli. The proposed approach resulted in a final model that presents a better fit to the in silico generated experimental data.

Conclusions

The presented MINLP-based optimization approach for nested-model selection and identification is a powerful methodology for model development in systems biology. This strategy can be used to perform model selection and parameter estimation in one single step, thus greatly reducing the number of experiments and computations of traditional modeling approaches.

Complexity analysis and parameter estimation of dynamic metabolic systems

A metabolic system consists of a number of reactions transforming molecules of one kind into another to provide the energy that living cells need. Based on the biochemical reaction principles, dynamic metabolic systems can be modeled by a group of coupled differential equations which consists of parameters, states (concentration of molecules involved), and reaction rates. Reaction rates are typically either polynomials or rational functions in states and constant parameters. As a result, dynamic metabolic systems are a group of differential equations nonlinear and coupled in both parameters and states. Therefore, it is challenging to estimate parameters in complex dynamic metabolic systems. In this paper, we propose a method to analyze the complexity of dynamic metabolic systems for parameter estimation. As a result, the estimation of parameters in dynamic metabolic systems is reduced to the estimation of parameters in a group of decoupled rational functions plus polynomials (which we call improper rational functions) or in polynomials. Furthermore, by taking its special structure of improper rational functions, we develop an efficient algorithm to estimate parameters in improper rational functions. The proposed method is applied to the estimation of parameters in a dynamic metabolic system. The simulation results show the superior performance of the proposed method.

Ensemble kinetic modeling of metabolic networks from dynamic metabolic profiles

Kinetic modeling of metabolic pathways has important applications in metabolic engineering, but significant challenges still remain. The difficulties faced vary from finding best-fit parameters in a highly multidimensional search space to incomplete parameter identifiability. To meet some of these challenges, an ensemble modeling method is developed for characterizing a subset of kinetic parameters that give statistically equivalent goodness-of-fit to time series concentration data. The method is based on the incremental identification approach, where the parameter estimation is done in a step-wise manner. Numerical efficacy is achieved by reducing the dimensionality of parameter space and using efficient random parameter exploration algorithms. The shift toward using model ensembles, instead of the traditional "best-fit" models, is necessary to directly account for model uncertainty during the application of such models. The performance of the ensemble modeling approach has been demonstrated in the modeling of a generic branched pathway and the trehalose pathway in Saccharomyces cerevisiae using generalized mass action (GMA) kinetics.

A quantitative model of the initiation of DNA replication in Saccharomyces cerevisiae predicts the effects of system perturbations

Background

Eukaryotic cell proliferation involves DNA replication, a tightly regulated process mediated by a multitude of protein factors. In budding yeast, the initiation of replication is facilitated by the heterohexameric origin recognition complex (ORC). ORC binds to specific origins of replication and then serves as a scaffold for the recruitment of other factors such as Cdt1, Cdc6, the Mcm2-7 complex, Cdc45 and the Dbf4-Cdc7 kinase complex. While many of the mechanisms controlling these associations are well documented, mathematical models are needed to explore the network’s dynamic behaviour. We have developed an ordinary differential equation-based model of the protein-protein interaction network describing replication initiation.

Results

The model was validated against quantified levels of protein factors over a range of cell cycle timepoints. Using chromatin extracts from synchronized Saccharomyces cerevisiae cell cultures, we were able to monitor the in vivo fluctuations of several of the aforementioned proteins, with additional data obtained from the literature. The model behaviour conforms to perturbation trials previously reported in the literature, and accurately predicts the results of our own knockdown experiments. Furthermore, we successfully incorporated our replication initiation model into an established model of the entire yeast cell cycle, thus providing a comprehensive description of these processes.

Conclusions

This study establishes a robust model of the processes driving DNA replication initiation. The model was validated against observed cell concentrations of the driving factors, and characterizes the interactions between factors implicated in eukaryotic DNA replication. Finally, this model can serve as a guide in efforts to generate a comprehensive model of the mammalian cell cycle in order to explore cancer-related phenotypes.

Parameter estimation method for improper fractional models and its application to molecular biological systems

Derived from biochemical principles, molecular biological systems can be described by a group of differential equations. Generally these differential equations contain fractional functions plus polynomials (which we call improper fractional model) as reaction rates. As a result, molecular biological systems are nonlinear in both parameters and states. It is well known that it is challenging to estimate parameters nonlinear in a model. However, in fractional functions both the denominator and numerator are linear in the parameters while polynomials are also linear in parameters. Based on this observation, we develop an iterative linear least squares method for estimating parameters in biological systems modeled by improper fractional functions. The basic idea is to transfer optimizing a nonlinear least squares objective function into iteratively solving a sequence of linear least squares problems. The developed method is applied to the estimation of parameters in a metabolism system. The simulation results show the superior performance of the proposed method for estimating parameters in such molecular biological systems.

Parameter estimation in systems biology models using spline approximation

Background

Mathematical models for revealing the dynamics and interactions properties of biological systems play an important role in computational systems biology. The inference of model parameter values from time-course data can be considered as a "reverse engineering" process and is still one of the most challenging tasks. Many parameter estimation methods have been developed but none of these methods is effective for all cases and can overwhelm all other approaches. Instead, various methods have their advantages and disadvantages. It is worth to develop parameter estimation methods which are robust against noise, efficient in computation and flexible enough to meet different constraints.

Results

Two parameter estimation methods of combining spline theory with Linear Programming (LP) and Nonlinear Programming (NLP) are developed. These methods remove the need for ODE solvers during the identification process. Our analysis shows that the augmented cost function surfaces used in the two proposed methods are smoother; which can ease the optima searching process and hence enhance the robustness and speed of the search algorithm. Moreover, the cores of our algorithms are LP and NLP based, which are flexible and consequently additional constraints can be embedded/removed easily. Eight system biology models are used for testing the proposed approaches. Our results confirm that the proposed methods are both efficient and robust.

Conclusions

The proposed approaches have general application to identify unknown parameter values of a wide range of systems biology models.

ON IDENTIFIABILITY OF NONLINEAR ODE MODELS AND APPLICATIONS IN VIRAL DYNAMICS

Ordinary differential equations (ODE) are a powerful tool for modeling dynamic processes with wide applications in a variety of scientific fields. Over the last 2 decades, ODEs have also emerged as a prevailing tool in various biomedical research fields, especially in infectious disease modeling. In practice, it is important and necessary to determine unknown parameters in ODE models based on experimental data. Identifiability analysis is the first step in determing unknown parameters in ODE models and such analysis techniques for nonlinear ODE models are still under development. In this article, we review identifiability analysis methodologies for nonlinear ODE models developed in the past one to two decades, including structural identifiability analysis, practical identifiability analysis and sensitivity-based identifiability analysis. Some advanced topics and ongoing research are also briefly reviewed. Finally, some examples from modeling viral dynamics of HIV, influenza and hepatitis viruses are given to illustrate how to apply these identifiability analysis methods in practice.

Sloppy models, parameter uncertainty, and the role of experimental design

Computational models are increasingly used to understand and predict complex biological phenomena. These models contain many unknown parameters, at least some of which are difficult to measure directly, and instead are estimated by fitting to time-course data. Previous work has suggested that even with precise data sets, many parameters are unknowable by trajectory measurements. We examined this question in the context of a pathway model of epidermal growth factor (EGF) and neuronal growth factor (NGF) signaling. Computationally, we examined a palette of experimental perturbations that included different doses of EGF and NGF as well as single and multiple gene knockdowns and overexpressions. While no single experiment could accurately estimate all of the parameters, experimental design methodology identified a set of five complementary experiments that could. These results suggest optimism for the prospects for calibrating even large models, that the success of parameter estimation is intimately linked to the experimental perturbations used, and that experimental design methodology is important for parameter fitting of biological models and likely for the accuracy that can be expected from them.

Global parameter estimation methods for stochastic biochemical systems

Background

The importance of stochasticity in cellular processes having low number of molecules has resulted in the development of stochastic models such as chemical master equation. As in other modelling frameworks, the accompanying rate constants are important for the end-applications like analyzing system properties (e.g. robustness) or predicting the effects of genetic perturbations. Prior knowledge of kinetic constants is usually limited and the model identification routine typically includes parameter estimation from experimental data. Although the subject of parameter estimation is well-established for deterministic models, it is not yet routine for the chemical master equation. In addition, recent advances in measurement technology have made the quantification of genetic substrates possible to single molecular levels. Thus, the purpose of this work is to develop practical and effective methods for estimating kinetic model parameters in the chemical master equation and other stochastic models from single cell and cell population experimental data.

Results

Three parameter estimation methods are proposed based on the maximum likelihood and density function distance, including probability and cumulative density functions. Since stochastic models such as chemical master equations are typically solved using a Monte Carlo approach in which only a finite number of Monte Carlo realizations are computationally practical, specific considerations are given to account for the effect of finite sampling in the histogram binning of the state density functions. Applications to three practical case studies showed that while maximum likelihood method can effectively handle low replicate measurements, the density function distance methods, particularly the cumulative density function distance estimation, are more robust in estimating the parameters with consistently higher accuracy, even for systems showing multimodality.

Conclusions

The parameter estimation methodologies described in this work have provided an effective and practical approach in the estimation of kinetic parameters of stochastic systems from either sparse or dense cell population data. Nevertheless, similar to kinetic parameter estimation in other modelling frameworks, not all parameters can be estimated accurately, which is a common problem arising from the lack of complete parameter identifiability from the available data.

Rule-based cell systems model of aging using feedback loop motifs mediated by stress responses

Investigating the complex systems dynamics of the aging process requires integration of a broad range of cellular processes describing damage and functional decline co-existing with adaptive and protective regulatory mechanisms. We evolve an integrated generic cell network to represent the connectivity of key cellular mechanisms structured into positive and negative feedback loop motifs centrally important for aging. The conceptual network is casted into a fuzzy-logic, hybrid-intelligent framework based on interaction rules assembled from a priori knowledge. Based upon a classical homeostatic representation of cellular energy metabolism, we first demonstrate how positive-feedback loops accelerate damage and decline consistent with a vicious cycle. This model is iteratively extended towards an adaptive response model by incorporating protective negative-feedback loop circuits. Time-lapse simulations of the adaptive response model uncover how transcriptional and translational changes, mediated by stress sensors NF-κB and mTOR, counteract accumulating damage and dysfunction by modulating mitochondrial respiration, metabolic fluxes, biosynthesis, and autophagy, crucial for cellular survival. The model allows consideration of lifespan optimization scenarios with respect to fitness criteria using a sensitivity analysis. Our work establishes a novel extendable and scalable computational approach capable to connect tractable molecular mechanisms with cellular network dynamics underlying the emerging aging phenotype.

The global process of aging disturbs a broad range of cellular mechanisms in a complex fashion and is not well understood. One important goal of computational approaches in aging is to develop integrated models in terms of a unifying aging theory, predicting progression of aging phenotypes grounded on molecular mechanisms. However, current experimental data incoherently reflects many isolated processes from a large diversity of approaches, biological model systems, and species, which makes such integration a challenging task. In an attempt to close this gap, we iteratively develop a fuzzy-logic cell systems model considering the interplay of damage, metabolism, and signaling by positive and negative feedback-loop motifs using relationships drawn from literature data. Because cellular biodynamics may be considered a complex control system, this approach seems particularly suitable. Here, we demonstrate that rule-based fuzzy-logic models provide semi-quantitative predictions that enhance our understanding of complex and interlocked molecular mechanisms and their implications on the aging physiome.

Parameter identifiability of power-law biochemical system models

Mathematical modeling has become an integral component in biotechnology, in which these models are frequently used to design and optimize bioprocesses. Canonical models, like power-laws within the Biochemical Systems Theory, offer numerous mathematical and numerical advantages, including built-in flexibility to simulate general nonlinear behavior. The construction of such models relies on the estimation of unknown case-specific model parameters by way of experimental data fitting, also known as inverse modeling. Despite the large number of publications on this topic, this task remains the bottleneck in canonical modeling of biochemical systems. The focus of this paper concerns with the question of identifiability of power-law models from dynamic data, that is, whether the parameter values can be uniquely and accurately identified from time-series data. Existing and newly developed parameter identifiability methods were applied to two power-law models of biochemical systems, and the results pointed to the lack of parametric identifiability as the root cause of the difficulty faced in the inverse modeling. Despite the focus on power-law models, the analyses and conclusions are extendable to other canonical models, and the issue of parameter identifiability is expected to be a common problem in biochemical system modeling.

An iterative identification procedure for dynamic modeling of biochemical networks

Background

Mathematical models provide abstract representations of the information gained from experimental observations on the structure and function of a particular biological system. Conferring a predictive character on a given mathematical formulation often relies on determining a number of non-measurable parameters that largely condition the model's response. These parameters can be identified by fitting the model to experimental data. However, this fit can only be accomplished when identifiability can be guaranteed.

Results

We propose a novel iterative identification procedure for detecting and dealing with the lack of identifiability. The procedure involves the following steps: 1) performing a structural identifiability analysis to detect identifiable parameters; 2) globally ranking the parameters to assist in the selection of the most relevant parameters; 3) calibrating the model using global optimization methods; 4) conducting a practical identifiability analysis consisting of two (a priori and a posteriori) phases aimed at evaluating the quality of given experimental designs and of the parameter estimates, respectively and 5) optimal experimental design so as to compute the scheme of experiments that maximizes the quality and quantity of information for fitting the model.

Conclusions

The presented procedure was used to iteratively identify a mathematical model that describes the NF-κB regulatory module involving several unknown parameters. We demonstrated the lack of identifiability of the model under typical experimental conditions and computed optimal dynamic experiments that largely improved identifiability properties.

The mathematics of tanning

BACKGROUND

The pigment melanin is produced by specialized cells, called melanocytes. In healthy skin, melanocytes are sparsely spread among the other cell types in the basal layer of the epidermis. Sun tanning results from an UV-induced increase in the release of melanin to neighbouring keratinocytes, the major cell type component of the epidermis as well as redistribution of melanin among these cells. Here we provide a mathematical conceptualization of our current knowledge of the tanning response, in terms of a dynamic model. The resolution level of the model is tuned to available data, and its primary focus is to describe the tanning response following UV exposure.

RESULTS

The model appears capable of accounting for available experimental data on the tanning response in different skin and photo types. It predicts that the thickness of the epidermal layer and how far the melanocyte dendrites grow out in the epidermal layers after UV exposure influence the tanning response substantially.

CONCLUSION

Despite the paucity of experimental validation data the model is constrained enough to serve as a foundation for the establishment of a theoretical-experimental research programme aimed at elucidating the more fine-grained regulatory anatomy underlying the tanning response.

Systematic identifiability testing for unambiguous mechanistic modeling--application to JAK-STAT, MAP kinase, and NF-kappaB signaling pathway models

Background

When creating mechanistic mathematical models for biological signaling processes it is tempting to include as many known biochemical interactions into one large model as possible. For the JAK-STAT, MAP kinase, and NF-κB pathways a lot of biological insight is available, and as a consequence, large mathematical models have emerged. For large models the question arises whether unknown model parameters can uniquely be determined by parameter estimation from measured data. Systematic approaches to answering this question are indispensable since the uniqueness of model parameter values is essential for predictive mechanistic modeling.

Results

We propose an eigenvalue based method for efficiently testing identifiability of large ordinary differential models and compare this approach to three existing ones. The methods are benchmarked by applying them to models of the signaling pathways mentioned above. In all cases the eigenvalue method proposed here and the orthogonal method find the largest set of identifiable parameters, thus clearly outperforming the other approaches. The identifiability analysis shows that the pathway models are not identifiable, even under the strong assumption that all system state variables are measurable. We demonstrate how the results of the identifiability analysis can be used for model simplification.

Conclusion

While it has undoubtedly contributed to recent advances in systems biology, mechanistic modeling by itself does not guarantee unambiguous descriptions of biological processes. We show that some recent signal transduction pathway models have reached a level of detail that is not warranted. Rigorous identifiability tests reveal that even if highly idealized experiments could be carried out to measure all state variables of these signaling pathways, some unknown parameters could still not be estimated. The identifiability tests therefore show that the level of detail of the investigated models is too high in principle, not just because too little experimental information is available. We demonstrate how the proposed method can be combined with biological insight, however, to simplify these models.