ABSTRACT In the field of optimisation, the accurate configuration of parameters in metaheuristic algorithms is a critical yet often arduous task that significantly impacts the efficiency and efficacy of the search process. This study was motivated by the need to address the inefficiencies and limitations associated with conventional methods of parameter configuration, which typically involve manual, trial‐and‐error approaches. These traditional methods can lead to suboptimal performance and increased computational overhead. To tackle these challenges, this study introduces a novel adaptive parameter control strategy for parallel island‐based metaheuristics, with a particular emphasis on the ant colony optimisation (ACO) algorithm. Our research process involved extensive experimentation to evaluate the effectiveness of this adaptive strategy. We conducted a series of tests to enable real‐time adjustment of key parameters based on the performance of ACO colonies, thereby enhancing both exploration and exploitation capabilities. The results indicate that the adaptive strategy consistently outperforms offline manual and automated tuning configurations, particularly in larger and more complex problem instances, providing a more efficient solution for parameter optimisation in metaheuristics. These findings highlight the potential of dynamic parameter control to reduce dependency on expert knowledge and manual tuning while improving algorithmic performance.

Microbial communities, complex ecological networks crucial for human and planetary health, remain poorly understood in terms of the quantitative principles governing their composition, assembly, and function. Dynamic modeling using ordinary differential equations (ODEs) is a powerful framework for understanding and predicting microbiome behaviors. However, developing reliable ODE models is severely hampered by their nonlinear nature and the presence of significant challenges, particularly critical issues related to identifiability. Here, we address the identification problem in dynamic microbial community models by proposing an integrated methodology to tackle key challenges. Focusing on nonlinear ODE-based models, we examine four critical pitfalls: identifiability issues (structural and practical), unstable dynamics (potentially leading to numerical blow-up), underfitting (convergence to suboptimal solutions), and overfitting (fitting noise rather than signal). These pitfalls yield unreliable parameter estimates, unrealistic model behavior, and poor generalization. Our study presents a comprehensive workflow incorporating structural and practical identifiability analysis, robust global optimization for calibration, stability checks, and rigorous predictive power assessment. The methodology's effectiveness and versatility in mitigating these pitfalls are demonstrated through case studies of increasing complexity, paving the way for more reliable and mechanistically insightful models of microbial communities.

Uncertainty quantification (UQ) is the process of systematically determining and characterizing the degree of confidence in computational model predictions. In systems biology, and particularly with dynamic models, UQ is critical due to the nonlinearities and parameter sensitivities that influence the behavior of complex biological systems. Addressing these issues through robust UQ enables a deeper understanding of system dynamics and more reliable extrapolation beyond observed conditions. Many state-of-the-art UQ approaches in this field are grounded in Bayesian statistical methods. While these frameworks naturally incorporate uncertainty quantification, they often require the specification of parameter distributions as priors and may impose parametric assumptions that do not always reflect biological reality. Additionally, Bayesian methods can be computationally expensive, posing significant challenges when dealing with large-scale models and seeking rapid, reliable uncertainty calibration. As an alternative, we propose using conformal predictions methods and introduce two novel algorithms designed for dynamic biological systems. These approaches can provide non-asymptotic guarantees, improving robustness and scalability across various applications, even when the predictive models are misspecified. Through several illustrative scenarios, we demonstrate that these conformal algorithms can serve as powerful complements-or even alternatives-to conventional Bayesian methods, delivering effective uncertainty quantification for predictive tasks in systems biology.

Abstract Binary combinatorial optimization plays a crucial role in various scientific and engineering fields. While deterministic approaches have traditionally been used to solve these problems, stochastic methods, particularly metaheuristics, have gained popularity in recent years for efficiently handling large problem instances. Ant Colony Optimization (ACO) is among the most successful metaheuristics and is frequently employed in non‐binary combinatorial problems due to its adaptability. Although for binary combinatorial problems ACO can suffer from issues such as rapid convergence to local minima, its eminently parallel structure means that it can be exploited to solve large and complex problems also in this field. In order to provide a versatile ACO implementation that achieves competitive results across a wide range of binary combinatorial optimization problems, we introduce a parallel multicolony strategy with an improved cooperation scheme to maintain search diversity. We evaluate our proposal (Binary Parallel Cooperative ACO, BiPCACO) using a comprehensive benchmark framework, showcasing its performance and, most importantly, its flexibility as a successful all‐purpose solver for binary combinatorial problems.

The mammalian cell cycle is regulated by a well-studied but complex biochemical reaction system. Computational models provide a particularly systematic and systemic description of the mechanisms governing mammalian cell cycle control. By combining both state-of-the-art multiplexed experimental methods and powerful computational tools, this work aims at improving on these models along four dimensions: model structure, validation data, validation methodology and model reusability. We developed a comprehensive model structure of the full cell cycle that qualitatively explains the behaviour of human retinal pigment epithelial-1 cells. To estimate the model parameters, time courses of eight cell cycle regulators in two compartments were reconstructed from single cell snapshot measurements. After optimisation with a parallel global optimisation metaheuristic we obtained excellent agreements between simulations and measurements. The PEtab specification of the optimisation problem facilitates reuse of model, data and/or optimisation results. Future perturbation experiments will improve parameter identifiability and allow for testing model predictive power. Such a predictive model may aid in drug discovery for cell cycle-related disorders.

