Inferring signaling pathway topologies from multiple perturbation measurements of specific biochemical species

FoxA1 knockdown promotes BMSC osteogenesis in part by activating the ERK1/2 signaling pathway and preventing ovariectomy-induced bone loss

The influence of deep learning in the medical and molecular biology sectors is swiftly growing and holds the potential to improve numerous crucial domains. Osteoporosis is a significant global health issue, and the current treatment options are highly restricted. Transplanting genetically engineered MSCs has been acknowledged as a highly promising therapy for osteoporosis. We utilized a random walk-based technique to discern genes associated with ossification. The osteogenic value of these genes was assessed on the basis of information found in published scientific literature. GO enrichment analysis of these genes was performed to determine if they were enriched in any certain function. Immunohistochemical and western blot techniques were used to identify and measure protein expression. The expression of genes involved in osteogenic differentiation was examined via qRT‒PCR. Lentiviral transfection was utilized to suppress the expression of the FOXA1 gene in hBMSCs. An in vivo mouse model of ovariectomy was created, and radiographic examination was conducted to confirm the impact of FOXA1 knockdown on osteoporosis. The osteogenic score of each gene was calculated by assessing its similarity to osteo-specific genes. The majority of the genes with the highest rankings were linked with osteogenic differentiation, indicating that our approach is useful for identifying genes associated with ossification. GO enrichment analysis revealed that these pathways are enriched primarily in bone-related processes. FOXA1 is a crucial transcription factor that controls the process of osteogenic differentiation, as indicated by similarity analysis. FOXA1 was significantly increased in those with osteoporosis. Downregulation of FOXA1 markedly augmented the expression of osteoblast-specific genes and proteins, activated the ERK1/2 signaling pathway, intensified ALP activity, and promoted mineral deposition. In addition, excessive expression of FOXA1 significantly reduced ALP activity and mineral deposits. Using a mouse model in which the ovaries were surgically removed, researchers reported that suppressing the FOXA1 gene in bone marrow stem cells (BMSCs) prevented the loss of bone density caused by ovariectomy. This finding was confirmed by analyzing the bone structure via micro-CT. Furthermore, our approach can distinguish genes that exhibit osteogenic differentiation characteristics. This ability can aid in the identification of novel genes associated with osteogenic differentiation, which can be utilized in the treatment of osteoporosis. Computational and laboratory evidence indicates that reducing the expression of FOXA1 enhances the process of bone formation in bone marrow-derived mesenchymal stem cells (BMSCs) and may serve as a promising approach to prevent osteoporosis.

Mathematical Modeling and Inference of Epidermal Growth Factor-Induced Mitogen-Activated Protein Kinase Cell Signaling Pathways

The mitogen-activated protein kinase (MAPK) pathway is an important intracellular signaling cascade that plays a key role in various cellular processes. Understanding the regulatory mechanisms of this pathway is essential for developing effective interventions and targeted therapies for related diseases. Recent advances in single-cell proteomic technologies have provided unprecedented opportunities to investigate the heterogeneity and noise within complex, multi-signaling networks across diverse cells and cell types. Mathematical modeling has become a powerful interdisciplinary tool that bridges mathematics and experimental biology, providing valuable insights into these intricate cellular processes. In addition, statistical methods have been developed to infer pathway topologies and estimate unknown parameters within dynamic models. This review presents a comprehensive analysis of how mathematical modeling of the MAPK pathway deepens our understanding of its regulatory mechanisms, enhances the prediction of system behavior, and informs experimental research, with a particular focus on recent advances in modeling and inference using single-cell proteomic data.

Posterior marginalization accelerates Bayesian inference for dynamical models of biological processes

Bayesian inference is an important method in the life and natural sciences for learning from data. It provides information about parameter and prediction uncertainties. Yet, generating representative samples from the posterior distribution is often computationally challenging. Here, we present an approach that lowers the computational complexity of sample generation for dynamical models with scaling, offset, and noise parameters. The proposed method is based on the marginalization of the posterior distribution. We provide analytical results for a broad class of problems with conjugate priors and show that the method is suitable for a large number of applications. Subsequently, we demonstrate the benefit of the approach for applications from the field of systems biology. We report an improvement up to 50 times in the effective sample size per unit of time. As the scheme is broadly applicable, it will facilitate Bayesian inference in different research fields.

Challenges and Perspectives in Target Identification and Mechanism Illustration for Chinese Medicine

Chinese medicine (CM) is an important resource for human life understanding and discovery of drugs. However, due to the unclear pharmacological mechanism caused by unclear target, research and international promotion of many active components have made little progress in the past decades of years. CM is mainly composed of multi-ingredients with multi-targets. The identification of targets of multiple active components and the weight analysis of multiple targets in a specific pathological environment, that is, the determination of the most important target is the main obstacle to the mechanism clarification and thus hinders its internationalization. In this review, the main approach to target identification and network pharmacology were summarized. And BIBm (Bayesian inference modeling), a powerful method for drug target identification and key pathway determination was introduced. We aim to provide a new scientific basis and ideas for the development and international promotion of new drugs based on CM.

New strategies for targeting kinase networks in cancer

Targeted strategies against specific driver molecules of cancer have brought about many advances in cancer treatment since the early success of the first small-molecule inhibitor Gleevec. Today, there are a multitude of targeted therapies approved by the Food and Drug Administration for the treatment of cancer. However, the initial efficacy of virtually every targeted treatment is often reversed by tumor resistance to the inhibitor through acquisition of new mutations in the target molecule, or reprogramming of the epigenome, transcriptome, or kinome of the tumor cells. At the core of this clinical problem lies the assumption that targeted treatments will only be efficacious if the inhibitors are used at their maximum tolerated doses. Such aggressive regimens create strong selective pressure on the evolutionary progression of the tumor, resulting in resistant cells. High-dose single agent treatments activate alternative mechanisms that bypass the inhibitor, while high-dose combinatorial treatments suffer from increased toxicity resulting in treatment cessation. Although there is an arsenal of targeted agents being tested clinically and preclinically, identifying the most effective combination treatment plan remains a challenge. In this review, we discuss novel targeted strategies with an emphasis on the recent cross-disciplinary studies demonstrating that it is possible to achieve antitumor efficacy without increasing toxicity by adopting low-dose multitarget approaches to treatment of cancer and metastasis.

Limited inhibition of multiple nodes in a driver network blocks metastasis

Metastasis suppression by high-dose, multi-drug targeting is unsuccessful due to network heterogeneity and compensatory network activation. Here, we show that targeting driver network signaling capacity by limited inhibition of core pathways is a more effective anti-metastatic strategy. This principle underlies the action of a physiological metastasis suppressor, Raf Kinase Inhibitory Protein (RKIP), that moderately decreases stress-regulated MAP kinase network activity, reducing output to transcription factors such as pro-metastastic BACH1 and motility-related target genes. We developed a low-dose four-drug mimic that blocks metastatic colonization in mouse breast cancer models and increases survival. Experiments and network flow modeling show limited inhibition of multiple pathways is required to overcome variation in MAPK network topology and suppress signaling output across heterogeneous tumor cells. Restricting inhibition of individual kinases dissipates surplus signal, preventing threshold activation of compensatory kinase networks. This low-dose multi-drug approach to decrease signaling capacity of driver networks represents a transformative, clinically relevant strategy for anti-metastatic treatment.

