*3.1. Transcriptome Data*

A total of six cDNA libraries were constructed in this study. In total, we obtained 23.75, 24.44, 25.65 million reads for male olfactory epithelium and 22.96, 21.70, 22.28 million reads for female olfactory epithelium, respectively. From these, 75.46–77.75% of the reads can be mapped to the annotated regions in zebrafish genome, indicating the good quality of reads (Supplementary Materials Table S1). Overall, we detected a total of 18,658 genes with FPKM > 1 in at least half of the six samples, which were defined as robustly expressed genes. The expression estimates for all annotated transcripts in each replicate were provided in Supplementary Materials Table S2.

We first evaluated the variation in gene expression among the three biologic replicate samples for each sex. We found that the correlation values were highly significant between them all (Figure 1, Pearson correlation coe fficients of at least 0.95, *p*-value < 2.2 × <sup>10</sup>−16). Only very small sets of genes are unusually variable among replicates (Figure 1). In general, the OE transcriptomes from the male and the female zebrafish were highly correlated (Figure 2A, Pearson correlation coe fficient, *r* = 0.89, *p* < 2.2 × <sup>10</sup>−16). We therefore averaged the FPKM values for each gene across each sex. In the olfactory epithelium from both male and female, a few gene are extremely highly expressed. For example, in male 312 most abundant genes account for almost 50% of the fragments and in female 333 most abundant genes account for almost 50% of the fragments obtained from the tissue (Supplementary Materials Table S3). Moreover, a total of 295 genes were found to be shared between male and female most abundant genes. These results were generally consistent with the patterns found in mouse olfactory transcriptomes [40].

#### *3.2. Sex Dimorphism in Zebrafish Olfactory Epithelium Expression Profiles*

To assess whether transcriptional di fferences in the zebrafish olfactory epithelium can account for sex-specific responses to olfactory cues, we examined their sexually dimorphic gene expression profiles. We found that the overall transcriptomes are highly similar between male and female zebrafish olfactory epithelia (Figure 2A,B). Briefly, we detected a total of 713 transcripts showing higher expression level in male olfactory epithelium and 605 transcripts showing higher expression level in female olfactory epithelium, respectively (Figure 2C). However, only 68 and 53 transcripts remained significantly di fferential expression between male and female olfactory epithelium after applying a false discovery rate (FDR) threshold of 0.05 (Figure 2D and Supplementary Materials Table S4). Among these, some genes identified to be di fferentially expressed between sexes are expected to be involved in response to olfactory cues (Table 1), such as gene *OR132-2* (odorant receptor, family H, subfamily 132, member 2), *OTX1* (orthodenticle homeobox 1) and *OR115-10* (odorant receptor, family F, subfamily 115, member 10). However, whether these di fferentially expressed olfactory receptor genes can detect sex relative pheromones needs to be verified in further functional experiments.

**Figure 1.** Correlations of gene expression profile between biologic replicates. The three male and three female samples are listed on the diagonal. Above the diagonal the rho value of the Pearson correlation coefficient is indicated. Below are pairwise comparisons between biologic replicates, shown as scatter plots of the abundances of transcripts(FPKM) expression values for all genes. \*\*\* means *p-*value > 0.95.

**Figure 2.** Sexual dimorphism in zebrafish olfactory system. (**A**) Heatmap of cross-correlations of all samples using differentially expressed transcripts; (**B**) pairwise comparison of gene expression abundances between male and female; (**C**) log2-fold change between male and female samples is plotted against their −log10 FDR; (**D**) heatmap of differentially expressed transcripts identified in zebrafish olfactory epithelium.



To gain further insights into the biologic processes that differed between male and female zebrafish olfactory epithelium, we performed a GO enrichment analysis with the upregulated genes detected in male and female zebrafish olfactory epithelium (Supplementary Materials Table S5). We found that gene ontology biologic process terms associated with response to external stimulus were highly overrepresented among the differentially expressed genes (DEGs) with male-biased expression (Figure 3 and Supplementary Materials Table S5), including response to other organism (GO:0051707), chemotaxis (GO:0006935) and taxis (GO:0042330). However, a significantly different pattern was found in the enriched GO terms for female olfactory epithelium upregulated genes. For example, female olfactory epithelium upregulated genes were enriched in ion import, such as potassium ion import (GO:0010107), inward rectifier potassium channel activity (GO:0005242) and ligand-gated channel activity (GO:0022834). Taken together, these results indicated that there are significant differences in the upregulated genes between male and female zebrafish olfactory epithelium.

**Figure 3.** Gene ontology (GO) terms associated with differentially expressed genes between male and female zebrafish olfactory epithelium.

