Morphological characteristics of diploid and triploid tea leaves
The leaf phenotypes of diploid and triploid annual shoots were determined under field growth conditions. As shown in Fig.1, the leaves of triploid tea trees have obvious growth advantages compared with diploid tea leaves, and the leaf length, leaf width and leaf area are larger than those of diploid leaves. Compared with diploid leaf length, leaf width and leaf area, they increased by 23.37%, 41.12% and 59.81%, respectively, with extremely significant differences. However, the diploid leaf thickness was significantly higher than that of the triploid leaf, which was 295.63 μm and 252.33 μm, respectively, and the difference was extremely significant. There was no significant difference in the length of the petiole of triploid and diploid tea, which were 0.7 cm and 0.667 cm, respectively. The diameter of the petiole of the triploid tea is larger than the diameter of the diploid petiole, with significant differences, 2.05 mm and 2.31 mm, respectively. The dry weight of the three-leaf leaves increased by 39.29% compared with the diploid, and the difference was extremely significant (Table 1).
Analysis of the number and size of stomata in diploid and triploid leaves
In order to investigate the difference in stomatal appearance on the epidermis of diploid and triploid tea leaves, we performed a statistical analysis of the size and density of the stomata of the epidermis in the same part of the tea leaves. The results showed that the stomata size (length and width) of triploid leaves were significantly larger than diploid, which increased by 57.20% and 84.44%, respectively, compared with diploid length and width. However, the stomata density of triploid leaves is significantly lower than that of diploids, and the stomata density of diploid tea leaves is twice that of triploid stomata (Fig. 2).
Paraffin section analysis of diploid and triploid leaves
In order to understand the difference in growth and development of diploid and triploid leaves, we observed paraffin sections. As a result, it can be seen from Fig.3 that the leaf vein portion of the diploid tea tree undergoes a significant change after triploidization. Among them, the xylem of the veins is the most obvious. The xylem of the triploid vein is more developed than the diploid, and the area is also larger. The area of the triploid xylem is three times that of the diploid, which is 0.476 mm2. The increase in the area of the triploid xylem is caused by two reasons: one is the increase in cell areas; the other is the increase in the number of xylem cells, and the number of diploid and triploid xylem cell layers averaged 19 and 25. The number of triploid xylem cell layers increased by 30.67% compared to diploid. There were no significant changes in the size of the diploid and triploid veins and the phloem and formation.
It can be seen from Fig. 4 that the shape and size of the diploid and triploid mesophyll cells are significantly different, and the epidermal cell thickness of the triploid mesophyll is 22.28% larger than that of the diploid, which is extremely significant. The cells of the leaf bark of diploid tea leaves are arranged closely, and the shape and size are relatively uniform. The cells of the triploid leaf mesophyll puncture are loosely arranged, the cell gap is large, and the cell shape and size are different. The average length of the diploid palisade tissue cells is larger than that of the triploid, and the length is 65 μm, which is 15.65% longer than the triploid. The width of the triploid palisade tissue cells was significantly larger than that of the diploid, which was 70% higher than that of the diploid. Diploid sponge tissue cells are small and dense, and the number is more than that of triploids, which is about twice the number of triploid cells.
Illumina sequencing and reads assembly
To investigate the molecular mechanisms of diploid and triploid leaf growth, and to understand the metabolic processes involved in leaf growth and development, we analyzed the gene expression profiles of diploid and triploid leaves. By de novo transcriptome sequencing, the test samples obtained an average of about 40 million original reads, with high quality reads reaching over 99%, six cDNA libraries (CaS419_1, CaS419_2, CaS419_3, CaS4_1, CaS4_2, and CaS4_3) were generated from diploid and triploid mRNAs, which were sequenced using Illumina deep-sequencing HiSeq ™ 2000.
The raw data obtained after sequencing on the machine is filtered, the sequencing error rate is checked, the GC content distribution is checked, and the GC content analysis is used to detect whether there is A/T or G/C separation phenomenon. Finally, 6 sample data (clean read) for subsequent analysis were obtained. The filtered data is summarized in the table below. Among all the raw reads, 96 % had Phred-like quality scores at the Q20 level (an error probability of 1 %).After removing adapters, low-quality sequences and ambiguous reads, we obtained approximately 45 million, 43 million, 60 million,58 million ,62 million and 48 million clean reads from the diploid samples (CaS419_1, CaS419_2 , and CaS419_3, ), and triploid samples (CaS4_1 , CaS4_2, and CaS4_3), respectively (Table 2). Raw reads are filtered and assembled by De novo assembly software Trinity. The assembled sequences are redundant and spliced by software TGICL to obtain the longest non-redundant unigene set and further statistics on unigene sets.
Functional annotation and Cluster Analysis
Due to the lack of a complete genome sequence in Camellia sinensis, only 27,031 unigenes were co-annotated into six databases (NR, NT, SwissProt, COG, GO, and KEGG), accounting for 26.13% of 103,448 unigenes. Among them, the most frequently cited genes in the NCBI NR and NT databases were 90,547 and 89,933 unigenes (87.53% and 86.93% of all annotated unigenes), while 35,298 (34.12%) and 61,318 (59.27%) unigenes could be annotated into the COG, Swiss-Prot database. We annotated 45,820 (44.29%) and 67.980 (65.71%) unigenes to GO and KEGG databases (Fig. 5).
The main GO terms included biological process (BP), cellular component (CC), and molecular function (MF). Based on sequence homology, 45,820 unigenes were mainly categorized into 55 functional groups (Fig. 6). In the category of BP, the two major groups of cellular processes and metabolic processes accounted for the highest proportion. Of these, approximately 24,750 genes have been annotated as metabolic process categories, which may allow the identification of novel genes involved in secondary metabolism pathways in triploid. As for the MF category, unigenes with binding and catalytic activity formed the largest groups. For CC, the top three largest categories were cell, cell part, and membrane. To further evaluate the reliability of our transcriptome results and the effectiveness of our annotation process, we searched the annotated sequences for genes with COG classifications (Fig. 7). Among the 26 COG categories, the cluster for “General function prediction only” (9,286) represented the largest group, followed by “Transcription” (5,108), “Posttranslational modification, protein turnover, chaperones” (4,438), and “Replication, recombination and repair” (4,369). The categories “Extracellular structures” (6) was the smallest group.
To analyze of the biological functions of the unigenes, we used the annotated sequences to comparisons against the KEGG database. In total, 67,980 annotated unigenes were assigned to 136 known pathways based on the KEGG BLAST analysis. The top 19 pathways with the largest numbers of unigenes are listed in Fig. 8. The majority of the unigenes (22,658; 31.95%) were involved in Global and overview maps pathways, followed by pathways in Carbohydrate metabolism (7,072 unigenes; 9.97%), Translation (6,418 unigenes; 9.05%), and Folding, sorting and degradation (4,408 unigenes; 6.22%).
Differentially Expressed Genes and qRT-PCR Validation Between Diploid and Triploid Camellia sinensis
To confirm the results of the Solexa/Illumina sequencing, twelve unigenes were selected for quantitative RT-PCR assays. The qRT-PCR analysis performed for ten upregulated and two downregulated DEGs growth-related genes confirmed the transcriptomic changes detected by RNA-seq (Fig. 9). Although the expression levels did not exactly match; however, Quantitative real-time PCR analysis showed that the patterns of gene expression were consistent with RNA-seq results. Thus, qRT-PCR results validated the reliability of RNA-seq data.
Analysis of key gene expression of stomatal development based on transcriptome results
In order to analyze the difference in stomatal development of the diploid and triploid leaves of tea plants, We identified 16 differentially expressed stomatal-related genes (P < 0.005), differences in expression of these genes lead to changes in triploid stomatal density and size. Nine genes belong to the Negative-regulatory factors, which are key regulators of stomatal in plants. In the negative regulatory family, SDD1, SERK1,2 and EPF1,2 have a negative regulatory effect on stomatal development, whileEPLF9/ stomagenhas a positive regulatory effect. Among them, the SDD1 and SERK1 genes are up-regulated and the EPF1 gene is down-regulated (Table 3). SERKs can interact with TMM in a non-ligand-dependent manner to form a multiprotein receptor complex and negatively regulate stomatal development through signal transduction. These genes are involved in biosynthesis and signal transduction of stomatal development by participating in different biological processes. For example, During stomatal development, cysteine-rich secretory peptides belonging to the EPF/EPFL family act as ligands to interact with the corresponding receptors to transmit developmental signals, ensuring proper stomatal density and distribution. The key transcription factor involved in stomatal development in plants is the bHLH type protein, which Family SPCH, FAMA and MUTE play important regulatory roles in stomatal development. There was no significant difference in expression between these three transcription factors in diploid and triploid.
The COP and HIC genes are stomatal development factors that respond to light and CO2 signaling, in which the COP gene is down-regulated and HIC is up-regulated.
Identification and expression analysis of candidate genes involved in leaves development based on transcriptome results
The morphological structure of the leaves seems simple, but the regulation mechanism of its development is very complicated. The final size of leaves, is tightly controlled by environmental and genetic factors that must spatially and temporally coordinate cell expansion and cell cycle activity. In this study, we identified 28 putative genes associated with leaf development that belong to different pathways, including cell division, photosynthesis, transcription factor, and auxin synthesis, which showed significant differential expression between diploid and triploid.
Transcripts involved in Photosynthesis-photoreaction phase were observed to be differentially expressed, including genes encoding Photosystem I (PS Ⅰ), photosystem II (PS Ⅱ), cytochrome b6/f complex, and ATP synthase. The center action of PS I is the pigment molecule P700, while Psa A and Psa B are key genes regulating the synthesis of P700 chlorophyll a apolipoprotein A1 and A2. The photoreaction center pigment of PSII is P680, and the PsbA and PsbE genes are involved in the synthesis of the P680 reaction center D1 protein. PetB is a key gene involved in the synthesis of the cytochrome b6f complex. F1B is a F-type H+ transport ATPase subunit β gene in ATP synthase. All these genes were down-regulated in leaves and are functional essential for carbon dioxide assimilation. Interestingly, it is contrary to the photoreaction phase of photosynthesis gene expression, all of the key enzymes in the carbon reaction phase is up-regulated. Among them, ribulose-1,5-bisphosphate carboxylase/oxygenase Rubisco gene, phosphoenolpyruvate carboxylase PPC and malate dehydrogenase MDH1 were significantly up-regulated (Table S1).
To explore the intracellular transcriptional activity of diploid and triploid plants, we analyzed the expression of genes involved in the regulation of RNA polymerase and transcription initiation factors during transcription. The results showed that RNA polymerase I, RNA polymerase II and transcription initiation factor were both up-regulated. The RPB2 gene and the TFIIA1 gene regulating RNA polymerase II and transcription initiation factor were significantly up-regulated, which was up-regulated by about 2.5-fold compared to diploid (Table 4).
The expression of cell cycle-associated genes that regulate cell division can alter the organ volume of a plant, showing an increase in the number of cells and expansion of the cell volume in time and space. Cyclin-dependent kinase (CDK), also known as the cell cycle engine molecule, plays a central role in the regulation of cell cycle function. The cyclin-dependent kinase 1 (CDK1) gene is up-regulated 6-fold compared to the diploid (Table 4). The CDK1 gene plays a key role in controlling the eukaryotic cell cycle by regulating centrosome circulation and mitosis initiation, promoting G2-M conversion, and regulating G1 and G1-S transformation by binding to multiple interphase cyclins. Cyclin-dependent kinase 7, CDK7, is a catalytic subunit of the CDK-activated kinase (CAK) complex that regulates cell cycle progression. The relative diploid expression of CDK7 gene was up-regulated by 22.66% in triploid leaves of tea (Table 5). The serine/threonine kinase BUB1 gene is involved in cell cycle control and RNA polymerase II-mediated RNA transcription, and is also up-regulated in triploids, which is about twice the expression of diploid genes.
Taken together, these results provide a framework for the regulatory network of leaf development response in diploid and triploid.