Hormetic response and transcriptome proling reveal the mechanism underlying metabolic and cuticular resistance to host-specic terpenoid defences in an emerging insect pest, Pagiophloeus tsushimanus (Coleoptera: Curculionidae)

The resistance mechanisms evolved by insects to overcome host-plant allelochemicals are a key consideration in pest management. Camphor oil (EO) and its main component (i.e., D-camphor) form a specic terpenoid-defensive system in camphor trees, Cinnamomum camphora. However, an emerging insect pest, Pagiophloeus tsushimanus, has recently caused serious damage to this intractable plant species and is largely elusive. Here, we used feeding bioassays and RNA-seq to investigate the mechanism underlying the resistance of the beetle to host-specic terpenoid defences. First, a hormetic response in both larval weight and developmental time, which is a highly generalized dose-response phenomenon in toxicology but occurs infrequently in the context of insect-plant interactions, was observed in terpenoid-feeding individuals. Then, comparative transcriptome analysis between terpenoid-feeding and control groups indicated that both CYP450-mediated metabolic resistance and CP-mediated cuticular resistance were jointly employed to cope with terpenoid-induced stress. In addition, a small portion of genes involved in the glucose transport pathway were upregulated at the low D-camphor dose, suggesting that an extra intake of glucose used for larval growth may contribute to a hormetic response. These ndings suggested that the dual terpenoid resistance mechanisms in this specialist are an essential precondition for the hormetic response in larval growth, ultimately contributing to the widespread successful colonization of host camphor trees. Overall, our study will open new avenues for understanding insect-plant coevolutionary adaptation and developing durable pest control strategies.


