Antibody-mediated targeting of TNFR2 activates CD8+ T cells in mice and promotes antitumor immunity.

Tumor necrosis factor receptor 2 (TNFR2) is the alternate receptor for TNF and can mediate both pro- and anti-inflammatory activities of T cells. Although TNFR2 has been linked to enhanced suppressive activity of regulatory T cells (Tregs) in autoimmune diseases, the viability of TNFR2 as a target for cancer immunotherapy has been underappreciated. Here, we show that new murine monoclonal anti-TNFR2 antibodies yield robust antitumor activity and durable protective memory in multiple mouse cancer cell line models. The antibodies mediate potent Fc-dependent T cell costimulation and do not result in significant depletion of Tregs Corresponding human agonistic monoclonal anti-TNFR2 antibodies were identified and also had antitumor effects in humanized mouse models. Anti-TNFR2 antibodies could be developed as a novel treatment option for patients with cancer.

Recipes for Analysis of Molecular Networks Using the Data2Dynamics Modeling Environment.

Mechanistic models of biomolecular processes are established research tools that enable to quantitatively investigate dynamic features of biological processes such as signal transduction cascades. Often, these models aim at describing a large number of states, for instance concentrations of proteins and small molecules, as well as their interactions. Each modeled interaction increases the number of potentially unknown parameters like reaction rate constants or initial amount of proteins. In order to calibrate these mechanistic models, the unknown model parameters have to be estimated based on experimental data. The complexity of parameter estimation raises several computational challenges that can be tackled within the Data2Dynamics modeling environment. The environment is a well-tested, high-performance software package that is tailored to the modeling of biological processes with ordinary differential equation models and using experimental biomolecular data.In this chapter, we introduce and provide "recipes" for the most frequent analyses and modeling tasks in the Data2Dynamics modeling environment. The presented protocols comprise model building, data handling, parameter estimation, calculation of confidence intervals, model selection and reduction, deriving prediction uncertainties, and designing informative novel experiments.

MM-131, a bispecific anti-Met/EpCAM mAb, inhibits HGF-dependent and HGF-independent Met signaling through concurrent binding to EpCAM.

Activation of the Met receptor tyrosine kinase, either by its ligand, hepatocyte growth factor (HGF), or via ligand-independent mechanisms, such as MET amplification or receptor overexpression, has been implicated in driving tumor proliferation, metastasis, and resistance to therapy. Clinical development of Met-targeted antibodies has been challenging, however, as bivalent antibodies exhibit agonistic properties, whereas monovalent antibodies lack potency and the capacity to down-regulate Met. Through computational modeling, we found that the potency of a monovalent antibody targeting Met could be dramatically improved by introducing a second binding site that recognizes an unrelated, highly expressed antigen on the tumor cell surface. Guided by this prediction, we engineered MM-131, a bispecific antibody that is monovalent for both Met and epithelial cell adhesion molecule (EpCAM). MM-131 is a purely antagonistic antibody that blocks ligand-dependent and ligand-independent Met signaling by inhibiting HGF binding to Met and inducing receptor down-regulation. Together, these mechanisms lead to inhibition of proliferation in Met-driven cancer cells, inhibition of HGF-mediated cancer cell migration, and inhibition of tumor growth in HGF-dependent and -independent mouse xenograft models. Consistent with its design, MM-131 is more potent in EpCAM-high cells than in EpCAM-low cells, and its potency decreases when EpCAM levels are reduced by RNAi. Evaluation of Met, EpCAM, and HGF levels in human tumor samples reveals that EpCAM is expressed at high levels in a wide range of Met-positive tumor types, suggesting a broad opportunity for clinical development of MM-131.

Analysis of Single-Cell RNA-Seq Identifies Cell-Cell Communication Associated with Tumor Characteristics.

