Genome-Wide Diversity Analysis of Valeriana o ﬃ cinalis L. Using DArT-seq Derived SNP Markers

: Common valerian ( Valeriana o ﬃ cinalis L.) is one of the most important medicinal plants, with a mild sedative, nervine, antispasmodic and relaxant e ﬀ ect. Despite a substantial number of studies on this species, the genetic diversity and population structure have not yet been analyzed. Here, we use a next-generation sequencing-based Diversity Array Technology sequencing (DArT-seq) technique to analyze Polish gene bank accessions that originated from wild populations and cultivars. The major and, also, the most astounding result of our work is the low level of observed heterozygosity of individual plants from natural populations, despite the fact that the species is widespread in the studied area. Inbreeding in naturally outcrossing species such as valerian decreases reproductive success. The analysis of the population structure showed the potential presence of a metapopulation in the central part of Poland and the formation of a distinct gene pool in the Bieszczady Mountains. The results also indicate the presence of the cultivated gene pool within wild populations in the region where the species is cultivated for the needs of the pharmaceutical industry, and this could lead to structural and genetic imbalances in wild populations. Bieszczady Mountains in the Podkarpackie Voivodeship. The second gene pool consisted of individuals from the rest of the country. The secondary population structure, for K = 3, indicated that the third gene pool was separated from the second one. Eight out of thirteen individuals from the populations that originated from the Mazowieckie Voivodeship had a share coefficient of over 80% of the third gene pool. In addition, its negligible presence was observed among individuals from the Kujawsko-Pomorskie and Ś wi ę tokrzyskie Voivodeships. The structure of the lower order was identified at the K = 5 level. The additional admixture, i.e., the fourth gene pool, has been shown for individuals of the accession PL 403183 (6) collected


Introduction
Common valerian (Valeriana officinalis L.) has a long tradition of use as a medicinal plant [1]. Valerian root is one of the best-tested plant materials in terms of pharmacology and clinical aspects. The raw material contains compounds from several chemical groups, including essential oil (according to European Pharmacopoeia [2], not less than 0.5%), valepotriates and sesquiterpenic acids (according to Eur. Ph. not less than 0.17%), with valerenic and acetoxyvalerenic acids considered as species-specific marker compounds. Valerian root and the extracts thereof are used in the states of anxiety or sleeping disorders. Their efficacy relates to an interplay between distinct groups of compounds rather than with individual compounds [3].
V. officinalis belongs to the genus Valeriana L. that is a part of the family Valerianaceae, now considered also as a part of the Caprifoliaceae Juss., in the order Dipsacales Juss. ex Bercht. & J.Presl [4]. Approximately 350 species are listed within this cosmopolitan genus, several of which are of commercial interest [5]. V. officinalis occurs natively throughout Eurasia, except for arctic and desert zones [6]. It occupies fresh  Table 1.   Table 1.

DNA Extraction and Genetic Analysis
For DNA analysis, 10 seedlings from wild populations and 20 from the cultivar were selected at random. Young, healthy leaves were collected separately from each plant and dried in the presence of active silica gel. Until DNA isolation, the dried material was stored in −18 • C in vacuum-sealed packages. The tissue was crushed into a fine powder in a bead mill (Mixer Mill MM 200, Retsch, Haan, Germany). Isolation of DNA was conducted using the Genomic Mini AX Plant kit (A&A Biotechnology, Gdynia, Poland) according to the manufacturer's protocol. DNA integrity was evaluated visually after electrophoresis in 2% agarose gel. Purity and concentration were evaluated spectrophotometrically (NanoDrop ND-1000, NanoDrop Technologies, Willmington, DA, USA). A set of 188 isolates with appropriate qualitative and quantitative parameters was sent for DArT-seq analysis to Diversity Arrays Technology lab, University of Canberra, Australia. Reduction of the genome complexity was performed using two restriction enzymes, i.e., MseI and PstI.

