Skip to main content

Comparative transcriptome analysis of high- and low-embryogenic Hevea brasiliensis genotypes reveals involvement of phytohormones in somatic embryogenesis



Rubber plant (Hevea brasiliensis) is one of the major sources of latex. Somatic embryogenesis (SE) is a promising alterative to its propagation by grafting and seed. Phytohormones have been shown to influence SE in different plant species. However, limited knowledge is available on the role of phytohormones in SE in Hevea. The anther cultures of two Hevea genotypes (Yunyan 73477-YT and Reken 628-RT) with contrasting SE rate were established and four stages i.e., anthers (h), anther induced callus (y), callus differentiation state (f), and somatic embryos (p) were studied. UPLC-ESI-MS/MS and transcriptome analyses were used to study phytohormone accumulation and related expression changes in biosynthesis and signaling genes.


YT showed higher callus induction rate than RT. Of the two genotypes, only YT exhibited successful SE. Auxins, cytokinins (CKs), abscisic acid (ABA), jasmonic acid (JA), salicylic acid (SA), gibberellins (GAs), and ethylene (ETH) were detected in the two genotypes. Indole-3-acetic acid (IAA), CKs, ABA, and ETH had notable differences in the studied stages of the two genotypes. The differentially expressed genes identified in treatment comparisons were majorly enriched in MAPK and phytohormone signaling, biosynthesis of secondary metabolites, and metabolic pathways. The expression changes in IAA, CK, ABA, and ETH biosynthesis and signaling genes confirmed the differential accumulation of respective phytohormones in the two genotypes.


These results suggest potential roles of phytohormones in SE in Hevea.

Peer Review reports


The Rubber tree (Hevea brasiliensis Muell. Arg. hereafter Hevea) is a tropical tree that was first discovered in the Amazonian basin of South America [1]. It is a member of the Euphorbiaceae family. It is frequently farmed in Southeast Asia to produce natural rubber. One of the major sources of latex is the Hevea, which serves as the basis for more than 40,000 products [2]. Hevea multiplication has been hampered by its lengthy growth cycles, sever inbreeding depression upon self-pollination, cross-pollination, recalcitrant seed, and poor seed germination [3, 4]. As a result, the tree is mostly propagated through grafting. However, the growth and natural rubber yield of a grafted tree can be influenced by the interaction between its scion and rootstock [4]. Somatic embryogenesis (SE) is a quick and efficient technique for plant regeneration and propagation [5]. SE has been successfully achieved in the case of rubber trees by employing a variety of explants including roots, inflorescences, internal integuments of immature fruits, and anthers [6, 7]. Self-rooted juvenile clones (SRJCs) are plants that have grown through the SE process and have their own roots and juvenile characteristics. In comparison to donor clones (DCs), SRJCs exhibit higher stress resistance, growth, and rubber yield, making them a promising option for future rubber tree planting material [5, 7, 8]. Considering urgent need for the development of natural rubber industry, it essential to develop somatic embryos propagated plants for the excellent varieties of rubber trees. However, rubber plant regeneration through SE is largely dependent on the genotype [9,10,11]. Consequently, the ability to obtain regenerated plants has been restricted to a narrow range of rubber tree genotypes [5, 12,13,14,15].

Phytohormones play a crucial role in somatic cell embryogenesis in plants. Auxins, cytokinins (CKs), and abscisic acid (ABA) are the most critical hormones involved in this process [16]. Several investigations have found that auxin and CK ratios are critical for determining the developmental fate of somatic cells during embryogenesis [17]. High levels of auxin and low levels of CK are required for the induction of SE, while a balance of both hormones is needed for the development of somatic embryos [18]. ABA has also been implicated in the regulation of embryo growth and maturation in SE [19]. Other hormones such as gibberellins (GAs) and ethylene (ETH) may also influence SE, but their roles are not yet fully understood. Overall, the complex interactions of various hormones are essential for the successful induction and development of somatic embryos in plants.

The ability to achieve SE in rubber tree is limited to certain genotypes, as documented in earlier studies [12,13,14]. Despite the significance of SE in rubber tree, there is limited understanding of the molecular regulatory mechanism driving this biological process. A pioneering study (two decade ago) identified five differentially expressed cDNAs in an embryogenically regenerating Hevea line that could potentially aid in early diagnosis of the embryogenic potential of friable rubber tree callus [12]. However, the functions of these cDNAs remain unknown. In another study [20], researchers used transcriptome analysis to investigate the molecular mechanisms underlying SE in Hevea. However, the analysis was limited to a single genotype. To investigate the molecular regulatory mechanisms of plant SE, RNA-seq analyses have been performed to identify SE-related genes in various plant species, including Arabidopsis [21], Gossypium hirsutum [22], maize [23], strawberry [24], rice [25], coconut palm [26], and Norway spruce [27]. These studies revealed regulatory mechanisms of SE at the molecular level and identified potential key genes such as somatic embryogenesis receptor-like kinase (SERK) [28], late embryogenesis abundant (LEA) protein [29], AGAMOUS-like 15 (AGL15), leafy cotyledon, WUSCHEL (WUS), WUS-homeobox 2, and Baby Boom(BBM) ([30] and references therein, [31]). Li, et al. [32] identified three differentially expressed MADS-box genes during early embryogenesis in rubber tree. Another study reported 11 AP2/ERF genes that may serve as expression markers for different stages of the SE process in rubber tree [14]. Other studies have also reported that AP2/ERF genes are involved in regulation of SE [33,34,35]. However, the molecular mechanisms regulating SE and the effects of phytohormones on genotypes with different SE capabilities remain poorly understood.

In this study, we explored the endogenous phytohormone content of two H. brasiliensis genotypes differing in SE potential. We used transcriptome sequencing analysis to explore the expression differences in the two genotypes regarding phytohormone biosynthesis and signaling.


Plant material

Two H. brasiliensis varieties, Yunyan 73477 (YT) and Reken 628 (RT), were used as plant material. The varieties YT and RT exhibit high and low SE, respectively. However, the latter is the main variety in rubber planting areas in China, therefore, it is an important resource to be improved. The plant materials were obtained from the Genebank of Yunnan Institute of Tropical Crops. They are conserved at the Genebank of Yunnan Institute of Tropical Crops under voucher numbers YITC4002 and YITC4009. No permission is needed to study these accessions. The formal identification of plant materials was undertaken by Prof Guoping Liang. Male flowers were used as starting material. Male flowers were soaked in 75% (v/v) ethanol for 30 s, then soaked in 0.1% HgCl2 (w/v) for 10 min followed by five rinses with distilled water. The YT immature anthers were isolated and inoculated on solid callus induction medium (MS medium containing 2,4-D 1 mg L−1, KT 1.0 mg L−1, NAA 1.0 mg L−1, sucrose 70 g L−1 and coconut water 50 ml L−1) for 50 days. Then transferred to a differentiation medium (solid MS medium with activated charcoal 1.0 g L−1, KT 2.0 mg L−1, NAA 0.1 mg L−1, GA3 0.5 mg L−1, ABA 0.2 mg L−1, sucrose 70 g L−1 and coconut water 50 ml L−1) to induce somatic embryos. After 25 days of differentiation, embryonic calli with spherical embryos were obtained. After differentiation and cultivation for 60 days, somatic embryos in cotyledonary stage were obtained. The RT anthers were inoculated for 50 days to induce calli and the non-embryonic calli were cultured for 25 days. Callus induction rate and SE rate were measured according to the following equations.

$$Callus\ induction\ rate=\frac{No.\ of\ anther\ explants\ forming\ callus}{number\ of\ anther\ explants\ inoculated}\times 100\%$$
$$Somatic\ embryogenesis\ rate=\frac{No.\ of\ calli\ with\ embryogenesis\ ability}{number\ of\ calli\ inoculated}\times 100\%$$