The Virtual Epileptic Patient (VEP) refers to a computer-based representation of a patient with epilepsy that combines personalized anatomical data with dynamical models of abnormal brain activities. It is capable of generating spatio-temporal seizure patterns that resemble those recorded with invasive methods such as stereoelectro EEG data, allowing for the evaluation of clinical hypotheses before planning surgery. This study highlights the effectiveness of calibrating VEP models using a global optimization approach. The approach utilizes SaCeSS, a cooperative metaheuristic algorithm capable of parallel computation, to yield high-quality solutions without requiring excessive computational time. Through extensive benchmarking on synthetic data, our proposal successfully solved a set of different configurations of VEP models, demonstrating better scalability and superior performance against other parallel solvers. These results were further enhanced using a Bayesian optimization framework for hyperparameter tuning, with significant gains in terms of both accuracy and computational cost. Additionally, we added a scalable uncertainty quantification phase after model calibration, and used it to assess the variability in estimated parameters across different problems. Overall, this study has the potential to improve the estimation of pathological brain areas in drug-resistant epilepsy, thereby to inform the clinical decision-making process.

MOTIVATION: Dynamic mechanistic modelling in systems biology has been hampered by the complexity and variability associated with the underlying interactions, and by uncertain and sparse experimental measurements. Ensemble modelling, a concept initially developed in statistical mechanics, has been introduced in biological applications with the aim of mitigating those issues. Ensemble modelling uses a collection of different models compatible with the observed data to describe the phenomena of interest. However, since systems biology models often suffer from a lack of identifiability and observability, ensembles of models are particularly unreliable when predicting non-observable states. RESULTS: We present a strategy to assess and improve the reliability of a class of model ensembles. In particular, we consider kinetic models described using ordinary differential equations with a fixed structure. Our approach builds an ensemble with a selection of the parameter vectors found when performing parameter estimation with a global optimization metaheuristic. This technique enforces diversity during the sampling of parameter space and it can quantify the uncertainty in the predictions of state trajectories. We couple this strategy with structural identifiability and observability analysis, and when these tests detect possible prediction issues we obtain model reparameterizations that surmount them. The end result is an ensemble of models with the ability to predict the internal dynamics of a biological process. We demonstrate our approach with models of glucose regulation, cell division, circadian oscillations and the JAK-STAT signalling pathway. AVAILABILITY AND IMPLEMENTATION: The code that implements the methodology and reproduces the results is available at https://doi.org/10.5281/zenodo.6782638. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

MOTIVATION: One of the main causes hampering predictability during the model identification and automated design of gene circuits in synthetic biology is the effect of molecular noise. Stochasticity may significantly impact the dynamics and function of gene circuits, specially in bacteria and yeast due to low mRNA copy numbers. Standard stochastic simulation methods are too computationally costly in realistic scenarios to be applied to optimization-based design or parameter estimation. RESULTS: In this work, we present IDESS (Identification and automated Design of Stochastic gene circuitS), a software toolbox for optimization-based design and model identification of gene regulatory circuits in the stochastic regime. This software incorporates an efficient approximation of the Chemical Master Equation as well as a stochastic simulation algorithm-both with GPU and CPU implementations-combined with global optimization algorithms capable of solving Mixed Integer Nonlinear Programming problems. The toolbox efficiently addresses two types of problems relevant in systems and synthetic biology: the automated design of stochastic synthetic gene circuits, and the parameter estimation for model identification of stochastic gene regulatory networks. AVAILABILITY AND IMPLEMENTATION: IDESS runs under the MATLAB environment and it is available under GPLv3 license at https://doi.org/10.5281/zenodo.7788692.

