Methylome and Transcriptome Analysis Reveals Differences in Callus Development and Plantlet Regeneration Capacity Between Two Eucalyptus Species

TL;DR

This study compares two Eucalyptus species to understand differences in their ability to regenerate plantlets through callus development using methylome and transcriptome data.

Contribution

The study identifies key genes and pathways involved in callus development and regeneration capacity differences between two Eucalyptus species.

Findings

Differentially expressed genes were enriched in plant hormone signal transduction and MAPK signaling pathways.

Genes like AHP, BAK1, and CRE1 are potentially involved in callus development and regeneration.

Observed differences provide a molecular basis for studying Eucalyptus regeneration mechanisms.

Abstract

Eucalyptus is a highly diverse genus of the Myrtaceae family that is planted worldwide. Many changes occur during callus development, an important process during in vitro plant regeneration. In this study, we conducted methylome and transcriptome analyses to reveal such changes. The results showed that differentially expressed genes between E. camaldulensis (voucher ID: c0009; high embryogenic potential) and E. grandis × urophylla (voucher ID: j0017; low embryogenic potential) during callus development were enriched in plant hormone signal transduction and MAPK (Mitogen-activated protein kinase) signaling pathways. qRT-PCR analysis showed AHP, BAK1, BSK, CRE1, GID1, MKS1, PR-1, PYL, RbohD, and TCH4 could be involved in the callus development and plantlet regeneration capacity. The differences observed in regenerative potential during callus maturation between the two species under study provide a reliable molecular basis for the study of Eucalyptus regeneration mechanisms.

Figures

Full text

1. Introduction

A highly diverse genus of the Myrtaceae family, Eucalyptus, is widely planted in China, Brazil, and India [1]. Eucalyptus trees are characterized by a high growth rate and large biomass production, an outstanding ability to grow in a wide range of environments and soils, and high wood quality for solid wood products and paper production [2]. Genus Eucalyptus consists of approximately 700 species, among which Eucalyptus camaldulensis Dehn (E. camaldulensis) and Eucalyptus urophylla × E. grandis have attracted considerable attention among researchers [3]. In particular, E. camaldulensis is widely cultivated for wood and timber, with approximately 20 million hectares planted globally [4]. These Eucalyptus species grow under rainfed conditions in the semiarid tropics and are known for their high degree of drought tolerance. Further, its widespread distribution and well-known ability to adapt to changing rainfall (250–1600 mm) and temperature regimes make it economically useful and industrially important [5]. In turn, Eucalyptus urophylla × E. grandis is a hybrid of Eucalyptus urophylla and Eucalyptus grandis that is characterized by high wood yield and carbon sequestration and is mainly planted in southern China [6]. The natural regeneration of Eucalyptus mainly relies on the seeds, and the breeding process is generally slow because of the length of the juvenile phase of the trees before flowering. Although natural hybridization of Eucalyptus has often been reported in association with the acquisition of certain advantageous characteristics (e.g., cold resistance), artificial hybridization is still a means of producing a limited number of seeds. Therefore, selected hybrids are commonly multiplied by vegetative means [7]. Vegetative propagation consists of both micropropagation (e.g., air layering, grafting, and rooting of cuttings) and micropropagation using in vitro tissue culture techniques, including adventitious budding, axillary shoot tips, and somatic embryogenesis (SE) [8]. Routine reproduction is a difficult and time-consuming process for Eucalyptus species because of their highly heterozygous genetic backgrounds [2].

Plants produce healing tissues when they experience stress such as wounding or pathogen infection. Further, the cells that make up healing tissues are totipotent and capable of regenerating the entire plant [9]. Thus, one of the most important aspects to consider during the establishment of an efficient Eucalyptus regeneration system is the regeneration of as many embryogenic calli as possible. Thus, callus development is a critical process during in vitro plant regeneration; indeed, it is a complex process that plays a key role while studying the early stages of the regulation of plant developmental and morphogenesis in a model system [10]. Further, many changes occur during callus development, including profound changes in RNA and DNA synthesis, pH, and oxygen uptake, as well as enzymatic activity enhancement, cytoskeleton transformation, and cytoplasmic inactivation, among others [11]. Recent studies showed that CIM-induced callus formation occurs from the pericycle or pericycle-like cells through a root developmental pathway. Compared to other tissues, callus tissue possesses a stronger regenerative capacity [12]. Our team’s preliminary research findings indicated that there were significant differences in the callus regeneration capacity between E. camaldulensis and E. grandis × urophylla, and initial investigations have been conducted into the molecular mechanisms [1].

