Improving replicability in single-cell RNA-Seq cell type discovery with Dune.

BACKGROUND: Single-cell transcriptome sequencing (scRNA-Seq) has allowed new types of investigations at unprecedented levels of resolution. Among the primary goals of scRNA-Seq is the classification of cells into distinct types. Many approaches build on existing clustering literature to develop tools specific to single-cell. However, almost all of these methods rely on heuristics or user-supplied parameters to control the number of clusters. This affects both the resolution of the clusters within the original dataset as well as their replicability across datasets. While many recommendations exist, in general, there is little assurance that any given set of parameters will represent an optimal choice in the trade-off between cluster resolution and replicability. For instance, another set of parameters may result in more clusters that are also more replicable. RESULTS: Here, we propose Dune, a new method for optimizing the trade-off between the resolution of the clusters and their replicability. Our method takes as input a set of clustering results-or partitions-on a single dataset and iteratively merges clusters within each partitions in order to maximize their concordance between partitions. As demonstrated on multiple datasets from different platforms, Dune outperforms existing techniques, that rely on hierarchical merging for reducing the number of clusters, in terms of replicability of the resultant merged clusters as well as concordance with ground truth. Dune is available as an R package on Bioconductor: https://www.bioconductor.org/packages/release/bioc/html/Dune.html . CONCLUSIONS: Cluster refinement by Dune helps improve the robustness of any clustering analysis and reduces the reliance on tuning parameters. This method provides an objective approach for borrowing information across multiple clusterings to generate replicable clusters most likely to represent common biological features across multiple datasets.

Successional adaptive strategies revealed by correlating arbuscular mycorrhizal fungal abundance with host plant gene expression.

The shifts in adaptive strategies revealed by ecological succession and the mechanisms that facilitate these shifts are fundamental to ecology. These adaptive strategies could be particularly important in communities of arbuscular mycorrhizal fungi (AMF) mutualistic with sorghum, where strong AMF succession replaces initially ruderal species with competitive ones and where the strongest plant response to drought is to manage these AMF. Although most studies of agriculturally important fungi focus on parasites, the mutualistic symbionts, AMF, constitute a research system of human-associated fungi whose relative simplicity and synchrony are conducive to experimental ecology. First, we hypothesize that, when irrigation is stopped to mimic drought, competitive AMF species should be replaced by AMF species tolerant to drought stress. We then, for the first time, correlate AMF abundance and host plant transcription to test two novel hypotheses about the mechanisms behind the shift from ruderal to competitive AMF. Surprisingly, despite imposing drought stress, we found no stress-tolerant AMF, probably due to our agricultural system having been irrigated for nearly six decades. Remarkably, we found strong and differential correlation between the successional shift from ruderal to competitive AMF and sorghum genes whose products (i) produce and release strigolactone signals, (ii) perceive mycorrhizal-lipochitinoligosaccharide (Myc-LCO) signals, (iii) provide plant lipid and sugar to AMF, and (iv) import minerals and water provided by AMF. These novel insights frame new hypotheses about AMF adaptive evolution and suggest a rationale for selecting AMF to reduce inputs and maximize yields in commercial agriculture.

mbkmeans: Fast clustering for single cell data using mini-batch k-means.

Single-cell RNA-Sequencing (scRNA-seq) is the most widely used high-throughput technology to measure genome-wide gene expression at the single-cell level. One of the most common analyses of scRNA-seq data detects distinct subpopulations of cells through the use of unsupervised clustering algorithms. However, recent advances in scRNA-seq technologies result in current datasets ranging from thousands to millions of cells. Popular clustering algorithms, such as k-means, typically require the data to be loaded entirely into memory and therefore can be slow or impossible to run with large datasets. To address this problem, we developed the mbkmeans R/Bioconductor package, an open-source implementation of the mini-batch k-means algorithm. Our package allows for on-disk data representations, such as the common HDF5 file format widely used for single-cell data, that do not require all the data to be loaded into memory at one time. We demonstrate the performance of the mbkmeans package using large datasets, including one with 1.3 million cells. We also highlight and compare the computing performance of mbkmeans against the standard implementation of k-means and other popular single-cell clustering methods. Our software package is available in Bioconductor at https://bioconductor.org/packages/mbkmeans.