Tumor ecosystems are composed of multiple cell types that communicate by ligand-receptor interactions. Targeting ligand-receptor interactions (for instance, with immune checkpoint inhibitors) can provide significant benefits for patients. However, our knowledge of which interactions occur in a tumor and how these interactions affect outcome is still limited. We present an approach to characterize communication by ligand-receptor interactions across all cell types in a microenvironment using single-cell RNA sequencing. We apply this approach to identify and compare the ligand-receptor interactions present in six syngeneic mouse tumor models. To identify interactions potentially associated with outcome, we regress interactions against phenotypic measurements of tumor growth rate. In addition, we quantify ligand-receptor interactions between T cell subsets and their relation to immune infiltration using a publicly available human melanoma dataset. Overall, this approach provides a tool for studying cell-cell interactions, their variability across tumors, and their relationship to outcome.

Modeling chemotherapy-induced stress to identify rational combination therapies in the DNA damage response pathway.

Cells respond to DNA damage by activating complex signaling networks that decide cell fate, promoting not only DNA damage repair and survival but also cell death. We have developed a multiscale computational model that quantitatively links chemotherapy-induced DNA damage response signaling to cell fate. The computational model was trained and calibrated on extensive data from U2OS osteosarcoma cells, including the cell cycle distribution of the initial cell population, signaling data measured by Western blotting, and cell fate data in response to chemotherapy treatment measured by time-lapse microscopy. The resulting mechanistic model predicted the cellular responses to chemotherapy alone and in combination with targeted inhibitors of the DNA damage response pathway, which we confirmed experimentally. Computational models such as the one presented here can be used to understand the molecular basis that defines the complex interplay between cell survival and cell death and to rationally identify chemotherapy-potentiating drug combinations.

Estimation of immune cell content in tumour tissue using single-cell RNA-seq data.

As interactions between the immune system and tumour cells are governed by a complex network of cell-cell interactions, knowing the specific immune cell composition of a solid tumour may be essential to predict a patient's response to immunotherapy. Here, we analyse in depth how to derive the cellular composition of a solid tumour from bulk gene expression data by mathematical deconvolution, using indication-specific and cell type-specific reference gene expression profiles (RGEPs) from tumour-derived single-cell RNA sequencing data. We demonstrate that tumour-derived RGEPs are essential for the successful deconvolution and that RGEPs from peripheral blood are insufficient. We distinguish nine major cell types, as well as three T cell subtypes. Using the tumour-derived RGEPs, we can estimate the content of many tumours associated immune and stromal cell types, their therapeutically relevant ratios, as well as an improved gene expression profile of the malignant cells.

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.

Predicting ligand-dependent tumors from multi-dimensional signaling features.

Targeted therapies have shown significant patient benefit in about 5-10% of solid tumors that are addicted to a single oncogene. Here, we explore the idea of ligand addiction as a driver of tumor growth. High ligand levels in tumors have been shown to be associated with impaired patient survival, but targeted therapies have not yet shown great benefit in unselected patient populations. Using an approach of applying Bagged Decision Trees (BDT) to high-dimensional signaling features derived from a computational model, we can predict ligand dependent proliferation across a set of 58 cell lines. This mechanistic, multi-pathway model that features receptor heterodimerization, was trained on seven cancer cell lines and can predict signaling across two independent cell lines by adjusting only the receptor expression levels for each cell line. Interestingly, for patient samples the predicted tumor growth response correlates with high growth factor expression in the tumor microenvironment, which argues for a co-evolution of both factors in vivo.

IL-1ß-induced and p38MAPK-dependent activation of the mitogen-activated protein kinase-activated protein kinase 2 (MK2) in hepatocytes: Signal transduction with robust and concentration-independent signal amplification.

