Genome-wide analysis of the R2R3-MYB transcription factor genes in Chinese cabbage (Brassica rapa ssp. pekinensis) reveals their stress and hormone responsive patterns
BMC Genomics volume 16, Article number: 17 (2015)
The MYB superfamily is one of the most abundant transcription factor (TF) families in plants. MYB proteins include highly conserved N-terminal MYB repeats (1R, R2R3, 3R, and atypical) and various C-terminal sequences that confer extensive functions. However, the functions of most MYB genes are unknown, and have been little studied in Chinese cabbage.
Here, we analyzed 256 (55.2% of total MYBs) R2R3-MYB genes from Chinese cabbage (Brassica rapa ssp. pekinensis) and anchored them onto the 10 chromosomes and three subgenomes. The R2R3-, 3R- and atypical MYB proteins in Chinese cabbage formed 45 subgroups based on domain similarity and phylogenetic topology. Organization and syntenic analysis revealed the genomic distribution and collinear relationships of the R2R3-BrMYBs. Synonymous nucleotide substitution (Ka/Ks) analysis showed that the Chinese cabbage MYB DNA-binding domain is under strong purifying selection. Moreover, RNA-seq data revealed tissue-specific and distinct R2R3-BrMYB expression profiles, and quantitative real-time PCR (qPCR) analysis in leaves showed stress responsive expression and crosstalk with ABA-auxin signaling cascades.
In this study, we identified the largest MYB gene family in plants to date. Our results indicate that members of this superfamily may be involved in plant development, stress responses and leaf senescence, highlighting their functional diversity.
Plant growth and development are regulated by the coordinated expression of thousands of genes at every moment throughout their lives. Transcription factors (TFs) play a key role in these processes by self-regulating or regulating the transcription of downstream target genes. They usually consist of at least four discrete domains, namely a DNA-binding domain (DBD), a nuclear localization signal, a transcription-activation domain, and an oligomerization site . These domains function together to mediate many physiological and biochemical processes, and to activate and/or repress transcription in response to endogenous and exogenous stimuli [2,3]. Additionally, most TFs are members of gene families, thereby making their regulation more complex, but also more orderly .
The MYB superfamily is one of the largest TF families in plants . MYB proteins are found in all eukaryotes  and are defined by a highly conserved MYB DBD at the N-terminus . The MYB domain is highly conserved among eukaryotes and forms 1–4 imperfect repeats (R0, R1, R2, and R3) with a consensus sequence of approximately 50 amino acid residues. Moreover, each repeat contains regularly spread triplet tryptophan (W) residues, forming a hydrophobic core structure . The higher structure of each repeat is composed of three α-helices. The latter two helices form the HTH (helix-turn-helix) structure and bind to the promoters of target genes . The third helix plays a crucial role in DNA recognition . In general, these DBDs are localized to the N-terminus of MYBs, while their C-termini function as trans-acting domains (TAD) and vary considerably, which leads to the wide range of regulatory roles for the MYB gene family . MYB transcription factors have been separated into four classes named 1R-, R2R3-, 3R- and 4R-MYB proteins according to the number of DBD repeats .
The first identified plant MYB gene was C1, isolated from Zea mays, and encodes a c-myb-like transcription factor that regulates anthocyanin biosynthesis . An increasing number of plant R2R3-MYB superfamily members have been identified subsequently and characterized in numerous plants, such as Arabidopsis, grape, maize, petunia and snapdragon [4,12-14]. Plant R2R3-MYB proteins play important roles in many biological processes including cell metabolism [12,15], cell fate, development  and stress responses . In addition, 3R-MYBs only account for a very small proportion; for example, Arabidopsis thaliana contains only five 3R-MYB genes, compared with up to 190 R2R3-MYB and MYB-related genes .
Recently, numerous studies have shown that MYB family transcription factors play roles in plant stress responses. AtMYB15 functions as a negative regulator in the CBF pathway in response to cold stress in Arabidopsis . OsMYB2, a rice MYB gene, has been shown to respond to salt, cold, and dehydration stresses . The wheat TaMYBsdu1 gene has been reported to act as a potentially important regulator in tolerance to salt and drought stresses . AmMYB1 from Avicennia marina regulates the response processes under salt stress and transgenic tobacco plants expressing it showed better tolerance to NaCl stress . Wang et al. have reported that transferring apple MdSIMYB1 to both tobacco and apple could increase tolerance to multiple stresses .
R2R3-MYB family transcription factors participate in multiple plant-specific processes, raising the hypothesis that their expansion may be responsible for the diversity of plant evolution . R2R3-MYB families from several sequenced plants such as Arabidopsis, rice, corn, wheat, barley and soybean have been identified [4,13,23,24]. However, studies on R2R3-MYB TFs from vegetable crops have been limited and unsystematic so far. Chinese cabbage (Brassica rapa ssp. pekinensis) is a vital Cruciferae Brassica vegetable, but the functions of only a few Chinese cabbage R2R3-MYB (R2R3-BrMYB) genes (MYBs) are known . Therefore, it is very important to characterize the roles of R2R3-BrMYBs and to achieve complete identification and classification of these genes. In this study, we first identified 256 MYB family members in Chinese cabbage and then systematically analyzed their organization, collinearity and stress-responsive expression patterns. Our results showed the functional diversity of the R2R3-BrMYB genes, which may be involved in plant development, stress responses and leaf senescence.
Results and discussion
Identification and conserved DBD analysis of MYB TFs in Chinese cabbage
To define the BrMYB gene family, we searched the entire B. rapa genome sequence for genes containing the MYB domain using the Pfam program with the MYB DBD model (PF00249) as a query. We identified more than 400 sequences containing MYB or MYB-like repeats (Additional file 1: Table S1). Firstly, 21 Golgi-associated retrograde protein (GAPRs) were excluded ; consequently, based on the identification numbers and chromosome locations, any redundant sequences were removed from the dataset. To verify the reliability of our results, we also performed SMART analysis to identify all of the putative MYB protein sequences in the Chinese cabbage genome. The results were consistent with the Pfam outcome. Finally, 191 MYB-related, 256 typical R2R3-MYB (2R-MYB) (including 3 AtCDC5 homologous genes) and 11 R1R2R3-MYB (3R-MYB) proteins were successfully identified in Chinese cabbage. Six atypical MYB proteins were also identified, including four 4R-like proteins and two 5R-MYB proteins . The resulting sequences were named according to the standard constructed by Stracke , and the corresponding relationships between the names we defined and their genomic IDs are shown in Additional file 1: Table S1. Our analysis revealed that the R2R3-MYB subfamily was the largest MYB subgroup, comprising 55.73% of Chinese cabbage MYB genes (Figure 1), which was consistent with previous studies in rice and Arabidopsis [4,24].
Wang et al. divided the Chinese cabbage genome into three subgenomes according to their fractionation degree, namely the least fractionated (LF), medium fractionated (MF1), and most fractionated (MF2) subgenomes, and the LF subgenome seemed to be fractionated later than the MF1 and MF2 subgenomes because that the earlier subgenomes evolutionally appeared, the more time they would have to proceed fractionation . In our study, the LF subgenome had the highest number of MYB genes (43.97%), and atypical MYB genes were distributed in all three subgenomes (Figure 1B), indicating that atypical MYB genes appeared before the MF subgenomes began to fractionate. In total, MYB genes represented approximately 1.1% of the 41,174 predicted Chinese cabbage protein-coding loci. We also counted MYB genes in plants ranging from algae to higher plants, except P. trichocarpa , G. max , A. thaliana , V. vinifera , Z. mays  and O. sativa  that had published MYB information, while there were few genome-wide studies of MYBs in other selected plants, thus MYB numbers in these plants were obtained through the strategy used in Chinese cabbage MYB identification in this study (Figure 2); among these species, land plants seemed to carry far more MYB genes than algae, indicating that a huge expansion of MYB family members occurred after the evolution of land plants. The R2R3-MYB family is the most abundant transcription factor family in most plants, with 130 members in Arabidopsis , 141 in rice [24,31], and 118 in grape . Moreover, species-specific members of this subgroup of the MYB gene superfamily have been identified.
To investigate the homologous domain sequence features, we performed multiple alignment analysis using the 130, 256 and 141 homologous domain amino acid sequences of R2R3 repeats from Arabidopsis, Chinese cabbage and rice, respectively (Figure 3). The basic regions of the MYB domains had around 103 amino acid residues, with rare deletions or insertions as previously reported . Figure 3 shows the distribution of amino acid residues at the corresponding positions of the R2 and R3 MYB repeats of each species. Generally speaking, the distribution of conserved amino acids among the MYB domains of Chinese cabbage was very similar to those of Arabidopsis and rice, suggesting evolutionary conservation of MYBs among plants. They all included highly conserved triplet tryptophan (Trp, W) residues in each DBD repeat, and the characteristic W residues were located at positions 3, 23, and 44 of the R2 repeat (Figure 3C) and 3, 24 and 44 of the R3 repeat in Chinese cabbage (Figure 3C,D); similar localization was observed in both Arabidopsis and rice (Figure 3A,B,E and F). Conserved W residues have also been found in MYB-related and 3R-MYB genes (Additional file 2: Figure S1), indicating the indispensable role of these residues in maintaining the helix-turn-helix structure of MYB domains . In the R3 repeat, the first tryptophan (Trp3) residue was generally replaced by phenylalanine (Phe, F). However, the second and the third tryptophan residues were apparent and showed high conservation. In each repeat, the major conserved residues in the MYB domain were mainly distributed at the second and third conserved Trp residues, suggesting that the first part of each repeat in the MYB domain was apparently less conserved . This was mainly because helix-3 is highly conserved in Chinese cabbage for its DNA recognition and direct contact functions. In addition to the highly conserved W residues, more than 90% of alternative residues were highly conserved in the Chinese cabbage R2R3-MYB domains, including E-7, D-8, L-11 and G-19 in R2 repeats and G-1, E-7, G-19, G-21, N-22 and R-35 in R3 repeats (Figure 3). However, the MYB domains in both repeats of the R2R3-MYB genes in Chinese cabbage and rice seemed to be larger than that in Arabidopsis ones; this was inferred from the space between neighboring W residues. An analogous phenomenon also existed in other types of MYB domains (Additional file 2: Figure S1). The largest insertions in MYB domains were observed in rice, while the size varied only slightly between Arabidopsis and Chinese cabbage.
Chromosomal distribution and collinearity analysis of duplicated R2R3-BrMYB genes
Genome chromosomal location analysis revealed that the Chinese cabbage MYBs were distributed on all 10 chromosomes and all three subgenomes (Figure 4 and Additional file 3: Figure S2). In total, 273 BrMYBs (256 MYB-type ones and 17 members contain MYB domains > 2) were separately mapped onto chromosomes A01–A10, except for three members (BrMYB254, BrMYB255 and BrMYB256) on the scaffolds. On average, one R2R3-MYB gene was present every 2.5 Mb relative to the whole genome. Relatively high densities of BrMYBs were observed in some chromosomal regions, including the top and bottom of chromosomes A01, A02, A03, A05, A06 and A09, and the bottom of chromosomes A04, A07, A08 and A10. In contrast, almost all central chromosome regions lacked R2R3-MYBs. Among the 10 chromosomes, chromosome A03 contained the most R2R3-MYB genes, while chromosome A04 possessed the least (~5%) (Figure 1A). Furthermore, the 273 BrMYB genes were also mapped onto the chromosomes in relation to the three subgenomes (LF, MF1, and MF2), including 116 in LF, 85 in MF1, and 72 in MF2 (Figure 1B and Figure 3). Therefore, the 273 putative R2R3 proteins (including R3- and atypical MYBs), could be divided into three groups accordingly. However, 3R-MYB and atypical MYB genes were seemingly not present on all chromosomes in Chinese cabbage; furthermore, chromosomes A06 and A07 only had MYB-related and R2R3-MYB genes (Figure 4 and Additional file 3: Figure S2).
It has been confirmed that gene duplication occurred during the process of plant evolution, thereby contributing to the establishment of new gene functions . The emergence of multigene families is attributed to gene duplication via region-specific duplication or genome-wide polyploidization. MCScanX was used to further analyze the collinear relationships of the R2R3-MYBs in Chinese cabbage . It had been well addressed that a genome duplication event in Chinese cabbage occurred approximately 5–9 million years ago (MYA) and resulted in a highly duplicated genome . The collinear relationships of the duplicated pairs in the R2R3-MYB gene family in Chinese cabbage are shown in Figure 5. In total, we identified 185 pairs (pairs and groups of three or more) of highly similar paralogs that shared a high degree of identity through their protein sequences (Table 1 and Figure 5). At least eight BrMYBs were located in duplicated segments on each chromosome. Interestingly, all of the R2R3-MYB genes in subgenome LF had one or more duplicates in the other subgenomes, suggesting that all R2R3-MYB genes were retained in the genomic triplication of Chinese cabbage; this might also have contributed to the expansion of R2R3-MYB gene family.
It is likely that after duplication, a series of synonymous and/or non-synonymous mutations in their ORFs generated new functions for the BrMYBs during evolution. Therefore, we calculated the synonymous (Ks) and nonsynonymous substitutions (Ka) per site between duplicated gene pairs to estimate the selection types and divergence timing. The calculation results for the 185 duplicated pairs are listed in Table 1. All duplicated R2R3-MYB gene pairs had a Ka/Ks ratio < 1, representing purifying selection (Table 1). All of the duplicated genes were found to be segmentally duplicated according to the classify method constructed previously , which are located on duplicated segments on 10 chromosomes and 3 subgenomes in Chinese cabbage. Among them, R2R3-MYB genes containing subgenome LF have one or more duplicated genes in other subgenomes, suggesting that all the BrMYB genes have been retained in Chinese cabbage after genome triplications. In previous reports, estimations of cruciferous plant evolutionary timescales were based on the synonymous substitution rate . The divergence times of the duplicated R2R3-MYBs were also calculated as described in the “Methods”. The divergence time ranged from 0.54 (for BrMYB60-BrMYB3R10 and BrMYB87-BrMYB255) to 12.83 (for BrMYB3-BrMYB163) million years (Table 1). This indicated that the duplication events for most R2R3-MYBs in Chinese cabbage occurred after the genome triplication event (i.e., 5–9 MYA) . Their duplication seemed to be inconsistent with the whole-genome duplication, it might be caused by the genome partly deletion, and the lacking degree of the three subgenomes were distinct, which further leaded to the delay of the calculated divergence times. In contrast, the divergence times of BrMYB3R11-BrMYB3R2, BrMYB177-BrMYB231, BrMYB59-BrMYB51, BrMYB39-BrMYB233, BrMYB169-BrMYB111, BrMYB169-BrMYB104, BrMYB222-BrMYB98 and BrMYB3-BrMYB163 were earlier than the triplication event.
Phylogenetic analysis and conserved motif identification of the R2R3-MYB family in Chinese cabbage
To evaluate the evolutionary relationships within the R2R3-MYB gene family, we performed a combined phylogenetic analysis of Arabidopsis and Chinese cabbage R2R3-MYB proteins (including 7 and 17 members with more than two MYB domains, respectively) to obtain a Maximum Likelihood (ML) tree using MEGA 5 (1000 bootstrap replicates, Figure 6A and Additional file 4: Figure S3). Because of the large number of taxa and relatively low support values for informative characters, we used NJ analysis to support our subgroup designations (Additional file 5: Figure S4). The tree topologies derived from the ML and NJ analyses were basically identical, which indicated that the two methods were in strong agreement. Five sequences did not belong to any of the subfamilies (Figure 6A and Additional file 6: Table S3). The sequence similarity and phylogenetic tree topology allowed us to divide the genes into 45 subfamilies, which ranged in size from 2 to 23 MYBs (Figure 6A and Additional file 6: Table S3 and Additional file 4: Figure S3). In our subfamily classification of MYB genes, we also referred to the classification model of Arabidopsis R2R3-MYB genes constructed by Stracke et al. and Dubos et al. [4,27]. In Arabidopsis, 90 of the 126 R2R3-MYBs had been divided into 25 subfamilies (S1–25), so we labeled the previously defined clades in the trees shown in Figure 6A and Additional file 4: Figure S3 to compare our results with these studies. Most of the large subgroups (e.g., C4, C10 and C14) were supported by previous studies while some small ones (e.g., C6–9) were not. The unequal distribution of R2R3-MYBs between Chinese cabbage and Arabidopsis further supported the existence of the B. rapa whole genome triplication event (Table 1 and Figure 6A). In most subgroups defined in our ML tree, there were more R2R3-MYBs in Chinese cabbage than Arabidopsis; by contrast, the C13 subgroup included an equal number of MYBs from Arabidopsis and Chinese cabbage. These findings indicated that the R2R3-MYBs in Chinese cabbage experienced duplications after the divergence of Chinese cabbage and Arabidopsis. Notably, subgroups C5 and C22 contained 22 and 5 R2R3-BrMYBs but no R2R3-AtMYBs, which suggested that the members of these subfamilies might have specialized roles that were either lost in Arabidopsis or acquired in the Chinese cabbage lineage after divergence from the last common ancestor with Arabidopsis (Figure 6A and Additional file 6: Table S3). To determine whether this Arabidopsis ortholog gene loss phenomenon was unique to dicots or also extended to monocots, we constructed a ML phylogenetic tree of R2R3-MYBs from Arabidopsis, Chinese cabbage and rice (Additional file 7: Figure S5). The tree topology showed there were also ancestral duplication and gene loss events in rice R2R3-MYBs. Taking previous studies on poplar into consideration , we suggest that the ancestral duplication of R2R3-MYB genes might extend to various types of land plant species.
Genes from the same subfamily sharing the same motifs are likely to share similar functions . Since our classification was based on the Arabidopsis model, which was grouped according to the functions of the AtMYBs, we could further explore the common motifs and potential functions of each R2R3-BrMYB group. Using MEME, we searched for conserved motifs outside of the MYB domains of each group. Thirty-one of the 45 classified subfamilies shared one or more motifs outside of the MYB domains, which provided further support for the subfamily definitions. We identified 45 conserved motifs in the C-terminal regions and two motifs located before the MYB domains (Figure 6B and Table 2), ranging in size from 8 to 110 amino acids. This was consistent with previous hypothesis that MYBs with similar protein structures were clustered into the same conserved subfamily . Most of these conserved motifs were novel, but some had been characterized with different functions . For instance, motif 8, which was conserved in C6 members, was characterized as involved in apoptotic signaling and has SPT5 protein-binding characteristics according to SMART analysis. These results suggest that these genes probably participate in cell death processes regulated by SPT5-mediated transcriptional elongation . In many plants, R2R3-MYB (such as AtMYB95), bHLH and WD-repeat proteins regulate the anthocyanin biosynthesis pathway . The companions of AtMYB95 in C43 contained motif 46, which is known to participate in WD-repeat interactions. In Arabidopsis, AtMYBCDC5 (AtMYB125) is distantly related to typical R2R3-MYB proteins [4,27]. AtMYBCDC5 contains an R3-repeat in its MYB domain that shows low homology to typical R3-repeats in R2R3-MYBs (Figure 5B, marked with a blue box), and a very long C-terminal region. In Chinese cabbage, four genes (BrMYB246, BrMYB286, BrMYB364 and BrMYB389) were highly homologous to AtMYBCDC5 (Figure 6A, subfamily C7). A MEME search identified a conserved 80 amino acid motif in the C-terminal region of this subfamily (Figure 6B and Additional file 6: Table S3) that might interact with histone deacetylase (HDAC) proteins, raising the possibility that these subfamily members might be involved in HDAC-mediated transcription inhibition .
Expression profiling of R2R3-MYB genes in Chinese cabbage
Previously developed RNA-seq web-based tools, including tissue-specific gene expression data, allowed us to analyze the transcriptome in Chinese cabbage . Then, different transcript patterns were identified for the 273 R2R3-MYB genes (including 3R and atypical MYBs) using BRAD data, and gene expression levels were calculated in RPKM units (Additional file 8: Table S2). Consequently, we obtained expression information for each subfamily and compared the expression profiles of MYB transcription factor subfamilies of Chinese cabbage in different tissues (root, stem and leaf). We subsequently summarized these expression profiles against the phylogenic tree (Figure 6C).
As with many genes encoding transcription factors, many of the R2R3-BrMYB genes had low transcript levels according to the RNA-seq analysis. However, different transcript abundance patterns were identified in the RNA-seq dataset for the R2R3-BrMYB genes. The RPKM values of the R2R3-MYB genes are shown in Additional file 8: Table S2. Among the 273 genes, 234 were expressed in at least one tissue, while the remaining 39 members either had no expression or their expression profiles could not be found in the RNA-seq database. Nearly 120 of the 234 genes (~51%) were expressed at relatively low levels in all three tissues. For example, the expression of 29 R2R3-BrMYBs (~12%) in the roots, 66 (~28%) in the stem, and 61 (~26%) in the leaves was downregulated (Figure 6C). However, 35 (~15%) R2R3-MYBs showed high transcript levels in all three tissues, indicating that they might be indispensable in maintaining normal growth and metabolic processes of Chinese cabbage. In contrast, 101 (~43%) of the 234 genes had marked peaks in transcript levels in only one tissue, including 73 in the roots, 18 in the stem and 10 in the leaves, which suggests that these R2R3-BrMYB proteins act as regulators limited to discrete tissues or organs. For instance, 12 of the R2R3-BrMYB members with the most abundant expression in roots encode proteins in subfamily C30 (Figure 6A). The C30 (S14) subfamily containing Arabidopsis MYBs (AtMYB37, 38, 68 and 84) has been reported to function in the regulation of root development and axillary meristems [43,44]. This suggests that R2R3-BrMYB members of this subfamily that are expressed in Chinese cabbage roots may have a similar role in determining root architecture. BrMYB246 had the highest transcript abundance in leaves, and its homologous gene AtMYB16 has been shown to regulate cuticle formation in trichomes and induce over-accumulation of waxy substances on leaves ; thus, we could deduce from the peak expression of BrMYB246 in leaves that it probably be involved in waxy substance formation to protect the leaves of Chinese cabbage. In addition, some of the R2R3-BrMYBs exhibited tissue-specific expression. For example, consistent with its Arabidopsis homolog AtMYB72, BrMYB191 was only expressed in the root, indicating a role in rhizobacteria-induced systemic resistance such as AtMYB72 performs in Arabidopsis roots . However, in contrast with the expression profile of AtMYB110, which was shown to function in seed size regulation , its Chinese cabbage homolog BrMYB211 was only expressed in the stem. Overall, although the functions of most R2R3-BrMYB genes are unknown, our phylogenetic and expression profiling analyses provide a foundation for further research on R2R3-BrMYB gene functions.
R2R3-BrMYBs involved in abiotic stresses and signal transduction
R2R3-MYB proteins that have been characterized mainly participate in plant-specific processes, such as primary and secondary metabolism, cell identity, developmental regulation and stress responses [4,24]. In nature, plants suffer various biotic and abiotic stresses throughout their growth and development. Some R2R3-MYBs, such as AtMYB2, AtMYB6 and AtMYB30 are involved in responses to these stresses . We selected forty-three R2R3-BrMYBs that had relatively remarkable expression in the expression profiles above for qPCR analysis of their responses to abiotic stresses (cold and osmotic stress) and signaling hormones (ABA and auxin), to explore whether these BrMYB genes had significant performance in response to exogenous stressors. The overall expression trends of these selected genes in response to cold stress were similar under osmotic stress, and more than half of the selected R2R3-BrMYB genes were differentially expressed under at least one stressor (cold and/or osmotic stress) (Figure 7A,B). Most of the R2R3-BrMYB genes up-regulated by cold or osmotic stress reached their peak expression at about 12 h after treatment, indicating that their stress response might be rapidly regulated. By contrast, some selected genes such as BrMYB80, BrMYB170, and BrMYB250 were continuously expressed (i.e., at 0, 12, 24, 48 and 96 h) under different abiotic stresses. BrMYB210, BrMYB137, BrMYB88, BrMYB154 and BrMYB222 were significantly upregulated by both cold and osmotic stress treatments respectively (>10 fold-change), suggesting that they have roles in abiotic stress responses, much like their Arabidopsis orthologs. Similarly, the C11 subfamily members (consisting of 3R-type BrMYBs; BrMYB3R5, and BrMYB3R9) were also up-regulated by both stresses (Figure 7). Previous studies have revealed that plant MYB3R factors participate in the transcriptional control of cyclins, especially in late G2 and M phase, and OsMYB3R-2 regulates a cyclin involved in the CBF pathway to increase tolerance to low temperatures and drought [5,49]. Our results were consistent with these findings, and indicate that the 3R-MYB factors of Chinese cabbage are probably involved in stress response regulation, and that some homologous genes (e.g., BrMYB3R5-BrMYB3R9) may be functionally redundant in these processes. However, BrMYB261 (an ortholog of AtMYB28) had no response to cold treatment, but was induced drastically by osmotic stress. Since AtMYB28 was identified as a regulator of glucosinolate biosynthesis , a process involved in leaf water balance in broccoli , our findings strengthen the suggestion that MYB28 homologs may be novel regulators in the plant water deficit response.
The plant hormones ABA and auxin control important cellular processes, including seed germination, leaf senescence, stomatal aperture and stress responses . The qPCR results showed that 26 of the 43 R2R3-MYB genes were up-regulated under ABA treatment, most of which showed similar patterns to their Arabidopsis orthologs, suggesting that they function in the same processes . For instance, Arabidopsis S1 subgroup members MYB60 and MYB96 act through the ABA signaling pathway to regulate stomatal movement and disease resistance ; likewise, the C36 (S1) members BrMYB137 and BrMYB210 had relatively high transcript levels in response to exogenous ABA, suggesting their roles in the ABA signaling cascade. Moreover, most of the 43 genes responded to auxin treatment. Among them, BrMYB140, BrMYB172, BrMYB229, BrMYB208, BrMYB137 and BrMYB210 were up-regulated by ABA but down-regulated by auxin, suggesting that these R2R3-BrMYB genes might act as regulators in ABA-auxin antagonistic regulation of senescence processes . In addition, AtMYB96 was shown to regulate lateral root meristem activation via ABA-auxin signaling crosstalk . Notably, BrMYB210 (an ortholog of AtMYB96) was down-regulated by auxin treatment in leaves; thus, we hypothesize that Chinese cabbage MYB96-homologous genes may also participate in ABA-auxin signaling crosstalk in the aerial parts of Chinese cabbage.
In total, 256 (~55% of total BrMYBs) R2R3-MYB TFs were identified in the whole Chinese cabbage genome, most of which were localized on the 10 chromosomes and three subgenomes. Duplicated gene pairs among the R2R3-BrMYB genes were detected by syntenic analysis, which supports the genome triplication event in Chinese cabbage. Phylogenetic analysis of the R2R3-MYB family in Chinese cabbage and Arabidopsis revealed the conserved organization of this family, which further indicates that R2R3-MYB family members from various plants underwent gene duplication events with a common origin and were retained over a long period by each genome. Additionally, the increased number of R2R3-MYBs that seemingly evolved independently in Chinese cabbage and rice may contribute to plant viability under adverse conditions and functional specialization of R2R3-MYB genes. In addition, the tissue-specific expression profiles of the R2R3-MYB genes suggest that some of them have important roles in developmental and metabolic processes. Moreover, qPCR analysis indicated that several genes might function in stress responses and ABA-auxin hormone signal-mediated morphogenesis and cell senescence, which further highlights the functional diversity and indispensability of the R2R3-MYB genes in the normal growth and development of Chinese cabbage. This study gives an overview of the R2R3-MYB genes in Chinese cabbage and enabled us to provide some insights into plant stress response mechanisms and how transcription factors act in complex signal transduction, but how R2R3-BrMYB genes participate in these processes will require further investigation.
Identification of MYB transcription factors in different plants
The whole-genome proteins of Chinese cabbage were downloaded from BRAD (http://brassicadb.org/brad/) and those of other species used were obtained from PlantGDB (http://www.plantgdb.org). Then, the Pfam program was employed to search for candidate MYB genes in the extracted full-length protein sequences (http://pfam.sanger.ac.uk/). Only hits with e-values < 1.0 were considered to be members of the MYB family . To confirm the obtained amino acid sequences, the putative MYB sequences were examined for the MYB domain using the hidden Markov model of the SMART tool (http://smart.embl-heidelberg.de/) and the ExPASy Proteomics Server (http://expasy.org/prosite/) . Manual inspection was performed to ensure that the putative MYB genes contained conserved Trp (W) residues. The sequences of all MYB members in the genomes of other species assessed were downloaded from the plant TFDB database (http://planttfdb.cbi.edu.cn/). However, gene identifiers for 132 Arabidopsis thaliana R2R3- and R1R2R3-MYB genes were obtained from TAIR (http://www.arabidopsis.org/).
Protein properties and conserved motif analysis
To investigate the protein properties of the putative BrMYB proteins, their molecular weights (MW) and isoelectric points (pI) were calculated using Pepstats (http://www.ebi.ac.uk/Tools/seqstats/emboss_pepstats/). The conserved motifs of the R2R3-MYB proteins were identified statistically with the MEME program (version 4.8.1) (http://meme.nbcr.net/meme/intro.html) . The following parameter settings were used: maximum number of motifs, 50; minimum width of motif, 6; maximum width of motif, 250. All putative motifs with expected values < 1e-10 were discarded. Subsequently, the MAST program (version 4.8.1) (http://meme.nbcr.net/meme/cgi-bin/mast.cgi) was used to align the conserved motifs of the proteins.
Multiple sequence alignment and phylogenetic analysis
Phylogenetic trees were produced individually using the full-length sequences of the R2R3-type MYB TFs. The DNA-binding domains (DBDs) of MYB proteins from Arabidopsis, Chinese cabbage and rice were subjected to multiple alignment analysis with ClustalW (http://www.ebi.ac.uk/Tools/msa/clustalw2/) and Weblogo analysis . Phylogenetic analyses were conducted using MEGA5 (http://www.megasoftware.net/) with the Maximum-Likelihood (ML) and Neighbor-Joining (NJ) methods; the bootstrap value was set to 1000.
Identification of orthologous and paralogous MYBs
The position of each BrMYB was marked on the chromosomes using a Perl script. The orthologous and paralogous MYB genes in Chinese cabbage and Arabidopsis were identified using OrthoMCl (http://orthomcl.org/orthomcl/). The relationships between the orthologous and paralogous genes among the three species were plotted using Circos (http://circos.ca/).
Syntenic analysis and Ka/Ks calculation
The duplicated R2R3-MYB genes were identified using MCScanX (http://chibba.pgml.uga.edu/mcscan2/) as previously described . The whole-genome protein sequences from Chinese cabbage were compared against each other using BLASTP, with a tabular output format and an e-value < 1e-20. The BLASTP results with simplified gene location files were used as an input for MCScanX to identify syntenic gene pairs and duplication types with default settings. We calculated the synonymous rate (Ks), non-synonymous rate (Ka) and evolutionary constraint (Ka/Ks) between the duplicated pairs of R2R3-BrMYBs (Table 1) based on their coding sequence alignments , and the divergence time was calculated according to the neutral substitution rate of 1.5 × 10−8 substitutions per site per year for chalcone synthase .
RNA-seq data analysis
To analyze the Chinese cabbage R2R3-MYB expression patterns, we used Illumina RNA-seq data reported previously . These data included three tissues (root, stem and leaf) of B. rapa. Gene expression levels were calculated as reads per kilobase of exon model per million mapped reads (RPKM) units (Additional file 8: Table S2). Heat maps were generated and hierarchical clustering was done using Cluster 3.0.
Plant materials, growth conditions and stress treatments
Seedlings of Chinese cabbage cultivar YANZA03 were germinated in plastic Petri dishes in darkness at 22°C for 2 days, and then transferred to pots containing soil growth medium under artificial growth conditions of 22°C, approximately 120 μmol photons m−2 s−1, a photoperiod of 16/8 h, and 60% relative humidity. Half-strength Murashige and Skoog liquid solution (pH 5.8) was added once every 3 days. Five-leaf-stage plants were subjected to various treatments under a continuous time course (0, 12, 24, 48, and 96 h). For cold treatment, the pots were exposed to low temperature (4°C) conditions; for osmotic stress treatment, the pots were irrigated with 15% (w/v) polyethylene glycol (PEG) and kept standing in the irrigation solution for 30 minutes under normal growth conditions; hormone treatments were performed with ABA (100 μM) and auxin (50 mg/L NAA). The seedlings were harvested under a continuous time course (0, 12, 24, 48, and 96 h) with three biological replicates for RNA preparation.
RNA isolation and quantitative real time-PCR (qPCR) analysis
Total RNA was isolated from treated leaves using Trizol (Invitrogen, San Diego, CA, USA) according to the manufacturer’s instructions. The total RNA was treated with DNase I (Invitrogen) and 1 μg treated RNA was reverse-transcribed using PrimeScript™ RT reagent Kit (Perfect Real Time) for qPCR (Takara, Dalian China). The GAPDH gene was used as an internal control . The qPCR assays were performed with three biological and technical replicates. The SYBR® select Master Mix (Invitrogen) was used to detect gene expression according to the manufacturer’s recommendations on the One-step Real-Time PCR System (Applied Biosystems). qPCR was carried out according to a previous report . Gene-specific primers that were used to detect transcripts are listed in Additional file 9: Table S4. The PCR conditions and relative gene expression calculation were as previously described .
Availability of supporting data
The supporting sequence data are available in the Additional file 10: Table S5, and were obtained from Brassica database (http://brassicadb.org/brad/index.php). The supporting expression profile data are available in the Additional file 8: Table S2, and were obtained from a public data set (http://brassicadb.org/brad/genomeDominanceData.php).
Ptashne M. How eukaryotic transcriptional activators work. Nature. 1988;335:683–9.
Riechmann J, Heard J, Martin G, Reuber L, Keddie J, Adam L, et al. Arabidopsis transcription factors: genome-wide comparative analysis among eukaryotes. Science. 2000;290(5499):2105–10.
Amoutzias G, Veron A, Weiner J, Robinson-Rechavi M, Bornberg-Bauer E, Oliver S, et al. One billion years of bZIP transcription factor evolution: conservation and change in dimerization and DNA-binding site specificity. Mol Biol Evol. 2007;24(3):827–35.
Dubos C, Stracke R, Grotewold E, Weisshaar B, Martin C, Lepiniec L. MYB transcription factors in Arabidopsis. Trends Plant Sci. 2010;15(10):573–81.
Feller A, Machemer K, Braun EL, Grotewold E. Evolutionary and comparative analysis of MYB and bHLH plant transcription factors. Plant J. 2011;66(1):94–116.
Lipsick JS. One billion years of Myb. Oncogene. 1996;13(2):223–35.
Kanei-Ishii C, Sarai A, Sawazaki T, Nakagoshi H, He D-N, Ogata K, et al. The tryptophan cluster: a hypothetical structure of the DNA-binding domain of the myb protooncogene product. J Biol Chem. 1990;265(32):19990–5.
Hanis C, Boerwinkle E, Chakraborty R, Ellsworth D, Concannon P, Stirling B, et al. A genome–wide search for human non–insulin–dependent (type 2) diabetes genes reveals a major susceptibility locus on chromosome 2. Nat Genet. 1996;13(2):161–6.
Ogata K, Morikawa S, Nakamura H, Hojo H, Yoshimura S, Zhang R, et al. Comparison of the free and DNA-complexed forms of the DMA-binding domain from c-Myb. Nat Struct Mol Biol. 1995;2(4):309–20.
Davidson C, Ray E, Lipsick J. Evolution of Myb proteins. In: Myb transcription factors: Their role in growth, differentiation and disease. Proteins Cell Regul. 2004:1–33.
Paz-Ares J, Ghosal D, Wienand U, Peterson P, Saedler H. The regulatory c1 locus of Zea mays encodes a protein with homology to myb proto-oncogene products and with structural similarities to transcriptional activators. EMBO J. 1987;6(12):3553.
Bedon F, Grima-Pettenati J, Mackay J. Conifer R2R3-MYB transcription factors: sequence analyses and gene expression in wood-forming tissues of white spruce (Picea glauca). BMC Plant Biol. 2007;7(1):17.
Du H, Feng B-R, Yang S-S, Huang Y-B, Tang Y-X. The R2R3-MYB transcription factor gene family in maize. PLoS One. 2012;7(6):e37463.
Matus JT, Aquea F, Arce-Johnson P. Analysis of the grape MYB R2R3 subfamily reveals expanded wine quality-related clades and conserved gene structure organization across Vitis and Arabidopsis genomes. BMC Plant Biol. 2008;8(1):83.
Stracke R, Ishihara H, Huep G, Barsch A, Mehrtens F, Niehaus K, et al. Differential regulation of closely related R2R3‐MYB transcription factors controls flavonol accumulation in different parts of the Arabidopsis thaliana seedling. Plant J. 2007;50(4):660–77.
Perez-Rodriguez M, Jaffe FW, Butelli E, Glover BJ, Martin C. Development of three different cell types is associated with the activity of a specific MYB transcription factor in the ventral petal of Antirrhinum majus flowers. Development. 2005;132(2):359–70.
Deluc L, Bogs J, Walker AR, Ferrier T, Decendit A, Merillon J-M, et al. The transcription factor VvMYB5b contributes to the regulation of anthocyanin and proanthocyanidin biosynthesis in developing grape berries. Plant Physiol. 2008;147(4):2041–53.
Ding Z, Li S, An X, Liu X, Qin H, Wang D. Transgenic expression of MYB15 confers enhanced sensitivity to abscisic acid and improved drought tolerance in Arabidopsis thaliana. J Genet Genomics. 2009;36(1):17–29.
Yang A, Dai X, Zhang W-H. A R2R3-type MYB gene, OsMYB2, is involved in salt, cold, and dehydration tolerance in rice. J Exp Bot. 2012;63(7):2541–56.
Liu H, Zhou X, Dong N, Liu X, Zhang H, Zhang Z. Expression of a wheat MYB gene in transgenic tobacco enhances resistance to Ralstonia solanacearum, and to drought and salt stresses. Funct Integr Genomics. 2011;11(3):431–43.
Ganesan G, Sankararamasubramanian H, Harikrishnan M, Ashwin G, Parida A. A MYB transcription factor from the grey mangrove is induced by stress and confers NaCl tolerance in tobacco. J Exp Bot. 2012;63(12):4549–61.
Wang RK, Cao ZH, Hao YJ. Overexpression of a R2R3 MYB gene MdSIMYB1 increases tolerance to multiple stresses in transgenic tobacco and apples. Physiol Plantarum. 2014;150(1):76–87.
Wilkins O, Nahal H, Foong J, Provart NJ, Campbell MM. Expansion and diversification of the Populus R2R3-MYB family of transcription factors. Plant Physiol. 2009;149(2):981–93.
Katiyar A, Smita S, Lenka SK, Rajwanshi R, Chinnusamy V, Bansal KC. Genome-wide classification and expression analysis of MYB transcription factor families in rice and Arabidopsis. BMC Genomics. 2012;13(1):544.
Yu S, Zhang F, Yu Y, Zhang D, Zhao X, Wang W. Transcriptome profiling of dehydration stress in the Chinese cabbage (Brassica rapa L. ssp. pekinensis) by tag sequencing. Plant Mol Biol Rep. 2012;30(1):17–28.
Hosoda K, Imamura A, Katoh E, Hatta T, Tachiki M, Yamada H, et al. Molecular structure of the GARP family of plant Myb-related DNA binding motifs of the Arabidopsis response regulators. Plant Cell Online. 2002;14(9):2015–29.
Stracke R, Werber M, Weisshaar B. The R2R3-MYB gene family in Arabidopsis thaliana. Curr Opin Plant Biol. 2001;4(5):447–56.
Wang X, Wang H, Wang J, Sun R, Wu J, Liu S, et al. The genome of the mesopolyploid crop species Brassica rapa. Nat Genet. 2011;43(10):1035–9.
Du H, Yang SS, Liang Z, Feng BR, Liu L, Huang YB, et al. Genome-wide analysis of the MYB transcription factor superfamily in soybean. Bmc Plant Biology. 2012;12(1):106.
Jiang CZ, Gu X, Peterson T. Identification of conserved gene structures and carboxy-terminal motifs in the Myb gene family of Arabidopsis and Oryza sativa L. ssp indica. Genome Biology. 2004;5(7):R46.
Streisfeld MA, Young WN, Sobel JM. Divergent selection drives genetic differentiation in an R2R3-MYB transcription factor that contributes to incipient speciation in mimulus aurantiacus. PLoS Genet. 2013;9(3):e1003385.
Durbarry A, Vizir I, Twell D. Male germ line development in Arabidopsis duo pollen mutants reveal gametophytic regulators of generative cell cycle progression. Plant Physiol. 2005;137(1):297–307.
Saikumar P, Murali R, Reddy EP. Role of tryptophan repeats and flanking amino acids in Myb-DNA interactions. Proc Natl Acad Sci. 1990;87(21):8452–6.
Tombuloglu H, Kekec G, Sakcali MS, Unver T. Transcriptome-wide identification of R2R3-MYB transcription factors in barley with their boron responsive expression analysis. Mol Genet Genomics. 2013;288(3–4):141–55.
Cannon SB, Mitra A, Baumgarten A, Young ND, May G. The roles of segmental and tandem gene duplication in the evolution of large gene families in Arabidopsis thaliana. BMC Plant Biol. 2004;4(1):10.
Tang J, Wang F, Hou X-L, Wang Z, Huang Z-N. Genome-wide fractionation and identification of WRKY transcription factors in Chinese Cabbage (Brassica rapa ssp. pekinensis) reveals collinearity and their expression patterns under abiotic and biotic stresses. Plant Mol Biol Rep. 2013;32(4):1–15.
Wang YP, Tang HB, DeBarry JD, Tan X, Li JP, Wang XY, et al. MCScanX: a toolkit for detection and evolutionary analysis of gene synteny and collinearity. Nucleic Acids Res. 2012;40(7):49.
Teichmann SA, Babu MM. Gene regulatory network growth by duplication. Nat Genet. 2004;36(5):492–6.
Cheng F, Wu J, Fang L, Sun S, Liu B, Lin K, et al. Biased gene fractionation and dominant gene expression among the subgenomes of Brassica rapa. PLoS One. 2012;7(5):e36442.
He X-J, Chen T, Zhu J-K. Regulation and function of DNA methylation in plants and animals. Cell Res. 2011;21(3):442–65.
Pourcel L, Irani NG, Koo AJ, Bohorquez‐Restrepo A, Howe GA, Grotewold E. A chemical complementation approach reveals genes and interactions of flavonoids with other pathways. Plant J. 2013;74(3):383–97.
De Ruijter A, Van Gennip A, Caron H, Kemp S, van Kuilenburg A. Histone deacetylases (HDACs): characterization of the classical HDAC family. Biochem J. 2003;370:737–49.
Keller T, Abbott J, Moritz T, Doerner P. Arabidopsis REGULATOR OF AXILLARY MERISTEMS1 controls a leaf axil stem cell niche and modulates vegetative development. The Plant Cell Online. 2006;18(3):598–611.
Müller D, Schmitz G, Theres K. Blind homologous R2R3 Myb genes control the pattern of lateral meristem initiation in Arabidopsis. The Plant Cell Online. 2006;18(3):586–97.
Oshima Y, Mitsuda N. The MIXTA-like Transcription factor MYB16 is a major regulator of cuticle formation in vegetative organs. Plant Signal Behav. 2013;8(11):e26826.
Ent S, Pozo MJ, Verhagen B, Bakker D, Van Loon L, Pieterse C. Transcription factors in roots and shoots of Arabidopsis involved in rhizobacteria-induced systemic resistance. IOBC/wprs Bulletin. 2006;29(2):157–61.
Zhang Y, Liang W, Shi J, Xu J, Zhang D. MYB56 encoding a R2R3 MYB transcription factor regulates seed size in Arabidopsis thaliana. J Integr Plant Biol. 2013;55(11):1166–78.
Raffaele S, Rivas S. Regulate and be regulated: integration of defense and other signals by the AtMYB30 transcription factor. Frontiers in plant science. 2013;4:98.
Ma Q, Dai X, Xu Y, Guo J, Liu Y, Chen N, et al. Enhanced tolerance to chilling stress in OsMYB3R-2 transgenic rice is mediated by alteration in cell cycle and ectopic expression of stress genes. Plant Physiol. 2009;150(1):244–56.
Gigolashvili T, Yatusevich R, Berger B, Müller C, Flügge UI. The R2R3‐MYB transcription factor HAG1/MYB28 is a regulator of methionine‐derived glucosinolate biosynthesis in Arabidopsis thaliana. Plant J. 2007;51(2):247–61.
López-Berenguer C, Martínez-Ballesta MC, García-Viguera C, Carvajal M. Leaf water balance mediated by aquaporins under salt stress and associated glucosinolate synthesis in broccoli. Plant Sci. 2008;174(3):321–8.
Bari R, Jones JD. Role of plant hormones in plant defence responses. Plant Mol Biol. 2009;69(4):473–88.
Abe H, Urao T, Ito T, Seki M, Shinozaki K, Yamaguchi-Shinozaki K. Arabidopsis AtMYC2 (bHLH) and AtMYB2 (MYB) function as transcriptional activators in abscisic acid signaling. The Plant Cell Online. 2003;15(1):63–78.
Nooden L. Abscisic acid, auxin, and other regulations of senescence. 1988.
Seo PJ, Park C-M. Auxin homeostasis during lateral root development under drought condition. Plant Signal Behav. 2009;4(10):1002–4.
Punta M, Coggill PC, Eberhardt RY, Mistry J, Tate J, Boursnell C, et al. The Pfam protein families database. Nucleic Acids Res. 2012;40(D1):D290–301.
Letunic I, Copley RR, Schmidt S, Ciccarelli FD, Doerks T, Schultz J, et al. SMART 4.0: towards genomic data integration. Nucleic Acids Res. 2004;32 suppl 1:D142–4.
Bailey TL, Williams N, Misleh C, Li WW. MEME: discovering and analyzing DNA and protein sequence motifs. Nucleic Acids Res. 2006;34 suppl 2:W369–73.
Crooks GE, Hon G, Chandonia J-M, Brenner SE. WebLogo: a sequence logo generator. Genome Res. 2004;14(6):1188–90.
Wang Y, Tang H, DeBarry JD, Tan X, Li J, Wang X, et al. MCScanX: a toolkit for detection and evolutionary analysis of gene synteny and collinearity. Nucleic Acids Res. 2012;40(7):e49–9.
Zhang Z, Li J, Zhao X-Q, Wang J, Wong GK-S, Yu J. KaKs_Calculator: calculating Ka and Ks through model selection and model averaging. Genomics Proteomics Bioinformatics. 2006;4(4):259–63.
Koch MA, Haubold B, Mitchell-Olds T. Comparative evolutionary analysis of chalcone synthase and alcohol dehydrogenase loci in Arabidopsis, Arabis, and related genera (Brassicaceae). Mol Biol Evol. 2000;17(10):1483–98.
Qi JN, Yu SC, Zhang FL, Shen XQ, Zhao XY, Yu YJ, et al. Reference Gene Selection for Real-Time Quantitative Polymerase Chain Reaction of mRNA Transcript Levels in Chinese Cabbage (Brassica rapa L. ssp pekinensis). Plant Mol Biol Rep. 2010;28(4):597–604.
Tang J, Wang F, Wang Z, Huang Z, Xiong A, Hou X. Characterization and co-expression analysis of WRKY orthologs involved in responses to multiple abiotic stresses in Pak-choi (Brassica campestris ssp. chinensis). BMC Plant Biology. 2013;13(1):188.
Wang F, Hou X, Tang J, Wang Z, Wang S, Jiang F, et al. A novel cold-inducible gene from Pak-choi (Brassica campestris ssp. chinensis), BcWRKY46, enhances the cold, salt and dehydration stress tolerance in transgenic tobacco. Mol Biol Rep. 2012;39:4553–64.
"This study is supported by National Natural Science Foundation of China (Key Program, No.31330067, Jiangsu Provincial Science and Technology Support Program of China (key program, No. BE2013429), A Project Funded by the Priority Academic Program Development of Jiangsu Higher Education Institutions, the National Natural ScienceFoundation of China (Grant No. 31201634, 41201241) and Jiangsu Province Natural Science Foundation (Grant No. BK2012074)
The authors declare that they have no competing interests.
ZW, JT and X-LH designed research; ZW, RH and PW performed research; A-SX and X-MS contributed new reagents/analytic tools; ZW, JT and RH analyzed data; and ZW and JT wrote the paper. All authors read and approved the final manuscript.
Listing of MYB transcription factor genes in Brassica rapa.
Comparisons of DNA-binding domain of MYB-related and 3R-MYB transcription factor proteins in Arabidopsis, Chinese cabbage and rice.
Distribution of MYB genes on 10 chromosomes and 3 subgenomes.
Phylogenetic relationships and subgroup designations in MYB proteins from Chinese cabbage and Arabidopsis.
The NJ tree of R2R3-MYBs from Arabidopsis and Chinese cabbage.
The statistic analysis of each group of ML tree.
ML phylogenetic tree of R2R3-MYBs from Arabidopsis, Chinese cabbage and rice.
The genomic distribution and RPKM values of R2R3 MYB gene family in Chinese cabbage.
Primers for quantitative PCR of R2R3-BrMYBs.
Listing of BrMYB gene sequences.
About this article
Cite this article
Wang, Z., Tang, J., Hu, R. et al. Genome-wide analysis of the R2R3-MYB transcription factor genes in Chinese cabbage (Brassica rapa ssp. pekinensis) reveals their stress and hormone responsive patterns. BMC Genomics 16, 17 (2015). https://doi.org/10.1186/s12864-015-1216-y