Molecular, cellular, and developmental organization of the mouse vomeronasal organ at single cell resolution

  1. Max Henry Hills
  2. Limei Ma
  3. Ai Fang
  4. Thelma Chiremba
  5. Seth Malloy
  6. Allison R Scott
  7. Anoja G Perera
  8. C Ron Yu  Is a corresponding author
  1. Stowers Institute for Medical Research, United States
  2. Department of Cell Biology and Physiology, University of Kansas Medical Center, United States
7 figures, 1 table and 1 additional file

Figures

Figure 1 with 3 supplements
Single cell transcriptomic profile of the whole vomeronasal organ.

UMAP visualization of integrated cell-type clusters for whole-VNO single-cell RNA-seq. (B) Cell-type marker-gene normalized expression across the cell clusters. C. UMAP of cell-type clusters split by age. (D) UMAP of cell-type clusters split by sex. (E) A representative image of transcript distribution for 9 genes in a VNO slice using the Molecular Cartography platform Resolve Biosciences. Insets (a and b) show the magnified image of areas identified in the main panel. Individual cell shapes can be determined from the transcript clouds. (F) Spatial location of individual VNO cells color-coded according to cell type prediction based on the spatial transcriptomic analysis. (G) Location of cell belonging to HBC, GBC, INP, and LP cell types, respectively. Heat indicates confidence of predicted values. BL: basal lamina; MZ: marginal zone.

Figure 1—figure supplement 1
Age differences in gene expression.

(A) Volcano plot of gene expression differences between P14 and P56 (Wilcoxon rank sum test, FDR ≤ 0.05). (B) 16 significantly differentially expressed genes with largest positive or negative log2-fold-change values (Wilcoxon rank sum test, FDR ≤ 0.05). (C) 50 significantly enriched GO terms with largest positive or negative log2-fold-change values (GSEA Permutation testing w/ FDR ≤ 0.05).

Figure 1—figure supplement 2
Sex differences in gene expression.

(A) Volcano plot of gene expression differences between male and female mice (Wilcoxon rank sum test, FDR ≤ 0.05). (B) 16 significantly differentially expressed chemosensory receptors with largest positive or negative log2-fold-change values (Wilcoxon rank sum test, FDR ≤ 0.05). (C) 35 significantly enriched GO terms (GSEA Permutation testing w/ FDR ≤ 0.05).

Figure 1—figure supplement 3
Zonal distribution of cell types in the VNO neuroepithelia.

(A) Schematic indicating quantification of cells in the marginal, intermediate, and main zones. (B) Stacked bar-plot of cell-type proportions by VNO zone (Wilcoxon rank sum test, FDR ≤ 0.05). (C) Box plots of GBC, INP, and immature VSN cell counts by zone, across 13 slides.

Figure 2 with 3 supplements
Novel neuronal lineage.

(A) UMAP visualization of cell-type clusters for the neuronal lineage. (B) Expression of Gnai2 and Gnao1 in the neuronal lineage. (C) Expression of Cnga2 and Gnal in the OSN lineage. (D–E) Location of mOSNs (D) and sVSNs (E) in a VNO slice. Color indicates prediction confidence. (F) Heatmap of normalized expression for a select set of mutually differentially expressed genes between sVSNs and mature V1Rs, V2Rs, and OSNs. (G) Enriched gene ontology (GO) terms in sVSNs when compared with V1R and V2R VSNs, respectively (GSEA Permutation testing w/ FDR ≤ 0.05). (H) Box plots of normalized expression for Muc2, Obp2a, Obp2b, and Lcn3 across mature sensory neurons (Wilcoxon rank sum test, FDR ≤ 0.05).

Figure 2—figure supplement 1
Neuronal lineage features and differentially expressed genes in sVSNs.

(A and B) Normalized Gap43 and Stmn2 expression in the neuronal lineage. (C) Normalized. Xist expression in the neuronal lineage, split by sex. (D) Violin plot of normalized Gnai2. expression in mature V2R neurons, split by sample. (E–G) Total number of genes, counts, and percent ribosomal gene expression detected in mature neurons, split by cell type (Wilcoxon rank sum test, FDR ≤0.05). (H) Heatmap of 503 significant (Wilcoxon rank sum test, FDR ≤0.001, fold-change ≥1.5) differentially expressed genes in sVSNs.

Figure 2—figure supplement 2
Significantly regulated genes in sVSNs.

(A-D) Feature plots showing odor binding protein gene expression (Muc2, Obp2a, Obp2b, and Lcn3). (E–L) Top upregulated genes in sVSNs. (M–T) Top downregulated genes in sVSNs.

Figure 2—figure supplement 3
sVSN Pseudotime, OSN markers and GO Terms.

