Genome-Wide Analysis of Citrus R2R3MYB Genes and Their Spatiotemporal Expression under Stresses and Hormone Treatments

The R2R3MYB proteins represent one of the largest families of transcription factors, which play important roles in plant growth and development. Although genome-wide analysis of this family has been conducted in many species, little is known about R2R3MYB genes in citrus, In this study, 101 R2R3MYB genes has been identified in the citrus (Citrus sinesis and Citrus clementina) genomes, which are almost equal to the number of rice. Phylogenetic analysis revealed that they could be subdivided into 21 subgroups. The evolutionary relationships and the intro-exon organizations were also analyzed, revealing strong gene conservation but also the expansions of particular functional genes during the plant evolution. Tissue-specific expression profiles showed that 95 citrus R2R3MYB genes were expressed in at least one tissue and the other 6 genes showed very low expression in all tissues tested, suggesting that citrus R2R3MYB genes play important roles in the development of all citrus organs. The transcript abundance level analysis during abiotic conditions (NaCl, abscisic acid, jasmonic acid, drought and low temperature) identified a group of R2R3MYB genes that responded to one or multiple treatments, which showed a promising for improving citrus adaptation to stresses. Our results provided an essential foundation for the future selection of the citrus R2R3MYB genes for cloning and functional dissection with an aim of uncovering their roles in citrus growth and development.