Triplicate samples of male flower (RT-h), and two callus stages: anther callus induction stage (RT-y) and differentiation culture stage (RT-f), were taken for RT, frozen in liquid nitrogen, and stored at -80 ℃ until further processing. Similarly, triplicate samples of male flower (YT-h, two callus stages: anther callus induction stage (YT-y) and differentiation culture stage (YT-f), and somatic embryo (YT-p) were collected for YT, frozen in liquid nitrogen, and stored at -80 ℃ until further processing.

Endogenous plant hormone analysis

Reagents were purchased from Merck (Darmstadt, Germany), Millipore (Bradford, USA), Olchemim Ltd. (Olomouc, Czech Republic), isoReag (Shanghai, China), and Sigma-Aldrich (St Louis, MO, USA). Stock solutions of the standards were prepared at a concentration of 1 mg/mL in methanol (MeOH) and stored at -20 ℃. Working solutions were prepared from the stock solutions by diluting with MeOH.

From samples stored at -80 ℃, 50 mg of each sample (each replicate) was ground to powder (30 Hz, 1 min), followed by extraction in 1 mL MeOH/water/formic acid (15:4:1, V/V/V). For quantification, we added 10 μL of internal standard mixed solution (100 ng/mL) to the extracts as internal standards (IS). The mixture was vortexed for 10 min, centrifuged at 12,000 r/min at low temperature (4 ℃) for 5 min, and the supernatants were transferred to clean plastic microtubes. The supernatants were dried by evaporation, dissolved in 100 μL 80% MeOH (v/v), and filtered through a 0.22 μm Durapor ® membrane filter (Merck KGaA, Darmstadt, Germany). The filtrates were then used for LC–MS/MS analyses using an UPLC-ESI–MS/MS system (UPLC, ExionLC™ AD,; MS, Applied Biosystems 6500 Triple Quadrupole, The analytical conditions were set as follows, LC: column, Waters ACQUITY UPLC HSS T3 C18 (100 mm × 2.1 mm, 1.8 µm); the solvent system included water with 0.04% acetic acid (A) and acetonitrile with 0.04% acetic acid (B); the gradient program started at 5% B for 0–1 min, increased to 95% B for 1–8 min, maintained at 95% B for 8–9 min, and finally decreased to 5% B for 9.1–12 min; the flow rate was 0.35 mL/min; the temperature was 40 °C, and the injection volume was 2 μL. Linear ion trap and triple quadrupole scans were acquired on a triple quadrupole-linear ion trap (QTRAP) mass spectrometer, QTRAP® 6500 + LC–MS/MS system. The system was equipped with an ESI Turbo Ion-Spray interface, operated in two ion modes: positive and negative, and controlled by Analyst 1.6.3 software (Sciex). The operating parameters of the ESI source were as reported earlier [36]. Phytohormones were analysed by scheduled multiple reaction monitoring (MRM) and data acquisitions were performed using Analyst 1.6.3 software. Multiquant 3.0.3 software (Sciex) was used for quantification. All the mass spectrometer parameters including the declustering potentials (DP) and collision energies (CE) for individual MRM transitions were performed with further DP and CE optimisation. A specific set of MRM transitions was monitored for each time period according to the metabolites eluted within that time period.

Statistical and bioinformatic analysis of phytohormone data

Hierarchical cluster analysis (HCA) was performed using the R package pheatmap and presented as heat maps of normalised signal intensities of metabolites with dendrograms. Significantly differentially accumulated metabolites (DAMs) between groups were determined by absolute Log2FC (fold change). The DAMs were annotated using the KEGG compound database ( and mapped to the KEGG pathway database ( Pathways with significantly regulated metabolites mapped were then fed into MSEA (metabolite sets enrichment analysis), their significance was determined by hypergeometric test’s p-values.

Transcriptome sequencing

RNA extraction, library preparation, and sequencing

Total RNAs were extracted from the triplicate samples of each tissue of the two genotypes using a Spin Column lant total RNA Purification Kit (Tiandz, Beijing, China) [37]. The RNA purity and integrity were checked, and sequencing libraries were prepared as previously reported [36]. The libraries were sequenced on the Illumina sequencer.

Sequencing data analysis

The quality of the sequencing libraries was checked using fastp [38] and removed reads with adapters, if the N content in the sequencing read exceeded 10% of the number of bases in the read, or if the number of low-quality bases (Q ≤ 20) in a sequencing read exceeded 50% of the number of bases in that read. This step was followed by checking the sequencing error rate distribution check and the GC content distribution. We used HISAT2 [39] to compare the clean reads to the reference genome [40]. Gene expression was calculated as Fragments Per Kilobase of transcript per Million fragments mapped (FPKM) and visualized as a boxplot. Pearson’s correlation coefficient (PCC) and principal component analysis (PCA) were calculated in R ( Differential gene expression was calculated in DESeq2 [41] followed by correction of the hypothesis testing probability using the Benjamini–Hochberg method [42] for multiple hypothesis testing and to obtain the false discovery rate (FDR). The differential genes were screened using the criteria |log 2-fold change|> = 1, and FDR < 0.05. Heat maps and Venn diagrams were generated in tBtools [43]. Differentially expressed genes (DEGs) were annotated in the KEGG [44] and GO [45] databases. The DEGs were enriched in KEGG pathways [46, 47].

Antioxidant enzyme activity and H2O2 content determination

The activities of the antioxidant enzymes, peroxidase (POD), superoxide dismutase (SOD), and catalase (CAT), and the content of H2O2 in the immature anther culture stages of YT and RT were studied according to established methods [48]. Whereas, the ascorbate peroxidase (APX) activity was determined according to the method used in an earlier study [49].


Genotype dependent SE in H. brasiliensis

The process of SE in two H. brasiliensis genotypes i.e., YT and RT was established following a previously described method [50] (Fig. 1). Both varieties of anthers (h) induced callus (y) in MS medium containing 2,4D, KT, and NAA. The callus induction rate of genotype YT (89.9%) was higher than that of genotype RT (82.8%). Twenty-five days after the callus was transferred to the differentiation medium containing GA3, ABA, high KT content and low NAA content (that is, somatic cell embryo induction culture), it was observed that only YT callus had SE ability. Its embryogenic callus had spherical or cardioid cell embryos (YT-f), and the somatic embryo induction rate was 83.04%. After 60 days of differentiation and culture, cotyledon type SEs (YT-p) were obtained. But RT callus was not embryogenic (RT-f). After differentiation and culture, no somatic cell embryogenesis occurred (Fig. 2a).

Fig. 1
figure 1

Immature anther culture of two H. brasiliensis varieties: Yunyan 73477 (YT) and Reken 628 (RT) differing in SE ability. SE: somatic embryogenesis, h: immature anther, y: anther induced callus, f: callus differentiation stage, p: cotyledonary embryo

Fig. 2
figure 2

a Callus induction rate and SE rate in Yunyan 73477 (YT) and Reken 628 (RT). b-n Relative metabolite intensities of phytohormones in anthers (h), anther induced callus (y), callus differentiation state (f), and somatic embryos (p) in YT and RT

UPLC-ESI–MS/MS based phytohormone detection in the two genotypes

The employed system detected seven major phytohormones: auxins, CKs, ABA, jasmonic acid (JA), salicylic acid (SA), GAs, and ETH and their conjugates/precursors (Supplementary Table S1). We observed that the relative content of ABA was lower in YT-y and YT-f compared to YT-h, while in the case of YT-p, the content was 43.78% higher than YT-h. In the case of RT, the relative ABA content was lower in RT-y compared to RT-h, however, it was higher in RT-f. Overall, we observed that the relative ABA content was higher in YT-h and YT-y compared to the RT-h and RT-y. However, in the case of differentiated callus (f), the ABA content was higher in RT than YT.

During SE, the levels of ABA glycosyl ester varied (Supplementary Table S1). The pattern of changes observed in ABA glycosyl ester levels was similar to that of ABA, except for in the case of YT-p, where ABA glycosyl ester levels were lower than ABA. These observations suggest potential conversion of ABA glycosyl ester into ABA during SE. The levels of IAA differed in YT-h and RT-h such that the latter had slightly higher content than the former. Interestingly, the IAA content decreased from YT-h to YT-y and then increased in YT-f and YT-p. However, despite the increase in latter stages, the content remained lower than in YT-h. On the contrary, we did not detect IAA in RT-y neither in RT-f. This is an important observation and may be a possible reason for the different SE capacity of the two genotypes.

There were 19 metabolites related to auxin (Supplementary Table S1, Supplementary Figure S1). Among them, L-tryptophan was highest in YT-h and RT-h stages (YT-h with higher content). The relative content of indole-3-acetyl-L-aspartic acid was highest in YT-f and YT-p. However, 1-O-indol-3-ylacetylglucose was the major IAA-component in YT-h, RT-h, and RT-y. For CKs 28 metabolites were identified. These were classified into precursors, storage forms, byproducts, degradation products, and intermediates (Supplementary Figure S2). Two major CKs i.e., cis-Zeatin (cZ) and trans-Zeatin (tZ) showed different accumulation trends in the two genotypes. We observed that the relative content of cZ decreased in both genotypes with the progression of the callus induction, redifferentiation, and SE (in YT), whereas, tZ relative content showed a similar trend in RT but not in YT (Fig. 2). The tZ content in RT-h was relatively lower than YT-h. In case of RT, tZ content decreased in the following stages whereas in case of YT, its content significantly increased from YT-h to YT-y. However, in the latter stages i.e., YT-f and YT-p, it decreased compared to YT-h. In YT, the precursor of CK (4-[[(9-beta-D-Glucopyranosyl-9H-purin-6-l) amino]methyl]phenol (also known as isopentenyladenine-adenine) followed an increasing trend from YT-h to YT-f and then decreased in YT-p. However, in the case of RT, it was detected from RT-y and RT-f (Supplementary Table S1). Various forms of active CK (kinetin, dihydrozeatin, cZ riboside, tZ riboside, meta-topolin riboside, 2-methylthio-cZ riboside, and dihydrozeatin ribonucleoside) were detected in different SE stages.

The relative content of ETH showed variable trends in both genotypes. In YT, it was detected only in “h” and “p” stages, whereas in case of RT, it was detected in “h” and “f” stages. The absence of ETH in YT-f and its high content in RT-f may suggest that it is negatively related to SE. Apart from these hormones, we also detected ten GAs including GA1, 3, 4, 7, 9, 15, 19, 20, 24, and 53 (Supplementary Table S1). Among the active GAs, GA1 and 7 had higher contents in YT-h and RT-h but were almost absent in other stages. The GA1 was also detected in YT-f and RT-f suggesting limited roles in SE. In the case of GA3, YT-f and RT-f showed accumulation but not in other stages. GA4 was detected only in YT-h and YT-p. These observations suggest that GAs have a limited role in the SE of Hevea.

Finally, we detected the differential accumulation of nine metabolites classified as JA or conjugates (Supplementary Table S1). Of these, we observed that YT-h and RT-h had higher levels of the major Jas: methyl jasmonate (MeJA), jasmonoyl-isoleucine (JA-Ile), and JA. The absence or detection of negligible amounts of JA and JA-Ile indicates limited role of JA in SE.

Transcriptome analysis of H. brasiliensis during SE

The RNA purity and integrity were evaluated, and overall high-quality RNA samples were used for transcriptome sequencing (Supplementary Table S2). The transcriptome sequencing analysis of 21 samples was completed, and a total of 186.64 Gb of clean data was obtained with an average of 7 Gb clean data per sample. A total of 1,244,105,756 clean reads were obtained after quality control. The percentage of Q30 bases was > 93% and the percentage of GC was 43% (Supplementary Table S2). Of the clean reads ~ 94% could be mapped to the H. brasiliensis genome. A total of 39,894 unigenes were identified and annotated using KEGG, NR, Swissprot, Trembl, KOG, GO, and Pfam databases. The overall gene expression levels of the different samples are shown as a box plot (Fig. 3A). PCA showed that all the replicates for each sample were close to each other. The PC1 and PC2 represented 20.95% and 32.71% variability, respectively (Fig. 3b). A heat map was generated to evaluate the uniformity of replicates for gene expression, and it was observed that the replicates were clustered together (Fig. 3c).

Fig. 3
figure 3

Transcriptome statistics. a FPKM distribution boxplot of each sample in this project. b Principal component analysis. c Heat map of differentially expressed transcripts. SE: Somatic embryogenesis, h: Immature Anther, y: Anther-induced callus, f: Callus differentiation stage, p: Cotyledonary embryo; YT: genotype- Y or Yunyan 73477; RT: RT or Reken 628

Differential gene expression in H. brasiliensis genotypes YT and RT during SE

A Venn diagram was used to illustrate overlapping or distinctly expressed genes at h, y, f and p stages of SE in both genotypes (Fig. 4). The top 10 most expressed genes in RT-h, RT-y, RT-f, YT-h, YT-y, YT-f, and YT-p are enlisted in Table 1. The Extensin-3-like gene was among the top 10 most expressed genes in all the samples. Therefore, it can be a good candidate for marker studies. Similarly, SRC1-like (LOC110638422) was uniformly expressed in all samples except for RT-f. The top-5 most expressed genes at YT-p (LOC110647147, LOC110659495, LOC110662110, LOC110635440) had very low or no expression in all stages, suggesting their essential role in the formation of cotyledonary embryos.

Fig. 4
figure 4

Summary of the DEGs during SE stages. a The Venn diagram of DEGs in SE stages in YT, b RT, c between both genotypes at relevant SE stage, d Number of differentially up-regulated and down-regulated genes in different stages of SE in both genotypes. SE: Somatic embryogenesis, h: Immature Anther, y: Anther induced callus, f: Callus differentiation stage, p: Cotyledonary embryo; Y: genotype- Y or Yunyan 73477; R: genotype-R or Reken 628

Table 1 List of top 10 highly expressed genes in the studied stages of YT and RT genotypes

In order to reveal the potential key factors and deeply understand the regulatory network of SE, the DEGs for each stage were identified (Fig. 4); only 1,947 and 2,471 genes were commonly expressed among all these stages in genotypes YT and RT, respectively. These might be the distinct set of genes in respective genotypes. A total of 5,208, 834, 762, 1,705, and 434 unique DEGs were identified for RT-h vs YT-h, RT-y vs YT-y, R-f vs YT-f, RT-f vs YT-p and YT-f vs YT-p, respectively. It is important to note that significant phenotypic variations (for embryogenesis) appear between both genotypes at primary embryo (y stage, Fig. 1). In YT, 2515 (up-regulated) and 3,303 (down-regulated) DEGs were identified for YT-y vs YT-f (Fig. 4d). Of these, 1,040 were exclusively expressed during YT-y vs YT-f (Fig. 4a). Similarly, 3,925 and 4,279 DEGs were respectively up- and down-regulated in YT-f vs YT-p. In contrast, a relatively higher number of DEGs (4,017 up- and 4,019 down-regulated) were identified in RT-y vs RT-f. Only 942 of these DEGs were uniquely expressed at this stage (Fig. 4b).

For further interpretation of gene function, KEGG pathway annotation analysis of DEGs was performed. The DEGs were assigned to 138 and 136 KEGG pathways in RT-h vs RT-y, and RT-y vs RT-f, respectively (Supplementary Figures S3 and S4). The most representative pathways include metabolic pathways (1890 unigenes), biosynthesis of secondary metabolites (1202 unigenes), plant hormones signal transduction (506 unigenes), and MAPK signaling pathway-plant (378 unigenes). On the other hand, the DEGs in YT-h vs YT-y, YT-y vs YT-f, and YT-f vs YT-p were assigned to > 130 KEGG pathways (Supplementary Figure S4), which were mostly the same as described above. Interestingly, the DEGs identified in the treatment comparisons between the two genotypes i.e., RT-h vs YT-h, RT-y vs YT-y, RT-f vs YT-f, RT-f vs YT-p, were also common to the individual list of KEGG pathways for each genotype. However, the representative number of genes/transcripts was different for each comparison.

Transcriptome analysis confirm the phytohormone accumulation trends

To further elucidate expression changes related to hormone biosynthesis and signaling regulation, we explored the differential expression of related genes (Supplementary Table S3). Tryptophan aminotransferase (TAA) and YUCCA (YUC) enzymes (like indole-3-pyruvate monooxygenase [EC:]) are the key genes involved in the auxin biosynthesis pathway. These genes showed relatively higher expression in YT as compared to RT in all four stages of SE (Fig. 5). From anther to callus induction, the differential expression of TAAs and YUCs was lower. However, compared to YT-h and YT-y, a significant increase in their expressions was observed in YT-f. This is the stage when the cells start to divide and undergo proliferation and differentiation. Moreover, the expression of these genes was higher in YT-f as compared to RT-f (Fig. 5). These changes in expression are consistent with the observed IAA accumulation in both genotypes.

Fig. 5
figure 5

Auxin biosynthesis and signaling pathway. a representation of auxin biosynthesis and signaling pathway genes and metabolites. b Heatmap representation of DEGs enriched in auxin biosynthesis and signaling pathway. Red: up-accumulation, Blue: down-accumulation, yellow: no variation of expression. SE: Somatic embryogenesis, h: immature anther, y: anther induced callus, f: callus differentiation stage, p: cotyledonary embryo

We further looked at the expression changes in the auxin signaling pathway (as per KEGG pathway map ko04075) and found eight auxin transporter-like (AUX1-like) proteins, two of which (LOC110647842, LOC110673399) were highly up-regulated in YT during callus induction (Fig. 5). On the other hand, LOC110633296, LOC110653578, and LOC110665570 were highly expressed in RT-f compared to YT-f. In addition, we observed upregulation of transport inhibitor response 1 (TIR1, LOC110673399, LOC110636674, LOC110640319, and LOC110643443) in YT-f and downregulation in RT-h. There were 34 DEGs annotated as AUX/IAA family proteins that showed variable expression trends. Generally, these transcripts were down-regulated during the “immature anther or h” stage, up-regulated during callus induction (YT-y and RT-y), and down-regulated again in later stages (Supplementary Table S3, Fig. 5).

Among 43 ARF DEGs, LOC110638638, LOC110640182, LOC110640184, LOC110650061, LOC110665593, LOC110667322 were down-regulated in transition from YT-h to YT-y and then up-regulated in next stage. Interestingly, these transcripts were mostly down-regulated in RT. In general, the FPKM values of these genes were higher in YT compared to RT. Apart from the ARFs, we also observed the expression changes in GH3 and SAUR family genes. The GH3s (LOC110631692, LOC110643427, LOC110656138, LOC110659801, LOC110664903) were up-regulated during YT-h vs YT-y and down-regulated in later stages. A similar trend was observed for in RT-h vs RT-y and the later stages. However, between both genotypes, these genes were down-regulated in YT compared to RT. Moreover, we identified 46 DEGs annotated as SAURs, of which LOC110633433, LOC110660228, LOC110661010 were the only DEGs not detected in other stages except for YT-f vs YT-p. Interestingly, the transcripts were also highly up-regulated in RT-f vs YT-p (Fig. 5), suggesting a potential involvement in transition from callus differentiation (f) to cotyledonary embryo formation (p) in YT. Important DEGs related to auxin pathway are resented in Fig. 5.

Together with endogenous auxin levels, it can be observed that in YT, IAA biosynthesis decreases in YT-y compared to YT-h, and then further increases in YT-f and YT-p stages, which can be attributed to the increased expression of the IAA biosynthesis genes TAA and YCC. Contrasting IAA accumulation and gene expression in RT, clearly suggest that IAA plays significant role in SE. Furthermore, the expression changes in IAA signaling genes confirm this proposition as detection of higher IAA changes the expression of downstream genes, thus suggesting that IAA could be a positive regulator of SE.

The biosynthesis of CKs starts with the production of isopentenyl diphosphate (IPP) and dimethylallyl diphosphate (DMAPP) through the mevalonic acid pathway or the non-mevalonic acid pathway [51]. The precursors IPP and DMAPP are then combined by an enzyme called isopentenyl transferase (IPT) to form isopentenyl adenine nucleotides, which are the primary forms of CKs in most plants. Our results did not show differential expression of the IPT gene (Fig. 6). However, there were two differentially up-regulated adenylate isopentenyl transferase (Aipt) transcripts (LOC110633565 and LOC110647477) in YT-y vs YT-f. The expression of LOC110633565 was higher in YT stages compared to those of RT. In addition to IPT and Aipt, other enzymes such as cytochrome P450 monooxygenase (also called Dwarf protein 11 or DWARF11; LOC110671515) were exclusive to RT (Fig. 6).

Fig. 6
figure 6

Cytokinin biosynthesis and signaling pathway. a representation of CK biosynthesis and signaling pathway genes and metabolites. b Heatmap representation of differentially expressed genes in CK biosynthesis and signaling pathway. Red: up-accumulation, Blue: down-accumulation, yellow: no variation of expression. SE: somatic embryogenesis, h: immature anther, y: anther-induced callus, f: callus differentiation stage, p: cotyledonary embryo

For the cytokinin signaling pathway (which is the part of main plant-hormone signaling pathway i.e., ko04075), eight DEGs annotated as CRE1 or histidine kinase 4 were identified. Only one DEG was found in the first two stages (Fig. 6). At this stage, the plant cells are still in a non-embryogenic state, and CK signalling is likely to be low. However, LOC110644203, LOC110657347, LOC110658259, and LOC110661394 were up-regulated during RT-y vs YT-y and RT-y vs YT-p, which is consistent with the higher CK content in YT-y (Fig. 2e). Eleven DEGs annotated as histidine-containing phosphotransfer protein 1-like (AHP) were identified. Of these, LOC110656915 and LOC110635512 showed down-regulation from YT-h to YT-p but did not differentially express in RT. However, both AHPs were up-regulated in YT-f as compared to RT-f. Our results showed the differential expression of 38 B-ARRs. Among these B-ARRs, the transcription repressor KAN1 (KANADI 1; LOC110652566, LOC110637165), MYB family transcription factor (TF) (PHL8-like; LOC110656263), Protein ABERRANT TESTA SHAPE (KANADI4; LOC110667027) represented the most up-regulated DEGs in YT with an increasing expression trend. However, these showed decreasing expression trends in RT. Like AHPs, the B-ARRs were up-regulated in YT-f as compared to RT-f (Fig. 6b). Of the 10 A-ARRs, ARR9 (LOC110666539) was down-regulated from YT-h to YT-p and up-regulated from RT-h to RT-y (Fig. 6). However, generally, its expression was higher in YT compared to RT. In summary, genes involved in CK signalling were globally up-regulated at callus differentiation stage (f) in YT compared to RT.

In case of ABA biosynthetic pathway (which is the part of carotenoid biosynthesis pathway in plants i.e., ko00906), a beta-carotene 3-hydroxylase (crtZ: LOC110643463) was down-regulated during callus induction stages (YT-y and RT-y), and up-regulated during callus differentiation stage (YT-f and RT-f) (Fig. 7). Zeaxanthin epoxidase (PA-ZE: LOC110662928) was up-regulated during the YT-f and RT-f stages. For 9-cis-epoxycarotenoid dioxygenase (NCED2), three transcripts (LOC110654913, LOC110641806, and LOC110633552) were up-regulated during YT-f and RT-f stages, but their expression was higher in RT-f compared to YT-f, which is consistent with the ABA content in these SE stages (Fig. 2).

Fig. 7
figure 7

ABA biosynthesis and signaling pathway. a Representation of ABA biosynthesis and signaling pathway genes and metabolites. b Heatmap representation of differentially expressed genes in ABA biosynthesis and signalling pathway. Red: up-accumulation, Blue: down-accumulation, yellow: no variation of expression

The ABA signaling pathway (which is the part of main plant-hormone signaling pathway i.e., ko04075) in plants involves a network of genes, including ABA receptors (PYR/PYL/RCARs), protein phosphatase 2Cs (PP2Cs), SNF1-related protein kinases 2 (SnRK2s), and ABA-responsive element binding factors (ABFs/AREBs) (Fig. 7). In the current study, 15 PYR/PYL genes were differentially expressed. The LOC110635302, LOC110633511, LOC110645564 represented the most significant DEGs. Here, LOC110635302 was the most up-regulated PYR/PYL in callus differentiation stages (f) of both genotypes (YT-y vs YT-f and RT-y vs YT-f). Whereas, of the 19 differentially expressed PP2Cs, five (LOC110636987, LOC110643225, LOC110644802, LOC110646778, and LOC110659164) genes were either down-regulated or weekly up-regulated in SE of YT. However, these were significantly up-regulated during RT-f. The SnRK2 transcripts (13 in total) demonstrated relatively lower differential expression. Though they had globally higher expression in YT compared to RT. Similarly, AREB TFs were significantly up-regulated in SE of YT as compared to RT (Fig. 7).

Differential expression of genes enriched in MAPK signaling pathways

We found 630 DEGs associated with 40 KO terms that could be enriched in the MAPK signaling pathways. Of these, 625, 602, 604, and 620 DEGs were expressed in YT-h, YT-y, YT-f, and YT-p, respectively. On the other hand, 627, 605, and 590 DEGs were expressed in RT-h, RT-y, and RT-f, respectively. Of the 630 DEGs, 254 were commonly enriched in the MAPK signaling and plant-hormone signaling pathways. These DEGs were associated with 13 KO terms related to ETH, JA, and ABA signaling. Since we observed increased ETH contents in RT-f compared to no detection in YT-f, therefore, we specially looked at the genes related to ETH signaling. In total, there were 73 transcripts associated with 13 KO terms enriched in ETH signaling. Of these, 38 genes, associated with 10 KO terms had higher FPKM values in RT-f compared to YT-f. These genes include ethylene response regulator (ERR), serine/threonine-protein kinase (SRK), ethylene insensitive 3-like (EIN3), EIN3-binding F-box protein 1-like (EBF1/2), ethylene-response factor 1B-like, copper-transporting ATPase (RANI1), MAPK3, chitinase-like protein 1 (CHiB1), and protein REVERSION-TO-ETHYLENE SENSITITY 1 (RTE1) (Fig. 8a-c). The higher expression of these genes is consistent with the increased ETH content detected in RT-f compared to YT-f. Additionally, we observed that several ETH signaling genes showed decreasing expression in successive SE stages in YT but not in RT. These include SRKs (STY8-like, LOC110662211, CTR1-like, LOC110662376, and HT1-like, LOC110662494), EIN3 (LOC110663488), ERF 1B-like (LOC110663726), RANI1 (HMA4-like, LOC110663743 and PAA2, LOC110664140), MAPK3 (LOC110664258), and CHiB2 (LOC110664259 and LOC110664469) (Supplementary Table S3).

Fig. 8
figure 8

Differential expression of genes involved in ethylene signaling (a, b, and c) and MAPK signaling (d, e, and f). Red: up-accumulation, Blue: down-accumulation, yellow: no variation of expression. SE: somatic embryogenesis, h: immature anther, y: anther induced callus, f: callus differentiation stage, p: cotyledonary embryo

In addition to ETH signaling, we also looked for the changes in the expression of genes related to oxidative stress as this is known to induce SE [52]. In the MAPK signaling pathway, the oxidative stress is perceived as H2O2, which then activates a signaling cascade of MAPKs, SRKs (OXI1), and downstream genes. We observed the changes in the expression of 45 MAPKs including MMK2, MAPK3, MAPKK (2, 3, 4, 9, 10), MAPKKK (1, 5, 17, 18, 20, and YODA), NTF3, 4, ANP1, and NPK1. Interestingly, the expression of ANP1 (LOC110668494), MAPKK2 (LOC110668244), MAPKK4 (LOC110633615), MAPKKK YODA (LOC110669673, LOC110669898, LOC110670175, LOC110670219), MAPKKK5 (LOC110669442, LOC110669547, and LOC110670173), MMK2 (LOC110667794 and LOC110667796), NTF3 (LOC110664204), and NTF4 (LOC110663028) increased from YT-h to YT-f but showed contrasting expressions in RT in the respective stages. However, these trends were not observed for all MAPKs. Particularly, we observed that MAPK3, MAPKKK1, MAPKKK18, and NPK1 had higher expressions in YT than RT in stage h but for other stages (y, f, and p) their expressions were higher in RT than YT. These expression trends indicate that oxidative stress plays important role in SE in H. brasiliensis, particularly the genes showing increasing expression trends may be positive regulators of SE (Fig. 8d-f).

Differential expression of somatic embryogenesis related genes

During SE, the WUS genes are expressed in embryogenic cells and are required to maintain their undifferentiated state. In comparative analysis of both genotypes, WUS gene expression was generally higher in later stages of YT compared to RT. The LEAFY COTYLEDON (LEC) TF has also been implicated in plant SE. Three transcripts of LEC genes (LOC110657040, LOC110636426, LOC110631489) were up-regulated during YT-y vs YT-f but were not differentially expressed in RT-y vs RT-f. This suggests a potentially important role for these genes in callus differentiation. Among the four differentially expressed LEA transcripts, LOC110631517 was only up-regulated during YT-f vs YT-p. Similarly, LOC110631895 was only up-regulated during RT-y vs RT-f (Fig. 9). Other than these genes, CLAVATA3 (CLV3) has been shown to play a role in regulating embryogenic competence and promoting the formation of somatic embryos. In the current study, CLV3 (LOC110632667) was not detected in YT and up-regulated in RT-y vs RT-f. Another CLV3 (LOC110649182) was highly up-regulated in RT-f vs YT-p (Fig. 9). Moreover, five SERK transcripts (LOC110651481, LOC110656809, LOC110657640, LOC110668667, LOC110670202) were generally down-regulated in SE in RT-y compared to RT-h. However, LOC110670202 (SERK2) was significantly up-regulated during YT-y vs YT-f and down-regulated during same stage (RT-y vs RT-f) in RT.

Fig. 9
figure 9

Expression profiling of SE related genes. a Gene expression analysis in YT and b) RT. c Comparative gene expression analysis during SE in both genotypes. Red: up-accumulation, Blue: down-accumulation, yellow: no variation of expression. SE: somatic embryogenesis; h: anther; y: embryogenic callus induction; f: callus differentiation for primary embryo; p: cotyledonary embryo

The GRDH (glucose and ribitol dehydrogenase, LOC110662067) was up-regulated only during YT-f vs YT-p. Three transcripts (LOC110645241, LOC110664260, LOC110665770) of LERKs (L-type lectin-domain containing receptors also called SE receptor kinase 1) were globally up-regulated in RT compared to YT. On the other hand, AGL15 transcripts were significantly up-regulated in RT-f vs YT-f. Furthermore, four BBM transcripts (LOC110642700, LOC110671423, LOC110635436, LOC110653362) were differentially up-regulated during the early SE stages and then their expressions decreased. Interestingly, BBM1 (LOC110642700) and BBM2 (LOC110635436) were up-regulated in RT-f vs YT-f. Similarly, LBD, CUP-SHAPED COTYLEDON, and FUSCA3 genes were globally significantly up-regulated in callus differentiation (YT-f) and somatic embryo (YT-p) stages. However, SD2-5 and IAA30 transcripts were up-regulated in RT-y and RT-f stages (Fig. 9). Overall, distinct gene expression profiles were observed for several SE-related genes in both genotypes, suggesting their potential roles in SE.

Changes in H2O2 content

Since the transcriptome expression suggested a possible role of oxidative stress in the differential SE of the two genotypes, we measured the H2O2 content. We observed that YT genotype had highest H2O2 content compared to RT. We observed that the content first decreased from “h” to “y” stages in both genotypes and then increased from “y” to “f” stages. However, the relative increase was higher in YT than in RT (Fig. 10). This validates the expression changes related to oxidative stress (Fig. 8d-f).

Fig. 10
figure 10

Relative H2O2 content in YT and HT genotypes during SE. Where h: anther; y: embryogenic callus induction; f: callus differentiation for primary embryo; p: cotyledonary embryo


The most commonly used method in Hevea spp. propagation is grafting onto seedling rootstocks, but the plants developed by this method show interclonal variability due to the inherent genetic variability of the rootstocks [53]. SE is an important tool used for the propagation of plant species of economic and agronomic significance. It has been considered as a key technology to shape the future of Hevea due to its advantages such as regaining ontogenetic juvenility and production of transgenic lines [5]. To efficiently utilize SE in Hevea, it is essential to explore the differences in SE potential of the available genotypes. A previous transcriptome-based study reported the involvement of phytohormone signaling, phenylpropanoid biosynthesis, and starch and sucrose metabolism pathways in Hevea SE [20]. However, limited or no data is available on the differential accumulation of phytohormones in different stages involved in SE and the corresponding changes in the expression of biosynthetic genes. We used the UPLC-ESI–MS/MS method to determine the differential accumulation of phytohormones in different stages of somatic embryo formation from anthers of two genotypes with different embryogenesis capacity. The genotype YT exhibited a higher callus induction compared to RT (Fig. 1), such differences could be due to the genetic background as reported in cotton [54], ryegrass [55], and Primula ssp. [56] etc. Further, the absence of SE in RT, lead us to hypothesize that this could be due to differences in endogenous phytohormone content at different SE stages in the two genotypes. This hypothesis was based on the known role of phytohormones in differential SE in many plant species such as Medicago truncatula [57], cotton [54], and Korean pine [58]. In Hevea SE, a three decade old study investigated that the presence of auxins such as 3,4-dichlorophenoxyacetic acid or CKs can suppress SE [59], which is contrary to recent studies highlighting that IAA signaling genes are activated in SEs [20]. Thus, we further explored the expression changes in genes related to biosynthesis and signaling of phytohormones including IAA, CKs, ET, and ABA.

Among the phytohormones, auxins are key regulators of plant cell division and cell differentiation, and are therefore widely used for callus and SE induction. Thus, increased endogenous auxin levels are related to regulation of SE process [60]. Our observations that IAA levels increased from YT-y to YT-f are consistent with significantly increased expression of TAA and YUC in YT-f (Fig. 2c), which are genes involved in IAA biosynthesis. Among these, YUC genes are transcriptionally regulated by BBM, which in turn increases IAA biosynthesis and auxin response, as recently reported in Arabidopsis [61]. In addition to BBM, YUCs are also direct targets of LECs. When LECs were activated in Arabidopsis seedlings they induce SE and auxin biosynthesis [62]. Since these genes showed higher expressions in YT compared to RT, it is possible that SE in Hevea is BBM and LEC mediated. These genes induce higher auxin biosynthesis during SE, which in turn triggers a signaling cascade that affects the expression of several SE-related genes and initiate vegetative-to-embryogenic transition. This auxin induced transition is WUS-mediated [63]. Our observations that WUS genes had higher expressions in YT-h compared to RT-h are consistent with these reports as well as other reports in M. truncatula [64], cotton [65], and coffee [66] (Fig. 9). These results indicate that auxin accumulates during SE progression in YT, but not in RT, and this accumulation is associated with changes in the expression of SE-related and auxin biosynthesis genes.

ETH plays a species-specific role in SE. Earlier studies on Daucus carota and Coffea canephora have shown that ETH reduces the SE potential [67, 68]. ETH overproducer Arabidopsis plants showed reduced free IAA levels [69]. This reduction was associated with the increased expression of AUX1 influx carriers that are responsible for increased auxin transport [70]. Our observations that ETH was not detected in YT-f, in contrast to relatively higher levels in RT-f, and the contrasting IAA levels (Fig. 2) in the two genotypes are consistent with the reports on D. carota and C. canephora. These observations suggest that ETH could be a possible cause of the absence of SE in the RT genotype. The relatively higher FPKM value of three AUX1/LAX genes and 38 ETH signaling related genes in RT-f compared to YT-f are also consistent with these findings (Supplementary Table S3) [70]. Thus, our results indicate that the genotype with no SE potential shows higher ETH content at the “f” stage compared to no content in YT at the corresponding stage. This higher ETH triggers downstream signaling genes. Hence, we propose that the absence of ETH in YT-f could be another possible cause of SE, as observed in black spruce, carrot, and coffee [67, 68, 71].

Studies have shown that ABA is associated with embryo development and maturation [72]. We also observed an increase in relative ABA levels from “y” to “f” stages of both genotypes (Fig. 2). In cotton, an appropriate concentration of ABA could promote SE [22]. Endogenous ABA levels are under the coordinated control of multiple factors, including its biosynthesis, transport, and stresses [73]. The increase in the expression of cZ from y to f stages together with higher expression of PA-ZE and NCED2 in YT-f and RT-f is consistent with the observed ABA levels (Figs. 2 and 7). PA-ZE converts zeaxanthin to violaxanthin, which is a key reaction in ABA biosynthesis [74], whereas NCED catalyzes the steps that convert neoxanthin to ABA [75]. The ABA accumulation triggers storage proteins and LEAs, suggesting that LEAs are components of ABA-inducible systems [76]. Since both genotypes had specific LEA transcripts that were differentially expressed (Fig. 9), it is premature to define and specify how ABA plays a role in SE in Hevea. However, the accumulation trends and expression changes clearly indicate that ABA plays a positive role in SE. Experiments specially targeting ABA in SE, and characterizing the key biosynthesis and signaling genes would significantly improve our understanding about SE in Heavea.

Auxins and CKs are known for their key roles as regulators of both cell division and differentiation in plants [77]. The significant increase in tZ levels in YT-y (Fig. 2) are consistent with earlier reports that CK is regulator of the shoot apical meristem [78]. This higher tZ content in YT could be due to the up-regulation of Aipt in YT compared to RT because this gene converts adenosine monophosphate into N6-(delta2-isopentenyl) adenosine-5'-monophosphate, which can be further phosphorylated and dephosphorylated to form various forms of CKs i.e., tZ [79]. Moreover, the differences in CK levels may also be due to the expression of different CK synthesis genes in YT and HT (Fig. 6) [80]. Nevertheless, the reducing levels of cZ and presence of higher quantities of tZ at the YT-y stage indicate their important role in SE. The higher CK content in earlier stages of SE observed in our results indicate that CK is required for cell proliferation, probably to increase mass [81]. Since from “y” to “f”, the callus starts dedifferentiation, and the fact that CKs, if produced constantly, they can produce irregular cell division. Hence the observed reduced content of CKs in latter SE stages is understandable [51]. The opposite auxin and CK accumulation trends suggest that CKs also regulate local auxin metabolism and modulate auxin pools in Hevea, as in other plants [82]. However, this needs further detailed exploration by modulating exogenous concentrations of the two hormones. The levels of CKs are perceived by CK signaling genes such as AHPs, which act as CK receptors [83], and the reduction in AHPs expression at YT stages is consistent with the decreasing trends in CK accumulation. After sensing, the CK response is mediated by A-ARRs (positive regulators) and B-ARRs (negative regulators) [84]. Thus the contrasting expression trends of genes annotated as A-ARRs and B-ARRs during SE are consistent with earlier reports [85]. Taken together, the accumulation trends of CKs and the expressions of biosynthesis and signalling genes indicate that CKs together with auxin, play significant role in SE. Particularly, our results that CKs and auxin had opposite accumulation trends are consistent with their antagonistic roles in plant development [51]. Since our results provide a broader understanding based on the expressed transcripts (genes) and accumulation trends, therefore, to underpin the specific role(s) of each of these phytohormones gene characterization studies would be needed.

Finally, in addition to stress stimuli, reactive oxygen species (ROS) are thought to be ubiquitous in normal cells and act as cellular messengers that can induce gene expression leading to SE in plants [52, 86]. In particular, hydrogen peroxide (H2O2) has been implicated in SE in recent years [87]. In the MAPK signaling pathway, H2O2 is perceived as an intracellular messenger and triggers a signal cascade. The up-regulation of MMK2, MAPK3, MAPKKs, MAPKKKs, NTF3, ANP1, and NPK1 from YT-h to YT-f indicates that H2O2 sensing is increased at this stage, inducing the expression of these signaling genes (Fig. 8). This was further confirmed by the observation that changes in H2O2 content were consistent to these expression changes (Fig. 9). Thus, our results provide preliminary understanding that an increased oxidative stress can also be a possible reason for SE in YT. Previous studies have shown a positive correlation between SE and H2O2 in Lycium barbarum [88]. Similarly, H2O2 has been reported to play a key role in conversion of embryogenic callus into somatic embryos [89]. The fact that these genes were significantly up-regulated in YT-f suggests that H2O2 may act as a messenger activating the expression of genes and the biosynthesis of proteins required to regulate the morphogenetic pathways leading to SE [90, 91].


We conducted a comparative transcriptome and phytohormone analysis using YT and RT Hevea genotypes with contrasting SE potential. The results indicate that YT has a higher callus induction rate (%) compared to RT and leads to SE, which is absent in RT. The UPLC-ESI–MS/MS system detected seven phytohormones of which IAA, CKs (tZ and cZ), ABA, and ETH showed differential accumulation trends. Phytohormone biosynthesis and signal transduction gene expression confirmed the UPLC-ESI–MS/MS results. Based on the two technologies applied, we conclude that the higher IAA content in YT compared to RT was associated with the higher expression of IAA biosynthesis, hence it could be a positive regulator of SE. Similarly, the detection of ETH in RT but not in YT at “f” stage and the corresponding gene expression changes may be related to the SE potential of the two genotypes. Finally, the accumulation trends of the CKs and the expressions of biosynthesis and signalling genes indicate that CKs play an important role in SE; most probably as a negative regulator. Our results also suggest the possible role of H2O2 in SE in Hevea.

Availability of data and materials

The raw sequencing data are available at NCBI SRA under the project number: PRJNA963065 (


  1. Priyadarshan P, de S. Goncalves P. Hevea gene pool for breeding. Genet Resour Crop Evol. 2003;50:101–14.

    CAS  Google Scholar 

  2. Mooibroek H, Cornish K. Alternative sources of natural rubber. Appl Microbiol Biotechnol. 2000;53:355–65.

    PubMed  CAS  Google Scholar 

  3. Clément-Demange A, Priyadarshan PM, Thuy Hoa TT, Venkatachalam P. Hevea Rubber Breeding and Genetics. In Plant Breeding Reviews, J. Janick (Ed.). 2007.

  4. Sobhana P, Gopalakrishnan J, Jacob J, Sethuraj M. Physiological and biochemical aspects of stock-scion interaction in Hevea brasiliensis. Indian J Nat Rubber Res. 2001;14:131–6.

    CAS  Google Scholar 

  5. Mignon E, Werbrouck S. Somatic embryogenesis as key technology for shaping the rubber tree of the future. Front Plant Sci. 1804;2018:9.

    Google Scholar 

  6. Wang Z, Zeng X, Chen C, Wu H, Li Q, Fan GA, Lu W. Induction of rubber plantlets from anther of Hevea brasiliensis Muell. Arg. in vitro. Chinese J Trop Crops. 1980;1:16–26.

    Google Scholar 

  7. Wu JL, Yang SQ, Hao BZ, Yuan XH, Xu LY. Characteristics related to higher rubber yield of Hevea brasiliensis Juvenile-type clone G11. J Rubber Res. 1998;125–32.

  8. Yang S, Mo Y. Some physiological propeptries of latex from somatic plants of Hevea brasiliensis. Chin J Trop Crops. 1994;15:13–20.

    Google Scholar 

  9. Wen L, Li W, Parris S, West M, Lawson J, Smathers M, Li Z, Jones D, Jin S, Saski CA. Transcriptomic profiles of non-embryogenic and embryogenic callus cells in a highly regenerative upland cotton line (Gossypium hirsutum L.). BMC Dev Biol. 2020;20:1–15.

    Google Scholar 

  10. Niskanen A-M, Lu J, Seitz S, Keinonen K, Von Weissenberg K, Pappinen A. Effect of parent genotype on somatic embryogenesis in Scots pine (Pinus sylvestris). Tree Physiol. 2004;24:1259–65.

    PubMed  Google Scholar 

  11. Trolinder NL, Xhixian C. Genotype specificity of the somatic embryogenesis response in cotton. Plant Cell Rep. 1989;8:133–6.

    PubMed  CAS  Google Scholar 

  12. Charbit E, Legavre T, Lardet L, Bourgeois E, Ferrière N, Carron M-P. Identification of differentially expressed cDNA sequences and histological characteristics of Hevea brasiliensis calli in relation to their embryogenic and regenerative capacities. Plant Cell Rep. 2004;22:539–48.

    PubMed  CAS  Google Scholar 

  13. Lardet L, Dessailly F, Carron M-P, Montoro P, Monteuuis O. Influences of aging and cloning methods on the capacity for somatic embryogenesis of a mature Hevea brasiliensis genotype. Tree Physiol. 2009;29:291–8.

    PubMed  Google Scholar 

  14. Piyatrakul P, Putranto R-A, Martin F, Rio M, Dessailly F, Leclercq J, Dufayard J-F, Lardet L, Montoro P. Some ethylene biosynthesis and AP2/ERF genes reveal a specific pattern of expression during somatic embryogenesis in Hevea brasiliensis. BMC Plant Biol. 2012;12:1–20.

    Google Scholar 

  15. Lardet L, Martin F, Dessailly F, Carron M-P, Montoro P. Effect of exogenous calcium on post-thaw growth recovery and subsequent plant regeneration of cryopreserved embryogenic calli of Hevea brasiliensis (Müll. Arg.). Plant Cell Rep. 2007;26:559–69.

    PubMed  CAS  Google Scholar 

  16. Su Y-H, Liu Y-B, Zhang X-S. Auxin–cytokinin interaction regulates meristem development. Mol Plant. 2011;4:616–25.

    PubMed  PubMed Central  CAS  Google Scholar 

  17. Radhakrishnan D, Kareem A, Durgaprasad K, Sreeraj E, Sugimoto K, Prasad K. Shoot regeneration: a journey from acquisition of competence to completion. Curr Opin Plant Biol. 2018;41:23–31.

    PubMed  CAS  Google Scholar 

  18. Kurepa J, Shull TE, Smalle JA. Antagonistic activity of auxin and cytokinin in shoot and root organs. Plant Direct. 2019;3:e00121.

    PubMed  PubMed Central  Google Scholar 

  19. Wang FX, Shang GD, Wang JW. Towards a hierarchical gene regulatory network underlying somatic embryogenesis. Trends Plant Sci. 2022;1209–17.

  20. Wang Y, Li H-L, Zhou Y-K, Guo D, Zhu J-H, Peng S-Q. Transcriptomes analysis reveals novel insight into the molecular mechanisms of somatic embryogenesis in Hevea brasiliensis. BMC Genomics. 2021;22:1–18.

    Google Scholar 

  21. Wickramasuriya AM, Dunwell JM. Global scale transcriptome analysis of Arabidopsis embryogenesis in vitro. BMC Genomics. 2015;16:1–23.

    CAS  Google Scholar 

  22. Jin F, Hu L, Yuan D, Xu J, Gao W, He L, Yang X, Zhang X. Comparative transcriptome analysis between somatic embryos (SE s) and zygotic embryos in cotton: evidence for stress response functions in SE development. Plant Biotechnol J. 2014;12:161–73.

    PubMed  CAS  Google Scholar 

  23. Salvo SA, Hirsch CN, Buell CR, Kaeppler SM, Kaeppler HF. Whole transcriptome profiling of maize during early somatic embryogenesis reveals altered expression of stress factors and embryogenesis-related genes. PLoS One. 2014;9:e111407.

    PubMed  PubMed Central  Google Scholar 

  24. Gao L, Zhang J, Hou Y, Yao Y, Ji Q. RNA-Seq screening of differentially-expressed genes during somatic embryogenesis in Fragaria× ananassa Duch.‘Benihopp.’ J Horticult Sci Biotechnol. 2015;90:671–81.

    Google Scholar 

  25. Indoliya Y, Tiwari P, Chauhan AS, Goel R, Shri M, Bag SK, Chakrabarty D. Decoding regulatory landscape of somatic embryogenesis reveals differential regulatory networks between japonica and indica rice subspecies. Sci Rep. 2016;6:1–15.

    Google Scholar 

  26. Rajesh M, Fayas T, Naganeeswaran S, Rachana K, Bhavyashree U, Sajini K, Karun A. De novo assembly and characterization of global transcriptome of coconut palm (Cocos nucifera L.) embryogenic calli using Illumina paired-end sequencing. Protoplasma. 2016;253:913–28.

    PubMed  CAS  Google Scholar 

  27. Yakovlev IA, Carneros E, Lee Y, Olsen JE, Fossdal CG. Transcriptional profiling of epigenetic regulators in somatic embryos during temperature induced formation of an epigenetic memory in Norway spruce. Planta. 2016;243:1237–49.

    PubMed  CAS  Google Scholar 

  28. Savona M, Mattioli R, Nigro S, Falasca G, Della Rovere F, Costantino P, De Vries S, Ruffoni B, Trovato M, Altamura M. Two SERK genes are markers of pluripotency in Cyclamen persicum Mill. J Exp Bot. 2012;63:471–88.

    PubMed  CAS  Google Scholar 

  29. Jamaluddin ND, Mohd Noor N, Goh H-H. Genome-wide transcriptome profiling of Carica papaya L. embryogenic callus. Physiol Mol Biol Plants. 2017;23:357–68.

    PubMed  PubMed Central  CAS  Google Scholar 

  30. Salaün C, Lepiniec L, Dubreucq B. Genetic and molecular control of somatic embryogenesis. Plants. 2021;10:1467.

    PubMed  PubMed Central  Google Scholar 

  31. Méndez-Hernández HA, Ledezma-Rodríguez M, Avilez-Montalvo RN, Juárez-Gómez YL, Skeete A, Avilez-Montalvo J, De-la-Peña C, Loyola-Vargas VM. Signaling overview of plant somatic embryogenesis. Front Plant Sci. 2019;10:77.

    PubMed  PubMed Central  Google Scholar 

  32. Li H-L, Wang Y, Guo D, Tian W-M, Peng S-Q. Three MADS-box genes of Hevea brasiliensis expressed during somatic embryogenesis and in the laticifer cells. Mol Biol Rep. 2011;38:4045–52.

    PubMed  CAS  Google Scholar 

  33. Srinivasan C, Liu Z, Heidmann I, Supena EDJ, Fukuoka H, Joosen R, Lambalk J, Angenent G, Scorza R, Custers JB. Heterologous expression of the BABY BOOM AP2/ERF transcription factor enhances the regeneration capacity of tobacco (Nicotiana tabacum L.). Planta. 2007;225:341–51.

    PubMed  CAS  Google Scholar 

  34. Mantiri FR, Kurdyukov S, Lohar DP, Sharopova N, Saeed NA, Wang X-D, VandenBosch KA, Rose RJ. The transcription factor MtSERF1 of the ERF subfamily identified by transcriptional profiling is required for somatic embryogenesis induced by auxin plus cytokinin in Medicago truncatula. Plant Physiol. 2008;146:1622–36.

    PubMed  PubMed Central  CAS  Google Scholar 

  35. Mantiri FR, Kurdyukov S, Chen S-K, Rose RJ. The transcription factor MtSERF1 may function as a nexus between stress and development in somatic embryogenesis in Medicago truncatula. Plant Signal Behav. 2008;3:498–500.

    PubMed  PubMed Central  Google Scholar 

  36. Chen L, Wu Q, He W, He T, Wu Q, Miao Y. Combined de novo transcriptome and metabolome analysis of common bean response to Fusarium oxysporum f. sp. phaseoli infection. Int J Mol Sci. 2019;20:6278.

    PubMed  PubMed Central  CAS  Google Scholar 

  37. Lee SI, Muthusamy M, Nawaz MA, Hong JK, Lim M-H, Kim JA, Jeong M-J. Genome-wide analysis of spatiotemporal gene expression patterns during floral organ development in Brassica rapa. Mol Genet Genomics. 2019;294:1403–20.

    PubMed  CAS  Google Scholar 

  38. Chen S, Zhou Y, Chen Y, Gu J. fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics. 2018;34:i884–90.

    PubMed  PubMed Central  Google Scholar 

  39. Kim D, Paggi JM, Park C, Bennett C, Salzberg SL. Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype. Nat Biotechnol. 2019;37:907–15.

    PubMed  PubMed Central  CAS  Google Scholar 

  40. Tang C, Yang M, Fang Y, Luo Y, Gao S, Xiao X, An Z, Zhou B, Zhang B, Tan X. The rubber tree genome reveals new insights into rubber production and species adaptation. Nat Plants. 2016;2:1–10.

    Google Scholar 

  41. Varet H, Brillet-Guéguen L, Coppée J-Y, Dillies M-A. SARTools: a DESeq2-and EdgeR-based R pipeline for comprehensive differential analysis of RNA-Seq data. PLoS One. 2016;11:e0157022.

    PubMed  PubMed Central  Google Scholar 

  42. Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J Roy Stat Soc Ser B Methodol. 1995;57:289–300.

    Google Scholar 

  43. Chen C, Chen H, Zhang Y, Thomas HR, Frank MH, He Y, Xia R. TBtools: An Integrative Toolkit Developed for Interactive Analyses of Big Biological Data. Mol Plant. 2020;13(8):1194–202.

  44. Kanehisa M, Goto S. KEGG: kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 2000;28:27–30.

    PubMed  PubMed Central  CAS  Google Scholar 

  45. Ashburner M, Ball CA, Blake JA, Botstein D, Butler H, Cherry JM, Davis AP, Dolinski K, Dwight SS, Eppig JT. Gene ontology: tool for the unification of biology. Nat Genet. 2000;25:25–9.

    PubMed  PubMed Central  CAS  Google Scholar 

  46. Kanehisa M, Sato Y. KEGG Mapper for inferring cellular functions from protein sequences. Protein Sci. 2020;29:28–35.

    PubMed  CAS  Google Scholar 

  47. Arakawa K, Kono N, Yamada Y, Mori H, Tomita M. KEGG-based pathway visualization tool for complex omics data. In Silico Biol. 2005;5:419–23.

    PubMed  CAS  Google Scholar 

  48. Konieczny R, Banaś AK, Surówka E, Michalec Ż, Miszalski Z, Libik-Konieczny M. Pattern of antioxidant enzyme activities and hydrogen peroxide content during developmental stages of rhizogenesis from hypocotyl explants of Mesembryanthemum crystallinum L. Plant Cell Rep. 2014;33:165–77.

    PubMed  CAS  Google Scholar 

  49. Chao J, Zhang S, Chen Y, Tian W-M. Cloning, heterologous expression and characterization of ascorbate peroxidase (APX) gene in laticifer cells of rubber tree (Hevea brasiliensis Muell. Arg.). Plant Physiol Biochem. 2015;97:331–8.

    PubMed  CAS  Google Scholar 

  50. Hua Y, Huang TA, Huang H. Micropropagation of self-rooting juvenile clones by secondary somatic embryogenesis in Hevea brasiliensis. Plant Breed. 2010;129:202–7.

    CAS  Google Scholar 

  51. Schaller GE, Bishopp A, Kieber JJ. The yin-yang of hormones: cytokinin and auxin interactions in plant development. Plant Cell. 2015;27:44–63.

    PubMed  PubMed Central  CAS  Google Scholar 

  52. Fehér A. Somatic embryogenesis—stress-induced remodeling of plant cell fate. Biochim Biophys Acta. 2015;1849:385–402.

    PubMed  Google Scholar 

  53. Priyadarshan P. Soil tillage, crop establishment and nutrition. In: Biology of Hevea rubber. Oxfordshire: CABI Wallingford UK; 2011. p. 126–133.

  54. Guo H, Guo H, Zhang L, Fan Y, Wu J, Tang Z, Zhang Y, Fan Y, Zeng F. Dynamic transcriptome analysis reveals uncharacterized complex regulatory pathway underlying genotype-recalcitrant somatic embryogenesis transdifferentiation in Cotton. Genes. 2020;11:519.

    PubMed  PubMed Central  CAS  Google Scholar 

  55. Grogg D, Rohner M, Yates S, Manzanares C, Bull SE, Dalton S, Bosch M, Studer B, Broggini GA. Callus induction from diverse explants and genotypes enables robust transformation of perennial ryegrass (Lolium perenne L.). Plants. 2022;11:2054.

    PubMed  PubMed Central  CAS  Google Scholar 

  56. Schween G, Schwenkel H-G. Effect of genotype on callus induction, shoot regeneration, and phenotypic stability of regenerated plants in the greenhouse of Primula ssp. Plant Cell Tissue Organ Cult. 2003;72:53–61.

    CAS  Google Scholar 

  57. Kępczyńska E, Orłowska A. Profiles of endogenous ABA, bioactive GAs, IAA and their metabolites in Medicago truncatula Gaertn. non-embryogenic and embryogenic tissues during induction phase in relation to somatic embryo formation. Planta. 2021;253:1–13.

    Google Scholar 

  58. Peng C, Gao F, Tretyakova IN, Nosov AM, Shen H, Yang L. Transcriptomic and metabolomic analysis of korean pine cell lines with different somatic embryogenic potential. Int J Mol Sci. 2022;23:13301.

    PubMed  PubMed Central  CAS  Google Scholar 

  59. El Hadrami I, Carron M-P, D’Auzac J. Influence of exogenous hormones on somatic embryogenesis in Hevea brasiliensis. Ann Bot. 1991;67:511–5.

    Google Scholar 

  60. Nic-Can GI, Loyola-Vargas VM. The Role of the Auxins During Somatic Embryogenesis. In: Loyola-Vargas, V., Ochoa-Alejo, N. editors. Somatic Embryogenesis: Fundamental Aspects and Applications. Cham: Springer; 2016.

  61. Li M, Wrobel-Marek J, Heidmann I, Horstman A, Chen B, Reis R, Angenent GC, Boutilier K. Auxin biosynthesis maintains embryo identity and growth during BABY BOOM-induced somatic embryogenesis. Plant Physiol. 2022;188:1095–110.

    PubMed  CAS  Google Scholar 

  62. Stone SL, Braybrook SA, Paula SL, Kwong LW, Meuser J, Pelletier J, Hsieh T-F, Fischer RL, Goldberg RB, Harada JJ. Arabidopsis LEAFY COTYLEDON2 induces maturation traits and auxin activity: implications for somatic embryogenesis. Proc Natl Acad Sci. 2008;105:3151–6.

    PubMed  PubMed Central  CAS  Google Scholar 

  63. Zuo J, Niu QW, Frugis G, Chua NH. The WUSCHEL gene promotes vegetative-to-embryonic transition in Arabidopsis. Plant J. 2002;30:349–59.

    PubMed  CAS  Google Scholar 

  64. Chen S-K, Kurdyukov S, Kereszt A, Wang X-D, Gresshoff P, Rose R. The association of homeobox gene expression with stem cell formation and morphogenesis in cultured Medicago truncatula. Planta. 2009;230:827–40.

    PubMed  PubMed Central  CAS  Google Scholar 

  65. Zheng W, Zhang X, Yang Z, Wu J, Li F, Duan L, Liu C, Lu L, Zhang C, Li F. AtWuschel promotes formation of the embryogenic callus in Gossypium hirsutum. PLoS One. 2014;9:e87502.

    PubMed  PubMed Central  Google Scholar 

  66. Arroyo-Herrera A, Ku Gonzalez A, Canche Moo R, Quiroz-Figueroa FR, Loyola-Vargas V, Rodriguez-Zapata L, Burgeff D′ Hondt C, Suárez-Solís V, Castaño E. Expression of WUSCHEL in Coffea canephora causes ectopic morphogenesis and increases somatic embryogenesis. Plant Cell Tissue Organ Cult. 2008;94:171–80.

    Google Scholar 

  67. Fuentes SR, Calheiros MB, Manetti-Filho J, Vieira LG. The effects of silver nitrate and different carbohydrate sources on somatic embryogenesis in Coffea canephora. Plant Cell Tissue Organ Cult. 2000;60:5–13.

    CAS  Google Scholar 

  68. Roustan J, Latche A, Fallot J. Stimulation of Daucus carota somatic embryogenesis by inhibitors of ethylene synthesis: cobalt and nickel. Plant Cell Rep. 1989;8:182–5.

    PubMed  CAS  Google Scholar 

  69. Lewis DR, Negi S, Sukumar P, Muday GK. Ethylene inhibits lateral root development, increases IAA transport and expression of PIN3 and PIN7 auxin efflux carriers. Development. 2011;138:3485–95.

    PubMed  CAS  Google Scholar 

  70. Vandenbussche F, Petrášek J, Žádníková P, Hoyerová K, Pešek B, Raz V, Swarup R, Bennett M, Zažímalová E, Benková E. The auxin influx carriers AUX1 and LAX3 are involved in auxin-ethylene interactions during apical hook development in Arabidopsis thaliana seedlings. Development. 2010;137:597–606.

    PubMed  CAS  Google Scholar 

  71. El Meskaoui A, Tremblay FM. Involvement of ethylene in the maturation of black spruce embryogenic cell lines with different maturation capacities. J Exp Bot. 2001;52:761–9.

    PubMed  Google Scholar 

  72. Walther M, Wagner I, Raschke J, et al. Abscisic acid induces somatic embryogenesis and enables the capture of high-value genotypes in Douglas fir (Pseudotsuga menziesii [MIRB.] Franco). Plant Cell Tiss Organ Cult. 2022;148:45–59.

  73. Dong T, Park Y, Hwang I. Abscisic acid: biosynthesis, inactivation, homoeostasis and signalling. Essays Biochem. 2015;58:29–48.

    PubMed  Google Scholar 

  74. Schwarz N, Armbruster U, Iven T, Brückle L, Melzer M, Feussner I, Jahns P. Tissue-specific accumulation and regulation of zeaxanthin epoxidase in Arabidopsis reflect the multiple functions of the enzyme in plastids. Plant Cell Physiol. 2015;56:346–57.

    PubMed  CAS  Google Scholar 

  75. Jiang C-C, Zhang Y-F, Lin Y-J, Chen Y, Lu X-K. Illumina® sequencing reveals candidate genes of carotenoid metabolism in three pummelo cultivars (Citrus maxima) with different pulp color. Int J Mol Sci. 2019;20:2246.

    PubMed  PubMed Central  CAS  Google Scholar 

  76. Chakraborty S, Roychoudhury A. Functional regulation of Responsive to abscisic acid (Rab) genes from representative plant species and their stress response. Plant Physiol Rep. 2022;27:653–64.

  77. Avilez-Montalvo JR, Quintana-Escobar AO, Méndez-Hernández HA, Aguilar-Hernández V, Brito-Argáez L, Galaz-Ávalos RM, Uc-Chuc MA, Loyola-Vargas VM. Auxin-cytokinin cross talk in somatic embryogenesis of Coffea canephora. Plants. 2013;2022:11.

    Google Scholar 

  78. Wu L-Y, Shang G-D, Wang F-X, Gao J, Wan M-C, Xu Z-G, Wang J-W. Dynamic chromatin state profiling reveals regulatory roles of auxin and cytokinin in shoot regeneration. Dev Cell. 2022;57:526-542 e527.

    PubMed  CAS  Google Scholar 

  79. Takei K, Sakakibara H, Sugiyama T. Identification of genes encoding adenylate isopentenyltransferase, a cytokinin biosynthesis enzyme, inArabidopsis thaliana. J Biol Chem. 2001;276:26405–10.

    PubMed  CAS  Google Scholar 

  80. Zhou Y, Tao Y, Zhu J, Miao J, Liu J, Liu Y, Yi C, Yang Z, Gong Z, Liang G. GNS4, a novel allele of DWARF11, regulates grain number and grain size in a high-yield rice variety. Rice. 2017;10:1–11.

    Google Scholar 

  81. Fehér A. Callus, dedifferentiation, totipotency, somatic embryogenesis: what these terms mean in the era of molecular plant biology? Front Plant Sci. 2019;10:536.

    PubMed  PubMed Central  Google Scholar 

  82. Chandler JW, Werr W. Cytokinin–auxin crosstalk in cell type specification. Trends Plant Sci. 2015;20:291–300.

    PubMed  CAS  Google Scholar 

  83. Sardesai N, Lee L-Y, Chen H, Yi H, Olbricht GR, Stirnberg A, Jeffries J, Xiong K, Doerge RW, Gelvin SB. Cytokinins secreted by Agrobacterium promote transformation by repressing a plant myb transcription factor. Sci Signal. 2013;6:ra100.

    PubMed  Google Scholar 

  84. Kim J. Phosphorylation of A-Type ARR to function as negative regulator of cytokinin signal transduction. Plant Signal Behav. 2008;3:348–50.

    PubMed  PubMed Central  Google Scholar 

  85. Yuan J, Chao Y, Han L. Uncovering a phenomenon of active hormone transcriptional regulation during early somatic embryogenesis in Medicago sativa. Int J Mol Sci. 2022;23:8633.

    PubMed  PubMed Central  CAS  Google Scholar 

  86. Blazquez S, Olmos E, Hernández JA, Fernández-García N, Fernández JA, Piqueras A. Somatic embryogenesis in saffron (Crocus sativus L.). Histological differentiation and implication of some components of the antioxidant enzymatic system. Plant Cell Tissue Organ Cult. 2009;97:49–57.

    CAS  Google Scholar 

  87. Prudente DDO, de Souza LB, Paiva R. Plant somatic embryogenesis: Modulatory role of oxidative stress. Proc Natl Acad Sci India Sect B Biol Sci. 2020;90:483–7.

    CAS  Google Scholar 

  88. Kairong C, Ji L, Gengmei X, Jianlong L, Lihong W, Yafu W. Effect of hydrogen peroxide on synthesis of proteins during somatic embryogenesis in Lycium barbarum. Plant Cell Tissue Organ Cult. 2002;68:187–93.

    Google Scholar 

  89. Cheng W-H, Wang F-L, Cheng X-Q, Zhu Q-H, Sun Y-Q, Zhu H-G, Sun J. Polyamine and its metabolite H2O2 play a key role in the conversion of embryogenic callus into somatic embryos in upland cotton (Gossypium hirsutum L.). Front Plant Sci. 2015;6:1063.

    PubMed  PubMed Central  Google Scholar 

  90. Gallego P, Martin L, Blazquez A, Guerra H, Villalobos N. Involvement of peroxidase activity in developing somatic embryos of Medicago arborea L. Identification of an isozyme peroxidase as biochemical marker of somatic embryogenesis. J Plant Physiol. 2014;171:78–84.

    PubMed  CAS  Google Scholar 

  91. Zhang W, Wang X-M, Rong F, Yin G-X, Ke W, Du L-P, Xiao L-L, Ye X-G. Effects of inter-culture, arabinogalactan proteins, and hydrogen peroxide on the plant regeneration of wheat immature embryos. J Integr Agricult. 2015;14:11–9.

    Google Scholar 

Download references


Not applicable.


This research was funded by the Sci-Tech Innovation System Construction for Tropical Crops Grant of Yunnan Province (RF2023-1), and The YITC Science fund for Distinguished Young Scholars (QNCZ2022-7).

Author information

Authors and Affiliations



Designed the experiments, writing—original draft preparation, review and editing, transcriptome and phytohormone test, data curation, analysis, and visualization, L.L.; performed experiments and carried out the phytohormone and transcriptomic data analysis, X.L.S. and W.C.Y.; Material preparation and morphological observation, M.C.G and Y.F.Q.; investigation, M.T and H.T.; project administration, G.P.L.

Corresponding author

Correspondence to Guoping Liang.

Ethics declarations

Ethics approval and consent to participate

All relevant institutional, national, and international guidelines and legislations were followed while conducting this experiment. The plant materials are obtained and conserved at the Genebank of Yunnan Institute of Tropical Crops under the number YITC4002 and YITC4009.

Consent for publication

Not applicable.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Supplementary Table 1.

Relative phytohormone based on UPLC-ESI-MS/MS in YT and HT genotypes at four stages of somatic embryogenesis of Hevea genotypes YT and RT.

Additional file 2: Supplementary Table 2.

Statistics of transcriptome sequencing of 21 libraries belonging to four stages of somatic embryogenesis of Hevea genotypes YT and RT.

Additional file 3: Supplementary Table 3.

Expression of genes enriched in MAPK signaling pathway in four stages of somatic embryogenesis of Hevea genotypes YT and RT.

Additional file 4: Supplementary Figure 1.

Details of metabolites identified related to IAA biosynthesis, storage and degradation.

Additional file 5: Supplementary Figure 2.

Metabolites related to CK biosynthesis and degradation.

Additional file 6: Supplementary Figure 3.

Classification of KEGG pathways to which the DEGs in RT-h vs RT-y and RT-y vs RT-f were enriched.

Additional file 7: Supplementary Figure 4.

Classification of KEGG pathways to which the DEGs in YT-h vs YT-y, RT-y vs RT-f, and YT-f vs YT-p were enriched.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Li, L., Sun, X., Yu, W. et al. Comparative transcriptome analysis of high- and low-embryogenic Hevea brasiliensis genotypes reveals involvement of phytohormones in somatic embryogenesis. BMC Plant Biol 23, 489 (2023).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: