Comparative Transcriptome Analysis of Two Populations of Dastarcus helophoroides (Fairmaire) under High Temperature Stress

: The differentially expressed genes (DEGs), key genes and metabolic pathways of the parasitic beetle, Dastarcus helophoroides (Fairmaire), were compared between the ﬁftieth commercially reared population and the ﬁrst natural population to reveal the adaptive mechanism in response to high temperature stress. The high-throughput sequencing technique was employed for transcriptome sequencing of two populations of D. helophoroides . In total, 47,763 non-redundant transcripts with the average length of 989.31 bp and the N50 of 1607 bp were obtained. Under high temperature stress, 1108 DEGs were found in the commercial population; while there were 3946 DEGs in the natural population, which were higher than those in the commercial population (3.56 times). High temperature stress of D. helophoroides promoted the expression of heat shock proteins (HSPs) and metabolism-related genes in both populations, but metabolism synthesis and hydrolysis of natural population was much higher, allowing them to produce more resistant substances (such as HSPs, superoxide dismutase (SOD), peroxiredoxin (Prx), etc.). Therefore, HSPs may play a major role in the high temperature adaptation of a commercial population, while the natural population probably respond to heat stress with more resistant substances (such as HSPs, SOD, Prx, etc.). These results provide a reference to select and domesticate a speciﬁc ecotype with stronger adaptability to the high temperature weather in the forest and further improve the efﬁciency of D. helophoroides as a bio-control factor. ME—membrane, MEL—


Introduction
Insects are ectotherms and unable to regulate their body temperature; consequently, their physiological and biochemical characteristics and life activities are significantly affected by high temperature stress [1][2][3]. The capability of a parasitic insect to attack its host depends largely on its environmental adaptability and in particular to that of increased temperature [4]. Currently, there are some successful examples in the mass rearing and commercial production of several parasitoid species, such as Trichogramma brassicae, T. dendrolimi, Diglyphus isaea, Fopius arisanus (Sonan) and Dichasmimorpha longicaudata (Ashmead), which are widely used to control some crops and forest pests [5][6][7]. However, the tolerance of commercial natural enemy insects to high temperature stress is reduced due to long-term rearing under constant temperature condition. For example, when mass rearing populations of predatory mites were released into the field, they were often unable to quickly adapt to natural high temperature stress and other variable environmental factors, leading to reduced control efficiency or even failure [8]. On the other hand, if the natural enemy insects have been exposed to extreme environmental conditions for a long time, and their adaptability has been strengthened, they can acquire high level of resistance under long term stress [9]. In particular, this paper explores the molecular differences between the commercially reared and natural populations of Dastarcus helophoroides (Fairmaire) (Coleoptera: Bothrideridae) in response to high temperature stress.
D. helophoroides is currently one of the dominant natural enemy species attacking old larvae or pupae of the Japanese pine sawyer beetle, Monochamus alternatus (Coleoptera: Cerambycidae) [10][11][12], which is a vector of pine wood nematodes (Bursaphelenchus xylophilus), causing millions of hectares of dead pine trees in China. Fortunately, a commercialized population of D. helophoroides has been widely used in southern China, effectively suppressing the population of longhorn beetles in the forests [13]. However, recent studies have revealed that their parasitic efficiency decreased in southern pine forests, especially in those forests containing dying Pinus massoniana Lamb caused by the long horned beetle [14,15]. Previous studies have shown that D. helophoroides adults roam free in the forests [16], and are often stressed by several ecological factors such as extreme temperature [17][18][19]. In particular, it has been shown that the ontogeny and reproductive characteristics of adults are adversely affected by high temperature conditions [20]. For example, heat shock protein genes such as HSP70 and sHSPs can respond to heat stress and play an important role under high temperature conditions (26-44 • C) [21], indicating that heat stress can affect and regulate the biological and physiological processes of D. helophoroides, which affects the population's stability and pest control effectiveness. In recent years, high temperatures have frequently occurred in southern China. For example, the average daily temperature from July to August at Fujian in 2020 was 35.0 • C, and the highest was 43.4 • C for 4 h, which continued for 8 days successively (https: //tianqi.2345.com/wea_history/58847.htm) (Accessed on 21 December 2021), and a new highest temperature for Canada was as high as 49.6 • C recorded on 29 June 2021 at Lytton, BC, Canada (https://www.currentresults.com/Weather-Extremes/Canada/hottest.php) (Accessed on 21 December 2021). Obviously, high temperature stress is one of the threshold factors which may determine the survival of commercial population once released in the forest. On the other hand, there are natural population of D. helophoroides in Fujian forests, which indicates that they have the better adaptability to high temperature stress than the commercial population of D. helophoroides.
In this study, we compared and analyzed the key genes and metabolic pathways of commercial and natural populations of D. helophoroides under high temperature (39 • C) stress by using transcriptome sequencing and bioinformatics methods, which will be helpful to reveal the adaptation mechanism to high temperature stress, and populations differentiation as well.

