Comparative analysis of latex MW and rubber particle size
The samples collected from 1-year-old virgin trees and 30-year-old regularly tapped trees of the clones Reyan7-33-97 and RRIM600 were designated Reyan73397y1, RRIM600y1, Reyan73397y30 and RRIM600y30, respectively. The number-average molecular weight (MN), MW and MWD (calculated as MW/MN) of the latex samples were determined by GPC, as shown in Table 1 and Fig. 1. As the age of Hevea trees increased from 1 year to 30 years, the MW of Reyan7-33-97 latex increased from 4.77×105 to 6.51×105, and the MWD also increased from 5.53 to 7.19, while the latex from RRIM600 trees showed that the MW and MN dramatically increase from 4.26×105 to 8.28×105 and 7.18×104 to 17.73×104, respectively, with a slight decrease in MWD (from 5.94 to 4.67). In trees of the same age, the Reyan73397y1 latex displayed a slightly higher MW than that of RRIM600y1, whereas the RRIM600y30 latex had a significantly higher MW than that of Reyan73397y30 in the same trial (Table 1). As shown in Fig. 1A, these latex samples have a MW range from 2.0×103 to 5.0×106 and comprise a high molecular weight fraction (> 106) and a low molecular weight fraction (< 105). Latex from Reyan7-33-97 trees aged 1 and 30 years showed a skewed unimodal MW distribution with a shoulder in the high molecular weight region, while the latex obtained from RRIM600 trees showed a typical bimodal distribution (Fig. 1A). With increasing tree age, the latex from both Reyan7-33-97 and RRIM600 trees showed a dramatic increase in the high molecular weight fraction. In 1-year-old trees, latex from Reyan73397y1 had more high MW rubbers than that from RRIM600y1. In 30-year-old trees, the latex of Reyan73397y30 comprises more low MW rubbers (< 2.0×104), and the position of the high MW peak is smaller than that of RRIM600y30 latex (Fig. 1A).
Table 1. The number-average molecular weight (MN), weight-average molecular weight (MW) and rubber particle size of latex from 1-year virgin trees and 30-year-old regularly tapped trees of Reyan7-33-97 and RRIM600 clones.
Latex Samples
|
MN (×104)
|
MW (×105)
|
MWD (MW/MN)
|
Geometric Mean Diameter (μm)
|
Reyan73397y1
|
8.62 ± 1.84
|
4.77 ± 0.65
|
5.53
|
0.72
|
RRIM600y1
|
7.18 ± 0.39
|
4.26 ± 0.22
|
5.94
|
0.69
|
Reyan73397y30
|
9.05 ± 2.16
|
6.51 ± 0.88
|
7.19
|
0.85
|
RRIM600y30
|
17.73 ± 4.78
|
8.28 ± 0.58
|
4.67
|
0.79
|
Fresh latex samples were analyzed by laser light scattering (LLS) to determine the particle size, a factor potentially affecting the MW of NR (Table 1 and Fig. 1B). The results showed that the geometric mean diameter of RPs varies with age and clone. Latex from 1-year-old trees showed a particle size distribution of predominately less than 1 μm, while those of 30-year-old trees contained particles having a larger size (predominately in between 0.5 μm and 2 μm). On average, RP sizes are smaller in latex from RRIM600 trees than in Reyan7-33-97 trees of the same age. Transmission electron microscopy (TEM) images of ultrathin cross sections showed that the rubber particles within bark laticifers are predominantly spherical in shape. As shown in Figs. 1 C-F, the size of the RPs correlates well with the LLS measurements: the latex from 30-year-old trees has larger RPs than that from 1-year-old trees, and the RPs in RRIM600 latex tend to be smaller than those in Reyan7-33-97 latex. In both laticifers of Reyan7-33-97 at ages 1 and 30, we observed that some elongated particles appeared to result from aggregation, while in the analyzed sections of RRIM600 at ages 1 and 30, the rubber particles are uniformly dispersed.
Transcriptome sequencing and analysis of differentially expressed genes (DEGs)
To understand the transcriptional variations underlying the clonal difference in rubber MW, we sequenced and analyzed the transcriptome files of these latex samples. Total RNA was isolated and sequenced using the Illumina paired-end sequencing approach. The adaptors, low-quality reads and contaminated reads were filtered out, and an average number of 6.45 Gb, 6.60 Gb, 6.35 Gb, and 6.39 Gb clean bases were generated for Reyan73397y1, RRIM600y1, Reyan73397y30 and RRIM600y30, respectively. These clean reads were then mapped to the Hevea reference genome, with mapping rates ranging from 95.53% - 96.51% (Table S1). The gene expression patterns for each latex sample were calculated as fragments per kilobase of transcript per million fragments mapped (FPKM). A total of 34907 genes with detectable expression levels were obtained from four latex samples. To analyze the sequencing results, the comparisons of RRIM600y1 vs. Reyan73397y1 (6y1vs7y1) and RRIM600y30 vs. Reyan73397y30 (6y30vs7y30) were performed to investigate the transcriptional regulation underlying the clonal MW difference, while RRIM600y30 (6y30) and Reyan73397y30 (7y30) were respectively compared with RRIM600y1 (6y1) and Reyan73397y1 (7y1) to study age-related latex variation. The numbers of DEGs for these four comparisons are shown in Table 2. A total of 2648 genes were upregulated and 2553 genes were downregulated in the latex of RRIM600y30 compared with that of Reyan73397y30. Interestingly, the comparison 7y30vs7y1 showed the most DEGs, at 8688, indicating distinct transcriptional regulation in Reyan7-33-97 latex as tree age increased. A Venn diagram was constructed to show the relationships of DEGs among different comparisons. As shown in Fig. 2A, 2124 DEGs overlapped among the comparisons 6y30vs6y1 and 7y30vs7y1, suggesting their roles in age-related differences in latex samples. The 653 overlapping DEGs between comparisons 6y1vs7y1 and 6y30vs7y30 are supposed to have clonally different expression patterns. Hierarchical cluster analysis of gene expression profiles indicated that DEGs obtained from the transcriptome sequencing analysis could be divided into 4 clusters, and the genes in the same subcluster had similar expression patterns among these latex samples (Fig. 2B and Fig. S1). The expression patterns of DEGs in Reyan73397y30 were less consistent with those in RRIM600y30, RRIM600y1 and Reyan73397y1 (Fig. 2B). To validate the results obtained from the transcriptome analyses, the expression levels of 17 key rubber biosynthesis genes and regulators were further studied by qRT-PCR. As shown in Fig. S2, the qRT-PCR results were significantly correlated with those of RNA-seq (R2=0.950), indicating that the expression data obtained by transcriptome sequencing were reliable.
Table 2. Numbers of the DEGs across four comparisons.
|
6y1vs7y1
|
6y30vs7y30
|
6y30vs6y1
|
7y30vs7y1
|
Upregulated
|
840
|
2648
|
1302
|
3797
|
Downregulated
|
1046
|
2553
|
2158
|
4891
|
Total
|
1886
|
5201
|
3460
|
8688
|
Expression levels of genes related to rubber biosynthesis
Previously, a total of 85 genes were identified as rubber biosynthesis genes, including 18 involved in the MVA pathway, 22 in MEP, in for APP initiator biosynthesis and 30 in rubber particle-associated rubber elongation [3, 6, 26]. Forty-nine genes among these rubber biosynthesis genes were identified as DEGs in the four comparisons. A simple mode diagram of the rubber biosynthesis pathway was drawn, and the FPKM values of the rubber biosynthesis genes in four latex samples are shown in the form of a heatmap (Fig. 3) and listed in Table S2. In Reyan73397y1 as compared with RRIM600y1, 2 DEGs (HMGR4 and SRPP6) were upregulated, while 3 DEGs (REF4, CPT9 and DXR1) were downregulated. In addition to HMGR4, the rest of the MVA pathway genes showed slightly higher expression levels in Reyan73397y1 than in RRIM600y1, although they did not reach the DEG thresholds. A total of 27 DEGs were upregulated and 3 DEGs (HMGR1, MCS2 and GGPS4) were downregulated in RRIM600y30 compared with Reyan73397y30, among which 13 genes in the MVA pathway and 9 genes in the MEP pathway were upregulated (Table S2 and Fig. 3). To initiate rubber biosynthesis, IPP was first converted to its highly electrophilic isomer, DMAPP, which is catalyzed by the isopentenyl diphosphate isomerase (IPPI) gene. Our results showed that both IPPI1 and IPPI2 were present at relatively low levels in latex. IPPI1 was comparably expressed among the four latex samples, while IPPI2 showed higher expression levels in latex from 1-year-old trees than in that from 30-year-old trees (Table S2 and Fig. 3). A group of trans-prenyltransferases, including GGPS, FPPS and GGPPS, catalyze the biosynthesis of APPs, which can serve as rubber molecule initiators and also the precursors of non-rubber isoprenoids. The expression levels of these trans-prenyltransferase genes were downregulated in Reyan73397y30 compared with those in RRIM600y30 but displayed no significant difference between those in RRIM600y1 and Reyan73397y1. In addition, RRIM600y30 latex accumulated significantly higher levels of CPT and CPTL transcripts than Reyan73397y30 latex.
KEGG pathway enrichment analysis and gene set enrichment analysis (GSEA) were performed to further investigate the expression levels of genes in pathways related to rubber biosynthesis. Enriched KEGG pathways in comparisons of 6y1vs7y1 and 6y30vs7y30 are respectively shown in Figs. 4A and 4B, and the enriched gene sets related to isoprenoid pathways from GSEA are listed in Table S3. Pyruvate metabolism (pop00620) and glycolysis/gluconeogenesis (pop00010) pathways, which are the pathways producing pyruvate and acetyl-CoA precursors for IPP biosynthesis, were enriched in RRIM600y30 compared with those in Reyan73397y30. The terpenoid backbone biosynthesis pathway (pop00900) was significantly enriched in RRIM600y30, which contributes to higher levels of rubber biosynthetic enzymes. In Reyan73397y30 as compared with RRIM600y30, the ribosome, proteasome and flavonoid pathways were significantly upregulated, while fatty acid biosynthesis and metabolism genes were downregulated (Fig. 4B). In the comparison of Reyan73397y1 vs. RRIM600y1, non-rubber isoprenoid pathways, including carotenoid biosynthesis (pop00906), zeatin biosynthesis (pop00908), ubiquinone and other terpenoid quinine biosynthesis (pop00130) and diterpenoid biosynthesis (pop00904), were enriched in the latex of Reyan73397y1, and the N-glycan biosynthesis (pop00510) pathway was enriched in that of RRIM600y1 (Table S3). In latex from 30-year-old trees, the non-rubber isoprenoid branches, including sesquiterpenoid and triterpenoid biosynthesis (pop00909), steroid biosynthesis (pop00100), zeatin biosynthesis (pop00908) and N-glycan biosynthesis (pop00510), showed higher expression levels in RRIM600 latex than in Reyan7-33-97 latex (Table S3). To clearly show the expression patterns of the DEGs in these enriched non-rubber isoprenoid pathways, their FPKM values are shown in Table S4 and Fig. 4C.
Expression levels of genes encoding rubber particle associated proteins
REF and SRPP are the major proteins associated with rubber particles. Comparative transcriptome analysis in our present study showed that the REF and SRPP genes are differentially expressed among these latex samples. As shown in Fig. 3 and Table S2, the majority of REF genes showed higher expression levels in latex from 30-year-old trees than in that from 1-year-old trees, except for REF5 and REF8, which showed lower expression levels in Reyan73397y30. Additionally, latex from 30-year-old trees displayed a higher abundance of transcripts for 4 SRPP isoforms (SRPP1, SRPP2, SRPP3 and SRPP4), while expressed the rest of SRPP genes at lower levels than those in 1-year-old trees (Fig. 3 and Table S2).
In addition to REF/SRPP genes, we identified 136 transcripts encoding rubber particle associated proteins based on the previously published proteome data of rubber particles [27]. Among these genes, 68 were differentially expressed among the four latex samples, and their expression levels are listed in Table S5 and Fig. 5A. Several genes that have been previously identified as rubber biosynthesis regulators were shown, including rubber biosynthesis stimulator protein (RBSP) and rubber biosynthesis inhibitor protein (RBIP). In trees of the same age, latex from RRIM600 showed higher expression levels of RP proteins and smaller RPs than latex from Reyan7-33-97 trees.
Analysis of differentially expressed TFs that are potentially related to IPP biosynthesis
A higher supply of IPP precursors is essential for increasing the yield and MW of rubber. However, current knowledge of the transcriptional regulation of IPP-generating pathways in Hevea is very limited. In our present study, a total of 3077 genes with detectable expression levels in Hevea latex were identified as transcriptional regulators. Among these genes, 466 were differentially expressed among the four latex samples, and their FPKM values are listed in Table S6. IPP is biosynthesized using sucrose as the raw material, which is first metabolized into pyruvate and then into acetyl-CoA [2, 35]. We thus investigated the TFs that are potentially involved in the regulation of pyruvate metabolism and IPP biosynthetic pathways. Based on the above findings, 32 genes in pyruvate metabolism and 22 genes in the MVA and MEP pathways showed significantly higher expression levels in RRIM600y30 than in Reyan73397y30. Therefore, the cis-regulatory elements on the promoters (−1 to −1000 bp from the transcription start site) of these upregulated genes were searched for enriched motifs using the MEME suite [36]. The top three statistically significant motifs (depending on the E-value) and their corresponding TF families are shown in Fig. 5B. The results showed that the recognition sites of four TF families, including BASIC PENTACYSTEINE (BPC), WRKY, NAC and AP2, were significantly enriched in the promoters of those upregulated DEGs in the pyruvate metabolism, MVA and MEP pathways in RRIM600y30. These TFs are thought to be involved in the regulation of IPP biosynthesis. As shown in Fig. 5C, among the differentially expressed TFs, we identified only one BPC gene (LOC110656288) whose expression level was relatively low in Hevea latex and appeared to be positively correlated with the expression levels of IPP biosynthetic pathways. Interestingly, a large group of WRKY, NAC and AP2 genes showed significantly higher expression levels in Reyan73397y30 than in RRIM600y30, revealing their potential negative roles in the regulation of IPP production.