Mechanistic dynamical models allow us to study the behavior of complex biological systems. They can provide an objective and quantitative understanding that would be difficult to achieve through other means. However, the systematic development of these models is a non-trivial exercise and an open problem in computational biology. Currently, many research efforts are focused on model discovery, i.e. automating the development of interpretable models from data. One of the main frameworks is sparse regression, where the sparse identification of nonlinear dynamics (SINDy) algorithm and its variants have enjoyed great success. SINDy-PI is an extension which allows the discovery of rational nonlinear terms, thus enabling the identification of kinetic functions common in biochemical networks, such as Michaelis-Menten. SINDy-PI also pays special attention to the recovery of parsimonious models (Occam's razor). Here we focus on biological models composed of sets of deterministic nonlinear ordinary differential equations. We present a methodology that, combined with SINDy-PI, allows the automatic discovery of structurally identifiable and observable models which are also mechanistically interpretable. The lack of structural identifiability and observability makes it impossible to uniquely infer parameter and state variables, which can compromise the usefulness of a model by distorting its mechanistic significance and hampering its ability to produce biological insights. We illustrate the performance of our method with six case studies. We find that, despite enforcing sparsity, SINDy-PI sometimes yields models that are unidentifiable. In these cases we show how our method transforms their equations in order to obtain a structurally identifiable and observable model which is also interpretable.

This work presents a methodology to evaluate the bioenergetic feasibility of alternative metabolic pathways for a given microbial conversion, optimising their energy yield and driving forces as a function of the concentration of metabolic intermediates. The tool, based on thermodynamic principles and multi-objective optimisation, accounts for pathway variants in terms of different electron carriers, as well as energy conservation (proton translocating) reactions within the pathway. The method also accommodates other constraints, some of them non-linear, such as the balance of conserved moieties. The approach involves the transformation of the maximum energy yield problem into a multi-objective mixed-integer linear optimisation problem which is then subsequently solved using the epsilon-constraint method, highlighting the trade-off between yield and rate in metabolic reactions. The methodology is applied to analyse several pathway alternatives occurring during propionate oxidation in anaerobic fermentation processes, as well as to the reverse TCA cycle pathway occurring during autotrophic microbial CO2 fixation. The results obtained using the developed methodology match previously reported literature and bring about insights into the studied pathways.

Mechanistic dynamic models have become an essential tool for understanding biomolecular networks and other biological systems. Biochemical stochasticity can be extremely important in some situations, e.g., at the single-cell level where there is a low copy number of the species involved. In these scenarios, deterministic models are not suitable to characterize the dynamics, so stochastic dynamic models should be considered. Here, we address the challenging problem of parameter estimation in stochastic dynamic models. Despite recent advances, this area is considerably less mature than its deterministic counterpart. We present a novel strategy based on two components: (i) global optimization via a hybrid stochastic-deterministic approach, and (ii) stochastic simulation techniques tailored to the sparsity of the available experimental data. Regarding the latter, for cases of dense population data we make use of a novel approach using a Partial Integro-Differential Equation (PIDE) model solved using a semilagrangian method. In order to further speed up the simulations, we also present efficient parallel implementations for multi-core CPUs and also for graphical processing units (GPUs). Importantly, whereas SDE and Fokker Planck approximations of the Chemical Master Equation (CME) apply when the reactant populations are sufficiently large, the PIDE approximation to the CME is valid for very low copy numbers, and therefore they enable us to tackle parameter estimation for systems with large intrinsic molecular noise, (highly stochastic regimes far from the thermodynamic limit). We test our strategy with four challenging problems: a Lotka-Volterra system, a polarization system in S. cerevisiae, a genetic toggle switch, and a genetic circadian oscillator. Our method could successfully solve these problems in very reasonable computation times (often a few minutes for the first two problems) using standard low-cost hardware, showing very significant speedups with respect to recent alternative methods. The code used to obtain the results reported here is available at https://doi.org/10.5281/zenodo.5195408.

Biological processes are often modelled using ordinary differential equations. The unknown parameters of these models are estimated by optimizing the fit of model simulation and experimental data. The resulting parameter estimates inevitably possess some degree of uncertainty. In practical applications it is important to quantify these parameter uncertainties as well as the resulting prediction uncertainty, which are uncertainties of potentially time-dependent model characteristics. Unfortunately, estimating prediction uncertainties accurately is nontrivial, due to the nonlinear dependence of model characteristics on parameters. While a number of numerical approaches have been proposed for this task, their strengths and weaknesses have not been systematically assessed yet. To fill this knowledge gap, we apply four state of the art methods for uncertainty quantification to four case studies of different computational complexities. This reveals the trade-offs between their applicability and their statistical interpretability. Our results provide guidelines for choosing the most appropriate technique for a given problem and applying it successfully.

Ordinary differential equation models are nowadays widely used for the mechanistic description of biological processes and their temporal evolution. These models typically have many unknown and nonmeasurable parameters, which have to be determined by fitting the model to experimental data. In order to perform this task, known as parameter estimation or model calibration, the modeller faces challenges such as poor parameter identifiability, lack of sufficiently informative experimental data and the existence of local minima in the objective function landscape. These issues tend to worsen with larger model sizes, increasing the computational complexity and the number of unknown parameters. An incorrectly calibrated model is problematic because it may result in inaccurate predictions and misleading conclusions. For nonexpert users, there are a large number of potential pitfalls. Here, we provide a protocol that guides the user through all the steps involved in the calibration of dynamic models. We illustrate the methodology with two models and provide all the code required to reproduce the results and perform the same analysis on new models. Our protocol provides practitioners and researchers in biological modelling with a one-stop guide that is at the same time compact and sufficiently comprehensive to cover all aspects of the problem.

Acquiring a functional comprehension of the deregulation of cell signaling networks in disease allows progress in the development of new therapies and drugs. Computational models are becoming increasingly popular as a systematic tool to analyze the functioning of complex biochemical networks, such as those involved in cell signaling. CellNOpt is a framework to build predictive logic-based models of signaling pathways by training a prior knowledge network to biochemical data obtained from perturbation experiments. This training can be formulated as an optimization problem that can be solved using metaheuristics. However, the genetic algorithm used so far in CellNOpt presents limitations in terms of execution time and quality of solutions when applied to large instances. Thus, in order to overcome those issues, in this paper we propose the use of a method based on ant colony optimization, adapted to the problem at hand and parallelized using a hybrid approach. The performance of this novel method is illustrated with several challenging benchmark problems in the study of new therapies for liver cancer.

Reproducibility and reusability of the results of data-based modeling studies are essential. Yet, there has been-so far-no broadly supported format for the specification of parameter estimation problems in systems biology. Here, we introduce PEtab, a format which facilitates the specification of parameter estimation problems using Systems Biology Markup Language (SBML) models and a set of tab-separated value files describing the observation model and experimental data as well as parameters to be estimated. We already implemented PEtab support into eight well-established model simulation and parameter estimation toolboxes with hundreds of users in total. We provide a Python library for validation and modification of a PEtab problem and currently 20 example parameter estimation problems based on recent studies.

BACKGROUND: Optimality principles have been used to explain the structure and behavior of living matter at different levels of organization, from basic phenomena at the molecular level, up to complex dynamics in whole populations. Most of these studies have assumed a single-criteria approach. Such optimality principles have been justified from an evolutionary perspective. In the context of the cell, previous studies have shown how dynamics of gene expression in small metabolic models can be explained assuming that cells have developed optimal adaptation strategies. Most of these works have considered rather simplified representations, such as small linear pathways, or reduced networks with a single branching point, and a single objective for the optimality criteria. RESULTS: Here we consider the extension of this approach to more realistic scenarios, i.e. biochemical pathways of arbitrary size and structure. We first show that exploiting optimality principles for these networks poses great challenges due to the complexity of the associated optimal control problems. Second, in order to surmount such challenges, we present a computational framework which has been designed with scalability and efficiency in mind, including mechanisms to avoid the most common pitfalls. Third, we illustrate its performance with several case studies considering the central carbon metabolism of S. cerevisiae and B. subtilis. In particular, we consider metabolic dynamics during nutrient shift experiments. CONCLUSIONS: We show how multi-objective optimal control can be used to predict temporal profiles of enzyme activation and metabolite concentrations in complex metabolic pathways. Further, we also show how to consider general cost/benefit trade-offs. In this study we have considered metabolic pathways, but this computational framework can also be applied to analyze the dynamics of other complex pathways, such as signal transduction or gene regulatory networks.

MOTIVATION: Kinetic models contain unknown parameters that are estimated by optimizing the fit to experimental data. This task can be computationally challenging due to the presence of local optima and ill-conditioning. While a variety of optimization methods have been suggested to surmount these issues, it is difficult to choose the best one for a given problem a priori. A systematic comparison of parameter estimation methods for problems with tens to hundreds of optimization variables is currently missing, and smaller studies provided contradictory findings. RESULTS: We use a collection of benchmarks to evaluate the performance of two families of optimization methods: (i) multi-starts of deterministic local searches and (ii) stochastic global optimization metaheuristics; the latter may be combined with deterministic local searches, leading to hybrid methods. A fair comparison is ensured through a collaborative evaluation and a consideration of multiple performance metrics. We discuss possible evaluation criteria to assess the trade-off between computational efficiency and robustness. Our results show that, thanks to recent advances in the calculation of parametric sensitivities, a multi-start of gradient-based local methods is often a successful strategy, but a better performance can be obtained with a hybrid metaheuristic. The best performer combines a global scatter search metaheuristic with an interior point local method, provided with gradients estimated with adjoint-based sensitivities. We provide an implementation of this method to render it available to the scientific community. AVAILABILITY AND IMPLEMENTATION: The code to reproduce the results is provided as Supplementary Material and is available at Zenodo https://doi.org/10.5281/zenodo.1304034. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

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).

