The emergence of transcriptional identity in somatosensory neurons

Abstract

More than twelve morphologically and physiologically distinct subtypes of primary somatosensory neuron report salient features of our internal and external environments1-4. It is unclear how specialized gene expression programs emerge during development to endow these subtypes with their unique properties. To assess the developmental progression of transcriptional maturation of each subtype of principal somatosensory neuron, we generated a transcriptomic atlas of cells traversing the primary somatosensory neuron lineage in mice. Here we show that somatosensory neurogenesis gives rise to neurons in a transcriptionally unspecialized state, characterized by co-expression of transcription factors that become restricted to select subtypes as development proceeds. Single-cell transcriptomic analyses of sensory neurons from mutant mice lacking transcription factors suggest that these broad-to-restricted transcription factors coordinate subtype-specific gene expression programs in subtypes in which their expression is maintained. We also show that neuronal targets are involved in this process; disruption of the prototypic target-derived neurotrophic factor NGF leads to aberrant subtype-restricted patterns of transcription factor expression. Our findings support a model in which cues that emanate from intermediate and final target fields promote neuronal diversification in part by transitioning cells from a transcriptionally unspecialized state to transcriptionally distinct subtypes by modulating the selection of subtype-restricted transcription factors.

Conflict of interest statement

Competing interests

The authors declare no competing financial interests

Figures

Extended Data Figure 1.. Quality control metrics for DRG sensory neuron scRNA-seq data and canonical correlation analysis.

a-e. Distribution of the number of genes discovered in each cell (individual points) in each population of sensory neuron (underlying violin plot) in Adult (A), Postnatal day 5 (B), Postnatal day 0 (C), Embryonic day 15.5 (D), Embryonic day 12.5 (E). Individual cells with <1000 genes (considered to be low quality) or >10000 genes (considered to be likely doublets) discovered were eliminated from subsequent analysis. Individual cells with <2000 UMIs (considered to be low quality) were excluded from subsequent analysis. UMI: unique molecular identifier f. Integration of Adult/P5 (1st plot), P5/P0 (2nd plot), P0/E15.5 (3rd plot), E15.5/E12.5 (4th plot) using canonical correlation analysis to find common sources of variation between timepoints. Single cells are labeled as individual points, with color representing identified cell types and gray representing cells in the preceding timepoint.

Extended Data Figure 2.. Somatosensory neuron subtype composition varies across axial levels.

a. (Left) Schematic representing which axial levels were quantified. (Right) Quantification of single molecule RNA-FISH to determine the percentage of C6/7, T7/8, and L4/5 DRG neurons that corresponds to each transcriptionally defined somatosensory neuron subtype. Black dotted lines are used to highlight the subtypes present at different percentages at different axial levels. b. Example images of single molecule RNA-FISH for transcriptionally distinct somatosensory neuron subtypes in C6/7 (top row), T7/8 (middle row) and L4/5 (bottom row) DRGs.

Extended Data Figure 3.. Dorsal root ganglia and trigeminal ganglia are constituted by similar somatosensory neuron subtypes.

a. tSNE visualization of trigeminal ganglia scRNA-seq data obtained from Adult (P28–42) mice. Colors denote principle cell types and dotted circles were added to aid in visualization of principal cell types. LTMR/proprioceptor specific gene expression overlaid onto t-SNE visualization of mature DRG sensory neurons. b. Distribution of the number of genes discovered in each population of sensory neuron in adult trigeminal ganglia displayed as violin plot. UMI: unique molecular identifier. c. Heatmap depicting expression of the genes enriched in somatosensory neuron subtypes resident in the dorsal root ganglia as well as their expression levels in cognate trigeminal ganglia subtype counterparts. d. Heatmap depicting expression of the genes enriched in somatosensory neuron subtypes resident in the trigeminal as well as their expression levels in cognate dorsal root ganglia subtype counterparts. For c,d boxes represent IQR, whiskers represent minimum and maximum values, and notches represent the 95% confidence interval of the median. TPT: tags per ten thousand. * denotes two sided Wilcoxon rank-sum test with Bonferroni corrected p < 0.01. TPT: tags per ten thousand.

Extended Data Figure 4.. Neural crest progenitors, sensory neuron progenitors and unspecialized sensory neurons express highly distinct gene programs.

