Genome-Wide Identiﬁcation and Expressional Proﬁling of the Metal Tolerance Protein Gene Family in Brassica napus

: The Cation Diffusion Facilitator ( CDF ) family, also named Metal Tolerance Protein ( MTP ), is one of the gene families involved in heavy metal transport in plants. However, a comprehensive study of MTPs in Brassica napus has not been reported yet. In the present study, we identiﬁed 33 BnMTP genes from the rapeseed genome using bioinformatic analyses. Subsequently, we analyzed the phylogenetic relationship, gene structure, chromosome distribution, conserved domains, and motifs of the BnMTP gene family. The 33 BnMTPs were phylogenetically divided into three major clusters (Zn-CDFs, Fe/Zn-CDFs, and Mn-CDFs) and seven groups (group 1, 5, 6, 7, 8, 9, and 12). The structural characteristics of the BnMTP members were similar in the same group, but different among groups. Evolutionary analysis indicated that the BnMTP gene family mainly expanded through whole-genome duplication (WGD) and segmental duplication events. Moreover, the prediction of cis-acting elements and microRNA target sites suggested that BnMTPs might be involved in plant growth, development, and stress responses. In addition, we found the expression of 24 BnMTPs in rapeseed leaves or roots could respond to heavy metal ion treatments. These results provided an important basis for clarifying the biological functions of BnMTPs , especially in heavy metal detoxiﬁcation, and will be helpful in the phytoremediation of heavy metal pollution in soil.


Introduction
With the development of industrial and agricultural modernization, and the urbanization of human society, the heavy metal contamination of agricultural land and water has become one of the main restrictive factors affecting plant growth and crop yield [1][2][3]. Metal cations, such as Zn, Mn, and Cu, are necessary trace elements that could serve as catalytic or structural cofactors of enzymes and regulatory factors in regulating the various life activities of plants [4][5][6]. However, excessive trace elements are toxic to plants [7,8]. Therefore, their concentration in plant cells must be strictly controlled to achieve optimal growth. On the other hand, some unnecessary heavy metals, such as Hg, Cd, Cr, and Pb, are toxic to plant growth and development even in small doses [9]. In addition, these toxic heavy metals could be easily transported and accumulated in different organs of crops (e.g., root, stem, leaf, and seed), thus causing severe damage to animal and human health [10,11]. Under long-term heavy metal stress, plants have gradually formed two adaptation mechanisms to heavy metal toxicity. The first is to restrict the metal influx across the membrane, increase the metal binding to the cell wall, and stimulate the efflux pumping of metals from cytosol, to maintain the low metal concentration in the cytoplasm [12]. The second is to chelate the absorbed heavy metals with organic molecules and transform them into the future, and provides a new idea for creating rapeseed germplasm with heavy metal resistance or high soil remediation abilities.

