The use of machine learning to discover regulatory networks controlling biological systems.

Biological systems are composed of a vast web of multiscale molecular interactors and interactions. High-throughput technologies, both bulk and single cell, now allow for investigation of the properties and quantities of these interactors. Computational algorithms and machine learning methods then provide the tools to derive meaningful insights from the resulting data sets. One such approach is graphical network modeling, which provides a computational framework to explicitly model the molecular interactions within and between the cells comprising biological systems. These graphical networks aim to describe a putative chain of cause and effect between interacting molecules. This feature allows for determination of key molecules in a biological process, accelerated generation of mechanistic hypotheses, and simulation of experimental outcomes. We review the computational concepts and applications of graphical network models across molecular scales for both intracellular and intercellular regulatory biology, examples of successful applications, and the future directions needed to overcome current limitations.

Leveraging multi-omics data to empower quantitative systems pharmacology in immuno-oncology.

Understanding the intricate interactions of cancer cells with the tumor microenvironment (TME) is a pre-requisite for the optimization of immunotherapy. Mechanistic models such as quantitative systems pharmacology (QSP) provide insights into the TME dynamics and predict the efficacy of immunotherapy in virtual patient populations/digital twins but require vast amounts of multimodal data for parameterization. Large-scale datasets characterizing the TME are available due to recent advances in bioinformatics for multi-omics data. Here, we discuss the perspectives of leveraging omics-derived bioinformatics estimates to inform QSP models and circumvent the challenges of model calibration and validation in immuno-oncology.

Multiomics Empowers Predictive Pancreatic Cancer Immunotherapy.

Advances in cancer immunotherapy, particularly immune checkpoint inhibitors, have dramatically improved the prognosis for patients with metastatic melanoma and other previously incurable cancers. However, patients with pancreatic ductal adenocarcinoma (PDAC) generally do not respond to these therapies. PDAC is exceptionally difficult to treat because of its often late stage at diagnosis, modest mutation burden, and notoriously complex and immunosuppressive tumor microenvironment. Simultaneously interrogating features of cancer, immune, and other cellular components of the PDAC tumor microenvironment is therefore crucial for identifying biomarkers of immunotherapeutic resistance and response. Notably, single-cell and multiomics technologies, along with the analytical tools for interpreting corresponding data, are facilitating discoveries of the systems-level cellular and molecular interactions contributing to the overall resistance of PDAC to immunotherapy. Thus, in this review, we will explore how multiomics and single-cell analyses provide the unprecedented opportunity to identify biomarkers of resistance and response to successfully sensitize PDAC to immunotherapy.

A field guide to cultivating computational biology.

Evolving in sync with the computation revolution over the past 30 years, computational biology has emerged as a mature scientific field. While the field has made major contributions toward improving scientific knowledge and human health, individual computational biology practitioners at various institutions often languish in career development. As optimistic biologists passionate about the future of our field, we propose solutions for both eager and reluctant individual scientists, institutions, publishers, funding agencies, and educators to fully embrace computational biology. We believe that in order to pave the way for the next generation of discoveries, we need to improve recognition for computational biologists and better align pathways of career success with pathways of scientific progress. With 10 outlined steps, we call on all adjacent fields to move away from the traditional individual, single-discipline investigator research model and embrace multidisciplinary, data-driven, team science.

Engineered bispecific antibodies targeting the interleukin-6 and -8 receptors potently inhibit cancer cell migration and tumor metastasis.

Simultaneous inhibition of interleukin-6 (IL-6) and interleukin-8 (IL-8) signaling diminishes cancer cell migration, and combination therapy has recently been shown to synergistically reduce metastatic burden in a preclinical model of triple-negative breast cancer. Here, we have engineered two novel bispecific antibodies that target the IL-6 and IL-8 receptors to concurrently block the signaling activity of both ligands. We demonstrate that a first-in-class bispecific antibody design has promising therapeutic potential, with enhanced selectivity and potency compared with monoclonal antibody and small-molecule drug combinations in both cellular and animal models of metastatic triple-negative breast cancer. Mechanistic characterization revealed that our engineered bispecific antibodies have no impact on cell viability, but profoundly reduce the migratory potential of cancer cells; hence they constitute a true anti-metastatic treatment. Moreover, we demonstrate that our antibodies can be readily combined with standard-of-care anti-proliferative drugs to develop effective anti-cancer regimens. Collectively, our work establishes an innovative metastasis-focused direction for cancer drug development.

