Genetic diversity and structure analysis of Croatian garlic collection assessed by SSR markers

Garlic (Allium sativum L.) is a valuable crop and one of the oldest-known horticultural species. It has been used as a condiment, food and medicinal plant since ancient times (Ipek et al., 2005). The production of garlic in 2018 reached 28,494,130 metric tons, with China being the largest producer and accounting for almost 80% of the world production (FAOSTAT, 2018). Its centre of origin is considered to be in Central Asia (Kyrgyzstan, Tajikistan, Turkmenistan and Uzbekistan) (Hong and Etoh, 1996) with Allium longicuspis proposed to be its wild ancestor (Egea et al., 2017). It belongs to the genus Allium whose recent classification is mostly based on nuclear ribosomal DNA sequences involving the internal transcribed spacer region and includes >800 species in 15 subgenera and >70 sections (Li et al., 2010; Friesen et al., 2020). Two regions with broad species diversity are from the Mediterranean basin to Central Asia and North-West of America (Fritsch and Friesen, 2002; Li et al., 2010).

Garlic is diploid (2n = 16) and asexually (vegetative) propagated species, primarily from cloves due to obligatory apomixis (Bradley et al., 1996), although a few fertile types had been reported (Etoh et al., 1988). Therefore, clonal selection is the method of choice for garlic breeding. Sexual reproduction and fertility are impaired among commercial cultivars (García-Lampasona et al., 2003), relying solely on mutations (spontaneous or induced), somaclonal variation and genetic transformation. Vegetative reproduction usually assures crop uniformity, but in garlic, there is a great degree of diversity in phenotypic traits, as well as in some agronomic traits such as drought tolerance and other stress responses (Bradley et al., 1996; Stavělíková et al., 2008). Large genetic diversity has been observed in garlic clones that are cultivated around the world and also in collected non-cultivated garlic accessions found in nature (Maass and Klaas, 1995; Bradley et al., 1996; Volk et al., 2004). Classification of cultivars based on morphological traits is difficult due to a wide variation in clove number, colour and size within the same cultivar, varying further phenotypically under different agro-climatic conditions, despite the lack of sexual reproduction.

Conservation of garlic genetic resources is important as such collections present a source of traits with potential economic value. A large number of garlic genetic resources are maintained in gene banks, which are a source of variability and present an opportunity for genetic research and breeding programs. Two main strategies for preserving garlic collections are field maintenance and cryopreservation, each of them having their benefits and drawbacks, but for long-term preservation cryopreservation generally becomes less costly than field culture after 8–13 years (Keller et al., 2013). Maintenance of field collections is expensive and laborious and demands equipment and consumables. Hence, characterisation and evaluation of accessions require a combination of agro-morphological traits and molecular markers to help decision-making in minimising the number of accessions and maximising the genetic diversity. Molecular studies of garlic are mainly focused on the analysis of genetic diversity of cultivated garlic clones, mostly studied through the use of several molecular marker techniques, such as random amplified polymorphic DNA (RAPD) (Maass and Klaas, 1995; Paredes et al., 2008), amplified fragment length polymorphism (AFLP) (García-Lampasona et al., 2003; Volk et al., 2004; Ipek et al., 2005, 2006) and isozymes (Ipek et al., 2003). These techniques are preferential in discriminating garlic taxonomy because they are not affected by environmental conditions, which is an issue with morphological traits.

Simple sequence repeats (SSRs) or microsatellites have been extensively used in the analysis of genetic diversity and population genetic studies on many plant species, including garlic (Ma et al., 2009; Wang et al., 2010; Zhao et al., 2011; Jo et al., 2012; Cunha et al., 2014; Ipek et al., 2015; Barboza et al., 2018; Kumar et al., 2019). Besides next-generation sequencing techniques and high-throughput genotyping-by-sequencing (GBS) approaches, SSR markers are still preferred and intensively used due to genome abundance, high polymorphism, co-dominance, high reproducibility, low cost and transferability in laboratories (Wang et al., 2010). Barboza et al. (2018) developed 110 SSR markers from expressed sequence tags (EST) and summarised in total 312 previously reported genomic and EST SSR markers in garlic species (Ma et al., 2009; Lee et al., 2011; Cunha et al., 2012; Chand et al., 2015; Ipek et al., 2015; Liu et al., 2015). However, their usage has been limited by the large (15.9 Gbp/1C) and duplicated genome that result, in some cases, in more than two alleles per SSR locus (Ipek et al., 2008; Ma et al., 2009; Egea et al., 2017; Peška et al., 2019). In the Adriatic region, garlic has been grown since ancient times (Dumičić et al., 2013). Six Croatian garlic landraces (traditional cultivars or varieties) – Brgudski ozimi sinonim Benkovački, Cerićki ozimi, Istarski crveni, Podravski ozimi, Slavonski ozimi and Šokac – are registered in the EU Plant Variety Database (v.3.2.1).

