Metabolic Proﬁling of Varronia curassavica Jacq. Terpenoids by Flow Modulated Two-Dimensional Gas Chromatography Coupled to Mass Spectrometry

: In this study, a metabolomic approach was used to investigate the e ﬀ ect of seasonality on the chemical composition and yield of anti-inﬂammatory active principle, α -humulene, in the essential oil of three genotypes of Varronia curassavica Jacq. (Syn. Cordia verbenaceae ). The essential oils were extracted by hydrodistillation and analyzed by comprehensive two-dimensional gas chromatography coupled to mass spectrometry (GC × GC-MS). The GC × GC approach a three-fold improvement in qualitative analysis (48 compounds were identiﬁed by GC-MS versus 135 by GC × GC-MS). The improved resolving power of GC × GC resolved important coelutions and enabled the detection of unusual substances in V. curassavica essential oil. The chromatographic data was analyzed by using peak table-based chemometrics, namely, principal component analysis (PCA) and hierarchical cluster analysis (HCA). The metabolic study showed that seasonality has a signiﬁcant e ﬀ ect on the chemical composition. The α -humulene content was a ﬀ ected by genotype and season. Spring and summer were the best harvest seasons for the yield of the active ingredient, found in higher concentrations in the VC2 genotype. The proposed metabolomic workﬂow was successfully applied to terpene analysis found in V. curassavica essential oil, and such results have broadened our understanding of the inﬂuence of seasonal factors on the specialized metabolism of the species.


Introduction
Essential oils (EOs) are sources of biologically active compounds that can be applied in food, agronomic, cosmetic, and pharmaceutical industries for the elaboration of perfumes, cosmetics, pharmaceuticals, and used as adjuvant in the formulation of medicines. An important example of this application is the essential oil of Varronia curassavica Jacq. (syn. Cordia verbenacea)-a Brazilian native species, which is mainly found in the Atlantic Rainforest biome, popularly known as 'erva-baleeira' [1].
V. curassavica is a medicinal and aromatic plant from the Boraginaceae family with significant economic importance. The EO produced from its leaves is used in the formulation of the topic anti-inflammatory medicine Acheflan ® , a phytotherapeutic product fully developed with Brazilian technology and currently exported to Japan, Chile, Mexico, Costa Rica, Equator, Peru, the United States, and Canada [2,3]. In Brazil, the product is the leading prescription drug for chronic tendonitis and muscle pain, being responsible for more than 25% of the market sales [3]. The active principle responsible for anti-inflammatory activity of the EO is the sesquiterpene α-humulene [4][5][6]. The substance α-humulene is considered a chemical marker in the pharmaceutical industry for the quality control of V. curassavica essential oil. Besides its anti-inflammatory activity, the species has been demonstrated to

Plant Material
Three genotypes of V. curassavica (VC1, VC2, and VC3) were evaluated, comprising a total of nine plants. The plants, obtained from commercial matrices and propagated by stables (VC1 and VC2) and tissue culture (VC3), were kept in 20-L plastic roofed pots at the Center of Plant Genetic Resources of the Agronomic Institute (IAC), located in the city of Campinas, State of São Paulo (22 • 52 south latitude, 47 • 04 west longitude, and 677 m altitude). The regional climate, according to the Köppen classification, is classified as Cwa-subtropical climate with hot summers and dry winters. The climatological data for the sampling period were obtained from the meteorological station of the Agrometric Information Center (Ciiagro) [50]. The aerial parts were harvested in winter (16 August 2018), spring (21 November 2018), summer (21 February 2019) and autumn (30 May 2019). The harvested leaves of V. curassavica were manually separated from the stems and dried at 40 • C in an air-circulation oven for 72 h.

Extraction of Essential Oil
Essential oils were extracted by hydrodistillation in a Clevenger apparatus for a period of two hours. The essential oil was stored in a freezer (−4 • C) for further characterization of the chemical composition.