Data Analysis
The results of genotyping were generated as a table listing SNP (codominant) markers that are single nucleotide polymorphisms detected in the sequenced fragments of genome representations. SNP markers were then scored as binary data, indicating the presence (1) or absence (0) of a marker in the genomic representation of each sample, as described by Cruz et al. [16] and, in the next step, were tested for reproducibility (%), call rate (%) and polymorphism information content (PIC). The following coefficients were calculated: observed heterozygosity (H o ), Nei unbiased genetic diversity (uH e ) and inbreeding coefficient (F IS ). The following equations were used: where n is the sample size, and p i is the frequency of the ith trait/marker.
Genetic distance was measured by using the Jaccard coefficient [21] according to the formula: where a is the number of variables present in both individuals being compared, b is the number of variables present in the first individual and absent from the second individual and c is the number of variables present in the second individual and absent from the first one. The symmetric genetic distance matrix was used to assess the amount of variation among the assigned regional groupings by Analysis of Molecular Variance (AMOVA) with tests including 999 permutations. In order to test for a correlation between genetic and geographical distances (in km) among populations, Mantel tests were performed (computing 10,000 permutations). Genetic structure was assessed using Principal Coordinate Analysis (PCoA) and model-based Bayesian clustering algorithms. Admixture and a correlated allele frequencies model was used to determine the number of clusters (K) in the range from one to ten. Ten independent runs were set for each number of clusters (K value), with a burn-in period of 100,000 and 500,000 Markov Chain Monte Carlo replications after burn-in with no prior information about the origin of individuals. The ∆K method was used to determine the most probable value of K [22]. All above mentioned analyses were performed using the Microsoft Excel 2016 (Microsoft, Redmond, WA, USA), R packages (vegan, dplyr and ape); XLSTAT Ecology (Addinsoft, Inc., Brooklyn, NY, USA); GenAlEx 6.501 [23]; STRUCTURE v2.3.4 [24] and CLUMPAK [25].

Marker Quality Analysis
A total of 65,510 SNP markers with an average reproducibility of 99.6% were obtained using the DArT-seq method for 188 individuals of V. officinalis. Over 98% of SNPs had reproducibility above 95%, of which 55,101 were completely reproducible (Figure 2a). The SNP reproducibility is defined as the proportion of successful replicates among all repetition attempts. The call rate was in the range from 20% to 100% (Figure 2b). Approximately 63% of SNPs had a call rate above 75%. The call rate is defined as the proportion of individuals in the study for which SNP information in particular loci is available. When the quality criteria were applied, i.e., reproducibility above 95% and call rate 100%, 23,507 SNP markers were assigned to a subsequent analysis. The PIC was in the range of 0.01-0.5, with an average of 0.181 and a median of 0.129. As much as 40% of SNP markers had a PIC below 0.1, and only 1% had a maximum value of 0.5 ( Figure 2c). A total of 65,510 SNP markers with an average reproducibility of 99.6% were obtained using the DArT-seq method for 188 individuals of V. officinalis. Over 98% of SNPs had reproducibility above 95%, of which 55,101 were completely reproducible ( Figure 2a). The SNP reproducibility is defined as the proportion of successful replicates among all repetition attempts. The call rate was in the range from 20% to 100% (Figure 2b). Approximately 63% of SNPs had a call rate above 75%. The call rate is defined as the proportion of individuals in the study for which SNP information in particular loci is available. When the quality criteria were applied, i.e., reproducibility above 95% and call rate 100%, 23,507 SNP markers were assigned to a subsequent analysis. The PIC was in the range of 0.01-0.5, with an average of 0.181 and a median of 0.129. As much as 40% of SNP markers had a PIC below 0.1, and only 1% had a maximum value of 0.5 ( Figure 2c).

