Deconstructing Retinal Organoids: Single Cell RNA‐Seq Reveals the Cellular Components of Human Pluripotent Stem Cell‐Derived Retina

Abstract The rapid improvements in single cell sequencing technologies and analyses afford greater scope for dissecting organoid cultures composed of multiple cell types and create an opportunity to interrogate these models to understand tissue biology, cellular behavior and interactions. To this end, retinal organoids generated from human embryonic stem cells (hESCs) were analyzed by single cell RNA‐sequencing (scRNA‐Seq) at three time points of differentiation. Combinatorial data from all time points revealed the presence of nine clusters, five of which corresponded to key retinal cell types: retinal pigment epithelium (RPE), retinal ganglion cells (RGCs), cone and rod photoreceptors, and Müller glia. The remaining four clusters expressed genes typical of mitotic cells, extracellular matrix components and those involved in homeostasis. The cell clustering analysis revealed the decreasing presence of mitotic cells and RGCs, formation of a distinct RPE cluster, the emergence of cone and rod photoreceptors from photoreceptor precursors, and an increasing number of Müller glia cells over time. Pseudo‐time analysis resembled the order of cell birth during retinal development, with the mitotic cluster commencing the trajectory and the large majority of Müller glia completing the time line. Together, these data demonstrate the feasibility and potential of scRNA‐Seq to dissect the inherent complexity of retinal organoids and the orderly birth of key retinal cell types. Stem Cells 2019;37:593–598


INTRODUCTION
It has been estimated that 285 million people are affected by visual impairment globally, with retinal diseases accounting for approximately 26% of blindness [1]. Many inherited and agerelated retinal dystrophies culminate in the loss of photoreceptors [2,3]. There are currently no treatments to reverse this degeneration, thus cell replacement has become a prerequisite on the path toward therapeutic transplantations. The derivation of human embryonic stem cells (hESCs) in 1998 [4] and induced pluripotent stem cells (hiPSCs) in 2007 [5] has provided the much needed breakthrough as both cell types can be expanded indefinitely in vitro as well as being able to generate photoreceptors and retinal pigment epithelium (RPE) [2,3]. The subsequent development of protocols to derive three-dimensional retinal organoids from hESCs/ hiPSCs demonstrated that structures akin to the developing eye and laminated retina arise and many retinal cell types are produced [6,7]. Despite these rapid advancements, methods to accurately define and characterize the cells that arise within these organoids over time have not been fully realized.
The improvements in next generation sequencing technologies and protocols for the application of these at a single cell level have broadened their application to multiple systems [8,9]. A few pioneering studies have applied single cell RNA sequencing (scRNA-Seq) in human retinal tissue and organoids [10][11][12]. These have largely focused on enriched specific cell types through the use of reporters or bait genes. Despite these important studies, we postulated that the full complexity of retinal organoids is yet to be fully resolved by scRNA-Seq. In this study we used high throughput Integrated Fluidic Circuits (IFC) for scRNA-Seq with a capacity of 800 capture sites to interrogate hESC-derived retinal organoids through a

MATERIALS AND METHODS
A detailed description of all experimental procedures is presented in the Supporting Information.