Matrix factorization and transfer learning uncover regulatory biology across multiple single-cell ATAC-seq data sets.

While the methods available for single-cell ATAC-seq analysis are well optimized for clustering cell types, the question of how to integrate multiple scATAC-seq data sets and/or sequencing modalities is still open. We present an analysis framework that enables such integration across scATAC-seq data sets by applying the CoGAPS Matrix Factorization algorithm and the projectR transfer learning program to identify common regulatory patterns across scATAC-seq data sets. We additionally integrate our analysis with scRNA-seq data to identify orthogonal evidence for transcriptional regulators predicted by scATAC-seq analysis. Using publicly available scATAC-seq data, we find patterns that accurately characterize cell types both within and across data sets. Furthermore, we demonstrate that these patterns are both consistent with current biological understanding and reflective of novel regulatory biology.

Enter the Matrix: Factorization Uncovers Knowledge from Omics.

Omics data contain signals from the molecular, physical, and kinetic inter- and intracellular interactions that control biological systems. Matrix factorization (MF) techniques can reveal low-dimensional structure from high-dimensional data that reflect these interactions. These techniques can uncover new biological knowledge from diverse high-throughput omics data in applications ranging from pathway discovery to timecourse analysis. We review exemplary applications of MF for systems-level analyses. We discuss appropriate applications of these methods, their limitations, and focus on the analysis of results to facilitate optimal biological interpretation. The inference of biologically relevant features with MF enables discovery from high-throughput data beyond the limits of current biological knowledge - answering questions from high-dimensional data that we have not yet thought to ask.

projectR: an R/Bioconductor package for transfer learning via PCA, NMF, correlation and clustering.

MOTIVATION: Dimension reduction techniques are widely used to interpret high-dimensional biological data. Features learned from these methods are used to discover both technical artifacts and novel biological phenomena. Such feature discovery is critically importent in analysis of large single-cell datasets, where lack of a ground truth limits validation and interpretation. Transfer learning (TL) can be used to relate the features learned from one source dataset to a new target dataset to perform biologically driven validation by evaluating their use in or association with additional sample annotations in that independent target dataset. RESULTS: We developed an R/Bioconductor package, projectR, to perform TL for analyses of genomics data via TL of clustering, correlation and factorization methods. We then demonstrate the utility TL for integrated data analysis with an example for spatial single-cell analysis. AVAILABILITY AND IMPLEMENTATION: projectR is available on Bioconductor and at https://github.com/genesofeve/projectR. CONTACT: gsteinobrien@jhmi.edu or ejfertig@jhmi.edu. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Digitizing omics profiles by divergence from a baseline.

Data collected from omics technologies have revealed pervasive heterogeneity and stochasticity of molecular states within and between phenotypes. A prominent example of such heterogeneity occurs between genome-wide mRNA, microRNA, and methylation profiles from one individual tumor to another, even within a cancer subtype. However, current methods in bioinformatics, such as detecting differentially expressed genes or CpG sites, are population-based and therefore do not effectively model intersample diversity. Here we introduce a unified theory to quantify sample-level heterogeneity that is applicable to a single omics profile. Specifically, we simplify an omics profile to a digital representation based on the omics profiles from a set of samples from a reference or baseline population (e.g., normal tissues). The state of any subprofile (e.g., expression vector for a subset of genes) is said to be "divergent" if it lies outside the estimated support of the baseline distribution and is consequently interpreted as "dysregulated" relative to that baseline. We focus on two cases: single features (e.g., individual genes) and distinguished subsets (e.g., regulatory pathways). Notably, since the divergence analysis is at the individual sample level, dysregulation can be analyzed probabilistically; for example, one can estimate the probability that a gene or pathway is divergent in some population. Finally, the reduction in complexity facilitates a more "personalized" and biologically interpretable analysis of variation, as illustrated by experiments involving tissue characterization, disease detection and progression, and disease-pathway associations.