a. Heatmap depicting cell cycle (S/G2/M) associated genes for the principal subtypes identified at E11.5. b. Heatmap depicting expression levels of unspecialized sensory neurons enriched genes in both in mature somatosensory neuron subtypes and unspecialized sensory neurons. c. (Left) Heatmap depicting expression of the genes enriched in unspecialized sensory neurons (USN) as well as their expression levels in neural crest progenitors (NCP) and sensory neuron progenitors (SNP). (Right) Violin and box plots depicting example genes enriched in USNs. d. Heatmap depicting expression of the genes enriched in neural crest progenitors (NCP) as well as their expression levels in sensory neuron progenitors (SNP) and unspecialized sensory neurons (USN). (Right) Violin and box plots depicting example genes enriched in NCPs. e. Heatmap depicting expression of the genes enriched in sensory neuron progenitors (SNP) as well as their expression levels in neural crest progenitors (NCP) and sensory neuron progenitors (SNP). (Right) Violin and box plots depicting example genes enriched in SNPs. For a-e boxes represent IQR, whiskers represent minimum and maximum values, and notches represent the 95% confidence interval of the median. TPT: tags per ten thousand. * denotes two sided Wilcoxon rank-sum test with Bonferroni corrected p < 0.01.

Extended Data Figure 5.. Force directed layout of putative subtype-restricted transcription factors.

a. Force directed layout representation of DRG displaying with expression patterns displayed for the remaining putative subtype-restricted transcription factors. b. tSNE visualization of Runx1, Runx3, Pou4f2 and Pou4f3 expression in the adult DRG. TPT: tags per ten thousand. c. (Left) Single molecule RNA FISH for Runx1 and Runx3 in E11.5, P0 or adult DRGs. For E11.5, the spinal cord and DRG are labeled as references. (Right) Single molecule RNA FISH for Pou4f2 and Pou4f3 in E11.5, P0 or adult DRGs. For E11.5, the spinal cord and DRG are labeled as references. (Bottom) Quantification of the RNA-FISH. * represents two-sided t-test p<0.01.

Extended Data Figure 6.. The expression of somatosensory neuron subtype specific genes during development.

A. Box plots representing subtype specific genes at E12.5, E15.5, P0, P5 and Adult (P28–42) for each identified somatosensory neuron subtype. Boxes represent IQR, whiskers represent minimum and maximum values, and notches represent the 95% confidence interval of the median. TPT: tags per ten thousand. * denotes two sided Wilcoxon rank-sum test with Bonferroni corrected p < 0.01. B. Normalized line plots displaying what percentage of adult-levels of subtype specific gene expression are detected at E12.5, E15.5, P0, and P5. The black line represents the median of each time point with Adult being defined as ‘100%’. Upper and lower bands represent a 95% confidence interval (defined as ± 1.87 * IQR/√n, where n=sample size, IQR: Interquartile range)

Extended Data Figure 7.. DRG counts and TF analysis in Pou4f2 and Pou4f3 mutants.

a. Representative images of Avil smRNA-FISH from T7/8 DRGs in Pou4f2KO(Cre)/WT (left) or Pou4f2KO(Cre)/KO(Cre) (right) littermate control DRGs. Quantification of estimated cell count per DRG presented to the right of the images. b. Representative images of Avil smRNA-FISH from T7/8 DRGs in Pou4f3WT/WT (left) or Pou4f3KO/KO (right) littermate control DRGs. Quantification of estimated cell count per DRG presented to the right of the images. c. Box plots displaying the expression levels of subtype-restricted TFs in each somatosensory neuron subtype in Pou4f2WT/WT (left column) or Pou4f2KO(Cre)/KO(Cre) (right column) littermates. d. Box plots displaying the expression levels of subtype-restricted TFs in each somatosensory neuron subtype in Pou4f3WT/WT (left column) or Pou4f3KO/KO (right column) littermates. For c,d boxes represent IQR, whiskers represent minimum and maximum values, and notches represent the 95% confidence interval of the median. TPT: tags per ten thousand.

Extended Data Figure 8.. Generation and validation of Bmpr1bT2a-Cre and Avpr1aT2a-Cre mouse lines.

