Identification and Expression Analyses of Olfactory Gene Families in the Rice Grasshopper, Oxya chinensis, From Antennal Transcriptomes

The rice grasshopper Oxya chinensis is an important agricultural pest of rice and other gramineous plants. Chemosensory genes are crucial factors in direct interactions with odorants in the olfactory process. Here we identified genes encoding 18 odorant-binding proteins (OBPs), 13 chemosensory proteins (CSPs), 94 olfactory receptors (ORs), 12 ionotropic receptors (IRs), and two sensory neuron membrane proteins (SNMPs) from O. chinensis using an transcriptomic approach. Semi-quantitative RT-PCR assays revealed that six OBP-encoding genes (OchiOBP4, 5, 8, 9, 10, and 14), one CSP gene (OchiOBP10) and two IR genes (OchiIR28 and 29) were exclusively expressed in antennae, suggesting their roles in olfaction. Real-time quantitative PCR analyses revealed that genes expressed exclusively or predominantly in antennae also displayed significant differences in expression levels between males and females. Among the differentially expressed genes, 17 OR-encoding genes, one CSP- and one SNMP-gene showed female-biased expression, suggesting that they may be involved in some female-specific behaviors such as seeking oviposition site; whereas the three remaining OR-encoding genes showed male-biased expression, indicating their possible roles in sensing female sex pheromones. Our results laid a solid foundation for future studies to reveal olfactory mechanisms as well as designing strategies for controlling this rice pest.

Oxya chinensis is an oligophagous grasshopper pest and primarily feeds on graminaceous grasses. Electroantennogram (EAG) and behavioral bioassay have showed that O. chinensis adults have significantly higher olfactory sensitivity to geraniol compared to other host volatiles (Lu et al., 2008). Identification of target genes that are involved in olfactory perception may lead to environmentally-friendly approaches for controlling this pest (Venthur and Zhou, 2018). The objectives of this study were to generate antennal transcriptomes to identify major olfactoryrelated gene families (OBPs, CSPs, ORs, IRs, and SNMPs) from O. chinensis, examine expression profiles of the genes in various tissues to identify the antenna-specific genes, and compare their sex-specific expression patterns to identify female-and malespecific olfactory-related genes.

Insects, Tissue Collections, and RNA Extraction
The colony of O. chinensis used in this study was derived from insects originally collected from rice fields in the environs of Leizhou,Guangdong,China (110 • 05 E,20 • 34 N). The colony was maintained in the laboratory on bristle grass in climatecontrolled chamber under 30 • C ± 1 • C with 60% relative humidity and a photoperiod of 16 h light versus 8 h.
For transcriptomic analyses, 100 pairs of antennae were dissected separately from O. chinensis adult females and males. For semi-quantitative RT-PCR analysis, antennae, maxillary palps, foreleg tarsus, wings, and genitals were separately collected from adults (male: female = 1:1, n = 25 each). Two replicates were included for each tissue. For RT-qPCR analyses, 50 pairs of antennae were obtained from both females and males, separately along with other body parts including 20 pairs of maxillary palps, 20 pairs of foreleg tarsus, 20 pairs of wings, and 10 genitals from both males and females. Three replicates were included for each tissue sample. Tissues were homogenized to powder immediately in liquid nitrogen and stored at −80 • C for further analyses.
Total RNA was extracted using Trizol reagent (Invitrogen, Carlsbad, CA, United States) and treated with RNase-free DNase I (Takara, Dalian, China) to remove potential genomic DNA contamination. RNA concentration, quality and quantity were analyzed on a NanoDrop ND-2000 Spectrophotometer (Nanodrop Technologies, United States) and a Qubit R 2.0 Fluorometer (Invitrogen, Life Technologies, United States). RNA integrity was assessed with an Agilent 2100 Bioanalyzer (Agilent Technologies, United States).
cDNA Library Construction, Illumina Sequencing and de novo Assembly Transcriptomes were generated through a commercial contract with Novogene Bioinformatics Technology, Co., Ltd. (Beijing, China). Libraries for sequencing were generated from 1.5 µg purified RNA per sample using a TruSeq RNA Sample Preparation Kit v2 (Illumina, San Diego, CA, United States) according to Illumina instructions and sequenced on an Illumina HiSeq 2500 platform (San Diego, CA, United States). Approximately 150 bp paired-end reads were generated.
Raw reads were firstly processed through in-house perl scripts. High-quality clean reads were obtained by removing reads containing adapter and unknown (poly-N) and low-quality reads. Clean reads were combined from both female and male antennae and were assembled into unigenes using Trinity (Version: r2013-11-10) with min_kmer_cov set to 2 under default settings (Grabherr et al., 2011). The clean reads from the O. chinensis female antennae and male antennae were deposited in the NCBI Sequence Read Archive (Female antennae: SAMN11484273; Male antennae: SAMN11484274).

