Three TSSM-resistant cassava cultivars (C1115, Miandian and SC9) and three TSSM-susceptible cassava cultivars (SC205, Bread and BAR900) that were identified in our previous study (Liang et al., 2022) were supplied by the National Cassava Germplasm Nursery of China, Chinese Academy of Tropical Agricultural Sciences (CATAS). Cassava stem of about 20 cm length were vertically planted with nutritive soil (equal quantity of soil, peat and perlite) in the pot and grown in a greenhouse for TSSM-resistance identification. The light/dark photoperiod was set as 14 h/10 h, temperature was maintained at 28 ± 1°C, and relative humidity was kept at 75 ± 5%.

TSSM rearing was conducted based on our previously study (Chen Q. et al., 2022). Three mature and healthy cassava leaves with identical growth status from the middle of the 3-month-old plants were selected, and 50 healthy 1-day old female adult mites were infested on abaxial leaves of either resistant or susceptible cassava cultivars, besides, petioles were coated with lanolin to avoid TSSM escape. Leaves without TSSM infestation, 1 day post infestation (dpi) and 4 dpi were sampled (make sure all the mites were eliminated from the sampled leaves). Each treatment repeated three times (Three plants for one treatment, and each plants sampled three leaves).

Total RNA was extracted from 0.1g of leaf sample using a RNAprep Pure Plant plus Kit (Polysaccharides & Polyphenolics-rich, Tiangen, China). RT EasyMix for qPCR (+2-Step gDNA Erase-out, Tolobio, China) was used for first-strand cDNA synthesis. The qPCR reaction system was prepared according to the 2×Q3 SYBR qPCR Master Mix kit (Tolobio, China). The qPCR reactions were performed in 10 μL volume in QuantStudio 6 Flex (Thermo Fisher, America) with three biological replicates as described previous. PAL, 4-coumarate–CoA ligase (4CL), shikimate O-hydroxycinnamoyltransferase (HCT), trans-cinnamate 4-monooxygenase (C4H), caffeoylshikimate esterase (CSE), caffeoyl-CoA O-methyltransferase (CCoAOMT), caffeic acid 3-O-methyltransferase (COMT), ferulate-5-hydroxylase (F5H), cinnamoyl-CoA reductase (CCR) and cinnamyl-alcohol dehydrogenase (CAD), Using the Primer5 to design specific primers according to the conservative domain of lignin biosynthesis gene. The qPCR primers were shown in Supplementary Table S1 , MeActin was used as reference gene. The RT-qPCR conditions were: an initial denaturation for 2 min at 95°C, followed by 40 cycles of denaturation at 95°C for 5 s and annealing at 60°C for 30 s, and a final elongation step at 72°C for 60 s. For the melting curve analysis, a dissociation step cycle (65°C for 5 s, and then an increase of 0.5°C every 10 s up to 95°C) was used. The relative gene transcription was calculated based on the comparative 2^-ΔΔCt method (Livak and Schmittgen, 2001).

The activities of lignin biosynthesis pathway enzymes, such as PAL, 4CL, CCR and HCT, and the lignin content were all analyzed using ELISA kits (Shanghai Enzyme-linked Biotechnology Co., Ltd, Shanghai, China) according to the manufacturer’s instructions (Guo et al., 2022). Each treatment repeated three times.

The whole cassava genome sequence (v.8.0) and annotation were derived from the phytozome database¹. The Hidden Markov Model profile of Aromatic amino acid lyase (PF00221) was retrieved from Pfam² (Finn et al., 2015). The TBtools program (Chen et al., 2020) was used to search for Aromatic amino acid lyase in cassava genome, and then the predicted putative domains were further confirmed by the NCBI Conserved Domain Database³. The MW (molecular weight) and pI (theoretical isoelectric point) of these identified MePAL genes were predicted by ExPASy⁴ (Artimo et al., 2012). Finally, each identified MePAL gene was mapped on the chromosome, and all the MePAL genes were named consecutively (start with MePAL1) according to the orders of the chromosomes.

The chromosomal mapping, gene structures and conserved motif analysis were visualized using the Tbtools. The locations of these MePAL genes were determined by querying the cassava genome. In addition, the conserved motifs among all MePAL genes were identified via online software MEME⁵ (Bailey et al., 2015).

