In this study, the tissue culture plants of Populus simonii × P. nigra were maintained in the biological clean room of Northeast Forestry University. The 20-day-old sterile seedlings were transferred to greenhouse for soil cultivation. After 30 days, the 72 seedlings with similar growth state were divided into 12 groups and those of 6 groups were selected for insect feeding. We conducted insect attack treatment by placing five third instar larvae of H. cunea on the leaves of poplar plants in mesh bag. After 48 hours, the damaged leaves of the treated plants and the leaves at the same position of the control plants were collected and stored at -80°C, respectively. The plant samples were transported with dry ice to Novogene (China), in which four biological replicates were used for RNA sequencing and proteome analysis, and six biological replicates for metabolome analysis.

The high-throughput output resulted in a large amount of high-quality raw data with scores of Q30 or higher for most bases. The low-quality reads, such as those containing connectors, were completely removed to obtain clean reads. And the resulting clean reads were mapped to Populus trichocarpa genome using HISAT 2 software (https://phytozome-next.jgi.doe.gov/). The mapped reads were then assembled and quantified using String Tie. The FPKM values of the genes were calculated using RSEM.

Quantitative proteome was performed using Tandem Mass Tags (TMT) technology as follows: proteins were extracted according to the previous method (Rao et al., 2021). Protein concentrations were assayed using BCA experimental kit. The peptides were solubilized with 0.5 M TEAB and labelled with the TMT kit (Thermo Fisher Scientific, USA). Liquid chromatography-electrospray ionization-tandem mass spectrometry (LC-ESI-TMS) analysis was then performed as described in the literature (Li et al., 2017). Based on Populus trichocarpa genome database, the raw files obtained from mass spectrometry detection were compared to identify proteins. Peptide, protein and parent ion mass tolerance distribution analyses were performed to assess the quality of the mass spectrometry data. Proteome Discoverer software was used for database searching and quantification of spectral peptide and protein.

Metabolites were analyzed based on high resolution mass spectrometry (HRMS) detection as follows: the samples were processed and detected according to the previous methods (Barri and Dragsted, 2013). The scan range was selected as m/z 100–1500. Offline raw data were analyzed using CD3.3 software. The data were firstly screened with parameters such as retention time and mass-to-charge ratio. The data were then compared with the secondary spectral databases mzCloud and mzVault and the primary database Mass List to identify metabolites. Metabolites with a coefficient of variance (CV) of less than 30% in the QC samples were retained as final identifications for subsequent analysis.

For DEGs, |log2(FC)|>1 and P-value<0.05 were used as thresholds. For DEPs, FC>1.2 or FC<0.833 and P-value<0.05 were used as thresholds. For DEMs, VIP>1.0, FC>1.2 or FC<0.833 and P-value<0.05 were used as thresholds. For metabolites, VIP>1.0, FC>1.2 or FC<0.833 and P-value<0.05 were used as thresholds. The obtained DEGs, DEPs and DEMs were analyzed for GO and KEGG enrichment respectively. Enrichment analysis and visualization were performed using TBtools software (Chen et al., 2023).

qRT-PCR was used to verify the accuracy of RNA-Seq results. Actin was used as the internal control, and three biological replicates were performed for each procedure (Regier and Frey, 2010). All primers used in this study are listed in Supplementary Table S5 . The normalization of the data was performed according to the instructions for the GoTaq 2-Step qRT-PCR System Kit and by the 2^−ΔΔCt method (Livak and Schmittgen, 2001).

Standard errors and standard deviations were calculated by t test. Statistically significant level was set to p-value ≤ 0.05. The data were presented as mean ± standard error (SE) with each SE being calculated from four or six biological repeats..

The differentially expressed genes (DEGs) of Populus simonii × P. nigra under H. cunea stress were analyzed by RNA-Seq. A total of 47.71 Gb high-quality clean data was generated from 8 sequencing libraries, with an average of 5.96 Gb per library. There were 86.18% of the reads mapped to the genome of Populus trichocarpa. Among them, the reads with score of Q30 accounted for 93.91%, and the average GC content was 43.94% ( Supplementary Table S1 ). Principal component analysis (PCA) showed that the eight samples were clustered into two groups: treatment group and control group ( Figure 1A ).

There were a total of 1039 DEGs in poplar under the H. cunea stress of which 992 DEGs were up-regulated and 47 were down-regulated ( Figure 1B ). To validate the accuracy of RNA-Seq, 12 DEGs were randomly selected to detect their relative expression level by qRT-PCR. The results showed that the expression trends of the 12 DEGs were consistent in RNA-Seq analysis and qRT-PCR results, indicating that the transcriptome sequencing data were reliable ( Figure 1D ). Cluster analysis showed that treatment group and control group were clustered into two different categories, and the gene expression pattern was similar within each category ( Supplementary Figure 1 ).

There were 98 transcription factor (TF) genes and 82 functional genes related to pest stress identified in Populus simonii × P. nigra transcriptome ( Supplementary Table S2 ). The functional genes related to adversity stress mainly include protease inhibitors (23), cytochrome P450 superfamily proteins (20), protein kinase superfamily proteins (11), ABC transporter family proteins (7), glutathione S transferases (7) and JA hormone-related genes (14) ( Figure 1C ). The 98 TFs contain ERF, WRKY, MYB, bHLH, C2H2, NAC and other TF family genes. Among them, ERF (25), WRKY (13), MYB (15), and bHLH (8) accounted for the majority ( Figure 1E ).

Based on GO enrichment analysis, the 1039 DEGs were significantly enriched into 146 GO terms, including 93 terms associated with biological processes, 48 molecular functions and 5 cellular components. Significant enrichment terms were involved in catalytic activity, plasma membrane, response to chemical and response to stress. In addition, multiple GO terms were identified to be related to plant stress resistance, including biotic stimulus, abiotic stimulus, transcription regulator activity, lipid metabolic process and enzyme regulator activity ( Figure 1F ). KEGG pathway enrichment analysis showed that many pathways related to plant stress resistance were significantly enriched, including α-linolenic acid metabolism, MAPK hormone signal transduction, plant hormone signal transduction, amino sugar and nucleotide sugar metabolism, phenylpropanoid biosynthesis, ABC transporter, linoleic acid metabolism pathway, glutathione metabolism and plant pathogen interaction ( Figure 1G ).

TMT-labeled quantitative proteomics was used to analyze the differentially expressed proteins of poplar under pest stress. A total of 7094 proteins were detected, and as many as 106 differentially expressed proteins (DEPs) were identified. Among them, 79 DEPs were up-regulated and 27 were down-regulated ( Figure 2A ).

Protein domain analysis showed that protease inhibitors, glycoside hydrolase, cytochrome P450 and nucleoside phosphorylase accounted for 12.62%, 7.77%, 6.80% and 4.85%, respectively, which were the four kinds of proteins with the highest frequency ( Figures 2B, C ). GO enrichment and KEGG functional annotation showed that these DEPs were significantly enriched into 35 GO terms (p ≤ 0.05), including 14 terms associated with biological processes, 20 molecular functions and 1 cellular component. Among the 14 terms associated with biological processes, the terms involved in response to stress, carbohydrate derivative metabolic process, nucleoside metabolic process and single metabolic organism process, which had the highest enrichment abundance. In addition, the terms associated with defense response and biotic stimulus were also significantly enriched. Among the 20 terms associated with molecular functions, the terms of ion binding, metal ion binding, enzyme inhibitor activity, endopeptidase inhibitor activity and transition metal ion binding were the most abundant. The one terms associated with cellular component was enriched in integral components of thylakoid membrane ( Figure 2D ). KEGG enrichment analysis revealed a total of 24 significant enrichment pathways, including α-linolenic acid metabolism, phenylpropanoid biosynthesis pathway, amino sugar and nucleotide sugar metabolism, and MAPK signaling pathway ( Figure 2E ).

The 106 DEPs accounted for 10.2% of the 1039 DEGs. The correlation between the expression level of DEPs and the corresponding gene expression level was analyzed. There were 50 DEPs displaying similar expression trends with their corresponding genes, both of which were significantly up-regulated ( Figure 3B ). The correlation coefficient between the two was 0.5791 ( Figure 3C ).

Among the 50 DEPs, there were 8 trypsin inhibitor family proteins, 7 phosphorylase superfamily proteins, 4 proteins involved in the α-linolenic acid pathway, and 3 cytochrome P450 family proteins identified. In particular, there was an EF-hand protein with a MYB domain ( Supplementary Table S3 ).

GO enrichment results showed that the significantly enriched terms included endopeptidase inhibitor activity, carbohydrate derivative metabolic process, ion binding, stress response, heme binding and single organism metabolism ( Figure 3A ). KEGG enrichment found 10 significantly enriched metabolic pathways, including α-linolenic acid metabolic pathway, phenylpropanoid biosynthesis pathway, cysteine and methionine metabolic pathway and MAKP signal transduction pathway. The α-linolenic acid metabolic pathway was significantly enriched at both transcriptome and proteome level ( Figure 3D ).

Differentially expressed metabolites (DEMs) of Populus simonii × P. nigra under pest stress were identified by liquid chromatography-mass spectrometry (LC-MS). Partial least squares discriminant analysis (PLS-DA) was used to distinguish the difference among the different treatments and different replicates. The variation between different treatments (PC1) and the variation between different replicates (PC2) accounted for 22.07% and 8.54% of the total variation, respectively ( Figure 4A ). With the standards of VIP > 1.0, FC > 1.2 or FC < 0.833 and P-value < 0.05, a total of 212 DEMs were identified, including 128 cationic DEMs, of which 83 DEMs were up-regulated and 45 were down-regulated. There were 84 anion DEMs including 61 up-regulated DEMs and 23 down-regulated DEMs ( Figure 4C ). And the DEMs can be divided into 12 categories, including lipid and lipid molecules, phenylpropanoids and polyketides, organic oxygenates and various free amino acids ( Figure 4B ).

The enriched cationic DEMs can be clustered into 36 metabolic pathways, including nicotinic acid and nicotinamide metabolic pathways, lysine biosynthesis metabolism, phenylpropanoid biosynthesis pathway, quinoline alkaloid biosynthesis, ubiquinone and other terpenoid-quinone biosynthesis pathways, and phenylalanine metabolic pathway ( Figure 4D ). The enriched anion DEMs were classified in 21 metabolic pathways, including purine metabolic pathway, α-linolenic acid metabolism, caffeine metabolism and phenylalanine, tyrosine and tryptophan biosynthesis pathways ( Figure 4E ). Both anion and cation DEMs were related to amino acid biosynthesis and energy metabolism pathways, such as the biosynthesis of cysteine, methionine, arginine, histidine, alanine, aspartic acid, and glutamic acid. At the same time, the DEMs associated with TCA cycle and glucose metabolism pathway also changed significantly.

Association analysis of transcriptome and metabolome showed that a total of 31 metabolic pathways were significantly enriched in DEGs and DEMs, including plant hormone signal transduction, α-linolenic acid metabolic pathway, amino sugar and nucleotide sugar metabolism, phenylpropanoid biosynthesis and other pathways ( Figures 5A, C ). And association analysis of proteomics and metabolomics revealed that there were 12 pathways enriched in DEPs and DEMs simultaneously, including metabolic pathways, biosynthesis of secondary metabolites, α-linolenic acid metabolism, etc. ( Figures 5B, D ).

Multi-omics analysis showed that the DEGs, DEPs and DEMs were significantly enriched in α-linolenic acid metabolic pathway ( Figure 6A ), cysteine and methionine metabolism ( Figure 6B ) and phenylpropanoid biosynthesis and metabolic pathway ( Figure 6C ), etc..

Among them, α-linolenic acid metabolism ( Figure 6A ) is an important way to produce JA, which plays a key role in plant insect resistance. Transcriptome sequencing revealed there were 20 DEGs closely related to α-linolenic acid metabolic pathway, all of which were up-regulated, encoding LOX2S, AOC, OPR, AOS, E2.1.1.141, ACAA1, ACX, DOX1, HPL1, OPCL1 and other 10 enzyme proteins, which were key enzymes in the α-linolenic acid metabolic pathway. The expression trends of 6 DEGs encoding LOX2S, AOC, OPR and AOS were consistent with those of 4 DEPs. In addition, five DEMs including 13(s) -HOTrE, 9(s)-HpOTrE, 12-OPDA, JA and MeJA were closely related to α-linolenic acid metabolic pathway. Among them, 9(s)-HpOTrE synthesis was closely related to the transcription and expression of LOX2S gene. 12-OPDA synthesis was closely related to the transcription and expression of two genes encoding oxidene cyclase (AOC). JA O-methyltransferase (E2.1.1.141) was a key enzyme in the conversion efficiency of JA to methyl jasmonate (MeJA) ( Figure 6A ).

Phenylpropanoid compounds are not only an important part of woods, but also widely involved in the process of insect resistance, disease resistance, UV radiation resistance and anti-aging in woody plants. RNA-seq revealed a total of 27 DEGs encoding 14 enzymes were involved in the biosynthesis and metabolism of phenylalanine, tyrosine and tryptophan, and 7 out of them encoded peroxidase (POD). Moreover, the DEGs encoding amidase (TAT) and prephenate dehydratase (PDT) in the phenylalanine metabolic pathway, the DEGs encoding polyphenol oxidase (PPO) in the tyrosine metabolic pathway, and the DEGs encoding tryptophan synthase (PALP) in the tryptophan synthesis pathway were all up-regulated. In addition, two POD DEPs were up-regulated ( Figure 6C ). Metabolomic analysis showed that the 7 DEMs including tyrosine, tryptophan, cinnamic acid, Caffeic acid, Isoeugenol, Protocatechuate and 3-Hydroxy-benzoate were up-regulated. Among them, the up-regulated expression of cinnamic acid, caffeic acid and isoeugenol in the metabolic group was significantly correlated with the up-regulated expression of the genes involved in phenylalanine ammonia lyase (PAL) and caffeoyl shikimate esterase (CSE) ( Figure 6C ). It shows that the biosynthesis and metabolic pathway of phenylpropanoid compounds is an important way in poplar under H. cunea stress.

Glutathione is a tripeptide containing γ-amide bond and sulfhydryl group, which is composed of glutamic acid, cysteine and glycine. And it has antioxidant effect in plant stress resistance. Transcriptome sequencing revealed that there were 8 DEGs encoding GST, DHAR and gpx in the glutathione metabolic pathway, of which 6 DEGs encoded glutathione S-transferase ( Figure 6B ). Correspondingly, the expression of glutathione sulfide GSSG was significantly up-regulated in the metabolome. In addition, there were 9 DEGs encoding 6 enzymes, including metE, ACHY, metK, E4.4.1.11, TAT and ACS1, which were related to cysteine and methionine metabolic pathways. Among them, three DEGs encoded aminocyclopropane-1-carboxylic acid synthase (ACS1), two of them were up-regulated and one (Potri.014G012600) was down-regulated. The DEGs encoding methionine-γ-lyase (E4.4.1.11) were significantly up-regulated, which showed same expression trend in the proteome, promoting methionine degradation. Metabolomics analysis showed that a total of four DEMs were related to cysteine and methionine metabolic pathways, including Methionine, L-Aspartate, 4-Methylthio-2-Oxobutanoate and O-phospho-L-serine, of which only L-Aspartate was up-regulated. And the down-regulated expression of 4-Methylthio-2-Oxobutanoate was related to the expression change of TAT genes ( Figure 6B ).

Plant tissue damage caused by caterpillars and beetles can activate JA-based signal transduction, which induces JA (Kandoth et al., 2007). JA precursors and JA derivatives play a key role in mediating plant response to pests (Xu et al., 2021). α-Linolenic acid is one of the precursors of JA biosynthesis, which can be converted into 12-oxo phytodienoic acid (12-OPDA) by lipoxygenase (LOX), allene oxide synthase (AOS) and allene oxide cyclase (AOC). And 12-OPDA can be converted into JA under the action of peroxisome-localized OPDA reductase (OPR), which can further derive MeJA (Wasternack and Feussner, 2018). The expression levels of JA biosynthesis genes, LmLOX2 and LmAOS were significantly increased in lupin seeds with MeJA treatment. In addition, a few JA -related defense genes such as LmTPS1, LmTPS4, LmPI2, LmMBL, LmLODC, LmCSD1, and LmPOD were also induced (Erazo-Garcia et al., 2021). In this study, the α-linolenic acid metabolic pathway was significantly enriched in the transcriptome, proteome and metabolome in Populus simonii × P. nigra under H. cunea stress. At transcriptome and proteome levels, the DEGs and DEPs encoding LOX, AOS, AOC, etc. were up-regulated in poplar plants. At the same time, the up-regulated DEMs detected in the metabolome included JA, MeJA and 12-OPDA, as well as JA derivatives such as propyl jasmonate and jasmone. In particular, the DEGs related to the transport of OPDA from chloroplasts to peroxisomes is up-regulated in the transcriptome of poplar under pest stress. It indicated that OPR genes play an important role in the biosynthesis of 12-OPDA, which was key enzymes in the metabolic pathway of α-linolenic acid and JA synthesis.

Studies have shown that MeJA can induce the production of protease inhibitors (PIs). Treatment of C.medinalis and MeJA in rice main stem significantly increased the level of PIs in the leaves of rice at early tillering stage (Tong et al., 2023). MeJA can also induce the accumulation of PI I and PI II in rice leaves (Farmer et al., 1992). In this study, we found that the DEGs encoding PIs was up-regulated by pest stress, which was the highest frequency in all DEGs, accounting for 2.2%. At the same time, the DEPs containing PI domain accounted for 12.62% in all DEPs, which was the highest frequency in all DEPs categories. GO enrichment analysis also indicated PIs were significantly enriched in poplar under pest stress. It shows that pest stress can induce the expression of PIs genes through JA pathway to resist pest feeding in poplar.

Previous study proved that glycoside-bound volatiles and various flavonoids were significantly increased in tea plant under green leafhopper stress, and the related DEGs encoding the above secondary metabolites were significantly up-regulated, while the glutathione content was significantly reduced (Zhao et al., 2020). Glutathione sulfide is an important oxidation product of glutathione after stress in plants (Dixon et al., 2005). In this study, the DEGs encoding glutathione peroxidase (gpx) and deoxyascorbate reductase (DHAR) were significantly up-regulated in poplar under pest stress, and the concentration of glutathione sulfide increased significantly. The contents of amino acids, alkaloids and phenols in the leaves of passion fruit were significantly increased after infection of A. pallida (Yang et al., 2023). Aspartic acid (Asp) is not only a component of protein, but also a downstream metabolite necessary for regulating plant growth and resisting environmental stress (Han et al., 2021). The content of Asp increased significantly when the tomato seedlings were infected by nematodes, indicating that Asp plays an important role in the response to nematode stress (Bali et al., 2020). In this study, metabolomics analysis showed that the content of Asp in poplar increased significantly under the stress of H. cunea. Through the combined analysis of transcriptome and metabolome, it was further found that amino acid metabolic pathways related to Asp were enriched, such as cysteine and methionine metabolism. The results of metabolomics also showed that the content of methionine and O-phospho-L-serine decreased in the metabolic pathway of cysteine and methionine. Among them, O-phospho-L-serine is the precursor of methionine, which participates in threonine pathway. And methionine can be converted into glucosinolates, which can accumulate nearly 30% of the total sulfur content in plant organs, displaying a repellent effect on pests (Hirai et al., 2007). In this study, the DEGs and DEPs encoding Methionine- γ- Lyase E4.4.1.11 were significantly up-regulated in poplar under pest stress. The results indicated that methionine was decomposed under the action of E4.4.1.11, which was conducive to further synthesis of glucosinolates, improving the resistance of plants to pests.

Aromatic amino acids (phenylalanine, tyrosine and tryptophan) in plants are not only essential components of protein synthesis (Maeda and Dudareva, 2012), but also precursors of various secondary metabolites. The secondary metabolites produced by aromatic amino acids are mainly phenylpropanoids, and their biosynthesis is initiated by the activity of phenylalanine ammonia lyase (PAL) (Tzin and Galili, 2010). In general, phenylpropanoid metabolism takes a small number of intermediates in the shikimic acid pathway as the core, producing a large number of secondary metabolites (Vogt, 2010). They are not only representative indicators of plant stress response to light or mineral treatment, but also key medium for pest resistance in plants. The expression of phenylalanine metabolic genes such as PAL and 4-coumarate-CoA ligase (4CL) can be induced by tea green leafhopper, as well as some genes encoding POD and shikimate-O-hydroxycinnamoyl transferase (HCT) in the phenylpropanoid pathway (Camera et al., 2004). In this study, two genes encoding PAL, one gene encoding HCT, three genes encoding 4CL, and seven genes encoding POD were up-regulated. And the proteins encoding POD were also up-regulated, indicating that these two genes play a more important role in pest response in poplar. Phenylpropanoid is traditionally synthesized by phenylalanine, however, the second synthetic method of phenylpropanoid is through tyrosine synthesis in several plants (Beaudoin-Eagan and Thorpe, 1985). In this study, the genes encoding PALP and TAT were significantly up-regulated, which promoted the accumulation of tryptophan and tyrosine. Studies have shown that MYB165 can inhibit the synthesis or accumulation of flavonoids in poplar suspension cells. Overexpression of MYB165 in transgenic poplars reduces the expression of Potri.011G024900, which was a key gene for β-chain of tryptophan synthase (Ma et al., 2018). In our study, one DEG (Potri.011G024900) encoding PALP was up-regulated in poplar under pest stress, and the accumulation of tryptophan was increased, indicating that the gene plays an important role in the synthesis and accumulation of tryptophan.

Lipids are main components of cell membranes, and some lipid molecules play an important role in signal transduction of plants under biotic and abiotic stresses (Xie et al., 2023). In this study, a large of DEMs were identified as lipids, accounting for 29.56% of all DEMs, which include lipid metabolites such as lysophospholipid (LPL) and cucurbitacin. LPL is a lipid-derived metabolic intermediate in the cell membrane, which regulates many important biological processes as an important biologically active signaling molecule (Cao et al., 2023). Previous studies have shown that LPL may be involved in the signal transduction of stress response, and its expression is induced by salt and osmotic stresses (Meijer et al., 2001). The accumulation of a variety of LPL may be more conducive to rapid transmission of stress signals in plants (Arisz and Munnik, 2011). Cucurbitacin (Cuc) is a highly oxidized triterpenoid compound, which is widely found in cucurbitaceae plants and some other plant families, which plays an important role in plant resistance to pathogens and pests. The bitterness produced by Cuc-B in melon is considered to be a defensive response to repel pests and herbivores (Dong et al., 2021). The addition of Cuc-B or Cuc-E to the diet of Drosophila larvae can inhibit its peeling and even cause its death (Toyofuku et al., 2021). Cuc-I, Cuc-B and Cuc-E all affect the feeding behavior of beetles on cucumber and zucchini (Yang et al., 2005). In this study, a variety of LPL, Cuc-B, Cuc-I and Cuc-E were induced to express in poplar under pest stress, which may enhance stress signal transduction and induce the repellent to pests, thereby improving pest resistance in poplar.