Metabolomics and DNA-Based Authentication of Two Traditional Asian Medicinal and Aromatic Species of Salvia subg. Perovskia

Subgenus Perovskia of the extended genus of Salvia comprises several Central Asian medicinal and aromatic species, of which S. yangii and S. abrotanoides are the most widespread. These plants are cultivated in Europe as robust ornamentals, and several cultivars are available. However, their medicinal potential remains underutilized because of limited information about their phytochemical and genetic diversity. Thus, we combined an ultra-high performance liquid chromatography quadrupole time of flight mass spectrometry (UHPLC-QTOF-MS) based metabolomics with DNA barcoding approach based on trnH-psbA and ITS2 barcodes to clarify the relationships between these two taxa. Metabolomic analysis demonstrated that aerial parts are more similar than roots and none of the major compounds stand out as distinct. Sugiol in S. yangii leaves and carnosic acid quinone in S. abrotanoides were mostly responsible for their chemical differentiation, whereas in roots the distinction was supported by the presence of five norditerpenoids in S. yangii and two flavonoids and one norditerpenoid in S. abrotanoides. To verify the metabolomics-based differentiation, we performed DNA authentication that revealed S. yangii and S. abrotanoides to be very closely related but separate species. We demonstrated that DNA barcoding coupled with parallel LC-MS profiling constitutes a powerful tool in identification of taxonomically close Salvia species.