(A) VNO neuronal lineage cell-types in UMAP space. (B) V1R, V2R, and sVSN lineage. pseudotimes in UMAP space. (C) Volcano plot of mOSN gene expression versus all other cells. from the neuronal lineage (Wilcoxon rank sum test, FDR ≤0.05). (D) 45 significant GO terms enriched in OSNs versus all other cells (GSEA Permutation testing w/ FDR ≤0.05).

Figure 3 with 1 supplement
Molecular cascades separating the neuronal lineages.

(A) PCA plot of V1R and V2R pseudotime principal curves across cell types. (B) Cell density plots across pseudotime for P14 and P56 mice (Two-sided Kolmogorov-Smirnov test). (C) Heatmap showing expression in pseudotime for genes that differentially expressed between the V1R and V2R lineages. Heat indicates Z-score values. (D) Zoomed in UMAP of cell types early in the neuronal lineage. (E–I) Feature plots for select genes expressed early in the neuronal lineage. (J) UMAP of INP, iVSN, iOSN, and mOSN cell types. (K–U) Normalized expression of candidate genes for V1R/V2R/OSN lineage determination. (V) A simplified model for lineage determination by transcription factors in VNO sensory neurons.

Figure 3—figure supplement 1
Genes expressed during INP to immature neuron transition.

(A and B) Violin plot of normalized expression for Neurog1 (A) and Neurod1 (B) split by cell type. (C–R) Feature plots of normalized gene expression for early differentially expressed genes between V1R, V2R, and OSN lineages. (S and T) Feature plots of normalized expression for Notch1 (S) and Dll4 (T).

Figure 4 with 2 supplements
Receptor expression in the VNO.

(A-D) Expression level of individual receptor (total raw counts) plot against the number of cells expressing the receptor for V1R (A), V2R (B), VNO-Olfr (C), and MOE-Olfr (D). (E-H) Ranked distribution of average expression per cell for the four receptor classes. Inset pie charts show the number of cells expressing a receptor at the specified range. (I and J) Heatmaps showing the Pearson correlation coefficient of transcription factor expression among the V1Rs (I) or V2Rs (J). (K) A simplified model of transcription factor selection in mVSNs.

Figure 4—figure supplement 1
Receptor distribution in the VNO.

(A) Log scale rank-frequency distribution plots for V1R, V2R, VSN-OR, and OSN-OR. (B–D) Number of cells expressing (nCells) vs. total raw counts for V1R, V2R, and VSN-OR. (E–F). Rank by average count distributions for V1R, V2R, and VSN-OR.

Figure 4—figure supplement 2
Co-expression between vomeronasal receptors and transcription factors.

(A) Distribution density plot showing relationship between similarity of transcription factor (TF) gene expression profiles and receptor sequence similarity. (B–C) Heatmaps showing the. V1R-TF (B) and V2R-TF (C) associations. Heat shows average expression level for TF genes for a given receptor type.

Figure 5 with 1 supplement
Co-expression among VNO receptor classes.

(A) Shannon Indices showing the specificity of receptor expression. High values indicate more co-expressions. (B–C) Prevalence and level of receptor co-expression by age and cell- type for the V1R (B) and V2R (C) lineages, respectively. (D–F) Circos plot of genomic loci for significantly co-expressed receptor pairs in the V1R (D) V2R, (E) and across-type populations (F). (G–K) Detection of receptor gene co-expression using Molecular Cartography. Individual dots represent single molecules. Colors represent different receptor genes. DAPI stain is shown as gray. Scale bar: 10 µm.

Figure 5—figure supplement 1
Receptor expression statistics.

(A-C) Percent of all read counts from a receptor gene (A), number of receptor species present. per cell (B), and total counts from vomeronasal receptors (C), separated by cell type (Wilcoxon rank sum test, FDR ≤0.05). (D–G) Proportions of first, second, and third most expressed receptor gene as percent of total receptor counts, separated by cell type.

Axon guidance molecules associated with receptor genes.

(A) Distribution density plot showing relationship between similarity of axon guidance (AG) gene expression profiles and receptor sequence similarity. The distribution of receptor similarity (x-mean and x-median) remains constant over the range of AG similarity. The AG similarity (y-mean and y-median) as a function of receptor similarity shows a strong correlation at the dense part of the curve. (B) Heatmap showing the Pearson correlation coefficient among VRs in their AG expression. (C and D) Heatmaps showing the V1R-AG (C) and V2R-AG (D) associations. Heat shows average expression level for AG genes for a given receptor type. (E) A simplified model of hierarchical distribution of AG in the mVSNs.

Transcriptional determinants of axon guidance molecules for individual receptor types.

