Photoreceptor-enriched full-length scRNA-seq of developing human retina.

A. Overview of sample collection and sequencing. B. UMAP plot and low-resolution cell type clusters. C. Expression of marker genes for RPC/MGs (LHX2), rods (NR2E3), S cones (OPN1SW), L/M cones (THRB). Insets: Gene expression violin plots (from left to right): RPC/MG (red), iPRP (brown), LM cone (green), S cone (teal), early rod (blue), late rod (pink). D. Expression of markers of rod maturation (PDE6G, RHO) and cone maturation (PDE6H, OPN1LW). Arrowheads: Late maturing RHO+ rods (top), late maturing OPN1LW+ cones (bottom). See Figure S1 for additional examples. UMAP and violin plots for any gene or Ensembl transcript isoform can be produced at https://docker.saban.chla.usc.edu/cobrinik/app/seuratApp/.

Regulon-defined RPC and photoreceptor precursor states

A. Ward-clustered heatmap of the highest scoring SCENIC regulons in each cluster, displaying Z-score normalized regulon activities. Late = late maturing L/M cones. B. Box plot of RAX regulon area under the curve (AUC) values for early and late L/M cones and rods. *, p<0.005; ***, p<0.0005, Dunn test. C,D. UMAP plots of regulon AUC values for (C) PAX6 (RPC/MG) and E2F2 (RPC), and (D) OTX2 (photoreceptors and photoreceptor-committed RPCs), NRL (rod) and THRB and ISL2 (L/M cone).

Differential expression of NRL isoforms in rod and cone precursor populations.

A. Expression of NRL gene and the most highly assigned Ensembl isoforms ENST00000397002 (FL-NRL) and ENST00000560550 (Tr-NRL). B. Mean NRL isoform assignments for clusters defined in Figure 1B, presented as total counts (top) and percentage of total counts (bottom). Significance for LM vs. ER fold change, colored by isoform. ****, p <0.0002; *****, <0.000001 (bootstrapped Welch’s t-test). C. Top: Mean read counts across Ensembl NRL exons for each cluster. Bottom: Transcript structures numbered according to amino acid positions. Minimal transactivation domain (MTD) in green. Arrowheads: Red/black: First exons where red is higher of two peaks. D. Relative difference box plot of raw reads mapping to truncated (Tr) and full length (F) transcript first exons in each cell, according to cluster. Relative difference is the difference in reads mapping to truncated and full-length NRL first exons (Tr-F) divided by the sum of both (Tr+F). Values >0 indicate more reads assigned to truncated isoform, values < 0 indicate more reads assigned to full-length isoform. ***, p <0.0001 (post-hoc Dunn test). E. NRL and RXRγ immunostaining and RNA FISH with probes specific to truncated Tr-NRL exon 1T (green puncta) and FL-NRL exons 1 and 2 (orange puncta) in FW16 retina. Boxed regions enlarged at right show a RXRlo, NRL+ rod with one Tr-NRL and six FL-NRL puncta (top) and an RXRγhi, NRL− cone with one Tr-NRL and no FL-NRL puncta (bottom), indicated with same-color arrows. Scale bar: 10 µm. F. Ratio of fluorescent puncta observed in experiment depicted in (E) for NRL+ or RXRγhi cells where Tr puncta >0. ***; p <0.0005 (Welch’s t-test). G. Top: PDE6B-luciferase reporter activity in NIH-3T3 cells transfected with indicated amounts (ng) of pcDNA4C-EF1α and derived FL-NRL and Tr-NRL. Bottom: Expression constructs. Error bars = standard deviation of triplicate measurements. *, p <0.05; **, <0.005 (Student’s t-test). Data representative of two experiments in NIH-3T3 and one in 293T. H. Long-read nanopore sequencing of pooled 5’ RACE reactions initiated with NRL exon 3 primers and performed on cDNA libraries from 23 ER cells (top) and 21 LM cells (bottom). Each schematic shows total exon coverage (above) and individual transcript structures (below), where expressed sequences are grey and introns light blue. Full-length (FL), alternatively spliced or internally initiated exon 2 (Δex2), and truncated (Tr) transcripts are indicated by brackets. Red arrow: Transcripts resembling DD10, with internal exon 2 transcription initiation and premature splicing to exon 3. Ensembl FL-NRL and Tr-NRL transcript isoforms and RACE primer positions are shown below.