Insects and Samples
The fiftieth generation of commercial population (Hunan population, HN)  Twenty healthy adults (containing 10 males and 10 females) were selected from both populations and reared under 39 • C (numbered HN-39 and YC-39, respectively) [22], and adults reared under 25 • C as the control (numbered HN-CK and YC-CK, respectively). All samples of D. helophoroides were quickly taken out and put into liquid nitrogen after 3 h, and then transferred to −80 • C refrigerator for RNA extraction. Three biological replications were performed with each treatment, and twenty adults were used in each biological replicate.

RNA Extraction and Illumina Sequencing
Total RNA was extracted from adults of D. helophoroides based on the instructions of Trizol reagent (Invitrogen, Carlsbad, CA, USA), and the quality of RNA was detected by 1% agarose gel electrophoresis. The purity, concentration and integrity were detected by Nanodrop (IMPLEN, Palo Alto, CA, USA), Qubit2.0 (Life Technologies, Carlsbad, CA, USA) and Agilent 2100 (Agilent Technologies, Palo Alto, CA, USA), respectively. These RNA samples were independently used as templates for first-strand cDNA synthesis with reverse transcriptase. The cDNA library construction and sequencing of total RNA were commissioned by Beijing Biomarker Technologies Co. (Beijing, China), and transcriptome sequencing was performed on the Illumina NovaSeq 6000 platform to generate 150 bp paired-end reads.

Sequence Assembly and Gene Annotation
The original sequencing data were filtered to remove the joint sequences, ploy-A and low-quality sequences to obtain high quality clean read sequences. Trinity software (Version: release-20131110) was used to complete sequence splicing and obtain the Unigenes sequence [23]. BLAST comparison was performed with protein databases Nr, Swiss-Prot, Pfam, COG, KOG, KEGG and GO to obtain information annotations. According to different functions, the annotated Unigenes sequences were divided into different GO functions and KEGG metabolic pathways [24].

Gene Expression Quantification
Bowtie2 (Version 2.2.5) software was used to map all the valid reading orders, and RSEM (Version 2.2.5) software was used to estimate the expression level [25]. FPKM value (Fragments Per KB Per Million Fragments) was used to evaluate the expression abundance of the corresponding Unigene. DEGseq method was used to compare the expression differences of transcripts between two different populations under high temperature stress and normal temperature respectively, and the DEGs sets between the two samples were obtained [26]. The error detection rate (FDR < 0.001) and |Log2(FC)| >1 were used as screening criteria to evaluate the key indicators of DEGs.

RNA Sequencing and Assembly
The total transcriptome sequencing data was 26.69 Gb, of which the GC Content was about 42.95-44.14%, and Q30 base percentage was at least 91.24%. The total and average length of these unigenes were 47,252,603 bp and 989.31 bp, respectively, and the N50 was 1607 bp (Table 1). Unigene length was mostly distributed in the range of 200 to 2000 bp ( Figure 1). In general, the results of sequencing and assembly of the transcriptome library are of a quality that meets the requirements of further bioinformatics analysis.

