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.

Data-driven structural analysis of small cell lung cancer transcription factor network suggests potential subtype regulators and transition pathways.

Small cell lung cancer (SCLC) is an aggressive disease and challenging to treat due to its mixture of transcriptional subtypes and subtype transitions. Transcription factor (TF) networks have been the focus of studies to identify SCLC subtype regulators via systems approaches. Yet, their structures, which can provide clues on subtype drivers and transitions, are barely investigated. Here, we analyze the structure of an SCLC TF network by using graph theory concepts and identify its structurally important components responsible for complex signal processing, called hubs. We show that the hubs of the network are regulators of different SCLC subtypes by analyzing first the unbiased network structure and then integrating RNA-seq data as weights assigned to each interaction. Data-driven analysis emphasizes MYC as a hub, consistent with recent reports. Furthermore, we hypothesize that the pathways connecting functionally distinct hubs may control subtype transitions and test this hypothesis via network simulations on a candidate pathway and observe subtype transition. Overall, structural analyses of complex networks can identify their functionally important components and pathways driving the network dynamics. Such analyses can be an initial step for generating hypotheses and can guide the discovery of target pathways whose perturbation may change the network dynamics phenotypically.

Automatic mechanistic inference from large families of Boolean models generated by Monte Carlo tree search.

Many important processes in biology, such as signaling and gene regulation, can be described using logic models. These logic models are typically built to behaviorally emulate experimentally observed phenotypes, which are assumed to be steady states of a biological system. Most models are built by hand and therefore researchers are only able to consider one or perhaps a few potential mechanisms. We present a method to automatically synthesize Boolean logic models with a specified set of steady states. Our method, called MC-Boomer, is based on Monte Carlo Tree Search an efficient, parallel search method using reinforcement learning. Our approach enables users to constrain the model search space using prior knowledge or biochemical interaction databases, thus leading to generation of biologically plausible mechanistic hypotheses. Our approach can generate very large numbers of data-consistent models. To help develop mechanistic insight from these models, we developed analytical tools for multi-model inference and model selection. These tools reveal the key sets of interactions that govern the behavior of the models. We demonstrate that MC-Boomer works well at reconstructing randomly generated models. Then, using single time point measurements and reasonable biological constraints, our method generates hundreds of thousands of candidate models that match experimentally validated in-vivo behaviors of the Drosophila segment polarity network. Finally we outline how our multi-model analysis procedures elucidate potentially novel biological mechanisms and provide opportunities for model-driven experimental validation.

Integrated Multi-omics Analysis of Early Lung Adenocarcinoma Links Tumor Biological Features with Predicted Indolence or Aggressiveness.

Lung adenocarcinoma (LUAD) is a heterogeneous group of tumors associated with different survival rates, even when detected at an early stage. Here, we aim to investigate the biological determinants of early LUAD indolence or aggressiveness using radiomics as a surrogate of behavior. We present a set of 92 patients with LUAD with data collected across different methodologies. Patients were risk-stratified using the CT-based Score Indicative of Lung cancer Aggression (SILA) tool (0 = least aggressive, 1 = most aggressive). We grouped the patients as indolent (x ≤ 0.4, n = 14), intermediate (0.4 > x ≤ 0.6, n = 27), and aggressive (0.6 > x ≤ 1, n = 52). Using Cytometry by time of flight (CyTOF), we identified subpopulations with high HLA-DR expression that were associated with indolent behavior. In the RNA sequencing (RNA-seq) dataset, pathways related to immune response were associated with indolent behavior, while pathways associated with cell cycle and proliferation were associated with aggressive behavior. We extracted quantitative radiomics features from the CT scans of the patients. Integrating these datasets, we identified four feature signatures and four patient clusters that were associated with survival. Using single-cell RNA-seq, we found that indolent tumors had significantly more T cells and less B cells than aggressive tumors, and that the latter had a higher abundance of regulatory T cells and Th cells. In conclusion, we were able to uncover a correspondence between radiomics and tumor biology, which could improve the discrimination between indolent and aggressive LUAD tumors, enhance our knowledge in the biology of these tumors, and offer novel and personalized avenues for intervention. Significance: This study provides a comprehensive profiling of LUAD indolence and aggressiveness at the biological bulk and single-cell levels, as well as at the clinical and radiomics levels. This hypothesis generating study uncovers several potential future research avenues. It also highlights the importance and power of data integration to improve our systemic understanding of LUAD and to help reduce the gap between basic science research and clinical practice.

