Identification of Immunity-Related Genes in Dialeurodes citri against Entomopathogenic Fungus Lecanicillium attenuatum by RNA-Seq Analysis

Dialeurodes citri is a major pest in citrus producing areas, and large-scale outbreaks have occurred increasingly often in recent years. Lecanicillium attenuatum is an important entomopathogenic fungus that can parasitize and kill D. citri. We separated the fungus from corpses of D. citri larvae. However, the sound immune defense system of pests makes infection by an entomopathogenic fungus difficult. Here we used RNA sequencing technology (RNA-Seq) to build a transcriptome database for D. citri and performed digital gene expression profiling to screen genes that act in the immune defense of D. citri larvae infected with a pathogenic fungus. De novo assembly generated 84,733 unigenes with mean length of 772 nt. All unigenes were searched against GO, Nr, Swiss-Prot, COG, and KEGG databases and a total of 28,190 (33.3%) unigenes were annotated. We identified 129 immunity-related unigenes in transcriptome database that were related to pattern recognition receptors, information transduction factors and response factors. From the digital gene expression profile, we identified 441 unigenes that were differentially expressed in D. citri infected with L. attenuatum. Through calculated Log2Ratio values, we identified genes for which fold changes in expression were obvious, including cuticle protein, vitellogenin, cathepsin, prophenoloxidase, clip-domain serine protease, lysozyme, and others. Subsequent quantitative real-time polymerase chain reaction analysis verified the results. The identified genes may serve as target genes for microbial control of D. citri.


Introduction
Dialeurodes citri Asahmaed (D. citri) belongs to the family Aleyrodidae, of the order Hemiptera. This important, widely spreading pest is found in the world's citrus producing areas. It originates from Southeast Asia and has been found in Asia, South and North America, and Europe [1,2]. D. citri has a piercing/sucking mouthpart, enabling them to suck the branches and leaves of citrus, which will cause leaf and fruit abscission. More seriously, D. citri can secrete honeydew, which can adsorb dust in the air and provide nutrition for some fungi. This leads to the occurrence of sooty mold and seriously affects photosynthesis and the quality of citrus fruits. In the state of Florida in the USA, D. citri was once the dominant pest for citrus trees [3]. In Oceania, D. citri was detected for the first time in New Zealand in 2000 [4], and it spread swiftly to other citrus-producing regions, leading to a reduction in output by about 90% [5]. The citrus-producing regions in Chongqing, China, such as Changshou, Tongnan, Tongliang, and Nanchuan, are at risk for large-scale outbreaks of D. citri, which bring much trouble to growers. At present, the main method for preventing D. citri infection is chemical control. However, it is difficult for a chemical pesticide to effectively control an outbreak of D. citri infection in a field. In addition, excessively frequent chemical spraying may lead to environmental pollution, pesticide residue, and drug resistance [6][7][8]. Therefore, an environmentally, friendly and effective control strategy for D. citri is needed.
However, it is difficult for pathogenic fungi to infect D. citri. Pests, similar to those higher animals, have a complete immune defense system, which decreases greatly the lethality of pathogenic fungi. It is well known that insects possess three major lines of immune defense: the body wall, cell-mediated immunity, and humoral immunity. The body wall of insects is their first barrier against invaders. The main functional parts of the body wall are chitin, protein, wax, and some recognition factors. The cellular immunity of the insect relies mainly on the ability of the hemocytes in the body to enwrap and devour the invaders and antigens. In addition, hemocytes are also involved in wound healing and blood coagulation, acting to prevent pathogenic microorganisms from entering the body through a wound [18]. Gillespie et al. found that, after treatment by the Metarhizium acridum, the desert locust shows a sharp rise in the number of hemocytes and begins, playing a rold in the immune response [11]. This means the hemocytes respond rapidly, playing an important role in cellular immunity. The insect humoral immune response includes melanization, lysozyme, antimicrobial peptides (AMPs), lectins, antiviral factors (AvFs), and proteinase inhibitors (Pis). In case of microorganism infection, the insect can activate the cellular and humoral immune responses through four steps [19]: 1) identify the invaders. The insect has special pattern recognition receptors (PRRs), such as peptidoglycan recognition proteins (PGRPs), thioester-containing proteins (TEPs), gramnegative binding proteins (GNBPs), scavenger receptors (SCRs), C-type lectins (CTLs), and galectin (GALE), that can identify the external pathogenic microorganisms, leading to a downstream immune response [20]; 2) an extracellular cascade reaction to activate serine proteases and remove serine protease inhibitor; 3) amplification of the signal of infection or removal of the false alarms; 4) stimulation of the transcription of the effecting factor through a signal transduction pathway to produce the immune response. The AMPs are mainly subjected to the Toll pathway, Imd pathway, and Jas/Stat signal transduction pathway.
At present, genes and proteins related to the immune response is mainly studied in a model insect. However, there have been few studies on the identification of genes in D. citri that are responsible for the immune response to L. attenuatum infection. Our analyses of the transcriptome of D. citri (no reference genome) and differential expression of genes in the digital gene expression (DGE) profile provide information on immunity-related genes of D. citri and a theoretical basis for research into the molecular mechanism of the immunity of D. citri against pathogenic fungi. With gene function annotation, it is possible to find the major genes responsible for the defense of D.citri against L. attenuatum, and thus, the limitations and blindness of research on a single gene are reduced. This research will benefit the biological control of insects using pathogenic fungi.