Methylome and transcriptome analyses effectively capture the above changes. Thus, for example, researchers have explored the genetic and epigenetic bases of callus development in grapevines by profiling the transcriptome, epigenome, and small RNAome of undifferentiated, embryogenic, and non-embryogenic callus tissues derived from two genotypes differing in regeneration capacity: Sangiovese and Cabernet Sauvignon [13]. Similarly, MET1 (Methyltransferase 1)-dependent DNA methylation repression of light signaling and its effects on Arabidopsis regeneration was studied by comparing genome-wide DNA methylomes and transcriptomes of wild-type and met1-3 calli [14]. Furthermore, transcriptome and epigenome analyses revealed that the expression of defense pathways was enhanced concomitant with the alteration of DNA methylation patterns in the embryonic calli of the maize inbred line A188 [15]. CMT3 (Chromomethylase 3), DRM1/2 (Domain Rearranged Methyltransferase 1/2), and MET1 were identified as DNA methyltransferase-regulated genes participating in plant regeneration [14,16]. DEMETER (DME), a DNA glycosylase that removes DNA methylation as part of the base excision repair (BER) pathway, also influences the regeneration of Arabidopsis [17]. E. grandis × urophylla was a globally significant forest tree species renowned for its rapid growth, high yield, and exceptional wood production efficiency. Gene annotation and evaluation revealed that the E. urograndis genome contains 32,151 genes, of which 93.50% were fully annotated using Benchmarking Universal Single-Copy Orthologs (BUSCOs) [18]. Integrated transcriptomic and gibberellin analyses revealed that the plant MAPK signaling pathway is related to Eucalyptus urophylla development [19]. It has been reported that genome-wide identification and expression analysis of plant growth and development-related gene families were identified in Eucalyptus grandis [20]. Researchers reported a haplotype-phased genome reference for E. grandis, which was highly contiguous and nearly complete. This improved reference provides an important resource for studies focusing on genome sequence of the species, compared to earlier studies [21]. Although the role of DNA methylation in regeneration has been reported in model plants such as Arabidopsis, the epigenetic basis for differences in regenerative capacity among Eucalyptus species—particularly the regulatory mechanisms involving DNA methylation—remains unclear.

In our research, combined methylome and transcriptome analysis indicated that the MAPK signaling pathway regulated plant hormone signal transduction, leading to differences in callus development and plantlet regeneration capacity between two Eucalyptus species. Prominently featuring genes such as AHP (histidine-containing phosphotransmitter), BAK1 (BRI1-associated receptor kinase 1), BSK (Brassinosteroid-signaling kinase1), CRE1 (Cytokinin response 1), GID1 (Gibberellin insensitive dwarf 1), MKS1 (MAP kinase 4 substrate 1), PR-1 (Pathogenesis-related protein 1), PYL (Pyrabactin resistance 1-like), RbohD (Respiratory burst oxidase homologues D), and TCH4 (Touch-induced protein 4) played essential roles in governing Eucalyptus callus formation and regeneration capacity. Our findings enhanced the comprehension of Eucalyptus callus tissue development and regeneration, providing a foundation for further investigation.

2. Results

2.1. Callus Induction and Incubation

In our preliminary studies, we conducted in vitro tissue culture experiments on stem samples from two Eucalyptus species (Figure 1A) and obtained callus tissues at various developmental stages (Figure 1A–E). Notably, there were no significant differences in the appearance of the primary callus (pri-callus) (Figure 1B), mature callus (mat-callus) (Figure 1C), shoot regeneration stage callus (SRS-callus) (Figure 1D), or senescence callus (sen-callus) (Figure 1E) between the two species.

Sample sections obtained at each stage of tissue culture were stained and examined. The results showed no significant differences between the two species in stem tissues before incubation on CIM (Figure 1F) or in the pri-callus tissues (Figure 1G). When the mat-callus tissues were transferred to SIM and incubated for five (Figure 1H) and 10 d (Figure 1I), the regenerated callus tissues showed sprout germination and elongation. Notably, at both stages, the number of regenerated sprouts in E. camaldulensis was greater than that in E. grandis × urophylla (Figure 1H,I). In contrast, sen-callus tissues incubated on CIM for approximately 35 d showed signs of aging and loss of regenerative capacity, with no significant differences in histological morphology between the two species (Figure 1J). These distinct histological phenotypes, particularly the species-specific divergence in regenerative capacity during the shoot induction phase, prompted us to investigate the underlying epigenetic and transcriptional regulatory mechanisms.

2.2. Physiological Differences Between the Two Eucalyptus Species Under Study

