Deep sequencing of the murine olfactory receptor neuron transcriptome.

The ability of animals to sense and differentiate among thousands of odorants relies on a large set of olfactory receptors (OR) and a multitude of accessory proteins within the olfactory epithelium (OE). ORs and related signaling mechanisms have been the subject of intensive studies over the past years, but our knowledge regarding olfactory processing remains limited. The recent development of next generation sequencing (NGS) techniques encouraged us to assess the transcriptome of the murine OE. We analyzed RNA from OEs of female and male adult mice and from fluorescence-activated cell sorting (FACS)-sorted olfactory receptor neurons (ORNs) obtained from transgenic OMP-GFP mice. The Illumina RNA-Seq protocol was utilized to generate up to 86 million reads per transcriptome. In OE samples, nearly all OR and trace amine-associated receptor (TAAR) genes involved in the perception of volatile amines were detectably expressed. Other genes known to participate in olfactory signaling pathways were among the 200 genes with the highest expression levels in the OE. To identify OE-specific genes, we compared olfactory neuron expression profiles with RNA-Seq transcriptome data from different murine tissues. By analyzing different transcript classes, we detected the expression of non-olfactory GPCRs in ORNs and established an expression ranking for GPCRs detected in the OE. We also identified other previously undescribed membrane proteins as potential new players in olfaction. The quantitative and comprehensive transcriptome data provide a virtually complete catalogue of genes expressed in the OE and present a useful tool to uncover candidate genes involved in, for example, olfactory signaling, OR trafficking and recycling, and proliferation.

Introduction

The sense of smell has been perfected to detect a tremendous range of different volatile chemical substances. Olfactory receptors (OR) constitute the largest superfamily of mammalian G protein-coupled receptors (GPCRs) and account for the vast discriminatory power of the olfactory system. In humans, this discriminatory power is achieved by the relatively small number of approximately 400 functional ORs, whereas mice possess a considerably larger number of OR genes with approximately 900 functional ORs [1-6]. Despite their lower number of functional genes, humans can discriminate more than 10,000 different odors. Presumably, this discrimination is achieved by a combinatorial code, in which single receptor types are able to respond to several different odorant molecules, and single odorant compounds are recognized by a number of different receptor types [7]. According to the classical olfactory signal transduction scheme, receptor-ligand interactions in the ciliary membranes of the olfactory receptor neurons (ORNs) lead to activation of olfactory-specific G-proteins (Gαolf) [8], stimulation of adenylyl cyclase type III (ACIII), production of cyclic adenosin monophosphate (cAMP) [9], and activation of cyclic nucleotide-gated (CNG) channels that are composed of CNGA2, CNGA4 and CNGB1b [10-12]. The initial depolarization produced by the CNG channels is enhanced by subsequent activation of the calcium activated chloride channel (CaCC) ANO2 [13-18]. Focused research on this model of olfactory perception has provided more detail, and the participation of many more components has been reported. The action of the Na+/K+/Cl--cotransporter Nkcc1, most likely together with Slc4a1 [17] and other transporters, provides the high intracellular chloride concentration necessary for signal amplification [19, 20]. The trafficking of ORs is accomplished by receptor transporting proteins (RTP) [21]. After OR activation, the signal termination and internalization processes are achieved through G protein-coupled receptor kinase (GRK) [22-24]and possibly β-arrestins [25-27]. ORs and accessory proteins bind to the scaffolding multi PDZ-domain proteins (MUPP) to form an assembled transduction complex [28]. The efficacy of the signal transduction cascade is enhanced by the nucleotide exchange factor Ric8b [29-31]. In contrast, signal adaptation is achieved through reduced ACIII activity mediated by calmodulin (CaM) kinases [32], degradation of cAMP by specific phosphodiesterases (PDEs) [33-37], and desensitization of CNG-channels by binding of CaM [38-41].

Stephan et al. (2012) showed that the principal Na+/Ca2+ exchanger Nckx4 is involved in rapid response termination and proper adaptation [42]. However, adaptation processes are complex, and there are many more proteins involved in the regulation of this mechanism.

Further, an inhibitory signal transduction cascade that involves phosphoinositide 3-kinase (PI3K) has been proposed [43-47]. Although the basic mechanisms involved in signal processing are known, many issues require further exploration before understanding of the complex olfactory perception process on molecular level will be achieved. Hence, it is desirable to identify the plenary repertoire of ORN-specific proteins to fully uncover the factors involved in olfactory signaling pathways.

In recent publications [48, 49], the transcriptome of the murine OE was characterized by microarray analysis. According to this DNA-array study, more than 10,000 genes are expressed in ORNs. ORN-specific genes were identified by the differential analysis of OMP-positive versus OMP-negative neurons using a transgenic OMP-GFP mouse [50], and the specific expression of 300 genes in ORNs was verified by in situ hybridization (ISH). Most recently, next generation sequencing (NGS)-based OE transcriptome data focused on sex-specific differences in OR genes were generated for BALB/c mice [51], and the study by Keydar et al. (2013) [52] generated a more general catalogue of genes based on NGS from data of C57BL/6J mice.

Additionally, transcriptome analyses are complemented by proteomic studies of OE plasma and ciliary membranes [53–56, 16, 15] that have led to the detection of approximately 2,346 proteins. Nevertheless, the detection of OR gene expression, which is fundamental for an olfactory tissue, was incomplete because membrane proteins are difficult to detect in such proteome studies, and in some of these studies, ORs failed to be detected at all.

Here, we present the first NGS based transcriptome of fluorescence-activated cell sorting- (FACS)-sorted ORNs in combination with a comprehensive transcriptome study of the murine OE. Our data provide a comprehensive list of transcripts for membrane proteins that included established and previously described proteins and new, previously unrecognized proteins. We present a detailed expression ranking for GPCRs and additionally analyze the expression patterns of several newly identified membrane proteins in the OE.

Compared with microarray [48, 49] and proteome data [53–56, 16, 15] from the OE, deep sequencing of murine the OE transcriptome detects more expressed genes and nearly all known ORs and allows better quantification of expression levels, as verified by quantitative PCR.

Our data provide an important new and sensitive tool to guide novel approaches and advance research on olfactory signaling mechanisms in the OE and, especially, ORNs.

Materials and Methods

Animals

CD1 and C57BL/6J mice were obtained from Charles River (Sulzfeld, Germany) and transgenic OMP-GFP mice [50] were kindly provided by Dr. Peter Mombaerts, Max Planck Institute of Biophysics, Frankfurt.

Mice were offered normal laboratory chow and water ad libitum in standard cages. All animal experiments were carried out in accordance with the European Union Community Council guidelines and approved by the competent state office of the Federal Land of Northrhine Westphalia (file number 87–51.04.2010.A180).

Preparation of mice OE and FACS of ORNs

OE of 4 weeks old male and female CD1 or C57BL/6J mice was prepared and subsequently RNA was isolated; each RNA sample (male and female CD1, female C57BL/6J) was prepared from an OE pool from eight different mice for CD1 mice; in case of C57BL/6J three mice were used for OE preparation. ORNs were obtained from homo- or heterozygous OMP-GFP mice [50]. OE was prepared and collected in 300 µl cold Ringer’s solution. Epithelia were minced and ORNs dissociated by adding papain followed by an incubation of 15 minutes at 37°C. Reaction was stopped by washing with Ringer’s solution; cell suspension was centrifuged 10 minutes at 121 g and cell pellet was resuspended in phosphate buffered saline. Enrichment of ORNs from cell suspension was done by FACS. Altogether, six (homozygous) or eight (heterozygous) adult mice were used to obtain about 100,000 sorted, fluorescent ORNs.

RNA isolation

Total RNA was isolated either out of the pooled OE or of the sorted ORN sample with the RNeasy Mini Kit (Qiagen, Hilden, Germany) according to the manufacturer’s protocol including the optional on-column DNaseI digestion.

NGS library preparation and Illumina sequencing

Libraries for NGS sequencing were prepared from total RNA and subjected to DSN normalization by standard Illumina protocols. Afterwards, Illumina sequencing was performed on a Genome Analyzer IIx using single end 36-bp (CD1 OE male and female) respectively 75-bp (ORNs homozygous) reads. RNA for heterozygous ORNs was sequenced on a HiSeq-2000 by standard Illumina protocols (101-bp, paired end). In case of C57BL/6J mice and additional probes of CD1 OE, mRNA was sequenced on HiSeq-2000 (101-bp, paired end).

Alignment of RNA-Seq reads using TopHat

We analyzed the raw sequence data in fastq format as previously described [57]. RNA-Seq reads were aligned to version mm9 of mouse reference genome and transcriptome using TopHat (v2.0.7) [58] which utilizes the ultra-fast short–read mapping program Bowtie to arrange the alignment [59]. TopHat output files in BAM format were sorted and indexed with SAMtools [60]. In order to reduce the alignment of repetitive reads a multiread-correction was used allowing up to 5 hits per read.

Gene expression calculation using Cufflinks

Aligned RNA-Seq reads for each sample were assembled into transcripts and their abundance was estimated by the program Cufflinks (v1.3.3) [61] using the RefSeq mm9 reference transcriptome in Gene Transfer Format (GTF) obtained from the UCSC Genome Bioinformatics database (University of California Santa Cruz). In order to estimate transcript expression, the GTF-file was supplied to Cufflinks. The parameter—compatible-hits-norm was set to ensure that FPKM normalization was performed based on reference transcriptome only.

Cufflinks was provided with a multifasta file (mm9.fa) to improve accuracy of the relatively transcript abundance estimation [62]. We further used a masked command –M and the mask GTF rmsk.gtf to mask all possible reads from RNA repeats (including tRNA, snRNA, scRNA, srpRNA) short and long interspersed nuclear elements (SINE, LINE) and other classes of repeats.

Cufflinks indicates and quantifies the relative abundances of transcripts in the unit FPKM [61]. For comparison to olfactory tissue, we reanalyzed already published raw RNA-Seq data from brain, liver, muscle and testes [63, 64] in the same manner as our own data. The data sets were visualized and investigated by the Integrative Genomic Viewer (www.broadinstitue.org/igv) for proving sequence alignements and correct mapping of reads for the top expressed genes. While the raw data analysis was performed on a Linux based computer further calculations were carried out with Microsoft Excel 2010. In this study, our intention was to monitor the expression of protein coding genes only. Therefore, we removed all entries for non-polyadenylated transcripts from our analyzed data including microRNA (miRNA) and small nucleolar RNA (snoRNAs). This also improves the comparability to other data sets such as mRNA-Seq versus total-RNA-Seq. DSN-normalization and different types of RNA preparation lead to differences especially in such small RNA species.

For a differential gene expression analysis, we used the program Cuffdiff, which identifies significant changes in transcript expression between two datasets [57].

Availability of raw data sets. The raw RNA-Seq data sets (FASTQ file format) for FACS-sorted ORNs and OE were deposited in Gene Expression Omnibus (GEO) repository (www.ncbi.nlm.nih.gov/geo/) under the following accession number: GSE53793.

Quantitative PCR

Mice were decapitated and OE was collected from septum und turbinates. Total RNA from single mouse samples was isolated with RNeasy Mini Kit (Qiagen, Hilden, Germany) including optional DNaseI digestion and cDNA was prepared using iScript cDNA Synthesis Kit (Bio-Rad, München, Germany). PCR reactions (1 min, 95°C; 1 min, 63°C; 1 min, 72°C; 35 cycles) were performed on an iQ5 thermal cycler using iQ SYBR Green Supermix (both Bio-Rad, München, Germany). At least three independent runs with three technical replicates each were performed and expression levels were calculated using the ddCt method.

Since OR genes contain no introns in their coding sequence, we used genomic DNA as efficiency control and normalized OR cDNA data to OR genomic DNA data.

Primers for housekeeping gene α1b-tubulin were used as described [65], others were designed with Primer-BLAST (sequence 5’ to 3’):

a1b-tubulin long rv (AGCAATGGCTGTGGTGTTGCTCA)

ACIII (ACGACCACAAGCGCTTTCAGG; ACTTGGAAGGCACCAGGGGCA)

mOR-EG short (TGACAGGTTTGTGGCCATCCGC; CAGTATCAAGGAGCATACCACCCCC)