Comprehensive Two-Dimensional Gas Chromatographic Analysis (GC×GC-FID/MS)
Analyses of the essential oils from V. curassavica were performed on a GC×GC system, which comprised of a TRACE 1300 GC equipped with a flame ionization detector (FID) and an ISQ single transmission quadrupole mass spectrometer (QMS) (Thermo Scientific-Waltham, MA, USA). Differential flow modulation using the reverse fill/flush configuration was attained using INSIGHT modulator (SepSolve Analytical-Frankfurt, Germany). Data digitalization using FID was attained at 120 Hz. A 40-450 m/z scanning range were used producing 31 scans s −1 . The transfer line and the ion source were operated at 300 • C and 220 • C, respectively. Helium was used as carrier gas and auxiliary gas at constant flow rates of 0.5 mL min −1 and 20.0 mL min −1 . The column set consisted of a primary 30 m × 0.25 mm-id × 0.25 µm (β of 250) MEGA-5HT and a secondary 5 m × 0.25 mm-id × 0.25 µm (β of 250) HP-50 + wall coated open tubular (WCOT) capillary columns (Agilent Technologies-Santa Clara, CA, USA). Deactivated fused silica capillaries were used as transfer lines for simultaneous MS (5 m × 0.18 mm-id) and FID (5 m × 0.32 mm-id) detection. This setup ensured a reproducible splitting of the 2 D effluent to FID (70%) and MS (30%). Sample introduction was performed using a split ratio of 1:20. The GC oven was programmed from 60 to 240 • C at 3 • C min −1 . Modulation period was set to 5.0 s with a re-injection (flush) pulse of 100 ms. The chromatographic analyses were performed in triplicate.
Xcalibur software (Thermo Scientific-Waltham, MA, USA) was used for data acquisition and processing was performed using GC Image software (Zoex-Houston, TX, USA). Identification of the terpenes was performed comparing the substance mass spectra with NIST database and filtered by the retention index (LTPRI) using the NIST WebBook and published literature [51]. Retention indices were determined by the injection of a homologous series of C 8 -C 20 n-alkanes (04070-1ML) (Merck-St. Louis, MO, USA) using the Van den Dool and Kratz equation [52]. The analytes were tentative identified through comparative analysis of mass spectra considering a minimum similarity correspondence of 80% and a deviation of ±25 LTPRI from the NIST and Adams [51].

Data Analysis
Analysis of variance (ANOVA) followed by Tukey test (p < 0.05) was applied to verify the differences between essential oil yields and α-humulene contents. The results of the chemical analyses were submitted to multivariate statistical analyses, namely, principal components analysis (PCA), and hierarchical clustering analysis (HCA). Such models were built using the software XLSTAT-2020 version (Addinsoft-Bordeaux, France). Heatmap were constructed using the on-line tool MetaboAnalist 4.0 [53].

Essential Oil Content
The average yield of the essential oils of the three genotypes harvested during the winter, spring, summer, and autumn seasons is shown in Table 1. The yield of essential oils varied from 0.22 to 0.60% (v/w). The yields are in agreement with those reported in other studies, which showed values ranging from 0.16 to 2.74%. [46][47][48]. The analysis of variance (ANOVA) showed that there was a significant difference in the yield of essential oils (p < 0.05) depending on the season. Spring was the season with the highest production (0.533%) and winter the season with the lowest production (0.291%). Between summer and autumn there was no significant difference (0.407 and 0.397%, respectively). The VC2 genotype showed the highest production of essential oil (0.60 ± 0.03%) in spring, a period that coincides with the flowering of the plant. Although this species blooms at any time of the year, it flowers with higher intensity during spring.
The correlation of the results presented in Table 1 with the meteorological data ( Figure 1) allows us to ascertain that the climatic conditions also constitute a strong factor. The rainfall exhibited a significant variation, since the accumulated rainfall indexes that varied from 75 to 484 mm. Matias and collaborators [46] also observed the influence of rainfall on the production of V. curassavica essential oil, which corroborates our findings. Analysis of variance (ANOVA) followed by Tukey test (p < 0.05) was applied to verify the differences between essential oil yields and α-humulene contents. The results of the chemical analyses were submitted to multivariate statistical analyses, namely, principal components analysis (PCA), and hierarchical clustering analysis (HCA). Such models were built using the software XLSTAT-2020 version (Addinsoft -Bordeaux, France). Heatmap were constructed using the on-line tool MetaboAnalist 4.0 [53].

Essential Oil Content
The average yield of the essential oils of the three genotypes harvested during the winter, spring, summer, and autumn seasons is shown in Table 1. The yield of essential oils varied from 0.22 to 0.60% (v/w). The yields are in agreement with those reported in other studies, which showed values ranging from 0.16 to 2.74%. [46][47][48]. The analysis of variance (ANOVA) showed that there was a significant difference in the yield of essential oils (p < 0.05) depending on the season. Spring was the season with the highest production (0.533%) and winter the season with the lowest production (0.291%). Between summer and autumn there was no significant difference (0.407 and 0.397%, respectively). The VC2 genotype showed the highest production of essential oil (0.60 ± 0.03%) in spring, a period that coincides with the flowering of the plant. Although this species blooms at any time of the year, it flowers with higher intensity during spring.
The correlation of the results presented in Table 1 with the meteorological data ( Figure 1) allows us to ascertain that the climatic conditions also constitute a strong factor. The rainfall exhibited a significant variation, since the accumulated rainfall indexes that varied from 75 to 484 mm. Matias and collaborators [46] also observed the influence of rainfall on the production of V. curassavica essential oil, which corroborates our findings.