The objective of this study was to assess the genetic diversity and structure of the Croatian garlic collection comprising 64 accessions from 3 Croatian regions (Istria, Dalmatia and continental Croatia) and 16 foreign landraces using SSR markers. To our knowledge, this is the first attempt to extensively evaluate garlic accessions from Croatia. Obtained results would give better insight on garlic genetic resources, and thus support the effective and sustainable management of Croatian collections.

Thirteen SSR markers were selected for this study (Table 1). We found 1–5 alleles per locus for each accession. Overall, 71 alleles were detected with a mean of 5.46 alleles per locus. The number of alleles per locus ranged from two at AS623 to eight at GB-ASM-053 (Table 1).

Among the 80 accessions evaluated, we identified 61 distinct MLGs, of which 51 represented unique genotypes. The remaining accessions were divided in 10 MLG groups, comprising potential duplicates or redundant genotypes. Group size ranged from two to seven accessions per group.

Our molecular diversity study of garlic accessions revealed a total of 48 distinct genotypes in the group of 64 Croatian accessions, while all of the 16 foreign accessions were distinct genotypes (Table 2). Croatian and foreign sets shared three genotypes (Figure 1; accessions P1S02/P4S75, P3S40/P4S70 and P4S78/P1S01/P3S47/P3S48/P3S52/P3S57/P3S59). Within the group of Croatian accessions, the highest number of distinct genotypes (19) was observed in the group originating from Dalmatia. The lowest number of distinct genotypes (13) was found in accessions from continental Croatia. However, the values of genotypic richness, i.e. the number of genotypes independent of sample size was relatively similar across regions, ranging from 12.66 to 13.61. The Croatian region of Istria had the highest share of private genotypes, as well as the highest private genotypic richness (Table 2), while the region of continental Croatia showed the highest Shannon’s information index.

Figure 1

A total of 48 distinct genotypes were identified in the Croatian set of accessions. The genotypes within groups (from one or more Croatian regions, but not among foreign accessions) were counted as one distinct genotype; e.g. accessions P1S02 (MLG G06) and P3S40 (MLG G09) are not considered as distinct genotypes (Table S1 in Supplementary Materials). The number of private genotypes attributed to the Croatian regions refers to all genotypes, except for foreign accessions (Table S1 in Supplementary Materials). The higher percentage of polymorphic loci and Shannon’s index of the foreign set of accessions reflected the geographically distant origins.

UPGMA tree based on Dice’s distance matrix was constructed to assess the relationships among accessions (Figure 1). Five main clusters (A–E) and two separate lineages (P2S30 and P4S77-SYR) were identified. Clusters B, C and E were supported with bootstrap values >50%.

The assessment of the genetic structure was performed on 61 distinct MLGs using STRUCTURE and BAPS softwares. The most probable number of clusters (K) was inferred by the ad hoc statistic ΔK based on the rate of change in the log probability of data between successive K values (Evanno et al., 2005). The highest ΔK value was detected at K = 3 (ΔK = 391.14) followed by K = 4 (ΔK = 21.18) and K = 5 (ΔK = 9.37), whereas the highest average estimate of the likelihood of the data, ln[Pr(X|K)], was obtained at K = 5 (−944.85) (Figure S1 in Supplementary Materials).

As expected, the results obtained by STRUCTURE (Figure 1B) were congruent with the UPGMA tree (Figure 1A). At K = 3 clusters A and B corresponded to the UPGMA clustering, whereas cluster C included the UPGMA clusters C–E. At K = 4, cluster C splitted into UPGMA clusters C and D, and finally, at K = 5, cluster D splitted into UPGMA clusters D and E. Therefore, STRUCTURE and UPGMA clustering were in accordance. Two separate lineages (P2S30 and P4S77-SYR) observed in the UPGMA tree represented accessions of mixed ancestry (Figure 1B).