mOR-EG long (TCTCCTCTCTCCTTTCACTTTCACGC; AGGCTTTGCGGCGTCCACTT)

Olfr2 (GGAGCGAAGGAACCACACTGGG; TCAGCACCAACACGTAGGCCA)

Olfr31 (TGCTCGTCTCACATGACTGTGGT; CCCAGTGCCCTCCTCAGACCT)

Olfr169 (TCCTTTTCCACCTGTTCCTTCCACA; TGCCAAAAATTTATCCTGCCCTGGA)

Olfr259 (TCAGAACACAAGTCTACGAGGCTGT; TGCCACATAGCGATCATAGGCCA)

Olfr309 (CTGACTTTCTGTGGGCCGAA; CTCCCACAGCCAGGGTGATA)

Olfr314 (ACCCCCATGCTGAACCCCCT; GGTCCAAGCCAGCCAAGGCA)

Olfr355 (TCGGCCTGAGGACCAAAAACCA; AGGGTCAGAACGGATGGCCAGA)

Olfr411 (CATGCGGCTGGGGTCAGTGG; ATGCCAAGTCTCGCTTTCCTGTGA)

Olfr525 (AGCTGCTGCTGCTCACTGCC; TGCACTGATGACCCACACACTGC)

Olfr533 (GTGACCCCATCTCTGAACCC; GCTTTTAGAGTTCATCTCCAAGCTC)

Olfr545 (GGCTCGCACTGCCTCCCAA; TCAGGGGCCTGCACTCAGGA)

Olfr632 (GGGCCCCGTGTGGCATTGAT; GCGGGTGTCAGAACAGGCCA)

Olfr705 (CCTGGCACTGACACTGGGTGG; TGTGGCCACCATAAGCCAGCA)

Olfr1126 (TGTGGTGAGCTGTGCCACACAA; AGGCTGCTGCCAGCTGAGTG)

Olfr1301 (TGCAACTGCCTCTCTGTGGTCAT; GCTCTGGCGAACAGTGAGAAGGA)

Olfr1347 (TTGACCAACATGACCAGAGTCCAAC; TCTCCAGGAGGGTCAGCAGGT)

Olfr1348 (ATGGACGGGGCCCTACCACA; TGTAGTGGGCGGCAGATGGC)

Olfr1349 (CACCAGCCTGGACTTCGGCC; GGCTGCAGCAGGGTCACCAG)

Olfr1395 (TCCTGGTGGCCCTCTCTGGC; GGAGCCACGGCCCAGGAGAA)

Olfr1507 (GGGCCATGTGGACAGCAGGG; AGGCCAGTTCGATCACCTGAGGT)

Olfr1508 (TGGCTTGCACTGACACCCACA; ACCACGGTGAGATGGGCTGC)

Adipor1 (CGGCAGACAAGAGCAGGAGTGT; AGAACCAGCCCATCTGGCCCA)

Gpr108 (TGAGTCCCGTGAGGAGGGTGC; TGCAGATGCCGAATGGACCAGA)

Gpr177 (AAAGGGGGTCGCAGAAATGGCTG; GGGTGCTGGAGCGATCAAGCCT)

Stomatin (CAGCAACCCGTCTTTTGGCAC; TCCCCCAGTCATCCGTGGCAT)

Stoml2 (TGGACTGCCCCGAAACACCG; CCGGATTCGGTCTAACACGGGG)

Stoml3 (GCCACGTTCCTGCTGGCTCAG; ATTCCCCAGAGCTCGGTGGCA)

RTP1 (GTACCTGGAACTTCATGCCTCAGGC; CACGCATGCGCACCGAACCAG)

RTP2 (GCAGTACCTGGAGCAACATGCCT; CGCATGCGCACTGAGCCAACA)

RTP3 (CGTTTTGAAGCCAGGATGGACGC; GGCCTTCTTCTCACACCACCGC)

RTP4 (GCACCACCTGGCTCGAACCC; GCATCACCCAGTCTGCCACAGT)

REEP1 (ATGGTGCTCCTGCTCCCTCGG; AGGGTTCTAGGCGGTGCCCG)

REEP2 (AACAGCGGCTCCTGCTGCCT; GGCTGGGTACAGGGTGCCAA)

REEP3 (CGCTGAGGCGATCCCAGAGCA; TCTTGGGAAGGGCAAGGTGTAGACA)

REEP4 (CTGAGTGAGCGTCCCTCCGC; ACGGCCTTGTAGGAAGCATACGC)

Tmem30a (CATCGTCGCTACGTGAAATCTCGA; CGCCACATGGCGCAATTGGT)

Tmem30b (GTTCCGCTCCGACGGTCTGC; GATCGAGCGGCACCTCGACG)

Tmem30c (TTTGCAGGAAGTTCGAAGCCTCTGC; AAAGGCAGCCGTCCGCATCC)

Tmem205 (GAGCTGCCACGGTTCTTGGTCC; CCAGGCACCAGACAAGACCAGC)

Tmem16b (GACGCCCAGATGCTGTGAGGA; AGGCGGGGGATAAAGTCGGAGG)

In situ hybridization

Digoxigenin-labeled antisense riboprobes, typically about 600–800 nucleotide length, were generated from cDNA of interest by in vitro transcription performed with the DIG RNA labeling mix (Roche, Palo Alto, CA) according to the manufacturer’s instructions. Coronal Cross sections of OE were obtained from wild type C57BL/6J mice of both sexes aged postnatal 21 days (P21). Mice were sacrificed using CO2 followed by decapitation. Subsequently the maxillary and anterior cranial region of the head was dissected free. The whole head was fixed in 4% paraformaldehyde in PBS at 4°C overnight, followed by a 24 hours decalcification step with 0.5 M ethylen-diamin-tetra-acetat (EDTA) (pH 8.0) and a final cryoprotecting incubation step with an increasing 10%–30% sucrose in PBS series (10% and 20% sucrose for 1 hour; 30% sucrose for at least 3 hours) at 4°C. Afterwards, coronal sections (12 µm) of quickly frozen heads embedded in tissue freezing medium OCT for supporting tissue during cryotomy (Leica Microsystems, Bensheim, Germany) were cut on a cryostat (Leica, CM 3050S) and mounted on Superfrost Plus Slides (Thermo Scientific, Menzel Gläser). After dehydration using an increasing ethanol series, slices were stored at -80°C until further use.

In situ hybridizations were performed as described with minor modifications [66]. Briefly, fixed cryosections were incubated in RIPA-buffer, followed by an acetylation step with acetic anhydride in TEA buffer. Next, a prehybridization step in 50% deionized formamide, 10% dextran sulfate, 5x Denhardt’s solution, 5x SSC, 10 mM DTT, 250 µg/ml yeast tRNA, 500 µg/ml sheared and denatured herring sperm, 50 µg/ml heparin, 2.5 mM EDTA, and 0.1% (v/v) Tween-20 was carried out for 1 h at 55°C to prevent unspecific binding of riboprobes. After each incubation step a wash step with SSC- or PBS—/T-buffers followed.

Finally, 50 ng antisense riboprobes were hybridized at 55–65°C on cyrosections mounted on slides overnight. The hybridized mRNA was visualized using an alkaline phosphatase-conjugated antibody to digoxigenin and hydrolysis of nitro-blue tetrazolium chloride/5-bromo-4-chloro-3-indolylphosphate p-toluidine. Antisense and a control sense probe were tested in parallel. The slides were covered with cover slips using polyvinyl-alcohol containing embedding medium (Mowiol, Immo-Mount). Digital images were obtained with Axiocam camera on an Axioscope2 microscope (Zeiss, Oberkochen, Germany). All images of sense and antisense samples were recorded under same conditions (brightness, contrast and light exposure time).

A signal was considered a positive when the antisense labeling was visually noticeably darker than the sense labeling.

Immunohistochemistry

Immunostainings of the OE were performed on 12 µm coronal cryosections of paraformaldehyde-fixed tissue of OMP-GFP mice [50] on Superfrost slides (Thermo Scientific, Menzel Gläser). After blocking with 1% gelatin or 10% goat serum in phosphate-buffered saline / 0.1% Triton X-100, sections were incubated with primary antibody (TRPC1: 5 µg/ml and TRPM7: 2 µg/ml; kindly provided by D.E. Clapham) and fluorescence-labeled secondary antibody (dilutions 1:500 or 1:1000) in phosphate-buffered saline / 1% gelatin / 0.1% Triton X-100. Stained sections were mounted in ProLong Antifade Gold medium (Molecular Probes). All fluorescence images were collected on a confocal laser scanning microscope (LSM510 Meta, Zeiss, Oberkochen, Germany). ImageJ (NIH) was used for image processing.

Literature search