A Probabilistic Approach to Explore Signal Execution Mechanisms With Limited Experimental Data

Mathematical models of biochemical reaction networks are central to the study of dynamic cellular processes and hypothesis generation that informs experimentation and validation. Unfortunately, model parameters are often not available and sparse experimental data leads to challenges in model calibration and parameter estimation. This can in turn lead to unreliable mechanistic interpretations of experimental data and the generation of poorly conceived hypotheses for experimental validation. To address this challenge, we evaluate whether a Bayesian-inspired probability-based approach, that relies on expected values for quantities of interest calculated from available information regarding the reaction network topology and parameters can be used to qualitatively explore hypothetical biochemical network execution mechanisms in the context of limited available data. We test our approach on a model of extrinsic apoptosis execution to identify preferred signal execution modes across varying conditions. Apoptosis signal processing can take place either through a mitochondria independent (Type I) mode or a mitochondria dependent (Type II) mode. We first show that in silico knockouts, represented by model subnetworks, successfully identify the most likely execution mode for specific concentrations of key molecular regulators. We then show that changes in molecular regulator concentrations alter the overall reaction flux through the network by shifting the primary route of signal flow between the direct caspase and mitochondrial pathways. Our work thus demonstrates that probabilistic approaches can be used to explore the qualitative dynamic behavior of model biochemical systems even with missing or sparse data.

Bayesian modeling reveals metabolite-dependent ultrasensitivity in the cyanobacterial circadian clock

Mathematical models can enable a predictive understanding of mechanism in cell biology by quantitatively describing complex networks of interactions, but such models are often poorly constrained by available data. Owing to its relative biochemical simplicity, the core circadian oscillator in Synechococcus elongatus has become a prototypical system for studying how collective dynamics emerge from molecular interactions. The oscillator consists of only three proteins, KaiA, KaiB, and KaiC, and near-24-h cycles of KaiC phosphorylation can be reconstituted in vitro. Here, we formulate a molecularly detailed but mechanistically naive model of the KaiA-KaiC subsystem and fit it directly to experimental data within a Bayesian parameter estimation framework. Analysis of the fits consistently reveals an ultrasensitive response for KaiC phosphorylation as a function of KaiA concentration, which we confirm experimentally. This ultrasensitivity primarily results from the differential affinity of KaiA for competing nucleotide-bound states of KaiC. We argue that the ultrasensitive stimulus-response relation likely plays an important role in metabolic compensation by suppressing premature phosphorylation at nighttime.

The road to ERK activation: Do neurons take alternate routes?

The ERK cascade is a central signaling pathway that regulates a wide variety of cellular processes including proliferation, differentiation, learning and memory, development, and synaptic plasticity. A wide range of inputs travel from the membrane through different signaling pathway routes to reach activation of one set of output kinases, ERK1&2. The classical ERK activation pathway beings with growth factor activation of receptor tyrosine kinases. Numerous G-protein coupled receptors and ionotropic receptors also lead to ERK through increases in the second messengers calcium and cAMP. Though both types of pathways are present in diverse cell types, a key difference is that most stimuli to neurons, e.g. synaptic inputs, are transient, on the order of milliseconds to seconds, whereas many stimuli acting on non-neural tissue, e.g. growth factors, are longer duration. The ability to consolidate these inputs to regulate the activation of ERK in response to diverse signals raises the question of which factors influence the difference in ERK activation pathways. This review presents both experimental studies and computational models aimed at understanding the control of ERK activation and whether there are fundamental differences between neurons and other cells. Our main conclusion is that differences between cell types are quite subtle, often related to differences in expression pattern and quantity of some molecules such as Raf isoforms. In addition, the spatial location of ERK is critical, with regulation by scaffolding proteins producing differences due to colocalization of upstream molecules that may differ between neurons and other cells.

Approximating multivariate posterior distribution functions from Monte Carlo samples for sequential Bayesian inference

An important feature of Bayesian statistics is the opportunity to do sequential inference: the posterior distribution obtained after seeing a dataset can be used as prior for a second inference. However, when Monte Carlo sampling methods are used for inference, we only have a set of samples from the posterior distribution. To do sequential inference, we then either have to evaluate the second posterior at only these locations and reweight the samples accordingly, or we can estimate a functional description of the posterior probability distribution from the samples and use that as prior for the second inference. Here, we investigated to what extent we can obtain an accurate joint posterior from two datasets if the inference is done sequentially rather than jointly, under the condition that each inference step is done using Monte Carlo sampling. To test this, we evaluated the accuracy of kernel density estimates, Gaussian mixtures, mixtures of factor analyzers, vine copulas and Gaussian processes in approximating posterior distributions, and then tested whether these approximations can be used in sequential inference. In low dimensionality, Gaussian processes are more accurate, whereas in higher dimensionality Gaussian mixtures, mixtures of factor analyzers or vine copulas perform better. In our test cases of sequential inference, using posterior approximations gives more accurate results than direct sample reweighting, but joint inference is still preferable over sequential inference whenever possible. Since the performance is case-specific, we provide an R package mvdens with a unified interface for the density approximation methods.

TopoFilter: a MATLAB package for mechanistic model identification in systems biology

Background

To develop mechanistic dynamic models in systems biology, one often needs to identify all (or minimal) representations of the biological processes that are consistent with experimental data, out of a potentially large set of hypothetical mechanisms. However, a simple enumeration of all alternatives becomes quickly intractable when the number of model parameters grows. Selecting appropriate dynamic models out of a large ensemble of models, taking the uncertainty in our biological knowledge and in the experimental data into account, is therefore a key current problem in systems biology.

Results

The TopoFilter package addresses this problem in a heuristic and automated fashion by implementing the previously described topological filtering method for Bayesian model selection. It includes a core heuristic for searching the space of submodels of a parametrized model, coupled with a sampling-based exploration of the parameter space. Recent developments of the method allow to balance exhaustiveness and speed of the model space search, to efficiently re-sample parameters, to parallelize the search, and to use custom scoring functions. We use a theoretical example to motivate these features and then demonstrate TopoFilter’s applicability for a yeast signaling network with more than 250’000 possible model structures.

Conclusions

TopoFilter is a flexible software framework that makes Bayesian model selection and reduction efficient and scalable to network models of a complexity that represents contemporary problems in, for example, cell signaling. TopoFilter is open-source, available under the GPL-3.0 license at https://gitlab.com/csb.ethz/TopoFilter. It includes installation instructions, a quickstart guide, a description of all package options, and multiple examples.

Exploiting network topology for large-scale inference of nonlinear reaction models

The development of chemical reaction models aids understanding and prediction in areas ranging from biology to electrochemistry and combustion. A systematic approach to building reaction network models uses observational data not only to estimate unknown parameters but also to learn model structure. Bayesian inference provides a natural approach to this data-driven construction of models. Yet traditional Bayesian model inference methodologies that numerically evaluate the evidence for each model are often infeasible for nonlinear reaction network inference, as the number of plausible models can be combinatorially large. Alternative approaches based on model-space sampling can enable large-scale network inference, but their realization presents many challenges. In this paper, we present new computational methods that make large-scale nonlinear network inference tractable. First, we exploit the topology of networks describing potential interactions among chemical species to design improved 'between-model' proposals for reversible-jump Markov chain Monte Carlo. Second, we introduce a sensitivity-based determination of move types which, when combined with network-aware proposals, yields significant additional gains in sampling performance. These algorithms are demonstrated on inference problems drawn from systems biology, with nonlinear differential equation models of species interactions.

