Transcriptome Response of Metallicolous and a Non-Metallicolous Ecotypes of Noccaea goesingensis to Nickel Excess

Root transcriptomic profile was comparatively studied in a serpentine (TM) and a non-metallicolous (NTM) population of Noccaea goesingensis in order to investigate possible features of Ni hyperaccumulation. Both populations were characterised by contrasting Ni tolerance and accumulation capacity. The growth of the TM population was unaffected by metal excess, while the shoot biomass production in the NTM population was significantly lower in the presence of Ni in the culture medium. Nickel concentration was nearly six- and two-fold higher in the shoots than in the roots of the TM and NTM population, respectively. The comparison of root transcriptomes using the RNA-seq method indicated distinct responses to Ni treatment between tested ecotypes. Among differentially expressed genes, the expression of IRT1 and IRT2, encoding metal transporters, was upregulated in the TM population and downregulated/unchanged in the NTM ecotype. Furthermore, differences were observed among ethylene metabolism and response related genes. In the TM population, the expression of genes including ACS7, ACO5, ERF104 and ERF105 was upregulated, while in the NTM population, expression of these genes remained unchanged, thus suggesting a possible regulatory role of this hormone in Ni hyperaccumulation. The present results could serve as a starting point for further studies concerning the plant mechanisms responsible for Ni tolerance and accumulation.


Introduction
Plants able to grow and reproduce on potentially toxic metal-rich soils without symptoms of toxicity are termed as "metallophytes" [1]. Tolerance to metal toxicity has evolved independently several times in members of distinct families [2]. Metallophytes have developed a range of cellular growing on a non-metalliferous soil, which may be expressed as a higher Ni accumulation capacity and also as more efficient Ni distribution (mainly into above-ground organs). This adaptation, as the effect of long-term, permanent plant exposure to Ni, is probably caused by the qualitative/quantitative changes in pathways involved in Ni uptake and compartmentalization. To determine the mechanisms responsible for Ni hyperaccumulation, we compared the Ni accumulation capacities of two N. goesingensis populations, one from a serpentine soil (TM) and one from a non-metalliferous soil (NTM), and we employed the RNA-seq approach to analyse the global transcriptome response in the roots of plants growing in pot cultures supplemented with Ni.

Plant Growth Response to Nickel Supplementation and Ni Accumulation
Plants from a non-metallicolous population (NTM) and from a serpentinite area (TM) were grown in control media and under Ni exposure. The fresh weight of roots of both tested populations remained unaffected by Ni in comparison to the respective untreated control. The fresh weight of shoots of the NTM population was significantly lower when grown in Ni in comparison with the untreated control, in contrast to the serpentine TM population, where there were no differences in the fresh weight of shoots in relation to control ( Figure 1).
Plants 2020, 9, x FOR PEER REVIEW 3 of 13 capacity and also as more efficient Ni distribution (mainly into above-ground organs). This adaptation, as the effect of long-term, permanent plant exposure to Ni, is probably caused by the qualitative/quantitative changes in pathways involved in Ni uptake and compartmentalization. To determine the mechanisms responsible for Ni hyperaccumulation, we compared the Ni accumulation capacities of two N. goesingensis populations, one from a serpentine soil (TM) and one from a nonmetalliferous soil (NTM), and we employed the RNA-seq approach to analyse the global transcriptome response in the roots of plants growing in pot cultures supplemented with Ni.

Plant Growth Response to Nickel Supplementation and Ni Accumulation
Plants from a non-metallicolous population (NTM) and from a serpentinite area (TM) were grown in control media and under Ni exposure. The fresh weight of roots of both tested populations remained unaffected by Ni in comparison to the respective untreated control. The fresh weight of shoots of the NTM population was significantly lower when grown in Ni in comparison with the untreated control, in contrast to the serpentine TM population, where there were no differences in the fresh weight of shoots in relation to control (Figure 1).
. Figure 1. Fresh biomass of roots, shoots and whole plants of both N. goesingensis ecotypes in control and Ni treatment conditions. The asterisk indicates the statistical significance of the t-test at p ≤ 0.05, error bars represent ± SD (N = 15).
The Ni concentration in roots of the NTM ecotype was significantly higher in comparison to the TM ecotype. On the other hand, the Ni concentration in shoots was higher in the TM ecotype ( Figure  2). The Ni concentration in roots of the NTM ecotype was significantly higher in comparison to the TM ecotype. On the other hand, the Ni concentration in shoots was higher in the TM ecotype ( Figure 2).

