Yan Lin†
Xiao Han
Huahong Huang
Junhong Zhang- State Key Laboratory of Subtropical Silviculture, School of Forestry & Bio-technology, Zhejiang A&F University, Hangzhou, China
Because of the immobility, plants encounter a series of stresses, such as varied nutrient concentrations in soil, which regulate plant growth, development, and phase transitions. Nitrogen (N) is one of the most limiting factors for plants, which was exemplified by the fact that low nitrogen (LN) has a great adverse effect on plant growth and development. In the present study, we explored the potential role of microRNAs (miRNAs) in response to LN stress in Betula luminifera. We identified 198 miRNAs using sRNA sequencing, including 155 known and 43 novel miRNAs. Among them, 98 known miRNAs and 31 novel miRNAs were differentially expressed after 0.5 h or 24 h of LN stress. Based on degradome data, 122 differential expressed miRNAs (DEmiRNAs) including 102 known miRNAs and 20 novel miRNAs targeted 203 genes, comprising 321 miRNA–target pairs. A big proportion of target genes were transcription factors and functional proteins, and most of the Gene Ontology terms were enriched in biological processes; moreover, one Kyoto Encyclopedia of Genes and Genomes term “ascorbate and aldarate metabolism” was significantly enriched. The expression patterns of six miRNAs and their corresponding target genes under LN stress were monitored. According to the potential function for targets of DEmiRNAs, a proposed regulatory network mediated by miRNA–target pairs under LN stress in B. luminifera was constructed. Taken together, these findings provide useful information to elucidate miRNA functions and establish a framework for exploring N signaling networks mediated by miRNAs in B. luminifera. It may provide new insights into the genetic engineering of the high use efficiency of N in forestry trees.
Introduction
Nitrogen (N) is one of the most indispensable elements for plant growth and development, which is considered the main limiting macronutrient (Perez-Tienda et al., 2014). N participates in a variety of biological processes, such as protein synthesis, nucleic acid synthesis, chlorophyll synthesis, ATP synthesis, and phytohormone synthesis, and acts as a signal (Zhao et al., 2015). Nevertheless, the majority of soil does not satisfy the plant’s needs (Perez-Tienda et al., 2014). Although a large amount of fertilizers are applied, relatively low efficiency for fertilizer is commonly observed, leading to environmental pollution and even harm to human health (Ahmed et al., 2017). Therefore, when a plant encounters N starvation, it is urgent to understand the molecular regulatory mechanism of adaption. Plants often grow in barren soil with limited N, especially forestry trees that are usually planted in infertile soil. Betula luminifera with rapid growth and excellent wood properties, as a pioneer tree, often grows in barren soil with limited N source, which restricts significantly the growth and productivity. It is urgent to elucidate the adaption mechanism of N starvation in this forestry tree species, which might help to make breeding strategies improve tolerance and adaptation to N starvation in B. luminifera and other forestry trees.
Amounting evidence showed that plants have evolved a complex but systematic adaptive mechanism to adapt the N starvation, such as gene transcription and regulation (Liu et al., 2017; Ye et al., 2019), morphological changes (Wang et al., 2019), and photosynthetic characteristics and metabolic pathways (Liu et al., 2019; Zhao et al., 2020). Among them, increasing studies suggest that miRNAs play an important role in response to N starvation, by regulating the expression of target genes, triggering the change of physiological and biochemistry characteristics in many species (Zuluaga and Sonnante, 2019; Li L. J. et al., 2021).
miRNA is a class of 18–24-nt endogenous non-coding RNA, which is processed from the MIRNA gene with stem-loop structure by Dicer-like (DCL) (Jones-Rhoades et al., 2006). Plant miRNAs are usually highly complementary with their target genes, triggering the target mRNA cleavage by RNA-induced silencing complex. However, the translational inhibition is occurred in the target genes, when the seed region of miRNAs harbors relatively low complementary with target mRNA (Brodersen et al., 2008). In 2002, miRNA was first found in plants, which was proved to exert various effects on plants (Llave et al., 2002). Previous studies showed that miRNAs are shown to take part in the regulation of growth and development, as well as acclimatization to adversity stress. For instance, miR396e and miR396f are two important regulators of grain size and plant architecture in rice (Miao et al., 2020). The upregulated miR160, miR167, and miR393 might be participated in the regulation of primary root length and lateral root number under N deficiency (Li L. J. et al., 2021). Cui et al. (2020) showed that miR1885 dynamically regulates both innate immunity and plant growth and responds to viral infection, by targeting both TIR-NBS-LRR and a photosynthesis-related gene for negative regulation through distinct modes of action (Cui et al., 2020). Amounting studies have demonstrated that miRNAs help plants to adapt to drought stress (Bhat et al., 2020), cold stress (Yang et al., 2017), heat stress (Pan et al., 2017), and heavy metal toxicity (Celik and Akdas, 2019). In recent years, several studies showed that miRNA–target modules played important roles in response to nutritive stress. For instance, the miR396e/f-GRF4/6/8 module takes part in the regulation of yield, N assimilation, and utilization in rice, and mir396e/f mutants have higher N content and yield in N-deficient conditions (Zhang et al., 2020). TamiR2275 was induced by N starvation stress, which regulated N acquisition and photosynthetic function when plants are challenged by N deprivation (Qiao et al., 2018). Our previous study demonstrated that a potential regulatory mechanism of the miR169c-NFYA10 module is involved in the low-nitrogen (LN) stress response of B. luminifera (Cheng et al., 2021). However, the expression profiles of miRNAs in response to LN stress are waiting to be elucidated in B. luminifera.
In recent years, the widespread high-throughput sequencing strategy was used for sRNAome (Wen et al., 2022). The ability to detect low-expressed miRNAs contributes to the tremendous progress of miRNAs, especially in non-model plant species (Dos Santos et al., 2019; Li H. et al., 2021). Furthermore, the combination of sRNA sequencing and degradome sequencing promotes the identification of miRNA target genes (Charles Addo-Quaye, 2008). In our previous study, we identified 114 known miRNAs and defined 49 targets for 26 miRNA families in B. luminifera (Zhang et al., 2016), and 44 miRNAs and 71 corresponding target genes were identified in response to heat stress (Pan et al., 2017). Therefore, in the present study, we sequenced and analyzed sRNA libraries and roots and shoots under control, 0.5 h N-starved, and 24 h N-starved treatments, respectively. We tested two hypotheses: 1) The abundance of miRNAs is affected by N starvation, and these miRNAs might have regulatory roles in N starvation adaptation in B. luminifera by regulating the abundance of their corresponding target genes. 2) The expression patterns of N starvation-responsive miRNAs change with stress duration.
Materials and Methods
Plant Cultivation and Low-Nitrogen Treatments
The B. luminifera G49# seedlings were propagated in vitro in 1/2 MS at the temperature of 25 ± 2°C with a 16-h photoperiod. After 2-month growth, seedlings were transferred to a mixed substance consisting of perlite, peat, and clay loam in proportion to 1:1:2, under 120 μmol m−2•s−1 light for 2 weeks. Then, 27 uniformly growing individuals were selected and randomly divided into two groups, namely, CK and LN. Seedlings in the CK group were cultivated using half-strength liquid MS media including 15 mm KNO3, whereas the seedlings in the LN group were cultivated using half-strength MS media containing 0.03 mm KNO3 (the only source of N), as described previously (Cheng et al., 2021). Then, the leaves and roots were sampled individually after 0.5 and 24 h after the CK and LN treatments, at 11:00 a.m., respectively. Three randomly selected seedlings were pooled together as one replicate, with three replicates per group, and samples were flash frozen in liquid N and stored at −80°C until RNA extraction.
The Expression Patterns of Low-Nitrogen-Responsive Marker Genes
Total RNA was extracted using TRIZOL (Invitrogen, Carlsbad, CA, United States) and was then reverse transcribed to cDNA using PrimeScript RT Reagent Kit (TaKaRa, Dalian, China). To verify the effectiveness of LN stress treatment, two marker genes NRT2.5 and NRT2.7 were quantified via qRT-PCR (Supplementary Table S1). The large ribosomal subunit 39 (RPL39) was used as a reference gene to normalize the expression levels of marker genes (Cheng et al., 2021).
sRNA Library Construction, Sequencing, and Bioinformatics
The construction and sequencing of the sRNA library were accomplished by LC Bio (Hangzhou, China) with Illumina Hiseq 2500. After the quality assessment of sequencing data, the retained data were processed using ACGT101-miR (LC Sciences, Houston, Texas, United States), including removing adaptor dimers, junk sequences, and low complexity sequences and filtering out mRNA, Rfam, and Repbase. Then, the remained sequences were blasted against miRBase21.0 and genome, and the hairpin structure of predicted miRNA precursor sequences was validated using RNAfold software (http://rna.tbi.univie.ac. at/cgi-bin/RNAfold. Cgi).
Analysis of Differentially Expressed miRNAs
The reads of miRNAs from the control, 0.5 h, and 24 h of LN treatment were normalized to tags per million. The combination of log2 ratio and Fisher’s exact test was used to identify the differential expressed miRNAs (DEmiRNAs). Only when the absolute value of log2 ratio ≥ 1 and p ≤ 0.05, the miRNA was significantly differentially expressed; however, in other cases, there was no significant difference in expression (Xie et al., 2015).
Target Prediction of DEmiRNAs and Functional Classification
Target prediction of DEmiRNAs both known and novel miRNAs was performed using Targetfinder (http://targetfinder.org/), and the mismatch score was less than four. Based on our published degradome sequencing data of B. luminifera (GEO: GSE80074), we used CleaveLand 3.0.1 pipeline to identify and classify the target genes of DEmiRNAs. To further elucidates the potential regulation of target genes, GO (Gene Ontology) analysis was conducted using agriGO2.0 (http://bioinfo.cau.edu.cn/agriGO/). In addition, Kyoto Encyclopedia of Genes and Genomes (KEGG) annotation was conducted using a database (http://www.genome.jp/kegg/kegg1.html).
The Expression Patterns of miRNAs and Their Target Genes
Total RNA was extracted as described above and was reversely transcribed with Mir-X miRNA First-Strand Synthesis Kit (TaKaRa, Dalian, China). qRT-PCR was conducted using SYBR Premix EX Taq Kit (TaKaRa). Gene U6 was used as a reference gene to normalize the expression levels of miRNAs. The RNA was reversely transcribed to cDNA using PrimeScript RT Reagent Kit (TaKaRa), which was used to quantify the target genes. The RPL39 was used as a reference gene, and the relative expression levels of genes were calculated via the 2−ΔΔCt method. The qRT-PCR primers are shown in Supplementary Table S1.
Results
The Expression Patterns of Low-Nitrogen-Responsive Marker Genes in B. luminifera
To confirm the physiological status of the B. luminifera seedlings (i.e., to confirm whether the plants were subjected to LN stress), the expression levels of the marker genes BlNRT2.5 and BlNRT2.7 were analyzed. BlNRT2.5 and BlNRT2.7 as HATS, responsible for nitrate transport and absorption under LN, are indicators of whether the plants significantly underwent LN stress, as reported previously (Wang et al., 2012). In both root and shoot, the expression level of BlNRT2.5 was increased significantly after 0.5 h of LN treatment and reached the peak value after 24 h of LN treatment (Figure 1A). For the BlNRT2.7, the expression level was induced 10-fold after 0.5 h of LN treatment and remained high level after 24 h of N starvation in the root, whereas the expression level was unchanged until 24 h of LN treatment in the shoot (Figure 1B).
FIGURE 1. The expression levels of marker genes under low nitrogen (LN) treatment. (A) The expression level of BlNRT2.5; (B) The expression level of BlNRT2.7 under LN treatment. The reference gene was BlRPL39, each sample had three replicates, and error bars were standard deviation.
sRNA Profiles in Response to Low-Nitrogen Stress
To investigate the response of miRNAs to N starvation, six mixed sRNA libraries (each library was a pool of RNAs from three biological samples) were generated from roots and shoots under LN stress for 0, 0.5, and 24 h, respectively, and then were sequenced via Solexa high-throughput sequencing technology. After the remove of adaptor and junk sequences, a total of 66,031,208 high-quality reads were obtained and then mapped to known RNAs, resulting in 2,305,703 (461,778), 5,012,446 (1,028,698), 2,663,870 (514,935), 3,475,390 (785,402), 4,953,286 (1,468,323), and 1,322,911 (319.385) redundant (unique) valid reads, which were used to identify miRNAs (Supplementary Table S2).
Length Distribution of sRNAs Under CK and Low-Nitrogen Conditions
Most of the sRNA sequences ranged from 18 to 24 nt, representing 97.87% (96.58%) of redundant (unique) valid reads, which was consistent with the characteristic of the DCL1 cleavage product. In the root, 24-nt sRNAs had the highest proportion of sRNAs in the unique valid reads (Figure 2A). Under LN stress, the length distribution patterns of sRNAs in the root were altered, which was exemplified by the fact that the ratio of 24-nt sRNAs increased from 26.31% (0 h) to 40.42% (0.5 h) and then recovered to 23.24% (24 h). A similar phenomenon was observed in the shoot, which was exemplified by the fact that 24-nt sRNAs accounted for 32.26% in 0 h and increased to 45.49% after 0.5 h of LN stress and then reduced to 12.20% (Figure 2B). However, the 18-nt sRNAs had the highest proportion, accounting for 21.66% in the shoot after 24-nt LN stress.
FIGURE 2. The length distribution of sRNAs under low-nitrogen stress. (A). Unique valid reads in root; (B). Unique valid reads in shoot; (C). Redundant valid reads in root; and (D). Redundant valid reads in shoot.
For redundant valid reads, 21-nt sRNAs were the most abundant class of sRNAs in the root, and the ratio of 24-nt sRNAs increased after 0.5 h of LN treatment (Figure 2C). In the shoot, 21-nt sRNAs were also the most abundant class, accounting for 26.76 and 25.69% in the 0 and 0.5 h of LN treatment, respectively. However, the ratio of 18-nt sRNAs reached the peak level after 24 h of LN stress (Figure 2D).
Identification of Known and Novel miRNAs in B. luminifera
A total of 155 known miRNAs from 55 families were identified, and 43 novel miRNAs from 32 families were discovered, resulting in 198 miRNAs from 87 families. Among them, 113 known miRNAs were identified with according precursors, whereas other 42 known miRNAs lacked the corresponding precursors, and miR166 was the most abundant miRNA family (Supplementary Table S3).
miRNA Expression Response to Low-Nitrogen Stress
In the shoot, 76 known miRNAs from 42 families and 22 novel miRNAs from 20 families were recognized as DEmiRNAs after 0.5 h and/or 24 h N-starved treatment, whereas 98 known miRNAs from 46 families and 31 novel miRNAs from 24 families were differentially expressed in the root (Supplementary Table S4). To be precise, 54 and 69 DEmiRNAs were identified in the shoot, whereas 105 and 72 DEmiRNAs were detected in the root after 0.5 and 24 h N-starved treatment, respectively (Figure 3A). Among them, 22 and 44 DEmiRNAs were specifically expressed in the shoot after 0.5 and 24 h treatment, respectively, whereas 48 and 38 DEmiRNAs were specifically identified in the root (Figure 3B, Supplementary Table S4).
FIGURE 3. The identification and validation of DEmiRNAs after 0.5 and 24 h N-starved treatment. (A) Total number of DEmiRNA in shoot and root; (B) The specific DEmiRNA in shoot and root; (C) The validation of DEmiRNAs after 0.5 h LN treatment; and (D) The validation of DEmiRNAs after 24 h of low nitrogen treatment.
To validate the results of the sRNA sequencing, eight DEmiRNAs were selected randomly for qRT-PCR analysis. Although the fold-changes of DEmiRNAs between qPCR and sRNA sequencing were different, similar expression tendencies were observed (Figures 3C,D).
Target Genes Regulated by DEmiRNA
A total of 167 miRNAs were predicted to target 836 genes, involving 1,094 miRNA–target pairs. Based on degradome data, 122 miRNAs including 102 known miRNAs and 20 novel miRNAs targeted 203 genes, comprising 321 miRNA–target pairs (Supplementary Table S5). These targets were divided into five categories: categories 0–4 (Figure 4). For known miRNAs, most of the target genes were annotated as transcription factors. For instance, SPL13A (Squamosa promoter-binding protein-like, SPL) was targeted by three miR156, NAC021 (NAM, ATAF, CUC, NAC) by miR164a, NFYA2/7 (Nuclear Transcription Factor Y Subunit A, NFYA) by five miR169, RAP2–7 (Ethylene-responsive transcription factor, RAP) by miR172a/b, and GRF4/6 (Growth-Regulating Factor, GRF) by miR396a/b. A big proportion of target genes for novel miRNAs were proteins or enzymes. For example, NSP2 (Nodulation-Signaling Pathway 2 protein) was targeted by blu-miR1-3p, CMTA2 (CalModulin-binding Transcription Activator 2-like) by blu-miR5-3p, OPT3 (OligoPeptide Transporter 3) by blu-miR7-3p, and FBX (F-box/kelch-repeat protein) by blu-miR8-3p.
FIGURE 4. T-plots of miRNA targets based on degradome sequencing. Asterisks were used to mark cleavage sites of target genes.
Gene Ontology and Kyoto Encyclopedia of Genes and Genomes Pathway Analyses of Target Genes
GO was used to predict the function of target genes. A total of 665 target genes were annotated with 1,115 GO terms (Supplementary Table S6). Among them, a large proportion (54.08%) of target genes were classified into biological processes, one-third of genes for molecular function (327 terms), and 169 items belonging to cellular components. Of these, 13 GO terms were extremely significantly enriched (p < 0.01), including leaf development, auxin-activated signaling pathway, cell differentiation, transcription factor activity, DNA binding, and nucleus (Figure 5A). One KEGG term “ascorbate and aldarate metabolism” (ko00053) was significantly enriched (p < 0.05) (Figure 5B).
FIGURE 5. Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses of genes targeted by DEmiRNAs. (A) GO enrichment analysis; (B) KEGG pathway enrichment analysis of genes targeted by DEmiRNAs.
Expression Patterns of miRNAs and Their Target Genes Under Low-Nitrogen Stress
The expression patterns of six miRNAs and their target genes were monitored after 0.5 and 24 h of LN treatment (Figure 6). The expression levels of miR156a and miR166b-5p were initially decreased after 0.5 h of N starvation but increased after 24 h of LN treatment, whereas the expression patterns of blu-miR1-3p and blu-miR8-3p showed an opposite trend. However, miR164a and miR167a exhibited continuous decreasing expression patterns under LN treatment. To determine the relationships between miRNAs and target genes, we examined the expression levels of target genes under LN treatment using qRT-PCR (Figure 6). The overall expression levels of target genes were unchanged or decreased slightly after 0.5 h of LN stress but were induced significantly by 24 h of LN stress.
FIGURE 6. The expression levels of six miRNAs and their predicted target genes under low-nitrogen stress. The left y axis represents the expression level of miRNA, and the right Y axis represents the expression level of the predicted target gene. All expression was normalized by BlRPL39, and three replicates were conducted. (A) miR156a-SPL13; (B) miR164a-NAC100; (C) miR166b-5p-HHP1; (D) miR167a-TuDP; (E) miR1-3p-NSP2; (F) miR8-3p-FBX.
Discussion
N is one of the most important restrictive elements in plant growth and development, and nitrate is the major form for most land plants (Glass et al., 2002). The limited nitrate in soil restricts plant growth; thus, it is of great importance to understand how plants respond to nitrate deficiency. Amounting studies showed that miRNAs played important roles in response to N starvation, such as Arabidopsis (Zhao et al., 2011), wheat (Liu et al., 2021), potato (Tiwari et al., 2020), and rice (Yu et al., 2018).
Length Distribution of sRNAs Varied Under Low-Nitrogen Stress
Several plant species, such as Arabidopsis, cotton, and Chinese white poplar, the most redundant sRNAs, were 21-nt sized. In the present study, 21-nt sRNAs were the most abundant class of sRNAs in both root and shoot of B. luminifera (Figures 2C,D). A similar phenomenon was observed in our previous study (Zhang et al., 2016), whereas the 24-nt class of sRNAs was the most dominant in heat-stressed B. luminifera (Pan et al., 2017). It might be attributed to developmental- and/or stress-specific differences. In the Chrysanthemum nankingense, 21- and 24-nt sRNAs were the most abundant classes, and the proportion of 24-nt sRNA in the root was higher than that of 21-nt sRNA, but it was the opposite in leaves (Song et al., 2015). In potatoes, 24-nt sRNAs accounted for the largest proportion in all classes, followed by 21-nt sRNAs (Tiwari et al., 2020). In maize, the proportion of 21-nt sRNA was the highest under all treatments of shoots, which was the same in the root under control treatment, but the sRNA of the root was mainly distributed in 24 and 22 nt under LN treatment (Zhao et al., 2012). Moreover, the ratio of 24-nt sRNAs increased in the root after 0.5 h of LN treatment (Figure 2C), and the ratio of 18-nt sRNAs reached the peak level after 24 h of LN stress in the shoot (Figure 2D). It was consistent with the major peak at 18 nt in B. luminifera roots under LP stress (Zhang et al., 2021). A previous study showed that the components of sRNAs with the size of 18–19 nt originated from tRNAs, and the accumulation of tRNA fragments in Pi-starved (Hsieh et al., 2009) and nitrate-starved Arabidopsis (Liang et al., 2012). Moreover, tRNA halves were detected in the phloem sap of Cucurbita maxima, which inhibited translational activity in vitro, indicating that these tRNA halves coordinate the metabolic status between source and sink tissues, acting as long-distance signals (Zhang et al., 2009). Loss-Morais et al. (2013) showed that tRNA fragments can be incorporated into argonaute (AGO) complexes, thus regulating gene expression post-transcriptionally, acting as miRNAs (Loss-Morais et al., 2013). Therefore, these 18-nt sRNAs induced by LN stress might play important regulatory roles in response to N starvation in B. luminifera, which need to be further studied.
miRNA Expression Response to Nitrogen Starvation
Previous studies have identified a series of LN-responsive miRNAs. For instance, 102 miRNAs from 42 families in maize (Zhao et al., 2012) and 404 and 628 putative miRNAs were identified in potato roots and shoots, respectively (Tiwari et al., 2020). Although three studies on miRNA of B. luminifera have been conducted (Zhang et al., 2016; Pan et al., 2017; Zhang et al., 2021), the LN-responsive miRNAs were waiting to be identified. In the present study, a total of 198 miRNAs from 87 families were identified, of which 113 known miRNAs were identified with precursors. Thus, it enriches the miRNA information of B. luminifera.
miR166 was the most abundant miRNA family, and miR166a accounted for the largest proportion, which was consistent with the result of Chinese cabbage (Ahmed et al., 2019), whereas miR166a-3p was the most abundant miRNA in our previous study (Pan et al., 2017). The inconsistency might relate to the alternative cleavage mechanism depending on the plant species, growth stage, and/or stress. For example, the miR535 family had the highest abundance, and miR156, miR166, and miR168 families also had high expression levels in bananas (Zhu et al., 2019). Several miRNAs were specifically expressed in some samples (Supplementary Table S4), such as 156 d-3p in a 0.5 h LN-treated shoot. A similar phenomenon was observed in potatoes, where five conserved miRNAs were specifically expressed in roots and nine in shoots. It indicated that the accumulation of some miRNAs is dependent on plant species, tissues, and/or stresses. In the present study, more downregulated DEmiRNAs were observed under LN stress (Figure 3). A similar phenomenon was observed in sorghum and durum wheat, which was exemplified by the fact that the number of downregulated DEmiRNAs was more than that of upregulated DEmiRNAs under LN stress (Liu et al., 2021; Zhu et al., 2021). The GO analysis of DEmiRNA targets showed that leaf development, auxin-activated signaling pathway, cell differentiation, and transcription factor activity were enriched (Figure 5A), and the KEGG term “ascorbate and aldarate metabolism” (ko00053) was significantly enriched (Figure 5B). It suggested that miRNAs play important roles in leaf development, hormone signaling pathway, cell differentiation, ascorbate and aldarate metabolism, and transcriptional regulations in response to LN stress in B. luminifera.
Expression Patterns of miRNAs and Targets in B. luminifera Under Low-Nitrogen Stress
Under LN stress, plants change their phenotypes and transcriptional activity (Curci et al., 2017). Amounting studies showed that miRNAs play essential roles in response to LN stress by cleaving target mRNAs or inhibiting translation of mRNAs. In our study, a series of miRNAs and validated target genes based on degradome sequencing might take part in response to environment starvation; thus, a putative model on miRNA–target interactions under LN stress was proposed (Figure 7). The biogenesis of miRNAs was processed by DCL1, which was regulated by miR162 (Xie et al., 2003). Then, the mature miRNA entered the AGO to regulate target genes. AGO1 was proved to be targeted by miR168, which was verified in this study. In addition, AGO2 was targeted by miR403, as described in our previous study (Pan et al., 2017). Therefore, the miRNA biogenesis and function were feedback regulated by miR162, miR168, and miR403, when the B. luminifera encounter the LN stress, thus helping plants to adapt to N starvation condition. In the present study, most of the target genes for known miRNAs were TFs, such as ARFs, HD-Zips, NFYAs, and SPLs (Supplementary Table S5), whereas some targets were protein-coding genes, such as Disease resistance protein (RGA), Threonine Synthase (TS1) and Laccase (LACs). These miRNA–target modules played important roles in transcriptional regulation, auxin response, energy metabolism, biosynthesis, and signal transduction.
FIGURE 7. Proposed model of the miRNA-mediated regulatory network in response to low nitrogen in B luminifera. The miRNAs negatively regulated their target genes, and these targets are involved in leaf, stem, root development, signal transduction, stress resistance, cell composition and nutrient acquisition, and substance synthesis.
The aging pathway in flowering regulation is controlled mainly by miR156. SPLs were targeted by miR156, and overexpression of SPL7 and SPL8 promotes flowering, and vice versa (Gou et al., 2019). miR156 was induced by LN stress, and the miR156-SPL3/SPL12 module directly activates MADS50 in the node to regulate crown root development in rice (Shao et al., 2019). In general, miR156 is induced by LN in Chrysanthemum nankingense (Song et al., 2015). In our study, miR156s could be divided into two groups based on expression patterns, which was exemplified by the fact that one group was downregulated continuously during the treatment period, whereas the other group was upregulated after 0.5 h treatment and then downregulated after 24 h, which was consistent with the miR156 expression in poplar (Bao et al., 2019). miR156a was repressed under LN stress, and SPL13 was downregulated after 0.5 h but was upregulated after 24 h of LN stress, which might promote the transition from the vegetative to the reproductive stage under LN stress in B. luminifera.
A previous study showed that miR164 played a negative role in adventitious root development by downregulating the expression of NAC1 (Li J. et al., 2021). A NAC gene EjNACLl47 had a positive effect on cell expansion and organ enlargement in Eriobotrya japonica (Chen Q. et al., 2021). Besides, NAC is related to leaf senescence; for instance, GmNAC81 was upregulated during the leaf senescence, and the overexpression of GmNAC81 accelerated leaf senescence by accumulating H2O2 (Ferreira et al., 2020). Leaf senescence is a process of recycling nutrients for other organs (Lim et al., 2007). In the present study, the expression level of miR164 was decreased under LN stress, whereas the target NAC100 was inhibited after 0.5 h and then was induced significantly after 24 h of LN stress (Figure 6). It suggested that miR164 might play an important role in leaf senescence mainly by targeting NAC100, to accelerate the N recycling under LN conditions in B. luminifera.
The Class III homeodomain Leu zipper (HD-Zip III) gene family plays essential roles in plant growth and development. For example, overexpression of HOX32 (an HD-Zip III family member) resulted in curled leaf, reduced leaf angle, semi-dwarf phenotype, and reduced photosynthetic capacity. OsHox32 may regulate plant architecture and leaf development by regulating YABBY genes in rice (Li et al., 2016). HD-Zip III including Hox32 targeted by miR166 and the OsmiR166b-OsHox32 module plays important roles in plant growth and development and plant architecture by regulating the cell wall-related gene expression (Chen H. et al., 2021). However, miR166b-5p was validated to target HHP1 (Heptahelical protein) in B. luminifera. Previous studies showed that HHP participated in the signal transduction pathway by cooperating with MYB to transduce hormone and environmental signals to cope with cold stress (Lee and Seo, 2015). Therefore, we speculated that the inhibited miR166b-5p to ensure the accumulated target HHP1 thus took part in the signal transduction pathway to cope with LN stress in B. luminifera.
blu-miR1-3p was validated to target NSP (nodulation-signaling pathway), an GRAS (GAI, RGA, SCR)-type transcription factor, which was a constituent of the Myc signaling pathway (Delaux et al., 2013). NSP is generally related to nodule symbiosis (Singh and Verma, 2021) and also participated in other stress responses, such as drought stress (Yao et al., 2019). Liu et al. (2011) showed that NSP1 and NSP2 are indispensable for strigolactone (SL) biosynthesis both in the legume and rice, which was exemplified by the fact that NSP2 is essential for the conversion of orobanchol into didehydro-orobanchol, dependent on DWARF27, a gene essential for SL biosynthesis (Liu et al., 2011). Therefore, we speculated that novel miRNA blu-miR1-3p might participate in SL biosynthesis by targeting the NSP gene, when B. luminifera encounters N starvation. However, the response mechanism of blu-miR1-3p under LN stress still needs to be further studied.
Conclusion
Integrated analysis of sRNA-Seq and degradome-Seq was conducted to elucidate the molecular mechanisms in response to LN stress in B. luminifera. sRNA sequencing identified 98 known miRNAs, and 31 novel miRNAs were differentially expressed after 0.5 or 24 h of LN stress. A total of 122 DEmiRNAs targeted 203 genes, comprising 321 miRNA–target pairs. Thus, a miRNA-mediated network was proposed. Taken together, these findings provide useful information to elucidate miRNA functions and establish a framework for exploring N signaling networks mediated by miRNAs in B. luminifera and other woody plants.
Data Availability Statement
The data presented in the study are deposited in the CNGB Sequence Archive (CNSA) of China National GeneBank DataBase (CNGBdb) repository, accession number CNP0003144.
Author Contributions
JZ and ZT conceived and designed the experiments. YL, JZ, and SC analyzed the sRNA sequencing data. HH, ZT, and JZ contributed reagents and materials. YL, XX, and XH performed the qRT-PCR experiment. YL, JZ, and SC wrote and revised the manuscript. All authors read and approved the final manuscript.
Funding
This work was financially supported by the National Natural Science Foundation of China (31971677), Zhejiang Science and Technology Major Program on Agricultural New Variety Breeding (2021C02070-10), and National Key Research and Development Program of China (2021YFD2200304-2).
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s Note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2022.957505/full#supplementary-material
References
Addo-Quaye, C., Eshoo, T. W., and Axtell, D. P. M. J. (2008). Endogenous siRNA and miRNA Targets Identified by Sequencing of the Arabidopsis Degradome. Curr. Biol. 18, 758–762. doi:10.1016/j.cub.2008.04.042
Ahmed, M., Rauf, M., Mukhtar, Z., and Saeed, N. A. (2017). Excessive Use of Nitrogenous Fertilizers: an Unawareness Causing Serious Threats to Environment and Human Health. Environ. Sci. Pollut. Res. 24, 26983–26987. doi:10.1007/s11356-017-0589-7
Ahmed, W., Xia, Y., Zhang, H., Li, R., Bai, G., Siddique, K. H. M., et al. (2019). Identification of Conserved and Novel miRNAs Responsive to Heat Stress in Flowering Chinese Cabbage Using High-Throughput Sequencing. Sci. Rep. 9, 14922. doi:10.1038/s41598-019-51443-y
Bao, H., Sun, F. S., Xu, Q. H., and Wang, Y. W. (2019). Differential Expression of 10 miRNAs in Poplar under Low Nitrogen Stress and the Identification and Analysis of Target Genes. Mol. Plant Breed. 17, 771. doi:10.13271/j.mpb.017.000771
Bhat, K. V., Mondal, T. K., Gaikwad, A. B., Kole, P. R., Chandel, G., and Mohapatra, T. (2020). Genome-wide Identification of Drought-Responsive miRNAs in Grass Pea (Lathyrus Sativus L.). Plant gene. 21, 100210. doi:10.1016/j.plgene.2019.100210
Brodersen, P., Sakvarelidze-Achard, L., Bruun-Rasmussen, M., Dunoyer, P., Yamamoto, Y. Y., Sieburth, L., et al. (2008). Widespread Translational Inhibition by Plant miRNAs and siRNAs. Science 320, 1185–1190. doi:10.1126/science.1159151
Çelik, Ö., and Akdaş, E. Y. (2019). Tissue-specific Transcriptional Regulation of Seven Heavy Metal Stress-Responsive miRNAs and Their Putative Targets in Nickel Indicator castor Bean (R. Communis L.) Plants. Ecotoxicol. Environ. Saf. 170, 682–690. doi:10.1016/j.ecoenv.2018.12.006
Chen, H., Fang, R., Deng, R., and Li, J. (2021a). The OsmiRNA166b-OsHox32 Pair Regulates Mechanical Strength of Rice Plants by Modulating Cell Wall Biosynthesis. Plant Biotechnol. J. 19, 1468–1480. doi:10.1111/pbi.13565
Chen, Q., Jing, D. L., Wang, S. M., Xu, F., Bao, C. Y., Luo, M., et al. (2021b). The Putative Role of the NAC Transcription Factor EjNACL47 in Cell Enlargement of Loquat (Eriobotrya Japonica Lindl.). Horticulturae 7, 90323. doi:10.3390/horticulturae7090323
Cheng, L. L., Wu, F. M., Lin, Y., Han, X., Xu, X. S., Zhang, Y. T., et al. (2021). A miR169c-NFYA10 Module Confers Tolerance to Low-Nitrogen Stress to Betula Luminifera. Industrial Crops Prod. 172, 113988. doi:10.1016/j.indcrop.2021.113988
Cui, C., Wang, J.-J., Zhao, J.-H., Fang, Y.-Y., He, X.-F., Guo, H.-S., et al. (2020). A Brassica miRNA Regulates Plant Growth and Immunity through Distinct Modes of Action. Mol. Plant 13, 231–245. doi:10.1016/j.molp.2019.11.010
Curci, P. L., Aiese Cigliano, R., Zuluaga, D. L., Janni, M., Sanseverino, W., and Sonnante, G. (2017). Transcriptomic Response of Durum Wheat to Nitrogen Starvation. Sci. Rep. 7, 1176. doi:10.1038/s41598-017-01377-0
Delaux, P. M., Bécard, G., and Combier, J. P. (2013). NSP1 Is a Component of the Myc Signaling Pathway. New Phytol. 199, 59–65. doi:10.1111/nph.12340
Dos Santos, T. B., Soares, J. D. M., Lima, J. E., Silva, J. C., Ivamoto, S. T., Baba, V. Y., et al. (2019). An Integrated Analysis of mRNA and sRNA Transcriptional Profiles in Coffea Arabica L. Roots: Insights on Nitrogen Starvation Responses. Funct. Integr. Genomics 19, 151–169. doi:10.1007/s10142-018-0634-8
Ferreira, D. O., Fraga, O. T., Pimenta, M. R., Caetano, H. D. N., Machado, J. P. B., Carpinetti, P. A., et al. (2020). GmNAC81 Inversely Modulates Leaf Senescence and Drought Tolerance. Front. Genet. 11, 601876. doi:10.3389/fgene.2020.601876
Glass, A. D. M., Britto, D. T., Kaiser, B. N., Kinghorn, J. R., Kronzucker, H. J., Kumar, A., et al. (2002). The Regulation of Nitrate and Ammonium Transport Systems in Plants. J. Exp. Bot. 53, 855–864. doi:10.1093/jexbot/53.370.855
Gou, J., Tang, C., Chen, N., Wang, H., Debnath, S., Sun, L., et al. (2019). SPL 7 and SPL 8 Represent a Novel Flowering Regulation Mechanism in Switchgrass. New Phytol. 222, 1610–1623. doi:10.1111/nph.15712
Hsieh, L.-C., Lin, S.-I., Shih, A. C.-C., Chen, J.-W., Lin, W.-Y., Tseng, C.-Y., et al. (2009). Uncovering Small RNA-Mediated Responses to Phosphate Deficiency in Arabidopsis by Deep Sequencing. Plant Physiol. 151, 2120–2132. doi:10.1104/pp.109.147280
Jones-Rhoades, M. W., Bartel, D. P., and Bartel, B. (2006). MicroRNAs and their Regulatory Roles in Plants. Annu. Rev. Plant Biol. 57, 19–53. doi:10.1146/annurev.arplant.57.032905.105218
Lee, H. G., and Seo, P. J. (2015). The MYB 96- HHP Module Integrates Cold and Abscisic Acid Signaling to Activate the CBF - COR Pathway in Arabidopsis. Plant J. 82, 962–977. doi:10.1111/tpj.12866
Li, H., Meng, H., Sun, X., Deng, J., Shi, T., Zhu, L., et al. (2021a). Integrated microRNA and Transcriptome Profiling Reveal Key miRNA-mRNA Interaction Pairs Associated with Seed Development in Tartary Buckwheat (Fagopyrum Tataricum). BMC Plant Biol. 21, 132. doi:10.1186/s12870-021-02914-w
Li, J., Zhang, H., Zhu, J., Shen, Y., Zeng, N., Liu, S., et al. (2021b). Role of miR164 in the Growth of Wheat New Adventitious Roots Exposed to Phenanthrene. Environ. Pollut. 284, 117204. doi:10.1016/j.envpol.2021.117204
Li, L. J., Li, Q., Davis, K. E., Patterson, C., Oo, S., Liu, W. Y., et al. (2021c). Response of Root Growth and Development to Nitrogen and Potassium Deficiency as Well as microRNA-Mediated Mechanism in Peanut (Arachis hypogaea L.). Front. Plant Sci. 12, 695234. doi:10.3389/fpls.2021.695234
Li, Y.-y., Shen, A., Xiong, W., Sun, Q.-l., Luo, Q., Song, T., et al. (2016). Overexpression of OsHox32 Results in Pleiotropic Effects on Plant Type Architecture and Leaf Development in Rice. Rice 9, 46. doi:10.1186/s12284-016-0118-1
Liang, G., He, H., and Yu, D. (2012). Identification of Nitrogen Starvation-Responsive microRNAs in Arabidopsis thaliana. PLoS One 7, e48951. doi:10.1371/journal.pone.0048951
Lim, P. O., Kim, H. J., and Gil Nam, H. (2007). Leaf Senescence. Annu. Rev. Plant Biol. 58, 115–136. doi:10.1146/annurev.arplant.57.032905.105316
Liu, H., Able, A. J., and Able, J. A. (2021). Nitrogen Starvation-Responsive MicroRNAs Are Affected by Transgenerational Stress in Durum Wheat Seedlings. Plants (Basel) 10, 50826. doi:10.3390/plants10050826
Liu, W., Kohlen, W., Lillo, A., Op Den Camp, R., Ivanov, S., Hartog, M., et al. (2011). Strigolactone Biosynthesis in Medicago Truncatula and Rice Requires the Symbiotic GRAS-type Transcription Factors NSP1 and NSP2. Plant Cell. 23, 3853–3865. doi:10.1105/tpc.111.089771
Liu, W., Sun, Q., Wang, K., Du, Q., and Li, W. X. (2017). Nitrogen Limitation Adaptation ( NLA ) Is Involved in Source-to-sink Remobilization of Nitrate by Mediating the Degradation of NRT 1.7 in Arabidopsis. New Phytol. 214, 734–744. doi:10.1111/nph.14396
Liu, Y., Li, M., Xu, J., Liu, X., Wang, S., and Shi, L. (2019). Physiological and Metabolomics Analyses of Young and Old Leaves from Wild and Cultivated Soybean Seedlings under Low-Nitrogen Conditions. BMC Plant Biol. 19, 389. doi:10.1186/s12870-019-2005-6
Llave, C., Xie, Z., Kasschau, K. D., and Carrington, J. C. (2002). Cleavage of Scarecrow-like mRNA Targets Directed by a Class of Arabidopsis miRNA. Science 297, 2053–2056. doi:10.1126/science.1076311
Loss-Morais, G., Waterhouse, P. M., and Margis, R. (2013). Description of Plant tRNA-Derived RNA Fragments (tRFs) Associated with Argonaute and Identification of Their Putative Targets. Biol. Direct 8, 6. doi:10.1186/1745-6150-8-6
Miao, C., Wang, D., He, R., Liu, S., and Zhu, J. K. (2020). Mutations in MIR 396e and MIR 396f Increase Grain Size and Modulate Shoot Architecture in Rice. Plant Biotechnol. J. 18, 491–501. doi:10.1111/pbi.13214
Pan, Y., Niu, M., Liang, J., Lin, E., Tong, Z., and Zhang, J. (2017). Identification of Heat-Responsive miRNAs to Reveal the miRNA-Mediated Regulatory Network of Heat Stress Response in Betula Luminifera. Trees 31, 1635–1652. doi:10.1007/s00468-017-1575-x
Pérez-Tienda, J., Corrêa, A., Azcón-Aguilar, C., and Ferrol, N. (2014). Transcriptional Regulation of Host NH4+ Transporters and GS/GOGAT Pathway in Arbuscular Mycorrhizal Rice Roots. Plant Physiology Biochem. 75, 1–8. doi:10.1016/j.plaphy.2013.11.029
Qiao, Q., Wang, X., Yang, M., Zhao, Y., Gu, J., and Xiao, K. (2018). Wheat miRNA Member TaMIR2275 Involves Plant Nitrogen Starvation Adaptation via Enhancement of the N Acquisition-Associated Process. Acta Physiol. Plant 40, 183. doi:10.1007/s11738-018-2758-9
Shao, Y., Zhou, H.-Z., Wu, Y., Zhang, H., Lin, J., Jiang, X., et al. (2019). OsSPL3, an SBP-Domain Protein, Regulates Crown Root Development in Rice. Plant Cell. 31, 1257–1275. doi:10.1105/tpc.19.00038
Singh, J., and Verma, P. K. (2021). NSP1 Allies with GSK3 to Inhibit Nodule Symbiosis. Trends Plant Sci. 26, 999–1001. doi:10.1016/j.tplants.2021.07.001
Song, A., Wang, L., Chen, S., Jiang, J., Guan, Z., Li, P., et al. (2015). Identification of Nitrogen Starvation-Responsive microRNAs in Chrysanthemum Nankingense. Plant Physiology Biochem. 91, 41–48. doi:10.1016/j.plaphy.2015.04.003
Tiwari, J. K., Buckseth, T., Zinta, R., Saraswati, A., Singh, R. K., Rawat, S., et al. (2020). Genome-wide Identification and Characterization of microRNAs by Small RNA Sequencing for Low Nitrogen Stress in Potato. PLoS One 15, e0233076. doi:10.1371/journal.pone.0233076
Wang, J., Song, K., Sun, L., Qin, Q., Sun, Y., Pan, J., et al. (2019). Morphological and Transcriptome Analysis of Wheat Seedlings Response to Low Nitrogen Stress. Plants (Basel) 8, 40098. doi:10.3390/plants8040098
Wang, Y.-Y., Hsu, P.-K., and Tsay, Y.-F. (2012). Uptake, Allocation and Signaling of Nitrate. Trends Plant Sci. 17, 458–467. doi:10.1016/j.tplants.2012.04.006
Wen, Z., Hong, Y., Qiu, Z. L., Yang, K., Hou, Q. D., Qiao, G., et al. (2022). Identification of miRNAs Mediating Shoot Growth of Grafted Sweet Cherry through Small RNA and Degradome Sequencing. Sci. Hortic. 291, 110557. doi:10.1016/j.scienta.2021.110557
Xie, F., Wang, Q., Sun, R., and Zhang, B. (2015). Deep Sequencing Reveals Important Roles of microRNAs in Response to Drought and Salinity Stress in Cotton. J. Exp. Bot. 66 (3), 789–804. doi:10.1093/jxb/eru437
Xie, Z., Kasschau, K. D., and Carrington, J. C. (2003). Negative Feedback Regulation of Dicer-Like1 in Arabidopsis by microRNA-Guided mRNA Degradation. Curr. Biol. 13, 784–789. doi:10.1016/s0960-9822(03)00281-1
Yang, Y., Zhang, X., Su, Y., Zou, J., Wang, Z., Xu, L., et al. (2017). miRNA Alteration Is an Important Mechanism in Sugarcane Response to Low-Temperature Environment. Bmc Genomics 18, 833. doi:10.1186/s12864-017-4231-3
Yao, Z., Li, X., and Asimuguli, A. (2019). Analysis of NSP Family Gene Sequence and Expression in Response to Drought Stress in Arabidopsis thaliana. Mol. Plant Breed. 17, 1169. doi:10.13271/j.mpb.017.001169
Ye, J. Y., Tian, W. H., and Jin, C. W. (2019). A Reevaluation of the Contribution of NRT 1.1 to Nitrate Uptake in Arabidopsis under Low-Nitrate Supply. FEBS Lett. 593, 2051–2059. doi:10.1002/1873-3468.13473
Yu, C., Chen, Y., Cao, Y., Chen, H., Wang, J., Tian, F., et al. (2018). Overexpression of miR169o, an Overlapping MicroRNA in Response to Both Nitrogen Limitation and Bacterial Infection, Promotes Nitrogen Use Efficiency and Susceptibility to Bacterial Blight in Rice. Plant & Cell. Physiology 59, 1234–1247. doi:10.1093/pcp/pcy060
Zhang, S., Sun, L., and Kragler, F. (2009). The Phloem-Delivered RNA Pool Contains Small Noncoding RNAs and Interferes with Translation. Plant Physiol. 150, 378–387. doi:10.1104/pp.108.134767
Zhang, J. H., Lin, Y., Wu, F. M., Zhang, Y. T., Cheng, L. J., Huang, M. H., et al. (2021). Profiling of MicroRNAs and Their Targets in Roots and Shoots Reveals a Potential MiRNA-Mediated Interaction Network in Response to Phosphate Deficiency in the Forestry Tree Betula Luminifera. Front. Genet. 12. doi:10.3389/fgene.2021.552454
Zhang, J., Huang, M., Liang, J., Pan, Y., Cheng, L., Wu, J., et al. (2016). Genome-wide Mining for microRNAs and Their Targets in Betula Luminifera Using High-Throughput Sequencing and Degradome Analyses. Tree Genet. Genomes 12, 99. doi:10.1007/s11295-016-1047-2
Zhang, J., Zhou, Z., Bai, J., Tao, X., Wang, L., Zhang, H., et al. (2020). Disruption of MIR396e and MIR396f Improves Rice Yield Under Nitrogen-Deficient Conditions. Natl. Sci. Rev. 7, 102–112. doi:10.1093/nsr/nwz142
Zhao, M., Ding, H., Zhu, J. K., Zhang, F., and Li, W. X. (2011). Involvement of miR169 in the Nitrogen-Starvation Responses in Arabidopsis. New Phytol. 190, 906–915. doi:10.1111/j.1469-8137.2011.03647.x
Zhao, M., Guo, R., Li, M., Liu, Y., Wang, X., Fu, H., et al. (2020). Physiological Characteristics and Metabolomics Reveal the Tolerance Mechanism to Low Nitrogen in Glycine Soja Leaves. Physiol. Plant. 168, 819–834. doi:10.1111/ppl.13022
Zhao, M., Tai, H., Sun, S., Zhang, F., Xu, Y., and Li, W.-X. (2012). Cloning and Characterization of Maize miRNAs Involved in Responses to Nitrogen Deficiency. PLoS One 7, e29669. doi:10.1371/journal.pone.0029669
Zhao, W., Yang, X., Yu, H., Jiang, W., Sun, N., Liu, X., et al. (2015). RNA-Seq-Based Transcriptome Profiling of Early Nitrogen Deficiency Response in Cucumber Seedlings Provides New Insight into the Putative Nitrogen Regulatory Network. Plant Cell. Physiology 56, 455–467. doi:10.1093/pcp/pcu172
Zhu, H., Zhang, Y., Tang, R., Qu, H., Duan, X., and Jiang, Y. (2019). Banana sRNAome and Degradome Identify microRNAs Functioning in Differential Responses to Temperature Stress. BMC Genomics 20, 33. doi:10.1186/s12864-018-5395-1
Zhu, Z., Li, D., Cong, L., and Lu, X. (2021). Identification of microRNAs Involved in Crosstalk between Nitrogen, Phosphorus and Potassium under Multiple Nutrient Deficiency in Sorghum. Crop J. 9, 465–475. doi:10.1016/j.cj.2020.07.005
Keywords: miRNA, nitrogen starvation, target genes, Betula luminifera, sRNA sequencing
Citation: Lin Y, Chu S, Xu X, Han X, Huang H, Tong Z and Zhang J (2022) Identification of Nitrogen Starvation-Responsive miRNAs to Reveal the miRNA-Mediated Regulatory Network in Betula luminifera. Front. Genet. 13:957505. doi: 10.3389/fgene.2022.957505
Received: 31 May 2022; Accepted: 24 June 2022;
Published: 17 August 2022.
Edited by:
Yuepeng Song, Beijing Forestry University, ChinaReviewed by:
Chunlong Li, Huazhong Agricultural University, ChinaJinhui Chen, Hainan University, China
Copyright © 2022 Lin, Chu, Xu, Han, Huang, Tong and Zhang. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Zaikang Tong, zktong@zafu.edu.cn; Junhong Zhang, zhangjunhong@zafu.edu.cn
†These authors have contributed equally to this work