a. Targeting strategy for inserting a T2a-Cre-TGASTOP codon; Frt-PGK:NeoR-pA-Frt cassette immediately upstream of the stop codon in the Bmpr1b gene. b. Single molecule RNA-FISH for both Bmpr1b and GFP in Bmpr1bT2aCre ; AAV-CAG:FLEX-GFPP14 I.V. mice to confirm the specificity and utility of the Bmpr1bT2a-Cre allele. c. Targeting strategy for inserting a T2a-Cre-TGASTOP codon; Frt-PGK:NeoR-pA-Frt cassette immediately upstream of the stop codon in the Avpr1a gene. d. Single molecule RNA-FISH for both Avpr1a and tdTomato in Avpr1aT2a-Cre(ΔNeo) ;Rosa26 LSL-tdTomato/WT mice to confirm the specificity and utility of the Avpr1aT2-aCre allele. e. (Top left) tSNE representation of transcriptionally mature DRG overlaying the expression pattern of Avpr1a, (Remaining images) Representative immunostaining images of tdTomato and CGRP in skin sections obtained from Avpr1aT2a-Cre; Rosa26LSL-tdTomato animals. TPT: tags per ten thousand. f. (Top left) tSNE representation of transcriptionally mature DRG overlaying the expression pattern of Bmpr1b. (Remaining images) Representative immunostaining images of tdTomato and CGRP in skin sections obtained from Bmpr1bT2a-Cre; AAV-CAG:FLEX-GFPP14 I.V. animals. TPT: tags per ten thousand. g. Representative immunostaining images of GFP in skin sections obtained from Pou4f2KO(Cre);AAV-CAG:FLEX-GFPP14 I.V animals. h. Quantification of ending morphology for CGRP-α (Avpr1aT2a-Cr);Rosa26LSL-tdTomato) and CGRP-η (Bmpr1bT2a-Cre;AAV-CAG:FLEX-GFPP14 I.V.) somatosensory neuron subtypes. i. Schematic representation of the skin with the distinct morphological ending types of CGRP-α and CGRP-η neurons displayed. j. Representative images of CGRP immunostaining in skin samples of 2–3 week old Pou4f3WT/WT (left) or Pou4f3KO/KO (right) littermate controls. Statistical comparisons were done using a two-tailed t-test. * represents p<0.01 k. Representative images of GFP immunostatining in skin samples of 3–4 week old Pou4f2KO(Cre)/WT (top left) or Pou4f2KO(Cre)/KO(Cre) (right) littermate controls and representative RNA-FISH for GFP in Pou4f2KO(Cre)/WT and Pou4f2KO(Cre)/KO(Cre) littermate controls are displayed below the skin immunostaining images. In h * represents two-way ANOVA with a Tukey’s HSD post-hoc analysis p < 0.01. In j,k * represents two-sided t-test p<0.01Bar graphs in h,j,k show mean +/- s.e.m.

Extended Data Figure 9.. Postnatal depletion of Pou4f3 results in a loss of subtype specific gene expression in CGRP-α and CGRP-η □eurons.

a. Quantitative RT-PCR analysis using cDNA generated from animals transduced with Luciferase or Pou4f3shRNA expressing AAVs. Error bars represent mean +/- standard error of the mean. Statistical comparisons were done using a paired two-sample t-test. * represents p < 0.01. b. Distribution of the number of UMIs discovered in each population of control sensory neuron. UMI: unique molecular identifier. c. Distribution of the number of UMIs counted in each each population of shRNA sensory neuron. UMI: unique molecular identifier d. t-SNE visualizations of scRNA-seq data for neurons generated from LuciferaseshRNA (Left) and Pou4f3shRNA littermates DRGs (Right). e. Boxplots showing the fold-change distribution of cell-type specific gene expression in sensory neuron subtypes expressing the highest and lowest levels of Pou4f3 when comparing LuciferaseshRNA and Pou4f3shRNA littermate control samples. f. Boxplots showing the fold-change distribution in a randomized and expression matched control gene set g. Control or Pou4f3-depleted mice were exposed to two surfaces with the indicated temperature (x-axis), and the percentage of time spent in the 30°C chamber over a 5-minute test period. h. Representative images of CGRP immunostaining in skin samples of 1–2 week old LuciferaseshRNA (left) or Pou4f3shRNA (right) littermate controls. In a,g,h mean +/- s.e.m displayed. In a * represents two-sided t-test p<0.01. In g * represents two-way ANOVA with a Tukey’s HSD post-hoc analysis p < 0.05.

Extended Data Figure 10.. Subtype restricted TF expression profiles in NGF−/−; Bax−/− cell clusters