The presence of distinct clusters among accessions was further evaluated using a Bayesian approach as implemented in BAPS. The partition with the highest log marginal likelihood (−1,305.93) consisted of five clusters that were highly congruent with both the UPGMA tree clustering and the results of STRUCTURE at K = 5 (Figure 1C).

Most accessions (38/80) were grouped into cluster A, followed by cluster B (24/80), and clusters C–E (with 8, 8 and 2 accessions, respectively) (Figure 1). Accessions from the Istria region prevailed in cluster A (18 of 38), followed by those from Dalmatia (9) and continental Croatia (7). Only four foreign accessions were grouped in this cluster. Accessions from Dalmatia grouped mostly into cluster B (15/24), together with four accessions from continental Croatia and five foreign accessions. Eight accessions grouped into cluster C, of which five were from continental Croatia, one from Dalmatia and two foreign accessions. Cluster D comprised three foreign and four Istrian accessions, as well as a Dalmatian one, whereas two foreign accessions grouped into cluster E. An accession from Istria and a foreign accession (Syria) assigned to cluster D showed highly admixed ancestry and formed separate lineages on the UPGMA tree.

Accessions from Slovenia and Serbia grouped in cluster A, including one accession from the USA and one from Brazil. Cluster B comprised accessions from Balkan countries (Bulgaria, Montenegro and Serbia), Spain and the USA. Accessions from the Czech Republic and Uzbekistan grouped in cluster D, whereas two Brazilian accessions formed the cluster E.

The AMOVA was used to partition the total molecular diversity among and within groups of accessions originating from the three Croatian regions and abroad (Table 3).

AMOVA showed that approximately 97.25% of the total genetic diversity was attributable to differences among accessions within the two origin groups resulting in a non-significant φ_ST value (Table 3(A)). However, the partition of the total genetic diversity among and within three Croatian regions revealed that a significant (p < 0.0001) proportion of total variation (12.29%) was due to genetic differences of accessions among the three Croatian regions (Table 3(B)). The same was true in AMOVA that included the accessions from three Croatian regions and a group of foreign accessions (Table 3(C)) although the proportion of total variation (8.89%) attributable to the among-group component was somehow lower, but still highly significant (p < 0.0001).

The pairwise comparisons between the four groups of accessions revealed that the partition between the group of foreign accessions (P4) and those from continental Croatia (P1) was non-significant as it was between the Dalmatian group (P3) and foreign accessions (P4) (Table 4). However, all the partitions including the Istrian group of accessions (P2) were significant, and the highest φ_ST was obtained when analysing the partition between the Istrian (P2) and Dalmatian (P3) groups.

The cultivation and use of garlic are a part of the Croatian tradition, and each region or village has its local landraces preserved for generations (Dumičić et al., 2014). In the last decades, significant efforts have been made to collect and characterise national genetic resources to ensure their long-term conservation in ex situ collections, and thus, their sustainable management. The Institute of Agriculture and Tourism (IPTPO) collection was established in 2014 and is continually expanded ever since. A molecular characterisation of the collection has not been carried out so far, and the available published data are sparse (Poljuha et al., 2018).

In this study, we aimed to molecularly characterise the current Croatian garlic collection and to assess its genetic diversity and structure. We genotyped 64 accessions originating from three Croatian regions, along with 16 introduced cultivars, using 13 polymorphic SSR markers. The study revealed a total of 48 distinct genotypes in the group of Croatian accessions, while all analysed foreign accessions were distinct genotypes. The diversity of the latter set (also confirmed by a higher percentage of polymorphic loci and Shannon’s index) was expected because the foreign garlic materials analysed in the present study were intentionally selected to broaden genetic variation in the sample set. They were obtained from the collection of the USDA (Pullman, WA, USA) and chosen to represent 7 out of 10 different phylogenetic groups identified in previous AFLP-based diversity analysis of garlic (Ipek et al., 2003). In addition, within this set of samples, we analysed three accessions from the Brazilian collection at Embrapa. The diversity of Croatian garlic accessions observed through different parameters was relatively similar across all three regions (Table 1). The most diverse group of Croatian accessions originated from the Dalmatia region, with 19 distinct genotypes, out of the 26 analysed. Although the observed diversity is probably related to the sample size, Dalmatia is a well-known location of garlic cultivation (Ožanić, 1955; Dumičić et al., 2013), with a few widely distributed and highly valued ecotypes, such as Ljubitovački šarac and Brgudski ozimi. However, the highest genotypic richness and the number of private genotypes were found in the set of accessions from the Istria region, indicating the potential of this area as a valuable source of genetic material for breeding. Genotypic richness could correlate with specific and diverse microclimate conditions in a relatively small area.