Unified tumor growth mechanisms from multimodel inference and dataset integration.

Mechanistic models of biological processes can explain observed phenomena and predict responses to a perturbation. A mathematical model is typically constructed using expert knowledge and informal reasoning to generate a mechanistic explanation for a given observation. Although this approach works well for simple systems with abundant data and well-established principles, quantitative biology is often faced with a dearth of both data and knowledge about a process, thus making it challenging to identify and validate all possible mechanistic hypothesis underlying a system behavior. To overcome these limitations, we introduce a Bayesian multimodel inference (Bayes-MMI) methodology, which quantifies how mechanistic hypotheses can explain a given experimental datasets, and concurrently, how each dataset informs a given model hypothesis, thus enabling hypothesis space exploration in the context of available data. We demonstrate this approach to probe standing questions about heterogeneity, lineage plasticity, and cell-cell interactions in tumor growth mechanisms of small cell lung cancer (SCLC). We integrate three datasets that each formulated different explanations for tumor growth mechanisms in SCLC, apply Bayes-MMI and find that the data supports model predictions for tumor evolution promoted by high lineage plasticity, rather than through expanding rare stem-like populations. In addition, the models predict that in the presence of cells associated with the SCLC-N or SCLC-A2 subtypes, the transition from the SCLC-A subtype to the SCLC-Y subtype through an intermediate is decelerated. Together, these predictions provide a testable hypothesis for observed juxtaposed results in SCLC growth and a mechanistic interpretation for tumor treatment resistance.

Data-driven structural analysis of Small Cell Lung Cancer transcription factor network suggests potential subtype regulators and transition pathways.

Small Cell Lung Cancer (SCLC) is an aggressive disease and challenging to treat due to its mixture of transcriptional subtypes and subtype transitions. Transcription factor (TF) networks have been the focus of studies to identify SCLC subtype regulators via systems approaches. Yet, their structures, which can provide clues on subtype drivers and transitions, are barely investigated. Here, we analyze the structure of an SCLC TF network by using graph theory concepts and identify its structurally important components responsible for complex signal processing, called hubs. We show that the hubs of the network are regulators of different SCLC subtypes by analyzing first the unbiased network structure and then integrating RNA-seq data as weights assigned to each interaction. Data-driven analysis emphasizes MYC as a hub, consistent with recent reports. Furthermore, we hypothesize that the pathways connecting functionally distinct hubs may control subtype transitions and test this hypothesis via network simulations on a candidate pathway and observe subtype transition. Overall, structural analyses of complex networks can identify their functionally important components and pathways driving the network dynamics. Such analyses can be an initial step for generating hypotheses and can guide the discovery of target pathways whose perturbation may change the network dynamics phenotypically.

Model certainty in cellular network-driven processes with missing data.

Mathematical models are often used to explore network-driven cellular processes from a systems perspective. However, a dearth of quantitative data suitable for model calibration leads to models with parameter unidentifiability and questionable predictive power. Here we introduce a combined Bayesian and Machine Learning Measurement Model approach to explore how quantitative and non-quantitative data constrain models of apoptosis execution within a missing data context. We find model prediction accuracy and certainty strongly depend on rigorous data-driven formulations of the measurement, and the size and make-up of the datasets. For instance, two orders of magnitude more ordinal (e.g., immunoblot) data are necessary to achieve accuracy comparable to quantitative (e.g., fluorescence) data for calibration of an apoptosis execution model. Notably, ordinal and nominal (e.g., cell fate observations) non-quantitative data synergize to reduce model uncertainty and improve accuracy. Finally, we demonstrate the potential of a data-driven Measurement Model approach to identify model features that could lead to informative experimental measurements and improve model predictive power.

A biochemical necroptosis model explains cell-type-specific responses to cell death cues.

