CoGAPS: an R/C++ package to identify patterns and biological process activity in transcriptomic data

Identification of shared gene expression programs activated in multiple modes of torpor across vertebrate clades

Torpor encompasses diverse adaptations to extreme environmental stressors such as hibernation, aestivation, brumation, and daily torpor. Here we introduce StrokeofGenus, an analytic pipeline that identifies distinct transcriptomic states and shared gene expression patterns across studies, tissues, and species. We use StrokeofGenus to study multiple and diverse forms of torpor from publicly-available RNA-seq datasets that span eight species and two classes. We identify three transcriptionally distinct states during the cycle of heterothermia: euthermia, torpor, and interbout arousal. We also identify torpor-specific gene expression patterns that are shared both across tissues and between species with over three hundred million years of evolutionary divergence. We further demonstrate the general sharing of gene expression patterns in multiple forms of torpor, implying a common evolutionary origin for this process. Although here we apply StrokeofGenus to analysis of torpor, it can be used to interrogate any other complex physiological processes defined by transient transcriptomic states.

Spatial multi-omics reveal intratumoral humoral immunity niches associated with tertiary lymphoid structures in pancreatic cancer immunotherapy pathologic responders

Pancreatic adenocarcinoma (PDAC) is a rapidly progressing cancer that responds poorly to immunotherapies. Intratumoral tertiary lymphoid structures (TLS) have been associated with rare long-term PDAC survivors, but the role of TLS in PDAC and their spatial relationships within the context of the broader tumor microenvironment remain unknown. We generated a spatial multi-omics atlas encompassing 26 PDAC tumors from patients treated with combination immunotherapies. Using machine learning-enabled H&E image classification models and unsupervised gene expression matrix factorization methods for spatial transcriptomics, we characterized cellular states within TLS niches spanning across distinct morphologies and immunotherapies. Unsupervised learning generated a TLS-specific spatial gene expression signature that significantly associates with improved survival in PDAC patients. These analyses demonstrate TLS-associated intratumoral B cell maturation in pathological responders, confirmed with spatial proteomics and BCR profiling. Our study also identifies spatial features of pathologic immune responses, revealing TLS maturation colocalizing with IgG/IgA distribution and extracellular matrix remodeling.

GRAPHICAL ABSTRACT

HIGHLIGHTS

Integrated multi-modal spatial profiling of human PDAC tumors from neoadjuvant immunotherapy clinical trials reveal diverse spatial niches enriched in TLS.TLS maturity is influenced by tumor location and the cellular neighborhoods in which TLS immune cells are recruited.Unsupervised machine learning of genome-wide signatures on spatial transcriptomics data characterizes the TLS-enriched TME and associates TLS transcriptomes with survival outcomes in PDAC.Interactions of spatially variable gene expression patterns showed TLS maturation is coupled with immunoglobulin distribution and ECM remodeling in pathologic responders.Intratumoral plasma cell and immunoglobin gene expression spatial dynamics demonstrate trafficking of TLS-driven humoral immunity in the PDAC TME.

Significance

We report a spatial multi-omics atlas of PDAC tumors from a series of immunotherapy neoadjuvant clinical trials. Intratumorally, pathologic responders exhibit mature TLS that propagate plasma cells into malignant niches. Our findings offer insights on the role of TLS-associated humoral immunity and stromal remodeling during immunotherapy treatment.

PanIN and CAF transitions in pancreatic carcinogenesis revealed with spatial data integration

This study introduces a new imaging, spatial transcriptomics (ST), and single-cell RNA-sequencing integration pipeline to characterize neoplastic cell state transitions during tumorigenesis. We applied a semi-supervised analysis pipeline to examine premalignant pancreatic intraepithelial neoplasias (PanINs) that can develop into pancreatic ductal adenocarcinoma (PDAC). Their strict diagnosis on formalin-fixed and paraffin-embedded (FFPE) samples limited the single-cell characterization of human PanINs within their microenvironment. We leverage whole transcriptome FFPE ST to enable the study of a rare cohort of matched low-grade (LG) and high-grade (HG) PanIN lesions to track progression and map cellular phenotypes relative to single-cell PDAC datasets. We demonstrate that cancer-associated fibroblasts (CAFs), including antigen-presenting CAFs, are located close to PanINs. We further observed a transition from CAF-related inflammatory signaling to cellular proliferation during PanIN progression. We validate these findings with single-cell high-dimensional imaging proteomics and transcriptomics technologies. Altogether, our semi-supervised learning framework for spatial multi-omics has broad applicability across cancer types to decipher the spatiotemporal dynamics of carcinogenesis.

Interface-guided phenotyping of coding variants in the transcription factor RUNX1

Single-gene missense mutations remain challenging to interpret. Here, we deploy scalable functional screening by sequencing (SEUSS), a Perturb-seq method, to generate mutations at protein interfaces of RUNX1 and quantify their effect on activities of downstream cellular programs. We evaluate single-cell RNA profiles of 115 mutations in myelogenous leukemia cells and categorize them into three functionally distinct groups, wild-type (WT)-like, loss-of-function (LoF)-like, and hypomorphic, that we validate in orthogonal assays. LoF-like variants dominate the DNA-binding site and are recurrent in cancer; however, recurrence alone does not predict functional impact. Hypomorphic variants share characteristics with LoF-like but favor protein interactions, promoting gene expression indicative of nerve growth factor (NGF) response and cytokine recruitment of neutrophils. Accessible DNA near differentially expressed genes frequently contains RUNX1-binding motifs. Finally, we reclassify 16 variants of uncertain significance and train a classifier to predict 103 more. Our work demonstrates the potential of targeting protein interactions to better define the landscape of phenotypes reachable by missense mutations.

Multi-sample non-negative spatial factorization

It is important to model biological variation when analyzing spatial transcriptomics data from multiple samples. One approach to multi-sample analysis is to spatially align samples, but this is a challenging problem. Here, we provide an alignment-free framework for generalizing a one-sample spatial factorization model to multi-sample data. Using this framework, we develop a method, called multi-sample non-negative spatial factorization (mNSF) that extends the one-sample non-negative spatial factorization (NSF) framework to a multi-sample dataset. Our model allows for a sample-specific model for the spatial correlation structure and extracts a low-dimensional representation of the data. We illustrate the performance of mNSF by simulation studies and real data. mNSF identifies true factors in simulated data, identifies shared anatomical regions across samples in real data and reveals region-specific biological functions. mNSFs performance is similar to alignment based methods when alignment is possible, but extends analysis to situations where spatial alignment is impossible. We expect multi-sample factorization methods to be a powerful class of methods for analyzing spatially resolved transcriptomics data.

Cytokines drive the formation of memory-like NK cell subsets via epigenetic rewiring and transcriptional regulation

