Abstract The use of single-cell transcriptomics has become a major approach to delineate cell subpopulations and the transitions between them. While various computational tools using different mathematical methods have been developed to infer clusters, marker genes, and cell lineage, none yet integrate these within a mathematical framework to perform multiple tasks coherently. Such coherence is critical for the inference of cell-cell communication, a major remaining challenge. Here we present similarity matrix-based optimization for single-cell data analysis (SoptSC), in which unsupervised clustering, pseudotemporal ordering, lineage inference, and marker gene identification are inferred via a structured cell-to-cell similarity matrix. SoptSC then predicts cell-cell communication networks, enabling reconstruction of complex cell lineages that include feedback or feedforward interactions. Application of SoptSC to early embryonic development, epidermal regeneration, and hematopoiesis demonstrates robust identification of subpopulations, lineage relationships, and pseudotime, and prediction of pathway-specific cell communication patterns regulating processes of development and differentiation.

ABSTRACT How stem cells give rise to human interfollicular epidermis is unclear despite the crucial role the epidermis plays in barrier and appendage formation. Here we use single cell-RNA sequencing to interrogate basal stem cell heterogeneity of human interfollicular epidermis and find at least four spatially distinct stem cell populations that decorate the top and bottom of rete ridge architecture and hold transitional positions between the basal and suprabasal epidermal layers. Cell-cell communication modeling through co-variance of cognate ligand-receptor pairs indicate that the basal cell populations distinctly serve as critical signaling hubs that maintain epidermal communication. Combining pseudotime, RNA velocity, and cellular entropy analyses point to a hierarchical differentiation lineage supporting multi-stem cell interfollicular epidermal homeostasis models and suggest the “transitional” basal stem cells are stable states essential for proper stratification. Finally, alterations in differentially expressed “transitional” basal stem cell genes result in severe thinning of human skin equivalents, validating their essential role in epidermal homeostasis and reinforcing the critical nature of basal stem cell heterogeneity.

Abstract Methods to model dynamic changes in gene expression at a genome-wide level are not currently sufficient for large (temporally rich or single-cell) datasets. Variational autoencoders offer means to characterize large datasets and have been used effectively to characterize features of single-cell datasets. Here we extend these methods for use with gene expression time series data. We present RVAgene: a recurrent variational autoencoder to model gene expression dynamics. RVAgene learns to accurately and efficiently reconstruct temporal gene profiles. It also learns a low dimensional representation of the data via a recurrent encoder network that can be used for biological feature discovery, and can generate new gene expression data by sampling from the latent space. We test RVAgene on simulated and real biological datasets, including embryonic stem cell differentiation and kidney injury response dynamics. In all cases, RVAgene accurately reconstructed complex gene expression temporal profiles. Via cross validation, we show that a low-error latent space representation can be learnt using only a fraction of the data. Through clustering and gene ontology term enrichment analysis on the latent space, we demonstrate the potential of RVAgene for unsupervised discovery. In particular, RVAgene identifies new programs of shared gene regulation of Lox family genes in response to kidney injury.

Clonal hematopoiesis becomes increasingly common with age, but its cause is enigmatic because driver mutations are often absent. Serial observations infer weak selection indicating variants are acquired much earlier in life with unexplained initial growth spurts. Here we use fluctuating CpG methylation as a lineage marker to track stem cell clonal dynamics of hematopoiesis. We show, via the shared prenatal circulation of monozygotic twins, that weak selection conferred by stem cell variation created before birth can reliably yield clonal hematopoiesis later in life. Theory indicates weak selection will lead to dominance given enough time and large enough population sizes. Human hematopoiesis satisfies both these conditions. Stochastic loss of weakly selected variants is naturally prevented by the expansion of stem cell lineages during development. The dominance of stem cell clones created before birth is supported by blood fluctuating CpG methylation patterns that exhibit low correlation between unrelated individuals but are highly correlated between many elderly monozygotic twins. Therefore, clonal hematopoiesis driven by weak selection in later life appears to reflect variation created before birth.

The monocytic/macrophage system is essential for skeletal muscle homeostasis, but its dysregulation contributes to the pathogenesis of muscle degenerative disorders. Despite our increasing knowledge of the role of macrophages in degenerative disease, it still remains unclear how macrophages contribute to muscle fibrosis. Here, we used single-cell transcriptomics to determine the molecular attributes of dystrophic and healthy muscle macrophages. We identified six novel clusters. Unexpectedly, none corresponded to traditional definitions of M1 or M2 macrophage activation. Rather, the predominant macrophage signature in dystrophic muscle was characterized by high expression of fibrotic factors, galectin-3 and spp1. Spatial transcriptomics and computational inferences of intercellular communication indicated that spp1 regulates stromal progenitor and macrophage interactions during muscular dystrophy. Galectin-3 + macrophages were chronically activated in dystrophic muscle and adoptive transfer assays showed that the galectin-3 + phenotype was the dominant molecular program induced within the dystrophic milieu. Histological examination of human muscle biopsies revealed that galectin-3 + macrophages were also elevated in multiple myopathies. These studies advance our understanding of macrophages in muscular dystrophy by defining the transcriptional programs induced in muscle macrophages, and reveal spp1 as a major regulator of macrophage and stromal progenitor interactions.