In order to identify novel genes in terms of OE, an intense Pubmed (http://www.ncbi.nlm.nih.gov/pubmed) literature search was required. Thereby, the query contained the name of gene interested in combined with the terms “olfactory” or “olfaction”.

Results

To assess gene expression, the Illumina RNA-Seq protocol was used to amplify and sequence up to 57 million fragments from OEs of male and female CD1 adult mice. Furthermore, 13 million reads were generated for FACS-sorted ORNs obtained from homozygous and 58 million reads from heterozygous transgenic OMP-GFP mice [50] with a C57BL/6 background (Table 1, S1 Table). These OMP-GFP mice are determined to be valid models for assessing ORN gene expression patterns in a study using microarrays [48, 49].

Table 1

Summary of RNA-Seq data.

Tissue/Cells	Mouse			mRNA-Seq	Total number of
	strain	genotype	Gender		Sequenced reads	Mapped reads	[%]
Sorted ORNs	C57BL6 transgene OMP-GFP	homozygous	mixed	Illumina/Genome Analyser II, 72 bp reads, single-end	13,574,176	8,857,058	65.25
Sorted ORNs	C57BL6 transgene OMP-GFP	heterozygous	mixed	Illumina/HiSeq2000, 101 bp reads, paired end	58,234,129	34,519,988	59.28
OE1	CD1	WT	male	RNA-Seq (DNS-Normalization) Illumina/Genome Analyser II, 36 bp reads, single-end	37,284,951	28,881,184	77.46
OE2	CD1	WT	male	Illumina/Genome Analyser II, 101 bp reads, paired-end	40,053,764	28,363,505	70.81
OE3	CD1	WT	male	Illumina/Genome Analyser II, 101 bp reads, paired-end	36,141,244	26,397,003	73.04
OE4	CD1	WT	male	Illumina/Genome Analyser II, 101 bp reads, paired-end	11,394,752	8,384,689	73.58
OE1	CD1	WT	female	RNA-Seq (DNS-Normalization)Illumina/Genome Analyser II, 36 bp reads, single-end	52,716,231	38,700,603	73.41
OE2	CD1	WT	female	Illumina/HiSeq 2000, 101 bp reads, paired-end	54,589,506	36,103,429	66.14
OE3	CD1	WT	female	Illumina/HiSeq 2000, 101 bp reads, paired-end	57,808,346	39,268,107	67.93
OE4	CD1	WT	female	Illumina/HiSeq 2000, 101 bp reads, paired-end	33,354,414	24,884,683	74.61
OE5	CD1	WT	female	Illumina/HiSeq 2000, 101 bp reads, paired-end	19,807,752	14,513,434	73.27
OE1	C57BL6	WT	female	Illumina/HiSeq 2000, 101 bp reads, paired-end	86,404,176	63,862,406	73.91
OE2	C57BL6	WT	female	Illumina/HiSeq 2000, 101 bp reads, paired-end	30,763,538	22,694,013	73.77
OE3	C57BL6	WT	female	Illumina/HiSeq 2000, 101 bp reads, paired-end	25,520,854	17,789,702	69,7
OE4	C57BL6	WT	female	Illumina/HiSeq 2000, 101 bp reads, paired-end	25,520,328	18,831,564	73.79

Additionally, we sequenced approximately up to 86 million 101 nt length fragments from the OEs of female C57BL/6 mice (Table 1, S1 Table).

We analyzed our sequencing results using TopHat and Cufflinks software [57], and reads were mapped onto the mouse reference genome mm9. Expression values were calculated from the number of reads per kilobase per million reads in each sample. The quantitative measurement of gene expression employed in the RNA-Seq experiments was the FPKM (fragments per kilobase of exon per million fragments mapped) value.

Our analysis enabled us to detect the expression of 13,955 protein-coding genes in male and 13,280 in female OEs and up to 13,103 genes in FACS-sorted ORNs (FPKM > 1) of the 21,550 genes interrogated with a Refseq-based gene model. Including genes with an extremely low expression level (FPKM < 1), we counted a total of up to 16,786 genes (S1 Table).

To quantify olfactory gene expression relative to other tissues, we reanalyzed the raw sequencing data of previously published murine RNA-Seq NGS transcriptome studies that focused on the brain, muscle, liver and testes [63, 64] in the same manner with which we analyzed our own sequencing data (S3 Table, S4 Table). To ensure the comparability of datasets, we studied expression patterns of different housekeeping genes. A general, rough scale regards FPKM values ~1 to indicate weakly expressed genes, ~10 indicates medium expression, and ~100 indicates highly expressed genes based on comparisons to housekeeping genes (S1 Fig).

Olfactory Receptor Genes

To analyze the transcriptome in more detail, we first focused on chemoreceptor gene families starting with OR genes. The high sensitivity of RNA-Seq was demonstrated by the detection of nearly all ORs. Of 1,125 annotated OR genes, we detected 1,075 in female and 1,001 in male CD1 mice and 1,044 in female C57BL/6J mice with FPKMs > 0.1 (Fig. 1, S5 Table, S2 Fig.).

Figure 1

OR detection using RNA-Seq: Bar chart showing the percentage of detected OR genes in female and male OE tissue and FACS-sorted ORNs.

Percentages were calculated based on the 1,125 OR genes annotated in the Refseq based gene model.

Regarding FACS-sorted ORNs, we detected the expression of 582 (homozygous) and 905 (heterozygous) OR genes with FPKM > 0.5. For a typical 1 kb OR coding sequence at a FPKM of 0.1, the expression was confirmed by 3 reads in the OE, and due to the lower sequencing depths for FACS-sorted ORNs at a FPKM of 0.5 (S3 Fig.). As Shiao et al. (2012), we set a similar threshold for the detection of ORs [51].

If very weakly expressed OR genes, whose expression was only supported by 1–2 fragments, were included, expression of an additional ~2% of ORs was detected in both male and female OEs. A total of 98 and 18 ORs had no detectable expression in male and female OEs, respectively.

In comparison, RNA-Seq results from mouse brains, livers and muscles revealed greatly reduced numbers of detected OR genes (3—10 ORs). In contrast, pronounced expression of 155 ORs was detected in the testes (S4 Fig.).

The strength of OR gene expression is strongly dependent on the receptor gene and reached up to 97 FPKM for Olfr1507 in the ORNs. However, most OR genes (up to 73%) were expressed at low levels (FPKMs < 3) due to their mosaic-like expression patterns (S5 Fig.) in the OE. The median expression was 2 in the OE and 2.7 FPKM in the FACS-sorted ORNs.

The expression levels of individual ORs was strongly conserved between the male and female OEs (Pearson correlation coefficient, r = 0.83), and the three most highly expressed receptors (Olfr533, Olfr1507 and Olfr309) were identical in the males and females. Furthermore, the three most highly expressed ORs in the sorted ORNs were found in the top ten most highly expressed genes in the male and female OEs (Fig. 2).

Figure 2

Correlations of expression levels plotted for each detected OR gene.

A. Correlation of the OR gene expression patterns between male and female CD1 mice. Only OR genes with detectable expression levels (FPKM>0.1) are shown. The FPKM values are logarithmically presented. The Pearson correlation coefficient of r = 0.83 confirmed the strong correlation of OR gene expression patterns between female and male CD1 mice. The three ORs (Olfr533, Olfr1507 and Olfr309) with the highest expression levels were also the most highly expressed in the RNA-Seq data from both sexes. B. Correlation of OR gene expression patterns between females of strain CD1 and C57BL6. The Pearson correlation coefficient of r = 0.75 confirmed the strong correlation between the expression patterns of OR genes between the different strains; however, these patterns exhibited greater divergence between strains than between the sexes of the CD1 strain. The most highly expressed ORs, Olfr533 and Olfr309, had the same expression ranking, and Olfr1507 was among the ten most highly expressed OR genes in both strains. C. Verification of RNA-Seq results for ORs by real time RT-PCR. Expression levels are relative to mOREG (Olfr73). Error bars represent the SEM.

These finding suggest that the expression levels of specific receptor types are not only conserved between sexes but also between different mouse strains, as OEs from CD1 mice and sorted ORNs from transgenic OMP-GFP mice with a C57BL/6J background were examined.

Furthermore, we detected the expression of 21 out of 24 annotated OR pseudo-genes in the OE. Among these pseudo-genes, Olfr1372-ps exhibited the highest expression level (S6 Table). Because previous studies have established that certain OR pseudo-genes are functional [67], Olfr1372-ps is an interesting candidate for future studies.

We additionally verified the RNA-Seq results for several highly expressed OR genes with real time RT-PCR. The expression levels were comparable to that of mOR-EG (Olfr73), a well-characterized mouse OR [68].

We detected similar expression rankings for the 20 tested ORs compared to our RNA-Seq data and confirmed that Olfr1507, Olfr533 and Olfr309 belong to the highly expressed OR genes (Fig. 2).

Correlation Analysis

Whole data set

We generated for each condition (OE from CD1 male, CD1 female and C57BL16) up to five biological replicates (Table 1, S1 Table). Correlation measurements for protein-coding gene expression patterns between biological replicates showed that Pearson correlation coefficient values ranged for CD1 male from r = 0.43 to 0.96 (rmean = 0.6), CD1 female from r = 0.35 to 0.91(rmean = 0.67) and C57BL/6 female from r = 0.37 to 0.91 (rmean = 0.61) (S6 Fig.).

Analyzing genes with the most diverging expression pattern between the biological replicates, interestingly, we found predominantly transcripts for odorant-binding proteins (OBPs) (Obp, Mup, Lcn and several other gene families) [51] (S1 Table, S7 Fig.). OBPs function as soluble transport proteins of the nasal mucus [69].

In addition, the expression pattern of several cytochrome P450 (Cyp) genes, which are known to be expressed in the olfactory mucosa [70], showed a high variability in expression level between the OE samples (S1 Table). The presence of these enzymes in the olfactory organ is probably necessary for the metabolic transformation of odorant molecules and they play a protective role by detoxifying inhaled chemicals [70]. It is therefore plausible that mice individually adapt and regulate their olfactory transport system to their appropriate surrounding olfactory environment.

In summary, the genes with the highest variability of expression between individual mice belong to proteins with a major olfactory function. After excluding these genes from the correlation analysis, the biological replicates strongly correlate up to r = 0.97. Even data between CD1 and C57BL6 strains show now a strong correlation (Pearson correlation coefficient values up to r = 0.95) (S6 Fig.). Therefore, these results bode for the reproducibility of our data and represent a good basis for comparative analysis.

OR subgenome

To ensure reproducibility of the conserved pattern for OR gene expression, we analyzed the additionally generated replicates for each condition (OE of CD1 male, CD1 and C57BL6 female) (S8 Fig., S9 Fig.).

By comparing the highly expressed ORs in each replicate, our analysis confirms that the ORs Olfr533, Olfr1507 and Olfr309 indeed belong to the top expressed OR genes in all OE samples (S10 Fig.).

Further, the expression levels for individual ORs within the biological replicates were strongly conserved and reached Pearson correlation coefficient values up to r = 0.95) (S9 Fig.). OR expression levels were comparable between female CD1 and C57BL/6J mice; the Pearson correlation coefficient were r = 0.65 to 0.86 between replicates (S9 Fig.)

The variance of OR gene expression pattern between replicates within one strain was lower than between strains (S9 Fig.). This shows that the OR gene expression is widely conserved in the OE, whereby the expression pattern of a small set of OR genes is strain specific.

Other Chemoreceptors

Next, we analyzed other classes of chemoreceptor genes. Our analysis emphasized the specific expression of the “olfactory” trace amine-associated receptors (TAARs) in the OE [71]. All TAARs, with the exception of TAAR1, known to be specifically expressed in the brain, were exclusively detected in OE and FACS-sorted ORNs (Fig. 3). We also detected several vomeronasal organ (VNO) receptors (V1R and V2R) in the OE. Out of the 207 annotated V1R receptors, we observed weak expression of 15 genes with FPKMs of 0.1–0.5 in the OE, and 24 genes with FPKM values of up to 2 in sorted ORNs. Regarding V2R receptors, of the 127 annotated V2Rs we detected 5 genes with low expression levels (0.1–0.7 FPKM) in the OE and up to 21 genes with FPKM values from 0.5 to 8 in ORNs (Fig. 3). Analyses of these expression patterns in C57BL/6J mice yielded basically the same results.

Figure 3

Expression pattern of chemoreceptor genes.

Heatmap showing the expression levels of the following chemoreceptor classes: TAARs, VNO receptors, GC-D, taste receptors and FPR in olfactory (male and female OE, and FACS-sorted ORNs) and non-olfactory tissue (brain, muscle, liver and testes). Higher FPKM values are indicated by deeper colors. Only genes with a FPKMs >1 are represented in this chart.

In contrast, the family of formylpeptide receptors (FPR), which function as chemoreceptors in the VNO [72, 73], were largely undetected. Only weak expression of FPR1 and FPR2 genes was noticed (Fig. 3). Likewise, guanylyl cyclase GC-D receptor [74] expression was low (FPKM 0.8) (Fig. 3).

Interestingly, we also detected the expression of taste receptors in the OE. Our data from sorted ORNs revealed an enrichment of Tas1r1 and Tas1r3, (FPKM 6.7), these receptors typically form heterodimers and function as the umami receptor [75]. Type 2 bitter taste receptors (Tas2r) were weakly detected in OE at a FPKM of 0.7 for Tas2r108 (Fig. 3).

Genes with Known OE-Specific Expression

Next, we analyzed whether our RNA-Seq results matched the reported expression patterns of molecules linked to olfactory signaling cascades. Due to the extreme mosaic-like pattern of OR gene expression, FPKM values for OR expression cannot be directly compared to the values for other genes. Therefore, we calculated the accumulated gene expression of all ORs and assumed that this represented the expression level of a single OR in a single ORN.

Under these assumptions, the FPKM value was approximately on average 4,000 (n = 13) which would make the ORs by far the most highly expressed protein coding genes. The other most highly expressed genes (Fig. 4) included OMP, which is one of the most abundant proteins in the OE, followed by the subunits of trimeric G-protein composed of Gαolf, Gβ1 and Gγ13, and, with a somewhat lower expression rate, Go and Gß2. Moreover, the guanine nucleotide exchange factor Ric8b and the receptor interacting proteins RTP1, RTP2 and Reep1 were also among the most highly expressed proteins (Fig. 4). Other classical signaling molecules enriched in ORNs include the following: ACIII (which is the only adenylyl cyclase that we found to be highly enriched in the ORNs), the CNG channel (formed by the CNGA2, -A4 and -B1 subunits [76]), and the calcium-activated chloride channel Ano2 [16, 77]. We also detected ßArr2, which is responsible for OR internalization [78, 27], and a weak expression of βArr1. Among the 1,000 most highly expressed genes that have been implicated in olfactory signaling were Nkcc1 (also called Slc12a2) [14, 19], 6 PKC genes (α, β, δ, ε, η, and ζ), of which only α and β are known to be expressed in mouse ORNs, 5 GRKs (Grk2, Grk3, Grk4, Grk5, Grk6) [79–81, 22–25], and voltage-gated sodium channels (S1 Table) of which Scn9a was the primary and most highly enriched isoform in ORNs[82].

Figure 4

Expression pattern of known genes of the olfactory signal transduction.

A. Heatmap showing the expression levels of genes known to be involved in olfactory signaling and other genes known to be highly expressed in ORNs as determined by NGS. The FPKM of OR total represents the accumulated gene expression of all ORs and shows that the OR is the most highly expressed gene in the OE. The main components of the signal transduction scheme were among the 200 most highly expressed genes. The FPKM values shown for OMP in sorted ORNs (homozygous) are rough estimations based on the calculation of reads located in the 3’-untranslated regions of OMP and are therefore only valid to a limited extent. Higher FPKM values are indicated by deeper colors.