The 2,000 bp sequences from upstream of the transcription start site ATG were extracted as the putative promoter region, and then PlantCare⁶ (Higo et al., 1999) database were used to screen and identify the cis-acting regulatory elements. The components predicted by the PlantCare online tool were first screened, and the unannotated components were eliminated while the annotated component were sorted out, after that, the components with consistent functions shared with the identical function annotations, and element visualization was conducted by using Tbtools.

Phylogenetic analysis of PAL genes from four plant species, i.e., M. esculenta, Ricinus communis, Hevea brasiliensis and A. thaliana, were conducted by MEGA (version 11.0). A phylogenetic tree was constructed with the neighbor-joining method and the parameters were Jones-Taylor-Thornton model, pairwise deletion, and 1,000 bootstrap replicates. Circos program (Krzywinski et al., 2009) was used to illustrate the relationships of duplicated genes. The collinearity analysis of PAL genes between M. esculenta and other three species were detected by One Step MCScanX. KaKs was used to calculate the nonsynonymous replacement rate (Ka) and synonymous replacement rate (Ks), where Ka/Ks1 means positive selection, Ka/Ks1 means purifying selection and Ka/Ks = 1 means neutral evolution.

The TSSM-resistant cassava cultivar SC9 was used as the transfected plant. For vector construction, 300 bp MePAL6 and 300 bp MeCHI (cassava chelatase subunit I gene that accounts for chlorophyll synthesis) DNA fragments were cloned and linked into CsCMV-NC vector as described by Tuo et al. (2021), and SGN-VIGS⁷ (Fernandez-Pozo et al., 2015) online tool was used to select the targeted regions. After silencing of MeCHI, the leaves were supposed to be chlorosis and whitening, which is used as a positive control (Tuo et al., 2021). CsCMV-NC (empty plasmid) was used as a negative control to detect the effect of no-load on cassava plants. All the plasmids (recombinant or empty plasmids) were transformed into Agrobacterium-competent cell GV3101 and cultured in a constant temperature incubator at 28°C for 2-3 days. Once the positive control present obvious chlorosis and whitening symptoms, the leaves of MePAL6-silenced lines were sampled and the silenced efficiency was measured.

After confirming the silencing of the MePAL6 gene, mites were inoculated on the MePAL6-silenced and negative control plants, according to the method described above (Section 2.2). The TSSM infestation symptoms were recorded on 0, 1 and 4 dpi, respectively. In addition, the TSSM-resistance levels of each tested plant were evaluated (conducted on 4 dpi) according to the mite damage index method ( Supplementary Figure S1 ). Furthermore, the expression of other MePAL members, the downstream lignin biosynthesis genes and the lignin content in MePAL6-silenced lines were also detected. The primer information was listed in Supplementary Table S1 .

All data analyses were performed using the SPSS (Version 26.0), and statistical analysis was conducted using one-way analysis of variance (ANOVA) with Tukey’s honestly significant difference (HSD) multiple comparison test. Pearson’s correlation test was used to calculate the correlation coefficients among mite damage indexes, transcription of genes, enzyme activities and lignin contents at 4 dpi. Significant and highly-significant difference were considered if P-values were < 0.05 and < 0.01, respectively.

After TSSM infestation, the transcription of lignin biosynthesis pathway genes ( Figure 1A ) in TSSM-resistant and TSSM-susceptible cassava cultivars presented distinct patterns (The mite damage indexes which were used to interpret the resistance level of each tested cultivars was presented in Figure 1B ). The transcriptions of MePAL, Me4CL, MeCCoAOMT, MeCCR and MeCAD in three TSSM-resistant cultivars (C1115, Miandian and SC9) were significantly increased over time (P < 0.05), while the transcriptions of MeC4H, MeHCT, MeCSE and MeF5H in certain TSSM-resistant cultivars (SC9 and Miandian) were first increased (1 dpi) and then decreased (4 dpi). In addition, the gene transcription in the three TSSM-susceptible cassava cultivars (SC205, Bread and BRA900) showed inconsistent change pattern with MePAL as an exception, in which the transcriptions were continuously increased (1.21- to 2.31-fold relative to those before infestation), in comparison, the transcriptions of MePAL in TSSM-resistant cultivars were much higher after TSSM infestation (3.83- to 6.49-fold relative to those before infestation) ( Figure 1C ).