Abstract Single-cell genomic technologies offer vast new resources with which to study cells, but their potential to inform parameter inference of cell dynamics has yet to be fully realized. Here we develop methods for Bayesian parameter inference with data that jointly measure gene expression and Ca 2+ dynamics in single cells. We propose to share information between cells via transfer learning: for a sequence of cells, the posterior distribution of one cell is used to inform the prior distribution of the next. In application to intracellular Ca 2+ signaling dynamics, we fit the parameters of a dynamical model for thousands of cells with variable single-cell responses. We show that transfer learning accelerates inference with sequences of cells regardless of how the cells are ordered. However, only by ordering cells based on their transcriptional similarity can we distinguish Ca 2+ dynamic profiles and associated marker genes from the posterior distributions. Inference results reveal complex and competing sources of cell heterogeneity: parameter covariation can diverge between the intracellular and intercellular contexts. Overall, we discuss the extent to which single-cell parameter inference informed by transcriptional similarity can quantify relationships between gene expression states and signaling dynamics in single cells.

Inference of gene regulatory networks (GRNs) can reveal cell state transitions from single-cell genomics data. However, obstacles to temporal inference from snapshot data are difficult to overcome. Single-nuclei multiomics data offer means to bridge this gap and derive temporal information from snapshot data using joint measurements of gene expression and chromatin accessibility in the same single cells. We developed popInfer to infer networks that characterize lineage-specific dynamic cell state transitions from joint gene expression and chromatin accessibility data. Benchmarking against alternative methods for GRN inference, we showed that popInfer achieves higher accuracy in the GRNs inferred. popInfer was applied to study single-cell multiomics data characterizing hematopoietic stem cells (HSCs) and the transition from HSC to a multipotent progenitor cell state during murine hematopoiesis across age and dietary conditions. From networks predicted by popInfer, we discovered gene interactions controlling entry to/exit from HSC quiescence that are perturbed in response to diet or aging.

Abstract Myeloid-derived suppressor cells (MDSCs) play a prominent and rising role in the tumor microenvironment. An understanding of the tumor-MDSC interactions that influence disease progression is critical, and currently lacking. To address this, we developed a mathematical model of metastatic growth and progression in immune-rich tumor microenvironments. We model the tumor-immune dynamics with stochastic delay differential equations, and study the impact of delays in MDSC activation/recruitment on tumor growth outcomes. We find when the circulating level of MDSCs is low, the MDSC delay has a pronounced impact on the probability of new metastatic establishment: blocking MDSC recruitment can reduce the probability of metastasis by as much as 50%. We also quantify the extent to which decreasing the immuno-suppressive capability of the MDSCs impacts the probability that a new metastasis will persist or grow. In order to quantify patient-specific MDSC dynamics under different conditions we fit individual tumors treated with immune checkpoint inhibitors to the tumor-MDSC model via Bayesian parameter inference. We reveal that control of the inhibition rate of natural killer cells by MDSCs has a larger influence on tumor outcomes than controlling the tumor growth rate directly. Posterior classification of tumor outcomes demonstrates that incorporating knowledge of the MDSC responses improves predictive accuracy from 63% to 82%. Our results illustrate the importance of MDSC dynamics in the tumor microenvironment and predict interventions that may shift environments towards a less immune-suppressed state. We argue that there is a pressing need to more often consider MDSCs in analyses of tumor microenvironments.

Mammalian organs exhibit distinct physiology, disease susceptibility and injury responses between the sexes. In the mouse kidney, sexually dimorphic gene activity maps predominantly to proximal tubule (PT) segments. Bulk RNA-seq data demonstrated sex differences were established from 4 and 8 weeks after birth under gonadal control. Hormone injection studies and genetic removal of androgen and estrogen receptors demonstrated androgen receptor (AR) mediated regulation of gene activity in PT cells as the regulatory mechanism. Interestingly, caloric restriction feminizes the male kidney. Single-nuclear multiomic analysis identified putative cis-regulatory regions and cooperating factors mediating PT responses to AR activity in the mouse kidney. In the human kidney, a limited set of genes showed conserved sex-linked regulation while analysis of the mouse liver underscored organ-specific differences in the regulation of sexually dimorphic gene expression. These findings raise interesting questions on the evolution, physiological significance, and disease and metabolic linkage, of sexually dimorphic gene activity.

Abstract In animal models, Nipbl -deficiency phenocopies gene expression changes and birth defects seen in Cornelia de Lange Syndrome (CdLS), the most common cause of which is Nipbl -haploinsufficiency. Previous studies in Nipbl +/- mice identified aberrant gene expression and heart defects as early as cardiac crescent (CC) stage. Here, we performed single-cell RNA-sequencing on wildtype (WT) and Nipbl +/- mouse embryos at CC- and earlier (gastrulation) stages. Nipbl +/- embryos had fewer mesoderm cells than WT and altered proportions of mesodermal cell subpopulations. These findings were associated with an underexpression of genes implicated in driving specific mesodermal lineages. Nipbl +/- embryos also misexpressed developmentally-critical genes, including the transcription factor, Nanog , and genes governing left-right and anterior-posterior patterning. These events of cell misallocation and transcriptional dysregulation foreshadowed defects in tissue composition and patterning that arise later in Nipbl +/- mice, offering insights into early developmental contributions to birth defects in CdLS. Teaser Gene expression changes during gastrulation of Nipbl -deficient mice shed light on early origins of structural birth defects.