Quantitative Analysis of Differentially Expressed Genes
There were 1108 and 3946 DEGs after heat shock treatment in commercial and natural populations, respectively, with a criterion of FDR < 0.001 and |Log2(FC)| >1. Among them, there were 800 up-regulated genes and 308 down-regulated genes in commercial population. The up-regulated and down-regulated genes in the natural population were 3.13 and 4.67 times as many as those in the commercial population ( Figure 3), indicating that high temperature increased the expression levels of most genes in D. helophoroides.
The Venn diagrams show that 672 unique DEGs and 3510 unique DEGs were detected after heat shock treatment in commercial and natural populations, respectively, and 436 common DEGs were detected within both populations ( Figure 4).  The Venn diagrams show that 672 unique DEGs and 3510 unique DEGs were detected after heat shock treatment in commercial and natural populations, respectively, and 436 common DEGs were detected within both populations ( Figure 4).

GO and KEGG Analysis of Differentially Expressed Genes
The GO data of the commercial population response to heat stress were mainly divided into three categories, namely, biological processes, cellular components, and molecular functions ( Figure 5A). For the up-regulated genes, cellular processes (47.40%), cell components (31.49%) and catalytic activity (52.94%) were enriched in each of three categories, respectively. For the down-regulated genes, metabolic processes (46.00%), membrane (27.00%), catalytic activity (57.00%) were enriched in each of three categories, respectively. In the KEGG analysis, endoplasmic reticulum protein synthesis (14.35%), phagosome (8.33%), longevity regulation pathway (6.02%) and cell cycle (2.78%) were most up-regulated components and processes; and the lysosome components accounted for the highest portion (19.64%) among the down-regulated genes ( Figure 6(A1,A2)).

GO and KEGG Analysis of Differentially Expressed Genes
The GO data of the commercial population response to heat stress were mainly divided into three categories, namely, biological processes, cellular components, and molecular functions ( Figure 5A). For the up-regulated genes, cellular processes (47.40%), cell components (31.49%) and catalytic activity (52.94%) were enriched in each of three categories, respectively. For the down-regulated genes, metabolic processes (46.00%), membrane (27.00%), catalytic activity (57.00%) were enriched in each of three categories, respectively. In the KEGG analysis, endoplasmic reticulum protein synthesis (14.35%), phagosome (8.33%), longevity regulation pathway (6.02%) and cell cycle (2.78%) were most up-regulated components and processes; and the lysosome components accounted for the highest portion (19.64%) among the down-regulated genes ( Figure 6(A1,A2)).

Analysis of Differentially Expressed Genes
The DEGs of commercial population to heat shock obtain 348 genes which functi was known (FDR < 0.001, |Log2(FC)| >1, the same below) by Nr database comparison, which 285 genes were up-regulated and 63 genes were down-regulated. Among the u

