1. Introduction
The vascular system provides the transport system for water, minerals, and photosynthetic products in plants. It also supports the upright growth of plants and stores most of the biomass accumulated from photosynthesis. The vascular system consists of the xylem, which transports water, and the phloem, which transports nutrients [1]. Xylogenesis is considered as a sequential and complex process involving cambial cell mitosis, differentiation, secondary wall biosynthesis and lignification [2]. Cellulose is generated from glucose during the formation of β-1,4 glycosidic bonds to produce linear glucans, and is catalysed by cellulose synthase [3]. Cellulose is the most abundant natural organic polymer in nature [4]. Hemicelluloses are a type of heteropolysaccharide consisting of two or more monosaccharides with variable contents and structures [5]. As an important component of the secondary cell wall of plants, lignin is distributed in the cell walls of transport and lignified tissues and plays an important role not only in improving the water barrier function and mechanical strength of the cell wall, but also in increasing the resistance of plants to disease and stress [6,7].
Noncoding RNAs (ncRNAs) are a newly defined class of regulatory RNAs that do not encode proteins and were discovered in last century. They can interact with protein-coding genes to regulate gene expression. The ncRNAs can be classified into long noncoding RNAs (long ncRNAs) and small noncoding RNAs (sRNAs) according to their length [8]. The former mainly include circular RNAs (circRNAs) and lncRNAs, which are longer than 200 nt and usually display low sequence conservation and vary greatly between species; sRNAs include microRNAs (miRNAs) and small interfering RNAs (siRNAs), which are usually shorter in length but more conserved [8]. Plant sRNAs can be divided into two major categories, siRNAs and miRNAs, according to their biosynthetic processes and functional differences [9]. Regulatory sRNAs, represented by miRNAs and siRNAs, have been hot topics of research in the field of plant molecular genetics for almost two decades. SiRNAs generally originate from double-stranded RNA (dsRNA), are derive from variable locations and can be subdivided into heterochromatic siRNAs (hc-siRNAs), secondary siRNAs, and natural antisense transcript siRNAs (NAT-siRNAs) [9]. MiRNAs are highly evolutionarily conserved sRNAs of 20–24 nt in length and generated via single-stranded RNA folding to form hairpin RNA (hpRNA); miRNAs are derived from relatively fixed locations and degrade the transcripts of target genes or inhibit their translation by complementary base pairing with a target site, thereby regulating the growth and development of plants [9]. As miRNA research continues, increasing numbers of miRNAs and their target genes have been reported in plants; a large number of studies have elucidated the biosynthesis and functional mechanisms of miRNAs, and their biological functions have been widely revealed [10].
The biosynthesis of miRNA occurs via a multistep, multilevel reaction that includes the transcription of miRNA precursors, processing and transport of miRNAs, and assembly of the RNA-induced silencing complex (RISC); the first two steps occur mainly in the cell nucleus, while the RISC complex is formed in the cytoplasm [11]. The Transcription of MIR genes is mostly dependent on RNA polymerase II (Pol II) and is rarely mediated by RNA polymerase Ⅲ (Pol III) [12,13,14]. MIR genes are transcribed into primary miRNAs (pri-miRNAs) via Pol II- or Pol III-mediated catalysis. Pri-miRNAs have a stem-loop structure with a cap structure at the 5’ end and a polyA tail at the 3’ end and vary in length from a few hundred to thousands of nucleotides. Under the combined action of Dicer-like1 (DCL1), HYPONASTIC LEAVES1 (HYL1), and SERRATE (SE), pri-miRNA is further processed to form a precursor miRNA (pre-miRNA) [15]. Subsequently, the miRNA-processing machinery, the DCL1-HYL1-SE complex, then processes the pre-miRNA to form an miRNA/miRNA* duplex of approximately 21–25 nt in length [16]. The duplex is transferred from the cell nucleus to the cytoplasm via HASTY (HST) after being methylated by the transmethylation enzyme HUA ENHANCER 1 (HEN1) [17,18,19]. SE, HYL1, and HST have also been reported to have novel functions in miRNA biosynthesis. In addition to assisting DCL1 processing, SE can promote pri-miRNA degradation by recruiting Nuclear Exosome Targeting (NEXT) complexes and can affect the miRNA processing efficiency by promoting D-body formation [20,21,22]. In addition to forming homodimers to promote the efficiency and accuracy of DCL1 processing, HYL1 protects pri-miRNA from nucleic acid exonuclease attacks and maintains the dynamic balance of miRNA biosynthesis [23,24,25]. HST can also promote MIR gene transcription and processing by interacting with DCL1 [26]. Furthermore, newly identified processing factors such as the TREX complex (transcription/export), MOS4-associated complex 5 (MAC5), DAWDLE (DDL), PROTEIN PHOSPHATASE4 (PP4), and SMALL1 (SMA1), are involved in miRNA biosynthesis and processing, and most of these proteins regulate miRNA biosynthesis by interacting with DCL1 and HYL1 [25,27,28,29]. In the cytoplasm, Argonaute proteins (AGOs) bind to miRNA/miRNA* duplexes and release miRNA*, which is usually degraded [30,31]. AGO consists of the following four structural domains: the N-terminal, PAZ, MID, and PIWI domains. The PIWI domain displays RNaseH activity and is an essential functional component in the degradation of mRNAs by miRNAs [32]. When an miRNA binds to its target gene through complementary base pairing, the PIWI structural domain cuts the target gene at the site corresponding to positions 10–11 at the 5’ end of the miRNA, and the unprotected target gene is rapidly degraded by nucleases in the cell, thus losing its function [32].
Poplar is an important industrial timber species in the mid-latitudes. Species of section Aigeiros have become some of the most widely cultivated species worldwide due to their rapid growth and wide adaptability [33]. MiRNAs have been reported to function in the xylem of P. tomentosa and P. trichocarpa, but not in P. deltoides [34,35]. MiRNAs that can regulate wood formation have also been reported. In P. trichocarpa, miR397a was found to target a range of laccase (LAC) genes [36]. In Arabidopsis thaliana, miR397b negatively regulates LAC4. An overexpression of miR397b reduces lignin deposition and decreases fiber cell wall thickness [37]. In maize, miR528 can regulate lignin biosynthesis by negatively regulating ZmLAC3 and ZmLAC5 [38]. Wood formation is regulated by a three-level regulatory network of transcription factors with NAC-MYB at its core [39]. Numerous studies have reported that miRNA can target MYB and NAC. In A. thaliana, miR159 and miR399 was able to target MYB; in Oryza sativa, miR159 was found to target MYB; in Gossypium hirsutum, miR447, miR828, miR858, and miR 5255 was able to target MYB; in A. thaliana and Z. mays, miR164 was found to target NAC; in Triticum aestivum, miR164 was found to target NAC1 [40,41,42,43,44]. However, it is unknown whether similar mechanisms exist in P. deltoides [39]. Here, we constructed six sRNA libraries for developing xylem and phloem derived from P. deltoides stems. We then systematically identified miRNAs involved in the developing xylem and phloem and further analysed miRNAs that may be involved in lignin and cellulose biosynthesis. This study will provide new insights into poplar wood formation.
2. Results
2.1. Construction of Small RNA Libraries and Sequencing
High quality RNA was extracted from developing xylem (X1, X2, and X3) and phloem (P1, P2, and P3) for sRNA sequencing. A total of 67,371,464 raw reads (3.37 G base) were obtained via small RNA sequencing, with mean Q20, Q30, and GC content values of 98.73%, 95.63%, and 50.65%, respectively. After filtering out low-quality reads (in which bases with a quality value of less than 20 account for more than 30% of the entire read), 275,357 reads in which N bases accounted for more than 10% of the read (unidentifiable bases accounted for more than 10% of the entire read), 106,181 reads with a 5′ adapter, 1,155,625 reads without a 3′ adapter, and 55,266 reads with poly A/T/C/G sequences, a total of 65,779,035 clean reads were obtained. Studies have shown that the length distribution of plant small RNAs ranges from 18 to 30 nt. A total of 52,233,007 reads of 18–30 nt in length were obtained using small RNA sequencing, including 23,684,340 and 28,548,667 reads specific to the developing xylem and phloem, respectively. Most of the sRNAs were 21 nt in length, and these sRNAs constituted 16.14% of the reads distributed within the six libraries (Figure S1). Finally, we used Bowtie to map the reads to the P. trichocarpa genome, and the X1, X2, X3, P1, P2, and P3 mapping rate values were 90.07%, 95.10%, 91.37%, 97.40%, 97.22%, and 97.39%, respectively. The above information indicates that the applied small RNA sequencing protocol met the requirements for the subsequent analysis. A brief summary and detailed information of the sRNA libraries are listed in Table 1 and Table S1.
2.2. Identification of MiRNAs Involved in the Developing Xylem and Phloem
A total of 194 miRNAs (97 mature miRNAs and 97 miRNAs*) were identified in P. deltoides according to the updated criteria for plant miRNAs, 96 of which were known miRNAs. We found that these known miRNAs belonged to 32 miRNA families. The largest miRNA family was the miR166 family, which had 13 members in the P. deltoides genome, followed by the miR319, miR396, miR160, miR167, miR395, and miR164 families, which had 10, 7, 6, 6, 6, and five members in the P. deltoides genome, respectively (Figure 1a, Table S2). After filtering out the miRNAs with low expression, we identified 134 miRNAs with high expression for the subsequent analysis, including 132 known and two novel miRNAs. Reads per ten million (RPTM) values were used to measure the expression levels of miRNAs, and the abundances of the miRNAs varied greatly. The RPTM values of pde-miR319d-mature, pde-miR166a-mature, pde-miR472a-mature, pde-miR396d-mature, and pde-482a-mature were greater than 10,000, and the values of pde-482a-star, pde-miR395f-star, pde-miR6427-star, pde-miR171b-star, and pde-miR476-star were lower than 100. The RPTM values also varied greatly between different members of the same miRNA family, as observed in the pma-miR159 family. The RPTM value of Pde-miR477c-mature in the phloem and developing xylem was greater than 500, but the RPTM value of Pde-miR477a-mature was lower than 300 (Table S3). In general, we found that miRNAs were highly expressed in the phloem relative to the developing xylem (Figure 1b). A total of 134 miRNAs ranging from 20–23 nt in length were identified; most of the miRNAs were 21 nt in length and accounted for more than 50% of the total miRNAs, which is consistent with the criteria established in previous studies (Figure 1c). Detailed information about the mature and precursor sequences of the 134 miRNAs is provided in Table S4.
2.3. Identification of Differentially Expressed miRNAs
A total of 58 differentially expressed miRNAs (DEmiRNAs) were identified in the developing xylem in contrast to the phloem. Among the 58 DEmiRNAs, 21 miRNAs (21 known miRNAs) were more highly expressed in developing xylem than in phloem, while 37 miRNAs (35 known miRNAs and two novel miRNAs) were less highly expressed in developing xylem than in phloem (Figure 2a). DEmiRNAs accounted for 43.28% (58 out of 198) of all miRNAs. The top five upregulated DEmiRNAs in developing xylem relative to phloem were miR171a-mature, miR171d-mature, miR171a-star, miR171d-star, and miR394-mature. The log2 FC values of the five miRNAs were all greater than five. The most changed miRNAs were miR171a-mature and miR171d-mature, whose expression values in developing xylem and phloem were 1064.13 and 50.09, respectively (Figure 2b, Table S5). Their expression in developing xylem was twenty times higher than that in phloem, suggesting that these miRNAs may play an important role in developing xylem growth and development. The top downregulated DEmiRNAs in developing xylem relative to phloem were from the miR166 family and included miR166a-mature, miR166b-mature, miR166c-mature, miR166e-mature, miR166f-mature, miR166g-mature, miR166i-mature, miR166j-mature, miR166l-mature, miR166m-mature and miR166n-mature. The log2 FC values of the five miRNAs were all greater than seven (Figure 2b, Table S5). The miR166 family is highly conserved, and miR166 family members have been identified in several species by the efforts of researchers [45]. MiR166 has been found to be involved in seed, vascular, and root development. MiR166 family members showed dramatic changes in developing xylem and phloem expression, suggesting that they also play an important role in the growth and development of poplar [45].
2.4. Identification of MiRNA Target Genes via Degradome Sequencing
A total of 2431 target genes of the 134 miRNAs were obtained via high-throughput sequencing. The miRNAs and their target genes could form 5180 regulatory pairs. The complete list of miRNA–mRNA regulatory pairs is shown in Table S6. Most miRNAs can target more than two target genes. For example, miR396 family members (miR396b-mature, miR396f-mature, miR396a-mature, miR396c-mature, miR396g-mature and miR396e-mature) could target approximately 30 different transcripts, while miR160 (miR160a-mature, miR160b-mature, miR160c-mature, miR160d-mature, miR160e-mature and miR160f-mature) and miR167 (miR167a-mature, miR167b-mature, miR167c-mature, miR167d-mature, miR167e-mature and miR167g-mature) family members could target approximately 20 different transcripts. Only 16 miRNAs appeared to target a single gene (Table S6). The target genes were enriched in nine KEGG pathways, and the top five terms were “chaperones and folding catalysts”, “ribosome biogenesis”, “glyoxylate and dicarboxylate metabolism”, “amino sugar and nucleotide sugar metabolism”, and “phenylalanine, tyrosine and tryptophan biosynthesis” (Figure 3). The top five terms among biological processes were “response to hormone”, “developmental process”, “DNA recombination”, “transcription, DNA-templated” and “histone modification”; the top two terms among cellular components were “COPI vesicle coat” and “spliceosomal complex”; and the top five terms among molecular function were “four-way junction helicase activity”, “lipid binding”, “unfolded protein binding”, “phosphatidylinositol phosphate kinase activity”, and “pyruvate kinase activity” (Figure S2). Table 2 lists several miRNAs target genes.
2.5. MicroRNAs Involved in the Regulation of Lignin and Cellulose Biosynthesis
A total of 13 miRNAs identified in this study may be involved in the lignin biosynthesis pathway, including 12 known miRNAs (miR482c-mature, miR530-star, miR319a/b/c/d/g/h/j-mature, miR403a/b-mature, and miR472a-star) and one novel miRNA (miRN8-mature). MiR482c-mature and miRN8-mature could target CAD genes, thereby inhibiting the conversion of ρ-coumaraldehyde, coniferaldetyde, and sinapaldehyde to ρ-coumarylacohol, coniferylalcohol, and sinapylalcohol, respectively. MiR530-star could target the LAC gene, and miR319a/b/c/d/g/h/j-mature, miR403a/b-mature, and miR472a-star could target POX genes, thereby inhibiting the polymerization of lignin monomers into polymeric lignin (Figure 4a). Among these 13 miRNAs, the top five highly expressed miRNAs were miR319d-mature, miR319h-mature, miR319a-mature, miR319b-mature, and miR319j-mature, all of which belong to the miR319 family. Additionally, miR319 family members could target MYB genes (Potri.001G036000.1, Potri.001G224500.2, Potri.003G189700.7) and the TCP gene (Potri.008G182300.4). MiR403 and miR472 family members could target AGO genes (Potri.012G118700.3, Potri.012G118700.4), which encode proteins that are essential for miRNAs function. MiR482 family members could target R genes (Potri.018G138500.2, Potri.003G201800.6), and an LRR gene (Potri.007G027500.4). MiR530 family members could target bHLH genes (Potri.014G099700.1, Potri.014G099700.2), a GRAS gene (Potri.003G110800.3), and a WRKY gene (Potri.018G107000.5). Zinc finger protein (Potri.005G115100.3) and F-box genes (Potri.018G048700.4) could be targeted by the novel miRN8 (Table S6).
A total of nine miRNAs may be involved in the cellulose biosynthesis pathway, all of which are known miRNAs. MiR319c-mature, miR477c-mature, and miR475a-mature could target SUSY genes, thereby inhibiting the mutual transformation of sucrose and UDP-glucose; UDP-glucose can also be transformed to cellulose when catalysed by CesA, and the CesA genes could be targeted by five miRNAs (miR396c-star, miR477a-star, and miR319e/f/i-mature). MiR166f-star could target SPS genes, thereby inhibiting UDP-glucose conversion to sucrose 6P (Figure 4b). Additionally, miR166 family members could target WRKY (Potri.010G160100.5), P450 reductase (Potri.005G153800.1), and NAC genes (Potri.002G037100.1). MiR475 family members could target citrate synthase genes (Potri.014G141900.1, and Potri.016G089300.1). MiR396 family members could target GRF genes (Potri.006G115200.2, Potri.014G007200.1, and Potri.007G007100.1), SCARECROW-like genes (Potri.016G009700.1, and Potri.016G009700.2), and an AUX gene (Potri.005G174000.1) (Table S6).
3. Discussion
3.1. Annotation of MiRNAs Involved in the Developing Xylem and Phloem of P. deltoides
A large number of studies have revealed the potential roles of miRNAs in plant growth and development and the responses to biotic and abiotic stresses, but no studies have yet revealed the functions of miRNAs in the developing xylem and phloem of P. deltoides [39]. In this study, we constructed six sRNA libraries (three for developing xylem and three for phloem) to systematically identify miRNAs involved in the developing xylem and phloem of P. deltoides. We obtained 67,371,464 raw reads via sRNA sequencing and 65,779,035 clean reads after quality control (Table S1). Studies have shown that the length distribution of plant sRNAs ranges from 18 to 30 nt. A total of 52,233,007 reads were obtained for miRNA prediction and analysis after the filtering of reads that were too long or too short (Table S1). Based on the most recent plant miRNA annotation criteria [46], a total of 134 miRNAs were identified as highly expressed in at least one library, including 132 known miRNAs and two novel miRNAs, among which six miRNAs (miR164c-star, miR166f-star, miR395f-star, miR399b-star, miR6427-star, and miR6434-mature) were only expressed in phloem, and their RPTM values were 32.58, 11.61, 7.79, 23.53, 10.79, and 93.30, respectively (Tables S2 and S3). We found substantial differences in the expression between developing xylem and phloem for 58 miRNAs, with 21 miRNAs significantly upregulated in developing xylem in contrast to the phloem and 37 miRNAs that were significantly down-regulated in developing xylem. We found that most miRNAs tend to be highly expressed in phloem (Figure 2b). We validated several miRNAs involved in lignin and cellulose biosynthesis via qRT-PCR, and the qRT-PCR and sRNA sequencing results of most of the selected miRNAs showed relatively high agreement (Figure 5). A total of 134 miRNAs were distributed among 32 different miRNA families. Most miRNA families presented more than one member in P. deltoides (miR160, miR166, miR167, miR319, and miR396), but some miRNA families were found to have only one member in P. deltoides, such as the miR172, miR394, miR475, miR476, miR530, miR6421, miR6427, miR6434, miR6459, miR6478, and miR7817 families (Table S2). The top five families were miR166, with 13 members, miR319 with 10 members, miR396 with seven members, miR160 with six members, and miR167 with six members. The length distribution of the 134 miRNAs ranged from 20 to 23 nt, and 21 nt miRNAs accounted for 79.10% of the total miRNAs, far exceeding the number of miRNAs of other lengths. A similar phenomenon is found in species such as Pinus massoniana (56.57%) and Cinnamomum camphora (more than 47%).
3.2. MiRNAs Families in Monocots, Eudicots, and Gymnosperms
MiRNAs belong to a class of non-coding genes and, similar to protein-coding genes, miRNAs are subject to constant evolution. Over the long period of plant evolution, many lineage-specific miRNA families have emerged, and these miRNA families have attracted much attention from researchers. For example, miR1515 is a species-specific miRNA present in Citrus sinensis [47]. Twelve plant species, including P. deltoides, were selected to investigate the presence of P. deltoides-specific, or genus populus-specific miRNAs, or eudicot-specific miRNA families. Among the 12 plants, the following five were monocots: Zea mays, Phoenix dactylifera, Oryza sativa, Musa acuminata, and Ananas comosus; the following five were eudicots: Vitis vinifera, P. trichocarpa, P. deltoides, Glycine max, and A. thaliana; and the last two, Picea abies and Ginkgo biloba, were gymnosperms.
A total of 386 miRNA families were found in 12 plant species. The eudicot V. vinifera had 152 miRNAs belonging to 43 different miRNA families; P. trichocarpa had 135 miRNAs belonging to 35 different miRNA families; G. max had 386 miRNAs belonging to 86 different miRNA families; A. thaliana had 152 miRNAs belonging to 82 different miRNA families; and P. deltoides exhibited 32 miRNA families, which were similar to those of P. trichocarpa, indirectly indicating that P. deltoides and P. trichocarpa are more closely related than the other plants (Table S2). The monocot Z. mays had 185 miRNAs belonging to 36 different miRNA families; O. sativa had 273 miRNAs belonging to 93 different miRNA families; M. acuminata had 210 miRNAs belonging to 31 different miRNA families; P. dactylifera had 128 miRNAs belonging to 28 different miRNA families, and A. comosus had 123 miRNAs belonging to 31 different miRNA families (Table S2). The gymnosperm P. abies had 547 miRNAs belonging to 171 different miRNA families, and G. biloba had 188 miRNAs belonging to 37 different miRNA families (Table S2). Eleven of the 386 families were present in all 12 plant species and are as follows: miR160, miR162, miR164, miR166, miR167, miR168, miR169, miR171, miR319, miR394, and miR396, all of which are highly conserved miRNA families (Tables S2 and S7). In addition, 18 miRNA families were present among gymnosperms, monocots, and eudicots. There were 70, 134, and 143 miRNA families that were specifically present in monocots, eudicots, and gymnosperms, respectively. Five miRNA families were present in both monocots and eudicots (miR530, miR5072, miR479, miR827, and miR158); five miRNA families were present in both gymnosperms and eudicots (miR4414, miR2111, miR2950, miR3627, and miR477); no miRNA families were exclusively present in gymnosperms and monocots; and 29 miRNA families were present across gymnosperms, monocots, and eudicots (Figure 6, Table S7). Most miRNA families (84.46%) were species specific, including 135 miRNA families present only in P. abies; 58 present only in O. sativa; 54 present only in A. thaliana; and 54 present only in G. max (Figure S3).
We found there were some differences in the miRNA family distribution and members among four species of Populus (P. deltoides, P. trichocarpa, P. alba, and P. tremula). P. trichocarpa had 135 miRNAs belonging to 35 different miRNA families; P. alba had 149 miRNAs belonging to 32 different miRNA families; P. tremula had 143 miRNAs belonging to 31 different miRNA families, and P. deltoides exhibited 32 miRNA families. There were 18 miRNA families present in all four species of Populus (P. deltoides, P. trichocarpa, P. alba, and P. tremula); one miRNA family (miR6459) present in P. deltoides, P. trichocarpa, and P. alba; three miRNA families (miR399, miR395, and miR6457) present in P. deltoides, P. alba, and P. tremula; six miRNA families (miR156, miR159, miR1448, miR397, miR390, and miR393) present in P. trichocarpa, P. alba, and P. tremula; three miRNA families (miR7841, miR828, and miR536) present in P. alba and P. tremula; one miRNA family (miR477) present in P. deltoides and P. alba; five miRNA families (miR476, miR475, miR530, miR6421, and miR2111) present in P. deltoides and P. trichocarpa; one miRNA family (miR1447) present in P. trichocarpa and P. tremula. Furthermore, four miRNA families (miR6434, miR6478, miR1446, and miR7817) were P. deltoides-specific; and four miRNA families (miR6456, miR1444, miR6425, and miR408) were P. trichocarpa-specific, (Figure S4). The top five miRNA families in P. deltoides were miR166, miR319, miR396, miR160, and miR164; in P. trichocarpa they were miR166, miR171, miR156, miR169, and miR164; in P. alba they were miR166, miR171, miR156, miR169, and miR167; and in P. tremula they were miR166, miR171, miR169, miR156, and miR167. Some miRNA families have similar members between the four species of Populus, such as miR166, miR160, miR168, miR390, miR394, miR396, and miR472, while some miRNA families are the opposite, such as miR1444, miR1446, and miR1447 (Figure S5, Table S8). The number of some miRNA families in P. deltoides was lower than in other species of Populus, such as miR156, miR390, and miR397. We think that this result is due to the fact that there were only two tissues in this study, whereas the miRNA data from other species of Populus originated from multiple tissues and differences in the depth of research on miRNAs in different plants.
3.3. Most Target Genes Are Transcription Factors
Two types of regulators play an important role in eukaryotic post-transcriptional gene expression, and are as follows: transcription factors and miRNAs, both of which perform regulatory functions at the transcriptional and posttranscriptional levels [48]. Transcription factors are a class of proteins with specific functions that initiate transcription by binding to the promoter region of a gene or recruiting other proteins to the promoter region [49].
MiRNAs are a class of small noncoding RNAs discovered in the last century that can directly degrade target mRNAs or inhibit the translation of target genes by binding to their 3’ UTRs [11]. Transcription factors and miRNAs have also been found to be mutually regulated, with transcription factors binding directly to the promoter regions of miRNAs and regulating miRNA biogenesis, while miRNAs can directly target the mRNAs of transcription factors, making it impossible for the proteins to perform their original roles [49]. In this study, the majority of identified miRNA target genes were found to encode transcription factors, including AP2, ARF, bHLH, bZIP, GRAS, GRF, MYB, NAC, TCP, and WRKY genes (Figure 7; Table S6). For example, miR172, targeting AP2, has been reported to be involved in floral development in A. thaliana and to be involved in nodule formation in G. max and Phaseolus vulgaris [50,51,52]. MiR160, targeting ARF, is reported to be involved in the vascular cambium initiation process in P. tomentosa and to be involved in root growth in A. thaliana [34,53]. In A. thaliana, miR159, targeting MYB, has been reported to be involved in seed germination [54]; in G. hirsutum, miR447 and miR5255, targeting MYB, have been reported to be involved in root and fibre development [43]. Furthermore, miR319, targeting TCP, is reported to be involved in regulating secondary growth in P. tomentosa and to be involved in floral and leaf development in A. thaliana [34,43].
3.4. Conserved MicroRNAs-mRNAs Modules Involved in the Biosynthesis of Lignin and Cellulose
A total of 13 miRNAs that may be involved in lignin biosynthesis in P. deltoides were identified in this study. MiR482c-mature and miRN8-mature, targeting CAD genes, inhibit the conversion of ρ-coumaraldehyde, coniferaldetyde, and sinapaldehyde to ρ-coumaraldehyde, coniferylalcohol, and sinapylalcohol, respectively. CAD is one of the key enzymes in lignin biosynthesis, producing precursors of lignin monomers (ρ-coumaryl alcohol, coniferyl alcohol, and sinapyl alcohol) [7]. A novel miRNA whose target gene was CAD was identified in Artemisia annua [55]; in addition, the level of Pde-miR5913 was shown to decrease after infestation with Marssonina brunnea, which led to an upregulation of its target gene CAD, assisting the plant in dealing with biotic stress [56]. MiR530-star, targeting LAC, inhibits the polymerization of lignin monomers into lignin. LAC has been reported to be regulated by miRNAs in Pinus massoniana, Z. mays, P. trichocarpa, Linum usitatissimum, and Pyrus bretschneideri, implying that the role of the miRNA-LAC module in lignin biosynthesis may be conserved in gymnosperms and angiosperms [36,38,39,57]. POX can be targeted by 10 miRNAs, most of which are members of the miR319 family. The target gene of miR3493b in T. aestivum encodes the peroxidase 52 protein, suggesting that this miRNA is also a conserved miRNA-mRNA regulatory module [58].
A total of nine miRNAs that may be involved in cellulose biosynthesis in P. deltoides were identified in this study. MiR319c-mature, miR477c-mature, and miR475a-mature can target SUSY by inhibiting sucrose conversion to UDP-glucose. SUSY is a key enzyme in plant sugar metabolism that catalyses the biosynthesis and breakdown of sucrose and plays an important regulatory role in starch and cellulose biosynthesis, fruit quality development, reproductive growth, and sugar signal transduction. In Panicum virgatum, miR156 has been reported to be a regulator of SUSY [59]. MiR396c-star, miR477a-star, miR319e-mature, miR319f-mature, and miR319i-mature can target CesA by inhibiting UDP-glucose conversion to cellulose. CesAs in plants are components of a cellulose synthase complex (CSC) that is localised on the cell membrane and contains six large subunits, each including at least three cellulose synthases [60]. The CSC contains a large number of proteins, including CesA proteins, cellulose-like synthase (CSL) proteins, and KORRIGAN proteins [61]. In A. thaliana, miR395a has been reported to regulate a transcript encoding cellulose synthase [62]. MiR166f-star can target SPS by inhibiting UDP-glucose conversion to sucrose 6P. SPS is widespread in plants and is a key rate-limiting enzyme in the biosynthesis of sucrose. The amino acid sequences of SPSs in plants are generally long, and these SPSs contain additional long amino acid sequences that can be subjected to posttranslational modifications that in turn regulate the enzyme activity of SPSs [63,64]. In Sorghum bicolor, a novel miRNA has been reported to target a transcript encoding sucrose phosphate synthase [65]. Furthermore, some miRNAs, such as miR166a-star and miR166e-mature, target homeobox genes, which have been linked to wood development. Overexpression of ATHB8 (a homeobox gene family member) in A. thaliana accelerates lignin accumulation [66].
Conserved miRNAs-mRNAs modules involved in wood formation are also reported in other species of Populus. A study on the potential role of miRNAs in wood formation in P. trichocarpa showed that there were 57 miRNAs enriched in xylem, of which 11 miRNAs were enriched in mechanically treated xylem (MTX). MiRX50A, an MTX-specific miRNA, has been predicted to target xyloglucan endotransglycosylase/hydrolase16 (XTH16) and is thought to play a critical role in G-fiber formation and function. The xylem-specific miRNAs, miRX50 and miRX87, may be involved in wood formation via targeting the NAC domain TFs NAC083 and NAC050, respectively. MiRX41 is an MTX-specific miRNA involved in wood formation via targeting CSLD4, which encodes cellulose synthase [35]. According to a study on the potential role of miRNAs in diverse phases of wood formation in P. trichocarpa, members of the miR858 family, miR858x and miR858y, directly regulate secondary wall biosynthesis by targeting MYB35, MYB52, MYB63, and MYB83. Novel-m0998-5p could target MYB5, which forms MYB–bHLH–WDR TF complexes to activate the PAL gene, which encodes the first enzyme of the lignin biosynthesis pathway. Novel-m0260-5p could target C4H; novel-m1190-3p could target LAC2; and miR397x could target LAC11. They are all related to lignin biosynthesis and deposition [67].
4. Materials and Methods
4.1. Plant Materials and Small RNA Sequencing
Phloem and developing xylem tissues were derived from three 10-year-old ramets of P. deltoides cv. ‘I-69′. The developing xylem was harvested separately from the vigorously growing stems of trees in mid-July of the current year. The bark was peeled to expose the youngest developing secondary xylem, which was gently scraped from the wood with a razor blade and used for RNA extraction. Three ramets of each clone represented three biological replicates. Total RNA was extracted from the developing xylem and phloem of P. deltoides using an RNAprep Pure kit (Tiangen, Beijing, China) according to the manufacturer’s instructions. In our previous publications, we described the methods applied for the construction and sequencing of small RNA libraries [39]. In brief, the concentration of total RNA extracted from the developing xylem and phloem of P. deltoides was measured using the Qubit® RNA Assay Kit (Life Technologies, Carlsbad, CA, USA), and RNA integrity was assessed according to RNA Integrity Number (RIN). Six small RNA libraries were constructed using the NEBNext® Multiplex Small RNA Library Prep Set for Illumina® kit (NEB, Ipswich, MA, USA), and the six sRNA libraries were then sequenced using the Illumina HiSeq 2500 platform for 1 × 50 bp sequencing (Novogene, Beijing, China).
4.2. MiRNA Annotation Pipeline
The raw data of six small RNA libraries were used to assess sequencing quality with FastQC [68]. We used Cutadapt to filter out the reads that did not meet the criteria for plant small RNA [35], which were less than 18 nt or more than 30 nt in length, after using dnapi.py to obtain the adaptor sequences [69]. To increase the quality of clean data, we also filtered out reads containing more than one N (ambiguous base) base and required the quality values of all bases in the reads to be greater than 20. Studies have reported that small RNA fragments derived from rRNA or tRNA may interfere with the accuracy of miRNA annotation [70]. We used the Bowtie to remove reads that were derived from rRNA, tRNA, or organelle genomes, with one base mismatch permitted [71]. Then, miRNA annotation was performed using sRNAminer according to the most recent plant miRNA annotation criteria [46,72,73]. Briefly, miRNA/miRNA* duplexes do not contain a large loop, which limits the precursor and mature miRNA lengths to 300 nt and 20–24 nt, respectively; and miRNA/miRNA* duplexes contain up to five mismatches, only three of which are in the loop region. After removing miRNAs that were represented by less than 10 read counts in all samples, a total of 134 miRNAs were finally obtained. Known and novel miRNAs were identified based on comparisons with the plant miRNAs stored in the miRBase database, revealing 132 known miRNAs and two novel miRNAs. In this study, we used DEseq2 to perform the differential expression analyses of miRNAs. Firstly, we filtered out the miRNAs that were not expressed in all samples and used the DESeqDataSetFromMatrix function to construct a DESeqDataSet object; secondly, we used the DESeq function to perform differential expression analyses and used the sizeFactors function to calculate the scaling factor; thirdly, we used the results function to extract the results of differential expression analyses and used the estimateSizeFactors function to calculate the normalized counts of miRNAs. MiRNAs with a |log2 (fold change)| value of more than 1.5 and an Padj value of less than 0.05 were identified as differentially expressed [74]. The normalised expression levels of miRNAs were measured as RPTM (reads per ten million), which were calculated as follows:
(1)
4.3. Identification of MiRNA Target Genes
The RNA isolated from the developing xylem and phloem was mixed together in equal amounts to construct the degradome library, which was subjected to 1 × 50 bp sequencing on the Illumina HiSeq 2500 platform (Lianchuan Bio, Hangzhou, China). The whole-genome sequence, coding sequence, and gene function annotations of P. trichocarpa v4.1 were retrieved from the Phytozome database [75]. The degradome sequencing data were then subjected to adaptor removal using fastp with the default parameters [76]. Finally, we used CleavelLand4 to address the sequencing data and obtain the target genes of miRNAs [77]. GO and KEGG enrichment analyses of miRNA target genes were performed via clusterProfiler [78].
4.4. qRT-PCR Validation
The reference genes selected for use in qRT-PCR validation in this study were 18S and EF1α, and the sequences of the primers used are shown in Table S9 [79]. The detailed methods and steps can be found in our previous publication [39].
5. Conclusions
We identified miRNAs in the phloem and developing xylem of P. deltoides via small RNA sequencing and target genes of miRNAs via degradome sequencing. We found that most target genes of miRNAs were transcription factors. Furthermore, 13 and nine miRNAs were found to be involved in lignin and cellulose biosynthesis, respectively. This study provides new insights into the wood formation of poplar.
M.X., X.S. and H.P. conceived the project. C.D., H.Z. and N.R. performed the experiments. C.D. and T.S. analysed the data and wrote the manuscript. All authors have read and agreed to the published version of the manuscript.
Please add: This work was funded by the National Key Research and Development Program of China [grant number 2021YFD2201205], the National Natural Science Foundation of China (31971679), the Natural Science Foundation of the Jiangsu Higher Education Institutions of China (19KJB180001) and the Priority Academic Program Development of Jiangsu Higher Education Institutions (PAPD).
Not applicable.
Not applicable.
All the data are shown in the main manuscript and in the
We are grateful to Xinhua Chen for his assistance in data analysis.
The authors declare no conflict of interest.
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Enlarge this image.Figure 1. The miRNAs identified in the phloem and developing xylem of P. deltoides in this study. (a) Total numbers of known miRNA members in 12 plants. Only families with more than 10 members in at least one species are shown. Detailed information is provided in Table S2. (b) Expression value distribution of 134 miRNAs in the phloem and developing xylem of P. deltoides; the expression levels were log10 (RPTM + 1) normalised. The details of the expression of the 134 miRNAs are shown in Table S3. (c) Length distribution of the 134 miRNAs. The detailed sequences of the 134 miRNAs and their precursors are shown in Table S4.
Enlarge this image.Figure 2. The DEmiRNAs found in the phloem and developing xylem of P. deltoides in this study. (a) Volcano plot of DEmiRNAs in the developing xylem relative to the phloem, with red and green dots representing up- and downregulated miRNAs, respectively. Detailed information on the log2 FC and Padj values of 134 miRNAs is shown in Table S5. (b) Clustering of the expression profiles of DEmiRNAs in the phloem and developing xylem; the expression levels were normalized by the log10 (RPTM + 1) normalised. The DEmiRNAs were divided into different clusters based on their expression profiles. Note: Due to some miRNAs in the same miRNA family having equal y-axis and x-axis, they overlap in Figure 2a. Based on the most recent plant miRNA annotation criteria, plant miRNA identification pipelines require the mapping of sRNA sequencing reads to the reference sequence and extending the reads that could be mapped to the reference sequence upstream and downstream by 150bp each, as putative miRNA precursors. In this process, there are a number of reads that may be mapped to multiple positions that share the same sequence, belong to the same miRNA family, perform similar functions, and are distinguished by appending a, b, and c to the nomenclature, respectively. Although they share the same sequence, in reality they are still different genes and, in theory, should be expressed differently, but at this stage, sRNA sequencing does not distinguish the extent to which each gene locus contributes to the overall expression, so the overall expression is used, leading to this phenomenon.
Enlarge this image.Figure 3. Statistics of KEGG pathway enrichment of miRNA target genes in the phloem and developing xylem of P. deltoides. The y-axis presents the KEGG pathway terms, and the x-axis presents the rich factors of the KEGG pathway terms. The calculation steps were as follows: the rich factor of certain KEGG pathways is equal to the number of genes in that pathway among miRNA target genes divided by the number of genes in the P. deltoides genome annotated in the pathway. The sizes and colours of the dots represent gene counts and p values, respectively. Larger dots indicate a greater number of genes, and deeper red coloration indicates a more significant p value.
Enlarge this image.Figure 4. Typical pathways of lignin (a) and cellulose (b) biosynthesis in P. deltoides. Intermediates are shown in black, and the enzymes involved in each step are shown in green. The “T” arrows indicate the enzymes that could be cleaved by miRNAs supported by the degradome; expression levels was log10 (RPTM + 1) normalized.
Enlarge this image.Figure 5. Tissue expression patterns of candidate miRNAs involved in lignin and cellulose biosynthesis validated by qRT-PCR.
Enlarge this image.Figure 7. MiRNA-mRNA regulatory networks. The red and green dots indicate target genes and miRNAs, respectively.
Summary of sRNA libraries.
| Sample | Raw Reads | Clean Reads | Q30 (%) | GC (%) | Mapped Rate (%) |
|---|---|---|---|---|---|
| P1 | 11,419,418 | 11,161,174 | 95.96 | 50.74 | 90.07 |
| P2 | 11,768,425 | 11,432,450 | 95.47 | 51.42 | 95.10 |
| P3 | 11,997,478 | 11,760,792 | 95.83 | 50.77 | 91.37 |
| X1 | 10,320,359 | 10,100,623 | 95.85 | 50.13 | 97.40 |
| X2 | 10,796,007 | 10,499,531 | 95.17 | 50.58 | 97.22 |
| X3 | 11,069,777 | 10,824,465 | 95.52 | 50.27 | 97.39 |
Target genes of miRNAs identified in the phloem and developing xylem of P. deltoides by degradome sequencing.
| miRNA | Transcript |
Degradome |
Degradome |
|---|---|---|---|
| miR171a-mature | GRAS | 0 | 0.000271151 |
| miR160f-mature | ARF17 | 0 | 0.000271151 |
| miR167d-mature | ARF8 | 0 | 0.000271151 |
| miR396b-mature | GRF4 | 0 | 0.000271151 |
| miR472a-star | AGO3 | 0 | 0.000271151 |
| miR319a-mature | MYB65 | 0 | 0.000542229 |
| miR164d-mature | NAC80 | 0 | 0.000542229 |
| miR172c-mature | AP2 | 0 | 0.000813233 |
| miR160d-mature | ARF10 | 0 | 0.002708204 |
| miR160e-mature | ARF16 | 0 | 0.002708204 |
| miR396a-mature | GRF9 | 0 | 0.002708204 |
| miR530-mature | bHLH | 0 | 0.017206273 |
| miR396d-mature | SCL13 | 0 | 0.019336135 |
| miR171d-star | WRKY | 1 | 0.027125399 |
Note: In this study, we used Cleaveland4 software to predict target gene of miRNAs via degradome sequencing, which classifies the target genes confidence (Category) into five categories, with smaller categories generating more reliable results. Category0 represents the highest abundance of degradome sequencing reads matched to this target site and the site is unique, as can be seen from the Target plot where a sharp peak is formed at Cateogry0; Category1 represents the highest abundance of degradome sequencing reads matched to the target site but the site is not unique; Category2 represents an above average abundance of degradome sequencing reads matched to the target site; Category3 represents the below average abundance of degradome sequencing reads matched to the target site; and Category4 means that only one of the degradome sequencing reads matched the target site. Category3 or 4 means that the results are not very reliable and often need to be filtered out.
Supplementary Materials
The following supporting information can be downloaded at:
References
1. Fukuda, H.; Ohashi-Ito, K. Chapter Six—Vascular tissue development in plants. Current Topics in Developmental Biology; Grossniklaus, U. Academic Press: Cambridge, MA, USA, 2019; Volume 131, pp. 141-160.
2. De Rybel, B.; Mähönen, A.P.; Helariutta, Y.; Weijers, D. Plant vascular development: From early specification to differentiation. Nat. Rev. Mol. Cell Biol.; 2016; 17, pp. 30-40. [DOI: https://dx.doi.org/10.1038/nrm.2015.6] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/26580717]
3. McFarlane, H.; Döring, A.; Persson, S. The Cell Biology of Cellulose Synthesis. Annu. Rev. Plant Biol.; 2014; 65, pp. 69-94. [DOI: https://dx.doi.org/10.1146/annurev-arplant-050213-040240] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/24579997]
4. Cosgrove, D. Loosening of plant cell walls by expansins. Nature; 2000; 407, pp. 321-326. [DOI: https://dx.doi.org/10.1038/35030000] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/11014181]
5. Scheller, H.V.; Ulvskov, P. Hemicelluloses. Annu. Rev. Plant Biol.; 2010; 61, pp. 263-289. [DOI: https://dx.doi.org/10.1146/annurev-arplant-042809-112315]
6. Vaahtera, L.; Schulz, J.; Hamann, T. Cell wall integrity maintenance during plant development and interaction with the environment. Nat. Plants; 2019; 5, pp. 924-932. [DOI: https://dx.doi.org/10.1038/s41477-019-0502-0]
7. Liu, Q.; Luo, L.; Zheng, L. Lignins: Biosynthesis and Biological Functions in Plants. Int. J. Mol. Sci.; 2018; 19, 335. [DOI: https://dx.doi.org/10.3390/ijms19020335]
8. Amin, N.; McGrath, A.; Chen, Y.-P.P. Evaluation of deep learning in non-coding RNA classification. Nat. Mach. Intell.; 2019; 1, pp. 246-256. [DOI: https://dx.doi.org/10.1038/s42256-019-0051-2]
9. Axtell, M. Classification and comparison of SMALL RNAs from plants. Annu. Rev. Plant Biol.; 2013; 64, pp. 137-159. [DOI: https://dx.doi.org/10.1146/annurev-arplant-050312-120043]
10. Song, X.; Li, Y.; Cao, X.; Qi, Y. MicroRNAs and Their Regulatory Roles in Plant–Environment Interactions. Annu. Rev. Plant Biol.; 2019; 70, pp. 489-525. [DOI: https://dx.doi.org/10.1146/annurev-arplant-050718-100334]
11. Yu, Y.; Jia, T.; Chen, X. The ‘how’ and ‘where’ of plant microRNAs. New Phytol.; 2017; 216, pp. 1002-1017. [DOI: https://dx.doi.org/10.1111/nph.14834]
12. Kim, Y.J.; Zheng, B.; Yu, Y.; Won, S.; Mo, B.; Chen, X. The role of Mediator in small and long noncoding RNA production in Arabidopsis thaliana. EMBO J.; 2011; 30, pp. 814-822. [DOI: https://dx.doi.org/10.1038/emboj.2011.3] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/21252857]
13. Xie, Z.; Allen, E.; Fahlgren, N.; Calamar, A.; Givan, S.; Carrington, J. Expression of Arabidopsis miRNA Genes. Plant Physiol.; 2005; 138, pp. 2145-2154. [DOI: https://dx.doi.org/10.1104/pp.105.062943] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/16040653]
14. Borchert, G.; Lanier, W.; Davidson, B. RNA polymerase III transcribes human microRNAs. Nat. Struct. Mol. Biol.; 2007; 13, pp. 1097-1101. [DOI: https://dx.doi.org/10.1038/nsmb1167] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/17099701]
15. Xie, Z.; Johansen, L.; Gustafson, A.; Kasschau, K.; Lellis, A.; Zilberman, D.; Jacobsen, S.; Carrington, J. Genetic and Functional Diversification of Small RNA Pathways in Plants. PLoS Biol.; 2004; 2, E104. [DOI: https://dx.doi.org/10.1371/journal.pbio.0020104]
16. Fang, Y.; Spector, D. Identification of Nuclear Dicing Bodies Containing Proteins for MicroRNA Biogenesis in Living Arabidopsis Plants. Curr. Biol. CB; 2007; 17, pp. 818-823. [DOI: https://dx.doi.org/10.1016/j.cub.2007.04.005]
17. Li, J.; Yang, Z.; Yu, B.; Liu, J.; Chen, X. Methylation Protects miRNAs and siRNAs from a 3′-End Uridylation Activity in Arabidopsis. Curr. Biol. CB; 2005; 15, pp. 1501-1507. [DOI: https://dx.doi.org/10.1016/j.cub.2005.07.029]
18. Park, W.; Li, J.; Song, R.; Messing, J.; Chen, X. CARPEL FACTORY, a dicer homolog, and HEN1, a novel protein, act in microRNA metabolism in Arabidopsis thaliana. Curr. Biol. CB; 2002; 12, pp. 1484-1495. [DOI: https://dx.doi.org/10.1016/S0960-9822(02)01017-5]
19. Bollman, K.; Aukerman, M.; Park, M.Y.; Hunter, C.; Berardini, T.; Poethig, R. HASTY, the Arabidopsis ortholog of exportin 5/MSN5, regulates phase change and morphogenesis. Development; 2003; 130, pp. 1493-1504. [DOI: https://dx.doi.org/10.1242/dev.00362]
20. Bajczyk, M.; Lange, H.; Bielewicz, D.; Szewc, L.; Bhat, S.; Dolata, J.; Kuhn, L.; Gagliardi, D.; Jarmolowski, A. SERRATE interacts with the nuclear exosome targeting (NEXT) complex to degrade primary miRNA precursors in Arabidopsis. Nucleic Acids Res.; 2020; 48, pp. 6839-6854. [DOI: https://dx.doi.org/10.1093/nar/gkaa373]
21. Xie, D.; Chen, M.; Niu, J.; Wang, L.; Li, Y.; Fang, X.; Li, P.; Qi, Y. Phase separation of SERRATE drives dicing body assembly and promotes miRNA processing in Arabidopsis. Nat. Cell Biol.; 2021; 23, pp. 1-8. [DOI: https://dx.doi.org/10.1038/s41556-020-00606-5]
22. Lee, S.; Martienssen, R. Phase separation in plant miRNA processing. Nat. Cell Biol.; 2020; 23, pp. 5-6. [DOI: https://dx.doi.org/10.1038/s41556-020-00614-5] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/33318660]
23. Gao, S.; Wang, J.; Jiang, N.; Zhang, S.; Wang, Y.; Zhang, J.; Li, N.; Fang, Y.; Yang, L.; Chen, S. et al. Hyponastic Leaves 1 protects pri-miRNAs from nuclear exosome attack. Proc. Natl. Acad. Sci. USA; 2020; 117, pp. 17429-17437. [DOI: https://dx.doi.org/10.1073/pnas.2007203117] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32636270]
24. Wang, Z.; Ma, Z.; Castillo-Gonzalez, C.; Sun, D.; Li, Y.; Yu, B.; Zhao, B.; Li, P.; Zhang, X. SWI2/SNF2 ATPase CHR2 remodels pri-miRNAs via Serrate to impede miRNA production. Nature; 2018; 557, pp. 516-521. [DOI: https://dx.doi.org/10.1038/s41586-018-0135-x] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/29769717]
25. Zhang, B.; You, C.; Zhang, Y.; Zeng, L.; Hu, J.; Zhao, M.; Chen, X. Linking key steps of microRNA biogenesis by TREX-2 and the nuclear pore complex in Arabidopsis. Nat. Plants; 2020; 6, pp. 957-969. [DOI: https://dx.doi.org/10.1038/s41477-020-0726-z] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32690891]
26. Cambiagno, D.; Giudicatti, A.; Arce, A.; Gagliardi, D.; Lei, L.; Yuan, W.; Lundberg, D.; Weigel, D.; Pablo, M. HASTY modulates miRNA biogenesis by linking pri-miRNA transcription and processing. Mol. Plant; 2020; 14, pp. 426-439. [DOI: https://dx.doi.org/10.1016/j.molp.2020.12.019]
27. Li, S.; Li, M.; Liu, K.; Zhang, H.; Zhang, S.; Zhang, C.; Yu, B. MAC5, an RNA-binding protein, protects pri-miRNAs from SERRATE-dependent exoribonuclease activities. Proc. Natl. Acad. Sci. USA; 2020; 117, pp. 23982-23990. [DOI: https://dx.doi.org/10.1073/pnas.2008283117]
28. Li, S.; Xu, R.; Li, A.; Liu, K.; Gu, L.; Li, M.; Zhang, H.; Zhang, Y.; Zhuang, S.; Wang, Q. et al. SMA1, a Homolog of the Splicing Factor Prp28, Has a Multifaceted Role in miRNA Biogenesis in Arabidopsis. Nucleic Acids Res.; 2018; 46, pp. 9148-9159. [DOI: https://dx.doi.org/10.1093/nar/gky591]
29. Wang, S.; Quan, L.; Li, S.; You, C.; Zhang, Y.; Gao, l.; Zeng, L.; Liu, L.; Qi, Y.; Mo, B. et al. The Protein Phosphatase 4 Complex Promotes Transcription and Processing of Primary microRNAs in Arabidopsis. Plant Cell; 2019; 31, pp. 486-501. [DOI: https://dx.doi.org/10.1105/tpc.18.00556]
30. Khvorova, A.; Reynolds, A.; Jayasena, S. Functional siRNAs and miRNAs exhibit strand bias. Cell; 2003; 115, pp. 209-216. [DOI: https://dx.doi.org/10.1016/S0092-8674(03)00801-8]
31. Schwarz, D.; Hutvagner, G.; Du, T.; Xu, Z.; Aronin, N.; Zamore, P. RISC—Asymmetry in the assembly of the RNAi enzyme complex. Cell; 2003; 115, pp. 199-208. [DOI: https://dx.doi.org/10.1016/S0092-8674(03)00759-1]
32. Gyorgy, H.; Simard, M.J. Argonaute proteins: Key players in RNA silencing. Mol. Cell Biol.; 2008; 9, pp. 22-32.
33. Liu, S.; Xu, M. Identification and characterization of long non-coding RNAs involved in the formation and development of poplar adventitious roots. Ind. Crops Prod.; 2018; 118, pp. 334-346. [DOI: https://dx.doi.org/10.1016/j.indcrop.2018.03.071]
34. Tang, F.; Wei, H.; Zhao, S.; Wang, L.; Zheng, H.; Lu, M. Identification of microRNAs Involved in Regeneration of the Secondary Vascular System in Populus tomentosa Carr. Front. Plant Sci.; 2016; 7, 724. [DOI: https://dx.doi.org/10.3389/fpls.2016.00724] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/27303419]
35. Puzey, J.R.; Karger, A.; Axtell, M.; Kramer, E.M. Deep Annotation of Populus trichocarpa microRNAs from Diverse Tissue Sets. PLoS ONE; 2012; 7, e33034. [DOI: https://dx.doi.org/10.1371/journal.pone.0033034]
36. Lu, S.; Li, Q.; Wei, H.; Chang, M.-J.; Tunlaya-Anukit, S.; Kim, H.; Liu, J.; Song, J.; Sun, Y.-H.; Yuan, L. et al. Ptr-miR397a is a negative regulator of laccase genes affecting lignin content in Populus trichocarpa. Proc. Natl. Acad. Sci. USA; 2013; 110, pp. 10848-10853. [DOI: https://dx.doi.org/10.1073/pnas.1308936110]
37. Wang, C.-Y.; Zhang, S.; Yu, Y.; Luo, Y.-C.; Liu, Q.; Ju, C.; Zhang, Y.-C.; Qu, L.-H.; Lucas, W.J.; Wang, X. et al. MiR397b regulates both lignin content and seed number in Arabidopsis via modulating a laccase involved in lignin biosynthesis. Plant Biotechnol. J.; 2014; 12, pp. 1132-1142. [DOI: https://dx.doi.org/10.1111/pbi.12222]
38. Sun, Q.; Liu, X.; Yang, J.; Liu, W.; Du, Q.; Wang, H.; Fu, C.; Li, W.-X. MicroRNA528 Affects Lodging Resistance of Maize by Regulating Lignin Biosynthesis under Nitrogen-Luxury Conditions. Mol. Plant; 2018; 11, pp. 806-814. [DOI: https://dx.doi.org/10.1016/j.molp.2018.03.013]
39. Shen, T.; Xu, M.; Qi, H.; Feng, Y.; Yang, Z.; Xu, M. Uncovering miRNA-mRNA Regulatory Modules in Developing Xylem of Pinus massoniana via Small RNA and Degradome Sequencing. Int. J. Mol. Sci.; 2021; 22, 10154. [DOI: https://dx.doi.org/10.3390/ijms221810154]
40. Reyes, J.L.; Chua, N.-H. ABA induction of miR159 controls transcript levels of two MYB factors during Arabidopsis seed germination. Plant J.; 2007; 49, pp. 592-606. [DOI: https://dx.doi.org/10.1111/j.1365-313X.2006.02980.x]
41. Baldoni, E.; Genga, A.; Cominelli, E. Plant MYB Transcription Factors: Their Role in Drought Response Mechanisms. Int. J. Mol. Sci.; 2015; 16, pp. 15811-15851. [DOI: https://dx.doi.org/10.3390/ijms160715811]
42. Wu, X.; Ding, D.; Shi, C.; Xue, Y.; Zhang, Z.; Tang, G.; Tang, J. microRNA-dependent gene regulatory networks in maize leaf senescence. BMC Plant Biol.; 2016; 16, 73. [DOI: https://dx.doi.org/10.1186/s12870-016-0755-y] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/27000050]
43. Xie, F.; Wang, Q.; Sun, R.; Zhang, B. Deep sequencing reveals important roles of microRNAs in response to drought and salinity stress in cotton. J. Exp. Bot.; 2015; 66, pp. 789-804. [DOI: https://dx.doi.org/10.1093/jxb/eru437] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/25371507]
44. Feng, H.; Duan, X.; Zhang, Q.; Li, X.; Wang, B.; Huang, L.; Wang, X.; Kang, Z. The target gene of tae-miR164, a novel NAC transcription factor from the NAM subfamily, negatively regulates resistance of wheat to stripe rust. Mol. Plant Pathol.; 2014; 15, pp. 284-296. [DOI: https://dx.doi.org/10.1111/mpp.12089] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/24128392]
45. Sakaguchi, J.; Watanabe, Y. miR165/166 and the development of land plants. Dev. Growth Differ.; 2012; 54, pp. 93-99. [DOI: https://dx.doi.org/10.1111/j.1440-169X.2011.01318.x] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/23741750]
46. Axtell, M.; Meyers, B. Revisiting criteria for plant miRNA annotation in the era of big data. Plant Cell; 2018; 30, pp. 272-284. [DOI: https://dx.doi.org/10.1105/tpc.17.00851] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/29343505]
47. Taylor, R.; Tarver, J.; Hiscock, S.; Donoghue, P. Evolutionary history of plant microRNAs. Trends Plant Sci.; 2014; 19, pp. 175-182. [DOI: https://dx.doi.org/10.1016/j.tplants.2013.11.008] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/24405820]
48. Jangra, S.; Chaudhary, V.; Yadav, N. Transcription Factors and MicroRNA Interplay: A New Strategy for Crop Improvement; Books on Demand: Norderstedt, Germany, 2018; pp. 103-123.
49. Mullany, L.; Herrick, J.; Wolff, R.; Stevens, J.; Samowitz, W.; Slattery, M. MicroRNA-Transcription Factor Interactions and Their Combined Effect on Target Gene Expression in Colon Cancer Cases. Genes Chromosom. Cancer; 2017; 57, pp. 192-202. [DOI: https://dx.doi.org/10.1002/gcc.22520]
50. Wu, G.; Park, M.Y.; Conway, S.; Wang, J.-W.; Weigel, D.; Poethig, R. The Sequential Action of miR156 and miR172 Regulates Developmental Timing in Arabidopsis. Cell; 2009; 138, pp. 750-759. [DOI: https://dx.doi.org/10.1016/j.cell.2009.06.031]
51. Yan, Z.; Hossain, M.S.; Wang, J.; Valdés-López, O.; Liang, Y.; Libault, M.; Qiu, L.; Stacey, G. miR172 Regulates Soybean Nodulation. Mol. Plant-Microbe Interact. MPMI; 2013; 26, pp. 1371-1377. [DOI: https://dx.doi.org/10.1094/MPMI-04-13-0111-R]
52. Nova-Franco, B.; Iñiguez, L.; Valdés-López, O.; Alvarado-Affantranger, X.; Leija Salas, A.; Fuentes, S.; Ramirez, M.; Paul, S.; Reyes, J.; Girard, L. et al. The miR172c-AP2-1 Node as a Key Regulator of the Common Bean—Rhizobia Nitrogen Fixation Symbiosis. Plant Physiol.; 2015; 168, pp. 273-291. [DOI: https://dx.doi.org/10.1104/pp.114.255547]
53. Wang, J.-W.; Wang, L.-J.; Mao, Y.-B.; Cai, W.-J.; Xue, H.-W.; Chen, X.-Y. Control of Root Cap Formation by MicroRNA-Targeted Auxin Response Factors in Arabidopsis. Plant Cell; 2005; 17, pp. 2204-2216. [DOI: https://dx.doi.org/10.1105/tpc.105.033076] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/16006581]
54. Roy, S. Function of MYB domain transcription factors in abiotic stress and epigenetic control of stress response in plant genome. Plant Signal. Behav.; 2015; 11, e1117723. [DOI: https://dx.doi.org/10.1080/15592324.2015.1117723] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/26636625]
55. Efroni, I.; Blum, E.; Goldshmidt, A.; Eshed, Y. A Protracted and Dynamic Maturation Schedule Underlies Arabidopsis Leaf Development. Plant Cell; 2008; 20, pp. 2293-2306. [DOI: https://dx.doi.org/10.1105/tpc.107.057521]
56. Khan, S.; Ali, A.; Saifi, M.; Saxena, P.; Ahlawat, S.; Abdin, M. Identification and the potential involvement of miRNAs in the regulation of artemisinin biosynthesis in A. annua. Sci. Rep.; 2020; 10, pp. 1-11. [DOI: https://dx.doi.org/10.1038/s41598-020-69707-3] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32788629]
57. Liao, Y.; Zhang, Q.; Cui, R.; Xu, X.; Zhu, F.-Y.; Cheng, Q.; Li, X. High-Throughput Sequencing Reveals the Regulatory Networks of Transcriptome and Small RNAs During the Defense Against Marssonina brunnea in Poplar. Front. Plant Sci.; 2021; 12, 719549. [DOI: https://dx.doi.org/10.3389/fpls.2021.719549] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34567031]
58. Le Roy, J.; Blervacq, A.-S.; Creach, A.; Huss, B.; Hawkins, S.; Neutelings, G. Spatial regulation of monolignol biosynthesis and laccase genes control developmental and stress-related lignin in flax. BMC Plant Biol.; 2017; 17, pp. 1-20. [DOI: https://dx.doi.org/10.1186/s12870-017-1072-9]
59. Xue, C.; Yao, J.L.; Qin, M.-F.; Zhang, M.; Allan, A.; Wang, D.-F.; Wu, J. PbrmiR397a regulates lignification during stone cell development in pear fruit. Plant Biotechnol. J.; 2018; 17, pp. 103-117. [DOI: https://dx.doi.org/10.1111/pbi.12950]
60. Cao, J.; Gulyas, Z.; Kalapos, B.; Boldizsar, A.; Liu, X.; Pal, M.; Yao, Y.; Galiba, G.; Kocsy, G. Identification of a redox-dependent regulatory network of miRNAs and their targets in wheat. J. Exp. Bot.; 2019; 70, pp. 85-99. [DOI: https://dx.doi.org/10.1093/jxb/ery339]
61. Stein, O.; Granot, D. An Overview of Sucrose Synthases in Plants. Front. Plant Sci.; 2019; 10, 95. [DOI: https://dx.doi.org/10.3389/fpls.2019.00095]
62. Doblin, M.; Kurek, I.; Jacob-Wilk, D.; Delmer, D. Cellulose biosynthesis in plants: From genes to rosettes. Plant Cell Physiol.; 2003; 43, pp. 1407-1420. [DOI: https://dx.doi.org/10.1093/pcp/pcf164]
63. Turner, S.; Kumar, M. Cellulose synthase complex organization and cellulose microfibril structure. Philos. Trans. A Math. Phys. Eng. Sci.; 2018; 376, 20170048. [DOI: https://dx.doi.org/10.1098/rsta.2017.0048] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/29277745]
64. Tiwari, B.; Habermann, K.; Arif, M.A.; Top, O. Identification of Small RNAs During High Light Acclimation in Arabidopsis thaliana. Front. Plant Sci.; 2021; 12, 1179. [DOI: https://dx.doi.org/10.3389/fpls.2021.656657] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34211484]
65. Langenkamper, G.; Fung, R.W.M.; Newcomb, R.; Atkinson, R.; Gardner, R.C.; Macrae, E. Sucrose phosphate synthase genes in plants belong to three different families. J. Mol. Evol.; 2002; 54, pp. 322-332. [DOI: https://dx.doi.org/10.1007/s00239-001-0047-4] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/11847558]
66. Chua, T.; Bujnicki, J.; Tan, T.-C.; Huynh, F.; Patel, B.; Sivaraman, J. The Structure of Sucrose Phosphate Synthase from Halothermothrix orenii Reveals Its Mechanism of Action and Binding Mode. Plant Cell; 2008; 20, pp. 1059-1072. [DOI: https://dx.doi.org/10.1105/tpc.107.051193] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/18424616]
67. Gyawali, B.; Barozai, M.Y.; Aziz, A. Comparative expression analysis of microRNAs and their targets in emerging bio-fuel crop sweet sorghum (Sorghum bicolor L.). Plant Gene; 2021; 26, 100274. [DOI: https://dx.doi.org/10.1016/j.plgene.2021.100274]
68. Ohashi-Ito, K.; Kubo, M.; Demura, T.; Fukuda, H. Class III Homeodomain Leucine-Zipper Proteins Regulate Xylem Cell Differentiation. Plant Cell Physiol.; 2005; 46, pp. 1646-1656. [DOI: https://dx.doi.org/10.1093/pcp/pci180]
69. Wang, R.; Reng, M.; Tian, S.; Liu, C.; Cheng, H.; Liu, Y.; Zhang, H.; Wei, H.; Wei, Z. Transcriptome-wide identification and characterization of microRNAs in diverse phases of wood formation in Populus trichocarpa. G3-Genes Genomes Genet.; 2021; 11, jkab195. [DOI: https://dx.doi.org/10.1093/g3journal/jkab195]
70. Andrews, S. FastQC A Quality Control Tool for High Throughput Sequence Data. Available online: http://www.bioinformatics.babraham.ac.uk/projects/fastqc/ (accessed on 8 January 2020).
71. Martin, M. CUTADAPT removes adapter sequences from high-throughput sequencing reads. EMBnet.J.; 2011; 17, pp. 10-12. [DOI: https://dx.doi.org/10.14806/ej.17.1.200]
72. Tsuji, J.; Deng, Z. DNApi: A De Novo Adapter Prediction Algorithm for Small RNA Sequencing Data. PLoS ONE; 2016; 11, e0164228. [DOI: https://dx.doi.org/10.1371/journal.pone.0164228]
73. Veneziano, D.; Di Bella, S.; Giovanni, N.; Laganà, A.; Ferro, A.; Croce, C. Non-coding RNA: Current Deep Sequencing Data Analysis Approaches and Challenges. Hum. Mutat.; 2016; 37, pp. 1283-1298. [DOI: https://dx.doi.org/10.1002/humu.23066]
74. Langmead, B.; Trapnell, C.; Pop, M.; Salzberg, S.L. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol.; 2009; 10, R25. [DOI: https://dx.doi.org/10.1186/gb-2009-10-3-r25] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/19261174]
75. Chen, C.; Li, J.; Feng, J.; Bo, L.; Feng, L.; Yu, X.; Li, G.; Zhai, J.; Meyers, B.; Xia, R. sRNAanno—a database repository of uniformly annotated small RNAs in plants. Hortic. Res.; 2021; 8, pp. 1-8. [DOI: https://dx.doi.org/10.1038/s41438-021-00480-8] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/33642576]
76. Feng, L.; Xia, R.; Liu, Y. Comprehensive Characterization of miRNA and PHAS Loci in the Diploid Strawberry (Fragaria vesca) Genome. Hortic. Plant J.; 2019; 5, pp. 255-267. [DOI: https://dx.doi.org/10.1016/j.hpj.2019.11.004]
77. Love, M.; Huber, W.; Anders, S. Moderated estimation of fold change and dispersion for RNA-Seq data with DESeq2. Genome Biol.; 2014; 15, 550. [DOI: https://dx.doi.org/10.1186/s13059-014-0550-8] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/25516281]
78. Goodstein, D.; Shu, S.; Howson, R.; Neupane, R.; Hayes, R.; Fazo, J.; Mitros, T.; Dirks, W.; Hellsten, U.; Putnam, N. et al. Phytozome: A comparative platform for green plant genomics. Nucleic Acids Res.; 2011; 40, pp. D1178-D1186. [DOI: https://dx.doi.org/10.1093/nar/gkr944] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/22110026]
79. Chen, S.; Zhou, Y.; Chen, Y.; Gu, J. fastp: An ultra-fast all-in-one FASTQ preprocessor. Bioinformatics; 2018; 34, pp. i884-i890. [DOI: https://dx.doi.org/10.1093/bioinformatics/bty560]
© 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/). Notwithstanding the ProQuest Terms and Conditions, you may use this content in accordance with the terms of the License.
Details
1 State Key Laboratory of Tree Genetics and Breeding, Key Laboratory of Tree Breeding and Cultivation of State Forestry Administration, Research Institute of Forestry, Chinese Academy of Forestry, Beijing 100091, China;
2 Co-Innovation Center for Sustainable Forestry in Southern China, Key Laboratory of Forest Genetics and Biotechnology of Ministry of Education, Nanjing Forestry University, Nanjing 210037, China;