Transcriptomic and Anatomic Profiling Reveal Etiolation Promotes Adventitious Rooting by Exogenous Application of 1-Naphthalene Acetic Acid in Robinia pseudoacacia L.

The process of etiolation contributes significantly to vegetative propagation and root formation of woody plants. However, the molecular interaction pattern of different factors for etiolated adventitious root development in woody plants remains unclear. In the present study, we explored the changes at different etiolation stages of adventitious root formation in Robinia pseudoacacia. Histological and transcriptomic analyses were performed for the etiolated lower portion of hypocotyls to ascertain the adventitious root responses. We found that the dark-treated hypocotyls formed roots earlier than the control. Exogenous application of NAA (0.3 mg/L) stimulated the expressions of about 310 genes. Among these, 155 were upregulated and 155 were downregulated. Moreover, differentially expressed genes (DEGs) were significantly enriched in multiple pathways, including the biosynthesis of secondary metabolites, metabolic pathway, plant hormone signal transduction, starch and sucrose metabolism, phenylpropanoid biosynthesis, and carbon metabolism. These pathways could play a significant role during adventitious root formation in etiolated hypocotyls. The findings of this study can provide novel insights and a foundation for further studies to elucidate the connection between etiolation and adventitious root formation in woody plants.


Introduction
R. pseudoacacia is native to North America and was first introduced to China in 1877 [1]. R. pseudoacacia possesses the most significant ecological and economic characteristics, including excellent coppicing, fast growth, high yield, and adaptability to a wide range of environments [2]. Several reports have shown that in vitro propagation of R. pseudoacacia is an effective method to produce large numbers of clonal plants [3] because woody species Under NAA treatment, the R. pseudoacacia-148 hypocotyl developed a small callus and exhibited root elongation at 72 h in donor-treated plants. In the control, R. pseudoacacia-148 did not exhibit any change during tested stages of adventitious root formation ( Figure 1). Under NAA treatment, the histology at 60 h (2.5 days), 66 h (2.75 days), and 72 h (3 days) showed obvious changes ( Figure 1), but under the control conditions, the hypocotyl tissues up to 72 h did not exhibit any change. At each time point (60 h, 66 h, and 72 h), there was a clear difference between the treatment and control. The small callus was formed on the wounded sites of the hypocotyl in culturing medium at 72 h. Furthermore, the gathering of cells and division was initiated at 60 h, a dome-shaped adventitious root primordium was visible at 66 h, and root emergence from the hypocotyl was recorded at 72 h ( Figure 1). Based on the histological observations, adventitious root formation indicated that dark pretreatment had boosted up the root formation as compared to the control. Anatomical visualization for transcriptome analysis was used to validate the accuracy of the sampling time point. These anatomical observations provide basic information about adventitious root formation in R. pseudoacacia.

Deep Sequencing, Functional Annotation of DEGs, and qRT-PCR Analysis
A total of six samples generated 38-60 million clean reads. On average, 7.72 Gb of data was produced (Table 1). For all sequenced data, the Q30 percentage was greater than 91% ( Table 1). The final mapping result of tested reads against a reference (R. pseudoacacia) was about 85.89% (Table 1). After reference genome alignment, a total of 45,164 new transcripts were found, out of which 29,472 belonged to a novel isoform, 4042 transcripts belonged to the new protein-coding genes, 33,514 were coding transcripts, while the remaining 11,650 were long-chain non-coding RNAs. On average, 24,488 genes were detected per Under NAA treatment, the R. pseudoacacia-148 hypocotyl developed a small callus and exhibited root elongation at 72 h in donor-treated plants. In the control, R. pseudoacacia-148 did not exhibit any change during tested stages of adventitious root formation ( Figure 1). Under NAA treatment, the histology at 60 h (2.5 days), 66 h (2.75 days), and 72 h (3 days) showed obvious changes ( Figure 1), but under the control conditions, the hypocotyl tissues up to 72 h did not exhibit any change. At each time point (60 h, 66 h, and 72 h), there was a clear difference between the treatment and control. The small callus was formed on the wounded sites of the hypocotyl in culturing medium at 72 h. Furthermore, the gathering of cells and division was initiated at 60 h, a dome-shaped adventitious root primordium was visible at 66 h, and root emergence from the hypocotyl was recorded at 72 h ( Figure 1). Based on the histological observations, adventitious root formation indicated that dark pretreatment had boosted up the root formation as compared to the control. Anatomical visualization for transcriptome analysis was used to validate the accuracy of the sampling time point. These anatomical observations provide basic information about adventitious root formation in R. pseudoacacia.

