Parameter estimation in ordinary differential equations for biochemical processes using the method of multiple shooting

Parameter Estimation for Dynamical Systems Using a Deep Neural Network

The deep neural network (DNN) was applied for estimating a set of unknown parameters of a dynamical system whose measured data are given for a set of discrete time points. We developed a new vectorized algorithm that takes the number of unknowns (state variables) and number of parameters into consideration. The algorithm, first, trains the network to determine weights and biases. Next, the algorithm solves the systems of algebraic equations to estimate the parameters of the system. If the right hand side function of the system is smooth and the system have equal numbers of unknowns and parameters, the algorithm solves the algebraic equation at the discrete point where absolute error between the neural network solutions and the measured data is minimum. This improves the accuracy and reduces computational time. Several tests were carried out in linear and non-linear dynamical systems. Last, we showed that the DNN approach is more successful in terms of computational time as the number of hidden layers increases.

Isotope-assisted metabolic flux analysis as an equality-constrained nonlinear program for improved scalability and robustness

Stable isotope-assisted metabolic flux analysis (MFA) is a powerful method to estimate carbon flow and partitioning in metabolic networks. At its core, MFA is a parameter estimation problem wherein the fluxes and metabolite pool sizes are model parameters that are estimated, via optimization, to account for measurements of steady-state or isotopically-nonstationary isotope labeling patterns. As MFA problems advance in scale, they require efficient computational methods for fast and robust convergence. The structure of the MFA problem enables it to be cast as an equality-constrained nonlinear program (NLP), where the equality constraints are constructed from the MFA model equations, and the objective function is defined as the sum of squared residuals (SSR) between the model predictions and a set of labeling measurements. This NLP can be solved by using an algebraic modeling language (AML) that offers state-of-the-art optimization solvers for robust parameter estimation and superior scalability to large networks. When implemented in this manner, the optimization is performed with no distinction between state variables and model parameters. During each iteration of such an optimization, the system state is updated instead of being calculated explicitly from scratch, and this occurs concurrently with improvement in the model parameter estimates. This optimization approach starkly contrasts with traditional “shooting” methods where the state variables and model parameters are kept distinct and the system state is computed afresh during each iteration of a stepwise optimization. Our NLP formulation uses the MFA modeling framework of Wiechert et al. [1], which is amenable to incorporation of the model equations into an NLP. The NLP constraints consist of balances on either elementary metabolite units (EMUs) or cumomers. In this formulation, both the steady-state and isotopically-nonstationary MFA (inst-MFA) problems may be solved as an NLP. For the inst-MFA case, the ordinary differential equation (ODE) system describing the labeling dynamics is transcribed into a system of algebraic constraints for the NLP using collocation. This large-scale NLP may be solved efficiently using an NLP solver implemented on an AML. In our implementation, we used the reduced gradient solver CONOPT, implemented in the General Algebraic Modeling System (GAMS). The NLP framework is particularly advantageous for inst-MFA, scaling well to large networks with many free parameters, and having more robust convergence properties compared to the shooting methods that compute the system state and sensitivities at each iteration. Additionally, this NLP approach supports the use of tandem-MS data for both steady-state and inst-MFA when the cumomer framework is used. We assembled a software, eiFlux, written in Python and GAMS that uses the NLP approach and supports both steady-state and inst-MFA. We demonstrate the effectiveness of the NLP formulation on several examples, including a genome-scale inst-MFA model, to highlight the scalability and robustness of this approach. In addition to typical inst-MFA applications, we expect that this framework and our associated software, eiFlux, will be particularly useful for applying inst-MFA to complex MFA models, such as those developed for eukaryotes (e.g. algae) and co-cultures with multiple cell types.

Isotope-assisted metabolic flux analysis (MFA) is a computationally intensive parameter estimation problem. Isotopically nonstationary MFA (inst-MFA) represents the most computationally burdensome MFA application. We present the formulation of the steady-state and inst-MFA problems as equality-constrained nonlinear programs (NLPs), solved by a state-of-the-art solver implemented in an algebraic modeling language. We show that this formulation leads to robust convergence properties compared to traditional approaches, particularly for inst-MFA. We developed a software, eiFlux that uses the NLP formulation to perform both steady-state and inst-MFA. We demonstrate the application of eiFlux on several examples, including a genome-scale inst-MFA model, and show that it has robust optimal convergence even when started from a very poor initial guess for the parameters. eiFlux is implemented using the Python programming language and the General Algebraic Modeling System (GAMS), using the CONOPT solver. eiFlux is available upon request, pending institutional approval, and is free for academic use.

Guidelines for benchmarking of optimization-based approaches for fitting mathematical models

Insufficient performance of optimization-based approaches for the fitting of mathematical models is still a major bottleneck in systems biology. In this article, the reasons and methodological challenges are summarized as well as their impact in benchmark studies. Important aspects for achieving an increased level of evidence for benchmark results are discussed. Based on general guidelines for benchmarking in computational biology, a collection of tailored guidelines is presented for performing informative and unbiased benchmarking of optimization-based fitting approaches. Comprehensive benchmark studies based on these recommendations are urgently required for the establishment of a robust and reliable methodology for the systems biology community.

