Transcriptomic analyses provide insight into adventitious root formation of Euryodendron excelsum H. T. Chang during ex vitro rooting

Euryodendron excelsum H. T. Chang, a critically endangered species endemic to China, is a source of valuable material for the furniture and construction industries. However, this species has some challenges associated with rooting during in vitro propagation that have yet to be resolved. In this study, we optimized rooting and conducted a transcriptomic analysis to appreciate its molecular mechanism, thereby promoting the practical application of in vitro propagation of E. excelsum, and providing technical support for the ecological protection of this rare and endangered species. Results showed that ex vitro rooting performed the highest rooting percentage with 98.33% at 25 days. During ex vitro rooting, there was a wide fluctuation of endogenous levels of indole-3-acetic acid (IAA) and hydrogen peroxide (H2O2) at the stage of root primordia formation. Transcriptome analysis revealed multiple differentially expressed genes (DEGs) involved in adventitious root (AR) development. DEGs involved in plant hormone signal transduction, such as genes encoding auxin-induced protein, auxin-responsive protein, and auxin transporter-like protein, and in response to H2O2, oxidative stress, abiotic and biotic stimuli were significantly up- or down-regulated by ex vitro treatment with 1 mM indole-3-butyric acid (IBA). Our results indicate that ex vitro rooting is an effective method to induce AR from E. excelsum plantlets during micropropagation. DEGs involved in the plant hormone signal transduction pathway played a crucial role in AR formation. H2O2, produced by environmental stimulation, might be related to AR induction as a result of the synergistic action with IBA, ultimately regulating the level of endogenous IAA. Under ex vitro rooting, a synergistic action between H2O2 produced by environmental stimulation and IBA played crucial role in the regulation of AR formation from E. excelsum plantlets during micropropagation.