Transcriptomic analysis of field-droughted sorghum from seedling to maturity reveals biotic and metabolic responses.

Drought is the most important environmental stress limiting crop yields. The C4 cereal sorghum [Sorghum bicolor (L.) Moench] is a critical food, forage, and emerging bioenergy crop that is notably drought-tolerant. We conducted a large-scale field experiment, imposing preflowering and postflowering drought stress on 2 genotypes of sorghum across a tightly resolved time series, from plant emergence to postanthesis, resulting in a dataset of nearly 400 transcriptomes. We observed a fast and global transcriptomic response in leaf and root tissues with clear temporal patterns, including modulation of well-known drought pathways. We also identified genotypic differences in core photosynthesis and reactive oxygen species scavenging pathways, highlighting possible mechanisms of drought tolerance and of the delayed senescence, characteristic of the stay-green phenotype. Finally, we discovered a large-scale depletion in the expression of genes critical to arbuscular mycorrhizal (AM) symbiosis, with a corresponding drop in AM fungal mass in the plants' roots.

A Tetrahymena Piwi bound to mature tRNA 3' fragments activates the exonuclease Xrn2 for RNA processing in the nucleus.

Emerging evidence suggests that Argonaute (Ago)/Piwi proteins have diverse functions in the nucleus and cytoplasm, but the molecular mechanisms employed in the nucleus remain poorly defined. The Tetrahymena thermophila Ago/Piwi protein Twi12 is essential for growth and functions in the nucleus. Twi12-bound small RNAs (sRNAs) are 3' tRNA fragments that contain modified bases and thus are attenuated for base pairing to targets. We show that Twi12 assembles an unexpected complex with the nuclear exonuclease Xrn2. Twi12 functions to stabilize and localize Xrn2, as well as to stimulate its exonuclease activity. Twi12 function depends on sRNA binding, which is required for its nuclear import. Depletion of Twi12 or Xrn2 induces a cellular ribosomal RNA processing defect known to result from limiting Xrn2 activity in other organisms. Our findings suggest a role for an Ago/Piwi protein and 3' tRNA fragments in nuclear RNA metabolism.

clusterExperiment and RSEC: A Bioconductor package and framework for clustering of single-cell and other large gene expression datasets.

Clustering of genes and/or samples is a common task in gene expression analysis. The goals in clustering can vary, but an important scenario is that of finding biologically meaningful subtypes within the samples. This is an application that is particularly appropriate when there are large numbers of samples, as in many human disease studies. With the increasing popularity of single-cell transcriptome sequencing (RNA-Seq), many more controlled experiments on model organisms are similarly creating large gene expression datasets with the goal of detecting previously unknown heterogeneity within cells. It is common in the detection of novel subtypes to run many clustering algorithms, as well as rely on subsampling and ensemble methods to improve robustness. We introduce a Bioconductor R package, clusterExperiment, that implements a general and flexible strategy we entitle Resampling-based Sequential Ensemble Clustering (RSEC). RSEC enables the user to easily create multiple, competing clusterings of the data based on different techniques and associated tuning parameters, including easy integration of resampling and sequential clustering, and then provides methods for consolidating the multiple clusterings into a final consensus clustering. The package is modular and allows the user to separately apply the individual components of the RSEC procedure, i.e., apply multiple clustering algorithms, create a consensus clustering or choose tuning parameters, and merge clusters. Additionally, clusterExperiment provides a variety of visualization tools for the clustering process, as well as methods for the identification of possible cluster signatures or biomarkers. The R package clusterExperiment is publicly available through the Bioconductor Project, with a detailed manual (vignette) as well as well documented help pages for each function.