Necroptosis is a form of regulated cell death associated with degenerative disorders, autoimmune and inflammatory diseases, and cancer. To better understand the biochemical mechanisms regulating necroptosis, we constructed a detailed computational model of tumor necrosis factor-induced necroptosis based on known molecular interactions from the literature. Intracellular protein levels, used as model inputs, were quantified using label-free mass spectrometry, and the model was calibrated using Bayesian parameter inference to experimental protein time course data from a well-established necroptosis-executing cell line. The calibrated model reproduced the dynamics of phosphorylated mixed lineage kinase domain-like protein, an established necroptosis reporter. A subsequent dynamical systems analysis identified four distinct modes of necroptosis signal execution, distinguished by rate constant values and the roles of the RIP1 deubiquitinating enzymes A20 and CYLD. In one case, A20 and CYLD both contribute to RIP1 deubiquitination, in another RIP1 deubiquitination is driven exclusively by CYLD, and in two modes either A20 or CYLD acts as the driver with the other enzyme, counterintuitively, inhibiting necroptosis. We also performed sensitivity analyses of initial protein concentrations and rate constants to identify potential targets for modulating necroptosis sensitivity within each mode. We conclude by associating numerous contrasting and, in some cases, counterintuitive experimental results reported in the literature with one or more of the model-predicted modes of necroptosis execution. In all, we demonstrate that a consensus pathway model of tumor necrosis factor-induced necroptosis can provide insights into unresolved controversies regarding the molecular mechanisms driving necroptosis execution in numerous cell types under different experimental conditions.

Processes in DNA damage response from a whole-cell multi-omics perspective.

Technological advances have made it feasible to collect multi-condition multi-omic time courses of cellular response to perturbation, but the complexity of these datasets impedes discovery due to challenges in data management, analysis, visualization, and interpretation. Here, we report a whole-cell mechanistic analysis of HL-60 cellular response to bendamustine. We integrate both enrichment and network analysis to show the progression of DNA damage and programmed cell death over time in molecular, pathway, and process-level detail using an interactive analysis framework for multi-omics data. Our framework, Mechanism of Action Generator Involving Network analysis (MAGINE), automates network construction and enrichment analysis across multiple samples and platforms, which can be integrated into our annotated gene-set network to combine the strengths of networks and ontology-driven analysis. Taken together, our work demonstrates how multi-omics integration can be used to explore signaling processes at various resolutions and demonstrates multi-pathway involvement beyond the canonical bendamustine mechanism.

Microbench: automated metadata management for systems biology benchmarking and reproducibility in Python.

MOTIVATION: Computational systems biology analyses typically make use of multiple software and their dependencies, which are often run across heterogeneous compute environments. This can introduce differences in performance and reproducibility. Capturing metadata (e.g. package versions, GPU model) currently requires repetitious code and is difficult to store centrally for analysis. Even where virtual environments and containers are used, updates over time mean that versioning metadata should still be captured within analysis pipelines to guarantee reproducibility. RESULTS: Microbench is a simple and extensible Python package to automate metadata capture to a file or Redis database. Captured metadata can include execution time, software package versions, environment variables, hardware information, Python version and more, with plugins. We present three case studies demonstrating Microbench usage to benchmark code execution and examine environment metadata for reproducibility purposes. AVAILABILITY AND IMPLEMENTATION: Install from the Python Package Index using pip install microbench. Source code is available from https://github.com/alubbock/microbench. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Programmatic modeling for biological systems.

Computational modeling has become an established technique to encode mathematical representations of cellular processes and gain mechanistic insights that drive testable predictions. These models are often constructed using graphical user interfaces or domain-specific languages, with community standards used for interchange. Models undergo steady state or dynamic analysis, which can include simulation and calibration within a single application, or transfer across various tools. Here, we describe a novel programmatic modeling paradigm, whereby modeling is augmented with software engineering best practices. We focus on Python - a popular programming language with a large scientific package ecosystem. Models can be encoded as programs, adding benefits such as modularity, testing, and automated documentation generators, while still being extensible and exportable to standardized formats for use with external tools if desired. Programmatic modeling is a key technology to enable collaborative model development and enhance dissemination, transparency, and reproducibility.

Selective Inhibition of JAK1 Primes STAT5-Driven Human Leukemia Cells for ATRA-Induced Differentiation.