Enzyme activities were analyzed to examine whether the transcription and post- transcription of the lignin pathway genes showed identical change trend. The four genes, i.e., MePAL, Me4CL, MeCCR and MeHCT that showed stable and significant elevation of transcription were selected. The result illustrated that only PAL demonstrate significantly higher enzyme activity in the TSSM-resistant cultivars compared with the susceptible ones, which was consistent with the transcription levels. However, the rest of three genes did not show consistency between gene transcription and enzyme activity, besides, the activities of 4CL, CCR and HCT in TSSM-resistant cultivars were not necessarily higher than those of susceptible ones ( Figure 1D ).

After TSSM infestation, the elevation of lignin contents in C1115-(R), Miandian-(R) and SC9-(R) posed a hysteretic manner, as compared with those before infestation, the lignin contents remained unchanged on 1 dpi but significantly increased on 4 dpi (P < 0.05). It is noteworthy that the lignin contents in SC205-(S) and BRA900-(S) did not showed statistically difference (P < 0.05), although Bread-(S) possessed constitutively low lignin content before mite infestation, TSSM infestation would significantly trigger the lignin accumulation. In comparison, during TSSM-cassava interaction, especially on long term exposure (4 dpi), the lignin contents in all the resistant cultivars were significantly higher than the susceptible ones ( Figure 1E ).

Correlation analysis was conducted among mite damage indexes, transcription of genes, enzyme activities and lignin contents. It was noticeable that the gene transcriptions and lignin contents were moderately- and highly-correlated with mite damage, while the enzyme activities only showed pretty low correlation with mite damage ( Table 1 ).

A total of 6 MePAL proteins were characterized from cassava and named them from MePAL1 to MePAL6. These 6 MePALs all contained the aromatic amino acid lyase domain based on Pfam analysis. MePAL protein lengths were ranged from 703 (MePAL2) to 790 (MePAL4) amino acids, MW were from 64.47 to 86.14 kDa, and pI were from 5.93 to 6.31. In addition, subcellular localization of MePALs was predicted by WoLF PSORT10⁸. Only MePAL1 were localized on the endoplasmic reticulum, while the remaining MePALs were localized on the chloroplast. More detailed information including MePAL gene annotation, gene accession, chromosome locus, protein length, MW, pI of all identified MePAL were shown in Supplementary Table S2 . To further investigate the chromosomal distribution of the MePAL genes, the DNA sequence of each MePAL was obtained using Blastn in cassava genome database, and these 6 MePAL genes were mapped on 6 chromosomes ( Figure 2A ).

The amino acid sequences of PALs deriving from three Euphorbiaceae plant species, i.e., M. esculenta, R. communis, H. brasiliensis together with a model plant A. thaliana were used to construct a phylogenetic tree by using neighbor-joining method. All the PALs from different plant species were clustered into 4 groups. Group 1 contained the largest PAL gene members (8 PALs), followed by Group 2 and Group 4 with 5 PALs. The 6 MePALs were evenly distributed in each group which contained at least one PAL gene ( Figure 2B ), In addition, it is noteworthy that MePAL6 were clustered with all the 4 AtPAL genes in Group 1, indicating that the evolution conservation between MePAL6 and AtPALs, moreover, as the PAL gene function in all the three Euphorbiaceae plants were not identified so far, while the function of AtPALs in dealing with environmental stress had been well-studied (Wanner et al., 1995), based on the above information, we assumed that specific gene in cassava (i.e., MePAL6) which was clustered together with the AtPALs maybe responsible for abiotic or biotic stresses.