Introduction
The subgenus Perovskia Kar. comprises several species with an ethnobotanical value and medicinal potential. These shrubby, aromatic plants are widespread throughout Central Asia, in an area that includes Iran, Pakistan, Afghanistan, Turkmenistan, Uzbekistan, and Tibet. The ranges of most of the accepted species are limited to the arid montane regions in Central Asia, but two of them, Salvia abrotanoides Karel. (previously known as Perovskia abrotanoides Kar.) and Salvia yangii B.T. Drew (formerly Perovskia atriplicifolia Benth.), are distributed over vast areas from Western Iran and Pakistan and as far as Tibet and Xinjiang, in China [1 -9].
There are a lot of discrepancies regarding the taxonomy of the subgenus Perovskia. According to primary classification based on flower morphology, species belonging to cultivars are variants of S. yangii or hybrids between S. abrotanoides and S. yangii [18,19]. This makes the seemingly simple determination of Perovskia species unexpectedly complicated. Both S. yangii and S. abrotanoides are used as a folk medicinal herb in the areas of natural occurrence. S. yangii is a folk medicinal plant used in Pakistan, Turkmenistan, and Afghanistan [21,22]. The plant has antibacterial activity and is also used as a cooling medicine in the treatment of fever [3,23]. It is known for its flavouring qualities and it is also smoked as a euphoriant and used as a decoction for chronic dysentery [21]. In Uzbekistan and Kyrgyzstan, the aboveground parts are used to heal wounds, and a decoction is used to treat scabies [17]. In traditional Tibetan and Chinese medicine, it has been claimed that S. yangii is a powerful analgesic and parasiticide agent [9,24].
S. abrotanoides, locally known as brazambol, is used to treat leishmaniasis in Iranian folk medicine. Villagers in the Isfahan Province apply a poultice, made of crushed roots of the plant, water, sesame oil, and wax, on lesions caused by cutaneous leishmaniasis [25,26]. The plant is also used locally for treatment of typhoid, fever, headache, gonorrhoea, vomiting, motion, toothache, cardiovascular diseases, liver fibrosis, painful urination, cough, and as a sedative, analgesic and antiseptic drug [21,[27][28][29][30]. Also, the essential oil of S. abrotanoides is effective as an antimicrobial agent [31].
Species of subgenus Perovskia contain a wide array of natural compounds belonging to different classes [22]. Most of the published data refer to the composition and activities of essential oils from aerial parts [32][33][34][35][36][37][38][39][40][41], indicating a significant variability associated with geographic origin and environmental factors Aerial parts contain also lignans, triterpenes, steroids, and their glycosides, as well as phenolic compounds, including considerable amounts of rosmarinic acid, which is a typical product of the Nepetoideae subfamily of Lamiaceae, known for its numerous health-promoting properties [2,3,7,15,42,43].
Roots of species belonging to the subgenus Perovskia are abundant in red-colored quinoid norabietanoids called tanshinones and are considered an alternative source of these valuable pharmacologically active compounds [25,44]. Tanshinones are the main constituents of the well-known Traditional Chinese Medicine herbal drug 'Danshen' (Salvia miltiorrhiza rhizome et radix a pharmacopoeial drug of China and Europe) [45][46][47][48].
The research performed to date on plants of subgenus Perovskia concentrated on either roots [25,44] or aerial parts [7,[49][50][51]. There are also reports on phytochemical analysis of the air-dried whole plants [2][3][4]9,24,52]. Plants were collected from their natural habitat or were cultivated. Various methods and solvents were used by researchers to extract natural products of interest which often makes the results hardly comparable. Therefore, we aimed at comprehensive analysis of phytochemical profiles of non-volatile compounds in both roots and leaves of S. yangii and S. abrotanoides cultivated in Europe under identical climatic and edaphic conditions (Wroclaw, Poland). The ethnopharmacological importance of S. yangii and S. abrotanoides has prompted us to reanalyse their chemical constituents with potential or proven bioactivity. Here, for the first time, we offer a detailed comparison of phytochemical non-volatile profile in roots and leaves of S. yangii and S. abrotanoides together with the analysis of a mutual genetic relationship between the two species. The proposed approach of a comparative analysis of closely related plants of a similar but not identical chemical profile may bring tangible results which might potentially pinpoint the outlook for follow-up work through genetic approaches, transcriptomics, proteomics, and other applications. Choosing non-model but potentially important medicinal plants allowed us to test this methodology and prove that the slight differences in morphological and metabolic phenotypes can be verified by molecular analyses, thus facilitating the proper material identification and accelerating the selection of desired traits, breeding, and herb quality control and improvement.

Plant Material
Salvia abrotanoides and S. yangii plants were grown at the experimental field in the Botanical Garden of Medicinal Plants at the Wroclaw Medical University (Poland). The geo-graphical location of the Garden is as follows: latitude 51 • 07 03" N, longitude 17 • 04 27" E, altitude 117 m a.s.l. S. abrotanoides was introduced from its original habitat in Iranian Karkas Mountains (coordinates: 33 • 24 57.3" N 51 • 46 45.4" E). Seeds of S. yangii were obtained from the collection of the Botanical Garden in Essen (Germany). Cultivars "Lacey Blue" and "Blue Spire" are available commercially as ornamental plants. Plants analyzed in this study were cultivated clonally via vegetative propagation since introduced into the Botanical Garden of Medicinal Plants at the Wroclaw Medical University. Plant voucher specimens were deposited in the herbarium of the Botanical Garden of Medicinal Plants at the Wroclaw Medical University (reference numbers: S. abrotanoides: P-122; S. yangii: P-123; cultivar "Blue Spire": P-124; cultivar "Lacey Blue": P-125).
Leaves and roots of S. abrotanoides and wild type S. yangii (Appendix A, Figure A1) were collected for phytochemical analysis and were immediately frozen in liquid nitrogen. Harvested leaf samples were finely ground using a mortar and pestle with liquid nitrogen. Root samples were finely ground using an automatic/electric grinder Analysette 3 Spartan and Pulverisette 0 Cryo-Box (Fritsch GmbH, Idar-Oberstein, Germany) filled with liquid nitrogen.
Three accessions of S. yangii (including wild type and cultivars "Lacey Blue" and "Blue Spire") and a single accession of vegetatively propagated S. abrotanoides were subject to the molecular study. For phylogenetic analysis, a single, fresh, mature leaf of each plant was harvested and ground to a fine powder in liquid nitrogen.
All samples were stored at −80 • C. The harvest of leaves and roots was done during the flowering season (August 2016).

Sample Preparation for Qualitative Analysis
For direct qualitative analysis, 100 mg of lyophilized plant material was extracted with 5 mL of MeOH: H 2 O (8:2 v/v) in an ultrasonic bath for 20 min at 25 • C. The supernatant was collected and evaporated to dryness. Extracts were suspended in 2 mL of 0.2% formic acid (FA) in H 2 O prior to preliminary purification and concentration using solid phase extraction (SPE). SPE was performed using an Oasis Max 1cc Vac Cartridge, 500 mg Sorbent per Cartridge, 30 µm Particle Size (Waters Corporation, Milford, MA, USA) conditioned with MeOH (4 mL) and H 2 O (4 mL). The sample was loaded and cartridges washed with 5 mL water in 0.1% FA. The phenolic compounds were eluted with 5 mL of 80% MeOH in 0.1% FA. Samples were re-evaporated to dryness and dissolved in 25% of MeOH to obtain 1 mg dry weight/mL baseline plant material concentration.

Liquid Chromatography-Mass Spectrometry Analysis
Ultra-high performance liquid chromatography quadrupole time of flight mass spectrometry (UHPLC-QTOF-MS) estimation of metabolite contents was carried out on a Thermo Ultimate 3000 RS (Thermo Fischer Scientific, Waltham, MA, USA) chromatographic system coupled to a Bruker Compact (Bruker, Billerica, MA, USA) quadrupole time of flight (QTOF) mass spectrometer, consisting of a binary pump, sample manager, column manager, and the PDA/diode array detector. Separations were performed on a Kinetex C18 column (2.1 × 100 mm, 2.6 µm, Phenomenex, Torrance, CA, USA), with mobile phase A consisting of 0.1% (v/v) FA in water and mobile phase B consisting of 0.1% (v/v) FA in acetonitrile. The samples (2 µL) were injected, and a linear gradient from 15% to 70% phase B in phase A was used for 27 min to separate all main compounds. The flow rate was 0.3 mL/min, and the column was held at 30 • C. Mass spectra were acquired in positive-ion mode over a mass range from m/z 100 to 1500 with 5 Hz frequency. Ultrapure nitrogen was used as drying and nebulizer gas and argon was used as collision gas. The collision energy was set automatically from 15 to 75 eV, depending on the m/z of the fragmented ion. Operating parameters of the ESI ion source were as follows: capillary voltage 3 kV, dry gas flow 6 L/min, dry gas temperature 200 • C, nebulizer pressure 0.7 bar, collision radio frequency 700.0 V, transfer time 100.0 µs, and pre pulse storage 7.0 µs. Acquired data were calibrated internally with sodium formate introduced to the ion source at the beginning and end of each separation via a 20 µL loop. Processing of spectra was performed with Bruker DataAnalysis 4.3 software. UV DAD data were collected at 190-450 nm, and chromatograms acquired at UV 280 nm were smoothed by use of the Savitzky−Golay algorithm (window width 5 points, one iteration). Then, the detected peaks at retention time and the corresponding MS/MS were integrated. All determinations were done in triplicate (biological replication meaning a different extract n = 3); each biological replication was measured at least twice as a technical replication.
ProfileAnalysis software (version 2.1, Bruker Daltonik GmbH, Bremen, Germany) was used to preprocess the raw UHPLC-QTOF-MS data. ProfileAnalysis parameters were set as follows: advanced bucket generation with retention time range of 1.0-24.0 min, mass range of 100-800 m/z, normalization to sum of peaks, background subtraction, and time alignment; then, analyses were processed with the Find Molecular Futures (FMF) function to create compounds with a signal-to-noise threshold of 3 for peak detection. The generated bucket table consisting of m/z-retention time pairs and respective compound intensity was exported and uploaded to Metabolists 5.0 (http://www.MetaboAnalyst.ca/) for multivariate statistical analysis. The initial and final retention times were set for data collection. Data were imported into Metabolists 5.0 (http://www.MetaboAnalyst.ca/) online software to estimate missing values and to filter and normalize data (normalization by the median). No transformation was generalized, and the data matrix was meancentered and divided by the square root of the standard deviation of each variable (Pareto scaling). The obtained data matrix was introduced into SIMCA-P+ 16.01 (Sartorius AG, Göttingen, Germany) software for in-review multivariate statistical analysis of principal component analysis (PCA). The PCA score plot was used to present a natural correlation between the observations. To identify differential compounds, the PLS-DA (Partial Least Squares Discriminant Analysis) model was used to explore differences in depth between the profile metabolome of Salvia yangii and S. abrotanoides The PLS-DA model with VIP values (VIP ≥ 1.0) and |p(corr)| ≥ 0.5 was used to select differential compounds (PLS-DA data not shown).

Genomic DNA Extraction
DNA was extracted with the modified Doyle & Doyle method [53]. Spectrometric DNA concentration measurements were conducted with a NanoDrop ND 2000C (Thermo Fisher Scientific, Waltham, MA, USA). If needed, DNA isolates were purified with the Syngen DNA clean-up Kit (Syngen Biotech, Wroclaw, Poland) and the DNA concentration was measured again.

DNA Barcoding
ITS rDNA regions were amplified using primers ITS1F [54] and ITS4 [55] while the trnH-psbA chloroplast intergenic spacer was amplified with trnHf_05 [56] and psbA3_f [57] primers (Table 1). PCR was performed in a T-100 thermal cycler (BIO-RAD Laboratories, Hercules, CA, USA) in 20 µL volumes containing template DNA and was set up according to the Q5™ High-Fidelity DNA Polymerase manufacturer's protocol (Q5™ High-Fidelity DNA Polymerase, New England BioLabs Ltd., Hitchin, UK). Thermal cycling conditions were as follows: 95 • C for 30 s; 35 cycles of 95 • C for 10 s; 60 • C (for ITS1F&ITS4) or 62 • C (for trnHf_05&psbA3_f) for 30 s; 72 • C for 45 s and 72 • C for 2 min for final extension. PCR products were visualized by gel electrophoresis in 2% agarose gel in TEA buffer with 50 bp DNA Ladder (New England Biolabs Ltd., Hitchin, UK) and stained with SimplySafe (EURx, Gdansk, Poland). Amplified DNA products were purified with the Syngen Gel/PCR Mini Kit (Syngen Biotech, Wroclaw, Poland). Table 1. Primers used for amplification and Sanger sequencing of the ITS rDNA region and trnH-psbA intergenic spacer.

Primer Sequence (5 to 3 ) Reference
CGCGCATGGTGGATTCACAATCC [56] Sanger sequencing was carried out with the BrilliantDye™ Terminator v3.1 Kit (Nimagen B.V., Nijmegen, The Netherlands). Reactions were set with a fourfold dilution of the reaction premix and the addition BrilliantDye ® Terminator 5X Sequencing Buffer (Nimagen B.V., Nijmegen, The Netherlands) according to the producer's instructions. Sequencing products were precipitated with ethanol, dissolved in TSR (Hi-Di Formamide) (Thermo Fisher Scientific, Waltham, MA, USA) and then separated by capillary electrophoresis on an Applied Biosystems™ 310 Genetic Analyzer (Thermo Fisher Scientific, Waltham, MA, USA). Two reads were collected for each sample. Eight sequences were generated for molecular and phylogenetic analyses performed along with additional sequences obtained from NCBI (https://www.ncbi.nlm.nih.gov).
Sequence quality was checked with the Sanger Quality Check App (Thermo Fisher Scientific, Waltham, MA, USA). Forward and reverse sequencing reads for each marker region were assembled into a contig using the BioEdit Sequence Alignment Editor [58]. The identity of all obtained sequences was confirmed as ITS or trnH-psbA through similarity to published sequences using the BLAST algorithm. The ITS2 fragment was extracted with ITSx [59].
Alignments were assembled separately for each marker with the ClustalW algorithm under default conditions in MEGA 7 [60] and were then entered into the program Gblocks v.091b to objectively eliminate poorly aligned positions and divergent regions [61]. Then, a final alignment was assembled by a combination of both alignments calculated for individual markers. The relationship among the taxa was estimated using the Neighbor-Joining (NJ) method [62] and presented as an -rooted phylogram. Due to the large number of gaps in the dataset, a Pairwise deletion Gaps/Missing data treatment was applied. The evolutionary distances were computed using the Maximum Composite Likelihood method [63] according to the general recommendation. To estimate reliability of the tree bootstrap, a 1000-replication test was used [64]. The analysis involved 29 nucleotide sequence combinations including L. maculatum and L. amplexicaule as an outgroup.

Statistical Analysis of UHPLC-QTOF-MS Data
In order to reveal group patterns in non-volatile metabolic profiles of analyzed leaves and roots, an unsupervised PCA was applied. The PC loadings revealed a clear trend to form two clusters among all tested samples, either in leaves or roots, which represent the data in the space formed by the principal components (Figures 3 and 4).
In leaves, the first two principal components (PCs) explained 47.3% and 17.1% of the variation in the spectral data, respectively. On the right side of the score plots (green dots) a strong cluster of S. abrotanoides samples formed, while on the left side (blue dots) S. yangii samples were well separated (Figure 3a). According to the loadings plot ( Figure 3b) these two groups were significantly positively correlated with particular compounds detected in leaves. The PC loadings analyses confirmed that the content of sugiol (38) in S. yangii leaves and of carnosic acid quinone (31) in leaves of S. abrotanoides were mostly responsible for their chemical differentiation and revealed a clear separation of samples derived from these two species. In roots, the first principal component explained 65.3% of the variation in the spectral data and the second, 10.2%. Samples of S. yangii, represented by blue dots, were closely grouped on the right side of the plot, while on the left side (green dots), S. abrotanoides samples occurred (Figure 4a). The loadings plot revealed that these two groups were significantly positively correlated with particular compounds found in roots (Figure 4b).

Genetic Relationship Analysis
Sanger sequencing resulted in four rDNA ITS2 complete sequences (MT599312, MT599313, MT599314, MT599315) and four trnH-psbA intergenic spacer sequences (MT8158 71, MT815872, MT815873, MT815874) obtained for the cultivars S.OBRL.01, S.OBRL.03, S.OBRL.11 and S.OBRL.12, respectively (Table 4). All sequences were subject to BLAST analysis. Due to the lack of publicly available molecular data of the trnH-psbA barcode in subgenus Perovskia, BLAST analysis of the sequenced trnH-psbA region was not able to indicate similarity to any species from the Perovskia subgenus. Instead, the closest homolog of S.OBRL.01 and S.OBRL.11 trnH-psbA regions indicated by BLAST analysis was the trnH-psbA sequence of Salvia miltiorrhiza, while for S.OBRL.03 and S.OBRL.12, Salvia californica and Salvia chionopeplica were indicated as the closest homologs, respectively. BLAST analysis of the abundantly available data in the ITS2 region revealed high similarity between Salvia abrotanoides and Salvia yangii (Table 4).
NJ analysis was performed to estimate the relationship among the studied taxa (S. yangii, including wild type and cultivars "Lacey Blue" and "Blue Spire", and S. abrotanoides) in comparison with 23 related Salvia species (Table 5). Lamium maculatum and Lamium amplexicaule were applied as an outgroup. Finally, 50 additional sequences used for the NJ analysis were downloaded from GenBank (25 sequences of ITS2 and 25 sequences of trnH-psbA). The manually adjusted alignment of 29 combined ITS2 and 29 trnH-psbA sequences was 540 characters long with 361/540 variable sites, 301/540 Parsimony-informative sites.
According to the phylogram reconstructed based on ITS2 and trnH-psbA regions, the subgenus Perovskia was separated from other sages, yet clearly embedded in Salvia, ( Figure 5). The topology of the phylogram exhibiteds low interspecies genetic variation between subgenus Perovskia representatives indicated by short branches of the tree.
Despite the close genetic relationship, S. abrotanoides and S. yangii are suggested to be separate species. Moreover, intraspecies genetic variation within the S. yangii taxon, between "wild type" S. yangii and its cultivars ('Lacey Blue' and 'Blue Spire'), was observed. Furthermore, the intermediate location of S. yangii 'Blue Spire' between the two Salvia species suggests that the 'Blue Spire' cultivar might have been created through hybridization between S. abrotanoides and S. yangii during the cultivation process.

Discussion
In this work, we have provided a comprehensive analysis of the phytochemical profile of non-volatile compounds separately for roots and leaves of S. yangii and S. abrotanoides cultivated in Europe (Wroclaw, Poland). A reliable comparison was possible due to ensuring the same cultivation conditions such as: latitude, weather conditions etc., the same sampling method (immediate freezing in liquid nitrogen), and time (middle of the growing season) as well as the same extraction procedure (MeOH: H2O; 8:2 v/v) and method for chemical analysis (UHPLC-QTOF-MS). Despite ensuring identical cultivation conditions and methodology, a detailed comparison of the phytochemical non-volatile profile in roots and leaves of S. yangii and S. abrotanoides resulted in detecting several metabolites present in one species and absent in the other. The PC loadings analyses showed that sugiol in S. yangii leaves and carnosic acid quinone in leaves of S. abrotanoides were mostly responsible for their chemical differentiation. In roots, isograndifoliol, acetyloxycryptotanshinone, ketoisograndofoliol, oxocryptotanshinone, and grandifoliol in S. yangii as well as chrysoeriol 7-rutinoside, hesperidin, and one unidentified cryptotanshinone derivative in S. abrotanoides allowed for chemical distinguishing of root samples from analyzed species.
As noticed before, most of the published data from aerial parts of species belonging to subgenus Perovskia refer to the composition and activities of essential oils, indicating a significant variation [67]. A thorough review of the chemical profile of essential oils from subgenus Perovskia is provided by Mohammadhosseini et al. [22]. There are very few reports, however, of non-volatile constituents from aerial parts of species belonging to subgenus Perovskia, as a vast majority of them present an analysis of a whole-plant material.
Moreover, the influence of geographic origin and environmental factors influenced the phytochemical profile in both the volatile and non-volatile subset. The extent of such influences in individual cases is hard to assess without extensive experimental research. For example, a light quality modified metabolic profile of cultivated S. abrotanoides [68].
Tarawneh et al. [51] isolated four flavoinoids from dried leaves of S. yangii, which were not present in our isolates, except for the hydroxy-trimethoxyflavone, detected in leaves of both S. yangii and S. abrotanoides. Additionally, we detected trihydroxy-dimethoxyflavone, not found by Tarawneh et al. [51].
Amongst ten compounds found in methanolic extracts of air-dried aerial parts of S. abrotanoides [7], we detected rosmarinic acid and hesperidin, both present in the leaves of S. yangii and S. abrotanoides. Interestingly, in our study, hesperidin was also found to be present in the roots of S. abrotanoides. Two new tertracyclic diterpens, abrotandiol and abrotanone, found by Khaliq et al. [7] were absent in our extracts. The dihydroxy-methoxymethylflavone, however, found in our leaf isolates, might be the same compound as cismaritin found by Khaliq et al. [7]. Confirming that would require additional NMR analysis.
In a previous study of S. yangii roots (air-dried crude material, extraction with nhexane), eight diterpenoids were detected [44]. Extraction with 80% methanol performed in this study allowed for the detection of 18 diterpenoids in S. yangii roots and 15 in roots of S. abrotanoides. 1βOH-cryptotanshinone, cryptotanshinone, tanshinone IIat, and arucadiol (a synonym of miltiodiol) found by Slusarczyk at al. [44] were also detected in our methanolic root extracts in both S. yangii and S.abrotanoides. Miltirone, however, was present only in the methanolic root extract of S. yangii and not in S. abrotanoides, which corroborates previous results [25,44]. Two dehydrogenation products, didehydrotanshinone and didehydromiltirone, which were previously detected in n-hexane root extract of S. yangii, were not present in our methanolic root extract, neither in S. yangii nor in S. abrotanoides. Interestingly, isograndifoliol, reported recently in n-hexane root extract of S. yangii [66] was also detected in the methanolic root extract of S. yangii and was proved to contribute most to chemical differentiation between S. yangii and S. abrotanoides.
The phytochemistry and chemotaxonomic significance of metabolites found in subgenus Perovskia representatives was reviewed by Mohammadhosseini et al. [22].
In this study, we also aimed at analysis of the mutual genetic relationship between S. abrotanoides and S. yangii rather than their taxonomical position within the Salvia subgenera. To our knowledge, no previous work has done that with using DNA barcoding. Usually, a single-entry representative of Perovskia subgenus was used for phylogenetic tree reconstructions of Salvia [11][12][13][14]70]. There are two reports of employing molecular markers (ISSR markers) to investigate genetic variability among populations from subgenus Perovskia growing wild in Iran [71,72]. Heshemifar and Rahimmalek [72] analysed 63 accessions of S. abrotanoides from 16 populations and four accessions of S. yangii. Thirteen ISSR markers were shown to have a relatively high efficiency in the assessment of population diversity, which was in most cases confirmed by morphological classifications. In the study of Pourhosseini at al. [71], twelve S. abrotanoides populations were analyzed and a cultivated population of S. yangii was used as an outgroup. The authors have shown that nine ISSR markers are able to genetically distinguish between all studied populations and that the S. yangii population was clearly separated from S. abrotanoides populations, proving the effectiveness of this DNA fingerprinting method in studying the genetic relationship within the Perovskia subgenus. Despite the suitability of the ISSR markers for studying genetic variation among species within subgenus Perovskia, we have chosen DNA barcoding as it is regarded as a more sophisticated, standardized, and worldwide acknowledged method, recommended for genetic relationship analysis of land plants [73][74][75].
Although the DNA barcoding has been endorsed as the method of choice for land plant species identification, no consensus has been reached on which plant DNA sequences should be regarded as barcodes. The generally low rate of nucleotide substitution in plant mitochondrial genomes precludes the use of cytochrome oxidase 1 (CO1), a potential equivalent for the standard animal barcode, as a universal plant barcode [76,77]. Instead, the search for a plant barcode has involved looking outside the mitochondrial genome and from the outset, many researchers have accepted that multiple markers will be required to obtain adequate species discrimination [78]. Consequently, the Consortium for the Barcode of Life (CBOL) Plant Working Group has recommended a core-barcode consisting of portions of two plastid-coding regions, rbcL+matK, to be supplemented with additional markers as required [75]. Limitations of the rbcL+matK barcodes and idiosyncratic performance of different markers in different taxonomic groups were well discussed by Fazekas et al. [77]. Using additional markers was also strongly suggested by these authors. Beyond the core rbcL+matK barcode, the most widely used plastid barcoding marker is the intergenic spacer trnH-psbA. This region is straightforward to amplify across land plants, and is one of the more variable intergenic spacers in plants [79]. Another suggestion is supplementing the phylogenetic analysis of plants with the internal transcribed spacer from nuclear ribosomal DNA (nrITS) [80,81]. A 20% gain in discriminatory power for full ITS and 10-15% gain even when just the ITS2 region is used is an observed and appreciable benefit [82].
Following these recommendations, we have chosen four most widely used and routinely included plant DNA barcodes out of seven leading and 15 known plant sequences considered as barcodes [78]. These were as follows: two plastid coding regions rbcL and matK, plastid intergenic spacer trnH-psbA, and nuclear internal transcribed spacer ITS2.
In our approach of an extended marker framework for DNA barcodes, a combination of trnH-psbA and ITS2 barcodes showed sufficient discriminatory power to identify analysed plant organisms to the species level. The phylogram tree constructed with trnH-psbA and ITS2 barcodes has shown the existence of a separate Perovskia clade, which is clearly embedded within Salvia, which corroborates the findings of Drew et al. [14], who formally transferred former non-Salvia genera, such as Perovskia, to Salvia and created a subgeneric designation within Salvia for each of these former genera.
In the present form, our molecular analysis may suggest that S. abrotanoides and S. yangii are very closely related but regarded as two separate species. To make this conclusion fully justified, multiple accessions of S. abrotanoides would have to be included into the analysis. Unfortunately, at present, we do not have access to other taxonomically undoubtful accessions of S. abrotanoides; thus, we decided to rely on what is available. Despite a suggested close genetic relationship, S. abrotanoides and S. yangii possess distinct phytochemical profiles including qualitative differences within the non-volatile chemicals, which make them a potentially useful model for further analysis. Interestingly, the LC-MS profiles of aerial parts are apparently more similar, and none of the major compounds stand out as distinct, but the PCA processing revealed an obvious separation between the species. Hence, presence or absence of the minor compounds or their varied proportions is also important with such examples as the abovementioned sugiol/trilobinol pair. In contrast, the roots of both species are phytochemically distinct at first sight, and chemometric analysis has verified this fact. Isograndifoliol (Figure 2, peak 27), a hydroxylated norditerpenoid present only in S. yangii, is a rarely detected compound that has a high pharmacological potential [66], and its presence in a relatively high amount points out S. yangii as an especially valuable medicinal plant. Sugiol, even if a minor compound only, also has significant pharmacological potential, reviewed recently by Bajpai et al. 2020 [83], that can, for example, attenuate STAT3-dependent proliferation and unlocks apoptosis pathways in cancer cells and can block virus entry to cells. Whether or not this compound can contribute to the overall activity of S. yangii leaf as a drug has not been explored so far. On the other hand, carnosic acid quinone, a phytochemical marker of S. abrotanoides leaves, was reported as a highly potent factor in antioxidant reactions in S. rosmarinus and S. officinalis [84].
Yet, the phytochemical profiling alone is not sufficient to prove the taxonomic differentiation, especially if easily accessible organs (leaves) would be used for sampling. On the other hand, the results from molecular studies are independent from the sample type and coupled with parallel LC-MS profiling constitute a powerful tool in monitoring of taxonomically close and probably easily crossing species.
Qualitative chemical differences between closely related plant species might be of a different etiology. They might result directly from genetic (and/or epigenetic) differences between species or be an outcome of various processes occurring at the transcriptional and post-transcriptional levels. Comparative analysis of S. abrotanoides and S. yangii transcriptomes might help discover particular transcripts responsible for detected chemical differences. Performing genetic crosses between S. abrotanoides and S. yangii would create a segregating population which, through parallel phenotyping and genotyping, might serve for tracing genetic loci potentially responsible for the observed chemical profiles. A similar approach in Papaver somniferum allowed for discovering that ten genes from the noscapine biosynthesis pathway are clustered in the poppy genome [85]. In fact, there are numerous reports published to date proving that many secondary metabolites are a product of gene-clusters [86]. Once a certain target genetic locus is found, it could be subject to different biotechnological approaches, such as gene cloning, gene silencing or heterologous expression of proteins, in order to perform a functional analysis of target genes. In the long run, this may lead to developing applied research aiming at optimization of the chemical profile of Salvia, which may in turn become an alternative source of bioactive chemicals such as tanshinones.