Improvement of the gene model for scRNA-seq of ferretes.

A Schematic representing cortical development and emergence of diverse neural progenitors in humans and ferrets. Progenitors sequentially generate DL- and subsequent UL-neurons, and finally glial cells, and form ependymal cells. RG and oRG have been morphologically and positionally identified in both humans and ferrets while tRG have been only reported in humans. VZ, ventricular zone; ISVZ, inner subventricular zone; OSVZ, outer subventricular zone; IZ, intermediate zone; DL, deep layer; UL, upper layer; CP, cortical plate; NEP cell, neuroepithelial cell; Mn, migrating neuron; and EP cell, ependymal cell.

B Schematic representing the experimental design and time points used to build the transcriptome atlas of developing somatosensory cortex of ferrets. Single cells were isolated using 10X Chromium (see Materials and Methods).

C Linkage of the 3’UTR sequence and the coding sequence of genes has been improved in this study. This linkage is necessary to assign the gene coding for a mRNA as far as 10x Genomics Chromium kit is used for scRNA-seq (see Materials and Methods).

D Table comparing quality-control metrics of an alignment with either MusPutFur 1.0 (UCSC gene models) or MusPutFur 2.60 (this study) using E34 samples. The total number of genes detected and median genes per cell were higher with MusPutFur 2.60.

E Violin plots showing the number of genes, mRNAs and the percentage of mitochondrial genes per cell in each sample and time point.

Single-cell RNA-seq reveals ferret transcriptome signatures and the cell types.

A UMAP visualization showing cells colored by Seurat clusters and annotated by cell types (as shown in C, D).

B Heatmap showing expression profiles of cluster marker genes based on log fold-change values. Cells were grouped by Seurat clustering (transverse). Cell types were assigned according to the expression of marker and differentially expressed genes in each cluster. Early RG was defined by the expression of HMGA2, LDHA, and LIX1; and late RG, by PTN, ALDOC, and FABP7. OPC, oligodendrocyte precursor cell; ITN, interneuron; and EN, endothelial cells. Other cell type abbreviations are shown in Fig. 1A and the main text. Color bar matches the Seurat clusters in (A). The ten most enriched representative genes in each cluster are shown, with typical marker genes noted in (CE).

C Normalized expression levels of representative marker genes of different cell types projected onto the UMAP plot in (A). MG, microglia.

D Expression pattern of CRYAB, a marker for tRG, only described in humans and other primates.

E The expression pattern of HOPX, a marker for oRG in humans and other primates in the UMAP plot.

Truncated tRG emerges during postnatal cortical development in ferrets.

A Representative images showing the cellular features of tRG and aRG in VZ at P0, and P5, respectively, stained for GFP (red) and CRYAB (green). RG cells were sparsely labelled with a GFP-expressing plasmid at E30 via IUE for P0 samples, and at P3 for P5 samples. MAX projection was performed on a 30 µm vibratome section with 5 µm or 2.5 µm interval for each z-image at P0 or P5, respectively. Scale bar, 50 µm.

B Expression patterns of CRYAB (green) and PAX6 (red) in ferret germinal zones during postnatal development (P1, P5, and P10) (cryosection thickness = 12 µm).

C Developmental profile of CRYAB expression in ferret cortices. Immunostaining for CRYAB (green) and DAPI (gray) on cryosections at E32, E38, P3, P5, and P10 (cryosection thickness = 12 µm).

D Quantification of CRYAB+ cells among all nuclei (DAPI) in VZ strips through immunostaining, indicating a great increase in CRYAB-expressing cells (strip width = 150 µm; n = 3 for P5; n = 2 for P10; Wilcoxon rank-sum test p-value = 6.574e-08).

E–G Expression of other markers in CRYAB+ cells. Representative images taken with a 100X- objective are shown with MAX projection. The border of VZ is shown with a white line. Below each double staining image, quantification of KI67-, TBR2-, or OLIG2-expressing cells among the CRYAB+ cell population in the VZ is shown (n = 2 for P5; n = 2 for P10 for each staining except n = 3 for TBR2 staining). Box and whisker plots indicate ranges (lines) and upper and lower quartiles (box) with the median. Scale bars, 20 µm.