The exon/intron organization and conserved motifs of all MePAL genes were analyzed. Eight conserved motifs with multiple repeats were identified among the 6 MePAL proteins, besides, the motifs number was ranged from 9 (MePAL4) to 14 (MePAL1, MePAL5), although MePAL4 contains only 9 motifs, it possessed the longest amino acid sequence among all MePALs ( Figure 3A ), in addition, sequence alignment indicated that the similarity of the 6 MePAL proteins was 77.01% ( Supplementary Figure S2 ). All MePAL genes had 2 exons and 2 introns, except for MePAL6, which had only one intron ( Figure 3B ).

PlantCARE database was used to characterize the cis-acting regulatory elements (CAREs) of each MePAL gene. Light responsiveness and MYB elements were presented in all 6 MePALs. Furthermore, CAREs related to hormone responses were also identified, such as Methyl jasmonate (MeJA)-responsiveness, abscisic acid responsiveness and gibberellin-responsive element. Interestingly, the CAREs of some MePALs also contained elements related to biotic and abiotic responses, including wound-responsive element, anaerobic induction, MYB binding site involved drought-inducibility, defense and stress responsiveness. In summary, diverse hormone and environmental factors might affect the expression of MePAL genes ( Figure 3C and Supplementary Table S3 ).

To reveal the expansion mechanism of the MePAL genes, CDS of all MePAL genes were subjected to Blastn within the cassava genome. Totally, 8 pairs (5 MePAL genes) of segmental duplications as well as 8 pairs of fragment duplications (MePAL1/MePAL2, MePAL1/MePAL3, MePAL1/MePAL4, MePAL1/MePAL5, MePAL2/MePAL3, MePAL2/MePAL4, MePAL2/MePAL5 and MePAL3/MePAL4) were identified ( Figure 3D ). In addition, the Ka/Ks values of all gene pairs were less than 1.0, which indicated that these genes evolved under purification selection ( Supplementary Table S4 ). In addition, to detect the synteny of PAL genes, a collinearity analysis between M. esculenta and other plant species using TBtools were performed. Finally, 7 paired collinearity relationships between 5 MePAL and 4 AtPAL genes were established, 15 paired collinearity relationships between 6 MePAL and 6 HbPAL genes, and 16 paired collinearity relationships between 6 MePAL and 4 RcPAL genes ( Figure 3E ).

In order to delicately speculate the transcription patterns of the 6 identified MePAL genes, primers were designed mapping to the specific regions of each gene ( Supplementary Table S1 ). After TSSM infestation., the transcriptions of most MePAL genes (MePAL1, MePAL2, MePAL4 and MePAL6) significantly increased over time, while the transcription of MePAL3 significantly decreased, and the transcription of MePAL5 did not show statistically difference (P < 0.05) ( Figure 4 ). In particular, the transcription of MePAL6 increased by approximately 6-fold (at 4 dpi), which was the most strongly induced gene ( Figure 4 ). Take into account the vigorous induction as well as the function similarity with model plant (Section 3.2), the MePAL6 gene, was finally selected to clarify the potential function associated with TSSM-resistance in cassava.

To study the roles of MePAL genes in lignin biosynthesis and TSSM-resistance, a MePAL6-silenced line CsCMV-MePAL6 (treatment) was constructed by using VIGS. The chlorosis and whitening phenotype of CsCMV-MeCHI line (positive control) indicated the effectiveness of viral inoculation on silenced plant ( Supplementary Figure S3 ). In addition, the transcription of MePAL6 in CsCMV-MePAL6 were drastically reduced compared with negative control CsCMV-NC ( Figure 5A ) (Silenced efficiency was approximately 78.35%). Furthermore, as the VIGS target region of MePAL6 showed high identity to the rest of 5 MePAL genes, thus, their transcriptions were also investigated. qPCR analysis showed that MePAL2 and MePAL5 were also suffered transcription depression, while the transcription of MePAL1, MePAL3 and MePAL4 were not suppressed but were increased at certain timepoints after TSSM infestation ( Figure 5A ).