Figure 2.
Nickel accumulation in shoots and roots of two Noccaea goesingensis ecotypes upon 21 days of exposure to 150 µM Ni. Error bars represent ± SD, asterisk indicates the statistical significance of the t-test at p ≤ 0.05 (N = 3).

RNA-seq Analysis of N. goesingensis Root Transcriptomes
The exposure to Ni altered the expression of 403 and 534 genes in the roots of the NTM and TM ecotypes, respectively. Only 19 genes were commonly regulated in the roots of both ecotypes ( Figure  3). Among them, we focused our attention on differentially expressed genes belonging to four categories possibly playing a key role in heavy metal stress responses: transport, stress-related, hormones and transcription factors (Figure 4a). Among all categories, a contrasting response to Ni treatment was observed between both ecotypes.

RNA-seq Analysis of N. goesingensis Root Transcriptomes
The exposure to Ni altered the expression of 403 and 534 genes in the roots of the NTM and TM ecotypes, respectively. Only 19 genes were commonly regulated in the roots of both ecotypes ( Figure 3). Among them, we focused our attention on differentially expressed genes belonging to four categories possibly playing a key role in heavy metal stress responses: transport, stress-related, hormones and transcription factors (Figure 4a). Among all categories, a contrasting response to Ni treatment was observed between both ecotypes.

Down-regulated Genes
In the TM ecotype, only one gene was classified into the metal transport category-COPT4 (COPPER ION TRANSMEMBRANE TRANSPORTER, At2g37925), which was down-regulated ( Figure 4b). In turn, in the NTM ecotype, two genes associated with metal transport were downregulated, i.e., IRT1 and F1K23.4 (HEAVY METAL TRANSPORT/DETOXIFICATION SUPERFAMILY PROTEIN, At1g29000) ( Figure 4b). Within the abiotic stress category, three genes were down-regulated in the TM ecotype, including T3K9.23 (HEAT SHOCK PROTEIN BINDING, At2g41000), RD29B (LOW-TEMPERATURE-INDUCED 65, At5g52300) and ATDI21 (ARABIDOPSIS THALIANA DROUGHT-INDUCED 21, At4g15910) ( Figure 4b). None of the genes from this category were down-regulated in the NTM ecotype. Two genes, At3g12900 and At2g25450, associated with response to ethylene, were down-regulated in the NTM ecotype ( Figure 4b). In turn, no negative gene regulation in this category was shown in the TM ecotype. Ten genes were down-regulated from the transcription factors category in TM ecotype, two of which belonged to the ethylene-responsive group: At1g06170 and CRF10 (CYTOKININ RESPONSE FACTOR 10, At1g68550) ( Figure 4b). In the NTM ecotype, the expression of eight genes was down-regulated in this category (Figure 4a). Among these genes, only two are associated with the response to ethylene: CRF7 (CYTOKININ RESPONSE FACTOR 7, At1g22985) and INTEGRASE-TYPE DNA-BINDING SUPERFAMILY PROTEIN (At3g60490) (Figure 4b).

Discussion
The comparison of the transcriptome of two N. goesingensis ecotypes, that differed in Ni accumulation capacity and tolerance, enabled us to provide an insight into mechanisms underlying Ni hyperaccumulation. Noccaea goesingensis is considered as a Ni hyperaccumulator and has frequently been studied in terms of Ni tolerance and accumulation [17][18][19][20]22]; however, data concerning the mechanisms involved in Ni hyperaccumulation in this plant species are still scarce. According to our results, plants from the TM population accumulated more Ni in shoots and simultaneously less in roots in comparison to the NTM ecotype, suggesting variation in metal uptake and distribution mechanisms. According to Himmelbauer et al. [22], the root system of N. goesingensis exhibits the potential for enhanced Ni uptake from the soil. Furthermore, it is considered that the constitutive overexpression of genes encoding transmembrane transporters in roots, such as members of ZIP and HMA, plays a determinant role in driving the translocation of metals to shoots and is a one of hallmarks distinguishing hyperaccumulators from non-hyperaccumulator plants (reviewed by Rascio et al. [37]). To determine the differences in mechanisms of Ni uptake and translocation between ecotypes, the analysis of the global transcriptome response was performed only in roots.
The two N. goesingensis populations responded to Ni by modulating the expression of a different gene set, as shown by the small overlap in the number of genes that were regulated by the metal with the same trend of expression (up-or down-regulated), indicating a more qualitative than quantitative variation in plant response ( Figure 5).

