Genome-Wide Identification of the HMA Gene Family and Expression Analysis under Cd Stress in Barley

In recent years, cadmium (Cd) pollution in soil has increased with increasing industrial activities, which has restricted crop growth and agricultural development. The heavy metal ATPase (HMA) gene family contributes to heavy metal stress resistance in plants. In this study, 21 HMA genes (HvHMAs) were identified in barley (Hordeum vulgare L., Hv) using bioinformatics methods. Based on phylogenetic analysis and domain distribution, barley HMA genes were divided into five groups (A–E), and complete analyses were performed in terms of physicochemical properties, structural characteristics, conserved domains, and chromosome localization. The expression pattern analysis showed that most HvHMA genes were expressed in barley and exhibited tissue specificity. According to the fragments per kilobase of exon per million fragments values in shoots from seedlings at the 10 cm shoot stage (LEA) and phylogenetic analysis, five HvHMA genes were selected for expression analysis under Cd stress. Among the five HvHMA genes, three (HvHMA1, HvHMA3, and HvHMA4) were upregulated and two (HvHMA2 and HvHMA6) were downregulated following Cd treatments. This study serves as a foundation for clarifying the functions of HvHMA proteins in the heavy metal stress resistance of barley.


Introduction
Cadmium (Cd) pollution is one of the negative consequences of industrialization. Cd is highly toxic to plants and is easily absorbed by the roots and accumulates in the tissues [1], which influences various processes including water and mineral uptake, respiration, and photosynthesis, and leads to the inhibition of growth and even death [2]. Cd ions with a lack of specificity enter the plant through other transporters (Fe 2+ /Fe 3+ , Zn 2+ , and Mn 2+ ) and compete with other nutrients for plant uptake, resulting in deficient nutrition [3,4]. In response to Cd poisoning, various defense mechanisms have evolved in plants, such as extrusion across plasma membrane, chelation in the cytosol, and vacuolar sequestration [5]. Previous studies have identified multiple proteins related to Cd transport, including heavy metal ATPase (HMA) [6], yellow stripe-like proteins (YSL) [7], and natural resistanceassociated macrophage proteins (NRAMP) [8], to name a few.
HMA, also known as P 1B -ATPase, is a type of protein combining ATP hydrolysis with metal ion transport across membranes [9,10] participating in absorbing and transporting heavy metal ions (Cu 2+ , Zn 2+ , Co 2+ , Cd 2+ , and Pb 2+ ) [6]. Typical HMA proteins contain the E1-E2 ATPase domain and haloacid dehalogenase-like hydrolase (Hydrolase) domain [11]. Additionally, both sides of the N-terminal and C-terminal metal-binding sites may possess one or more soluble metal-binding domains (MBDs) that interact with or bind to specific metal ions [6]. At present, a number of HMA genes have been identified in plants, including 8 in Arabidopsis thaliana [12], 9 in rice (Oryza sativa L.) [12] 11 in maize (Zea mays L.) [13], 11 in sorghum (Sorghum bicolor L.) [13], 17 in Populus trichocarpa [11], and 20 in soybean (Glycine max L.) [14]. Studies have demonstrated that AtHMA2 and AtHMA4 are two essential genes mediating Cd translocation in A. thaliana [15]. The translocation of Cd from the roots to shoots was near-completely abolished in the hma2 hma4 double mutant. TcHMA3, a tonoplast-localized transporter highly specific for Cd, is responsible for sequestering Cd into the leaf vacuoles so as to detoxify Cd in Thlaspi caerulescens [16]. OsHMA3, which localizes to vacuolar membranes, was identified as the gene that controls root-to-shoot Cd translocation rates in rice [17]. These results indicate that the HMA gene family plays diverse roles in plant resistance to Cd stress.
Barley (Hordeum vulgare L., Hv), an important cereal crop, is widely used in numerous industries, including animal feed, brewing, and food [18]. The exploration of vital genes related to heavy metal stress resistance is beneficial for cultivating Cd-tolerant barley varieties. Recent barley genome sequencing accomplishments have facilitated further studies on barley genomics. Studies on the SBP-box [19], WRKY [20], ABC [21], F-box [22], and SOD [23] gene family in barley have been successfully completed. Although the functions of several HMA proteins in barley have been reported [24,25], genome-wide analysis of the HvHMA family is lacking. In this study, the HvHMA gene family was genome-widely identified in barley, and the phylogenetic relationships, structural characteristics, physicochemical properties, chromosomal location, as well as the tissue expression of identified members were analyzed. Moreover, the expression of some members following Cd treatment was investigated using quantitative real-time polymerase chain reaction (qRT-PCR). These combined analyses of the biological characteristics and expression changes of the HvHMA gene family provide helpful information for studying the function of HvHMA genes and improve the Cd tolerance of barley varieties.