Activation of natural killer (NK) cells with the cytokines interleukin-12 (IL-12), IL-15, and IL-18 induces their differentiation into memory-like (ML) NK cells; however, the underlying epigenetic and transcriptional mechanisms are unclear. By combining ATAC-seq, CITE-seq, and functional analyses, we discovered that IL-12/15/18 activation results in two main human NK fates: reprogramming into enriched memory-like (eML) NK cells or priming into effector conventional NK (effcNK) cells. eML NK cells had distinct transcriptional and epigenetic profiles and enhanced function, whereas effcNK cells resembled cytokine-primed cNK cells. Two transcriptionally discrete subsets of eML NK cells were also identified, eML-1 and eML-2, primarily arising from CD56bright or CD56dim mature NK cell subsets, respectively. Furthermore, these eML subsets were evident weeks after transfer of IL-12/15/18-activated NK cells into patients with cancer. Our findings demonstrate that NK cell activation with IL-12/15/18 results in previously unappreciated diverse cellular fates and identifies new strategies to enhance NK therapies.

Gene expression signatures in blood from a West African sepsis cohort define host response phenotypes

Our limited understanding of the pathophysiological mechanisms that operate during sepsis is an obstacle to rational treatment and clinical trial design. There is a critical lack of data from low- and middle-income countries where the sepsis burden is increased which inhibits generalized strategies for therapeutic intervention. Here we perform RNA sequencing of whole blood to investigate longitudinal host response to sepsis in a Ghanaian cohort. Data dimensional reduction reveals dynamic gene expression patterns that describe cell type-specific molecular phenotypes including a dysregulated myeloid compartment shared between sepsis and COVID-19. The gene expression signatures reported here define a landscape of host response to sepsis that supports interventions via targeting immunophenotypes to improve outcomes.

Transfer Learning Reveals Cancer-Associated Fibroblasts Are Associated with Epithelial-Mesenchymal Transition and Inflammation in Cancer Cells in Pancreatic Ductal Adenocarcinoma

Pancreatic ductal adenocarcinoma (PDAC) is an aggressive malignancy characterized by an immunosuppressive tumor microenvironment enriched with cancer-associated fibroblasts (CAF). This study used a convergence approach to identify tumor cell and CAF interactions through the integration of single-cell data from human tumors with human organoid coculture experiments. Analysis of a comprehensive atlas of PDAC single-cell RNA sequencing data indicated that CAF density is associated with increased inflammation and epithelial-mesenchymal transition (EMT) in epithelial cells. Transfer learning using transcriptional data from patient-derived organoid and CAF cocultures provided in silico validation of CAF induction of inflammatory and EMT epithelial cell states. Further experimental validation in cocultures demonstrated integrin beta 1 (ITGB1) and vascular endothelial factor A (VEGFA) interactions with neuropilin-1 mediating CAF-epithelial cell cross-talk. Together, this study introduces transfer learning from human single-cell data to organoid coculture analyses for experimental validation of discoveries of cell-cell cross-talk and identifies fibroblast-mediated regulation of EMT and inflammation.

SIGNIFICANCE

Adaptation of transfer learning to relate human single-cell RNA sequencing data to organoid-CAF cocultures facilitates discovery of human pancreatic cancer intercellular interactions and uncovers cross-talk between CAFs and tumor cells through VEGFA and ITGB1.

An integrated single-nucleus and spatial transcriptomics atlas reveals the molecular landscape of the human hippocampus

The hippocampus contains many unique cell types, which serve the structure's specialized functions, including learning, memory and cognition. These cells have distinct spatial topography, morphology, physiology, and connectivity, highlighting the need for transcriptome-wide profiling strategies that retain cytoarchitectural organization. Here, we generated spatially-resolved transcriptomics (SRT) and single-nucleus RNA-sequencing (snRNA-seq) data from adjacent tissue sections of the anterior human hippocampus across ten adult neurotypical donors. We defined molecular profiles for hippocampal cell types and spatial domains. Using non-negative matrix factorization and transfer learning, we integrated these data to define gene expression patterns within the snRNA-seq data and infer the expression of these patterns in the SRT data. With this approach, we leveraged existing rodent datasets that feature information on circuit connectivity and neural activity induction to make predictions about axonal projection targets and likelihood of ensemble recruitment in spatially-defined cellular populations of the human hippocampus. Finally, we integrated genome-wide association studies with transcriptomic data to identify enrichment of genetic components for neurodevelopmental, neuropsychiatric, and neurodegenerative disorders across cell types, spatial domains, and gene expression patterns of the human hippocampus. To make this comprehensive molecular atlas accessible to the scientific community, both raw and processed data are freely available, including through interactive web applications.

Inferring cellular and molecular processes in single-cell data with non-negative matrix factorization using Python, R and GenePattern Notebook implementations of CoGAPS

Non-negative matrix factorization (NMF) is an unsupervised learning method well suited to high-throughput biology. However, inferring biological processes from an NMF result still requires additional post hoc statistics and annotation for interpretation of learned features. Here, we introduce a suite of computational tools that implement NMF and provide methods for accurate and clear biological interpretation and analysis. A generalized discussion of NMF covering its benefits, limitations and open questions is followed by four procedures for the Bayesian NMF algorithm Coordinated Gene Activity across Pattern Subsets (CoGAPS). Each procedure will demonstrate NMF analysis to quantify cell state transitions in a public domain single-cell RNA-sequencing dataset. The first demonstrates PyCoGAPS, our new Python implementation that enhances runtime for large datasets, and the second allows its deployment in Docker. The third procedure steps through the same single-cell NMF analysis using our R CoGAPS interface. The fourth introduces a beginner-friendly CoGAPS platform using GenePattern Notebook, aimed at users with a working conceptual knowledge of data analysis but without a basic proficiency in the R or Python programming language. We also constructed a user-facing website to serve as a central repository for information and instructional materials about CoGAPS and its application programming interfaces. The expected timing to setup the packages and conduct a test run is around 15 min, and an additional 30 min to conduct analyses on a precomputed result. The expected runtime on the user's desired dataset can vary from hours to days depending on factors such as dataset size or input parameters.

Digitize your Biology! Modeling multicellular systems through interpretable cell behavior

Cells are fundamental units of life, constantly interacting and evolving as dynamical systems. While recent spatial multi-omics can quantitate individual cells' characteristics and regulatory programs, forecasting their evolution ultimately requires mathematical modeling. We develop a conceptual framework-a cell behavior hypothesis grammar-that uses natural language statements (cell rules) to create mathematical models. This allows us to systematically integrate biological knowledge and multi-omics data to make them computable. We can then perform virtual "thought experiments" that challenge and extend our understanding of multicellular systems, and ultimately generate new testable hypotheses. In this paper, we motivate and describe the grammar, provide a reference implementation, and demonstrate its potential through a series of examples in tumor biology and immunotherapy. Altogether, this approach provides a bridge between biological, clinical, and systems biology researchers for mathematical modeling of biological systems at scale, allowing the community to extrapolate from single-cell characterization to emergent multicellular behavior.

