Division of labor by dual feedback regulators controls JAK2/STAT5 signaling over broad ligand range

Mechanistic modeling of cell viability assays with in silico lineage tracing

Data from cell viability assays, which measure cumulative division and death events in a population and reflect substantial cellular heterogeneity, are widely available. However, interpreting such data with mechanistic computational models is hindered because direct model/data comparison is often muddled. We developed an algorithm that tracks simulated division and death events in mechanistically detailed single-cell lineages to enable such a model/data comparison and suggest causes of cell-cell drug response variability. Using our previously developed model of mammalian single-cell proliferation and death signaling, we simulated drug dose response experiments for four targeted anti-cancer drugs (alpelisib, neratinib, trametinib and palbociclib) and compared them to experimental data. Simulations are consistent with data for strong growth inhibition by trametinib (MEK inhibitor) and overall lack of efficacy for alpelisib (PI-3K inhibitor), but are inconsistent with data for palbociclib (CDK4/6 inhibitor) and neratinib (EGFR inhibitor). Model/data inconsistencies suggest (i) the importance of CDK4/6 for driving the cell cycle may be overestimated, and (ii) that the cellular balance between basal (tonic) and ligand-induced signaling is a critical determinant of receptor inhibitor response. Simulations show subpopulations of rapidly and slowly dividing cells in both control and drug-treated conditions. Variations in mother cells prior to drug treatment all impinging on ERK pathway activity are associated with the rapidly dividing phenotype and trametinib resistance. This work lays a foundation for the application of mechanistic modeling to large-scale cell viability assay datasets and better understanding determinants of cellular heterogeneity in drug response.

Canonical and Non-Canonical Functions of Erythropoietin and Its Receptor in Mature Nucleated Erythrocytes of Western Clawed Frog, Xenopus tropicalis

Enucleated erythrocytes are characteristic of adult mammals. In contrast, fish, amphibians, reptiles, birds, and fetal mammals possess nucleated erythrocytes in their circulation. Erythroid maturation is regulated by erythropoietin (EPO) and its receptor (EPOR), which are conserved among vertebrates. In mammals, EPOR on the erythroid progenitor membrane disappears after terminal differentiation. However, in western clawed frog, Xenopus tropicalis, mature erythrocytes maintain EPOR expression, suggesting that they have non-canonical functions of the EPO-EPOR axis rather than proliferation and differentiation. In this study, we investigated the non-canonical functions of EPOR in Xenopus mature erythrocytes. EPO stimulation of peripheral erythrocytes did not induce proliferation but induced phosphorylation of intracellular proteins, including signal transducer and activator of transcription 5 (STAT5). RNA-Seq analysis of EPO-stimulated peripheral erythrocytes identified 45 differentially expressed genes (DEGs), including cytokine inducible SH2 containing protein gene (cish) and suppressor of cytokine signaling 3 gene (socs3), negative regulators of the EPOR-Janus kinase (JAK)-STAT pathway. These phosphorylation studies and pathway analysis demonstrated the activation of the JAK-STAT pathway through EPO-EPOR signaling in erythrocytes. Through comparison with EPO-responsive genes in mouse erythroid progenitors obtained from a public database, we identified 31 novel EPO-responsive genes indicating non-canonical functions. Among these, we focused on ornithine decarboxylase 1 gene (odc1), which is the rate-limiting enzyme in polyamine synthesis and affects hematopoietic progenitor differentiation and the endothelial cell response to hypoxic stress. An EPO-supplemented culture of erythrocytes showed increased odc1 expression followed by a decrease in polyamine-rich erythrocytes, suggesting EPO-responsive polyamine excretion. These findings will advance our knowledge of the unknown regulatory systems under the EPO-EPOR axis and functional differences between vertebrates' nucleated and enucleated erythrocytes.

Signal execution modes emerge in biochemical reaction networks calibrated to experimental data

Mathematical models of biomolecular networks are commonly used to study cellular processes; however, their usefulness to explain and predict dynamic behaviors is often questioned due to the unclear relationship between parameter uncertainty and network dynamics. In this work, we introduce PyDyNo (Python dynamic analysis of biochemical networks), a non-equilibrium reaction-flux based analysis to identify dominant reaction paths within a biochemical reaction network calibrated to experimental data. We first show, in a simplified apoptosis execution model, that despite the thousands of parameter vectors with equally good fits to experimental data, our framework identifies the dynamic differences between these parameter sets and outputs three dominant execution modes, which exhibit varying sensitivity to perturbations. We then apply our methodology to JAK2/STAT5 network in colony-forming unit-erythroid (CFU-E) cells and provide previously unrecognized mechanistic explanation for the survival responses of CFU-E cell population that would have been impossible to deduce with traditional protein-concentration based analyses.

Controllability and accessibility analysis of nonlinear biosystems

BACKGROUND

We address the problem of determining the controllability and accessibility of nonlinear biosystems. We consider models described by affine-in-inputs ordinary differential equations, which are adequate for a wide array of biological processes. Roughly speaking, the controllability of a dynamical system determines the possibility of steering it from an initial state to any point in its neighbourhood; accessibility is a weaker form of controllability.

METHODS

While the methodology for analysing the controllability of linear systems is well established, its generalization to the nonlinear case has proven elusive. Thus, a number of related but different properties - including different versions of accessibility, reachability or weak local controllability - have been defined to approach its study, and several partial results exist in lieu of a general test. Here, leveraging the applicable results from differential geometric control theory, we source sufficient conditions to assess nonlinear controllability, as well as a necessary and sufficient condition for accessibility.

RESULTS

We develop an algorithmic procedure to evaluate these conditions efficiently, and we provide its open source implementation. Using this software tool, we analyse the accessibility and controllability of a number of models of biomedical interest. While some of them are fully controllable, we find others that are not, as is the case of some models of EGF and NFκB signalling networks.

CONCLUSIONS

The contributions in this paper facilitate the accessibility and controllability analysis of nonlinear models, not only in biomedicine but also in other areas in which they have been rarely performed to date.

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.

Uncovering specific mechanisms across cell types in dynamical models

Ordinary differential equations are frequently employed for mathematical modeling of biological systems. The identification of mechanisms that are specific to certain cell types is crucial for building useful models and to gain insights into the underlying biological processes. Regularization techniques have been proposed and applied to identify mechanisms specific to two cell types, e.g., healthy and cancer cells, including the LASSO (least absolute shrinkage and selection operator). However, when analyzing more than two cell types, these approaches are not consistent, and require the selection of a reference cell type, which can affect the results. To make the regularization approach applicable to identifying cell-type specific mechanisms in any number of cell types, we propose to incorporate the clustered LASSO into the framework of ordinary differential equation modeling by penalizing the pairwise differences of the logarithmized fold-change parameters encoding a specific mechanism in different cell types. The symmetry introduced by this approach renders the results independent of the reference cell type. We discuss the necessary adaptations of state-of-the-art numerical optimization techniques and the process of model selection for this method. We assess the performance with realistic biological models and synthetic data, and demonstrate that it outperforms existing approaches. Finally, we also exemplify its application to published biological models including experimental data, and link the results to independent biological measurements.

Likelihood-ratio test statistic for the finite-sample case in nonlinear ordinary differential equation models