Both species successfully produced primary and mature calli with similar appearance and induction rates. However, during the shoot regeneration stage, bud regeneration induction rate was significantly higher in E. camaldulensis than in E. grandis × urophylla, suggesting the likelihood of certain differences in gene expression and/or physiological indices between the two. The activation of cellular reprogramming of molecular pathways during CIM incubation is critical for the acquisition of callus cell totipotency [22], and endogenous and exogenous plant hormones can act as triggers of such cellular reprogramming [23]. We hypothesized that the species-specific shoot regeneration efficiency might be driven by distinct endogenous hormonal balances. To test this, we quantified the levels of in-dole-3-acetic acid (IAA), zeatin riboside (ZR), gibberellin A_3_ (GA_3_), and abscisic acid (ABA)—hormones central to auxin signaling, stress adaptation, and growth promotion—in mature calli prior to regeneration induction (Figure 1K).

HPLC analysis confirmed excellent linearity for all standard curves of hormones with R^2^ values exceeding 0.999 (Table S1). The retention times for IAA, ZR, ABA, and GA_3_ were 15.100 min, 5.170 min, 23.600 min, and 11.300 min, respectively. The results of the hormone content revealed that the levels of IAA, ABA, and GA_3_ in E. camaldulensis were significantly higher than those in E. grandis × urophylla (Figure 1K), suggesting that the differences in hormone detected in this study were likely related to the differences observed in the regenerative potential between the two species. Specifically, the higher IAA may promote cellular competency for reprogramming, elevated GA_3_ could enhance bud elongation, and increased ABA might modulate stress-responsive pathways that facilitate differentiation, together fostering a more efficient regeneration process in E. camaldulensis. Given that hormonal signaling can orchestrate cellular reprogramming through epigenetic modifications, we next investigated whether the observed hormonal and phenotypic disparities were associated with genome-wide DNA methylation changes.

2.3. Global DNA Methylation in Two Eucalyptus Species During Callus Maturation

To further explore whether methylation is involved in the regulation of callus maturation, we conducted WGBS in callus tissues from the two species at various stages of callus development (Table S2). The average unique mapping rates were 45.61% for E. grandis × urophylla and 26.72% for E. camaldulensis (Table S2). For the cytosine sites that were successfully mapped, the sequencing coverage was high (Table S3). The mean coverage depth across all samples was 13.32×, indicating that a substantial proportion of sites achieved coverage sufficient for reliable methylation calling (Table S3).

Cluster analysis of methylation levels across gene features (promoter, 5′ UTR, exon, 3′ UTR) revealed distinct methylation modification patterns between the two species during callus maturation (Figure 2). Analysis by sequence context confirmed methylation in CG, CHG, and CHH contexts in Eucalyptus (Figure 2A–C). A notable species-specific difference emerged in promoter methylation: E. grandis × urophylla exhibited significantly higher promoter methylation levels compared to E. camaldulensis (Figure 2). This hypermethylation in E. grandis × urophylla occurred predominantly in the CG context. In contrast, methylation levels in the 5′ UTR were relatively low and conserved across all samples.

To elucidate the functional implications of these global methylation differences, particularly the elevated promoter methylation in E. grandis × urophylla, we systematically identified differentially methylated regions (DMRs) and the associated genes (DMGs) between the two species across three key developmental stages (pri-, mat-, and SRS-callus), focusing on the CG context (Figure 3). Pathway enrichment analysis revealed distinct functional themes associated with differential methylation, dependent on genomic location. DMGs located within gene bodies were consistently enriched in pathways related to phenylpropanoid biosynthesis, flavonoid biosynthesis, flavone and flavonol biosynthesis, isoquinoline alkaloid biosynthesis, and tropane, piperidine, and pyridine alkaloid biosynthesis (Figure 3; Tables S4–S6, false discovery rate (FDR) < 0.05). Conversely, DMGs linked to promoter DMRs were enriched in plant–pathogen interaction and amino acid metabolism pathways (Figure 3; Tables S7–S9, false discovery rate (FDR) < 0.05). Notably, the DMR-related genes across all three stages were primarily enriched in sesquiterpenoid and triterpenoid biosynthesis pathways (Figure 3; Tables S4–S9).

2.4. Combined Analysis of Methylome and Transcriptome

The interplay between epigenetic regulation and gene expression has been confirmed by numerous studies. To further explore the regulatory relationship between DNA methylation and gene expression, the genes with both hypermethylation and downregulated expression, and genes with both hypomethylation and upregulated expression in E. grandis × urophylla compared to E. camaldulensis were subjected to KEGG pathway enrichment analysis (Figure 4). Notably, in SRS-callus, genes that were hypomethylated and upregulated showed significant enrichment in the Plant hormone signal transduction pathway (Figure 4). During the primary and mature callus stages, the plant–pathogen interaction pathway was concurrently enriched among two contrasting groups: genes that were hypermethylated and downregulated, as well as those that were hypomethylated and upregulated (Figure 4). In SRS-callus, the plant–pathogen interaction pathway remained significantly enriched specifically in the group of genes that were hypermethylated and downregulated (Figure 4).