Materials and Methods
Fungus culture and conidia suspension preparation L. attenuatum strain TL001 was cultured on potato dextrose agar (PDA) plates at 25°C and 80% humidity. Conidia (spores) used for infection were harvested from 3-4 weeks old cultures by scraping the surface of the mycelia with sterile cell scrapers into sterile deionized water containing 0.1% Tween-80. Conidia were separated from other mycelial structures over a sterile funnel packed with autoclaved glass wool, washed twice with ddH 2 O by centrifugation at 4,000 rpm, counted, and diluted to 1×10 8 spores/ml. Freshly prepared conidia were used for all experiments.
Collection of D. citri and RNA extraction D. citri eggs, nymphs of all ages, and adults were collected in the net house of D. citri and transferred to a 1.5-ml EP tube in equal proportions. The EP tube was numbered #1, frozen in liquid nitrogen, and stored at -80°C.
The citrus leaves to which the D. citri larvae were attached were collected and washed with clear water. These leaves were separated into two groups. One group was treated with the spore suspension for 3 days, whereas the other group, as the control group, was sprayed with sterile water. The successful infection of D. citri larvae with L. attenuatum was confirmed under a dissecting microscope by observation that the attached spores on the cuticle of larvae had germinated and enwrapped pests' bodies. D. citri larvae in the control and treatment groups were picked and transferred to two 1.5-ml EP tubes, numbered #2 and #3, respectively. The contents of the EP tubes were frozen with liquid nitrogen and stored at -80°C.
Using the RNA Isolater Total RNA Extraction Reagent (Vazyme, China) method, total RNA was extracted from samples #1, #2 and #3. The genomic DNA in the total RNA was removed using DNase I (Invitrogen, USA). The quality of the extracted RNA was evaluated through 1% agarose gel electrophoresis, and the concentration was determined with an Nanodrop 2000N spectrophotometer (Thermo Fisher Scientific, USA), recorded as optical density ratios OD260/OD280 and OD260/OD230.