Deep Sequencing, Functional Annotation of DEGs, and qRT-PCR Analysis
A total of six samples generated 38-60 million clean reads. On average, 7.72 Gb of data was produced (Table 1). For all sequenced data, the Q30 percentage was greater than 91% ( Table 1). The final mapping result of tested reads against a reference (R. pseudoacacia) was about 85.89% (Table 1). After reference genome alignment, a total of 45,164 new transcripts were found, out of which 29,472 belonged to a novel isoform, 4042 transcripts belonged to the new protein-coding genes, 33,514 were coding transcripts, while the remaining 11,650 were long-chain non-coding RNAs. On average, 24,488 genes were detected per sample by comparing reads to genes, using the "bowtie2" software. For differential gene analysis, an average of 310 genes were detected, in which 155 were upregulated and 155 were downregulated. The distribution intervals of gene expression level in RNA-seq data of each sample are shown in Figure S1. In the present study, DEGs were analyzed among different comparisons of treatment and with their respective control at 60 h (CK60). There were 12, 51, and 383 DEGs in a total and 4, 51, and 367 specific DEGs between the control and treatment, respectively (Figure 2a showed that the BL-60 group was in the same quadrant, and two replications of each sample were in the same quadrant in BL-60, BL-72, and BL-CK60 ( Figure S2b). Almost all groups were closely associated with each other, which shows a strong correlation ( Figure S2a). A heat map of all DEGs on a selected time point is presented in Figure 2d.
Gene ontology (GO) assignment was performed to classify the functions of DEGs at different time points (60 h, 66 h, and 72 h) ( Figure S3, Table S1). The DEGs between the control and treatment were categorized into different comparisons (BL-CK60-VS-BL-60: 18,  13, Figure S3d-f), divided into the biological process, cellular component, and molecular function, respectively. Most DEGs were assigned to metabolic process, cellular process, biological regulation, and response to stimulus in the category of biological process; cell part, cell, and membrane in the category of cellular component; and binding, catalytic activity, and transporter activity in the category of molecular function ( Figure S3). Within the biological category, GO terms of DEGs including metabolic process (GO:0008152), cellular process (GO:0009987), single-organism process (GO:0044699), and biological regulation (GO:0065007) were significantly enriched ( Figure 3).
To further elaborate the biological interpretation, all DEGs were mapped to the KEGG database ( Figure 4, Table S2). According to KEGG pathway enrichment analysis, 10

Effect of Exogenous NAA on Plant Hormone Signal Transduction Pathway
Transcriptome analysis revealed differentially expressed genes related to phytohormones, cytokinin, auxin, jasmonic acid, brassinosteroids, gibberellin, and salicylic acid. Phytohormones play important roles in the regulation of adventitious root formation. In the present study, the auxin signaling pathway genes (AUX/IAA, SAUR, and GH3), cytokinin signaling genes (B-ARRs, APL, and EFM), abscisic acid-related genes (PP2C, SNRK2-7 and ABI5), salicylic acid signaling genes, jasmonic acid, brassinosteroids (MYC2, TCH4, BRI1, SPL12, TGA9, and PR-1-LIKE), and gibberellin-related SCL15 gene were identified. These genes exhibited changes in expression levels at various stages of adventitious root formation at 60 h, 66 h, and 72 h. A heat map of all DEGs related to hormones on FPKM log 2 fold values and their relative qRT-PCR is presented ( Figure 6, Table S3).

Auxin
The heat map of DEGs is presented in Figure 6a. In the present study, exogenous NAA treatment upregulated the relative expression of auxins, AUX/IAA; indole-3-acetic acid-inducible 19 (gene26601), which decreased from BL-60 h to BL-72 h (Figure 6b). The expression of the SAUR (gene3728) gene decreased from BL-60 h to BL-66 h and then increased from BL-66 h to BL-72 h (Table S3). In contrast, GH3.3 (gene829) was upregulated from BL-60 h to BL-66 h and decreased from BL-66 h to BL-72 h (Figure 6b). The qRT-PCR data showed that an exogenous NAA application induced auxin-related gene expression compared to the control ( Figure 6).

Cytokinin
The relative expressions of cytokinin and two pseudo-response regulators 2(PRR2) were observed under the NAA treatment. The expression of PRR2 (gene957) first decreased then increased from BL-60 h to BL-72 h (Figure 6b), while the expression of PRR2 (gene30838) continuously increased from BL-60 h to BL-72 h (Table S3). The expression of

Auxin
The heat map of DEGs is presented in Figure 6a. In the present study, exogenous NAA treatment upregulated the relative expression of auxins, AUX/IAA; indole-3-acetic acid-inducible 19 (gene26601), which decreased from BL-60 h to BL-72 h (Figure 6b). The expression of the SAUR (gene3728) gene decreased from BL-60 h to BL-66 h and then increased from BL-66 h to BL-72 h (Table S3). In contrast, GH3.3 (gene829) was upregulated from BL-60 h to BL-66 h and decreased from BL-66 h to BL-72 h (Figure 6b). The qRT-PCR data showed that an exogenous NAA application induced auxin-related gene expression compared to the control ( Figure 6).

Cytokinin
The relative expressions of cytokinin and two pseudo-response regulators 2(PRR2) were observed under the NAA treatment. The expression of PRR2 (gene957) first decreased then increased from BL-60 h to BL-72 h (Figure 6b), while the expression of PRR2 (gene30838) continuously increased from BL-60 h to BL-72 h (Table S3). The expression of another two genes, altered phloem development (APL) (gene25497) and early flowering myb protein (EFM) (gene36128), first decreased then increased from BL-60 h to BL-72 h (Table S3).

Gibberellic Acid
There are relatively few studies on GA-regulated root development in woody plants. DELLA, scarecrow-like 15 (SCL15) (gene1989) has been found with decreasing expression to increasing expression from BL-60 h to BL-72 h (Table S3).

Starch and Sucrose Metabolism Pathway
The initiation and development of the adventitious root primordium require more energy triggered by starch, which is a significant source of carbohydrates in woody plants [20]. The heat map of differentially expressed genes (DEGs) and qRT-PCR are presented in Figure 7a. Under NAA treatment, the relative expression level of SUS3 (gene9499) slightly increased from BL-60 h to BL-72 h. In contrast, AIR3 (gene30616) expression continuously decreased from BL-60 h to BL-72 h (Figure 7b). Three genes of the plant invertase/pectin methyl esterase inhibitor superfamily were found. Two genes showed opposite trends. Gene19401 continuously increased from BL-60 h to BL-72 h, and gene1053 continuously decreased from BL-60 h to BL-72 h (Table S3). Gene198 first decreased from BL-60 h to BL-66 h and then increased from BL-66 h to BL-72 h (Table S3). Gene17148 showed slightly decreasing expression trends, and while the gene1118 downregulated from BL-60 h to BL-66 h and then upregulated from BL-66 h to BL-72 h (Table S3), respectively. In addition, BMY5 (gene15899) showed an increasing trend from BL-60 h to BL-66 h and then decreased from BL-66 h to BL-72 h (Figure 7b). RT-qPCR also validated the expression pattern of selected genes. Under NAA treatment, most of the expression profiles for the aforementioned genes were consistent with the profiles obtained by RT-qPCR ( Figure 7).

Phenylpropanoid Biosynthesis Pathway
Heat maps of all DEGs expressed during phenylpropanoid biosynthesis were examined during the three tested time points (60 h, 66 h, 72 h) (Figure 8a) of adventitious root formation, which were all upregulated. Among phenylpropanoid DEGs, Arabidopsis thaliana peroxidase CB (gene5) expression increased continuously from BL-60 h to BL-72 h (Figure 8b). In addition, peroxidase superfamily protein (gene12931) (Figure 8b) and peroxidase 2 (gene12930) showed a similar trend, first increasing from BL-60 h to BL-66 h and then decreasing from BL-66 h to BL-72 h (Table S3). Contrarily, the expression of two genes, peroxidase 40 (gene22473) (Figure 8b) and telomeric DNA-binding protein 2 (gene22265), increased continuously from BL-60 h to BL-72 h (Table S3). Most of the expression profiles for the aforementioned genes were consistent with the profiles obtained by RT-qPCR ( Figure 8).

Phenylpropanoid Biosynthesis Pathway
Heat maps of all DEGs expressed during phenylpropanoid biosynthesis were examined during the three tested time points (60 h, 66 h, 72 h) (Figure 8a) of adventitious root formation, which were all upregulated. Among phenylpropanoid DEGs, Arabidopsis thaliana peroxidase CB (gene5) expression increased continuously from BL-60 h to BL-72 h (Figure 8b). In addition, peroxidase superfamily protein (gene12931) (Figure 8b) and peroxidase 2 (gene12930) showed a similar trend, first increasing from BL-60 h to BL-66 h and then decreasing from BL-66 h to BL-72 h (Table S3). Contrarily, the expression of two genes, peroxidase 40 (gene22473) (Figure 8b) and telomeric DNA-binding protein 2 (gene22265), increased continuously from BL-60 h to BL-72 h (Table S3). Most of the expression profiles for the aforementioned genes were consistent with the profiles obtained by RT-qPCR ( Figure 8).

Carbon Metabolism Pathway
The heat map of all DEGs involved in carbon metabolism is presented in Figure 9a, and carbon metabolism is an essential pathway in plants for the production of structural components and energy sources. DEGs expressed in carbon metabolism might be playing a plausible role during adventitious root formation. These genes included PSAT1 (phosphoserine aminotransferase 1), SERAT3;2 (serine acetyltransferase 3;2), GOX1 (glycolate oxidase 1), FBA1 (fructose-bisphosphate aldolase 1), PFK7 (phosphofructokinase 7), PRK (phosphoribulokinase) and SHM1 (serine hydroxy methyltransferase 1). With NAA treatment, the relative expression level of GOX1 (gene25254), FBA1 (gene16722), and PFK7 (gene23954) first decreased from stage BL-60 h to BL-66 h and then increased from stage BL-66 h to BL-72 h, and the expression of the first two genes was observed at BL-72 h and for the last one at BL-60 h, respectively (Figure 9b). The relative expression continuously increased from BL-60 h to BL-72 h in the case of SERAT3;2 (gene2215) (Figure 9b) and PRK (gene30172) (Table S3), with the highest expression at BL-72 h. PSAT1 (gene30347) exhibited an increasing trend from BL-60 h to BL-66 h and then decreased from BL-66 h to BL-72 h (Table S3). Most of the above-mentioned genes showed positive expression, but SHM1 (gene30059) showed negative expression at 60 h and 66 h and then changed to positive expression at 72 h (Figure 9a). Collectively, the RT-qPCR results indicated that the relative expression levels of carbon metabolism-related genes were affected by the exogenous NAA application compared to the control, while results were consistent with RNA-sequencing ( Figure 9).

Carbon Metabolism Pathway
The heat map of all DEGs involved in carbon metabolism is presented in Figure 9a, and carbon metabolism is an essential pathway in plants for the production of structural components and energy sources. DEGs expressed in carbon metabolism might be playing a plausible role during adventitious root formation. These genes included PSAT1 (phosphoserine aminotransferase 1), SERAT3;2 (serine acetyltransferase 3;2), GOX1 (glycolate oxidase 1), FBA1 (fructose-bisphosphate aldolase 1), PFK7 (phosphofructokinase 7), PRK (phosphoribulokinase) and SHM1 (serine hydroxy methyltransferase 1). With NAA treatment, the relative expression level of GOX1 (gene25254), FBA1 (gene16722), and PFK7 (gene23954) first decreased from stage BL-60 h to BL-66 h and then increased from stage BL-66 h to BL-72 h, and the expression of the first two genes was observed at BL-72 h and for the last one at BL-60 h, respectively (Figure 9b). The relative expression continuously increased from BL-60 h to BL-72 h in the case of SERAT3;2 (gene2215) (Figure 9 b) and PRK (gene30172) (Table S3), with the highest expression at BL-72 h. PSAT1 (gene30347) exhibited an increasing trend from BL-60 h to BL-66 h and then decreased from BL-66 h to BL-72 h (Table S3). Most of the above-mentioned genes showed positive expression, but SHM1 (gene30059) showed negative expression at 60 h and 66 h and then changed to positive expression at 72 h (Figure 9a). Collectively, the RT-qPCR results indicated that the relative expression levels of carbon metabolism-related genes were affected by the exogenous NAA application compared to the control, while results were consistent with RNAsequencing ( Figure 9).

Discussion
Since antiquity, root formation has been improved by propagating avocado stock plants in darkness [21,22], signifying the role of etiolation in the adventitious root for-

Discussion
Since antiquity, root formation has been improved by propagating avocado stock plants in darkness [21,22], signifying the role of etiolation in the adventitious root formation of stem cuttings. Etiolation maintains a balance among phytohormones (auxin and cytokinin) in different conditions to aid developmental processes and organogenesis [23]. Etiolation declines the lignification of shoots which results in less development of sclerenchyma tissues. Therefore, the penetration of phytohormones into plants is easy because of less accumulation of cell wall deposits and thin cell walls [24].
The present study focused on identifying the key DEGs with their functions in the process of etiolated hypocotyls in R. pseudoacacia during the formation of the adventitious root. We examined the global gene expression of the control and lower segments of etiolated hypocotyl treated with NAA at three time points (60 h, 66 h, and 72 h). The RNA-seq data depicted that Glycine max had the greatest number of hits from R. pseudoacacia during transcriptome assembly, reflecting the high sequence resemblance.
Based on enriched pathways (metabolic pathway, biosynthesis of secondary metabolites, phenylpropanoid biosynthesis, protein processing in the endoplasmic reticulum, carbon metabolism, starch and sucrose metabolism, and plant hormone transduction) analysis, we identified key functional genes, playing a significant role in adventitious root formation with exogenous NAA.

Histological Changes for the Process of Adventitious Root Formation in R. pseudoacacia
Generally, three stages were observed during adventitious root formation, including induction, initiation, and elongation [21]. The origin of the dome-shaped primordium was near the vascular cambium and secondary phloem tissues. Root formation is a continuous process and, therefore, cell division, cell elongation, and cell differentiation occur continuously to form root primordia [22] and differentiated adventitious root [23]. In the present study, we observed that the root primordia emerged out from the hypocotyl's outer layer at 72 h of hormone treatment (Figure 1). The histology at 60 h, 66 h, and 72 h showed obvious changes ( Figure 1); however, in control conditions, the hypocotyl tissues up to 72 h did not exhibit any growth. At each time point (60 h, 66 h, and 72 h), there was a clear difference between the treatment and control. The small callus was formed on the wounded sites of the hypocotyl, culturing in the medium at 72 h. Anatomical changes also showed that both new and old primordia were present in which the organization and differentiation of cells had started, which contributed to the formation of adventitious roots in poplar "NL895". The results are consistent with previous findings [24,25]. Therefore, it is suggested that there were no clear boundaries between stages.

Functional Annotation of DEGs and Expression Analysis
RNA-sequencing was performed for understanding the gene expression changes during the formation of adventitious roots in R. pseudoacacia. The data depicted that the significant number of DEGs during the early stages of adventitious root formation are not similar with the DNA microarray analysis.
Significantly enriched GO terms for the DEGs in all comparisons included biological processes, including the metabolic process, catalytic activity, and cellular process (Figure 3). According to the KEGG enrichment analysis, pathways associated with metabolic processes and the biosynthesis of secondary metabolites were the most influenced by etiolation. For example, the 40S ribosomal protein S5 A, classified in the ribosome pathway, was expressed only in etiolated tissues. The 40S ribosomal protein S5 A encodes the ribosomal protein S5. Cell division processes were delayed or distributed in a heterozygous Arabidopsis thaliana mutant, and development was arrested completely at an early embryonic stage in a homozygous mutant [26]. This indicates that cell division and the development of etiolated tissues may be postponed or restricted.

Auxin
In the present study, a significant number of DEGs were expressed in the auxin signaling pathways (AUX/IAA, SAUR, and GH3), cytokinin signaling (B-ARRs, APL, and EFM), abscisic acid signaling (PP2C, SNRK2-7 and ABI5), salicylic acid signaling, jasmonic acid, brassinosteroids (MYC2, TCH4, BRI1, SPL12, TGA9, and PR-1-LIKE), and gibberellin. These genes were expressed at various stages of adventitious root formation at 60 h, 66 h, and 72 h. The adventitious root formation produced by wounding in cuttings is dependent on polar auxin transport and the early gathering of IAA in the rooting zone [27]. These results are inconsistent with previous findings reported by Ren,Hu [28] that during adventitious root formation in R. pseudoacacia stem cuttings, the most auxin-related genes (AUX/IAA(IAA19), GH3, and SAUR) were upregulated with exogenous auxin treatment [29]. In addition, the GH3.3 gene was upregulated during adventitious root formation through hypocotyl tissue culture for hypocotyls of R. pseudoacacia, which contradicts with the previous findings in other species [30]. Previously, SAUR and GH3. The type B cytokinin response regulators greatly influence the number of adventitious roots by increasing or decreasing the expression levels. The type B cytokinin response regulators greatly influence the number of adventitious roots by increasing or decreasing the expression levels were upregulated, and IAA was downregulated in etiolated walnut stems [31]. However, in rice, the OsGH3 gene was upregulated on auxin treatment, and GH3.3 was positively correlated with the regulation of adventitious root formation [30] The above mentioned genes may play a key role during the adventitious rooting of R. pseudoacacia.

Cytokinin
Cytokinin is antagonistic to auxin and mainly acts as a suppressor of adventitious root formation in many species [32]. There are two main types of cytokinin (type A and type B), which positively regulate cytokinin signaling. A-ARRs negatively regulate the response of cytokinin and are transcriptionally regulated by B-ARRs [33]. We found that genes B-ARRs (gene957, gene30838), APL (25497), and EFM (gene36128) are involved and encode type B-ARR, which showed upregulation on three tested time points (60 h, 66 h, 72 h). The relative expressions of cytokinin and two pseudo-response regulators 2(PRR2) were observed under the NAA treatment. The expression of PRR2 (gene957) first de-creased and then increased from BL-60 h to BL-72 h (Figure 6b), while the expression of PRR2 (gene30838) continuously increased from BL-60 h to BL-72 h (Table S3). The type B cytokinin response regulators greatly influence the number of adventitious roots by increasing or decreasing the expression levels [32,34], whereas APL is essential for maintaining phloem cell identity [35]. AtMYB96 acts as an important molecular link between salicylic acid and abscisic acid crosstalk, through which it can enhance the resistance against pathogens in Arabidopsis thaliana [36]. These results indicate that the differential expression of type B-ARRs may play a role in R. pseudoacacia adventitious root formation induced by auxin. These proposed interactions and the mechanisms underlying them need further experimental proof.

Abscisic Acid
It is well known that abscisic acid acts as a repressor in root initiation and negatively affects adventitious root formation [37]. Our data indicated that under NAA treatment, the differentially expressed genes associated with abscisic acid are PP2C (gene22236), SNRK2-7 (gene4461), and ABI5 (gene9262). The relative expression of abscisic acid, protein phosphatase 2CA (PP2C, gene22236) continuously increased from BL-60 h to BL-72 h under NAA treatment. ABI5 is a primary leucine zipper transcription factor that functions in the core abscisic acid signaling pathway comprising PYR/PYL/RCAR receptors, the negative regulator of PP2C phosphatases, and the positive regulator of SnRK2, through the regulation of gene expression. Moreover, during partial submergence, the abscisic acid concentration in the petioles of Rumex palustris and stems of S. dulcamara sharply decreased due to a fast downregulation of abscisic acid biosynthesis and upregulation of abscisic acid breakdown [38]. In addition, the abscisic acid downstream signaling gene ABF is upregulated in Red Fife [39]. Abscisic acid insensitive 5 (ABI5, gene9262) expression decreased from BL-60 h to BL-66 h and increased from BL-66 h. Additionally, abscisic acid induces ABI5 expression in the lateral root tips. These results indicate that ABI5 acts as a negative regulator of lateral root development in the presence of stress [40], which may play a role during adventitious root formation of R. pseudoacacia. These results indicate that differential expression of abscisic acid genes may play a vital role during the adventitious rooting of R. pseudoacacia.

Gibberellic Acid
The BnSCL1 gene is regulated by auxin and is functionally associated with root development [41]. Scarecrow-like genes are induced in rooting-competent cells at the earliest stages of adventitious root formation in the presence of exogenous auxin [42]. In this study, DELLA, scarecrow-like 15 (SCL15) (gene1989) has been found with decreasing expression followed by increasing expression pattern from BL-60 h to BL-72 h. The exogenous application of GA inhibits adventitious root formation in rice, and rice mutants deficient in GA biosynthesis develop more adventitious roots [43]. GA promotes the proteasome-mediated degradation of DELLA proteins [44]. Therefore, it is concluded that GA, through DELLA, repressed the expression of other genes [45]. Thus, GA-related genes may mediate the interaction of GA and auxin to regulate adventitious root formation.
Jasmonic acid is a positive or a negative regulator of adventitious root formation [48]. MYC2 and bHLH are very important for the entire jasmonic acid signaling pathway [49], and help in the regulation of root growth and development, cellular proliferation, and differentiation [50]. The expression of TCH4 is upregulated by darkness, increased auxin, and brassinosteroids [51]. The expression of the TCH4-GUS reporter gene in transgenic Arabidopsis thaliana is significantly higher in growing parts, such as the elongating hypocotyls of seedlings grown under low light intensity, auxin and brassinosteroids. Increased TCH4 expression may be an important step in hormone-induced alterations for cell expansion [52]. Moreover, the expressions of SPL genes may be regulated by light and phytohormones [53]. SPLs are also involved in the appropriate timing of the root developmental process [54]. The relative expression of salicylic acid, motif-binding protein 9 (TGA9) (gene24193) expression continuously decreased from BL-60 h to BL-72 h (Table S3). Salicylic acid production, which usually occurs as a response to various stresses, has also been reported to trigger plant root development. PR1 (gene13913) and TGA9 (gene24193) are important genes in salicylic acid signaling. The present findings suggested that genes involved in auxin, abscisic acid, brassinosteroid, jasmonic acid, cytokinin, gibberellin, and salicylic acid signaling may interact with auxin to mediate adventitious root formation. However, the interaction of hormones with auxin during adventitious root formation in R. pseudoacacia is not clear.

Starch and Sucrose Metabolism Pathway
Auxins help in the mobilization and hydrolysis of starch and sugars to the cutting base [55]. The DEGs encoded by the starch and sucrose metabolism include SUS3 (sucrose synthase 3), AIR3 (auxin-induced in root cultures 3), ADG2 (ADP glucose pyro phosphorylase 2), plant invertase/pectin methyl esterase inhibitor superfamily, PMEI9, PMEI3 (pectin methylesterase inhibitor 19, 13) and BMY5, BAM6 (beta-amylase 5, beta-amylase 6). In the present study, the relative expression level of SUS3 (gene9499) slightly increased from BL-60 h to BL-72 h (Figure 7b). β-Amylase helps in the degradation of starch [55], maltose, sucrose, and cell wall components, but decreases starch content [56], which means etiolation requires more energy as compared to the control. In addition, SUSs, as a biocatalyst, encode sucrose synthase and improve adventitious root formation during lateral stages [57]. Additionally, AIR3 encodes subtilisin-like proteases [58] that help to facilitate easy root emergence by weakening the connection among cells in Arabidopsis thaliana during root development [59]. In contrast, AIR3 (gene30616) expression continuously decreased from BL-60 h to BL-72 h (Figure 7b). Three genes of the plant invertase/pectin methyl esterase inhibitor superfamily were found. The DEGs were upregulated in the starch and sucrose metabolism pathway [60]. Moreover, glucose pyro-phosphorylase is a key metabolite and precursor of carbohydrate formation [61]. It plays a role in the formation of the pericycle during starch/sucrose production, which can be used to fine-tune roots in Arabidopsis by the accessibility of carbon sources [62]. Our study shows that adventitious root formation and elongated etiolated seedlings may require more sucrose as an energy source. These expression changes and results concluded that the mechanism in which etiolation promotes adventitious root formation might be linked with more energy supply, hormone synthesis, or different signaling pathways.

Phenylpropanoid Biosynthesis Pathway
The phenylpropanoids and related flavonoids of the phenylpropanoid metabolism pathway might play their role in NAA-induced adventitious root formation. Phenylpropanoids convert phenylalanine to an activated cinnamic acid derivative, which leads to the production of flavonoids or iso-flavonoids [63]. Among phenylpropanoid DEGs, Arabidopsis thaliana peroxidase CB (gene5) expression increased continuously from BL-60 h to BL-72 h (Figure 8b). In addition, peroxidase superfamily protein (gene12931) (Figure 8b) and peroxidase 2 (gene12930) showed a similar trend, first increasing from BL-60 h to BL-66 h and then decreasing from BL-66 h to BL-72 h (Table S3). The phenylpropanoids from the secondary metabolic pathway play an important role in auxin-induced adventitious rooting [64]. Phenolic acids and flavonoids from the phenylpropanoid pathway [65] majorly contribute to cell division and differentiation [66] during in vitro rooting [67]. Therefore, prior changes in phenylpropanoids contribute to Arabidopsis thaliana root development. In this study, the expression of peroxidase 40 (gene22473) (Figure 8b) and telomeric DNAbinding protein 2 (gene22265) increased continuously from BL-60 h to BL-72 h (Table S3). Moreover, specific peroxidases help to maintain the cell division and the differentiation in root meristematic cells [68] during root emergence [69]. In stem cuttings, plant growth, differentiation, and development are controlled by total peroxidase activity [70]. Moreover, peroxidases may promote cell elongation by producing more free radicals that may loosen the cell wall [71], which helps in the penetration of phytohormones for the formation of adventitious roots. Therefore, DEGs involved in the phenylpropanoid pathway, which may have a key role in adventitious root formation in R. pseudoacacia.

Carbon Metabolism Pathway
Carbon metabolism shares sucrose metabolites, which are converted to hexose by SUS and INV (invertase) for the root primordium and root elongation [72]. PSAT1 promotes cell proliferation and chemoresistance in vitro [73]. Therefore, the phosphorylated route for serine synthesis is reportedly associated with meristematic or rapidly proliferating plant tissues [72]. In this study, DEGs in carbon metabolism might be playing a plausible role during adventitious root formation. These genes included PSAT1 (phosphoserine aminotransferase 1), SERAT3;2 (serine acetyltransferase 3;2), GOX1 (glycolate oxidase 1), FBA1 (fructose-bisphosphate aldolase 1), PFK7 (phosphofructokinase 7), PRK (phosphoribulokinase) and SHM1 (serine hydroxy methyltransferase 1). Moreover, SERAT3;2 has low substrate affinities and low gene expression levels [74] but plays roles in later developmental stages of root formation and metabolism [75]. In addition, fructose-bisphosphate aldolase catalyzes an aldol cleavage of fructose-1,6-bisphosphate to dihydroxyacetonephosphate [76], which plays a significant role in regulating signal transduction [77]. It has been reported that fructose-bisphosphate aldolase (FBA) plays an important physiological role in accelerating cell elongation to promote root elongation in rice [72]. The seed germination and root elongation of AtFBA1 knock-out plants exhibited sensitivity to abscisic acid, indicating that AtFBA1 expression and aldolase activity are important for changing abscisic acid expression in plants [78]. Under NAA treatment, the relative expression level of GOX1 (gene25254), FBA1 (gene16722), and PFK7 (gene23954) first decreased at BL-60 h to BL-66 h, followed by an increase from stage BL-66 h to BL-72 h, and the highest expression of the first two genes was observed at BL-72 h and for the last one at BL-60 h, respectively ( Figure 9b). Phosphofructokinase (PFK) is an important enzyme in central carbon metabolism that promotes carbon flux into the glycolytic pathway [79]. PFK was upregulated only in JM-262. In addition, when plants were grown under nitrogen-free conditions, the genes involved in nitrogen compound metabolism, carbon metabolism, amino acid metabolism, and photosynthesis were downregulated in roots and leaves [80]. Besides, only three genes (GOX1, GOX2, and GOX3) out of five from Arabidopsis were expressed in the elongation zone. GOX1 and GOX2 are the major glycolate oxidase enzymes [81], which catalyze the oxidation of glycolate to glyoxylate [82]. Therefore, the expression suggests that these GOX family members are associated with distinct functions in roots [81]. All the DEGs in carbon metabolism may play a role during adventitious root formation.

Plant Material
R. pseudoacacia-148 is a half-sib, and their seeds were selected based on national improved variety, from Qingshui County, Gansu Province, China. To investigate the kinetics of adventitious root formation, preparation of RNA extraction for high-throughput sequencing was performed, using tissue culture plantlets of R. pseudoacacia-148 as plant material. Sterilized seeds were germinated on Murashige and Skoog medium without hormones. Dark and light treatments were performed according to Masomi et al. [83], with a little modification-the hypocotyl was sliced from the primary root and hypocotyl junction, but also hypocotyls were wounded for easy uptake of hormones. Seeds were given pretreatment of darkness for 4 days continuously followed by 4 days' photoperiod 16/8 h, and controlled samples were kept under 16/8 h photoperiod for 8 days. After completing 8 days, the 1 cm lower part of elongated (treated) and non-elongated (control) hypocotyls was sliced from the lower portion of the seedlings. These parts of treated and control were evaluated for adventitious root response by growing into Murashige and Skoog medium containing hormone (NAA 0.3 mg/L) under 16/8 h photoperiod. Samples for RNA-sequencing were collected at 60 h, 66 h, and 72 h from control and treated samples. Dark pretreated samples were given the names BL-60 h, BL-66 h, and BL-72 h, while control samples were given the names BL-CK60 h, BL-CK66 h, and BL-CK72 h.

Anatomical and Microscopic Analysis
For histological analysis, the same method and growth condition as RNA-Seq were employed for collecting the 1 cm lower portion of the hypocotyl. At least 20 hypocotyls were collected for each time point with three biological replicates. The collected tissues were fixed in formaldehyde acetic acid FAA (50% ethanol, 38% formaldehyde, acetic acid 5 mL/100 mL, and glycerol 5 mL/100 mL) and then treated with the protocol as described by [84]. A lower portion of the hypocotyl (8 µm) was prepared with a rotary microtome (LEICA RM2235) for histological study and then stained with safranin (1%) and fast-green (0.5%). All sections were examined and photographed with LEICA DMI40008 microscope for different tested time points.

RNA Extraction, Deep Sequencing, Functional Annotation, Pathway Enrichment, and qRT-PCR Analysis
Total RNA was extracted from all 12 samples (control and treatment), and RIN value (28S/18S) and RNA purity (OD260/280) were detected using ultraviolet spectrophotometer Nanodrop. The first-strand cDNAs were synthesized by using the Takara Reverse Transcription System (TaKaRa, Tokyo, Japan). The specific PCR primers for the selected candidate gene were designed and primer sequences are listed in Table S4. The qRT-PCR reactions were performed in a Mini Opticon Real-Time PCR System (Bio-Rad, Hercules, CA, USA). SYBR Green Master Mix Reagent (TaKaRa, Shuzo, Otsu, Tokyo, Japan) was used per protocol. The comparative 2-∆∆CTmethod was used to calculate the relative expression levels [47]. The heat map for the gene expression profiles was generated with Mev 4.0 software (http://www.tm4.org/) [48]. Actin was used as a housekeeping gene.
Constructed libraries were tested using the Agilent 2100 bioanalyzer and concentration was quantified using the ABI Step One Plus real-time PCR system. After quality inspection, sequencing was performed by Illumina Hi-Seq sequencer. Raw image data obtained by high-throughput Illumina Hi-Seq-2500 were converted to sequences by CASAVA bas call analysis and stored in FASTQ file format, which contained sequence information. The sequencing data and mapping results are shown in Table 1. All samples of treatment and control were grouped. We  [85]. HISAT combined with String-Tie [86] gives more accurate gene reconstruction and better prediction of expression level than programs like cufflinks. DEG analysis was performed with DESeq2 = (Fold Change) ≥ 2 with adjusted value ≤ 0.05 and RSEM (v1.2.12) to calculate the expression levels of genes and transcripts (http://deweylab. biostat.wisc.edu/rsem). For functional annotation, the data were searched against the NR ( http://www.ncbi.nlm.nih.gov/), Kyoto Encyclopedia of Genes and Genomes (KEGG) (http: //www.genome.jp/kegg/) and GO (http://www.ncbi.nlm.nih.gov/COG/) databases with a typical cut-off value of E-value < 1 × 10 −5 . Blast2GO v2.5 (E-value < 1 × 10 −6 ) software was used to assign Gene Ontology (GO) annotations. After the GO and KEGG annotation, enrichment analysis was carried out using the Phyper function in R software. The p-values and functional divergence ratio (FDR) were evaluated by considering an FDR ≤ 0.01 value as significant enrichment.

Conclusions
The present study describes the important histological events of early root formation and biological pathways occurring during adventitious root formation in R. pseudoacacia-148 under etiolation as a pretreatment. The histological analysis provided the origin and timing of different stages of adventitious root formation. Etiolation increased the number and length of adventitious roots and promoted root formation earlier than the control at 72 h after 8 days of pretreatment. The present study provides novel information for understanding the relationship and involvement of candidate genes (AUX/IAA, SAUR, GH3.3, PRR2, SCL15, MYC2, SUS3, AIR3, and SHM1) in the various important pathways (plant hormone signal pathway, starch and sucrose metabolism, phenylpropanoids, carbon metabolism) that regulate adventitious root formation. The current study provides an overview of the transcriptomic changes that occur during adventitious root formation in R. pseudoacacia. Based on annotation and enrichment analysis, the response of R. pseudoacacia during adventitious root formation has been demonstrated. This information will be helpful for future studies on how etiolation responds during the formation of adventitious roots in R. pseudoacacia and other woody species. The whole-genome sequence of R. pseudoacacia has not yet been released; our study will provide reference data for molecular research. However, many details of the molecular mechanism by which etiolation promotes adventitious root formation in R. pseudoacacia remain unclear. These will be our research focus in future studies. Adventitious root formation is a complex biological process regulated by different pathways. Our work provides a transcriptomic overview of NAA-induced adventitious root formation. However, owing to technical limitations, our study did not offer any post-transcriptional or (post)-translational evidence for adventitious root formation, which is also an important regulatory pathway for adventitious root formation. Whether these pathways are also involved in NAA-induced adventitious root formation still needs further analysis.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/f12060789/s1, Figure S1: The distribution interval of gene expression in RNA-Seq data of each sample, Figure S2: Correlation and principal component analysis, Figure S3: Gene ontology classification of DEGs among control and treatment and among treatment three time points, Table S1: Gene ontology classification of different comparisons, Table S2: Kyoto Encyclopedia of Genes and Genomes (KEGG) of different comparisons, Table S3: DEGs identified in different pathways at three time points with respective control (CK60, 60h, 66h, 72h) with expression values in log2(FPKM), Table S4: Sequence of primers used in the qRT-PCR analysis to validate gene expression.