Agronomic and Metabolomic Side-Effects of a Divergent Selection for Indol-3-Ylmethylglucosinolate Content in Kale (Brassica oleracea var. acephala)

Brassica oleracea var. acephala (kale) is a cruciferous vegetable widely cultivated for its leaves and flower buds in Europe and a food of global interest as a “superfood”. Brassica crops accumulate phytochemicals called glucosinolates (GSLs) which play an important role in plant defense against biotic stresses. Studies carried out to date suggest that GSLs may have a role in the adaptation of plants to different environments, but direct evidence is lacking. We grew two kale populations divergently selected for high and low indol-3-ylmethylGSL (IM) content (H-IM and L-IM, respectively) in different environments and analyzed agronomic parameters, GSL profiles and metabolomic profile. We found a significant increase in fresh and dry foliar weight in H-IM kale populations compared to L-IM in addition to a greater accumulation of total GSLs, indole GSLs and, specifically, IM and 1-methoxyindol-3-ylmethylGSL (1MeOIM). Metabolomic analysis revealed a significant different concentration of 44 metabolites in H-IM kale populations compared to L-IM. According to tentative peak identification from MS interpretation, 80% were phenolics, including flavonoids (kaempferol, quercetin and anthocyanin derivates, including acyl flavonoids), chlorogenic acids (esters of hydroxycinnamic acids and quinic acid), hydroxycinnamic acids (ferulic acid and p-coumaric acid) and coumarins. H-IM kale populations could be more tolerant to diverse environmental conditions, possibly due to GSLs and the associated metabolites with predicted antioxidant potential.