Genetic Diversity
The percentage of polymorphic loci ranged from 31% to 61% (Figure 3a). The comparison of wild populations and cultivar regarding the occurrence of unique alleles showed that their number ranged from 0 to 1059 (Figure 3b). No unique alleles were found in two populations (PL 403186 (8) and PL 403187 (9) ) originating from the Świętokrzyskie Voivodeship. Twenty-two percent of alleles that occurred in wild populations were not represented in the cultivar. The contribution of SNP markers unique for the regions is shown in Figure 3c. The observed heterozygosity (Ho) ranged from 0.196 (PL 401941 (12) ) to 0.354 PL 401951 (20) ) ( Figure 3d). The mean heterozygosity of the wild populations was 0.250. Individual specimens in accessions PL 406738 (18) and PL 403192 (19) showed significantly higher Ho levels. It was comparable to the Ho observed to the cultivar, which is tetraploid. The expected heterozygosity (uHe) varied between 0.274 (PL 403192 (19) ) and 0.348 (PL 403183 (6) ), with an average for the wild populations of 0.315 ( Figure 3d). The lowest value of the inbreeding coefficient (FIS) was found for cultivar PL 401951 (20) (-0.104), while the highest for accession PL 401941 (12) (0.354) (Figure 3d). In general, in wild populations, this coefficient reached values above zero, which

Genetic Diversity
The percentage of polymorphic loci ranged from 31% to 61% (Figure 3a). The comparison of wild populations and cultivar regarding the occurrence of unique alleles showed that their number ranged from 0 to 1059 (Figure 3b). No unique alleles were found in two populations (PL 403186 (8) and PL 403187 (9) ) originating from theŚwiętokrzyskie Voivodeship. Twenty-two percent of alleles that occurred in wild populations were not represented in the cultivar. The contribution of SNP markers unique for the regions is shown in Figure 3c. The observed heterozygosity (H o ) ranged from 0.196 (PL 401941 (12) ) to 0.354 PL 401951 (20) ) ( Figure 3d). The mean heterozygosity of the wild populations was 0.250. Individual specimens in accessions PL 406738 (18) and PL 403192 (19) showed significantly higher H o levels. It was comparable to the H o observed to the cultivar, which is tetraploid. The expected heterozygosity (uH e ) varied between 0.274 (PL 403192 (19) ) and 0.348 (PL 403183 (6) ), with an average for the wild populations of 0.315 (Figure 3d). The lowest value of the inbreeding coefficient (F IS ) was found for cultivar PL 401951 (20) (-0.104), while the highest for accession PL 401941 (12) (0.354) (Figure 3d). In general, in wild populations, this coefficient reached values above zero, which indicated inbreeding. The accession PL 401945 (15) from theŚwiętokrzyskie Voivodeship was the closest to the state of Hardy-Weinberg equilibrium. A negative F-value for the cultivar indicates excessive heterozygosity. The genetic distance, according to Jaccard's coefficient, between the examined individuals ranged from 0.08 (individuals from accessions PL 406738 (18) and PL 403192 (19) ) to 0.380 (individuals from accessions PL 401941 (12) and PL 406738 (18) ). The mean distance between individuals from natural populations was 0.239.
Agronomy 2020, 10, x FOR PEER REVIEW 6 of 15 indicated inbreeding. The accession PL 401945 (15) from the Świętokrzyskie Voivodeship was the closest to the state of Hardy-Weinberg equilibrium. A negative F-value for the cultivar indicates excessive heterozygosity. The genetic distance, according to Jaccard's coefficient, between the examined individuals ranged from 0.08 (individuals from accessions PL 406738 (18) and PL 403192 (19) ) to 0.380 (individuals from accessions PL 401941 (12) and PL 406738 (18) ). The mean distance between individuals from natural populations was 0.239.

Mantel Test
The Mantel test, in which genetic distance matrices between populations with a matrix of geographical distance between them, were collated and showed a moderate positive relationship (0.304; p < 0.0001).

