Genome-wide analysis and expression profiles of glyoxalase gene families in Chinese cabbage (Brassica rapa L)

The glyoxalase pathway is composed of glyoxalase I (GLYI) and glyoxalase II (GLYII) and is responsible for the detoxification of a cytotoxic metabolite methylglyoxal (MG) into the nontoxic S-D-lactoylglutathione. The two glyoxalase enzymes play a crucial role in stress tolerance in various plant species. Recently, the GLY gene families have well been analyzed in Arabidopsis, rice and soybean, however, little is known about them in Chinese cabbage (Brassica rapa). Here, 16 BrGLYI and 15 BrGLYII genes were identified in the B. rapa genome, and the BrGLYI and BrGLYII proteins were both clustered into five subfamilies. The classifications, chromosomal distributions, gene duplications, exon–intron structures, localizations, conserved motifs and promoter cis-elements were also predicted and analyzed. In addition, the expression pattern of these genes in different tissues and their response to biotic and abiotic stresses were analyzed using publicly available data and a quantitative real-time PCR analysis (RT-qPCR). The results indicated that the expression profiles of BrGLY genes varied among different tissues. Notably, a number of BrGLY genes showed responses to biotic and abiotic stress treatments, including Plasmodiophora brassicae infection and various heavy metal stresses. Taken together, this study identifies BrGLYI and BrGLYII gene families in B. rapa and offers insight into their roles in plant development and stress resistance, especially in heavy metal stress tolerance and pathogen resistance.