Gene Identification
To identify genes coding for OBPs, CSPs, ORs, IRs, and SNMPs, known protein sequences (OBPs, CSPs, ORs, IRs, and SNMPs) from other Orthopteran species were selected as queries to search the O. chinensis antennal transcriptomes. L. migratoria, S. gregaria, O. asiaticus, O. infernalis, and Ceracris kiangsu were used for OBPs; L. migratoria, S. gregaria, and O. asiaticus for CSPs; L. migratoria, S. gregaria, and O. asiaticus for ORs; L. migratoria and O. asiaticus for IRs; S. gregaria and O. asiaticus for SNMPs (Supplementary Table S1). tBLASTn was used to search and identify candidate chemosensory genes against the O. chinensis antennal transcriptomes, with a cutof E-value of 10 −5 . Putative O. chinensis chemosensory genes were in turn used as queries to identify additional genes (tBLASTx and BLASTp). Repetitions were completed until no new candidates were identified. Candidate unigenes from the initial search were further manually examined using BLASTx against Genbank. Open reading frames (ORFs) were predicted with the ORF Finder in NCBI 4 . Amino acid sequences were aligned with MAFFT (version 7.308) (E-INS-I parameter set) (Katoh and Standley, 2013) and visualized with Geneious (version 9.1.3) (Kearse et al., 2012). Transmembrane domains, signal peptides, conserved cysteine locations in candidates were analyzed using the InterProScan tool plug-in in Geneious (Version: 9.1.3.) (Quevillon et al., 2005). Candidate genes coding for OBPs, CSPs, ORs, IRs, and SNMPs were listed in Supplementary Table S2, together with genetic characteristics, best matches in NCBI-nr database, protein domains and estimated expression levels.