Introduction
Adventitious root (AR) system is one of plant roots systems that arises from parts of the plant rather than from the roots of the embryo (Barlow 1986). ARs derived from non-root tissues, are the main path by which new plantlets root in vegetative propagation, and usually generated during normal development or stress conditions (Steffens and Rasmussen 2016). In vitro propagation via tissue culture has become an important technology in plant conservation strategies given its advantages, such as high propagation coefficient, Communicated by Paloma Moncaleán.
Yuping Xiong and Shuangyan Chen have contributed equally to this work. freedom from restrictions imposed by season, especially for rare and endangered species (Bhardwaj et al. 2018;Khater and Benbouza 2019;Rameshkumar et al. 2017). However, in some in vitro propagation systems, plants may display rooting-recalcitrance problems, for example by Juniperus thurifera L. (Khater and Benbouza 2019), Zeyheria montana Mart. (Cardoso and Teixeira da Silva 2013), Elegia capensis (Burm. f.) Schelpe (Verstraeten and Geelen 2015), and Cariniana legalis (Lerin et al. 2021). Rooting-related problems limit the application of in vitro propagation for plant breeding and conservation efforts. Therefore, AR formation during plant in vitro culture is a top research objective for plant asexual propagation breeders.
Plant growth regulators (PGRs) are commonly AR inducers used in in vitro culture, such as 1-naphthaleneacetic acid (NAA), indole-3-acetic acid (IAA) and indole-3-butyric acid (IBA), but these tend to show species-and concentrationdependent AR induction efficiency. For Laburnum anagyroides Medic., a low concentration of IAA, IBA or NAA induced AR normally, but high concentrations induced callus formation in shoot tips and subsequently, plant death (Timofeeva et al. 2014). 0.25 or 0.5 mg/L of NAA promoted rooting during in vitro culture of Cornus alba L., while IBA had an adverse effect on root growth and even inhibited AR induction at 1.0 mg/L (Ilczuk and Jacygrad 2016). Moreover, in vitro rooting and ex vitro rooting also display some differences in AR formation. In vitro rooting to induce AR is always performed under aseptic conditions (Barpete et al. 2014;Guo et al. 2019;Nourissier and Monteuuis 2008). In contrast, ex vitro rooting employs unrooted plantlets that are removed from aseptic in vitro conditions culture to induce AR in an open environment (Revathi et al. 2018;Shekhawat and Manokari 2016). Although the two culture methods are affected by various factors, ex vitro rooting can enhance rooting percentage and survival during plant acclimatization, and reducing limiting factors in micropropagation (Benmahioul et al. 2012;Yan et al. 2010). For example, Ceratonia siliqua L. plantlets treated with 4.8 μM IBA displayed a 46.3% rooting response, forming a fragile root system when rooted in vitro whereas the induction of AR from ex vitro shoots treated with 14.4 μM IBA showed significantly higher rooting percentage (91.7%), and a normal morphological appearance, and were successfully acclimatized, showing more than 90% survival (Lozzi et al. 2019).
The mechanism of AR induction involves various key genes, proteins, and pathways (Chen et al. 2020a;Qi et al. 2020;Stevens et al. 2018). Most PGRs promote AR development by regulating the level of endogenous IAA, thus genes and pathways related to the biosynthesis and transport of IAA are considered to play a significant role in AR formation . Transcriptome sequencing revealed that candidate genes involved in AR formation of Mangifera indica L. cv. Zihua cotyledon segments were predicted to encode polar auxin transport carriers, auxin-regulated proteins and cell wall remodeling enzymes (Li et al. 2017). In Arabidopsis thaliana, IBA induced AR formation in thin cell layers by conversion into IAA involving nitric oxide activity, and by positive action on IAA transport and biosynthesis (Fattorini et al. 2017). Genes related to the synthesis, transport, metabolism and recognition of plant hormone were involved in the in vitro induction and elongation of ARs in Populus euramericana (Zhang et al. 2019b). However, knowledge of the molecular aspects of adventitious rooting in plants, especially in woody species that are recalcitrant to rooting, remains scanty. Understanding the mechanism of AR formation is of great importance to strategize plant breeding and conservation efforts to maximize the marketable yield and research value, especially of rare and endangered plants.
Euryodendron excelsum H. T. Chang, a monotypic genus endemic to China, is fine-textured and colorful, making it a source of valuable material for the furniture and construction industries (Chang 1963). However, mainly as a result of habitat destruction and deforestation caused by human activity, only a single population of E. excelsum can now be found at Bajia Zhen, Yangchun County, Guangdong Province, in southern China (Shen et al. 2008;Ye et al. 2002). E. excelsum is naturally propagated only by seeds, but seed germination and seedling growth toward adulthood are fragile stages that limit natural recruitment and regeneration (Shen et al. 2009. Based on the categorization of the International Union for Conservation of Nature and Natural Resources (IUCN), E. excelsum has been listed as a critically endangered plant since 1998, and continues to maintain this status and faces a high risk of extinction, implying that some strategies were put forward for the conservation of E. excelsum populations (Barstow 2020).
In our previous study, a micropropagation system for E. excelsum was established by in vitro culture. When treated with either IBA or NAA, in vitro E. excelsum showed lower rooting percentage in agarized woody plant medium (WPM) (Lloyd and McCown 1980) than in agar-free vermiculitebased WPM after culture for 2 months, and callus formed at the base of stems in these media, hampering the successful transplantation of plantlets (Chen et al. 2020b). We inferred from that there may be other factors that can stimulate AR formation in E. excelsum when cultured in vitro, or that enhance the induction efficiency of PGRs. Thus, the objectives of this study were to improve the micropropagation system of E. excelsum by optimizing the AR induction conditions, and to reveal the key influencing factors and related genes and processes underlying AR formation by transcriptomic analysis. By better elucidating the mechanism of AR formation of E. excelsum in vitro, research on biological conservation and genetic engineering of E. excelsum can be promoted and advanced.

Culture of plantlets
The basic micropropagation system for E. excelsum that was previously established (Chen et al. 2020b), was employed in this study. In vitro plantlets were maintained and propagated on WPM supplemented with 4.44 μM BA (Solarbio, Beijing, China) and 0.53 μM NAA (Macklin, Shanghai, China). Single shoots with more than four leaves and two nodes cut from multiple shoots were inoculated on PGR-free WPM for 30 days. During this period, AR was not induced.

Adventitious root induction
For the in vitro rooting treatment, 2 mm of the base of single shoots cultured on PGR-free WPM for 30 days were cut and trimmed shoots were inoculated on WPM supplemented with 0, 0.05, 0.5, and 5 μM NAA, IBA or IAA (Macklin). Shoots in the control group were inoculated on PGR/auxin-free WPM. Ten shoots were placed in each jar, and four jars were prepared for each treatment. Three replicates were performed for each treatment (n = 12 jars; 120 shoots in total).
For the ex vitro rooting treatment, single shoots cultured on PGR-free WPM for 30 days were removed from culture jars, and about 2 mm was trimmed from the base. Trimmed shoots were treated with 0, 1, 2, and 3 mM NAA, IBA or IAA for 10 min, then transferred to plates (5 cm in height; 27 cm in width; 47 cm in length) for raising seedlings supplemented with vermiculite and perlite (v/v, 1:1). Trimmed shoots cultured on PGR/auxin-free WPM served as the control. Forty shoots were planted in each plate, and three replicate plates were prepared for each treatment (n = 3 plates; 120 shoots in total).
Rooting percentage as well as average root number and root length were calculated for each treatment. After one-way analysis of variance (ANOVA), treatment means were compared by Duncan's multiple range test (DMRT) in SPSS Statistics version 20.0 (IBM, New York, USA) and were considered to be significant difference between the designated treatments at P < 0.05.

Histological analysis
The base of shoots (0.5-1.0 cm) was collected at 0, 2, 4, 6, 8, 10 and 12 days after the optimum treatment method under ex vitro rooting and fixed for 24 h in formalin/acetic acid/alcohol at 25 ± 2 °C. At least 15 bases were collected for each time point. Fixed material was dehydrated in a 70-100% alcohol dehydration series followed by infiltration with molten paraffin (Mackin), and embedded in paraffin wax. Sections (8-10 μm thick) were made with a rotary microtome (KEDEE, Zhejiang, China) and stained in 0.02-0.05% toluidine blue (Mackin). Sections were viewed with a Nikon Eclipse E200 microscope (Nikon, Tokyo, Japan) and micrographs were captured using a HQimage C630 digital camera (Hengqiao, Hangzhou, China).

Determination of endogenous IAA and hydrogen peroxide (H 2 O 2 ) content
To analyze IAA and H 2 O 2 content, the same method and growth conditions were employed as for histological analysis. Material was stored at − 80 ℃. Three biological replicates of 10 cut stem bases were harvested as 0.1 g fresh weight (FW) to assess endogenous IAA and H 2 O 2 content, according to the instructions of an IAA Enzyme Linked Immunosorbent Assay kit (Dogesce, Beijing, China)  and Hydrogen Peroxide Assay kit (Solarbio, Beijing, China) . After one-way ANOVA, treatment means were assessed by DMRT in SPSS Statistics version 20.0 and were considered to be significantly different between the designated treatments at P < 0.05.

Isolation of RNA and cDNA library construction
The same samples used to analyze IAA and H 2 O 2 content were employed for RNA-seq analysis. Samples collected from 0, 2, 4, 6, 8, 10, and 12 days were marked as ER0, ER2, ER4, ER6, ER8, ER10, ER12, respectively, and stored at − 80 ℃. The Column Plant RNA OUT Extraction kit (Tiandz, Beijing, China) was used to isolate total RNA from each sample, using the methods suggested by the manufacturer. The concentration and quality of all RNA samples was examined by agarose electrophoresis on an Agilent 2100 bioanalyzer (Agilent Technologies, Palo Alto, CA, USA). Sequencing libraries were generated using the TruSeq RNA Sample Preparation Kit (Illumina, San Diego, CA, USA). Magnetic beads with oligo (dT) were used to purify mRNA, which was fragmented into short fragments (200-300 bp). Cleaved mRNA fragments were primed with a random hexamer primer for first-strand and second-strand cDNA synthesis. After purification, end repair, and ligation to sequencing adapters, 21 cDNA libraries of three biological replicates for each treatment were prepared and sequenced using the Illumina Novaseq 6000 platform by Personal Biotechnology Co., Ltd. (Shanghai, China).

Quantitative reverse transcription-polymerase chain reaction (qRT-PCR) analysis
For qRT-PCR analysis, ten candidate DEGs (Six D E G s , T R I N I T Y _ D N 1 1 7 3 6 _ c 0 _ g 1 , T R I N -ITY_DN14475_c0_g1, TRINITY_DN299_c1_g1, TRINITY_DN4858_c0_g3, TRINITY_DN5423_c0_ g2 and TRINITY_DN5557_c0_g2 were identified as log2(FC) > 5; three DEGs, TRINITY_DN5677_c0_g1, TRINITY_DN2748_c0_g1 and TRINITY_DN5677_c0_g2, were enriched in "plant hormone signal transduction"; One DEGs, TRINITY_DN7281_c1_g1, was enriched in "tryptophan metabolism" pathway) were randomly chosen to validate the transcriptomic data. qRT-PCR was performed with the LightCycler 480 System (Roche Diagnostics, Mannheim, Germany) using PerfectStart Green qPCR Supermix (TransGen Biotech, Beijing, China). E. excelsum actin was used as the internal control and the 2 −ΔΔCt method (Livak and Schmittgen 2001) was used to analyze the differential expression of candidate DEGs. Gene-specific primers were designed by Primer Premier 5.0 and are listed in Table S1. Three biological replicates and three technical replicates were performed for each candidate gene.

Adventitious root formation during in vitro and ex vitro rooting
During in vitro rooting, compared with the control group, NAA did not induce ARs whereas IBA or IAA could. High concentrations of IBA and IAA inhibited rooting percentage (Fig. 1a), root number (Fig. 1c) and root length (Fig. 1e). Highest rooting percentage (72.50%) was obtained by 0.5 μM IAA at 60 d after treatment. During ex vitro rooting, NAA, IBA and IAA significantly increased rooting percentage (Fig. 1b), root number (Fig. 1d) and root length (Fig. 1f) compared with the control group. A low concentration of IBA (1 mM) most effectively induced rooting, resulting in the highest rooting percentage (98.33%) and root length (2.72 cm) at 25 d after treatment.
Ex vitro rooting induced the highest rooting percentage (98.33%) at 25 days, while highest rooting percentage during in vitro rooting was 72.50% at 60 days. Thus, ex vitro rooting induced AR from E. excelsum plantlets earlier (faster) than in vitro rooting. The samples collected from ex vitro rooting were used for the next analysis. 8 days after the 1 mM IBA ex vitro rooting treatment, AR primordia were evident, and ARs emerged from the epidermis after 10 days (Fig. 2a, b). ARs elongated, rooting percentage was almost 100% by 25 days, and plantlet survival reached 100%.

IAA and H 2 O 2 content analysis
In the 1 mM IBA treatment during ex vitro rooting, IAA content in stem bases increased gradually from 0 to 8 days, then dropped at 10 and 12 days. A sharp increase in IAA content was observed at 8 days (Fig. 3a). The trend of H 2 O 2 content was different from that of IAA content (Fig. 3b). H 2 O 2 accumulated rapidly after treatment, peaked at 2 days, then sharply decreased at 8 days. The highest content of IAA and lowest content of H 2 O 2 at 8 days corresponded to the timing of AR primordia formation.

De novo assembly and sequence analysis
To identify genes involved in AR induction of E. excelsum plantlets during ex vitro rooting, 21 cDNA libraries were prepared from three repeat mRNA samples collected from 0 (ER0), 2 (ER2), 4 (ER4), 6 (ER6), 8 (ER8), 10 (ER10) and 12 (ER12) days after 1 mM IBA treatment (Table 1). The total number of raw reads produced for each library ranged from 42,845,192 to 52,575,518 with Q20 > 97.48% and Q30 > 93.42%. After filtering, the clean reads per library ranged from 39,721,158 to 49,142,018 with the percentage of clean reads > 91.07% (Table 1). Trinity software v2.5.1 was used to assemble clean reads and obtain transcripts and unigenes for subsequent analysis. The quality and length distribution of transcripts and unigenes are shown in Table 2 and Fig. S1, respectively.
The unigenes were processed in six databases to perform best hits by Blast with E values < 10 -5 , and inferred putative functions of the sequences were assigned. A total of 52,188 (40.15%) unigenes were matched to known genes in the NR database, 23,159 (17.82%) sequences to Pfam and 37,827 (29.10%) sequences in the Swiss-Prot database ( Table 3). The NR database queries revealed that the annotated unigenes were assigned with a best score to sequences from Adventitious root induction of Euryodendron excelsum shoots during in vitro rooting at 60 days (a, rooting percentage; c, root number; e, root length) and ex vitro rooting at 25 d (b, rooting percentage; d, root number; f, root length) after treatment. Bars indicate means ± SE. Different letters indicate statistically significant differences based on Duncan's multiple range test (P < 0.05) between the designated treatments. Forty shoots were prepared for each treatment, and three replicates were performed for each treatment the top seven species (Fig. 4): Vitis vinifera (21.76%), Theobroma cacao (4.34%), Coffea canephora (4.01%), Nelumbo nucifera (3.88%), Sesamum indicum (3.13%), Ziziphus jujuba (2.67%) and Manihot esculenta (2.25%).
The annotation of GO terms revealed that 24,939 unigenes (19.19%) were assigned to biological processes, molecular functions, and cellular components (Fig. S2). Most annotated unigenes in biological processes were involved in "cellular process", "metabolic process", and "single-organism process". In the cellular component category, most annotated unigenes were annotated as "cell", "cell part" and "membrane". In the molecular functions, most annotated unigenes were categorized as "binding", "catalytic activity" and "transporter activity".
A total of 22,160 unigenes (17.05%) and 33 pathways were assigned based on metabolism, genetic information processing, environmental information processing, cellular processes and organismal systems pathway (Fig. S3). On the basis of KEGG analysis, most unigenes were annotated into "carbohydrate metabolism" of metabolism, "translation" of genetic information processing, "signal transduction" of environmental information processing, "transport and catabolism" of cellular processes, and "endocrine system" of organismal system.
The possible functions of unigenes were predicted and classified by alignment to the eggNOG database. A total of 50,632 unigenes (38.95%) were distributed into 25 categories (Fig. S4). Among them, the NOG category "general function prediction only" represented the largest group, followed by "function unknown", "signal transduction mechanisms", and "posttranslational modification, protein turnover, chaperones".

DEGs in response to IBA-induced ex vitro rooting
Hierarchical clustering was used to analyze the expression patterns of DEGs in ER0, ER2, ER4, ER6, ER8, ER10 and ER12 libraries with three biological replicates. These DEGs were divided into nine main clusters (Fig. S5). DEGs in cluster 1, 2, 3 and 4 always showed high expression in the ER0 library with different trends in the other five libraries. The remaining five clusters represented DEGs with high expression levels induced by IBA treatment. The highest number of up-regulated genes was observed in the ER2 library (7364) and fewest in the ER8 library (5649) (Fig. 5a). Upset plot diagram analysis showed that 4635 unigenes maintained differential expression after IBA-induced treatment from 2 to 12 days (Fig. 5b).

GO enrichment analysis
According to GO enrichment analysis, the degree of enrichment was measured based on the rich factor (higher rich factor represents greater enrichment), the FDR value (range from 0 to 1; a score close to 0 indicates more significant enrichment) and the number of genes enriched to a GO term. The significant enrichment GO terms of DEGs showed a few differences in the six libraries (Fig. S6). In the ER2 library, the significantly enriched terms were "monooxygenase activity", "response to auxin", "oxidoreductase activity, acting on paired donors, with incorporation or reduction of molecular oxygen" and "oxidoreductase activity, acting on paired donors, with incorporation or reduction of molecular oxygen, NAD(P)H as one donor, and incorporation of one atom of oxygen". More enrichment terms were categorized into biological processes (BP) and molecular functions (MF). In the ER4 library, more enrichment terms were categorized into cellular component (CC) and MF, and the significantly enriched terms were the same as the ER2 library, but the term "response to auxin" was replaced by "photosynthesis, light reaction". In the ER6 and ER8 libraries, more enrichment terms were categorized into CC and MF, and "oxidoreductase activity, acting on paired donors, with incorporation or reduction of molecular oxygen, NAD(P)H as one donor, and incorporation of one atom of oxygen", "monooxygenase activity" and "photosystem" were listed as significantly enriched terms. In the ER10 library, more enrichment terms were categorized into CC and MF, and "photosynthesis, light harvesting", "chlorophyll binding" and "photosystem I" represented the main significantly enriched terms. In the ER12 library, more enrichment terms were categorized into BP and MF, and the significantly enriched terms were "hydrogen peroxide metabolic process", "hydrogen peroxide catabolic process" and "phenylpropanoid metabolic process".
Besides "hydrogen peroxide metabolic process" and "hydrogen peroxide catabolic process", several terms related to adversity stress were also identified in the GO enrichment analysis (Fig. 6). The terms "hydrogen peroxide metabolic process" and "hydrogen peroxide catabolic process" shared the same number and type of DEGs, most of which, including DEGs for the "response to oxidative stress" term, were up-regulated at 8 days after ex vitro treatment with IBA (Table S2). Most DEGs were associated with the term "response to abiotic stimulus", followed by "response to oxidative stress" and "response to biotic stimulus", while Different letters indicate statistically significant differences based on Duncan's multiple range test (P < 0.05) between the designated treatments. Three biological replicates of 10 cut stem bases were harvested as 0.1 g fresh weight (FW) to assess endogenous IAA and H 2 O 2 content the fewest DEGs were associated with the term "response to hydrogen peroxide" (Fig. 6a). Most of the DEGs in the terms "response to abiotic stimulus" and "response to biotic stimulus" were up-regulated throughout the entire process of AR formation (Table S2). The terms "hydrogen peroxide metabolic process", "hydrogen peroxide catabolic process" and "response to oxidative stress" encompassed 43 DEGs simultaneously (Fig. 6b), and these were mainly identified Reads (No.)-total number of reads, Clean reads (No.)-number of high quality sequence reads, Bases (bp)-total number of bases, Clean data-number of high quality sequence bases, N (%)-proportion of unknown nucleotides in clean reads. Q20 (%) is the percentage of base recognition accuracy of more than 99.0% of bases. Q30 (%) is the percentage of base recognition accuracy of more than 99.9% of bases  Table 2 The quality of transcripts and unigenes in Euryodendron excelsum N50 (bp)-All sequences are arranged from long to short, and the sequence length are added in this order. When the added length reaches 50% of the total length of the sequence, the length of the last sequence is N50, N90 (bp)-All sequences are arranged from long to short, and the sequence length are added in this order. When the added length reaches 90% of the total length of the sequence, the length of the last sequence is N90, N50 Sequence No.-The total number of sequences that the length is greater than N50, N90 Sequence No.-The total number of sequences that the length is greater than N90, GC%-The proportion of guanidine and cytosine nucleotides among total nucleotides as genes related to cationic peroxidase and peroxidase (Table S2).
The 12 up-regulated genes related to auxin-induced proteins, identified as AX6B, AX10A, A10A5, AX15A, AUX22, AUX22D, and AUX28, were sharply enriched in the "plant hormone signal transduction" pathway (Fig. 8a). Most of those genes were also differentially expressed at an early stage (2, 4 days) after IBA treatment. Among these genes, TRINITY_DN39834_c0_g1 and TRIN-ITY_DN31248_c0_g1 were extremely highly differentially expressed with log2(FC) > 8 in the ER0 vs ER2 and ER0 Fig. 4 Species distribution of the top BLAST hits against the NR database in Euryodendron excelsum. Pie charts were generated by R software vs ER4 comparisons. Four DEGs were up-regulated with log2(FC) > 1 by IBA treatment under ex vitro rooting in the ER0 vs ER2, ER0 vs ER4, ER0 vs ER6 and ER0 vs ER8 comparisons while TRINITY_DN5677_c0_g1 maintained up-regulated expression at all stages.
In addition, five up-regulated LAX genes (Fig. 8b) were significantly enriched in the "plant hormone signal transduction" pathway. Only TRINITY_DN9654_c0_g1 was up-regulated at all stages while the other four LAX genes were up-regulated at ER10 and ER12, except for TRIN-ITY_DN380_c4_g1, which was up-regulated at ER6.

QRT-PCR analysis of gene expression
To further validate the results from the RNA-seq data, 10 candidate DEGs related to adventitious root formation were selected for qRT-PCR analysis of E. excelsum samples that were collected 0, 2, 4, 6, 8, 10 and 12 days after 1 mM IBA treatment during ex vitro rooting. In the seven time points, the expression trend of the unigenes from qRT-PCR and RNA-seq analysis were largely consistent (Fig. 9). These results demonstrate that the transcriptome data accurately reflects the ex vitro response of IBA-induced AR formation of E. excelsum plantlets. Upset plot diagram was generated by R software. a Number of upand down-regulated DEGs in 2 (ER2), 4 (ER4), 6 (ER6), 8 (ER8), 10 (ER10) and 12 (ER12) days treatments compared with 0 days (ER0). The x-axis represents the comparison group, while the y-axis represents the number of DEGs; b Upset plot diagram analysis of the DEGs. The connections on the x-axis between points in vertical lines represent the intersection between corresponding data sets, while the y-axis represents the number of DEGs in each intersection

Discussion
AR development is a vital step in plant vegetative propagation, such as in vitro propagation and cuttings. Rooting recalcitrance is a critical factor limiting the application and further development of vegetative propagation (Diaz-Sala 2020; Stevens et al. 2018). IBA is the most frequently used plant hormone for clonal propagation in horticulture and forestry. Although IAA is a primarily native auxin in plants, IBA is more stable and effective in promoting ARs (Ludwig-Müller et al. 2005;Quan et al. 2017;Rout 2006). It is necessary to screen PGRs to find the optimal speciesconcentration ratio for AR induction during in vitro culture. In this study, IBA and IAA treatment significantly promote AR formation of E. excelsum, especially during ex vitro rooting. Furthermore, ex vitro rooting was more suitable for E. excelsum plantlets, with a higher rooting percentage and earlier rooting than in vitro rooting. Ex vitro rooting of shoots has also been applied to many difficult-to-root woody plant species, such as pistachio (Benmahioul et al. 2012), Dalbergia sissoo Roxb. (Vibha et al. 2014) and Bauhinia racemosa Lam (Sharma et al. 2017). In general, the chances of root damage during transplantation to substrates are less during ex vitro rooting, and plantlets tend to be more vigorous, allowing them to cope with environmental stresses during hardening (Arya et al. 2003;Vengadesan and Pijut 2009). Thus, ex vitro rooting is an obvious choice for AR induction during the micropropagation of woody species with further improvement in the choice of PGRs, substrates and other factors.
E. excelsum plantlets experienced a radical environmental change from in vitro aseptic conditions to open ex vitro rooting conditions, which may constitute an adversity stress. In the annotation and enrichment analysis of GO terms, we identified multiple DEGs involved in H 2 O 2 -related biological activities, oxidative stress, abiotic and biotic stimulus. AR formation is also a stress response of plants under adversity stress, and plays a key function in the adaptation of plants to abiotic and biotic stresses (Bellini et al. 2014;Steffens and Rasmussen 2016). The external environmental change may stimulate oxidative damage and increase the production of reactive oxygen species in plants (You and Chan 2015). H 2 O 2 is viewed mainly as a type of reactive oxygen species and a signaling messenger of many biological processes in plants, such as fruit growth and development (Khandaker et al. 2012), leaf senescence (Lin et al. 2019), stomatal closure (Zhang et al. 2019a), and root growth (Xiong et al. 2015). H 2 O 2 and IBA may also act synergistically to regulate adventitious rooting, dependent on the auxin pathway, in marigold explants (Liao et al. 2011). The exogenous application of H 2 O 2 to cucumber plants significantly increased the emergence of ARs (Li et al. 2016b). In this study, the wide fluctuation of endogenous IAA and H 2 O 2 content in E. excelsum plantlets were observed at the stage of root primordia formation. And most DEGs involved in significantly enriched pathway of "plant hormone signal transduction" were up-regulated at the stage corresponded to the timing of H 2 O 2 accumulation. Those results indicate that adversity stress may have a positive effect on AR induction of E. excelsum plantlets under the synergistic action of AR formation involves a series of responses by genes, proteins and metabolites. Multiple biological activities and pathways have specific roles during AR development (de Almeida et al. 2020;Wei et al. 2014). In mungbean seedlings, KEGG pathway enrichment during transcriptomic analysis showed that ribosome biogenesis, plant hormone signal transduction, pentose and glucuronate interconversions, photosynthesis, phenylpropanoid biosynthesis, sesquiterpenoid and flavonoid biosynthesis, and phenylalanine metabolism were the pathways most highly regulated by IBA-induced AR formation, indicating their potential contribution to adventitious rooting (Li et al. 2016a). For apple rootstocks, the most heavily enriched KEGG pathways involved in AR formation were metabolic, biosynthesis of secondary metabolites, plant hormone signal transduction, phenylpropanoid biosynthesis and phenylalanine metabolism pathways, etc. (Li et al. 2018). In sugarcane shoots, DEGs associated with plant hormone signaling, flavonoid and phenylpropanoid biosynthesis, cell cycle, and cell wall modification, and transcription factors were involved in AR formation . During AR development of E. excelsum, we found that more DEGs were enriched in "plant hormone signal transduction" and "phenylpropanoid biosynthesis" pathways, similar to a number of previous studies. Therefore, we conclude that these two pathways have a vital influence on AR formation in E. excelsum plants.
IAA is the most abundant natural auxin, and endogenous IAA is closely related to the development of ARs in plants. The conversion of exogenous hormones to endogenous auxin and the synthesis of auxin are key factors regulating AR development (Olatunji et al. 2017). Tissue that produces ARs requires high levels of auxin, and the enrichment of high concentrations of auxin depends on polar auxin transport (Ahkami et al. 2013;Garrido et al. 2002). Thus, genes related to the synthesis, signaling and polar transport of auxin, like AUX, LAX, and PIN, are closely related to plant adventitious rooting (Druege et al. 2016). For example, auxin influx carriers MiAUX3 and MiAUX4 might play important roles during AR formation in mango cotyledon segments, and the expression levels of MiAUX3 and MiAUX4 resulted in a significant promotive effect of IBA on adventitious rooting (Li et al. 2012). Papaya plantlets not exposed to IBA could not form ARs and displayed a low expression of all auxin transporter genes in stem base tissues whereas IBA-treated plants were able to produce ARs and showed significantly increased expression of most auxin transporter genes, especially CpLAX3 and CpPIN2 (Estrella-Maldonado et al. 2016). In E. excelsum, DEGs for AUX and LAX, which were significantly enriched in plant hormone signal transduction, showed a high fold change during AR development. This implies that the expression patterns of those genes were linked to AR induction from E. excelsum plantlets. AUX/IAA protein is an early auxin response protein that always participates in the auxin signaling pathway by interacting with auxin response factor (ARF) protein or other genes (Salehin et al. 2015). During AR formation in petunia cuttings, the expression of genes of the Aux/IAA family showed strong temporal variation, supporting their important role in the induction and transition to subsequent root formation phases (Druege et al. 2014). The auxin receptor (TRANSPORT INHIBITOR1) TIR1 homolog gene, PagFBL1, interacted strongly with both PagIAA28.1 and PagIAA28.2 in the presence of NAA to regulate AR induction in poplar stem segments (Shu et al. 2019). In Arabidopsis thaliana, Aux/IAA proteins, IAA6, IAA9, and IAA17, interacted with ARF6 and/or ARF8 and likely repressed their activity in AR development, and complexed with TIR1 and (AUXIN-SIGNALLING F-BOX) AFB2 to form specific sensing to modulate jasmonic acid homeostasis and control AR initiation . In this study, 13 IAA genes were significantly enriched in the plant hormone signal transduction pathway, suggesting a significant relationship between AUX/IAA and AR formation in E. excelsum. The mechanism and interaction with other IAA genes would need to be revealed in future research.
SAUR , the largest family of early auxin response genes in plants, mediate the regulation of several aspects of plant growth and development (Ren and Gray 2015). SAUR proteins showed positive or negative effects on primary, lateral and adventitious root development. In A. thaliana, plants overexpressing SAUR41 exhibited increased primary root growth and a higher number of lateral roots (Kong et al. 2013). AtSAUR15 acts downstream of the auxin response factors ARF6,8 and ARF7,19 to regulate auxin signaling-mediated lateral root and AR formation, and plants overexpressing AtSAUR15 exhibit more lateral roots and ARs (Yin et al. 2020). In contrast to AtSAUR41 and AtSAUR15, overexpression of OsSAUR39 in rice resulted in reduced root elongation and lateral root development (Kant et al. 2009). SAUR proteins may display a species-or type-dependent positive function in AR formation. In E. excelsum, three SAUR genes maintained up-regulated expression after IBA-induced treatment, indicating a close association with AR formation.
We also found several highly up-regulated GH3 genes at all stages of AR formation in E. excelsum. GH3 proteins are also an early auxin response protein, play a crucial role in conjugating IAA to amino acids, and are critical in maintaining auxin homeostasis (Brunoni et al. 2020). Three GH3 genes, GH3.3, GH3.5, and GH3.6, were required for fine-tuning AR initiation in A. thaliana hypocotyls (Gutierrez et al. 2012). In cucumber hypocotyls, salicylic acid plays an inducible role in AR formation through competitive inhibition of the auxin conjugation enzyme CsGH3.5, and salicylic acid-induced IAA accumulation was also associated with the enhanced expression of CsGH3.5 (Dong et al. 2020). In apple plants, overexpression of MsGH3.5 significantly reduced the content of free IAA and increased the content of some IAA-amino acid conjugates, and MsGH3.5-overexpressing lines produced fewer ARs than the control (Zhao et al. 2020). Those results demonstrated that GH3 proteins were intricately involved in AR development, but did not only perform a positive role.

Conclusion
Here, we confirmed that ex vitro rooting was an obvious choice for AR formation during the micropropagation of E. excelsum plantlets. DEGs enriched in the pathway of plant hormone signal transduction played a crucial role in AR formation. H 2 O 2 produced by environmental stimulation might be related to AR induction in E. excelsum ex vitro by the synergistic action with IBA, ultimately regulating the level of endogenous IAA. The knowledge gained from this study will help researchers understand the molecular traits of IBA-based regulation of adventitious rooting of E. excelsum plantlets. These results will provide technical support for the ecological protection of this rare and endangered species and are important for research and commercial applications aimed at overcoming rooting recalcitrance in plant species of economic value, in difficult-to-root woody plants, or in rare or endangered plants. Fig. 9 Analysis of the fold changes of 10 candidate genes in Euryodendron excelsum determined by RNA-seq and qRT-PCR. Bars indicate means ± SD. The x-axis represents the time point after IBA treatment under ex vitro rooting while the y-axis represents the log2(fold change). Graphs were generated by Microsoft Excel and Adobe Photoshop CC 2018 ◂