Introduction
Plants and herbivores are continuously engaged in a struggle for survival and employ different strategies to prevent the evolution of resistance or adaptation to defensive compounds (Stahl et al., 2018). The preexisting constitutive defences in plants (e.g., toxic secondary metabolites and proteins) are everincreasingly investigated for the development of novel biopesticides (Mossa, 2016;Ntalli et al., 2019;Magalhaes et al., 2020;Maurya et al., 2020). E cient plant resistance exerts strong selection pressure on insects, which in turn develop ways to cope with the deleterious effect of phytotoxin in a stepwise process (Heckel, 2014; Heidel-Fischer and Vogel, 2015; Petschenka and Agrawal, 2016). In recent decades, most reports have documented a suite of factors that confer insect resistance to xenobiotics, including metabolic resistance, cuticular resistance (i.e., penetration resistance), target-site resistance, and behaviour resistance (Tangtrakulwanich and Reddy, 2014). There is, however, a lack of understanding of how these different mechanisms interact (additively, synergistically, or even antagonistically) and to what extent they affect the resistance phenotype.
Evidence is accumulating that cytochrome P450 monooxygenases (CYP450s) are pivotal to metabolic resistance (Hrycay and Bandiera, 2015;Lu et al., 2021;Nauen et al., 2021), and very recently, a comprehensive review discussing the roles of CYP450s in response to phytochemicals was published elsewhere (Vandenhole et al., 2021). Notably, despite considerable progress in the identi cation of shortterm inducible P450s, the relationship between upregulated expression and metabolism remains underexplored, and the extent to which this response is biologically functional should be further determined. Moreover, few studies have examined whether genes induced by initial exposure become constitutively overexpressed after long-term adaptation (Wybouw et al., 2015;Müller et al., 2017;Salehipourshirazi et al., 2020). CYP450-mediated detoxi cation in insect adaptation to xenobiotics has been studied extensively, while other mechanisms exist outside of this paradigm. Cuticular resistance refers to modi cations in the cuticle (i.e., thickening of the cuticle or altered cuticle composition) that ultimately minimize penetration of xenobiotics into insect bodies (Balabanidou et al., 2018). It is commonly believed that structural cuticular proteins (CPs), as basic components of the cuticle, are frequently overexpressed to accomplish cuticular thickening. To date, a small number of CP genes associated with cuticular resistance to insecticides have been identi ed, but a few studies have proven that CPs also confer resistance to host allelochemicals ( . These studies have provided a battery of candidate genes (e.g., the CYP6 family and a suite of CPs) potentially involved in adaptation to plant defensive chemicals. Nevertheless, only a handful of models for insect-plant interactions offered unbiased transcriptional pro les related to host plant adaptations in herbivorous insects, and few studies have revealed these two resistance mechanisms simultaneously in the same experimental model.
The insect pest Pagiophloeus tsushimanus (Morimoto) is a trunk-boring specialist herbivore with only camphor trees as its preferred host in natural environments (Chen et al., 2020;Li et al., 2021). Since this weevil was reported in 2014, it has rapidly reached epidemic levels in Shanghai, China. The camphor tree Cinnamomum camphora (L.) Presl (Laurales: Lauraceae), has been widely exploited as a shade tree or as industrial material in several countries, and it is less likely to be vulnerable to various herbivorous insects because of its speci c terpenoid defences (e.g. a class of monoterpene-derived metabolites such as Dcamphor) (Li et al., 2021). Thus, it would be interesting to explore how this pest could mount an e cient counter-defence response. In this work, we rst evaluated larval performance with different dosages of terpenoid-infused arti cial diets (i.e., treatment with camphor oil and D-camphor in the diet) by measuring a set of developmental parameters (i.e., lethal concentration, larval viability, body weight and developmental time). To determine whether P. tsushimanus exhibits both metabolic and cuticular resistance to host-speci c terpenoid metabolites, the transcriptional responses of P. tsushimanus larvae to different dosages of camphor oil and D-camphor were investigated by using RNA-seq analysis. Our ndings not only provide a new example of the molecular mechanisms underlying insect adaptation to plant secondary metabolites but also establish a foundation for future pest management.

Plant materials, insects and reagents
The camphor tree Cinnamomum camphora (L.) Presl, used as an arti cial diet substrate, was provided by Nanjing Forestry University (32°4′43″N/118°48′30″E), Jiangsu Province, China, and never treated with pesticides. Twigs aged 1-3 years were collected and crushed and then passed through a 100 mesh screen. The resulting powder was stored at 4 ℃ for subsequent experiments.
Pagiophloeus tsushimanus was obtained from a culture established at the Entomology Biological Laboratory, Nanjing Forestry University (EBL, NJFU), Jiangsu Province, China, in 2019 and originally collected from C. camphora plantations in Songjiang District, Shanghai Province, China (30°56′6.15″N, 121°12′32.76″E). The breeding system applied to P. tsushimanus eggs, larvae, pupae and adults was described in our previous study (Li et al., 2021). In brief, eggs were deposited in a Petri dish (5 cm i.d. × 1 cm) with a water-soaked cotton pad for incubation, and neonate larvae were then individually reared in a six-well Petri dish (3.5 cm i.d. × 1.5 cm each well) with arti cial diets. The formulation of the arti cial diets was 50% (w/w) distilled water, 40% (w/w) host-plant powder, 5% (w/w) agar, 3% (w/w) sucrose, 1% (w/w) yeast, 0.5% (w/w) sorbic acid and 0.5% (w/w) sodium benzoate. Both eggs and larvae were kept in an environmental incubator (MTI-201B, Tokyo Rikakikai, Japan) at 28 ± 1°C with a 70 ± 5% relative humidity and a photoperiod of 16:8 h (L:D). Feeding bioassays: larval performance with terpenoidinfused arti cial diets To evaluate their dose response to host terpenoids, neonate larvae were fed arti cial diets mixed with different doses of EO (0, 2, 4, 6, 10, 12, 30, 50, 100, 150, and 200 mg g − 1 (w/w)) and D-camphor (0,5,10,15,17,20,30,50,70, 100, 120, 150, 170, and 200 mg g − 1 (w/w)). Since pure D-camphor is solid at room temperature, it was added directly to the host-plant powder and dissolved through vortexing to prepare a homogenous mixture for formulating the diets. EO was pre-dissolved as a liquid in a hydrosolvent. No fewer than sixty individuals were used for each dose treatment. The mortality at 72 h was recorded, and follow-up observations were performed 24 h later to con rm mortality. Data generated from this trial were used to determine the chronic working doses of EO and D-camphor in accordance with the maximum nonlethal concentration (LC 00 ) and lethal concentration killing 15% of the beetle population (LC 15 ).
According to the results of the dose-response assessment, a total of 800 neonate larvae were continuously reared on arti cial diets containing two doses of EO and D-camphor (LC 00 and LC 15 ), respectively, until pupation (200 larvae per treatment). The control group of larvae (200 larvae) were fed a terpenoid-free diet. These larvae were segregated and reared on fresh diets, which were replaced every three days, and their weight, cumulative mortality and developmental time were simultaneously recorded. The duration of the pupal stage and pupation rates were also determined.

RNA sample collection, library construction, and Illumina sequencing
The treatment procedures were carried out as described above. Several vigorous larvae were randomly sampled for RNA sequencing (RNA-seq) after continuous exposure for 36 days (approximately 4th-instar larvae). Three replicates with four larvae per replicate were used for each treatment (total of 12 larvae per treatment). All RNA-seq samples were processed for total RNA extraction using TRIzol Reagent (TaKaRa, Japan) according to the manufacturer's protocol. RNA quantity was checked using a Nanodrop 2000 (Thermo Scienti c, Waltham, MA, United States). RNA integrity was monitored on 1% agarose gels. The Transcriptome assembly and functional annotation Raw reads were ltered by removing reads containing adapters and > 10% ambiguous 'N' nucleotides and low-quality reads (de ned as reads with > 20% of the bases having quality scores < 10) in SeqPrep software (https://github.com/jstjohn/SeqPrep) before transcriptome assembly. The Q20, Q30 and GC content of the cleaned datasets were calculated to assess sequencing quality in fastx_toolkit_0.0.14 software (http://hannonlab.cshl.edu/fastx_toolkit/). These clean reads of all samples were pooled and then de novo assembled without a reference genome in the TRINITY program with default parameters ( statistical routines for determining differentially expressed genes (DEGs) between the control group (CK) and terpenoid-exposed groups (i.e., EO and D-camphor treatments) (Anders and Huber, 2010). Fold change (FC) values for gene expression were considered signi cant if P-value < 0.05 (false discovery rate ≤ 0.01) and |log2FC| ≥ 1. GO enrichment analysis was performed for the identi ed DEGs of interest using the Goatools program (https://github.com/tanghaibao/GOatools), and an adjusted P-value < 0.05 was the threshold value for signi cant enrichment results.

Weighted gene co-expression network analysis (WGCNA)
Using the WGCNA package in R (https://horvath.genetics.ucla.edu/html/CoexpressionNetwork/Rpackages/WGCNA/), a Pearson's correlation matrix was constructed for all pairs of genes, and a weighted adjacency matrix was constructed with a power function: (where α mn is the adjacency between gene m and gene n, S mn is Pearson's correlation between gene m and gene n, and β is the soft thresholding parameter allowing for a maximal correlation coe cient). The power of β = 4 (model tting index R 2 = 0.89) was set in accordance with the scale-free topology criterion ( Supplementary Fig. S2). The adjacency matrix was further transformed into a topological overlap matrix using the WGCNA package to determine the network connectivity of the genes. The genes with similar expression pro les were classi ed into gene modules by using the cut-tree Hybrid function in Dynamic Tree Cut package version 1.62 with default parameters (minModuleSize = 30 and mergeCutHeight = 0.25). To identify the modules that were signi cantly associated with terpenoid exposure treatments, module signi cance (MS) was set as Spearman's correlation coe cient, and plotted in the WGCNA package. A module with a high MS value was considered to be a key module. Hub genes were de ned based on the degree of network connectivity using the WGCNA package (Langfelder and Horvath, 2008). The top 30 target genes (i.e., those with the highest connectivity) were selected as hub genes associated with terpenoid exposure treatments and visualized using Cytoscape version 3.4 (Shannon et al., 2003).

Real-time quantitative PCR (RT-qPCR) validation
To validate the RNA-seq results, 14 genes of interest were selected for RT-qPCR experiments. Insect samples were prepared as described above for RNA-Seq library preparation. Three biological replicates were included for each treatment, and each biological replicate contained three technical replicates. The cDNA of each sample was synthesized using a HiScript 1st Strand cDNA Synthesis Kit according to the manufacturer's protocol (Vazyme Biotech Co., Ltd, Nanjing, China). RT-qPCR experiments were performed on an Applied Biosystem 7500 Real-Time PCR System (Thermo Fisher, Massachusetts, United States) using SYBR Premix Ex Taq II (TaKaRa, Japan) according to the manufacturer's protocol. The reaction conditions were as follows: 5 min at 95 ℃,followed by 40 cycles of 10 s at 95 ℃ and 40 s at 60 ℃. TATA-binding protein (TBP) was selected as the housekeeping gene based on the geNorm calculation results. Speci c primers for each gene were designed using Primer Premier version 5.0 (Supplementary  Table S5) (Singh et al., 1998) and synthesized at Bgi Gene Technology Co., Ltd. (Beijing, China). Relative expression levels were determined by the 2 −△△Ct method (Schmittgen and Livak, 2008).

Data analyses
The mortality rates of the neonate larvae reared on terpenoid-infused diets after 72 h were regressed against the concentrations using a log-logistic regression model, and the lethal concentrations (i.e., LC 15 and LC 50 ) was calculated using the tted function: See function in the supplementary les section, where y and x are response and dose variables, respectively. This model is generally applied to describe logistic dose responses in pharmacology and chemistry. To combine the effects of terpenoid exposure on developmental time (DT), body weight (BW), and larval viability (LV) in feeding bioassays, we calculated a relative performance index (RPI) for each developmental stage using (Soto et al., 2014). Furthermore, we used the chi-square contingency table test (P < 0.05) to determine the maximum nonlethal concentration (LC 00 ) in the feeding bioassays at low concentrations, and Yates's correction for continuity was performed when expected frequencies were below ve. The developmental durations and larval and pupal weights at different stages were analysed using the Kruskal-Wallis ANOVA test (P < 0.05), and a post hoc test was conducted to evaluate the differences among the treatments. The age-stage speci c survival rates of the larvae reared on terpenoid-infused arti cial diets were compared using Kaplan-Meier survival analysis with the log-rank test. Statistical signi cance of RT-qPCR data was determined using one-way ANOVA followed by Tukey's HSD test (P < 0.05). To verify the reliability of RNA-Seq data in gene expression analysis, Pearson's correlation between fold changes (Log2 transformed) in RT-qPCR and RNA-Seq results was analysed.
All analyses were conducted in SPSS version 20.0 software (IBM SPSS Statistics, Chicago, IL, United States), and plots were generated with OriginPro version 9.0 software (OriginLab Inc., Northampton, United Kingdom).

Effect of EO and D-camphor on larval performance
When neonate larvae were reared on an arti cial diet with different low concentrations of EO and Dcamphor for 3 d, signi cant deaths initially occurred at 12 mg·g − 1 and 20 mg·g − 1 , respectively (P < 0.05, Fig. 1a). Therefore, the LC 00 values of the EO and D-camphor treatments were 10 mg·g − 1 and 17 mg·g − 1 , respectively. As shown in Fig. 1b There were variations in developmental time, body weight, and larval viability between the terpenoidexposed and control groups (Fig. 2). The developmental time of larvae feeding on the terpenoid-infused diets was longer than that of the control group when they reached the 2nd instar. However, a signi cant increase in the growth rate of larvae exposed to terpenoids at LC 00 (EO_10mg_g & D_camphor_17mg_g) was observed afterwards (Fig. 2a). Body weight was either unaffected or slightly inhibited at LC 15 (EO_30mg_g & D_camphor_41mg_g), and we even found that the highest weight gain was observed at LC 00 , with a gain 10% greater than when terpenoid was absent in terms of pupal weight (Fig. 2b).
Terpenoid exposure signi cantly affected the age-speci c survival rate (log-rank test, chi-square = 247.772, P < 0.001). The survivorship curves showed that the larvae reared on terpenoid-infused diets had rapidly declining survivorship from 9 to 24 days after egg hatching, and increasing doses of EO and Dcamphor markedly affected larval viability (Fig. 2c). Our nal step was to combine the developmental parameters described above into a single index (i.e., the RPI) that serves as a proxy of overall performance. According to the heatmap, the populations reared on terpenoid-infused diets at lower dosages had greater z-scores than the control group, whereas the lowest value occurred at higher dosages (Fig. 2d).
Transcriptome responses to EO and D-camphor exposure RNA-seq of P. tsushimanus whole larvae fed terpenoid-infused or terpenoid-free diets yielded a total of 794,599,342 high-quality, clean reads, with a 42.94% − 43.62% GC content and over 96% and 91% Q20 and Q30 values, respectively (Table S1). De novo assembly of transcriptome data resulted in 34,486 unigenes and 44,513 transcripts. The total and average lengths of these unigenes were 41,545,312 bp and 1204.7 bp, respectively, with an N50 (length N for which 50% of all bases in the assembly were located in a unigene of length < N) of 2,245 bp and a mapped fragment percentage of 59.75% (Table S2).
A total of 34,452 unigenes were annotated by comparing them to the COG, GO, KEGG, NR, Swiss-Prot or Pfam databases (Table S3). As shown in Figure S1, NR, GO and KEGG annotation provided comprehensive information on gene function. The species distribution of the top BLAST hits in the NR database showed that 6303 unigenes (34.1%) had the highest scoring blast matches to Dendroctonus ponderosae, followed by Sitopjilus oryzae (3434 unigenes, 18.6%) (Fig. S1a). The observed GO terms for each domain (biological process, cellular component, and molecular function) and KEGG pathways for each category (metabolism, human diseases, organismal systems, cellular processes, environmental information processing, and genetic information processing) are detailed in Figure S1b and c.
To annotate the putative functions of the DEGs identi ed in each treatment, GO enrichment analysis was performed (Fig. 4). In brief, terpenoid exposure induced similar molecular responses in all groups: most upregulated unigenes were signi cantly enriched in monooxygenase activity related to CYP450s and structural constituent of cuticle related to CPs, while downregulated unigenes signi cantly enriched in protein folding and de novo posttranslational protein folding were only observed in the EO_10mg_g treatment (Fig. 4a). Subsequently, we placed emphasis on two GO pathways: monooxygenase activity (GO: 0004497) and structural constituent of cuticle (GO: 0042302). CYP450 genes accounted for 87.6% of the monooxygenase activity (Fig. 5a). Among these genes, several CYP450 genes were signi cantly upregulated after terpenoid exposure, notably the genes called 'TRINITY_DN986_C0_g1', 'TRINITY_DN32261_C0_g1', and 'TRINITY_DN1 2993_C0_g1' (Fig. 5b-e). Analogous to the results for the monooxygenase activity pathway, among the genes associated with structural constituent of cuticle pathway, 88.9% belonged to the CP gene family (Fig. 6a). We also found that most of the CP genes were commonly upregulated in all treatments, while a few genes were found to be downregulated after EO treatments, including the genes called 'TRINITY_DN3873_C0_g1' and 'TRINITY_DN17856_C0_g1' (Fig. 6be).
WGCNA was used to further identify the target genes involved in countering host-speci c terpenoid defence. As shown in the gene dendrogram, the genes with similar expression patterns were globally divided into 13 unique module eigengenes using average linkage clustering ( Supplementary Fig. S3). The relevance between each module and terpenoid-exposure treatment showed that the MEblack module had a higher correlation with D-camphor treatment than other modules (MS-value = 0.441), whereas the MEtan module was the most relevant to EO treatment (MS-value = 0.787) (Fig. 7a). Next, the top 30 genes with the highest connectivity (i.e., hub genes) in the MEblack and MEtan modules were used to analyse gene expression and co-expression networks (Fig. 7b-g). We found that the hub genes in the MEtan module exhibited higher expression levels in EO treatments than in D-camphor treatments (Fig. 7b & c). With respect to hub genes in the MEblack module, there was no obvious difference in expression level between the EO and D-camphor groups (Fig. 7e & f). Using co-expression networks, ve core hub genes in the MEtan and MEblack modules were identi ed ( Fig. 7d & g). Remarkably, the core hub genes in the MEtan module, namely TRINITY_DN986_C0_g1, TRINITY_DN32261_C0_g1, TRINITY_DN12993_C0_g1, and TRINITY_DN1904_ C0_g1, belong to the CYP450 family, whereas the gene TRINITY_DN1820_C0_g1 is of unknow function (Fig. 7d). In the MEblack modules, the core hub genes, namely, TRINITY_DN39845_C0_g1, TRINITY_DN1984_C0_g1, TRINITY_DN624_C0_g1, and TRINITY_DN5193_C0_g1, were annotated to the cuticular protein family, whereas the gene TRINITY_DN129_C0_g1 was uncharacterized (Fig. 7g).

Validation of common DEGs by RT-qPCR
We inspected the expression of fourteen interesting unigenes encoding CYP450s and cuticular proteins across all terpenoid treatments by qRT-PCR (Supplementary Table S4). All the unigenes encoding cytochrome P450s were upregulated to varying degrees, which was highly consistent with the RNA-seq results (Pearson's r = 0.951, P < 0.001) (Fig. 8). Similar trends of upregulation of the unigenes encoding cuticular proteins between RT-qPCR and RNA-seq data were also observed (Pearson's r = 0.604, P < 0.001) (Fig. 9). These results supported the accuracy and reliability of the RNA-Seq data.

A hormetic response in larval growth
Dose responses, known as the backbone of toxicology, are crucial to understanding the capacity of organisms (e.g., insect pests) to tolerate toxic exposures to chemical, biological, and physical agents (e.g., insecticides or phytotoxins) (Agathokleous and Calabrese, 2019). In this context, hormesis, a term that was coined by Southam and Ehrlich (1943), describes biphasic dose-response relationships separated into two dose zones, one (i.e., a relatively low dose) with stimulatory responses, and the other one (i.e., a relatively high dose) with detrimental effects. The results herein showed that hormesis seems to be a notable feature of the biological response of P. tsushimanus to host-speci c terpenoid defences, as demonstrated by the rapid accumulation of body weight within a shorter developmental time when the individuals were exposed to a terpenoid-infused diet at LC 00 (i.e., EO_10mg_g & D_camphor_17mg_g).
This phenomenon was also exempli ed by studies on the response of cotton pests to gossypol, whereby a similar stimulatory effect was triggered by low gossypol concentrations (Stipanovic et (Stebbing, 1998;Calabreseci, 1999), while the latter posits that the stress-exposed individuals shift the balance between energy-con icting physiological trade-offs (e.g., reproduction vs. longevity), favouring one at the expense of another (Forbes, 2000;Jager et al., 2013). Although in this regard we have no data on the relevant parameters of adult reproduction and longevity to further support these hypotheses, the combination of these hypotheses suggests an evolutionary expectation that some emerging specialists that e ciently utilize host-plant chemicals as a feeding stimulus could build up a large population with damaging numbers in a short period, ultimately leading to an outbreak of insect pests.
In addition, the concept of phytotoxin-induced hormesis has not received su cient attention, because its molecular basis either is completely unpredictable, or falls into a general pattern of activating stressresponse mechanisms that are generally bene cial to the organism (e.g., hormones, heat shock proteins and antioxidant enzymes) (Vaiserman, 2010;Yu et al., 2010). In evaluating whether transcriptional responses can offer a potential mechanism for hormesis, the results of this study indicated that low-dose terpenoid treatments generally elicited molecular responses (i.e., CYP450 and CP gene families) similar to those in treatments with terpenoids at high doses. However, one unanticipated nding was the upregulation of the glucose transport pathway after D-camphor exposure at a relatively low dose (Fig. 4c), suggesting a massive intake of glucose (i.e., carbohydrate biotransformed possibly from terpenoid-derived products) applied to larval growth (pictured in Fig. 10). Overall, hormesis is a highly generalized dose-response phenomenon and occurs frequently in various experimental designs (Calabrese, 2005;Calabrese, 2017;Kim et al., 2018), and in the context of insect-plant interactions, it deserves to be further investigated.

CYP450s function in exogenous terpenoid metabolism
Prior studies have noted that the presence of phytochemicals in host plants is considered as a major factor driving the evolution of CYP450-mediated metabolic resistance in both specialist and generalist herbivores (Alyokhin and Chen, 2017). Specialists deploy more e cient CYP450s and superior detoxi cation mechanisms against host-speci c chemicals than generalists (Mao et al., 2006;Zhu et al., 2016). Consistent with the literature, our research provided the rst evidence that several CYP450 genes belonging to the CYP6 and CYP9 subfamilies are jointly responsive to host terpenoids (i.e., EO and Dcamphor) in an emerging specialist herbivore, P. tsushimanus. It can be speculated that a complete metabolic pathway takes place in this beetle. Speci cally, once the terpenoid molecules enter the body cell, they are immediately captured by nuclear receptors to form a complex before migrating into the nucleus. Then the complex increases the transcriptional capacity of downstream genes and achieves the massive synthesis of CYP450 enzymes engaging in the degradation of terpenoid molecules (e.g., oxidation, epoxidation, dehydrogenation, hydrolysis or reduction) and their elimination from the cell (pictured in Fig. 10). Although the transcriptional responses are expected to be directed towards the detoxi cation function of CYP450s, an essential next step in con rming the metabolic fate of terpenoids in P. tsushimanus larvae needs to be carried out. For example, using an arti cial diet containing 14 Clabelled gossypol for Heliothis virescens larvae, Rojas et al. (1992) identi ed that 25% of the labelled gossypol was metabolized by conjugation with six sugar molecules per mole of gossypol and excreted in the frass. Another classical case illustrated that after topically treatment with R-(+)/-S-(-)-α-pinene, the following compounds were detected using GC-MS in the hindgut extracts of the spruce bark beetle, Ips typographus: 4S-(-)-cis-verbenol, 4R-(-)-/4S-(+)-trans-verbenol, R-(+)-/S-(-)-verbenone, and 1R(-)-/1S-(+)myrtenol (Fang et al., 2021). Meanwhile, the metabolic α-pinene products were found to be sensed as aggregation pheromones by conspeci cs of I. typographus. Such sequestration of putative toxic metabolites has been observed in more than 250 herbivorous insect species (Opitz and Müller, 2009). Therefore, exactly as we inferred from the issue of hormesis, in our current study, P. tsushimanus larvae have the capacity to sequester certain terpenoid metabolites as nutrients (e.g., carbohydrates), but further evidence is needed.

CPs constitute a physical barrier against terpenoids
The insect cuticle is composed of the polysaccharide chitin, proteins and lipids that are distributed in distinct horizontal layers, acting mainly as a physical barrier against external compounds (Moussian, 2010;Balabanidou et al., 2018). As the inner chitinous layer, the procuticle that accounts for most of the cuticular mass is thought to be closely associated with penetration resistance (i.e., cuticular resistance). Early works employing genomic analyses revealed that a variety of genes encoding structural CPs are constitutively deposited in the procuticle (Futahashi et al., 2008;Dittmer et al., 2015). Subsequently, there has been growing evidence that cuticular thickening attributed to the overexpression of CPs is one of the mechanisms underlying cuticular resistance (Vontas et al., 2007;Fang et al., 2015;Yahouedo et al., 2017b). In this study, transcriptome analysis in P. tsushimanus indicated that a suite of CP genes were expressed at higher levels in terpenoid treatments. Consequently, in addition to metabolic resistance, CPmediated cuticular resistance may be concomitantly recruited to decrease xenobiotic penetration (pictured in Fig. 10). It, however, remains to be further elucidated whether the elevated expression of CP proteins is consistent with the cuticle-thickening resistance mechanism. In combination with phenotypic analysis (e.g., measurement of cuticle thickness using scanning or transmission electron microscopy), a gene silencing system (e.g., RNA interference) will be a powerful tool for reversely verifying the physiological function of one or groups of CP genes. Using the above methods, Huang et al. (2018b) found that silencing of a single gene, CpCPLCG5, decreased cuticle thickness and increased pyrethroid susceptibility in the resistant strain Culex pipiens pallens.
Overall, cuticular resistance was reported in our present work in combination with metabolic resistance, suggesting an interaction between these two mechanisms. It has been widely accepted that cuticular resistance acts as a modulator that allows metabolic enzymes to work in su cient time and consequently strengthens resistance conferred by metabolic mechanisms. Considering all these ndings, we further predict that the dual terpenoid resistance mechanisms existing in this specialist are an essential precondition for the hormetic response in larval growth, ultimately contributing to the widespread successful colonization of host camphor trees. Concentration-response curves of Pagiophloeus tsushimanus larvae to camphor oil (EO) and D-camphor. a is the maximum nonlethal concentration (LC00) determined in the feeding bioassays at low concentrations. b is the lethal concentrations (i.e., LC15 and LC50) calculated using quadratic logistic regression. Asterisks indicate that the survival rates at different concentrations are signi cantly different according to the chi-square test. *P < 0.05, **P < 0.01, ***P < 0.001, NS, not signi cant.   Red, white, and green sectors indicate that the gene-set expression levels of terpenoid treatments were higher than, similar to, and lower than those of the control group, respectively.

Figure 5
Identi cation of key functional genes based on DEGs and GO enrichment analysis -Cytochrome P450 monooxygenases (CYP450s). a is a Venn diagram depicting the assignment of CYP450 genes to the 'monooxygenase activity' pathway (GO: 0004497). b, c & d are heatmaps and hierarchical clustering graphs of the CYP450 genes that were differentially expressed between the terpenoid treatment and control groups. Red bars indicate upregulated genes, while blue bars indicate downregulated genes.

Figure 6
Identi cation of key functional genes based on DEGs and GO enrichment analysis -Cuticular proteins (CPs). a is a Venn diagram depicting the assignment of CP genes to the 'structural constituent of cuticle' pathway (GO: 0042302). b, c & d are heatmaps and hierarchical clustering graphs of the CP genes that were differentially expressed between the terpenoid treatment and control groups. Red bars indicate upregulated genes, while blue bars indicate downregulated genes.  Relative expression levels of seven DEGs encoding CYP450s in P. tsushimanus treated with EO and Dcamphor. a is comparative analysis of RNA-Seq and RT-qPCR results. b is the correlation between RNA-Seq and RT-qPCR results determined using the Pearson correlation coe cient.

Figure 9
Relative expression levels of seven DEGs encoding CPs in P. tsushimanus treated with EO and Dcamphor. a is comparative analysis of RNA-Seq and RT-qPCR results. Different lowercase letters indicate signi cant differences among treatment groups according to one-way ANOVA followed by Tukey's HSD test at P < 0.05. b is the correlation between RNA-Seq and RT-qPCR results obtained using Pearson correlation coe cient.

Figure 10
Overview of the proposed mechanisms underlying resistance to host-speci c terpenoid defences in P. tsushimanus larvae. Upon exposure to host terpenoids by means of feeding ingestion and penetration, the larvae deploy cuticular resistance [including ( ) activation of the transcriptional capacity of CP genes; ( ) massive production of CPs; ( ) translocation of CPs; ( ) thickening of the procuticle, contributing to an enhanced physical barrier] and metabolic resistance [including ( ) capture of terpenoid molecules; ( ) activation of the transcriptional capacity of CYP450 genes; ( ) massive production of CYP450 enzyme; ( ) high degradation of terpenoid molecules; and ( ) excretion of degradation products]. As these responses occur, certain products may be converted into carbohydrates used for energy, resulting in a stimulatory effect on the larval growth phenotype (i.e., hormesis).