Analysis of Differentially Expressed Genes
The DEGs of commercial population to heat shock obtain 348 genes which function was known (FDR < 0.001, |Log2(FC)| >1, the same below) by Nr database comparison, of which 285 genes were up-regulated and 63 genes were down-regulated. Among the up-regulated genes, the number of genes related to metabolic expression and heat shock proteins (HSPs) were the largest (3.51%). Among the metabolic expressions, functional genes related to protein and nucleic acid metabolism were dominant (2.11%), followed by genes related to lipid metabolism and glycolysis (both 0.71%). Carbohydrate metabolism, toxic factors, signal transduction and energy transport were increased in down-regulated expression ( Figure 7A). 892 genes with known functions were obtained from the natural population, which was 2.56 times as those genes from the commercial population. Among them, 565 genes were up-regulated and 327 genes were down-regulated, which were 1.98 and 5.19 times of data from the commercial population, respectively. The most up-regulated genes were signal transduction and energy transport, with 146 genes (25.84%). The following were 129 metabolic expressions (22.83%) and 43 toxic factors (7.61%). Among the metabolismrelated genes, functional genes related to protein and nucleic acid metabolism were mainly expressed (5.84%), followed by genes related to lipid metabolism and glycolysis (both 3.01%) ( Figure 7B). The results showed that the two populations had similar differential expression gene in metabolic expression in response to the same high temperature stress. Genes related to signal transduction and energy transport are differentially expressed.
The DEGs of two populations exposed to 39 • C were compared by Nr database, and 893 genes with known functions were obtained. Among them, 287 genes were upregulated. The most up-regulated genes were metabolic expression (62 genes, accounting for 21.60%), signal transduction and energy transport (22 genes, accounting for 7.67%), and reproduction related genes (16 genes, accounting for 5.57%). At the same time, the number of down-regulated genes was higher than that of up-regulated genes, among which metabolic expression was the most abundant (29.54%), and the metabolic expression was mainly related to protein and nucleic acid metabolism (11.72%). This was followed by lipid metabolism (3.30%), other carbohydrate metabolism (2.81%), and glycolysis (1.65%) ( Figure 7C). regulated genes, the number of genes related to metabolic expression and heat shock proteins (HSPs) were the largest (3.51%). Among the metabolic expressions, functional genes related to protein and nucleic acid metabolism were dominant (2.11%), followed by genes related to lipid metabolism and glycolysis (both 0.71%). Carbohydrate metabolism, toxic factors, signal transduction and energy transport were increased in down-regulated expression ( Figure 7A). 892 genes with known functions were obtained from the natural population, which was 2.56 times as those genes from the commercial population. Among them, 565 genes were up-regulated and 327 genes were down-regulated, which were 1.98 and 5.19 times of data from the commercial population, respectively. The most up-regulated genes were signal transduction and energy transport, with 146 genes (25.84%). The following were 129 metabolic expressions (22.83%) and 43 toxic factors (7.61%). Among the metabolismrelated genes, functional genes related to protein and nucleic acid metabolism were mainly expressed (5.84%), followed by genes related to lipid metabolism and glycolysis (both 3.01%) ( Figure 7B). The results showed that the two populations had similar differential expression gene in metabolic expression in response to the same high temperature stress. Genes related to signal transduction and energy transport are differentially expressed.
The DEGs of two populations exposed to 39 °C were compared by Nr database, and 893 genes with known functions were obtained. Among them, 287 genes were up-regulated. The most up-regulated genes were metabolic expression (62 genes, accounting for 21.60%), signal transduction and energy transport (22 genes, accounting for 7.67%), and reproduction related genes (16 genes, accounting for 5.57%). At the same time, the number of down-regulated genes was higher than that of up-regulated genes, among which metabolic expression was the most abundant (29.54%), and the metabolic expression was mainly related to protein and nucleic acid metabolism (11.72%). This was followed by lipid metabolism (3.30%), other carbohydrate metabolism (2.81%), and glycolysis (1.65%) ( Figure 7C).

Discussion
Temperature is one of the most important factors affecting the distribution of insects worldwide [27]. It is of great significance to study the tolerance mechanism of natural enemy to high temperature stress for evaluation of stability and colonization. A previous study indicated that transcriptome sequencing of D. helophoroides, obtained 2.48 Gb of nucleotides, which laid a foundation for revealing the mechanism of longevity and sustainable reproduction at the molecular level [28]. In comparison, our work has benefited from the update of sequencing technology, and obtained more abundant data (24.21 Gb). In addition, this study was designed to reveal the differential adaptation to higher temperature stress, which is conducive for scientific evaluation of the survival and pest control efficiency of commercially reared populations in the south China during summer.