BACKGROUND: All-trans retinoic acid (ATRA), a derivate of vitamin A, has been successfully used as a therapy to induce differentiation in M3 acute promyelocytic leukemia (APML), and has led to marked improvement in outcomes. Previously, attempts to use ATRA in non-APML in the clinic, however, have been underwhelming, likely due to persistent signaling through other oncogenic drivers. Dysregulated JAK/STAT signaling is known to drive several hematologic malignancies, and targeting JAK1 and JAK2 with the JAK1/JAK2 inhibitor ruxolitinib has led to improvement in survival in primary myelofibrosis and alleviation of vasomotor symptoms and splenomegaly in polycythemia vera and myelofibrosis. OBJECTIVE: While dose-dependent anemia and thrombocytopenia limit the use of JAK2 inhibition, selectively targeting JAK1 has been explored as a means to suppress inflammation and STAT-associated pathologies related to neoplastogenesis. The objective of this study is to employ JAK1 inhibition (JAK1i) in the presence of ATRA as a potential therapy in non-M3 acute myeloid leukemia (AML). METHODS: Efficacy of JAK1i using INCB52793 was assessed by changes in cell cycle and apoptosis in treated AML cell lines. Transcriptomic and proteomic analysis evaluated effects of JAK1i. Synergy between JAK1i+ ATRA was assessed in cell lines in vitro while efficacy in vivo was assessed by tumor reduction in MV-4-11 cell line-derived xenografts. RESULTS: Here we describe novel synergistic activity between JAK1i inhibition and ATRA in non-M3 leukemia. Transcriptomic and proteomic analysis confirmed structural and functional changes related to maturation while in vivo combinatory studies revealed significant decreases in leukemic expansion. CONCLUSIONS: JAK1i+ ATRA lead to decreases in cell cycle followed by myeloid differentiation and cell death in human leukemias. These findings highlight potential uses of ATRA-based differentiation therapy of non-M3 human leukemia.

MuSyC is a consensus framework that unifies multi-drug synergy metrics for combinatorial drug discovery.

Drug combination discovery depends on reliable synergy metrics but no consensus exists on the correct synergy criterion to characterize combined interactions. The fragmented state of the field confounds analysis, impedes reproducibility, and delays clinical translation of potential combination treatments. Here we present a mass-action based formalism to quantify synergy. With this formalism, we clarify the relationship between the dominant drug synergy principles, and present a mapping of commonly used frameworks onto a unified synergy landscape. From this, we show how biases emerge due to intrinsic assumptions which hinder their broad applicability and impact the interpretation of synergy in discovery efforts. Specifically, we describe how traditional metrics mask consequential synergistic interactions, and contain biases dependent on the Hill-slope and maximal effect of single-drugs. We show how these biases systematically impact synergy classification in large combination screens, potentially misleading discovery efforts. Thus the proposed formalism can provide a consistent, unbiased interpretation of drug synergy, and accelerate the translatability of synergy studies.

Unsupervised logic-based mechanism inference for network-driven biological processes.

Modern analytical techniques enable researchers to collect data about cellular states, before and after perturbations. These states can be characterized using analytical techniques, but the inference of regulatory interactions that explain and predict changes in these states remains a challenge. Here we present a generalizable, unsupervised approach to generate parameter-free, logic-based models of cellular processes, described by multiple discrete states. Our algorithm employs a Hamming-distance based approach to formulate, test, and identify optimized logic rules that link two states. Our approach comprises two steps. First, a model with no prior knowledge except for the mapping between initial and attractor states is built. We then employ biological constraints to improve model fidelity. Our algorithm automatically recovers the relevant dynamics for the explored models and recapitulates key aspects of the biochemical species concentration dynamics in the original model. We present the advantages and limitations of our work and discuss how our approach could be used to infer logic-based mechanisms of signaling, gene-regulatory, or other input-output processes describable by the Boolean formalism.

Thunor: visualization and analysis of high-throughput dose-response datasets.

High-throughput cell proliferation assays to quantify drug-response are becoming increasingly common and powerful with the emergence of improved automation and multi-time point analysis methods. However, pipelines for analysis of these datasets that provide reproducible, efficient, and interactive visualization and interpretation are sorely lacking. To address this need, we introduce Thunor, an open-source software platform to manage, analyze, and visualize large, dose-dependent cell proliferation datasets. Thunor supports both end-point and time-based proliferation assays as input. It provides a simple, user-friendly interface with interactive plots and publication-quality images of cell proliferation time courses, dose-response curves, and derived dose-response metrics, e.g. IC50, including across datasets or grouped by tags. Tags are categorical labels for cell lines and drugs, used for aggregation, visualization and statistical analysis, e.g. cell line mutation or drug class/target pathway. A graphical plate map tool is included to facilitate plate annotation with cell lines, drugs and concentrations upon data upload. Datasets can be shared with other users via point-and-click access control. We demonstrate the utility of Thunor to examine and gain insight from two large drug response datasets: a large, publicly available cell viability database and an in-house, high-throughput proliferation rate dataset. Thunor is available from www.thunor.net.

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

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