B. RT-PCR verification of the highly expressed genes in the OE. Gene expression was normalized to the level of adenylyl cyclase type III (ACIII) RNA. The Investigated genes were RTP1–4, REEP1–4, Stom, Stoml2–3 and transmembrane proteins, including ANO2 (Tmem16b). Error bars represent the SEM.

The expression patterns determined by NGS were in good agreement with previous reports about ORN gene expression and confirmed that the key proteins in olfactory signal transductions are strongly expressed in ORNs.

We also verified the RNA-Seq results for several highly expressed genes by real time RT-PCR using ACIII as a standard (Fig. 4). As expected, RTP genes showed the highest expression levels for RTP1 and RTP2. Stomatin-like protein 3 (Stoml3) is known to be expressed in ORNs [83] and showed the highest expression level of the stomatin family genes tested. The expression levels of the tested genes correlated with the FPKM values determined by RNA-Seq.

The top expressed genes in FACS-sorted ORNs from homo- and heterozygous OMP-GFP mice have similar expression levels. Therefore, we assume that a lack of OMP expression has no general influence on the expression of other genes. We calculate a Pearson coefficient of r = 0.9 for gene expression patterns between these two groups, which suggests a strong correlation of gene expression pattern. Additionally, using Cuffdiff analysis, we identified only few statistical significant changes in gene expression between both datasets. Our analysis revealed that only 13 genes were differentially expressed (S2 Table). Therefore, we assume that the data of FACS-sorted homozygous ORNs represent a nearly normal transcriptome of ORNs. These results are in accordance with the study of Sammeta et al. (2007)[48], in which no statistically significant differences in mRNA abundance between these two genotypes were detected by microarray technique. However, the limited sequencing depth of our homozygous ORNs (~13 million) complicates the detection especially of lower differentially expressed genes with a statistical significance, so that more regulated genes will be probably detected at a higher sequencing depth. The OMP expression level in FACS-sorted ORNs with a heterozygous genotype was comparable to the level in the OE datasets.

Differences in Gene Expression Patterns between OE and ORNS

The OE is composed of several different cell types. Next to the ORNs are sustentacular cells, basal cells, including globose and horizontal cells, microvillar cells, and cells lining the Bowman’s glands and duct are found in the OE [84, 85]. Accordingly, the differential expression pattern of the OE compared to FACS-sorted ORNs revealed a catalog of genes expressed in ORNs and/or other cell types of the OE. To avoid strain specific differences, we only compared the ORN data with C57BL/6J transcriptome data. We found that, in nearly all instances, genes known to be expressed in ORNs had 2–8 times higher FPKM values in sorted ORNs than in the OE, while genes known to be expressed in other cell types of the OE had 2–119-fold higher FPKMs in the OE compared to the sorted ORNs (Fig. 5).

Figure 5

Differences in the gene expression patterns between ORNs and the OE.

Comparison of FACS-sorted ORNs and the OE (CD1 OE both sexes; C57BL6 female OE) revealed that genes that are known to be expressed in mature ORNs were expressed in ORNs at levels that were about two to three-fold higher those of the OE. Genes specific for non-neuronal cell types were expressed at levels that were at least two to 119-fold greater in the OE.

To examine the expression patterns of non-neuronal cell types, we analyzed the expression of the marker genes reported by Sammeta et al. (2007) [48]. FPKMs of the sustentacular cell markers Cbr2 and Pax6 were higher in the OE but were also detected in ORNs [66], which suggests a small proportion of these cells were present in our FACS-sorted ORN probe. Expression of Krt14 and Krt5, markers of horizontal basal cells, was either weak or absent in the FACS-sorted ORN probe, and stronger expression of these markers was detected in the OE.

Reg3g, a marker of respiratory epithelium, was strongly expressed in the OE and expressed at lower levels in the ORNs. Moreover, Ascl1 and Neurog1, markers for transiently amplifying progenitor cells, were found in the OE and only weakly detected in the ORN probe set. Further, GAP43, a marker of immature neurons, was present in the OE and, to a lesser extent, in the ORN probe.

These expression pattern analyses confirm the purity of the FACS-sorted ORN sample and support the reliability of our NGS-analysis.

Non-Olfactory GPCRS (nGPCRS)

Our analysis allowed us to detect the transcripts of known and previously undescribed membrane proteins in the OE. GPCRs form the largest family of transmembrane proteins and function in various signaling pathways. Out of 407 annotated non-olfactory GPCRs (nGPCRs) in our data set (S7 Table), we detected the expression of 114 (sorted ORNs) and 159 (male/female OE) nGPCRs, which displayed FPKM values >1 (Table 2).

Table 2

The total number of GPCRs detected in mouse OE (RefSeq gene model).

	ORNs (heterozygous)	OE1 CD1 male	OE1 CD1 female	Number of genes annotated in RefSeq gene model
Chemosensory GPCRs
ORs	905	1,001	1,075	1,225
TAARs	13	11	13	15
Vomeronasal receptors 1	24	5	4	218
Vomeronasal receptors 2	20	4	5	130
Taste receptors	6	3	5	38
Non-chemosensory GPCRs	238	264	250	407
Total	1,206	1,288	1,352	2,033

We classified the nGPCRs into five main GPCR families: glutamate, rhodopsin, adhesion, frizzled and secretin. This classification reveals that the nGPCRs (including genes with FPKMs below 1) that were expressed in ORNs mainly belong to the rhodopsin family (64%). Ten percent belong to the adhesion family, approximately 4% to the frizzled family, 6% to the glutamate family, and 0.2% to the secretin family. An additional set of nGPCRs with unknown functions, no known ligands and atypical structure could not be classified (approximately 12%) nGPCRs that are highly expressed may hypothetically function as co-receptors in ORNs [86, 87]. Therefore, we created a ranking of the expression of all nGPCRs to identify the highly expressed candidates (Fig. 6). Interestingly, only six nGPCRs were ranked among the 1,000 most highly expressed genes in ORNs: Adipor1, Gpr178, Gabbr1, Gprc5c, Drd2, and Lphn3 (Fig. 6). An additional 30 nGPCRs exhibited medium expression profiles (FPKMs between 52–10), and 65 nGPCRs exhibited low expression profiles (FPKMs between 1–10) (S11 Fig.). Some of these GPCRs have previously been reported to be expressed in ORNs [48, 49]. Regarding the most highly expressed 30 nGPCRs, our study describes the expression of an additional 17 genes: Gpr178, Lphn1, Gpr137, Gpr162, Gpr155, Gpr63, Paqr9, Gpr108, Lgr4, Gpr89, Wls, A030009H04Rik, Gpr107, Gpr87, Gpr137b, Ptger1 and P2ry14 (Fig. 6). A few of these genes have also been tabularly presented in a recent transcriptome wide study of the total OE [52].

Figure 6

Expression ranking of non-olfactory GPCRs (nGPCR).

A. Heatmap showing the ranking of the 30 most highly expressed nGPCRs in the FACS-sorted ORNs. Only six nGPCRs were among the 1,000 most highly expressed genes (Adipor1, Gpr178, Gabbr1, Gprc5c, Drd2, and Lphn3). Of the 30 most highly expressed genes, the expression of 17 in the OE was no previously known. B. Heat map showing the ranking of nGPCRs that are specifically enriched in ORNs according the criteria that the nGPCRs’ FPKMs were greater than 1 and 5 times greater in the ORNs than in non-olfactory tissue (brain, liver, muscle and testes). A total of 18 specifically enriched nGPCRs were found in ORNs. Excluding the specifically enriched candidates that were already presented in the list of the 30 most highly expressed genes, an additional 10 genes were found to be specifically enriched and are new in terms of olfaction. Regarding the genes that were among the 30 most highly expressed and were specifically enriched in ORNs, 60% had neither been shown to be expressed in the OE or been ascribed any function in the OE in any previous study with the exception of for several candidates in a tabular form in a recent transcriptome-wide study of the total OE [52].

In a subsequent analysis, we classified the nGPCRs according to their expression levels and expression patterns. To identify nGPCRs that were specifically expressed in ORNs, we compared expression levels between ORNs and non-chemosensory tissue. We defined a gene as being enriched in ORNs, if the expression level of that gene was at least five times higher (as determined by RNA-Seq data) in the ORNs than in any of the control tissues (brain, liver, muscle and testes) (S8 Table). We found that 14 (FPKM > 1) nGPCRs were specifically enriched in ORNs, and 100 nGPCRs had broader tissue distributions (Fig. 6 and S7 Table).

Focusing on nGPCRs specifically enriched in ORNs and excluding genes already presented under the 30 top expressed candidates, we additionally detected the expression of 10 genes: Gpr171, Sfrp1, Rho, Gpr18, Gpr183, F2rl1, Gpr152, Gpr126, Gpr35 and Gpr82.

Finally, we constructed a quantitative ranking of the expression of all GPCRs detected in the OE and ORNs (S7 Table). We selected some of these new nGPCRs and verified their expression in ORNs with in situ hybridization (Fig. 7). All selected probes produced signals in the mature ORN layer as predicted by the FPKM values from ORNs. Additionally, we tested the most highly expressed GPCR, Adipor1 and two other identified GPCRs with RT-PCR and determined their relative expression levels (Fig. 8). As expected, we found that Adipor1 was most highly expressed followed by a lesser extent for Gpr177. These results match our RNA-Seq and in situ hybridization results.

Figure 7

In situ hybridization for nGPCRs mRNA.

A: Adipor1 (adiponectin receptor 1) B: Gpr178; C: Gpr155; D: Gpr177 (aka Wls); E: Ptgdr (prostaglandin D receptor). All transcripts were detected in the mature ORN cell layer as predicted by the expression levels observed in sorted ORNs. Scale bar = 30 µm.

Figure 8

Expression levels of GPCRs by real time RT-PCR.

Bar chart shows the relative expression levels of selected GPCRs relative to adenylyl cyclase type III (ACIII) as determined by RT-PCR. In accordance with our RNA-Seq data, the most highly expressed GPCR was Adipor1. Error bars represent the SEM.

Non-GPCR Membrane Proteins

The detection of unknown transmembrane proteins in several recent proteome studies has led to the discovery of proteins with essential functions in olfactory signaling processes. For example, Ano2 was identified as the olfactory CaCC [16, 15].

To identify other potential candidates, we searched for new non-GPCR membrane proteins that were highly and/or specifically expressed in ORNs. We assembled a catalogue of integral membrane proteins by manual inspection and GO terms (integral to membrane: GO: 0016021) using Ontologizer [88]. This process led to the detection of the expression (FPKMs > 1) of up to 2,339 genes in ORNs (S8 Table), 2,706 in male OE, and 2,575 in female OE. We next ranked these genes based on their expression levels and expression patterns.

This ranking showed that, of the 1,000 most highly expressed genes, 20% were genes that code for membrane proteins. We found the following known components of basic signaling among the 30 most highly expressed membrane protein genes: ACIII, Cnga2, Cnga4 Ano2, and Rtp1, which is involved in OR trafficking [21]. In addition to these genes that are known to have roles in basic signaling, we found Stoml3, Stom, Umodl1, Plekhb1, Atp1a1, Nsg1, Tmbim6, Atp1b1, Aplp2, Aplp1, Tmem66, Sgpl1, Sec14l3, Kcnc4, Clstn1, Rtn1, Ormdl3, Flrt1, Mslnl, Igsf8, Ncam1, Olfm1, Acsl6 and Faim2 (Fig. 9).

Figure 9

Expression patterns and ranking of genes coding for non-GPCR membrane proteins.

A. Heatmap showing the ranking of the 30 most highly expressed genes in the FACS-sorted ORNs. B. Heatmap showing the ranking of the 30 most highly expressed genes that were specifically enriched in ORNs according to criteria that their FPKMs > 1 and their expression level in ORNs was at least 5x greater than that in non-olfactory tissue (brain, liver, muscle and testes).