Characterization of the Essential Oil
The chemical compositions of the EOs of the three genotypes of V. curassavica, harvested in four different seasons, were analyzed by GC-MS and GC×GC-MS, as shown in Figure 2. Many substances that were not apparent in 1D-GC, due to overlapping peaks or because they are found in low intensities, are clearly revealed in the GC×GC chromatogram. An example of peaks that were co-eluted in the first dimension ( 1 D) and were resolved in the second dimension ( 2 D), as shown in Figure 3. The composite system resolved the coelution between germacrene D, (E)-β-ionone and tau-elemene ( Figure 3A) and between δ-cadinene, cadin-1(2),4-diene, and two unidentified compounds ( Figure 3B).  The GC×GC analysis enabled the tentative identification of 135 peaks compared to 48 peaks with GC-MS resulting in a three-fold improvement in qualitative analysis. Previous reports corroborated this improved feature for qualitative analysis [37,38,40,54]. For instance, Santos et al. [24] evaluated the essential oil of two species of Piperaceas (Manekia obtusa and Piper cubataonum). In M. obtuse oil, 80 compounds were identified by GC×GC-MS, while 22 only were assigned by GC-MS. For the leaf and branch oils of P. cubataonum, 57 and 66 compounds were identified by GC×GC compared with only 14 and 20 compounds by GC/MS, respectively.
It is important to note that although 168 compounds were detected in the GC×GC analysis, only 135 (~80%) were tentatively identified. The presence of very similar mass spectra imposed restraints to fully elucidate the compounds present in the EO.  Table 2. The chemical classes of the analytes found in the EO of V. curassavica were monoterpenes and sesquiterpenes (~96% of the oil).
Among the identified compounds, neryl acetate was detected only in genotypes VC2 and VC3 in the winter samples. β-copaene was detected in genotype VC1 (winter and autumn) and genotype VC2 (summer only). Terpene cis-muurol-5-en-4-ol was detected only in genotype VC2 (winter and spring), while cis-sesquisabinene hydrate was found only in VC3 (winter and spring). This variation in EO composition with seasonality is supported by previous reports [55,56]. Furthermore, different genotypes cultivated under the same conditions may exhibit differential expressions resulting in distinct chemical profiles, including the occurrence or absence of metabolites [57].  Table 2. The chemical classes of the analytes found in the EO of V. curassavica were monoterpenes and sesquiterpenes (~96% of the oil).
Among the identified compounds, neryl acetate was detected only in genotypes VC2 and VC3 in the winter samples. β-copaene was detected in genotype VC1 (winter and autumn) and genotype VC2 (summer only). Terpene cis-muurol-5-en-4-ol was detected only in genotype VC2 (winter and spring), while cis-sesquisabinene hydrate was found only in VC3 (winter and spring). This variation in EO composition with seasonality is supported by previous reports [55,56]. Furthermore, different genotypes cultivated under the same conditions may exhibit differential expressions resulting in distinct chemical profiles, including the occurrence or absence of metabolites [57]. Table 2. Chemical composition (relative %) of Varronia curassavica essential oils harvested in different seasons and analyzed by GC×GC-FID/MS. LTPRI -Linear temperature programmed retention indices; Exp.-LTPRI experimental obtained by the injection of a homologous series of C8-C20 n-alkanes using the Van den Dool and Kratz equation [52]; Lit. -LTPRI obtained from literature [51] and NIST WebBook. (*)-LTPRI not found; (**)-Substance identified by comparison with Sciarrone and collaborators [44]; (ta)-Trace amounts (≤0.003); (-)-Absence of the metabolite.      The substance α-humulene is the chemical marker in the pharmaceutical industry for the quality control of V. curassavica essential oil, which varied from 1.47% to 3.93%. Variance analysis (ANOVA) demonstrated significant differences (p < 0.05) between the seasons and between the genotypes evaluated (Table 3). Spring and summer were the seasons with the highest average production of the active principle α-humulene (3.82% and 3.93%, respectively), while autumn was the season with the lowest average production (1.38%). The highest average relative proportion was observed in the VC2 genotype during summer (3.93 ± 0.10) ( Table 3). Besides environmental factors, genetic characteristics can also determine and modify the production rate of essential oils [58][59][60]. Table 3. Contents of α-humulene (% Relative and SE Average) in the EO obtained from V. curassavica genotypes harvested at different seasons and analyzed by GC×GC-MS. Distinct uppercase letters (A, B and AB) for the rows and lowercase (a, b and ab) for columns represent statistically significant differences (p < 0.05). The heatmap plot (Figure 4), generated from hierarchical clustering, provides an overview of V. curassavica metabolic profile in response to the seasonal cycle and plant genetics. It is possible to observe that the genotype was responsible for most of the variance observed in the samples, highlighting two main blocks. The first block (upper eighteen rows) comprises a group of metabolites that are most abundant and present/absent in the VC3 genotype. For instance, α-pinene (2) is most abundant in the samples of VC3 genotype. The (E)-α-bergamotenal (53), (E)-α-santalal (54), e α-santalene (27) are present only in the VC3 genotype. The second block showed similar metabolic profiles, which is the case for the VC1 and VC2 genotypes, which is ascertained by the results attained with PCA analysis, as presented in the following paragraphs.