Differential expression of THRB isoforms in rod and cone precursors.

A. Expression of THRB and highly assigned isoforms ENST00000280696 (encoding TRβ2) and ENST00000396671 (TRβ1). B. Mean THRB isoform assignments for each cluster presented as counts (top) and percentage of counts (bottom); β1 PTT isoforms indicated. C. Top: Mean read counts across Ensembl THRB exons. Bottom: Transcript structures for TRβ1, TRβ2, and two TRβ1 truncations. Green arrowhead: First TRβ2 exon. ENST00000396671 exon numbers are indicated above and protein domains (AF1, DNA-binding (DBD), and ligand binding (LBD)) below. D. Read coverage for LR cells across THRB exons 4 and 6 splice donor sites. E. Percentage of exon splice donor reads that are spliced or readthrough to the subsequent intron. F. Long-read sequencing of pooled 3’ RACE reactions initiated with exon 4 (left) or TRβ2 exon 1 (right) performed on cDNA libraries from 21 LM cells (top) and 5 LR cells (bottom). Schematics show total coverage (above) and individual transcripts (below). TRβ1 and TRβ2 first exons (green boxes) are enlarged at right. Red arrowheads: intronic premature transcription termination (PTT).

Two post-mitotic photoreceptor precursor populations rapidly expressing rod or cone markers.

A. UMAP plot colored by high resolution clusters. B. RNA Velocity plots with cell clusters as in A and cell velocity vectors. C. Enlarged view highlighting RPC and MG clusters (left), RPC-localized iPRP and TR clusters (middle), and RNA velocity (right). Black line: limit of the RPC-localized region. Arrows depict inferred trajectories. D. Violin plots depicting expression of selected genes in RPC, MG, and RPC-localized iPRP and TR cells. Colored asterisks compare cluster of same color to cluster at right of line. E. SCENIC regulon violin and box plots for RPC-localized cells in each cluster, selected from most specific regulons for MG, RPC, ER and LM clusters. ^, p <0.1; *, <0.05; **, <0.005; ***, <0.0005 (post-hoc Dunn test).

An iPRP population with cone-rod properties.

A-C. UMAP ‘bridge’ region cells colored by ER and iPRP cluster and RNA velocity (A), rod and cone marker gene expression (B), and NRL and THRB regulon activity (C, arrows indicate cells with both regulon signals). D. Box plot of Z-score-normalized NRL and THRB regulon AUCs for each cluster. E,F, Tiled composite fluorescence image of FW14 retinal section after immunofluorescence staining of RXRγ and NR2E3 and RNA FISH of GNAT2 and NR2E3 (E) and diagram of same retina with the section’s most peripheral (i.e., youngest and least mature) and most central (i.e., oldest and most mature) rods and cones labeled (F). Boxes indicate regions further evaluated as in Supplementary Figure S11 and quantitated in G and H. G,H. Quantitation of (G) outer and (H) inner layer photoreceptor precursors expressing combinations of RXRγ and NR2E3 proteins and GNAT2 and NR2E3 RNAs.

CHRNA1 marks early-born cone and rod precursors.

A. UMAP plots of iPRP bridge region marker genes. B-D Combined RXRγ/NRL immunohistochemical staining and CHRNA1 FISH of FW12 retina. B. Tiled images of retina section with nuclei stained with DAPI. White boxes: fields used for quantitative fluorescent imaging. Distances along apical edge of tissue marked in mm from midpoint of central image (0 mm, C). *: Imaged region shown in C. Scale bar = 500 µm. C. Top: Retinal nuclear and cellular segmentation and identification of cells as RXRγ+ (green outline) or NRL+ (red outline). Yellow box: Field shown below. Bottom: RXRγ or NRL immunofluorescence staining with CHRNA1 FISH. Arrows: RXRγ+CHRNA1+ (green), NRL+CHRNA1+ (yellow). Scale bars = 15 µm. H. Quantitation of fluorescent puncta in RXRγ+ and NRL+ cells by image field. X-axis: Distance from the midpoint of each image to retina center (0 mm, C). **, p <0.005 (Wald test, images from 0-6 mm).