Bayesian inference for parameter estimation in lactoferrin-mediated iron transport across blood-brain barrier

BACKGROUND

In neurodegenerative diseases such as Alzheimer's and Parkinson's, excessive irons as well as lactoferrin (Lf), but not transferrin (Tf), have been found in and around the affected regions of the brain. These evidences suggest that lactoferrin plays a critical role during neurodegenerative diseases, although Lf-mediated iron transport across blood-brain barrier (BBB) is negligible compared to that of transferrin in normal condition. However, the kinetics of lactoferrins and lactoferrin-mediated iron transport are still unknown.

METHOD

To determine the kinetic rate constants of lactoferrin-mediated iron transport through BBB, a mass-action based ordinary differential equation model has been presented. A Bayesian framework is developed to estimate the kinetic rate parameters from posterior probability density functions. The iron transport across BBB is studied by considering both Lf- and Tf-mediated pathways for both normal and pathologic conditions.

RESULTS

Using the point estimates of kinetic parameters, our model can effectively reproduce the experimental data of iron transport through BBB endothelial cells. The robustness of the model and parameter estimation process are further verified by perturbation of kinetic parameters. Our results show that surge in high-affinity receptor density increases lactoferrin as well as iron in the brain.

CONCLUSIONS

Due to the lack of a feedback loop such as iron regulatory proteins (IRPs) for lactoferrin, iron can transport to the brain continuously, which might increase brain iron to pathological levels and can contribute to neurodegeneration.

GENERAL SIGNIFICANCE

This study provides an improved understanding of presence of lactoferrin and iron in the brain during neurodegenerative diseases.

A unified approach for sparse dynamical system inference from temporal measurements

Motivation

Temporal variations in biological systems and more generally in natural sciences are typically modeled as a set of ordinary, partial or stochastic differential or difference equations. Algorithms for learning the structure and the parameters of a dynamical system are distinguished based on whether time is discrete or continuous, observations are time-series or time-course and whether the system is deterministic or stochastic, however, there is no approach able to handle the various types of dynamical systems simultaneously.

Results

In this paper, we present a unified approach to infer both the structure and the parameters of non-linear dynamical systems of any type under the restriction of being linear with respect to the unknown parameters. Our approach, which is named Unified Sparse Dynamics Learning (USDL), constitutes of two steps. First, an atemporal system of equations is derived through the application of the weak formulation. Then, assuming a sparse representation for the dynamical system, we show that the inference problem can be expressed as a sparse signal recovery problem, allowing the application of an extensive body of algorithms and theoretical results. Results on simulated data demonstrate the efficacy and superiority of the USDL algorithm under multiple interventions and/or stochasticity. Additionally, USDL’s accuracy significantly correlates with theoretical metrics such as the exact recovery coefficient. On real single-cell data, the proposed approach is able to induce high-confidence subgraphs of the signaling pathway.

Availability and implementation

