Single-cell analysis of olfactory neurogenesis and differentiation in adult humans

Abstract

The presence of active neurogenic niches in adult humans is controversial. We focused attention to the human olfactory neuroepithelium, an extracranial site supplying input to the olfactory bulbs of the brain. Using single-cell RNA sequencing analyzing 28,726 cells, we identified neural stem cell and neural progenitor cell pools and neurons. Additionally, we detailed the expression of 140 olfactory receptors. These data from the olfactory neuroepithelium niche provide evidence that neuron production may continue for decades in humans.

Conflict of interest statement

Competing Interests

J.W.H. is the inventor of intellectual property related to prognostic testing for uveal melanoma. He is a paid consultant for Castle Biosciences, licensee of this intellectual property, and he receives royalties from its commercialization. No other authors declare a potential competing interest.

Figures

Extended Data Fig. 1. DotPlot visualization listing scRNA-seq clusters.

a, Cell phenotypes listed on y-axis, showing unbiased gene expression for the top 8 genes per cluster identified by log Fold Change; genes (features) are listed along the x-axis. Dot size reflects percentage of cells in a cluster expressing each gene; dot color reflects expression level (as indicated on legend). The plot depicts clusters from 28,726 combined olfactory and respiratory mucosal cells, n=4 patients. b, DotPlot visualization of the Heatmap shown in Fig. 1d. The plot depicts clusters from 28,726 combined olfactory and respiratory mucosal cells, n=4 patients. c, Histogram showing the neuronal lineage cell types captured by scRNA-seq, as a percentage of the total cells analyzed per patient; see Fig. 1c.

Extended Data Fig. 2. Additional human immunohistochemistry of basal cell populations.

Co-staining for SOX2, Ki67 or the HBC marker TP63. Proliferative activity has been used as a hallmark of the globose basal cell (GBC) phenotype. We reasoned that, although some proliferating cells in the olfactory epithelium (OE) might be immune or inflammatory cells, proliferative cells within the GBC layers of the OE that are SOX2+/Ki67+/TP63− would be categorized as GBCs. a, Sustentacular cell nuclei at the top of the OE are SOX2-bright; and horizontal basal cells (HBCs) and a subset of GBCs are SOX2+, although less intensely. Arrow marks a SOX2+/KI67+ cell among the proliferative KI67+ basal region, consistent with the GBC phenotype. b, SOX2 co-localizes with TP63 in HBCs; arrows mark Sox2+/TP63− GBCs. c, TP63+ HBCs are mitotically quiescent, while many GBCs are actively proliferating, often in scattered cell clusters. Arrow marks a cluster of KI67+/TP63− basal cells. Dashed line indicates basal lamina. Immunostaining for a-c was conducted in triplicate with similar results. Bar, 50 μm.

Extended Data Fig. 3. Additional human immunohistochemistry (IHC) of immature and mature olfactory sensory neurons (OSN) populations.

a, Co-staining for LHX2 and OMP demonstrates many LHX2+/OMP− neurons, distributed in deeper layers of the OE, which are immature OSNs; OMP is a marker for fully differentiated OSNs, while LHX2 expression in differentiating OSNs orchestrates OR receptor expression. b, DCX was identified by scRNA-seq here as enriched in immature OSNs (see Fig. 2a). IHC confirms scattered DCX+ neuronal somata and dendrites in the OE. c, Similarly, the bHLH transcription factor OLIG2 was identified to be enriched in immature OSNs; IHC confirms nuclear expression in the deeper OSN layers of OE tissue. Dashed line indicates basal lamina. Immunostaining for a-c was conducted in triplicate with similar results. Bar, 25 μm.

Extended Data Fig. 4. Analysis of immune cell populations.

Feature plots indicate expression of inflammatory cell markers in human nasal biopsy samples. UMAP clustering identifies lymphocyte populations; a, cytotoxic T cell markers co-localize including CD8A, PRF1, GZMA and GZMB. The plot depicts clusters from 28,726 combined olfactory and respiratory mucosal cells, n=4 patients. b, Within the monocyte/macrophage populations (CD14+, CD68+ cells), markers for activated M2 macrophages, such as CD163 and IL10, are indicated. The plot depicts clusters from 28,726 combined olfactory and respiratory mucosal cells, n=4 patients. c, DotPlot visualization of additional immune cell gene expression from combined aggregate samples. Cell cluster identity is listed on the y-axis, genes (Features) are listed on x-axis. The plot depicts clusters from 28,726 combined olfactory and respiratory mucosal cells, n=4 patients. d, Immune cell DotPlots showing the contribution of cell types and gene expression patterns by individual patient sample. The plot depicts clusters from 28,726 combined olfactory and respiratory mucosal cells, n=4 patients.

Extended Data Fig. 5. Focused UMAP plot of OE neuronal lineage populations, with cell phenotype assignments indicated.

Compare with gene expression feature plots in Fig. 1e and 2f. The plot depicts clusters from 694 GBCs, immature olfactory neurons and mature olfactory neurons, n=4 patients.

Extended Data Fig. 6. Gene set enrichment analysis on differential expression data from selected cell clusters.

The top 50 Reactome pathways ranked by adjusted p-value were plotted in the visualization. a, mOSNs versus GBCs; many top terms involve neuronal, transduction and synapse function. The differential expression was calculated from 222 mOSNs and 115 GBCs, n=4 patients. The default two-sided non-parametric Wilcoxon rank sum test was utilized with bonferroni correction using all genes in the dataset. b, GBCs versus olfactory HBCs; top terms include cell cycle or neurogenesis functions. The differential expression was calculated from 115 GBCs and 2,182 olfactory HBCs, n=4 patients. The default Seurat two-sided non-parametric Wilcoxon rank sum test was utilized with bonferroni correction using all genes in the dataset.