CoGAPS 3: Bayesian non-negative matrix factorization for single-cell analysis with asynchronous updates and sparse data structures.

BACKGROUND: Bayesian factorization methods, including Coordinated Gene Activity in Pattern Sets (CoGAPS), are emerging as powerful analysis tools for single cell data. However, these methods have greater computational costs than their gradient-based counterparts. These costs are often prohibitive for analysis of large single-cell datasets. Many such methods can be run in parallel which enables this limitation to be overcome by running on more powerful hardware. However, the constraints imposed by the prior distributions in CoGAPS limit the applicability of parallelization methods to enhance computational efficiency for single-cell analysis. RESULTS: We developed a new software framework for parallel matrix factorization in Version 3 of the CoGAPS R/Bioconductor package to overcome the computational limitations of Bayesian matrix factorization for single cell data analysis. This parallelization framework provides asynchronous updates for sequential updating steps of the algorithm to enhance computational efficiency. These algorithmic advances were coupled with new software architecture and sparse data structures to reduce the memory overhead for single-cell data. CONCLUSIONS: Altogether our new software enhance the efficiency of the CoGAPS Bayesian matrix factorization algorithm so that it can analyze 1000 times more cells, enabling factorization of large single-cell data sets.

Integrated single-cell and bulk gene expression and ATAC-seq reveals heterogeneity and early changes in pathways associated with resistance to cetuximab in HNSCC-sensitive cell lines.

BACKGROUND: Identifying potential resistance mechanisms while tumour cells still respond to therapy is critical to delay acquired resistance. METHODS: We generated the first comprehensive multi-omics, bulk and single-cell data in sensitive head and neck squamous cell carcinoma (HNSCC) cells to identify immediate responses to cetuximab. Two pathways potentially associated with resistance were focus of the study: regulation of receptor tyrosine kinases by TFAP2A transcription factor, and epithelial-to-mesenchymal transition (EMT). RESULTS: Single-cell RNA-seq demonstrates heterogeneity, with cell-specific TFAP2A and VIM expression profiles in response to treatment and also with global changes to various signalling pathways. RNA-seq and ATAC-seq reveal global changes within 5 days of therapy, suggesting early onset of mechanisms of resistance; and corroborates cell line heterogeneity, with different TFAP2A targets or EMT markers affected by therapy. Lack of TFAP2A expression is associated with HNSCC decreased growth, with cetuximab and JQ1 increasing the inhibitory effect. Regarding the EMT process, short-term cetuximab therapy has the strongest effect on inhibiting migration. TFAP2A silencing does not affect cell migration, supporting an independent role for both mechanisms in resistance. CONCLUSION: Overall, we show that immediate adaptive transcriptional and epigenetic changes induced by cetuximab are heterogeneous and cell type dependent; and independent mechanisms of resistance arise while tumour cells are still sensitive to therapy.

Correction: Integrated single-cell and bulk gene expression and ATAC-seq reveals heterogeneity and early changes in pathways associated with resistance to cetuximab in HNSCC-sensitive cell lines.

An amendment to this paper has been published and can be accessed via a link at the top of the paper.

Inferring causal molecular networks: empirical assessment through a community-based effort.

It remains unclear whether causal, rather than merely correlational, relationships in molecular networks can be inferred in complex biological settings. Here we describe the HPN-DREAM network inference challenge, which focused on learning causal influences in signaling networks. We used phosphoprotein data from cancer cell lines as well as in silico data from a nonlinear dynamical model. Using the phosphoprotein data, we scored more than 2,000 networks submitted by challenge participants. The networks spanned 32 biological contexts and were scored in terms of causal validity with respect to unseen interventional data. A number of approaches were effective, and incorporating known biology was generally advantageous. Additional sub-challenges considered time-course prediction and visualization. Our results suggest that learning causal relationships may be feasible in complex settings such as disease states. Furthermore, our scoring approach provides a practical way to empirically assess inferred molecular networks in a causal sense.

Splice Expression Variation Analysis (SEVA) for inter-tumor heterogeneity of gene isoform usage in cancer.

