In-silico genome wide identication, expressional proling and regulation of NHX (Sodium/ Hydrogen Antiporter) gene family in Camellia sinensis

Salt stress affects the plant growth and productivity worldwide and NHX is one of those genes that are well known to improve salt tolerance in transgenic plants. It is well characterized in several plants such as Arabidopsis and cotton however not much is known about NHXs in tea plant. In the present study, NHX genes of tea were obtained through a genome wide search using Arabidopsis thaliana as reference genome. Out of the 9 NHX genes in tea, 7 genes were localized in vacuole while the remaining 2 genes were localized in the endoplasmic reticulum (ER) (TEA014468.1) and plasma membrane (PM) (TEA006997.1) respectively. Furthermore, phylogenetic relationships along with structural analysis which includes gene structure, location as well as protein conserved motifs and domains, were systematically examined and further, predictions were validated by the expression analysis. The dN/dS values show that the majority of tea NHX genes are subjected to strong purifying selection under the course of evolution. Also, functional interaction was carried out in C. sinensis based on the orthologous genes in Arabidopsis. The expression proles linked to various stress treatments revealed wide involvement of NHX genes from tea in response to various abiotic factors. study provides the targets for further comprehensive identication, functional study, also


Result
In the present study, NHX genes of tea were obtained through a genome wide search using Arabidopsis thaliana as reference genome. Out of the 9 NHX genes in tea, 7 genes were localized in vacuole while the remaining 2 genes were localized in the endoplasmic reticulum (ER) (TEA014468.1) and plasma membrane (PM) (TEA006997.1) respectively. Furthermore, phylogenetic relationships along with structural analysis which includes gene structure, location as well as protein conserved motifs and domains, were systematically examined and further, predictions were validated by the expression analysis. The dN/dS values show that the majority of tea NHX genes are subjected to strong purifying selection under the course of evolution. Also, functional interaction was carried out in C. sinensis based on the orthologous genes in Arabidopsis. The expression pro les linked to various stress treatments revealed wide involvement of NHX genes from tea in response to various abiotic factors.

Conclusion
This study provides the targets for further comprehensive identi cation, functional study, and also contributed for a better understanding of the NHX regulatory network in C. sinensis.

Background
Excessive use of inorganic fertilizers is making the land infertile and unavailable for agriculture due to over accumulation of salts in it [1]. Moreover, abiotic stresses such as drought and heat stress etc. cause an additive effect and overall decrease the crop yield and quality [2]. Therefore, to keep up with growing demands of population, there is a pressing need to identify and characterize more salt tolerant genes from different plant species and use them for improvement of salt tolerance in crop plants.
Sodium chloride is one of major salts present in the soil and most of the salt tolerance mechanism focuses on the transport and compartmentalization of sodium ions. Na + in ux is controlled by either NHX family of cation/H + transporters [3] or nonselective cation channels (NSCCs), or high-a nity K + transporters (HKTs) [4]. HKT can regulate the long-distance transport of Na + [5] while Na + / H + antiporter (NHX) are involved in the transport of Na + ions from cytoplasm to vacuole or outside of the cell. To achieve this, it utilizes the H + electrochemical gradient formed by two proton pumps i.e. H + -ATPase and H + -PPase thereby avoiding the cell from the toxic effects of sodium ions [3].
NHX proteins belong to the cation/proton antiporter 1 (CPA1) superfamily and most of NHX proteins possess 10 transmembrane helices [6,7,8,9]. Localization of NHX proteins is mainly restricted to plasma membranes, vacuoles and endosomes [10,11]. The rst plant NHX gene was recognized in barley root tips [12] followed by its identi cation and characterization in Arabidopsis thaliana (At) [13], and a total of 8 NHX genes have been reported in Arabidopsis till date. Out of 8 NHXs in Arabidopsis, 2 genes (AtNHX7 and AtNHX8) belong to PM-class (plasma membranes), 2 genes (AtNHX5 and AtNHX6) belong to Endoclass (endosomes), and 4 genes (AtNHX1-4) belong to Vac-class (vacuoles). This classi cation is done on the basis of their subcellular localization [6, 10, 14, & 15]. Apart from the involvement of these genes in salt tolerance, NHX antiporters are involved in the regulation of wide variety of physiological processes such as vesicle tra cking, pH regulation, K + homeostasis, protein transport, growth/ development [16, 17, 18, & 19] C. sinensis is native to East Asia, the Indian Subcontinent, and Southeast Asia, but it is today cultivated across the world in tropical and subtropical regions. Tea plant (Camellia sinensis L.) is an important economic crop, leaves of which are an important source of non-alcoholic beverage. As a leaf-harvested crop, tea plant is unavoidably threatened with various adverse environment stresses throughout the whole life cycle, such as drought [20], salt [21] and cold [22] stresses, which critically hinders the development of the tea industry. With drastic environmental changes leading to a decline in the cultivated land area, like many other economic crops, tea planting elds are moving to salinity and drought-affected areas. In this study, we performed a genome-wide analysis of NHX genes in C. sinensis including the phylogenetic relationships, a motif analysis, promoter analysis, gene expression pattern and the gene structures. Through a systematic analysis of all the members of the NHX gene, we can understand the gene regulation, expression pattern, and eventually their biological functions in tea.