Slingshot: cell lineage and pseudotime inference for single-cell transcriptomics.

BACKGROUND: Single-cell transcriptomics allows researchers to investigate complex communities of heterogeneous cells. It can be applied to stem cells and their descendants in order to chart the progression from multipotent progenitors to fully differentiated cells. While a variety of statistical and computational methods have been proposed for inferring cell lineages, the problem of accurately characterizing multiple branching lineages remains difficult to solve. RESULTS: We introduce Slingshot, a novel method for inferring cell lineages and pseudotimes from single-cell gene expression data. In previously published datasets, Slingshot correctly identifies the biological signal for one to three branching trajectories. Additionally, our simulation study shows that Slingshot infers more accurate pseudotimes than other leading methods. CONCLUSIONS: Slingshot is a uniquely robust and flexible tool which combines the highly stable techniques necessary for noisy single-cell data with the ability to identify multiple trajectories. Accurate lineage inference is a critical step in the identification of dynamic temporal gene expression.

Clustering of mRNA-Seq data based on alternative splicing patterns.

Sequencing of messenger RNA (mRNA) can provide estimates of the levels of individual isoforms within the cell. It remains to adapt many standard statistical methods commonly used for analyzing gene expression levels to take advantage of this additional information. One novel question is whether we can find clusters of samples that are distinguished not by their gene expression but by their isoform usage. We propose a novel approach for clustering mRNA-Seq data that identifies such clusters. We show via simulation that our methods are more sensitive to finding clusters based on isoform usage than standard clustering techniques. We demonstrate its performance by finding a technical artifact that resulted in different batches having different isoform usage patterns, and illustrate its usage on several The Cancer Genome Atlas datasets.

Phylogenetic analyses of melanoma reveal complex patterns of metastatic dissemination.

Melanoma is difficult to treat once it becomes metastatic. However, the precise ancestral relationship between primary tumors and their metastases is not well understood. We performed whole-exome sequencing of primary melanomas and multiple matched metastases from eight patients to elucidate their phylogenetic relationships. In six of eight patients, we found that genetically distinct cell populations in the primary tumor metastasized in parallel to different anatomic sites, rather than sequentially from one site to the next. In five of these six patients, the metastasizing cells had themselves arisen from a common parental subpopulation in the primary, indicating that the ability to establish metastases is a late-evolving trait. Interestingly, we discovered that individual metastases were sometimes founded by multiple cell populations of the primary that were genetically distinct. Such establishment of metastases by multiple tumor subpopulations could help explain why identical resistance variants are identified in different sites after initial response to systemic therapy. One primary tumor harbored two subclones with different oncogenic mutations in CTNNB1, which were both propagated to the same metastasis, raising the possibility that activation of wingless-type mouse mammary tumor virus integration site (WNT) signaling may be involved, as has been suggested by experimental models.

Subtype and pathway specific responses to anticancer compounds in breast cancer.

Breast cancers are comprised of molecularly distinct subtypes that may respond differently to pathway-targeted therapies now under development. Collections of breast cancer cell lines mirror many of the molecular subtypes and pathways found in tumors, suggesting that treatment of cell lines with candidate therapeutic compounds can guide identification of associations between molecular subtypes, pathways, and drug response. In a test of 77 therapeutic compounds, nearly all drugs showed differential responses across these cell lines, and approximately one third showed subtype-, pathway-, and/or genomic aberration-specific responses. These observations suggest mechanisms of response and resistance and may inform efforts to develop molecular assays that predict clinical response.

Methods and challenges in timing chromosomal abnormalities within cancer samples.