#### *3.3. Alternative Splicing between the Two Sexes of Zebrafish Olfactory Epithelium*

To examine the influence of sex on alternative splicing regulation in zebrafish olfactory epithelium, alternative splicing events (ASEs) including skipped exon (SE), alternative 5 splice site (A5SS), alternative 3 splice site (A3SS), mutually exclusive exons (MXE) and retained intron (RI) (Figure 4A) between male and female zebrafish olfactory epithelium were identified using rMATS [36].

**Figure 4.** Alternative splicing events (ASEs) between male and female zebrafish olfactory epithelium. (**A**) Schematic representation of five basic types of alternative splicing events. Alternative exons are shown as orange boxes and flanking constitutive exons are shown as blue boxes; (**B**) distribution of isoform numbers for genes in zebrafish genome; (**C**) pie chart showing the percentage distribution of ASEs.

Roughly 46% of genes have 2 detectable isoforms in zebrafish genome (Figure 4B). Among them, we detected a total of 15,193 ASEs, which were distributed in 7585 genes. The most abundant ASEs were skipped exon, accounting for 76.4% of all ASEs, followed by retained intron (7.0%), alternative 3 splice site (6.5%), mutually exclusive exons (6.3%) and alternative 5 splice site (3.9%) (Figure 4C). Using *p* value cutoff of 0.05, differential ASEs were identified between the two sexes of zebrafish (Supplementary Materials Table S6). A total of 86 significantly differential ASEs were identified between the male and female zebrafish olfactory epithelium, which included 50 SE, 4 A5SS, 17 A3SS, 6 MXE, and, 9 RI (Table 2). For example, *maptb*, which is expressed in the developing central nervous system [41], was found to generate different isoforms by sex. *stxbp5l*, which is involved in syntaxin binding, also showed differential isoforms between sexes, which was identified as sexually dimorphic gene between cattle and rat species [42].


**Table 2.** Of selected differential alternative splicing events.

#### *3.4. The Chemosensory Receptor Repertoires*

As in mouse olfactory sensory neurons (OSNs), each of the OSNs in mature zebrafish randomly expressed only one olfactory receptor [43–45]. Thus, the expression level of any given olfactory receptor within the olfactory epithelium will be low. Consistent with this prediction, among the 170 *or*, 126 *taar*, 5*ora*/*V1r* and 57 *olfC*/*V2r* genes in the zebrafish genome, we detected 147 out of 170 *or*, 104 out of 126 *taar*, 4 out of 5*ora*/*V1r* and 57 out of 57 *olfC*/*V2r* genes expressed in either male or female olfactory epithelium with FPKM > 1, confirming that the chemosensory receptors are mainly expressed in the olfactory epithelium (Supplementary Materials Table S7). However, almost all (90%) of these chemosensory receptor genes were found to be very lowly expressed in the zebrafish olfactory epithelium transcriptomes (FPKM < 15), which is consistent with the results from mice [46].

To further assess whether there are differences in chemosensory receptor gene expression levels between the two sexes in zebrafish, we compared the expression levels of chemosensory receptor genes between male and female zebrafish olfactory epithelium transcriptome. From this analysis, we made the following interesting observations. First, the relative receptor abundance level vary greatly between each member, suggesting that each member of the chemosensory receptor genes was expressed asymmetrically. Some members have high expression levels, whereas other members expressed with very low levels (Figure 5). Second, almost all of the chemosensory receptor genes expressed with no differences between male and female zebrafish olfactory epithelium, with only six *or* genes showing differentially expressed between the two sexes (Supplementary Materials Table S7). Third, all of the 6 *or* genes (Table 1) that showed sexual preference were expressed at a higher level in the male olfactory epithelium than in females, which is generally consistent with the results from mice [40,46]. In order to verify the expression results from RNA-seq, we further employed qRT-PCRto measure relative mRNA levels for a total of 12 candidate genes. Our results demonstrated that the expression patterns for these genes were highly consistent between RNA-seq and qRT-PCR (Figure 5E), suggesting the reliability of our RNA-seq datasets.

**Figure 5.** Expression pattern of the chemosensory receptor genes in the zebrafish olfactory epithelium. Mean FPKM expression values across the three samples between male (blue) and female (red) for each of the *or* (**A**), *taar* (**B**), *ora*/*V1r* (**C**) and *olfC*/*V2r* (**D**) genes. Phylogenetic trees were reconstructed using RAxML (version 8.1.17) under the GTRGAMMMAI model with bootstrap support values determined using 1000 replicates. \* means differentially expressed genes. Bootstrap values for basal nodes are provided; (**E**) validation of RNA-seq data using quantitative real-time PCR (qRT-PCR).