The IL-1ß induced activation of the p38MAPK/MAPK-activated protein kinase 2 (MK2) pathway in hepatocytes is important for control of the acute phase response and regulation of liver regeneration. Many aspects of the regulatory relevance of this pathway have been investigated in immune cells in the context of inflammation. However, very little is known about concentration-dependent activation kinetics and signal propagation in hepatocytes and the role of MK2. We established a mathematical model for IL-1ß-induced activation of the p38MAPK/MK2 pathway in hepatocytes that was calibrated to quantitative data on time- and IL-1ß concentration-dependent phosphorylation of p38MAPK and MK2 in primary mouse hepatocytes. This analysis showed that, in hepatocytes, signal transduction from IL-1ß via p38MAPK to MK2 is characterized by strong signal amplification. Quantification of p38MAPK and MK2 revealed that, in hepatocytes, at maximum, 11.3% of p38MAPK molecules and 36.5% of MK2 molecules are activated in response to IL-1ß. The mathematical model was experimentally validated by employing phosphatase inhibitors and the p38MAPK inhibitor SB203580. Model simulations predicted an IC50 of 1-1.2 µm for SB203580 in hepatocytes. In silico analyses and experimental validation demonstrated that the kinase activity of p38MAPK determines signal amplitude, whereas phosphatase activity affects both signal amplitude and duration. p38MAPK and MK2 concentrations and responsiveness toward IL-1ß were quantitatively compared between hepatocytes and macrophages. In macrophages, the absolute p38MAPK and MK2 concentration was significantly higher. Finally, in line with experimental observations, the mathematical model predicted a significantly higher half-maximal effective concentration for IL-1ß-induced pathway activation in macrophages compared with hepatocytes, underscoring the importance of cell type-specific differences in pathway regulation.

Driving the Model to Its Limit: Profile Likelihood Based Model Reduction.

In systems biology, one of the major tasks is to tailor model complexity to information content of the data. A useful model should describe the data and produce well-determined parameter estimates and predictions. Too small of a model will not be able to describe the data whereas a model which is too large tends to overfit measurement errors and does not provide precise predictions. Typically, the model is modified and tuned to fit the data, which often results in an oversized model. To restore the balance between model complexity and available measurements, either new data has to be gathered or the model has to be reduced. In this manuscript, we present a data-based method for reducing non-linear models. The profile likelihood is utilised to assess parameter identifiability and designate likely candidates for reduction. Parameter dependencies are analysed along profiles, providing context-dependent suggestions for the type of reduction. We discriminate four distinct scenarios, each associated with a specific model reduction strategy. Iterating the presented procedure eventually results in an identifiable model, which is capable of generating precise and testable predictions. Source code for all toy examples is provided within the freely available, open-source modelling environment Data2Dynamics based on MATLAB available at http://www.data2dynamics.org/, as well as the R packages dMod/cOde available at https://github.com/dkaschek/. Moreover, the concept is generally applicable and can readily be used with any software capable of calculating the profile likelihood.

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

Model calibration and uncertainty analysis in signaling networks.

For a long time the biggest challenges in modeling cellular signal transduction networks has been the inference of crucial pathway components and the qualitative description of their interactions. As a result of the emergence of powerful high-throughput experiments, it is now possible to measure data of high temporal and spatial resolution and to analyze signaling dynamics quantitatively. In addition, this increase of high-quality data is the basis for a better understanding of model limitations and their influence on the predictive power of models. We review established approaches in signal transduction network modeling with a focus on ordinary differential equation models as well as related developments in model calibration. As central aspects of the calibration process we discuss possibilities of model adaptation based on data-driven parameter optimization and the concomitant objective of reducing model uncertainties.

CERENA: ChEmical REaction Network Analyzer--A Toolbox for the Simulation and Analysis of Stochastic Chemical Kinetics.

Gene expression, signal transduction and many other cellular processes are subject to stochastic fluctuations. The analysis of these stochastic chemical kinetics is important for understanding cell-to-cell variability and its functional implications, but it is also challenging. A multitude of exact and approximate descriptions of stochastic chemical kinetics have been developed, however, tools to automatically generate the descriptions and compare their accuracy and computational efficiency are missing. In this manuscript we introduced CERENA, a toolbox for the analysis of stochastic chemical kinetics using Approximations of the Chemical Master Equation solution statistics. CERENA implements stochastic simulation algorithms and the finite state projection for microscopic descriptions of processes, the system size expansion and moment equations for meso- and macroscopic descriptions, as well as the novel conditional moment equations for a hybrid description. This unique collection of descriptions in a single toolbox facilitates the selection of appropriate modeling approaches. Unlike other software packages, the implementation of CERENA is completely general and allows, e.g., for time-dependent propensities and non-mass action kinetics. By providing SBML import, symbolic model generation and simulation using MEX-files, CERENA is user-friendly and computationally efficient. The availability of forward and adjoint sensitivity analyses allows for further studies such as parameter estimation and uncertainty analysis. The MATLAB code implementing CERENA is freely available from http://cerenadevelopers.github.io/CERENA/.