MOTIVATION: Tumors acquire many chromosomal amplifications, and those acquired early in the lifespan of the tumor may be not only important for tumor growth but also can be used for diagnostic purposes. Many methods infer the order of the accumulation of abnormalities based on their occurrence in a large cohort of patients. Recently, Durinck et al. (2011) and Greenman et al. (2012) developed methods to order a single tumor's chromosomal amplifications based on the patterns of mutations accumulated within those regions. This method offers an unprecedented opportunity to assess the etiology of a single tumor sample, but has not been widely evaluated. RESULTS: We show that the model for timing chromosomal amplifications is limited in scope, particularly for regions with high levels of amplification. We also show that the estimation of the order of events can be sensitive for events that occur early in the progression of the tumor and that the partial maximum likelihood method of Greenman et al. (2012) can give biased estimates, particularly for moderate read coverage or normal contamination. We propose a maximum-likelihood estimation procedure that fully accounts for sequencing variability and show that it outperforms the partial maximum-likelihood estimation method. We also propose a Bayesian estimation procedure that stabilizes the estimates in certain settings. We implement these methods on a small number of ovarian tumors, and the results suggest possible differences in how the tumors acquired amplifications. AVAILABILITY AND IMPLEMENTATION: We provide implementation of these methods in an R package cancerTiming, which is available from the Comprehensive R Archive Network (CRAN) at http://CRAN.R-project.org/.

A single-cell atlas of IL-23 inhibition in cutaneous psoriasis distinguishes clinical response.

Psoriasis vulgaris and other chronic inflammatory diseases improve markedly with therapeutic blockade of interleukin-23 (IL-23) signaling, but the genetic mechanisms underlying clinical responses remain poorly understood. Using single-cell transcriptomics, we profiled immune cells isolated from lesional psoriatic skin before and during IL-23 blockade. In clinically responsive patients, a psoriatic transcriptional signature in skin-resident memory T cells was strongly attenuated. In contrast, poorly responsive patients were distinguished by persistent activation of IL-17-producing T (T17) cells, a mechanism distinct from alternative cytokine signaling or resistance isolated to epidermal keratinocytes. Even in IL-23 blockade-responsive patients, we detected a recurring set of recalcitrant, disease-specific transcriptional abnormalities. This irreversible immunological state may necessitate ongoing IL-23 inhibition. Spatial transcriptomic analyses also suggested that successful IL-23 blockade requires dampening of >90% of IL-17-induced response in lymphocyte-adjacent keratinocytes, an unexpectedly high threshold. Collectively, our data establish a patient-level paradigm for dissecting responses to immunomodulatory treatments.

Visualizing scRNA-Seq Data at Population Scale with GloScope.

Increasingly scRNA-Seq studies explore the heterogeneity of cell populations across different samples and its effect on an organism's phenotype. However, relatively few bioinformatic methods have been developed which adequately address the variation between samples for such population-level analyses. We propose a framework for representing the entire single-cell profile of a sample, which we call its GloScope representation. We implement GloScope on scRNA-Seq datasets from study designs ranging from 12 to over 300 samples. These examples demonstrate how GloScope allows researchers to perform essential bioinformatic tasks at the sample-level, in particular visualization and quality control assessment.

A Latent Activated Olfactory Stem Cell State Revealed by Single Cell Transcriptomic and Epigenomic Profiling.

The olfactory epithelium is one of the few regions of the nervous system that sustains neurogenesis throughout life. Its experimental accessibility makes it especially tractable for studying molecular mechanisms that drive neural regeneration after injury-induced cell death. In this study, we used single cell sequencing to identify major regulatory players in determining olfactory epithelial stem cell fate after acute injury. We combined gene expression and accessible chromatin profiles of individual lineage traced olfactory stem cells to predict transcription factor activity specific to different lineages and stages of recovery. We further identified a discrete stem cell state that appears poised for activation, characterized by accessible chromatin around wound response and lineage specific genes prior to their later expression in response to injury. Together these results provide evidence that a subset of quiescent olfactory epithelial stem cells are epigenetically primed to support injury-induced regeneration.