L/M cone subcluster marker genes, regulons, and pseudotemporal trajectory with successive lncRNA gene expression.

A. Violin plots of high-resolution cone cluster marker genes with increasing maturation-associated expression. All significant differences between adjacent clusters are indicated. B. Violin plots of regulons with highest LM1-4 cone cluster specificity scores. *, p <0.05; **, <0.005; ***, <0.0005 (post-hoc Dunn test). C. Pseudotemporal trajectory through the L/M cone population derived with Monocle 3. *: Root cell used to define endpoint. D. Trendlines of relative count expression (left) and UMAP plots for lncRNAs correlating with early or late-upregulating modules. Line color matched to label boxes. E, F. Combined RXRγ immunohistochemical staining and FISH of lncRNAs on FW16 retina. E. Tiled images of retina with nuclei stained with DAPI. White boxes: fields used for quantitative fluorescent imaging. Distances along apical edge of tissue marked in mm from fovea to ciliary margins. Scale bar = 500 µm. Asterisks identify fields shown in F. F. Combined RXRγ immunostaining and multiplex FISH for four lncRNAs, of which two are shown in peripheral and central retina regions. Arrows: White: RXRγ/CTC-378H22.2+. Magenta: RXRγ/HOTAIRM1+. Blue: RXRγ/CTD-2034I21.2+. Green: RXRγ/RP13-143G15.4+. Scale bar = 15 µm. G. Quantitation of lncRNA fluorescent puncta assigned to RXRγ+ cells after segmentation. Colored bars mark lncRNA expression regions as described in the text. Arrows indicate significant count differences between images (p <0.05, Wald test).

Cone intrinsic MYCN and SYK contribute to the proliferative response to pRB loss.

A. Immunofluorescent staining of MYCN in ARR3+ cones but not in NRL+ rods in FW18 retina. B-D. UMAP plots of MYCN expression (B), MYCN regulon activity (C), and SYK expression (D). E. SYK and MYCN gene expression violin plots by cluster. *, p <0.05; ns = not significant (t-test). F. Immunohistochemical staining of SYK and cone arrestin (ARR3, LUMIF-hCAR antibody) in FW18 and FW16 retinae. Green arrow: ARR3+, SYK+. White arrow: ARR3+, SYK-. Scale bar = 25 µm. G. Effect of SYK inhibitor GS-9876 on KI67+ expression in shRB1- and YFP transduced RXRγ+ cone cells from dissociated FW16.5 fetal retina. Values represent means of three immunocytochemical analyses from two treatment replicates. Error bars: standard deviation. *, p <0.05; **, <0.005 (Student’s T-test with equal variance, 2-tailed). RXRγ+ cells: Experiment 1, n=1,340. Experiment 2, n=804. Range 107 – 366 cells per condition. H. Model of SYK expression and function in cone maturation and retinoblastoma development.

scRNA-seq sample and sequencing summary.

A. Retinae and cell numbers examined at each age. B. Histogram of total read counts for each cell, colored by sequencing run. Seq 1 cells were isolated by C1, Seq 2 cells by FACS or C1 as indicated, and all others by FACS. Dotted line: 100,000 read cutoff for cell exclusion. C. Box plots of read counts, genes detected, and Ensembl transcript isoforms detected per cell ordered by fetal age and specimen number. D-F. UMAP plots colored by sequencing run (D), isolation method (E), or retina ID (F).

Expression of marker genes of RPCs and Müller glia (A), photoreceptors (B), rods (C), and cones (D).

Insets: Gene expression violin plots (from left to right): RPC/MG (red), iPRP (brown), LM cone (green), S cone (teal), early rod (blue), late rod (pink). UMAP and violin plots for any gene or Ensembl transcript can be produced at https://docker.saban.chla.usc.edu/cobrinik/app/seuratApp/.

Differential expression between rod and cone maturation states.