Motivation: Mathematical modeling using ordinary differential equations is used in systems biology to improve the understanding of dynamic biological processes. The parameters of ordinary differential equation models are usually estimated from experimental data. To analyze a priori the uniqueness of the solution of the estimation problem, structural identifiability analysis methods have been developed. Results: We introduce GenSSI 2.0, an advancement of the software toolbox GenSSI (Generating Series for testing Structural Identifiability). GenSSI 2.0 is the first toolbox for structural identifiability analysis to implement Systems Biology Markup Language import, state/parameter transformations and multi-experiment structural identifiability analysis. In addition, GenSSI 2.0 supports a range of MATLAB versions and is computationally more efficient than its previous version, enabling the analysis of more complex models. Availability and implementation: GenSSI 2.0 is an open-source MATLAB toolbox and available at https://github.com/genssi-developer/GenSSI. Contact: [email protected] or [email protected]. Supplementary information: Supplementary data are available at Bioinformatics online.

Motivation: Optimality principles have been used to explain many biological processes and systems. However, the functions being optimized are in general unknown a priori. Here we present an inverse optimal control framework for modeling dynamics in systems biology. The objective is to identify the underlying optimality principle from observed time-series data and simultaneously estimate unmeasured time-dependent inputs and time-invariant model parameters. As a special case, we also consider the problem of optimal simultaneous estimation of inputs and parameters from noisy data. After presenting a general statement of the inverse optimal control problem, and discussing special cases of interest, we outline numerical strategies which are scalable and robust. Results: We discuss the existence, relevance and implications of identifiability issues in the above problems. We present a robust computational approach based on regularized cost functions and the use of suitable direct numerical methods based on the control-vector parameterization approach. To avoid convergence to local solutions, we make use of hybrid global-local methods. We illustrate the performance and capabilities of this approach with several challenging case studies, including simulated and real data. We pay particular attention to the computational scalability of our approach (with the objective of considering large numbers of inputs and states). We provide a software implementation of both the methods and the case studies. Availability and implementation: The code used to obtain the results reported here is available at https://zenodo.org/record/1009541. Supplementary information: Supplementary data are available at Bioinformatics online.