Charting the Fragmented Landscape of Drug Synergy.

Even as the clinical impact of drug combinations continues to accelerate, no consensus on how to quantify drug synergy has emerged. Rather, surveying the landscape of drug synergy reveals the persistence of historical fissures regarding the appropriate domains of conflicting synergy models - fissures impacting all aspects of combination therapy discovery and deployment. Herein we chronicle the impact of these divisions on: (i) the design, interpretation, and reproducibility of high-throughput combination screens; (ii) the performance of algorithms to predict synergistic mixtures; and (iii) the search for higher-order synergistic interactions. Further progress in each of these subfields hinges on reaching a consensus regarding the long-standing rifts in the field.

ACDC: Automated Cell Detection and Counting for Time-Lapse Fluorescence Microscopy.

Advances in microscopy imaging technologies have enabled the visualization of live-cell dynamic processes using time-lapse microscopy imaging. However, modern methods exhibit several limitations related to the training phases and to time constraints, hindering their application in the laboratory practice. In this work, we present a novel method, named Automated Cell Detection and Counting (ACDC), designed for activity detection of fluorescent labeled cell nuclei in time-lapse microscopy. ACDC overcomes the limitations of the literature methods, by first applying bilateral filtering on the original image to smooth the input cell images while preserving edge sharpness, and then by exploiting the watershed transform and morphological filtering. Moreover, ACDC represents a feasible solution for the laboratory practice, as it can leverage multi-core architectures in computer clusters to efficiently handle large-scale imaging datasets. Indeed, our Parent-Workers implementation of ACDC allows to obtain up to a 3.7× speed-up compared to the sequential counterpart. ACDC was tested on two distinct cell imaging datasets to assess its accuracy and effectiveness on images with different characteristics. We achieved an accurate cell-count and nuclei segmentation without relying on large-scale annotated datasets, a result confirmed by the average Dice Similarity Coefficients of 76.84 and 88.64 and the Pearson coefficients of 0.99 and 0.96, calculated against the manual cell counting, on the two tested datasets.

Interactive Multiresolution Visualization of Cellular Network Processes.

Visualization plays a central role in the analysis of biochemical network models to identify patterns that arise from reaction dynamics and perform model exploratory analysis. To facilitate these analyses, we developed PyViPR, a visualization tool that generates static and dynamic representations of biochemical network processes within a Python-based environment. PyViPR embeds network visualizations within Jupyter notebooks, thus enabling integration with modeling, simulation, and analysis workflows. To present the capabilities of PyViPR, we explore execution mechanisms of extrinsic apoptosis in HeLa cells. We show that community-detection algorithms identify groups of molecular species that capture key biological functions and ease exploration of the apoptosis network. We then show how different kinetic parameter sets that fit the experimental data equally well exhibit significantly different signal-execution dynamics as the system progresses toward mitochondrial outer-membrane permeabilization. Therefore, PyViPR aids the conceptual understanding of dynamic network processes and accelerates hypothesis generation for further testing and validation.

Systems-level network modeling of Small Cell Lung Cancer subtypes identifies master regulators and destabilizers.

Adopting a systems approach, we devise a general workflow to define actionable subtypes in human cancers. Applied to small cell lung cancer (SCLC), the workflow identifies four subtypes based on global gene expression patterns and ontologies. Three correspond to known subtypes (SCLC-A, SCLC-N, and SCLC-Y), while the fourth is a previously undescribed ASCL1+ neuroendocrine variant (NEv2, or SCLC-A2). Tumor deconvolution with subtype gene signatures shows that all of the subtypes are detectable in varying proportions in human and mouse tumors. To understand how multiple stable subtypes can arise within a tumor, we infer a network of transcription factors and develop BooleaBayes, a minimally-constrained Boolean rule-fitting approach. In silico perturbations of the network identify master regulators and destabilizers of its attractors. Specific to NEv2, BooleaBayes predicts ELF3 and NR0B1 as master regulators of the subtype, and TCF3 as a master destabilizer. Since the four subtypes exhibit differential drug sensitivity, with NEv2 consistently least sensitive, these findings may lead to actionable therapeutic strategies that consider SCLC intratumoral heterogeneity. Our systems-level approach should generalize to other cancer types.