Co-occurrence networks reveal more complexity than community composition in resistance and resilience of microbial communities.

Plant response to drought stress involves fungi and bacteria that live on and in plants and in the rhizosphere, yet the stability of these myco- and micro-biomes remains poorly understood. We investigate the resistance and resilience of fungi and bacteria to drought in an agricultural system using both community composition and microbial associations. Here we show that tests of the fundamental hypotheses that fungi, as compared to bacteria, are (i) more resistant to drought stress but (ii) less resilient when rewetting relieves the stress, found robust support at the level of community composition. Results were more complex using all-correlations and co-occurrence networks. In general, drought disrupts microbial networks based on significant positive correlations among bacteria, among fungi, and between bacteria and fungi. Surprisingly, co-occurrence networks among functional guilds of rhizosphere fungi and leaf bacteria were strengthened by drought, and the same was seen for networks involving arbuscular mycorrhizal fungi in the rhizosphere. We also found support for the stress gradient hypothesis because drought increased the relative frequency of positive correlations.

Defining Patient-Level Molecular Heterogeneity in Psoriasis Vulgaris Based on Single-Cell Transcriptomics.

Identifying genetic variation underlying human diseases establishes targets for therapeutic development and helps tailor treatments to individual patients. Large-scale transcriptomic profiling has extended the study of such molecular heterogeneity between patients to somatic tissues. However, the lower resolution of bulk RNA profiling, especially in a complex, composite tissue such as the skin, has limited its success. Here we demonstrate approaches to interrogate patient-level molecular variance in a chronic skin inflammatory disease, psoriasis vulgaris, leveraging single-cell RNA-sequencing of CD45+ cells isolated from active lesions. Highly psoriasis-specific transcriptional abnormalities display greater than average inter-individual variance, nominating them as potential sources of clinical heterogeneity. We find that one of these chemokines, CXCL13, demonstrates significant correlation with severity of lesions within our patient series. Our analyses also establish that genes elevated in psoriatic skin-resident memory T cells are enriched for programs orchestrating chromatin and CDC42-dependent cytoskeleton remodeling, specific components of which are distinctly correlated with and against Th17 identity on a single-cell level. Collectively, these analyses describe systematic means to dissect cell type- and patient-level differences in cutaneous psoriasis using high-resolution transcriptional profiles of human inflammatory disease.

Classification of human chronic inflammatory skin disease based on single-cell immune profiling.

Inflammatory conditions represent the largest class of chronic skin disease, but the molecular dysregulation underlying many individual cases remains unclear. Single-cell RNA sequencing (scRNA-seq) has increased precision in dissecting the complex mixture of immune and stromal cell perturbations in inflammatory skin disease states. We single-cell-profiled CD45+ immune cell transcriptomes from skin samples of 31 patients (7 atopic dermatitis, 8 psoriasis vulgaris, 2 lichen planus (LP), 1 bullous pemphigoid (BP), 6 clinical/histopathologically indeterminate rashes, and 7 healthy controls). Our data revealed active proliferative expansion of the Treg and Trm components and universal T cell exhaustion in human rashes, with a relative attenuation of antigen-presenting cells. Skin-resident memory T cells showed the greatest transcriptional dysregulation in both atopic dermatitis and psoriasis, whereas atopic dermatitis also demonstrated recurrent abnormalities in ILC and CD8+ cytotoxic lymphocytes. Transcript signatures differentiating these rash types included genes previously implicated in T helper cell (TH2)/TH17 diatheses, segregated in unbiased functional networks, and accurately identified disease class in untrained validation data sets. These gene signatures were able to classify clinicopathologically ambiguous rashes with diagnoses consistent with therapeutic response. Thus, we have defined major classes of human inflammatory skin disease at the molecular level and described a quantitative method to classify indeterminate instances of pathologic inflammation. To make this approach accessible to the scientific community, we created a proof-of-principle web interface (RashX), where scientists and clinicians can visualize their patient-level rash scRNA-seq-derived data in the context of our TH2/TH17 transcriptional framework.