Phylogenetic Analyses
Multiple amino acid sequence alignments were made with MAFFT (E-INS-I parameter) (Katoh and Standley, 2013). Phylogenetic trees were constructed using maximum likelihood analyses as implemented in FastTree2 [Jones-Taylor-Thornton (JTT) amino acid substitution model, 1,000 bootstrap replications] (Price et al., 2009(Price et al., , 2010. Dendrograms were created and colored in FigTree 5 . Sequences used for phylogenetic analysis are listed in Supplementary Table S3.

Analyses of Expression Levels and Differentially Expressed Genes
Gene expression levels were estimated using RSEM (Version: 1.2.15) (Li and Dewey, 2011) for each sample. The clean reads were mapped back to the assembled transcriptome. The read counts normalized using the Trimmed Mean of M-value normalization method (Robinson and Oshlack, 2010) for each mapped gene was used to calculate gene expression levels following the FPKM (fragments per kilobase per million read) method (Cock et al., 2010;Trapnell et al., 2010), and were used to identify differentially expressed genes between females and males using DEGseq edgeR package (Version: 3.4.2) . Criteria for estimating significantly differentially expressed genes were set as q-value < 0.005 and the absolute value of log2 Fold Change > 1.

Semi-Quantitative RT-PCR and Real-Time Quantitative PCR
Semi-quantitative RT-PCR was carried out to compare expression levels of genes coding for OBPs, CSPs, and IRs in various tissues including antennae, maxillary palps, foreleg tarsus, wings and genitals. cDNA template was synthesized from total RNA using a PrimeScript RT reagent Kit (Takara, China). Two housekeeping genes, β-actin (Ochiβ-actin) and ribosomal protein 49 (Ochirp49) were used as controls. PCR reactions were conducted in a thermal cycler from Bio-Rad (CA, United States). PCR conditions were as follows: 95 • C for 2 min, followed by 35 cycles of 95 • C for 30 s, 56 • C for 30 s, 72 • C for 1 min, and a final extension for 10 min at 72 • C. PCR products were separated on 1.5% agarose gels. Individual PCR reactions were repeated twice with independently isolated RNA samples.
Expression profiles of the antenna-predominant candidate genes for OBPs, CSPs, and IRs as well as all ORs and SNMPs were analyzed on a LightCycler 480 system (Roche Applied Science). RT-qPCR cycling parameters were set at 95 • C for 5 min, followed by 45 cycles of 95 • C for 10 s, and 60 • C for 20 s. A melting curve analysis was then performed at 95 • C for 20 s, 60 • C for 30 s, and 95 • C for 30 s in order to determine the specificity of primers. Negative controls without cDNA template (nuclease-free water) were included for each reaction. Three independent biological replications were performed with each biological replication measured in three technique replications. The results were analyzed using the LightCycler 480 Gene Scanning Software. Relative quantification was calculated following the 2 − CT method (Livak and Schmittgen, 2001), and normalized against the two reference genes (Ochiβ-actin and Ochirp49). Genespecific primers were designed using Primer3 6 and listed in Supplementary Table S4.

Statistical Analysis
RT-qPCR data were analyzed and plotted using Prism 6.0 (GraphPad Software, CA, United States). The comparative analyses of each target gene among various tissues were determined using a one-way nested analysis of variance (ANOVA) followed by Duncan's new multiple range test (α = 0.05) using SPSS 22.0 (SPSS Inc., Chicago, IL, United States). Values are presented as mean ± SE.

Illumina Sequencing and de novo Assembly
A total of 53.4 and 44.8 million raw reads were obtained from male and female antennae cDNA libraries of O. chinensis, respectively. After removing adaptor sequences, low quality reads and contaminant sequences, 52.2 and 44.7 million of clean reads were generated. The reads were assembled into 120,803 unigenes with a mean length of 1,213 bp, and an N50 of 2,458 bp. The number of unigenes longer than 500 bp was 33,079, which accounted for 35.21% of all unigenes (Supplementary Table S5).

Candidate Genes Coding for Odorant Binding Proteins
Eighteen OBP unigenes were identified from the O. chinensis antennal transcriptomes (Supplementary Table S2: Sheet 1). All identified unigenes except OchiOBP18 had a full-length ORF with sizes ranged from 128 to 266 amino acid residues. Except three OBP genes (OchiOBP12, 14, and 17), all predicted proteins had a predicted signal peptide. Sequence identities of FIGURE 1 | Phylogenetic tree of putative odorant-binding proteins (OBPs). Branch support (circles at the branch nodes) was estimated using an approximate likelihood ratio test based on the scale indicated on the top left. Bars indicate branch lengths in proportion to amino acid substitutions per site. Classic OBPs are in red; plus-C OBPs type-A in green; plus-C OBPs type-B in blue; and atypical OBPs in purple.
Frontiers in Physiology | www.frontiersin.org predicted OBPs with those from other Orthopteran insects in the NCBI-nr database ranged from 45.9 to 88.6%, with an average of 73.8%. Consistent with previous Orthopteran OBP classification by Jiang et al. (2017), twelve OBPs (OchiOBP1, 3-8, 10, 12-14, 16) were classic OBPs, with the typical six conserved C-residues (Supplementary Figure S2). Two (OchiOBP9 and 17) were atypical OBPs, with extraordinary long stretches between conserved C1 and C2 (Supplementary Figure S3). Four OBPs were plus-C OBPs, which were further divided into plus-C OBP type-A (OchiOBP2, 11, and 15) (with the extra C-residue both in C4-C5 and C5-C6) (Supplementary Figure S4) and plus-C OBP type-B (OchiOBP18) (with the extra C-residue in front of C1 and behind C6) (Supplementary Figure S5). Phylogenetic analyses of OBPs from O. chinensis and other five locust species (O. infernalis, O. asiaticus, L. migratoria, S. gregaria, and C. kiangsu) showed that all the locust OBPs formed distinct clades based on the number and position of cysteine residues, and segregated into the classic OBP, atypical OBP, and plus-C OBP sub-families as well as other locust OBPs (Figure 1).

Candidate Genes Coding for Chemosensory Proteins
Thirteen unigenes encoding CSPs were identified from the O. chinensis antennal transcriptomes (Supplementary Table S2: Sheet 2). All but two (OchiCSP11 and 13) had full-length ORFs encoding 120-297 amino acid residues. All predicted FIGURE 2 | Phylogenetic tree of putative chemosensory proteins (CSPs). Branch support (circles at the branch nodes) was estimated using an approximate likelihood ratio test based on the scale indicated at the top left. Bars indicate branch lengths in proportion to amino acid substitutions per site.
Frontiers in Physiology | www.frontiersin.org OchiCSPs contained four highly conserved four-cysteine profiles (Supplementary Figure S6). Except OchiCSP1, all predicted proteins possessed a signal peptide. A phylogenetic tree was built using all predicted OchiCSPs together with those from other locust species. Most OchiCSPs had orthologs with CSPs from other locust species, and no O. chinensis-specific CSP lineage was evident (Figure 2).

Candidate Genes Coding for Odorant Receptors
Ninety-four OR-encoding unigenes were identified, including one Orco (OchiOR1) and 93 conventional OR genes (OchiOR2-94), which were classified to the 7-transmembrane receptor superfamily (Supplementary Table S2: Sheet 1). Among these OR unigenes, 54 had full-length ORFs encoding proteins with 367 to 487 amino acid residues. The highly conserved coreceptor OchiOR1 shared 94.25% identity to a co-receptor from L. migratoria (ALD51504.1), while other OchiORs shared 33.56 to 91.13% identity with average 68.7% with the respective L. migratoria ORs. A phylogenetic tree was generated using our identified ORs along with a data set containing representative ORs from L. migratoria and S. gregaria (Figure 3). The locust ORs formed three distinct clades, with the O. chinensis OR co-receptor (OchiOR1) formed a clade with the protein Orco from two other related species. The vast majority of OchiORs were clustered with orthologs from other locust species. Only a few speciesspecific clades and sister pairs were observed. Three OchiORs (OchiOR27, 53 and 56) were segregated into unique clades, while OchiOR36/60 and OchiOR28/76 formed the sister pairs. FIGURE 3 | Phylogenetic tree of putative odorant receptors (ORs). The distance tree was rooted by the conservative ORco gene orthologs (red). The speciesspecific clades and sister pairs are labeled with red dots. Branch support (circles at the branch nodes) was estimated using an approximate likelihood ratio test based on the scale indicated at the top left. Bars indicate branch lengths in proportion to amino acid substitutions per site. The clade I is in orange; clade II in green; and clade III in blue.

Candidate Genes Coding for Sensory Neuron Membrane Proteins
Two SNMP unigenes were identified, which matched the CD36 family. The genes contained complete ORFs that encoded proteins with two transmembrane domains (Supplementary Table S2: Sheet 3). Phylogenetic analyses revealed that all SNMPs were classified into two distinct subgroups, SNMP1 and SNMP2 (Figure 5). In the SNMP1 clade, OchiSNMP1 and OchiSNMP2 were clustered into two subclades along with those from other locust species.

Tissue-and Sex-Specific Expression
Based on semi-quantitative RT-PCR analyses, OchiOBP4,5,8,9,10, and 14 were almost exclusively expressed in antennae, while OchiOBP3 and 17 were abundant in antennae and maxillary palps ( Figure 6A). The remaining OBP-encoding genes were abundant in multiple body parts. For CSP-encoding genes, OchiCSP10 was exclusively expressed in antennae, while other OchiCSPs were present multiple body parts ( Figure 6B). For iGluR/IR-encoding genes, OchiIR28 and 29 were expressed predominately in antennae ( Figure 6C).

DISCUSSION
Identification and characterization of chemosensory genes are important steps toward understanding the evolution and primary functions of the insect olfactory system. In this study, we identified a set of candidate chemosensory genes from O. chinensis via analyzing antennal transcriptomes. Genetic and phylogenetic analyses of candidate chemosensory genes in O. chinensis was carried out to examine the similarities and differences of molecular components in chemosensory pathways. We further analyzed the expression profiles of chemosensory genes to identify olfaction-specific genes for future functional studies.
Odorant-binding proteins play vital roles in carrying odorants through the hemolymph to OR neurons, and transducing the resultant signals to downstream effectors in the olfactory system (Hallem et al., 2006). The number of OBPs identified in our antennal transcriptome was comparable to its related Locust species, namely 18 from O. infernalis, 16 from L. migratoria, 15 from O. asiaticus, and 14 from S. gregaria) (Ban et al., 2003;Xu et al., 2009;  Yu et al., 2009;Zhang et al., 2015Zhang et al., , 2018Jiang et al., 2017). This may reflect physiological and evolutionary consistency between O. chinensis and other Locust species. OBPs from O. chinensis can be classified into four types: classic, atypical, plus-C (type-A and B), as well as Locust-specific OBPs . Expression analysis revealed that five classic OBP (OchiOBP4,5,8,10,and 14) and one atypical OBP (OchiOBP9) were expressed exclusively in antennae, suggesting that these genes are likely to have a role in antennal chemicalrecognition processes. CSPs and OBPs can act as carriers for odorant molecules (Liu et al., 2012;Yi et al., 2013;Zhang T. et al., 2013;Zhang et al., 2014). However, in our RT-PCR analysis, only OchiCSP10 showed significantly antennae-biased expression, suggesting that this CSP could be involved in olfactory function.
Olfactory receptors, which connect binding proteins and OSNs to transduce olfactory signals, are the best known group of insect chemoreceptors. Like other locust species, O. chinensis possesses an expanded OR family (94) compared with Dipterans, Lepidopterans, and Hemipterans (Vieira and Rozas, 2011). A large number of OR-encoding genes may be associated with its ability for diverse host-odor perception to feed on a wider range of host plants. Great variation in the number of OR genes was found in different Locust species. Specifically, 94 OR genes were identified in O. chinensis, less than 120 from S. gregaria , and 141 from L. migratoria , but significantly more than 60 from O. asiaticus (Zhou et al., 2019). The variation in gene numbers may reflect that O. chinensis inhabits a different ecological niche than other locust species. Generally, sexually dimorphic expression of ORs in antennae indicates possible pheromone receptors contributing sexual behaviors. Typically, Lepidopteran sex pheromones are released via females to attract males. Several moth sex pheromone for ORs have been functionally characterized, and most of them are expressed at higher levels in the male antennae (Krieger et al., 2004;Wanner et al., 2010;Zhang and Löfstedt, 2013). ORs expressed predominantly in female antennae are predicted to function in the detection of egg-laying-related odorant (Pelletier et al., 2010) or male released pheromones (or signals) (Anderson et al., 2009); ORs expressed evenly in the male and female antennae are predicted to function in general odorant perception (Yan et al., 2015). Therefore, we hypothesize that some or all of female-predominant OchiORs, including OchiOR2,6,8,9,11,15,33,38,39,41,43,44,48,50,53,56,and 62), are likely involved in female specific behaviors such as finding plant hosts for oviposition. 46,and 86) may be associated with detecting female sex pheromones. The remaining ORs with roughly equal expression levels in female and male antennae might be dedicated to general odorant detection.
We also identified two SNMP unigenes. SNMPs are conserved throughout holometabolous insects and the SNMP1 subfamily are usually expressed in pheromonesensitive OSNs, and mediate responses to lipid pheromones (Jin et al., 2008;Nichols and Vogt, 2008;Vogt et al., 2009;Gomez-Diaz et al., 2016). All SNMPs from locust species were clustered into the SNMP1 subclade. Both SNMPs were predominantly expressed in antennae in O. chinensis, supporting that these two SNMPs may perform functions in pheromone perception. Interestingly, OchiSNMP1 was predominantly expressed in female antennae (Figure 7B), suggesting that OchiSNMP1 may participate in sex pheromone reception.

DATA AVAILABILITY STATEMENT
The datasets generated for this study can be found in the clean reads from the O. chinensis female antennae and male antennae were deposited in the NCBI Sequence Read Archive (Female antennae: SAMN11484273; Male antennae: SAMN11484274).

AUTHOR CONTRIBUTIONS
YC and CK performed the experiments. ZW analyzed the data. ZW and JL wrote and revised the manuscript.

ACKNOWLEDGMENTS
We thank Dr. Mingshun Chen (Kansas State University, United States) for comments and editorial assistance on the manuscript.