AMOVA
The AMOVA analysis for SNP markers of all samples showed that 23% of the variance in the allele frequency could be attributed to differences between populations. In the cases when only wild populations were analyzed, the value of the population differentiation (FPT) coefficient was 0.199 (p < 0.001). The hierarchical AMOVA of wild populations revealed that 78% of the total variance was found among individuals and 16% was found among populations, while the rest of the variation (6%) was among regions. The genetic variation between the cultivar and wild populations was 29% ( Figure 4).

Mantel Test
The Mantel test, in which genetic distance matrices between populations with a matrix of geographical distance between them, were collated and showed a moderate positive relationship (0.304; p < 0.0001).

AMOVA
The AMOVA analysis for SNP markers of all samples showed that 23% of the variance in the allele frequency could be attributed to differences between populations. In the cases when only wild populations were analyzed, the value of the population differentiation (ɸPT) coefficient was 0.199 (p < 0.001). The hierarchical AMOVA of wild populations revealed that 78% of the total variance was found among individuals and 16% was found among populations, while the rest of the variation (6%) was among regions. The genetic variation between the cultivar and wild populations was 29% ( Figure  4).

Population Structure
The first three main coordinates of PCoA explained only 24.5% of the total variability ( Figure 5). In the first and second coordinate systems, four groups-i.e., the cultivar, the population originating from the Podkarpackie region, with a separate population of PL 403183 (6) , and a large group comprising populations from the other three Voivodeships-were clearly visible. The third coordinate allowed to isolate, as well, a group of individuals from the Mazowieckie Voivodeship. In the group formed by the cultivar, there were also five specimens from the Mazowieckie Voivodeship. These were specimens with a higher level of Ho, which seems to confirm that they were tetraploid.
In order to further elaborate the genetic structure of the populations, a model-based Bayesian clustering was conducted using Structure 2.3.4 software [24]. The analysis was run separately for diand tetraploids. To find the suitable value of K, the number of clusters (K) was tested in the range from 1 to 10 and was plotted against ΔK. Genotypes that the estimated proportion of membership (Q) scored > 0.80 were considered as pure, while < 0.80 as admixed. A sharp peak for K = 2 ( Figure  6a) indicated that two distinct pools were found to contribute significant genetic information across tested diploid accessions. However, lower-order structures were also recorded for K = 3 and K = 5. For tetraploid individuals, a clear peak was present for K = 2. (Figure 6b). The proportion of membership of individuals in each pool was illustrated in the bar plot of the population assignment test in the structure analysis (Figure 6c,d). Among the diploids, two main pools formed 30.4% and 66.7% of the individuals, respectively. The remaining 2.9% of the individuals were classified as admixed. To the first gene pool were assigned only individuals of accessions that originated from the