Out of these 30 most highly expressed genes in ORNs, we detected 11 genes (Fig. 9), whose expression was not reported in the OE. That these genes ranked among the 200 most highly expressed genes indicates their potential importance in olfaction.

Among the top 30 specifically enriched genes, we also found candidates that were among the 30 most highly expressed genes overall: Stom, Stoml3, Umodl1, Cnga2, Ano2, Rtp1, ACIII and Mslnl. Further known components of olfactory signal transduction that specifically were enriched in ORNs were the other subunits of the CNG channel (Cngb1, Cnga4), Scn9a, Scn5a and Trpm5, which is involved in the detection of semiochemicals [89]. Additionally, the following genes exhibited specifically enriched expression patterns: Sec14l3, Gcnt7, Tmem211, B4galnt3, Pcdhb1, Cldn9, Abca13, Pirt, Kcnmb3, Mcoln3, Svopl, Atp2c2, Flrt1, Manea, Gldn and Kcnh4.

We next focused on genes involved in the transport or flux of ions across the membrane.

We detected 132 ion channels, 199 members of the solute carrier (SLC) superfamily and 122 active transporters in ORNs with FPKMs >1 (S10 Table).

Our analysis revealed that the CNG channel subunits, Ano2, Scn9a and Scn4b were among the 10 most highly expressed channel genes in ORNs. Several potassium channel genes showed relatively high and specific expression patterns in ORNs. The most highly expressed was Kcnc4, the expression of which in ORNs has also been reported by Sammeta et al. (2007) [48].

In addition to Kcnc4, we detected two members of the voltage-gated KCNH channel family, Kcnh3 and Kcnh4 and several members of the KCNA family: Kcna1, Kcna2, Kcna5, Kcna6 and Kcna10. Interestingly, the expression of these genes in ORNs has not previously been reported (S10 Table).

As mentioned above, the expression levels of Ano2 were high in ORNs. We also detected other members of the anoctamin family, namely Ano1, 3, 6, 7, 8, 9 and 10 (S10 Table). A small subset of the OE RNA-Seq data for anoctamines 1–10 and other genes (Ttyh1–3, Trpa1, Trpm8 and Trpv1) were previously published by Schöbel and colleagues (2013) [90] in a Heat Map figure.

Members of the transmembrane channel-like protein (TMCs) family are evolutionarily related to the anoctamines [91]. We detected expression of Tmc4, 5, and 7 in ORNs (S1 Table). We further identified the expression of up to 7 aquaporin genes; Aqp3 was strongly expressed in ORNs (S10 Table). Hcn2, a channel protein that is abundantly expressed in the olfactory bulb [92], was also detected.

We next focused on SLC members that were specifically expressed in ORNs (FPKMs >1 and 5x greater expressing in ORNs than in non-olfactory tissue) and have reported function in ion homeostasis. We detected Nkcc1 (Slca12a2), which is important for Cl- accumulation in ORNs [18, 19], Nckx4 (Slc24a4), a Na+/Ca2+-exchanger involved in rapid response termination and adaptation of the olfactory response [42], and Nckx2 (Slc24a2). Our data also revealed the expression of Slc8a1, a Na+/Ca2+ exchanger, in ORNs [93]; this exchanger may contribute to the regulation of Ca2+ flux in ORNs (S10 Table).

Interestingly, we, for the first time, report the detection of the specific expression of additional membrane proteins in ORNs/OE; these novel proteins ABC-transporters, ATPases, tetraspanines, TMEM-proteins, WD-repeat domain proteins, Gram domain-containing proteins and several other proteins with unknown functions (Table 3).

Table 3

Further genes for non-GPCR membrane proteins that were among the 1,000 most highly expressed genes.

ABC transporter	Abca13, Abca5, Abcc4, Abcc1, Abcg1, Abcc10, Abca7
further ATPases	Atp6v1b1, Atp2c2, Atp6v1c2, Atp10d, Atp6v0a4, Atp13a5, Atp11a, Atp2c1
other subunits of pumps or transporter	Slco15a2, Sec14l3, Slc44a2, Slc15a2, Slc27a2, Slc22a20, Slc6a6
tetraspanines	Tspan3, Tspan13, Tspan7
transmembrane proteins	Tmem66, Tmem64, Tmem59, Tmem30a, Tmem63b, Tmem151b, Tmem50a, Tmem205, Tmem213
WD-repeat domain proteins	Wdr17
GRAM domain containing proteins	Gramd1c
other proteins	Tusc5, Dnajb14, Efr3b, Gpm6a, Homer2, Tm9sf3 Ttc9, Nehrf1, Nherf2

We selected several genes and confirmed their expression patterns with ISH (Fig. 10). All probes produced positive signals in the mature ORN layer.

Figure 10

In situ hybridization for mRNA of non-GPCR membrane proteins.