Introduction
The glyoxalase system is a ubiquitous pathway in all organisms that consists of the following two enzymes: glyoxalase I (GLYI) and glyoxalase II (GLYII). The major function of this pathway is the detoxification of the potent cytotoxin methylglyoxal (MG) into D-lactate through two sequential reactions [1]. GLYI catalyzes the conversion of MG into S-D-lactoylglutathione (S-LG) with glutathione (GSH). GLYII catalyzes S-LG to yield D-lactate and replenishes the GSH that was consumed in the GLYI reaction step. The functions of the glyoxalases have been studied in animals and microbial systems (Thornalley, 1990 and 12 h dark. The roots were sampled from young seedlings. Fresh flower buds were obtained, and the other tissues were sampled approximately 25 day after flowering, including stems, leaves, siliques and seeds. Three biological replicates of each tested tissue were prepared by harvesting samples from three different individuals. The samples were quickly frozen in liquid nitrogen and stored at -80˚C until RNA isolation. A B. rapa landrace (Wuxianzangcaizi) was used for the heavy metal treatment. Healthy seeds of similar sizes were surface-sterilized, dried and then germinated in sterilized moist filter paper. The seeds were treated with fresh medium supplemented with 20 mL 15 mg/L CdCl 2 or 50 mg/L Pb (NO3) 2 [36]. Seeds treated with an equal amount of distilled water served as controls. Three replicates of 50 seeds were used for each treatment. The treatment and control seeds were cultured in darkness for 24 h at 22˚C and then cultured during the photoperiod (16 h light /8 h dark cycle) for seven days. The shoots and roots from seedlings of similar sizes were harvested separately and washed three times with deionized water. The samples were frozen in liquid nitrogen until the RNA extraction.

Identification and analysis of glyoxalase proteins in B. rapa
Pfam (http://pfam.sanger.ac.uk/) accessions PF00903 for GLYI and PF00753 for GLYII were used for a Hidden Markov Model (HMM) search [34]. The whole genomic sequence of B. rapa was obtained from the Brassica database (BRAD, http://brassicadb.org/brad/) [35]. The GLYI and GLYII genes were extracted from the whole genomic sequence according to the descriptions provided by Wang et al. [37].

Sequences analysis and construction of the phylogenetic tree
Clustal X software (ftp://ftp-igbrmc.u-strasbg.fr/pub/clustalX/) was used for amino acid (aa) alignments. Phylogenetic analysis was constructed with the MEGA 5.05 software using the neighbor-joining (NJ) method and 1,000 bootstrap tests [40].

Sub-cellular localization of the predicted GLYI proteins
The sub-cellular localizations of all predicted BrGLYI and BrGLYII proteins were analyzed using different online tools, i.e., Wolf pSORT [41], TargeP, and ChloroP [42].
To reveal the response of the BrGLY genes to biotic stress in Chinese cabbage, the expression of all BrGLYI and BrGLYII genes in response to pathogen infection was analyzed using the reported RNA-seq data (http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE74044) [46].
The raw data obtained using the tag-based transcriptome sequencing approach were used to confirm the response of the BrGLY genes to the heavy metal stress, which was accessible through the GEO database (http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE55264) [47].

RT-qPCR analyses
Total RNA was isolated using an isolation kit (BioTeke, RP3201). The cDNA was synthesized using a synthesis kit (TransGen Biotech), and the RT-qPCR was carried out as descried by Li et al. [39]. The relative expression of BrGLY was analyzed with the Actin as a housekeeping gene using a previously described method [48]. The specific primers designed are listed in S1 Table.

Identification of the GLYI/GLYII genes in B. rapa
According to B. rapa genome sequence, we identified the GLY proteins in B. rapa. Proteins that contained the glyoxalase domain (Pfam databases, PF00903) and had a putative lactoylglutathione lyase function were classified as BrGLYI proteins. Likewise, proteins that contained the metallo-beta-lactamase domain (Pfam databases, PF00753) and had a putative hydroxyacyl glutathione hydrolase function were classified as BrGLYII proteins. In B. rapa, 16 BrGLYI and 15 BrGLYII genes were identified. The coding sequences and amino acid sequences of BrGLY genes were shown in S2 and S3 Tables.

Detailed information for the identified BrGLYI and BrGLYII genes
We analyzed all the identified BrGLYI and BrGLYII genes in detail. The chromosomal locations, orientation, DNA length, exons and introns, coding DNA sequence (CDS) length, polypeptide (PP) length and isoelectric point (pI) of each BrGLY gene are shown in Table 2. The full DNA sequence length of BrGLYI varied from 555 bp (BrGLYI8) to 8430 bp (BrGLYI6), and their CDS length varied from 414 bp (BrGLYI2 and BrGLYI13) to 3393 bp (BrGLYI6). Accordingly, BrGLYI6 encodes the largest protein of the family (1131 aa, 125.1 kDa), and BrGLYI2 and BrGLYI13 encode the smallest protein (137 aa, 15.24 kDa) ( Table 2). In addition, the proteins showed a large variation in pI value from 4.78 (BrGLYI8) to 8.59 (BrGLYI6). Most of the BrGLYI proteins were acidic, and only four proteins (i.e., BrGLYI2, BrGLYI4, BrGLYI5 and BrGLYI6) showed a basic pI value ( Table 2). Most of the BrGLYI proteins were localized in the chloroplast, followed by the cytosol, nucleus and mitochondria (Table 2). Similarly, the full DNA sequence length of BrGLYII varied from 1433 bp (BrGLYII8) to 4835 bp (BrGLYII12), and the CDS length of BrGLYII varies from 729 bp (BrGLYII13) to 2724 bp (BrGLYII11) ( Table 2). The largest protein (908 aa, 100.1 KDa) of the BrGLYII family is encoded by BrGLYII11, and the smallest protein is BrGLYII13 (243 aa, 26.4 kDa) ( Table 2). BrGLYII proteins also showed a deviation in pI values, which varied from 5.14 (BrGLYII12) to 8.86 (BrGLYII4). Overall, 10 of the 15 BrGLYII proteins showed an acidic pI value, while only five showed a basic pI value. These results are similar to those obtained for the BrGLYI proteins. The localization analysis indicated that BrGLYII proteins localized more in the mitochondria than at the other sites, such as the chloroplast, cytosol and nucleus ( Table 2). chromosomes, which is highly uneven (Fig 1A). Chromosome 6 harbored the most BrGLYI genes (four BrGLYI genes). Chromosome 9 contained three BrGLYI genes, which is ranked second. Two BrGLYI genes are located on chromosomes 7, 8, and 10, and one BrGLYI gene is located on chromosomes 2, 3, and 5 ( Fig 1A and Table 2). No BrGLYI genes were found on chromosomes 1 and 4. Regarding the BrGLYII genes, there were three genes on chromosomes 5 and 6. Two genes were identified on chromosomes 1, 3 and 9, while chromosomes 2, 4, and 8 harbored one BrGLYII gene each (Fig 1B and Table 2). No BrGLYII genes were present on chromosomes 7 and 10.

Chromosomal distribution of the BrGLYI and BrGLYII genes
Duplication events have been previously studied in most plant species. Among the BrGLYI proteins, nine duplicated genes, which shared relatively high sequence similarities (aa identity >90%), were identified in the B. rapa genome (Fig 1A, S4 Table). Three of the duplicated genes were categorized into one group (BrGLYI10/BrGLYI12/BrGLYI7) that exhibited a high sequence similarity (>95%). The other six duplicated genes were divided into three groups (BrGLYI15/ BrGLYI14, BrGLYI4/BrGLYI13, BrGLYI2), each of which contained only two duplicated genes. Three duplicated genes were located on chromosome 9, and two of these duplicated genes were distributed on chromosomes 6 and 10. Chromosomes 3 and 8 had one duplicated gene (Fig 1A).
There are seven duplicated BrGLYII genes in the B. rapa genome. The duplicated genes were divided into three groups, and the aa similarity of the genes in a group was above 95%. One gene group contained three genes (BrGLYII7/BrGLYII6/BrGLYII 5), and the other two groups contained two genes (BrGLYII13/BrGLYII10 , BrGLYII8/BrGLYII4). Chromosome 3 had two duplicated genes, and chromosomes 4, 5, 6 and 8 contained one duplicated gene (Fig 1B, S5 Table).
Additionally, by comparing the GLYI and GLYII genes between Arabidopsis and the B. rapa subgenomes, we found that there are seven BrGLYI genes in least fractionated blocks (LF), four BrGLYI genes are located in the medium fractionated blocks (MF1), and two BrGLYI genes are located in the most fractionated blocks (MF2) ( Table 3). Two of the BrGLYII genes Table 3. Identification of homologous GLYI and GLYII genes between A. thaliana and subgenomes in B. rapa.

A. thaliana
Locus identifier are located in the LF and MF2 blocks, and one gene is distributed in the MF1 blocks ( Table 3).
In addition, only one gene, the AtGLYI4 gene, is triplicated; three genes, including AtGLYI1, AtGLYI11 and AtGLYII, are duplicated in the subgenome of B. rapa ( Table 3). There are no Arabidopsis genes that are homologous to BrGLYI3, BrGLYI6 and BrGLYI14.

Phylogenetic and structure analyses of the BrGLYI and BrGLYII gene families
The gene structure analysis of the BrGLYI and BrGLYII indicated that the BrGLY genes had one to 23 introns except BrGLYI8 (Fig 2 and Table 2). BrGLYI13 and BrGLYI14 only contained one intron. BrGLYI2 gene had the largest number of introns. The BrGLYII genes also contained varied numbers of introns; for example, eighteen introns were identified in BrGLYII1, and four introns were predicted in BrGLYII2. As shown in Fig 2, the GLY proteins that clustered together possess a similar structure. To examine the evolutionary relationships of the GLY genes among the predicted GLY proteins in Chinese cabbage, Arabidopsis and rice, a phylogenetic tree was drawn using their aa sequences. The results indicated that the GLYI and GLYII proteins were divided into five subfamilies (Fig 3). Among the GLYI proteins, the largest clade (Clade I) contained 15 members, whereas the smallest group (Clade IV) contained only two members from Arabidopsis ( Fig  3A). The results indicated that the homology between BrGLYI and OsGLYI was much lower than that between BrGLYI and AtGLYI (Fig 3A). Clade I included six members of B. rapa, whereas four proteins were from Arabidopsis and five proteins were from rice. In this group, AtGLYI4 transcription can be induced by osmotic, extreme temperature and wounding stress. Furthermore, AtGLYI7 is highly up-regulated under salt, osmotic, extreme temperature and wounding stress [34]. Three BrGLYI proteins (BrGLYI7, BrGLYI10 and BrGLYI12) in Chinese cabbage had a high sequence similarity with AtGLYI4. BrGLYI1 showed a high similarity to AtGLYI7. We hypothesized that the similar BrGLYI proteins may play similar roles in the stress response. Group II contained one GLYI protein each in rice and Arabidopsis and two proteins in Chinese cabbage. Group III contained four BrGLYI proteins in Chinese cabbage. The functions of the proteins in this group may be related to salt stress because the OsGLYI 11 protein in this group improved the transgenic tobacco adaptation to lower Na + /K + ratio stress [7]. Group IV only included two Arabidopsis proteins, i.e., the AtGLY5 and AtGLY10 proteins. Three BrGLYI proteins from Chinese cabbage, two GLYI proteins from rice and one protein from Arabidopsis belonged to Group V. OsGLYI3 in this group was found to be stress responsive (salinity stress, oxidative stress, and exogenous MG) in rice, which indicated its possible function in stress tolerance [34].
Similarly, the GLYII proteins formed five distinct clades (Fig 3B). Two BrGLYII proteins were clustered in groups III and IV, three BrGLYII proteins were classified in group V, and only one BrGLYII protein, BrGLYII9, was in Group II, whereas group I included six BrGLYII proteins (Fig 3B).
A phylogenetic relationship analysis revealed that GLYI and GLYII shared a closer relationship at the interspecific level, such as BrGLYI1, BjGLYI and BnGLYI. In addition, the proteins in Chinese cabbage showed a much closer evolutionary distance to Arabidopsis than rice; for example, BrGLYII9 displayed a closer relationship to AtGLYII5 than to OsGLYII (Fig 3).
To further analyze the protein sequence features of BrGLYI and BrGLYII, the conserved motifs of each protein were also identified using MEME (S1 Fig). We found that most proteins in the same group had similar motifs, and the LOGOs of these protein motifs were obtained by MEME (S2 Fig).

Expression profiles of BrGLYI and BrGLYII in different tissues
The transcription level of the BrGLYI and BrGLYII genes was analyzed using genome-wide transcription profiling data of Chinese cabbage (B. rapa). The expression data in roots, stems, leaves, flowers, siliques and callus were supplied. The Fragments Per Kilo base of exon sper Million fragments mapped (FPKM) values of the BrGLY gene are shown in Fig 4 and S6 Table. The expression clustering clearly reveals that the BrGLYI and BrGLYII genes were classified into different groups (Fig 4). By comparing the expression clustering data with the phylogeny analysis, we found that there was no direct correlation between the gene evolution and expression profiles. BrGLYI3, BrGLYI5, BrGLYI9, BrGLYI11, BrGLYI14 and BrGLYI15 showed a high level of ubiquitous expression during all developmental stages. Of the 16 BrGLYI genes, only the expression of BrGLYI4 was undetectable during the six stages. BrGLYI13 showed a very faint expression in siliques, whereas BrGLYI8 showed a very faint expression in flowers and siliques. BrGLYI6 showed a weak expression during all developmental stages, except for callus. Certain genes showed tissue-specific expression; for example, BrGLYI12 was a root-specific gene. The expression of the other BrGLYI genes showed variable expression levels across different tissues. In contrast, the BrGLYII genes were all expressed in the six organs with variable expression levels; However, BrGLYII3, BrGLYII14 and BrGLYII15 were weakly expressed, and BrGLYII1, BrGLYII4 and BrGLYII7 showed a low expression level during all development stages (Fig 4).
To determine the expression patterns of the BrGLYI and BrGLYII genes obtained from the GEO data, we performed a RT-qPCR analysis of several genes from seven different organs (roots, stems, leaves, flower bud, siliques, silique wall and seeds) of B. rapa. After verifying the specificity for each primer pair, suitable RT-qPCR primer pairs for a total of 11 BrGLYI genes and 5 BrGLYII genes were selected (S6 Table). The expression of the other genes was not detected due to the unspecific primer design. The PCR products amplified ranged from 80 to 250 bp (S6 Table). According to the data, the expression pattern of the different BrGLYI and BrGLYII genes varied among the tissues (Fig 5). The expression of BrGLYI4 was undetected, and the BrGLYI8 and BrGLYI6 genes were faintly expressed, which was consistent with the GEO data (Fig 5). BrGLYI9 appeared to be expressed only in the root and stem, and BrGLYI15 showed a lower expression level in the root. In addition, the expression of several genes could not be detected in certain tissues, e.g., the expression of BrGLYI1, BrGLYI2 and BrGLYI6 was not detected in the root (Fig 5, S6 Table). The above-mentioned results were consistent with the GEO data. The expression patterns of BrGLYI2, BrGLYI3, BrGLYI9, BrGLYI11 and BrGLYI14 were similar to the GEO data (Fig 4 and Fig 5). However, several genes showed a lower or higher expression level in specific tissues, which was inconsistent with the GEO data. For example, BrGLYI6 and BrGLYI11 showed particularly high expression levels in the siliques, and BrGLYI9 was strongly expressed in the silique wall; BrGLYI5 and BrGLYII2 were highly expressed in the root. BrGLYI8 and BrGLYI11 were much more highly expressed in the flower buds (Fig 4 and Fig 5). Furthermore, the expression of BrGLYI12 was obviously inconsistent with the public data due to its constitutive expression (Fig 4 and Fig 5).

Expression analysis of the BrGLY genes under stress conditions
To reveal the response of the glyoxalase genes to biotic and abiotic stresses in Chinese cabbage, the expression of all BrGLYI and BrGLYII genes in response to stress conditions (including P. brassicae infection and FeD, ZnD, ZnE and CdE stress) were analyzed using the publicly available data regarding GSE74044 and GSE55264 in the GEO database. Among all BrGLY genes, 14 BrGLYI and BrGLYII genes were analyzed after P. brassicae infection, and 10 BrGLYI and 11 BrGLYII genes were analyzed under heavy metal stress. Different BrGLY genes showed diverse expression levels under these stresses (Fig 6 and Fig 7).
The BrGLYI genes in the two different lines showed relatively similar expression patterns (Fig 6A and 6B). In the clubroot-resistant line, four, five and seven BrGLYI genes showed an up-regulation at 12 hai, 72 hai and 96 hai after the P. brassicae infection, respectively. Among the up-regulated genes, the expression of BrGLYI1, BrGLYI11 and BrGLYI16 was up-regulated by more than 1.5-fold compared with their corresponding expression under the control condition. In the clubroot-susceptible line, four, five and seven BrGLYI genes showed an up-regulation at 12 hai, 72 hai and 96 hai after the P. brassicae infection, respectively. BrGLYI1, BrGLYI6, BrGLYI11 and BrGLYI16 were up-regulated by more than 1.5-fold at different infection times. Most interestingly, the BrGLYI2 gene was down-regulated in the clubroot-resistant line; however, it was up-regulated in the clubroot-susceptible line.
In the case of BrGLYII, ten, eight and nine BrGLYII genes were induced in the clubroot-resistant line, while six, seven and seven BrGLYII genes were induced in the clubroot-susceptible line at 12 hai, 72 hai and 96 hai, respectively (Fig 6C and 6D). BrGLYII8 and BrGLYII10 were highly expressed after P. brassicae infection at 72 hai in both lines. BrGLYII13 was more highly induced in the resistant line than susceptible line at 12 hai. Moreover, the expression of BrGLYII15 was upregulated in the resistant line; however, it was down-regulated in the susceptible line.
The expression of several BrGLYI genes was induced under the heavy metal stress conditions. FeD causes seven BrGLYI genes to become up-regulated and two BrGLYI genes to become down-regulated. Two BrGLYI genes were induced and four BrGLYI genes were downregulated under the ZnD condition. Six BrGLYI genes were induced and four BrGLYI genes were down-regulated under the ZnE condition. The expression level of two BrGLYI genes increased and that of eight BrGLYI genes decreased under the CdE condition. Among these genes, BrGLYI1 was induced by more than 2-fold under the ZnD condition, whereas BrGLYI13 was induced by more than 1.5-fold under the ZnE condition. BrGLYI15 was significantly induced under the FeD condition (Fig 7A).
By analyzing the response of the BrGLYII genes to the heavy metal stress, we found that the expression of BrGLYII5, BrGLYII10 and BrGLYII13 was significantly up-regulated (over 1.5-fold) under the FeD condition. BrGLYII4, BrGLYII5, BrGLYII9, BrGLYII10 and BrGLYII13 showed an up-regulation under the ZnD stress condition, whereas BrGLYII5 showed a significant up-regulation under the ZnE stress condition. BrGLYII11 was induced by approximately 1.5-fold under the CdE stress condition. These results illustrate the diverse responses of different BrGLY genes in the stress regulatory pathways in Chinese cabbage. Among the BrGLYII gene members, BrGLYII5 was induced under the FeD, ZnD and ZnE stress conditions, which suggested that it may play a crucial role in heavy metal stress and its function requires further validation (Fig 7B). To verify the response of the glyoxalase genes to heavy metals, RT-qPCR was performed to validate the nine candidate BrGLYI genes (BrGLYI-1, 2, 3, 4, 5, 6, 8, 11, 15 and 16) under the Pb and Cd treatment conditions (Fig 8). In the shoot, BrGLYI8 were significantly up-regulated under the Cd condition and were approximately 1.8-fold higher than the expression under the control condition. The expression of BrGLYI11 had no significant change under the Cd stress condition. Moreover, the expression of BrGLYI3 and BrGLYI6 showed significant increase under the Pb treatment compared with that in the control. Most interestingly, although the expression of BrGLYI1 was almost undetected under the control conditions, its expression level was significantly induced under the Pb treatment conditions. In addition, the expression of BrGLYI16 showed no significant difference under Cd and Pb treatments compared with control. In the root, BrGLYI1, BrGLYI6 and BrGLYI8 showed a significant up-regulation in response to Pb stress, and the expression level of BrGLYI11 did not show any change under the Pb stress condition. However, the other BrGLYI genes were clearly suppressed under the stress conditions (Fig 8). These results indicated that BrGLYI1, BrGLYI6, BrGLYI8, BrGLYI11 and BrGLYI16 may play an important role in heavy metal resistance.

Analysis of the regulatory elements in the BrGLY promoter
The cis-acting elements in promoter regions are known as regulation of gene transcription and their response to stress. Therefore, an analysis of 19 stress-responsive cis-acting elements in each BrGLY gene promoter was performed using PlantCARE database [40,43], including ABRE, ACE, AE-box, AuxRR-core etc. (Fig 9). All these elements played a critical role in regulating gene transcription induced by various biological processes, such as biotic and abiotic stress responses, developmental processes, etc. Thus, the preliminary analyses of these elements will be helpful for understanding the gene , responses to different stresses [49,50]. These elements are distributed randomly in the BrGLY promoter sequences (including both positive and negative strands) without following a particular rule (Fig 9). Among the BrGLYI genes, the BrGLYI12 promoter only has eight elements, while BrGLYI2 and BrGLYI13 have the maximum number of cis-elements (21 elements). Among the BrGLYII members, BrGLYII12 has a maximum of 28 elements, while BrGLYII2 has only eight cis elements. Almost every promoter region in the BrGLY genes contained ARE, Skn-1_motif and TGACG motif. Although the relationship between these elements and the responses of the genes under stress conditions requires further experimental investigation, our analysis results suggested that the BrGLY genes had a certain stress-responsive characteristic.

Discussion
The genus Brassica is one of the most significant genera and is grown because its seeds, oil and vegetables have high nutritional value and include nutrients such as iron, vitamins, phytosterols and fiber [51]. The genus Brassica comprises six crop species. Among them, B. rapa (AA), B. nigra (BB), and B. oleracea (CC) were the ancestors of the three amphidiploid species B. napus (AC), B. juncea (AB) and B. carinata (BC) [52]. Chinese cabbage (B. rapa subsp. pekinensis), which is a type of B. rapa, is one of the most important vegetable crops in the world. The 'A' genome of B. rapa is valuable for gaining a better understanding of the genetic evolution of Brassica and expediting the genetic improvement of Brassica crops. Recently, many genomes of crop species, including B. rapa, have been sequenced, and the data have been released. Furthermore, bioinformatics analyses have developed rapidly. Therefore, we have the ability to identify large gene families in these species systematically.
MG is a cytotoxic metabolite generated from carbohydrate and lipid metabolism [1]. Previous reports have indicated that the level of MG increases when plants encounter various abiotic stresses [53]. The glyoxalase system, which contains GLYI and GLYII, can detoxify MG into D-lactate. The glyoxalase activity can be up-regulated under stress conditions, which reduces MG accumulation and protects plants from MG damage to a certain extent [10,16,26,[53][54][55]. Therefore, the two genes are suggested to be important candidate genes for improving plant tolerance by gene engineering. Recently, a genome-wide identification of the GLY gene has been performed preliminary in Arabidopsis and rice [34]. The analysis was also completed in soybean (Glycine max) and the results illustrate their developmental and stress specific responses [40]. However, the two gene families have not been analyzed in any other plant, including Brassica plants. In our study, GLY gene families in Chinese cabbage were identified at the genome level. The chromosomal location, gene structure, protein localization, protein motifs and expression patterns were then analyzed. In this study, we found 16 BrGLYI and 15 BrGLYII genes in Chinese cabbage. In a previous report, Arabidopsis and rice were shown to contain 11 GLYI genes and five and three GLYII genes [34], respectively. The number of GLY genes was lager in B. rapa than that in Arabidopsis, particularly the GLYII gene. Moreover, AtGLYI10 did not have homologous genes in the three B. rapa subgenomes, and one BrGLYI and nine BrGLYII genes did not show homology to the AtGLYI genes. Therefore, the processes of polyploid evolution are likely accompanied by gene mutations and losses in addition to duplications or triplications.
The glyoxalase system is located in cellular organelles and cytoplasm. The widespread distribution of the GLY protein in living organisms indicates that it fulfills a function that is important to biological life. Previous studies suggested that the glyoxalase enzymes play a crucial role in tissue proliferation, cell division and malignancy [56][57][58]. In higher plants, GLYI activity was reported to be related to the cell division in pea, a Datura callus suspension and Brassica [16,17]. Subsequently, the effects of GLYI on cell division and hormone levels were confirmed in soybean cell-suspension cultures [59]. In our study, we found that the expression of the five BrGLYI genes (BrGLYI1, BrGLYI7, BrGLYI10, BrGLYI14 and BrGLYI15) in callus was much higher than that in the other tissues (FPKM > 200); however, BrGLYII did not show a similar expression pattern in callus. These results indicated that BrGLYI may play an important regulatory role in cell division as previously reported; however, its precise regulatory mechanism in cell division remains unclear and requires further study.
To investigate the response of the glyoxalase genes to various abiotic stress factors at the transcription level, the expression patterns of the BrGLYI and BrGLYII genes were analyzed using publicly available expression data and RT-qPCR. The expression of BrGLYI4 was undetected. It may be that BrGLYI4 had no expression or had spatial and temporal expression patterns. BrGLYI8 and BrGLYI6 expressed faintly; however, they were up-regulated under the Pb and Cd treatments. The two genes were selected to further study their functions. Several genes showed a high expression in specific tissues, such as the expression of BrGLYI6 in siliques and that of BrGLYI9 in silique walls. The abundant transcription of a gene in a specific organ usually suggests that the gene may play an important role in the development of the corresponding tissue. Many genes were highly expressed in more than one tissue and some genes were constitutively expressed in all the seven tissues, such as BrGLYII12 and BrGLYII14. These genes may be required for development throughout the whole life. Moreover, the expression patterns of several BrGLY genes were inconsistent between the RT-qPCR and GEO data, such as the pattern for BrGLYII2. The possible reasons may be as follows: first, the plant materials were not sampled at precisely the same time, and some genes showed spatial and temporal expression patterns, and, second, the GEO data may not be specific to a gene because highly homologous genes might be difficult to distinguish.
The expression analysis of the GLY genes under the biotic and abiotic stress conditions showed that several GLY genes were stress responsive. BrGLYI1, BrGLYI2, BrGLYI6, BrGLYI11 and BrGLYI16 were up-regulated by more than 1.5-fold at different times when infected by P. brassicae in both the clubroot-resistant and clubroot-susceptible lines. Moreover, BrGLYII8 and BrGLYII10 were expressed at a high level after P. brassicae infection in both lines. Previous studies also showed that the GLYI genes were induced by pathogenic microorganism [9,60,61]. Thus, glyoxalases may play a crucial role in defending plants against infection by pathogens [21,62,63], and their function in plant disease resistance requires further investigation. In addition, BrGLYI1 was significantly up-regulated under the ZnD condition, BrGLYI13 was induced under the ZnE condition, and BrGLYI15 was significantly induced under the FeD condition. Moreover, BrGLYII5 was the most stress-inducible gene and was induced under the FeD, ZnD and ZnE stress conditions. The RT-qPCR analysis indicated that BrGLYI1, BrGLYI3, BrGLYI6 and BrGLYI8 were up-regulated under the Cd and Pb treatment conditions. In summary, using different Chinese cabbage varieties, we found that BrGLYI6 and BrGLYI1 may play an important role in tolerance to clubroot disease and heavy metal stress. The results will facilitate further functional exploration of these candidate genes in stress tolerance.
Moreover, many studies have confirmed that the glyoxalase pathway plays an important role in stress tolerance. In plants, previous reports have shown that transgenic plants overexpressing the GLYI genes have an improved tolerance to stress. Transgenic tobacco and V. mungo overexpressing GLYI from B. juncea had a high salt tolerance [3,23]. Tobacco transgenically overexpressing GLYI and GLYII showed an enhanced tolerance to salinity and MG stress compared to that in wild type plants. Furthermore, when GLYI from rice, wheat, and sugar beet was expressed in tobacco, the transgenic tobacco showed an increased tolerance to salinity, heavy metal and MG stress [7,9,28]. Recently, we found that BnGLYI-3 transgenic yeast cells enhanced their tolerance to extreme temperature stress [2]. Jain et al. have found that the overexpression of AtGLYI2, AtGLYI3 and AtGLYI6 in Escherichia coli provides multistress tolerance (including salinity, exogenous MG, oxidative, mannitol and heat stress) [64]. Thus, the glyoxalase pathway is directly related to stress resistance in plant. In our study, certain BrGLY genes shared a high similarity with previously reported genes, and we speculated that these genes may have a similar function in Arabidopsis; for example, BrGLYI5 shared approximately 87% identity with ATGLYI2 and BnGLYI-3. Similarly, BrGLYI11 and BrGLYI9 showed 93% and 87% identity with ATGLYI3 and ATGLYI6. In addition, further investigations should explore the mechanism of the response of the glyoxalase pathway to stress tolerance in plants to generate more stress-tolerant varieties using molecular approaches.

Conclusion
We conducted a comprehensive analysis of glyoxalase gene families (BrGLYI and BrGLYII) in Chinese cabbage and then characterized 16 BrGLYI and 15 BrGLYII genes based on a genome wide sequence analysis. Detailed information, including chromosomal distribution, gene structure, duplication, phylogenetic relationships, conserved motifs, promoter cis-elements and the expression profiling in different organs and under biotic and abiotic stress conditions, was predicted and analyzed. Based on the phylogenetic analysis, the presence of conserved motifs and their corresponding expression, we provided insight into the possible function of these gene families in plant development and responses to specific stresses (pathogen infection and heavy metal stress). Our data shed light on the selection of candidate genes for stress tolerance and lay the foundation for further functional investigation on the Glyoxalase genes.