Identification of Olfactory Genes in Monochamus saltuarius and Effects of Bursaphelenchus xylophilus Infestation on Their Expression

The pine wood nematode (PWN) Bursaphelenchus xylophilus has caused disastrous losses of pine forests in many countries, and the success of PWN depends strongly on interactions with its insect vectors. Monochamus saltuarius is a newly recorded vector in Northeast China. Feeding (i.e., immature) and egg-laying (i.e., mature) Monochamus spp. target different host plants, and olfactory cues play important roles regarding host choice. Whether infestation with PWN affects olfactory mechanisms in M. saltuarius related to feeding and oviposition is of interest as this may affect the spread of nematodes to new healthy hosts. However, little is known about molecular mechanisms of the olfactory system of M. saltuarius. We identified chemosensory-related genes in adult M. saltuarius and examined the influence of B. xylophilus on the respective expression patterns. Fifty-three odorant-binding proteins (OBPs), 15 chemosensory proteins, 15 olfactory receptors (ORs), 10 gustatory receptors, 22 ionotropic receptors (IRs), and two sensory neuron membrane proteins were identified, and sex bias among non-infested beetles was mainly found with respect to expression of OBPs. Interestingly, OBPs and ORs were markedly down-regulated in male M. saltuarius infested with B. xylophilus, which may reduce olfactory sensitivity of male M. saltuarius and affect the spreading of B. xylophilus to new hosts. Our results will help understand the interactions between B. xylophilus and M. saltuarius, which may lead to the identification of new control targets in the olfactory system of M. saltuarius.


Introduction
Pine is the common name of Pinus species (family Pinaceae), and it is the largest genus of conifers, which is widely distributed throughout the world [1]. Pines are not only pioneer species in afforestation, but also have very high economic value owing to their by-products such as pine resin, rosin, and turpentine, which are important raw materials for producing medicinal compounds, organic solvents, and fuels. Currently, pine forests are threatened by pine wood nematodes (PWNs; Bursaphelenchus xylophilus) (Steiner and Buhrer, 1934) Nickle (1970) [2]. More than 60 pine species can be infected with B. xylophilus, and vast areas of pine forests in several countries have been devastated as a result of pine wilt disease. In China alone, an estimated 50 million trees were killed in an area of 80,000 ha between 1982 and 2000 [3].
In 1982, PWN was first found in Nanjing, China, and in only a few decades, it spread rapidly from Nanjing to vast areas of China (Li and Zhang, 2018). According to the M. saltuarius at the mating stage were collected on the Dahuofang Experimental Forest Farm, Fushun City, Shenyang Province, China, in May 2019. Samples (n = 40) were collected in the forest. We searched M. saltuarius carefully in the area where some of the host trees were weak, and once mating couples were found (a male and a female M. saltuarius were mating, usually located on the trunk), they would be put it in a breathable plastic tube and taken back to the laboratory. Antennae were collected as follows: both antennae of each M. saltuarius were cut off at the base, placed in a marked centrifuge tube, and were immediately frozen in liquid nitrogen. Then, the respective donor individual was examined for nematodes using the Baermann funnel method. Four categories were used, i.e., males and females without nematodes (here termed 'healthy') and with nematodes ('infested'), respectively. Five biological replicates were used per category, and each of the five biological replicates consisted of antennae coming from one beetle.

RNA-Seq Library Preparation and Sequencing
Total RNA was isolated from M. saltuarius antennae using TRIzol reagent (Invitrogen, Carlsbad, CA, USA) [39,40], and genomic DNA was removed using DNase I (TaKara, Dalian, Liaoning, China). Integrity and purity of the extracted RNA were determined using a 2100 Bioanalyzer (Agilent Technologies, Inc., Santa Clara, CA, USA), and RNA was quantified using a ND-2000 device (NanoDrop Thermo Scientific, Wilmington, DE, USA).