As for plant hormone signal transduction, genes involved in the biosynthesis of zeatin, diterpenoid, carotenoids, and brassinosteroids, and those in phenylalanine metabolism, differed significantly in promoter methylation levels and gene expression between the two species (Figure 5). The gene CRE1 was downregulated in both methylation and expression levels in E. camaldulensis during the A2 and A3 stages, compared to E. grandis × urophylla, whereas the opposite occurred at the A4 stage (Figure 5). The AHP gene showed downregulated methylation and expression levels, specifically at the A2 stage, with no significant differences at the other stages (Figure 5). The BAK1 gene showed downregulated methylation and expression levels during the A3 stage. The BSK gene showed upregulated methylation and expression at A2, while the TCH4 gene showed downregulation in methylation with upregulation in expression during A4 (Figure 5). The MKS1 gene showed downregulation in both methylation and expression levels during A2, whereas the RbohD gene showed upregulated methylation and downregulated expression levels at A2, but upregulated methylation and expression levels at the A4 stage (Figure 5).

To independently verify the expression patterns of key regulatory genes identified by our integrated methylome–transcriptome analysis (particularly those within the plant hormone signal transduction and MAPK signaling pathways showing coordinated methylation and expression differences), qRT-PCR was performed on ten candidate genes (AHP, BAK1, BSK, CRE1, GID1, MKS1, PR-1, PYL, RbohD, and TCH4) across callus developmental stages 2, 3, and 4 in E. camaldulensis and E. grandis × urophylla. Using E. grandis × urophylla as control, Figure 6 showed the ratio of relative expression level between E. camaldulensis and E. grandis × urophylla. The qRT-PCR results were largely concordant with the RNA-seq-derived expression trends (Figure 5), thereby providing independent experimental confirmation of the stage-specific, species-dependent differential expression of these epigenetically regulated candidate genes.

3. Discussion

As part of a fundamental technique in genetic engineering, callus development has great potential for epigenetic variation, for which DNA methylation for genome activity is paramount. In this study, we performed WGBS of callus tissues at various stages of callus development in the two Eucalyptus species under study. Revealing the potential epigenetic mechanisms involved in callus development could improve our understanding of the regeneration process of plant cells, which will be beneficial for overcoming regeneration recalcitrance in Eucalyptus. Consistent with previous reports on Brassica napus [24], Beta vulgaris L. [25], and Oryza Sativum [26], our results showed that genome-wide hypomethylation was predominant during callus development (Figure 2), underlining the importance of regulating DNA methylation during this process. Particularly, low methylation levels during callus development may be induced at the tissue culture phase [26]. In vitro culture of callus was accompanied by DNA hypomethylation and H3K27me3 (histone H3 at lysine 27) demethylation, which could activate auxin- and cytokinin-related regulators to induce callus development. The DNA methylation inhibitor, 5-azacytidine, could significantly increase callus development [27]. Moreover, the promoter of LEAFY COTYLEDON1 (LEC1), a master regulator under the control of multiple epigenetic layers, undergoes hypomethylation before somatic embryo formation, suggesting that DNA methylation contributes to plantlet regeneration [28]. Furthermore, the coordination of multiple chromatin modifiers during somatic embryogenesis is also closely associated with the methylation modification [29].

In this study, integrated methylome–transcriptome identified hormone signaling as a key regulator of callus development and plantlet regeneration capacity between two Eucalyptus species. We constructed a comprehensive regulatory model linking differences in hormone content to differences in promoter methylation, gene expression, and regenerative capacity between two Eucalyptus species, E. camaldulensis, and E. grandis × urophylla. We found that sesquiterpenoid and triterpenoid biosynthesis signaling pathways may be involved in callus development through the regulation of methylation (Figure 3). Further, sesquiterpenoids and triterpenoids have been recognized as precursors of endogenous hormones such as ABA, JA, and IAA; moreover, growth and differentiation in callus cultures are controlled by interactions between exogenous and endogenous hormones [30]. Hormones and differentiation of callus are relevant to methylation levels during cullus culture, which has been confirmed in the in vitro culture of Brassica napus [24] and garlic tissues [31]. High auxin induced callus through a histone methylation-dependent mechanism. This epigenetic regulation, in turn, was required for the transcriptional activation of these genes during callus formation [32]. The massive transcriptional reprogramming for cell fate transition by auxin during callus formation requires epigenetic modifications [33]. Indeed, endogenous hormone levels have large effects on physiological processes, plant architecture, and initiation of proliferation centers in explants [31]. Additionally, they are frequently used to control cell division [34]. Specifically, ABA and JA are predominantly associated with callus induction [9], while a correlation between IAA levels and callus propagation and maintenance has been proved [35]. Plant–pathogen interactions help plants defend themselves against pathogens and other external adverse stimuli [36]. In our study, the DMR-related genes in the promoter regions of the two species were primarily enriched in plant–pathogen interactions across the three callus development stages analyzed herein. This finding indicates that the plant–pathogen interaction signaling pathway regulated callus regeneration via methylation, in agreement with results reported for Moringa oleifera Lam [37], Taxus media [38], Manihot esculenta [39], and Picea balfouriana [40].