Likelihood ratios are frequently utilized as basis for statistical tests, for model selection criteria and for assessing parameter and prediction uncertainties, e.g. using the profile likelihood. However, translating these likelihood ratios into p-values or confidence intervals requires the exact form of the test statistic's distribution. The lack of knowledge about this distribution for nonlinear ordinary differential equation (ODE) models requires an approximation which assumes the so-called asymptotic setting, i.e. a sufficiently large amount of data. Since the amount of data from quantitative molecular biology is typically limited in applications, this finite-sample case regularly occurs for mechanistic models of dynamical systems, e.g. biochemical reaction networks or infectious disease models. Thus, it is unclear whether the standard approach of using statistical thresholds derived for the asymptotic large-sample setting in realistic applications results in valid conclusions. In this study, empirical likelihood ratios for parameters from 19 published nonlinear ODE benchmark models are investigated using a resampling approach for the original data designs. Their distributions are compared to the asymptotic approximation and statistical thresholds are checked for conservativeness. It turns out, that corrections of the likelihood ratios in such finite-sample applications are required in order to avoid anti-conservative results.

Realistic simulation of time-course measurements in systems biology

In systems biology, the analysis of complex nonlinear systems faces many methodological challenges. For the evaluation and comparison of the performances of novel and competing computational methods, one major bottleneck is the availability of realistic test problems. We present an approach for performing realistic simulation studies for analyses of time course data as they are typically measured in systems biology. Since the design of experiments in practice depends on the process of interest, our approach considers the size and the dynamics of the mathematical model which is intended to be used for the simulation study. To this end, we used 19 published systems biology models with experimental data and evaluated the relationship between model features (e.g., the size and the dynamics) and features of the measurements such as the number and type of observed quantities, the number and the selection of measurement times, and the magnitude of measurement errors. Based on these typical relationships, our novel approach enables suggestions of realistic simulation study designs in the systems biology context and the realistic generation of simulated data for any dynamic model. The approach is demonstrated on three models in detail and its performance is validated on nine models by comparing ODE integration, parameter optimization, and parameter identifiability. The presented approach enables more realistic and less biased benchmark studies and thereby constitutes an important tool for the development of novel methods for dynamic modeling.

Benchmarking tools for a priori identifiability analysis

MOTIVATION

The theoretical possibility of determining the state and parameters of a dynamic model by measuring its outputs is given by its structural identifiability and its observability. These properties should be analysed before attempting to calibrate a model, but their a priori analysis can be challenging, requiring symbolic calculations that often have a high computational cost. In recent years, a number of software tools have been developed for this task, mostly in the systems biology community. These tools have vastly different features and capabilities, and a critical assessment of their performance is still lacking.

RESULTS

Here, we present a comprehensive study of the computational resources available for analysing structural identifiability. We consider 13 software tools developed in 7 programming languages and evaluate their performance using a set of 25 case studies created from 21 models. Our results reveal their strengths and weaknesses, provide guidelines for choosing the most appropriate tool for a given problem and highlight opportunities for future developments.

AVAILABILITY AND IMPLEMENTATION

https://github.com/Xabo-RB/Benchmarking_files.

Stochastic dynamics of Type-I interferon responses

Interferon (IFN) activates the transcription of several hundred of IFN stimulated genes (ISGs) that constitute a highly effective antiviral defense program. Cell-to-cell variability in the induction of ISGs is well documented, but its source and effects are not completely understood. The molecular mechanisms behind this heterogeneity have been related to randomness in molecular events taking place during the JAK-STAT signaling pathway. Here, we study the sources of variability in the induction of the IFN-alpha response by using MxA and IFIT1 activation as read-out. To this end, we integrate time-resolved flow cytometry data and stochastic modeling of the JAK-STAT signaling pathway. The complexity of the IFN response was matched by fitting probability distributions to time-course flow cytometry snapshots. Both, experimental data and simulations confirmed that the MxA and IFIT1 induction circuits generate graded responses rather than all-or-none responses. Subsequently, we quantify the size of the intrinsic variability at different steps in the pathway. We found that stochastic effects are transiently strong during the ligand-receptor activation steps and the formation of the ISGF3 complex, but negligible for the final induction of the studied ISGs. We conclude that the JAK-STAT signaling pathway is a robust biological circuit that efficiently transmits information under stochastic environments.

We investigate the impact of intrinsic and extrinsic noise on the reliability of interferon signaling. Information must be transduced robustly despite existing biochemical variability and at the same time the system has to allow for cellular variability to tune it against changing environments. Getting insights into stochasticity in signaling networks is crucial to understand cellular dynamics and decision-making processes. To this end, we developed a detailed stochastic computational model based on single cell data. We are able to show that reliability is achieved despite high noise at the receptor level.

Erythropoietin-driven dynamic proteome adaptations during erythropoiesis prevent iron overload in the developing embryo

Erythropoietin (Epo) ensures survival and proliferation of colony-forming unit erythroid (CFU-E) progenitor cells and their differentiation to hemoglobin-containing mature erythrocytes. A lack of Epo-induced responses causes embryonic lethality, but mechanisms regulating the dynamic communication of cellular alterations to the organismal level remain unresolved. By time-resolved transcriptomics and proteomics, we show that Epo induces in CFU-E cells a gradual transition from proliferation signature proteins to proteins indicative for differentiation, including heme-synthesis enzymes. In the absence of the Epo receptor (EpoR) in embryos, we observe a lack of hemoglobin in CFU-E cells and massive iron overload of the fetal liver pointing to a miscommunication between liver and placenta. A reduction of iron-sulfur cluster-containing proteins involved in oxidative phosphorylation in these embryos leads to a metabolic shift toward glycolysis. This link connecting erythropoiesis with the regulation of iron homeostasis and metabolic reprogramming suggests that balancing these interactions is crucial for protection from iron intoxication and for survival.

BlotIt—Optimal alignment of Western blot and qPCR experiments

Biological systems are frequently analyzed by means of mechanistic mathematical models. In order to infer model parameters and provide a useful model that can be employed for systems understanding and hypothesis testing, the model is often calibrated on quantitative, time-resolved data. To do so, it is typically important to compare experimental measurements over broad time ranges and various experimental conditions, e.g. perturbations of the biological system. However, most of the established experimental techniques such as Western blot, or quantitative real-time polymerase chain reaction only provide measurements on a relative scale, since different sample volumes, experimental adjustments or varying development times of a gel lead to systematic shifts in the data. In turn, the number of measurements corresponding to the same scale enabling comparability is limited. Here, we present a new flexible method to align measurement data that obeys different scaling factors and compare it to existing normalization approaches. We propose an alignment model to estimate these scaling factors and provide the possibility to adapt this model depending on the measurement technique of interest. In addition, an error model can be specified to adequately weight the different data points and obtain scaling-model based confidence intervals of the finally scaled data points. Our approach is applicable to all sorts of relative measurements and does not need a particular experimental condition that has been measured over all available scales. An implementation of the method is provided with the R package blotIt including refined ways of visualization.

Modelling of glucose repression signalling in yeast Saccharomyces cerevisiae

Saccharomyces cerevisiae has a sophisticated signalling system that plays a crucial role in cellular adaptation to changing environments. The SNF1 pathway regulates energy homeostasis upon glucose derepression; hence, it plays an important role in various processes, such as metabolism, cell cycle and autophagy. To unravel its behaviour, SNF1 signalling has been extensively studied. However, the pathway components are strongly interconnected and inconstant; therefore, elucidating its dynamic behaviour based on experimental data only is challenging. To tackle this complexity, systems biology approaches have been successfully employed. This review summarizes the progress, advantages and disadvantages of the available mathematical modelling frameworks covering Boolean, dynamic kinetic, single-cell models, which have been used to study processes and phenomena ranging from crosstalks to sources of cell-to-cell variability in the context of SNF1 signalling. Based on the lessons from existing models, we further discuss how to develop a consensus dynamic mechanistic model of the entire SNF1 pathway that can provide novel insights into the dynamics of nutrient signalling.

A Quantitative Systems Pharmacology Perspective on the Importance of Parameter Identifiability

There is an inherent tension in Quantitative Systems Pharmacology (QSP) between the need to incorporate mathematical descriptions of complex physiology and drug targets with the necessity of developing robust, predictive and well-constrained models. In addition to this, there is no “gold standard” for model development and assessment in QSP. Moreover, there can be confusion over terminology such as model and parameter identifiability; complex and simple models; virtual populations; and other concepts, which leads to potential miscommunication and misapplication of methodologies within modeling communities, both the QSP community and related disciplines. This perspective article highlights the pros and cons of using simple (often identifiable) vs. complex (more physiologically detailed but often non-identifiable) models, as well as aspects of parameter identifiability, sensitivity and inference methodologies for model development and analysis. The paper distills the central themes of the issue of identifiability and optimal model size and discusses open challenges.

Mini-batch optimization enables training of ODE models on large-scale datasets

Quantitative dynamic models are widely used to study cellular signal processing. A critical step in modelling is the estimation of unknown model parameters from experimental data. As model sizes and datasets are steadily growing, established parameter optimization approaches for mechanistic models become computationally extremely challenging. Mini-batch optimization methods, as employed in deep learning, have better scaling properties. In this work, we adapt, apply, and benchmark mini-batch optimization for ordinary differential equation (ODE) models, thereby establishing a direct link between dynamic modelling and machine learning. On our main application example, a large-scale model of cancer signaling, we benchmark mini-batch optimization against established methods, achieving better optimization results and reducing computation by more than an order of magnitude. We expect that our work will serve as a first step towards mini-batch optimization tailored to ODE models and enable modelling of even larger and more complex systems than what is currently possible.

Relationship Between Dimensionality and Convergence of Optimization Algorithms: A Comparison Between Data-Driven Normalization and Scaling Factor-Based Methods Using PEPSSBI

Ordinary differential equation models are used to represent intracellular signaling pathways in silico, aiding and guiding biological experiments to elucidate intracellular regulation. To construct such quantitative and predictive models of intracellular interactions, unknown parameters need to be estimated. Most of biological data are expressed in relative or arbitrary units, raising the question of how to compare model simulations with data. It has recently been shown that for models with large number of unknown parameters, fitting algorithms using a data-driven normalization of the simulations (DNS) performs best in terms of the convergence time and parameter identifiability. DNS approach compares model simulations and corresponding data both normalized by the same normalization procedure, without requiring additional parameters to be estimated, as necessary for widely used scaling factor-based methods. However, currently there is no parameter estimation software that directly supports DNS. In this chapter, we show how to apply DNS to dynamic models of systems and synthetic biology using PEPSSBI (Parameter Estimation Pipeline for Systems and Synthetic Biology). PEPSSBI is the first software that supports DNS, through algorithmically supported data normalization and objective function construction. PEPSSBI also supports model import using SBML and repeated parameter estimation runs executed in parallel either on a personal computer or a multi-CPU cluster.

Efficient gradient-based parameter estimation for dynamic models using qualitative data

MOTIVATION

Unknown parameters of dynamical models are commonly estimated from experimental data. However, while various efficient optimization and uncertainty analysis methods have been proposed for quantitative data, methods for qualitative data are rare and suffer from bad scaling and convergence.

RESULTS

Here, we propose an efficient and reliable framework for estimating the parameters of ordinary differential equation models from qualitative data. In this framework, we derive a semi-analytical algorithm for gradient calculation of the optimal scaling method developed for qualitative data. This enables the use of efficient gradient-based optimization algorithms. We demonstrate that the use of gradient information improves performance of optimization and uncertainty quantification on several application examples. On average, we achieve a speedup of more than one order of magnitude compared to gradient-free optimization. In addition, in some examples, the gradient-based approach yields substantially improved objective function values and quality of the fits. Accordingly, the proposed framework substantially improves the parameterization of models from qualitative data.

AVAILABILITY AND IMPLEMENTATION

The proposed approach is implemented in the open-source Python Parameter EStimation TOolbox (pyPESTO). pyPESTO is available at https://github.com/ICB-DCM/pyPESTO. All application examples and code to reproduce this study are available at https://doi.org/10.5281/zenodo.4507613.

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

Cell-to-cell variability in JAK2/STAT5 pathway components and cytoplasmic volumes defines survival threshold in erythroid progenitor cells

Survival or apoptosis is a binary decision in individual cells. However, at the cell-population level, a graded increase in survival of colony-forming unit-erythroid (CFU-E) cells is observed upon stimulation with erythropoietin (Epo). To identify components of Janus kinase 2/signal transducer and activator of transcription 5 (JAK2/STAT5) signal transduction that contribute to the graded population response, we extended a cell-population-level model calibrated with experimental data to study the behavior in single cells. The single-cell model shows that the high cell-to-cell variability in nuclear phosphorylated STAT5 is caused by variability in the amount of Epo receptor (EpoR):JAK2 complexes and of SHP1, as well as the extent of nuclear import because of the large variance in the cytoplasmic volume of CFU-E cells. 24-118 pSTAT5 molecules in the nucleus for 120 min are sufficient to ensure cell survival. Thus, variability in membrane-associated processes is sufficient to convert a switch-like behavior at the single-cell level to a graded population-level response.

IFN-γ inhibits ovarian cancer progression via SOCS1/JAK/STAT signaling pathway

PURPOSE

Ovarian cancer (OC) is a common malignancy, and IFN-γ, a multifunctional cytokine, is unveiled to impede the multiplication and enhance apoptosis in diverse tumor cells in previous research. Nonetheless, its function and mechanism in OC are blurred.

METHODS

OC cell lines SKOV3 and OVCAR3 were dealt with different concentrations (0-40 ng/ml) of IFN-γ. CCK-8 experiment was utilized to examine cell multiplication; Flow cytometry was executed to detect apoptosis and cell cycle; Wound healing assay was utilized to detect cell migration; and Transwell experiment was implemented to examine cell invasion. qRT-PCR analysis was applied to detect STAT5, STAT3, JAK2 and JAK3 mRNA expression in OC cell lines. Western blot experiment was applied to detect the protein and phosphorylation levels of SOCS1, STAT5 and STAT3.

RESULTS

IFN-γ suppressed OC cell multiplication in a concentration-dependent manner. Relative to the control group, IFN-γ restrained OC cell migration, invasion, enhanced apoptosis and prevented cell transformation from G0/G1 to S phase. Further analysis revealed that IFN-γ up-modulated SOCS1 expression and impeded STAT5 and STAT3 protein phosphorylation levels, and knockdown of SOCS1 partially counteracted the inhibitory effect of IFN-γ on STAT5 and STAT3 protein phosphorylation levels.

CONCLUSION

IFN-γ represses OC progression by facilitating SOCS1 to suppress STAT3 and STAT5 protein phosphorylation.

Efficient gradient-based parameter estimation for dynamic models using qualitative data

Motivation

Unknown parameters of dynamical models are commonly estimated from experimental data. However, while various efficient optimization and uncertainty analysis methods have been proposed for quantitative data, methods for qualitative data are rare and suffer from bad scaling and convergence.

Results

Here, we propose an efficient and reliable framework for estimating the parameters of ordinary differential equation models from qualitative data. In this framework, we derive a semi-analytical algorithm for gradient calculation of the optimal scaling method developed for qualitative data. This enables the use of efficient gradient-based optimization algorithms. We demonstrate that the use of gradient information improves performance of optimization and uncertainty quantification on several application examples. On average, we achieve a speedup of more than one order of magnitude compared to gradient-free optimization. In addition, in some examples, the gradient-based approach yields substantially improved objective function values and quality of the fits. Accordingly, the proposed framework substantially improves the parameterization of models from qualitative data.

Availability and implementation

The proposed approach is implemented in the open-source Python Parameter EStimation TOolbox (pyPESTO). pyPESTO is available at https://github.com/ICB-DCM/pyPESTO. All application examples and code to reproduce this study are available at https://doi.org/10.5281/zenodo.4507613.

Supplementary information

Supplementary data are available at Bioinformatics online.

STAT5 as a Key Protein of Erythropoietin Signalization

Erythropoietin (EPO) acts on multiple tissues through its receptor EPOR, a member of a cytokine class I receptor superfamily with pleiotropic effects. The interaction of EPO and EPOR triggers the activation of several signaling pathways that induce erythropoiesis, including JAK2/STAT5, PI3K/AKT, and MAPK. The canonical EPOR/JAK2/STAT5 pathway is a known regulator of differentiation, proliferation, and cell survival of erythroid progenitors. In addition, its role in the protection of other cells, including cancer cells, is under intense investigation. The involvement of EPOR/JAK2/STAT5 in other processes such as mRNA splicing, cytoskeleton reorganization, and cell metabolism has been recently described. The transcriptomics, proteomics, and epigenetic studies reviewed in this article provide a detailed understanding of EPO signalization. Advances in this area of research may be useful for improving the efficacy of EPO therapy in hematologic disorders, as well as in cancer treatment.

Increased expression of IFN-γ in preeclampsia impairs human trophoblast invasion via a SOCS1/JAK/STAT1 feedback loop

The weakening of extravillous trophoblast (EVT) invasion results in shallow placenta implantation. In HTR8/SVneo cells, IFN-γ can activate STAT1 and reduce cell invasion, and suppressor of cytokine signaling (SOCS) is an important negative regulatory protein in the Janus kinase (JAK)/STAT activator pathway and has a negative feedback function on JAK/STAT1. The aim of the present study was to elucidate how SOCS1 feedback regulates JAK/STAT1 and affects EVT cell invasion, which in turn affects the development of preeclampsia (PE). MTT and Annexin V/phosphatidylserine (PS) assays were performed to evaluate the viability and apoptosis of HTR8/SVneo cells treated with IFN-γ, respectively. Wound healing and invasion assays were also conducted to measure the migratory and invasive abilities of IFN-γ-treated HTR8/SVneo cells. The mRNA and protein expression levels of genes were detected using reverse transcription-quantitative PCR and western blot analysis. Small interfering RNA knockdown of SOCS1 was used to verify the role of feedback regulation in the IFN-γ-activated JAK/STAT1 signaling pathway. IFN-γ can inhibit HTR8/SVneo migration and invasion, and promote apoptosis by increasing the expression of phosphorylated (p)-JAK, p-STAT1 and caspase3, and reducing the expression of platelet-derived growth factor receptor A and Ezrin. Furthermore, SOCS1 may negatively regulate JAK/STAT1 and affect HTR-8/SVneo invasiveness. Evaluation of clinical samples demonstrated that the expression levels of SOCS1 and IFN-γ were higher in patients with PE compared with the healthy group. Collectively, the present results indicated that IFN-γ reduced the invasion of HTR-8/SVneo cells by activating JAK/STAT1, concurrently leading to an increase in SOCS1, which negatively regulates JAK/STAT1 and eliminates the pro-inflammatory effects of IFN-γ, thus forming a feedback loop.

Nonlinear Observability Algorithms with Known and Unknown Inputs: Analysis and Implementation

The observability of a dynamical system is affected by the presence of external inputs, either known (such as control actions) or unknown (disturbances). Inputs of unknown magnitude are especially detrimental for observability, and they also complicate its analysis. Hence, the availability of computational tools capable of analysing the observability of nonlinear systems with unknown inputs has been limited until lately. Two symbolic algorithms based on differential geometry, ORC-DF and FISPO, have been recently proposed for this task, but their critical analysis and comparison is still lacking. Here we perform an analytical comparison of both algorithms and evaluate their performance on a set of problems, while discussing their strengths and limitations. Additionally, we use these analyses to provide insights about certain aspects of the relationship between inputs and observability. We found that, while ORC-DF and FISPO follow a similar approach, they differ in key aspects that can have a substantial influence on their applicability and computational cost. The FISPO algorithm is more generally applicable, since it can analyse any nonlinear ODE model. The ORC-DF algorithm analyses models that are affine in the inputs, and if those models have known inputs it is sometimes more efficient. Thus, the optimal choice of a method depends on the characteristics of the problem under consideration. To facilitate the use of both algorithms, we implemented the ORC-DF condition in a new version of STRIKE-GOLDD, a MATLAB toolbox for structural identifiability and observability analysis. Since this software tool already had an implementation of the FISPO algorithm, the new release allows modellers and model users the convenience of choosing between different algorithms in a single tool, without changing the coding of their model.

Parameterization of mechanistic models from qualitative data using an efficient optimal scaling approach

Quantitative dynamical models facilitate the understanding of biological processes and the prediction of their dynamics. These models usually comprise unknown parameters, which have to be inferred from experimental data. For quantitative experimental data, there are several methods and software tools available. However, for qualitative data the available approaches are limited and computationally demanding. Here, we consider the optimal scaling method which has been developed in statistics for categorical data and has been applied to dynamical systems. This approach turns qualitative variables into quantitative ones, accounting for constraints on their relation. We derive a reduced formulation for the optimization problem defining the optimal scaling. The reduced formulation possesses the same optimal points as the established formulation but requires less degrees of freedom. Parameter estimation for dynamical models of cellular pathways revealed that the reduced formulation improves the robustness and convergence of optimizers. This resulted in substantially reduced computation times. We implemented the proposed approach in the open-source Python Parameter EStimation TOolbox (pyPESTO) to facilitate reuse and extension. The proposed approach enables efficient parameterization of quantitative dynamical models using qualitative data.

Disentangling molecular mechanisms regulating sensitization of interferon alpha signal transduction

Tightly interlinked feedback regulators control the dynamics of intracellular responses elicited by the activation of signal transduction pathways. Interferon alpha (IFNα) orchestrates antiviral responses in hepatocytes, yet mechanisms that define pathway sensitization in response to prestimulation with different IFNα doses remained unresolved. We establish, based on quantitative measurements obtained for the hepatoma cell line Huh7.5, an ordinary differential equation model for IFNα signal transduction that comprises the feedback regulators STAT1, STAT2, IRF9, USP18, SOCS1, SOCS3, and IRF2. The model-based analysis shows that, mediated by the signaling proteins STAT2 and IRF9, prestimulation with a low IFNα dose hypersensitizes the pathway. In contrast, prestimulation with a high dose of IFNα leads to a dose-dependent desensitization, mediated by the negative regulators USP18 and SOCS1 that act at the receptor. The analysis of basal protein abundance in primary human hepatocytes reveals high heterogeneity in patient-specific amounts of STAT1, STAT2, IRF9, and USP18. The mathematical modeling approach shows that the basal amount of USP18 determines patient-specific pathway desensitization, while the abundance of STAT2 predicts the patient-specific IFNα signal response.

Benchmark problems for dynamic modeling of intracellular processes

Motivation

Dynamic models are used in systems biology to study and understand cellular processes like gene regulation or signal transduction. Frequently, ordinary differential equation (ODE) models are used to model the time and dose dependency of the abundances of molecular compounds as well as interactions and translocations. A multitude of computational approaches, e.g. for parameter estimation or uncertainty analysis have been developed within recent years. However, many of these approaches lack proper testing in application settings because a comprehensive set of benchmark problems is yet missing.

Results

We present a collection of 20 benchmark problems in order to evaluate new and existing methodologies, where an ODE model with corresponding experimental data is referred to as problem. In addition to the equations of the dynamical system, the benchmark collection provides observation functions as well as assumptions about measurement noise distributions and parameters. The presented benchmark models comprise problems of different size, complexity and numerical demands. Important characteristics of the models and methodological requirements are summarized, estimated parameters are provided, and some example studies were performed for illustrating the capabilities of the presented benchmark collection.

Availability and implementation

The models are provided in several standardized formats, including an easy-to-use human readable form and machine-readable SBML files. The data is provided as Excel sheets. All files are available at https://github.com/Benchmarking-Initiative/Benchmark-Models, including step-by-step explanations and MATLAB code to process and simulate the models.

Supplementary information

Supplementary data are available at Bioinformatics online.

Cytokine-mediated communication: a quantitative appraisal of immune complexity

Intercellular communication mediated by cytokines is the main mechanism by which cells of the immune system talk to each other. Many aspects of cytokine signalling in the immune system have been explored in great detail at the structural, biophysical, biochemical and cellular levels. However, a systematic understanding of the quantitative rules that govern cytokine-mediated cell-to-cell communication is still lacking. Here, we discuss recent efforts in the field of systems immunology to bring about a quantitative understanding of cytokine-mediated communication between leukocytes and to provide novel insights into the orchestration of immune responses and inflammation.

A Dual-Sensing DNA Nanostructure with an Ultrabroad Detection Range

Despite considerable interest in the development of biosensors that can measure analyte concentrations with a dynamic range spanning many orders of magnitude, this goal has proven difficult to achieve. We describe here a modular biosensor architecture that integrates two different readout mechanisms into a single-molecule construct that can achieve target detection across an extraordinarily broad dynamic range. Our dual-mode readout DNA biosensor combines an aptamer and a DNAzyme to quantify adenosine triphosphate (ATP) with two different mechanisms, which respond to low (micromolar) and high (millimolar) concentrations by generating distinct readouts based on changes in fluorescence and absorbance, respectively. Importantly, we have also devised regulatory strategies to fine-tune the target detection range of each sensor module by controlling the target-sensitivity of each readout mechanism. Using this strategy, we report the detection of ATP at a dynamic range spanning 1-500 000 μM, more than 5 orders of magnitude, representing the largest dynamic range reported to date with a single biosensor construct.

Hierarchical optimization for the efficient parametrization of ODE models

Motivation

Mathematical models are nowadays important tools for analyzing dynamics of cellular processes. The unknown model parameters are usually estimated from experimental data. These data often only provide information about the relative changes between conditions, hence, the observables contain scaling parameters. The unknown scaling parameters and corresponding noise parameters have to be inferred along with the dynamic parameters. The nuisance parameters often increase the dimensionality of the estimation problem substantially and cause convergence problems.

Results

In this manuscript, we propose a hierarchical optimization approach for estimating the parameters for ordinary differential equation (ODE) models from relative data. Our approach restructures the optimization problem into an inner and outer subproblem. These subproblems possess lower dimensions than the original optimization problem, and the inner problem can be solved analytically. We evaluated accuracy, robustness and computational efficiency of the hierarchical approach by studying three signaling pathways. The proposed approach achieved better convergence than the standard approach and required a lower computation time. As the hierarchical optimization approach is widely applicable, it provides a powerful alternative to established approaches.

Availability and implementation

The code is included in the MATLAB toolbox PESTO which is available at http://github.com/ICB-DCM/PESTO.

Supplementary information

Supplementary data are available at Bioinformatics online.

Extensions of ℓ1 regularization increase detection specificity for cell-type specific parameters in dynamic models

BACKGROUND

Ordinary differential equation systems are frequently utilized to model biological systems and to infer knowledge about underlying properties. For instance, the development of drugs requires the knowledge to which extent malign cells differ from healthy ones to provide a specific treatment with least side effects. As these cell-type specific properties may stem from any part of biochemical cell processes, systematic quantitative approaches are necessary to identify the relevant potential drug targets. An ℓ₁ regularization for the maximum likelihood parameter estimation proved to be successful, but falsely predicted cell-type dependent behaviour had to be corrected manually by using a Profile Likelihood approach.

RESULTS

The choice of extended ℓ₁ penalty functions significantly decreased the number of falsely detected cell-type specific parameters. Thus, the total accuracy of the prediction could be increased. This was tested on a realistic dynamical benchmark model used for the DREAM6 challenge. Among Elastic Net, Adaptive Lasso and a non-convex ℓ_q penalty, the latter one showed the best predictions whilst also requiring least computation time. All extended methods include a hyper-parameter in the regularization function. For an Erythropoietin (EPO) induced signalling pathway, the extended methods ℓ_q and Adaptive Lasso revealed an unpublished alternative parsimonious model when varying the respective hyper-parameters.

CONCLUSIONS

Using ℓ_q or Adaptive Lasso with an a-priori choice for the hyper-parameter can lead to a more specific and accurate result than ℓ₁. Scanning different hyper-parameters can yield additional pieces of information about the system.

Antiviral interferon response at single-cell resolution

Experimental studies of the innate immune response of mammalian cells to viruses reveal pervasive heterogeneity at the level of single cells. Interferons are induced only in a fraction of virus-infected cells; subsequently a fraction of cells exposed to interferons upregulate interferon-stimulated genes. Nevertheless, quantitative experiments and linked mathematical models show that the interferon response can be effective in curbing viral spread through two distinct mechanisms. First, paracrine interferon signals from scattered source cells can protect many uninfected cells, and the self-amplification of interferon production might serve to calibrate response amplitude to strength of viral infection. Second, models of the tug-of-war between viral replication and the innate interferon response imply a pivotal role of interferon action on already infected cells in curbing viral spread, through effectively lowering virus replication rate. This finding is in line with the observation that several pathogenic viruses selectively abrogate interferon action on infected cells. Thus, interferons may delay viral spread in acute infections by acting as sentinels, warning uninfected cells of imminent danger, or as negative feedback regulators of virus replication in infected cells. The timing of the interferon response relative to the onset of viral replication is critical for its effectiveness in curbing viral spread.