Temporal fates of RG cells predicted by pseudo-time trajectory analysis.

A Pseudo-time trajectory tree of progenitor cells in the ferret developing cortex (assembled with Monocle v2). Cellular distribution at each state is the same as in the UMAP plot shown in Fig. 1C. Cell types representing each state and their gene markers are shown next to each state (see below).

B Trajectory trees split by collection stages (AG samples are shown).

C Distribution of cells expressing marker genes for stem cell states (HES1 and HMGA2) and the neurogenic state (EOMES) along trajectories. Color densities indicate the log-normalized unique molecular identifier count for each gene.

D Normalized expression of three genes: SPARCL1, CRYAB, and FOXJ1 in the state 6 and state 7 cells in the pseudotime trajectory tree (A). The combination of expression levels of the four markers gene discriminates state 6 (ependymal) and state 7 (astrogenic). See the text for details.

E Distribution of tRG cells along the tree and tRG-focused UMAP visualization. The tRG cells (n = 162) were found on the three states, NPC3 (9.3%; 15 cells), ependymal-tRG (61.1%; 99 cells), and astrogenic-tRG (29.6%; 48 cells).

F Composition of the three types of tRG in (D) by the collection stage (AG and T samples are combined). No tRG cells were collected from E25 and E34.

tRG adopt both ependymal cells and astrogenic fates in ferrets.

A Immunostaining of 12 µm cortical cryosections with CRYAB (green) and FOXJ1 (red), focusing on the VZ at E40, P3, P5, P10, P14, and P35. Scale bar, 50 µm. CRYAB and FOXJ1 double positive nuclear rows are shown with arrowheads at P5 and P10 when this nuclear alignment is visible.

B Percentage of FOXJ1 expression in CRYAB-expressing cells (n = 3 for P5; n = 2 for P10; Wilcoxon rank sum test p-value = 1.749e-05). Images with merged channels in A are shown with the same color codes, antibodies and scale bars as A.

C Number of cells expressing both CRYAB and FOXJ1 in the mid and late RG, or tRG clusters shown in Fig. 1C (with unique molecular identifier counts higher than 0.25). Colored bars indicate the stages of sample collection. CRYAB- and FOXJ1-expressing cells increased over time only in the tRG cluster.

D Normalized expression levels of FOXJ1 in each cell in the indicated ferret clusters in Fig. 1C.

E Cortical origin and shape of FOXJ1-expressing ependymal cells indicated with white arrows. Staining for FOXJ1 at P14 after labeling cortical progenitors with an mCherry-expressing vector via IUE at E30. The maximum projection images with 1 µm z-step size are shown. Cryosection thickness = 12 µm; scale bar = 20 µm.

F Normalized expression of SPARCL1 and FOXJ1 in each cell in ferret tRG clusters separated by the pseudo-time trajectory analysis (see Table S3).

G Relationship between SPARCL1 and FOXJ1 transcripts in individual cells within each tRG cluster. Within the 15 cells that were classified in state 5 tRG, only one expressed FOXJ1 mRNA. Pearson relationship analysis was performed using cells from state 6 and 7 (R=0.19 p value=0.018).

H Immunostaining of a P10 cortical tissue for FOXJ1 protein and SPARCL1 mRNA.

I Cells within the VZ in (H) were clustered into four classes by FOXJ1 protein± and SPARCL1 mRNA±. Cluster A: SPARCL1/FOXJ1-/-, B: SPARCL1/FOXJ1 +/-, C: SPARCL1/FOXJ1-/+, D: SPARCL1/ FOXJ1 +/+ (Fig. S5C, S5D, Table S4). See Materials and Methods for clustering. Scale bar = 20 μm.

Comparison of molecular identity of RG subtypes between ferret and human.

A UMAP visualization of integrated human (n = 2,672; left) and ferret (n=30,234; right) single-cell datasets colored according to the different clusters. The names of clusters from human and ferret cells begin with “H” and “F,” respectively.

B Homologous temporal pattern of transcriptomic characters of progenitors and their progenies between ferrets and humans. Corresponding neurogenic stages between ferrets and humans can be assigned by cell distribution at homologous positions in UMAP plots. Gliogenic RG cells (a subtype of the “late_RG” group and OPC in Fig. 2A) were first distinguished transcriptionally at ferret E40 and at human GW14, as pointed out by arrows.