A, C, E. Volcano plots of differential expression (pAdj <0.05, log2FC >|0.5|) between ER and LR (A), LM and S (C), and spatially separated late maturing L/M cones and remaining LM cluster (E). Labels indicate genes with highest significance and fold change. B,D. Overrepresentation of gene ontologies for genes upregulated in LR over ER (GO-Molecular Function instead of Hallmark50) (B) or upregulated in LM over S (D). F. UMAP plots of upregulated genes in the late maturing L/M cone group.

Differential NRL transcript splicing in early rod and cone precursors.

Top: Sashimi plots of NRL transcript splicing based on nanopore long read sequencing of 5’ RACE reactions from 23 ER cells and 21 LM cells. The minimum number of reads required to display transcript isoforms was set at 20 for ER cells and at 10 for LM cells, proportional to the total NRL reads in each sample. Bottom: Exon structure of FL and Tr NRL isoforms and RACE primer location. The NRL gene is oriented relative to chromosome 14 coordinates. Hash marks indicate sequences preceding a rarely used far upstream first exon. Note increased intra-exon 2 splicing and ex1T use in the LM population.

Cell type-specific NRL protein expression.

A. NRL antibody validation by immunoblot (IB) of ectopic FL-NRL and Tr-NRL in HEK293T. Left: anti-total NRL. Right: anti-N-terminal NRL specific to FL-NRL. B. Average mapped NRL isoform expression in scRNA-seq of RB31 retinoblastoma cell line. C. Immunoblot of endogenous total NRL in RB31 cells treated with or without 10 μM retinoic acid (RA) and 10 μM proteasome inhibitor MG132 and in FL-NRL transfected HEK293T. D. pULTRA-EGFP-P2A-Tr-NRL lentiviral vector used for explanted retina transduction. E. NRL Immunofluorescent staining of HEK293T cells 48 h after lentiviral transduction of Tr-NRL-GFP. Scale bar = 50 μm. F,G. Immunostaining of total NRL (F) or N-terminal-NRL (G) in explanted FW16.5 fetal retina 7 d after lentiviral transduction of Tr-NRL-GFP. Arrows: Yellow = GFP+, NRL-. Green = GFP+,RXRγ+. White = GFP-,FL-NRL+. Scale bar, 25 μm.

Differential THRB isoform assignments in rod and cone precursors.

A. Mean THRB isoform assignments for each cluster presented as total transcript counts (top) and percentage of total counts (bottom) as in Figure 4B, with isoforms grouped according to capacity to encode full-length canonical TRβ1 or TRβ2 or to have PTT following ex4, ex5, or ex6 (exons numbered as for ENST00000396671) regardless of differences in 5’ noncoding exons and 3’ poly(A) sites. B. THRB Ensemble isoforms identified in panel A with ENST00000396671 exon numbers indicated, accessed June 22, 2024, at: http://useast.ensembl.org/Homo_sapiens/Gene/Summary?g=ENSG00000151090;r=3:24117153-24495756.

THRB isoforms and premature transcription termination in cone and rod precursors.

Top: Sashimi plots of THRB transcript splicing based on long read nanopore sequencing of 3’ RACE reactions from 23 LM cells and 5 LR cells. RACE reactions were initiated from primers specific to TRβ1 exon 4 (left) or specific to the unique TRβ2 exon (right). Gray arrows: premature transcription termination following TRβ1 exon 4, the TRβ2 first exon, the shared exon 6, and a novel exon used solely in TRβ2 transcripts. Bottom: TRβ1 and TRβ2 exon structures and RACE primer positions. The THRB gene is oriented relative to chromosome 3 coordinates.

Similar RXRG isoform expression in rod and cone precursors.