Source code is available at Bioinformatics online. USDL algorithm has been also integrated in SCENERY (http://scenery.csd.uoc.gr/); an online tool for single-cell mass cytometry analytics.

Supplementary information

Supplementary data are available at Bioinformatics online.

Mathematical Analysis of a Transformed ODE from a PDE Multiscale Model of Hepatitis C Virus Infection

Mathematical modeling has revealed the quantitative dynamics during the process of viral infection and evolved into an important tool in modern virology. Coupled with analyses of clinical and experimental data, the widely used basic model of viral dynamics described by ordinary differential equations (ODEs) has been well parameterized. In recent years, age-structured models, called "multiscale model," formulated by partial differential equations (PDEs) have also been developed and become useful tools for more detailed data analysis. However, in general, PDE models are considerably more difficult to subject to mathematical and numerical analyses. In our recently reported study, we successfully derived a mathematically identical ODE model from a PDE model, which helps to overcome the limitations of the PDE model with regard to clinical data analysis. Here, we derive the basic reproduction number from the identical ODE model and provide insight into the global stability of all possible steady states of the ODE model.

Finite-time parametric identification for the model representing the metabolic and genetic regulatory effects of sequential aerobic respiration and anaerobic fermentation processes in Escherichia coli

Mathematical modelling applied to biological systems allows for the inferring of changes in the dynamic behaviour of organisms associated with variations in the environment. Models based on ordinary differential equations are most commonly used because of their ability to describe the mechanisms of biological systems such as transcription. The disadvantage of using this approach is that there is a large number of parameters involved and that it is difficult to obtain them experimentally. This study presents an algorithm to obtain a finite-time parameter characterization of the model used to describe changes in the metabolic behaviour of Escherichia coli associated with environmental changes. In this scheme, super-twisting algorithm was proposed to recover the derivative of all the proteins and mRNA of E. coli associated to changes in the concentration of oxygen available in the growth media. The 75 identified parameters in this study maintain the biological coherence of the system and they were estimated with no more than 20% error with respect to the real ones included in the proposed model.

Dissecting Cell-Fate Determination Through Integrated Mathematical Modeling of the ERK/MAPK Signaling Pathway

The past three decades have witnessed an enormous progress in the elucidation of the ERK/MAPK signaling pathway and its involvement in various cellular processes. Because of its importance and complex wiring, the ERK pathway has been an intensive subject for mathematical modeling, which facilitates the unraveling of key dynamic properties and behaviors of the pathway. Recently, however, it became evident that the pathway does not act in isolation but closely interacts with many other pathways to coordinate various cellular outcomes under different pathophysiological contexts. This has led to an increasing number of integrated, large-scale models that link the ERK pathway to other functionally important pathways. In this chapter, we first discuss the essential steps in model development and notable models of the ERK pathway. We then use three examples of integrated, multipathway models to investigate how crosstalk of ERK signaling with other pathways regulates cell-fate decision-making in various physiological and disease contexts. Specifically, we focus on ERK interactions with the phosphoinositide-3 kinase (PI3K), c-Jun N-terminal kinase (JNK), and β-adrenergic receptor (β-AR) signaling pathways. We conclude that integrated modeling in combination with wet-lab experimentation have been and will be instrumental in gaining an in-depth understanding of ERK signaling in multiple biological contexts.

An Improved Cuckoo Search Optimization Algorithm for the Problem of Chaotic Systems Parameter Estimation

This paper proposes an improved cuckoo search (ICS) algorithm to establish the parameters of chaotic systems. In order to improve the optimization capability of the basic cuckoo search (CS) algorithm, the orthogonal design and simulated annealing operation are incorporated in the CS algorithm to enhance the exploitation search ability. Then the proposed algorithm is used to establish parameters of the Lorenz chaotic system and Chen chaotic system under the noiseless and noise condition, respectively. The numerical results demonstrate that the algorithm can estimate parameters with high accuracy and reliability. Finally, the results are compared with the CS algorithm, genetic algorithm, and particle swarm optimization algorithm, and the compared results demonstrate the method is energy-efficient and superior.

A Parameter Estimation Method for Biological Systems modelled by ODE/DDE Models Using Spline Approximation and Differential Evolution Algorithm

The inverse problem of identifying unknown parameters of known structure dynamical biological systems, which are modelled by ordinary differential equations or delay differential equations, from experimental data is treated in this paper. A two stage approach is adopted: first, combine spline theory and Nonlinear Programming (NLP), the parameter estimation problem is formulated as an optimization problem with only algebraic constraints; then, a new differential evolution (DE) algorithm is proposed to find a feasible solution. The approach is designed to handle problem of realistic size with noisy observation data. Three cases are studied to evaluate the performance of the proposed algorithm: two are based on benchmark models with priori-determined structure and parameters; the other one is a particular biological system with unknown model structure. In the last case, only a set of observation data available and in this case a nominal model is adopted for the identification. All the test systems were successfully identified by using a reasonable amount of experimental data within an acceptable computation time. Experimental evaluation reveals that the proposed method is capable of fast estimation on the unknown parameters with good precision.

Biochemical systems identification by a random drift particle swarm optimization approach

Background

Finding an efficient method to solve the parameter estimation problem (inverse problem) for nonlinear biochemical dynamical systems could help promote the functional understanding at the system level for signalling pathways. The problem is stated as a data-driven nonlinear regression problem, which is converted into a nonlinear programming problem with many nonlinear differential and algebraic constraints. Due to the typical ill conditioning and multimodality nature of the problem, it is in general difficult for gradient-based local optimization methods to obtain satisfactory solutions. To surmount this limitation, many stochastic optimization methods have been employed to find the global solution of the problem.

Results

This paper presents an effective search strategy for a particle swarm optimization (PSO) algorithm that enhances the ability of the algorithm for estimating the parameters of complex dynamic biochemical pathways. The proposed algorithm is a new variant of random drift particle swarm optimization (RDPSO), which is used to solve the above mentioned inverse problem and compared with other well known stochastic optimization methods. Two case studies on estimating the parameters of two nonlinear biochemical dynamic models have been taken as benchmarks, under both the noise-free and noisy simulation data scenarios.

Conclusions

The experimental results show that the novel variant of RDPSO algorithm is able to successfully solve the problem and obtain solutions of better quality than other global optimization methods used for finding the solution to the inverse problems in this study.

Stochastic differential equations as a tool to regularize the parameter estimation problem for continuous time dynamical systems given discrete time measurements

In this paper we consider the problem of estimating parameters in ordinary differential equations given discrete time experimental data. The impact of going from an ordinary to a stochastic differential equation setting is investigated as a tool to overcome the problem of local minima in the objective function. Using two different models, it is demonstrated that by allowing noise in the underlying model itself, the objective functions to be minimized in the parameter estimation procedures are regularized in the sense that the number of local minima is reduced and better convergence is achieved. The advantage of using stochastic differential equations is that the actual states in the model are predicted from data and this will allow the prediction to stay close to data even when the parameters in the model is incorrect. The extended Kalman filter is used as a state estimator and sensitivity equations are provided to give an accurate calculation of the gradient of the objective function. The method is illustrated using in silico data from the FitzHugh-Nagumo model for excitable media and the Lotka-Volterra predator-prey system. The proposed method performs well on the models considered, and is able to regularize the objective function in both models. This leads to parameter estimation problems with fewer local minima which can be solved by efficient gradient-based methods.

Bayesian model comparison and parameter inference in systems biology using nested sampling

Inferring parameters for models of biological processes is a current challenge in systems biology, as is the related problem of comparing competing models that explain the data. In this work we apply Skilling's nested sampling to address both of these problems. Nested sampling is a Bayesian method for exploring parameter space that transforms a multi-dimensional integral to a 1D integration over likelihood space. This approach focuses on the computation of the marginal likelihood or evidence. The ratio of evidences of different models leads to the Bayes factor, which can be used for model comparison. We demonstrate how nested sampling can be used to reverse-engineer a system's behaviour whilst accounting for the uncertainty in the results. The effect of missing initial conditions of the variables as well as unknown parameters is investigated. We show how the evidence and the model ranking can change as a function of the available data. Furthermore, the addition of data from extra variables of the system can deliver more information for model comparison than increasing the data from one variable, thus providing a basis for experimental design.

Prediction stability in a data-based, mechanistic model of σF regulation during sporulation in Bacillus subtilis

Mathematical modeling of biological networks can help to integrate a large body of information into a consistent framework, which can then be used to arrive at novel mechanistic insight and predictions. We have previously developed a detailed, mechanistic model for the regulation of σ ^F during sporulation in Bacillus subtilis. The model was based on a wide range of quantitative data, and once fitted to the data, the model made predictions that could be confirmed in experiments. However, the analysis was based on a single optimal parameter set. We wondered whether the predictions of the model would be stable for all optimal parameter sets. To that end we conducted a global parameter screen within the physiological parameter ranges. The screening approach allowed us to identify sensitive and sloppy parameters, and highlighted further required datasets during the optimization. Eventually, all parameter sets that reproduced all available data predicted the physiological situation correctly.

Identification of parameter correlations for parameter estimation in dynamic biological models

Background

One of the challenging tasks in systems biology is parameter estimation in nonlinear dynamic models. A biological model usually contains a large number of correlated parameters leading to non-identifiability problems. Although many approaches have been developed to address both structural and practical non-identifiability problems, very few studies have been made to systematically investigate parameter correlations.

Results

In this study we present an approach that is able to identify both pairwise parameter correlations and higher order interrelationships among parameters in nonlinear dynamic models. Correlations are interpreted as surfaces in the subspaces of correlated parameters. Based on the correlation information obtained in this way both structural and practical non-identifiability can be clarified. Moreover, it can be concluded from the correlation analysis that a minimum number of data sets with different inputs for experimental design are needed to relieve the parameter correlations, which corresponds to the maximum number of correlated parameters among the correlation groups.

Conclusions

The information of pairwise and higher order interrelationships among parameters in biological models gives a deeper insight into the cause of non-identifiability problems. The result of our correlation analysis provides a necessary condition for experimental design in order to acquire suitable measurement data for unique parameter estimation.

Systematic parameter estimation in data-rich environments for cell signalling dynamics

Motivation: Computational models of biological signalling networks, based on ordinary differential equations (ODEs), have generated many insights into cellular dynamics, but the model-building process typically requires estimating rate parameters based on experimentally observed concentrations. New proteomic methods can measure concentrations for all molecular species in a pathway; this creates a new opportunity to decompose the optimization of rate parameters.

Results: In contrast with conventional parameter estimation methods that minimize the disagreement between simulated and observed concentrations, the SPEDRE method fits spline curves through observed concentration points, estimates derivatives and then matches the derivatives to the production and consumption of each species. This reformulation of the problem permits an extreme decomposition of the high-dimensional optimization into a product of low-dimensional factors, each factor enforcing the equality of one ODE at one time slice. Coarsely discretized solutions to the factors can be computed systematically. Then the discrete solutions are combined using loopy belief propagation, and refined using local optimization. SPEDRE has unique asymptotic behaviour with runtime polynomial in the number of molecules and timepoints, but exponential in the degree of the biochemical network. SPEDRE performance is comparatively evaluated on a novel model of Akt activation dynamics including redox-mediated inactivation of PTEN (phosphatase and tensin homologue).

Availability and implementation: Web service, software and supplementary information are available at www.LtkLab.org/SPEDRE

Supplementary information:Supplementary data are available at Bioinformatics online.

Contact: LisaTK@nus.edu.sg

Parameter estimation using meta-heuristics in systems biology: a comprehensive review

This paper gives a comprehensive review of the application of meta-heuristics to optimization problems in systems biology, mainly focussing on the parameter estimation problem (also called the inverse problem or model calibration). It is intended for either the system biologist who wishes to learn more about the various optimization techniques available and/or the meta-heuristic optimizer who is interested in applying such techniques to problems in systems biology. First, the parameter estimation problems emerging from different areas of systems biology are described from the point of view of machine learning. Brief descriptions of various meta-heuristics developed for these problems follow, along with outlines of their advantages and disadvantages. Several important issues in applying meta-heuristics to the systems biology modelling problem are addressed, including the reliability and identifiability of model parameters, optimal design of experiments, and so on. Finally, we highlight some possible future research directions in this field.

Analyzing and constraining signaling networks: parameter estimation for the user

The behavior of most dynamical models not only depends on the wiring but also on the kind and strength of interactions which are reflected in the parameter values of the model. The predictive value of mathematical models therefore critically hinges on the quality of the parameter estimates. Constraining a dynamical model by an appropriate parameterization follows a 3-step process. In an initial step, it is important to evaluate the sensitivity of the parameters of the model with respect to the model output of interest. This analysis points at the identifiability of model parameters and can guide the design of experiments. In the second step, the actual fitting needs to be carried out. This step requires special care as, on the one hand, noisy as well as partial observations can corrupt the identification of system parameters. On the other hand, the solution of the dynamical system usually depends in a highly nonlinear fashion on its parameters and, as a consequence, parameter estimation procedures get easily trapped in local optima. Therefore any useful parameter estimation procedure has to be robust and efficient with respect to both challenges. In the final step, it is important to access the validity of the optimized model. A number of reviews have been published on the subject. A good, nontechnical overview is provided by Jaqaman and Danuser (Nat Rev Mol Cell Biol 7(11):813-819, 2006) and a classical introduction, focussing on the algorithmic side, is given in Press (Numerical recipes: The art of scientific computing, Cambridge University Press, 3rd edn., 2007, Chapters 10 and 15). We will focus on the practical issues related to parameter estimation and use a model of the TGFβ-signaling pathway as an educative example. Corresponding parameter estimation software and models based on MATLAB code can be downloaded from the authors's web page ( http://www.bsse.ethz.ch/cobi ).

A parameter optimization method to determine ski stiffness properties from ski deformation data

The deformation of skis and the contact pressure between skis and snow are crucial factors for carved turns in alpine skiing. The purpose of the current study was to develop and to evaluate an optimization method to determine the bending and torsional stiffness that lead to a given bending and torsional deflection of the ski. Euler-Bernoulli beam theory and classical torsion theory were applied to model the deformation of the ski. Bending and torsional stiffness were approximated as linear combinations of B-splines. To compute the unknown coefficients, a parameter optimization problem was formulated and successfully solved by multiple shooting and least squares data fitting. The proposed optimization method was evaluated based on ski stiffness data and ski deformation data taken from a recently published simulation study. The ski deformation data were used as input data to the optimization method. The optimization method was capable of successfully reproducing the shape of the original bending and torsional stiffness data of the ski with a root mean square error below 1 N m2. In conclusion, the proposed computational method offers the possibility to calculate ski stiffness properties with respect to a given ski deformation.

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.

Fixed point characterization of biological networks with complex graph topology

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.

Physiological modeling of isoprene dynamics in exhaled breath

Human breath contains a myriad of endogenous volatile organic compounds (VOCs) which are reflective of ongoing metabolic or physiological processes. While research into the diagnostic potential and general medical relevance of these trace gases is conducted on a considerable scale, little focus has been given so far to a sound analysis of the quantitative relationships between breath levels and the underlying systemic concentrations. This paper is devoted to a thorough modeling study of the end-tidal breath dynamics associated with isoprene, which serves as a paradigmatic example for the class of low-soluble, blood-borne VOCs. Real-time measurements of exhaled breath under an ergometer challenge reveal characteristic changes of isoprene output in response to variations in ventilation and perfusion. Here, a valid compartmental description of these profiles is developed. By comparison with experimental data it is inferred that the major part of breath isoprene variability during exercise conditions can be attributed to an increased fractional perfusion of potential storage and production sites, leading to higher levels of mixed venous blood concentrations at the onset of physical activity. In this context, various lines of supportive evidence for an extrahepatic tissue source of isoprene are presented. Our model is a first step towards new guidelines for the breath gas analysis of isoprene and is expected to aid further investigations regarding the exhalation, storage, transport and biotransformation processes associated with this important compound.

Parameter inference and model selection in signaling pathway models

To support and guide an extensive experimental research into systems biology of signaling pathways, increasingly more mechanistic models are being developed with hopes of gaining further insight into biological processes. In order to analyze these models, computational and statistical techniques are needed to estimate the unknown kinetic parameters. This chapter reviews methods from frequentist and Bayesian statistics for estimation of parameters and for choosing which model is best for modeling the underlying system. Approximate Bayesian computation techniques are introduced and employed to explore different hypothesis about the JAK-STAT signaling pathway.

Identifikation biochemischer Reaktionsnetzwerke: Ein beobachterbasierter Ansatz (Identification of Biochemical Reaction Networks: An Observer based Approach)

Die Basis für viele Untersuchungen im Bereich der Systembiologie bilden mathematische Modelle. Diese hängen jedoch oftmals von einer Vielzahl unbekannter oder nur ungenau bekannter Reaktionsparameter ab. Die direkte Messung dieser Parameter basierend auf Einzelreaktionen ist aufgrund der damit verbundenen hohen Kosten, des notwendigen Zeitaufwands zur Durchführung der Experimente oder der sich ergebenden Komplexität durch die Vielzahl an Einzelreaktionen nicht möglich. Deshalb müssen diese Parameter aus indirekten Messungen an der realen Zelle, zum Beispiel aus Zeitreiheninformationen, gewonnen werden. Existierende Parameteridentifikationsverfahren berücksichtigen meistens nicht die speziellen Strukturen biochemischer Reaktionsnetzwerke. Oftmals liefern sie keine zuverlässigen Parameterschätzungen oder Garantien über die Eindeutigkeit der gefundenen Parameter. Im Rahmen dieser Arbeit beschreiben wir ein Parameterschätzverfahren, das einen Zustandsbeobachter verwendet und es erlaubt, die auftretenden Strukturen direkt zu berücksichtigen. Unter der Annahme, dass die auftretenden Reaktionsmechanismen durch polynomiale oder rationale Funktionen beschrieben werden können, stellen wir ein dreistufiges Schätzverfahren vor. In einem ersten Schritt wird ein erweitertes Modell hergeleitet, mit dem die Parameterschätzung von der notwendigen Schätzung der unbekannten Zustände getrennt werden kann. In einem zweiten Schritt wird der erweiterterte Zustand mittels eines geeigneten Beobachters rekonstruktiert. Basierend auf den hieraus erhaltenen Zuständen werden die Reaktionsparameter im letzten Schritt geschätzt. Die Vor- und Nachteile des vorgeschlagenen Verfahrens werden an einem vereinfachten Modell für den zirkadischen Rythmus verdeutlicht.

Estimating parameters for generalized mass action models with connectivity information

BACKGROUND

Determining the parameters of a mathematical model from quantitative measurements is the main bottleneck of modelling biological systems. Parameter values can be estimated from steady-state data or from dynamic data. The nature of suitable data for these two types of estimation is rather different. For instance, estimations of parameter values in pathway models, such as kinetic orders, rate constants, flux control coefficients or elasticities, from steady-state data are generally based on experiments that measure how a biochemical system responds to small perturbations around the steady state. In contrast, parameter estimation from dynamic data requires time series measurements for all dependent variables. Almost no literature has so far discussed the combined use of both steady-state and transient data for estimating parameter values of biochemical systems.

RESULTS

In this study we introduce a constrained optimization method for estimating parameter values of biochemical pathway models using steady-state information and transient measurements. The constraints are derived from the flux connectivity relationships of the system at the steady state. Two case studies demonstrate the estimation results with and without flux connectivity constraints. The unconstrained optimal estimates from dynamic data may fit the experiments well, but they do not necessarily maintain the connectivity relationships. As a consequence, individual fluxes may be misrepresented, which may cause problems in later extrapolations. By contrast, the constrained estimation accounting for flux connectivity information reduces this misrepresentation and thereby yields improved model parameters.

CONCLUSION

The method combines transient metabolic profiles and steady-state information and leads to the formulation of an inverse parameter estimation task as a constrained optimization problem. Parameter estimation and model selection are simultaneously carried out on the constrained optimization problem and yield realistic model parameters that are more likely to hold up in extrapolations with the model.

Identification of neutral biochemical network models from time series data

BACKGROUND

The major difficulty in modeling biological systems from multivariate time series is the identification of parameter sets that endow a model with dynamical behaviors sufficiently similar to the experimental data. Directly related to this parameter estimation issue is the task of identifying the structure and regulation of ill-characterized systems. Both tasks are simplified if the mathematical model is canonical, i.e., if it is constructed according to strict guidelines.

RESULTS

In this report, we propose a method for the identification of admissible parameter sets of canonical S-systems from biological time series. The method is based on a Monte Carlo process that is combined with an improved version of our previous parameter optimization algorithm. The method maps the parameter space into the network space, which characterizes the connectivity among components, by creating an ensemble of decoupled S-system models that imitate the dynamical behavior of the time series with sufficient accuracy. The concept of sloppiness is revisited in the context of these S-system models with an exploration not only of different parameter sets that produce similar dynamical behaviors but also different network topologies that yield dynamical similarity.

CONCLUSION

The proposed parameter estimation methodology was applied to actual time series data from the glycolytic pathway of the bacterium Lactococcus lactis and led to ensembles of models with different network topologies. In parallel, the parameter optimization algorithm was applied to the same dynamical data upon imposing a pre-specified network topology derived from prior biological knowledge, and the results from both strategies were compared. The results suggest that the proposed method may serve as a powerful exploration tool for testing hypotheses and the design of new experiments.

Identification of neutral biochemical network models from time series data

Background

The major difficulty in modeling biological systems from multivariate time series is the identification of parameter sets that endow a model with dynamical behaviors sufficiently similar to the experimental data. Directly related to this parameter estimation issue is the task of identifying the structure and regulation of ill-characterized systems. Both tasks are simplified if the mathematical model is canonical, i.e., if it is constructed according to strict guidelines.

Results

In this report, we propose a method for the identification of admissible parameter sets of canonical S-systems from biological time series. The method is based on a Monte Carlo process that is combined with an improved version of our previous parameter optimization algorithm. The method maps the parameter space into the network space, which characterizes the connectivity among components, by creating an ensemble of decoupled S-system models that imitate the dynamical behavior of the time series with sufficient accuracy. The concept of sloppiness is revisited in the context of these S-system models with an exploration not only of different parameter sets that produce similar dynamical behaviors but also different network topologies that yield dynamical similarity.

Conclusion

The proposed parameter estimation methodology was applied to actual time series data from the glycolytic pathway of the bacterium Lactococcus lactis and led to ensembles of models with different network topologies. In parallel, the parameter optimization algorithm was applied to the same dynamical data upon imposing a pre-specified network topology derived from prior biological knowledge, and the results from both strategies were compared. The results suggest that the proposed method may serve as a powerful exploration tool for testing hypotheses and the design of new experiments.

A perturbation-based estimate algorithm for parameters of coupled ordinary differential equations, applications from chemical reactions to metabolic dynamics

Conversion of complex phenomena in medicine, pharmaceutical and systems biology fields to a system of ordinary differential equations (ODEs) and identification of parameters from experimental data and theoretical model equations can be treated as a computational engine to arrive at the best solution for chemical reactions, biochemical metabolic and intracellular pathways. Particularly, to gain insight into the pathophysiology of diabetes's metabolism in our current clinical studies, glucose kinetics and insulin secretion can be assessed by the ODE model. Parameter estimation is usually performed by minimizing a cost function which quantifies the difference between theoretical model predictions and experimental measurements. This paper explores how the numerical method and iteration program are developed to search ODE's parameters using the perturbation method, instead of the Gauss-Newton or Levenberg-Marquardt method. Several interesting applications, including Lotka-Volterra chemical reaction system, Lorenz chaos, dynamics of tetracycline hydrochloride concentration, and Bergman's Minimal Model for glucose kinetics are illustrated.

Positive- and negative-feedback regulations coordinate the dynamic behavior of the Ras-Raf-MEK-ERK signal transduction pathway

The Ras-Raf-MEK-ERK pathway (or ERK pathway) is an important signal transduction system involved in the control of cell proliferation, survival and differentiation. However, the dynamic regulation of the pathway by positive- and negative-feedback mechanisms, in particular the functional role of Raf kinase inhibitor protein (RKIP) are still incompletely understood. RKIP is a physiological endogenous inhibitor of MEK phosphorylation by Raf kinases, but also participates in a positive-feedback loop in which ERK can inactivate RKIP. The aim of this study was to elucidate the hidden dynamics of these feedback mechanisms and to identify the functional role of RKIP through combined efforts of biochemical experiments and in silico simulations based on an experimentally validated mathematical model. We show that the negative-feedback loop from ERK to SOS plays a crucial role in generating an oscillatory behavior of ERK activity. The positive-feedback loop in which ERK functionally inactivates RKIP also enhances the oscillatory activation pattern of ERK. However, RKIP itself has an important role in inducing a switch-like behavior of MEK activity. When overexpressed, RKIP also causes delayed and reduced responses of ERK. Thus, positive- and negative-feedback loops and RKIP work together to shape the response pattern and dynamical characteristics of the ERK pathway.

Parameter estimation and optimal experimental design

Mathematical models are central in systems biology and provide new ways to understand the function of biological systems, helping in the generation of novel and testable hypotheses, and supporting a rational framework for possible ways of intervention, like in e.g. genetic engineering, drug development or treatment of diseases. Since the amount and quality of experimental 'omics' data continue to increase rapidly, there is great need for methods for proper model building which can handle this complexity. In the present chapter we review two key steps of the model building process, namely parameter estimation (model calibration) and optimal experimental design. Parameter estimation aims to find the unknown parameters of the model which give the best fit to a set of experimental data. Optimal experimental design aims to devise the dynamic experiments which provide the maximum information content for subsequent non-linear model identification, estimation and/or discrimination. We place emphasis on the need for robust global optimization methods for proper solution of these problems, and we present a motivating example considering a cell signalling model.

The impact of time delays on the robustness of biological oscillators and the effect of bifurcations on the inverse problem

Differential equation models for biological oscillators are often not robust with respect to parameter variations. They are based on chemical reaction kinetics, and solutions typically converge to a fixed point. This behavior is in contrast to real biological oscillators, which work reliably under varying conditions. Moreover, it complicates network inference from time series data. This paper investigates differential equation models for biological oscillators from two perspectives. First, we investigate the effect of time delays on the robustness of these oscillator models. In particular, we provide sufficient conditions for a time delay to cause oscillations by destabilizing a fixed point in two-dimensional systems. Moreover, we show that the inclusion of a time delay also stabilizes oscillating behavior in this way in larger networks. The second part focuses on the inverse problem of estimating model parameters from time series data. Bifurcations are related to nonsmoothness and multiple local minima of the objective function.

Computational procedures for optimal experimental design in biological systems

Mathematical models of complex biological systems, such as metabolic or cell-signalling pathways, usually consist of sets of nonlinear ordinary differential equations which depend on several non-measurable parameters that can be hopefully estimated by fitting the model to experimental data. However, the success of this fitting is largely conditioned by the quantity and quality of data. Optimal experimental design (OED) aims to design the scheme of actuations and measurements which will result in data sets with the maximum amount and/or quality of information for the subsequent model calibration. New methods and computational procedures for OED in the context of biological systems are presented. The OED problem is formulated as a general dynamic optimisation problem where the time-dependent stimuli profiles, the location of sampling times, the duration of the experiments and the initial conditions are regarded as design variables. Its solution is approached using the control vector parameterisation method. Since the resultant nonlinear optimisation problem is in most of the cases non-convex, the use of a robust global nonlinear programming solver is proposed. For the sake of comparing among different experimental schemes, a Monte-Carlo-based identifiability analysis is then suggested. The applicability and advantages of the proposed techniques are illustrated by considering an example related to a cell-signalling pathway.

Hybrid optimization method with general switching strategy for parameter estimation

BACKGROUND

Modeling and simulation of cellular signaling and metabolic pathways as networks of biochemical reactions yields sets of non-linear ordinary differential equations. These models usually depend on several parameters and initial conditions. If these parameters are unknown, results from simulation studies can be misleading. Such a scenario can be avoided by fitting the model to experimental data before analyzing the system. This involves parameter estimation which is usually performed by minimizing a cost function which quantifies the difference between model predictions and measurements. Mathematically, this is formulated as a non-linear optimization problem which often results to be multi-modal (non-convex), rendering local optimization methods detrimental.

RESULTS

In this work we propose a new hybrid global method, based on the combination of an evolutionary search strategy with a local multiple-shooting approach, which offers a reliable and efficient alternative for the solution of large scale parameter estimation problems.

CONCLUSION

The presented new hybrid strategy offers two main advantages over previous approaches: First, it is equipped with a switching strategy which allows the systematic determination of the transition from the local to global search. This avoids computationally expensive tests in advance. Second, using multiple-shooting as the local search procedure reduces the multi-modality of the non-linear optimization problem significantly. Because multiple-shooting avoids possible spurious solutions in the vicinity of the global optimum it often outperforms the frequently used initial value approach (single-shooting). Thereby, the use of multiple-shooting yields an enhanced robustness of the hybrid approach.

Data-based mathematical modeling of vectorial transport across double-transfected polarized cells

Vectorial transport of endogenous small molecules, toxins, and drugs across polarized epithelial cells contributes to their half-life in the organism and to detoxification. To study vectorial transport in a quantitative manner, an in vitro model was used that includes polarized MDCKII cells stably expressing the recombinant human uptake transporter OATP1B3 in their basolateral membrane and the recombinant ATP-driven efflux pump ABCC2 in their apical membrane. These double-transfected cells enabled mathematical modeling of the vectorial transport of the anionic prototype substance bromosulfophthalein (BSP) that has frequently been used to examine hepatobiliary transport. Time-dependent analyses of (3)H-labeled BSP in the basolateral, intracellular, and apical compartments of cells cultured on filter membranes and efflux experiments in cells preloaded with BSP were performed. A mathematical model was fitted to the experimental data. Data-based modeling was optimized by including endogenous transport processes in addition to the recombinant transport proteins. The predominant contributions to the overall vectorial transport of BSP were mediated by OATP1B3 (44%) and ABCC2 (28%). Model comparison predicted a previously unrecognized endogenous basolateral efflux process as a negative contribution to total vectorial transport, amounting to 19%, which is in line with the detection of the basolateral efflux pump Abcc4 in MDCKII cells. Rate-determining steps in the vectorial transport were identified by calculating control coefficients. Data-based mathematical modeling of vectorial transport of BSP as a model substance resulted in a quantitative description of this process and its components. The same systems biology approach may be applied to other cellular systems and to different substances.