The 11-Kilodalton Nonstructural Protein of Human Parvovirus B19 Facilitates Viral DNA Replication by Interacting with Grb2 through Its Proline-Rich Motifs

Lytic infection of human parvovirus B19 (B19V) takes place exclusively in human erythroid progenitor cells of bone marrow and fetal liver, which disrupts erythropoiesis. During infection, B19V expresses three nonstructural proteins (NS1, 11-kDa, and 7.5-kDa) and two structural proteins (VP1 and VP2). While NS1 is essential for B19V DNA replication, 11-kDa enhances viral DNA replication significantly. In this study, we confirmed the enhancement role of 11-kDa in viral DNA replication and elucidated the underlying mechanism. We found that 11-kDa specially interacts with cellular growth factor receptor-bound protein 2 (Grb2) during virus infection and in vitro We determined a high affinity interaction between 11-kDa and Grb2 that has an equilibrium dissociation constant (KD ) value of 18.13 nM. In vitro, one proline-rich motif was sufficient for 11-kDa to sustain a strong interaction with Grb2. In consistence, in vivo during infection, one proline-rich motif was enough for 11-kDa to significantly reduce phosphorylation of extracellular signal-regulated kinase (ERK). Mutations of all three proline-rich motifs of 11-kDa abolished its capability to reduce ERK activity and, accordingly, decreased viral DNA replication. Transduction of a lentiviral vector encoding a short hairpin RNA (shRNA) targeting Grb2 decreased the expression of Grb2 as well as the level of ERK phosphorylation, which resulted in an increase of B19V replication. These results, in concert, indicate that the B19V 11-kDa protein interacts with cellular Grb2 to downregulate ERK activity, which upregulates viral DNA replication.IMPORTANCE Human parvovirus B19 (B19V) infection causes hematological disorders and is the leading cause of nonimmunological fetal hydrops during pregnancy. During infection, B19V expresses two structural proteins, VP1 and VP2, and three nonstructural proteins, NS1, 11-kDa, and 7.5-kDa. While NS1 is essential, 11-kDa plays an enhancing role in viral DNA replication. Here, we elucidated a mechanism underlying 11-kDa protein-regulated B19V DNA replication. 11-kDa is tightly associated with cellular growth factor receptor-bound protein 2 (Grb2) during infection. In vitro, 11-kDa interacts with Grb2 with high affinity through three proline-rich motifs, of which at least one is indispensable for the regulation of viral DNA replication. 11-kDa and Grb2 interaction disrupts extracellular signal-regulated kinase (ERK) signaling, which mediates upregulation of B19V replication. Thus, our study reveals a novel mechanism of how a parvoviral small nonstructural protein regulates viral DNA replication by interacting with a host protein that is predominately expressed in the cytoplasm.

Statistical Network Inference for Time-Varying Molecular Data with Dynamic Bayesian Networks

In this chapter, we review the problem of network inference from time-course data, focusing on a class of graphical models known as dynamic Bayesian networks (DBNs). We discuss the relationship of DBNs to models based on ordinary differential equations, and consider extensions to nonlinear time dynamics. We provide an introduction to time-varying DBN models, which allow for changes to the network structure and parameters over time. We also discuss causal perspectives on network inference, including issues around model semantics that can arise due to missing variables. We present a case study of applying time-varying DBNs to gene expression measurements over the life cycle of Drosophila melanogaster. We finish with a discussion of future perspectives, including possible applications of time-varying network inference to single-cell gene expression data.

Integrating -omics data into genome-scale metabolic network models: principles and challenges

At genome scale, it is not yet possible to devise detailed kinetic models for metabolism because data on the in vivo biochemistry are too sparse. Predictive large-scale models for metabolism most commonly use the constraint-based framework, in which network structures constrain possible metabolic phenotypes at steady state. However, these models commonly leave many possibilities open, making them less predictive than desired. With increasingly available -omics data, it is appealing to increase the predictive power of constraint-based models (CBMs) through data integration. Many corresponding methods have been developed, but data integration is still a challenge and existing methods perform less well than expected. Here, we review main approaches for the integration of different types of -omics data into CBMs focussing on the methods' assumptions and limitations. We argue that key assumptions - often derived from single-enzyme kinetics - do not generally apply in the context of networks, thereby explaining current limitations. Emerging methods bridging CBMs and biochemical kinetics may allow for -omics data integration in a common framework to provide more accurate predictions.

Bayesian parameter estimation for biochemical reaction networks using region-based adaptive parallel tempering

Motivation

Mathematical models have become standard tools for the investigation of cellular processes and the unraveling of signal processing mechanisms. The parameters of these models are usually derived from the available data using optimization and sampling methods. However, the efficiency of these methods is limited by the properties of the mathematical model, e.g. non-identifiabilities, and the resulting posterior distribution. In particular, multi-modal distributions with long valleys or pronounced tails are difficult to optimize and sample. Thus, the developement or improvement of optimization and sampling methods is subject to ongoing research.

Results

We suggest a region-based adaptive parallel tempering algorithm which adapts to the problem-specific posterior distributions, i.e. modes and valleys. The algorithm combines several established algorithms to overcome their individual shortcomings and to improve sampling efficiency. We assessed its properties for established benchmark problems and two ordinary differential equation models of biochemical reaction networks. The proposed algorithm outperformed state-of-the-art methods in terms of calculation efficiency and mixing. Since the algorithm does not rely on a specific problem structure, but adapts to the posterior distribution, it is suitable for a variety of model classes.

Availability and implementation

The code is available both as Supplementary Material and in a Git repository written in MATLAB.

Supplementary information

Supplementary data are available at Bioinformatics online.

C-RAF function at the genome-wide transcriptome level: A systematic view

C-RAF was the first member of the RAF kinase family to be discovered. Since its discovery, C-RAF has been found to regulate many fundamental cell processes, such as cell proliferation, cell death, and metabolism. However, the majority of these functions are achieved through interactions with different proteins; the genes regulated by C-RAF in its active or inactive state remain unclear. In the work, we used RNA-seq analysis to study the global transcriptomes of C-RAF bearing or C-RAF knockout cells in quiescent or EGF activated states. We identified 3353 genes that are promoted or suppressed by C-RAF. Gene ontology and Kyoto Encyclopedia of Genes and Genomes analyses revealed that these genes are involved in drug addiction, cardiomyopathy, autoimmunity, and regulation of cell metabolism. Our results provide a panoramic view of C-RAF function, including known and novel functions, and have revealed potential targets for elucidating the role of C-RAF.

A Bayesian approach to estimating hidden variables as well as missing and wrong molecular interactions in ordinary differential equation-based mathematical models

Ordinary differential equations (ODEs) are a popular approach to quantitatively model molecular networks based on biological knowledge. However, such knowledge is typically restricted. Wrongly modelled biological mechanisms as well as relevant external influence factors that are not included into the model are likely to manifest in major discrepancies between model predictions and experimental data. Finding the exact reasons for such observed discrepancies can be quite challenging in practice. In order to address this issue, we suggest a Bayesian approach to estimate hidden influences in ODE-based models. The method can distinguish between exogenous and endogenous hidden influences. Thus, we can detect wrongly specified as well as missed molecular interactions in the model. We demonstrate the performance of our Bayesian dynamic elastic-net with several ordinary differential equation models from the literature, such as human JAK-STAT signalling, information processing at the erythropoietin receptor, isomerization of liquid α-Pinene, G protein cycling in yeast and UV-B triggered signalling in plants. Moreover, we investigate a set of commonly known network motifs and a gene-regulatory network. Altogether our method supports the modeller in an algorithmic manner to identify possible sources of errors in ODE-based models on the basis of experimental data.

