Different segments of stem bark exhibit distinct morphological features
Ramie fibers continuously develop along the stem during plant’s growth, while the internodes of the stem show obvious elongation only until the plant is fully elongated. To analyze the developmental stages of fiber formation and the gene expression profiles, three parts of ramie’s shoot were harvested, including top bud (TB), internode elongating region (ER) of stem and internode fully elongated region (FER) of stem (Fig.1A). The top buds and the barks from both ER and FER regions peeled off from the woody center pillar were used for histological analysis and RNA extraction. The strategy for the subsequent RNA-seq data analysis is illustrated in Fig. 1B. The cross and longitudinal sections of TB, ER and FER were analyzed (Fig.1C and D). In the TB sample, ramie has amphicribral vascular bundle, which is different from flax or hemp plants but is similar to woody plant with continuous cambia within and outside the vascular bundles (Fig. 1C). The vascular structure in TB is characteristic of multiple layers of primary phloem without obvious boundary between vascular bundles. In ER and FER, clear differences were observed between these two regions (Fig. 1D). Firstly, FER has thicker bark than ER, and FER barks consist of more enlarged cells and more layers of phloem tissues. Secondly, fiber cells show thicker cell wall in FER without an increase in cell size; the thickness of the fiber cell wall is about 5.38 µm in FER vs. 1.87 µm in ER (Supplemental table 1). Thirdly, the cell wall of the fibers from FER phloem contains more lignin than that from ER, which is indicated by stronger red color resulting from the staining of safranin dye (Fig. 1b and 1d). The differences among these three samples indicate different developing stages of phloem fiber cells. Therefore, we used these samples for gene expression profiling attempting to identify genes important for fiber development in ramie.
Assembly of de novo transcriptome and identification of unigenes
Thirty-three RNA samples were collected and subjected to the next generation sequencing (NGS), and the RNA-seq data, including the 9 submitted SRA files (SRR9112644-SRR9112651), were analyzed. More than 5G sequences with clean bases from each sample was obtained, and thus the total analyzed clean bases were about 1.7E+11. The genome size of Zhongzhu No. 1 is approximately 340 Mb [3, 4]. Therefore, the depth of the RNA-seq data used in this study was expected to be enough for a high quality de novo assembly of transcriptome for the expressed genes from the top bud and stem bark tissues. The 10 species with the most matching reads to our RNA-seq data were shown in Fig. 2A. Among all the reads generated, 3048 reads match with those in Boehmeria nivea, and the highest matching ratio (28%) was found to be with Morus notabils. The whole genome sequencing analysis supported a closest evolutionary relationship between ramie and Morus notabils [4]. Overall, there were 59486 unigenes assembled with the length longer than 300 bp, 47016 unigenes longer than 500 bp, and 31395 unigenes longer than 1000 bp. The GC content distribution of all unigenes was shown in Fig. 2B, and two peaks appeared between the range of 30% and 45%. The detailed size distribution of all unigenes was illustrated in Fig. 2C. The sequence of each unigene was subsequently processed by blast to NR, SWISSPROT and KOG databases, respectively, and the annotations were obtained according to the most similar protein or gene with e<1e-5.
Identification of differentially expressed genes (DEGs) and expression patterns among TB, ER and FER
DEGs among the three tissues were identified following the scheme shown in Fig. 1B. When compared with TB, there were 4138 unigenes up-regulated and 6638 unigenes down-regulated in the ER, and 3853 unigenes up-regulated and 5075 unigenes down-regulated in the FER (Fig. 3A). The VEN diagram showed that the DEGs among these 3 samples were grouped in 6 distinct clusters (Fig. 3B). The heatmaps of the expression of these clustered genes were shown in Fig. 4A, and the schematic map of the expression patterns and GO analysis were illustrated in Fig. 4B.
The cluster 1 and 2 contain the most DEGs with 4354 up- and 2046 down-regulated unigenes only in TB (Figure 4A and 4B). These two clusters may contain tissue-specific DEGs resulting in the differences between the top bud and stem bark samples. The cluster 1 DEGs consist of the unigenes with higher expression level in TB but lower expression level in both bark regions. GO analysis showed that these DEGs are involved in meiotic chromosome segregation and cell division, which could be explained by the active gene expression required for vigorous cell division in the SAM region of the top bud. Other DEGs involved in stomatal or leaf development in the cluster 1 could be due to the remaining young emerging leaves in the TB samples, while no leaves were present in the ER and FER samples (Fig.4B). GO analysis of the DEGs in cluster 2 showed that up-regulated transcription factors or transcription processes and the plant-type secondary cell wall biogenesis are among the top categories. Fifty-five unigene contigs of cell wall components or cell wall biogenesis involving factors were identified in the cluster 2 (Table S2). These factors include Cellulose Synthase A Catalytic Subunit 3 and 8 (CesA 3 and 8), Fasciclin-like Arabinogalactan Protein (FLA), beta-galactosidase (BGAL), several pectinesterase/pectinesterase inhibitors (PEMs) and the enzymes for the synthesis of other cell wall components, such as glucuronoxylan glucuronosyltransferase, galacturonosyltransferase, endochitinase, callose synthase, xyloglucan glycosyltransferase, xyloglucan endotransglucosylase, etc. (Tab. S2).
The cluster 3 and 4 show the unigenes up- or down-regulated only in FER (Fig. 4A and 4B). In these two clusters, there were 93 unigenes down-regulated and 476 unigenes up-regulated only in the FER. In cluster 3, a small amount of unigenes for membrane construction were down-regulated in FER. In cluster 4, relatively more unigenes were up-regulated in FER comparing with the down-regulated unigenes in cluster 3. Among these up-regulated unigenes, ethylene signaling pathway genes were the most enriched unigenes. There were totally 39 transcription factors up-regulated only in FER, and 18 out of the 39 were ethylene activating unigenes (Tab. S3). The DEGs only in the ER were clustered in cluster 5 and 6. Interestingly, some phloem development related unigenes were found to be down-regulated only in ER when compared with those in both TB and CER (Fig.4B). According to the GO analysis of DEGs in cluster 4 and 5, phytohormons such as ethylene and gibberellin might play distinct roles in these two bark stages, as some of the gibberellin-responsive genes were only up-regulated in ER while some of the ethylene signaling pathway related genes were only activated in FER. Furthermore, higher expression of several DEGs (FER vs. ER, cluster 6) involved in the phloem development suggested a more vigorous secondary phloem development in the relatively mature stem or bark from FER.
In addition to the expression patterns analyzed among TB, ER and FER, DEGs between TB and ER or FER were also analyzed and GO analyses were performed. The top 10 items of three GO terms were shown in Fig. S1 and S2. When compared with TB sample, barks from FER showed distinct features in gene expression patterns comparing with those from ER. Therefore, DEG identification between ER and FER is essential to uncover the differences between the barks in different fiber developmental stages.
GO analysis of DEGs between ER and FER
There were 1628 up-regulated unigenes and 757 down-regulated unigenes identified in ramie’s bark of FER when compared with ER (Fig.5). GO analysis shown in Fig. 6 revealed the top 10 up-regulated biological processes including phloem development, response to chitin, ethylene-activated signaling pathway, DNA replication, salicylic acid mediated signaling, defense response, protein transmembrane transport and vasculature development, and the top 10 down-regulated biological processes including cytoplasmic translation, tricarboxylic acid cycle, indole glucosinolate metabolic process, plant-type secondary cell wall biogenesis, etc.. Processes such as phloem development, vasculature development and DNA replication were up-regulated in the FER, which suggests that the secondary phloem formation is activated in the bark of the relative mature FER rather than in ER. The activation of ethylene signaling pathway was evidenced by the up-regulation of the genes in this pathway in FER (Tab. S3). Overall 21 unigenes or contigs of 14 Ethylene Respond Factors (ERFs) were up-regulated in FER, which include ERF1, ERF1A, ERF1B, ERF2, ERF3, ERF5, ERF17, ERF22, ERF53, ERF61, ERF71, ERF109, PAR2-13 and PAP2-4 (Tab. S5).
KEGG analysis of DEGs between ER and FER
The KEGG analysis of total DEGs from FER vs. ER revealed additional information to the GO analysis. The KEGG analysis indicated that these DEGs are involved in the pathways of starch and sucrose metabolism, citrate cycle, nitrogen metabolism, cysteine and methionine metabolism, ribosome, diterpenoid biosynthesis, phenylpropanoid biosynthesis, DNA replication, cell cycle, etc. (Fig. 7 and Tab. S7).
The secondary cell wall synthesis is a very complex biological process involving the biosynthesis of multiple components or species-specific secondary metabolisms. However, starch and sucrose metabolisms are the important pathways linking with the secondary cell wall synthesis. From the KEGG analysis, we found that the expression of 23 unigenes encoding 11 enzymes in the starch and sucrose metabolisms differed between ER and FER. These enzymes include sucrose synthase (EC2.4.1.13), sucrose-phosphate synthase (EC2.4.1.14), bata-amylase EC3.2.1.2, endoglucanase (EC3.2.1.4), bata-glucosidase (EC3.2.1.21), glucan endo-1, 3-beta-glucosidase (EC3.2.1.39), glucose-6-phosphate isomerase (EC5.3.1.9), phosphoglucomutase (EC5.4.2.2), UTP-glucose-1-phosphate uridylyltransferase (EC2.7.7.9), trehalose phosphatase (EC3.1.3.12) and trehalase (EC3.2.1.28) (Fig. 8). Most of these enzyme-encoding unigenes were up-regulated in FER, which suggests that multiple pathways for free D-glucose production might be enhanced in FER. In addition, other sugar producing processes such as sucrose-6P, maltose and dextrin might also be enhanced in FER. The increase in these sugar precursors could be important in providing building materials for the secondary cell wall biogenesis in ramie.
More lignin accumulation in FER was observed by the staining with safranin dye (Fig. 1d), which might be due to the up-regulation of enzymes responsible for phenylpropanoid biosynthesis in this region according to our KEGG analysis (Fig.7 and 9). In the pheylpropanoid biosynthesis pathway (KO00940), 16 up-regulated unigenes encode enzymes including peroxidase (EC1.11.1.7), trans-cinnamate 4-monooxygenase (EC1.14.13.11), anthranilate N-methyltransferase (EC2.1.1.68), flavonoid 3',5'-methyltransferase (EC2.1.1.104), beta-glucosidase (EC3.2.1.21), caffeylshikimate esterase and vinorine synthase, some of which were among the DEGs contributing to secondary cell wall synthesis (Fig. 9 and Tab. S6). Most enzyme-encoding genes involved in the biosynthesis of lignin were up-regulated in FER, supporting the observation of more lignin accumulation in FER. Similarly, up-regulation of several PEMs may account for differential accumulation of pectin in FER compared with ER (Tab. S6).
Interestingly, in the diterpenoid biosynthesis pathway, unigenes encoding enzymes such as Ent-kaurenoic acid oxidase 2 (EC1.14.13.79) and gibberellins 20 oxidase (EC1.14.11.12) for converting the precursors to active GA isoforms were up-regulated, while the transcript level of the enzyme gibberellins 2-beta-dioxygenese 8 (EC1.14.11.13) for inactivation of GAs was decreased in FER (Fig.10). These results suggest that a higher concentration of active GA might be needed in FER than in ER, and GA could be another phytohormone regulating the development of bast fiber in ramie.