Plant Materials and Treatment
The 'ZJU3' barley variety was used in this study. Seeds uniform in size and with a full shape were selected and sterilized in 2.5% NaClO for 10 min, rinsed with distilled water four times, and then germinated at 28 • C under dark conditions. After 48 h, seedlings with a root length of approximately 0.5 cm were moved to hydroponic culture boxes (day/night temperatures of 26 • C/24 • C, light/dark photoperiod of 14 h/10 h, and light intensity of 18000 Lx). At the one-leaf stage, the seedlings were treated with 1/4 Hoagland's nutrient solution. At the two-leaf stage, Cd stress experiments were performed. The CdCl 2 solutions (50 µmol/L and 100 µmol/L) prepared with Hoagland's nutrient solution were used to simulate Cd stress, and Hoagland's nutrient solution without CdCl 2 was used as the control. After 120 h of treatment, more than 10 barley seedlings were selected for each sample, and quickly stored at −80 • C until analysis. The experiment was performed in triplicate. The Hoagland's nutrient solution formula was as described by Zhang et al. [26].

Phylogenetic Analysis of the Barley HMA Family
The HMA protein sequences of A. thaliana, rice, and barley were imported into the MEGA v7.0 program (Sudhir Kumar, Temple University, Philadelphia, PA, USA) and multiple sequence alignments were performed using ClustalW. The alignment file was then used to construct a neighbor-joining (NJ) phylogenetic tree, with the following parameters: p-distance model, 1000 bootstrap replications, and other default parameters [33]. The tree was displayed and modified using the iTOL website (https://itol.embl.de/; accessed on 12 August 2020) [34].

Quantitative RT-PCR Analysis of Barley HMA Genes
Five pairs of primers related to specific genes were designed using Primer Premier v5.0 (PREMIER Biosoft, San Francisco, CA, USA) for qRT-PCR (Table 1). The barley actin gene HvActin (HORVU1Hr1G002840) was used as an internal control. The qRT-PCR analysis was performed on the CFX96 Real-Time PCR Detection System (Bio-Rad, Hercules, CA, USA), and the data were analyzed using the 2 −∆∆Ct method with three biological replicates [36]. IBM SPSS Statistics v20 (IBM, Armonk, NY, USA) statistical software was then used to analyze significance (*, **, and *** indicates p < 0.05, p < 0.01, and p < 0.001 respectively). Histograms were drawn with SigmaPlot v10.0 (SYSTAT, San Jose, CA, USA).

Identification and Molecular Characteristics of Barley HMA Proteins
Through multiple bioinformatics analyses, a total of 21 barley HMA proteins were screened by removing redundant sequences and validating domains, which were named HvHMA1-21. The basic information on HvHMA1-21, including the number of amino acid residues, molecular weight, theoretical isoelectric point, grand average of hydropathicity, and subcellular localization, which indicated molecular characteristics of barley HMA proteins, is listed in Table 2. The results showed that the 21 HvHMA proteins contained 672 (HvHMA5) to 1083 (HvHMA21) amino acid residues with molecular weights ranging from 73112.29 (HvHMA5) to 118116.93 (HvHMA21) Da. The theoretical isoelectric points of the 21 HvHMA proteins ranged from 5.00 (HvHMA8) to 7.82 (HvHMA17), with the majority constituting acidic proteins. The GRAVY numeric values of the 21 HvHMA proteins varied from −0.090 (HvHMA1) to 0.422 (HvHMA3), indicating that these proteins were likely amphoteric proteins. Additionally, the subcellular localization results showed that the 18 HvHMA proteins were localized in the plasma membrane, whereas three proteins (HvHMA9, HvHMA14, and HvHMA19) were localized in the endoplasmic reticulum.

Phylogenetic Analysis and Classification of Barley HMA Genes
To explore the evolutionary characteristics of HvHMA genes and the evolutionary relationships between the AtHMA, OsHMA, and HvHMA genes, HMA sequences from A. thaliana, rice, and barley, including 8 AtHMA proteins, 8 OsHMA proteins, and 21 HvHMA proteins, were subjected to phylogenetic analysis. As shown in Figure 1, the HvHMA genes were divided into five groups (A-E). Among the 21 HvHMA genes, 2 belong to group A, 6 to group B, 6 to group C, 3 to group D, and 4 to group E. The phylogenetic tree indicated that members of groups A, D, and E were homologous to the AtHMA and OsHMA proteins. Moreover, compared to dicotyledonous A. thaliana, monocotyledonous barley and rice were more closely related. With the exception of the E1-E2 ATPase and Hydrolase domains, some HvHMA proteins contained other domains including an HMA domain, Cation_ATPase_N domain, Cation_ATPase_C domain, Cation_ATPase domain, Hydrolase_3 domain, and CaATP_NAI domain, which revealed that the HvHMA proteins contained more abundant domains than the AtHMA and OsHMA proteins. Thus, it is inferred that members of the barley HMA gene family are more functionally diverse and therefore worth exploring. Additionally, there were some differences among groups in the domain distribution of the HvHMA proteins. The HMA domains were concentrated in groups D and E, the Cation_ATPase_N domains and Cation_ATPase domains were distributed in groups B and C, and Cation_ATPase_C domains were all distributed in group B. In relation to other groups, there were more types of domains in group B, indicating that members of group B might be complex in terms of function.

Chromosomal Location of Barley HMA Genes
According to the genome annotations, 18 of the 21 HvHMA genes were distributed on the six barley chromosomes (Figure 2), with the largest number of genes located on Chr4 (5), followed by Chr5 (4), Chr7 (4), Chr1 (2), Chr2 (2), and Chr6 (1). However, HvHMA7, HvHMA12, and HvHMA16, without clear localization information, were not OsHMAs, and the red letters represent AtHMAs. Different-colored rectangles represent different structural domains. The three green rectangles from light to dark represent Cation_ATPase, E1-E2 ATPase, and Hydrolase_3, respectively; the yellow rectangle represents Hydrolase; the pink rectangle represents HMA; the red rectangle represents Cation_ATPase_C; and the two purple rectangles from light to dark represent Cation_ATPase_N and CaATP_NAI, respectively.

Chromosomal Location of Barley HMA Genes
According to the genome annotations, 18 of the 21 HvHMA genes were distributed on the six barley chromosomes (Figure 2), with the largest number of genes located on Chr4 (5), followed by Chr5 (4), Chr7 (4), Chr1 (2), Chr2 (2), and Chr6 (1). However, HvHMA7, HvHMA12, and HvHMA16, without clear localization information, were not positioned onto barley chromosomes. In addition, most of the HvHMA genes were concentrated on or near the end of the barley chromosomes. These results suggested that the distribution of HvHMA genes was uneven.

Motif Composition of the Barley HMA Proteins
Conserved motif analysis of the HvHMA proteins helped elucidate the conservation as well as the diversification of this family, and a total of 10 distinct conserved motifs were revealed. As exhibited in Figure 3, all HvHMA proteins contained common motifs including motif 1 and motif 10, which suggested that the two motifs might be the characteristic motifs of HvHMA proteins. With the exception of HvHMA21, all HvHMA proteins contained motif 3. Additionally, HvHMA proteins within the same group were generally found to show a similar motif composition. For example, motif 2 was distributed in all groups except group A, whereas motif 4, motif 5, motif 6, motif 8, and motif 9 were unique to group C. Moreover, the differences in motif composition among groups combined with the phylogenetic analysis results supported the reliability of the group classifications and indicated that HvHMA genes in distinct groups might be functionally divergent.

Intron-Exon Structure of Barley HMA Genes
The introns disrupted the coding sequences of most HvHMA genes. As shown in Figure 4, there were some differences among the HvHMA genes in terms of the number and size of the introns, which might be caused by intron deletion and insertion events. With the exception of HvHMA16 without introns, all HvHMA genes contained 2-33 introns (13 with 2-8 introns, 6 with 12-20 introns, and 1 with 33 introns). Additionally, several non-coding regions were distributed in the 18 HvHMA genes, with the exception of HvHMA16, HvHMA18, and HvHMA20.

Motif Composition of the Barley HMA Proteins
Conserved motif analysis of the HvHMA proteins helped elucidate the conservation as well as the diversification of this family, and a total of 10 distinct conserved motifs were revealed. As exhibited in Figure 3, all HvHMA proteins contained common motifs including motif 1 and motif 10, which suggested that the two motifs might be the characteristic motifs of HvHMA proteins. With the exception of HvHMA21, all HvHMA proteins contained motif 3. Additionally, HvHMA proteins within the same group were generally found to show a similar motif composition. For example, motif 2 was distributed in all groups except group A, whereas motif 4, motif 5, motif 6, motif 8, and motif 9 were unique to group C. Moreover, the differences in motif composition among groups combined with the phylogenetic analysis results supported the reliability of the group classifications and indicated that HvHMA genes in distinct groups might be functionally divergent.

Intron-Exon Structure of Barley HMA Genes
The introns disrupted the coding sequences of most HvHMA genes. As shown in Figure 4, there were some differences among the HvHMA genes in terms of the number and size of the introns, which might be caused by intron deletion and insertion events. With the exception of HvHMA16 without introns, all HvHMA genes contained 2-33 introns (13 with 2-8 introns, 6 with 12-20 introns, and 1 with 33 introns). Additionally, several noncoding regions were distributed in the 18 HvHMA genes, with the exception of HvHMA16, HvHMA18, and HvHMA20.

Expression Pattern Analysis of Barley HMA Genes and Target Gene Screening
The analysis of gene expression patterns contributed to gene function prediction. The expression profiles of the HvHMA genes ( Figure 5) revealed that the expression of some genes was tissue-specific. For example, the expression levels of HvHMA5, HvHMA9, and HvHMA15 were high during grain development; the expression levels of HvHMA3, HvHMA6, HvHMA7, HvHMA8, HvHMA13, and HvHMA18 were high in the leaves; HvHMA2, HvHMA12, and HvHMA19 were specifically expressed in the inflorescences; and HvHMA4, HvHMA10, and HvHMA11 were specifically expressed in the tillers, roots, and epidermal strips, respectively. These results indicated that these genes might play specific roles in the relevant tissues. Moreover, the clustering results of the expression data unclearly correspond to the groupings based on phylogenetic analysis, implying that the expression pattern similarity incompletely depended on the sequence similarity. Based on a comprehensive consideration of the FPKM values at the 10 cm shoot stage and phylogenetic analysis, a total of five HvHMA genes were screened for expression analysis under Cd stress.

Expression Pattern Analysis of Barley HMA Genes and Target Gene Screening
The analysis of gene expression patterns contributed to gene function prediction. The expression profiles of the HvHMA genes ( Figure 5) revealed that the expression of some genes was tissue-specific. For example, the expression levels of HvHMA5, HvHMA9, and HvHMA15 were high during grain development; the expression levels of HvHMA3 ,  HvHMA6, HvHMA7, HvHMA8, HvHMA13, and HvHMA18 were high in the leaves; HvHMA2, HvHMA12, and HvHMA19 were specifically expressed in the inflorescences; and HvHMA4, HvHMA10, and HvHMA11 were specifically expressed in the tillers, roots, and epidermal strips, respectively. These results indicated that these genes might play specific roles in the relevant tissues. Moreover, the clustering results of the expression data unclearly correspond to the groupings based on phylogenetic analysis, implying that the expression pattern similarity incompletely depended on the sequence similarity. Based on a comprehensive consideration of the FPKM values at the 10 cm shoot stage and phylogenetic analysis, a total of five HvHMA genes were screened for expression analysis under Cd stress.

Expression Analysis of Barley HMA Genes in Response to Cd Treatment
HMA proteins participate in the distribution of non-essential heavy metal ions

Expression Analysis of Barley HMA Genes in Response to Cd Treatment
HMA proteins participate in the distribution of non-essential heavy metal ions including Cd 2+ , which greatly affect the plant response to heavy metal stress. To analyze the expression of HvHMA genes under Cd stress, five members (Table 1) were carefully selected from 21 HvHMA genes, and qRT-PCR experiments were further performed at the seedling stage. The results ( Figure 6) revealed that three genes (HvHMA1, HvHMA3, and HvHMA4) were upregulated and two genes (HvHMA2 and HvHMA6) were downregulated. Compared with the control, the expression levels of all five genes were significantly different under the 100 µmol/L CdCl 2 treatment, whereas four genes, except for HvHMA4, were significantly different under the 50 µmol/L CdCl 2 treatment. With the increase of Cd concentration, the expression of HvHMA1 and HvHMA4 was significantly higher under the high-concentration stress than under the low-concentration stress. Moreover, the expression changes of HvHMA2, HvHMA3, and HvHMA6 were similar after the two Cd treatments: the expression of HvHMA2, HvHMA3, and HvHMA6 was slightly higher under the high-concentration stress than under the low-concentration stress.
lants 2021, 10, x FOR PEER REVIEW increase of Cd concentration, the expression of HvHMA1 and HvHMA4 wa higher under the high-concentration stress than under the low-concen Moreover, the expression changes of HvHMA2, HvHMA3, and HvHMA6 after the two Cd treatments: the expression of HvHMA2, HvHMA3, and H slightly higher under the high-concentration stress than under the lowstress.

Discussion
The HMA family, which plays a significant role in heavy metal tra widely in plants. In this study, 21 HMA genes were identified in barley. phylogenetic analysis, the 21 HvHMA genes could be divided into five g Compared to groups B and C, members belonging to groups A, D, and E pos homology to the proteins of A. thaliana and rice. Except for the E1-E2 Hydrolase domains, members of groups B and C contained the Catio domains. In addition, members of group B contained the Cation_ATPase_C Cation_ATPase_N domains and the Cation_ATPase_C domains, which are domains, participate in metal ion transport in plants. As indicated by the co analysis, motifs 4-6 and motifs 8-9 were only distributed in group C. T

Discussion
The HMA family, which plays a significant role in heavy metal transport, exists widely in plants. In this study, 21 HMA genes were identified in barley. According to phylogenetic analysis, the 21 HvHMA genes could be divided into five groups (A-E). Compared to groups B and C, members belonging to groups A, D, and E possessed higher homology to the proteins of A. thaliana and rice. Except for the E1-E2 ATPase and Hydrolase domains, members of groups B and C contained the Cation_ATPase_N domains. In addition, members of group B contained the Cation_ATPase_C domains. The Cation_ATPase_N domains and the Cation_ATPase_C domains, which are metal-binding domains, participate in metal ion transport in plants. As indicated by the conserved motif analysis, motifs 4-6 and motifs 8-9 were only distributed in group C. Therefore, it is speculated that the domains as well as motifs unique to groups B and C resulted in the separation of groups B and C from the other groups in the phylogenetic tree. Furthermore, the characteristics of groups B and C illustrated the differences in evolution between barley and other species.
The subcellular localization results revealed that most of the HvHMA members were predicted as plasma membrane proteins, with the exception of HvHMA9, HvHMA14, and HvHMA19, which were all located in the endoplasmic reticulum and were placed in the same group (group B). These results suggested that there were certain corresponding relationships between the phylogenetic groupings based on sequence similarity and subcellular localization. Therefore, homologous genes may be similar in gene function and signal transduction process.
According to the presence or absence of introns, eukaryotic genes can be divided into intron-containing and intronless genes. Most eukaryotic genes belong to the former, but some belong to the latter. Previous studies identified that there are 5846 (21.7%) intronless genes in A. thaliana and 11,109 (19.9%) in rice [37]. Among the 21 HvHMA genes, 20 members with 2-33 introns are intron-containing genes, whereas HvHMA16 is an intronless gene. Some plausible explanations may account for the origin of intronless genes. It has been suggested that intronless genes evolved owing to a loss of introns [38]. Another probability is that intronless genes formed as a result of reverse transcription [39]. During the process of retroposition, mRNAs are reverse-transcribed into cDNAs and inserted into new genomic positions that lack introns [40]. Therefore, it can be inferred that intron loss or retrotransposition events impacted on the intron-exon structures of HvHMA genes, leading to the presence of an intronless gene (HvHMA16) in the barley HMA gene family. By comprehensively analyzing the results of the evolutionary tree and the expression values at the 10 cm shoot stage, five genes were screened that might be related to stress responses. Among them, the expressions of HvHMA1, HvHMA3, and HvHMA4 were significantly upregulated under Cd stress. HvHMA1 was highly homologous to OsHMA2, which participates in the root-to-shoot translocation of Cd [41][42][43]. Compared to the wild-type (WT), the Cd concentration in the grains of OsHMA2-overexpressing rice was decreased by approximately half [44]. HvHMA3 was homologous to OsHMA3, which sequesters Cd into the vacuoles of root cells in rice, thereby controlling the rate of Cd translocation from the roots to shoots [45]. Additionally, HvHMA4 was homologous to OsHMA9, whose expression was induced by a high concentration of Cd [46]. The above results indicated that the changes in expression after Cd treatment were in line with theoretical expectations and that the phylogenetic analysis results were credible.
Cd stress negatively effects plant growth and development. The qRT-PCR results suggested that Cd stress can promote or inhibit the expression of HvHMA genes, which indicated that different HvHMAs exhibit diverse mechanisms to protect barley from stress damage. The expression of HvHMA1 gradually increased with the growth of Cd concentration. Mills [47] found that the corresponding HvHMA gene conferred Cd sensitivity to wild-type yeast due to transport activity. It was speculated that HvHMA1 changed the transport activity of Cd in response to Cd stress in barley. Lei [48] found that HvHMA3 played a crucial role in grain Cd accumulation. In this study, HvHMA3 were significantly upregulated after Cd treatment. Therefore, it was speculated that HvHMA3 might be involved in Cd distribution in stress condition. Following Cd stress, the expression levels of HvHMA2 and HvHMA6 decreased significantly, indicating that Cd stress negatively regulated the expression of the two genes. According to this, it is inferred that the expression regulation pathways related to HvHMA2 and HvHMA6 may be similar. Among experimental groups, HvHMA1 and HvHMA4 were significantly upregulated with increasing Cd concentration, which indicated that the two genes were sensitive to the change of Cd stress. However, the molecular mechanisms of these HvHMA genes in response to Cd stress need further exploration.

Conclusions
The barley HMA gene family was explored herein using bioinformatics analysis. The results revealed the characteristics of the barley HMA gene family in terms of physicochemical properties, phylogenetic relationships, domain distribution, chromosomal location, motif composition, intron-exon structure, as well as tissue expression. Moreover, five HvHMA genes were selected for expression analysis which indicated that the five genes responded differently to Cd stress. HvHMA1, HvHMA3, and HvHMA4 were strongly activated by Cd stress, whereas HvHMA2 and HvHMA6 were significantly restrained. This study preliminarily confirms that HvHMA1, HvHMA3, and HvHMA4 play vital roles in Cd tolerance, providing a theoretical basis for further research on the functions of related genes and the improvement of barley varieties.