Our research group has been working on Eucalyptus calli for several years. Previously, we reported key miRNAs [41], genes [18], and proteins [42] involved in plantlet regeneration in the two Eucalyptus species: E. camaldulensis (with high embryogenic potential) and E. grandis × urophylla (with low embryogenic potential). Transcriptome and proteome integration analysis has revealed that DEGs/DEPs are predominantly enriched in the plant hormone signal transduction and MAPK signaling pathways [43]. This finding is consistent with the results reported herein. Combined analysis of the methylome and transcriptome has revealed that DEGs are mostly enriched in plant hormone signal transduction and MAPK signaling pathways, which could be a critical factor regulating the regenerative capacity of Eucalyptus callus tissue. Indeed, plant hormone signal transduction is indispensable for callus development [44], and it is widely accepted that plant hormone signaling-related genes and methylation levels are critical for cell differentiation in different plants [45]. Thus, for example, transcriptome and endogenous hormone analyses showed that DEGs were significantly enriched in plant hormone signal transduction pathways during callus maturation in Osmanthus fragrans [46], similar to what has been reported for Zea mays L. [47], Picea spruce [39], and Agapanthus praecox [48].

Mitogen-activated protein kinase (MAPK) cascades (MAPK-MAPKK-MAPKKK) were important signaling modules in eukaryotes that function downstream of sensors/receptors to coordinate cellular responses to achieve normal growth and development of organisms and their adaptation to the ever-changing environment [49]. MAPK cascade was conserved in all eukaryotic organisms. The MAPK cascades could be regulating N assimilation, which enhanced plant growth and development [50]. The genetic expressions of MAPKs and the MAPKK were slightly regulated by N. However, the genetic expressions of MAPKKKs RAF14 and RAF79 showed a very strong repression by ammonium, which suggests that they may have a key role in the regulation of N assimilation [51]. Further, the MAPK signaling pathway is involved in the biosynthesis and/or signaling of ethylene, JA, ABA, and other plant hormones [1]; it is an upstream regulator that controls hormone biosynthesis and transport, and MAPKs play a critical role in hormone signaling as downstream regulators [52]. Indeed, the MAPK signaling pathway plays an essential role in plant growth and development through its control over phytohormones [1]. In our study, we found that plant hormone signal transduction and MAPK signaling pathways consistently co-occurred in plants (Figure 4), likely due to the regulation of plant hormone signal transduction by the MAPK signaling pathway, which is consistent with results obtained for Setaria italica [53] and Oryza sativa [54].