An easy and efficient approach for testing identifiability

Motivation

The feasibility of uniquely estimating parameters of dynamical systems from observations is a widely discussed aspect of mathematical modelling. Several approaches have been published for analyzing this so-called identifiability of model parameters. However, they are typically computationally demanding, difficult to perform and/or not applicable in many application settings.

Results

Here, an approach is presented which enables quickly testing of parameter identifiability. Numerical optimization with a penalty in radial direction enforcing displacement of the parameters is used to check whether estimated parameters are unique, or whether the parameters can be altered without loss of agreement with the data indicating non-identifiability. This Identifiability-Test by Radial Penalization (ITRP) can be employed for every model where optimization-based parameter estimation like least-squares or maximum likelihood is feasible and is therefore applicable for all typical systems biology models. The approach is illustrated and tested using 11 ordinary differential equation (ODE) models.

Availability and implementation

The presented approach can be implemented without great efforts in any modelling framework. It is available within the free Matlab-based modelling toolbox Data2Dynamics. Source code is available at https://github.com/Data2Dynamics.

Contact

ckreutz@fdm.uni-freiburg.de.

Supplementary information

Supplementary data are available at Bioinformatics online.

Leukemia inhibitory factor inhibits erythropoietin-induced myelin gene expression in oligodendrocytes

Background

The pro-myelinating effects of leukemia inhibitory factor (LIF) and other cytokines of the gp130 family, including oncostatin M (OSM) and ciliary neurotrophic factor (CNTF), have long been known, but controversial results have also been reported. We recently overexpressed erythropoietin receptor (EPOR) in rat central glia-4 (CG4) oligodendrocyte progenitor cells (OPCs) to study the mechanisms mediating the pro-myelinating effects of erythropoietin (EPO). In this study, we investigated the effect of co-treatment with EPO and LIF.

Methods

Gene expression in undifferentiated and differentiating CG4 cells in response to EPO and LIF was analysed by DNA microarrays and by RT-qPCR. Experiments were performed in biological replicates of N ≥ 4. Functional annotation and biological term enrichment was performed using DAVID (Database for Annotation, Visualization and Integrated Discovery). The gene-gene interaction network was visualised using STRING (Search Tool for the Retrieval of Interacting Genes).

Results

In CG4 cells treated with 10 ng/ml of EPO and 10 ng/ml of LIF, EPO-induced myelin oligodendrocyte glycoprotein (MOG) expression, measured at day 3 of differentiation, was inhibited ≥4-fold (N = 5, P < 0.001). Inhibition of EPO-induced MOG was also observed with OSM and CNTF. Analysis of the gene expression profile of CG4 differentiating cells treated for 20 h with EPO and LIF revealed LIF inhibition of EPO-induced genes involved in lipid transport and metabolism, previously identified as positive regulators of myelination in this system. In addition, among the genes induced by LIF, and not by differentiation or by EPO, the role of suppressor of cytokine signaling 3 (SOCS3) and toll like receptor 2 (TLR2) as negative regulators of myelination was further explored. LIF-induced SOCS3 was associated with MOG inhibition; Pam3, an agonist of TLR2, inhibited EPO-induced MOG expression, suggesting that TLR2 is functional and its activation decreases myelination.

Conclusions

Cytokines of the gp130 family may have negative effects on myelination, depending on the cytokine environment.

Electronic supplementary material

The online version of this article (10.1186/s10020-018-0052-3) contains supplementary material, which is available to authorized users.

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.

A scalable moment-closure approximation for large-scale biochemical reaction networks

Motivation

Stochastic molecular processes are a leading cause of cell-to-cell variability. Their dynamics are often described by continuous-time discrete-state Markov chains and simulated using stochastic simulation algorithms. As these stochastic simulations are computationally demanding, ordinary differential equation models for the dynamics of the statistical moments have been developed. The number of state variables of these approximating models, however, grows at least quadratically with the number of biochemical species. This limits their application to small- and medium-sized processes.

Results

In this article, we present a scalable moment-closure approximation (sMA) for the simulation of statistical moments of large-scale stochastic processes. The sMA exploits the structure of the biochemical reaction network to reduce the covariance matrix. We prove that sMA yields approximating models whose number of state variables depends predominantly on local properties, i.e. the average node degree of the reaction network, instead of the overall network size. The resulting complexity reduction is assessed by studying a range of medium- and large-scale biochemical reaction networks. To evaluate the approximation accuracy and the improvement in computational efficiency, we study models for JAK2/STAT5 signalling and NF κ B signalling. Our method is applicable to generic biochemical reaction networks and we provide an implementation, including an SBML interface, which renders the sMA easily accessible.

Availability and implementation

The sMA is implemented in the open-source MATLAB toolbox CERENA and is available from https://github.com/CERENADevelopers/CERENA .

Contact

jan.hasenauer@helmholtz-muenchen.de or atefeh.kazeroonian@tum.de.

Supplementary information

Supplementary data are available at Bioinformatics online.

The topological requirements for robust perfect adaptation in networks of any size

Robustness, and the ability to function and thrive amid changing and unfavorable environments, is a fundamental requirement for living systems. Until now it has been an open question how large and complex biological networks can exhibit robust behaviors, such as perfect adaptation to a variable stimulus, since complexity is generally associated with fragility. Here we report that all networks that exhibit robust perfect adaptation (RPA) to a persistent change in stimulus are decomposable into well-defined modules, of which there exist two distinct classes. These two modular classes represent a topological basis for all RPA-capable networks, and generate the full set of topological realizations of the internal model principle for RPA in complex, self-organizing, evolvable bionetworks. This unexpected result supports the notion that evolutionary processes are empowered by simple and scalable modular design principles that promote robust performance no matter how large or complex the underlying networks become.

Robust perfect adaptation (RPA), the ability of a system to return to its pre-stimulus state in the presence of a new signal, enables organisms to respond to further changes in stimuli. Here, the authors identify the modular structure of the full set of network topologies that can confer RPA on complex networks.

Profile likelihood-based analyses of infectious disease models

Ordinary differential equation models are frequently applied to describe the temporal evolution of epidemics. However, ordinary differential equation models are also utilized in other scientific fields. We summarize and transfer state-of-the art approaches from other fields like Systems Biology to infectious disease models. For this purpose, we use a simple SIR model with data from an influenza outbreak at an English boarding school in 1978 and a more complex model of a vector-borne disease with data from the Zika virus outbreak in Colombia in 2015-2016. Besides parameter estimation using a deterministic multistart optimization approach, a multitude of analyses based on the profile likelihood are presented comprising identifiability analysis and model reduction. The analyses were performed using the freely available modeling framework Data2Dynamics (data2dynamics.org) which has been awarded as best performing within the DREAM6 parameter estimation challenge and in the DREAM7 network reconstruction challenge.