C, D Correlation coefficient (C) and significance (D) between indicated clusters of ferrets and humans, calculated by marker gene scores. See the text and Materials and Methods for details.

E Normalized expression levels of the indicated genes in humans and ferrets from each progenitor cluster. See the text and Materials and Methods for details.

F tRG scores of the indicated clusters are presented as box and whisker plots for humans and ferrets. As for the definition of cluster scores, see the text and Material and Methods. These represent the range without outliers (lines) and upper and lower quartiles with the median (box). Outliers are represented as points outside the box. Outliers are 1.5-fold larger or smaller than interquartile range from the third or first quartile, respectively.

- Identification of tRG subtypes in humans and oRG in ferrets.

AUMAP visualization of human brain cells at GW25 (and removing neurons and other cell types). Cells are colored by cell type and identified by marker genes (data from Bhaduri et al, 2021).

B Schematic of the integration strategy of human and ferret subtypes. We merged the two datasets by pairing mutual nearest neighbor cells (MNN) following canonical correlation analysis (CCA) across species (Stuart and Butler et al, 2019).

C UMAP visualization of integrated ferret and human datasets colored and numbered by different clusters.

D, E Identification of three tRG subtypes in ferrets (D) and humans (E). The red dots highlight the indicated tRG subtypes named after pseudotime trajectory analysis in Fig. 4: presumptively “uncommitted”, “astrogenic” and “ependymal”.

F Distribution of humans (left) and ferrets (right) tRG in the integrated dataset. Human and ferret tRG were identified from the separated dataset. The numbering corresponds to cluster numbers in the integrated clustering (C) of human and ferret datasets.

G Distribution of ferret tRG subtypes (identified by the pseudo-time analysis) in the integrated dataset. A major cell population of each of the 3 ferret tRG subtypes (classified by the pseudo-time analysis), more or less, belongs to a single cluster in ferret-human integrated clustering; astrogenic tRG corresponds to cluster 21; ependymal tRG to cluster 20; uncommitted tRG to cluster 7.

H Identification of oRG-like cells in ferret (left) and oRG cells in human (right), by the same way that we have done for tRG. All cells highlighted as red dots.

I oRG scores for oRG-like cells and other NPCs in ferrets. The result indicates that oRG-like cells share more transcriptomic features than other NPCs, whereas ordinary clustering among ferret NPCs failed to distinguish oRG-like cells from the rest. Box and whisker plots indicate ranges (lines) and upper and lower quartiles with the median (box).

J The expression levels of CRYM, CLU, and HOPX in oRG-like cells and other NPCs in ferrets.

Temporal pattern of neurogenesis and gliogenesis in the cerebral cortex of ferrets.

A Coronal section of the dorsal cortices immunostained for RG (PAX6, red), IPC (TBR2, green), and OPC (OLIG2, cyan) from early to late neurogenesis at E25, E32, E36, E40, P5, and P10. Scale bars = 100 µm. The position of the image strip at each stage is shown in the corresponding dorsal hemisphere image above the strip. The approximative boundaries of dorsal cortex area used for single-cell RNA sequencing are highlighted with dotty line segments in the dorsal cortex hemisphere above each strip.

B Immunostaining for CTIP2 (green) and SATB2 (red) at the same developmental stages as in (A), showing the onset of DL- (E25 or earlier) and UL-neurogenesis (E32–E34). Scale bars = 100 µm.

C Immunostaining for GFAP, a marker showing gliogenic potential, in cortical germinal layers at E32, E36, E40, P5, and P10. Gliogenic progenitors emerge around E40, while mature astrocytes with a typical astrocytic morphology appear later at approximately P10 onward in the outside of the germinal layers (Reillo & Borrell, 2012).

D, E mCherry labeling of the dorsal cortex at E30 showing that a population of OLIG2+ oligodendrocyte precursors are generated in the cerebral cortex by P10.