(A) Correlation heatmap between receptor types calculated from co-expressed AG and TF genes. Receptor types are color-coded. (B) 3-D heatmap showing the Jaccard Indices between AG and TF genes for each of the VRs in the dataset. (C–G) Heatmaps showing Jaccard Indices of TF-AG associations for V1R (C and D) and V2R types (E–G). The lists of TF and AG genes here are abridged from the full list to enhance visualization. (C) These receptors share Meis2 expression but different AG genes. Note that Vmn1r185 and vmn1r69 both recognize female identify pheromones. (D) Similarity and distinction of AG/TF expression for three V1R types that are located in the same genomic location and with high sequence homology. (E and F) Shared TFs and AG genes for broadly (E) and uniquely (F) expressed V2R types. (G), distinct TFs and AGGs for uniquely expressed V2R types.

Tables

Key resources table
Reagent type (species) or resourceDesignationSource or referenceIdentifiersAdditional information
Strain, strain background (Mus musculus; 2 females, 2 males)C57BL/6In-house breeding
Strain, strain background (Mus musculus; 4 males, 4 females)CD-1In-house breeding
Commercial assay or kitChromium Next GEM Single Cell 3’ GEM, Library and Gel Bead Kit v3.110 X Genomics1000120
Commercial assay or kitChromiumSingle Cell 3’ GEM, Library & Gel Bead Kit v3.010 X Genomics1000075
Commercial assay or kitDNase I (RNase free)NEBM0303
Commercial assay or kitPapainCalbiochem5125
Commercial assay or kitDAPIThermo Fisher Scientific62247
Commercial assay or kitL-cysteineCalbiochem243005
Commercial assay or kitBSASigma-AldrichA8806
Commercial assay or kitFrozen section mediaLeica3801481
Commercial assay or kitDRAQ5Invitrogen65-0880-96
Commercial assay or kitHBSSVWRVWRL0121-0500
Commercial assay or kitPBSGibco10010023
Commercial assay or kitUrethaneSigma-AldrichU2500
Commercial assay or kitNovaSeq S1Illumina20012865
software, algorithmFiji ImageJGoldstein et al., 2018https://imagej.net/software/fiji/
Software, algorithmQuPath v0.4.3Bankhead et al., 2017https://qupath.github.io
Software, algorithmSeuratHao et al., 2021https://satijalab.org/seurat/
Software, algorithmkallisto | bustoolsMelsted et al., 2019https://www.kallistobus.tools
Software, algorithmDropletUtilsLun et al., 2019https://bioconductor.org/packages/release/bioc/html/DropletUtils.html
Software, algorithmIllustratorAdobehttps://www.adobe.com/illustrator
Software, algorithmSoupXYoung and Behjati, 2020https://cran.r-project.org/web/packages/SoupX/index.html
Software, algorithmclustreeZappia and Oshlack, 2018https://cran.r-project.org/web/packages/clustree/index.html
Software, algorithmggplot2Wickham et al., 2016https://cran.r-project.org/web/packages/ggplot2/index.html
Software, algorithmglmGamPoiAhlmann-Eltze and Huber, 2021https://bioconductor.org/packages/release/bioc/html/glmGamPoi.html
Software, algorithmveganOksanen et al., 2019https://cran.r-project.org/web/packages/vegan/index.html
Software, algorithmScrublet v0.2.3Wolock et al., 2019https://github.com/swolock/scrublet
Software, algorithmreticulateUshey et al., 2017https://cran.r-project.org/web/packages/reticulate/index.html
Software, algorithmGeneOverlapShen, 2019https://bioconductor.org/packages/release/bioc/html/GeneOverlap.html
Software, algorithmcirclizeGu et al., 2014https://cran.r-project.org/web/packages/circlize/index.html
Software, algorithmSlingshotStreet et al., 2018https://www.bioconductor.org/packages/release/bioc/html/slingshot.html
Software, algorithmtradeSeqVan den Berge et al., 2020https://www.bioconductor.org/packages/release/bioc/html/tradeSeq.html
Software, algorithmmsigdbrDolgalev, 2020https://cran.r-project.org/web/packages/msigdbr/vignettes/msigdbr-intro.html
Software, algorithmfgseaKorotkevich et al., 2021https://bioconductor.org/packages/release/bioc/html/fgsea.html
Software, algorithmMolecular Cartography Resolve Bioscienceshttps://resolvebiosciences.com/

Additional files

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

  1. Max Henry Hills
  2. Limei Ma
  3. Ai Fang
  4. Thelma Chiremba
  5. Seth Malloy
  6. Allison R Scott
  7. Anoja G Perera
  8. C Ron Yu
(2024)
Molecular, cellular, and developmental organization of the mouse vomeronasal organ at single cell resolution
eLife 13:RP97356.
https://doi.org/10.7554/eLife.97356.3