Extended Data Fig. 7. OR gene expression in human olfactory neurons.

a, Range of gene expression in our datasets (binned to 0.01). We identified 4.80E+07 observations (gene expression measurements) expressing >0. Genes with no expression (=0, n=532063621) were excluded in this plot. The distribution plot shows that choosing a cutoff of 0.5 (red vertical dotted line). b, Doublet analysis. Box plots depicting the number of UMIs (“nCount_RNA”, left plot), and genes (“nFeature_RNA”, right plot) in immature and mature neurons expressing 1 or 2 ORs.

Extended Data Fig. 8. Principal component determination analysis.

The 4-patient combined data set was analyzed in Seurat to explore p principal components (PCs) contributing to heterogeneity, and to determine an appropriate PC selection. a, Using the JackStraw approach, approximately 100 PCs had low a p-value. The plot depicts PC calculated from from 28,726 combined olfactory and respiratory mucosal cells, n=4 patients. The jackstraw test implemented in Seurat was used to calculate p-values of PCs. b, To select a suitable number of PCs for downstream analysis, we used the elbow plot heuristic approach, indicating that beyond 20–30 PCs, very little additional variation is explained. Therefore, for downstream analysis we chose to include 30 PCs.

Extended Data Fig. 9. scRNA-seq quality control plots.

a, number of genes (features) per cluster. The plot depicts clusters from 28,726 combined olfactory and respiratory mucosal cells, n=4 patients. b, number of UMIs (nCount) per cluster. Cluster cell type identities are listed along the x-axis. Violin plot widths are proportional to the density of the distribution.

Fig. 1 |. Aggregate analysis of 28,726 single cells from human olfactory cleft mucosa.

a, Schematic diagram of the respiratory epithelium versus olfactory epithelium. Abbreviations: globose basal cells (GBCs), Horizontal basal cells (HBCs), Bowman’s ducts (BD), Bowman’s glands (BG), Vascular smooth muscle (VSM), Endothelial cell (EC), Pericyte (PC) White blood cell (WBC), Macrophage (MP), Olfactory ensheathing cell (OEC), Cranial Nerve I (CN1), basal lamina (BL) . b, UMAP dimensionality reduction plot of 28,726 combined olfactory and respiratory mucosal cells, n=4 patients. Cell cluster phenotype is noted on color key legend/labels. c, Plots of individual patient samples, n=4 patients. Patient 1: 5,683 cells; Patient 2: 11,184 cells; Patient 3: 5,538 cells; Patient 4: 6,321 cells; see also Extended Data Fig. 1c. d, Heatmap depicting selected gene expression among olfactory cell clusters. e, UMAP depicting GNG8 and GNG13 expression in 694 GBCs, immature olfactory neurons, and mature olfactory neurons, n=4 patients. f, GNG13 immunostaining (red) in adult human and mouse OE; dashed line marks basal lamina; nuclei are stained with DAPI (blue). Immunostaining was conducted in triplicate with similar results. Scale bar, 25 μm.

Fig. 2 |. Gene expression analysis of human OE.

a, DotPlot visualization of neuron lineage cell populations from adult human OE, n=4 patients; iOSN, immature olfactory sensory neuron; mOSN, mature olfactory sensory neuron. b-e, Cell type-specific marker validation in human versus mouse olfactory mucosal sections. TUJ1 labels somata of iOSNs, which are more abundant in our adult human samples (n=3, two-sided Welch’s t-test, p=0.015). SOX2 marks basal and sustentacular cells; KRT5 labels HBCs, which in many areas in human samples have a reactive rather than flat morphology (boxed region, enlarged) and are abundant (n=3, two-sided Welch’s t-test, p=0.04). Measure of center and error bars for b-e are mean ± standard deviation. KI67 marks proliferative GBCs. LHX2 was highly expressed in iOSNs (DotPlot in a); IHC confirms widespread expression in human and variable expression in mouse; nuclei are stained with DAPI (blue); dashed line marks basal lamina. Scale bar, 50μm; 10 μm in inset in c. f, Focused UMAP plots visualizing gene expression of HES6 and NEUROD1 in GBCs, n=4 patients.

Fig. 3 |. Analysis of OR expression in human olfactory epithelium.

a, Expression level of all ORs (n = 545 total receptors) in immature (GNG8+) and mature (GNG13+) OSNs. Expression levels of a VN1R1 receptor in 1 cell is indicated in red. Y-axis represents normalized expression values, x axis individually expressed receptors. b, ORs expressed in individual immature and mature neurons. * p<0.05; *** p<0.001, two-sided χ test without Yates’ correction, n=3 patients. Y-axis represents percentage of cell population, x-axis the number of unique ORs per cell. c, Most commonly identified ORs in OSNs. The top 8 ORs in this list were detected statistically more than expected. d, Expression of OR families in immature and mature neurons; blue=Class II, orange=Class I ORs. e, Co-expression matrix of ORs; X- and Y- axis contain every OR found to be co-expressed with at least one other OR (n = 141, including VN1R1). Blue squares show expression of intersecting ORs co expressed in 1 cell, orange indicates 2 cells, and red 3 cells. f, List of most co-expressed ORs from Fig. 3e. “Total ORs” indicates sum of all ORs found co-expressed with indicated member, and “Unique ORs” the sum of unique OR genes co-expressed, in the indicated number of neuronal cells. The top six ORs in this list are statistically more co-expressed than expected.