Results
Genome wide identi cation of NHX genes in C. sinensis In order to retrieve the members of the NHX gene family in tea, the published NHX protein sequences of Arabidopsis (8) and rice (7) were retrieved from TAIR database (https://www.arabidopsis.org/) and Rice Genome Annotation Project database (http://rice.plantbiology.msu.edu/) respectively. These peptide sequences were then used as queries to search against the genome database of tea, Tea Plant Information Archive (TPIA) (http://tpia.teaplant.org/) by making use of the BLASTp algorithm with e value and identity percentages set to 1e-5 and 50% as threshold respectively (Additional File 1: Table S1). The tea NHX peptide sequences identi ed were further screened using the Hidden Markov Model (PF00999) to con rm the presence of the Na + /H + _Exchanger domain. Based on the results, 9 putative tea NHX genes were incorporated into the nal dataset.
The physicochemical properties of the identi ed tea NHX protein sequences were evaluated and analysed by the ExPASy ProtParam tool ( Table 1). The length of the NHX peptide sequences ranged from 201 (TEA011468.1) to 1204 (TEA021179.1) amino acid residues while the molecular weights varied from 21764.56 (TEA011468.1) to 134630.87 (TEA021179.1) kDa. The predicted isoelectric points (pI) values ranged from 5.82 (TEA012245.1) to 8.79 (TEA012286.1). 5 out of the 9 NHX peptide sequences had more positive residues than negative ones, 3 had more of negative residues and remaining one (TEA006997.1) had equal number of positive and negative residues. All the 9 NHX peptide sequences had positive grand average of hydropathy (GRAVY index) values, ranging from 0.209 (TEA021179.1) to 0.695 (TEA011468.1). This indicated that all the 9 NHX peptides identi ed are hydrophobic in nature. The instability index scores revealed that 2 out of 9 NHX peptides (TEA012286.1 and TEA006997.1) were above 40 while the rest 7 had scores below the given level, indicating that most of the screened peptides had a stable nature [23]. The aliphatic index of the peptides ranged from 102.02 (TEA006997.1) to 114.38 (TEA011468.1). The subcellular localization revealed that most of the NHX genes in tea were localized in vacuole (7 out of 9), while the remaining 2 genes were localized in the endoplasmic reticulum (ER) (TEA014468.1) and plasma membrane (PM) (TEA006997.1) respectively. Additionally the presence of transmembrane helices was also analysed and it revealed that all the NHX peptides had a considerable number of transmembrane helices, ranging from a minimum of 6 in TEA011468.1 to a maximum of 12 in TEA006997.1 (Additional File 2: Fig. S1).

Phylogenetic analysis of tea NHXs
To explore the evolutionary relationships of the NHX genes among the different plant species, a phylogenetic analysis was conducted comparing the identi ed tea NHX genes along with NHXs from 10 other plants. For this study, we retrieved the NHX peptide sequences from Arabidopsis thaliana (At), Oryza sativa (Os), Solanum lycopersicum (Sl), Solanum tuberosum (St), Medicago truncatula (Mt), Populus trichocarpa (Pt), Gossypium hirsutum (Gh), Sorghum bicolor (Sb), Zea mays (Zm) and Glycine max (Gm) from their respective genome databases. The sizes of all the NHX gene family from the 11 members ranged from a minimum of 5 in Solanum tuberosum to a maximum of 23 in Gossypium hirsutum ( Table  2).
The phylogenetic tree was then constructed using all the 100 NHX peptide sequences from the 11 species. MEGA 7.0.14 was used to generate the phylogenetic trees, using the Neighbor-Joining (NJ) algorithm, at default parameters and 1000 bootstrap replicates. The phylogenetic tree shows a direct relation with the subcellular localization as all the NHX peptides clustered into 3 different clades based on their localizations (Fig. 1). The 3 different clades were the Vac-class (Vacuole), Endo-class (Endosomal) and PM-class (Plasma membrane). Among these 3 classes, the Vac-class was the most abundantly present class of NHXs in all the 11 species with 71 genes, followed by the Endo-class and PM-class with 18 and 11 genes respectively. Motif composition of tea NHXs To evaluate the structural characteristics and diversity of the tea NHXs, a correlative study of the conserved motifs from the NHX peptides of Arabidopsis thaliana, Camellia sinensis and Oryza sativa was conducted using the MEME suite (

Gene structure analysis of tea NHXs
To identify the structural characteristics of the tea NHXs, the intron/exon architecture of the genes were analysed using Gene Structure Display Server v2.0. Study of the intron/exon patterns revealed some signi cant differences concerning the number of introns and exons, which further contributes to the variation in gene lengths. Abundant presence of non-coding sequences within a genome is regarded to be an indicator of genome complexity [24,25,26]. Analysing these intron arrangements thereby provides signi cant information regarding the evolution, regulation and function of the NHXs [27,28,29,30]. The analysis of the tea NHX gene structures indicated considerable differences with respect to the number of introns and exons across the 3 classes (Fig. 3A). Among the 9 tea NHXs, only TEA012938. among the 3 classes with only 6 exons and 5 intron segments. However, the PM-class NHX (TEA006997.1) had the most share of IEs with 25 exons and 24 introns. It was observed that the genes belonging to the same clade had a similar distribution of introns and exons. The intron segments and exon lengths were relatively conserved among the genes of the same class. Additionally, analysing the amino acid sequence identity also supported the sequence conservation among the tea NHXs (Fig. 3B).

Retrieval of tea NHX promoter regions and analysis of CAREs
Cis-acting regulatory elements (CAREs) plays a key role in determining gene regulation, function, transcription and gene expression [31,32]. Analysis of these regulatory elements help in de ning the plant responses to various environmental stimuli, stress factors, thereby affecting the growth regulation [33]. To explore the transcriptional potential of the tea NHX genes, the promoter sequences of 2000 bp upstream of the transcriptional start codon "ATG" were retrieved from the TPIA database. These promoter sequences were then used to predict and analyse the CAREs using the PlantCARE database. 41 total CAREs were identi ed randomly distributed across the promoter regions of the 9 tea NHXs (Additional File 1: Table S2). Based on the speci c biological functions of the identi ed CAREs, they were grouped together into a pie chart under 20 different sections (Fig. 4A). Most of the CAREs had sequence lengths of 6 and 9 bp, while the others ranged between 5 to 13 bp (Fig. 4B). Analysing the 41 CAREs, it was observed that 18 elements were involved in light responsiveness, 9 elements in phyto-hormonal as well as plant growth and regulation each and 5 elements in stress response. The light responsive elements had the largest share of CAREs and were present in all the 9 tea NHXs. Among these 18 light responsive elements, the Box-4 and G-box elements were abundantly present in 8 and 6 NHXs respectively. Few of the other light responsive elements were TCCC-motif, AE-box, AT1-motif, Box-II, TCT-motif, chs-CMA1a and chs-CMA2a in 4, 2, 2, 1, 1, 3 and 3 tea NHXs respectively. NHXs are mostly involved in response to various environmental stresses and regulation [33]. The stress responsiveness elements comprised of elements responding to drought stress (MBS), low temperature (LTR), defense and stress (TC-rich repeats) and anaerobic induction (ARE) in 1, 2, 3 and 7 tea NHXs correspondingly. Another element was involved in Genomic distribution map and evolutionary pressures on tea NHXs In an attempt to understand the genome distribution pattern of the tea NHXs, the genes were mapped onto their genomic scaffolds. Due to the lack of chromosome-level assembly data in the TPIA database, the genes had to be mapped onto their scaffolds instead of the chromosomes. The 9 tea NHXs were distributed evenly across 9 different scaffolds (Fig. 5). The genes were positioned such that a single scaffold housed individual genes. Additionally, the Ka/Ks or dN/dS (non-synonymous substitution rate/ synonymous substitution rate) ratios were calculated in-order to understand the evolutionary pressures and gene divergence mechanisms (Additional File 1: Table S3). The dN/dS ratio helps determine whether Darwinian selection pressures were involved in the duplication events [26,34]. If the value of the dN/dS ratio is > 1, it implies a positive or Darwinian selection. If the ratio is equal to 1, it implies a neutral selection and if the ratio is < 1, it determines a purifying selection [35,36]. Pairwise comparisons of the 9 tea NHXs revealed 13 gene pairs having their dN/dS ratios > 1, indicating a positive selection. The rest 23 gene pairs had their ratios < 1, indicating a negative or purifying selection. Additionally, a cumulative graph of the tea NHXs was also generated (Additional File: 2 Supplementary Fig. S3). The results from the gene distribution pattern and dN/dS ratios showed that the NHXs were extensively distributed across the C. sinensis genome. Tandem duplication events were however absent across the tea NHXs. The dN/dS ratios are conclusive proof that strong purifying selection pressures had occurred during the evolution thereby enabling a number of different factors to regulate the NHXs in tea genome.

Functional interaction network of tea NHX proteins
To understand and explore the interaction pattern of NHX genes in tea, a protein interaction network was constructed using the STRING server based on an Arabidopsis association model (Fig. 7). The Arabidopsis model had to be employed due to the absence of tea database in the STRING server. NHXs are largely involved in a variety of biological processes of which response to salt stress is one of the prominent ones [34]. The interaction network was therefore constructed based on 5 tea NHXs (TEA012938.1, TEA012286.1, TEA021179.1, TEA012245.1 and TEA000661.1) involved in response to salt stress according to GO ontology. The Arabidopsis homolog for these 5 tea NHXs, AtNHX2 (AT3G05030) was used as the central node to build the full network. The tea proteins, homologs to the Arabidopsis proteins participating in the network were also added. These homologous proteins were designated as STRING proteins and were selected based on high bit scores in BLAST results [26].
Similarity search program such as BLAST are frequently used to produce accurate statistical estimates that help ensuring protein sequences with signi cant similarity to have similar structures [37]. In addition, proteins sharing higher degree of sequence and structural similarities often tend to have similar functions as well [38]. AtNHX2 (AT3G05030) is involved in active K + uptake at the tonoplast and involved in regulating stomatal closure [39]. AtNHX1 (AT5G27150) encodes a vacuolar sodium/proton antiporter involved in salt tolerance, ion homeostasis, and leaf development. Two of the tea NHXs (TEA012938.1 and TEA012286.1) which are homologous to AtNHX1 are massively involved in response to salt stress. AVP1 (AT1G15690) is involved in regulation of apoplastic pH and auxin transport [39]. CLC-C (AT5G49890) is a chloride channel protein and is involved in the Cl − transmembrane transporter activity.
CLC-F (AT1G55620) is another chloride channel protein, localized in chloroplast and golgi apparatus and is involved in voltage-gated Cl − channel activity [39]. HKT1 (AT4G10310) encodes for a sodium transporter expressed in xylem parenchyma cells and is involved in response to osmotic stress and salt stress. SOS2 (AT5G35410) encodes a member of the CBL-interacting protein kinase family and is a regulatory component controlling plant K + nutrition [39]. TEA006066.1, which is homologous to SOS2, also shows response to nutrient and water deprivation according to GO ontology. SOS3 (AT5G24270) encodes for a calcium sensor that is essential for K + nutrition, K + /Na + selectivity, and salt tolerance. CHX18 (AT5G41610), CHX17 (AT4G23700) and CHX15 (AT2G13620) are all involved in regulation of pH and are members of the putative Na + /H + antiporter family [39].
Tissue speci c gene expression of tea NHXs The tissue speci c expression levels of the 9 tea NHXs in 8 different tissues were retrieved from the TPIA database wherein the levels of expressions were evaluated in transcripts per million (TPM). The database has the expression pro le data of all the C. sinensis genes, which have been experimentally validated [40]. The plant tissues that have been assessed in the study involved apical bud, ower, fruit, young leaf, mature leaf, old leaf, root and stem (Additional File 1: Table S5). All the 9 tea NHXs exhibited varying levels of expression in these 8 different tissues. Few of the genes had high levels of expression while the rest had negligible transcript levels (Fig. 8). TEA012938.1 was expressed the most in apical bud, closely followed by TEA000661.1. This similar pattern was observed when the expression levels were checked in ower, young leaf, root and stem. The highest expression level was recorded by TEA012938.1 in fruit, followed by TEA012286.1 and TEA000661.1. In mature leaf and old leaf, TEA012286.1 was expressed the most, followed by TEA012938.1 and TEA000661.1. These results suggested that 3 out of the 9 tea NHXs (TEA012938.1, TEA012286.1 and TEA000661.1) were signi cantly expressed in all the 8 tissues. The rest of the 6 NHXs were minimally expressed in these 8 tissues with TEA023041.1 and TEA012245.1 being the least.

Expression pro les of tea NHXs under cold and drought stress
In-order to check how the 9 tea NHXs respond to varying levels of cold and drought stress, their expression data was retrieved from the TPIA database. The TPIA database has experimentally veri ed expression data for all the C. sinensis genes under cold [41] and drought stress [42]. The cold acclimated data comprised of 5 stages of expression: 1. 25 ~ 20°C (CK), 2. Fully acclimated at 10°C for 6 h (CA1-6h), 3. 10 ~ 4°C for 7 days (CA1-7d), 4. Cold response at 4 ~ 0°C for 7 days (CA2-7d) and 5. Recovering under 25 ~ 20°C for 7 days (DA-7d) [43], where CK is the control (Additional File 1: Table S6). Expression levels for CA1-6h, showed that 7 out of 9 tea NHXs were upregulated while the rest 2 were downregulated. Out of these 7 upregulated genes, TEA012938.1, TEA000661.1 and TEA012286.1 were upregulated the most. When the cold stress was increased to the next stage (CA1-7d), again 7 genes showed upregulation with TEA012286.1, TEA012938.1 and TEA000661.1 being the highest. TEA006997.1, which was initially upregulated in the rst condition, was downregulated in this present condition. Further increasing the cold stress levels at CA2-7d, expression data revealed 5 NHXs being upregulated. TEA025916.1 was slightly upregulated at CA1-7d but was downregulated at CA2-7d. TEA012286.1 was also upregulated for the previous two levels of cold stress but at CA2-7d it was downregulated. TEA000661.1 was expressed the most followed by TEA12938.1. Expression levels under the recovery phase (DA-7d) showed that only 3 NHXs were upregulated (Fig. 9). Throughout the cold stress conditions, 2 genes namely TEA012938.1 and TEA000661.1 consistently maintained high levels of expression, followed by TEA012286.1 and TEA021179.1. These results indicated the active participation of these 4 tea NHX genes in response to cold stress. The expression levels were further checked under drought stress. The expression data in the TPIA database with respect to 25% polyethylene glycol (PEG) treatment includes four stages: 1. 0h, 2. 24h, 3. 48h and 4. 72 h [42], where 0h was taken as the control (Additional File 1: Table S7). Under the rst drought stress period of 24h, 5 of the 9 tea NHXs were upregulated with TEA012938.1 being expressed the most. The same set of genes were upregulated when the drought stress condition was extended to a period of 48 hours and then for 72 hours (Fig. 10). These 5 genes showed upregulated levels of expression throughout the experimental conditions and thereby suggest their roles in response to drought stress.

Expression pro les of tea NHXs under salt stress
The primary role of the NHX genes is response to salt stress [34]. To understand the potential role of the 9 tea NHXs in response to high levels of salinity, the expression data was analysed. the salt stress data in TPIA database is recorded based on treatment with 200 mM NaCl under 4 stages: 1. 0h, 2. 24h, 3. 48h, and 4. 72h [42] where 0h was taken as the control (Additional File: Table S8). Expression data under the 24 hours salt stress condition revealed 3 genes being upregulated. Among these 3 tea NHXs, TEA012938.1 was expressed the most. A similar pattern was observed when the salt stress conditions were extended for periods of 48 hours and 72 hours (Fig. 11). TEA023041.1 and TEA011468.1 were upregulated to a fair extent while TEA012938.1 maintained very high levels of expression throughout the experimental condition, with increasing transcript levels at each stage. GO ontology data too suggested the involvement of TEA012938.1 in response to salt stress. These results clearly indicate the active role of these tea NHXs in response to prolonged levels of salt stress.

Response of tea NHXs to MeJA treatment
The analysis of the cis-acting elements in the promoter regions of the 9 tea NHXs had revealed the presence of 2 MeJA (methyl-jasmonate) responsive elements (CGTCA-motif and TGACG-motif) (Additional File:1 Table S2). To further understand the effect of MeJA on the 9 tea NHXs, their expression data was retrieved from the TPIA database and analysed. This data is recorded based on the results of exposing the plant parts to aqueous solution of MeJA, under 4 stages: 1. 0h, 2. 12h, 3. 24h and 4. 48h [44] where, 0h was used as the control (Additional File 1: Table S9). 7 out of 9 tea NHXs showed upregulation in expression levels when exposed to the MeJA treatment for a period of 12 hours. Extending the duration of the experiment to 24 hours showed a few minor changes in the genes showing upregulation. TEA012286.1, which was showing upregulation initially now was slightly downregulated. On the other hand, TEA025916.1 was downregulated in the initial phase but showed upregulated levels of expression in the present condition. 7 genes were upregulated in total at 24 hours of MeJA treatment.
Further extending the experiment to 48 hours revealed that 5 genes were upregulated while the rest 4 were downregulated (Fig. 12). TEA012938.1 and TEA000661.1 consistently maintained high levels of expression throughout the 48 hours of the exposure to MeJA. These results suggested that the transcription levels of the tea NHXs might have a close relation to the regulation of MeJA.

Discussion
NHX gene families have already been identi ed and functionally characterized for several plants, including Arabidopsis, Rice, wheat, sweet beet, cotton and other [9,45,46,47,48]. However, the NHX genes in C. sinensis has not been studied yet. In this study, the gene structure, phylogenetic relationship, genomic distribution and expression of NHX genes in C. sinensis were all analysed at the genomic level. A total of 9 NHX genes have been identi ed in C. sinensis based on the Na + /H + exchanger domain (Table 1). However, there are 7 NHX genes in Oryza sativa (Os), Solanum lycopersicum (Sl), Medicago truncatula (Mt), Sorghum bicolor (Sb), Zea mays (Zm). Gossypium hirsutum (Gh), Glycine max (Gm) have 23 and 12 NHX genes and Solanum tuberosum (St) has 5. Gene duplication and loss speci c to different subfamilies of NHX over the course of evolution could explain these differences in the number of NHX genes in plants.
In-silico studies based on subcellular localizations showed that NHXs are grouped into three classes (Vac-, Endo-, and PM-class). In Arabidopsis, both NHX7 and NHX8 are localized in the plasma membrane [49] whereas in tea, TEA006997.1 localized in the plasma membrane, TEA011468.1 is localized in endosome and the others in the vacuole (Table 1). Members in each of the classes from algae to higher plants, showed that the NHX families were fairly similar, indicating that NHXs had conserved functions throughout the evolutionary process [7]. The function of NHX transporters may be in uenced by their subcellular localization. Members of the NHX family, which are found on both the plasma membrane and tonoplast, help to maintain ionic homeostasis by excluding and compartmentalizing excess Na + . Furthermore, endomembrane-bound NHX members have been discovered to be important for cellular cargo tra cking, growth development, and protein processing regulation [15,50]. The exon/intron structural diversity, which plays an important role in the evolution of gene families, brings to the evidence for phylogenetic groupings. In C. sinensis, TEA021179.1 possesses a greater number of introns (18) and exons (19)  The putative amiloride-binding site and membrane-spanning pore in the NHX gene families, which contain the amino acid sequence "FFIYLLPPI" [18], have been found to be highly conserved [6,16,18]. In the presence of the drug amiloride and/or its derivatives, this domain inhibits the cation/H + exchange [52]. In the motif study, the amiloride-binding site is located in the N-terminus of motif 3 and it is found in 6 NHX genes of C. sinensis out of the 9. However, motif 3 is not conserved in TEA023041.1, TEA011468.1 and TE006997.1 (Fig. 2). The C-terminus of NHX proteins was diverse in contrast to the conserved Nterminus. Studies have shown that the deletion of the C-terminal hydrophilic region results in increased Na + /H + transport activity, implying that the C-terminus is important not only for subcellular localization but also for transport activity regulation [8,53]. Additionally, conserved motif analysis also showed that motif 8 and motif 14 conserved in all the NHX genes (Fig. 2). Similar results are observed in Arabidopsis [10] and Oryza sativa [46]. The phylogenetic analysis indicated that, the NHXs in P. trichocarpa [34], S. bicolor [45], and B. vulgaris [9] showed three phylogenetic clusters based on their location in the cell; we found the same results for tea NHX transporters. So according to these ndings, the NHX family genes have remained relatively conserved throughout evolution.
Cis-acting regulatory elements function as key molecular switches in transcriptional regulation of gene activities that control a variety of biological processes such as hormonal response, abiotic stress response and development [54,55]. Hormones including ABA, ethylene, SA and IAA play signi cant roles in plants development and stress response [56,57,58,59]. In this study, cis-acting regulatory elements related to transcription factors were identi ed to be randomly distributed across the promotor region of the 9 tea NHXs (Additional File 1: Table S2). One ABA-responsive element (ABRE) has been discovered in 6 NHXs (TEA012938.1, TEA012286.1, TEA012245.1, TEA000661.1, TEA025916.1, TEA011468.1) of C. sinensis (Additional File 1: Table S2). Whereas in poplar, one or two ABREs were observed [34]. This analysis showed that NHX genes may play a role in the ABA signaling pathway. Furthermore, ARE (anaerobic induction), DRE (drought-responsive cis-acting element), LTR (low-temperature responsive element), MBS (drought response,) and STRE (stress-response) were identi ed as stress responsive regulatory elements in tea. Similarly, in PtNHXs from poplar and SbNHXs from S. bicolor are also found to contain similar elements [45]. The results indicated that the identi ed regulatory elements in this study aid in understanding their roles in various abiotic and biotic stress-related pathways.
Further to understand the distribution pattern of the tea NHXs, the genomic distribution mapping was performed. Tandem duplication events were absent across the tea NHXs (Fig. 5). The duplication of genes increases the functional divergence, which is an essential factor in adaptability under changing environmental conditions [60]. The dN/dS ratio indicates the selection pressure on amino acid substitutions, with a ratio less than 1 indicating purifying selection and a ratio greater than 1 indicating positive selection. Wang et al [61] found that positive selection of a gene during evolution increases its potential and transcription levels under stress conditions in T. aestivum and TaBT1. Whereas in tea, the dN/dS ratios provided conclusive evidence that strong purifying selection pressure existed during evolution, allowing a variety of factors to regulate the genes (Additional File 1: Table S3 [70], G. max [71], O. sativa [72,73] and S. bicolor [45]. The expression data for various tissues and stress conditions showed that the tea NHXs may be involved in developmental processes and abiotic stress responses. Our study revealed that in C. sinensis, the NHX genes express differentially in 8 different tissues where TEA012938.1, TEA012286.1 and TEA000661.1 genes showed the highest level of expression in all the tissues and belonged to Vac-class (Additional File 1: Table S5). The different expression patterns in various tissues (Fig. 8), indicated that the NHX gene family provide opportunities to break the functional constraint from the original gene during the course of evolution.
Based on data from other species, functional annotation and interaction analysis of NHX proteins can help us predict their potential regulatory roles. An interaction network was built using Arabidopsis as the model plant and 5 tea NHX proteins described as responsive to salt stress in GO annotations as well as in expression analysis. The electrochemical gradient of protons across tonoplasts, generated by two vacuolar H + -pumps, H + -APTase and H + -PPase, has been shown to drive the Vac-class NHXs [6, 74,75] In this analysis, all the tea genes considered for building the interaction network, belongs to the Vac-class. By increasing cation accumulation, co-expression of ZxNHX and ZxVP1 genes can improve salt tolerance in transgenic plant species such as Lotus corniculatus [76], Alfalfa [77], and sugar beet [78]. These nding suggested that when plants were exposed to salt stress, Vac-class NHXs might work together to transport Na + across tonoplasts. Calcineurin B-like (CBL) is well known for its ability to interact and modulate CBL interacting protein kinases (CIPK), which then mediate Ca + signal transduction [34,79].
During the salinity response, CBL regulates NHX7 (SOS1) and CIPK mediates the Ca 2+ signalling pathway [80]. A salt-stress elicited Ca 2+ signal activates a protein kinase complex consisting of CBL4 (SOS3) and CIPK24 (SOS2), and the complex then phosphorylates and activates the SOS1 protein to extrude Na + out of the cell in Arabidopsis under salt stress [81]. In transgenic tobacco, overexpression of SOS1 gene increased salt tolerance by maintaining a higher K + /Na + ratio [82]. In the current study, CLBs are hypothesized to interact with TEA006066.1, TEA012938.1, TEA012286.1, TEA21179.1, and TEA012245.1 but not with CIPK (Fig. 7). Similarly, NHX7 (SOS1) interactions with CBLs were predicted in poplar [34] and S. bicolor [50]. However, in the future, yeast two hybrid research will need to con rm these proteins interactions.
Ion transporters are important in many biological processes, including ion uptake and sequestration, energy provision and cell expansion [83]. Previous studies in plants found that Na + /H + antiporters as important members in transporters, mediate the coupled exchange of Na + or K + for H + in all cellular compartments [83,84]. The NHX genes primarily use two proton pumps, the H + -ATP enzyme and H + -PPase, to produce H + electrochemical gradients that transport Na + from the cytoplasm to vacuoles or outside the cell, there by maintain Na + ion stability and avoiding the toxic effect of Na + accumulation in cells [85,86].
Stress response analysis showed that each tea NHX genes were responsive to abiotic stresses of drought, cold and salt. Under PEG treatment, the expression of TEA012938.1, TEA012245.1,TEA00066.1, TEA023041.1, TEA011468.1 reached the highest level at 12h (Fig. 10), and TEA012938.1, TEA023041.1 and TEA011468.1 also responded to salt stress in varying degree, demonstrating these genes may be associated with salt and drought stress. MeJA was found to be linked to salt tolerance in few studies [87,88]. MeJA expression analysis revealed that TEA012938.1 and TEA000661.1, expressed high levels throughout the MeJA treatment (Fig. 12), demonstrating these genes may respond to MeJA hormone regulation. Further in the study, the expression levels of TEA012938.1, TEA023041.1 and TEA011468.1 were signi cantly up-regulated by various concentrations of NaCl over a 48-h period and 72-h period (Fig. 11), and their expression levels under high-salt stress were relatively higher than those under either mild or moderate-salt stress. In Reaumuria trigyna, the expression levels of RtNHX1 in leaves showed an increase and reached a high level at 3 hours, and then reduced after 6 hours when exposed to high salt stress (200 mM NaCl) [89]. A similar expression pattern was found in sweet potato, where IbNHX2 was signi cantly up-regulated at 4 hours after treatment of 200 mM NaCl [90]. Another study [91] found that the transcription level of TaNHX3 in both leaves and roots sharply increased at 24 hours and then gradually decreased after 48 hours over a 96 hours period in different wheat cultivars subjected to salt stress. Moreover TEA012938.1 belonging to the Vac-class NHX, showed the highest level of expression for all the salt stress condition. The study showed that the expression levels of Vac-class NHXs are signi cantly higher than other class genes thereby con rming that Vac-class NHXs might play critical roles in salt tolerance. The study also notices that TEA012938.1 and TEA023041.1 showed signi cant expression levels under all abiotic stress conditions thereby providing a comprehensive understanding of the functions of NHXs in C. sinensis.

Conclusion
Among the numerous transporters in monovalent cation/proton antiporter (CAP1) family, the Na + /H + antiporters (NHXs) are secondary ion transporters to exchange H + and transfer the Na + or K + across membrane. Although much progress has been made in identifying the functions of NHX genes in many organisms, these genes had not been studied in C. sinensis. The objective of this study was to identify, characterize and determine the role of NHX genes in tea at the genomic stage. In total, 9 NHX genes were found and classi ed. Using phylogenetic relationship, the tea NHXs were grouped into three major classes (Vac-, Endo-, PM-classes). Localization within the cell supported the pattern of grouping. The Amiloride-binding site (FFIYLLPPI) was found in the N-terminus of motif 3. The 9 NHXs genes taken under comparative study clearly revealed that tandem duplication events were missing due to their close distribution pattern in tea. These genes were then mapped onto the scaffold by retrieving publicly available information from the TPIA database. The functional interacting network of NHX genes was studied to reveal the interaction patterns of the NHXs. ABRE responsive element was found in 6 genes, implying NHX gene might be involved in ABA signalling pathway. A comparative study of the gene expression pro les in different plant tissues was conducted to demonstrate the possible role of the NHXs in developmental processes. Furthermore, responses of tea NHX to drought, cold, salinity indicated that the genes were involved in single or multiple stress responses. Adding to the tally, the responses to MeJA treatment on the tea NHXs were also analysed. Tea is a commercial crop, grown all over the world and abiotic stresses are one the major factors that limit the crop productivity worldwide. This work will therefore serve as a basis to provide valuable information for future studies and exploration of the role of NHX genes in various developmental process, as well as the elucidation of other potential functions in C. sinensis.

Identi cation of NHX genes of tea plant
The tea plant genome sequence was recovered from the Tea Plant Information Archive, TPIA [43] (http://tpia.teaplant.org/). The NHX genes from Arabidopsis and rice were retrieved from TAIR database [92] (https://www.arabidopsis.org/) and Rice Genome Annotation Project database [93] (http://rice.plantbiology.msu.edu/) respectively. These sequences were then used as a query sequences to scan the tea genome database using the BLASTp algorithm with an e-value of 1e-5 and an identity match of 50% as the threshold. To further con rm the presence of Na+/H+_Exchanger domain, the NHX genes were submitted to SMART [94]  GO ontology annotation and functional interaction network QuickGO (QuickGO(ebi.ac.uk)) was used to perform GO Ontology (GO) analysis for all the 9 tea NHX genes. Furthermore, the network of functionally interacting homologous genes between tea and Arabidopsis was identi ed and constructed using STRING online tool [106] (https://string-db.org/) with default parameters.

Expression pro le of tea NHXs
The tissue speci c expression pro les in 8 plant tissues, which include expression levels in apical bud, ower, fruit, young leaf, mature leaf, old leaf, root, and stem were retrieved from TPIA database and analysed [43]. Furthermore, gene expression data under cold, drought and salt stresses were analysed to understand the potential role of tea NHXs in response to the abiotic stress factors. Additionally, to check the effects of MeJA (methyl-jasmonate) treatment, its expression data was retrieved from TPIA database and analysed for the 9 tea NHXs. Respective graphs for the gene expression for all the tea NHX genes were generated. Heat maps for the same were generated using heatmapper online server [107] (Heatmapper). Phylogenetic tree of NHX genes from Arabidopsis thaliana, Camellia sinensis, Oryza sativa, Solanum lycopersicum, Solanum tuberosum, Medicago truncatula, Populus trichocarpa, Gossypium hirsutum, Sorghum bicolor, Zea mays, Glycine max. The full-length NHX protein sequences were aligned using MUSCLE, and the phylogenetic tree was constructed using MEGA 7.0.14 by the Neighbor-Joining (NJ) method with default parameters and 1000 bootstrap replicates. The tree is divided into three major classes of NHX genes, consisting of the Vac-, Endo-and PM-classes. The motif analysis of NHX genes in Arabidopsis thaliana, C. sinensis and Oryza sativa. The motif gures were generated by MEME suite. A total of 15 motifs were identi ed and are marked individually.   GO analysis of all the 9 NHX genes in C. sinensis. The results have been grouped into three main categories: Biological Process, Cellular Component and Molecular function. The y-axis represents the frequency of genes while the x-axis represents the potential functions.

Figure 7
Functional interaction networks of NHX proteins in C. sinensis. The interaction network was formed according to homologs in Arabidopsis thaliana.