Summary Optimization problems appear in many different applications in science and engineering. A large number of different algorithms have been proposed for solving them; however, there is no unique general optimization method that performs efficiently across a diverse set of problems. Thus, a multimethod optimization, in which different algorithms cooperate to outperform the results obtained by any of them in isolation, is a very appealing alternative. Besides, as real‐life optimization problems are becoming more and more challenging, the use of HPC techniques to implement these algorithms represents an effective strategy to speed up the time‐to‐solution. In addition, a parallel multimethod approach can benefit from the effortless access to q large number of distributed resources facilitated by cloud computing. In this paper, we propose a self‐adaptive cooperative parallel multimethod for global optimization. This proposal aims to perform a thorough exploration of the solution space by means of multiple concurrent executions of a broad range of search strategies. For its evaluation, we consider an extremely challenging case‐study from the field of computational systems biology. We also assess the performance of the proposal on a public cloud, demonstrating both the potential of the multimethod approach and the opportunity that the cloud provides for these problems.

Inferring the structure of unknown cellular networks is a main challenge in computational biology. Data-driven approaches based on information theory can determine the existence of interactions among network nodes automatically. However, the elucidation of certain features-such as distinguishing between direct and indirect interactions or determining the direction of a causal link-requires estimating information-theoretic quantities in a multidimensional space. This can be a computationally demanding task, which acts as a bottleneck for the application of elaborate algorithms to large-scale network inference problems. The computational cost of such calculations can be alleviated by the use of compiled programs and parallelization. To this end, we have developed PREMER (Parallel Reverse Engineering with Mutual information & Entropy Reduction), a software toolbox that can run in parallel and sequential environments. It uses information theoretic criteria to recover network topology and determine the strength and causality of interactions, and allows incorporating prior knowledge, imputing missing data, and correcting outliers. PREMER is a free, open source software tool that does not require any commercial software. Its core algorithms are programmed in FORTRAN 90 and implement OpenMP directives. It has user interfaces in Python and MATLAB/Octave, and runs on Windows, Linux, and OSX (https://sites.google.com/site/premertoolbox/).

Systems biology is an emerging approach focused in generating new knowledge about complex biological systems by combining experimental data with mathematical modeling and advanced computational techniques. Many problems in this field are extremely challenging and require substantial supercomputing resources to be solved. This is the case of parameter estimation in large-scale nonlinear dynamic systems biology models. Recently, Cloud Computing has emerged as a new paradigm for on-demand delivery of computing resources. However, scientific computing community has been quite hesitant in using the Cloud, simply because traditional programming models do not fit well with the new paradigm, and the earliest cloud programming models do not allow most scientific computations being efficiently run in the Cloud. In this paper we explore and compare two distributed computing models: the MPI (message-passing interface) model, that is high-performance oriented, and the Spark model, which is throughput oriented but outperforms other cloud programming solutions adding improved support for iterative algorithms through in-memory computing. The performance of a very well known metaheuristic, the Differential Evolution algorithm, has been thoroughly assessed using a challenging parameter estimation problem from the domain of computational systems biology. The experiments have been carried out both in a local cluster and in the Microsoft Azure public cloud, allowing performance and cost evaluation for both infrastructures.

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.

Despite significant efforts and remarkable progress, the inference of signaling networks from experimental data remains very challenging. The problem is particularly difficult when the objective is to obtain a dynamic model capable of predicting the effect of novel perturbations not considered during model training. The problem is ill-posed due to the nonlinear nature of these systems, the fact that only a fraction of the involved proteins and their post-translational modifications can be measured, and limitations on the technologies used for growing cells in vitro, perturbing them, and measuring their variations. As a consequence, there is a pervasive lack of identifiability. To overcome these issues, we present a methodology called SELDOM (enSEmbLe of Dynamic lOgic-based Models), which builds an ensemble of logic-based dynamic models, trains them to experimental data, and combines their individual simulations into an ensemble prediction. It also includes a model reduction step to prune spurious interactions and mitigate overfitting. SELDOM is a data-driven method, in the sense that it does not require any prior knowledge of the system: the interaction networks that act as scaffolds for the dynamic models are inferred from data using mutual information. We have tested SELDOM on a number of experimental and in silico signal transduction case-studies, including the recent HPN-DREAM breast cancer challenge. We found that its performance is highly competitive compared to state-of-the-art methods for the purpose of recovering network topology. More importantly, the utility of SELDOM goes beyond basic network inference (i.e. uncovering static interaction networks): it builds dynamic (based on ordinary differential equation) models, which can be used for mechanistic interpretations and reliable dynamic predictions in new experimental conditions (i.e. not used in the training). For this task, SELDOM's ensemble prediction is not only consistently better than predictions from individual models, but also often outperforms the state of the art represented by the methods used in the HPN-DREAM challenge.

The concept of dynamical compensation has been recently introduced to describe the ability of a biological system to keep its output dynamics unchanged in the face of varying parameters. However, the original definition of dynamical compensation amounts to lack of structural identifiability. This is relevant if model parameters need to be estimated, as is often the case in biological modelling. Care should we taken when using an unidentifiable model to extract biological insight: the estimated values of structurally unidentifiable parameters are meaningless, and model predictions about unmeasured state variables can be wrong. Taking this into account, we explore alternative definitions of dynamical compensation that do not necessarily imply structural unidentifiability. Accordingly, we show different ways in which a model can be made identifiable while exhibiting dynamical compensation. Our analyses enable the use of the new concept of dynamical compensation in the context of parameter identification, and reconcile it with the desirable property of structural identifiability.

MOTIVATION: Many problems of interest in dynamic modeling and control of biological systems can be posed as non-linear optimization problems subject to algebraic and dynamic constraints. In the context of modeling, this is the case of, e.g. parameter estimation, optimal experimental design and dynamic flux balance analysis. In the context of control, model-based metabolic engineering or drug dose optimization problems can be formulated as (multi-objective) optimal control problems. Finding a solution to those problems is a very challenging task which requires advanced numerical methods. RESULTS: This work presents the AMIGO2 toolbox: the first multiplatform software tool that automatizes the solution of all those problems, offering a suite of state-of-the-art (multi-objective) global optimizers and advanced simulation approaches. AVAILABILITY AND IMPLEMENTATION: The toolbox and its documentation are available at: sites.google.com/site/amigo2toolbox CONTACT: [email protected] information: Supplementary data are available at Bioinformatics online.

MOTIVATION: The design of de novo circuits with predefined performance specifications is a challenging problem in Synthetic Biology. Computational models and tools have proved to be crucial for a successful wet lab implementation. Natural gene circuits are complex, subject to evolutionary tradeoffs and playing multiple roles. However, most synthetic designs implemented to date are simple and perform a single task. As the field progresses, advanced computational tools are needed in order to handle greater levels of circuit complexity in a more flexible way and considering multiple design criteria. RESULTS: This works presents SYNBADm (SYNthetic Biology Automated optimal Design in Matlab), a software toolbox for the automatic optimal design of gene circuits with targeted functions from libraries of components. This tool makes use of global optimization to simultaneously search the space of structures and kinetic parameters. SYNBADm can efficiently handle high levels of circuit complexity and can consider multiple design criteria through multiobjective optimization. Further, it provides flexible design capabilities, i.e. the user can define the modeling framework, library of components and target performance function(s). AVAILABILITY AND IMPLEMENTATION: SYNBADm runs under the popular MATLAB computational environment, and is available under GPLv3 license at https://sites.google.com/site/synbadm CONTACT: [email protected] or [email protected].

MOTIVATION: Systems biology models can be used to test new hypotheses formulated on the basis of previous knowledge or new experimental data, contradictory with a previously existing model. New hypotheses often come in the shape of a set of possible regulatory mechanisms. This search is usually not limited to finding a single regulation link, but rather a combination of links subject to great uncertainty or no information about the kinetic parameters. RESULTS: In this work, we combine a logic-based formalism, to describe all the possible regulatory structures for a given dynamic model of a pathway, with mixed-integer dynamic optimization (MIDO). This framework aims to simultaneously identify the regulatory structure (represented by binary parameters) and the real-valued parameters that are consistent with the available experimental data, resulting in a logic-based differential equation model. The alternative to this would be to perform real-valued parameter estimation for each possible model structure, which is not tractable for models of the size presented in this work. The performance of the method presented here is illustrated with several case studies: a synthetic pathway problem of signaling regulation, a two-component signal transduction pathway in bacterial homeostasis, and a signaling network in liver cancer cells. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online. CONTACT: [email protected] or [email protected].

BACKGROUND: The inference of complex networks from data is a challenging problem in biological sciences, as well as in a wide range of disciplines such as chemistry, technology, economics, or sociology. The quantity and quality of the data greatly affect the results. While many methodologies have been developed for this task, they seldom take into account issues such as missing data or outlier detection and correction, which need to be properly addressed before network inference. RESULTS: Here we present an approach to (i) handle missing data and (ii) detect and correct outliers based on multivariate projection to latent structures. The method, called trimmed scores regression (TSR), enables network inference methods to analyse incomplete datasets by imputing the missing values coherently with the latent data structure. Furthermore, it substitutes the faulty values in a dataset by proper estimations. We provide an implementation of this approach, and show how it can be integrated with any network inference method as a preliminary data curation step. This functionality is demonstrated with a state of the art network inference method based on mutual information distance and entropy reduction, MIDER. CONCLUSION: The methodology presented here enables network inference methods to analyse a large number of incomplete and faulty datasets that could not be reliably analysed so far. Our comparative studies show the superiority of TSR over other missing data approaches used by practitioners. Furthermore, the method allows for outlier detection and correction.

BACKGROUND: Adjusting the capacity of metabolic pathways in response to rapidly changing environmental conditions is an important component of microbial adaptation strategies to stochastic environments. In this work, we use advanced dynamic optimization techniques combined with theoretical models to study which reactions in pathways are optimally targeted by regulatory interactions in order to minimize the regulatory effort that is required to adjust the flux through a complex metabolic network. Moreover, we analyze how constraints in the speed at which an organism can respond on a proteomic level influences these optimal targets of pathway control. RESULTS: We find that limitations in protein biosynthetic rates have a strong influence. With increasing protein biosynthetic rates the regulatory effort targeting the initial enzyme in a pathway is reduced while the regulatory effort in the terminal enzyme is increased. Studying the impact of allosteric regulation for different pathway topologies, we find that the presence of feedback inhibition by products of metabolic pathways allows organisms to reduce the regulatory effort that is required to control a metabolic pathway in all cases. In a linear pathway this even leads to the case where the sole transcriptional regulatory control of the terminal enzyme is sufficient to control flux through the entire pathway. We confirm the utilization of these pathway regulation strategies through the large-scale analysis of transcriptional regulation in several hundred prokaryotes. CONCLUSIONS: This work expands our knowledge about optimal programs of pathway control. Optimal targets of pathway control strongly depend on the speed at which proteins can be synthesized. Moreover, post-translational regulation such as allosteric regulation allows to strongly reduce the number of transcriptional regulatory interactions required to control a metabolic pathway across different pathway topologies.

BACKGROUND: Optimization is the key to solving many problems in computational biology. Global optimization methods, which provide a robust methodology, and metaheuristics in particular have proven to be the most efficient methods for many applications. Despite their utility, there is a limited availability of metaheuristic tools. RESULTS: We present MEIGO, an R and Matlab optimization toolbox (also available in Python via a wrapper of the R version), that implements metaheuristics capable of solving diverse problems arising in systems biology and bioinformatics. The toolbox includes the enhanced scatter search method (eSS) for continuous nonlinear programming (cNLP) and mixed-integer programming (MINLP) problems, and variable neighborhood search (VNS) for Integer Programming (IP) problems. Additionally, the R version includes BayesFit for parameter estimation by Bayesian inference. The eSS and VNS methods can be run on a single-thread or in parallel using a cooperative strategy. The code is supplied under GPLv3 and is available at http://www.iim.csic.es/~gingproc/meigo.html. Documentation and examples are included. The R package has been submitted to BioConductor. We evaluate MEIGO against optimization benchmarks, and illustrate its applicability to a series of case studies in bioinformatics and systems biology where it outperforms other state-of-the-art methods. CONCLUSIONS: MEIGO provides a free, open-source platform for optimization that can be applied to multiple domains of systems biology and bioinformatics. It includes efficient state of the art metaheuristics, and its open and modular structure allows the addition of further methods.

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.

UNLABELLED: Chemical reaction network theory is widely used in modeling and analyzing complex biochemical systems such as metabolic networks and cell signalling pathways. Being able to produce all the biologically and chemically important qualitative dynamical features, chemical reaction networks (CRNs) have attracted significant attention in the systems biology community. It is well-known that the reliable inference of CRN models generally requires thorough identifiability and distinguishability analysis together with carefully selected prior modeling assumptions. Here, we present a software toolbox CRNreals that supports the distinguishability and identifiability analysis of CRN models using recently published optimization-based procedures. AVAILABILITY AND IMPLEMENTATION: The CRNreals toolbox and the associated documentation are available at http://www.iim.csic.es/~gingproc/CRNreals/. The toolbox runs under the popular MATLAB computational environment and supports several free and commercial linear programming and mixed integer linear programming solvers.

MOTIVATION: Mathematical models of complex biological systems usually consist of sets of differential equations which depend on several parameters which are not accessible to experimentation. These parameters must be estimated by fitting the model to experimental data. This estimation problem is very challenging due to the non-linear character of the dynamics, the large number of parameters and the frequently poor information content of the experimental data (poor practical identifiability). The design of optimal (more informative) experiments is an associated problem of the highest interest. RESULTS: This work presents AMIGO, a toolbox which facilitates parametric identification by means of advanced numerical techniques which cover the full iterative identification procedure putting especial emphasis on robust methods for parameter estimation and practical identifiability analyses, plus flexible capabilities for optimal experimental design. AVAILABILITY: The toolbox and the corresponding documentation may be downloaded from: http://www.iim.csic.es/~amigo CONTACT: [email protected].

SUMMARY: Mathematical modeling has a key role in systems biology. Model building is often regarded as an iterative loop involving several tasks, among which the estimation of unknown parameters of the model from a certain set of experimental data is of central importance. This problem of parameter estimation has many possible pitfalls, and modelers should be very careful to avoid them. Many of such difficulties arise from a fundamental (yet often overlooked) property: the so-called structural (or a priori) identifiability, which considers the uniqueness of the estimated parameters. Obviously, the structural identifiability of any tentative model should be checked at the beginning of the model building loop. However, checking this property for arbitrary non-linear dynamic models is not an easy task. Here we present a software toolbox, GenSSI (Generating Series for testing Structural Identifiability), which enables non-expert users to carry out such analysis. The toolbox runs under the popular MATLAB environment and is accompanied by detailed documentation and relevant examples. AVAILABILITY: The GenSSI toolbox and the related documentation are available at http://www.iim.csic.es/%7Egenssi. CONTACT: [email protected].

SUMMARY: SensSB (Sensitivity Analysis for Systems Biology) is an easy to use, MATLAB-based software toolbox, which integrates several local and global sensitivity methods that can be applied to a wide variety of biological models. In addition to addressing the sensitivity analysis problem, SensSB aims to cover all the steps involved during the modeling process. The main features of SensSB are: (i) derivative and variance-based global sensitivity analysis, (ii) pseudo-global identifiability analysis, (iii) optimal experimental design (OED) based on global sensitivities, (iv) robust parameter estimation, (v) local sensitivity and identifiability analysis, (vi) confidence intervals of the estimated parameters and (vii) OED based on the Fisher Information Matrix (FIM). SensSB is also able to import models in the Systems Biology Mark-up Language (SBML) format. Several examples from simple analytical functions to more complex biological pathways have been implemented and can be downloaded together with the toolbox. The importance of using sensitivity analysis techniques for identifying unessential parameters and designing new experiments is quantified by increased identifiability metrics of the models and decreased confidence intervals of the estimated parameters. AVAILABILITY: SensSB is a software toolbox freely downloadable from http://www.iim.csic.es/ approximately gingproc/SensSB.html. The web site also contains several examples and an extensive documentation. CONTACT: [email protected] SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

BACKGROUND: Mathematical optimization aims to make a system or design as effective or functional as possible, computing the quality of the different alternatives using a mathematical model. Most models in systems biology have a dynamic nature, usually described by sets of differential equations. Dynamic optimization addresses this class of systems, seeking the computation of the optimal time-varying conditions (control variables) to minimize or maximize a certain performance index. Dynamic optimization can solve many important problems in systems biology, including optimal control for obtaining a desired biological performance, the analysis of network designs and computer aided design of biological units. RESULTS: Here, we present a software toolbox, DOTcvpSB, which uses a rich ensemble of state-of-the-art numerical methods for solving continuous and mixed-integer dynamic optimization (MIDO) problems. The toolbox has been written in MATLAB and provides an easy and user friendly environment, including a graphical user interface, while ensuring a good numerical performance. Problems are easily stated thanks to the compact input definition. The toolbox also offers the possibility of importing SBML models, thus enabling it as a powerful optimization companion to modelling packages in systems biology. It serves as a means of handling generic black-box models as well. CONCLUSION: Here we illustrate the capabilities and performance of DOTcvpSB by solving several challenging optimization problems related with bioreactor optimization, optimal drug infusion to a patient and the minimization of intracellular oscillations. The results illustrate how the suite of solvers available allows the efficient solution of a wide class of dynamic optimization problems, including challenging multimodal ones. The toolbox is freely available for academic use.

BACKGROUND: We consider the problem of parameter estimation (model calibration) in nonlinear dynamic models of biological systems. Due to the frequent ill-conditioning and multi-modality of many of these problems, traditional local methods usually fail (unless initialized with very good guesses of the parameter vector). In order to surmount these difficulties, global optimization (GO) methods have been suggested as robust alternatives. Currently, deterministic GO methods can not solve problems of realistic size within this class in reasonable computation times. In contrast, certain types of stochastic GO methods have shown promising results, although the computational cost remains large. Rodriguez-Fernandez and coworkers have presented hybrid stochastic-deterministic GO methods which could reduce computation time by one order of magnitude while guaranteeing robustness. Our goal here was to further reduce the computational effort without loosing robustness. RESULTS: We have developed a new procedure based on the scatter search methodology for nonlinear optimization of dynamic models of arbitrary (or even unknown) structure (i.e. black-box models). In this contribution, we describe and apply this novel metaheuristic, inspired by recent developments in the field of operations research, to a set of complex identification problems and we make a critical comparison with respect to the previous (above mentioned) successful methods. CONCLUSION: Robust and efficient methods for parameter estimation are of key importance in systems biology and related areas. The new metaheuristic presented in this paper aims to ensure the proper solution of these problems by adopting a global optimization approach, while keeping the computational effort under reasonable values. This new metaheuristic was applied to a set of three challenging parameter estimation problems of nonlinear dynamic biological systems, outperforming very significantly all the methods previously used for these benchmark problems.