qRT-PCR analysis showed the relative expression of CRE1, BAK1, BSK, TCH4, MKS1, and RbohD. The differential expression patterns of the above genes across Stage 2 to Stage 4 suggest that they could be involved in the callus development and regeneration capacity of Eucalyptus callus. CRE1 (Cytokinin Response 1), which was identified as a kind of cytokinin receptor [55], played an important role in de novo organ formation and growth regulation [56]. The cytokinin signaling cascade is initiated upon ligand binding to the CRE1’s N-terminal cyclase/histidine kinase-associated sensor extracellular chase domain, leading to a conformational change, allowing for ATP binding in the histidine kinase domain [57]. Cytokinin promoted early development by upregulating MAPK6, a member of the MAPK signaling pathway [58]. Our research found that the MAPK signaling pathway and BAK1 participated in callus development and plantlet regeneration of Eucalyptus. Plant growth hormone triggered heterodimerization with the BRI1-Associated Receptor-Like Kinase 1 (BAK1) (also known as Somatic Embryogenesis Receptor Kinase 3, SERK3) and its closest family members SERK1 and SERK4, which further relayed the signaling to a MAPK cascade, to regulate plant development [59]. Brassinosteroid Signaling Kinase (BSK) was detected as a phosphorylation product of BAK1, which activated the mitogen-activated protein kinase (MAPK) pathway upon embryogenesis of Arabidopsis [60]. It was consistent with the results we compared between E. camaldulensis and E. grandis × urophylla callus. BAK1 activates the brassinosteroid signaling pathway through phosphorylation, promoting callus regeneration [61]. TCH4 was a xyloglucan endotransglucosylase/hydrolase (XTH) family member. Extensive studies have shown that XTHs were very important in cell wall homeostasis for plant growth and development [62]. The relative high expression of TCH4 in the fourth stage of callus tissue may be related to the active synthesis of cell wall components. Cell wall functioned as a foundational barrier and scaffold supporting plant cells, which was closely related to the high regeneration rate of callus tissue [63]. In our research, the expression levels of CRE1 and TCH4 increased during callus developmental stages 2–4. BAK1 expression levels remained largely unchanged during the first and fourth stages, while showing a slight decrease during the third stage (Figure 6). CRE1, BAK1, and TCH4 were simultaneously enriched in the plant hormone signaling pathway, promoting metabolic and developmental processes [64]. In contrast, BSK and GID1 exhibited a decreasing expression pattern from Stage 2 to Stage 4 (Figure 6), potentially participating in the regulation of callus regeneration by enriching the redox pathways and plant hormone signaling pathways [65]. Furthermore, our reserch showed MKS1 and RbohD also played a critical role in regulating the callus development and regeneration capacity of Eucalyptus. MAP kinase 4 nuclear substrate 1 MKS1 gene, which was involved in the MAPK signaling pathway, showed downregulation in both methylation and expression levels during A2. Mitogen-activated protein kinase MAPK4 was predominantly localized in nuclei, where it interacts with a substrate protein—MKS1. As may be expected for a kinase substrate, MKS1 under-expression partially suppresses these MAPK4 mutant phenotypes, affecting plant development [66]. RbohD (Respiratory burst oxidase homologues D) gene exhibited distinct expression patterns during the A2 and A4 periods. RbohD promoted reactive oxygen species (ROS) burst to further enhance the long-term MAPK activation [67]. This also explained our experimental findings: differentially expressed genes during callus development were enriched in the MAPK signaling pathway. The expression levels of RbohD increased during callus developmental stages 2–4. PYL exhibited an expression pattern that first increased and then decreased. During Eucalyptus callus development, RbohD and PYL respond to the regeneration process through plant hormone signaling pathways and MAPK signaling pathway [68]. The synergistic regulation among genes was a critical factor underlying the differences in callus regeneration capacity between E. camaldulensis and E. grandis × urophylla.

This study revealed that methylation-mediated regulation of plant hormone signaling and the MAPK pathway served as a key mechanism underlying the variation in regeneration capacity in Eucalyptus. Integrated methylome and transcriptome analyses demonstrated that genes such as AHP, BAK1, BSK, and CRE1 were predominantly enriched in plant hormone signal transduction and the MAPK signaling pathway, which played critical roles in callus development and regeneration ability. These findings not only advanced the understanding of callus development in Eucalyptus but also provided a strategic direction for future regeneration improvement—specifically, through targeted modulation of methylation to enhance the efficiency of hormone and MAPK signaling. Furthermore, functional validation of these identified genes will be an essential next step to clarify their precise roles in regeneration and to support the development of targeted biotechnological strategies for optimizing Eucalyptus regeneration.

4. Materials and Methods

4.1. Plant Material and Culture Conditions

The original seeds of E. camaldulensis and E. grandis × urophylla were obtained from the wild in 1984 and were planted in the experimental fields of Guangxi Forestry Research Institute. Bud strips from E. camaldulensis individual (voucher ID: c0009) and E. grandis × urophylla individual (voucher ID: j0017) were collected for in vitro cultivation. Then, plantlets of two Eucalyptus species obtained via in vitro tissue culture were used in this study. To investigate the molecular basis of their contrasting regenerative potentials, we sampled tissues across a developmental series designed to capture the acquisition, manifestation, and loss of embryogenic competence: the initial explant (Stage 1), early and mature callus stages under induction conditions (stages 2–3), the regenerative stage upon hormonal trigger (Stage 4), and a senescent negative control (Stage 5). Samples were systematically labeled with a species prefix (‘A’ for E. camaldulensis, ‘B’ for E. grandis × urophylla) and a stage number suffix indicating stages (e.g., A1/B1 = Stage 1). Specifically, stem tissues (Stage 1) of E. camaldulensis (A1) and E. grandis × urophylla (B1) were extracted from in vitro tissue culture-induced seedlings of both species. Primary calli (pri-callus/Stage 2) were generated by incubating Stage 1 tissues (A1/B1) for 10 d in callus induction medium (CIM), MS medium supplemented with 20 mg/L calcium nitrate, 1 mg/L kinetin (KT), and 0.5 mg/L 2,4-dichlorophenoxyacetic acid (2,4-D). These were further incubated in CIM for an additional 11 d to obtain mature calli (mat-callus/Stage 3) samples A3/B3. Upon transferring to shoot-induction medium (SIM), which is MS medium supplemented with 20 mg/L calcium nitrate, 2.0 mg/L benzylaminopurine (6-BA), and 0.2 mg/L naphthaleneacetic acid (NAA), and incubating for 5–10 d, the tissues developed into the SRS-callus (SRS-callus/Stage 4; samples A4/B4), characterized by the emergence of buds, which is indicative of successful propagation. Conversely, the mat–calli cultured on CIM for an additional 14 d resulted in senescent calli (sen-callus/Stage 5; samples A5/B5), which completely lost their embryogenic capacity.

