Author response: Variation in olfactory neuron repertoires is genetically controlled and environmentally modulated Article Swipe

Full text Figures and data Side by side Abstract eLife digest Introduction Results Discussion Materials and methods Data availability References Decision letter Author response Article and author information Metrics Abstract The mouse olfactory sensory neuron (OSN) repertoire is composed of 10 million cells and each expresses one olfactory receptor (OR) gene from a pool of over 1000. Thus, the nose is sub-stratified into more than a thousand OSN subtypes. Here, we employ and validate an RNA-sequencing-based method to quantify the abundance of all OSN subtypes in parallel, and investigate the genetic and environmental factors that contribute to neuronal diversity. We find that the OSN subtype distribution is stereotyped in genetically identical mice, but varies extensively between different strains. Further, we identify cis-acting genetic variation as the greatest component influencing OSN composition and demonstrate independence from OR function. However, we show that olfactory stimulation with particular odorants results in modulation of dozens of OSN subtypes in a subtle but reproducible, specific and time-dependent manner. Together, these mechanisms generate a highly individualized olfactory sensory system by promoting neuronal diversity. https://doi.org/10.7554/eLife.21476.001 eLife digest Smells are simply chemicals in the air that are recognized by nerves in our nose. Each nerve has a receptor that can identify a limited number of chemicals, and the nerve then relays this information to the brain. Animals have hundreds to thousands of different types of these nerves meaning that they can detect a wide array of smells. Smell receptors are proteins, and the genes that encode these proteins can be very different in two unrelated people. This could partly explain, for example, why some people find certain odors intense and unpleasant while others do not. However, having different genes for smell receptors does not by itself completely explain why some people are more sensitive than others to particular smells. The amounts of each nerve type in the nose might also differ between people and have an effect, but to date it has not been possible to accurately count them all. Ibarra-Soria et al. have now devised a new method to essentially count the number of each nerve type in the noses of mice from different breeds. The method makes use of a technique called RNA-sequencing, which can reveal which genes are active at any one time, and thus show how many nerves are producing each type of smell receptor. Ibarra-Soria et al. learned that different breeds of mice had remarkably different compositions of nerves in their noses. Further analysis revealed that this was due to changes to the DNA code near to the genes that encode the smell receptor. Next, Ibarra-Soria et al. sought to find out how the amount of each nerve type is controlled by giving mice water with different smells for weeks and looking how this affected their noses. These experiments revealed that a small number of the nerve types became more or less common after exposure to a smell. The altered nerves were directly involved in recognizing the smells, proving that the very act of smelling can change the make-up of nerves in a mouse's nose. These results confirm that the diversity in the nose of each individual is not only dictated by the types of receptors found in there, but also by the number of each nerve type. The next challenge is to understand better how these differences change the way people perceive smells. https://doi.org/10.7554/eLife.21476.002 Introduction Mapping the neuronal diversity within a brain remains a fundamental challenge of neuroscience. Quantifying variance in a population of neurons within and between individuals first requires precise discrimination of cellular subtypes, followed by an accurate method of counting them. While this has been achieved in a simple invertebrate model containing hundreds of neurons (White et al., 1986), applying the same approach to mammalian brains that encompass many millions of neurons represents a significant challenge (Wichterle et al., 2013). The main olfactory epithelium (MOE) is an essential component of the olfactory sensory system. It contains olfactory sensory neurons (OSNs) that express olfactory receptors (ORs), the proteins that bind odorants (Buck and Axel, 1991; Zhao et al., 1998). The mouse genome codes for over a thousand functional OR genes, but each mature OSN expresses only one abundantly, in a monoallelic fashion (Hanchate et al., 2015; Saraiva et al., 2016; Tan et al., 2015; Chess et al., 1994). This results in a highly heterogeneous repertoire of approximately 10 million OSNs (Kawagishi et al., 2014) within the nose of a mouse, stratified into more than a thousand functionally distinct subpopulations, each one characterized by the particular OR it expresses. The monogenic nature of OR expression serves as a molecular barcode for OSN subtype identity. Thus, the MOE offers a unique opportunity to generate a comprehensive neuronal map of a complex mammalian sensory organ, and investigate the mechanisms that influence its composition and maintenance. To date only a few studies have quantified the number of OSNs that express a given OR (Bressel et al., 2016; Fuss et al., 2007; Khan et al., 2011; Rodriguez-Gil et al., 2010; Royal and Key, 1999; Young et al., 2003). For the scarce data available (<10% of the full repertoire) reproducible differences in abundance have been observed between OSNs expressing different ORs (Bressel et al., 2016; Fuss et al., 2007; Khan et al., 2011; Young et al., 2003). This suggests variance in the representation of OSN subtypes exists within an individual, but the extent of variation between individuals is unknown. Moreover, the mechanisms that dictate the abundance of OSN subtypes are poorly understood. Most promoters of OR genes contain binding sites for Olf1/Ebf1 (O/E) and homeodomain (HD) transcription factors (Young et al., 2011), and these are involved in determining the probability with which the OR genes are chosen for expression (Rothman et al., 2005; Vassalli et al., 2011). Enhancer elements also regulate the gene choice frequencies of nearby, but not distally located, ORs (Khan et al., 2011). To date, these studies have focused only on a handful of OSN subtypes. In addition to OR gene choice regulation exerted by genetic elements, it is conceivable that the olfactory system adapts to the environment. The MOE is continually replacing its OSN pool and the birth of every neuron presents an opportunity to shape the proportion of different subpopulations. It is also possible that relative OSN abundances could be altered by regulating the lifespan of each OSN subtype. Indeed activation extends a sensory neuron's life-span (Santoro and Dulac, 2012), suggesting that persistent exposure to particular odorants may, over time, increase the relative proportions of the OSNs responsive to them. Some OSN subtypes do reportedly increase in number in response to odor activation, but others do not (Cadiou et al., 2014; Cavallin et al., 2010; Watt et al., 2004). Whether this variation reflects differences in the biology of OSN subtypes or experimental procedures is unclear. Here, we fully map OSN diversity in the MOE and characterize the influence of both genetic and environmental factors on its regulation. We show that RNA sequencing (RNAseq) is an accurate proxy for measuring the number of OSNs expressing a particular OR type, and use this approach to quantify, in parallel, the composition of 1115 OSN subtypes in the MOE. We report that, while the repertoire of OSN subtypes is stable across individuals from the same strain, it reproducibly and extensively differs between genetically divergent strains of laboratory mice. We show that under controlled environmental conditions, these stereotypic differences in OSN abundance are directed by genetic variation within regulatory elements of OR genes that predominantly act in cis and are independent of the function of the OR protein. However, we find that persistent, but not continuous, exposure to specific odorants can also subtly alter abundance of the OSN subtypes that are responsive to such stimuli. Taken together, these results show that the OSN repertoire is shaped by both genetic and environmental influences to generate a unique nose for each individual. Results Olfactory sensory neuron diversity measured by RNAseq Previously, we characterized the transcriptional profile of the whole olfactory mucosa (WOM) in adult C57BL/6J animals (hereafter termed B6) to generate hundreds of new, extended OR gene annotations (Ibarra-Soria et al., 2014). As each OR gene is expressed in only a small fraction of cells within WOM, differences in their abundance are difficult to distinguish from sampling bias. We hypothesized that mapping RNAseq data to significantly extended OR transcripts should increase detection sensitivity. With these models, OR gene mRNA level estimates in adult WOM increase, on average, 2.3-fold, but some increase almost 20-fold (Figure 1—figure supplement 1A). Despite this improvement, most OR mRNAs still have relatively low-expression values (Figure 1A). Nevertheless, they show a dynamic range of abundance levels (Figure 1A, inset) that are consistent between biological replicates, as indicated by their very high correlation values (median rho = 0.89, p<2.2 × 10−16). Figure 1 with 1 supplement see all Download asset Open asset RNAseq is highly sensitive for OR mRNA detection and provides a measurement of OSN diversity. (A) Barplot of the mean normalized expression of 1249 OR genes from six biological replicates, accounting for gene length. Genes are ordered by decreasing abundance. The horizontal line is the median expression (32.06) and all the genes below it are shown in the inset. (B) Mean normalized mRNA expression values for the OR genes in chromosome 9 of the Olfr7 cluster deletion mouse line (green; n = 3). The corresponding abundances in wild-type animals (orange) are shown as a mirror image (n = 3). The break on the x-axis separates the two OR clusters. The dotted box encloses the deleted ORs. (C) Unequal RNAseq expression levels for different OR genes can be explained by two scenarios: (left) an OR gene with high RNAseq levels is expressed by a larger number of OSNs than a gene with low RNAseq abundance; and/or (right) an OR with high RNAseq values is expressed in the same number of OSNs as one with low RNAseq values, but at higher levels per OSN. (D) Comparison of the number of OSNs that express nine OR genes assessed by in situ hybridization (ISH; x-axis) to the corresponding RNAseq values (y-axis). Error bars are the standard error of the mean (ISH n = 4, RNAseq n = 6). The line is the linear regression and the Spearman's correlation coefficient (rho) indicates a very strong correlation. Representative ISH images of two OR genes (in red) are shown. (E) In single-cell RNAseq experiments, 63 OSNs were randomly collected from the MOE. The distribution of OR mRNA expression in WOM samples is plotted (left), alongside the equivalent values for the ORs that were present in the picked single-OSNs (right). There is a significant enrichment (p<6.44 × 10−9) toward collecting OSNs that express OR genes with high RNAseq counts in WOM. (F) Comparison of the normalized expression value for the highest OR detected in each of the 63 single-OSNs (y-axis) to the corresponding mean value in WOM (x-axis, n = 3). The line is the linear regression and the Spearman's correlation coefficient (rho) indicates there is no correlation. See also Figure 1—figure supplement 1. https://doi.org/10.7554/eLife.21476.003 To assess whether these low OR mRNA expression values are biologically meaningful or if they represent low-level technical artifacts of RNAseq analysis, we sequenced RNA from WOM of a mouse strain that has a targeted homozygous deletion of the Olfr7 OR gene cluster on chromosome 9 (Xie et al., 2000; Khan et al., 2011), and compared their gene expression profile to control mice. From the 94 OR genes of the cluster that have been deleted, 83 (88.3%) have no counts in any of the three biological replicates. The 11 remaining genes have just one or two fragments mapped in only one of the replicates (Figure 1B), resulting in normalized counts of less than 0.4. In contrast, the control mice have from 14.2 to 498.1 normalized counts for the same genes. Together these experiments demonstrate that the use of extended gene models significantly increases the sensitivity to detect OR mRNA expression in WOM, and that the full dynamic range of abundances reflects true measures of OR gene expression. The wide range of stereotypic OR gene expression can be explained by two scenarios, acting alone or in combination (Figure 1C): either (1) OR genes with high-expression values are monogenically expressed in more OSNs than those with low-expression values; and/or (2) OR genes are consistently expressed at different levels per OSN. To differentiate between these possibilities, we performed in situ hybridization (ISH) of probes specific to nine OR genes with expression values distributed across the dynamic range. We then counted the number of OSNs in which each OR is expressed (Figure 1D). We find very strong correlation between OSN number and RNAseq expression value (rho = 0.98, p=5 × 10−5). We additionally compared OR gene RNAseq expression levels with three independent measures of the number of OSNs expressing the same ORs (Bressel et al., 2016; Fuss et al., 2007; Khan et al., 2011). In all three cases, we find high correlations (Figure 1—figure supplement 1B–D). We next collected 63 single mature OSNs from WOM, and determined the OR gene most abundantly expressed in each using a single-cell RNAseq approach (Saraiva et al., 2016, unpublished data). If OR expression levels in WOM reflect the proportion of OSNs that express each receptor (Figure 1C), the probability of isolating each OSN type is not equal. Indeed, we find a strong selection bias towards OSNs that express OR genes with high RNAseq levels in WOM (hypergeometric test, p=6.44 × 10−9; Figure 1E), suggesting those OSN types are more numerous in the olfactory epithelium. Thus, consistent with a recent analysis in zebra fish (Saraiva et al., 2015), OR RNAseq values are an accurate measure of the number of each OSN subtype in the mouse WOM (scenario 1). But do consistent differences in OR mRNA levels per cell also contribute (scenario 2)? To test this, we quantified the mRNA levels of the most abundant OR gene in each of the 63 single, mature OSNs, normalized to three stably expressed OSN marker genes (Khan et al., 2013). We find OR mRNA levels vary within the single cells, but this does not correlate with expression levels across the WOM (rho = −0.04, p=0.7518) (Figure 1F). Analysis of ERCC spike-ins confirmed that the levels of OR mRNAs in single OSNs are reliable. Moreover, the single OSN transcript levels also positively correlate with transcript levels in pools of millions of OSNs (Saraiva et al., 2016). Together, these data demonstrate that OR mRNA levels obtained by RNAseq are an accurate proxy for quantifying the diversity of OSN subtypes that express each receptor. The OSN repertoire differs between strains of mouse The relative proportion of each OSN subtype is stable between genetically identical animals. We have previously reported the expression of OR genes in B6 male and female mice (Ibarra-Soria et al., 2014). By applying full gene models to these data, here we confirm their OSN distribution profiles are remarkably similar (Figure 2A); only 1.2% of the OR gene repertoire is significantly differentially expressed (Figure 2B). To investigate whether this OSN distribution is a stereotypic feature of the species, we next reconstructed the WOM transcriptome of a different laboratory strain, 129S5SvEv (hereafter termed 129). The 129 genome has 4.4 million single nucleotide polymorphisms (SNPs) and 0.81 million small indels compared to B6 (Keane et al., 2011), of which we find 13,484 SNPs and 1,936 indels within our extended OR gene transcripts. As OR genes are particularly variable in coding sequence between strains of mice (Logan, 2014), mapping RNAseq data from other strains to a B6 reference genome results in biases in OR gene expression estimates (Figure 2—figure supplement 1A). We therefore generated a pseudo-129 genome on which to map the RNAseq data, by editing the reference genome at all polymorphic sites. We confirmed that the RNAseq expression estimates correlate with the number of OSNs that express the corresponding receptor genes in 129 animals, as judged by in situ hybridization (rho = 1, p=5.5 × 10−6; Figure 2—figure supplement 1B). From the 1,249 OR genes, we find 462 are significantly differentially expressed (DE) compared to B6 (false discovery rate (FDR) < 5%), representing 37% of the whole repertoire (Figure 2C,D). Figure 2 with 1 supplement see all Download asset Open asset OSN diversity varies between mouse strains. (A) Mirrored barplot of the mean normalized RNAseq expression values for the OR genes in male (dark blue, top) and female (light blue, bottom) B6 animals (n = 3). (B) Scatter plot for the same data, with the Spearman's correlation (rho) indicating a strong correlation. The red line is the 1:1 diagonal. Significantly differentially expressed (DE) OR genes are represented in blue, non-DE genes are in black. (C) Same as in (A) but with the average B6 expression in blue (both males and females, n = 6) compared to the corresponding 129 expression values in yellow (n = 3). (D) Corresponding scatter plot, with the significant DE genes in green. (E) Same as in (A) but comparing the B6 expression in blue (n = 6) to the CAST abundances in red (n = 3). (F) Corresponding scatter plot, with DE genes in purple. (G) Venn diagram illustrating the intersection of DE OR genes between the pairwise comparisons of the three strains. (H) An example of an OR gene, Olfr6, that is DE in all strain comparisons. (I) Representative in situ hybridizations (ISH) on coronal slices of B6 and 129 MOEs for two OR genes, Olfr31 and Olfr736, that are DE between these strains. (J) The quantification of OSNs expressing each OR gene in B6 (blue) and 129 (yellow) are plotted alongside the corresponding RNAseq counts. The log2 fold-changes between the strains are indicated. (K) Fold-change between the strains obtained from ISH data (x-axis) or RNAseq counts (y-axis) for four DE OR genes expressed higher in B6 (blue), four expressed higher in 129 (yellow) and one expressed at equivalent levels in both strains (grey); these include Olfr31 and Olfr736. The line is the linear regression and the Spearman's coefficient (rho) indicates a strong correlation between OSN and RNAseq counts. See also Figure 2—figure supplement 1. https://doi.org/10.7554/eLife.21476.005 Figure 2—source data 1 OR expression data in three strains of mice. Excel workbook containing the normalized expression data for all OR genes in B6, 129 and CAST, along with the results of the differential expression Download Figure 2—source data 2 OR in the CAST with the of the in the CAST Download To whether genetic diversity influences the variance in OSN we this using a strain from the This strain has more than million SNPs and million indels relative to B6 (Keane et al., of we counted that SNPs and indels are found within our extended OR transcripts. mapping to a genome (Figure 2—figure supplement 1C), OR genes are significantly differentially expressed < compared to B6, of the whole OR repertoire (Figure The changes in expression for some OR genes are genes have differences of at all pairwise comparisons into 129 CAST, Figure 2—figure supplement OR genes are DE between at two strains. and of these are DE in all three pairwise comparisons (Figure for example, there are consistently different of OSNs in each strain (Figure To if the DE OR genes between strains reflect differences in the proportions of OSN subtypes, we performed ISH with probes specific to OR genes with significantly different expression values between B6 and 129 (Figure We then counted the number of OSNs that express nine different OR in each strain, and compared this with their RNAseq expression values (Figure We find a high correlation between the in OSN number and the in RNAseq expression values between B6 and 129 (rho = 0.98, p=5 × Figure our approach accurately measures the in OSN between strains. OR gene are in number between individual et al., 2007; Young et al., and mouse strain et al., Thus, it is possible that variance in OSN subtype are a of different of highly similar OR genes between strains. To assess this, we CAST genome sequence data (Keane et al., for SNPs within OR genes. We ORs that contain 10 or more an of genome sequencing data from these genes, we or OR genes. We the CAST RNAseq data to a genome these new OR and the expression of the OR The abundance profile remains for genes (Figure 2—figure supplement 1F). To assess whether this for the observed differential expression between we compared these estimates to of OR genes their DE and OR genes now DE (Figure 2—figure supplement Thus, while differences in OR gene number contribute to the diversity in OSN repertoire between three strains of mice, other mechanisms are for most of the OSN diversity independent of odor divergent mouse strains different in their et al., et al., and et al., each strain of mouse, in is to a unique and olfactory environment. As odor exposure the life-span of OSNs in an et al., and Dulac, Watt et al., genetic variation could regulate OSN population either directly or devised an to test and differentiate the influence of the olfactory from the genetic We four to B6 and 129 to to they an identical in environment. B6 were to B6 and 129 to 129 In B6 a single 129 and 129 a single B6 each a olfactory but one had a different genetic from the others (Figure 10 weeks of we quantified the OSN of six animals and six using We found that the OSN cluster in two by genetic (Figure The correlation coefficient for any two B6 samples was on average with no significant between the In contrast, the correlation for any B6 with a 129 had a mean of 0.89, which is significantly × and OR genes, other genes are DE between these mice by strain (Figure supplement 1A). In contrast, across the whole transcriptome we find only mRNA from two genes that show differences in expression to odor both of which encode ORs (Figure Figure supplement 1B). These data demonstrate that the diversity in OSN repertoire we between strains is almost dictated by genetic In a controlled the influence of on the and of the MOE is to only a few OSN subtypes. Figure with 1 supplement see all Download asset Open asset OSN diversity is determined by the genetic and not by the olfactory environment. (A) to differentiate genetic from environmental influences on OSN diversity. B6 (blue) and 129 (yellow) as were into the were to B6 and 129 Each B6 one 129 and 10 the WOM was collected for RNAseq from three from each strain, and one (B) of the expression of the OR genes in all sequenced animals cluster by the genetic of the animals. The strain and of each mouse is indicated (right). (C) expression revealed mRNA from only two genes, and that are significantly altered on the olfactory environment. values are shown for each and yellow B6 or 129 animals and the indicates the olfactory environment. See also Figure supplement 1. OSN diversity profiles are independent of OR function and are controlled in cis The of the OSN repertoire to the olfactory suggests its and is not by the specific of OR proteins by their coding To test this, we the OSN repertoire of We identify the of OSN subtypes across a dynamic range of abundance (Figure The differential proportions of OSNs expressing particular OR genes are therefore present the of the suggesting that it is not on the of the OSNs on differences in OSN Figure with 1 supplement see all Download asset Open asset OSN diversity is independent of OR and is controlled in (A) Mean normalized expression of the OR mRNA in the WOM of B6 animals, from most to abundant (n = 3). (B) Mean normalized mRNA expression of OR in the B6 adult WOM (n = 6). (C) genetically mouse line was that contains the coding sequence of in of The the use of to on either side of the and a DNA containing the along with for (D) Mirrored barplot of the mean normalized mRNA expression values for the OR repertoire in B6 animals (light blue, n = and in homozygous (dark blue, n = the most abundant OR and is no expressed in the genetically (E) Scatter plot of the mean normalized counts (x-axis) of OR genes the log2 between and n = OR genes that are significantly differentially expressed are represented in and are the of the repertoire is equivalent or very (F) Comparison of the of the CAST B6 OR expression (y-axis) to the between the CAST and B6 in the The genes on the 1:1 indicating the mRNA expression observed in the is in the and thus OR abundance is controlled in The correlation coefficient is which the correlation between the two estimates while for on the (G) of the normalized mRNA expression in the strains of an OR gene that is more abundant in CAST or in B6 The corresponding mRNA abundance of each in the is The log2 is indicated for each Error bars are the standard error of the See also Figure supplement 1. Figure data 1 OR expression data in the Excel workbook containing the normalized expression data for all OR genes in the mouse line and B6 along with the results of the