The JWS online simulation database

Summary

JWS Online is a web-based platform for construction, simulation and exchange of models in standard formats. We have extended the platform with a database for curated simulation experiments that can be accessed directly via a URL, allowing one-click reproduction of published results. Users can modify the simulation experiments and export them in standard formats. The Simulation database thus lowers the bar on exploring computational models, helps users create valid simulation descriptions and improves the reproducibility of published simulation experiments.

Availability and Implementation

The Simulation Database is available on line at https://jjj.bio.vu.nl/models/experiments/ .

Contact

jls@sun.ac.za .

Model Based Targeting of IL-6-Induced Inflammatory Responses in Cultured Primary Hepatocytes to Improve Application of the JAK Inhibitor Ruxolitinib

IL-6 is a central mediator of the immediate induction of hepatic acute phase proteins (APP) in the liver during infection and after injury, but increased IL-6 activity has been associated with multiple pathological conditions. In hepatocytes, IL-6 activates JAK1-STAT3 signaling that induces the negative feedback regulator SOCS3 and expression of APPs. While different inhibitors of IL-6-induced JAK1-STAT3-signaling have been developed, understanding their precise impact on signaling dynamics requires a systems biology approach. Here we present a mathematical model of IL-6-induced JAK1-STAT3 signaling that quantitatively links physiological IL-6 concentrations to the dynamics of IL-6-induced signal transduction and expression of target genes in hepatocytes. The mathematical model consists of coupled ordinary differential equations (ODE) and the model parameters were estimated by a maximum likelihood approach, whereas identifiability of the dynamic model parameters was ensured by the Profile Likelihood. Using model simulations coupled with experimental validation we could optimize the long-term impact of the JAK-inhibitor Ruxolitinib, a therapeutic compound that is quickly metabolized. Model-predicted doses and timing of treatments helps to improve the reduction of inflammatory APP gene expression in primary mouse hepatocytes close to levels observed during regenerative conditions. The concept of improved efficacy of the inhibitor through multiple treatments at optimized time intervals was confirmed in primary human hepatocytes. Thus, combining quantitative data generation with mathematical modeling suggests that repetitive treatment with Ruxolitinib is required to effectively target excessive inflammatory responses without exceeding doses recommended by the clinical guidelines.

L1 regularization facilitates detection of cell type-specific parameters in dynamical systems

Motivation: A major goal of drug development is to selectively target certain cell types. Cellular decisions influenced by drugs are often dependent on the dynamic processing of information. Selective responses can be achieved by differences between the involved cell types at levels of receptor, signaling, gene regulation or further downstream. Therefore, a systematic approach to detect and quantify cell type-specific parameters in dynamical systems becomes necessary.

Results: Here, we demonstrate that a combination of nonlinear modeling with L₁ regularization is capable of detecting cell type-specific parameters. To adapt the least-squares numerical optimization routine to L₁ regularization, sub-gradient strategies as well as truncation of proposed optimization steps were implemented. Likelihood-ratio tests were used to determine the optimal regularization strength resulting in a sparse solution in terms of a minimal number of cell type-specific parameters that is in agreement with the data. By applying our implementation to a realistic dynamical benchmark model of the DREAM6 challenge we were able to recover parameter differences with an accuracy of 78%. Within the subset of detected differences, 91% were in agreement with their true value. Furthermore, we found that the results could be improved using the profile likelihood. In conclusion, the approach constitutes a general method to infer an overarching model with a minimum number of individual parameters for the particular models.

Availability and Implementation: A MATLAB implementation is provided within the freely available, open-source modeling environment Data2Dynamics. Source code for all examples is provided online at http://www.data2dynamics.org/.

Contact:bernhard.steiert@fdm.uni-freiburg.de

Direct targets of pSTAT5 signalling in erythropoiesis

Erythropoietin (EPO) acts through the dimeric erythropoietin receptor to stimulate proliferation, survival, differentiation and enucleation of erythroid progenitor cells. We undertook two complimentary approaches to find EPO-dependent pSTAT5 target genes in murine erythroid cells: RNA-seq of newly transcribed (4sU-labelled) RNA, and ChIP-seq for pSTAT5 30 minutes after EPO stimulation. We found 302 pSTAT5-occupied sites: ~15% of these reside in promoters while the rest reside within intronic enhancers or intergenic regions, some >100kb from the nearest TSS. The majority of pSTAT5 peaks contain a central palindromic GAS element, TTCYXRGAA. There was significant enrichment for GATA motifs and CACCC-box motifs within the neighbourhood of pSTAT5-bound peaks, and GATA1 and/or KLF1 co-occupancy at many sites. Using 4sU-RNA-seq we determined the EPO-induced transcriptome and validated differentially expressed genes using dynamic CAGE data and qRT-PCR. We identified known direct pSTAT5 target genes such as Bcl2l1, Pim1 and Cish, and many new targets likely to be involved in driving erythroid cell differentiation including those involved in mRNA splicing (Rbm25), epigenetic regulation (Suv420h2), and EpoR turnover (Clint1/EpsinR). Some of these new EpoR-JAK2-pSTAT5 target genes could be used as biomarkers for monitoring disease activity in polycythaemia vera, and for monitoring responses to JAK inhibitors.

Modeling cytokine regulatory network dynamics driving neuroinflammation in central nervous system disorders

A central goal of pharmacological efforts to treat central nervous system (CNS) diseases is to develop systemic therapeutics that can restore CNS homeostasis. Achieving this goal requires a fundamental understanding of CNS function within the organismal context so as to leverage the mechanistic insights on the molecular basis of cellular and tissue functions towards novel drug target identification. The immune system constitutes a key link between the periphery and CNS, and many neurological disorders and neurodegenerative diseases are characterized by immune dysfunction. We review the salient opportunities for applying computational models to CNS disease research, and summarize relevant approaches from studies of immune function and neuroinflammation. While the accurate prediction of disease-related phenomena is often considered the central goal of modeling studies, we highlight the utility of computational modeling applications beyond making predictions, particularly for drawing counterintuitive insights from model-based analysis of multi-parametric and time series data sets.

STAT5-Interacting Proteins: A Synopsis of Proteins that Regulate STAT5 Activity

Signal Transducers and Activators of Transcription (STATs) are key components of the JAK/STAT pathway. Of the seven STATs, STAT5A and STAT5B are of particular interest for their critical roles in cellular differentiation, adipogenesis, oncogenesis, and immune function. The interactions of STAT5A and STAT5B with cytokine/hormone receptors, nuclear receptors, transcriptional regulators, proto-oncogenes, kinases, and phosphatases all contribute to modulating STAT5 activity. Among these STAT5 interacting proteins, some serve as coactivators or corepressors to regulate STAT5 transcriptional activity and some proteins can interact with STAT5 to enhance or repress STAT5 signaling. In addition, a few STAT5 interacting proteins have been identified as positive regulators of STAT5 that alter serine and tyrosine phosphorylation of STAT5 while other proteins have been identified as negative regulators of STAT5 via dephosphorylation. This review article will discuss how STAT5 activity is modulated by proteins that physically interact with STAT5.

Protein abundance of AKT and ERK pathway components governs cell type-specific regulation of proliferation

