Functional characterisation of Dof gene family and expression analysis under abiotic stresses and melatonin-mediated tolerance in pitaya (Selenicereus undatus)

Characterisation of the Dof gene family members in pitaya

To identify the Dof gene family members, we used the pitaya genome database (http://www.pitayagenomic.com, accessed on 13 March 2023) (Chen et al. 2022a). Dof family members of maize, rice, and Arabidopsis were retrieved from Phytozome (https://phytozome-next.jgi.doe.gov, accessed on 17 October 2023) (Goodstein et al. 2012) by blasting Dof sequence. To verify the domains of Dof proteins, we used the InterPro tool (https://www.ebi.ac.uk/interpro, accessed on 13 March 2023). Additionally, for the physiological and chemical characteristics of Dof proteins, we calculated their molecular weight, isoelectric point (pI), and GRAVY score (grand average of hydropathicity) using the Expasy ProtParam database (https://web.expasy.org/protparam, accessed on 13 March 2023). Tomato (Solanum lycopersicum), maize (Zea mays), and Arabidopsis thaliana Dof gene family members were also considered in our study.

Identification of domains and phylogenetic analysis of the HuDof genes

Domains of HuDof protein sequences were identified by using the NCBI conserved domain tool (https://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi, accessed on 13 March 2023). We obtained the protein sequences for Dof protein family members of maize (Z. mays – Refgen V4), tomatoes (S. lycopersicum – ITAG4.0), and Arabidopsis (A. thaliana – TAIR-7) from Phytozome-13 (https://phytozome-next.jgi.doe.gov, accessed on 17 March 2023) (Goodstein et al. 2012). Comparative phylogenetic analysis of these plants were carried out by using MEGA X, through maximum likelihood (ML) (Kumar et al. 2018) and the tree was designed through the iTOL database (Letunic and Bork 2019). The boundary was designed through Adobe Illustrator.

Characterisation of gene structure, conserved motif, and sequence logos of HuDof genes

The gene structures of HuDof genes were determined by downloading the pitaya genome from (http://www.pitayagenomic.com, accessed on 25 March 2023). Motif analysis of 26 HuDof genes was carried out by using MEME suite database (https://meme-suite.org/meme/tools/meme, accessed on 13 March 2023), and the motifs’ structures were drawn through TBTools software (Chen et al. 2020). Sequence logos of Dof motifs were created by the WebLogo tool (Crooks et al. 2004).

Chromosomal mapping, duplication, and synteny analysis of HuDof genes on the pitaya chromosomes

Mapping of all 26 HuDof genes on chromosomes was caried out by downloading the pitaya genome (http://www.pitayagenomic.com, accessed on 13 March 2023) (Chen et al. 2022a). The length of the chromosomes was identified through TBTools by using the ‘FASTA stats’ feature (Chen et al. 2020). Its ideogram was generated through PhenoGram (http://visualization.ritchielab.org/phenograms/plot, accessed on 19 March 2023). The output image from this ideogram was used to pinpoint the mapping of all 26 HuDof genes on the chromosomes. Duplication analysis of HuDof genes was carried out through the multiple collinearity scan toolbox (MCScanX) (Chen et al. 2020), with an e-value threshold of 1 × 10⁻⁵. Sequences were subjected to one-step MCScanX in the TBTools software (Chen et al. 2020) to create the collinearity data, in order to determine the syntenic connection between the Dof genes of Hylocereus undatus, O. sativa, S. lycopersicum, and A. thaliana. The identified HuDof gene pairs were visualised through CIRCUS analysis.

Sub-cellular localisation, selection pressure analysis, and orthologous gene clustering of HuDof genes

The localisation of HuDof genes was predicted by using Wolf PSORT (https://www.genscript.com/wolf-psort.html, accessed on 25 March 2023) and CELLO (ver. 2.5, http://cello.life.nctu.edu.tw, accessed on 25 March 2023). For selection pressure analysis, stop codons (UAA, UAG, UGA) were removed from HuDof protein sequences and then subjected to analysis in Selecton Server (http://selecton.tau.ac.il/index.html, accessed on 26 March 2023). Using OrthoVenn2, (https://orthovenn2.bioinfotoolkits.net/home, accessed on 28 March 2023) the orthologous genes in H. undatus, Z. mays, S. lycopersicum, and A. thaliana have been detected.

Correlation networks of HuDof genes

We have identified the interactions of all 26 HuDof genes by obtaining data from the pitaya database (http://www.pitayagenomic.com/coexpression.php, accessed on 27 March 2023) (Chen et al. 2022a). Interacting genes were visualised through the Cytoscape program (Shannon et al. 2003).

Identification of promoter and cis-acting elements in HuDof promoter regions

We extracted the data of promoter regions from the pitaya genomic file (2000 bp upstream of the start codon) through TBTool (Chen et al. 2020). For analysis of cis-elements in promoter regions, we used the PlantCARE database (http://bioinformatics.psb.ugent.be/webtools/plantcare/html, accessed on 28 March 2023).

HuDof gene expression profiling by RNA-seq

The expression data of the HuDof genes were retrieved from pitaya database (http://www.pitayagenomic.com, accessed on 26 March 2023) and the circular heat map was constructed through TBTools (Chen et al. 2020). All 26 HuDof genes’ expression levels were evaluated in various plant tissues representing five flower phases (F1, F2, F3, F4, F5), four flowering bud stages (FB1, FB2, FB3, FB4), three pericarp phases (PeriC-45d, PeriC-65d, and PeriC-85d) and the main phases of fruit pulp (Pulp-29d, Pulp-35d, Pulp-49d).

Plant materials and stress treatments

The seeds of pitaya (Selenicereus undatus L.) were collected from Hainan-Shengda Modern Agriculture Development Co., Qionghai, Hainan, China, and grown in a controlled environment at the Nanfan Research Institute in Yazhou Bay Seed Science and Technology City, Hainan University, China. The seeds were sown in trays filled with peat-moss soil and kept in a growth room with a 16-h day (24°C) and an 8-h night (70% relative humidity). After 40 days, seedlings were chosen for heavy metal treatments, coinciding with the onset of stem differentiation into cladodes. Separate experiments were conducted to get the optimised values for each treatment including vanadium (35 mg L⁻¹), 30% 5-polyethylene glycol (PEG-6000), melatonin (100 μM), and 42°C heat by following the protocol (Zaman et al. 2023b).

RNA isolation and RT-qPCR evaluation

Total RNA was extracted from diverse plant tissues, such as stems at different stages (1-month, 1-year, and 2-year), flower buds, pericarp, and pulp by using the RNAprep Pure Plant Kit (TIANGEN, Beijing, China). The NanoDrop 2000C spectrophotometer measured RNA concentrations, and DNase I was applied to eliminate genomic DNA. Subsequently, 1 g of RNA from each sample was reverse transcribed into complementary DNA (cDNA) using the QuantiTect Reverse Transcription Kit. Quantitative reverse transcription PCR (RT-qPCR) was employed to analyse gene expression, with actin as the internal control. Three biological and three technical replications were conducted for each sample. The double-delta CT method was utilised for data analysis, and the SYBER Green Master Mix, along with the LightCycler 480 real-time PCR system, were used in the experiment (Zaman et al. 2023a).

Dof gene family expression studies in response to exogenous melatonin and several stresses

We utilised the FPKM values (fragments per kilobase of transcript per million mapped reads) for the 26 HuDof genes to analyse their expression profiles under both single and multifactorial stress combinations. Plant samples from eight treatments and two replications of single and multifactorial stresses, as described in a previous study (Zaman et al. 2023c), were collected, quickly frozen in liquid nitrogen, and stored at −80°C. The cDNA synthesis process adhered to established protocols from an earlier study (Zaman et al. 2023c). Using RNA-seq data, eight genes were chosen to assess their expression patterns in pitaya tissues across all treatments. To confirm Dof gene expression, RT-qPCR was conducted with primers listed in Supplementary file 1, following methods outlined in a previous study (Zaman et al. 2023c).

Genome-wide identification and analysis of Dof genes in pitaya genome

HuDof-1 to HuDof-26 are the 26 gene family members that were identified by study in the pitaya genome. Physio-chemical characteristics of 26 gene family members were measured, including number of amino acids, molecular weight (kDa), isoelectric point (pI), grand average of hydropathicity (GRAVY), each gene’s chromosomal number, position on the chromosome, coding sequence (CDS) length, protein length, and protein molecular weight. Protein length ranged from 48 amino acids (HU01G02303.1) to 366 amino acids (HU03G01783.1), and their molecular weights ranged from 5 to 40 kDa. Due to the fact that all 26 pitaya HuDof proteins had typical GRAVY values that ranged from −1.586 (HU07G00122.1) to −0.536 (HU02G01516.1), hydrophilicity was predicted for all 26 of them (Table 1). pI net electrical charge at a specific pH value ranged from 9.7 (HU03G01783.1) to 10.5 (HU10G00371.1), with an average value of 9.95 for all HuDof proteins (Supplementary file S1).

Identification of domains

The analysis of domain for all 26 HuDof proteins was retrieved by using the NCBI domain database. The acquired information and TBTools were used to illustrate the domain structure. This domain-based investigation found that the Dof domains were included in all 26 proteins selected for sequencing (Fig. 1).

Fig. 1.

All 26 HuDof sequences comprised of Dof protein domains. Yellow domain shows Dof super-family and black line shows no conserved sequences.

Motif analysis, gene structure, and phylogenetic tree of HuDof genes

Each Dof protein was examined for homologous motifs using the MEME program. There were ten recognised conserved motifs. Highly conserved domains were present in the C-terminus regions of the HuDof proteins. The majority of the Dof proteins had motifs that showed similarities within the group. Our understanding of the intimate evolutionary ties among the same protein group was also aided by these indications. We clustered them into seven subgroups with appropriate motif sequences. Group A contained four genes, Group B and F contained two genes, Group C, D, and G contained one gene, and Group E contained 15 genes and share similar motif pattern witin the group of HuDof family (Fig. 2).

Fig. 2.

HuDof-1 to HuDof-26 were analysed for homologous motifs using the MEME program. Motif 1 was found in almost all the HuDof genes. Motif 3, 10, 6, 4, and 9 were present in two different types of genes, and Motif 5 was present in five genes.

We obtained the Dof proteins from the genomes of each species using the Phytozome database, along with their respective Dof protein sequences from pitaya, maize, tomato, and Arabidopsis (Supplementary file S2). A total of 115 protein sequences were retrieved, including 26 HuDof from pitaya, 37 OsDof from rice, 26 ZmDof from Z. mays, and 26 AtDof from A. thaliana. These were used to construct a phylogenetic tree. The Dof genes were clustered into seven subclasses in the phylogenetic tree on the basis of conserved motifs and their phylogenetic tree (Fig. 3: subfamily A, orange; B, pink; C purple; D, light pink; E, lilac; F, blue; and G, green). Among all, Group E contained the 15 most highly conserved genes and shared similar motif pattern within the gene family, Group C contained 4 genes, Group F contained 1 gene, Group G contained 8 genes.

Fig. 3.

Comparative phylogenetic analysis of HuDof family proteins of four species, with a bootstrap value of 1000 repetitions. Stars show Zea mays, triangles shows Solanum lycopersicum, circles show Arabidopsis thaliana and squares Selenicereus undatus. Data alignment was carried out using MEGA X software, and the visualisation was crafted with iTool software and finalised with Adobe Illustrator.

It was possible to identify the gene structures by aligning the intron and exon structures of each HuDof gene with a rectangular phylogenetic tree (Fig. 4). Exons and introns can be found in the majority of HuDof gene sequences.

Fig. 4.

Phylogenetic tree and gene structure of HuDof sequences ranging from HuDof-1 to HuDof-26; the blue colour shows CDS, red is UTR and the black line is exon.

Distribution of HuDof genes amongst pitaya chromosomes, duplication, and synteny analysis

All Dof genes were mapped on specific pitaya chromosomes and illustrated through PhenoGram plot software (Fig. 5). The HuDof genes were irregularly distributed across 11 chromosomes. Chromosome 5 (Chr 5) contained five HuDof genes. There are three HuDof genes in Chr 1, 3, 4, and 8. Two genes on Chr 2, 6, 7, and 10. Chr 11 contains only one gene. Their classification was based on chromosome length and common ancestors in the rectangular phylogenetic tree. The pitaya genome’s chromosome length and each gene’s location are presented in Supplementary file S3.

Fig. 5.

Dof genes were distributed throughout 11 pitaya chromosomes. The gene names are written in black. Dof genes were divided into six groups on the basis of their domains.

The 26 HuDof genes exhibited an uneven distribution across pitaya’s 11 chromosomes, with five genes on Chr 5, two gene family members on Chr 2, 6, 7, and 10, three on Chr 1, 3, 4, and 8, and one gene on Chr 11. The pitaya genome in the TBTool was examined using the MCScanX software to see if gene duplication events occurred.

To understand divergence and duplication events among the HuDof gene family members, K_a/K_s ratios for gene pairs were calculated. Divergence time of the gene pairs was computed using the synonymous mutation rate of substitution per synonymous site per year (T = K_s/2x), which was set at x = 6.5 × 10⁻⁹ (Zaman et al. 2023b). This helps reveal the gene duplication timing in the HuDof family. Our study identified three pairs of duplication events in HuDofs, which were diverged by segmental duplications in this gene family (Table 2). A circos plot visually displays chromosome connections among HuDof proteins (Fig. 6).

Fig. 6.

Pitaya has duplicated pairs of the HuDof genes, as seen in this circos plot. The genes that have undergone gene duplication events are highlighted in various colours. All pairs show segmental duplications including HuDof23–HuDof26, HuDof3–HuDof7, and HuDof6–HuDof20. However, the rest of the pairs shows duplications, but time of divergence was not predictable. Purple digits indicate Chr numbers.

The nature of selection pressure is identified by using the K_a/K_s ratio. A K_a/K_s value below 1 indicates purifying selection, equal to 1 shows neutral selection and exceeding 1 shows positive selection. The amino acids listed in the positive and purifying selection outcomes are presented using a colour scale (Fig. 7), illustrating that purifying selection has been applied to all proteins.

Fig. 7.

Purifying selection of HuDof genes. Dark red highlights indicate purifying or negative selection on codons, grey and white highlights indicate neutral selection. Yellow and orange highlights indicate positive selection.

Synteny in S. undatus was assessed in relation to S. lycopersicum, O. sativa, and A. thaliana. Several conserved regions/genes were revealed on a number of chromosomes of the pitaya plant (Fig. 8). The evolutionary relationship of pitaya genes HuDof-6, HuDof-5, HuDof-14, HuDof-12, HuDof-17, HuDof-11, HuDof-16, HuDof-20, HuDof-26, and HuDof-25 was witnessed in relation to Arabidopsis genes AtDof-4, AtDof-10, AtDof-28, AtDof-30, and AtDof-31, respectively.

Fig. 8.

Visualisation of the syntenic analysis. Synteny analysis of S. undatus with (a) S. lycopersicum (b) O. sativa, and (c) A. thaliana were identified by applying TBtools program to do an MCScanX examination. Grey lines that link the genomes indicate synteny blocks, while red lines indicate the synteny between the genes inside these genomes specifically.

Sub-cellular localisation, selection pressure, and orthologous gene clustering of HuDof genes

The OrthoVenn2 database was used to study orthologous relationships among Dof genes in four plant species: A. thaliana, Z. mays, S. undatus, and S. lycopersicum. Seventy-three Dof proteins were clustered into seven orthologous groups (Fig. 9).

Fig. 9.

Orthologous Dof proteins clustering analysis. Clusters of the orthologous protein between the Dof gene family of Z. mays, S. undatus, S. lycopersicum, and A. thaliana have been discovered and visualised by utilising the OrthoVenn2 online tool. For analysis, 1e⁻¹⁰ was used as the e-value cut-off.

We predicted the role of HuDof genes in several cell organelles, such as the mitochondria, chloroplasts, peroxisomes, and nuclei, using gene localisation (Fig. 10). For example, HuDof-22 was highly localised in chloroplasts, with intermediate levels in the nucleus and lower in the mitochondria. HuDof-4 was also highly localised in chloroplasts, and observed at intermediate levels in nucleus, plastid and mitochondria. Thus, most of the genes were engaged in the nucleus and chloroplasts.

Fig. 10.

A heat map visually representing HuDof genes’ sub-cellular localisation using a colour scheme (red for high, yellow for intermediate, blue for low protein numbers) in various cellular compartments; i.e. vacuole (Vac), cytoplasm (Cyto), mitochondria (Mito), nucleus (Nucl), plastid (Plas) and peroxisome (Pero). Extra indicates that these genes are localised other than the given organelles in the figure.

HuDof protein network analysis

The correlation network of 26 HuDof genes revealed strong interactions among them. HuDof-2, HuDof-5, HuDof-11, HuDof-14, HuDof-18, HuDof-21, and HuDof-22 were more interactive, as shown in Fig. 11 in red.

Fig. 11.

Protein–protein interactions of HuDof family members in the pitaya. Red shows highest levels of interaction, identified as hub genes. Blue shows intermediate interaction levels and yellow shows weakest interaction. Blue lines show interaction. Green arrows represent the direction of protein interaction.

Cis-acting elements in the promoter region of HuDof genes

We examined the sequences of the 26 HuDof genes in pitaya (Figs 12 and 13), paying particular attention to the 2000 base pairs upstream of the start codon, to ascertain their biological functions. To locate cis-elements for this investigation, the PlantCARE web tool was used. Nine genes from the HuDof gene family showed elevated expression and a variety of cis-acting regulatory mechanisms (Supplementary file S4). These cis-acting elements fell into eight distinct categories (Fig. 14): AURE (auxin responsive elements), LIRE (light responsive elements), AARE (abscisic acid responsive elements), MejARE (methyl jasmonate responsive elements), DSRE (defence and stress response), GIRE (gibberellin responsive elements), LTRE (low-temperature elements), and SARE (salicylic acid responsive elements).

Fig. 12.

Sequence logos of HuDof genes showing that most of the amino acids are conserved.

Fig. 13.

Twenty-six HuDof genes have cis-acting elements in their promoter regions, which are situated 2000 bp before the start codon. This graphic was built using TBTool.

Fig. 14.

A pie chart showing the different types of elements found in the promoter regions of HuDof genes. These elements are related to various functions, such as responding to light, stress, and hormones like auxin and abscisic acid.

The seven cis-elements linked to the Dof function were examined using the PLACE database. These seven cis-elements were revealed in nearly all the HuDof: the CCAATBOX1, LTRE1HVBLT49, LTRECOREATCOR15, GT1GMSCAM4, MYBCOR, MYCCONSENSUSAT, and DOFCOREZM. DOFCOREZM has also been found in the promoters of genes that are expressed in other parts of the maize plant, such as the roots and leaves, and is involved in various types of functions including development, metabolism, and stress response. MYCCONSENSUSAT plays a role in regulating flowering in soybeans. It is also involved in the promoters of genes that are expressed in other parts of the maize plant, such as the roots and leaves (Supplementary file S5).

Expression analysis of HuDof genes

To explore the patterns of HuDof genes, RNA-seq data from 15 tissues was obtained from the pitaya genome database through TBTool. We examined the HuDof gene expression in 15 distinct tissues of the pitaya plant, encompassing pericarp (three stages: 45 days, 65 days, and 85 days), flowers (five stages: F1–F5), flower buds (four stages: FB1–FB4), and fruit pulp (three stages: 29 days, 35 days, and 49 days). To validate the projected gene expressions using the transcriptome data, we conducted RT-qPCR for HuDof-1, HuDof-2, HuDof-3, HuDof-5, HuDof-6, HuDof-8, HuDof-9, HuDof-14, HuDof-23, and HuDof-24 on the basis of high expression in the pitaya plant. Then we designed primers for eight candidate genes. HuDof-1, HuDof-2, and HuDof-8 all had greater levels of expression across the tissues. HuDof-1 exhibited high expression levels across all tissues, including stems, flower buds, and pericarp. However, fruit pulp displayed comparatively lower expression. In contrast, HuDof-2 demonstrated high expression during flowering and pericarp development, with intermediate levels in flower buds. HuDof-8 showed active expression in pulp and flower buds, with weaker expression in flowers.

The circular heat map (Fig. 15) shows that how HuDof family genes react to various stress combinations, like drought, heat, vanadium, and others. In the colour bar, red shows highly active genes, with blue indicating less activity. Some genes showed differential expression when the plants were treated with a single stress compared to double abiotic stresses, and even further under multifactorial stresses. Applying melatonin externally boosted the expression of these genes, helping to increase the impact of combined abiotic stresses, which included vanadium, heat, and drought stresses. From the results, HuDof genes play an important role in reducing both single and multiple stresses. Three genes, HuDof-1, HuDof-2, and HuDof-8, showed higher expression in response to specific stresses (D, H, and V) compared to the control (Figs 16 and 17).

Fig. 15.

The heat map showing the expression of 26 HuDof genes across various pitaya tissues and developmental stages. The 45 columns represent individual genes, and the 15 rows correspond to different tissues and stages, including flower bud, flower, pericarp, and pulp phases. Colour variations within the heat map signify gene expression levels, transitioning from light blue (low or no expression) to dark blue, and from light yellow to dark red (indicating increased expression).

Fig. 16.

An examination of the HuDof gene family involved a transcriptome analysis under single and multifactorial stress combinations. The expression levels of each gene were depicted using a colour spectrum from dark blue to dark red, measured in log10-based FPKM. The HuDof genes’ expression profiles were assessed under specific stress conditions: CK (control), D (drought), H (heat), V (vanadium), VD (vanadium + drought), VH (vanadium + heat), VHD (vanadium + heat + drought), and MVHD (melatonin + vanadium + heat + drought).

Fig. 17.

Relative expression analysis of Dof genes was conducted through RT-qPCR. The study involved the validation of the expression of eight HuDof genes in pitaya seedlings subjected to both single and multiple stress conditions. Error bars were included to indicate the standard deviation across replicates. The experimental treatments encompassed eight variations: control (CK), vanadium (V), drought (D), heat (H), vanadium + heat (VH), vanadium + drought (VD), vanadium + heat + drought (VHD), and melatonin + vanadium + heat + drought (MVHD).