E. As reported for both mice and humans (Zheng et al, 2018; Rash et al, 2019b; Huang et al, 2020; Kessaris et al, 2006). Arrowheads in the cropped images of the left panels (E area in Fig. S1D) indicate cells that co-express mCherry, PAX6, and OLIG2 in the OSVZ. Scale bars = 20 µm.

– Quality and reproducibility of ferret single cell transcriptome dataset.

A The brain region used for the cellular isolation for a scRNA-seq are shown. Dotty lines indicate the approximative boundaries where the cortical slices were cut.

B Separate UMAP plots for each cell collection series that was prepared by different methods, after merged clustering. Left is the “AG” sample: FACS-based sorting of the neural stem cell fraction labeled with an AzamiGreen (AG)-driven by HES5 promoter (Ohtsuka et al, 2006). The “T” samples are shown in the right: cells forming the VZ, SVZ, and intermediate zones (IZ) of cerebral cortices after discarding the cortical plate. E25 cells, which were mostly RG, were not collected by cell sorting and put between AG and T in UMAP plots.

C Effects of cell filtration according to the mitochondrial content. We set the threshold to 10%, resulting in 28,686 cells in our dataset, and then performed the workflow from the normalization step with the same settings that had been applied to our original ferret dataset (Methods). Left is clusters after the 10%-filtered subset on UMAP, and the original clusters in Fig. 2A is right. Both conditions gave 26 clusters, and both major cell types and subtypes were conserved after filtering.

Cell types in developing ferret cerebral cortex.

A. Dotplot showing average expression and percentage of expression of representative marker genes for individual cell types in the merged dataset. The largest cluster was composed of UL neurons, and their four clusters (Fig. 2A) were combined as a single “UL” group. P.E., percentage of expressing cells per cell type; A.E., average expression per cell type.

B. Percentage of sample collection stages in each cell type (top half) and percentage of each cell type in the merged dataset (bottom half).

C. Heatmap showing the expression level of cell cycle marker genes based on log fold-change values. PCNA, DUT, SIVA1, and TYMS are markers of the S-phase; UBE2C, CDCA4, CENPF, and CCNB1 are markers of the G2M-phase. Three subclusters among early RG, late RG, and IPC clusters expressed different markers for cell-cycle states. The color bar at the top indicates Seurat clusters matching those in Fig. 2A.

D. Immunostaining for CRYAB+ tRG (green) and HOPX+ RG (red) together with DAPI on P1 and P10 germinal layers.

Development of ferret tRG during postnatal cortical development.

APie charts showing the percentage of CRYAB expression in the ferret clusters of the merged transcriptome dataset of all stages (Fig. 2). Cells were counted as positive if unique molecular identifier counts for CRYAB were > 1 using Seurat “WhichCells” function (slot = “counts”).

B Number of cells in the tRG clustered by age annotated in the ferret transcriptome dataset (Fig. 2).

C. Number of CRYAB-expressing cells in RG clusters of the transcriptome dataset (Fig. 2). CRYAB expression was undetectable at E25. CRYAB-expressing cells in the midRG and lateRG subtypes did not change in all stages (cell numbers: 58, 46, 68, 55, and 61 at E34, E40, P1, P5, and P10, respectively), while the number of CRYAB+ tRG cells increased after birth (cell numbers: 0, 21, 50, 91, and 118 at E34, E40, P1, P5, and P10, respectively).

D–F. Violin plots indicating the normalized expression of CRYAB (D), TBR2 or EOMES (E), and OLIG2 (F) in the ferret clusters shown in Fig. 2.

G. Cellular features of a tRG cell in cortex labelled with GFP at embryonic stages via IUE. Vibratome sections at P0 and P4 (thickness = 200 µm) were immunostained

Temporal fates of RG cells identified by pseudo-time trajectory analysis

A. Ferret cells used for pseudo-time trajectory analysis were lined up along the stretched pseudo-time axis. Expression levels of markers for typical cell types are also shown along the pseudo-time axis: LIX1 and HMGA2 for early RG, EOMES and NEUROG2 for IP, RGS20, and PTN for late RG, PDGFR and OLIG2 for glial cells, and FOXJ1 and FAM216B for ependymal cells.

B. Branching trees labelled by pseudo-time score. Color scale indicates log10 (pseudo-time +0.1) values. Each state is named after the major cell type.