A primer on correlation-based dimension reduction methods for multi-omics analysis

The continuing advances of omic technologies mean that it is now more tangible to measure the numerous features collectively reflecting the molecular properties of a sample. When multiple omic methods are used, statistical and computational approaches can exploit these large, connected profiles. Multi-omics is the integration of different omic data sources from the same biological sample. In this review, we focus on correlation-based dimension reduction approaches for single omic datasets, followed by methods for pairs of omics datasets, before detailing further techniques for three or more omic datasets. We also briefly detail network methods when three or more omic datasets are available and which complement correlation-oriented tools. To aid readers new to this area, these are all linked to relevant R packages that can implement these procedures. Finally, we discuss scenarios of experimental design and present road maps that simplify the selection of appropriate analysis methods. This review will help researchers navigate emerging methods for multi-omics and integrating diverse omic datasets appropriately. This raises the opportunity of implementing population multi-omics with large sample sizes as omics technologies and our understanding improve.

iPSC-based modeling of helicase deficiency reveals impaired cell proliferation and increased apoptosis after NK cell lineage commitment

Cell proliferation is a ubiquitous process required for organismal development and homeostasis. However, individuals with partial loss-of-function variants in DNA replicative helicase components often present with immunodeficiency due to specific loss of natural killer (NK) cells. Such lineage-specific disease phenotypes raise questions on how the proliferation is regulated in cell type-specific manner. We aimed to understand NK cell-specific proliferative dynamics and vulnerability to impaired helicase function using iPSCs from individuals with NK cell deficiency (NKD) due to hereditary compound heterozygous GINS4 variants. We observed and characterized heterogeneous cell populations that arise during the iPSC differentiation along with NK cells. While overall cell proliferation decreased with differentiation, early NK cell precursors showed a short burst of cell proliferation. GINS4 deficiency induced replication stress in these early NK cell precursors, which are poised for apoptosis, and ultimately recapitulate the NKD phenotype.

Activity-regulated gene expression across cell types of the mouse hippocampus

Activity-regulated gene (ARG) expression patterns in the hippocampus (HPC) regulate synaptic plasticity, learning, and memory, and are linked to both risk and treatment responses for many neuropsychiatric disorders. The HPC contains discrete classes of neurons with specialized functions, but cell type-specific activity-regulated transcriptional programs are not well characterized. Here, we used single-nucleus RNA-sequencing (snRNA-seq) in a mouse model of acute electroconvulsive seizures (ECS) to identify cell type-specific molecular signatures associated with induced activity in HPC neurons. We used unsupervised clustering and a priori marker genes to computationally annotate 15,990 high-quality HPC neuronal nuclei from N = 4 mice across all major HPC subregions and neuron types. Activity-induced transcriptomic responses were divergent across neuron populations, with dentate granule cells being particularly responsive to activity. Differential expression analysis identified both upregulated and downregulated cell type-specific gene sets in neurons following ECS. Within these gene sets, we identified enrichment of pathways associated with varying biological processes such as synapse organization, cellular signaling, and transcriptional regulation. Finally, we used matrix factorization to reveal continuous gene expression patterns differentially associated with cell type, ECS, and biological processes. This work provides a rich resource for interrogating activity-regulated transcriptional responses in HPC neurons at single-nuclei resolution in the context of ECS, which can provide biological insight into the roles of defined neuronal subtypes in HPC function.

Interface-guided phenotyping of coding variants in the transcription factor RUNX1 with SEUSS

Understanding the consequences of single amino acid substitutions in cancer driver genes remains an unmet need. Perturb-seq provides a tool to investigate the effects of individual mutations on cellular programs. Here we deploy SEUSS, a Perturb-seq like approach, to generate and assay mutations at physical interfaces of the RUNX1 Runt domain. We measured the impact of 115 mutations on RNA profiles in single myelogenous leukemia cells and used the profiles to categorize mutations into three functionally distinct groups: wild-type (WT)-like, loss-of-function (LOF)-like and hypomorphic. Notably, the largest concentration of functional mutations (non-WT-like) clustered at the DNA binding site and contained many of the more frequently observed mutations in human cancers. Hypomorphic variants shared characteristics with loss of function variants but had gene expression profiles indicative of response to neural growth factor and cytokine recruitment of neutrophils. Additionally, DNA accessibility changes upon perturbations were enriched for RUNX1 binding motifs, particularly near differentially expressed genes. Overall, our work demonstrates the potential of targeting protein interaction interfaces to better define the landscape of prospective phenotypes reachable by amino acid substitutions.

Machine learning in rare disease

High-throughput profiling methods (such as genomics or imaging) have accelerated basic research and made deep molecular characterization of patient samples routine. These approaches provide a rich portrait of genes, molecular pathways and cell types involved in disease phenotypes. Machine learning (ML) can be a useful tool for extracting disease-relevant patterns from high-dimensional datasets. However, depending upon the complexity of the biological question, machine learning often requires many samples to identify recurrent and biologically meaningful patterns. Rare diseases are inherently limited in clinical cases, leading to few samples to study. In this Perspective, we outline the challenges and emerging solutions for using ML for small sample sets, specifically in rare diseases. Advances in ML methods for rare diseases are likely to be informative for applications beyond rare diseases for which few samples exist with high-dimensional data. We propose that the method community prioritize the development of ML techniques for rare disease research.

Uncovering the spatial landscape of molecular interactions within the tumor microenvironment through latent spaces

Recent advances in spatial transcriptomics (STs) enable gene expression measurements from a tissue sample while retaining its spatial context. This technology enables unprecedented in situ resolution of the regulatory pathways that underlie the heterogeneity in the tumor as well as the tumor microenvironment (TME). The direct characterization of cellular co-localization with spatial technologies facilities quantification of the molecular changes resulting from direct cell-cell interaction, as it occurs in tumor-immune interactions. We present SpaceMarkers, a bioinformatics algorithm to infer molecular changes from cell-cell interactions from latent space analysis of ST data. We apply this approach to infer the molecular changes from tumor-immune interactions in Visium spatial transcriptomics data of metastasis, invasive and precursor lesions, and immunotherapy treatment. Further transfer learning in matched scRNA-seq data enabled further quantification of the specific cell types in which SpaceMarkers are enriched. Altogether, SpaceMarkers can identify the location and context-specific molecular interactions within the TME from ST data.

Molecular and cellular evolution of the amygdala across species analyzed by single-nucleus transcriptome profiling