Summary of the DREAM8 Parameter Estimation Challenge: Toward Parameter Identification for Whole-Cell Models.

Whole-cell models that explicitly represent all cellular components at the molecular level have the potential to predict phenotype from genotype. However, even for simple bacteria, whole-cell models will contain thousands of parameters, many of which are poorly characterized or unknown. New algorithms are needed to estimate these parameters and enable researchers to build increasingly comprehensive models. We organized the Dialogue for Reverse Engineering Assessments and Methods (DREAM) 8 Whole-Cell Parameter Estimation Challenge to develop new parameter estimation algorithms for whole-cell models. We asked participants to identify a subset of parameters of a whole-cell model given the model's structure and in silico "experimental" data. Here we describe the challenge, the best performing methods, and new insights into the identifiability of whole-cell models. We also describe several valuable lessons we learned toward improving future challenges. Going forward, we believe that collaborative efforts supported by inexpensive cloud computing have the potential to solve whole-cell model parameter estimation.

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.

Deciphering the mechanism behind Fibroblast Growth Factor (FGF) induced biphasic signal-response profiles.

BACKGROUND: The Fibroblast Growth Factor (FGF) pathway is driving various aspects of cellular responses in both normal and malignant cells. One interesting characteristic of this pathway is the biphasic nature of the cellular response to some FGF ligands like FGF2. Specifically, it has been shown that phenotypic behaviors controlled by FGF signaling, like migration and growth, reach maximal levels in response to intermediate concentrations, while high levels of FGF2 elicit weak responses. The mechanisms leading to the observed biphasic response remains unexplained. RESULTS: A combination of experiments and computational modeling was used to understand the mechanism behind the observed biphasic signaling responses. FGF signaling involves a tertiary surface interaction that we captured with a computational model based on Ordinary Differential Equations (ODEs). It accounts for FGF2 binding to FGF receptors (FGFRs) and heparan sulfate glycosaminoglycans (HSGAGs), followed by receptor-phosphorylation, activation of the FRS2 adapter protein and the Ras-Raf signaling cascade. Quantitative protein assays were used to measure the dynamics of phosphorylated ERK (pERK) in response to a wide range of FGF2 ligand concentrations on a fine-grained time scale for the squamous cell lung cancer cell line H1703. We developed a novel approach combining Particle Swarm Optimization (PSO) and feature-based constraints in the objective function to calibrate the computational model to the experimental data. The model is validated using a series of extracellular and intracellular perturbation experiments. We demonstrate that in silico model predictions are in accordance with the observed in vitro results. CONCLUSIONS: Using a combined approach of computational modeling and experiments we found that competition between binding of the ligand FGF2 to HSGAG and FGF receptor leads to the biphasic response. At low to intermediate concentrations of FGF2 there are sufficient free FGF receptors available for the FGF2-HSGAG complex to enable the formation of the trimeric signaling unit. At high ligand concentrations the ligand binding sites of the receptor become saturated and the trimeric signaling unit cannot be formed. This insight into the pathway is an important consideration for the pharmacological inhibition of this pathway.

Bridging the gaps in systems biology.

Systems biology aims at creating mathematical models, i.e., computational reconstructions of biological systems and processes that will result in a new level of understanding-the elucidation of the basic and presumably conserved "design" and "engineering" principles of biomolecular systems. Thus, systems biology will move biology from a phenomenological to a predictive science. Mathematical modeling of biological networks and processes has already greatly improved our understanding of many cellular processes. However, given the massive amount of qualitative and quantitative data currently produced and number of burning questions in health care and biotechnology needed to be solved is still in its early phases. The field requires novel approaches for abstraction, for modeling bioprocesses that follow different biochemical and biophysical rules, and for combining different modules into larger models that still allow realistic simulation with the computational power available today. We have identified and discussed currently most prominent problems in systems biology: (1) how to bridge different scales of modeling abstraction, (2) how to bridge the gap between topological and mechanistic modeling, and (3) how to bridge the wet and dry laboratory gap. The future success of systems biology largely depends on bridging the recognized gaps.