C. The composition of each state is shown by clusters defined in the UMAP plot (Fig. 2A). Percentages for clusters were calculated from all cells in each state.

D. Branching trees split by UMAP clusters (Fig. 2).

Ferret tRG adopt both ependymal and astrogenic fates in ferrets

A. Onset of ciliogenesis in the developing cortex of ferrets, shown via staining for ADENINE CYCLASE III (cyan) along the VZ surface at E40, P5, P10, and P35. Scale bars = 10 µm.

B. Normalized expression of SPARCL1/cell in the major cell types shown in Fig. 2A.

C. Plots showing SPARCL1 mRNA count (left) and normalized FOXJ1 protein levels (right) in each cell clustered into individual classes (A–D) within the VZ shown in Fig. 4H (Table S4). The box and whisker plots indicate the range (lines) and the upper and lower quartiles with the median (box).

D. Distances of each cell (its center of gravity) from the apical surface in the individual clusters (A–D), as indicated by the box and whisker plots in (C). Cluster C and D cells are significantly located on the apical surface compared to cluster B cells. Data normality was analyzed by Anderson-Darling test and statical significance was calculated by Kruskal-Wallis test with Dunn’s multiple comparisons test.

Comparison of molecular identity of RG subtypes between ferret and human

A. UMAP visualization of integrated ferret (n = 30,234) and human (n = 2,672) single-cell datasets colored according to the different clusters. The names of clusters from human and ferret cells begin with “H” and “F,” respectively.

B. Expression patterns of the indicated genes in each cluster. Circle color indicates the expression level of each gene. Circles size indicates the percentage of cells expressing the gene in the indicated cluster.

Comparison of the tRG and oRG subtypes between ferrets and humans

A, B Matrices of correlations between two arbitrary RG subtypes from humans and ferrets by calculating marker gene scores as described for Fig. 5. (A) Correlations and (B) p-values between pairs in human and ferret clusters are shown. The early RG in ferrets showed the highest similarity with “OLIG1” in human data. This might be because of two reasons; (1) cells only at GW25 without early RG were used for comparison, and (2) cells in “OLIG1” highly expressed genes related with cell proliferation, as the early RG also did (Table S6). We did not observe any cell types in human data corresponding to the ferret “mid RG” with the same reason described above.

C Distribution of human and ferret tRGs in the integrated dataset. The number in the X-axis indicates the cluster number identified by the UMAP in Fig. 7C.

D Expression levels of the subtype maker genes in the three tRG-enriched clusters (7, 21, and 28). Marker enrichment is consistent with the presumptive identification of tRG subtypes according to pseudotime trajectory analysis.

E Expression patterns of CRYM, CLU, and HOPX in oRG cells and other NPCs in humans.

F In situ hybridization of CLU at the ferret cortex showed its expression in all VZ, ISVZ, and OSVZ at P5.

The absence of a tRG cluster in human organoid datasets.

A. UMAP visualization of human organoid-driven cells colored by the original annotations for the cell line (left), age (middle) or state (right). Human organoid dataset from Bhaduri et al., 2020 was used. Bhaduri dataset contained organoids generated from 4 different lines, which showed a variability in terms of cell distribution on UMAP while overall temporal and differentiation axes were recapitulated.

B. UMAP visualization of cells by clusters (left), original cluster annotations (middle), or by the normalized expression level of CRYAB (right). Red circles on UMAP highlight the group of cells expressing CRYAB. CRYAB expression was restricted to two out of four cell lines. CRYAB was expressed in clusters 26 and 30 enriched in the YH10 line, and cluster 29 enriched in the13234 cell line.

C. UMAP visualization of integrated primary and organoid datasets from Bhaduri et al., 2020. Cells are colored by the originated dataset (left), primary tRG cluster (middle) or organoid dataset (right). Black circles were used to indicate the overlap of the region where human primary tRG was found. Human tRG did not overlap with organoid cells. We found that CRYAB-expressing organoid clusters 26 and 30 overlapped with “oRG/astrocyte” clusters of primary tissues.

D. UMAP visualization of cells derived from Herring et al., 2022 human organoid dataset, colored by clusters (left) or by the normalized expression level of CRYAB (right).