Introduction
MYB gene family is large, functionally diverse and present in all eukaryotes, the proteins encoding by which usually function as transcription factors with MYB binding domain conferring the ability to bind DNA. The MYB domain is consisting of one to four imperfect repeats (R), and each repeat has about 52 amino acid residues [1]. Based on the number of repeat(s) in the MYB domain, MYB proteins are divided into four types: 4RMYB contains four repeats, 3RMYB (R1R2R3MYB) has three consecutive repeats, R2R3MYB possesses two repeats, and the MYB-related type usually, but not always, has a single repeat [2]. Among these four types, R2R3MYB is specific to higher plants and quantitatively predominant in most plants, which is characterized by the presence of a conserved MYB domain and a highly variable C-terminal region [3,4].
Based on their well conserved DNA-binding domains, genome-wide identification of R2R3MYB members has been conducted in various plants, such as Arabidopsis (126 members) [5], Oryza sativa (102 members) [5], Vitis vinifera (117 members) [6], Populus trichocarpa (192 members) [7], Zea mays (157 members) [4], Glycine max (over 200 members) [3] and Cucumis sativus (55) [2]. In Arabidopsis, the members of the R2R3MYB family were classified into 25 subgroups [1]. By comparative phylogenetic analysis, Wilkins et al. identified new R2R3MYB subgroups from Pupulus trichocarpa [7] that had no representatives in Arabidopsis, the same case to some other plant species, suggesting these proteins possess specialized biological functions that have obtained after divergence from the last common ancestor or were lost in Arabidopsis or both. The expansion of the R2R3MYB transcription factors in plants is well favor of the observation that numerous R2R3MYB proteins play central roles in plant-specific processes [7].
The functions of R2R3MYB genes have been extensively studied in various plant species, which provided us a better understanding of this gene subfamily. However, very little is known about this gene subfamily in citrus. To the best of our knowledge, all available data, thus far, about R2R3MYB genes in citrus are related to anthocyanin biosynthesis [27,28]. Citrus as one of the most important economic crops for its high healthy value is widely grown all over the world. However, its growth and production are severely affected by numerous biotic and abiotic stresses including drought, temperature extremes, salinity and pathogens. Therefore, identification and functional analysis of citrus defense-and stressrelated genes may help to elucidate the molecular mechanisms underlining the plant defense and to improve plant stress tolerance.
Recently, two citrus genome sequences including sweet orange (Citrus sinesis) and clementine (Citrus clementina) were released (http://www.phytozome.net), thereby enabling genome-wide identification and analysis of citrus R2R3MYB genes to be conducted. In the present research, R2R3MYB genes were isolated based on genomic information available at (http://www.phytozome.net). Phylogenetic and structural analysis was conducted using the citrus R2R3MYB genes according to sequence data. Organ specific transcription profiles of R2R3MYB genes were constructed for various organs from Citrus junos cv. 'Ziyang'. Furthermore, the temporal expression analysis of R2R3MYB genes in response to stresses and hormones was also performed, which helped us to identified the potentially genes that participate in the stress signal transduction pathway in citrus. Additionally, these results, for the first time, provide information upon the relationship between functional divergence and evolution in citrus R2R3MYB subfamily.

Identification of Citrus R2R3MYB Genes
To identify the citrus R2R3MYB genes from citrus (sweet orange and clementine) genome (http://www.phytozome.net), a BLASTP search has been performed at the Join Genome Institute (JGI) (http://www.phytozome.net) using the amino acid sequences of Arabidopsis R2R3MYBs. All of the sequences having MYB domain were obtained from the citrus genome database. To further confirm the reliability of our results, the functional and structural domains were predicted by PROSITE profiling (http://www.expasy.org/tools/scanprosite/) [29] and SMART analysis (http://smart.embl-heidelberg.de/) [30], respectively. Only the sequences having two repeats (R2 and R3) were selected as the candidates.

Sequence Analysis
In order to analyze the sequence features of the 101 predicted Citrus R2R3MYB proteins, multiple sequence alignment of the MYB domains was performed using ClustalX [32] with default parameters. The deduced amino acid sequences were adjusted manually using BioEdit (version 7.0.0) with default parameters (The National Resource for Biomedical Supercomputing (NRBSC): http://www.nrbsc. org/) [33] with the aim to obtain optimized alignment.
The intron pattern can serve as an independent criterion in support of subgroup designations of phylogenetic analysis. Therefore, intron pattern analysis of CitMYB genes was performed. The genomic and cDNA sequences corresponding to each predicted CitMYB genes were uploaded based on the results of BLASTP searches in the citrus genome database, and their intron distribution patterns, phases and intro-exon boundaries were analyzed using the GSDS web-based bioinformatics tool (http://gsds.cbi.pku.edu.cn/).
The chromosomal position of genes was provided by the Citrus Genome Database. The distribution of CitRR3MYB genes throughout the sweet orange and clementine genomes was drawn manually.

Phylogenetic Analysis
Based on the aligned sequences of the citrus R2R3MYB proteins, a Neighbor Joining (NJ) tree was constructed using MEGA version 4.0 with a bootstrap of 1000 replicates, aiming to investigate the evolutionary history of the R2R3MYB genes in citrus. In order to predict the functions of the citrus R2R3MYB genes, a combined CitMYB (101 members) and AtMYB (126 members) phylogentic tree was created, also using MEGA 4.0 with NJ method and a bootstrap of 1000 replicates.

Expression Profiling of Citrus R2R3MYB Genes
To investigate the expression profiles of CitMYB genes in response to abiotic stress and plant hormone, citrus (Citrus junos Sieb. cv. 'Ziyang') seeds were peeled, and germinated on moist filter paper in a dark chamber with 28˚C and 100% relative humidity for 6 days. The germinated seeds were sown into nutritive soil and then placed in an illuminated chamber (28˚C, 80% relative humidity and 350 mmol m 22 s 21 light intensity) throughout the experiment, which were irrigated with water every 3 days. When seedlings were at two-true-leaf stage, six treatments were treated, respectively: 200 mM NaCl, dehydrate, 0˚C low temperature, 150 mM abscisic acid (ABA) and 200 mM jasmonic acid (MeJA). Roots and leaves used for RNA extraction were harvested at 0, 1 and 6 h after six treatments, of which the materials collected at 0 h served as the control. The flower, fruitlet, root and leaves of mature were collected for tissue specific expression analysis. All the samples were stored at 280˚C until used.
Total RNA was isolated from different tissues using RNApre pure plant Kit (TIANGEN, China) according to the manufacturer's instructions. Two mg total DNA-free RNA was used to synthesized first strand cDNA with PrimeScript 1 st Strand cDNA Synthesis Kit (TaKaRa, Japan).
Real-time PCR using SAND gene [34] as normalize was performed according to the manufacturer's specifications (SYBR PrimeScrip RT-PCR Kit; TaKaRa, Dalian, Liaoning, P. R. China). SYBR Green PCR was carried out using the iCycler iQ5 real-time detection system (Bio-Rad) for 30 s at 95˚C, followed by 40 cycles of 10 s at 95˚C, 30 s at 60˚C, and 30 s at 72˚C, with a final step at 72˚C for 1 min. Each expression profile was independently verified in 3 replicate experiments performed under indentical conditions. Each relative level of gene expression was calculated by the 2 -DDCt method [35]. The PCR primers were designed outside the conserved region to produce amplification products with 130-200 bp. All primer sequences were detailedly listed in Table S1.
The data obtained were statistically analysed using DPS Version 7.55 (http:// www.chinadps.net; Zhejiang University, Hangzhou, P. R. China). One-way ANOVA and Duncan's new multiple-range test were used to determine significant differences in mean values among materials at 0 h, 1 h and 6 h and P#0.05 was regarded as significant.

Identification of the Citrus R2R3MYB Genes in Citrus Genome
One hundred and twenty-six Arabidopsis R2R3MYB proteins were used as a query to search against the citrus (sweet orange and clementine) genomes at the Join Genome Institute (JGI) (http://www.phytozome.net) with BLASTP program. A total of 128 MYB related sequences with MYB domain were identified. To confirm putative R2R3MYB genes, PROSITE and SMART were employed to search for the amino acid sequences of all 128 proteins, and 101 typical R2R3MYB genes (named CitMYB001 to CitRMYB101) were confirmed. These 101 CitMYB genes were used to further analysis (Table 1), of which 10 genes including CitMYB013, 014, 028, 030, 033, 060, 075, 092, 093 and 097 were specifically present in the clementine genome, one (CitMYB099) specifically in the sweet orange genome and 90 in both genomes ( Fig. 1). In addition, comparative analysis showed that the genome distribution of R2R3MYB genes was highly conserved between the sweet orange and clementine.

Conserved Residues in the MYB Domain
To gain insight into the citrus R2R3MYB binding domains, amino acid sequence alignment was conducted to examine how well conserved the R2 and R3 repeats were in the R2R3MYB proteins within each residue position. As shown in Fig. 2, the basic regions of citrus R2R3MYB domains contained, on average, approximate 100 basic residues, with rare frequency of insertion or deletion. By contrast, the region outside the MYB binding domain was poorly conserved in terms of length as well as amino acid composition. Based on previous reports, the R2 and R3 repeats possessed characteristic amino acids, including a series of evenly distributed and highly conserved Trp (W) residues [2,4]. Within the 101 citrus R2R3MYB proteins, 97 of their R2 repeat sequences contained three tryptophan residues, which located at 4, 25 and 47, forming a hydrophobic core and serve as landmarks in plant MYB binding domain. However, in the R3 repeat, the first tryptophan residue (located at 4) of most members was replaced by phenylalanine. The second (located at 23) and third tryptophan residues (located at 42) were well conserved in almost all citrus R2R3MYB proteins, especially the second one which exist in all members. Phylogenetic Analysis of the Citrus R2R3MYB Family The phylogenetic relationship between the citrus R2R3MYB proteins has been examined by multiple sequence alignment of their whole protein sequences using the NJ method with bootstrap analysis (1,000 replicates). The 101 members of the citrus R2R3MYB family were subdivided into 21 subgroups, designated C1 to C21, according to clades with at least 50% bootstrap support (Fig. 3). Additionally, our results also showed that the phylogenetic trees established with MYB binding domains and whole protein sequences, respectively, were composed of nearly identical subgroups, despite the classification of only a few member varied ( Fig. 3; Fig. S1). This result indicated that the phylogenetic relationship between citrus R2R3MYBs based on the whole protein sequence was mainly decided by MYB binding domains, and those citrus R2R3MYBs within the same subgroup may bind to the same MYB recognition sequence, while the regulatory functions of which probably were divergent because of the dramatic divergence of their C-terminal regions that is the main transcriptional activation domain responsible for functional activity or/and specificity [36]. Thus, the phylogenetic tree built, in this study, with C-terminal regions of citrus R2R3MYBs seems more appropriate for revealing the similarity and divergence of regulatory function of the corresponding proteins. To better understand the functional clades with the citrus R2R3MYB genes, an unrooted NJ phylogenetic tree using bootstrap analysis (1000 replicates) was established by alignments of the whole protein sequences of R2R3MYBs from citrus (101), Arabidopsis (126), apple (214), grape (126), peach (110) and populus (192) as well as 25 well characterized R2R3MYBs of other plant species such as pear, Chrysanthemum, wheat, tobacco, rice and Leucaena leucocephala (Fig. 4). The resulting tree generated 68 subgroups (sequentially termed as S1 to S68) with at least 50% bootstrap support, which was similar to the results previously reported [2,4]. As shown in Fig. 4, 41 out of 68 subgroups were shared with citrus and other plant species. This indicated that most R2R3MYB genes in these species highly conserved during plant evolution. Meanwhile, ten species-specific subgroups such as S8, 10, 11, 15 and 32 were observed, indicating that these R2R3MYB genes may have evolved or been lost in a plant species following divergence. Interestingly, of these ten species-specific subgroups, none contained MYB members from citrus, which suggested that these genes may possess specialized roles in other plant species, while are probably dispensable in citrus. As expected, two CitMYBs (CitMYB007, CitMYB008 and CitMYB022) were not contained in any one of 68 subgroups, the functions of which are worth to detailedly elucidate in future.

Intro-exon Structure of the Citrus R2R3MYB Family
Among 101 citrus R2R3MYBs, up to 97 of them possessed at least 1 intron in the R2 and R3 domains. According to their relative positions and phases, all genes could be grouped into 12 patterns (P1-12) (Fig. 5). By contrast, outside the MYB domain, all but 20 of the 101 citrus R2R3MYBs lacked introns.
Pattern P1-3, composed of one or two intron (s) distributed at two highly conserved specific positions, accounting for approximately 73% of CitMYBs. Patterns P5-P12 contained 1 to 5 introns at varying positions in the R2 or R3 domain, were observed in less than 20% of the 101 CitMYBs. In addition, approximately 7% of CitMYBs have no intron in their R2 or R3 domain, forming the third main intron pattern. Intron phases in regard to codons were also investigated in this study. Fig. 5 showed that in the major splicing patterns P1, P2  and P3, the introns phases were 1 and/or 2, respectively, where the phase at the same position of the R2 domain was 1, and that at the R3 was 2.

Organ-Specific Expression Analysis
The expression profiles of citrus R2R3MYB gene family were also analyzed using root, leaf, flower and fruitlet. The results revealed that most of CitMYBs could express in at least one organ. However, few CitMYB genes, including CitMYB024, 035, 059, 074, 077, 096, did not show expression signals (Fig. 3), suggesting that these genes may be pseudogenes, or may be expressed at specific development or under special conditions. The rest showed remarkably variation in transcript abundance, characterized by high level of transcript abundance in one or some organs and low transcript abundance in others. The wide expression of CitMYB genes in different organs indicated that they may play important roles in the development of all citrus organs.

Expression Profiles of the Citrus R2R3MYB Genes in Response to Abiotic Stresses and Hormones
In this study, the transcript abundances of the 101 citrus R2R3MYB genes in leaf and root of citrus seedling at two-true-leaf stage were investigated under cold (0˚C), drought, NaCl (200 mM), ABA (150 mM) and MeJA (200 mM) treatments. The results showed that almost all genes responded to at least one treatment in root and/or leaf (Fig. 6). Numerous CitMYB genes could be positively regulated by one treatment, while being negatively regulated by others. For instance, the of 101 CitMYBs and the bootstrap values less than 50 are not shown in NJ tree. The tree shows the 21 phylogenetic subgroups (C1-C21) marked with colored backgrounds. Twelve proteins did not fit well into clusters. B, The gene structure is presented by exon (black boxes) and intron (black line between the black boxes). The number indicated the phases of corresponding introns. C, The expression patterns of 101 CitMYB genes in different organs. The letter R above the column of expression data refers to root, L refers to leaf, F refers to flower and G refers to fruitlet. doi:10.1371/journal.pone.0113971.g003 CitMYB038 gene was up-regulated by NaCl, ABA and draught, while downregulated by MeJA. There were a few genes, for example CitMYB085, which could be induced by all treatments, suggesting it is a pleiotropic regulator. Interestingly, none of 101 CitMYB genes has been found to be repressed by all treatments. Several genes only responded to a single treatment in leaf and/or root. The expression of ten CitMYB genes including CitMYB007, 009-011, 061-062, 069, 094, 096-097, was very low under all treatments, especially given that CitMYB096 Analysis of the Citrus R2R3MYB Genes gene also showed no expression signals in different organs, further indicating it is a pseudogene. As shown in Fig. 6, the expression patterns of most CitMYB genes in leaf were significantly different from that in root. For instance, CitMYB037 in leaf was repressed by NaCl and drought and induced by ABA, MeJA and cold, while in root was induced by all treatments, suggesting the mechanism of this gene in response to NaCl and drought was different in these two organs.

Characterization of the Citrus R2R3MYB Family
In this study, 101 citrus R2R3MYB genes were identified and detailedly characterized. The size of the R2R3family in citrus was smaller than that of Arabidopsis (126) [5], populus (197) [7], grape (117) [6], maize (157) [4] and soybean (244) [3], was almost equal to that of rice (102) [5], and was larger than that of cucumber (55) [2], suggesting the R2R3MYB gene family in citrus had shrink compared to Arabidopsis, polar, grape and soybean, but expanded compared to cucumber. the number of R2R3MYB genes in our study may be particularly true, considering that the total gene number predicted in Citrus Sinensis or Citrus clementina (24533) is even lower than that in Arabidopsis (26819). However, our data showed that species-specific R2R3MYB members were present in sweet orange or clementina, which let us undoubtedly believe that some new R2R3MYB members were contained in other citrus species such as mandarin (Citrus reticulata), sour orange (Citrus aurantium), pomelo (Citrus grandis), lemon (Citrus limon), citron (Citrus medica). By this reasoning, the number of the R2R3MYB genes in citrus probably is far more than 101, which remain to be further validated.
The MYB binding domains were hightly conserved in 101 citrus R2R3MYB proteins and most of them possessed characteristic amino acids, which were in line with those from Cucumber [2], Populus [7], Arabidopsis [31] and Triticum [37]. In addition, our data showed that the splicing phases and patterns of citrus R2R3MYB genes were highly conserved during the evolution, which were also observed in the MYB gene families of cucumber [2], soybean [3], maize [4], rice and Arabidopsis [5].
Phylogenic Relationship and Function of Citrus R2R3MYB Family Genes. It is well known that R2R3MYB proteins are involved in a range of different physiological processes, including the response to various stress conditions, secondary metabolism, cell shape and hormone responses. In spite of their importance and large number, very little information has been presented on citrus. Generally, the functions of a gene could be preliminarily predicated through phylogenetic analysis, because the genes grouping into a clade had similar functions and gene structure, and were considered to evolve from a recent common ancestor. By phylogenetic analysis, our results revealed that most subgroups contained the R2R3MYB members from Arabidopsis, citrus, apple, grape, peach, populus and other plant species, indicating their functions were highly conserved during plant evolution. Thus, it is entirely feasible that we could predict the functions of the citrus R2R3MYBs according to their phylogenetic relationship with well characterized R2R3MYBs from other plant species.
It is well known that plant R2R3MYBs widely take part in the control of plant development. Subgroup 43 contained 20 R2R3MYBs, 2 members from citrus, 3 from Arabidopsis, 6 from apple, 2 from peach, 4 from populus and 3 from grape, of which AtMYB93, to date, was well elucidated to negatively regulate lateral root development as an interaction partner of the lateral-root-promoting ARABIDILLO proteins [47]. Subgroup 23 consisted of the members that involved in axillary meristems [1]. In subgroup 21, two members from citrus was grouped with AtMYB35 which was involve in pollen development [48]. The CitMYB077, 078 and 079 in subgroup 60 also seem likely to regulate pollen development [11]. However, whether the abovementioned CitMYBs also have similar functions in the control of plant development remain to further demonstrate. Another wellknown role of R2R3MYBs is the regulation of cell fate. For example, subgroup4 and 61 clustered with several R2R3MYB proteins which potentially function in the determination of sperm cell differentiation, cell shape and trichome branching [49][50][51].
Recently, accumulating data demonstrated that numerous R2R3MYBs were widely involved into plant adaptation and tolerance to biotic and abiotic stresses. For example, Subgroup 45 consists of 24 R2R3MYBs, including AtMYB41 and AtMYB102, implicated in regulating the resistance to draught and insect [52,53]. Another example was provided by subgroup 63, which consisted of 56 members, several of them were involved into abiotic and biotic stresses such as disease, salt, draught [54][55][56][57]. In addition, the citrus R2R3MYB members in subgroup 5, 36, 38, 49 may possess the functions in stress resistance and those in subgroup 22,34,37,46,47,48,55,58, 62 and 65 maybe play roles in morphogenesis, lignin biosynthesis, cuticle development, secondary wall biosynthesis and stomatal development. Surely, the functions of many subgroups such subgroup 3, 6, 7 couldn't be predicted due to the absence of well characterized MYB members in them.

Expression Analysis of CitMYB Genes in Response to Abiotic Conditions
A large number of R2R3MYB proteins from different plant species have been characterized by genetic analysis and have been found to play important roles in various abiotic and biotic stresses [1,16,[53][54][55][56][57]. However, no information, to date, is available about citrus R2R3MYB gene involved into abiotic and biotic stresses. In general, we could preliminarily predict the biological functions of a gene by its expression patterns. For this reason, the expression patterns of the 101 CitMYB genes were investigated under cold, drought, NaCl, ABA and MeJA treatments. The results indicated that most CitMYB genes could be induced by at least one treatment, some of them responded to multiple treatments such as CitMYB022, CitMYB080 and CitMYB085. These CitMYB genes show a promise for improving citrus adaptation to stresses, especially the CitMYB genes that responded to multiple treatments, since plants often undergo multiple stresses concurrently. Additionally, some genes showed opposing expression patterns under different stress conditions, such as CitMYB016, CitMYB017 and CitMYB030, which indicated that they played a major role in the plant response to abiotic conditions and involved in communication between different signal transduction pathways [2].