4.2. Morphological and Histological Observations

Paraffin sections were obtained from stem tissue, primary callus, and sen-callus at 0, 10, and 35 d after inoculation to CIM. Additionally, for the SRS-callus, mature calluses inoculated on SIM and cultured for 5 and 10 d were sectioned. Observations were performed on biological triplicates for each species at each stage (i.e., three independent tissue samples).

For histological analysis, tissues from key developmental stages (stem tissue, primary callus, and SRS-callus) were fixed in FAA fixative (70% ethanol:glacial acetic acid:formaldehyde, 18:1:1) at 4 °C. Following fixation, samples were processed through a standard dehydration and clearing series using ethanol and xylene, embedded in paraffin, and sectioned at 10 μm thickness using a rotary microtome. Sections were stained with an aniline blue solution [69] and imaged under a Motic BA410 optical microscope (Motic Instruments Inc., Xiamen, China).

4.3. Determination of Endogenous Hormone Contents

Callus samples from three biological replicates per species and stage (mature callus, Stage 3) were used for hormone quantification. Briefly, callus samples (0.1 g) were thoroughly ground under liquid nitrogen, mixed with 1 mL of pre-cooled methanol (80% v/v) and left to extract overnight at 4 °C in the dark. The extraction mixture was centrifuged at 8000 g for 10 min, and the supernatant was collected. Following extraction, 0.5 mL of 80% (v/v) methanol was added to the residue, and another 2 h extraction at 4 °C was conducted. The mixture was then centrifuged and the resulting supernatants were combined. The supernatants were dried under nitrogen at 40 °C until no methanol remained. The remaining aqueous phase was extracted three times with petroleum ether (0.5 mL) and the upper ether phase was discarded. The aqueous phase was adjusted to pH 2.8 and then extracted three times with 1 mL ethyl acetate each time. The combined organic phase was dried at 40 °C under a nitrogen stream. The residue was dissolved with 0.5 mL methanol and filtered with a 0.45 μm pinhead filter for measurement. Chromatographic detection was performed using an ACCHROM high-performance liquid chromatograph with an Alphasil VC-C18 (4.6 mm × 250 mm, 5 μm) column (Dalian ACCHROM Science & Technology Co., Ltd., Dalian, China). The mobile phase was methanol: water contained 1% acetic acid (v/v) in a 40:60 ratio, and the flow rate was 0.8 mL min^−1^. The injection volume was 10 μL, the column temperature was maintained at 35 °C, and UV detection wavelength was set at 254 nm. Standard compounds GA_3_, ABA, ZR, and IAA were dissolved in methanol at a concentration of 100 mg/L to yield stock solutions. Working standard solutions at gradient concentrations were prepared through serial dilution of the stock solution. Subsequently, working standard solutions were analyzed under the chromatographic conditions described above to establish the standard curve, from which the regression equation and correlation coefficient (R^2^) were calculated. Levels of the hormones, including GA3, ABA, ZR, and IAA, were compared between the two Eucalyptus species based on Student’s t-test.

4.4. High-Throughput DNA Methylation Sequencing

We extracted whole-genome DNA from callus tissues at various developmental stages using the super plant genomic DNA kit (polysaccharides and polyphenolics-rich) (Tiangen Biotech (Beijing) Co., Ltd., Beijing, China), followed by quality control, including agarose gel electrophoresis to analyze the degree of DNA degradation and the presence of RNA contamination, Nanodrop measurement of DNA purity based on the A260/280 ratio, and Qubit for precise quantification of DNA concentration. After quality checks, a certain proportion of negative control lambda DNA was added to the DNA samples. Genomic DNA was then randomly sheared to 200–300 bp using a Covaris S220 (Covaris, Inc., Woburn, MA, USA). Sheared DNA fragments were subjected to end repair, adenylation, and ligation of sequencing adapters. Subsequently, DNA was subjected to bisulfite treatment using the EZ DNA Methylation Gold Kit (Zymo Research Corporation, Irvine, CA, USA), followed by PCR amplification to produce the final DNA library and quantification using Qubit 2.0 (Thermo Fisher Scientific, Waltham, MA, USA). The library was diluted to a concentration of 1 ng μL^−1^, and the insert size was then assessed using an Agilent 2100 (Agilent Technologies, Inc., Santa Clara, CA, USA). Qualified libraries were sequenced on an Illumina HiSeq 4000 platform (Illumina, Inc., San Diego, CA, USA) with paired-end 2 × 150 bp reads. Each callus tissue from each species at each developmental stage was subjected to whole-genome bisulfite sequencing (WGBS); three biological replicates were included.