a. Heatmap depicting expression of the subtype-restricted TFs in P0 somatosensory subtypes (left) or the clusters from NGF−/−; Bax−/− mutants (Right). b. Single molecule RNA FISH for pairs of subtype restricted TFs in NGF−/−;Bax−/− (Top) or littermate NGF−/−; Bax−/− mutants (Bottom). c. Quantification of the RNA-FISH data showing the number of Pou4f3/Shox2 double-positive, Pou4f3/Hopx double-positive, Bcl11a/Hopx double-positive, Neurod1 single-positive or Neurod6 single-positive neurons. d. Schematized model of gene expression programs as cells traverse development milestones. Transcriptionally unspecialized sensory neurons that emerge from Sox10+ and Neurogenin+ progenitors co-express multiple TFs which become restricted to select subtype as neurons mature. These TFs are responsible for establishing the transcriptional specializations found in each neuronal subtype. In c * represents wo-sided t-test p<0.01. TPT: tags per ten thousand.

Figure 1.. scRNA-seq of developing and mature DRG sensory neurons.

a. t-SNE visualizations DRG scRNA-seq data. b. UMAP visualization of DRG scRNA-seq data from E11.5 with developmental trajectory and gene expression information overlaid. TPT: tags per ten thousand. c. Quantification of tdTomato+ neurons and representative image. Mean +/- s.e.m. is indicated. d. Heatmap and quantification of genes enriched in each somatosensory neuron subtype as well as their expression levels in unspecialized sensory neurons. USN: unspecialized sensory neuron. Boxes represent IQR, whiskers represent minimum and maximum values, and notches represent the 95% confidence interval of the median. TPT: tags per ten thousand. * denotes two-sided Wilcoxon rank-sum test with Bonferroni corrected p < 0.0001.

Figure 2.. Transcriptional development of DRG neuron subtypes.

a. Force-directed layout of DRG sensory neurons overlaid with time point or cell type information. b. Force-directed layout of DRG development overlaid with expression of indicated genes. c. Heatmap of subtype-restricted TFs in each somatosensory neuron subtype of adult ganglia. TPT: tags per ten thousand.

Figure 3.. Refinement of TF expression in developing somatosensory neurons.

a. Developmental trajectory of sensory neurons (left) and tSNE visualization with TF expression overlaid. TPT: tags per ten thousand. b. Schematic description of the AvilCreERT2; Rosa26LSL-tdTomato labeling strategy. c. tSNE visualization of AvilCreERT2; Rosa26LSL-tdTomato scRNA-seq with cell type identity or tdTomato expression overlaid. TPT: tags per ten thousand. d. Schematic representing strategy for labeling neurons with Pou4f2Cre/WT mice. e. smRNA-FISH and quantification for the indicated transcripts. Mean +/- s.e.m is indicated. * represents two-tailed t-test. p<0.01.

Figure 4.. Pou4f2 and Pou4f3 regulate select somatosensory neuron subtype identities.

a. t-SNE visualizations of scRNA-seq data for neurons generated from Pou4f3WT/WT and Pou4f3KO/KO littermates. b. t-SNE visualizations of scRNA-seq data for neurons generated from Pou4f2WT/WT and Pou4f2KO/KO littermates. c. Fold-change distribution of cell-type specific genes when comparing Pou4f3KO/KO and Pou4f3WT/WT littermates control samples. d. Fold-change distribution of cell-type specific genes when comparing Pou4f2KO/KO and Pou4f2WT/WT littermates. In c,d Boxes represent IQR, whiskers represent minimum and maximum values, and notches represent the 95% confidence interval of the median and * denotes two sided Wilcoxon rank-sum test with Bonferroni corrected p < 0.01. TPT: tags per ten thousand.

Figure 5.. The extrinsic cue NGF is required for subtype specific gene expression and TF expression patterns.

a. t-SNE visualizations of scRNA-seq data for neurons generated from P0 Bax−/− and NGF−/−; Bax−/− littermates. b. Cell-type specific gene expression in proprioceptor and A-fiber mechanoreceptor sensory neuron subtypes in Bax−/− and NGF−/−; Bax−/− littermates. c. Cell-type specific gene expression in all other sensory neuron subtypes in Bax−/− and NGF−/−; Bax−/− littermates, Boxes represent IQR, whiskers represent minimum and maximum values, and notches represent the 95% confidence interval of the median. In b,c all clusters are different from controls by two sided Wilcoxon rank-sum test with Bonferroni corrected p < 0.01. TPT: tags per ten thousand.