Up-Regulated Genes
Among the metal transport category, four genes were upregulated in the TM ecotype, including IRT1 (IRON-REGULATED TRANSPORTER 1, At419690), IRT2 (IRON ION TRANSMEMBRANE

Down-Regulated Genes
In the TM ecotype, only one gene was classified into the metal transport category-COPT4 (COPPER ION TRANSMEMBRANE TRANSPORTER, At2g37925), which was down-regulated ( Figure 4b).
In turn, in the NTM ecotype, two genes associated with metal transport were down-regulated, i.e., IRT1 and F1K23.4 (HEAVY METAL TRANSPORT/DETOXIFICATION SUPERFAMILY PROTEIN, At1g29000) (Figure 4b). Within the abiotic stress category, three genes were down-regulated in the TM ecotype, including T3K9.23 (HEAT SHOCK PROTEIN BINDING, At2g41000), RD29B (LOW-TEMPERATURE-INDUCED 65, At5g52300) and ATDI21 (ARABIDOPSIS THALIANA DROUGHT-INDUCED 21, At4g15910) (Figure 4b). None of the genes from this category were down-regulated in the NTM ecotype. Two genes, At3g12900 and At2g25450, associated with response to ethylene, were down-regulated in the NTM ecotype (Figure 4b). In turn, no negative gene regulation in this category was shown in the TM ecotype. Ten genes were down-regulated from the transcription factors category in TM ecotype, two of which belonged to the ethylene-responsive group: At1g06170 and CRF10 (CYTOKININ RESPONSE FACTOR 10, At1g68550) (Figure 4b). In the NTM ecotype, the expression of eight genes was down-regulated in this category (Figure 4a). Among these genes, only two are associated with the response to ethylene: CRF7 (CYTOKININ RESPONSE FACTOR 7, At1g22985) and INTEGRASE-TYPE DNA-BINDING SUPERFAMILY PROTEIN (At3g60490) (Figure 4b).

Discussion
The comparison of the transcriptome of two N. goesingensis ecotypes, that differed in Ni accumulation capacity and tolerance, enabled us to provide an insight into mechanisms underlying Ni hyperaccumulation. Noccaea goesingensis is considered as a Ni hyperaccumulator and has frequently been studied in terms of Ni tolerance and accumulation [17][18][19][20]22]; however, data concerning the mechanisms involved in Ni hyperaccumulation in this plant species are still scarce. According to our results, plants from the TM population accumulated more Ni in shoots and simultaneously less in roots in comparison to the NTM ecotype, suggesting variation in metal uptake and distribution mechanisms. According to Himmelbauer et al. [22], the root system of N. goesingensis exhibits the potential for enhanced Ni uptake from the soil. Furthermore, it is considered that the constitutive overexpression of genes encoding transmembrane transporters in roots, such as members of ZIP and HMA, plays a determinant role in driving the translocation of metals to shoots and is a one of hallmarks distinguishing hyperaccumulators from non-hyperaccumulator plants (reviewed by Rascio et al. [37]). To determine the differences in mechanisms of Ni uptake and translocation between ecotypes, the analysis of the global transcriptome response was performed only in roots.
The two N. goesingensis populations responded to Ni by modulating the expression of a different gene set, as shown by the small overlap in the number of genes that were regulated by the metal with the same trend of expression (up-or down-regulated), indicating a more qualitative than quantitative variation in plant response ( Figure 5). Among the differentially regulated genes, genes classified into abiotic stress-related category were mainly up-regulated in the roots of the NTM ecotype, whereas in the TM ecotype, genes classified into the same group were down-regulated. This indicates that, although the root biomass of both ecotypes was unaffected by Ni treatment, the NTM ecotype was more prone to metal toxicity, activating mechanisms associated with the stress response. Differences in gene expression were also observed among genes encoding metal transporters, genes associated with ethylene metabolism and genes encoding ethylene-responsive transcription factors. The expression pattern of metal transporter encoding genes suggests that IRT1, IRT2, ZIP10, and MTP8 may contribute to enhanced Ni accumulation in the TM ecotype. Similar results were described for N. caerulescens metallicolous accession [24], where IRT1, IRT2 and ZIP10 encoding genes were suggested to be associated with the hyperaccumulation of several metals, including Ni. The analogous group of genes, including IRT1 and IRT2, was also up-regulated in Arabidopsis halleri Zn hyperaccumulating population [25]. IRT1 as well as IRT2 have been shown to possess low substrate specificity, transporting iron but also Zn, Cd and Ni [38][39][40][41]. The low specificity of metal transport applies also to ZIP10, which has been shown to transport Fe and Zn [42]. Even slight variations in the transporter sequence may influence the substrate specificity [38], and it remains to be investigated if differentially regulated transporters (including IRT1, IRT2, ZIP10, MTP8) transport Ni more preferably over other metals in the TM N. goesingensis ecotype.
The concentration of Ni in shoots of the TM ecotype was higher in comparison to the NTM ecotype, which may indicate changes in the mechanisms associated with Ni distribution in TM ecotype. Up to date, there is a lack of evidence concerning transporters involved in the Ni xylem loading process; however, existing studies suggested the involvement of Zn or/and Fe transport system [43,44]. Among upregulated genes encoding metal transporters, the IRT1 may be considered in the context of Ni distribution, since it was found to translocate other metals, such as Mn in Hordeum vulgare [45]. The IRT1 encoding gene was found to be expressed mainly in root epidermal cells in A. thaliana, where it participates in metal uptake from the soil [46]. The IRT1 localization in roots of plants from the Noccaea genus has not been shown to date, and therefore we cannot exclude a possible role of this transporter in Ni uptake and distribution, which needs further study.
Toxic metal uptake and distribution within plant organisms are under strict control of regulatory pathways mainly controlled by phytohormones, with particular emphasis on ethylene (reviewed by Khan et al., 2017). According to our results, the genes associated with ethylene biosynthesis and response were up-regulated only in the TM ecotype. Additionally, with a view to the abovementioned upregulation of genes involved in Fe transport, we can hypothesize that ethylene may serve here as a main regulating molecule implicated in the N. goesingensis adaptation to Ni excess. We observed the up-regulation of genes involved in both ethylene biosynthesis and ethylene response. However, this consideration is based solely on RNA-seq analysis-the ethylene levels were not determined in the tested ecotypes, and thus this hypothesis requires further verification. There Among the differentially regulated genes, genes classified into abiotic stress-related category were mainly up-regulated in the roots of the NTM ecotype, whereas in the TM ecotype, genes classified into the same group were down-regulated. This indicates that, although the root biomass of both ecotypes was unaffected by Ni treatment, the NTM ecotype was more prone to metal toxicity, activating mechanisms associated with the stress response. Differences in gene expression were also observed among genes encoding metal transporters, genes associated with ethylene metabolism and genes encoding ethylene-responsive transcription factors. The expression pattern of metal transporter encoding genes suggests that IRT1, IRT2, ZIP10, and MTP8 may contribute to enhanced Ni accumulation in the TM ecotype. Similar results were described for N. caerulescens metallicolous accession [24], where IRT1, IRT2 and ZIP10 encoding genes were suggested to be associated with the hyperaccumulation of several metals, including Ni. The analogous group of genes, including IRT1 and IRT2, was also up-regulated in Arabidopsis halleri Zn hyperaccumulating population [25]. IRT1 as well as IRT2 have been shown to possess low substrate specificity, transporting iron but also Zn, Cd and Ni [38][39][40][41]. The low specificity of metal transport applies also to ZIP10, which has been shown to transport Fe and Zn [42]. Even slight variations in the transporter sequence may influence the substrate specificity [38], and it remains to be investigated if differentially regulated transporters (including IRT1, IRT2, ZIP10, MTP8) transport Ni more preferably over other metals in the TM N. goesingensis ecotype.
The concentration of Ni in shoots of the TM ecotype was higher in comparison to the NTM ecotype, which may indicate changes in the mechanisms associated with Ni distribution in TM ecotype. Up to date, there is a lack of evidence concerning transporters involved in the Ni xylem loading process; however, existing studies suggested the involvement of Zn or/and Fe transport system [43,44]. Among upregulated genes encoding metal transporters, the IRT1 may be considered in the context of Ni distribution, since it was found to translocate other metals, such as Mn in Hordeum vulgare [45]. The IRT1 encoding gene was found to be expressed mainly in root epidermal cells in A. thaliana, where it participates in metal uptake from the soil [46]. The IRT1 localization in roots of plants from the Noccaea genus has not been shown to date, and therefore we cannot exclude a possible role of this transporter in Ni uptake and distribution, which needs further study.
Toxic metal uptake and distribution within plant organisms are under strict control of regulatory pathways mainly controlled by phytohormones, with particular emphasis on ethylene (reviewed by Khan et al., 2017). According to our results, the genes associated with ethylene biosynthesis and response were up-regulated only in the TM ecotype. Additionally, with a view to the above-mentioned upregulation of genes involved in Fe transport, we can hypothesize that ethylene may serve here as a main regulating molecule implicated in the N. goesingensis adaptation to Ni excess. We observed the up-regulation of genes involved in both ethylene biosynthesis and ethylene response. However, this consideration is based solely on RNA-seq analysis-the ethylene levels were not determined in the tested ecotypes, and thus this hypothesis requires further verification. There are studies that revealed increased levels of this phytohormone in response to metal presence, including Cd, Cu [47], Ni and Zn [29]. Interestingly, in this latter study, the ethylene-regulated antioxidant metabolism was shown to maintain a higher level of reduced glutathione, alleviating photosynthetic inhibition in metal-exposed mustard plants. Even though recent advances are constantly showing that ethylene plays a pivotal role in heavy metal tolerance (see, for example, reference [48]), more future studies are required to clarify the involvement of ethylene in toxic metal accumulation leading to enhanced tolerance.
Considering that this comparative study was mainly performed at the transcriptomic level, it is worth paying attention to changes observed in transcription factors playing a major role in changes in gene regulatory networks at transcriptional level. The TM population showed a higher number of upregulated transcription factors than NTM. Among them, four (WRKY70, MYB31, TEM1, CRF1) might be key players in the modulation of a signalling cascade of downstream response genes involved in heavy metal hyperaccumulation. WRKY and MYB transcriptional factors have been indicated as master regulators of several stress-related genes [49]. WRKY70, already shown to be upregulated by metal nanoparticles in Arabidopsis [50], is known to be involved in brassinosteroid-regulated plant growth [33], and such hormones can be in turn involved in the heavy metal response [51]. MYB31 is thought to act as a repressor of lignin biosynthesis [52] and was found to be responsive to auxin [53]; nonetheless, to our knowledge, it has never been linked to plant metal response. Together with the abovementioned induction of TEM1, involved in the ethylene-mediated response, the up-regulation of these transcription factors could suggest the involvement of a different hormonal scenario in the TM population with respect to the NTM one. In addition to such transcription factors, our results could suggest CRF1 as another interesting candidate, worthy of further investigated in the context of Ni hyperaccumulation, as it is known to be involved in heavy metal response in eukaryotic organisms, mediating copper and cadmium resistance in the yeast Yarrowia lipolytica [54].
In conclusion, the performed studies revealed clear differences in the response to Ni in TM and NTM N. goesingensis ecotypes at the level of accumulation capacity and in the whole transcriptome response in roots. The study allowed us to indicate increased expression of major genes encoding metal transporters, IRT1 and IRT2, and involved in ethylene biosynthesis and signalling, ACS7, ACO5, ERF104, and ERF105, probably associated with Ni hyperaccumulation in contrasting populations of N. goesingensis.

