Genome-Wide Identi cation and Expression Analysis of Metal Tolerance Protein (MTP) Gene Family in Medicago Truncatula Under a Broad Range of Heavy Metals Stress

Metal tolerance proteins (MTP) encompass plant membrane divalent cation transporters to specically participate in heavy metal stress resistance and minerals acquisition. However, the molecular behaviors and biological functions of this family in M. truncatula are scarcely known. We identied 12 potential MTP candidate genes and analyzed for a phylogenetic relationship, chromosomal distributions, gene structures, protein structures, gene ontology, and previous RNA-seq data. MtMTPs were classied into three major cation diffusion facilitator (CDFs) groups; Mn-CDFs, Zn-CDFs, and Fe/Zn-CDFs. Structural analysis of SlMTPs displayed high gene similarity within the same group where all of them have cation_eux domain or ZT_dimer. RNA-seq and gene ontology analysis revealed a signicant role of MTP genes during M. truncatula growth and development.MTP genes showed tissue-specic and variable expression levels under the stress of the following ve divalent heavy metals (Cd 2+ , Co 2+ , Mn 2+ , Zn 2+, and Fe 2+ ). Expression levels of Fe 2+ /MtMTP11 and Mn 2+ /MtMTP4 were upregulated, while Mn 2+ /MtMTP5 was downregulated. In conclusion, MtMTP1.1, MtMTP1.2, and MtMTP4 play a key role under heat and heavy metal stress in M. truncatula.