A single-cell transcriptional gradient in human cutaneous memory T cells restricts Th17/Tc17 identity.

The homeostatic mechanisms that fail to restrain chronic tissue inflammation in diseases, such as psoriasis vulgaris, remain incompletely understood. We profiled transcriptomes and epitopes of single psoriatic and normal skin-resident T cells, revealing a gradated transcriptional program of coordinately regulated inflammation-suppressive genes. This program, which is sharply suppressed in lesional skin, strikingly restricts Th17/Tc17 cytokine and other inflammatory mediators on the single-cell level. CRISPR-based deactivation of two core components of this inflammation-suppressive program, ZFP36L2 and ZFP36, replicates the interleukin-17A (IL-17A), granulocyte macrophage-colony-stimulating factor (GM-CSF), and interferon gamma (IFNγ) elevation in psoriatic memory T cells deficient in these transcripts, functionally validating their influence. Combinatoric expression analysis indicates the suppression of specific inflammatory mediators by individual program members. Finally, we find that therapeutic IL-23 blockade reduces Th17/Tc17 cell frequency in lesional skin but fails to normalize this inflammatory-suppressive program, suggesting how treated lesions may be primed for recurrence after withdrawal of treatment.

Cobolt: integrative analysis of multimodal single-cell sequencing data.

A growing number of single-cell sequencing platforms enable joint profiling of multiple omics from the same cells. We present Cobolt, a novel method that not only allows for analyzing the data from joint-modality platforms, but provides a coherent framework for the integration of multiple datasets measured on different modalities. We demonstrate its performance on multi-modality data of gene expression and chromatin accessibility and illustrate the integration abilities of Cobolt by jointly analyzing this multi-modality data with single-cell RNA-seq and ATAC-seq datasets.

Cell Wall Compositions of Sorghum bicolor Leaves and Roots Remain Relatively Constant Under Drought Conditions.

Renewable fuels are needed to replace fossil fuels in the immediate future. Lignocellulosic bioenergy crops provide a renewable alternative that sequesters atmospheric carbon. To prevent displacement of food crops, it would be advantageous to grow biofuel crops on marginal lands. These lands will likely face more frequent and extreme drought conditions than conventional agricultural land, so it is crucial to see how proposed bioenergy crops fare under these conditions and how that may affect lignocellulosic biomass composition and saccharification properties. We found that while drought impacts the plant cell wall of Sorghum bicolor differently according to tissue and timing of drought induction, drought-induced cell wall compositional modifications are relatively minor and produce no negative effect on biomass conversion. This contrasts with the cell wall-related transcriptome, which had a varied range of highly variable genes (HVGs) within four cell wall-related GO categories, depending on the tissues surveyed and time of drought induction. Further, many HVGs had expression changes in which putative impacts were not seen in the physical cell wall or which were in opposition to their putative impacts. Interestingly, most pre-flowering drought-induced cell wall changes occurred in the leaf, with matrix and lignin compositional changes that did not persist after recovery from drought. Most measurable physical post-flowering cell wall changes occurred in the root, affecting mainly polysaccharide composition and cross-linking. This study couples transcriptomics to cell wall chemical analyses of a C4 grass experiencing progressive and differing drought stresses in the field. As such, we can analyze the cell wall-specific response to agriculturally relevant drought stresses on the transcriptomic level and see whether those changes translate to compositional or biomass conversion differences. Our results bolster the conclusion that drought stress does not substantially affect the cell wall composition of specific aerial and subterranean biomass nor impede enzymatic hydrolysis of leaf biomass, a positive result for biorefinery processes. Coupled with previously reported results on the root microbiome and rhizosphere and whole transcriptome analyses of this study, we can formulate and test hypotheses on individual gene candidates' function in mediating drought stress in the grass cell wall, as demonstrated in sorghum.