Population Structure
The first three main coordinates of PCoA explained only 24.5% of the total variability ( Figure 5). In the first and second coordinate systems, four groups-i.e., the cultivar, the population originating from the Podkarpackie region, with a separate population of PL 403183 (6) , and a large group comprising populations from the other three Voivodeships-were clearly visible. The third coordinate allowed to isolate, as well, a group of individuals from the Mazowieckie Voivodeship. In the group formed by the cultivar, there were also five specimens from the Mazowieckie Voivodeship. These were specimens with a higher level of H o , which seems to confirm that they were tetraploid.
In order to further elaborate the genetic structure of the populations, a model-based Bayesian clustering was conducted using Structure 2.3.4 software [24]. The analysis was run separately for diand tetraploids. To find the suitable value of K, the number of clusters (K) was tested in the range from 1 to 10 and was plotted against ∆K. Genotypes that the estimated proportion of membership (Q) scored > 0.80 were considered as pure, while < 0.80 as admixed. A sharp peak for K = 2 (Figure 6a) indicated that two distinct pools were found to contribute significant genetic information across tested diploid accessions. However, lower-order structures were also recorded for K = 3 and K = 5. For tetraploid individuals, a clear peak was present for K = 2. (Figure 6b). The proportion of membership of individuals in each pool was illustrated in the bar plot of the population assignment test in the structure analysis (Figure 6c,d). Among the diploids, two main pools formed 30.4% and 66.7% of the individuals, respectively. The remaining 2.9% of the individuals were classified as admixed. To the first gene pool were assigned only individuals of accessions that originated from the Bieszczady Mountains in the Podkarpackie Voivodeship. The second gene pool consisted of individuals from the rest of the country. The secondary population structure, for K = 3, indicated that the third gene pool was separated from the second one. Eight out of thirteen individuals from the populations that originated from the Mazowieckie Voivodeship had a share coefficient of over 80% of the third gene pool. In addition, its negligible presence was observed among individuals from the Kujawsko-Pomorskie andŚwiętokrzyskie Voivodeships. The structure of the lower order was identified at the K = 5 level. The additional admixture, i.e., the fourth gene pool, has been shown for individuals of the accession PL 403183 (6) collected in the Bieszczady Mountains. The presence of the fifth gene pool was recorded among the accessions that originated from theŚwiętokrzyskie Voivodeship. It dominated in three accessions, i.e., PL 403186 (8) , PL 403187 (9) and PL 401938 (10) , while, in the remaining ones, it was an admixture. For tetraploids, two gene pools associated with biological status, i.e., cultivar and individuals from wild populations, were clearly visible (Figure 6d). Bieszczady Mountains in the Podkarpackie Voivodeship. The second gene pool consisted of individuals from the rest of the country. The secondary population structure, for K = 3, indicated that the third gene pool was separated from the second one. Eight out of thirteen individuals from the populations that originated from the Mazowieckie Voivodeship had a share coefficient of over 80% of the third gene pool. In addition, its negligible presence was observed among individuals from the Kujawsko-Pomorskie and Świętokrzyskie Voivodeships. The structure of the lower order was identified at the K = 5 level. The additional admixture, i.e., the fourth gene pool, has been shown for individuals of the accession PL 403183 (6) collected in the Bieszczady Mountains. The presence of the fifth gene pool was recorded among the accessions that originated from the Świętokrzyskie Voivodeship. It dominated in three accessions, i.e., PL 403186 (8) , PL 403187 (9) and PL 401938 (10) , while, in the remaining ones, it was an admixture. For tetraploids, two gene pools associated with biological status, i.e., cultivar and individuals from wild populations, were clearly visible (Figure 6d).   (c) the results of 100,000 iterations of Bayesian clustering algorithms using an admixture and correlated allele frequencies model by STRUCTURE software [24] with K = 2, K = 3 and K = 5, where K is the number of groups assumed for diploids, and (d) the results of 100,000 iterations of Bayesian clustering algorithms using an admixture and correlated allele frequencies model by STRUCTURE software with K = 2, where K is the number of groups assumed for tetraploids. Each vertical bar represents a single plant. The length of the colored segment shows the estimated proportion of membership of that sample to each group. Accession numbering is in accordance with Table 1.