Annotation
All assembled transcripts were searched against six public databases (NR, Swiss-Prot, Pfam, COG, GO, and KEGG) using BLAST to identify proteins with the highest sequence similarity with the produced transcripts to produce functional annotations (cutoff e-value < 1 −5 ). BLAST2GO (http://www.blast2go.com/b2ghome, accessed on 25 November 2019) [45,46]. Software was used to produce gene ontology (GO) annotations of unique assembled transcripts to describe biological processes, molecular functions, and cellular components [47,48]. Metabolic pathway analysis was performed using the Kyoto Encyclopedia of Genes and Genomes (KEGG; http://www.genome.jp/kegg/, accessed on 26 November 2019).
Chemosensory genes were identified in two ways: first, contig tBLASTx searches were performed using an M. saltuarius assembled unigene database with query sequences as chemosensory genes from other insects, such as Drosophila melanogaster, Bombyx mori, and other long-horned beetles. Second, annotation information of the unigenes was used. Open reading frames of the candidate genes from the two ways were identified and used to as confirmation through BLAST (http://blast.ncbi.nlm.nih.gov/Blast.cgi, accessed on 6 January 2022).

Gene Expression Quantification
The expression level of each transcript was calculated according to the fragments per kilobase of exon per million mapped reads (FRKM) method [49]. The R software package EdgeR vision 2.12 (http://www.bioconductor.org/packages/2.12/bioc/html/ edgeR.html, accessed on 10 October 2019) was used for differential expression analysis [50]. A compatible-hits-norm model was used to normalize unigene expression levels and DEGs differentially expressed genes [51].

Phylogenetic Analyses
IQ-TREE web service was used to construct phylogenetic trees of the predicted protein sequences and orthologous genes from other insects [52]. The sequences were first aligned using BioEdit software with the ClustalW multiple alignment tool. Bootstrap phylogenetic trees were generated using 1000 replicates, and the best models were selected using the model fit result of IQ-TREE. iTOL was used to annotate the phylogenetic trees [53].

M. saltuarius Antennae Transcriptomes
More than 7 Gbp clean sequence data were produced from most samples (except for one sample from an infested male, which produced 6.98 Gb data), and the Q30 value was > 94%. De-novo transcriptome assembly resulted in 38,198 contigs with an N50 of 1798 bp, and the largest unigene was 32,017 bp long. The length distribution of the transcriptome assembly is shown in Figure S1. The mapping rates of clean data were >78.16% for each sample. Gene annotation based on the NR, Swiss-Prot, Pfam, COG, GO, and KEGG databases yielded annotations for 20,931transcripts.

OBPs of M. saltuarius with and without B. xylophilus
We analyzed sex-specific expression patterns of OBPs in non-B. xylophilus-infested M. saltuarius, and 15 OBPs were expressed differentially between sexes. Most differentially expressed OBPs (14/15) showed higher expression in females ( Figure 1A). Among OBPs with female-biased expression in healthy beetles, OBP28, OBP52, OBP38, and OBP39 had high expression levels, and the other female-biased OBPs (OBP32, OBP18, OBP19, OBP31, OBP43, OBP3, OBP6, OBP7, OBP40, and OBP51) showed low expression. OBP14 expression was male-biased among healthy beetles, and, the expression of OBP14 was the highest in males and very low in females ( Figure 1A). Thus, OBP14 may be related to sex pheromone recognition in male M. saltuarius. Further studies are needed to understand the reason of its differential expression between males and females. Effects of B. xylophilus on the expression patterns of OBPs in male and female beetles was examined, and none of the OBPs showed altered expression levels in females after one month of infestation with B. xylophilus ( Figure 1B). However, expressions of five OBPs (including OBP14, OBP30, OBP52, OBP44, and OBP35) were significantly lower in infested male beetles than in noninfested ones ( Figure 1C). Thus, OBP expression in male M. saltuarius was affected by infestation with B. xylophilus.
ined, and none of the OBPs showed altered expression levels in females after one month of infestation with B. xylophilus ( Figure 1B). However, expressions of five OBPs (including OBP14, OBP30, OBP52, OBP44, and OBP35) were significantly lower in infested male beetles than in non-infested ones ( Figure 1C). Thus, OBP expression in male M. saltuarius was affected by infestation with B. xylophilus. A phylogenetic tree was constructed to analyze the characteristics of OBPs from M. saltuarius, based on orthologs to OBP sequences from several other Coleoptera insects A phylogenetic tree was constructed to analyze the characteristics of OBPs from M. saltuarius, based on orthologs to OBP sequences from several other Coleoptera insects Semanotus bifasciatus [54], Tribolium castaneum [55], Anomala corpulenta [56], Agrilus planipennis, and Dendroctonus ponderosae [57] (Figure 2). OBPs from M. saltuarius seems to be homologous to OBPs of Coleopteran insects whose genomic information are available, and indicates the relative complete identification of OBPs from M. saltuarius. MsalOBP14, which expressed at an especially high level in male adults, seems not to be a typical PBPs. We did not get a special PBP clade in Figure 2  Semanotus bifasciatus [54], Tribolium castaneum [55], Anomala corpulenta [56], Agrilus planipennis, and Dendroctonus ponderosae [57] (Figure 2). OBPs from M. saltuarius seems to be homologous to OBPs of Coleopteran insects whose genomic information are available, and indicates the relative complete identification of OBPs from M. saltuarius. MsalOBP14, which expressed at an especially high level in male adults, seems not to be a typical PBPs. We did not get a special PBP clade in Figure 2 (MsalOBP14 was clustered together with other M. saltuarius OBPs) and MsalOBP14 was not crusted in PBP clade in another phylogenetic tree containing more model insect species ( Figure S2). However, Pfam search hits PBP/GOBP family (PF01395.24, E-value is 4.8 × 10 -18 ) (http://pfam.xfam.org, 2 February 2019). The function of MsalOBP14 deserves further study.

CSPs of M. saltuarius with and without B. xylophilus
Sex-specific expression of CSPs was assessed in healthy M. saltuarius ( Figure 3A), and three CPSs were expressed differently between sexes; two of these (CSP1 and CSP13) showed higher expression in females, and one (CSP5) showed extremely high expression in males. Infestation with B. xylophilus slightly affected CSP expression, only CSP13 in female ( Figure 3B), and CSP9 and CSP13 male ( Figure 3C) were down-regulated in the infested beetles (please see supplementary Table S2). Sex-specific expression of CSPs was assessed in healthy M. saltuarius ( Figure 3A), and three CPSs were expressed differently between sexes; two of these (CSP1 and CSP13) showed higher expression in females, and one (CSP5) showed extremely high expression in males. Infestation with B. xylophilus slightly affected CSP expression, only CSP13 in female ( Figure 3B), and CSP9 and CSP13 male ( Figure 3C) were down-regulated in the infested beetles (please see supplementary Table S2).  A phylogenetic tree of CSPs was constructed based on orthologs to CPS sequences from M. saltuarius, S. bifasciatus, T. castaneum, A. corpulenta, A. planipennis, and D. ponderosae (Figure 4). The identified CSPs from M. saltuarius were homologous to CSPs of Coleopteran insects whose genomic information were available, indicating relative complete identification of CSPs from M. saltuarius. Transcript levels are expressed as reads per kilobase of exon per million mapped reads (RPKM), and the expression data were shown as mean±SE. The stars above each bar of figure indicate significant differences (p < 0.05 * and p < 0.05 **).
A phylogenetic tree of CSPs was constructed based on orthologs to CPS sequences from M. saltuarius, S. bifasciatus, T. castaneum, A. corpulenta, A. planipennis, and D. ponderosae ( Figure 4). The identified CSPs from M. saltuarius were homologous to CSPs of Coleopteran insects whose genomic information were available, indicating relative complete identification of CSPs from M. saltuarius.

ORs of M. saltuarius with and without B. xylophilus
No sex bias was observed regarding OR expression in healthy M. saltuarius ( Figure 5A). Infestation with B. xylophilus did not affect OR expression in female M. saltuarius ( Figure 5B). Marked changes in OR expression were detected in male M. saltuarius carrying B. xylophilus compared with non-infested individuals ( Figure 5C). OR1, OR5, OR10, OR12, OR14, and OR16 showed higher expression in healthy M. saltuarius than those that were infested. Thus, B. xylophilus can reduce olfactory sensitivity of male M. saltuarius.
No sex bias was observed regarding OR expression in healthy M. saltuarius ( Figure  5A). Infestation with B. xylophilus did not affect OR expression in female M. saltuarius ( Figure 5B). Marked changes in OR expression were detected in male M. saltuarius carrying B. xylophilus compared with non-infested individuals ( Figure 5C). OR1, OR5, OR10, OR12, OR14, and OR16 showed higher expression in healthy M. saltuarius than those that were infested. Thus, B. xylophilus can reduce olfactory sensitivity of male M. saltuarius. A phylogenetic analysis of ORs was performed based on the OR sequence orthologs from M. saltuarius, S. bifasciatus, A. corpulenta, A. planipennis, and D. ponderosae ( Figure 6).
Generally, almost all M. saltuarius ORs found orthology to other Coleoptera ORs, including the insect OR obligatory co-receptor, Orco. However, the total number of M. saltuarius ORs may have not been recovered due to its integrity in this transcriptomic study. This is obvious after comparison with genomic data from A. planipennis and D. ponderosae.
arius. (A) Differences between sexes of healthy M. saltuarius; (B) Differences between healthy and infested female M. saltuarius; (C) Differences between healthy and infested male M. saltuarius. Transcript levels are expressed as reads per kilobase of exon per million mapped reads (RPKM), and the expression data were shown as mean±SE. The stars above each bar of figure indicate significant differences (p < 0.05 * and p < 0.05 **).
A phylogenetic analysis of ORs was performed based on the OR sequence orthologs from M. saltuarius, S. bifasciatus, A. corpulenta, A. planipennis, and D. ponderosae ( Figure 6). Generally, almost all M. saltuarius ORs found orthology to other Coleoptera ORs, including the insect OR obligatory co-receptor, Orco. However, the total number of M. saltuarius ORs may have not been recovered due to its integrity in this transcriptomic study. This is obvious after comparison with genomic data from A. planipennis and D. ponderosae.

IRs, GRs, and SNMPs in Healthy and Infested M. saltuarius
Regarding IRs, IR10 and IR12 were expressed higher in healthy females compared to healthy males, IR4 expressed higher in healthy than infested females, while IR25 and IR21 were expressed higher in healthy than infested males (Figure 7(A1-A3)). Regarding GRs, no sex bias was observed in healthy adults, and only GR8 showed higher expression in female M. saltuarius infested with B. xylophilus (Figure 7(B1-B3)). GR8 was the only up-regulated chemosensory gene in M. saltuarius infested with B. xylophilus. Regarding SNMP, no sex bias was observed in healthy M. saltuarius and SNMP2 showed lower expression in infested males than in healthy ones (Figure 7(C1-C3)).

Discussion
Olfactory cues play important roles in the host choice of immature and mature M. saltuarius. The spread of nematodes to new healthy hosts mainly relies on the feeding activity of immature M. saltuarius. If infestation with B. xylophilus for approximately one month can affect olfactory recognition in longicorn beetles it could also affect the nematode's spreading interval duration. In the present study, we identified 53 OBPs, 15CSPs, 15 ORs, 10 GRs, 22 IRs, and 2 SNMPs in M. saltuarius, a new vector of B. xylophilus in northern China. Sex-specific expression patterns in healthy beetles occurred mainly regarding OBPs. Interestingly, we found that OBPs and ORs were down-regulated in male M. saltuarius infested with B. xylophilus, which may reduce olfactory sensitivity and delay their maturation; consequently, this may prolong the time period, during which B. xylophilus can be transmitted through infested beetles.
The numbers of OBPs and CSPs in M. saltuarius are consistent with other Coleopteran insects whose genomic information are available [55,57], indicating relative complete identification of these two chemosensory families in M. saltuarius. However, the number of OR genes in M. saltuarius was lower than the number of other Coleopteran insects, this may be due to some missing genes down-regulated in the adults, especially in sexually mature adults. Another reason maybe that these insects have mated and some of the chemosensory genes were down-regulated. However, S. bifasciatus ORs were also identified from antennae transcriptome data, and its numbers were much higher than M. saltuarius [54]. To elucidate the underlying causation, we reviewed olfactory gene identifications of other longicorn beetles. To date, olfactory-related gene identification and respective studies have only been performed on a few longicorn beetles. Fifty-two OBPs and 132 ORs were identified in the genome of Anoplophora glabripennis [58], and 57 ORs were identified in Megacyllene caryae [59]. However, 31 OBPs and only nine ORs were identified in the antenna transcriptome of M. alternatus [3]. Both Monochamus species possess considerably fewer ORs than other longicorn beetles, which may be a particular characteristic of this genus. However, further studies on the olfactory system of Monochamus spp. based on genome data are needed to confirm this prediction.
The present study is the first to confirm that infestation with B. xylophilus can influence olfactory gene expression in the vector M. saltuarius. Interactions between B. xylophilus and their vector insects are complex throughout all stages of development. Before accumulation in its vector, B. xylophilus can secrete ascarosides to promote beetle pupation. After developing into the adult stage, the beetle secretes ascarosides to attract the B. xylophilus 4th instar larvae [8]. B. xylophilus 4th instar larvae development relies on C16 and C18 fatty acid ethyl esters, which are abundantly produced by beetle pupae at emergence and newly emerged adults [60]. Thus, B. xylophilus can promote and use signals from vector beetles to facilitate its accumulation on its vector before the newly emerged beetles leave their old host plants. When the beetles disperse for foraging, the time span for B. xylophilus to spread on healthy trees is critical; a longer foraging time increases the probability that B. xylophilus will spread to healthy trees. Our results indicate that B. xylophilus can reduce the olfactory sensitivity of male M. saltuarius, which may delay sexual maturity and further prolong the beetle's feeding activity, thus increasing the likelihood of spreading PWNs to new healthy hosts. Conversely, the olfactory system of female M. saltuarius was not affected by longtime B. xylophilus infestation, which may be a compromise for B. xylophilus as it requires female M. saltuarius to lay eggs, thereby initiating a new spreading cycle. However, further research is required to test whether GR8 upregulation in female M. saltuarius infested with B. xylophilus (Figure 7(B1-B3)) can facilitate feeding activity in female beetles and therefore enhance egg maturity and efficient oviposition. Further studies should be undertaken to test whether GR up regulation in infested mated females may help discriminate adequate hosts for oviposition and thereby facilitate the spread of nematode to damaged trees.
On the other hand, the blocking of chemosensation in infested M. saltuarius may be for the self-preservation of the insect. Sexually mature males need to find dead logs and trees where females are waiting for mating. Chemosensory gene down-regulation in males would reduce their sensitivity and influence mating success, and thus transfer less nematodes to its host trees. Anyway, whether these hypotheses are correct requires further physiological and gene function researches, and clarification of the complex relationships between nematodes-vectors-trees will imagine new pest control strategies.
Genes affected by B. xylophilus infestation may be used as new molecular targets for M. saltuarius control in the future. However, the functions of the regulated ORs and OBPs genes in male M. saltuarius are unclear, and further research is needed to elucidate the respective interactions between nematodes and their beetle vectors. This work extends our understanding of the complex relationship between nematodes-vectors-trees; this is critical for us to grasp new control targets and develop new control methods to protect large areas of pine forest.
Conclusion: firstly, we identified 53 OBPs, 15 CSPs, 15 ORs, 10 GRs, 22 IRs, and two SNMPs in M. saltuarius. Secondly, sex expression bias of the chemosensory genes was analyzed. Finally, we found that B. xylophilus infestation mainly down-regulated the OBPs and ORs expression levels in male M. saltuarius at the mating stage.