Influence on the expression of downstream genes involved in lignin biosynthesis pathway were also examined. The transcription of Me4CL and MeCCR in the CsCMV-MePAL6 line significantly increased on 4 dpi, nevertheless, it was still significantly lower than that of negative control CsCMV-NC. In addition, MeHCT was also induced after TSSM-infestation, but there was no significant difference between CsCMV-MePAL6 line and CsCMV-NC ( Figure 5B ). Moreover, enzyme activities were analyzed to examine whether the transcription and post-transcription of the downstream genes showed identical change trend. Results indicated that compared with control, the PAL activity in MePAL6-silenced line was significantly decreased (P < 0.05), nevertheless, it still can be induced under TSSM-infestation ( Figure 5C ). By contrast, the activities of 4CL, CCR and HCT between control and MePAL6-silenced line basically did not show significant differences (P < 0.05) ( Figure 5C ).

Distinct TSSM infestation symptoms were observed between MePAL6-silenced line and control. During short-term exposure to TSSM, there were only very slight white and yellow spots on both two lines without difference in symptom. However, on 4 dpi, the MePAL6-silenced line suffered serious TSSM infestation symptoms, where dense TSSM damage spots covering the whole leaf, by comparison, there was no symptom deterioration in the control plant ( Figure 6A ). Moreover, the mite damage index of MePAL6-silenced line was 68.51% (the TSSM-resistant level was identified as ‘susceptible’), which was higher than that in CsCMV-NC (21.35%, the TSSM-resistant level was identified as ‘resistance’), indicating the TSSM-resistance level of SC9 was shift from resistance to susceptible, after silencing of lignin biosynthesis gene MePAL6 ( Figure 6B ). In addition, before TSSM infestation, the lignin content in MePAL6-silenced line was significantly lower than that in control, suggesting silencing of MePAL gene would decrease the lignin production. Moreover, after TSSM infestation, the lignin content in both two cassava lines significantly increased, but the CsCMV-NC presented higher lignin accumulation, especially when exposed to TSSM for a long term (4 dpi) ( Figure 6C ). We speculated that the decrease of lignin content might be the reason of TSSM-resistance deterioration.

Lignin and its intermediate products are crucial defensive substances coping with pest infestation. By serving as physical and chemical barrier, effective accumulation of lignin can significantly enhance plant resistance (Lee et al., 2019). A battery of genes participates in lignin biosynthesis ( Figure 1A ), and those genes could be significantly induced during insect pest invasion. The expression levels of the PAL, 4CL, COMT and CAD genes in the lignin synthesis pathway increased significantly after Panax notoginseng was inoculate with fungal (Yang et al., 2022). Transcriptomic analysis depicted that when Chinese chestnut (Castanea mollissima) in response to the chestnut gall wasp (Dryocosmus kuriphilus) infestation, the majority of genes associated with the lignin biosynthesis pathway, including PAL, CAD, CCOAOMT, COMT and HCT were significantly upregulated (Zhu et al., 2019), besides, identical genes with similar upregulation was found in the study regarding soybean-bean pyralid larvae interaction (Zeng et al., 2017). Similar to former studies, here we noticed that almost all the lignin biosynthesis pathway genes were significantly induced while cassava plants encountered TSSM-feeding stress, probably imply the universality of lignin-based defense response in different plants.

The pest-resistant and pest-susceptible plants usually display distinct expression pattern of lignin synthesis genes. For instances, the expression level of PAL increased rapidly after insect-resistant rice cultivars was infested by small brown planthopper (Nilaparvata lugens), and was significantly higher than the susceptible varieties (Duan et al., 2014), in another study conducted by Panatda (Jannoey et al., 2017), higher expression of C4H in resistant rice variety than the susceptible one also can be observed. In addition, when infested by root-knot nematode, both resistant and susceptible tomato cultivars showed expression elevation of lignin biosynthesis genes (PAL, C4H, HCT and F5H) at early times (2–4 dpi), and the induction was faster and greater in resistant cultivars after infection (Veronico et al., 2018). Similarily, in present study, stronger gene induction on TSSM-resistant cassava cultivars was also observed, which probably indicated the significance of lignin pathway contributing to cassava resistance to TSSM. Furthermore, as the MePAL gene demonstrated the most robust induction (in both transcription and enzymatic level) among all the lignin biosynthesis genes, we therefore focused on revealing the potential function of MePAL in conferring cassava resistance to TSSM.