Data quality control was performed using fastqc (v0.11.5, https://www.bioinformatics.babraham.ac.uk/projects/fastqc/, accessed on 10 June 2023) and Trimmomatic v0.36 (Usadellab, Aachen, Germany) [70] to obtain clean data, which were then aligned to the Eucalyptus reference genome (v2.0, https://plantgenie.org, accessed on 10 June 2023) using Bismark v0.16.3 [71]. Additionally, Bismark (Babraham Institute, Cambridge, UK) was used for de-duplication and methylation (mC) calling. The lambda genome was used for calculations, and the reliability of the methylation sites was assessed in subsequent analyses. A binomial test, B (n, p), was conducted for each C site to identify reliable methylation sites. Assuming x methylated Cs at a given site with n read coverages and a bisulfite non-conversion rate p, we verified the reliability of x methylated Cs occurring under the conditions of probability p, and sequencing depth, n. A set of thresholds [72,73] were applied during the analysis to accurately detect methylation sites: (1) sequencing depth ≥ 5; (2) q-value ≤ 0.01 [74]. The level of methylation of the identified methylation sites was calculated using the formula: ML = mC/(mC + umC), where ML represents the level of methylation, and mC and umC denote the counts of methylated and unmethylated Cs, respectively.

Differentially methylated regions (DMRs) were identified using DSS v2.12.0 software (Emory University, Atlanta, GA, USA) with the Bayesian hierarchical model [75,76,77]. Annotations for various functional genomic regions, such as promoters, exons, introns, CpG Islands (CGIs), CGI shores, and repeats, were applied to the DMRs, with the promoter region defined as an area 2kb upstream of the transcription start site (TSS) [78]. We focused on the CG context because extensive remodeling of CG methylation in gene regions occurs after callus regeneration and shoot formation [79]. Based on the genomic distribution of DMRs, genes with overlaps between the gene body (from TSS to the transcription end site, TES) and DMRs in the CG context were defined as DMR-related genes in the gene body, whereas those with overlaps in the promoter region were defined as DMR-related genes in the promoter. Enrichment analysis of gene ontology (GO) functions and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways were conducted for the aforementioned genes separately.

4.5. Integrated Analysis of WGBS and RNA-Seq Data

In preliminary studies, we conducted transcriptome sequencing and differential expression analysis across various stages of callus tissue development [80]. The raw data can be accessed through the NCBI Sequence Read Archive (SRA) under accession number PRJNA761197. The original transcriptome sequencing was also conducted with three biological replicates per sample. Based on the canonical model where promoter hypermethylation often leads to transcriptional repression [81], whereas hypomethylation is permissive for expression, we hypothesized that genes showing concordant changes (hyper/hypomethylation with corresponding low/high expression) are more likely to be directly regulated by DNA methylation. Therefore, we categorized these genes into four groups based on their levels of methylation and expression: hypermethylated coupled with low expression, hypermethylated coupled with high expression, hypomethylated coupled with low expression, and hypomethylated coupled with high expression. Subsequently, we performed GO functional and KEGG pathway enrichment analyses for these genes.

4.6. Functional Annotation Enrichment Analysis

Eucalyptus genes were annotated by mapping to the GO and KEGG databases [82]. Subsequently, we identified enriched GO terms and KEGG pathways associated with DEGs and/or genes related to DMRs using Fishers’ exact test. GO terms and KEGG pathways with false discovery rate (FDR) < 0.05 (Benjamini–Hochberg corrected p-value) were defined as significantly enriched.

4.7. qRT-PCR Analysis

Quantitative real-time PCR (qRT-PCR) was performed to validate expression changes in the genes in key pathways, following established protocols [83]. Briefly, ten target genes (AHP, BAK1, BSK, CRE1, GID1, MKS1, PR-1, PYL, RbohD, and TCH4) and the internal control H2B were analyzed using gene-specific primers designed with Primer3 (Table S10). Three biological replicates were measured for each gene, and cycle threshold (Ct) values were normalized to H2B (ΔCt). The relative expression differences in callus tissues between E. camaldulensis and E. grandis × urophylla were calculated as ΔΔCt, and the relative normalized expression was determined using the formula 2^−ΔΔCt^.