α-Humulene
The most different chemical composition was found for the VC3 genotype (Figure 4), suggesting the contribution of other regulatory levels, besides environmental and genetic, in the chemical composition of the essential oils, e.g., physiological cycle of the plant and interaction between genotype and environment. The observations on the intra-class variance may contribute to V. curassavica breeding programs aiming at increasing the production of substances with biological activity.
To investigate the effects of seasonality on the volatile compounds of V. curassavica, a peak table-based multivariate approach, namely, principal component (PCA) was adopted for pattern recognition of the three evaluated genotypes. A two-component PCA model expressed 71.8% of the total variance with PC1 being responsible for describing 61.8% and PC2 10% of the variance. The scores graph illustrated in Figure 5A indicated that the information contained in PC1 may be used to distinguish the samples into two groups (VC3 and VC1 + VC2). The loadings plot ( Figure 5B) describes the most important variables that contributed to such pattern.
The substances responsible for the observed clustering between the samples were α-pinene (2), found in a higher relative proportion in the VC3 genotype samples, and (E)-α-bisabolene  To investigate the effects of seasonality on the volatile compounds of V. curassavica, a peak table-based multivariate approach, namely, principal component (PCA) was adopted for pattern recognition of the three evaluated genotypes. A two-component PCA model expressed 71.8% of the total variance with PC1 being responsible for describing 61.8% and PC2 10% of the variance. The scores graph illustrated in Figure 5A indicated that the information contained in PC1 may be used to distinguish the samples into two groups (VC3 and VC1 + VC2). The loadings plot ( Figure 5B) describes the most important variables that contributed to such pattern.
The substances responsible for the observed clustering between the samples were α-pinene (2), found in a higher relative proportion in the VC3 genotype samples, and (E)-α-bisabolene (42)  In the current study, our results suggested that the genetic make-up is likely the most important factor contributing to the chemical composition of V. curassavica essential oils, whereas seasonal contributes to modulate the metabolite profile. Several factors have been demonstrated to contribute to the metabolic profile of plant essential oils [61,62], such as abiotic stress [63][64][65], genetics [58,59], and plant-pathogen interactions [60,[66][67][68]. Furthermore, environmental biotic and abiotic conditions modulate several regulatory levels determined by the plant genotype affecting the profile of specialized metabolites of essential oils [57,69]. The determination of the contribution of distinct factors to the chemical profile of plant essential oils may help breeding and agricultural strategies to In the current study, our results suggested that the genetic make-up is likely the most important factor contributing to the chemical composition of V. curassavica essential oils, whereas seasonal contributes to modulate the metabolite profile. Several factors have been demonstrated to contribute to the metabolic profile of plant essential oils [61,62], such as abiotic stress [63][64][65], genetics [58,59], and plant-pathogen interactions [60,[66][67][68]. Furthermore, environmental biotic and abiotic conditions modulate several regulatory levels determined by the plant genotype affecting the profile of specialized metabolites of essential oils [57,69]. The determination of the contribution of distinct factors to the chemical profile of plant essential oils may help breeding and agricultural strategies to increase the yield of metabolites of interest, as well as in the selection of genotypes for the development of breeding programs.

Conclusions
Two-dimensional gas chromatography has proved to be a versatile tool for analyzing complex mixtures such as V. curassavica essential oil. GC×GC-QMS allowed the identification of 135 constituents in the essential oil, three times more when compared to 1D-GC. The improved resolving power of GC×GC has resolved co-eluting peaks facilitating the detection of trace substances in the essential oil of V. curassavica. The metabolic study of volatile constituents showed that the metabolic profile and the content of the anti-inflammatory active principle, α-humulene, were significantly influenced by the genotype and seasonality. The content of the active principle was higher during spring and summer, being recommended for leaf harvest. Among the investigated genotypes, VC2 displayed higher amounts of the active principle.