Bayesian data integration for quantifying the contribution of diverse measurements to parameter estimates

Motivation

Computational models in biology are frequently underdetermined, due to limits in our capacity to measure biological systems. In particular, mechanistic models often contain parameters whose values are not constrained by a single type of measurement. It may be possible to achieve better model determination by combining the information contained in different types of measurements. Bayesian statistics provides a convenient framework for this, allowing a quantification of the reduction in uncertainty with each additional measurement type. We wished to explore whether such integration is feasible and whether it can allow computational models to be more accurately determined.

Results

We created an ordinary differential equation model of cell cycle regulation in budding yeast and integrated data from 13 different studies covering different experimental techniques. We found that for some parameters, a single type of measurement, relative time course mRNA expression, is sufficient to constrain them. Other parameters, however, were only constrained when two types of measurements were combined, namely relative time course and absolute transcript concentration. Comparing the estimates to measurements from three additional, independent studies, we found that the degradation and transcription rates indeed matched the model predictions in order of magnitude. The predicted translation rate was incorrect however, thus revealing a deficiency in the model. Since this parameter was not constrained by any of the measurement types separately, it was only possible to falsify the model when integrating multiple types of measurements. In conclusion, this study shows that integrating multiple measurement types can allow models to be more accurately determined.

Availability and implementation

The models and files required for running the inference are included in the Supplementary information.

Contact

l.wessels@nki.nl.

Supplementary information

Supplementary data are available at Bioinformatics online.

Comprehensive benchmarking of Markov chain Monte Carlo methods for dynamical systems

BACKGROUND

In quantitative biology, mathematical models are used to describe and analyze biological processes. The parameters of these models are usually unknown and need to be estimated from experimental data using statistical methods. In particular, Markov chain Monte Carlo (MCMC) methods have become increasingly popular as they allow for a rigorous analysis of parameter and prediction uncertainties without the need for assuming parameter identifiability or removing non-identifiable parameters. A broad spectrum of MCMC algorithms have been proposed, including single- and multi-chain approaches. However, selecting and tuning sampling algorithms suited for a given problem remains challenging and a comprehensive comparison of different methods is so far not available.

RESULTS

We present the results of a thorough benchmarking of state-of-the-art single- and multi-chain sampling methods, including Adaptive Metropolis, Delayed Rejection Adaptive Metropolis, Metropolis adjusted Langevin algorithm, Parallel Tempering and Parallel Hierarchical Sampling. Different initialization and adaptation schemes are considered. To ensure a comprehensive and fair comparison, we consider problems with a range of features such as bifurcations, periodical orbits, multistability of steady-state solutions and chaotic regimes. These problem properties give rise to various posterior distributions including uni- and multi-modal distributions and non-normally distributed mode tails. For an objective comparison, we developed a pipeline for the semi-automatic comparison of sampling results.

CONCLUSION

The comparison of MCMC algorithms, initialization and adaptation schemes revealed that overall multi-chain algorithms perform better than single-chain algorithms. In some cases this performance can be further increased by using a preceding multi-start local optimization scheme. These results can inform the selection of sampling methods and the benchmark collection can serve for the evaluation of new algorithms. Furthermore, our results confirm the need to address exploration quality of MCMC chains before applying the commonly used quality measure of effective sample size to prevent false analysis conclusions.

Scalable Parameter Estimation for Genome-Scale Biochemical Reaction Networks

Mechanistic mathematical modeling of biochemical reaction networks using ordinary differential equation (ODE) models has improved our understanding of small- and medium-scale biological processes. While the same should in principle hold for large- and genome-scale processes, the computational methods for the analysis of ODE models which describe hundreds or thousands of biochemical species and reactions are missing so far. While individual simulations are feasible, the inference of the model parameters from experimental data is computationally too intensive. In this manuscript, we evaluate adjoint sensitivity analysis for parameter estimation in large scale biochemical reaction networks. We present the approach for time-discrete measurement and compare it to state-of-the-art methods used in systems and computational biology. Our comparison reveals a significantly improved computational efficiency and a superior scalability of adjoint sensitivity analysis. The computational complexity is effectively independent of the number of parameters, enabling the analysis of large- and genome-scale models. Our study of a comprehensive kinetic model of ErbB signaling shows that parameter estimation using adjoint sensitivity analysis requires a fraction of the computation time of established methods. The proposed method will facilitate mechanistic modeling of genome-scale cellular processes, as required in the age of omics.

In this manuscript, we introduce a scalable method for parameter estimation for genome-scale biochemical reaction networks. Mechanistic models for genome-scale biochemical reaction networks describe the behavior of thousands of chemical species using thousands of parameters. Standard methods for parameter estimation are usually computationally intractable at these scales. Adjoint sensitivity based approaches have been suggested to have superior scalability but any rigorous evaluation is lacking. We implement a toolbox for adjoint sensitivity analysis for biochemical reaction network which also supports the import of SBML models. We show by means of a set of benchmark models that adjoint sensitivity based approaches unequivocally outperform standard approaches for large-scale models and that the achieved speedup increases with respect to both the number of parameters and the number of chemical species in the model. This demonstrates the applicability of adjoint sensitivity based approaches to parameter estimation for genome-scale mechanistic model. The MATLAB toolbox implementing the developed methods is available from http://ICB-DCM.github.io/AMICI/.

Respectful Modeling: Addressing Uncertainty in Dynamic System Models for Molecular Biology

Although there is still some skepticism in the biological community regarding the value and significance of quantitative computational modeling, important steps are continually being taken to enhance its accessibility and predictive power. We view these developments as essential components of an emerging 'respectful modeling' framework which has two key aims: (i) respecting the models themselves and facilitating the reproduction and update of modeling results by other scientists, and (ii) respecting the predictions of the models and rigorously quantifying the confidence associated with the modeling results. This respectful attitude will guide the design of higher-quality models and facilitate the use of models in modern applications such as engineering and manipulating microbial metabolism by synthetic biology.

BCM: toolkit for Bayesian analysis of Computational Models using samplers

BACKGROUND

Computational models in biology are characterized by a large degree of uncertainty. This uncertainty can be analyzed with Bayesian statistics, however, the sampling algorithms that are frequently used for calculating Bayesian statistical estimates are computationally demanding, and each algorithm has unique advantages and disadvantages. It is typically unclear, before starting an analysis, which algorithm will perform well on a given computational model.

RESULTS

We present BCM, a toolkit for the Bayesian analysis of Computational Models using samplers. It provides efficient, multithreaded implementations of eleven algorithms for sampling from posterior probability distributions and for calculating marginal likelihoods. BCM includes tools to simplify the process of model specification and scripts for visualizing the results. The flexible architecture allows it to be used on diverse types of biological computational models. In an example inference task using a model of the cell cycle based on ordinary differential equations, BCM is significantly more efficient than existing software packages, allowing more challenging inference problems to be solved.

CONCLUSIONS

BCM represents an efficient one-stop-shop for computational modelers wishing to use sampler-based Bayesian statistics.

