Differences in growth and wood properties between Pd and Ps
Wood is the secondary xylem of trees, mainly composed of cellulose, hemicellulose and lignin. All xylem cell types first undergo secondary cell wall (SCW) thickening and undergo programmed cell death [22, 23]. The plant height and ground diameter of annual Pd were significantly higher than that of Ps (Fig. 1a and b). Further material property determination found that the content of the three major elements of Pd is higher (Fig. 1c-i). SAM is an important regulatory site of plant height growth, phloem and developing xylem are important parts of plant radial growth. The anatomical structure showed that Pd had obvious shoot apex growth point (central zone), while Ps is not obvious (Fig. 1j). In addition, Pd exhibited wider phloem fiber and cambium, and larger vessel area and proportion, which may be beneficial to transports nutrient (Fig. 1k). The small basic density of Pd may be related to the large area occupied by the vessels (Fig. 1l).
Identification of lncRNAs from SAM, phloem, and developing xylem RNA-seq datasets
As an important fast-growing tree species, it is very important to understand the molecular pathways of growth and development of poplar. After trimming adapters and removing low quality and contaminated reads, in total, 246.24 Gb clean data was obtained from 18 libraries, with an average Q30 of 93.00% (Additional file 6: Table S1). Finally, we identified 6, 355 lncRNAs,with protections of 2,455 sense_overlapping lncRNAs, 2,005 of lincRNA, 1897 of antisense lncRNAs (Additional file 6: Table S2).
In order to characterize the characteristics of these lncRNAs, we evaluated the distribution of chromosome location, transcript length, exon number and expression level of lncRNAs. In general, lincRNA, sense_overlapping and antisense LncRNAs were evenly distributed on 19 chromosomes, although they had different emphases (Fig. 2a).The average length of lncRNAs was 990 bp and about 63.4% contained two exons (Fig. 2b and c). Antisense lncRNAs ranged in length from 201 to 9830 bp and average was 940 bp. LincRNAs ranged between 201and 4734 bp (average=783 bp). And average length of sense_overlapping lncRNAs is 1197 bp. The GC content of antisense lncRNA was 41.36%, which was significantly higher than that of lincRNAs and sense_overlapping lncRNAs (Fig. 2d). For expression levels, the lncRNAs expression levels were different and showed fewer average counts (FPKM=4.44) than the coding transcripts (FPKM=18.96) (Fig. 2e and f).
A principal component analysis (PCA) plot of the whole data set revealed a sequential order of the different samples. The results showed that SAM, xylem and phloem of the two species were clustered into one group, and the similarity between the two tissues was greater than that between genotypes (Fig. 3a).
Differentially expressed anaylsis between Pd and Ps
In order to further analyze whether these genes were differentially expressed between the two genotypes and different tissues, nine comparative combinations were carried out. Finally, 3, 573 differentially expressed (DE) lncRNAs and 27, 582 DEmRNAs were obtained. among them, the DElncRNAs of Pd_S vs. Ps_S is 1957 (Fig. 3b and c, Additional file 1:Figure S1). DEmRNAs participate in molecular functions such as “ADP binding” and “catalytic activity” (Additional file 1:Figure S1C).There are 980 DElncRNAs in two genotypes of the corresponding tissue comparison combination (Pd_ X vs. Ps_ X, Pd_ P vs. Ps_ P and Pd_ S vs. Ps_ S) (Fig. 3d). There are 322, 333 and 420 specifically expressed lncRNAs in SAM, xylem and phloem, respectively. Since lncRNAs play important roles in regulating gene expression, identification and analysis of their target genes may help us to explore their potential functions. Computational prediction identified a set of 12, 875 mRNAs corresponding to 8, 931 lncRNA–target pairs, including 6,565 colocation-regulated mRNAs for 3,551 lncRNAs and 7, 421 coexpression- regulated mRNAs for 3,285 lncRNAs.
To further analyze the function of these lncRNAs, we performed GO and KEGG analysis on their target genes. The colocation target genes of DElncRNAs were mainly enriched in 60 GO terms such as “ADP binding” and “nucleoside binding” (Additional file 2: Figure S2). Some target genes were enriched in photosynthesis pathway, including 39 lncRNAs and 68 mRNAs. And TCONS_00135489 showed the same trend as its target, and the expression level of related genes were high in SAM (Additional file 3: Figure S3).
The coexpression genes of Pd_ X vs. Ps_ X, Pd_ P vs. Ps_ P, and Pd_ S vs. Ps_ S were mainly related to “ADP binding”, “heme binding” and other biological functions (Fig. 4a). And they significantly enriched in the “phenylpropanoid biosynthesis” and “cyanoamino acid metabolism” pathway (Fig. 4b). These possible target genes provide new insight into the role of lncRNAs in poplar development.
Regulation of lncRNAs and transcription factors in the phenylpropanoid biosynthesis
Based on the predicted GO conditions of growth differential lncRNAs and the pathways associated with target genes, we speculated that lncRNAs might play an important role in phenylalanine biosynthesis in poplars. Phenylpropanoids are a group of plant secondary metabolites derived from phenylalanine, which has a variety of structural and signal molecular functions [24]. It is the starting compounds for biosynthesis of lignin, flavonoids, anthocyanins and so on, serves as the core mediators of crosstalk between developmental and defenserelated pathways [25]. Further analysis co-expression network of these lncRNAs and lignin, flavonoid biosynthesis, found that TCONS_00128372, incRNA, located in Chr12 and interact with MYB46, SECONDARY WALL- ASSOCIATED NAC DOMAIN2 (SND2), cinnamate-4-hydroxylase (C4H), caffeoyl-CoA 3-O-methyltransferase (CcoAMT), and laccase (LAC). And sense_overlapping lncRNA TCONS_00079190 co-expressed with MYB83, MYB46, NAC SECONDARY WALL THICKENING PROMOTIBNG FACTOR1 (NST1), and LAC (Fig. 5a). TCONS_00174042, TCONS_00101258, and TCONS_00136338 regulated the structural genes of flavonoid biosynthesis. TCONS_00174042 and TCONS_00101258 co-expressed with MYB3,chalcone and stilbene synthase(CHS), leucoanthocyanidin reductase (LAR), leucoanthocyanidin dioxygenase (LDOX), and naringenin 3- dioxygenase (F3H). In addition, TCONS_00136338 can coexpress with LDOX, CHS, UFGT, LAR, and dihydroflavonol 4-reductase (DFR) (Fig. 5b). Lignin related lncRNAs and structural genes were highly expressed in xylem of Pd and Ps, while flavonoid related genes were highly expressed in SAM (Additional file 4: Figure S4).
LncRNAs involved in the plant hormone biosynthesis signal transduction pathway
Hormones, such as auxin, cytokinin, and gibberellin, play an important role in primay and secondary growth [26, 27]. In order to further analyze the relationship between lncRNA and plant hormones, we constructed a co-expression network including auxin, cytokinin, and gibberellin related PCgenes and their regulatory lncRNAs (Fig. 6a). TCONS_00134627 could be co expressed with GA2OX8 (Potri.011G134000), GASA10 (Potri.009G092600), and SAUR94 (Potri.009G127300). The expression levels of DEGs and DElncRNAs related to auxin, cytokinin and gibberellin biosynthesis were shown in the heatmap (Fig. 6b).
miRNA involved in hormones and phenylalanine pathway
miRNAs, a major class of small RNAs with 20–24 nucleotides, create various aspects of plant development and stress responses through post transcriptionally regulation gene expression [28]. A total of 658 miRNA-lncRNA pairs consisting of 188 miRNAs and 200 lncRNAs were identified, including 19 plant hormone related pairs and 28 phenylalanine related pairs (Additional file 6: Table S3). TCONS_00066905, hormone related lncRNA, was predicted to be a target mimic of miR396a and miR396b. In addition, phenylalanine related TCONS_00023606 and TCONS_ 00093325 is the target mimic of flavonoid related regulatory genes miR156h and miR828a (Fig. 7a). miR396-GRFs is an important regulatory module of plant growth and development. We found that 12 GRFs were differentially expressed and highly expressed in SAM of Pd, which may be an important reason for the rapid growth of Pd (Fig. 7b).
Validation of lncRNA and gene expression by qRT-PCR
To examine whether these differentially expressed lncRNAs had a role in xylem development and hormone signal transduction, related lncRNAs and mRNAs were selected based on their expression levels in RNA-seq. The results proved the reliability of transcriptome data (Additional file 5: Figure S5).