PAL is encoded by a small gene family, and different PAL genes usually have multiple functions (Shine et al., 2016). For examples, there were 4 PAL genes identified in A. thaliana, and three of them (AtPAL1, AtPAL2, and AtPAL4) were highly expressed in lignifying cells-rich inflorescent stems, whereas AtPAL3 was expressed at a very low level, indicating only the former three genes involved in lignin biosynthesis (Raes et al., 2003). There are 14 PALs identified in potato (Solanum tuberosum), among them StPAL1, StPAL6, StPAL8, StPAL12, and StPAL13 functioned in the stress defense against high temperature and drought, while StPAL1, StPAL2, and StPAL6 participated in chemical defense mechanisms (Mo et al., 2022). In pepper (Capsicum annuum) only CaPAL1 was found to be responsible for defense against microbial pathogens (Kim and Hwang, 2014). Although the genome of cassava had been sequenced, the gene family identification of PAL has not been carried out yet, the function of PAL in cassava is largely unknown.

In this study, six MePAL genes were identified in cassava genome, and these six MePAL genes were distributed on six chromosomes, which were located at the top or bottom of chromosomes. The PAL gene family is represented as a mini-family in most plants, for examples, there are four genes in A.thaliana (Cochrane et al., 2004), five genes in Populus trichocarpa (Shi et al., 2013), nine genes in Oryza sativa (Gho et al., 2020), fifteen genes in Vitis vinifera (Zhao et al., 2021) and seven genes in Camellia sinensis (Chen X et al., 2022). The relative MW and pI values of MePAL family genes were similar, indicating that the evolutionary conservation of MePAL family genes in cassava. Phylogenetic analysis showed that the MePAL6 gene was clustered with all the 4 AtPALs, as different AtPAL genes were validated to cope with abiotic or biotic stress such as pathogen infection, trauma, nutrient depletion, ultraviolet radiation, and extreme temperatures (Huang et al., 2010; Joshi et al., 2022), we assumed that MePAL6 may also possessed similar function in dealing with environmental stress.

Cis-acting elements analysis of the critical genes will help in the elucidating the molecular mechanisms that associated with plant stress responses (Moghadam et al., 2012; Moghadam et al., 2013). This study indicated that considerable hormone response elements (i.e., MeJA-, salicylic acid-, gibberellin- and abscisic acid responsivenesses) as well as numerous stress response elements (i.e., anaerobic induction, MYB binding site involved drought-inducibility, wound-responsive element and low-temperature responsiveness) probably accounted for the induction of MePALs when cassava plant suffered TSSM infestation. Valifard et al. (2014) conducted promoter analysis of PAL genes from three plants (Mimulus guttatus, Zea mays and A. thaliana), and several common elements were predicted to be responsible for abscisic acid, alicylic acid, anaerobic induction, heat stress responses, light responses, MeJA-responses and wound responses. In agreement with these findings, studies have also shown the existence of TC-rich repeats (Zhao et al., 2012);, CAAT box, G box, CGTCA motif, TCA-element (Jiang et al., 2013) in PAL promoter in accordance with stress hormones may dealing with environmental stress. From present and previous studies, we presumably speculated that the cis-acting elements in different plant species shared similar biological function, and the induction of PAL gene could be considered as a universal stress response while plant was under abiotic or biotic stress.

Manipulation the gene expression involved in lignin biosynthesis will significantly alter the lignin accumulation. To date, most genes that participate in lignin biosynthesis had been subjected to establish genetic modified plants, of which the lignin content and disposition were further analyzed. For examples, overexpression of CCR, F5H, CSE, CCoAOMT will increase the lignin content in Paspalum dilatatum (Giordano et al., 2014), O. sativa (Takeda et al., 2017), while silenced the 4CL, CAD and HCT will decreased the lignin content in O. sativa (Liu, H. et al., 2017), Populus trichocarpa (Van Acker et al., 2017), Populus nigra (Vanholme et al., 2013). In addition, modification the expression of CCR, HCT and CAD can result in reconfiguration of different lignin units, i.e., S-unit, G-units and H-unit. Besides, some antioxidant enzymes, like POD (Herrero et al., 2013) SOD and APX (Shafi et al., 2015), and certain transcription factors such as MYB (Guo et al., 2017), bHLH (Gao et al., 2019) and WRKY (Guillaumie et al., 2009) also can regulate the lignin biosynthesis. Among those genes, PAL attracts great concern. In one hand, molecular genetics methods have been used to silence or disrupt PAL genes in numerous plants, like A. thaliana (Huang et al., 2010), tobacco (Korth et al., 2001), Scutellaria baicalensis (Park et al., 2012) and rice (Fang et al., 2013), not surprisingly, lignin modified plants were successfully constructed. In the other hand, delicate investigations were also performed to make out the specific biological function of certain PAL gene member (Kim and Hwang, 2014; Van Acker et al., 2017; Mo et al., 2022)