Tailored parameter optimization methods for ordinary differential equation models with steady-state constraints

BACKGROUND

Ordinary differential equation (ODE) models are widely used to describe (bio-)chemical and biological processes. To enhance the predictive power of these models, their unknown parameters are estimated from experimental data. These experimental data are mostly collected in perturbation experiments, in which the processes are pushed out of steady state by applying a stimulus. The information that the initial condition is a steady state of the unperturbed process provides valuable information, as it restricts the dynamics of the process and thereby the parameters. However, implementing steady-state constraints in the optimization often results in convergence problems.

RESULTS

In this manuscript, we propose two new methods for solving optimization problems with steady-state constraints. The first method exploits ideas from optimization algorithms on manifolds and introduces a retraction operator, essentially reducing the dimension of the optimization problem. The second method is based on the continuous analogue of the optimization problem. This continuous analogue is an ODE whose equilibrium points are the optima of the constrained optimization problem. This equivalence enables the use of adaptive numerical methods for solving optimization problems with steady-state constraints. Both methods are tailored to the problem structure and exploit the local geometry of the steady-state manifold and its stability properties. A parameterization of the steady-state manifold is not required. The efficiency and reliability of the proposed methods is evaluated using one toy example and two applications. The first application example uses published data while the second uses a novel dataset for Raf/MEK/ERK signaling. The proposed methods demonstrated better convergence properties than state-of-the-art methods employed in systems and computational biology. Furthermore, the average computation time per converged start is significantly lower. In addition to the theoretical results, the analysis of the dataset for Raf/MEK/ERK signaling provides novel biological insights regarding the existence of feedback regulation.

CONCLUSION

Many optimization problems considered in systems and computational biology are subject to steady-state constraints. While most optimization methods have convergence problems if these steady-state constraints are highly nonlinear, the methods presented recover the convergence properties of optimizers which can exploit an analytical expression for the parameter-dependent steady state. This renders them an excellent alternative to methods which are currently employed in systems and computational biology.

Construction of feasible and accurate kinetic models of metabolism: A Bayesian approach

Kinetic models are essential to quantitatively understand and predict the behaviour of metabolic networks. Detailed and thermodynamically feasible kinetic models of metabolism are inherently difficult to formulate and fit. They have a large number of heterogeneous parameters, are non-linear and have complex interactions. Many powerful fitting strategies are ruled out by the intractability of the likelihood function. Here, we have developed a computational framework capable of fitting feasible and accurate kinetic models using Approximate Bayesian Computation. This framework readily supports advanced modelling features such as model selection and model-based experimental design. We illustrate this approach on the tightly-regulated mammalian methionine cycle. Sampling from the posterior distribution, the proposed framework generated thermodynamically feasible parameter samples that converged on the true values, and displayed remarkable prediction accuracy in several validation tests. Furthermore, a posteriori analysis of the parameter distributions enabled appraisal of the systems properties of the network (e.g., control structure) and key metabolic regulations. Finally, the framework was used to predict missing allosteric interactions.

Synergistic combination of DT-13 and topotecan inhibits human gastric cancer via myosin IIA-induced endocytosis of EGF receptor in vitro and in vivo

Combination therapy has a higher success rate for many cancers compared to mono-therapy. The treatment of Topotecan (TPT) on gastric cancer (GC) is limited by its toxicity and the potential drug resistance. We found that the combination of the saponin monomer 13 from the dwarf lilyturf tuber (DT-13), performing anti-metastasis and anti-angiogenesis effects, with TPT synergistically induced apoptotic cytotoxicity in GCs with high EGF receptor (EGFR) expression, which was dependent on DT-13-induced endocytosis of EGFR. With TPT, DT-13 promoted EGFR ubiquitin--mediated degradation through myosin IIA-induced and Src/ caveolin-1 (Cav-1)-induced endocytosis of EGFR; inhibited EGFR downstream signalling and then increased the pro-apoptotic effects. Moreover, the synergistic pro-apoptotic efficacy of DT-13 and TPT in GCs with high EGFR expression was eliminated by both the NM II inhibitor (-)-blebbistatin and MYH-9 shRNA. The combination therapy of DT-13 with TPT showed stronger anti-tumour effects in vivo compared with their individual effects. Moreover, the results of combination therapy revealed selective upregulation of pro-apoptotic activity in TUNEL assays and cleaved caspase-3 and NM IIA in immunohischemical analysis; while specific downregulation of p-extracellular regulated kinase 1/2 (p-ERK1/2), EGFR and Cav-1 in immunohischemical analysis. Collectively, these findings have significant clinical implications for patients with tumours harbouring high EGFR expression due to the possible high sensitivity of this regimen.

Analyzing Th17 cell differentiation dynamics using a novel integrative modeling framework for time-course RNA sequencing data

Background

The differentiation of naive CD 4⁺ helper T (Th) cells into effector Th17 cells is steered by extracellular cytokines that activate and control the lineage specific transcriptional program. While the inducing cytokine signals and core transcription factors driving the differentiation towards Th17 lineage are well known, detailed mechanistic interactions between the key components are poorly understood.

Results

We develop an integrative modeling framework which combines RNA sequencing data with mathematical modeling and enables us to construct a mechanistic model for the core Th17 regulatory network in a data-driven manner.

Conclusions

Our results show significant evidence, for instance, for inhibitory mechanisms between the transcription factors and reveal a previously unknown dependency between the dosage of the inducing cytokine TGF β and the expression of the master regulator of competing (induced) regulatory T cell lineage. Further, our experimental validation approves this dependency in Th17 polarizing conditions.

Electronic supplementary material

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

Non-canonical Activation of Akt in Serum-Stimulated Fibroblasts, Revealed by Comparative Modeling of Pathway Dynamics

The dynamic behaviors of signaling pathways can provide clues to pathway mechanisms. In cancer cells, excessive phosphorylation and activation of the Akt pathway is responsible for cell survival advantages. In normal cells, serum stimulation causes brief peaks of extremely high Akt phosphorylation before reaching a moderate steady-state. Previous modeling assumed this peak and decline behavior (i.e., “overshoot”) was due to receptor internalization. In this work, we modeled the dynamics of the overshoot as a tool for gaining insight into Akt pathway function. We built an ordinary differential equation (ODE) model describing pathway activation immediately upstream of Akt phosphorylation at Thr³⁰⁸ (Aktp³⁰⁸). The model was fit to experimental measurements of Aktp³⁰⁸, total Akt, and phosphatidylinositol (3,4,5)-trisphosphate (PIP3), from mouse embryonic fibroblasts with serum stimulation. The canonical Akt activation model (the null hypothesis) was unable to recapitulate the observed delay between the peak of PIP3 (at 2 minutes), and the peak of Aktp³⁰⁸ (at 30–60 minutes). From this we conclude that the peak and decline behavior of Aktp³⁰⁸ is not caused by PIP3 dynamics. Models for alternative hypotheses were constructed by allowing an arbitrary dynamic curve to perturb each of 5 steps of the pathway. All 5 of the alternative models could reproduce the observed delay. To distinguish among the alternatives, simulations suggested which species and timepoints would show strong differences. Time-series experiments with membrane fractionation and PI3K inhibition were performed, and incompatible hypotheses were excluded. We conclude that the peak and decline behavior of Aktp³⁰⁸ is caused by a non-canonical effect that retains Akt at the membrane, and not by receptor internalization. Furthermore, we provide a novel spline-based method for simulating the network implications of an unknown effect, and we demonstrate a process of hypothesis management for guiding efficient experiments.