Construction and sequencing of cDNA library
The RNA (#1) with a poly(A) tail was purified with magnetic beads containing polynucleotide T and ruptured with buffer solution. Enriched poly(A) RNA of each sample was fragmented into 200-700 nt pieces with RNA Fragmentation Reagents. Then, the first-strand cDNA was synthesized with reverse transcriptase and arbitrary primers. Next, the second-strand cDNA was synthesized with DNA polymerase I and RNA enzyme. The double-stranded DNA was modified with Klenow fragment and T4 polynucleotide kinase successively. A base A was added on the 3'!5' exonuclease of the Klenow fragment. Then the base was connected to the corresponding adaptor using T4 ligase. Finally, fragments with around 200bp length were purified with QiaQuick GelPurify Kit (Qiagen, Hilden, Germany), and used as templates for PCR amplification to create the cDNA library. The library was paired-end sequenced using PE100 strategy on Illumina HiSeqTM 2500 (Illumina, San Diego, CA, USA) in Biomarker Technologies (Beijing, China).

Assembly and annotation of transcriptomes
Trinity (http://trinityrnaseq.sourceforge.net/) software was applied to perform the de novo assembly for the filtered high-quality data. The used parameters were as follow: min_glue = 2, V = 10, edge-thr = 0.05, min_kmer_cov = 2, kmer size = 25, path_reinforcement_distance = 80, and group_pairs_distance = 250. The other parameters were set as the default. The data samples were merged and assembled, and contigs were obtained through overlap of the assembled sequences. Then the contigs were clustered according to the paired-end information of sequences and similarity of contigs. Local assembly was conducted to generate transcripts. The longest transcript in each local region was selected for use as a unigene. The following parameters were used to ensure a high quality of assembly: a minimum of 95% identity, a minimum of 35 overlapping bases, a minimum of 35 scores and a maximum of 25 unmatched overhanging bases at sequence ends. The consensus cluster sequences and singletons make up the final unigene dataset. For functional annotations, we first searched all unigene sequences against various protein databases such as Nr, SwissProt, COG, and KEGG using BLASTX, and then searched nucleotide database Nt using BLASTN, with an E-value cut-off of 10 −5 . For inconsistencies between the unigene alignment results in different databases, priority was given to the nr, Swiss-Prot, KEGG, and COG data sequentially. Unigenes with no matched data in any of these databases were analyzed with ESTScan software to predict the coding region and direction of the sequence.

Analysis on the digital gene expression profile
According to data in the constructed cDNA library, differences in mRNA expression between samples #2 and #3 were identified using the digital gene expression (DGE) profile. Gene expression was quantified as reads per kb per million reads (RPKM) [21]. According on Audic et al.'s detection method for differentially expressed genes based on sequencing [22], a strict algorithm was adopted to screen the differentially expressed genes. Those with a false discovery rate (FDR) 0.001 and |fold change| !2 were considered to be differentially expressed. Combined with the function annotation of differentially expressed genes, the pathway (KEGG) and gene ontology (GO) enrichment analysis as well as the pattern clustering of differentially expressed genes were performed. Fisher's exact test was applied in the enrichment analysis. The results were adjusted using the Bonferroni correction method. Thus, the pathways with obviously enriched differentially expressed genes and GO functional categories were obtained for further analysis.
qRT-PCR veritification on DGE D. citri larvae of the same age were collected and sprayed with the spore suspension of L. attenuatum at the concentration of 1×10 8 spores/mL. The infection period lasted 5 days. A proper amount of D. citri was picked out every day, frozen with liquid nitrogen, and preserved at -80°C. D. citri larvae treated with sterile water were used as the control. The primers for differentially expressed genes and reference genes were designed using Primer Premier 5 software (S1 Table). The total RNA from the six treatments above was extracted separately, and 1 μg was used in the reverse transcription with HiScript II Q RT SuperMix for qPCR with gDNA wiper (Vazyme, China), followed by real-time fluorescent quantitative PCR with ChamQ SYBR qPCR Master Mix (Vazyme). The reaction was carried out on an ABI 7500 Fast Real-Time PCR System (Applied Biosystems, USA). Relative gene expression was calculated using the Pfaffl method [23]. Both D. citri α-tubulin and β-actin were used as reference genes. The analytic software for qPCR was 7500 softwear v2.0.6. SPSS Statistics v19.0.0 software was used to perform independent sample t-tests (P<0.05).

Transcriptome sequencing and assembly
Through transcriptome sequencing (runs accession number: SRR2980521), we found a total of 24,071,734 reads. The mean CycleQ20 value of the samples reached 100.00%, and the base Q30 was 81.41%, suggesting that the sequencing was reliable. All high-quality reads were assembled de novo into 2,318,371 contigs (Table 1) with a mean length of 82 nt. Contigs were clustered according to the paired-end information of sequences and similarity of contigs and then assembled into 119,428 transcripts with the mean length of 952 nt. Among them, 31,748 transcripts (26.59%) had a length greater than 1,000 nt. These transcripts were finally assembled into 84,733 unigenes, with a mean length of 772 nt. Among them, 16,602 unigenes (19.6%) had a length greater than 1,000 nt. The N50 lengths of the contigs, transcripts, and unigenes were 98 nt, 1,809 nt and 1,320 nt, respectively. Regarding the length distribution of transcripts and unigenes of D. citri, most were 200-300 nt, followed by 300-500 nt, with the smallest group being sequences longer than 2,000 nt.

Function annotation of unigenes
Among the 84,733 assembled unigenes, 28,190 unigenes had been annotated, as found through the sequence alignment on the NCBI (National Center for Biotechnology Information database) website. Specifically, 11,516 unigenes were annotated in the COG database, 17,081 in the GO database, 19,537 in the Swiss-Prot database, 9,387 in the KEGG database, and 27,746 in the Nr database (Table 2).
It was observed in GO analysis that 17,081 unigenes were named as 148,314 GO terms, of which most genes included more than one GO term. GO analysis was used mainly to predict the functions of D. citri proteins. Generally, GO terms are classified into three categories ( Fig  1): biological process, molecular function, and cellular component. In our analysis, most of the GO terms participated in biological process (90,520, accounting for 61.03% of the total), followed by molecular function (32,827,22.13%), and cellular component (24,967,16.83%). The three largest sub-categories of GO terms were cellular process (10,798 GO terms) and metabolic process (10,357 GO terms) in biological process and binding (8,959 GO terms) in the molecular function. We also used COG classifications to analyze the putative protein functions. In total, 11,516 unigenes were functionally classified into 25 COG categories (Fig 2). The largest category was "General function prediction only" (3,976, 35%), followed by "Translation, ribosomal structure and biogenesis" (1,134, 9.9%), "Carbohydrate transport and metabolism" (1,097, 9.5%), "Replication, recombination and repair" (1,067, 9.3%), "Amino acid transport, metabolism" (1,057, 9.2%). There were only a few of unigenes taking part in "Nuclear structure" (7 unigenes), "Cell motility" (20 unigenes), "RNA processing and modification" (81 unigenes), and "Chromatin structure and dynamics" (94 unigenes). No unigenes took part in Extracellular structures.
There were 9,387 unigenes that had been annotated and classified in KEGG database, distributed over 208 pathways (Fig 3). The numbers of unigenes involved in different pathways varied. Those pathways with less than 100 unigenes were classified into one category (others), and the remaining 3,416 pieces of unigenes were classified into 19 categories. The ribosome pathway (ko03010) was the largest category with 456 pieces of unigenes. The pathways of protein processing in the endoplasmic reticulum (ko04141, 398 pieces), spliceosome (ko03040, 358 pieces), and RNA transportation (ko03013, 336 pieces) also involved many unigenes.  Immunity-Related Genes Expressed in D. citri Infected by Lecanicillium attenuatum

Identification of immunity-related genes
By searching the transcriptome, we preliminarily identified 129 pieces of immunity-related unigenes, involving pattern recognition receptors, signal transduction factors, response factors, and so on (Table 3). These genes play important roles in endogenous and exogenous immunoreactions of D. citri. For example, prophenoloxidase (PPO) can participate in melanization; lysozyme can be capable of dissolving cell walls of fungi and Gram-positive bacteria; and the antimicrobial peptide and clip domain serine protease also perform functions related to humoral immunity when a pathogen invasion occurs. In terms of the pattern recognition receptors, we totally identified a total of 49 unigenes including 1 PGRP, 6 GNBPs, 19 βGRPs, 2 TEPs, 8 SCRs, 3 CTLs, 5 GALEs, 3 hemocytins, and 2 integrins. These genes play an important role for the insect to recognize external pathogens, triggering the downstream reaction. PGRPs mainly identify the special ingredient in the bacterial cell walls-peptidoglycan-and then trigger the transcription of antibacterial peptide or PPO activation cascade. The protein was originally found in the silkworm, and until now, 12 PGRP genes have been identified [24]. However, in the transcriptome of D. citri, we identified only one PGRP gene (PGRP-S2). In dealing with the infection of gram negative bacteria and fungi, GNBP/βGRP can identify and combine with the beta-1,3-glucan and trigger PPO cascade. We performed phylogenetic analyses on the complete six genes, DcβGRP- 5,7,8,16,17,19,in the sequence of D. citri, and found they have a close genetic relationship with AtβGRP-4a and MhGNBP-2 in Hemiptera (Fig 4). Four signal transduction pathways, Toll, Imd, JNK, and JAK/STAT are known to be involved ininsect immunity [25]. However, in D. citri, we only identified unigenes related to the Toll pathway, including 16 Tolls, 1 ECSIT, 6 Pelle, 1 pellino, and 3 NF-κB genes, and did not find genetic information in the other three pathways. We perform phylogenetic analysis on the Toll protein and found DcToll-7 has high similarity with other species (Fig 5). In the cascade mediated by serine protease, D. citri has 19 clip-domain serine proteases (clip-domain SP), 2 serine protease homologs (SPHs), and 6 serpins. Among them, DcSPH-1 and DcSPH-2 are prophenoloxidase-activating factors, which can catalyze conversion of PPO into PO (phenoloxidase). Serpin can negatively control the PPO activation, competing for PPO with SP.
As immune factors of D. citri, we identified 18 PPO unigenes, 6 lysozyme unigenes, and 3 AMPs. Most of the insects contained two kinds of PPOs, in the form of dimmers. Through clustering analysis, we found DcPPO-6 and DcPPO-7 with PPO of Hemiptera, and DcPPO1 and DcPPO-17 had a more distant relationship with other species (Fig 6). This may suggest that D. citri has more than two kinds of PPO genes. AMP includes cecropin, defensin, attacin, gloverin, and so on. D. citri has 3 AMP genes with sequences that are too short to permit identification of the type of AMP.

Expression profile analysis
We obtained the RPKM values for unigenes, and through further calculation of Log 2 fold changes [Log 2 Ratio(treatment group RPKM/control group RPKM)], we found 441 differentially expressed unigenes in control group and treatment group, among which 313 unigenes had been annotated (Data Accessibility) in the transcriptome database. The clustering analysis of differentially expressed genes suggested that Log 2 FC values of most genes in the treatment group were positive (Fig 7A). Namely, these genes were upregulated in comparison to the control group. It was found in statistical analysis of differentially expressed genes that the genes with Log 2 FC values in the range of -2 to 2 were in the majority (210 pieces, 47.62%; Fig 7B). Genes with Log 2 Ratio values in the range of 2 to 4, ranked second with 148 pieces (33.56%).
Generally, only small numbers of differentially expressed genes were downregulated. There were 38 genes (8.62%) with Log 2 FC values of -4 to -2, and 9 genes (2.04%) with Log 2 FC values of -6 to -4. In our two experimental groups, only one gene was upregulated by more than 10-fold, and its gene ID was unigene78104_c1. The gene ID of the most downregulated gene (by -6.38-fold) was unigene80382_c2.
We analyzed 24 genes that may participate in the response against infection by L. attenuatum (Table 4). Among them, function annotation identified that the cuticle protein was the highest upregulated gene with a Log 2 Ratio value of 10.04, followed by the vitellogenin gene. This implies that the two proteins may play roles in the response of D. citri against infection by L. attenuatum directly or indirectly. Regarding the melanization, we found three PPO genes and three clip-domain SP genes. Among them, the expression difference of DcSP-3 was 2.51-fold, and those of the two PPO genes were 1-2-fold. In addition, we found two lysozyme unigenes with fold changes of 2.47 and 2.60. Regarding the metabolic reaction process of insects, we examined the expression of cathepsin B.
Pathway analysis of upregulated genes of interest suggested that the differentially expressed genes were mainly involved in two pathways, the lysosome (ko04142) and MAPK signaling pathways (ko04010). In the lysosome pathway, 9 pieces of unigenes exhibited differential expression, with 8 upregulated and 1 downregulated. All of these genes are crucial to lysosome phagocytosis and immunoreactions. The MAPK signaling pathway plays the role of signal transduction in the processes of stress adaptation and inflammatory response. For this pathway, there were two pieces of differentially expressed genes, both of which were related to heat shock cognate 70. Immunity-Related Genes Expressed in D. citri Infected by Lecanicillium attenuatum

qRT-PCR analysis of differentially expressed genes
To verify the accuracy of fold changes in the digital gene expression profile and investigate the dynamic tendency of key gene expression with infection by L. attenuatum, we used qRT-PCR to analyze key genes in the response of D. citri to infection for 1-5 days. qRT-PCR analysis was conducted for the two pieces of cuticle protein unigenes most upregulated. After 1 day of infection, no significant change (P<0.05) in their expression was observed compared with that in the control group (CG) treated with clear water. However, by the second day of infection, DcCP-1 exhibited increased expression (Fig 8). The relative expression change after 2 days was 2.20-fold, and that after 3 days was 3.60-fold. Unexpectedly, the gene expression increased sharply on the fourth day, with 267.83-fold upregulation and 865.54-fold upregulation on the fifth day. The gene DcCP-2 exhibited obvious downregulation (P<0.05) on the second day and rebounded on the third day. On the fourth day, expression of the gene increased sharply as for DcCP-1, to 92.87-fold, and further increased to 202.63-fold on the fifth day.
We verified the two upregulated vitellogenin genes, DcVg-1 and DcVg-2, and found that they showed the same expression tendency (Fig 8). During the first 2 days of infection, their expression was inhibited compared to that in the control group (P<0.05). On the third day, their expression levels were comparable to those in the control group, and after infection for 4 days, obvious upregulation of the two genes occurred (215.90-fold and 199.98-fold respectively). On Immunity-Related Genes Expressed in D. citri Infected by Lecanicillium attenuatum the fifth day, these genes were upregulated by 327.19-fold and 287.84-fold, respectively. Thus, expression of the two vitellogenin genes was positively correlated with the infection time.
Among the genes participating in the lysosome pathway, two cathepsin B genes, DcCatB-1 and DcCatB-2, were analyzed through qPCR. DcCatB-1 was upregulated significantly (P<0.05) by 2.17-fold by the first day of infection and continued to be upregulated on days 2-5, with a peak fold change value of 3.69-fold on the second day (Fig 8). Expression of DcCatB-2 was inhibited on the first day, with a 0.43-fold difference from that in the control group. However, it was upregulated obviously on the second day and reached maximum upregulation by 3.93-fold on the third day. On the fourth and fifth days, its expression decreased, which might be caused by depletion of the enzyme in the latter stage.
For the melanization pathway, we tested the expression levels of DcPPO-7 and DcSP-3 ( Fig  8). The results indicated that these two genes were significantly upregulated (P<0.05) on the first day of infection by 4.40-fold and 4.12-fold, respectively, and by 40.88-fold and 35.71-fold, respectively, by the fifth day. During melanization, serine protease can activate PPO and convert it to PO, which participates in the generation of melanin. This is an immediate immune response to the invasion of fungi. The two genes exhibited obvious upregulation on the first day of infection, which supports this theory.
We measured the expression level of the lysozyme gene DcLys-4 (Fig 8) and found that, like PPO, it was significantly upregulated (P<0.05) by 2.81-fold on the first day of infection, 17.37-fold on the third day, and 27.00-fold on the fifth day.
Heat shock cognate 70 (DcH70-1, Fig 8) acts in the MAPK signaling pathway and responds to exogenous stimuli. Expression of this gene did not change significantly in the first 2 days of Immunity-Related Genes Expressed in D. citri Infected by Lecanicillium attenuatum infection. But during the last 3 days, it was significantly upregulated (P<0.05) with a 7.34-fold change on the fifth day.
In the expression profile, the hypothetical protein gene DcHP-1 also showed a large fold change of 8.13-fold on the fourth day and by 22.82-fold (Fig 8) on the fifth day. Further investigation of the importance of this gene is warranted.

Discussion
We obtained a transcriptome database containing 84,733 unigenes of D. citri using RNA-Seq technology. This amount of data was greater than that obtained previously [26]. By searching and screening the transcriptome database, we identified 129 immunity-related unigenes. These unigenes are related to pattern recognition receptors, information transduction factors and response factors. Using the digital gene expression profile, we identified 441 differentially expressed genes in D. citri infected with L. attenuatum. Among these genes, we found the response factors such as the PPO, lysozyme, and clip-domain SP participated in cascade. However, we did not find genes related to the pattern recognition receptors, suggesting that as the final effector, response factors show fluctuating expression The cuticle is the first barrier for insects to defend against pathogen infection, in addition to being indispensible for maintaining the shape and mobility of insects [27,28]. The major Table 4. Differentially expressed unigenes likely involved in the antifungal response of D. citri after treatment with L. attenuatum. Immunity-Related Genes Expressed in D. citri Infected by Lecanicillium attenuatum components of the insect cuticle are chitin and cuticle protein [29]. Cuticle protein contributes much to the stress resistance, drug resistance, and immunity of insects. When an insect is suffering from adverse environment conditions, cuticle protein genes are induced to strengthen or stabilize the cuticular structure, resist the effects of adverse factors, and maintain the insect's survival [30][31][32]. In researching the aphid's insecticide resistance mechanism, Silva et al. found that two RR2-type cuticle protein genes are upregulated [33]. Asano et al. found that when B. mori larvae are subjected to bacterial infection, the cuticle protein gene BmCb10 is significantly upregulated. They speculated that the gene could transmit the exogenous adverse stimulation to activate melanization [34]. He et al. proposed that cuticle protein may play a role in wound healing in Anopheles gambiae larvae and adults [35]. In this study, two cuticle protein genes showed significant upregulation upon D. citri infection by L. attenuatum. This finding provides some evidence for the role of cuticle protein in the immune defense of D. citri. Melanization, an important aspect of the insect defense system, involves the regulation of the melanin cascade mediated by PPO [36]. Upon pathogen invasion, it activates PPO and transforms it into PO, which can transform phenolic substances into quinone intermediates, that aggregate and form melanin before enclosing, isolating, and killing pathogens. In addition, PO takes part in the processes of wound healing and skin hardening. The activation and transformation of PPO into PO is generally considered to be completed through the cascade of clip-domain SPs [37]. Gillespie et al. found that when Schistocerca gregaria was infected by Metarhizium anisoplia, the PPO level in the body increased while PO activity decreased, and the lysozyme level exhibited a significant decrease in comparison with the control group [11]. By knocking out the PPO gene in the mosquito Armigeres subalbatus, it was found that the melanization function was influenced greatly [38]. When fungal spores were injected into the insect body, the PO expression level was increased significantly [39]. We measured the expression levels of clip-domain SP and PPO in D. citri infected by L. attenuatum. On the first day of infection, both were significantly upregulated, and the fold changes in expression exceeded 35-fold by the fifth day. This indicates that these enzymes contribute much to the response against exogenous pathogenic fungi.

Gene Name Unigene ID CG-RPKM TG-RPKM Log 2 Ratio (TG/CG) Up-Down FDR Nr_annotation
In the Digital Gene Expression Profiling, vitellogenin and cathepsin B with high expression levels received our attention. vitellogenin is an important source of energy for D. citri, and cathepsin B is an indispensable metabolic enzyme of lysosome. However, research on their roles in the immune response remains incomplete. Guo et al. showed that when Bemisia tabaci were fed a virus-infected plant, vitellogenin expression increased significantly [40]. Soderhall et al. demonstrated that a clotting protein belonging to the vitellogenin superfamily participates in the autologous immune defense in freshwater crayfish [41]. Raikhel et al. found that the upstream regulatory region of the vitellogenin gene of Aedes aegypti participates in the immune defense against pathogens [42]. Shi et al. demonstrated that vitellogenin can agglutinate erythrocytes of toad and chicken and has an inhibitory effect on various bacteria [43]. Amdma et al. showed that vitellogenin participates in the regulation of immune function and life in bees [44]. Cathepsin B enzymes can part in the processes of immune evasion [45,46]. Futahashi et al. found upregulated expression of the cathepsin gene in many tissues of Burkholderiainfected Riptortus pedestris [47]. Zhang et al. cloned the cathepsin O gene of B. mori and detected its expression in the hemolymph of B. mori treated with Escherichia coli, they speculated that the enzyme participates in the immune response in the body [48]. Kocks et al. found that after immune stimulation, cysteine protease L was activated and highly expressed in lysosomes [49]. Wang et al. investigated two Cathepsin L genes of channel catfish and proved that the two genes acted in mucosal immunity [50].
In future studies, we may take advantage of the RNAi technique to verify the effect of these genes in the response of D. citri to infection of entomopathogenic fungi, and improve the fungal toxicity and field application effectiveness.