Motivation: Current bioinformatics methods to detect changes in gene isoform usage in distinct phenotypes compare the relative expected isoform usage in phenotypes. These statistics model differences in isoform usage in normal tissues, which have stable regulation of gene splicing. Pathological conditions, such as cancer, can have broken regulation of splicing that increases the heterogeneity of the expression of splice variants. Inferring events with such differential heterogeneity in gene isoform usage requires new statistical approaches. Results: We introduce Splice Expression Variability Analysis (SEVA) to model increased heterogeneity of splice variant usage between conditions (e.g. tumor and normal samples). SEVA uses a rank-based multivariate statistic that compares the variability of junction expression profiles within one condition to the variability within another. Simulated data show that SEVA is unique in modeling heterogeneity of gene isoform usage, and benchmark SEVA's performance against EBSeq, DiffSplice and rMATS that model differential isoform usage instead of heterogeneity. We confirm the accuracy of SEVA in identifying known splice variants in head and neck cancer and perform cross-study validation of novel splice variants. A novel comparison of splice variant heterogeneity between subtypes of head and neck cancer demonstrated unanticipated similarity between the heterogeneity of gene isoform usage in HPV-positive and HPV-negative subtypes and anticipated increased heterogeneity among HPV-negative samples with mutations in genes that regulate the splice variant machinery. These results show that SEVA accurately models differential heterogeneity of gene isoform usage from RNA-seq data. Availability and implementation: SEVA is implemented in the R/Bioconductor package GSReg. Contact: bahman@jhu.edu or favorov@sensi.org or ejfertig@jhmi.edu. Supplementary information: Supplementary data are available at Bioinformatics online.

CancerInSilico: An R/Bioconductor package for combining mathematical and statistical modeling to simulate time course bulk and single cell gene expression data in cancer.

Bioinformatics techniques to analyze time course bulk and single cell omics data are advancing. The absence of a known ground truth of the dynamics of molecular changes challenges benchmarking their performance on real data. Realistic simulated time-course datasets are essential to assess the performance of time course bioinformatics algorithms. We develop an R/Bioconductor package, CancerInSilico, to simulate bulk and single cell transcriptional data from a known ground truth obtained from mathematical models of cellular systems. This package contains a general R infrastructure for running cell-based models and simulating gene expression data based on the model states. We show how to use this package to simulate a gene expression data set and consequently benchmark analysis methods on this data set with a known ground truth. The package is freely available via Bioconductor: http://bioconductor.org/packages/CancerInSilico/.

StereoGene: rapid estimation of genome-wide correlation of continuous or interval feature data.

MOTIVATION: Genomics features with similar genome-wide distributions are generally hypothesized to be functionally related, for example, colocalization of histones and transcription start sites indicate chromatin regulation of transcription factor activity. Therefore, statistical algorithms to perform spatial, genome-wide correlation among genomic features are required. RESULTS: Here, we propose a method, StereoGene, that rapidly estimates genome-wide correlation among pairs of genomic features. These features may represent high-throughput data mapped to reference genome or sets of genomic annotations in that reference genome. StereoGene enables correlation of continuous data directly, avoiding the data binarization and subsequent data loss. Correlations are computed among neighboring genomic positions using kernel correlation. Representing the correlation as a function of the genome position, StereoGene outputs the local correlation track as part of the analysis. StereoGene also accounts for confounders such as input DNA by partial correlation. We apply our method to numerous comparisons of ChIP-Seq datasets from the Human Epigenome Atlas and FANTOM CAGE to demonstrate its wide applicability. We observe the changes in the correlation between epigenomic features across developmental trajectories of several tissue types consistent with known biology and find a novel spatial correlation of CAGE clusters with donor splice sites and with poly(A) sites. These analyses provide examples for the broad applicability of StereoGene for regulatory genomics. AVAILABILITY AND IMPLEMENTATION: The StereoGene C ++ source code, program documentation, Galaxy integration scripts and examples are available from the project homepage http://stereogene.bioinf.fbb.msu.ru/. CONTACT: favorov@sensi.org. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