Signaling through the AKT and ERK pathways controls cell proliferation. However, the integrated regulation of this multistep process, involving signal processing, cell growth and cell cycle progression, is poorly understood. Here, we study different hematopoietic cell types, in which AKT and ERK signaling is triggered by erythropoietin (Epo). Although these cell types share the molecular network topology for pro‐proliferative Epo signaling, they exhibit distinct proliferative responses. Iterating quantitative experiments and mathematical modeling, we identify two molecular sources for cell type‐specific proliferation. First, cell type‐specific protein abundance patterns cause differential signal flow along the AKT and ERK pathways. Second, downstream regulators of both pathways have differential effects on proliferation, suggesting that protein synthesis is rate‐limiting for faster cycling cells while slower cell cycles are controlled at the G1‐S progression. The integrated mathematical model of Epo‐driven proliferation explains cell type‐specific effects of targeted AKT and ERK inhibitors and faithfully predicts, based on the protein abundance, anti‐proliferative effects of inhibitors in primary human erythroid progenitor cells. Our findings suggest that the effectiveness of targeted cancer therapy might become predictable from protein abundance.

Defective Signaling in the JAK-STAT Pathway Tracks with Chronic Inflammation and Cardiovascular Risk in Aging Humans

Chronic inflammation, a decline in immune responsiveness, and reduced cardiovascular function are all associated with aging, but the relationships among these phenomena remain unclear. Here, we longitudinally profiled a total of 84 signaling conditions in 91 young and older adults and observed an age-related reduction in cytokine responsiveness within four immune cell lineages, most prominently T cells. The phenotype can be partially explained by elevated baseline levels of phosphorylated STAT (pSTAT) proteins and a different response capacity of naive versus memory T cell subsets to interleukin 6 (IL-6), interferon α (IFN-α), and, to a lesser extent, IL-21 and IFN-γ. Baseline pSTAT levels tracked with circulating levels of C-reactive protein (CRP), and we derived a cytokine response score that negatively correlates with measures of cardiovascular disease, specifically diastolic dysfunction and atherosclerotic burden, outperforming CRP. Thus, we identified an immunological link between inflammation, decreased cell responsiveness in the JAK-STAT pathway, and cardiovascular aging. Targeting chronic inflammation may ameliorate this deficiency in cellular responsiveness and improve cardiovascular function.

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.

Identification of Cell Type-Specific Differences in Erythropoietin Receptor Signaling in Primary Erythroid and Lung Cancer Cells

Lung cancer, with its most prevalent form non-small-cell lung carcinoma (NSCLC), is one of the leading causes of cancer-related deaths worldwide, and is commonly treated with chemotherapeutic drugs such as cisplatin. Lung cancer patients frequently suffer from chemotherapy-induced anemia, which can be treated with erythropoietin (EPO). However, studies have indicated that EPO not only promotes erythropoiesis in hematopoietic cells, but may also enhance survival of NSCLC cells. Here, we verified that the NSCLC cell line H838 expresses functional erythropoietin receptors (EPOR) and that treatment with EPO reduces cisplatin-induced apoptosis. To pinpoint differences in EPO-induced survival signaling in erythroid progenitor cells (CFU-E, colony forming unit-erythroid) and H838 cells, we combined mathematical modeling with a method for feature selection, the L₁ regularization. Utilizing an example model and simulated data, we demonstrated that this approach enables the accurate identification and quantification of cell type-specific parameters. We applied our strategy to quantitative time-resolved data of EPO-induced JAK/STAT signaling generated by quantitative immunoblotting, mass spectrometry and quantitative real-time PCR (qRT-PCR) in CFU-E and H838 cells as well as H838 cells overexpressing human EPOR (H838-HA-hEPOR). The established parsimonious mathematical model was able to simultaneously describe the data sets of CFU-E, H838 and H838-HA-hEPOR cells. Seven cell type-specific parameters were identified that included for example parameters for nuclear translocation of STAT5 and target gene induction. Cell type-specific differences in target gene induction were experimentally validated by qRT-PCR experiments. The systematic identification of pathway differences and sensitivities of EPOR signaling in CFU-E and H838 cells revealed potential targets for intervention to selectively inhibit EPO-induced signaling in the tumor cells but leave the responses in erythroid progenitor cells unaffected. Thus, the proposed modeling strategy can be employed as a general procedure to identify cell type-specific parameters and to recommend treatment strategies for the selective targeting of specific cell types.

A major challenge in the development of therapeutic interventions is the selective inhibition of a signal transduction pathway in one cell type such as a cancer cell leaving the other cell type such as a healthy cell as unaffected as possible. Here, we propose a new approach that combines mathematical modeling based on quantitative experimental data with statistical methods. We demonstrate based on simulated data that our approach can determine which parameters are the same and which parameters differ in two exemplary cell types. We compare a lung cancer cell line to the precursor cells of red blood cells. We show that the same signal transduction network induced by erythropoietin (EPO), a hormone that is frequently employed to treat anemia in cancer patients, regulates survival of both cell types. Based on our experimental data in combination with our computational approach, we identify seven cell type-specific differences in this signaling pathway. Our strategy allows predicting therapeutic targets that could be inhibited to interfere with survival of lung cancer cells while leaving production of red blood cells unaffected.

Inference for Stochastic Chemical Kinetics Using Moment Equations and System Size Expansion

Quantitative mechanistic models are valuable tools for disentangling biochemical pathways and for achieving a comprehensive understanding of biological systems. However, to be quantitative the parameters of these models have to be estimated from experimental data. In the presence of significant stochastic fluctuations this is a challenging task as stochastic simulations are usually too time-consuming and a macroscopic description using reaction rate equations (RREs) is no longer accurate. In this manuscript, we therefore consider moment-closure approximation (MA) and the system size expansion (SSE), which approximate the statistical moments of stochastic processes and tend to be more precise than macroscopic descriptions. We introduce gradient-based parameter optimization methods and uncertainty analysis methods for MA and SSE. Efficiency and reliability of the methods are assessed using simulation examples as well as by an application to data for Epo-induced JAK/STAT signaling. The application revealed that even if merely population-average data are available, MA and SSE improve parameter identifiability in comparison to RRE. Furthermore, the simulation examples revealed that the resulting estimates are more reliable for an intermediate volume regime. In this regime the estimation error is reduced and we propose methods to determine the regime boundaries. These results illustrate that inference using MA and SSE is feasible and possesses a high sensitivity.

In this manuscript, we introduce efficient methods for parameter estimation for stochastic processes. The stochasticity of chemical reactions can influence the average behavior of the considered system. For some biological systems, a microscopic, stochastic description is computationally intractable but a macroscopic, deterministic description too inaccurate. This inaccuracy manifests itself in an error in parameter estimates, which impede the predictive power of the proposed model. Until now, no rigorous analysis on the magnitude of the estimation error exists. We show by means of two simulation examples that using mesoscopic descriptions based on the system size expansions and moment-closure approximations can reduce this estimation error compared to inference using a macroscopic description. This reduction is most pronounced in an intermediate volume regime where the influence of stochasticity on the average behavior is moderately strong. For the JAK/STAT pathway where experimental data is available, we show that one parameter that was not structurally identifiable when using a macroscopic description becomes structurally identifiable when using a mesoscopic description for parameter estimation.