The amygdala, or an amygdala-like structure, is found in the brains of all vertebrates and plays a critical role in survival and reproduction. However, the cellular architecture of the amygdala and how it has evolved remain elusive. Here, we generated single-nucleus RNA-sequencing data for more than 200,000 cells in the amygdala of humans, macaques, mice, and chickens. Abundant neuronal cell types from different amygdala subnuclei were identified in all datasets. Cross-species analysis revealed that inhibitory neurons and inhibitory neuron-enriched subnuclei of the amygdala were well-conserved in cellular composition and marker gene expression, whereas excitatory neuron-enriched subnuclei were relatively divergent. Furthermore, LAMP5⁺ interneurons were much more abundant in primates, while DRD2⁺ inhibitory neurons and LAMP5⁺SATB2⁺ excitatory neurons were dominant in the human central amygdalar nucleus (CEA) and basolateral amygdalar complex (BLA), respectively. We also identified CEA-like neurons and their species-specific distribution patterns in chickens. This study highlights the extreme cell-type diversity in the amygdala and reveals the conservation and divergence of cell types and gene expression patterns across species that may contribute to species-specific adaptations.

Integrated T cell cytometry metrics for immune-monitoring applications in immunotherapy clinical trials

Mass cytometry, or cytometry by TOF (CyTOF), provides a robust means of determining protein-level measurements of more than 40 markers simultaneously. While the functional states of immune cells occur along continuous phenotypic transitions, cytometric studies surveying cell phenotypes often rely on static metrics, such as discrete cell-type abundances, based on canonical markers and/or restrictive gating strategies. To overcome this limitation, we applied single-cell trajectory inference and nonnegative matrix factorization methods to CyTOF data to trace the dynamics of T cell states. In the setting of cancer immunotherapy, we showed that patient-specific summaries of continuous phenotypic shifts in T cells could be inferred from peripheral blood-derived CyTOF mass cytometry data. We further illustrated that transfer learning enabled these T cell continuous metrics to be used to estimate patient-specific cell states in new sample cohorts from a reference patient data set. Our work establishes the utility of continuous metrics for CyTOF analysis as tools for translational discovery.

Non-negative Independent Factor Analysis disentangles discrete and continuous sources of variation in scRNA-seq data

MOTIVATION

Single-cell RNA-seq analysis has emerged as a powerful tool for understanding inter-cellular heterogeneity. Due to the inherent noise of the data, computational techniques often rely on dimensionality reduction (DR) as both a pre-processing step and an analysis tool. Ideally, DR should preserve the biological information while discarding the noise. However, if the DR is to be used directly to gain biological insight it must also be interpretable-that is the individual dimensions of the reduction should correspond to specific biological variables such as cell-type identity or pathway activity. Maximizing biological interpretability necessitates making assumption about the data structures and the choice of the model is critical.

RESULTS

We present a new probabilistic single-cell factor analysis model, Non-negative Independent Factor Analysis (NIFA), that incorporates different interpretability inducing assumptions into a single modeling framework. The key advantage of our NIFA model is that it simultaneously models uni- and multi-modal latent factors, and thus isolates discrete cell-type identity and continuous pathway activity into separate components. We apply our approach to a range of datasets where cell-type identity is known, and we show that NIFA-derived factors outperform results from ICA, PCA, NMF and scCoGAPS (an NMF method designed for single-cell data) in terms of disentangling biological sources of variation. Studying an immunotherapy dataset in detail, we show that NIFA is able to reproduce and refine previous findings in a single analysis framework and enables the discovery of new clinically relevant cell states.

AVAILABILITY AND IMPLEMENTATION