Influential pathways of cell signalling (such as PI3K/Akt) are routinely communicated using simple textbook-like diagrams that show only the most widely-accepted steps of the pathway. At the same time, there are countless other molecular influences relevant to each pathway, documented in the published literature, and more are being published every week. It should perhaps come as little surprise that during a routine observation of the Akt activation pathway, a simulation of the canonical model was mathematically incompatible with our observed dynamics. To progress beyond the standard, simplified model without testing an unreasonable number of molecular candidates individually, we employed computational modeling to analyze the dynamics of pathway activation. We asked when and where a non-canonical deviation could occur, relative to the canonical pathway. We used the timing of downstream activation to solve for the possible times of upstream initiation. By categorizing unknown effects by their dynamics, we were able to prune away implausible hypotheses using an efficient number of in vitro experiments. At the end we had a single plausible explanation for non-canonical Akt activation in our cells, and we confirmed experimentally that Akt is retained at the membrane after PIP3 is no longer present.

In silico model-based inference: an emerging approach for inverse problems in engineering better medicines

Identifying the network of biochemical interactions that underpin disease pathophysiology is a key hurdle in drug discovery. While many components involved in these biological processes are identified, how components organize differently in health and disease remains unclear. In chemical engineering, mechanistic modeling provides a quantitative framework to capture our understanding of a reactive system and test this knowledge against data. Here, we describe an emerging approach to test this knowledge against data that leverages concepts from probability, Bayesian statistics, and chemical kinetics by focusing on two related inverse problems. The first problem is to identify the causal structure of the reaction network, given uncertainty as to how the reactive components interact. The second problem is to identify the values of the model parameters, when a network is known a priori.

Efficient Characterization of Parametric Uncertainty of Complex (Bio)chemical Networks

Parametric uncertainty is a particularly challenging and relevant aspect of systems analysis in domains such as systems biology where, both for inference and for assessing prediction uncertainties, it is essential to characterize the system behavior globally in the parameter space. However, current methods based on local approximations or on Monte-Carlo sampling cope only insufficiently with high-dimensional parameter spaces associated with complex network models. Here, we propose an alternative deterministic methodology that relies on sparse polynomial approximations. We propose a deterministic computational interpolation scheme which identifies most significant expansion coefficients adaptively. We present its performance in kinetic model equations from computational systems biology with several hundred parameters and state variables, leading to numerical approximations of the parametric solution on the entire parameter space. The scheme is based on adaptive Smolyak interpolation of the parametric solution at judiciously and adaptively chosen points in parameter space. As Monte-Carlo sampling, it is “non-intrusive” and well-suited for massively parallel implementation, but affords higher convergence rates. This opens up new avenues for large-scale dynamic network analysis by enabling scaling for many applications, including parameter estimation, uncertainty quantification, and systems design.

In various scientific domains, in particular in systems biology, dynamic mathematical models of increasing complexity are being developed and analyzed to study biochemical reaction networks. A major challenge in dealing with such models is the uncertainty in parameters such as kinetic constants; how to efficiently and precisely quantify the effects of parametric uncertainties on systems behavior remains a question. Addressing this computational challenge for large systems, with good scaling up to hundreds of species and kinetic parameters, is important for many forward (e.g., uncertainty quantification) and inverse (e.g., system identification) problems. Here, we propose a sparse, deterministic adaptive interpolation method tailored to high-dimensional parametric problems that allows for fast, deterministic computational analysis of large biochemical reaction networks. The method is based on adaptive Smolyak interpolation of the parametric solution at judiciously chosen points in high-dimensional parameter space, combined with adaptive time-stepping for the actual numerical simulation of the network dynamics. It is “non-intrusive” and well-suited both for massively parallel implementation and for use in standard (systems biology) toolboxes.

Biased signalling: the instinctive skill of the cell in the selection of appropriate signalling pathways

GPCRs (G-protein-coupled receptors) are members of a family of proteins which are generally regarded as the largest group of therapeutic drug targets. Ligands of GPCRs do not usually activate all cellular signalling pathways linked to a particular seven-transmembrane receptor in a uniform manner. The fundamental idea behind this concept is that each ligand has its own ability, while interacting with the receptor, to activate different signalling pathways (or a particular set of signalling pathways) and it is this concept which is known as biased signalling. The importance of biased signalling is that it may selectively activate biological responses to favour therapeutically beneficial signalling pathways and to avoid adverse effects. There are two levels of biased signalling. First, bias can arise from the ability of GPCRs to couple to a subset of the available G-protein subtypes: Gαs, Gαq/11, Gαi/o or Gα12/13. These subtypes produce the diverse effects of GPCRs by targeting different effectors. Secondly, biased GPCRs may differentially activate G-proteins or β-arrestins. β-Arrestins are ubiquitously expressed and function to terminate or inhibit classic G-protein signalling and initiate distinct β-arrestin-mediated signalling processes. The interplay of G-protein and β-arrestin signalling largely determines the cellular consequences of the administration of GPCR-targeted drugs. In the present review, we highlight the particular functionalities of biased signalling and discuss its biological effects subsequent to GPCR activation. We consider that biased signalling is potentially allowing a choice between signalling through ‘beneficial’ pathways and the avoidance of ‘harmful’ ones.

Noise-driven causal inference in biomolecular networks

Single-cell RNA and protein concentrations dynamically fluctuate because of stochastic ("noisy") regulation. Consequently, biological signaling and genetic networks not only translate stimuli with functional response but also random fluctuations. Intuitively, this feature manifests as the accumulation of fluctuations from the network source to the target. Taking advantage of the fact that noise propagates directionally, we developed a method for causation prediction that does not require time-lagged observations and therefore can be applied to data generated by destructive assays such as immunohistochemistry. Our method for causation prediction, "Inference of Network Directionality Using Covariance Elements (INDUCE)," exploits the theoretical relationship between a change in the strength of a causal interaction and the associated changes in the single cell measured entries of the covariance matrix of protein concentrations. We validated our method for causation prediction in two experimental systems where causation is well established: in an E. coli synthetic gene network, and in MEK to ERK signaling in mammalian cells. We report the first analysis of covariance elements documenting noise propagation from a kinase to a phosphorylated substrate in an endogenous mammalian signaling network.

Disentangling the Complexity of HGF Signaling by Combining Qualitative and Quantitative Modeling

Signaling pathways are characterized by crosstalk, feedback and feedforward mechanisms giving rise to highly complex and cell-context specific signaling networks. Dissecting the underlying relations is crucial to predict the impact of targeted perturbations. However, a major challenge in identifying cell-context specific signaling networks is the enormous number of potentially possible interactions. Here, we report a novel hybrid mathematical modeling strategy to systematically unravel hepatocyte growth factor (HGF) stimulated phosphoinositide-3-kinase (PI3K) and mitogen activated protein kinase (MAPK) signaling, which critically contribute to liver regeneration. By combining time-resolved quantitative experimental data generated in primary mouse hepatocytes with interaction graph and ordinary differential equation modeling, we identify and experimentally validate a network structure that represents the experimental data best and indicates specific crosstalk mechanisms. Whereas the identified network is robust against single perturbations, combinatorial inhibition strategies are predicted that result in strong reduction of Akt and ERK activation. Thus, by capitalizing on the advantages of the two modeling approaches, we reduce the high combinatorial complexity and identify cell-context specific signaling networks.