Plant Nickel Concentration
Approximately 50 mg of root and shoot samples (N = 3) were weighted to analytical accuracy, transferred into teflon autoclave (Speedwave ® Entry, Berghof, Germany) and pre-digested in 5.00 mL of 65% HNO 3 (Argenta, Poland) at room temperature for 1 h. The digestion was carried out for 30 min after the addition of 30% hydrogen peroxide (Sigma Aldrich, Saint Louis, MO, USA). Microwave digestion was carried out for 35 min (temp profile: step 1-ramp 5 • C/min, time-5 min, temp: 145 • C; step 2-ramp 3 • C/min, time-10 min, temp: 190 • C; step 3-ramp 10 • C/min, time-1 min, temp: 75 • C). Subsequently, the solution was cooled to room temperature and quantitatively transferred into a 25 mL volumetric flask and made up with deionized water. The blank samples were processed simultaneously according to the same analytical procedure. Nickel was measured by using atomic absorption spectrometry (graphite furnace [GF-AAS], equipped with an auto-sampler [Thermo Scientific, iC3000, US]). The external standard calibration method was applied using AAS standard solutions (Sigma Aldrich, USA). All chemicals were trace metal grade.

RNA Sequencing and Functional Data Mining
Total RNA Mini Kit (Bio-Rad, Hercules, CA, USA) was used for total RNA extraction from roots (from 5 plants per sample, N = 3). RNA digestion was performed with DNase (DNA free kit, Ambion Bioscience, Austin, TX USA) and RNA purity and quantity were analyzed by Biospec-Nano (SHIMADZU, Japan). The integrity of RNA was assessed with the Agilent 2100 Bioanalyzer (Agilent, Santa Clara, CA, USA) and RNA 6000 Nano Kit (Agilent, Santa Clara, CA, USA). Ion Total RNA-seq Kit v2 (ThermoFisher Scientific, Waltham, MA, USA) was utilized for whole transcriptome libraries preparation. Poly(A) RNA selection was performed with 1000 ng of total RNA with the use of Dynabeads mRNA DIRECT Micro Kit (ThermoFisher Scientific, Waltham, MA, USA) following the manufacturer's protocol. Next, the poly(A) RNA was fragmented with RNase III, purified and checked for the quality with Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA) and the RNA 6000 Pico Kit (Agilent Technologies, Santa Clara, CA, USA). In the next step, the RNA was hybridized and ligated. Reverse transcription was carried out using the Ion Total RNA-seq Kit v2 and the cDNA was purified, amplified and barcoded with the Ion Xpress RNA-Seq Barcode 1-16 Kit (ThermoFisher Scientific, Waltham, MA, USA) and purified. The High Sensitivity DNA Kit with the use of the 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA) was utilized for yield and size distribution of amplified DNA determination. All of the prepared libraries were diluted to equimolar concentrations (100 pM) and merged into sets of 6 samples. Template-positive ion (PI) ion sphere particles (ISPs) with 200 base-pair average insert libraries for sequencing were prepared and enriched using Ion PI Hi-Q OT2 200 Kit (ThermoFisher Scientific, Waltham, MA, USA) and Ion OneTouch 2 System (ThermoFisher Scientific, USA) according to the manufacturer's protocol. Sequencing was conducted with the use of Ion PI Hi-Q Sequencing 200 Kit (ThermoFisher Scientific, Waltham, MA, USA) and Ion PI Chip Kit v3 (ThermoFisher Scientific, Waltham, MA, USA) on the Ion Proton System (ThermoFisher Scientific, Waltham, MA, USA) according to the standard protocol provided by the manufacturer. The raw data were processed by trimming low quality bases and the adaptor sequence was removed to obtain high-quality clean reads using cutadapt version 1.8.1. The pre-processed high-quality reads were mapped to A. thaliana genome with HISAT2 version 2.1.0 [55] using the default parameters. The alignment statistics for the same are given in Table S1. The resulting output of HISAT2 was employed to determine the differentially expressed genes using Cuffdiff tool in Cufflinks version 2.2.1 pipeline with default parameters. We considered up-and down-regulated genes with a p-value < 0.05 and an absolute value of log2FC > 0.05, and they were used for downstream functional analysis. The Differentially Expressed Genes (DEGs) selected were annotated using the A. thaliana genome.
The differentially expressed genes were called upon by comparing Ni with non-Ni-treated libraries separately for the TM and NTM ecotypes. Using the p.adjust function of R, all the statistical tests were corrected for multiple comparisons using the Benjamini-Hochberg false discovery rate [56]. Differences among the selected studies were adjusted using the sample normalization. To remove systematic variation between different species, the normalization procedure served as a crucial pre-processing step to adjust for the different sample sequencing depths and other confounding technical effects. We used the geometric normalization method where Fragments Per Kilobase Million (FPKMs) and fragment counts are scaled via the median of the geometric means of fragment counts across all libraries. All the differentially regulated genes IDs of N. goesingensis were mapped to closely related A. thaliana genome alignment (an overall sequence identity of ca. 86%, based on ITS1 and ITS2 [57]) and the corresponding best hit TAIR ID (hypothetical ortholog) was determined using the annotation file downloaded from Phytozome. MapMan [58] was used with the A. thaliana mapping file to map and visualize the metabolic overview, hormone regulation, transcription factors, and transport-related categories. Approximately 65% of N. goesingensis reads were mapped to A. thaliana genome (TAIR). The comparison was performed within each genotype (TM and NTM) between nickel treatment and the control condition. Up-and down-regulated genes were, respectively, more and less expressed in response to nickel presence.

Data Availability
The RNA-seq data reported in this article were deposited in NCBI's Gene Expression Omnibus (GEO) and are accessible through GEO Series accession number GSE155054.

Statistical Analysis
Statistical analysis was performed using Statistica ver. 13.0 (StatSoft) for all the obtained data. Normality of all data was evaluated by the Shapiro-Wilk test. Variance homogeneity was tested by Levene's test. For pairwise comparisons, the t-Student test was used.
Supplementary Materials: The following are available online at http://www.mdpi.com/2223-7747/9/8/951/s1, Table S1: The number of reads for each ecotype/treatment combination and the number of reads successfully mapped to the reference genome vs. the number of reads not mapped.