Author response: Human embryoid bodies as a novel system for genomic studies of functionally diverse cell types
Year: 2021
Type: peer-review
Abstract: Article Figures and data Abstract Editor's evaluation eLife digest Introduction Results Discussion Materials and methods Data availability References Decision letter Author response Article and author information Metrics Abstract Practically all studies of gene expression in humans to date have been performed in a relatively small number of adult tissues. Gene regulation is highly dynamic and context-dependent. In order to better understand the connection between gene regulation and complex phenotypes, including disease, we need to be able to study gene expression in more cell types, tissues, and states that are relevant to human phenotypes. In particular, we need to characterize gene expression in early development cell types, as mutations that affect developmental processes may be of particular relevance to complex traits. To address this challenge, we propose to use embryoid bodies (EBs), which are organoids that contain a multitude of cell types in dynamic states. EBs provide a system in which one can study dynamic regulatory processes at an unprecedentedly high resolution. To explore the utility of EBs, we systematically explored cellular and gene expression heterogeneity in EBs from multiple individuals. We characterized the various cell types that arise from EBs, the extent to which they recapitulate gene expression in vivo, and the relative contribution of technical and biological factors to variability in gene expression, cell composition, and differentiation efficiency. Our results highlight the utility of EBs as a new model system for mapping dynamic inter-individual regulatory differences in a large variety of cell types. Editor's evaluation The authors generated embryoid bodies (EBs) from induced pluripotent stem cells (iPSCs) using a strong mixed-pool study design and performed scRNA-seq profiling. From this data, they identify dozens of cell types and infer differentiation trajectories that align well with known developmental gene expression dynamics. This system is likely to be a good platform for larger eQTL studies that interrogate new cell states. https://doi.org/10.7554/eLife.71361.sa0 Decision letter eLife's review process eLife digest One major goal of human genetics is to understand how changes in the way genes are regulated affect human traits, including disease susceptibility. To date, most studies of gene regulation have been performed in adult tissues, such as liver or kidney tissue, that were collected at a single time point. Yet, gene regulation is highly dynamic and context-dependent, meaning that it is important to gather data from a greater variety of cell types at different stages of their development. Additionally, observing which genes switch on and off in response to external treatments can shed light on how genetic variation can drive errors in gene regulation and cause diseases. Stem cells can produce more cells like themselves or differentiate – acquire the characteristics – of many cell types. These cells have been used in the laboratory to research gene regulation. Unfortunately, these studies often fail to capture the complex spatial and temporal dynamics of stem cell differentiation; in particular, these studies are unable to observe gene regulation in the transient cell types that appear early in embryonic development. To overcome these limitations, scientists developed systems such as embryoid bodies: three-dimensional aggregates of stem cells that, when grown under certain conditions, spontaneously develop into a variety of cell types. Rhodes, Barr et al. wanted to assess the utility of embryoid bodies as a model to study how genes are dynamically regulated in different cell types, by different individuals who have distinct genetic makeups. To do this, they grew embryoid bodies made from human stem cells from different individuals to examine which genes switched on and off as the stem cells that formed the embryoid bodies differentiated into different types of cells. The results showed that it was possible to grow embryoid bodies derived from genetically distinct individuals that consistently produce diverse cell types, similar to those found during human fetal development. Rhodes, Barr et al.’s findings suggest that embryoid bodies are a useful model to study gene regulation across individuals with different genetic backgrounds. This could accelerate research into how genetics are associated with disease by capturing gene regulatory dynamics at an unprecedentedly high spatial and temporal resolution. Additionally, embryoid bodies could be used to explore how exposure to different environmental factors during early development affect disease-related outcomes in adulthood in different individuals. Introduction Genome-wide association studies (GWAS) have identified thousands of genetic variants associated with human traits and diseases, many of which are located in noncoding regions of the genome and are putatively regulatory in function (Albert and Kruglyak, 2015). To understand regulatory and functional effects of trait-associated variants, it is necessary to perform molecular assays in the relevant cell types at the relevant stages of development, and potentially to also model different environmental exposures (Umans et al., 2020). However, most efforts to identify genetic variants that regulate gene expression (expression quantitative trait loci, or eQTLs) have relied on adult tissue samples collected at a single time point. While such efforts have mapped millions of static, steady-state eQTLs across dozens of tissues and cell types (GTEx Consortium, 2020), most disease-associated variants were not found to also be classified as eQTLs (Yao et al., 2020; Aguet et al., 2017). It is possible that dynamic and variable regulatory genetic effects, including those that are specific to a given cell type, time point, or environment, may underlie the mechanisms for many unexplained phenotypic associations. For example, recent efforts to characterize gene regulatory dynamics in human induced pluripotent stem cells (iPSCs) and their derived cell types have identified dynamic eQTLs that are associated with disease risk, supporting the intuitive notion that changes in gene regulation during development may play a role in shaping human adult phenotypes, including disease (Strober et al., 2019; Cuomo et al., 2020). Still, iPSCs are limited in their potential for identifying dynamic regulatory effects. The number of cell types that can be obtained from iPSCs using directed differentiation protocols is quite modest, and time-course experiments, although useful for studying gene regulation at discrete points along a continuum, are inefficient, expensive, and laborious to perform. With these challenges in mind, we wanted to develop and characterize a new in vitro model for studying gene regulation – a model capable of measuring gene expression continuously along the developmental trajectories towards multiple cell types, and to be able to do so in multiple individuals. To do this, we used iPSCs to form embryoid bodies (EBs), which are three-dimensional aggregates of spontaneously and asynchronously differentiating cells. EB formation has been used to verify stem cell pluripotency for decades; yet, until recently, the complexity of EB cellular composition has precluded their use in genomic studies. With single-cell RNA-sequencing (scRNA-seq), it is now possible to characterize the numerous spatially and developmentally distinct cell types within EBs, including transient cell types that would otherwise be inaccessible. Indeed, recent scRNA-seq studies of human EB differentiation have revealed the diversity of cell types composing these structures and the transcriptional dynamics governing early fate decisions (Han et al., 2018; Guo et al., 2019). To date, however, the only studies that have sequenced EB cells have relied on a small sample of cells from a single individual, leaving a gap in our understanding of technical, biological, and inter-individual variation present in this system (Han et al., 2018; Guo et al., 2019). Understanding the sources of variation that affect cell composition and scRNA-seq data from EBs is crucial for evaluating the utility of EBs as a novel system for population-level studies of gene regulation. To this end, we used a batch-controlled study design to generate and sequence EBs from multiple individuals using multiple replicates. This allowed us to measure the degree of technical and biological variability in cell identity and gene expression levels associated with repeated independent EB differentiations. We evaluated the consistency in cell type composition across replicates and individuals, characterized the structure of variation in gene expression across the entire data set, and finally, captured patterns of dynamic gene expression along distinct developmental trajectories. Our results indicate that scRNA-seq of differentiating EBs has the potential to be a powerful model system for the study of inter-individual variation in gene regulation across an array of functionally and temporally diverse cell types. Results We have performed a pilot study to establish and characterize the EB system. Toward the ultimate goal of performing dynamic eQTL studies using EBs, we designed a study that allowed us to effectively estimate different sources of variation in single-cell data from EBs. In this pilot study, we focused on consistency of the non-directed differentiation process, and the proportion of gene expression variability that can be explained by technical or biological factors. Study Design, data collection, and preprocessing To characterize sources of variation in gene expression in human EBs, we initially differentiated EBs from three human iPSC lines (18511, 18858, and 19160) in three replicates (see Materials and methods). We performed the experiment in three batches, where each batch includes one replicate from each of the three individuals. EBs differentiate quickly, with cell types representing endoderm, mesoderm, and ectoderm present after 8 days (Han et al., 2018). In this study, we maintained EBs for 3 weeks after formation, allowing cells to continue to differentiate and mature. After 21 days, we collected scRNA-seq data, targeting equal numbers of cells from each individual in each replicate. After filtering and quality control (Materials and methods), we retained high-quality data from a sample of 42,488 cells (an average of 4721 cells per individual/replicate). For these cells, we obtained a median of 16,712 UMI counts per cell, which allowed us to measure the expression of a median of 4274 genes per cell (Figure 1—figure supplement 1). We integrated data from all cells using Harmony, which anchors the data sets by cell type (Korsunsky et al., 2019). After these initial collections, we found that one individual, 18858, had lower differentiation efficiency than the other two lines (Figure 1G, see Differentiation efficiency across individuals). Too assess the robustness of differentiation efficiency and cell type composition among a larger sample of individuals from our YRI iPSC panel, we differentiated EBs in one replicate from each of five additional randomly chosen lines (18856, 18912, 19140, 19159, and 19210). After filtering and quality control, we retained an average of 5243 cells per individual in this new set of data, a median of 5983 UMI counts per cell, and a median of 2775 genes per cell (Figure 1—figure supplement 4). Throughout the text, when we use data from the five lines that were subsequently collected using only a single replicate (18856, 18912, 19140, 19150, and 19210), we refer to them as the ‘additional lines’. The initial replicated data collected from individuals 18511, 18858, and 19160 are used in all analyses throughout this study, while the additional lines are used only to demonstrate consistency of cell type composition across individuals. Figure 1 with 4 supplements see all Download asset Open asset Characterization of EB cell type composition using marker gene expression and clustering. (A–F) Visualization of EB cells with UMAP. (A) Cells from lines 18511, 18858, and 19160 colored by expression of pluripotent marker gene POU5F1, (B) Cells from lines 18511, 18858, and 19160 colored by expression of endoderm marker gene SOX17, (C) Cells from lines 18511, 18858, and 19160 colored by expression of mesoderm marker gene HAND1, (D) Cells from lines 18511, 18858, and 19160 colored by expression of early ectoderm marker gene PAX6. In A-D cells are colored by normalized counts. (E) Cells from lines 18511, 18858, and 19160 colored by Seurat cluster assignment at clustering resolution 0.1. (F) Cells from lines 18511, 18858, and 19160 colored by Seurat cluster assignment at clustering resolution 1. (G) Proportions of cells from replicates of lines 18511, 18858, and 19160 assigned to Seurat clusters at clustering resolution 0.1. (H) Proportions of cells from additional lines assigned to broad cell types present in EBs. Cell type composition To validate our expectation that EBs should contain cells from each germ layer, we first characterized the expression of early developmental marker genes. We found cells expressing markers for endoderm (SOX17, FOXA2), mesoderm (HAND1), and ectoderm (PAX6), in addition to cells still expressing pluripotency markers (POU5F1, MYC, NANOG). We visualized the data with uniform manifold approximation and projection (UMAP) and observed that cells expressing each of these germ layer markers occupied distinct groups in UMAP space (Figure 1A–D; Becht et al., 2018). Moreover, we found that every replicate in our experiment, regardless of the individual, includes cells from all three germ layers (Figure 1E and G). We next sought to further explore the heterogeneous cell types present in these EBs. In studies of scRNA-seq data from tissues and samples with well-characterized cell type composition, clustering is often applied to demarcate populations of pure cell types within heterogeneous samples. In these studies, clustering resolution, which determines the number of clusters identified by the algorithm, is typically chosen to recapitulate the expected number of known cell types. The identified clusters can be annotated based on the expression of known marker genes. In our case, however, we had no a priori knowledge of the exact number or types of cells that would result from the spontaneous differentiation of the EBs. Hence, we used three complementary approaches to annotate cells, capturing various perspectives on what might define a cell type in this data set. First, we identified cell types by clustering cells and annotating the cell types based on the genes that are highly expressed in each cluster. Second, we annotated cell types by considering the correlation of gene expression in our data with a reference data set of known primary cell types. For our third approach, we used a different perspective, and applied topic modeling to consider a less discrete definition of cell type. For the first approach, we used a standard clustering analysis, the Louvain algorithm in Seurat, to identify groups of cells with similar transcriptomes (Blondel et al., 2008). To avoid making assumptions about the true number of cell types present, we repeated this analysis across different clustering resolutions (resolution 0.1, 0.5, 0.8, and 1). As expected, the number of clusters we identified varied greatly depending on the resolution (Figure 1E–F, Figure 1—figure supplement 2). We performed each subsequent analysis using clusters defined at multiple resolutions, to ensure that our qualitative conclusions are robust with respect to the number of clusters identified. For each clustering resolution, we calculated pseudobulk gene expression levels using cells from the same cluster, individual, and replicate. To identify marker genes expressed in each cluster, we used Limma and voom to perform differential expression analysis (Materials and methods) using the pseudobulk estimates. For example, considering the gene expression data of the seven clusters identified at resolution 0.1, we found that the most significantly upregulated genes in each cluster included known marker genes for pluripotent cells (cluster 0), early ectoderm (cluster 1), mesoderm (cluster 2), neural crest (cluster 3), endoderm (cluster 4), neurons (cluster 5), and endothelial cells (cluster 6) (Figure 1E, Table 1). Using this approach provides a confident set of broad cell type categories present in these data. At higher resolutions, DE analysis between clusters enabled annotation of some more specific cell types; for example, at clustering resolution 1, cluster 19 is characterized by higher expression of hepatocyte marker genes FGB, TTR, and AFP. More generally, however, confident cell type classification of Seurat clusters at higher resolution based on DE alone proved difficult (Supplementary files 4-7). Table 1 Classification of Seurat cluster identity (clustering resolution 0.1) based on differential expression of marker genes. Cluster (Res. 0.1)# cells in clusterTop 10 marker genes by adj. PTop 10 marker genes by logFCAnnotation017,693TERF1, PHC1, SEPHS1, UGP2, DPPA4, TBC1D23, JARID2, USO1, ZNF398, LRRC47DPPA5, DPPA3, GDF3, NANOG, FGF4, POU5F1, CBR3, PRDM14, DPPA2, TRIML2Pluripotent Cells114,383TPBG, FGFBP3, FZD3, LIX1, SDK2, BTBD17, DACH1, PLAGL1, DEK, ZNF219FEZF2, EMX2, LHX2, SOX3, PAX6, WNT7, BARX, SOX1, ZIC1, SIX3Early Ectoderm23086TNNI1, COL6A3, COL5A1, RGS4, ACTA2, TMEM88, DOK4, SLC40A1, HAND2, COL3A1RGS13, LUM, TECRL, DCN, HAND1, PITX1, COL3A1, SLN, IGF2, FIBINMesoderm32673NR2F1, CNP, S100B, EDNRA, FGFBP3, ATP1A2, DNAJC1, ZEB2, PHACTR3, METRNMPZ, PRSS56, ROPN1, SOX10, S100B, SCRG1, NPR3, MOXD1, TFAP2B, PHACTR3Neural Crest42368S100A16, LGALS3, GATA3, CST3, KRT19, FN1, EPSTI1, DYNLT3, HDHD3, PKP2APOA2, CST1, APOA1, APOC3, FGB, RBP4, S100A14, TTR, FGA, APOBEndoderm51990TAGLN3, RTN1, NHLH1, STMN2, ELAVL2, FNDC5, PCBP4, ELAVL4, DCX, MLLT11NEUROD1, NHLH1, STMN2, NEUROD4, TBR1, STMN4, NEUROG1, SST, ELAVL3, SLC17A6Neurons6295EGFL7, GNG11, RAMP2, IGFBP4, PPM1F, RASGRP3, RCSD1, MAP4K2, PLVAP, DOCK6PLVAP, CD34, CD93, CDH5, DIPK2B, PECAM1, EMCN, CRHBP, ESAM, ECSCREndothelial Cells To pursue the second approach, we annotated cells by comparing our gene expression data to available reference sets of scRNA-seq data from primary fetal tissues, human embryonic stem cells (hESCs), and hESC-derived EBs (Cao et al., 2020; Han et al., 2020). To do this, we first integrated our data set with the reference data sets and visualized cells with UMAP (Figure 2A, Figure 2—figure supplement 1, Materials and Methods). We observed that reference hESCs cluster closely with pluripotent EB cells. We also observed that the hESC-derived EBs and our iPSC-derived EBs tend to occupy the same areas in UMAP space, implying high overall similarity in cell type composition despite differing protocols for EB differentiation (and despite the fact that the experiments were performed in different labs). EB cells also show overlap with many primary fetal cell types (Figure 2B–C, Figure 2—figure supplement 2). For example, EB cells annotated as neural crest based on our gene expression analysis, overlap with primary fetal cell types derived from neural crest, such as Schwann cells and enteric nervous system (ENS) glia (Figure 2B). EB cells annotated as neurons based on our gene expression analysis overlap with fetal neuronal subtypes, including and EB cells also show overlap with populations of cells that are in the fetal data set, including cells cells, and cells (Figure Figure with 3 supplements see all Download asset Open asset and cell type annotation with lines 18511, 18858, and (A) UMAP of EB cells from this study and cells from reference data sets of fetal cell types, EBs, and hESCs after Cells are colored by data set. (B) UMAP of EB cells from this study and cells from the fetal reference after Cells are colored by Seurat cluster identity at clustering resolution 0.1, with points representing cells from the fetal reference set. (C) UMAP of EB cells from this study and data from the fetal reference after Cells are colored by cell types present in the fetal reference data set, with points representing EB cells. (D) UMAP of EB cells from this data set with from the fetal and reference by these we the annotation of our EB cells to this were based on the expression of known marker by using the known of the reference primary fetal cell type data set. we cell to EB cells based on the reference cells in space (Figure Using this approach, we found EB cells representing of the primary cell types present in the reference fetal data set (Supplementary 1). The most annotation was of EB this can be to the high proportion of pluripotent cells in our EB data set, also to the fact that the reference fetal data set not many early developmental cell types. Indeed, many cells annotated as are likely to differentiating cells which are no pluripotent transcriptional more closely hESCs than the more highly differentiated fetal cell types present in the reference data set. In this EB data sets may capture transient developmental cell types that are difficult or to study in fetal primary samples. the many fetal cell types are only by small populations of EB cells. For example, only one EB cell is annotated as a and only one cell is annotated as a These indicate that, in the we can from a of EB single cells in order to explore their true cell type annotation based on the reference set revealed the of dozens of diverse cell types in EBs. Differentiation efficiency across individuals To assess the differentiation efficiency of each individual in each we calculated the proportion of cells assigned to each cluster as resolution (Figure While EBs from two of the individuals in our study differentiated across all we observed that of cells from individual were assigned to cluster the cluster annotated as pluripotent cells based on differential expression of marker genes 1, The EBs from this do high quality cells assigned to clusters representing each germ layer (Figure 1G, Figure 1—figure supplement 3), these EB have overall lower differentiation efficiency than EBs from individuals and To individual is an and more estimate how often EB differentiation is less we differentiated an additional five human iPSC lines from individuals 18912, 19140, 19159, and (Materials and methods). We were to that EBs from additional lines differentiated with cell type composition similar to and 19160 (Figure These results suggest that differentiation efficiency is expected to be among the YRI iPSC We further explored the robustness of cell type composition by the additional lines with the fetal and reference data sets using the same as was used for the lines (Figure Figure supplement 1, Figure supplement 2). We that assigned to cells from these additional lines fetal cell types; this set of included cell types that were in the three also that were in the three lines (Supplementary 2). we observed that most EB cells from the additional lines are annotated as of although many no pluripotency markers and do markers of various germ layers as we these results our conclusions that EBs contain many diverse cell types, many of which likely capture stages of development than are captured in fetal data. Figure 3 with supplements see all Download asset Open asset and cell type annotation with additional (A) UMAP of EB cells from this study and cells from reference data sets of fetal cell types, EBs, and hESCs after Cells are colored by data set. (B) UMAP of EB cells from this study and cells from the fetal reference after Cells are colored by broad cell type assigned using clustering and marker gene expression, with points representing cells from the fetal reference set. (C) UMAP of EB cells from this study and data from the fetal reference after Cells are colored by cell types present in the fetal reference data set, with points representing EB cells. (D) UMAP of EB cells from this data set with from the fetal and reference modeling of the single-cell gene expression data of the approaches we and to a reference data that are discrete each cell has a single true and cell type categories are to be and This definition of a cell type is intuitive and it to consider results from single-cell analysis in the of the of knowledge from However, cells into discrete groups is to capture the degree of heterogeneity in gene expression of single cells. particular cluster or cell may multiple cell functionally distinct such as cells in different stages of the cell This more in data sets that differentiating cells, which are expected to show of similarity to a cell type. In an cell type can be as with the expression of each cell representing of to multiple categories et al., 2017). One used to capture cell identity in this is topic which major patterns in gene expression within the data set, or and each cell as a of these We applied topic modeling using at a of topic resolutions, identifying and in our data. closely to Seurat on cells of a given cluster not on For example, in the topic analysis, topic 1 is on cells assigned to Seurat cluster 4 (cluster resolution 0.1) which we annotated as endoderm Table 1). to other topic 1 an in expression of and which are known markers of (Figure Table 2). Seurat clustering at higher resolution (resolution results in further of this large endoderm of cells into endoderm and (Figure modeling revealed that these cells variable of in topic 1 topic this a temporal of Figure 4 with 3 supplements see all Download asset Open asset modeling of EB cells. (A) the results of topic modeling at includes a of EB cells by Seurat cluster at resolution 0.1. (B) UMAP projection of cells colored by of topic 1. (C) the of topic 1 from the topic analysis on each Seurat cluster at clustering resolution 0.1. (D) the of topic 1 from the topic analysis on each Seurat cluster at clustering resolution 1. (E) genes expressed between topic 1 and all other from the topic are colored by the average on the Table genes of each topic from the topic model based on FN1, CST3, KRT19, LGALS3, DPPA4, KRT19, FN1, APOA1, are across cells assigned to different Seurat clusters (Figure Figure supplement For example, topic from the topic analysis is across all Seurat to all other topic expression of many
Source: eLife
Institutions: University of Chicago, Johns Hopkins University
Language: English
Subfield: Molecular Biology
Domain: Life Sciences
Open Access status: gold
Estimated APC paid: $2,500
OpenAlex ID: W4211168531