The silicon trypanosome: a test case of iterative model extension in systems biology

The African trypanosome, Trypanosoma brucei, is a unicellular parasite causing African Trypanosomiasis (sleeping sickness in humans and nagana in animals). Due to some of its unique properties, it has emerged as a popular model organism in systems biology. A predictive quantitative model of glycolysis in the bloodstream form of the parasite has been constructed and updated several times. The Silicon Trypanosome is a project that brings together modellers and experimentalists to improve and extend this core model with new pathways and additional levels of regulation. These new extensions and analyses use computational methods that explicitly take different levels of uncertainty into account. During this project, numerous tools and techniques have been developed for this purpose, which can now be used for a wide range of different studies in systems biology.

Raf-interactome in tuning the complexity and diversity of Raf function

Raf kinases have been intensely studied subsequent to their discovery 30 years ago. The Ras-Raf-mitogen-activated protein kinase/extracellular signal-regulated kinase kinase-extracellular signal-regulated kinase/mitogen-activated protein kinase (Ras-Raf-MEK-ERK/MAPK) signaling pathway is at the heart of the signaling networks that control many fundamental cellular processes and Raf kinases takes centre stage in the MAPK pathway, which is now appreciated to be one of the most common sources of the oncogenic mutations in cancer. The dependency of tumors on this pathway has been clearly demonstrated by targeting its key nodes; however, blockade of the central components of the MAPK pathway may have some unexpected side effects. Over recent years, the Raf-interactome or Raf-interacting proteins have emerged as promising targets for protein-directed cancer therapy. This review focuses on the diversity of Raf-interacting proteins and discusses the mechanisms by which these proteins regulate Raf function, as well as the implications of targeting Raf-interacting proteins in the treatment of human cancer.

Cutting the wires: modularization of cellular networks for experimental design

Understanding naturally evolved cellular networks requires the consecutive identification and revision of the interactions between relevant molecular species. In this process, initially often simplified and incomplete networks are extended by integrating new reactions or whole subnetworks to increase consistency between model predictions and new measurement data. However, increased consistency with experimental data alone is not sufficient to show the existence of biomolecular interactions, because the interplay of different potential extensions might lead to overall similar dynamics. Here, we present a graph-based modularization approach to facilitate the design of experiments targeted at independently validating the existence of several potential network extensions. Our method is based on selecting the outputs to measure during an experiment, such that each potential network extension becomes virtually insulated from all others during data analysis. Each output defines a module that only depends on one hypothetical network extension, and all other outputs act as virtual inputs to achieve insulation. Given appropriate experimental time-series measurements of the outputs, our modules can be analyzed, simulated, and compared to the experimental data separately. Our approach exemplifies the close relationship between structural systems identification and modularization, an interplay that promises development of related approaches in the future.

Causal network inference using biochemical kinetics

MOTIVATION

Networks are widely used as structural summaries of biochemical systems. Statistical estimation of networks is usually based on linear or discrete models. However, the dynamics of biochemical systems are generally non-linear, suggesting that suitable non-linear formulations may offer gains with respect to causal network inference and aid in associated prediction problems.

RESULTS

We present a general framework for network inference and dynamical prediction using time course data that is rooted in non-linear biochemical kinetics. This is achieved by considering a dynamical system based on a chemical reaction graph with associated kinetic parameters. Both the graph and kinetic parameters are treated as unknown; inference is carried out within a Bayesian framework. This allows prediction of dynamical behavior even when the underlying reaction graph itself is unknown or uncertain. Results, based on (i) data simulated from a mechanistic model of mitogen-activated protein kinase signaling and (ii) phosphoproteomic data from cancer cell lines, demonstrate that non-linear formulations can yield gains in causal network inference and permit dynamical prediction and uncertainty quantification in the challenging setting where the reaction graph is unknown.

AVAILABILITY AND IMPLEMENTATION

MATLAB R2014a software is available to download from warwick.ac.uk/chrisoates.

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

Optimal experiment design for model selection in biochemical networks

Background

Mathematical modeling is often used to formalize hypotheses on how a biochemical network operates by discriminating between competing models. Bayesian model selection offers a way to determine the amount of evidence that data provides to support one model over the other while favoring simple models. In practice, the amount of experimental data is often insufficient to make a clear distinction between competing models. Often one would like to perform a new experiment which would discriminate between competing hypotheses.

Results

We developed a novel method to perform Optimal Experiment Design to predict which experiments would most effectively allow model selection. A Bayesian approach is applied to infer model parameter distributions. These distributions are sampled and used to simulate from multivariate predictive densities. The method is based on a k-Nearest Neighbor estimate of the Jensen Shannon divergence between the multivariate predictive densities of competing models.

Conclusions

We show that the method successfully uses predictive differences to enable model selection by applying it to several test cases. Because the design criterion is based on predictive distributions, which can be computed for a wide range of model quantities, the approach is very flexible. The method reveals specific combinations of experiments which improve discriminability even in cases where data is scarce. The proposed approach can be used in conjunction with existing Bayesian methodologies where (approximate) posteriors have been determined, making use of relations that exist within the inferred posteriors.

Reverse engineering and identification in systems biology: strategies, perspectives and challenges

The interplay of mathematical modelling with experiments is one of the central elements in systems biology. The aim of reverse engineering is to infer, analyse and understand, through this interplay, the functional and regulatory mechanisms of biological systems. Reverse engineering is not exclusive of systems biology and has been studied in different areas, such as inverse problem theory, machine learning, nonlinear physics, (bio)chemical kinetics, control theory and optimization, among others. However, it seems that many of these areas have been relatively closed to outsiders. In this contribution, we aim to compare and highlight the different perspectives and contributions from these fields, with emphasis on two key questions: (i) why are reverse engineering problems so hard to solve, and (ii) what methods are available for the particular problems arising from systems biology?

A framework for parameter estimation and model selection from experimental data in systems biology using approximate Bayesian computation

As modeling becomes a more widespread practice in the life sciences and biomedical sciences, researchers need reliable tools to calibrate models against ever more complex and detailed data. Here we present an approximate Bayesian computation (ABC) framework and software environment, ABC-SysBio, which is a Python package that runs on Linux and Mac OS X systems and that enables parameter estimation and model selection in the Bayesian formalism by using sequential Monte Carlo (SMC) approaches. We outline the underlying rationale, discuss the computational and practical issues and provide detailed guidance as to how the important tasks of parameter inference and model selection can be performed in practice. Unlike other available packages, ABC-SysBio is highly suited for investigating, in particular, the challenging problem of fitting stochastic models to data. In order to demonstrate the use of ABC-SysBio, in this protocol we postulate the existence of an imaginary reaction network composed of seven interrelated biological reactions (involving a specific mRNA, the protein it encodes and a post-translationally modified version of the protein), a network that is defined by two files containing 'observed' data that we provide as supplementary information. In the first part of the PROCEDURE, ABC-SysBio is used to infer the parameters of this system, whereas in the second part we use ABC-SysBio's relevant functionality to discriminate between two different reaction network models, one of them being the 'true' one. Although computationally expensive, the additional insights gained in the Bayesian formalism more than make up for this cost, especially in complex problems.

Topological augmentation to infer hidden processes in biological systems

MOTIVATION