A. Expression of RXRG and the two most highly assigned RXRG isoforms (ENST00000359842 and ENST00000619224). B. Mean RXRG isoform assignments for each cluster presented as total transcript counts (top) and as a percentage of the total counts (bottom). Lines indicate significant LM vs. ER fold change, colored by isoform. *****, p <0.000001 (bootstrapped Welch’s t-test). C. Top: Mean read counts across Ensembl RXRG exons. Bottom: Transcript structures for ENST00000359842 and ENST00000619224. Red arrowheads: Extended (Ext) RXRG 5’ UTR. D. Long-read nanopore sequencing of pooled 5’ RACE initiated with RXRG exon 4 primers from 23 ER cells (top) and 21 LM cells (bottom). Schematic show total exon coverage (above) and individual transcript (below). Panels at left show RXRG locus with exon 1 in green box enlarged at right, revealing a similar range of RNA 5’ ends in ER and LM populations. E. Immunohistochemical analysis of RXRγ and NRL in FW16 retina. Arrows: Yellow = RXRγ-high cone. White = NRL+ rod with weak RXRγ. Scale bar = 25 μm. Right: Boxplot of RXRγ mean grey values in NRL+,RXRγlo and NRL−,RXRγhi cells. ***, p<0.0001 (Welch’s T-test).

High resolution cluster marker genes and differential RPC versus MG gene expression.

A. Marker gene dot plot for high resolution Louvain clusters. Dot size: Percent of cells expressing gene. B. Volcano plot of differential expression (pAdj <0.05, log2FC >|0.5|) between MG and RPC clusters with genes labeled when p <10-4 except for RPL34. C. Overrepresentation analysis for genes upregulated in RPC above pAdj and log2FC cutoff.

High-resolution regulon-defined RPC and photoreceptor precursor states.

Ward-clustered heatmap of the highest scoring SCENIC regulons in each high-resolution cluster, displaying Z-score normalized regulon activities. Subcluster labels: RPC-l = RPC-localized, Br = Bridge-localized, O = Other/remainder of original cluster, Late = late maturing LM4 cones.

NR2E3 and GNAT2 co-expression in early-born photoreceptor precursors.

A. Composite images of combined RXRγ and NR2E3 immunofluorescence staining and GNAT2 and NR2E3 RNA FISH in boxed regions 2-13 of Figure 6E. B. Enlarged images of boxed regions in images 2, 3, 5, and 8 in panel A. Dotted yellow lines demarcate the outer and inner NBL cells analyzed separately in Figure 6G,H. White arrows = GNAT2 puncta. Green arrows = NR2E3 puncta. Boxed regions in column 3 (GNAT2 RNA plus NR2E3 protein) and column 4 (NR2E3 RNA) are enlarged in column 6 (GNAT2 RNA plus NR2E3 RNA) and illustrate region 2) initial expression of GNAT2 RNA without NR2E3 RNA in outer NBL nascent cones; region 3) initial co-expression of GNAT2 and NR2E3 RNA in early outer NBL cones; region 5) initial expression of NR2E3 RNA with or without co-expressed NR2E3 protein; and region 8) extensive co-expression of GNAT2 and NR2E3 RNA without NR2E3 protein in RXRγhi outer NBL cones.

iPRP bridge region marker genes and regulon activities.

A. UMAP plots of iPRP marker S100A6 and photoreceptor transition state marker DLL3 gene expression. B. UMAP plots of OLIG2, LHX9, and THRB regulon activities (left) and OLIG2, LHX9, and ONECUT1 gene expression (right), illustrating high OLIG2 regulon activity with minimal OLIG2 gene expression, concurrent LHX9 regulon activity and LHX9 gene expression, and ONECUT1 gene expression in cone-directed iPRP cells preceding increased THRB regulon activity in early LM cones. Boxes indicate bridge region.

Heatmaps of coregulated gene modules across cone precursor pseudotime.

A. Heatmaps of individual genes in each module across the pseudotime trajectory shown in Figure 6C. lncRNA genes of interest are labeled. B. Heatmap of averaged gene module expression across pseudotime.

Differential gene expression in early cone and rod precursors.

A. Volcano plot of differential expression between early rod cluster ER and cone LM, excluding late-maturing OPN1LW+ population. pAdj cutoff =0.05, log2FC cutoff = |0.4|. Labeled genes: pAdj <10e-32, except for SYK and MYCN. B. Overrepresentation analysis of cone-enriched genes pAdj <0.05, log2FC ≥ |0.4|.