By silencing one of the PAL gene members, the MePAL6, we got a cassava line with 48% reduction of lignin content. Interestingly, the expressions of other PAL gene members (MePAL2 and MePAL5) were also suppressed to some extent. Similarily, some studies also speculated that manipulation the expression of single gene member would also affect the rest of family gene members. There are seven Rac genes in rice (O. sativa), by using highly conserved regions of the two members (OsRac1 and OsRac5) of the whole OsRac gene family, Miki et al. (2005) established transgenic lines by specifically RNA silencing of these two genes, meanwhile, suppression of all members of the gene family with variable efficiencies were also observed. Similarily, silencing GhOPR9(12-oxo-phytodienoic acid reductases) gene in cotton (Gossypium hirsutum) also cause the suppression of GhOPR3 expression (Liu et al., 2020). In plants, RNAi- and VIGS-based technologies successfully silenced the specific gene members without influencing the transcriptions of the rest of closely related family members or simultaneously silenced a few family members to overcome functional redundancy (Burch-Smith et al., 2004; Hwang and Gelvin, 2004). However, in certain rigid situation, 100 percent of none off-target was inevitable. In present study, the VIGS primers were specifically designed that only target MePAL6 ( Supplementary Figure S4 ), we assumed that off-target effect can be minimize, in fact, the target gene MePAL6 showed sufficient suppression here. In the other hand, silencing MePAL6 may also change the expression of some downstream synthesis genes, i.e, CCR and 4CL. This cascade effect of gene expression in specific metabolites biosynthesis pathway is common in several studies. In the GhENODL6 (Early nodulin-like protein 6) silenced cotton, the transcriptions of both PAL and 4CL genes significantly decreased, resulting in the reduction of SA content (Zhang et al., 2022). In addition, in the GhOPR9 silenced cotton, expression alterations of other jasmonic acid pathway related genes, like lipoxygenase (LOX), allene oxide cyclase (AOC) and allene oxide synthase (AOS), were also measured. In present study, silencing of MePAL6 cause systematic impact on the gene expression of lignin biosynthesis pathway, which resulted in decrease of lignin content. The gene expression and metabolite changes lead to the deterioration of TSSM-resistance, as the cassava cultivar SC9 would shift from TSSM-resistant plant into a TSSM-susceptible one. In addition, PAL genes not only participate in lignin biosynthesis, but also involved in biosynthesis of several important plant defense activators, i.e., salicylic acid, jasmonic acid and ethylene in several plant species (Kim and Hwang, 2014; Jiang and Yan, 2018; Fujimoto et al., 2015). Manipulating the expression of PAL gene might simultaneously affect other defensive responses. Thus, we assumed that apart from regulating lignin accumulation, MePAL6 would probably attribute to several other secondary metabolite-related defensive pathways, which collectively shaped cassava resistance to TSSM.

However, more delicate efforts are still needed to fill the knowledge gap of MePAL function in cassava. The PAL gene is associates with abiotic and biotic stress such as pathogen infection, mechanical damage, ultraviolet radiation, chemical treatment and extreme temperatures. This study is the first attempt to know about MePAL6 involved in insect pest resistance. However, gene functions of other MePALs have not been sufficiently illustrated. To fully demonstrated how every single MePAL gene work independently or collaboratively in specific biological process, will largely benefit the molecular design breeding of novel cassava varieties that are adaptive to abiotic and biotic stress, which can effectively reduce the adverse environmental impacts and resulting cost, ultimately, promoting the cassava production.