In this research, we chose 13 polymorphic SSR markers, developed in previous studies (Ma et al., 2009; Cunha et al., 2012; Ipek et al., 2015) and applied in the analyses of garlic genetic resource collections and germplasms (Zhao et al., 2011; Jo et al., 2012; Cunha et al., 2014; Barboza et al., 2020). The SSR markers revealed considerable genetic diversity among 61 MLGs identified, with Dice’s distance ranging from 0.019 to 0.773, with an average value of 0.346. Usually, SSR-based methods rely on identifying a maximum of two peaks in an electropherogram, measuring the fluorescence of correctly matched tags for either homozygous or heterozygous individuals. This becomes challenging in species where genome duplication is common. In this study, we observed from one to five peaks on a single SSR marker per accession (Table S2 in Supplemenatry Materials). Moreover, the observed variations in garlic SSR loci could also occur in part as a result of replication slippage (Hosseinzadeh-Colagar et al., 2016). Since we were able to observe clear peaks and other authors also reported more than two alleles per SSR locus in diploid Allium species (Fischer and Bachmann, 2000; Ipek et al., 2008, 2015; Ma et al., 2009; Jo et al., 2012), we disregarded slippage in favour of suggested genome duplication. Even if the number of peaks could be matched, there is no objective method to determine a cut-off value that would accurately describe the actual phylogenetic relationship of the accessions, as confirmed by other molecular methods, breeding programme records or even phenotypic observation. This does not devalue the use of SSR for phylogenetic conclusions if a stringent methodology for recognising peaks is observed. A few authors reported a significant diversity within garlic collections, using both morphological and molecular data (Zhao et al., 2011; García-Lampasona et al., 2012; Chen et al., 2014; Cunha et al., 2014; Barboza et al., 2020). In this study, AMOVA analysis revealed 87.71% of the total molecular diversity within and 12.29% among Croatian groups. The same pattern was observed when evaluating the three Croatian regions and a group of foreign accessions. These results support regional structuring and the existence of significant diversity within local populations as a consequence of probable ancestral sexual reproduction followed by independent domestication history (Cunha et al., 2014). These values are in accordance with the results of Zhao et al. (2011) and Kumar et al. (2019), who observed similar percentages of variance (84.4% and 84.0%, respectively) attributable to within-population differences.

Clustering analysis revealed five main clusters among the 64 Croatian and 16 foreign accessions (clusters A–E, Figure 1). The distinction of clusters based on geographical origin was not clear, mainly due to unknown historical origin of accessions, frequent exchange of propagule material among farmers from different Croatian regions and the vegetative propagation of the species. Istrian garlic accessions were mostly grouped in cluster A and Dalmatian prevailed in cluster B, whereas the accessions from continental Croatia were evenly distributed in clusters A–C (Figure 1A). Since our primary intention was to catalogue the collection and detect duplicate or redundant genotypes, we were especially interested in the closely related samples. With 13 markers used, we detected only 10 MLG groups in 80 analysed samples. All identical and closely related accessions were observed in clusters A–C, and their grouping was based on geographical origin (e.g. P2S19, P2S20 and P2S21; P2S32, P2S33 and P2S34; P3S42, P3S45 and P3S49; P1S09, P1S10 and P1S11; P1S15 and P1S16; etc.) or traditional ways of plant material exchange in these areas (e.g. P3S40 and P4S70; P3S47, P3S48, P3S52, P3S57, P3S59 and P1S01; etc.).

We also noticed a great diversity among accessions from some narrow geographical areas, such as the Central Istrian Pazin area, where we identified diverse groupings of collected samples. Closely related accessions P2S23, P2S24 and P2S25 grouped in cluster A, together with the Istrian genotype P2S27 with a peculiar morphological and aromatic characteristics. In contrast, two accessions from the same area, P2S18 and P2S37, grouped in cluster D, together with geographically closed P2S26 accession. Such areas with pronounced diversity are a valuable source of plant material for breeding purposes. In this study, we also identified the relationship between some widely distributed and highly valued cultivars, such as ‘Ljubitovečki šarac’ (accession P3S53), which displayed admixture of two populations with membership shared between clusters A and B (0.63:0.37, respectively).

