Genome-wide analysis of MYB transcription factors in Cymbidium species reveals their roles in anthocyanin synthesis

MYB transcription factors (TFs) play diverse roles in regulating plant growth and development, including secondary metabolism, signal transduction, and resistance to biotic and abiotic stresses^[1−4]. According to the repeat number, MYB proteins are generally classified into four subfamilies: 1R-MYB (MYB-related), 2R-MYB (R2R3-MYB), 3R-MYB (R1R2R3-MYB), and 4R-MYB^[1,5−7]. The first plant MYB was identified in Zea mays associated with anthocyanin biosynthesis in aleurone tissues^[8]. Since then, a large number of MYB superfamily genes has been identified in more than 70 plant genomes^[4,5,9].

MYB proteins occupy a dominant role in regulating the biosynthesis and accumulation of anthocyanin. In Arabidopsis, MYB members from the S4, S5, S6, S7, and S9 subfamilies are highly correlated with flavonoid biosynthesis^[1,10,11]. The early biosynthesis genes (EBGs), such as CHS, CHI, F3H, and F3'H, are directly bound with and regulated by the MYB TFs in the anthocyanin biosynthesis pathway (ABP)^[12]. In contrast, the late biosynthesis genes (LBGs), including DFR, ANS, and UFGT, are regulated by the MBW (MYB-bHLH-WD40) protein complex that controls the downstream accumulation of anthocyanins^[13−15]. Overexpression of the photoreactive transcription factor SmMYB35 promotes anthocyanin accumulation in the stems and petals of transgenic eggplant (Solanum melongena)^[16]. In addition to transcriptional activators, there are also R2R3-MYB transcriptional inhibitors. For example, LvMYB5 promotes anthocyanin synthesis by activating the ANS gene promoter, whereas LvMYB1 inhibits anthocyanin synthesis^[17].

The Orchidaceae family comprises over 30,000 species, characterized by their diverse floral organ morphology, and highly abundant corolla-color diversity^[18,19]. The Cymbidium species have a long cultivation history of more than 2,000 years, and are still dominating the world floriculture markets due to their excellent ornamental value^[20,21]. C. sinense and C. haematodes are terrestrial orchids with various flower colors and strong floral fragrance, making them valuable resources for Cymbidium cultivar improvement. C. haematodes was initially classified as a subspecies in C. ensifolium, and later as a variety of C. sinense (C. sinense var. haematodes). It is now widely recognized as a distinct species closely related to C. sinense^[22]. The regulation mechanisms in the anthocyanin biosynthesis pathway have been studied in several Orchidaceae species, mainly focused on tropical orchids, such as Phalaenopsis^[23], C. hybridum^[24], and Catilan^[25]. In Catilan, three R2R3-MYB transcription factors (RcPCP1, RcPAP1, and RcPAP2) in the 'KOVA' variety participate in the biosynthesis and accumulation of carotenoids and anthocyanins^[25]. However, there are few reports on the regulation of anthocyanin biosynthesis in Cymbidium, and none have been reported in C. haematodes. Identification and functional characterization of the MYB transcription factors in Cymbidium could shed light on the mechanisms underlying anthocyanin biosynthesis and flower color formation in orchids.

In the present study, a genome-wide identification of MYB TFs in C. sinense was performed, and 170 high-quality CysMYBs identified. A comprehensive analysis of the molecular characteristics, phylogenetic relationships, and gene/protein structures of the CysMYBs was conducted. Combined with flower color phenotypes, transcriptome sequencing, and qRT-PCR, functional divergence of CysMYBs were revealed in anthocyanin biosynthesis regulation between C. sinense and C. haematodes. Taken together, this study provides a valuable resource for future investigations on the role of MYB proteins in regulating anthocyanin biosynthesis in orchids.

The tepals (including sepals and petals) were collected from C. sinense variety 'Qihei', and six natural C. haematodes accessions with distinct flower colors: dark red (DR), red (R), light red (LR), yellow-green (YG), white-green (WG), and white (W), were used as the experimental materials. All plants were cultivated in the greenhouses at the Orchid Germplasm Resources Paddy of Zhejiang Province (Hangzhou, China). Healthy individuals from each C. sinense 'Qihei' and C. haematodes accessions were selected. The tepals were sampled 10 d after flowering under the same blooming condition. For each sample, three biological replicates were performed. All samples were immediately frozen in liquid nitrogen, and then stored at −80 °C for subsequent analysis.

Identification of CysMYBs