A common problem in understanding a biochemical system is to infer its correct structure or topology. This topology consists of all relevant state variables-usually molecules and their interactions. Here we present a method called topological augmentation to infer this structure in a statistically rigorous and systematic way from prior knowledge and experimental data.

RESULTS

Topological augmentation starts from a simple model that is unable to explain the experimental data and augments its topology by adding new terms that capture the experimental behavior. This process is guided by representing the uncertainty in the model topology through stochastic differential equations whose trajectories contain information about missing model parts. We first apply this semiautomatic procedure to a pharmacokinetic model. This example illustrates that a global sampling of the parameter space is critical for inferring a correct model structure. We also use our method to improve our understanding of glutamine transport in yeast. This analysis shows that transport dynamics is determined by glutamine permeases with two different kinds of kinetics. Topological augmentation can not only be applied to biochemical systems, but also to any system that can be described by ordinary differential equations.

AVAILABILITY AND IMPLEMENTATION

Matlab code and examples are available at: http://www.csb.ethz.ch/tools/index

Deregulated expression of TANK in glioblastomas triggers pro-tumorigenic ERK1/2 and AKT signaling pathways

Signal transmission by the noncanonical IkappaB kinases (IKKs), TANK-binding kinase 1 (TBK1) and IKKɛ, requires interaction with adapter proteins such as TRAF associated NF-κB activator (TANK). Although increased expression or dysregulation of both kinases has been described for a variety of human cancers, this study shows that deregulated expression of the TANK protein is frequently occurring in glioblastomas (GBMs). The functional relevance of TANK was analyzed in a panel of GBM-derived cell lines and revealed that knockdown of TANK arrests cells in the S-phase and prohibits tumor cell migration. Deregulated TANK expression affects several signaling pathways controlling cell proliferation and the inflammatory response. Interference with stoichiometrically assembled signaling complexes by overexpression or silencing of TANK prevented constitutive interferon-regulatory factor 3 (IRF3) phosphorylation. Knockdown of TANK frequently prevents constitutive activation of extracellular signal-regulated kinases 1 and 2 (ERK1/2). TANK-mediated ERK1/2 activation is independent from the canonical MAP kinase or ERK kinase (MEK) 1/2-mediated pathway and utilizes an alternative pathway that uses a TBK1/IKKɛ/Akt signaling axis, thus identifying a novel pathway suitable to block constitutive ERK1/2 activity.

Near-optimal experimental design for model selection in systems biology

MOTIVATION

Biological systems are understood through iterations of modeling and experimentation. Not all experiments, however, are equally valuable for predictive modeling. This study introduces an efficient method for experimental design aimed at selecting dynamical models from data. Motivated by biological applications, the method enables the design of crucial experiments: it determines a highly informative selection of measurement readouts and time points.

RESULTS

We demonstrate formal guarantees of design efficiency on the basis of previous results. By reducing our task to the setting of graphical models, we prove that the method finds a near-optimal design selection with a polynomial number of evaluations. Moreover, the method exhibits the best polynomial-complexity constant approximation factor, unless P = NP. We measure the performance of the method in comparison with established alternatives, such as ensemble non-centrality, on example models of different complexity. Efficient design accelerates the loop between modeling and experimentation: it enables the inference of complex mechanisms, such as those controlling central metabolic operation.

AVAILABILITY

Toolbox 'NearOED' available with source code under GPL on the Machine Learning Open Source Software Web site (mloss.org).

Programming biological models in Python using PySB

PySB is a framework for creating biological models as Python programs using a high-level, action-oriented vocabulary that promotes transparency, extensibility and reusability. PySB interoperates with many existing modeling tools and supports distributed model development.

PySB models are programs and leverage existing programming tools for documentation, testing, and collaborative development.

Reusable functions can encode common low-level biochemical processes as well as high-level modules, making models transparent and concise.

Modeling workflow is accelerated through close integration with Python numerical tools and interoperability with existing modeling software.

We demonstrate the use of PySB to encode 15 alternative hypotheses for the mitochondrial regulation of apoptosis, including a new ‘Embedded Together' model based on recent biochemical findings.

Mathematical equations are fundamental to modeling biological networks, but as networks get large and revisions frequent, it becomes difficult to manage equations directly or to combine previously developed models. Multiple simultaneous efforts to create graphical standards, rule-based languages, and integrated software workbenches aim to simplify biological modeling but none fully meets the need for transparent, extensible, and reusable models. In this paper we describe PySB, an approach in which models are not only created using programs, they are programs. PySB draws on programmatic modeling concepts from little b and ProMot, the rule-based languages BioNetGen and Kappa and the growing library of Python numerical tools. Central to PySB is a library of macros encoding familiar biochemical actions such as binding, catalysis, and polymerization, making it possible to use a high-level, action-oriented vocabulary to construct detailed models. As Python programs, PySB models leverage tools and practices from the open-source software community, substantially advancing our ability to distribute and manage the work of testing biochemical hypotheses. We illustrate these ideas using new and previously published models of apoptosis.

In silico modeling of Itk activation kinetics in thymocytes suggests competing positive and negative IP4 mediated feedbacks increase robustness

The inositol-phosphate messenger inositol(1,3,4,5)tetrakisphosphate (IP₄) is essential for thymocyte positive selection by regulating plasma-membrane association of the protein tyrosine kinase Itk downstream of the T cell receptor (TCR). IP₄ can act as a soluble analog of the phosphoinositide 3-kinase (PI3K) membrane lipid product phosphatidylinositol(3,4,5)trisphosphate (PIP₃). PIP₃ recruits signaling proteins such as Itk to cellular membranes by binding to PH and other domains. In thymocytes, low-dose IP₄ binding to the Itk PH domain surprisingly promoted and high-dose IP₄ inhibited PIP₃ binding of Itk PH domains. However, the mechanisms that underlie the regulation of membrane recruitment of Itk by IP₄ and PIP₃ remain unclear. The distinct Itk PH domain ability to oligomerize is consistent with a cooperative-allosteric mode of IP₄ action. However, other possibilities cannot be ruled out due to difficulties in quantitatively measuring the interactions between Itk, IP₄ and PIP₃, and in generating non-oligomerizing Itk PH domain mutants. This has hindered a full mechanistic understanding of how IP₄ controls Itk function. By combining experimentally measured kinetics of PLCγ1 phosphorylation by Itk with in silico modeling of multiple Itk signaling circuits and a maximum entropy (MaxEnt) based computational approach, we show that those in silico models which are most robust against variations of protein and lipid expression levels and kinetic rates at the single cell level share a cooperative-allosteric mode of Itk regulation by IP₄ involving oligomeric Itk PH domains at the plasma membrane. This identifies MaxEnt as an excellent tool for quantifying robustness for complex TCR signaling circuits and provides testable predictions to further elucidate a controversial mechanism of PIP₃ signaling.

Model selection in systems and synthetic biology

Developing mechanistic models has become an integral aspect of systems biology, as has the need to differentiate between alternative models. Parameterizing mathematical models has been widely perceived as a formidable challenge, which has spurred the development of statistical and optimisation routines for parameter inference. But now focus is increasingly shifting to problems that require us to choose from among a set of different models to determine which one offers the best description of a given biological system. We will here provide an overview of recent developments in the area of model selection. We will focus on approaches that are both practical as well as build on solid statistical principles and outline the conceptual foundations and the scope for application of such methods in systems biology.