Discussion
V. officinalis is a species commonly used in herbal medicine. It is widespread in Europe and Asia and is currently not considered being endangered. Previous studies concerned mainly the content and specification of active compounds and their use for therapeutic purposes [1,[26][27][28][29]. Several studies using molecular markers, such as amplified fragment length polymorphism (AFLP) or random amplification of polymorphic DNA (RAPD), were performed for this species, but none of them concerned the intraspecies diversity analysis [30][31][32]. In this context, Valeriana wallichii DC (syn. Valeriana jatamansi Jones) is much better recognized [33][34][35][36][37][38]. The NGS high-throughput sequencing technology, i.e., RADseq, has already been used in Valerianaceae phylogenetic research, but in the research presented in this paper, only the species native to South America were concerned, and V. officinalis was not included in it [39].

Utility of DArT-seq Markers
To the best of our knowledge, in the research presented in this paper, the NGS technology was used for the first time for analyzing intraspecies variations within the Valeriana genus. Until now, the used here DArT-seq technique has been mainly applied for the analyses of crop and animal species [40][41][42][43][44][45]. For 188 individuals of V. officinalis representing 19 accessions that originated from natural populations and one cultivar using the DArT-seq technique, over 65,000 SNPs with very high reproducibility coefficients were obtained. For comparison, in Triticum durum Desf., a DArT-seq analysis allowed to identify over than 20,000 SNPs [46]. In phylogenetic studies of the genus Secale, this technique allowed to identify about 14,000 SNPs, and in studies of rye inbreeding lines, almost 5000 SNP markers were obtained [47,48]. However, it should be noted that the genome of V. officinalis is over two point five times smaller than the genome of Secale cerale L. and more than four times smaller than the genome of T. durum [30,49,50]. Thus, a considerably larger part of the genome was covered here by the SNP analysis. In the population genomics of Eucalyptus species, 54,000 SNP markers were recorded [17]. In earlier studies of V. officinalis, using AFLP markers, only 311 fragments were generated [30]. The PIC value obtained in the presented study was quite low and amounted to 0.181. It was a derivative of a significant fraction of monomorphic loci. In the above-mentioned cereal studies, the PIC value was significantly higher and amounted to 0.302 in wheat and 0.37 in rye for the SNPs [17,48]. On the other hand, similar (~0.18) PIC values were obtained in the study of Ctenophorus caudicinctus, Litoria ewingii and Litoria paraewingi (Australian lizard and frogs) [45]. Based on the above data, therefore, it can be concluded that the DArT-seq technique is efficient and useful in the genetic analysis of V. officinalis.

Diversity
Due to the lack of earlier studies on genetic variability in the V. officinalis species, it is difficult to refer the level of parameters obtained in the presented paper to species-specific diversity patterns. A substantial fraction of monomorphic loci or such with a very low level of polymorphism and a negligible proportion of unique alleles/loci shows a very homogenous genetic background of the examined natural populations. The low level of observed heterozygosity and positive values of the inbreeding coefficients indicate a prevalence of mating between closely related individuals. Considering that V. officinalis produces rhizome and is capable of clonal propagation, there is a high probability that adjacent individuals at natural sites are, in fact, ramets clusters, i.e., a group of clones formed as a result of the vegetative proliferation of a genet (mother plant). They can still be linked to the mother plant or exist separately as a result of the division of the genet [51]. Moreover, according to the research of Konon and Novikova [52], in valerian, the mechanism of self-incompatibility appears to be absent. Thus, pollination with pollen of the same genotype is possible, although fertilization with foreign pollen is promoted by the delayed development of the stigma in relation to the stamens [32]. This is also confirmed by the lack of correlation between the level of variation and the abundance of V. officinalis in the natural sites. Moreover, in the cultivar, the level of genetic variability is significantly higher than in natural populations, and the inbreeding coefficient indicates the excess of heterozygotes, even though the cultivar was developed as a result of targeted selection. However, this selection was conducted in terms of the content of active compounds. Therefore, individuals representing a similar phenotype and a different genotype were probably selected. Their free crossbreeding led to an increase in heterozygosity. The higher variability of the cultivar results also from the farming practice and pollination biology of this species. Valerian has a high outcrossing rate of 76.5% to 97.7% under field conditions, as shown by studies of Penzkofer et al. [32]. Under the cultivation conditions for herbal raw material, plantations are established from the seedlings, and the plants are unearthed in the second year of growing. The result is that no clonal propagation takes place on plantations, i.e., a single plant is a single genet, and neighboring plants will probably have a different genetic makeup because of mixing of seeds before sowing and the random planting of seedlings. It is also important that the cultivar is tetraploid in contrast to accessions originated from the natural sites. Thus, the observed heterozygosity may result from the presence of different alleles on homologous chromosomes.

Gene Pool
Despite, to a large extent, the homogeneous genetic background of natural populations, based on the analysis of the population structure, the differentiation of gene pools typical for particular regions of the country was clearly visible. It is most noticeable in the case of Podkarpackie Voivodeship, where there occurs a spatial isolation due to the presence of the Bieszczady Mountains. However, the structure of the population also indicates that there is no complete isolation of the natural populations, since, in each of the regions, some admixture from the adjacent gene pool(s) appears, which is also indicated by the fact that the correlation between the genetic and geographical distance was rather limited. The presence of the additional gene pool in different degrees in populations originating from Podkarpackie andŚwiętokrzyskie Voivodeships from which the wild populations of V. officinalis originated may indicate that an additional gene pool may be present in Southern Poland. This may be facilitated by the topography of that part of the country, i.e., the presence of highlands and high mountains (Tatra Mountains) with a more severe climate. The presence of V. officinalis on natural sites in that region has been proven [6,51]. From the point of view of the genetic structure of the population, the similarity of genotypes occurring in theŚwiętokrzyskie and Kujawsko-Pomorskie Voivodeships also seems to be interesting. This may indicate that a metapopulation with free gene flow occurs in a significant part of Poland. This hypothesis is supported by the widespread occurrence of valerian across the country; however, its verification will be possible only after investigating a larger number of natural populations, especially from regions not included in this study [6]. In the populations inhabiting the Mazowieckie Voivodeship, the presence of individuals that are potential tetraploids were detected. Based on the morphological features, which show high variability, additionally affected by environmental factors, it is impossible to distinguish between particular cytotypes [10]. The analysis of the population structure revealed the presence of a gene pool typical for the cultivar within individuals representing a wild population. This is probably due to the presence of valerian plantations in the area. In other words, the uncontrolled spread of the cultivar beyond the plantation area and crossbreeding of the wild and cultivated specimens could be possible. However, the tetraploid individuals could also represent a separate cytotype that occurs natively in the area, although with a much lower abundance. In order to verify which of the above hypotheses is true, a more detailed study of the populations occurring in this region would be necessary. These studies, in addition to the genetic fingerprint and cytotype, should also determine the chemotype, because the cultivated forms have a distinct chemical profile and contain more active compounds. Valerian cultivation is very widespread and common in the Lubelskie Voivodeship, so, in that area, the probability of contamination of natural populations by the cultivars is extremely high. Unfortunately, no environmental monitoring is carried out in this context.

Vulnerabilities
Although V. officinalis is still considered to be a commonly occurring species not only in Poland but, also, throughout its entire range of distribution, due to the loss and degradation of habitats, the abundance of this species may start to decline rapidly. Climate change and increasing periods of drought are not without significance here [53]. Considering the low level of heterozygosity and definitely lower than assumed level of genetic variation, in a very short time, this species may face the threshold of extinction, as inbreeding decreases the reproductive capacity of outbreeding species [54]. An additional risk that is not considered is the uncontrolled release into the environment of cultivars. Since valerian is a native species, no environmental monitoring is carried out in this area, especially since, without the use of molecular techniques, it is basically not possible. However, it should be expected that cultivars developed by selection and not representing a complete genetic makeup of the species will have an additional negative impact on the size of the gene pool. The cultivation of medicinal plant species, on the one hand, enhances the protection of the wild populations that are not drained to obtain raw herbal resources and, on the other hand, may cause a disturbance of the structure and genetic balance of the wild populations through the uncontrolled spread of individuals of cultivars beyond the plantation area. Based on the results presented here, it is impossible to assess to what extent this threat is real and what long-term effects it will have. However, there is no doubt about the existence of this problem, and this issue should be explored further in the future.

Conclusions
This study showed that DArT-seq technology can be applied effectively in genetic studies of V. officinalis. The analysis of genetic diversity showed that the wild populations had highly similar genetic makeups, and the individuals displayed high levels of homozygosity. These may threaten the species, especially in the context of climate change. Individuals representing a gene pool typical for the cultivar were found in the natural environment. Further detailed research is necessary to determine their origin, and their impact on the wild populations should be monitored.