The genome assembly and annotation files of C. sinense were downloaded from NCBI (National Center for Biotechnology Information), BioSample ID: SAMN20059972. A total of 197 MYB protein sequences of Arabidopsis thaliana (AtMYB) were obtained from The Arabidopsis Information Resource (TAIR) (www.arabidopsis.org). The Hidden Markov Models (HMMs) of MYB conserved domains (PF00249, PF13921, PF14379) were retrieved from the Pfam database (http://pfam.xfam.org). Candidate MYB proteins were initially identified by performing a local BLASTP alignment against the C. sinense protein database using the 197 AtMYB protein sequences as query sequences, with an e-value cutoff of < 1e-5. Subsequently, the conserved MYB domain was used as a model to search and identify all candidate MYB proteins in the C. sinense genome, using the HMMER 3.0 software based on the HMMs. The parameters of HMMER were set as -E 1e-5 --domE 1e-5. The candidate sequences from both approaches were combined, and the domains were validated using InterProScan (www.ebi.ac.uk/interpro) and the NCBI Conserved Domain Database (CDD) (www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi). Redundant and structurally incomplete proteins, and those lacking typical domains were removed, resulting in the final set of C. sinense MYB (CysMYB) proteins. The physicochemical properties of the CysMYB proteins were analyzed using ProtParam online tool (https://web.expasy.org/protparam), and the subcellular localization of CysMYB family members was performed by Wolf PSORT (https://wolfpsort.hgc.jp) Prediction.

Phylogenetic analysis

A phylogenetic tree was constructed using the MYB proteins from five species, including A. thaliana, Oryza sativa, C. sinense, C. ensifolium, and C. goeringii. The MYB protein sequences of O. sativa were retrieved from the PFAM database (https://phytozome.jgi.doe.gov). The MYB protein sequences of C. ensifolium and C. goeringii were obtained by BLAST against the genome data^[26,27]. In addition, three separate phylogenetic trees were also constructed for R2R3-MYB, 1R-MYB, and the other MYBs from A. thaliana and C. sinense. All phylogenetic trees were constructed using MEGA X (V10.0, Tokyo Metropolitan University, Tokyo, Japan), with the Maximum Likelihood (ML) method. The bootstrap value was set to 1,000, and other parameters were retained as defaults. The C. sinense MYBs were classified into subgroups, based on the established A. thaliana classification system^[1,6] and the branching topology of the phylogenetic tree. Finally, the evolutionary trees were refined using MEGA X (V10.0), and the online tool iTOL (https://itol.embl.de).

Multiple sequence alignment and WebLogo

The R2R3-MYB protein sequences from CysMYB were subjected to multiple sequence alignment using the DNAMAN software. After that, the conserved amino acid sequences within the R2 and R3 repeats were extracted and analyzed. These trimmed and aligned sequences were then submitted to WebLogo (https://weblogo.threeplusone.com/create.cgi) to generate a sequence logo, which visually represents the consensus and conservation of the R2 and R3 domains in C. sinense MYB proteins.

Conserved motif, gene structure, and chromosomal distribution analysis

The conserved motifs of the CysMYB protein were identified using the online tool MEME (http://meme-suite.org/tools/meme), with the maximum number of motifs set to 15, and other parameters set as default. The phylogenetic tree of the CysMYB protein sequences was constructed using MEGA X software (V10.0). The genome annotation file of CysMYB genes (GFF format) was obtained from the C. sinense genome assembly^[28]. Subsequently, by integrating gene annotation, the conserved motifs and gene structures were visualized using TBtools software (v1.120)^[29].

The chromosomal location of CysMYB gene family members was extracted from the C. sinense GFF file, and visualized using TBtools software (v1.120). To analyze segmental and tandem duplication events caused by gene duplication in the C. sinense genome, CysMYB protein sequences were analyzed using MCscanX software with default parameters, and TBtools (v1.120). This analysis generated a synteny map illustrating the collinear relationships.

Flavonoid content detection

The flavonoid contents were analyzed by MetWare (www.metware.cn) using an ultra-performance liquid chromatography-tandem mass spectrometry (UPLC-MS/MS) platform. The petal samples of C. sinense 'Qihei', and five C. haematodes accessions (DR, R, LR, YG, and WG) at the full blooming stage were collected from storage at −80 °C. The samples were ground into a fine powder using a ball mill, and 50 mg of the powder was weighed and transferred into a 2 mL centrifuge tube. Subsequently, 500 μL of extraction solution (50% methanol aqueous solution containing 0.1% hydrochloric acid) was added to the tube. The mixture was vortexed for 5 min using a multi-tube vortex mixer, followed by ultrasonication for 5 min. After overnight incubation at 4 °C, the samples were centrifuged at 12,000 rpm, and 4 °C for 3 min. The supernatant was collected, and the extraction procedure was repeated once. The two supernatants were combined and filtered through a 0.22 μm syringe filter. The resulting filtrate was transferred into a 2 mL sample vial as the final anthocyanin extract for subsequent measurement and analysis. Three biological replicates, and three technical replicates were performed for each sample.

Gene expression analysis of CysMYB genes

RNA sequencing was performed on tepals of C. sinense 'Qihei', and five C. haematodes accessions (DR, R, LR, YG, and WG). Total RNA was extracted from the tepals, and subjected to RNA-seq analysis, following a protocol adapted from Wang et al.^[30]. Three biological replicates of each sample were conducted. Hisat2 (http://ccb.jhu.edu/software/hisat2) was used to compare reads with the C. sinense genome (version 4.03, http://solanaceae.plantbiology.msu.edu/pgsc_download.shtml). Transcripts were assembled and quantified using StringTie (http://ccb.jhu.edu/software/stringtie), where the reference guide model was used to estimate transcript profiling. Expression levels of the CysMYB genes were calculated as FPKM (Fragments Per Kilobase of exon per million fragments Mapped) value using Cufflinks software (http://cole-trapnell-lab.github.io/cufflinks, v2.2.1). Differential expression analysis between sample groups was performed using DESeq2 to obtain differentially expressed gene sets under different biological conditions. TBtools software (v1.120) was used to draw the heatmaps of the CysMYB genes according to the FPKM values (Supplementary Table S1).

qRT-PCR analysis

qRT-PCR analysis was performed to analyze the expression of candidate genes in the tepals of C. sinense 'Qihei' (Mol), and six C. haematodes accessions with distinct flower colors (DR, R, LR, YG, WG, and W). Total RNA was extracted using the RNA extraction kit according to Wang et al.^[30], and qualified RNA was reversely transcribed into cDNA using a cDNA synthesis system (Invitrogen, Shanghai, China). qRT-PCR primers (Supplementary Table S2) were designed by Primer Premier 5 software, and synthesized by Genepioneer Biotechnologies Co., Ltd. (Nanjing, China). The gene Mol017122, whose stable expression was experimentally validated, was used as the internal standard gene. Three biological replicates, and three technical replicates were performed for each sample. The relative expression levels (Supplementary Table S3) were calculated using the 2^−ΔΔCᴛ method, and are presented as mean ± standard deviation (SD) values of three biological replicates.

Subcellular localization

The full-length coding sequences of CysMYB131 cDNA were cloned into the pCAMBIA1300 vector. An artificial green fluorescent protein (GFP) tag was fused to the C-terminus of CysMYB131, under the control of the 35S promoter. To verify the subcellular localization, the construct was co-transformed with a nuclear marker-DsRed into epidermal cells of Nicotiana benthamiana leaves via transient transformation. The primer sequences for vector construction were designed using CE Design software (Supplementary Table S4). Fluorescence signals from the fusion protein, and the DsRed marker were observed using an Olympus FV3000 confocal laser microscope (Hachioji City, Tokyo, Japan).

Statistical analysis

Statistical analysis and plotting were conducted using Graphpad Prism software (v. 8.4.3) using one-way ANOVA (p < 0.05). * indicates significant differences in values: * p < 0.05, ** p < 0.01, *** p < 0.001, and **** p < 0.0001.

A total of 180 MYB transcription factors were identified in C. sinense, and designated CysMYB1 to CysMYB180. The CysMYBs were analyzed with InterPro and CDD database to remove redundant and structurally incomplete sequences. Finally, 170 high-quality CysMYB proteins were retained for subsequent analysis (Supplementary Table S5). In total, we 104 R2R3-MYBs, 43 1R-MYBs, four 3R-MYBs, two atypical MYBs (CysMYB5R and CysMYBCDC), and 17 MYB-CC members were identified containing the conserved LHEQLE motif. Notably, no typical 4R-MYB was identified, whereas CysMYB5R with five MYB repeats was found instead. The CysMYB proteins varied in length from 122 aa (CysMYB105) to 1106 aa (CysMYB122), with molecular weights (MW) ranging from 13,528.72 to 124,716.93 kD. The theoretical isoelectric point (pI) ranged from 4.49 (CysMYB152) to 10.32 (CysMYB175). Subcellular localization prediction indicated that 164 MYB proteins (96.47%) were located in the nucleus. In addition to nuclear proteins, the CysMYBs were also predicted with mitochondrial, cytoplasmic, or chloroplast localizations. Detailed information, including gene/protein name, protein length, number of introns/exons, MYB-domain type, chromosome localization, isoelectric point, molecular weight, and subcellular localization prediction of the 170 CysMYBs wis presented in Supplementary Table S5.

Phylogenetic tree analysis of the MYB transcription factors

To investigate the phylogenetic relationships of MYB proteins, a multiple sequence alignment of 881 protein sequences from five different species, including A. thaliana, rice, and three Cymbidium species (C. sinense, C. geogerii, and C. ensifolium) were performed. A Maximum-Likelihood (ML) phylogenetic tree was constructed, which classified the proteins into six groups (Fig. 1a). Consistent with findings in other plant species, R2R3-MYBs constituted the predominant account of all proteins. To further elucidate the evolution and potential functions of the CysMYB, three phylogenetic trees were constructed for Cys1R-MYB (Supplementary Fig. S1), CysR2R3-MYB (Fig. 1b), and the remaining CysMYB proteins (Supplementary Fig. S2), respectively. In the present study, the 42 Cys1R-MYB in C. sinense and 52 At1R-MYB in Arabidopsis are clustered into 16 subgroups (Supplementary Fig. S1). Based on the established classification of Arabidopsis R2R3-MYB^[1], the C. sinense R2R3-MYBs were divided into 33 subfamilies (Fig. 1b). These subfamilies are named continuously with the 25 AtMYBs subgroups (S1−S25), while the remaining CysR2R3-MYB subgroups designated S26−S33. Given that genes within the same phylogenetic clade or group often share similar functions, S4, S5, S6, S7 and S9 were focused on, which are known in Arabidopsis to positively or negatively regulate anthocyanin synthesis.

Figure 1. 

Phylogenetic analysis of MYB proteins. The phylogenetic trees were constructed using MEGA X with the Maximum-likelihood (ML) method and 1,000 bootstrap replicates. (a) Phylogenetic tree of MYB gene families from five plant species, including A. thaliana, O. sativa, C. sinense, C. goeringii, and C. ensifolium. The tree was divided into six major groups, highlighted in different colors. Members from R2R3-S4, S5, S6, S7, and S9 subfamilies of the five species are marked with different colors. (b) Phylogenetic tree of R2R3-MYBs proteins of C. sinense and A. thaliana. The tree was classified into 33 phylogenetic subgroups (S1–S33), each indicated by a different color.

The S4 subfamily contains seven CysMYB genes and six AtMYB genes, suggesting a potential role in the negative regulation of anthocyanin synthesis. Meanwhile, three proteins (CysMYB6, CysMYB131, and CysMYB92) were grouped into the S5 subfamily, two (CysMYB168 and CysMYB107) in the S6 subfamily, and two (CysMYB68 and CysMYB74) in the S7 subfamily. The S9 subfamily comprised five proteins (CysMYB88, CysMYB71, CysMYB72, CysMYB114, and CysMYB134), which were grouped with known AtMYBs that are involved in controlling flower color. In addition to R2R3-S4 subfamily proteins, the four 3R-CysMYB transcription factors clustered with Arabidopsis MYB-LIKE2 (MYBL2, AT1G71030)^[31,32] were also identified as potential negative regulators of anthocyanin biosynthesis. These 23 CysMYB genes were regarded as candidates with similar functions to the known AtMYBs in regulating anthocyanin biosynthesis and accumulation.

Conserved domain and structural characteristics of CysMYBs

The amino acid frequency within the repeat domains was investigated using 104 R2R3-MYB protein sequences. Multiple sequence alignment revealed that the conserved amino acids of R2 and R3 domains are [-W-(X20)-W-(X20)-W-] and [-F/I-(X19)W-(X19)-W-], respectively, where W represents tryptophane residues, F is phenylalanine, I is isoleucine, and X denotes any amino acid. Furthermore, the R3 domain contains a specific binding site, '(D/E)Lx2(R/K)x3Lx6Lx3R', which interacts with the bHLH protein (Fig. 2).

Figure 2. 

Consensus sequences, and the level of conservation of R2R3-type and R3-MYB domains from C. sinense. (a) R2 and R3 domains, based on nine R2R3-MYB alignments. The conserved tryptophan residues (Trp, W) in the MYB domain and the replaced amino acid residues in the R3 repeat were indicated with red boxes. (b) Sequence logos of the R2 and R3 MYB repeats were obtained based on full-length alignments of 104 R2R3-CysMYB. Black asterisks indicate the W residues in R2 and R3 domains, while the blue asterisk indicates the replaced amino acid residues in the R3 repeat.

Gene structure and motif analysis of CysMYBs

The diversity in gene structure and motif composition reflects functional divergence among gene family members. The exon/intron pattern of the CysMYB family genes showed that most genes within the same subfamily share similar exon/intron structures, which are closely related to subfamily classification and functional roles. Furthermore, exon/intron structure analysis revealed considerable variations in the exon numbers among CysMYB subfamilies, ranging from 2 to 14, indicating that the exons underwent loss or gain during gene evolution of the MYB family (Fig. 3; Supplementary Fig. S3).

Figure 3. 

Conserved motif and gene structure analysis of 78 selected CysMYBs, including 24 1R-MYBs (red bar), 31 R2R3-MYBs (blue bar), four 3R-MYBs (yellow bar), one 5R-MYB (cyan bar), one MYBCDC (purple bar), and 17 MYB-CC (green bar). (a) Protein names and subgroups. (b) Schematic representation of conserved motifs in CysMYB proteins. Different motifs are indicated with of different color boxes, and non-conserved sequences are represented by black lines. (c) Exon/intron structures of the CysMYB genes. Yellow boxes indicate exons, green boxes indicate UTR, and the lines between boxes indicate introns. The scale bars below (b) and (c) illustrate the length of the genes and motifs, respectively.

The most common exon numbers in CysMYBs were two or three. Specifically, 71 MYB genes were identified, with three exons and two introns, accounting for 41.76% of all MYBs. CysMYB91 and CysMYB135 possessed the highest number of exons and introns. To further investigate the sequence features, 15 conserved motifs in the CysMYB proteins (Motif1 to Motif15) were identified using MEME (Supplementary Fig. S3). Motifs 1−4 were present in the vast majority of the CysMYB proteins, and represent the characteristic DNA-binding domains of R2R3-CysMYB proteins. For example, motifs 3, 10, 12, and 13 are specific patterns to R2R3-S32, while motif 19 is specific to R2R3-S27 (Fig. 3; Supplementary Fig. S3). Meanwhile, the motif6 was unique to 1R-MYB, and motif12 was unique to R2R3-MYB. All seven CysMYBs in the S4 subfamily contain the 'EAR repression motif' (PDLNL(D/E)L) in motif 5 (PDLNLEL) (Supplementary Fig. S4). Except for motif1, 2, and 3, CysMYB88 also share motif8 with AtMYB106 and AtMYB16 (Supplementary Fig. S4). Whereas CysMYB114 and CysMYB134 do not contain motif10 or motif8, which may be due to the comparatively shorter length resulting from genome assembly errors (Supplementary Fig. S4).

Genomic distribution and syntenic analysis of the CysMYBs

Chromosomal localization analysis showed that the CysMYB genes were unevenly distributed across the 19 chromosomes of C. sinense. No MYB genes were on chromosome 11, while several chromosomes showed notable gene clustering (Fig. 4a). Chromosome 2 contained the largest number of CysMYBs (20 genes), followed by chromosomes 10 (17 genes), and 8 (15 genes). In contrast, chromosomes 15 and 20 contained the fewest, each with only two CysMYB members. The positions of several genes on the same chromosome are relatively close, such as CysMYB103 and CysMYB104 on chromosome 3, CysMYB26 and CysMYB27 on chromosome 9, and CysMYB110, CysMYB111, and CysMYB112 on chromosome 10. This indicates that these genes may be homologous genes arising from segmental or tandem duplication. Intra-genomic collinearity analysis further identified 23 syntenic gene pairs within the CysMYB family genes, such as CysMYB55/CysMYB73, CysMYB167/CysMYB124, and CysMYB118/CysMYB65 (Fig. 4b). These results indicate that the gene duplication has contributed to the expansion of the CysMYB family, potentially accompanied by functional divergence among its members.

Figure 4. 

Chromosomal distribution and gene duplications of the C. sinense MYBs genes. (a) Chromosomal locations of the CysMYBs in C. sinense that are mapped based on the of the genome annotation of C. sinense. The chromosome/contig ID is indicated on the left side of the chromosomes/contigs. The scale is indicated in megabase (Mb). (b) Collinearity analysis of the C. sinense MYBs genes. The syntenic CysMYB gene pairs are indicated by different color lines and marked on the out layer of the chromosomes.

Analysis of flavonoid components

The composition and content of flavonoids were analyzed in the tepals of C. sinense 'Qihei' (Mol), and five C. haematodes accessions with distinct flower colors including dark red (DR), red (R), light red (LR), yellow green (YG), and white green (WG) using LC-MS/MS (Supplementary Table S6). Five classes of flavonoid metabolites were detected: proanthocyanidins, anthocyanins, flavonols, dihydroflavonols, and flavanones (Supplementary Fig. S5). The flavonoid composition in C. sinense 'Qihei', and five C. haematodes accessions exhibited remarkable differences. Flavonols were the most abundant class in the five C. haematodes accessions, accounting for 86.2%−96.6% of total flavonoids (Supplementary Fig. S5). Anthocyanins were the second most abundant, with significantly higher content in DR (1,557.10 ug/g), R (2,096.26 ug/g), and LR (1,766.88 ug/g) than in the other C. haematodes accessions and C. sinense (Supplementary Table S6; Supplementary Fig. S5).

Cyanidin-based anthocyanins were identified as the major pigments responsible for tepal coloration in both C. haematodes and C. sinense (Fig. 5). Total contents of anthocyanidins in the DR, R, LR, YG, and WG accessions were 1,483.44 μg/g (96.8%), 1,985.37 μg/g (95.9%), 1,622.87 μg/g (92.3%), 297.51 μg/g (67.3%), and 438.43 μg/g (93.9%), respectively (Fig. 5b). Peonidin ranked the third most abundant anthocyanidin in the five C. haematodes accessions, whereas delphinidin content was relatively higher in C. sinense (13.91 ug/g) than in C. haematodes (1.55−5.14 ug/g). Cyanidin levels were significantly elevated in the red-colored accessions (DR, R, and LR) compared to light-colored accessions (YG and WG). These results demonstrate that cyanidin derivatives are the key anthocyanin contributing to color formation in the tepals of C. haematodes and C. sinense orchids, despite the variations in the composition of other anthocyanins.

Figure 5. 

The major anthocyanin composition and floral organ structure of C. sinense (Mol) and C. haematodes-DR, R, LR, YG, and WG. (a) Foral organ structure of the Cymbidium species for anthocyanin content detection. DR: dark-red tepals; R: red tepals; LR: light-red tepals; YG: yellow-green tepals; WG: white-green tepals. (b) Anthocyanin composition and content among different samples. The x-axis represents samples with different flower colors, C. sinense (Mol) and C. haematodes accessions with different flower colors (DR, R, LR, YG, and WG). The y-axis represents the content of pigments (ug/g FW). The pigments are indicated in different colors of dots, lines, and shades. Cyanidin: red; Delphinidin: blue; Malvidin: brown; Pelargonidin: pink; Peonidin: yellow; Petunidin: purple.

Diverse expression patterns of the CysMYB genes

To investigate the roles of CysMYB genes in anthocyanin synthesis, the expression patterns in C. sinense and the closely related species C. haematodes were analyzed. Transcriptome sequencing was performed using tepals of C. sinense 'Qihei', and five C. haematodes accessions (designated as C. haematodes-DR, R, LR, YG, WG). A total of 148 CysMYBs were expressed across all 18 samples, exhibiting species-specific expression patterns (Supplementary Fig. S6; Supplementary Table S1). Genes within the same phylogenetic subfamily generally exhibit similar expression patterns. Based on the classification of CysMYB and AtMYB, 23 R2R3-MYB genes in subfamilies S4, S5, S6, S7, and S9 may be involved in anthocyanin synthesis. The expression levels of 31 selected CysMYB genes are presented in the heatmap in Fig. 6. Among the positively regulated genes, cymsin_Mol011408 (CysMYB92) and cymsin_Mol016709 (CysMYB131) were highly expressed in R and DR tepals of C. haematodes, but showed relatively low expression in other samples. Considering that cyanidin-based derivatives were the predominant anthocyanins (Fig. 5), these genes may participate in the regulation of cyanidin biosynthesis in C. haematodes and C. sinense tepals. CysMYB79 that belongs to MYB-CC also has relative high expression levels in R and DR.

Figure 6. 

Expression patterns of 31 CysMYBs in tepals of C. sinense (Mol) and C. haematodes accessions with different flower colors (DR, R, LR, YG, and WG), according to the FPKM values. The CysMYB groups are indicated with different color blocks: 1R-MYB (black), R2R3-MYB (yellow), 3R-MYB (red), 5R-MYB (green), MYB-CC (blue), MYBCDC (purple). DR: dark-red tepals; R: red tepals; LR: light-red tepals; YG: yellow-green tepals; WG: white-green tepals.

In addition, several genes showed distinct sample-specific expression. For instance, CysMYB87 was highly expressed in R tepals. In contrast, CysMYB114 (cymsin_Mol013646) was strongly expressed in C. sinense tepals, indicating a species-specific function in anthocyanin synthesis. Similarly, CysMYB30, CysMYB72, CysMYB99, and CysMYB157 were all highly expressed in C. sinense, but lowly expressed across C. haematodes accessions, suggesting that they may be C. sinense-specific regulators. Moreover, the expression of these genes arepositively correlated with the relatively high level of pelargonidins and delphinidins specifically detected in C. sinense tepals. Six genes, including CysMYB107, CysMYB102, CysMYB39, CysMYB70, CysMYB77, and CysMYB6, were exclusively highly expressed in C. haematodes-YG, and low in other samples, suggesting a potential role in flavonoid synthesis and accumulation in YG. In conclusion, CysMYB genes exhibit functional divergence, and species-specific regulatory roles in anthocyanin synthesis between C. sinense, and the close relative C. haematodes.

qPCR analyses of CysMYB genes in different colors of tepals

To verify the function of CysMYB genes in anthocyanin biosynthesis, qPCR was performed to detect the expression levels in the tepals of C. sinense, and six C. haematodes accessions with different flower colors (DR, R, LR, YG, WG, and W). Nine candidate CysMYBs were selected and verified using qPCR according to phylogenetic analysis and gene expression abundance in RNA-seq data (Fig. 7; Supplementary Table S3). The results revealed significant differences in the expression levels of the nine genes across the tepals of different colors (Fig. 7). Among them, CysMYB72 showed a significantly high expression level in C. sinense, but low or nearly undetectable expression in C. haematodes, suggesting that CysMYB72 may be a species-specific regulator of anthocyanin synthesis in C. sinense. Similarly, CysMYB114 was also highly expressed in C. sinense. Different from the RNA-seq expression, CysMYB88 exhibits high expression levels in red spectrum samples (Mol, DR, R and LR), but low expression level in light-colored samples YG, WG, and W. These results suggest that CysMYB72 and CysMYB114 may play a function in flower color determination in C. sinense, whereas CysMYB131 and CysMYB88 are potential regulators of red coloration in C. haematodes. However, their specific roles in anthocyanin synthesis require further experimental validation.

Figure 7. 

The expression of nine CysMYBs genes among different flower colors of C. sinense and C. haematodes using qRCR. Values are mean ± standard deviation (SD) of three biological replicates. Statistics were analyzed using one-way ANOVA (p < 0.05), and * indicates significant differences in values. * p < 0.05, ** p < 0.01, *** p < 0.001, and **** p < 0.0001. Mol, C. sinense 'Qihei'; DR, R, LR, YG, WG, and W indicate C. haematodes accessions with different colors. DR: dark-red tepals; R: red tepals; LR: light-red tepals; YG: yellow-green tepals; WG: white-green tepals; W: white tepals.

Subcellular localization analysis

According to online prediction, the majority of the CysMYB proteins (96.47%) were localized to the nucleus (Supplementary Table S5). To verify this prediction, CysMYB131 was selected for experimental subcellular localization analysis. The CysMYB131-GFP protein was transiently expressed in the epidermal cells of Nicotiana benthamiana leaves, which showed that CysMYB131 was localized to the nucleus (Fig. 8). Moreover, the nucleus marker DsRed was used to co-transfect with the CysMYB131 protein. The results confirmed that CysMYB131 is a nucleus-localized transcription factor (Fig. 8).

Figure 8. 

Subcellular localization of the CysMYB131 protein. Gene-GFP fusion construct was transiently expressed in Nicotiana benthamiana. From left to right: fluorescence of the gene-GFP fusion construct (green), fluorescence of the nuclear marker DsRed (red), bright-field, and merged images.

The MYB transcription factor family is among the largest transcription factor superfamily in all eukaryotes^[1,3], which has been identified in numerous plant species. C. sinense and its close relative C. haematodes are valuable resources for orchid breeding due to their unique floral structures and rich coloration. However, research on the molecular mechanisms underlying floral pigmentation in Cymbidium species remains limited. The recent release of the whole-genome sequence of C. sinense^[28] has enabled genome-wide analysis of the MYB genes. In the present study, 170 high-quality CysMYBs were identified in the C. sinense genome, which is larger than that of Solanum lycopersicum (127 SlMYBs)^[33], but lower than that in Arabidopsis (197 AtMYBs)^[34], O. sativa (239 OsMYBs)^[35], and S. tuberosum (253 StMYBs)^[36]. R2R3-MYB proteins represent the most predominant and prevalent MYB in higher plants, with functions in responses to stresses^[37−39], specialized metabolism^[1,10,40], cell differentiation^[41], and development^[42]. One hundred and four R2R3-CysMYBs were identified, accounting for 61.18% of all CysMYBs (Fig. 1). The number was similar to that in other Orchidaceae species, such as C. goeringii^[43] (104 R2R3-MYBs), C. ensifolium^[44] (102 R2R3-MYBs), Dendrobium officinale^[45] (101 R2R3-MYBs), and P. equestris (96 R2R3-MYBs)^[46]. The predominance of R2R3-MYBs aligns with findings in other plant species, reinforcing their critical role in regulating various biological processes. Moreover, no 4R-CysMYB members were identified, while a unique 5R-MYB (CysMYB5R) was discovered in C. sinense, suggesting evolutionary divergence in the MYB family. Conserved tryptophan residues in R2 and R3 domains (Fig. 2) are critical for DNA binding and interaction with bHLH proteins, forming the MBW complex that regulates anthocyanin biosynthesis. Structural analysis of the CysMYBs, including conserved motifs and gene architectures (Fig. 3), support functional divergence and specialization among MYB subfamilies. Meanwhile, the variability in exon numbers and motif distributions among subfamilies underscores the evolutionary plasticity of MYB genes, which may drive their functional diversification.

Anthocyanins comprise various anthocyanidin pigments that impart floral organs with hues ranging from red, orange, to blue and purple^[47]. Previous studies have shown that cyanidin, pelargonidin, and peonidin are the predominant anthocyanidins in Orchidaceae plants^[21,48,49]. Cyanidins are identified as the main contributors to the coloration of perianths in certain cultivars of Pleione spp.^[49], Oncidium^[50], and Cattleya hybrid^[25]. Similarly, in C. hybrids, cyanidin and peonidin have been reported as the basis for flower color variation^[51]. In agreement with these findings, it was found that the anthocyanins in the tepals of C. sinense and C. haematodes are primarily cyanidin derivatives (Fig. 5; Supplementary Fig. S5). Interestingly, the anthocyanin content in C. sinense exhibited a significant reduction compared to C. haematodes accessions. Moreover, although DR accession exhibited a darker floral hue than R, the contents of total flavonoid, anthocyanin, and cyanidin all peaked in R, followed by LR and DR. These results suggested the participation of other pigments (i.e., flavonols and carotenoids) in determining the flower color of C. sinense and C. haematodes.

Studies in many plant species have established that the S4, S6, and S7 subgroups of the R2R3-MYB subfamily play key roles in regulating anthocyanin biosynthesis. The mechanisms underlying anthocyanin biosynthesis pathway have been studied in several Orchidaceae species, mainly focused on the tropical orchids. For example, in Phalaenopsis, PeMYB2, PeMYB11, and PeMYB12 regulate the formation of red background coloration, spot patterns, and red veins in sepals and petals^[23]. The coloration of the base color and spots in the lip flap is determined by PeMYB12 and PeMYB11, respectively^[23]. In the present study, among the 104 R2R3-CysMYBs identified in C. sinense, seven R2R3-CysMYBs were from the S4 subgroup, two from the S6 subgroup, and two from the S7 subgroup. Previous studies in related orchids have implicated S6 subfamily members in flower color variation. In C. ensifolium, CeMYB52 and CeMYB104 (S6 subfamily) are proposed as the key genes controlling flower color variation^[44]. Similarly, the CgMYB91 of the S6 subfamily was also associated with anthocyanin formation in C. goeringii, as well as CgMYB32 of S4 subfamily^[43]. However, the present expression profile revealed that most CysMYBs associated with anthocyanin biosynthesis from subfamilies S4, S6, and S7 exhibited relatively low transcript abundance across red-flowered accessions of both C. haematodes and C. sinense. In contrast, two S5 subfamily members, CysMYB131 and CysMYB92, were identified as positive regulators of anthocyanin accumulation in red-flowered C. haematodes accessions (DR and R). Their expression patterns closely mirrored the variation in cyanidin content across different color accessions of C. haematodes.

Phylogenetic analysis classified CysMYB131 as a TT2-type regulatory gene. Previous studies reported that AtMYB123 (TT2), a member of the S5 subfamily in A. thaliana, promotes proanthocyanidin accumulation in seed coats^[52]. However, heterologous overexpression of the TT2 ortholog FcMYB123 from Ficus carica significantly increased anthocyanin levels in apple fruits and calli^[53]. Similarly, the peach TT2-like gene influences floral coloration by promoting red speckle formation in white petal regions^[54]. Collectively, these findings suggest that CysMYB131 may also function in anthocyanin biosynthesis, and accumulation in C. haematodes. Given that cyanidin-based derivatives are the major anthocyanins in C. haematodes petals (Fig. 5), CysMYB131 is proposed as the potential gene in regulating cyanidin biosynthesis in C. haematodes. Conversely, the species-specific expression of CysMYB72 and CysMYB114 (both from subfamily S9) in C. sinense suggests divergent regulatory mechanisms between closely related species. The relative high levels of pelargonidin and delphinidin in C. sinense tepals imply that these species-specific expressed genes may influence flower coloration by regulating the biosynthesis of pelargonidin and delphinidin. Recent studies have revealed that AtMYB106 modulates floral coloration by suppressing petal epidermal cell morphogenesis^[55], and also frequently forms the MBW complex to promote anthocyanin biosynthesis in plants^[56]. Similarly, CysMYB88 that also from the S9 subfamily is hypothesized to be a potential MYB transcription factor regulating anthocyanin synthesis in C. sinense and C. haematodes (Fig. 7).

Further research is needed to investigate the regulatory effects of these candidate MYB transcription factors on the expression of upstream structural genes, and to validate the functions of the MYB TFs along with those of their target genes. These findings highlight the functional diversification of MYB TFs in Orchidaceae, contributing to the vast floral color diversity observed in this family.

This study presents a comprehensive genome-wide identification and functional characterization of MYB transcription factors in C. sinense, and its close relative C. haematodes, revealing their critical roles in anthocyanin biosynthesis and floral coloration. Transcriptome and qRT-PCR analyses highlighted species-specific expression patterns, such as the high expression of CysMYB72 and CysMYB114 in C. sinense. On the other hand, CysMYB131 and CysMYB88 were strongly associated with cyanidin accumulation in C. haematodes. These findings underscore the evolutionary diversification of MYB TFs in orchids and their pivotal role in regulating the biosynthesis of distinct anthocyanins. The findings provide valuable genetic resources for future studies on floral pigmentation. Future research should focus on functional validation of these candidates and their interactions within the anthocyanin regulatory network.