Phylogenetic Analysis and Characterization of BnMTPs
The multiple sequence alignment of MTPs was performed using ClustalX (http:// www.clustal.org/clustal2/, accessed on 10 November 2021). The phylogenetic tree was constructed using MEGA v. 6.0 with the neighbor-joining (NJ) method using the p-distance and pairwise deletion option, and with a bootstrap analysis of 1000 replicates.
The number of amino acids, molecular weight (MW), and theoretical isoelectric point (pI) of the BnMTPs were calculated using the protein isoelectric point calculator (http: //isoelectric.org/, accessed on 10 November 2021). The subcellular localization of BnMTPs were predicted with ProtComp v.9.0 in Softberry (http://linux1.softberry.com/, accessed on 10 November 2021). The TMHMM Server v.2.0 was used to predict the transmembrane region of BnMTPs. The exon-intron structures of BnMTPs were extracted from the genome annotation file of B. napus. The conserved motif and domain in BnMTPs were predicted by MEME (http://meme-suite.org/tools/meme, accessed on 10 November 2021) and Pfam, respectively. The gene structure, conserved motifs, and domains were visualized by TBtools [53].

Gene Location and Duplication Analysis of BnMTPs
MG2C (http://mg2c.iask.in/mg2c_v2.0/, accessed on 12 November 2021) was used to generate the location image of MTPs on the chromosomes of the B. napus genome. MCScanX was used to analyze the collinear blocks and gene replication events among B. napus, A. thaliana, B. rapa, and B. oleracea with default parameter settings. The collinear maps among B. napus and other species were visualized by TBtools. Nonsynonymous (Ka) and synonymous (Ks) values, and the Ka/Ks ratios were calculated by ParaAT 2.0 (https://bigd.big.ac.cn/tools/paraat, accessed on 12 November 2021) and KaKs_Calculator 2.0 [54].

Temporospatial Expression Analysis of BnMTPs Based on RNA-seq Data
To examine the tissue expression pattern of the BnMTP genes, three replicates of 15 samples representing the major developmental tissues and organs of B. napus were used for RNA-seq analysis, including leaf, cotyledon, hypocotyl, root, shoot apical meristem (SAM), stem, bud, flower, endosperm, silique, and five seed samples at 3, 4, 5, 6 weeks after flowering (WAF) and mature stage. Using the RNA-seq data, a heat map of BnMTPs in different developmental stages was generated based on the log 10 transformed values of fragments per kilobase of transcript per million fragments mapped (FPKM) values. If FPKM = 0, log 10 FPKM was artificially defined as −3.

Plant Materials, Growth Conditions and Treatments
To acknowledge the expression level of BnMTPs under different abiotic stresses, the seeds of rapeseed were germinated and transferred to 1/2 MS medium containing 15% PEG, 150 mM NaCl, or 150 mM mannitol according to previous studies with several modifications [56,57]. The seedlings were grown in a climate chamber under a photoperiod of 16 h light/8 h dark, at 22 • C for 14 days. For cold treatment, the seedlings were grown at 22 • C for 12 days and 4 • C for 2 days. Three biological replicates of 10 seedlings were pooled for each treatment. In order to examine the expression pattern of BnMTPs under hormone stimuli, we conducted hormone treatments as described by previous studies with several modifications [56,[58][59][60]. Briefly, five-week-old B. napus seedlings cultured in nutrition soil in a greenhouse were sprayed with different hormone solutions, including 100 µM abscisic acid (ABA), 500 µM gibberellin (GA), 50 µM indoleacetic acid (IAA), 100 µM kinetin (KT), and 10 µM strigolactone (SL). Leaf samples were collected at 0 h, 1 h, and 3 h after treatment. Three biological replicates of six leaves from three seedlings were pooled for each treatment. The abovementioned samples were ground into powder and stored at −80 • C for RNA-seq. The heatmap of BnMTP expression under abiotic stresses and hormone treatments was generated based on the log 2 transformed FPKM ratios, using the control group as CK. If the FPKM values were lower than 1, it was considered that there was no significant difference between the CK and test samples, and the log 2 transformed ratio was defined as 0. The BnMTPs were clustered according to hierarchical clustering.
To investigate the response of BnMTPs under different heavy metal treatments and selenium, rapeseed seeds were germinated on wet filter papers for two days and moved into an aperture disk containing vermiculite. The seedlings were watered with Hoagland nutrient solution, and grown with a 16 h light/8 h dark photocycle. The three-week-old plants were irrigated with Hoagland nutrient solution containing different metal ions, including 100 µM CuCl 2 , 100 µM ZnSO 4 , 100 µM MnCl 2 , 180 uM HgCl 2 , 100 µM Pb(NO 3 ) 2 , 400 µM K 2 Cr 2 O 7 and 1 µM Na 2 SeO 3 . After 8 days of treatment, the roots and leaves from five seedlings were pooled separately for RNA isolation and qRT-PCR analysis. The plants irrigated with Hoagland nutrient solution were used as a control.

RNA Isolation and qRT-PCR Analysis
The total RNA was isolated using RNA isolator Total RNA Extraction Reagent (Vazyme, Nanjing, China) according to the manufacturer's instructions. A total of 3 µg of RNA was reverse transcribed using the HiScript III RT SuperMix for qPCR (Vazyme, Nanjing, China) to generate the cDNA. PowerUp SYBR Green Master Mixes (Thermo, Waltham, MA, USA) and a StepOnePlus Real-Time PCR System (Thermo, Waltham, MA, USA) were used to perform qRT-PCR analysis. B. napus Actin-7 (NC_027775.2) was used as an internal control. The reactions were carried out according to the following program: 95 • C for 15 min, 40 cycles followed by 95 • C for 15 s, 60 • C for 1 min. The 2 − Ct method was used to calculate the relative expression level of BnMTPs. Each experiment was technically repeated three times. All of the primers were synthesized by TSINKE Biotech and are listed in Table S1.

Statistical Analysis
Statistical analysis was performed using SPSS 19.0. One-way ANOVA or independentsamples t-test were used to analyze significant differences among multiple samples or between each pair of samples at a 0.05 level, respectively.

Identification and Phylogenetic Analysis of BnMTP Proteins
A total of 33 BnMTPs, 17 BrMTPs, and 17 BoMTPs were identified in B. napus, B. rapa, and B. oleracea, respectively. The BnMTPs were named as BnMTP1.1 to BnMTP12.2 according to the similarity and phylogenetic relationship to AtMTPs (Table S2). For each AtMTP gene, at least two BnMTP orthologs were identified on the B. napus genome except for AtMTP7. A phylogenetic tree of 79 MTPs from A. thaliana, B. rapa, B. oleracea, and B. napus classified the MTPs into three major clusters: Zn-CDFs, Zn/Fe-CDFs, and Mn-CDFs, and the three clusters were further divided into seven groups (i.e., group 1, 5, 6, 7, 8, 9, and 12) (Figure 1). Group 9 was the largest group, containing 11 BnMTPs, but no BnMTP was included in group 7. There were ten BnMTPs in group 7 and six BnMTPs in group 8, and only two BnMTPs were classified in each of the group 5, 6, and 12.

Structural and Chromosomal Localization Analysis of BnMTP Genes
The genome annotation file from the Genoscope database of B. napus was used to analyze the exon-intron structure of the BnMTPs. The closely related BnMTP members showed similar exon-intron structure and gene length, which was consistent with the phylogenetic relationship mentioned above. We found the number of introns in BnMTPs varied dramatically, ranging from 0 to 11, but seven BnMTPs in the Zn-CDF clade had no intron (Figure 2A). Other Zn-CDFs in B. napus contained one to three introns, except for BnMTP5.1 and BnMTP5.2 in group 5, which contained nine introns. The two Zn/Fe-CDFs in B. napus were identified with 11 introns. The intron number of Mn-CDFs in B. napus ranged from two to six, of which members in group 8 had six introns, and BnMTPs in  50 kDa. The pI of most BnMTPs was lower than 7, except for BnMTP6.1, indicating that most BnMTPs may function in an acidic environment. All of the 33 BnMTPs were predicted with localization in the vacuole, and four of them were also located in the cell membrane. We found 24 BnMTPs containing 4-6 TMDs, whereas other BnMTPs harbored 3 (BnMTP2.1, BnMTP11.3, and BnMTP11.4), 7 (BnMTP5.1 and BnMTP10.4), and 16 (BnMTP12.1 and BnMTP12.2) TMDs. However, no TMD was identified in BnMTP6.1 or BnMTP6.2.

Structural and Chromosomal Localization Analysis of BnMTP Genes
The genome annotation file from the Genoscope database of B. napus was used to analyze the exon-intron structure of the BnMTPs. The closely related BnMTP members showed similar exon-intron structure and gene length, which was consistent with the phylogenetic relationship mentioned above. We found the number of introns in BnMTPs varied dramatically, ranging from 0 to 11, but seven BnMTPs in the Zn-CDF clade had no intron ( Figure 2A). Other Zn-CDFs in B. napus contained one to three introns, except for BnMTP5.1 and BnMTP5.2 in group 5, which contained nine introns. The two Zn/Fe-CDFs in B. napus were identified with 11 introns. The intron number of Mn-CDFs in B. napus ranged from two to six, of which members in group 8 had six introns, and BnMTPs in group 9 contained four to six introns (except for BnMTP10.2).  Among the 33 BnMTPs (15 on A-subgenome and 18 on C-subgenome), we found that 25 BnMTPs were unevenly distributed on 13 chromosomes of B. napus, and 8 BnMTPs were not assigned to specific chromosomes ( Figure 3). Chromosome C04 contained the largest number of BnMTPs (seven genes), chromosome A04 contained three BnMTPs, chromosomes A03, A05, C05, and C08 contained two BnMTP genes. While chromosomes A06, A07, A09, C02, C03, C06, and C09 were assigned with only one BnMTP gene, no BnMTP gene was located on chromosomes A01, A02, A08, A10, C01, or C07. In addition, we found no significant correlation between the chromosome length and number of BnMTPs.

Conserved Domain and Motif of BnMTP Proteins
The multiple sequence alignment of AtMTPs and BnMTPs revealed that all three clusters contained a conserved CDF signature (44 amino acids) at the N terminus. Similar to previous studies, two conserved HxxxD (x = any amino acid) residues were identified in the Zn-CDFs of B. napus, except for BnMTP2.1 and BnMTP2.2. One conserved HxxxD residue was contained in Zn/Fe-CDFs, and two DxxxD residues were found in the Mn-

Conserved Domain and Motif of BnMTP Proteins
The multiple sequence alignment of AtMTPs and BnMTPs revealed that all three clusters contained a conserved CDF signature (44 amino acids) at the N terminus. Similar to previous studies, two conserved HxxxD (x = any amino acid) residues were identified in the Zn-CDFs of B. napus, except for BnMTP2.1 and BnMTP2.2. One conserved HxxxD residue was contained in Zn/Fe-CDFs, and two DxxxD residues were found in the Mn-CDFs of B. napus. Additionally, the BnMTPs in group 1 and 12 contained a His-rich region ( Figure S1). We used HMMER to analyze the domains of BnMTPs, and found that all of the members contained the cation_efflux domain. The BnMTPs in group 6, 8, and 9 also contained a dimerization domain of Zinc Transporter (ZT_dimer), except for BnMTP11.3 and BnMTP11.4. In addition, a coiled-coil structure was identified in the BnMTPs of group 6, 8, and 9 (except for BnMTP8.3, BnMTP8.4, BnMTP9.1, BnMTP9.2, BnMTP11.2, BnMTP11.3, and BnMTP11.4) ( Figure 4). The motif analysis revealed 15 conserved motifs in BnMTPs, with a length ranging from 6 to 50 amino acids ( Figure 2B, Table S3). The motif sequences were submitted to Pfam and five conserved motifs might encode functional domains. Motif 1, 3, and 7 encoded Cation_efflux, motif 2 encoded ZT_dimer, while motif 10 might encode vacuolar ATP synthase subunit S1 (ATP-synt_S). Motif 6 was the most conserved motif in the BnMTPs, except for BnMTP10.2. Motif 5 and 13 existed in 28 and 26 BnMTP proteins, respectively. In general, BnMTPs in the same cluster or group were identified with a similar type and distribution of conserved motifs. For instance, most Mn-CDFs contained motif 1, 4, 9, and 12, except for BnMTP10.2, BnMTP11.3, and BnMTP11.4. Motif 2 and 3 were only present in Mn-CDFs and Zn/Fe-CDFs. Motif 10 and 11 were specific to the Mn-CDFs in group 9 and 8, respectively. Nearly all of the Zn-CDFs in B. napus contained motif 7 and 15, except for BnMTP2.1. Motif 8 only existed in Zn-CDFs in group 1, and motif 14 was specific to Zn-CDFs in group 1 and 12. Among the three clusters, Mn-CDFs contained the largest number of motifs (group 8 contained 10 motifs, and group 9 contained 9-10 motifs), except for BnMTP10.2, BnMTP11.3, and BnMTP11.4 that contained 5-6 motifs. In addition, the number, type, structure, and position of the motifs in Mn-CDFs were more similar than in the other two clusters. These unique and different conserved motif structures may lead to functional differences in BnMTP proteins. group 1, and motif 14 was specific to Zn-CDFs in group 1 and 12. Among the three clusters, Mn-CDFs contained the largest number of motifs (group 8 contained 10 motifs, and group 9 contained 9-10 motifs), except for BnMTP10.2, BnMTP11.3, and BnMTP11.4 that contained 5-6 motifs. In addition, the number, type, structure, and position of the motifs in Mn-CDFs were more similar than in the other two clusters. These unique and different conserved motif structures may lead to functional differences in BnMTP proteins.

Synteny Analysis of MTP Genes
The BLASTP and MCScanX were used to identify the homologous genes and the duplication events in the BnMTP gene family. A total of 29 pairs of BnMTP paralogs were found in B. napus ( Figure 5, Table S4). We found that all of the BnMTPs resulted from duplication events, of which 22 BnMTPs were derived from whole-genome duplication (WGD) or segmental duplication, 9 BnMTPs from dispersed duplication, and 2 BnMTPs from proximal duplication (Table S5). These results suggest that gene duplication events, especially WGD and segmental duplication, played essential roles in the generation and evolution of the BnMTP gene family.
To further understand the evolutionary relationship of MTP family members, collinearity analyses were conducted between B. napus and A. thaliana, B. napus and B. rapa, and B. napus and B. oleracea. A total of 22 BnMTPs (66.67%) were identified with collinearity relationship to MTPs in other Brassicaceae species, of which 19, 42, and 38 orthologous gene pairs were identified from B. napus-A. thaliana, B. napus-B. rapa, and B. napus-B. oleracea, respectively ( Figure 6). Moreover, we screened the BnMTPs in the rapeseed pangenome and found 36 BnMTPs in ZS11, Shengli, and Zheyou7, 35 MTPs in Gangan, and 32 genes in No2172, Quinta, Tapidor, and Westar, which were similar to the MTPs in Damor-bzh ( Figure 7). We found five, eight, five, four, one, three, three, and four accession-specific MTPs in Damor-bzh, Shengli, Gangan, No2172, Quinta, Tapidor, Westar, and Zheyou7, respectively. In general, these results show that gene duplication events were the main forces promoting the expansion of the MTP gene family in B. napus, and the MTPs were conserved in different B. napus accessions.

Cis-acting Elements in the Promoter Regions of BnMTPs
The cis-acting element is a non-coding DNA region upstream of the gene co region, which regulates the transcription of adjacent genes through the combinatio some regulatory molecules. Using PlantCARE, a total of 2702 annotated cis-a elements were identified in the promoter region of the BnMTPs, which were classified 10 types, including gene transcription (1781 elements), light responsiveness elements), phytohormone responsive (323 elements), abiotic stress responsiveness elements), biotic stress responsiveness (3 elements), site-binding (20 elements), t expression (10 elements), secondary metabolism (18 elements), circadian contr elements), and cell cycle (4 elements) ( Figure 8, Table S6). We found that transcription elements were the most abundant, including 1358 TATA-boxes and CAAT-boxes. Light responsive elements were also common in the promoter regio the BnMTPs, including 27 types of elements. The number of these common elemen BnMTPs ranged from 6 (BnMTP 2.1 and BnMTP 11.3) to 23 (BnMTP 11.1), of which G box 4, GT1 motif, and TCT motif were the most abundant elements. Moreover, 10 t of phytohormone responsive elements were identified in all BnMTPs, such as T element, TGA-box, and AuxRR-core (auxin responsiveness), TATC-box, GARE-moti P-box (GA responsiveness), TCA-element (salicylic acid responsiveness), CGTCA-m and TGACG-motif (MeJA responsiveness), and ABRE (ABA responsiveness). We f the abscisic acid responsive elements were the most abundant, followed by the M responsive elements. Additionally, abiotic stress elements including ARE and GC-m for anaerobic induction, LTR for low temperature responsiveness, MBS for dro inducibility, and TC-rich repeats for defense/stress responsiveness were also identifi the promoters of BnMTP genes. The cis-element analysis indicated that BnMTPs mig regulated by a variety of stimuli and be involved in plant growth and developm response to various stresses and hormones.

Cis-acting Elements in the Promoter Regions of BnMTPs
The cis-acting element is a non-coding DNA region upstream of the gene coding region, which regulates the transcription of adjacent genes through the combination of some regulatory molecules. Using PlantCARE, a total of 2702 annotated cis-acting elements were identified in the promoter region of the BnMTPs, which were classified into 10 types, including gene transcription (1781 elements), light responsiveness (385 elements), phytohormone responsive (323 elements), abiotic stress responsiveness (151 elements), biotic stress responsiveness (3 elements), site-binding (20 elements), tissue expression (10 elements), secondary metabolism (18 elements), circadian control (7 elements), and cell cycle (4 elements) (Figure 8, Table S6). We found that gene transcription elements were the most abundant, including 1358 TATA-boxes and 423 CAAT-boxes. Light responsive elements were also common in the promoter regions of the BnMTPs, including 27 types of elements. The number of these common elements in BnMTPs ranged from 6 (BnMTP 2.1 and Bn-MTP 11.3) to 23 (BnMTP 11.1), of which G-box, box 4, GT1 motif, and TCT motif were the most abundant elements. Moreover, 10 types of phytohormone responsive elements were identified in all BnMTPs, such as TGA-element, TGA-box, and AuxRR-core (auxin responsiveness), TATC-box, GARE-motif and P-box (GA responsiveness), TCA-element (salicylic acid responsiveness), CGTCA-motif and TGACG-motif (MeJA responsiveness), and ABRE (ABA responsiveness). We found the abscisic acid responsive elements were the most abundant, followed by the MeJA-responsive elements. Additionally, abiotic stress elements including ARE and GC-motif for anaerobic induction, LTR for low tempera-ture responsiveness, MBS for drought inducibility, and TC-rich repeats for defense/stress responsiveness were also identified in the promoters of BnMTP genes. The cis-element analysis indicated that BnMTPs might be regulated by a variety of stimuli and be involved in plant growth and development, response to various stresses and hormones.

Potential miRNA Target Sites in BnMTP Genes
miRNAs are small, non-coding RNAs that regulate gene expression by interfering with mRNA transcription, translation, or epigenetic processes [61]. We analyzed the potential miRNA target sites in the BnMTP genes by psRNATarget (expectation score < 5.0), and found that 12 miRNAs targeting eight BnMTPs with the value of target accessibility-maximum energy to unpair the target site (UPE) of the miRNA/BnMTP varied from 4.093 to 23.236 (Table 1)  3.6. Potential miRNA Target Sites in BnMTP Genes miRNAs are small, non-coding RNAs that regulate gene expression by interfering with mRNA transcription, translation, or epigenetic processes [61]. We analyzed the potential miRNA target sites in the BnMTP genes by psRNATarget (expectation score < 5.0), and found that 12 miRNAs targeting eight BnMTPs with the value of target accessibilitymaximum energy to unpair the target site (UPE) of the miRNA/BnMTP varied from 4.093 to 23.236 (Table 1)

Expression Pattern of BnMTPs under Abiotic Stress and Hormone Treatment
As indicated above, cis-acting elements related to abiotic stress and hormone response were abundant in the promoter of BnMTP genes. We examined the expression pattern of the BnMTPs under different abiotic stresses and hormone treatments using

Expression Pattern of BnMTPs under Abiotic Stress and Hormone Treatment
As indicated above, cis-acting elements related to abiotic stress and hormone response were abundant in the promoter of BnMTP genes. We examined the expression pattern of the BnMTPs under different abiotic stresses and hormone treatments using transcriptome sequencing ( Figure 10,

Expression of BnMTPs under Heavy Metal and Selenium Treatment
To further explore the potential biological function of BnMTP genes in heavy metal transportation, we analyzed the expression pattern under different metal treatments, including microelements (e.g., Zn, Cu, and Mn) and non-essential elements (e.g., Cr, Pb, and Hg). We also included selenium (Se) treatment, since rapeseed is known as an ideal crop with strong selenium enrichment ability ( Figure 11). We found 24 BnMTPs expressed in the roots and/or leaves, and 10 BnMTPs with low expression in the roots and leaves under normal conditions or heavy metal treatments. Under normal conditions, 13 BnMTPs were more highly expressed in the leaves than in the roots, and 8 BnMTPs were more highly expressed in the roots than in the leaves. BnMTP6.2, BnMTP12.1, and BnMTP12.2 were identified with similar expression in roots and leaves. Under metal treatments, the expression of 24 BnMTPs was significantly changed, and each BnMTP gene responded to one or more metal ion treatments. We summarized the fold change of BnMTP genes under heavy metal treatments ( Table 2). In leaf, nine BnMTP genes (BnMTP1.1, BnMTP4.2, BnMTP8.3, BnMTP8.4, BnMTP8.5, BnMTP9.1, BnMTP10.4

Expression of BnMTPs under Heavy Metal and Selenium Treatment
To further explore the potential biological function of BnMTP genes in heavy metal transportation, we analyzed the expression pattern under different metal treatments, including microelements (e.g., Zn, Cu, and Mn) and non-essential elements (e.g., Cr, Pb, and Hg). We also included selenium (Se) treatment, since rapeseed is known as an ideal crop with strong selenium enrichment ability ( Figure 11). We found 24 BnMTPs expressed in the roots and/or leaves, and 10 BnMTPs with low expression in the roots and leaves under normal conditions or heavy metal treatments. Under normal conditions, 13 BnMTPs were more highly expressed in the leaves than in the roots, and 8 BnMTPs were more highly expressed in the roots than in the leaves. BnMTP6.2, BnMTP12.1, and BnMTP12.2 were identified with similar expression in roots and leaves. Under metal treatments, the expression of 24 BnMTPs was significantly changed, and each BnMTP gene responded to one or more metal ion treatments. We summarized the fold change of BnMTP genes under heavy metal treatments ( BnMTP8.3, BnMTP8.4, and BnMTP8.6 expression in roots, but increased BnMTP8.2 in roots and leaves, and BnMTP8.5 expression in roots. Figure 11. The expression level of BnMTPs under different heavy metal treatments. CK represents the control samples under normal growth conditions. Data represent means ± SD. Different letters (a and b) indicate significant differences between root and leaf under normal conditions. Asterisks indicate significant differences between the treatment samples and CK (n = 3, p < 0.05). Figure 11. The expression level of BnMTPs under different heavy metal treatments. CK represents the control samples under normal growth conditions. Data represent means ± SD. Different letters (a and b) indicate significant differences between root and leaf under normal conditions. Asterisks indicate significant differences between the treatment samples and CK (n = 3, p < 0.05). "+" and "-" means: 2 < fold change < 4; "++" and "--" means: 4 < fold change < 8; "+++" and "---" means: fold change > 8.

Discussion
Heavy metals such as Cu and Zn are necessary for normal plant growth, but a high concentration of essential and non-essential metals would lead to growth inhibition and toxic symptoms [12]. Plants have gradually evolved a series of cellular mechanisms to improve their tolerance to heavy metals [62]. Several ion transport families have been identified with functions in response to heavy metals, including the plasma membrane transporters such as ZIP proteins, the HMA family, YSL transporter proteins, the NRAMP family, and IRT proteins, as well as tonoplast-localized transporters such as the CDF family transporter (also called MTP) and vacuolar iron transporter family [63]. In this study, a total of 33 BnMTP genes were identified in the B. napus genome. These BnMTP genes were unevenly distributed on 13 chromosomes, with seven genes located on chromosome C04.
The phylogenetic analysis of MTP proteins among B. napus, B. rapa, B. oleracea, and A. thaliana indicated that the MTPs could be divided into three clusters. We found 14, 2, and 17 BnMTPs included in Zn-CDF, Fe/Zn-CDF, and Mn-CDF, respectively. Since the homologous to AtMTP7 was absent in B. napus and only six groups of MTPs were identified in B. napus, we speculated that a small portion of BnMTPs may have undergone a gene-loss event during evolution. The protein length and intron number of BnMTPs varied significantly among the different groups, ranging from 151 to 769 amino acids and 0 to 12 introns. This might indicate that the BnMTPs has diverse functions. In particular, the protein size and molecular weight of BnMTP12.1 (769 amino acids and 86.67 kDa) and BnMTP12.2 (757 amino acids and 85.32 kDa) were significantly larger than other BnMTPs, indicating that they may have unique functions and specific evolutionary processes. Similar to previous reports [18], the modified CDF signature and cation_efflux domain were identified in all BnMTPs. Moreover, the zinc transporter dimer domain was observed in the BnMTPs of group 6, 8, and 9 (except BnMTP11.3 and BnMTP11.4). The ZT_dimer was reported as the dimerization region of the MTPs [64]; these BnMTPs with ZT_dimer structures might form homodimers or heterodimers to transport metal ions. The coiled-coil structure was identified in most BnMTPs of group 6, 8, and 9. Coiled-coils were involved in various processes, ranging from providing structural stiffness to the transduction of conformational changes [65,66], but whether this structure is related to the functional divergence of BnMTPs is unclear. The histidine-rich loop in MTP was considered to determine metal selectivity [67]. In the present study, typical histidine-rich regions were found in the BnMTPs of group 1 and 12. The length difference of these regions may be related to the transport ability of BnMTPs to specific metal ions. The consensus sequence HxxxD and DxxxD were differently distributed in three clusters of BnMTPs. Previous studies have shown that the different amino acid residues may be related to the functional differentiation and metal specificity of different CDF groups [18]. Interestingly, BnMTP6.1 and BnMTP6.2 did not possess any TMDs, suggesting that they might also play novel roles. Most BnMTPs in the same clade were identified with similar exon/intron structures, and BnMTP proteins in the same clade had similar conserved domains and motifs; this agreed with the phylogenetic tree constructed with the multi-sequence alignment. In general, the structure characteristics of BnMTPs were similar in the same group, but distinct among different groups, indicating they might have conserved yet diverse functions in B. napus.
Compared with the MTPs in three ancestral species, including 12 AtMTPs in Arabidopsis, 17 BrMTPs in B. rapa, and 17 BoMTPs in B. oleracea, we found that gene family expansion occurred in the BnMTP gene family. This might be due to the multiple polyploidization events during the evolution of B. napus [68]. In the present study, 66.7% (22/33) of the BnMTPs were derived from WGD or segmental duplication, which might be the main forces driving the expansion of the BnMTP gene family. Other duplication events were also found in the BnMTP gene family. Dispersed duplication produces two gene copies that are neither adjacent nor collinear, which is common in different plant genomes [69,70]. We identified nine BnMTPs derived from dispersed duplication, which were dispersed on different chromosomes. BnMTP11.3 and BnMTP11.4 were the only gene pair derived from proximal duplication, with a distance close to 6 kb. These results indicated that subgenomic duplication events such as dispersed and proximal duplication also play important roles in the expansion of the BnMTP gene family. Since B. napus is a young polyploid formed about 75 million years ago, the rapid genome expansion accompanied a large number of gene loss and recombination events [68]. The absence of AtMTP6 orthologs might be caused by gene loss. Since the diploid parents (B. rapa and B. oleracea) of B. napus have experienced triploidization events, the number of gene family members in B. napus should be six times of that in A. thaliana. However, we found that only 45.8% (33/72) of BnMTPs were retained, indicating that extensive gene loss happened during B. napus polyploidization. After WGD or polyploidization, positive selection is important in the early stage of duplicate gene retention [71,72]. In this study, the Ka/Ks ratios of all duplicate MTP gene pairs were less than 0.5, indicating that they underwent strong purification selection after polyploidization.
MTP genes have been confirmed as having functions in the transport and tolerance of different heavy metals, including plant trace elements and non-essential elements [73,74], which might play important roles in improving heavy metal resistance or enrichment in plants. In Brassica, B. juncea and B. napus have been reported to have strong capacity in the uptake of trace elements and heavy metals [41,42]. Therefore, it is of great significance to study the functional characteristics of the MTP gene family in B. napus. Cis-acting regulatory elements play essential roles in gene transcription. In the present study, 323 and 154 cis-acting elements involved in hormone and stress response were identified. Three elements associated with stress responsiveness (ARE, MBS, and LTR) and five hormonal responsive elements (TGACG-motif, CGTCA-motif, ABRE, TGA, and TCA) are enriched in the promoters of BnMTP genes (Table S6), suggesting that BnMTPs could be regulated by multiple environmental and hormonal stimuli. To date, except for the transport and tolerance of heavy metals [73,74], there is no functional report of MTP genes in response to abiotic stresses and hormones. To analyze the potential function of BnMTPs in abiotic stress and hormone response, we screened the expression of BnMTPs in rapeseed under different treatments, including abiotic stresses such as cold, PEG, NaCl, mannitol, and hormone treatments such as ABA, GA, IAA, KT, and SL. Based on the transcriptome sequencing data, we found that multiple BnMTPs were regulated by hormones and abiotic stresses. Three BnMTP genes, BnMTP1.4, BnMTP10.1, and BnMTP10.3, were induced by cold, mannitol, PEG, and salt stress, while BnMTP8.1 and BnMTP8.6 were strongly induced by PEG treatment. In addition, several BnMTP genes, including BnMTP9.1, BnMTP10.1, BnMTP10.2, BnMTP10.3, and BnMTP11.1, were upregulated by ABA treatment, which was in accordance with the enrichment of cis-elements of ABA responsive elements in their promoters. Therefore, these genes might be potential targets or regulators in response to abiotic stress or hormone treatment in B. napus. Abiotic stress seriously affects the yield of agricultural crops around the world [75]. The species in Brassica are also sensitive to abiotic stresses; drought and salinity stresses greatly hinder the yield and adaptation of Brassicas across the world [76,77]. Therefore, it is very important to explore the stress response mechanism to improve rapeseed production under adverse environmental conditions. With the publication of the rapeseed genome sequence [68], a few studies focused on candidate genes and regulatory factors were reported to improve the stress tolerance (e.g., drought, salt tolerance) of rapeseed, including BnPtdIns-PLC2, BnLAS, BnTTG2, BnALA, BnLEA3, BnVOC, BnaA6.RGA, and BnMRD107 [58,[78][79][80][81][82][83]. There should be more stress-responsive genes yet to be identified in B. napus. In our study, several BnMTP genes (e.g., BnMTP1.4, BnMTP8.1, BnMTP8.6, BnMTP10.1, and BnMTP10.3) were identified with putative functions in the stress response, and these genes would be valuable in the genetic modification and breeding of rapeseed with stress tolerance. Phytohormones are key regulators of plant growth and development, and can also enhance plant adaptation to various abiotic stresses [84]. ABA is a phytohormone with important roles in plant responses to adverse environmental stimuli, and is involved in regulating physiological processes ranging from stomatal opening to protein storage, as well as increasing plant resistance to many abiotic stresses such as drought, salt, and cold stress [84,85]. In this study, several BnMTPs, such as BnMTP9.1, BnMTP10.1, and BnMTP10.3, were induced by both ABA and abiotic stresses, suggesting that these BnMTP genes might regulate plant response and tolerance to abiotic stresses via the ABA signaling pathway. miRNAs and their target genes have been reported with functions in various physiological and biochemical processes of plants under heavy metal stress, including metal uptake and transport, metal chelation, reactive oxygen species clearance, and hormone signal transduction [86,87]. A total of 12 miRNAs were identified as targeting BnMTP genes, and some of these miRNAs have been reported in response to different stress conditions. For instance, miR390 and miR156 mediated lateral root growth under salt stress, and were involved in cadmium tolerance and accumulation [88][89][90][91]. In addition, bna-miR6031 expression was suppressed under various stress treatments [92]. Thus, the functional study of BnMTP genes' response to stresses would be valuable in the future.
We also identified the specific tissue expression pattern of different BnMTP genes, which might be helpful to understand the potential function of BnMTPs. For example, Bn-MTP8.1/8.6 and BnMTP10.1/10.3/10.4 were highly expressed in flower and bud, whereas BnMTP9.1 and BnMTP11.1/11.5 were highly expressed in cotyledon. These BnMTPs might play roles in floral and cotyledon development. In addition, we found that the homologous BnMTPs in the same group exhibited similar tissue expression patterns, indicating that they might have similar and redundant functions, and their functional diversification might be related to the polyploidization and sequence expansion [70]. During allopolyploidization, gene duplication was usually accompanied by epigenetic-induced gene silencing, to maintain normal plant growth and development [93,94]. In B. napus, BnMTP2.1/2.2, BnMTP3.1, and BnMTP8.2 exhibited no expression in any of the tested tissues, and this might be important to maintain the intracellular metabolic system and adaptation to environmental changes.
To investigate the potential roles of BnMTPs in response to heavy metal stress, we examined the BnMTP expression under treatment of six heavy metal ions and found the expression of several BnMTP genes increased or decreased significantly after Hg and Cr treatment. Hg and Cr are highly toxic to plants, with strong impacts on plant growth and crop yield [95,96]. To date, there is no report of MTP gene response to Hg and Cr toxicity; thus, further functional studies of these BnMTPs are necessary. Previous studies have shown that the expression of some MTPs was relatively constant at the transcription and translation levels, similar to housekeeper genes, which were not affected by the potential metal substrates [7,26,27,97]. For instance, AtMTP1 was involved in transferring excess zinc from the cytoplasm to vacuoles to maintain zinc homeostasis, but its expression level was relatively stable and was not affected by the Zn concentration [7]. In addition, the response of some MTP genes (e.g., CsMTP1 in cucumber) to metal substrates might occur at the post-transcriptional level rather than the transcriptional level [97]. Recent studies have found that MTP12 could form a functional complex with MTP5 to transport Zn to the Golgi apparatus, but the expression level was not affected by an excess or deficiency of zinc [26,27]. Similarly, we found the gene expression level of BnMTPs in Zn-CDFs and Mn-CDFs did not change significantly (fold change > 2) in the presence of excess metals, except for BnMTP8.6 that was downregulated by excess Mn. Therefore, future studies focusing on the expressional changes of BnMTP proteins under excessive heavy metal treatment are necessary to fully elucidate the function of BnMTP genes and potential BnMTP complexes.

Conclusions
In this study, 33 MTP members in B. napus were identified and divided into three main clusters and seven groups. The BnMTP gene family experienced gene expansion and loss during polyploidization, and the homologous gene of AtMTP7 was lost in the evolutionary history. The temporospatial expression pattern and response to different heavy metal stresses suggest that BnMTPs played important roles in the growth, development, and stress response of B. napus, especially in heavy metal transport, detoxification, tolerance, and enrichment. In particular, the expression levels of several BnMTPs were significantly increased or decreased after Hg or Cr treatment, indicating that these BnMTPs may be involved in plant response to Hg or Cr. In general, these results will be helpful for the functional study of BnMTPs in the future.
Supplementary Materials: The following supporting information can be downloaded at: https://www. mdpi.com/article/10.3390/genes13050761/s1. Figure S1. Multiple sequence alignment of BnMTP and AtMTP proteins, Table S1. The primers used in this study, Table S2. The information of BnMTP genes in B. napus , Table S3. The sequence and Pfam annotation of conserved motifs in BnMTP proteins, Table S4 Table S5. Duplication type of BnMTP genes, Table S6. Cis-acting regulatory elements identified in the promoter regions of BnMTP genes, Table S7. The RNA-seq data (FPKM) of BnMTP genes in different tissues of B. napus, Table S8. The expression pattern of BnMTP genes under abiotic stresses and hormone treatments.

Conflicts of Interest:
The authors declare no conflict of interest.