PatternMarkers & GWCoGAPS for novel data-driven biomarkers via whole transcriptome NMF.

SUMMARY: Non-negative Matrix Factorization (NMF) algorithms associate gene expression with biological processes (e.g. time-course dynamics or disease subtypes). Compared with univariate associations, the relative weights of NMF solutions can obscure biomarkers. Therefore, we developed a novel patternMarkers statistic to extract genes for biological validation and enhanced visualization of NMF results. Finding novel and unbiased gene markers with patternMarkers requires whole-genome data. Therefore, we also developed Genome-Wide CoGAPS Analysis in Parallel Sets (GWCoGAPS), the first robust whole genome Bayesian NMF using the sparse, MCMC algorithm, CoGAPS. Additionally, a manual version of the GWCoGAPS algorithm contains analytic and visualization tools including patternMatcher, a Shiny web application. The decomposition in the manual pipeline can be replaced with any NMF algorithm, for further generalization of the software. Using these tools, we find granular brain-region and cell-type specific signatures with corresponding biomarkers in GTEx data, illustrating GWCoGAPS and patternMarkers ascertainment of data-driven biomarkers from whole-genome data. AVAILABILITY AND IMPLEMENTATION: PatternMarkers & GWCoGAPS are in the CoGAPS Bioconductor package (3.5) under the GPL license. CONTACT: gsteinobrien@jhmi.edu or ccolantu@jhmi.edu or ejfertig@jhmi.edu. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Comparative mutational landscape analysis of patient-derived tumour xenografts.

BACKGROUND: Screening of patients for cancer-driving mutations is now used for cancer prognosis, remission scoring and treatment selection. Although recently emerged targeted next-generation sequencing-based approaches offer promising diagnostic capabilities, there are still limitations. There is a pressing clinical need for a well-validated, rapid, cost-effective mutation profiling system in patient specimens. Given their speed and cost-effectiveness, quantitative PCR mutation detection techniques are well suited for the clinical environment. The qBiomarker mutation PCR array has high sensitivity and shorter turnaround times compared with other methods. However, a direct comparison with existing viable alternatives are required to assess its true potential and limitations. METHODS: In this study, we evaluated a panel of 117 patient-derived tumour xenografts by the qBiomarker array and compared with other methods for mutation detection, including Ion AmpliSeq sequencing, whole-exome sequencing and droplet digital PCR. RESULTS: Our broad analysis demonstrates that the qBiomarker's performance is on par with that of other labour-intensive and expensive methods of cancer mutation detection of frequently altered cancer-associated genes, and provides a foundation for supporting its consideration as an option for molecular diagnostics. CONCLUSIONS: This large-scale direct comparison and validation of currently available mutation detection approaches is extremely relevant for the current scenario of precision medicine and will lead to informed choice of screening methodologies, especially in lower budget conditions or time frame limitations.

An unbiased in vivo functional genomics screening approach in mice identifies novel tumor cell-based regulators of immune rejection.

The clinical successes of immune checkpoint therapies for cancer make it important to identify mechanisms of resistance to anti-tumor immune responses. Numerous resistance mechanisms have been identified employing studies of single genes or pathways, thereby parsing the tumor microenvironment complexity into tractable pieces. However, this limits the potential for novel gene discovery to in vivo immune attack. To address this challenge, we developed an unbiased in vivo genome-wide RNAi screening platform that leverages host immune selection in strains of immune-competent and immunodeficient mice to select for tumor cell-based genes that regulate in vivo sensitivity to immune attack. Utilizing this approach in a syngeneic triple-negative breast cancer (TNBC) model, we identified 709 genes that selectively regulated adaptive anti-tumor immunity and focused on five genes (CD47, TGFß1, Sgpl1, Tex9 and Pex14) with the greatest impact. We validated the mechanisms that underlie the immune-related effects of expression of these genes in different TNBC lines, as well as tandem synergistic interactions. Furthermore, we demonstrate the impact of different genes with previously unknown immune functions (Tex9 and Pex14) on anti-tumor immunity. Thus, this innovative approach has utility in identifying unknown tumor-specific regulators of immune recognition in multiple settings to reveal novel targets for future immunotherapies.