Phenotype characteristics of peanut dwarf mutant
The main stem height of sdm2 is significantly shorter than that of FH1 (P<0.01) (Fig. 1a). We measured the main stem height of mutant and FH1 grown in the experimental farm every two weeks during the vegetative growing stage. Results showed that during the first two weeks, main stem height between sdm2 and FH1 had no obviously difference. At 28 DAP (days after planting), main stem height of sdm2 was significantly shorter than that of FH1 (P<0.01). In the following developmental stages, the height difference between sdm2 and FH1 became more and more obvious. At 84 DAP, the sdm2 was 15.9 cm in height, only about 60% of FH1 (26.8 cm) (Fig. 1b). The leaf of sdm2 was smaller, thicker and dark green in color, and showed delayed senescence compare to that of FH1. At harvest time, the number of internodes of sdm2 was 3~4 less than FH1. The number of branches (P<0.05) and pods (P>0.05) per plant was obviously more, while the 100-pod weight and 100-kernel weight were significantly decreased in sdm2 (P<0.01).
To further accurately measure the change of mutant main stem height, we conducted hydroponic culture of sdm2 and FH1 in Hogland solution. The length of each internode and hypocotyl of 14 DAP seedlings were measured. We found that it was different from the results in experimental farm. For 14 DAP seedlings, the sdm2 height was already significantly shorter than that of FH1 (P<0.05) (Fig. 1c, 1d). There were total five internodes in FH1 and four in sdm2 seedlings. The hypocotyl length of sdm2 was slightly shorter than FH1 and there was no significant difference (P>0.05) for the 14 DAP seedlings. The length of most internodes of sdm2 was significantly shorter than FH1 (Fig. 1d). These results indicated that the height difference between sdm2 and FH1 was caused by both reduced total number of internode and length of each internode.
Sequencing data analysis
To clarify global gene expression changes in sdm2, 12 cDNA libraries were constructed with stem and leaf from mutant and FH1 seedlings. The libraries were sequenced using BGISEQ-500 platform. A total of 0.26 billion raw reads were generated, and the average output of each sample was 21.7 million. After removing adaptor sequences, low-quality and N-containing reads, an average of 21.58 million clean reads were obtained from each library with the clean read ratio of 99.47% (Table 1). Approximately 19.27 million and 16.00 million clean reads in each library matched the reference genome and gene set perfectly with the average mapping ratio of 89.31% and 74.18%, respectively. To investigate gene expression correlation among samples, the Pearson correlation coefficient of all gene expression levels between each two samples were calculated (Additional file 1: Figure S1). It showed a higher correlation of similar tissue between sdm2 and FH1. The principal component analysis (PCA) could clearly divide all samples into two clusters according to two kinds of tissues (Additional file 2: Figure S2).
DEGs between sdm2 and FH1
From all the samples, a total of 65,435 genes were identified. The gene expression of stem and leaf between sdm2 and FH1 was analyzed, and 3,733 and 3,715 differentially expressed genes (DEGs) were identified, respectively (Additional file 3 and 4: Table S1 and S2). In sdm2 stem, 1,786 genes were up-regulated and 1,947 genes were down-regulated to compare with FH1 (Fig. 2a and Additional file 3: Table S1). In sdm2 leaf, the up- and down-regulated genes were 1,644 and 2,071, respectively (Fig. 2a and Additional file 4: Table S2). The majority of DEGs had different expression pattern in stem and leaf. Among 3,733 DEGs in stem, only 873 were also differentially expressed in leaf. There were 473 down-regulated and 163 up-regulated DEGs with the same change trend in stem and leaf. It is worth noting that the number of common down-regulated DEGs was obviously higher than that of up-regulated ones. Interestingly, a large number of DEGs showed opposite expression pattern in stem and leaf. For example, there were 149 genes down-regulated in mutant stem while up-regulated in leaf and 88 genes up-regulated in mutant stem but down-regulated in leaf (Fig. 2b).
Functional analysis of DEGs
GO classification showed that 1,724 DEGs from stem and 1,716 DEGs from leaf were classified into three categories: biological process, cellular component, and molecular function (Fig. 3 and Additional file 5: Table S3). In molecular function category, catalytic activity and binding were the most abundant terms. For cellular component category, cell, membrane, membrane part and organelle were the main terms. The top two terms of biological process were metabolic process and cellular process. The top 20 terms in stem and leaf were shown in Additional file 6: Fig. S3 and Additional file 7: Table S4.
KEGG pathway analysis was also carried out. According to the enrichment results, the top 20 pathways in stem and leaf were shown in Fig. 4 and Additional file 8: Table S5. In stem, many pathways involved in hormone biosynthesis and signal transduction including terpenoid backbone biosynthesis, monoterpenoid biosynthesis, diterpenoid biosynthesis, sesquiterpenoid and triterpenoid biosynthesis, zeatin biosynthesis, and brassinosteroid biosynthesis, indole alkaloid biosynthesis as well as plant hormone signal transduction were all enriched in the mutant (Fig. 4 and Additional file 8: Table S5). In addition, MAPK signaling pathway was also enriched. These results suggested that plant hormone play major roles in regulation of the mutant phenotype. There were other pathways including photosynthesis, phenylpropanoid biosynthesis, flavonoid biosynthesis, isoflavonoid biosynthesis, and flavone and flavonol biosynthesis pathways were found to be enriched.
In leaf, brassinosteroid biosynthesis, indole alkaloid biosynthesis, plant hormone signal transduction and MAPK signaling pathways were enriched in accordance with those in stem (Fig. 4 and Additional file 8: Table S5). Phenylpropanoid biosynthesis, flavonoid biosynthesis, isoflavonoid biosynthesis and plant-pathogen interaction were also in the top 20 enriched pathways. Different from that in stem, some glucide and glucolipid metabolic pathways were enriched in leaf.
Cell wall related genes
RNA-seq results showed that several CesA genes and all CSL genes were down-regulated in stem. The expression of most expansin genes was down-regulated in sdm2 stem compared with those in FH1 (Fig. 5 and Additional file 3: Table S1). Interestingly, the expression of CesA, CSL and expansin genes was up-regulated in leaf (Fig. 5 and Additional file 4: Table S2). These results were coincided with the phenotype of sdm2 leaf, dark green in color, thicker, smaller than FH1 leaf.
Some XTHs were significantly up-regulated and some were down-regulated in sdm2 stem (Fig. 5 and Additional file 3: Table S1). In leaf, the expressions of all XTHs were up-regulated, and some XTHs were 16-fold higher than those in FH1 (Fig. 5 and Additional file 4: Table S2). Some genes encoding endo-1,4-β-glucanases (EGases) in stem and leaf were down-regulated.
Expression changes of hormone biosynthesis and signaling genes
DWF4 was designated as CYP90B1 and proposed to catalyze steroid 22α-hydroxylation in BR biosynthesis [65]. A gene annotated as cytochrome P450 90B1 showed 73% identity with Arabidopsis DWF4 gene, which play key roles in BR biosynthesis. In stem, the expression of this gene was detected in FH1, while not detected in sdm2. In leaf, the expression of this gene in sdm2 was also lower than that in FH1. Previous studies confirmed that cytochrome P450 monooxygenase CYP90A1 (CPD) acted as BR C-3 oxidase [37]. The expression level of CYP90A1 (CPD) in sdm2 stem decreased by 1.2-fold than that in FH1. While in leaf, its expression was slightly higher in sdm2 than in FH1. The cytochrome P450 family members, CYP90C1/D1, CYP85 (Dwarf) and CYP71A1, mediate rate-limiting reactions in BR biosynthesis [66-68]. Our results showed that the expression of CYP90C1/D1 gene decreased by 2.64-fold and 1.41-fold in stem and leaf respectively in sdm2. The Dwarf genes which encoding cytochrome P450 85A-like in stem and cytochrome P450 85A in leaf were all up-regulated. CYP71A1 which catalyze the generation of deoxocastasterone (6-deoxoCS) and castasterone (CS) was down-regulated both in stem and leaf in sdm2. A cytochrome P450 monooxygenase (CYP734A1) which degrades BRs was up-regulated in sdm2 leaf (Table 2, Additional file 3 and 4: Table S1 and S2).
In both stem and leaf, the expression levels of BRI1 and BAK1 were down-regulated in sdm2. Noticeably, the expression levels of BES1/BZR1 gene were extremely low in stem and leaf in sdm2, 170-fold and 79-fold lower than in FH1, respectively. The expression levels of IWS1 in stem were up-regulated in sdm2, while in leaf, it was not differentially expressed. The expression of cyclin-D3-3 encoding gene was down-regulated in sdm2 stem, which was coincident with the expression trends of BR-related genes. While the expression levels of cyclin-D1 and cyclin-D5 increased in sdm2 leaf (Table 2, Additional file 3 and 4: Table S1 and S2).
Several GA biosynthesis and signal transduction genes were differentially expressed between sdm2 and FH1. In sdm2, most of GA20ox and 2-oxoglutarate-dependent dioxygenase encoding genes, which promote GA biosynthesis, were up-regulated in stem, while most 2-oxoglutarate-dependent dioxygenase decreased in leaf. While two genes encoding GA3ox were down-regulated in stem and its differential expression was not detected in leaf. GA2ox convert excess GAs to inactive forms through 2β-hydroxylation [29]. In mutant stem and leaf, genes encoding GA2ox1 and GA2ox2 were decreased in expression levels while the expression of three genes encoding GA2ox8 were up-regulated in stem (Table 2, Additional file 3 and 4: Table S1 and S2).
Two GID1 genes in stem showed opposite change trend and one GID1 gene was down-regulated in leaf. SCL3 (Scarecrow-like 3) seemed to attenuate DELLA protein and acted as a positive regulator in GA pathway [69]. In mutant stem, SCL3 gene was up-regulated while two SCL14 and one SCL21 genes were down-regulated. Interestingly, all genes encoding gibberellin-regulated protein (GRP) were reduced in sdm2 stem while up-regulated in leaf, which were consistent with the change trends of cell wall related genes including CesA, CSL, expansin and XTH genes (Table 2, Additional file 3 and 4: Table S1 and S2).
Four AUX1 genes in stem and two in leaf were all significantly up-regulated, while the expressional levels of two genes encoding auxin efflux carrier (PIN) in leaf decreased. In addition, the expression change trends of MDR/ABCB genes in stem and leaf of sdm2 were irregular, most of which were down-regulated and some were up-regulated. In sdm2 stem, two genes encoding AUX22 were down-regulated and almost all ARF members were up-regulated, while in leaf, the expressions of all Aux/IAA and ARF genes were increased. In addition, most SAUR genes in stem and almost all SAUR genes in leaf were up-regulated in sdm2 (Table 2, Additional file 3 and 4: Table S1 and S2).
Transcription factors involved in plant growth and development
Our results showed that several transcription factor genes involved in hormone signaling pathways and their downstream changed significantly. Gene encoding SHORT INTERNODES-like protein was down-regulated both in stem and leaf of sdm2. Whereas, SHI-RELATED SEQUENCE 3-like (SRS) gene increased in stem and decreased in leaf. All WRKY30/46/70 genes were down-regulated in sdm2 stem and leaf. The expression levels of four PRE6 genes decreased in stem and three PRE6 genes were up-regulated in leaf (Table 2, Additional file 3 and 4: Table S1 and S2).
Verification of DEGs using qRT-PCR
In order to validate the RNA-Seq data, the expression levels of DEGs were detected using Quantitative Real-time PCR (qRT-PCR). A total of 10 genes related to hormone signal transduction and cell wall organization were selected (Fig. 6 and Additional file 9: Table S6). The relative expression levels (log2 sdm2/FH1) of these genes estimated by qRT-PCR were generally consistent with those by RNA-seq. The overall correlation coefficient of a liner regression analysis in stem and leaf was 0.7737 and 0.8328, respectively (Fig. 6).