Network topology and parameter estimation: from experimental design methods to gene regulatory network kinetics using a community based approach.

BACKGROUND: Accurate estimation of parameters of biochemical models is required to characterize the dynamics of molecular processes. This problem is intimately linked to identifying the most informative experiments for accomplishing such tasks. While significant progress has been made, effective experimental strategies for parameter identification and for distinguishing among alternative network topologies remain unclear. We approached these questions in an unbiased manner using a unique community-based approach in the context of the DREAM initiative (Dialogue for Reverse Engineering Assessment of Methods). We created an in silico test framework under which participants could probe a network with hidden parameters by requesting a range of experimental assays; results of these experiments were simulated according to a model of network dynamics only partially revealed to participants. RESULTS: We proposed two challenges; in the first, participants were given the topology and underlying biochemical structure of a 9-gene regulatory network and were asked to determine its parameter values. In the second challenge, participants were given an incomplete topology with 11 genes and asked to find three missing links in the model. In both challenges, a budget was provided to buy experimental data generated in silico with the model and mimicking the features of different common experimental techniques, such as microarrays and fluorescence microscopy. Data could be bought at any stage, allowing participants to implement an iterative loop of experiments and computation. CONCLUSIONS: A total of 19 teams participated in this competition. The results suggest that the combination of state-of-the-art parameter estimation and a varied set of experimental methods using a few datasets, mostly fluorescence imaging data, can accurately determine parameters of biochemical models of gene regulation. However, the task is considerably more difficult if the gene network topology is not completely defined, as in challenge 2. Importantly, we found that aggregating independent parameter predictions and network topology across submissions creates a solution that can be better than the one from the best-performing submission.

Comparison of approaches for parameter identifiability analysis of biological systems.

MOTIVATION: Modeling of dynamical systems using ordinary differential equations is a popular approach in the field of Systems Biology. The amount of experimental data that are used to build and calibrate these models is often limited. In this setting, the model parameters may not be uniquely determinable. Structural or a priori identifiability is a property of the system equations that indicates whether, in principle, the unknown model parameters can be determined from the available data. RESULTS: We performed a case study using three current approaches for structural identifiability analysis for an application from cell biology. The approaches are conceptually different and are developed independently. The results of the three approaches are in agreement. We discuss strength and weaknesses of each of them and illustrate how they can be applied to real world problems. AVAILABILITY AND IMPLEMENTATION: For application of the approaches to further applications, code representations (DAISY, Mathematica and MATLAB) for benchmark model and data are provided on the authors webpage. CONTACT: andreas.raue@fdm.uni-freiburg.de.

Lessons learned from quantitative dynamical modeling in systems biology.

Due to the high complexity of biological data it is difficult to disentangle cellular processes relying only on intuitive interpretation of measurements. A Systems Biology approach that combines quantitative experimental data with dynamic mathematical modeling promises to yield deeper insights into these processes. Nevertheless, with growing complexity and increasing amount of quantitative experimental data, building realistic and reliable mathematical models can become a challenging task: the quality of experimental data has to be assessed objectively, unknown model parameters need to be estimated from the experimental data, and numerical calculations need to be precise and efficient. Here, we discuss, compare and characterize the performance of computational methods throughout the process of quantitative dynamic modeling using two previously established examples, for which quantitative, dose- and time-resolved experimental data are available. In particular, we present an approach that allows to determine the quality of experimental data in an efficient, objective and automated manner. Using this approach data generated by different measurement techniques and even in single replicates can be reliably used for mathematical modeling. For the estimation of unknown model parameters, the performance of different optimization algorithms was compared systematically. Our results show that deterministic derivative-based optimization employing the sensitivity equations in combination with a multi-start strategy based on latin hypercube sampling outperforms the other methods by orders of magnitude in accuracy and speed. Finally, we investigated transformations that yield a more efficient parameterization of the model and therefore lead to a further enhancement in optimization performance. We provide a freely available open source software package that implements the algorithms and examples compared here.