A: Stom (stomatin), B: Stoml3: Stomatin-like protein 3, C: Gpm6a (glycoprotein m6a), D: Gpm6b (glycoprotein m6b), E: Unc45a (protein unc-45 homolog A), F: Efr3b (EFR3 homolog B), G: Ttc9 (tetratricopeptide repeat domain 9), H: Homer2 (homer protein homolog 2), I: Tspan7 (tetraspanin 7), J: Wdr17 (WD repeat domain 17), K: Tm9sf3 (transmembrane 9 superfamily member 3), L: Tusc5 (tumor suppressor candidate 5), M: Gramd1c (GRAM domain-containing protein 1), N: Mslnl (mesothelin-like protein precursor), O: Pcdhb1 (protocadherin beta 1), P: Manea (glycoprotein endo-alpha-1,2-mannosidase), Q: Tmc5 (transmembrane channel-like gene family 5), R: Tmc4 (transmembrane channel-like gene family 4), S: Slc15a2 (solute carrier family 15 (H+/peptide transporter), member 2), T: Slc27a2 (solute carrier family 27 (fatty acid transporter), member 2), U: Slco1a5 (solute carrier organic anion transporter family, member 1a5), V: Slc44a2 (solute carrier family 44, member 2 / choline transporter-like protein 2), W: Slc22a20 (solute carrier family 22 member 20), X: Slc6a6 (sodium- and chloride-dependent taurine transporter), Y: Tmem205 (transmembrane protein 205), Z: Tmem66 (transmembrane protein 66), AA: Tmem213 (transmembrane protein 213), AB: Slc9r3a1 (solute carrier family 9 (sodium/hydrogen exchanger), member 3 regulator 1, AC: Sec14l3 (Sec14-like protein 3). All transcripts were detected in the mature ORN cell layer as predicted by the expression in sorted ORNs. Scale bar = 30 µm

Classes of Transcripts Enriched in ORNs

To obtain an overview of the transcript classes that were enriched in ORNs, we analyzed the gene expression patterns according to gene ontology (GO) categories using the Ontologizer software [88].

Previous microarray and proteome studies have described several transcript classes that are enriched in ORNs [48, 49, 54, 55, 53, 56, 16].

Genes associated with the sensory perception of smell (including ORs, OMP, Gαolf, and ACIII) represented the predominant class/GO term (GO: 0007608) in our list (followed by the GO categories of ion transport (GO:0006811) and cilia morphogenesis (GO:0005929) (Fig. 11). Proliferating basal cells continuously replace dying or aging ORNs in the OE. Thus, it is not surprising to find that transcripts for neuronal differentiation (GO:0030182) are enriched in the OE [84]. The Ca2+-dependent odorant induced excitation and adaptation process [94] is the origin for the class/GO-term (GO:005509) that includes genes that encode proteins with calcium binding properties. Biological processes, such as protein transport (GO:0015031), RNA-processing (GO:0006396), and chromatin modification (GO:016568), are ubiquitous functions in cells. The corresponding categories contain the genes that support basic cellular mechanisms [95].

Figure 11

Classification of GO terms enriched in ORNs.

Based on results from Ontologizer [88], specifically enriched GO terms are represented in this bar chart. For simplification, related classes were combined.

Genes classifications that include ubiquitously expressed genes but are also likely to include genes that are more restricted to the machinery of ORNs [48] include the following: regulation of gene expression (GO: 10468), zinc ion binding (GO:0008270), kinase (GO:0016301) or phosphatase activity (GO:0016791), cell adhesion (GO:0007155), calcium binding (GO:0005509), TPR (GO:0030911), and WD repeat (GO:0005515).

Genes Related to cAMP-Signaling

In addition to the highly expressed Gαolf, ACIII and CNG channel subunits, we analyzed the expression patterns of other proteins involved in cAMP signaling. We detected strong expression of phosphodiesterases (PDEs), which are possible candidates for rapid termination of olfactory signal transduction due to their action in the degradation of the second messenger cAMP [22, 96]. Previously, three PDEs, Pde1C [97, 34], Pde4A [35] and Pde2 [98] had been identified in mammalian ORNs. Our transcriptome data revealed the expression of new PDEs in ORNs and the OE (Fig. 12). Our dataset demonstrated that all 11 PDE gene families are represented by the expression of at least one gene in ORNs and/or the OE. We report for the first time that Pde7b is expressed in ORNs and that it is one of the three most highly expressed PDEs in ORNs. Pde7b is a cAMP-specific phosphodiesterase that may have a modulating effect in olfactory signal transduction via its ability to hydrolyze cAMP. Moreover, we discovered that the Pde6d subunit is also highly expressed. Pde6 is localized in rod and cone photoreceptors where it regulates cytoplasmic cGMP concentrations [99, 100]. The high level of expression of the cGMP-specific Pde6d raises the question of what impact cGMP has on olfactory signal transduction. The IBMX-insensitive PDEs Pde8 (Pde8a, Pde8b) and Pde9 (Pde9a) were detected at low expression levels in ORNs and the OE.

Figure 12

Expression patterns of genes involved in cAMP-dependent signaling.

A. Expression patterns of cyclic nucleotide phosphodiesterase (PDE) genes. Heatmap showing the expression levels of PDEs in olfactory tissue and non-olfactory tissues (brain, liver, muscle and testes). B. Expression patterns of cAMP-dependent protein kinases and associated regulatory proteins in olfactory tissue and non-olfactory tissues (brain, liver, muscle and testes). Higher FPKM values are indicated by deeper colors.

Additionally, we detected high levels of expression of cAMP-dependent protein kinases (Prk), which form tetramers that consist of 2 regulatory and 2 catalytic subunits [101]. Regarding the catalytic subunits, we detected the alpha, beta and Prkx subunits, and we detected 4 isoforms the regulatory subunit (alpha I/II and beta I/II). Moreover, we detected weak expression of 3 cAMP-dependent PKA inhibitor genes (Pkia, Pkib and Pkig) that regulate the action of PKAs (Fig. 12).

Genes Related to IP3 / PI Signaling

The physiological function of phosphatidylinositol (PI) signaling mechanisms in ORNs have been discussed for several years [45]. Recent findings suggest that PI3K-dependent signaling mediates the inhibition of odorant responses in ORNs, which express GPCR-activated isoforms of PI3K and exhibit odorant-induced PI3K activity [43, 102, 103, 47]. Several PI-kinases (S11 Table) were found in ORNs. The expression of the phosphatidylinositol 3-kinase type catalytic subunits PI3K-β and γ were of special interest because these subunits can be activated by the β/γ subunits of G-proteins, and both have been detected in the OE [103]. Our sequencing data revealed that ORNs mainly express PI3kα and PI3kβ and, to a lesser extent, PIk3γ. Nearly all subtypes were expressed more highly in ORNs than in non-OE control tissues (S11 Table); however, the expression levels of these kinases were only moderate compared to the expression of elements of cAMP signaling. Nevertheless, all four types of PI3-kinases were detected in ORNs by in situ hybridization (S12 Fig.).

Phospholipase A2 (XIIA), phospholipase C, (b3, 4) phospholipase C-like 2 and phospholipase D3 are the most highly expressed phospholipases, but none of these are specific to ORNs (S11 Table). The same is true for the most highly expressed PKCs. The expression of IP3-receptors was previously shown by Restrepo et al. (1990) [104], Fadool and Ache (1992) [105], Kalinoski et al. (1992) [106], Restrepo et al. (1992) [107], Cunningham et al. (1993) [108] and Munger et al. (2000) [39]. Among the three subtypes, IPTr1 is most highly expressed, although its expression levels are about 100-fold lower level than those of the CNG channels. In contrast to genes involved in the cAMP-mediated signaling pathway, none of the genes involved in PI signaling were among the 500 most highly expressed genes.

Expression Patterns of TRP Channels

TRP channels constitute a family of proteins that respond to a variety of stimuli [109]. The expression of several TRP channels has been reported in the OE [110, 89]. Our sequencing data provide a comprehensive overview of the expression levels of TRP channel transcripts in both the OE and FACS-sorted ORNs (Fig. 13). Within the TRPC subfamily, TRPC1, TRPC2 and TRPC4 were detected in ORNs, and TRPC1 transcripts were the most enriched. While the transcripts of seven TRPM subfamily members were amplified from pooled cDNAs of whole OE, only four members showed somewhat higher expression levels in ORNs (Fig. 13). Transcript levels of TRPM7, a channel fused to a protein kinase [111], were the highest. Among TRPV subfamily members, only TRPV2 transcripts were present in ORNs. We also detected transcripts for two of the intracellular mucolipin TRP proteins, TRPML1 (Mcoln1) and TRPML3 (Mcoln3). Finally, and consistent with previous studies [112], we found both PKD1 and PKD2 transcripts in ciliated ORNs.

Figure 13

Expression pattern of TRP channel genes.

Heatmap showing the expression levels of annotated TRP channels in olfactory (male and female OE, and FACS-sorted ORNs) and non-olfactory tissues (brain, muscle, liver and testes). Higher FPKM values are indicated by deeper colors.

TRPC1 and TRPM7 transcripts were among the most enriched TRP channel transcripts in ORNs, yet they constituted only ~ 7% of the total number of CNGA2 channel subunit transcripts. We examined the expression by immunohistochemistry and detected TRPC1 and TRPM7 immunoreactivity in ORNs (Fig. 14).

Figure 14

Expression of TRPC1 and TRPM7 in the OE.

Protein expression was assessed staining coronal sections of OMP-GFP mice with anti-TRPC1 and TRPM7 antibodies. TRPC1 immunoreactvity was observed in a few ORNs in 1 of 3 stained sections, while TRPM7 immunoreactivity was detected in a larger number of ORNs (1 stained section). Omission of the primary antibody served as control. Scale bar, 20 μm.

Discussion

A comprehensive understanding of the transcriptome is fundamental for the study of the functionality of ORNs and the machinery of their crucial elements. Several previous studies, including microarray [66, 48, 113, 49] and proteomics studies [53, 55, 56, 16, 15, 54], have revealed the expression of a large number of known and novel genes in the OE. While proteomic studies have facilitated progress, for example, the identification of the Ano2 gene that encodes the olfactory CaCC, most of these studies lack proper quantification of genes and remain incomplete, especially regarding the detection of OR gene expression. Recent advances in next generation sequencing provided the means to obtain a “complete” list of protein-coding genes expressed in the OE in a high-throughput fashion.

Here, we present a comprehensive analysis of the murine olfactory transcriptome that was generated by the high-throughput NGS technique. For the first time, we analyzed FACS-sorted ORNs to gain knowledge about the differences in expression profiles between ORNs and other cell types in the OE with the NGS technique.

We compared our data to transcriptome data from non-olfactory tissues (brain, liver, muscle and testes) to identify ORN-specific genes. We further described the expression of known and novel membrane proteins and generated a list of all non-olfactory GPCRs that were expressed in the OE.

During the preparation of our manuscript, a murine olfactory transcriptome study was published that focused mainly on sex-specific differences [51] in the OE. Using NGS, these authors detected nearly all classified OR genes. A second study by Keydar et al. (2013) [52] focused on a more general catalogue of genes expressed in the OE, but did not provide any data on specificity of expressed genes for ORNs.

Our analysis, and particularly our comprehensive list of genes expressed in ORNs, complements previous studies and provides a basis for the discovery and study of novel genes expressed in the OE, especially ORNs.

Chemoreceptors

In the mouse, the OR family comprises ~1209 OR genes, of which 913 have been identified as functional genes and 296 have been identified as OR-pseudogenes [1]. Using NGS, we successfully detected up to ~97% (FPKM > 0.1) of all OR genes based on the Refseq gene model in our probes. A previous study by Shiao et al. (2012) [51] detected 99% of all OR genes, and this result was possibly due to the even greater sequencing depth of that study. The fact that virtually all ORs are highly expressed in the OE and virtually none are expressed in other tissues (except for the testes) is astonishing for a gene class with ~1,200 members and suggests an extremely thorough regulation of gene expression.

Sequencing of the FACS-sorted ORNs detected the expression of fewer OR genes, and this was likely due to the lower number of sequences generated for homozygous ORNs (~13 million) and the limited number of sorted neurons. As we used OMP-GFP transgenic animals to obtain ORN samples for sequencing, we cannot exclude the possibility that some receptors were underrepresented. OMP is a marker of mature olfactory neurons [114]. The expression level of a particular OR gene is strongly dependent on the receptor type. In CD1 mice, the expression patterns of OR genes were strongly correlated between female and male OEs (Pearson coefficient, r = 0.83). The expression rankings for the OR genes Olfr533, Olfr1507 and Olfr309 were similar in male and female mice, and the three most highly expressed OR genes found in the sorted ORNs were among the ten most highly expressed OR genes in the OE. That similar percentages of OR genes were detected in our and previous studies suggests that the OR subgenome does not only strongly correlate between sexes, but also between mouse strains because OE samples were derived from female CD1 and C57BL/6J mice (S8 Fig., S9 Fig.). We calculate a Pearson coefficient of r = 0.75 for between these two groups, which suggest that only a small percentage of OR genes is strain-specific and differentially expressed. However, further studies are needed for detailed analyses of strain-specific differences in the OR subgenome.

One of the most highly expressed OR genes, Olfr1507, was also the most transcribed OR gene in the FACS-sorted ORNs sample in our study (FPKM ~ 97) [115]. The number of ORNs expressing an individual OR varies and is coupled to the strict regulation of gene choice [2]. Olfr1507 (also called MOR28) is located in a gene cluster consisting of seven ORs and an upstream H regulatory element on chromosome 14. The H element interacts with the most proximal OR gene (Olfr1507) and may facilitate its expression [116]. The expression levels of the remaining 6 genes in the cluster (Olfr1508, Olfr1509, Olfr1510, Olfr1511, Olfr1512) are considerably lower, suggesting that the most proximal position to the regulatory H element indeed defines the expression strength of these ORs in the OE.

In addition to ORs, we detected the expression of other types of chemoreceptors in the OE. Mice have 15 TAAR genes, 14 of which are expressed in the OE [71]. We detected all 14 TAAR genes in the OE, which underlines the quality of our data set. Interestingly, we also detected the expression of vomeronasal receptors in the OE. We identified the expression of 15 V1R and 5 V2R genes in the OE and 24 V1R and 21 V2R genes in FACS-sorted ORNs. The FPKM values were often small (< 1) and comparable to the minimally expressed ORs, suggesting a possible mosaic-like gene expression pattern for these chemoreceptors in only a few cells in the OE. The expression of the V1- and V2R genes in ORNs implies that the main OE contributes to pheromone detection. Interestingly, in previous studies, pheromone-induced behavior was not altered in mice after removal of the VNO [117], and electrophysiological data additionally support the detection of several pheromones in the OE [118]. Although we also detected the expression of TRPC2 and PLC (both of which are part of the transduction cascade in VNO neurons) [119-121] in our transcriptome analysis of the OE, the low number of detected vomeronasal receptors likely limits their contribution to OE-mediated pheromone sensing. Thus, it seems reasonable that ORs or other receptor classes are also involved in pheromone perception in the OE.

We also detected an enrichment of two taste receptors in ORNs, Tas1r1 and Tas1r3, which are known to heterodimerize and subsequently form the umami-receptor. The detection of taste receptors in the OE could likewise confirm the presence of solitary chemosensory neurons (SCN) in the OE, which are known to express TasRs [122]. Chemosensory information derived from taste and olfaction is used by organism to value the quality of food. The detection of taste receptors in the OE may underline the fact that the sense of olfaction and taste are linked modalities.

Signal Transduction in ORNs

The commonly known components of the classical olfactory signal transduction pathway were among the 200 most highly expressed genes in our murine OE transcriptome analysis.

It is reasonable to assume that other global players in olfactory signaling should be as highly expressed. In this respect, we detected several GPCRs and other membrane proteins with high expression levels. Their relevance to olfactory signaling remains elusive. Future studies may uncover the function of the presented candidate genes in ORNs.

Among the 200 most highly expressed genes were the previously unrecognized phosphodiesterases Pde6d and Pde7b. Previous studies suggest that hydrolyzation of cAMP by phosphodiesterases is involved in the termination of the olfactory signal transduction. The Ca2+/calmodulin-stimulated Pde1c appears to be enriched in the olfactory cilia, whereas the cAMP-specific Pde4a occurs throughout the cell but not in the cilia [34, 98]. Simultaneous disruption of the Pde1c and Pde4a genes in mice leads to prolonged response termination in electro-olfactogram measurements [97]. The potential regulatory function of this novel expression of Pde7b and Pde6d in the OE in olfactory signaling will be the subject of future studies. While Pde6 is primarily known to be localized to photoreceptors, where it regulates the cytoplasmatic cGMP concentration [99], less is known about the subunit Pde6d. Previous studies have reported a regulatory function of Pde6d in the membrane association of Ras and Rap GTPases [123] and a potential contribution to cell proliferation [124].

We additionally focused on molecules supporting alternative signaling pathways in ORNs. The involvement of odorant-stimulated PI signaling in ORNs has previously been shown [45, 43, 103]. This pathway involves PLC- and PI3K-dependent signaling. Our data provide support for a potential PI-mediated signaling of odorants as we detected weak expression of PLC, PI3K and IP3 receptors.

Interestingly, our data revealed the expression of several TRP channel transcripts. TRP channels are known downstream targets of PI-mediated signaling in chemosensory cells [125, 126]. The detection of semiochemicals is thought to involve Ca2+-activated TRPM5 channels that are expressed in a subset of ORNs in which both CNGA2 and PLC are co-expressed [89]. TRP channels form non-selective cation channels that are permeable to calcium, which may suggest a contribution to processes in the OE such as transduction, transcriptional control or proliferation. Among the detected TRP channel transcripts, TRPC1 and TRPM7 were the most enriched. Immunohistochemistry indicated expression of these proteins in ORNs, but their function in ORNs is unknown. A recent study proposed a regulatory function of TRPC1 in the fine-tuning of neuronal migration [127]. Another study by Kerstein et al. (2013) [128] addressed a mechanosensitive role of TRPC1 in spinal neuron growth cones. Additionally, TRPC1 has been proposed to be involved in store-operated Ca2+ entry, which regulates cell proliferation [129]. TRPM7 is ubiquitously expressed and has been implicated in a variety of cellular functions including magnesium homeostasis, cell cycle progression and control of the production of ROS [130]. Knock-out studies revealed a crucial role of TRPM7 in embryonic development [131, 132]. Future studies are required to examine potential functions of TRPC1 and TRPM7 in the OE.

In this study, we present a comprehensive expression pattern analysis of GPCRs in the murine OE. In FACS-sorted ORNs, we detected the expression of 114 nGPCRs (FPKM>1). The expression of most of these proteins had previously been documented in at least one study that employed proteome or transcriptome analytic approaches [48, 49, 53, 55, 56, 16, 15, 54].

Our analysis revealed the expression of several unrecognized GPCRs for which no expression or functions in ORNs have been described. In summary, we detected the expression of 27 GPCRs, whose expression in ORNs has not been reported before, and 14 of these exhibited specific enrichment (FPKM >1 and an expression level that was at least 5 times greater than that in control tissues) in ORNs: Gpr178, Lphn1, Gpr63, Gpr137, Gpr162, Gpr155, Gpr87, Paqr9, Gpr108, Lgr4, Gpr89, Wls, A030009H04Rik, Gpr171, Gpr107, Gpr137b, Ptger1, P2ry14, Sfrp1, Rho, Gpr18, Gpr183, F2rl1, Gpr152, Gpr126, Gpr35 and Gpr82.

Interestingly, among 1,000 most highly expressed genes in ORNs, only 6 non-olfactory GPCRs were found: Adipor1, Gpr 178, Gabbr1, Gprc5c, Drd2, and Lphn3.

GPCRs involved in food intake

The most highly expressed GPCR is Adipor1, which was first detected by Hass et al. (2008) [133] in mature ORNs. Another highly expressed adiponectin receptor is Paqr3 [134]. It is possible that adiponectin regulates food intake by acting as an appetite stimulator and conveying a starvation signal to the brain via AdipoR1 and AdipoR2 [135]. In animals, the sense of smell plays a crucial role in finding food resources. Thus, an interaction of adiponectin and AdipoR1 could modulate the function of ORNs depending on the nutritional status of the body. In this context, it is interesting to note that the modulatory effects of two anorectic peptides, insulin and leptin, on ORNs have already been demonstrated [136]. However, receptors for insulin, leptin, ghrelin and other peptides that are involved in the regulation of hunger and satiety were only weakly expressed in the OE (S13 Fig.); these findings lead to the assumption that adiponectin is the most important regulatory hormone in the OE. Therefore, future studies should concentrate on investigating the functions of AdipoR1 and Paqr3 and the modulatory effect of adiponectin on olfaction.

Furthermore, two orphan GPCRs, Gpr162, and Gpr82, are both possibly involved in food intake and glucose homeostasis [137, 138].

mGluR family

Gprc5c was discovered in the mouse OE through proteomics [16], is an orphan receptor of the mGluR family [139], and is also expressed in goldfish OE [140]. Like Gprc5c, Gpr158 belongs to the mGluR family and is most similar to the metabotropic GABA receptors [141, 142].

Inhibitory effect on olfactory signaling

Gabbr1, and the less-strongly expressed subunit Gabbr2, form the heterodimeric GABAB receptor [143], which has been shown to inhibit ORN axonal outgrowth [144] and is possibly located presynaptically. Activation of GABAB receptors, which couple to Gαi/o, stimulates increases in cAMP through βγ-mediated activation of adenylyl cyclase 2 and simultaneously inhibits Gαs-mediated activation of other adenylyl cyclases in the rat olfactory bulb [140, 145]. Supporting this GABAergic effect, Drd2 dopamine receptors are inhibitory on the input from ORNs and provide lateral inhibition of mitral cells, which provides olfactory discrimination in rodents [146]. Moreover Gpr63, a sphingosine 1-phosphate receptor, is expressed in mitral cells of the OB, but its function in this region in unknown [147].

Furthermore, many neurons that are classified as GABAergic might express GPR155 [148]. These findings imply that Gpr155 has an important role in GABAergic neurotransmission. GABAergic input has an inhibitory effect on olfaction [144]; thus, Gpr155 could also have a supportive function for the GABAergic effect on ORNs.

Cell Architecture / Cell development

Additionally, we detected a few GPCRs that have possible roles in anatomical structural development.

The protocadherin receptors Celsr2 and Celsr3 are key regulators of the correct positioning of cilia and, consequently, cilia function [149]. The latrophilin receptors (Lphn1, Lphn3) are involved in signaling of tissue polarity and morphogenesis [150].

Fzd3 is required for neurogenesis and target innervation during sympathetic nervous system development [151].

Lgr4 is part of the Wnt signaling pathway that is involved in the cell proliferation of the intestinal epithelium [152].

Orphan GPCRs with unknown functions

Furthermore, we detected several GPCRs, the functions of which are presently unknown, that were specific to ORNs.

Gpr178 (also named Tmem181a) is an orphan GPCR with unknown function that was not only among the 1,000 most highly expressed genes but was also specifically enriched in ORNs. This gene shares homologies with the newly identified Gpr177 (aka Wls) gene, which has been reported to be involved in the Wnt signaling pathway [153].

No endogenous ligands have been identified for Gpr137, Gpr152, Gpr158, Gpr89, Gpr171, Gpr108 or A030009H04Rik, the exact functions of these receptors are unknown.

We selected a few candidates and confirmed their expression in ORNs with ISH experiments: these experiments revealed that the mRNA transcripts of these candidates were prominent in olfactory neurons. Because these candidates were specifically expressed in mature ORNs and because of the described functions of these genes in other tissues, we suggest that they might play important roles in ORNs.

Previous studies have reported that ORs can co-assemble with other GPCRs such as beta adrenergic and muscarinergic receptors [86, 87]. However, according to our analysis, both of these GPCRs are only weakly expressed in ORNs. The six GPCRs we found to be highly expressed could form such co-receptors in principle, but our expression data did not reveal any candidates that were expressed at levels comparable those of the ORs. Consequently, the expression pattern of GPCRs contradicts the hypothesis of a GPCR co-receptor for ORs. In contrast, our data confirm that RTP1 and RTP2 are highly expressed; in addition to their known chaperone function, RTP1 and RTP2 have been reported to be required as co-receptors for ORs [21]. Our data support the hypothesis that RTP1 or RTP2 co-assemble with ORs in a stoichiometric manner because the expression levels of RTP1 and RTP2 are comparable that of the ORs.

Altogether, in our analysis, we provided a detailed expression ranking of all GPCRs detected in the OE and noted GPCRs new in terms of olfaction. Hence, our data directly leads to new perspectives to focus on so far unknown GPCRs. Due to a strong ORNs specific expression pattern; these GPCRs have a supposable important role in ORNs.

As membrane proteins form key nodes in the olfactory signaling process, we focused on finding further membrane proteins that were not previously known to be expressed in, or have a specific function in, the OE. We detected the expression of 2,339 (FPKM > 1) genes for non-GPCR membrane proteins in ORNs. We ranked the expression of these genes and specifically marked those that were enriched in the OE compared to non-olfactory tissues to highlight genes that possible have important roles in the function of the ORNs. The main known components of olfaction (ACIII, CNG channels, Ano2, Rtp1) were found among the 30 most highly expressed genes.

Additionally, we identified 22 novel genes (11 that were specifically enriched in ORNs) that because of their high expression levels and their described functions in other tissues, are promising and important candidates for future research in olfaction: Abca13, Aplp1, Atp1b1, Atp2c2, B4galnt3, Faim2, Flrt1, Gcnt7, Gramd1c, Kcnh4, Kcnmb3, Olfm1, Ormdl3, Pcdhb1, Pirt, Plekhb1, Rtn1, Sec14l3, Svopl, Tmem211 and Tmem66.

Transport

We newly detected the expression of a Na+/K+-transporting ATPase (Atp1a1, Atp1b1) in the OE; this ATPase modulates membrane potential [154]. We also detected Nsg1 in the OE; this molecule is involved in regulating receptor recycling [155]. The expression of both genes is important for olfactory signaling, as the adjustment of ion homeostasis is the basis for the depolarization of neurons and the recycling of receptors.

Further, the following transporter proteins were present: Abca13, an ABC transporter [156], Svopl, a putative transporter [157] and Atp2c2, a Ca2+ ATPase [158].

Aplp1 and Tmem66 are proteins that regulate the Ca2+ homeostasis [159, 160]. Rtn1 and Ormdl3 are involved in membrane vesicle trafficking and protein folding [161, 162], and these two processes are necessary for proper membrane expression of GPCRs and other membrane proteins.

Pirt, a phosphoinositide-binding protein, has been reported to function as a regulatory subunit of TRPV1 [163].

Ion channels

It has been reported that 80–90% of the receptor currents that result from odorant exposure and OR interaction are mainly generated by Cl- ions, which pass through the CaCC Ano2. Billig et al. (2011) showed that EOG recordings from Ano2-/- deficient mice were only reduced by up to ~40% [164]; nevertheless, these mice were able to smell. Therefore, it is highly probable that other CaCC channels are expressed in ORNs and that these channels, together with Ano2, mediate the major part of the receptor current. Although Pifferi et al. (2006) showed that the CaCC bestrophin-2 (Best2) is expressed in ORNs, wild type, and mice lacking Best2, exhibited no significant differences in olfactory ability [165]. Thus, Best2 is not a CaCC that makes a primary contribution to the odorant-induced chloride current. Therefore, we searched for new CaCC candidates that are important and account for olfactory transduction. In our analysis, we detected ORN-specific expression of genes coding for members of the TMC (transmembrane channel-like) protein family, which has an evolutionary relationship with anoctamines [91]. Based on expression patterns and the homology to anoctamines, we suggest that Tmc5, together with Tmc4 and Tmc7, are suitable candidates for the other CaCCs and could be important in olfaction. Future studies should examine the involvement of these proteins in olfaction.

Regarding voltage-gated sodium channels, we detected high levels of expression of Scn9a. These expression levels confirm the vital role of this channel in olfaction, as it has been shown that the loss of this protein leads to anosmia [82].

Among the variety of other K+ channels (S10 Table), the Kcnc4 was the most highly expressed voltage-dependent potassium channel subunit in ORNs. In the OE, the expression of Kcnc4 was already shown by Sammeta et al. (2010) [95], where they postulate a Kcnc4 expression sensitive to neuronal damage. Additionally, for the first time, we detected high and specific levels of expression of members of the voltage-gated Kcnh channel family, which is involved in the regulation of neuronal excitability [166]. Kcnh3 and Kcnh4, which are potentially members of the EAG-like (ELK) K+ subfamily, exhibited strong expression patterns in ORNs, indicating that these channels could be primarily responsible for determining and raising action potential thresholds by acting at voltages around the firing threshold to suppress excitability [166, 167]. Hagendorf et al. (2009) have previously demonstrated the expression of Kcnh channels in the sensory neurons of the VNO and their key roles in determining neuronal excitability [168]. Thus, it is conceivable that these candidate genes in ORNs have similar roles. As it has been shown that members of the ELK subfamily are able to co-assemble with each other [169], we further suggest that Kcnh3 and Kcnh4 could form functional heteromultimers and contribute to the hyperpolarizing effect in ORNs.

Similarly, strong and specific patterns of expression were observed for Kcna5, Kcna2, Kcna6, and a weaker pattern of expression was observed for Kcna1; these genes all code for members of the shaker-related voltage-gated K+-channel family. Interestingly, Eldstrom et al. (2002) reported that PDZ domains are able to bind to Kcna5 and other members of this family and affect potassium currents by regulating the assembly of the subunits [170]. Thus, it is of great interest to determine what impact this channel region has and if it interacts with Mupp1, which was recently identified as a PDZ-protein that is expressed in ORNs [28].

Moreover, we, for the first time, detected a medium-level expression of Kcnk12 and low-level expressions of Kcnk4, Kcnk5 and Kcnk1 in ORNs. In neurons, members of the two-pore domain potassium channel family Kcnk have been reported to determine membrane potentials and membrane input resistances, which influence the magnitudes and kinetics of responses to synaptic inputs [171]. Further focus should also be placed on Kcnk10, which was highly expressed, and Kcnk2, although this protein was expressed at lower levels. The described mechano- and thermosensitivity of these proteins [172, 173], is remarkable; these proteins are novel targets for the study of the physical impact of mechano- and thermosensation during the detection of odorants.

Others

Of the remaining membrane proteins, we discovered that Mslnl is highly and specifically expressed in ORNs. No studies of the function of Mslnl in ORNs or any other tissue exist; Mslnl is a mesothelin-like protein. The related mesothelin protein is known to have cell adhesive properties [174]. Additionally, our data revealed that, in addition to Gpm6b, Gpm6a was highly expressed; Gpm6a is a four-transmembrane protein that is abundantly expressed in the nervous system [175]. The expression and function of Gpm6a in the OE have not been described previously. In the murine retina, Gpm6a is known to regulate retinal development by mediating cell-cell interactions that are involved in axon fasciculation [176]. As the visual and olfactory signal transduction systems are nearly identical on the molecular level [177], we suggest that Gpm6a and Gpm6b have similar roles during the development of ORN neurites in the OE. Of the many chaperones that are necessary for proteostasis, we detected a specific enrichment of Dnajb14 (Hsp40) in ORNs. Dnajb14 acts as a co-chaperone for the Hsp70 protein folding machinery through its role in determining substrate specificity [178]. It has previously been reported that another member of the Hsp40 family is able to bind to, and assist in the folding of, a GPCR (progesterone receptor) [179]. Neuhaus et al. (2006) showed that the specific enrichment of Hsc70t, a variant of the Hsp70 family of heat shock proteins, in the OE, assists in the folding and trafficking of particular ORs to the plasma membrane [180]. Thus, we propose that this co-chaperone may contribute to the proteostasis of ORs or other types of GPCRs. Wdr17 is a retina-specific transcript that is thought to be involved in signaling events within and between cells [181]. Homer2, whose function in the OE has not previously been described, is highly expressed in ORNs. Homer proteins are known to promote the targeting and expression processes of GPCRs by interacting with them [182, 183]. Furthermore, it has been shown that Homer proteins are expressed in the VNO and form complexes with TRPC2 channels and IP3-receptors [184]. Due to the strong expression level of Homer2, we suggest that this protein is strongly involved in the chaperoning and guiding of olfactory and non-olfactory GPCRs to different membrane sites in the ORNs. Additionally, Homer2 may also be capable of interacting with and modulating, among others channels, TRP channels and thus influence the signal cascades of ORNs. In addition to Homer2, we, newly detected the expression of Slc9a3r1 and Slc9a3r2, which code for the PDZ scaffolding proteins Nherf1 and Nherf2 (Na+/H+ exchanger regulatory factor), respectively. Both candidates have been reported to interact via their PDZ domains with several GPCRs [185]. Nisar et al. (2012) demonstrated an important role for Nherf1 in potentiating GPCR internalization [186]. Here, after receptor stimulation, Nherf1 interacts with GPCRs via the scaffolding protein arrestin. Additionally, Nherf1, which is known to interact with β2-AR [187], regulates receptor-mediated Na+/H+-exchange and down-regulates the receptor by increasing the recycling of the receptor [188, 189]. Nherf2 is known to specifically couple LPA2-receptors and PLCβ3 and to regulate activity by this process [190]. Additionally, both proteins demonstrably increase the Gαq-mediated signaling [191, 192]. It remains to be examined whether Nherf1 and Nherf2 are also capable of regulating GPCR localization and diversifying the signal cascade of ORs in ORNs. Of the several tetraspanins that are expressed in ORNs, we found high levels of expression of Tspan7. Tetraspanins are known to regulate the signaling, trafficking and biosynthetic processing of associated proteins. Specifically, Tspan7 is involved in synaptic maturation and function and promotes the formation of filopodia and dendritic spines [193]. We assume that, due to the high levels of expression and reported functions of Tspan7, this protein is highly involved in the formation of the neurites of ORNs. Furthermore, we detected an enrichment of Ttc9 in ORNs. Ttc9 belongs to a family of tetratricopeptide repeat (TPR)-containing proteins [194] that are involved with, among other things, protein transport and folding [195, 196] and cell cycle control [197] and transcription and splicing events [198]. However, the exact function of Ttc9 in ORNs has yet to be elucidated. Next, we detected high levels of expression of two members of the Tm9sf family. Tm9sf2 and Tm9sf3 are expressed throughout the epithelium. A member of the Tm9sf family has been shown to have functional ligand properties and has consequently been suggested to function as channels or small molecule transporter or receptors [199]. Unc45a is highly and specifically expressed in ORNs, and the functions of this protein have not been well described in olfaction. It has been reported that Unc45a is involved in cytokinesis and motility via chaperoning myosin and further cooperates with Hsp90 to chaperon progesterone receptors [200, 201]. Additionally, Tusc5, an adipocyte-specific transcript, was abundant in ORNs. The reported co-expression of Tusc5 in adipocytes and peripheral somatosensory neurons [202] indicates a possible connection between the energy status of the body and distinct sensory systems. Accordingly, we suppose that Tusc5 also has a regulatory function in olfactory perception that is similar to that of Adipor1 or Paqr3 and depends on nutritional status.

In this study, we have presented several membrane proteins that have not previously been identified in OE. As the expression levels of these proteins were similar to, or even exceeded, those of the known major players in olfactory signal transduction, these proteins potentially have important roles in olfactory processes. ISH experiments confirmed the pronounced expressions of these selected genes in ORNs. Therefore, we suggest that the proteins encoded by these genes are indeed involved in the function of olfactory neurons. Future studies should concentrate on uncovering the role of these proteins in the machinery of ORNs.

Our data provide a nearly complete catalogue of the genes expressed in, and involved in, the function and maintenance of the OE, especially ORNs. The molecular portrait of the OE revealed by this quantitative and comprehensive analysis of the murine transcriptome has uncovered new and valuable approaches that will be beneficial for the advancement of knowledge regarding the molecular mechanisms underlying olfaction and the functionality of ORNs.

Conclusion

The unmatched power of RNA-Seq in terms of quantitative and differential transcriptome analysis and the simplicity of the practical usage of this technique clearly prove that is this technique is a useful and important tool for OE transcriptomics. In this study, we were able to identify new potential players in olfaction. Furthermore, we demonstrated that these data provide a valuable framework for the interpretation and understanding of the function of ORNs. Finally, this technique currently enables the most comprehensive analyses and the easiest integration of the vast knowledge gained by previous studies.

Distribution of different housekeeping genes in murine tissues.

Heatmap showing the expression levels of different housekeeping genes in olfactory and non-olfactory tissue. Higher FPKM values are indicated by deeper colors. Gapdh: glyceraldehyde-3-phosphate dehydrogenase, Actb: actin, cytoplasmic 1, Ldha: L-lactate dehydrogenase A chain isoform 2, Ubc: polyubiquitin-C, Tubb3: tubulin beta-3 chain, Hprt: hypoxanthine-guanine phosphoribosyltransferase.

(TIF)

Click here for additional data file.

OR detection using RNA-Seq in all replicates of the OE: Bar chart showing the percentage of detected OR genes in all replicates of OE tissue (n = 13) and FACS-sorted ORNs.

Percentages were calculated based on the 1,125 OR genes annotated in the Refseq based gene model. Bars in light blue: FACS-sorted ORNs (homo- and heterozygous); blue: OE replicates of CD1 male mice (n = 4), dark blue: OE replicates of CD1 female mice (n = 5), black: OE replicates of C57BL6 female mice (n = 4).

(TIF)

Click here for additional data file.

Detection of weakly and highly expressed ORs with RNA-Seq.

Sample representation of read coverage of ORs with different expression strength can be visualized by the Integrative Genomics Viewer. Shown are exemplary cufflinks data for OE of CD1 female mice; A. Olfr1507, FPKM = 70; B. Olfr730, FPKM = 5.4; C. Olfr1265, FPKM = 0.1. The exons are indicated by blue bars and introns by thin lines. The grey segments indicate reads that were mapped onto reference genome and red bridges exon spanning reads.. Above, the read coverage is shown (detected and mapped counts/bases at each respective position). In highly expressed ORs, 5’ UTRs can be identified by exon-spanning reads. For medium or low expressed ORs, this is not possible due to the lower number of mapped reads.

(TIF)

Click here for additional data file.

Ectopic expression of ORs in testes, brain, liver and muscle.

Heatmap showing the expression of OR genes expressed in non-olfactory tissues (testes, brain, liver, muscle). Higher FPKM values are indicated by deeper colors.

(TIF)

Click here for additional data file.

Classification of FPKM values.

A: Classification of FPKM values of OR genes. Graph showing the distribution of OR genes according to FPKM value classes in olfactory tissues (FACS sorted ORNs, OE CD1 male and female, OE C57BL/6J female. B: Classification of FPKM values of genes. Graph showing the distribution of genes according to FPKM value classes in olfactory tissues. Distribution confirms that the RNA-Seq data comprise the similar number of genes classified into the same range of FPKM values.

(TIF)

Click here for additional data file.

Correlation matrix of the whole data set.

Chart showing the Pearson correlation coefficient values for protein-coding gene expression pattern between all replicates of the OE (upper matrix). In the lower matrix, OBP and Cpy genes were excluded from the analysis. Higher correlation between replicates is indicated by a color scale from blue to orange.

(TIF)

Click here for additional data file.

Correlations of expression levels plotted for each detected gene.

Shown is the correlation of the protein-coding gene expression pattern between OE of two (exemplary chosen) individual mice for each condition (CD1 male, CD1 female and C57BL6 female). Only genes with detectable expression levels (FPKM>0.1) are shown. The FPKM values are logarithmically presented. Genes with the most diverging expression pattern belong to OBPs genes; which are marked in red.

(TIF)

Click here for additional data file.

Correlations of expression levels plotted for each detected OR gene in exemplary two biological replicates of OE from female CD1 mice.

For investigation of the reproducibility of the expression pattern, two biological replicates of the transcriptome of the female OE of CD1 mice were prepared. These new datasets were based on RNA of 8 pooled OE analyzed by mRNA Illumina sequencing on a HiSeq 2000 platform which generated 54–57 million reads (101 bp, paired end). Correlation of the OR gene expression between two biological replicates of female CD1 mice is shown. A detailed analysis of these data will be given elsewhere. Only OR genes with detectable expression levels (FPKM>0.1) are shown. The FPKM values are logarithmically presented. The Pearson correlation coefficient of r = 0.9 confirmed the strong correlation of OR gene expression patterns between biological replicates.

(TIF)

Click here for additional data file.

Correlation matrix of the OR subgenome.

Chart showing the Pearson correlation coefficient values for OR gene expression between all replicates of the OE (n = 13). Higher correlation between replicates is indicated by a color scale from blue to orange.

(TIF)

Click here for additional data file.

Top expressed OR genes.

Chart showing the expression ranking of the top 20 OR genes in OE replicates (CD1 male, CD1 female and C57BL6 mice) and ORNs. The OR genes Olfr1507, Olfr533 and Olfr309 are highly expressed. These receptors can be detected among the 20 most highly expressed OR genes in each replicate of the OE and ORNs.

(TIF)

Click here for additional data file.

Distribution of FPKM values of detected nGPCRs.

Bar chart showing the distribution of FPKM classes in ORNs, OE of CD1 mice (both sexes) and OE of female C57BL/6J.

(TIF)

Click here for additional data file.

In situ hybridization for mRNA of catalytic PI3K subunits.

Expression of transcripts for p110α (1), p110β (2), p110γ (3) and p110δ (4) were detected in the mature ORN cell layer as predicted by the expression in sorted ORNs RNA-Seq data. Sense (a) and antisense (b) RNA probes were tested in parallel and show the antisense specific staining

(TIF)

Click here for additional data file.

Expression level of receptor genes regulating food intake.

Bar chart showing the expression level of receptor genes regulating food intake in the OE. Adipor1 is by far the most highly expressed gene. Receptors for insulin, leptin or ghrelin, are weakly expressed. FPKM values are presented exemplary from CD1 male OE. Adipor1: adiponectin receptor1, Insr: insulin receptor, Paqr3: progestin and adipoQ receptor family member III, Ghrl: ghrelin receptor, Lepr: leptin receptor, Npy1r: neuropeptide Y receptor, Glp1r/Glp2r: Glucagon-like peptide receptor, Hctr1/Hctr2: orexin receptors, Mc3r/Mc4r: melanocortin receptors.

(TIF)

Click here for additional data file.

RNASeq data OE (CD1 male, CD1 and C57BL/6J female) and FACS sorted ORNs.

(XLSX)

Click here for additional data file.

Cuffdiff analysis for differential gene expression between homozygous and heterozygous ORNs datasets.

(XLSX)

Click here for additional data file.

RNASeq data reference tissue: brain, liver, muscle and testis.

(XLSX)

Click here for additional data file.

Raw data sets of RNASeq used in this study.

(DOCX)

Click here for additional data file.

OR genes.

(XLSX)

Click here for additional data file.

OR pseudo-genes.

(XLSX)

Click here for additional data file.

Non-olfactory GPCRs.

(XLSX)

Click here for additional data file.

Non-GPCR membrane proteins.

(XLSX)

Click here for additional data file.

Transporter: channels, SLCs and active transporters.

(XLSX)

Click here for additional data file.

PI3 kinases.

(XLSX)

Click here for additional data file.

Enriched genes in ORNs.

(XLSX)

Click here for additional data file.