Discussion
Temperature is one of the most important factors affecting the distribution of insects worldwide [27]. It is of great significance to study the tolerance mechanism of natural enemy to high temperature stress for evaluation of stability and colonization. A previous study indicated that transcriptome sequencing of D. helophoroides, obtained 2.48 Gb of nucleotides, which laid a foundation for revealing the mechanism of longevity and sustainable reproduction at the molecular level [28]. In comparison, our work has benefited from the update of sequencing technology, and obtained more abundant data (24.21 Gb). In addition, this study was designed to reveal the differential adaptation to higher temperature stress, which is conducive for scientific evaluation of the survival and pest control efficiency of commercially reared populations in the south China during summer.
Under high temperature stress, 1108 DEGs (800 up-regulated and 308 down-regulated) were found in the commercial population of D. helophoroides. In the natural population, there were 3946 DEGs (2507 up-regulated and 1439 down-regulated). Obviously, the number of DEGs from the commercial population was lower than that from the natural population, which might be due to the stable living environment temperature, relatively simple food resources and reduced gene exchange during long-term indoor rearing. This would likely lead to population degradation and the decrease of adaptive ability at the level of molecular expression. The complex and variable wild environment enables the natural population to be continuously trained and domesticated, and the gene exchange is frequent, leading to a high level of differential expression [29][30][31]. GO enrichment analysis showed that both commercial and natural populations were enriched in catalytic activity under high temperature stress. We assumed that the catalytic activities of molecular functions are involved in the regulation of various resistance genes and promote their adaptation to high temperature environment, which is similar to that of the beetle, Galeruca daurica [32].
KEGG analysis showed that there were more up-regulated expressions of endoplasmic reticulum protein synthesis, phagosome, longevity regulation pathway and cell cycle in commercial population. Among them, endoplasmic reticulum processing is mainly involved in assisting the folding, assembly and transportation of heat shock proteins and other resistant proteins [33]. Therefore, the synthesis and metabolism of resistant proteins such as heat shock proteins may play an important role in the high temperature adaptation of commercial population. However, in natural population, glycolysis and glycosylation, carbon metabolism, biosynthesis of amino acids, pyruvate metabolism, glycosphingolipid biosynthesis-ganglio series, glycosphingolipid biosynthesis-globo series and lysosome were more up-regulated. Glycolysis and glycosylation are pathways for energy generation in insects, which is required for various life activities and survival of natural population under adverse environment [34]. In addition, the more up-regulated genes related to ribosomes and amino acid biosynthesis indicate a more vigorous metabolism [35]. Compared with commercial population, pyruvate metabolism has a higher enrichment degree in natural population. Pyruvate kinase catalyzes the production of pyruvate from phosphoenolpyruvate and the release of ATP in Bombyx mori, thereby enhancing the functional level of glycolysis [36,37]. Therefore, we assumed that natural population may regulate pyruvate metabolism levels to meet the body's energy demand. Glycosphingolipid ensure the stability and strength of the membrane structure, and involved in the signal transduction of cell membrane and activate signal transmission [38]. Lysosomes degrade the denatured proteins produced by high temperature stress, thereby reducing the damage of denatured proteins to cells, which is very important for maintaining cell homeostasis under high temperature environment [39]. Therefore, it is speculated that the adaptation of natural population to temperature stress involves a network system of metabolic processes such as energy production and transport, signal transduction, ion osmosis along with multi-gene expression and complex regulation [40]. It has been shown in honey bees that the ability to adapt and the level of ecological plasticity of populations to the habitat depend on the value of genetic diversity. Natural populations are characterized by increased genetic diversity. Natural populations are characterized by increased genetic diversity compared to commercial population [41]. The adaptation of the natural population of D. helophoroides to temperature stress can be explained by the increased level of genetic diversity that has evolved as a result of long-term adaptive evolution and balancing selection in favor of heterozygotes.
Heat shock proteins (HSPs) are a class of proteins synthesized in vivo in response to environmental stress, which act as molecular chaperones to prevent misfolding and aggregation of other proteins, and are also involved in protein transmembrane transport [42]. In this study, 30 differentially expressed HSPs were extracted from the transcriptome annotation results. Among them, 14 genes were detected in the commercial population of D. helophoroides, and 24 genes were detected in the natural population (Tables A1 and A2 in Appendix A). In both populations, sHSPs and HSP70 were up-regulated, indicating that these two HSPs play an important role in the heat resistance of D. helophoroides, which is consistent with the previous research [21]. In contrast, three new heat shock proteins were present in D. helophoroides under high temperature stress [21], namely, HSP40, HSP60 and HSP90. Among them, the expression of HSP40 was up-regulated in both populations. This result is similar to those found in the endoparasitic wasps, Pteromalus puparum and Cotesia vestalis under heat stress [42,43]. Previous studies had also found that HSP40 generally acts as a molecular chaperone to protect cells, preventing the adhesion of abnormal proteins to each other, repairing denatured proteins, and assisting HSP70 or heat shock homologous protein 70 to play a physiological role together [44]. It is speculated that HSP40 plays an important role in the adaptation to high temperature in both populations of D. helophoroides. HSP60 is involved not only in the folding and assembly of proteins encoded by nuclear genes into mitochondria, but also in the folding and assembly of proteins encoded by mitochondrial genes [45]. However, we only found HSP60 within the natural population of D. helophoroides, which had only down-regulated expression. Up to date, little is known about insect HSP60, which may be related to its poor resistance to adverse temperatures and slow action time [46]. HSP90 can activate protein kinases and steroid hormone receptors by binding them to form binding complexes [47]. This protein was identified in both populations and its expression was upregulated, suggesting that D. helophoroides could regulate the development process by regulating the hormone content and metabolism level in vivo, and reducing the energy consumption of its own development to cope with high temperature stress [48].
Temperature stress induces cells to produce a large amount of reactive oxygen species (ROS), which destroys cell membrane fluidity and causes oxidative damage [49]. Superoxide dismutase (SOD) is the most important antioxidant enzyme in living organisms, which plays an important role in maintaining the balance of ROS in living cells [50]. In eukaryotes, it exists in two forms: Mn-SOD in mitochondrial matrix, Cu-ZnSOD in cytoplasm and extracellular. Mn-SOD plays an important role in scavenging intracellular free radicals [51]. We found that there was no upregulation of the SOD gene in the commercial population of D. helophoroides, while the expression of SOD in both forms was found to be upregulated in the natural population (Table A2), and Mn-SOD was upregulated ten times, indicating that SOD (especially Mn-SOD) played an important role in the adaptability to high temperature in the natural population. Additionally, Peroxiredoxin (Prx) plays an important role in ROS clearance and signal transduction [52,53], but we specifically detected the upregulation of two Prx genes within the natural population (Table A2), indicating that the natural population has much stronger ability to clear ROS than commercial population, or receive signal transduction. Additionally, nicotinamide adenine dinucleotide dehydrogenase (NADH dehydrogenase), acyl-coA dehydrogenase and ATPase are the sources responsible for providing energy for cell activities [40]. In the work, only 4 ATPase genes were found to be up-regulated within the commercial population, while 2 NADH dehydrogenase genes, 1 acyl-coA dehydrogenase gene and 3 ATPase genes were found within the natural population. At the same time, there were a large number of hydrolase genes within the natural population that were not clearly involved with specific metabolic functions. Hydrolase, as group of hydrolytic proteins and metal complexes, is supposed to provide more energy and metal ions (mainly Ca 2+ ) which is helpful to promote signal transduction and energy transport. Therefore, the natural population of D. helophoroides showed higher ROS scavenging mechanism and energy supply capacity than the commercial population.

Conclusions
The transcriptome revealed that these two populations had different ways to adapt the high temperature stress, the commercial population mainly regulated the synthesis of heat shock proteins under the high temperature stress, while the natural population mainly changed the metabolism mode and level to produce a large amount of energy molecules and regulated the neurobiochemical reaction, that produces more resistant substances (such as heat shock protein, SOD, Prx, etc.), so as to protect the cells and other physiological structures in the body. Therefore, the natural populations have more flexible and diverse molecular mechanisms to respond to high temperature stress, and thus have stronger heat resistance than the commercial populations. This result laid a foundation for further revealing the differences in the molecular adaptation mechanisms to high temperature between the two populations, and provided a reference for evaluation of the survival and pest control effectiveness during utilization of commercial population in the forests.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.