RESULTS AND DISCUSSION
A hESC (H9) cell line was differentiated to retinal organoids. Samples were collected at 60, 90, and 200 days, dissociated, partitioned into single cells using the Fluidigm C1 Single-Cell mRNA-Seq HT IFC and processed for scRNA-Seq. Following quality control and filtering (Supporting Information Fig. S1), data from scRNA-Seq of each time point were normalized (Supporting Information Fig. S2) and then merged using the Seurat package to allow analysis of a higher cell number  Table S1 and Fig. S3). The findMarkers function found the marker genes for each cluster [13]. The top 10 marker genes were used to manually identify each cluster. Five of the clusters could be assigned to a certain cell type, with cluster 0 identified with genes commonly expressed in Müller glia, cluster 1 to cone photoreceptors, cluster 2 to retinal ganglion cells (RGCs), cluster 3 to rod photoreceptors and cluster 4 to RPE cells (Fig. 1). Cluster 5 expressed a number of extracellular matrix (ECM) genes shown to be expressed in retinal cells [14,15], whereas cluster 6 expressed genes associated with mitosis, indicating a progenitor population. The remaining minor clusters 7 and 8 were comprised largely of ECM components, signaling molecules and metabolites, the majority of which have been proposed to play a role in retinal homeostasis (for example [16][17][18][19]).
Furthermore, clustering analyses were then performed on each time point (Fig. 2). This was compared with the combined data by superimposing the original clusters (identified in Fig. 1 and Supporting Information Fig. S3) to those identified at each individual time point (Fig. 2). This analysis revealed changes over the timeline of differentiation. For example, mitotic cells and RGCs decreased over time, photoreceptors resolved into cones and rods, a distinct RPE cluster formed and Müller glia cells increased toward day 200, recapitulating retinal development where Müller glia are the last retinal cell type to become postmitotic [20][21][22]. To confirm this increase in complexity the percentage of the cell types at each time point were compared. Mitotic cells decreased from 7.7% at day 60 to 2.4% at day 200, RGCs decreased from 21% at day 60 to 6.6% at 200 and Müller glia increased from 43.7% at day 60 to 53% at day 200. This was corroborated by immunohistochemical analysis, which indicated the decreased presence of mitotic cells and RGCs from day 60 to day 200, stepwise photoreceptor differentiation (shown by the presence of Recoverin photoreceptor precursors, NRL rod precursors and RxRy cone precursors at days 60 and 90 and separation of rod and cone photoreceptors at day 200) and increased Müller glia presence from day 60 to 200, shown by Vimentin immunostaining from day 60 to 200 and the expression of CRALBP at day 200 (Fig. 3A, 3C). It is interesting to note that although horizontal (PROX1), rod bipolar (PKC-α), and amacrine (AP2α)-like cells were detected by immunohistochemical staining these were not depicted by scRNA-Seq analysis. Therefore, to be able to resolve such clusters it may be necessary to increase the cell numbers analyzed. Although the decreased presence of mitotic cells during the differentiation time course was expected, that of RGCs was not and could be due to either culture conditions, which are not able to maintain this cell type for prolonged periods, or programmed cell death as seen during development [20,23].
Interestingly, this analysis also depicted cell types with transcriptional profiles shared by several clusters (Fig. 2). For example, at day 60, RPE shared transcriptional similarity with cluster 7 and ECM with cluster 8. At day 90, cell clusters with transcriptional profiles shared by Müller glia and photoreceptors, corroborating findings reported in adult murine retina [24] but not transcriptionally documented in human retinal development as yet, and Müller glia and cluster 8 were identified. Over time, the shared transcriptome diminished and distinct clusters emerged, with the exception of Müller glia cells, which retained a shared transcriptional profile with clusters 7 and 8.
As the complexity of the organoids increased and distinct cell types could be resolved over time additional pseudo-time analysis was conducted. Monocle [25] was used to analyze the highly variable genes from the mitotic, RPE, RGC, Müller glia,  (Fig. 4A). This analysis identified the mitotic cluster as the proliferating population from which the rest of the cells emerged. A Müller glia subpopulation resides next along this branch. This is potentially due to the expression of genes commonly expressed in retinal progenitor cells (RPCs) within the Müller glia cluster (Supporting Information Fig. S4), corroborating published data for murine Müller glia cells [24,26]. Continuing down this branch, a few cone-like photoreceptors begin to arise, here the pseudo-time trajectory forked with an upper branch bearing the majority and the remaining Müller glia cells along with minor cluster 7, which was characterized by genes that have been shown to be involved in retinal homeostasis, in which Müller glia play an important role. The lower branch is initially populated with RGCs and then terminates with two closely related clusters containing the rod and cone photoreceptors. This pseudo-time trajectory resembles the order of retinal cell development with a progenitor population giving rise firstly to RGCs, followed by cones, then rods and finally the majority of Müller glia, which are situated on a separate branch, potentially indicating their unique role in retinal development and homeostasis, and thus their transcriptome. The RPE and ECM clusters reside outside of this trajectory. RPE are known to differentiate relatively early in development during optic cup formation [27] and thus have a more distinct transcriptome. The genes from the ECM cluster are not exclusive to the retina, being common ECM components expressed throughout many tissue types and thus by ontology are unlikely to be associated with the retina. This is likely to explain why these two clusters are not associated with the time line (Fig. 4A). The cell ordering plot (Fig. 4B) further corroborates the order of cell type emergence within the organoids with mitotic cells, some Müller glia and cluster 8 present initially at the first time point. Cone photoreceptors and RGCs are next, arising at day 60, followed with diminishing RGC birth and the majority of cones appearing on day 90, alongside the emergence of rods and more Müller glia. A large number of Müller glia cells (along with some rods) were present at day 200 as the later cell types to mature.

CONCLUSION
Our data demonstrate the feasibility and potential of scRNA-Seq to dissect the inherent complexity of retinal organoids and the orderly birth of key retinal cell types therein, which recapitulates the order of retinal development.