Introduction
Crops belonging to the genus Brassica are among the top ten most agronomically and economically important vegetable species in the world. These crops show high morphological and agronomic diversity and are cultivated mainly in temperate regions of the Northern Hemisphere [1]. Kale (Brassica oleracea var. acephala) is a leafy vegetable crop that is becoming popular as a "superfood", due to its nutritional value (rich in Ca 2+ , folate, riboflavin, vitamins C, K and A), phytochemical composition (including polyphenols, glucosinolates, terpenoids and carotenoids), and its high anticarcinogenic and antioxidant potential [2]. Among Brassica phytochemicals, glucosinolates (GSLs)-sulfur compounds Error bars represent ± standard deviation (SD). Within each p cate significant differences between populations (ANOVA, p-value ≤ 0.05 results is presented as supplementary material (Table S1).

GSL Profiles
Our results indicated that the IM content was significantly hi lation than in the L-IM population ( Figure 2). The analysis of to divergent populations indicated significantly larger amounts in dry weight) population than in the L-IM population (24.3 μmol/g This higher levels of GSLs are mainly due to the accumulation of i the H-IM population (21.4 μmol/g dry weight) since we did not  (Table S1).

GSL Profiles
Our results indicated that the IM content was significantly higher in the H-IM population than in the L-IM population ( Figure 2). The analysis of total GSL content of the divergent populations indicated significantly larger amounts in the H-IM (34.1 µmol/g dry weight) population than in the L-IM population (24.3 µmol/g dry weight) (p < 0.01). This higher levels of GSLs are mainly due to the accumulation of indole GSLs (p < 0.01) in the H-IM population (21.4 µmol/g dry weight) since we did not observe differences between both populations for total aliphatic GSLs content (p = 0.0989). Aside from IM, we only observed a significant increase on the accumulation of 1MeOIM (p < 0.01) in the H-IM population (6.6 µmol/g dry weight) and a decrease on the accumulation of 3mSOp (6.4 µmol/g dry weight) (p = 0.0008) (Figure 2), indicating that the selection method was quite specific in increasing IM content.

Metabolomic Profiles
In order to identify the metabolomic changes that could take place in the kale p lation due to divergent selection on IM content, we performed a non-targeted metab ics analysis. Statistical univariate analyses reported 109 features that were differen accumulated in the H-IM and L-IM populations ( Figure 3). Data was then hand-fil taking into account retention time and correlation coefficients to remove feature were most likely due to in-source fragmentation of metabolites. Ultimately, 67 featu these were considered to be true metabolites (30 detected in negative and 37 in po ionization modes) ( Table 1). Forty-four of these metabolites had increased concentr in H-IM populations compared to L-IM populations, while 23 metabolites had decr concentrations. IM was not among the selected metabolites since though showed si cantly higher levels in the H-IM population vs. L-IM (p value = 0.001) in the metabol analysis, did not fit the fold change condition. Error bars represent ± standard error (SE). Within each glucosinolate category different letters indicate significant differences between divergent populations (ANOVA, p-value ≤ 0.05). Complete ANOVA table results is presented as Supplementary Material (Table S2). Abbreviations: GSLs: Glucosinolates, 3mSOp: 3-(methylsulfinyl)propylGSLs, 2-propenyl: 2-propenylGSL, IM: indol-3-ylmethylGSL, 4HOIM: 4-hydroxyindol-3-ylmethylGSL, 4MeOIM: 4-methoxyindol-3-ylmethylGSL and 1MeOIM: 1-methoxyindol-3-ylmethylGSL.

Metabolomic Profiles
In order to identify the metabolomic changes that could take place in the kale population due to divergent selection on IM content, we performed a non-targeted metabolomics analysis. Statistical univariate analyses reported 109 features that were differentially accumulated in the H-IM and L-IM populations ( Figure 3). Data was then hand-filtered, taking into account retention time and correlation coefficients to remove features that were most likely due to in-source fragmentation of metabolites. Ultimately, 67 features of these were considered to be true metabolites (30 detected in negative and 37 in positive ionization modes) ( Table 1). Forty-four of these metabolites had increased concentration in H-IM populations compared to L-IM populations, while 23 metabolites had decreased concentrations. IM was not among the selected metabolites since though showed significantly higher levels in the H-IM population vs. L-IM (p value = 0.001) in the metabolomics analysis, did not fit the fold change condition.

Metabolomic Profiles
In order to identify the metabolomic changes that could take place in the kale population due to divergent selection on IM content, we performed a non-targeted metabolomics analysis. Statistical univariate analyses reported 109 features that were differentially accumulated in the H-IM and L-IM populations ( Figure 3). Data was then hand-filtered, taking into account retention time and correlation coefficients to remove features that were most likely due to in-source fragmentation of metabolites. Ultimately, 67 features of these were considered to be true metabolites (30 detected in negative and 37 in positive ionization modes) ( Table 1). Forty-four of these metabolites had increased concentration in H-IM populations compared to L-IM populations, while 23 metabolites had decreased concentrations. IM was not among the selected metabolites since though showed significantly higher levels in the H-IM population vs. L-IM (p value = 0.001) in the metabolomics analysis, did not fit the fold change condition.  When possible, a molecular formula was assigned to each metabolite based on the exact mass and the isotopic pattern. Tentative identification was performed based on the molecular formula and MS/MS fragmentation pattern. We were able to tentatively assign compound names to 52 out of 67 metabolites. The majority of these compounds are classified as phenolics (80%). Among them, 70% are classified as flavonoids (kaempferol, quercetin and anthocyanins derivates, including acyl flavonoids), 18% as chlorogenic acids (esters of hydroxycinnamic acids and quinic acid), 6% as hydroxycinnamic acids (ferulic acid and p-coumaric acid) and 6% as coumarins. Flavonoids were identified on the bases of the aglycone fragment (Figure 4a). Deviation of the aglycone m/z can be observed in some signals on Table 1 (i.e., 285.04 or 284.03 on kaempferol glycosides) due to homolytic or heterolytic fragmentation. The homolytic fragmentation of flavonoid glycosides produces a radical aglycone ion [Y 0 −H] −• (m/z 284.03 for kaempferol), whereas the heterolytic fragmentation produces an aglycone fragment ion [Y 0 ] − (m/z 285.04 for kaempferol).
We propose the tentatively identification of the ion at m/z 980.26 (RT: 11.2 min) as kaempferol-sophoroside-(dihydroxymethoxy)-sophoroside ( Table 1). The most abundant fragment of this compound corresponds to a kaempferol-sophoroside (m/z 609.15) (Figure 4b). The kaempferol aglycone was also confirmed by the presence of a peak at m/z 284.03. Identification of a kaempferol-3-O-sophoroside-7-O-sophoroside (m/z 934.25) was previously reported in B. oleracea [38]. The neutral loss of m/z 371.11 (980.26 → 609.15) may indicate the loss of an anhydrohexose attached to a dihydroxymethoxy cyclohexane. Finally, the neutral loss of m/z 47.01 (980.26 → 933.25), supports the hypothesis of the presence of a dihydroxymethoxy group. This is, however, just a proposed structure and a conclusive elucidation will require further analysis.
Significant metabolites are evenly distributed between the groups of compounds with higher and lower concentration in the H-IM population compared to the L-IM population. The only exception is the group of coumarins that showed a higher concentration in the H-IM population. Aside from phenolics, we identified two compounds: kynurenic acid, a product of the kynurenine branch of tryptophan metabolism and an indolylacetyl dihexoside, a carbohydrate derivative.

Discussion
Changes on the metabolome allow plants to adapt to fluctuations in the environmental conditions. The accumulation of specific metabolites, especially those with antioxidant properties, act as a metabolic buffer under stressful conditions. It has been demonstrated that different environmental and cultivation conditions modify the profile and content of GSLs in Brassica crops [39,40]. In general terms, abiotic stresses tend to increase the content of GSLs in these plants, suggesting that GSLs may play a role on plant adaptation to different environments. However, the direct role of these compounds in plant adaptation has

Discussion
Changes on the metabolome allow plants to adapt to fluctuations in the environmental conditions. The accumulation of specific metabolites, especially those with antioxidant properties, act as a metabolic buffer under stressful conditions. It has been demonstrated that different environmental and cultivation conditions modify the profile and content of GSLs in Brassica crops [39,40]. In general terms, abiotic stresses tend to increase the content of GSLs in these plants, suggesting that GSLs may play a role on plant adaptation to different environments. However, the direct role of these compounds in plant adaptation has yet to be addressed. In this work we used two divergently selection populations (H-IM and L-IM) to study that possible role.
Our results may indicate that the leaf productivity of Brassica crops could be directly or indirectly affected by IM content and that this effect is stable across different environments. We observed a higher yield in the H-IM population in the various experimental locations from southern Spain to northern Norway. To the best of our knowledge this is the first time that the possible role of the IM has been studied directly in relation to local adaptation. Interestingly, both populations barely differ in the amount of other GSLs, so differences in yield could be attributed to a great degree to the accumulation of IM. However, it is plausible that a higher accumulation of IM after three cycles of divergent selection may also produce a reorganization of the plant metabolome, that could contribute to increased plant adaptability.
To study the extent of potential metabolome reorganization we performed an untargeted metabolomics analysis. Extraction conditions (80% MeOH), chromatographic setup (reverse-phase UPLC) and ionization interface (ESI) used in our analysis allowed for detection of a wide range of polar and mid-polar metabolites, but with a lack of information about highly polar (elute with the dead volume of the chromatography system) or apolar compounds. With this limitation in mind, our analysis indicate that divergent selection mainly affected phenolic biosynthetic pathways. More than 80% of the metabolites identified were phenols. Previous studies have reported a simultaneous increase of IM and total phenols in various crucifers (Isatis canescens, B. oleracea var. italica, or B. rapa ssp. rapa) [41][42][43]. However, our study represents the first example of a possible relationship between higher IM content and higher phenolic compound content in kale.
Phenolic compounds constitute a complex group of secondary metabolites that are widespread in the plant kingdom. They have allelopathic, antimicrobial and antioxidant activity in plants [44][45][46] and can be precursors of other secondary metabolites (e.g., hydroxycinnamic acids are precursor of lignin) [47]. It is not surprising that most of the phenolics we identified were flavonoids given that they are the most prominent phenolics in Brassica species [24]. Severe stress conditions activate the biosynthesis of flavonoids, which in turn act as an antioxidant system preventing cellular damage. It may be hypothesized that the high levels of IM, a stress-promoted molecule, are perceived by the plant as an indicator of stressful conditions, resulting in the activation of flavonoid biosynthesis. Based on our agronomic results, this is unlikely since plants with an imbalanced defensive response show lower growth rates [48,49]. Some authors suggest a direct or indirect role of flavonoids as growth regulators. Grandmaison et al. [50] reported that flavonoids can regulate cell development by interaction with nuclear proteins. Supporting this idea, Saslowsky et al. [51] demonstrate that the end products of flavonoid biosynthesis are located in the cytoplasm and the nuclei of the tip cells of Arabidopsis roots where they can interact with auxin biosynthesis. In vitro analysis shows that in both subcellular compartments, cytoplasmic and nuclear, flavonoids can also interact with actin, regulating its polymerization [52]. This interaction is structure dependent, with flavonols acting as inhibitors and flavanes as stimulators of actin polymerization [52].
We tentatively identified several derivatives of flavonol (quercetin, kaempferol, (iso) rhamnetin) and anthocyanin (cyanidin, petunidin and delphinidin) that accumulated differentially in the H-IM and L-IM populations. Only cyanidin glycosides, along with the group of coumarins, accumulate in a higher extent in the H-IM than in the L-IM population. Curiously, these two groups of compounds have been reported to inhibit plant development [53,54], so further studies will be necessary to elucidate the role of these compounds in kale growth.

Plant Populations
Two divergently selected kale populations, one with high (H-IM) and one with low (L-IM) IM content, were used in this study. These two populations were selected from a local Spanish population (MBG-BRS0062), kept at the Brassica germplasm bank at Misión Biológica de Galicia (MBG-CSIC) (Pontevedra, Spain). These populations had been subjected to three selection cycles (details explained in Sotelo et al.) [21].

Growing Conditions and Locations
Kale seeds (H-IM and L-IM) were sown in multipot-trays in a greenhouse. At the 5-6 leaf stage, plants (

Agronomic Parameters
Fresh weight was quantified using twenty-five fully developed leaves (7th-8th leaf from the apex) from each plot, harvested randomly. The same leaves were subsequently dried at 70 • C until a constant weight was reached to record the dry mass. Plant height was measured from the soil surface to the base of the upper leaf in 10 plants from each plot.

Biochemical Analysis
For GSLs and non-targeted metabolomics analyses, the 4th leaf from the apex of 15 plants/plot were collected in liquid nitrogen and stored at −80 • C until freeze-dried in a lyophilizer (GAMMA 2-16 LSC plus; Christ, Osterode am Harz, Germany). Samples were mechanically milled to a fine powder in a grinder (Janke and Kunkel A10 mill; IKA-LabortechnikStaufen, Staufen, Germany) before metabolite extraction.

Metabolomic Analysis
Freeze-dried powder (50 mg) was dissolved in 500 mL of 80% aqueous methanol and then sonicated for 15 min. After centrifugation for 10 min (16,000× g, at room temperature), the extract was filtered through a 0.20 µm micropore PTFE membrane and placed in vials for further analysis. For metabolomic composition analysis we used ultra-performance liquid chromatography coupled with electrospray ionization quadrupole (Thermo Dionex Ultimate 3000 LC; Thermo Fisher Scientific, Waltham, MA, USA) time-of-flight mass spectrometry (UPLC-Q-TOF-MS/MS) (Bruker Compact™) with a heated electrospray ionization (ESI) source. Chromatographic separation was performed in an Intensity Solo 2 C18 column (2.1× 100 mm 1.7 µm pore size; Bruker Daltonics, Billerica, MA, USA) using a binary gradient solvent mode consisting of 0.1% formic acid in water (solvent A) and acetonitrile (solvent B). The following gradient was used: 3 % B (0-4 min), from 3% to 25 % B (4-16 min), from 25 to 80% B (16-25min), from 80 to 100% B (25-30 min), hold 100% B until 32 min, from 100% to 3% B (32-33 min), hold 3% B until 36 min. The injection volume was 5 µL, the flow rate was established at 0.4 mL/min and column temperature was controlled at 35 • C. MS analysis was operated in spectra acquisition range from 50 to 1200 m/z. Both polarities (±) of ESI mode were used under the following specific conditions: gas flow 9 L/min, nebulizer pressure 38 psi, dry gas 9 L/min, and dry temperature 220 • C. Capillary and end plate offset were set to 4500 and 500 V, respectively. The instrument was calibrated externally with a calibration solution of 1mM sodium formate/acetate in iPrOH/H 2 O 50/50 with 0.2% formic acid directly infused to the source. Before sample injections, LC-qTOF system stability was tested by three consecutive injections of chloramphenicol (ESI-mode; ∆RT= 0.02 min; ∆m/z = 0.002) and triphenyl phosphate (ESI + mode; ∆RT = 0.02 min; ∆m/z = 0.001). The calibration solution was injected at the beginning of each run and all the spectra were calibrated prior to statistical analysis. MS/MS analysis was performed based on the previously determined accurate mass and RT and fragmented by using different collision energy ramps to cover a range from 15 to 50 eV. The algorithm T-Rex 3D from the MetaboScape 4.0 software (Bruker Daltonics, Billerica, MA, USA) was used for peak alignment and detection.

Statistical Analysis
Parametric statistical analysis was performed using the GLM procedure of SAS 9.4 (SAS Institute Inc., Cary, NC, USA) for agronomic traits and GSLs content. Populations were considered fixed effects and locations were considered random effects. A post hoc ANOVA analysis was performed using the Fisher's protected least significant difference (LSD) at p ≤ 0.05.
Statistical analysis of metabolomic data was performed using the web-based software Metaboanalyst [56]. In order to remove non-informative variables, data were filtered using the interquantile range filter (IQR). Moreover, Pareto variance scaling was used to remove the offsets and adjust the importance of high-and low-abundance ions to an equal level. The resulting three-dimensional matrix (peak indices, samples and variables) was further subjected to statistical analysis. Univariate analysis (one-way ANOVA) with a p value ≤ 0.05 was carried out to find differentially expressed metabolites. Using the Volcano Plot (VP) approach, which measure differentially accumulated metabolites based on t-statistics and fold changes simultaneously, we also highlighted the metabolites with a |log 2 (FC)| ≥ 1 and statistically significant difference (FDR ≤ 0.05) between populations.

Tentative Metabolite Identification
Tentative compound identification was performed using accurate mass metabolites reported in different publicly available databases such as METLIN, KEGG, Pubchem, HMDB and Plant Metabolic Network. Additionally, further partial identification of the most significant metabolites was made by comparison of MS/MS fragmentation patterns against reference compounds found in previously mentioned databases and bibliography on plants of the Brassicaceae family.

Conclusions
We reported a higher yield of H-IM populations of kale across different environments compared to L-IM populations indicating a potentially greater adaptive capacity of the H-IM populations to varied contexts, as measured by a higher production of foliar biomass. The GSL profiles analysis showed a higher content in indole GSLs in H-IM populations, previously described secondary metabolites which are thought to impart higher tolerance to abiotic stresses such as salinity [57] or drought [58]. The H-IM populations of kale had higher concentrations of compounds, which tentatively can be predicted to have antioxidant potential that may contribute to tolerance of abiotic stresses by reducing the generation of reactive oxygen species [59]. The high indole GSL content and the accumulation of other secondary metabolites may give the H-IM populations of kale an improved adaptive capacity under varied environmental conditions, which may be responsible for an observed higher yield of the high indole GSL population.