NFIA is a R package which is freely available at GitHub (https://github.com/wgmao/NIFA). The test dataset is archived at https://zenodo.org/record/6286646.

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

Entinostat Decreases Immune Suppression to Promote Antitumor Responses in a HER2+ Breast Tumor Microenvironment

Therapeutic combinations to alter immunosuppressive, solid tumor microenvironments (TME), such as in breast cancer, are essential to improve responses to immune checkpoint inhibitors (ICI). Entinostat, an oral histone deacetylase inhibitor, has been shown to improve responses to ICIs in various tumor models with immunosuppressive TMEs. The precise and comprehensive alterations to the TME induced by entinostat remain unknown. Here, we employed single-cell RNA sequencing on HER2-overexpressing breast tumors from mice treated with entinostat and ICIs to fully characterize changes across multiple cell types within the TME. This analysis demonstrates that treatment with entinostat induced a shift from a protumor to an antitumor TME signature, characterized predominantly by changes in myeloid cells. We confirmed myeloid-derived suppressor cells (MDSC) within entinostat-treated tumors associated with a less suppressive granulocytic (G)-MDSC phenotype and exhibited altered suppressive signaling that involved the NFκB and STAT3 pathways. In addition to MDSCs, tumor-associated macrophages were epigenetically reprogrammed from a protumor M2-like phenotype toward an antitumor M1-like phenotype, which may be contributing to a more sensitized TME. Overall, our in-depth analysis suggests that entinostat-induced changes on multiple myeloid cell types reduce immunosuppression and increase antitumor responses, which, in turn, improve sensitivity to ICIs. Sensitization of the TME by entinostat could ultimately broaden the population of patients with breast cancer who could benefit from ICIs.

Learning biologically-interpretable latent representations for gene expression data: Pathway Activity Score Learning Algorithm

Molecular gene-expression datasets consist of samples with tens of thousands of measured quantities (i.e., high dimensional data). However, lower-dimensional representations that retain the useful biological information do exist. We present a novel algorithm for such dimensionality reduction called Pathway Activity Score Learning (PASL). The major novelty of PASL is that the constructed features directly correspond to known molecular pathways (genesets in general) and can be interpreted as pathway activity scores. Hence, unlike PCA and similar methods, PASL's latent space has a fairly straightforward biological interpretation. PASL is shown to outperform in predictive performance the state-of-the-art method (PLIER) on two collections of breast cancer and leukemia gene expression datasets. PASL is also trained on a large corpus of 50000 gene expression samples to construct a universal dictionary of features across different tissues and pathologies. The dictionary validated on 35643 held-out samples for reconstruction error. It is then applied on 165 held-out datasets spanning a diverse range of diseases. The AutoML tool JADBio is employed to show that the predictive information in the PASL-created feature space is retained after the transformation. The code is available at https://github.com/mensxmachina/PASL.

DAMPs/PAMPs induce monocytic TLR activation and tolerance in COVID-19 patients; nucleic acid binding scavengers can counteract such TLR agonists

Millions of COVID-19 patients have succumbed to respiratory and systemic inflammation. Hyperstimulation of toll-like receptor (TLR) signaling is a key driver of immunopathology following infection by viruses. We found that severely ill COVID-19 patients in the Intensive Care Unit (ICU) display hallmarks of such hyper-stimulation with abundant agonists of nucleic acid-sensing TLRs present in their blood and lungs. These nucleic acid-containing Damage and Pathogen Associated Molecular Patterns (DAMPs/PAMPs) can be depleted using nucleic acid-binding microfibers to limit the patient samples' ability to hyperactivate such innate immune receptors. Single-cell RNA-sequencing revealed that CD16+ monocytes from deceased but not recovered ICU patients exhibit a TLR-tolerant phenotype and a deficient anti-viral response after ex vivo TLR stimulation. Plasma proteomics confirmed such myeloid hyperactivation and revealed DAMP/PAMP carrier consumption in deceased patients. Treatment of these COVID-19 patient samples with MnO nanoparticles effectively neutralizes TLR activation by the abundant nucleic acid-containing DAMPs/PAMPs present in their lungs and blood. Finally, MnO nanoscavenger treatment limits the ability of DAMPs/PAMPs to induce TLR tolerance in monocytes. Thus, treatment with microfiber- or nanoparticle-based DAMP/PAMP scavengers may prove useful for limiting SARS-CoV-2 induced hyperinflammation, preventing monocytic TLR tolerance, and improving outcomes in severely ill COVID-19 patients.

A phenotypically supervised single-cell analysis protocol to study within-cell-type heterogeneity of cultured mammalian cells

Here, we describe a protocol combining functional metrics with genomic data to elucidate drivers of within-cell-type heterogeneity via the phenotype-to-genotype link. This technique involves using fluorescence tagging to label and isolate cells grown in 3D culture, enabling high-throughput enrichment of phenotypically defined cell subpopulations by fluorescence-activated cell sorting. We then perform a validated phenotypically supervised single-cell analysis pipeline to reveal unique functional cell states, including genes and pathways that contribute to cellular heterogeneity and were undetectable by unsupervised analysis. For complete details on the use and execution of this protocol, please refer to Chen et al. (2020).

Variation of Human Neural Stem Cells Generating Organizer States In Vitro before Committing to Cortical Excitatory or Inhibitory Neuronal Fates

Better understanding of the progression of neural stem cells (NSCs) in the developing cerebral cortex is important for modeling neurogenesis and defining the pathogenesis of neuropsychiatric disorders. Here, we use RNA sequencing, cell imaging, and lineage tracing of mouse and human in vitro NSCs and monkey brain sections to model the generation of cortical neuronal fates. We show that conserved signaling mechanisms regulate the acute transition from proliferative NSCs to committed glutamatergic excitatory neurons. As human telencephalic NSCs develop from pluripotency in vitro, they transition through organizer states that spatially pattern the cortex before generating glutamatergic precursor fates. NSCs derived from multiple human pluripotent lines vary in these early patterning states, leading differentially to dorsal or ventral telencephalic fates. This work furthers systematic analyses of the earliest patterning events that generate the major neuronal trajectories of the human telencephalon.

Micali et al. report that human telencephalic NSCs in vitro transition through the organizer states that pattern the neocortex. Human pluripotent lines vary in organizer formation, generating divergent neuronal differentiation trajectories biased toward dorsal or ventral telencephalic fates and opening further analysis of the earliest cortical specification events.

GeneWalk identifies relevant gene functions for a biological context using network representation learning

A bottleneck in high-throughput functional genomics experiments is identifying the most important genes and their relevant functions from a list of gene hits. Gene Ontology (GO) enrichment methods provide insight at the gene set level. Here, we introduce GeneWalk ( github.com/churchmanlab/genewalk ) that identifies individual genes and their relevant functions critical for the experimental setting under examination. After the automatic assembly of an experiment-specific gene regulatory network, GeneWalk uses representation learning to quantify the similarity between vector representations of each gene and its GO annotations, yielding annotation significance scores that reflect the experimental context. By performing gene- and condition-specific functional analysis, GeneWalk converts a list of genes into data-driven hypotheses.

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.

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.

Extracellular Vesicles Released by Tumor Endothelial Cells Spread Immunosuppressive and Transforming Signals Through Various Recipient Cells

Head and neck squamous cell carcinoma (HNSCC) has a high recurrence and metastatic rate with an unknown mechanism of cancer spread. Tumor inflammation is the most critical processes of cancer onset, growth, and metastasis. We hypothesize that the release of extracellular vesicles (EVs) by tumor endothelial cells (TECs) induce reprogramming of immune cells as well as stromal cells to create an immunosuppressive microenvironment that favor tumor spread. We call this mechanism as non-metastatic contagious carcinogenesis. Extracellular vesicles were collected from primary HNSCC-derived endothelial cells (TEC-EV) and were used for stimulation of peripheral blood mononuclear cells (PBMCs) and primary adipose mesenchymal stem cells (ASCs). Regulation of ASC gene expression was investigated by RNA sequencing and protein array. PBMC, stimulated with TEC-EV, were analyzed by enzyme-linked immunosorbent assay and fluorescence-activated cell sorting. We validated in vitro the effects of TEC-EV on ASCs or PBMC by measuring invasion, adhesion, and proliferation. We found and confirmed that TEC-EV were able to change ASC inflammatory gene expression signature within 24-48 h. TEC-EV were also able to enhance the secretion of TGF-β1 and IL-10 by PBMC and to increase T regulatory cell (Treg) expansion. TEC-EV carry specific proteins and RNAs that are responsible for Treg differentiation and immune suppression. ASCs and PBMC, treated with TEC-EV, enhanced proliferation, adhesion of tumor cells, and their invasion. These data indicate that TEC-EV exhibit a mechanism of non-metastatic contagious carcinogenesis that regulates tumor microenvironment and reprograms immune cells to sustain tumor growth and progression.

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.

Decomposing Cell Identity for Transfer Learning across Cellular Measurements, Platforms, Tissues, and Species

Analysis of gene expression in single cells allows for decomposition of cellular states as low-dimensional latent spaces. However, the interpretation and validation of these spaces remains a challenge. Here, we present scCoGAPS, which defines latent spaces from a source single-cell RNA-sequencing (scRNA-seq) dataset, and projectR, which evaluates these latent spaces in independent target datasets via transfer learning. Application of developing mouse retina to scRNA-Seq reveals intrinsic relationships across biological contexts and assays while avoiding batch effects and other technical features. We compare the dimensions learned in this source dataset to adult mouse retina, a time-course of human retinal development, select scRNA-seq datasets from developing brain, chromatin accessibility data, and a murine-cell type atlas to identify shared biological features. These tools lay the groundwork for exploratory analysis of scRNA-seq data via latent space representations, enabling a shift in how we compare and identify cells beyond reliance on marker genes or ensemble molecular identity.

Large-scale Analyses of Disease Biomarkers and Apremilast Pharmacodynamic Effects

Finding biomarkers that provide shared link between disease severity, drug-induced pharmacodynamic effects and response status in human trials can provide number of values for patient benefits: elucidating current therapeutic mechanism-of-action, and, back-translating to fast-track development of next-generation therapeutics. Both opportunities are predicated on proactive generation of human molecular profiles that capture longitudinal trajectories before and after pharmacological intervention. Here, we present the largest plasma proteomic biomarker dataset available to-date and the corresponding analyses from placebo-controlled Phase III clinical trials of the phosphodiesterase type 4 inhibitor apremilast in psoriasis (PSOR), psoriatic arthritis (PsA), and ankylosing spondylitis (AS) from 526 subjects overall. Using approximately 150 plasma analytes tracked across three time points, we identified IL-17A and KLK-7 as biomarkers for disease severity and apremilast pharmacodynamic effect in psoriasis patients. Combined decline rate of KLK-7, PEDF, MDC and ANGPTL4 by Week 16 represented biomarkers for the responder subgroup, shedding insights into therapeutic mechanisms. In ankylosing spondylitis patients, IL-6 and LRG-1 were identified as biomarkers with concordance to disease severity. Apremilast-induced LRG-1 increase was consistent with the overall lack of efficacy in ankylosing spondylitis. Taken together, these findings expanded the mechanistic knowledge base of apremilast and provided translational foundations to accelerate future efforts including compound differentiation, combination, and repurposing.

Pathway Activity Score Learning for Dimensionality Reduction of Gene Expression Data

Molecular gene-expression datasets consist of samples with tens of thousands of measured quantities (e.g., high dimensional data). However, there exist lower-dimensional representations that retain the useful information. We present a novel algorithm for such dimensionality reduction called Pathway Activity Score Learning (PASL). The major novelty of PASL is that the constructed features directly correspond to known molecular pathways and can be interpreted as pathway activity scores. Hence, unlike PCA and similar methods, PASL’s latent space has a relatively straight-forward biological interpretation. As a use-case, PASL is applied on two collections of breast cancer and leukemia gene expression datasets. We show that PASL does retain the predictive information for disease classification on new, unseen datasets, as well as outperforming PLIER, a recently proposed competitive method. We also show that differential activation pathway analysis provides complementary information to standard gene set enrichment analysis. The code is available at https://github.com/mensxmachina/PASL.

DeepMF: deciphering the latent patterns in omics profiles with a deep learning method

BACKGROUND

With recent advances in high-throughput technologies, matrix factorization techniques are increasingly being utilized for mapping quantitative omics profiling matrix data into low-dimensional embedding space, in the hope of uncovering insights in the underlying biological processes. Nevertheless, current matrix factorization tools fall short in handling noisy data and missing entries, both deficiencies that are often found in real-life data.

RESULTS

Here, we propose DeepMF, a deep neural network-based factorization model. DeepMF disentangles the association between molecular feature-associated and sample-associated latent matrices, and is tolerant to noisy and missing values. It exhibited feasible cancer subtype discovery efficacy on mRNA, miRNA, and protein profiles of medulloblastoma cancer, leukemia cancer, breast cancer, and small-blue-round-cell cancer, achieving the highest clustering accuracy of 76%, 100%, 92%, and 100% respectively. When analyzing data sets with 70% missing entries, DeepMF gave the best recovery capacity with silhouette values of 0.47, 0.6, 0.28, and 0.44, outperforming other state-of-the-art MF tools on the cancer data sets Medulloblastoma, Leukemia, TCGA BRCA, and SRBCT. Its embedding strength as measured by clustering accuracy is 88%, 100%, 84%, and 96% on these data sets, which improves on the current best methods 76%, 100%, 78%, and 87%.

CONCLUSION

DeepMF demonstrated robust denoising, imputation, and embedding ability. It offers insights to uncover the underlying biological processes such as cancer subtype discovery. Our implementation of DeepMF can be found at https://github.com/paprikachan/DeepMF.

Glucocorticoid exposure during hippocampal neurogenesis primes future stress response by inducing changes in DNA methylation

Prenatal stress exposure is associated with risk for psychiatric disorders later in life. This may be mediated in part via enhanced exposure to glucocorticoids (GCs), which are known to impact neurogenesis. We aimed to identify molecular mediators of these effects, focusing on long-lasting epigenetic changes. In a human hippocampal progenitor cell (HPC) line, we assessed the short- and long-term effects of GC exposure during neurogenesis on messenger RNA (mRNA) expression and DNA methylation (DNAm) profiles. GC exposure induced changes in DNAm at 27,812 CpG dinucleotides and in the expression of 3,857 transcripts (false discovery rate [FDR] ≤ 0.1 and absolute fold change [FC] expression ≥ 1.15). HPC expression and GC-affected DNAm profiles were enriched for changes observed during human fetal brain development. Differentially methylated sites (DMSs) with GC exposure clustered into 4 trajectories over HPC differentiation, with transient as well as long-lasting DNAm changes. Lasting DMSs mapped to distinct functional pathways and were selectively enriched for poised and bivalent enhancer marks. Lasting DMSs had little correlation with lasting expression changes but were associated with a significantly enhanced transcriptional response to a second acute GC challenge. A significant subset of lasting DMSs was also responsive to an acute GC challenge in peripheral blood. These tissue-overlapping DMSs were used to compute a polyepigenetic score that predicted exposure to conditions associated with altered prenatal GCs in newborn's cord blood DNA. Overall, our data suggest that early exposure to GCs can change the set point of future transcriptional responses to stress by inducing lasting DNAm changes. Such altered set points may relate to differential vulnerability to stress exposure later in life.

Discovering Pathway and Cell Type Signatures in Transcriptomic Compendia with Machine Learning

Pathway and cell type signatures are patterns present in transcriptome data that are associated with biological processes or phenotypic consequences. These signatures result from specific cell type and pathway expression but can require large transcriptomic compendia to detect. Machine learning techniques can be powerful tools for signature discovery through their ability to provide accurate and interpretable results. In this review, we discuss various machine learning applications to extract pathway and cell type signatures from transcriptomic compendia. We focus on the biological motivations and interpretation for both supervised and unsupervised learning approaches in this setting. We consider recent advances, including deep learning, and their applications to expanding bulk and single-cell RNA data. As data and computational resources increase, there will be more opportunities for machine learning to aid in revealing biological signatures.

Hypoxia tolerance in the Norrin-deficient retina and the chronically hypoxic brain studied at single-cell resolution

The mammalian CNS is capable of tolerating chronic hypoxia, but cell type-specific responses to this stress have not been systematically characterized. In the Norrin KO (Ndp ^KO ) mouse, a model of familial exudative vitreoretinopathy (FEVR), developmental hypovascularization of the retina produces chronic hypoxia of inner nuclear-layer (INL) neurons and Muller glia. We used single-cell RNA sequencing, untargeted metabolomics, and metabolite labeling from ¹³C-glucose to compare WT and Ndp ^KO retinas. In Ndp ^KO retinas, we observe gene expression responses consistent with hypoxia in Muller glia and retinal neurons, and we find a metabolic shift that combines reduced flux through the TCA cycle with increased synthesis of serine, glycine, and glutathione. We also used single-cell RNA sequencing to compare the responses of individual cell types in Ndp ^KO retinas with those in the hypoxic cerebral cortex of mice that were housed for 1 week in a reduced oxygen environment (7.5% oxygen). In the hypoxic cerebral cortex, glial transcriptome responses most closely resemble the response of Muller glia in the Ndp ^KO retina. In both retina and brain, vascular endothelial cells activate a previously dormant tip cell gene expression program, which likely underlies the adaptive neoangiogenic response to chronic hypoxia. These analyses of retina and brain transcriptomes at single-cell resolution reveal both shared and cell type-specific changes in gene expression in response to chronic hypoxia, implying both shared and distinct cell type-specific physiologic responses.

Dynamic modeling and network approaches for omics time course data: overview of computational approaches and applications

Inferring networks and dynamics of genes, proteins, cells and other biological entities from high-throughput biological omics data is a central and challenging issue in computational and systems biology. This is essential for understanding the complexity of human health, disease susceptibility and pathogenesis for Predictive, Preventive, Personalized and Participatory (P4) system and precision medicine. The delineation of the possible interactions of all genes/proteins in a genome/proteome is a task for which conventional experimental techniques are ill suited. Urgently needed are rapid and inexpensive computational and statistical methods that can identify interacting candidate disease genes or drug targets out of thousands that can be further investigated or validated by experimentations. Moreover, identifying biological dynamic systems, and simultaneously estimating the important kinetic structural and functional parameters, which may not be experimentally accessible could be important directions for drug-disease-gene network studies. In this article, we present an overview and comparison of recent developments of dynamic modeling and network approaches for time-course omics data, and their applications to various biological systems, health conditions and disease statuses. Moreover, various data reduction and analytical schemes ranging from mathematical to computational to statistical methods are compared including their merits, drawbacks and limitations. The most recent software, associated web resources and other potentials for the compared methods are also presented and discussed in detail.

A review on machine learning principles for multi-view biological data integration

Driven by high-throughput sequencing techniques, modern genomic and clinical studies are in a strong need of integrative machine learning models for better use of vast volumes of heterogeneous information in the deep understanding of biological systems and the development of predictive models. How data from multiple sources (called multi-view data) are incorporated in a learning system is a key step for successful analysis. In this article, we provide a comprehensive review on omics and clinical data integration techniques, from a machine learning perspective, for various analyses such as prediction, clustering, dimension reduction and association. We shall show that Bayesian models are able to use prior information and model measurements with various distributions; tree-based methods can either build a tree with all features or collectively make a final decision based on trees learned from each view; kernel methods fuse the similarity matrices learned from individual views together for a final similarity matrix or learning model; network-based fusion methods are capable of inferring direct and indirect associations in a heterogeneous network; matrix factorization models have potential to learn interactions among features from different views; and a range of deep neural networks can be integrated in multi-modal learning for capturing the complex mechanism of biological systems.

Aging-like Spontaneous Epigenetic Silencing Facilitates Wnt Activation, Stemness, and BrafV600E-Induced Tumorigenesis

We addressed the precursor role of aging-like spontaneous promoter DNA hypermethylation in initiating tumorigenesis. Using mouse colon-derived organoids, we show that promoter hypermethylation spontaneously arises in cells mimicking the human aging-like phenotype. The silenced genes activate the Wnt pathway, causing a stem-like state and differentiation defects. These changes render aged organoids profoundly more sensitive than young ones to transformation by BrafV600E, producing the typical human proximal BRAFV600E-driven colon adenocarcinomas characterized by extensive, abnormal gene-promoter CpG-island methylation, or the methylator phenotype (CIMP). Conversely, CRISPR-mediated simultaneous inactivation of a panel of the silenced genes markedly sensitizes to BrafV600E-induced transformation. Our studies tightly link aging-like epigenetic abnormalities to intestinal cell fate changes and predisposition to oncogene-driven colon tumorigenesis.

A Mechanism of Resistance to Antibody-Targeted Immune Attack

Targeted monoclonal antibody therapy is a promising therapeutic strategy for cancer, and antibody-dependent cell-mediated cytotoxicity (ADCC) represents a crucial mechanism underlying these approaches. The majority of patients have limited responses to monoclonal antibody therapy due to the development of resistance. Models of ADCC provide a system for uncovering immune-resistance mechanisms. We continuously exposed epidermal growth factor receptor (EGFR+) A431 cells to KIR-deficient NK92-CD16V effector cells and the anti-EGFR cetuximab. Persistent ADCC exposure yielded ADCC-resistant cells (ADCCR1) that, compared with control ADCC-sensitive cells (ADCCS1), exhibited reduced EGFR expression, overexpression of histone- and interferon-related genes, and a failure to activate NK cells, without evidence of epithelial-to-mesenchymal transition. These properties gradually reversed following withdrawal of ADCC selection pressure. The development of resistance was associated with lower expression of multiple cell-surface molecules that contribute to cell-cell interactions and immune synapse formation. Classic immune checkpoints did not modulate ADCC in this unique model system of immune resistance. We showed that the induction of ADCC resistance involves genetic and epigenetic changes that lead to a general loss of target cell adhesion properties that are required for the establishment of an immune synapse, killer cell activation, and target cell cytotoxicity.

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.

Integrated time course omics analysis distinguishes immediate therapeutic response from acquired resistance

BACKGROUND

Targeted therapies specifically act by blocking the activity of proteins that are encoded by genes critical for tumorigenesis. However, most cancers acquire resistance and long-term disease remission is rarely observed. Understanding the time course of molecular changes responsible for the development of acquired resistance could enable optimization of patients' treatment options. Clinically, acquired therapeutic resistance can only be studied at a single time point in resistant tumors.

METHODS

To determine the dynamics of these molecular changes, we obtained high throughput omics data (RNA-sequencing and DNA methylation) weekly during the development of cetuximab resistance in a head and neck cancer in vitro model. The CoGAPS unsupervised algorithm was used to determine the dynamics of the molecular changes associated with resistance during the time course of resistance development.

RESULTS

CoGAPS was used to quantify the evolving transcriptional and epigenetic changes. Applying a PatternMarker statistic to the results from CoGAPS enabled novel heatmap-based visualization of the dynamics in these time course omics data. We demonstrate that transcriptional changes result from immediate therapeutic response or resistance, whereas epigenetic alterations only occur with resistance. Integrated analysis demonstrates delayed onset of changes in DNA methylation relative to transcription, suggesting that resistance is stabilized epigenetically.

CONCLUSIONS

Genes with epigenetic alterations associated with resistance that have concordant expression changes are hypothesized to stabilize the resistant phenotype. These genes include FGFR1, which was associated with EGFR inhibitors resistance previously. Thus, integrated omics analysis distinguishes the timing of molecular drivers of resistance. This understanding of the time course progression of molecular changes in acquired resistance is important for the development of alternative treatment strategies that would introduce appropriate selection of new drugs to treat cancer before the resistant phenotype develops.

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.

CoGAPS matrix factorization algorithm identifies transcriptional changes in AP-2alpha target genes in feedback from therapeutic inhibition of the EGFR network

Patients with oncogene driven tumors are treated with targeted therapeutics including EGFR inhibitors. Genomic data from The Cancer Genome Atlas (TCGA) demonstrates molecular alterations to EGFR, MAPK, and PI3K pathways in previously untreated tumors. Therefore, this study uses bioinformatics algorithms to delineate interactions resulting from EGFR inhibitor use in cancer cells with these genetic alterations. We modify the HaCaT keratinocyte cell line model to simulate cancer cells with constitutive activation of EGFR, HRAS, and PI3K in a controlled genetic background. We then measure gene expression after treating modified HaCaT cells with gefitinib, afatinib, and cetuximab. The CoGAPS algorithm distinguishes a gene expression signature associated with the anticipated silencing of the EGFR network. It also infers a feedback signature with EGFR gene expression itself increasing in cells that are responsive to EGFR inhibitors. This feedback signature has increased expression of several growth factor receptors regulated by the AP-2 family of transcription factors. The gene expression signatures for AP-2alpha are further correlated with sensitivity to cetuximab treatment in HNSCC cell lines and changes in EGFR expression in HNSCC tumors with low CDKN2A gene expression. In addition, the AP-2alpha gene expression signatures are also associated with inhibition of MEK, PI3K, and mTOR pathways in the Library of Integrated Network-Based Cellular Signatures (LINCS) data. These results suggest that AP-2 transcription factors are activated as feedback from EGFR network inhibition and may mediate EGFR inhibitor resistance.

Computational dynamic approaches for temporal omics data with applications to systems medicine

Modeling and predicting biological dynamic systems and simultaneously estimating the kinetic structural and functional parameters are extremely important in systems and computational biology. This is key for understanding the complexity of the human health, drug response, disease susceptibility and pathogenesis for systems medicine. Temporal omics data used to measure the dynamic biological systems are essentials to discover complex biological interactions and clinical mechanism and causations. However, the delineation of the possible associations and causalities of genes, proteins, metabolites, cells and other biological entities from high throughput time course omics data is challenging for which conventional experimental techniques are not suited in the big omics era. In this paper, we present various recently developed dynamic trajectory and causal network approaches for temporal omics data, which are extremely useful for those researchers who want to start working in this challenging research area. Moreover, applications to various biological systems, health conditions and disease status, and examples that summarize the state-of-the art performances depending on different specific mining tasks are presented. We critically discuss the merits, drawbacks and limitations of the approaches, and the associated main challenges for the years ahead. The most recent computing tools and software to analyze specific problem type, associated platform resources, and other potentials for the dynamic trajectory and interaction methods are also presented and discussed in detail.

Toward Signaling-Driven Biomarkers Immune to Normal Tissue Contamination

The goal of this study was to discover a minimally invasive pathway-specific biomarker that is immune to normal cell mRNA contamination for diagnosing head and neck squamous cell carcinoma (HNSCC). Using Elsevier's MedScan natural language processing component of the Pathway Studio software and the TRANSFAC database, we produced a curated set of genes regulated by the signaling networks driving the development of HNSCC. The network and its gene targets provided prior probabilities for gene expression, which guided our CoGAPS matrix factorization algorithm to isolate patterns related to HNSCC signaling activity from a microarray-based study. Using patterns that distinguished normal from tumor samples, we identified a reduced set of genes to analyze with Top Scoring Pair in order to produce a potential biomarker for HNSCC. Our proposed biomarker comprises targets of the transcription factor (TF) HIF1A and the FOXO family of TFs coupled with genes that show remarkable stability across all normal tissues. Based on validation with novel data from The Cancer Genome Atlas (TCGA), measured by RNAseq, and bootstrap sampling, the biomarker for normal vs. tumor has an accuracy of 0.77, a Matthews correlation coefficient of 0.54, and an area under the curve (AUC) of 0.82.

Preserving biological heterogeneity with a permuted surrogate variable analysis for genomics batch correction

MOTIVATION

Sample source, procurement process and other technical variations introduce batch effects into genomics data. Algorithms to remove these artifacts enhance differences between known biological covariates, but also carry potential concern of removing intragroup biological heterogeneity and thus any personalized genomic signatures. As a result, accurate identification of novel subtypes from batch-corrected genomics data is challenging using standard algorithms designed to remove batch effects for class comparison analyses. Nor can batch effects be corrected reliably in future applications of genomics-based clinical tests, in which the biological groups are by definition unknown a priori.

RESULTS

Therefore, we assess the extent to which various batch correction algorithms remove true biological heterogeneity. We also introduce an algorithm, permuted-SVA (pSVA), using a new statistical model that is blind to biological covariates to correct for technical artifacts while retaining biological heterogeneity in genomic data. This algorithm facilitated accurate subtype identification in head and neck cancer from gene expression data in both formalin-fixed and frozen samples. When applied to predict Human Papillomavirus (HPV) status, pSVA improved cross-study validation even if the sample batches were highly confounded with HPV status in the training set.

AVAILABILITY AND IMPLEMENTATION

All analyses were performed using R version 2.15.0. The code and data used to generate the results of this manuscript is available from https://sourceforge.net/projects/psva.

Pattern identification in time-course gene expression data with the CoGAPS matrix factorization

Patterns in time-course gene expression data can represent the biological processes that are active over the measured time period. However, the orthogonality constraint in standard pattern-finding algorithms, including notably principal components analysis (PCA), confounds expression changes resulting from simultaneous, non-orthogonal biological processes. Previously, we have shown that Markov chain Monte Carlo nonnegative matrix factorization algorithms are particularly adept at distinguishing such concurrent patterns. One such matrix factorization is implemented in the software package CoGAPS. We describe the application of this software and several technical considerations for identification of age-related patterns in a public, prefrontal cortex gene expression dataset.