Introduction
M. truncatula is a diploid plant species (2n = 2x = 16) with a comparatively small and its genome has been sequenced, that is being recursively used as a model plant for legume genetic research (Benedito et al. 2008, Tang et al. 2014, Young et al. 2011. It is one of the best essential forage crops widely cultivated across the world, same as alfalfa (Young &Udvardi 2009). M. truncatulais the best species to explore functional genomics of metal tolerance, although it does not accumulate metals by itself (Zhou et al. 2012).
Generally, metals act as co-factor, which has essential implications inactivating enzymes in plant cells to perform the speci c biological reaction. In contrast, the higher accumulation of these metals in the cytosol is toxic for plant cells and resulted in necrosis (Thomine &Vert 2013).
Plant response to heavy metal stress activates a complex signal transduction network. Heavy metal acts as external stimuli to activate in-vivo biosynthesis of stress-related proteins and signaling molecules which subsequently activate transcription of speci c metal-responsive genes to counter metal stress  (Schutzendubel &Polle 2002). In response to counter ROS, plants harbor antioxidants i.e., superoxide dismutase, peroxidase, catalase, glutathione reductase, and low non-enzymatic molecule antioxidants i.e., proline, tocopheroles, carotenoids, glutathione, ascorbic acid (Apel &Hirt 2004).
In nature, many transporters are belonging to different gene families play an essential role in metal regulatory processes which the plant divalent cation transporters or the plant metal tolerance proteins (MTPs) are one of them (Migeon et al. 2010). The cation Diffusion Facilitator (CDF) family genes, known in the plant as MTP, are integral membrane divalent cation transporters that play essential roles in transportation the divalent metals ions between the cell and the vacuole or the extracellular space beside their role in homeostasis and provide tolerance of cells to divalent metal ions, such as Cd + 2 , Zn + 2 , and Co + 2 , Nies &Silver 1995. The plant MTP proteins are homodimers with six transmembrane domains (TMDs), in which hydrophobic C-terminal is embedded in cell membrane while hydrophilic N-terminal faces towards cytosol (Gao et al. 2020, Kolaj-Robin et al. 2015. Based on Arabidopsis MTP sequencing annotation, MTPs are further classi ed into the following seven groups; 1, 5, 6, 7, 8, 10, and 12. A total of 10 MTPs in the Arabidopsis genome and the rst identi ed MTP gene was AtMTP1 ( All retrieved MTP protein sequences were analyzed at E-value < 10 − 5 to identify the MTP domain via SMART (http://smart.embl-heidelberg.de/) tools (Letunic et al. 2004). All detailed genetic information of the putative MTP gene family, such as chromosomal location and CDS, was obtained from the phytozome database (https://phytozome.jgi.doe.gov/). Furthermore, MTP family proteins were analyzed for their molecular weight, a number of atoms, amino acids, isoelectric point, and instability index using EXPASY PROTOPARAM (http://www.expasy.org/tools/protparam.html) (Gasteiger et al. 2003). Finally, theoretical PI and molecular weight were obtained using ProtParam Tool (http://web.expasy.org/portparam).

Chromosomal locations, and synteny analysis
M. truncatula genetic database (https://phytozome.jgi.doe.gov/) was analyzed to retrieve information about chromosomal localization of MTP genes, which were subsequently used to construct genetic map by MapChart software (https://www.wur.nl/en/show/Mapchart.htm). Two genes of the same species located in the same clade were de ned as coparalogs to identify tandem or segmental duplication events. Simultaneously, the Phytozome database (https://phytozome.jgi.doe.gov/) was analyzed to identify segmental duplication in M. truncatula genes. The paralogs were identi ed, which were considered results of tandem duplication due to splicing of two genes into ve or more genes within a 100 kb region . Similarly, coparalogs located within duplicated chromosomal regions were considered segmental duplications (Wei et al. 2007). Smith-Waterman algorithm (http://www.ebi.ac.uk/Tools/psa/) was employed for the calculation of local alignments of two protein sequences. The synteny analysis of the MtMTP gene family was performed using circos (http://circos.ca/) to localize different alleles distributed on different chromosomes (Krzywinski et al. 2009).

Gene structures and motif analyses
Structural analysis of all MTP gene family members was performed to explore the organization of intron/exon by using both gDNA and CDS sequences by deploying an online tool, namely Genes Structure Display Server program (GSDS)(http://gsds.cbi.pku.edu.cn/index.php) (Hu et al. 2015). The conserved gene family motifs were detected using a Multiple EM for motif elicitation (MEME) (http://meme.nbcr.net/meme3/meme.html) by adjusting the following parameters; maximum 20 motifs and 6-200 amino acids per motif (Bailey et al. 2006).

Growth conditions and heavy metal treatments
In this study, M. truncatula (cv. JemalongA17) line was cultivated during the Autom of 2020 at the experimental greenhouse of Yibin University (China). First, the seeds were washed with 10% hypochlorous acid and distilled water. The seeds have been germinated using watersaturated lter paper, and then transferred to fertilized pittmoss soil with germination conditions of 16 hours light (27ºC) and 8.0 dark (18ºC) with a relative humidity of 70%. Four seeds were planted in each plastic pot. After emergence, thinning was performed to maintain two uniform seedlings per pot.Thirty-day-old M. truncatula was placed in 1/2 Hoagland solutions (pH 6.0) with different heavy metal concentrations 0.1 mM CdCl 2 , 0.1 mM CoCl 2 , 0.5 mM FeSO 4 ,1 mM MnSO 4 and 0.5 mM ZnSO 4 respectively, while normal 1/2 Hoagland solutions as the control (CK) (Desoky et al. 2020, Gao et al. 2020). The experimental pots were positioned in a complete randomized block design. The experiment was composed of 6 treatments, as shown above, and each treatment was repeated with 3 pots. Then, 24 h later, the leaves and roots of tube plantlets were collected and used as RNA extraction materials. Three biological replicates of expression analyses have been performed for each treatment.

RNA extraction and qRT-PCR analysis
Trizol reagent (Invitrogen, USA) was used for RNA extraction from all plant samples (leaf, stem, and root), and subsequently reverse transcribed to cDNA using SuperMix Kit (Transgen, Beijing). Primer 5.0 tool was used to design speci c primers of all selected genes, includingβ-actin as a housekeeping gene (Table S1). Total, 20µL reaction mixture was used to perform real-time PCR containing the following reagents; 10 µL 2×SYBR premix Taq, 1 µL cDNA,0.5µL of each primer, and 8 µL ddH 2 O. Real-time PCR reaction conditions were adjusted as; 95°C for 10 min,95°C for 15 sec,60°C for 60 sec, and 40 cycles in total. The relative expression level was calculated by applying Livak Eq. 2 − ΔΔCT with three replications for each sample (Livak &Schmittgen 2001).

Statistical analysis
The data is presented in the form of mean ± standard deviation (SD). After testing the homogeneity of the experimental errors by Bartlett's test (Steel 1997), the data was also subjected to one-way analysis of variance (ANOVA)and Duncan's Multiple Range Test to show signi cant variations among means that were compared at p ≤ 0.05.COSTAT computer software (CoHort Software version 6.303, Berkeley, CA, USA) was used for the statistical analysis.

Identi cation of MTP genes in M. truncatula
In total, 27 genes were identi ed via blast analysis; subsequently, genes with incomplete functional domain were excluded for the next study, and nally, 12 candidate genes were selected for further analysis. Every gene was assigned with a speci c name i.e., MtMTP1.1, MtMTP1.2, MtMTP2, MtMTP4, MtMTP5, MtMTP7, MtMTP8.1, MtMTP8.2, MtMTP9, MtMTP10.1, MtMTP10.2and MtMTP11. Characteristics of all 12 genes such as gene locus, molecular weight, number of amino acids, grand average of hydropathicity, and isoelectric points (Table 1). Except chromosome 6, all other 7 chromosomes of M. truncatula were the locus of MTP genes. The molecular weight of all MTP protein molecules varied from 39491.59 to 53259.81 kDa. The total number of inter and intra protein ionic residues were variable i.e., the highest anionic residues were in MTP1.1, and lowest in MTP5. Similarly, the highest cationic residues were in MTP10.2 and lowest in MTP4 and MTP11. All MtMTP members harbor a variable number of introns, but MtMTP1.1, MtMTP1.2 and MtMTP4 were without any intron.

Phylogenetic analysis of MTP gene families
In order to unreveal the evolutionary footprints of the MTP gene family in M. truncatuala, comparison among MTP gene families in different species was performed. We retrieved 12 AtMTP genes of Arabidopsis thaliana, 9 CsMTP genes of Cucumis sativus, 21 PtMTP genes of Populus trichocarpa, 10 OsMTP genes of Oryza sativa, and 8 TaMTP genes of Triticum aestivum were aligned against 12 MtMTP genes of M. truncatulaand phylogenetic tree was constructed for comparison of evolutionary relationship. All MTPgene families were divided into seven groups i.e., Group1, 5, 6, 7, 8, 10, and 12 (Fig. 1).

Chromosomal locations and synteny analysis of MtMTP gene family
Synteny analyses were performed to unreveal the distribution of genes on different chromosomes. We observed that MtMTP genes are distributed among all seven chromosomes. Furthermore, for evaluation of gene family expansion and novel functions, we also investigated gene duplication and divergence with the help of circos. We observed only segmental gene pair duplication from PGDD (Plant Genome Duplication Database). Collinearity due to excision of segmental duplication was observed in many gene pairs with 70-100% identity percentage (Table S2). Segment duplication resulted in many homologies of MTP genes between M. truncatula chromosome pairs, such as what occurs with the genes, MtMTP1.1/MtMTP1.2,MtMTP4/MtMTP5,and MtMTP8.2/MtMTP10.2 (Fig.2). Except for MtMTP7, all rest of MtMTP genes in M. truncatula displayed single and multiple genetic duplications. Noticeably, we did not observe any obvious tandem duplication among all MtMTPs.

Gene structures and motif analyses
All MTP family genes were further divided into six subfamilies (A, B, C, D, E, and F) (Fig. 3a). Subfamilies A and F were the largest among all subfamilies with six members in each, followed by subfamily B with 2 members, whereas subfamilies C, D, and Eeach contain only one gene (Fig.3a). Intron and exon analysis of all MTP genes revealed that each retrieved sequence of the MtMTP gene family is a correct and true member of six subfamilies (Fig.3c). Although there was variation in size and location of intron and exon of all MtMTP genes, the similarity index was higher among all subfamilies, which proved a close evolutionary relationship among MtMTP gene family members. All MTP family genes contain a variable number of introns, but all members of subfamily F were without any intron. Amino acid sequence-based conserved motifs of MTP were analyzed using MEME (Fig.3b and Table S3). All conserved motifs comprised 50 amino acids except motif 10which was comprised of only 41 amino acids. The largest motifs were 3 and 6, which were observed in all subfamilies, followed by motif 10, 83.3% of aforementioned. Noticeably, the number, type, and order of motifs were similar in intrasubfamily than inter subfamilies.

Protein modeling, sub-cellular localization, and GO enrichment analysis
Phyre 2 web portal (http://www.sbg.bio.ic.ac.uk/phyre2/html/page.cgi?id=index) was employed for protein modeling using all MTP aminoacid sequences (Fig.4, and Table S4). All twelve predicted models for MTP proteins were 100% based on c6xpdB, c3j1zP, c2q B, and d2q a2 templates. Similarly, sub-cellular localization, molecular function, and biological process were predicted by GO enrichment analysis ( Fig.5 and Table S5). In sub-cellular localization analysis, the predicted distribution scores of MTP proteins were as following;12/60% in all membranes, 3/15% in the plasma membrane and vacuole, and 1/5% in Golgi apparatus and root hair. Noticeably, the MtMTP1.2 gene was localized in 12 sub-cellular compartments out of all 14, which underlined the signi cant role of MtMTP1.2in metal stress resistance. Collective scores of MTP protein molecules during biological processes were as following; trans-membrane transport of Zn + and Mn + ions was 3/43%, while trans-membrane transport of cations was 1/14%. More precisely, MtMTP1.1, MtMTP1.2, and MtMTP4 play a key role in transmembrane transport of Zn + , while MtMTP8.1, MtMTP8.2, and MtMTP11 play a crucial role in transmembrane transport of Mn + . Molecular function analysis revealed signi cant roles of MtMTP8.2 and MtMTP11 in heavy metal processes.

Gene expression analysis by RNA-seq data
MtMTP family genes expression pro ling was performed by analyzing previously sequenced RNA-seq data (https://mtgea.noble.org/v3/) of the following tissues of M. truncatula; leaf, bud, shoot, hypocotyl, stem, ower, seed coat, pod, root, and root tip (Table S6). A heatmap diagram was constructed to show the differential expression level of each MtMTP gene in all tissues (Fig. 6). Comparatively, MtMTP1.1 displayed the highest expression level in a pod, while MtMTP1.2 displayed the highest expression level in the root tip. Similarly, the highest expression level of MtMTP2 was observed in plant shoots, while the expression level of MtMTP5 was mild in hypocotyl, root, and root tip. The expression level of MtMTP7 was also mild in buds while higher in the shoot. MtMTP11 displayed the highest expression level in approximately all tissues, while MtMTP10.1only in hypocotyl and seed coat, MtMTP8.1 only in the root, and MtMTP9 only in ower tissue. Noticeably, MtMTP8.2 and MtMTP10.2 displayed the lowest expression level in all tissues but the highest expression in root.

qRT-PCR analysis of MtMTPs under the effect of heavy metals
All MtMTPgenes displayed differential gene expression levels under treatment of different types of heavy metals investigated in the following tissues; root, stem, and leaf (Fig. 7). In roots, MtMTP1. Based on phylogenetic and functional domain analysis, MTP proteins in M. truncatula were further divided into the following six subfamilies (A, B, C, D, E, and F ( Fig. 3a and b). MtMTP gene structure analysis revealed that introns, exons, and motif sequences share great similarities with previously explored MTP gene subfamilies (Liu et al. 2018). For example, all members of the A-subfamily harbor ve introns. At the same time, F-subfamily did not have even a single intron, which shows structural evolutionary changes in the MTP gene family in M. truncatula. The absence of introns in the F-subfamily also revealed a lower selection ability of introns gain or loss rate due to higher selection pressure of exons sequences (Harrow et al. 2006 Annotation and expansion analysis of the MTP gene family in M. truncatula revealed the existence of segmental duplication ( Fig. 2 and Table  S2), similar to other species (Schlueter et al. 2007). In total, we observed 28 segmental duplications in following gene pairs; MtMTP2/MtMTP8.1, MtMTP2/MtMTP10.2 ,and MtMTP8.2/MtMTP9. The existence of genetic duplication con rms that MtMTP gene family members are predominantly involved in secondary metabolism because genes involved in secondary metabolism often go under different types of gene duplication (Ober 2005). Amino acid sequences of each member of the MtMTP gene family were analyzed to predict the 3D structure of proteins, which is responsible for protein function (Fig. 4 and Table S4) (Büyükköroğlu et al. 2018). The c3j1zP template was employed to model MtMTP8.1 and MtMTP8.2 proteins, and all rest of the MtMTP proteins were analyzed with c6xpdB template at 100% con dence level. Noticeably, c6xpdB template speci cally used for transport proteins and in PDB database entitled as; cryo-em structure of human znt8 double mutant-d110n and d224n, 2 which determine outward-facing conformation. Similarly, the c3j1zP template is speci c to metal transport and in the PDB database entitled as; the inward-facing conformation of the zinc transporter yip revealed by cryo-electron microscopy. Besides major templates c6xpdB (for all protein) and c3j1zP (SlMTP2 and SlMTP9), other templates used in protein modeling were c2q B_ (transport protein) and d2q a2 (Cation e ux protein transmembrane domain-like). All four templates predicted the possible role of MtMTP proteins as heavy metal-uptake proteins, which provide a path for transport of essential metals towards cytoplasm and toxic heavy metals export out of the cytoplasm. Similarly, metal-e ux proteins play a key role in detoxifying cells by eliminating excessive toxic heavy metals (Mani &Sankaranarayanan 2018).
Gene ontology (GO) term enrichment analysis revealed vital roles of MtMTPs in heavy metal management such as; trans-membrane transporter activity, cation transmembrane transporter activity, transporter activity, and ion transmembrane transporter activity (Fig. 5 and Table S5). Similarly, RNA-seq is a robust technique that is used for the detection and quanti cation of mRNA pro ling (Zambounis et al. 2020). RNA-seq analysis revealed a signi cant change in the expression pattern of MTP genes in different tissues and at different developmental stages of M. trunactula under various metal ions stress, which was similar to P. trichocarpa, and these ndings endorsed GO enrichment analysis (Gao et al. 2020). Noticeably, the elevated expression level of the following six genes, MtMTP1.2, MtMTP2, MtMTP4, MtMTP7, MtMTP9, and MtMTP11, in plant root tips indicated a key role of them in metal ions uptake from the soil. Moreover, MtMTP10.1 was abundantly expressed in the seed coat, which suggests this gene might have a role in the deposition of metal ions in seed coat, but further studies are invited for authentication ( Fig. 6 and Table S6). Similarly, a higher expression level of MtMTP11 in buds shows its roles in shoots and leaves development similar to OsMTP11 in rice (Tsunemitsu et al., 2018). Regulation of gene expression is a fundamental phenomenon to regulate metabolism associated with duplicate genes (Qian et al. 2010), so down regulation in the expression of MtMTP8.2 and MtMTP10.2 was in order to maintain equilibrium during the metabolic pathway. qRT-PCR was performed to authenticate RNA-seq data analysis (Fig. 7); however, the slight variations between the results of both techniques are due to non-uniform growth conditions and varietal response. We examined the behaviour of these cation family with ve divalent metals where numerous studies in other plants indicated the major roles of them with plant cell to enhance their tolerance against divalent metals such as Mn + 2 , Cd + 2 , Co + 2 and Fe + 2 (Gao et al. 2020). The variable transcriptional pattern of MtMTPs in response to various heavy metals was complicated, which reveal their key role in heavy metal stress management, just like a tonoplast-localized Zn + transporter AtMTP1 in Arabidopsis (Dräger et al. 2004, Kobae et al. 2004, and steady CsMTP1 expression under high concentration of Zn 2+ in cucumber (Migocka et al. 2015).
As mentioned before, upregulation of expression of AtMTP12 was not due to Zn + treatment, but it was due to the formation of a heterodimeric complex with AtMTP5 for Zn + transport (Fujiwara et al. 2015), same in tobacco (Liu et al. (2019). Moreover, variable concentrations of Mn 2+ slightly affect expression of Mn-CDFs (AtMTP8, AtMTP9, AtMTP10, and AtMTP11) (

Conclusion
Genome-wide identi cation analysis revealed 12 MTP genes in M. truncatula, which were further phylogenetically and comprehensively analyzed. The MtMTPs were divided into three major substrate-speci c clusters (Zn/Fe-CDFs, Zn-CDFs, and Mn-CDFs), and six groups that were subjected to polyploidization under segmental duplication. All MtMTPs are predicted to harbor cation_e ux domain and/or ZT_dimerdomain, during each MTP share the same structural characteristics within the same group. The expression patterns of each MtMTPs gene in response to various heavy metals in different tissues indicated that these genes play a vital role in the growth and development of M. truncatula. Furthermore, our gene expression analysis revealed signi cant MtMTP1.1, MtMTP1.2, and MtMTP4 in heavy metal stress tolerance in plants. (-), (+), MW, aa, GRAVY and PI refer to, the total number of negatively charged residues (Asp + Glu), the total number of positively charged residues (Arg + Lys), molecular weight, amino acid number, Grand average of hydropathicity and isoelectric points, respectively    Predicted 3D models of M. truncatula MtMTP proteins. Models have been generated by using the Phyre 2 server in intensive mode. Models were visualized by rainbow colour from N to C terminus.