A few authors have reported a moderate level of correlation between genetic clustering based on SSR markers with geographical origin (Ma et al., 2009; Zhao et al., 2011; Jo et al., 2012; Kumar et al., 2019), whereas other studies reported none (Ipek et al., 2003, 2008; Volk et al., 2004; Paredes et al., 2008; Barboza et al., 2018). Clustering by phenotypic characteristics observed by different authors is also variable (Bradley et al., 1996; García-Lampasona et al. 2010; Ovesná et al., 2014; Ovesná et al., 2015; Barboza et al., 2020). In our study, accessions were not clustered according to their flowering behaviour (Table S1 in Supplementary Materials). The lack of correlation between morphological variation and marker pattern may originate from the non-association of the markers used for genotyping with phenotypic characteristics but also can be rooted in epigenetic variation in garlic lines, propagated in some regions for centuries. Barboza et al. (2020) observed and explained the discrepancy between the level of genetic and phenotypic variation found in the Argentine garlic collection by studies of Gimenez and García-Lampasona (2018) and Gimenez et al. (2016). They indicated frequent and dynamic epigenetic changes in non-coding and coding regions of the garlic genome, even under normal field-cultivation conditions, in comparison to less frequent genetic changes, as possible causes for the observed discrepancies between molecular and phenotypic variation.

In our analysis, foreign accessions were present in all clusters (Figure 1). Brazilian cultivars ‘Chinês Real’ (P4S65) and ‘Ugarte’ (P4S67) formed the genetically most distant cluster E. However, the third Brazilian accession P4S66 (cv. ‘Peruano’) grouped in cluster A, in a subcluster with accession P2S28 from the Istrian region. All the above-mentioned Brazilian cultivars were previously analysed by SSRs in the study of Cunha et al. (2014) using a set of 130 accessions. In their study, cultivar ‘Peruano’ formed a separate lineage as the only accession with admixed ancestry. Two Polish accessions (P4S71 and P4S72) grouped in cluster C, whereas accession P4S69 from North Moravia in Czech Republic and P4S80 originated from Uzbekistan grouped in cluster D. Accession P4S80 (UO85) was classified as A. longicuspis L. in the GRIN (USDA) database but in our analysis did not separate from other accessions of A. sativum L. This is in accordance with the study of Ipek et al. (2015) in which accessions of both species were clustered together based on their flowering ability but not on their geographic origin. Volk et al. (2004), who examined the diversity of 211 commercial and federally held accessions maintained by the US National Plant Germplasm System (NPGS), including the clone U085, also report the indistinguishability of A. longicuspis and A. sativum.

Our analysis revealed two accessions that formed separate lineage having admixture ancestry – Syrian P4S77 and P2S30 (Figure 1). Accession P4S77 is classified as cv. ‘Kisswani’ in the GRIN database (USDA) and the cultivar originates from Syria. Accession P2S30 is recorded as cv. ‘Music’ in the IPTPO collection. It belongs to A. sativum ssp. ophioscorodon, and the hardneck garlic type that produces scapes or flower stalks originates from Italy. This cultivar is characterised by high overwintering ability, increased vigour and high marketable yields (Walters, 2008).

This study confirms the usefulness of SSR genotyping for examination of diversity within garlic collections and provides important information for further management of Croatian garlic collection, as well as its relationship to known international garlic cultivars.

Garlic is a staple crop present in the Mediterranean culture for hundreds of years. As highly valued vegetative species, it was enjoyed and traded wherever people moved, and it is evident that farmers went through great efforts to select unique varieties. Assessment of genetic diversity and structure of the Croatian garlic germplasm supported regional structuring and the existence of significant diversity within local populations. Genetic differentiation of the accessions based on origin provides insights into the movement of garlic genotypes, as well as ethnological and historical data. This study is the first attempt of a molecular characterisation of collected genetic resources of garlic in Croatia with the aim of setting the course for future preservation strategies with particular emphasis on the value of diversity in the context of climate change both on macro and micro levels.