The Development of New SSR Markers and an Assay for Genotyping Sweet Cherry (Prunus avium L.) in One Reaction

TL;DR

This paper introduces new SSR markers and a one-step genotyping method for sweet cherry that improves accuracy and efficiency.

Contribution

A new set of SSR markers and a one-tube genotyping assay with significantly higher resolution than existing methods.

Findings

Thirty new SSR markers were developed from 299 sweet cherry genotypes.

A one-PCR assay using 16 markers achieved high resolution for genotyping.

The new system reduces the chance of random genetic identity by seven orders of magnitude compared to existing markers.

Abstract

Sweet cherry (Prunus avium L.) exhibits relatively low genetic diversity because of the self-compatibility of some varieties and repeated crossings of the same genotypes. High-quality markers are therefore needed for their reliable discrimination. However, the most currently used simple sequence repeat (SSR) markers offer only limited resolution for genotyping purposes. Here, thirty new highly polymorphic SSR markers were extracted from whole-genome sequences of 299 sweet cherry genotypes. Then, 16 highly polymorphic SSR markers were selected, multiplexed into one PCR, and successfully verified on a collection containing 294 unique genotypes. Compared with the set of 16 SSR markers recommended by the European Cooperative Programme for Plant Genetic Resources (ECPGR) for sweet cherry genotyping, our newly developed system has a seven orders of magnitude lower probability of the random identity of two genetically distinct samples than the ECPGR set (10−19 vs. 10−12). This higher resolution not only enables more precise genotyping but can also be successfully used for parentage or population analyses. This new and unique one-tube approach for sweet cherry genotyping will substantially simplify genotyping workflows, minimize errors, and save labor, time, and cost.

Figures

Full text

1. Introduction

Sweet cherry (Prunus avium L.) produces one of the earliest ripening and most tasty fruits in temperate climates, making it a popular choice among consumers. It is therefore not surprising that worldwide production is constantly increasing, reaching nearly 3 million tons in 2023 (www.fao.org (accessed on 12 September 2025)). The three main producers are Turkey, Chile, and the USA, followed by Uzbekistan, Iran, Greece, and Spain. It is estimated that there are thousands of sweet cherry varieties in the world, and breeding programs are ongoing to develop new and better varieties with extended ripening periods. For instance, the U.S. Department of Agriculture alone has more than 800 sweet cherry entries on file (www.usda.gov (accessed on 12 September 2025)), and Schuster et al. [1] described 1700 different varieties.

For cherries in particular, it is extremely important to distinguish varieties with certainty, as they differ in critical parameters for successful cultivation. Some are self-compatible, enabling them to grow in monocultures, but most varieties are self-incompatible [1] and need a pollination partner. To cross-pollinate each other, they have to possess at least one different so-called S-allele, a variant of the locus responsible for controlling pollination success [1,2]. Furthermore, for successful harvesting, it is necessary to grow a variety in a given area that is capable of meeting its chilling hour requirements to end dormancy, ensuring ample blooming [2,3]. When dormancy is not released naturally, dormancy breakers such as Dormex^®^ have to be used, increasing the cost of production [4]. The time of flowering is another critical trait, as growing early-flowering varieties increases the risk of crop loss in areas with frequent spring frosts, and expensive anti-freezing measures must be used when temperatures drop below zero [5]. Finally, varieties differ in such parameters as fruit size, firmness, and color; trees may be more or less susceptible to biotic and abiotic stresses, or ripen in different periods. Selection of the right sweet cherry variety is, therefore, a critical step for a farmer’s success. Hence, it is not surprising that growers want to be sure about the variety they plant in newly established orchards that will be managed for many years.

The identification/verification of sweet cherry varieties is usually done by phenological and phenotypical evaluation. However, assessing fruit-bearing trees requires at least two years to minimize inter-seasonal variation, and the evaluation must be performed by highly experienced specialists. On the other hand, molecular genetic approaches offer much quicker and more reliable techniques for identifying varieties, use a minimal amount of biological material, and can be performed throughout the year, even in young seedlings or trees without fruits. This approach is called genotyping or fingerprinting, and uses molecular markers with sufficient discrimination power to distinguish individual varieties.

Nowadays, genotyping is commonly performed with SSR (Simple Sequence Repeat) or SNP (Single Nucleotide Polymorphism) markers. Both types of markers have their advantages and disadvantages. SSR markers, suitable for genotyping, should be highly polymorphic, independently inherited (no genetic linkage between them), and should exhibit low stuttering and consistent reproducibility [6]. Twelve or fewer SSR markers are usually sufficient to distinguish fruit crop varieties, including those capable of self-pollination. The number of SSR markers recommended for reliable variety differentiation depends on the marker informativeness, the genetic diversity of the fruit species [6], and ploidy (higher species ploidy may lower the number of SSR markers needed for genotyping [7,8,9]). On the other hand, SNPs are suitable for extensive genotyping projects and for the construction of genetic maps, lowering the cost per analysis. As SNPs are usually bi-allelic, it is necessary to analyze substantially more SNPs than SSR markers, which requires advanced biostatistic methods to interpret results [10]. Despite advances in whole-genome sequencing techniques, SSR markers remain the gold standard for genotyping not only in fruit crops but also in humans, as SNPs are less suitable and too expensive for routine genotyping of low numbers of samples [6].

Many SSR markers have been developed in sweet cherries [11,12,13], and SSRs from other closely related diploid Prunus species—such as peach [14,15,16], Japanese plum [17], and almond [18]—have also been used for their genotyping. In 2009, the European Cooperative Programme for Plant Genetic Resources (ECPGR) recommended 16 SSR markers for genotyping sweet cherries by fragment analysis [19]. These markers have been used for genotyping several collections, e.g., 383 unique genotypes from the German Fruit Genebank [20]; however, some authors have revised the ECPGR SSR marker set and added new markers as some varieties could not be distinguished [21,22].

For SSR marker analysis, a capillary genetic analyzer is usually used to separate PCR fragments with one-nucleotide accuracy. This method uses several different fluorescent dyes in parallel with amplicons ranging from 80–600 bp in length, allowing the analysis of several SSR markers in one reaction. This approach, familiar from human genetics, has also been adapted for fruit crop genotyping. For example, Chambers et al. [7] amplified and analyzed eight SSRs in a single reaction to genotype up to octaploid strawberries; Akin et al. [23] multiplexed 14 SSRs in one reaction to genotype diploid hazelnuts; Zurn et al. [24] designed primers to co-amplify 10 SSRs in one reaction to distinguish accessions in the United States National Pyrus collection (diploids and triploids); Bassil et al. [8] published a 10-SSR multiplex for fingerprinting up to octaploid blueberries; Čmejlová et al. [25] optimized a single reaction amplifying 17 SSRs for genotyping diploid and triploid apples; and Čmejlová et al. [9] developed a set of eight SSR markers for the identification of hexaploid plum varieties in one reaction.

Genotyping by SSR markers may not only be useful for confirming the identity of varieties, but can also be used for further analyses. If a sufficient number of SSR markers is analyzed in a sufficient number of varieties (the number of used varieties depends on SSR/population heterogeneity), it is possible to use SSRs for parentage, pedigree, and/or population analyses. For these purposes, highly polymorphic markers distributed across the genome are needed. Information about parents is highly valued by breeders, while pedigree and population analyses are important mainly for scientists exploring relationships and diversity, especially in gene pools. The usefulness of SSR markers for these analyses has been confirmed, for example, by Barreneche et al. [21], Marchese et al. [26], or Liu et al. [27] in population and pedigree analyses of sweet cherries; Liu et al. [28] developed SSR markers for paternal identification in Camellia sinensis; Zhang et al. [29] identified parents in walnut F1 progeny by SSRs; and Vuletin Selak et al. [30] used SSR markers to study paternity in seeds from mixed olive orchards.

The objectives of this work were (i) to identify new high-quality SSR markers using whole-genome sequences of 299 sweet cherry varieties or hybrid material; (ii) to develop an SSR marker set for genotyping sweet cherries with high discrimination power; (iii) to provide an easy and quick genotyping assay in one reaction; and (iv) to use the results obtained by the assay for advanced analyses (pedigree, population, and parentage analyses). This one-reaction approach should make the determination of sweet cherry varieties accessible not only to the scientific community but also to national regulatory authorities and the general public for routine analyses at relatively low cost. The final aim of this study was to compare the newly developed identification system with the previously used genotyping system recommended by the ECPGR.

2. Results

This study was conducted according to the workflow depicted in Figure 1 and is described in detail in the following sections. The list of varieties/genotypes used in this study is presented in Supplementary Table S1.

2.1. Newly Identified SSR Markers

Based on 299 whole-genome sequences, 30 new SSR markers were extracted from the sweet cherry genome of the ‘Tieton’ variety (Table 1) that passed the following criteria: minimally, three alleles with a maximum allele frequency of 0.5; relatively short repetition to avoid high stuttering, and polymorphic information contents based on allele sequences (PIC_HipSTR_) greater than 0.8 for SSR markers with dinucleotide repetition, 0.7 for markers with trinucleotide motif, and 0.6 for SSRs with longer repetition unit. Most of the identified SSRs possessed a dinucleotide repetition (19 SSRs), while five SSRs were characterized by a repetition composed of three or four nucleotides, and one SSR marker contained a hexanucleotide repetition. Some of these SSR markers were located in a gene, but in non-coding regions (introns, untranslated regions); most were situated in intergenic regions (Table 1).

All 30 SSRs were successfully verified by fragment analysis on 49 samples to allow statistical evaluation (Table 2; primers used for this analysis are shown in Supplementary Table S2). The number of alleles in individual markers ranged from 3 to 11, with the number of effective alleles (Ae) from 2.553 (Pav_chr7_863) to 6.651 (Pav_chr3_002). Although PIC_HipSTR_ was one of the criteria for selection of high-quality SSR markers, the observed PIC values based on fragment analysis in the chosen subpopulation of 49 varieties were lower for some SSR markers, ranging from 0.540 (Pav_chr7_863) to 0.833 (Pav_chr3_002). Lower PIC values from fragment analysis compared to PIC_HipSTR_ values were expected since the HipSTR algorithm works with the sequence of alleles, while fragment analysis discriminates only the length of alleles, irrespective of their sequence. The existence of a one-length allele with two different sequences was confirmed by the presence of SNPs in samples homozygous for the particular SSR marker in the fragment analysis (for example, Pav_chr2_995, Pav_chr6_387, Pav_6_989, or Pav_chr7_863). Expected heterozygosity varied between 0.608 (Pav_chr7_863) and 0.850 (Pav_chr3_002); the minimal observed heterozygosity of 0.571 was found in Pav_chr7_863, and the maximal observed heterozygosity of 0.918 was detected in Pav_chr3_706.

2.2. Multiplexing of 16 SSR Markers into One Reaction

Two SSR markers per chromosome with the highest Ho, He, PIC, and Ae values were chosen for multiplexing into one PCR (Table 2). Positions of the SSR markers were also taken into account since two SSR markers identified on chromosome 2 were separated by only a 200 kbp region, and on chromosome 5, Pav_chr5_144 and Pav_chr5_903 lay only 70 kbp from each other. However, SSR markers Pav_chr1_028, Pav_chr5_059, Pav_chr6_641, Pav_chr7_286, and Pav_chr8_224 failed to amplify in a multiplex setup, even if very high concentrations of primers were used. Therefore, alternative SSR markers extracted from previously known genetic maps were also tested, and the best ones from each chromosome were chosen: EMPA005 (chr1), UCD-CH11 (chr2), BPPCT037 (chr5; ECPGR recommended), and CPPCT006 (chr8; ECPGR recommended). The main criteria for their selection were high polymorphism, specificity under the used PCR conditions, and low stuttering. Thirty-two primers for PCR amplification in one reaction were successfully multiplexed, and concentrations of primers were fine-tuned to obtain comparable signal heights for all markers. Alleles with a very low signal were observed for UCD-CH11. Whole-genome sequences in the respective samples were checked, and an A/C polymorphism in the forward primer-annealing site was found. To improve annealing, the primer was degenerated in this position, and while the signal of allele 131 remained somewhat lower, it was sufficient for reliable analysis. Primers for the amplification of the final set of 16 selected SSR markers, hereafter termed the 16in1 system, are summarized in Table 3. Figure 2 demonstrates a typical output from fragment analysis on a capillary genetic analyzer.

2.3. Statistical Evaluation

The 16in1 genotyping system was verified on 315 samples (most of them analyzed in multiple biological replicates); allele combinations for each marker in each sample are summarized in Supplementary Table S3. During testing, clones of one variety (Supplementary Table S3, differing alleles are background colored), and several misidentifications and/or synonyms that showed the same genotypes but different names were identified (Supplementary Table S4). These discrepancies were also confirmed using the ECPGR marker set (Supplementary Table S5). On the other hand, samples identified as identical by the ECPGR system were separated by the 16in1 kit into two or even three varieties, differing by only one or two alleles (Supplementary Table S6). This means that these samples had identical genotypes or differed by only one or two alleles in 30 analyzed SSR markers (two SSR markers are common for both 16 SSR marker sets); therefore, they were labeled herein as varieties/clones. To statistically evaluate the quality of individual SSR markers and the whole set, only unique 16in1 genotypes/varieties were used (Supplementary Table S3).

The results of 294 unique genotypes using the 16in1 genotyping system were statistically evaluated in POLYGENE software (version 1.7) (Table 4). The highest polymorphism was observed for SSR marker Pav_chr3_002 (PIC 0.846; Ae 7.223), and the lowest for EMPA005 (PIC 0.606; Ae 2.832). The whole 16in1 genotyping system has an average PIC of 0.767 with Ae of 5.227. Additional statistical parameters, including the allelic richness, number of private alleles, and Shannon’s information index, can be found in Supplementary Table S7. The probability of the random identity of two genetically different samples for the whole set of SSR markers (PI_t_) was calculated to be 1.27 × 10^−19^. All markers were amplified in all samples tested; no null alleles were recorded (the probability of null alleles ranged from 0.011 to -0.035). Statistical parameters of the collection of 49 cultivars were compared with those of the full set of 294 cultivars, and the PIC values for individual markers differed by no more than 5%. The cultivar collection used for the primary analyses therefore provided a reasonably good representation of all analyzed cultivars). An overview of all alleles with their frequencies is presented in Supplementary Table S8.

For comparison, the same set of unique genotypes was also analyzed with the ECPGR SSR marker set. The statistical evaluation of results obtained by ECPGR markers (Table 5) showed that EMPAS06 was the best SSR marker with PIC 0.806 and Ae 5.813; however, EMPA017 had the lowest PIC value of 0.235 and Ae of only 1.329. Other statistical characteristics are summarized in Supplementary Table S9. Varieties with homozygous null alleles (i.e., no amplification of this marker in some samples) were observed in SSR markers EMPaS10, EMPaS06, and EMPA026 (Supplementary Table S10; also contains a list of all identified alleles and their frequencies). The whole ECPGR genotyping system showed an average PIC of only 0.564 and Ae of 2.994, and the PI_t_ value for the analyzed samples was calculated to be 1.00 × 10^−12^.

2.4. Population Analysis

The results obtained with both SSR marker sets were further used for population analysis in the STRUCTURE software (version 2.2.4). Dendrograms for both 16in1 and ECPGR marker sets were created in POLYGENE to establish the rank of samples for their input to STRUCTURE (Supplementary Table S11). The analysis divided the 16in1 analyzed samples into 3 and 9 groups with the highest delta Ks of 189.02 and 37.11, respectively, while the ECPGR analyzed samples were separated into 2 and 10 groups with the highest delta Ks of 58.82 and 45.39, respectively (Supplementary Table S12). The composition of individual groups was evaluated from a pedigree perspective and considering geographical origin (K3 and K2 groups). ‘Bing’, ‘Germersdorfer’, ‘Kordia’, ‘Lambert’, ‘Rube’, ‘Schmidt’, ‘Stella’, and ‘Van’ were repeatedly observed in pedigrees of varieties with high Score–Inferred clusters (S-Ic), and the frequency of their incidence in groups was therefore evaluated. For the finer subdivision into groups (K9 and K10), only those representatives with an S-Ic > 0.9 were highlighted, as some groups contained only a small number of varieties. Figure 3 depicts the results of STRUCTURE analysis for both 16in1 and ECPGR marker sets sorted by Q.

For further analyses, the STRUCTURE groups resulting from the 16in1 system were labeled K3.1 to K3.3, and K9.1 to K9.9 (Supplementary Tables S13 and S14). The K3.1 group contained 58 samples in total; 29 varieties had an S-Ic higher than or equal to 0.8. In this group, modern cultivars prevailed—25 out of 29 varieties with the highest S-Ic had ‘Van’ (18/29), ‘Stella’ (10/29), ‘Kordia’ (7/29), or ‘Bing’ (5/29) in its pedigree, while 4/29 had an unknown origin. From a geographical point of view, 13 out of 29 varieties originated in Canada, and 9/29 were bred in the Czech Republic. The K3.2 group was composed of 117 highly heterogeneous samples, many of which have an unknown pedigree and may be considered landraces. Sixty-five of them had an S-Ic higher than 0.8. Six varieties out of twenty with known parents, and an S-Ic ≥ 0.8 had ‘Germersdorfer’ in their pedigree, while the ancestor of four varieties was ‘Schmidt’. Geographically, 17/65 varieties originated from Germany, and 10/65 from the USA. The remaining 119 varieties belonged to the K3.3 group. This group was also very heterogeneous, and contained many varieties with unknown pedigree and landraces. Sixty varieties had an S-Ic higher than 0.8; ‘Rube’ was one of the parents in 7/23 varieties with a known pedigree in this subgroup. Varieties coming from Europe clearly prevailed—only 4 varieties out of 60 with an S-Ic ≥ 0.8 had their origin outside of Europe (USA); the highest number of varieties from one country had their origin in Germany (17/60 varieties).

At a finer K9 distribution, the highest S-Ics (>0.9) in individual groups were found for ‘Droganova’, ‘Germersdorfer’, ‘Boambe de Cotnari’ (a clone of ‘Napoleonova’), ‘Hedelfingenska’, ‘Rebekka’, ‘Lyonska rana’, ‘Van’, and ‘Nero II c1’. ‘Reidler’ had the best score in group K9.8, but only 0.848. In the 16in1 genotyping system, K9.4 (‘Hedelfingenska’ group) and K9.7 (‘Van’ group) prevailed, with 54 and 51 members, respectively. Many identified representatives of individual groups corresponded to varieties frequently used in crossings, e.g., ‘Germersdorfer’, ‘Van’, or ‘Nero II c1’ (Supplementary Table S2).

Similarly, groups analyzed by the ECPGR system were labeled K2.1 to K2.2, and K10.1 to K10.10 (Supplementary Tables S15 and S16). K2.1 contained 152 varieties; 110 of them showed an S-Ic higher than or equal to 0.8; modern cultivars with a known pedigree predominated. As for geographical origin, 27 of them came from Canada, 21 originated from Germany, 19 were registered in the Czech Republic, and 13 of them were found/bred in the USA. K2.2 was represented by 142 varieties; 93 of them had an S-Ic ≥ 0.8. European landraces represented the main part of this group. Thirty of them came from Germany, and twelve from the Czech Republic. The representatives of individual groups for K10 with S-Ic scores >0.9 were ‘Van’, ‘Lyonska rana’, ‘Nero II c1’, ‘Alfa’, ‘Cristimar’, ‘Rebekka’, ‘Napoleonova’, ‘Zukunft’, and ‘Mednanska’. In group K10.4, two varieties (‘Germersdorfer’ and ‘Schneiders Spaete Knorpelkirsche’) sharing the same ECPGR genotype were identified (score 0.895). The largest ECPGR groups were K10.1, comprising 46 varieties (the ‘Van’ group), and K10.10 (‘Mednanska’) with 40 genotypes. Identified representatives of individual groups corresponded to varieties frequently used in crossings (Supplementary Table S2).

To summarize the results of the K3 (16in1) and K2 (ECPGR) group analyses, the 16in1 genotyping system resulted in a clearer distribution of genotypes. Most of the modern cultivars were grouped into K3.1; K3.2 and K3.3 were heterogeneous with a high content of landraces. Non-Europe genotypes concentrated in K3.2, whereas K3.3 contained only a small number of them. Varieties with ‘Germersdorfer’ and ‘Schmidt’ in their pedigree were assigned to K3.2; known relatives of ‘Rube’ were only found in group K3.3. In contrast, the ECPGR groups were less clearly defined from both pedigree and geographical origin points of view; however, modern cultivars were preferentially found in the K2.1 group.

Comparing representatives of individual groups K9 and K10 obtained by the 16in1 and ECPGR systems, respectively, six groups exhibited the same best representative (‘Germersdorfer’, ‘Boambe de Cotnari’ (clone of ‘Napoleonova’), ‘Rebekka’, ‘Lyonska rana’, ‘Van’, and ‘Nero II c1’). Moreover, according to the STRUCTURE results and dendrograms, ‘Hedelfingenska’ as a representative of K9.4 is closely related to ‘Mednanska’ (K10.10). The size of related groups was similar, e.g., K9.3 and K10.8, both represented by ‘Napoleonova’, contained 23 and 24 varieties, respectively; but some groups differed substantially, for example, group K9.5 of ‘Rebekka’ consisted of 41 varieties, while corresponding K10.7 comprised 20 varieties.

Next, a tanglegram comparing dendrograms generated from the results of both the 16in1 and ECPGR systems was constructed to further evaluate similarities between the genotyping systems (Supplementary Table S17). Dendrograms were rearranged with respect to each other, and although some similarity between sample groups may be observed, only several varieties were distributed in exactly the same subtrees (branches) in both dendrograms (highlighted by colored lines between individuals in the two dendrograms). The tanglegram also shows that several varieties genotyped by ECPGR-recommended SSR markers having the same genotype were distinguished by 16in1 SSR markers, but by only one or two differing alleles. The tanglegram Mantel’s statistic value of r = 0.5552 indicated a moderately strong positive correlation between both dendrograms, and the genetic distances between individuals were similar in both sets of markers (p-value 0.001).

2.5. 16in1 Kit Validation: Parentage Analysis

Parentage analysis in POLYGENE software was used to look for parent pairs for all samples, using a very strict criterion of null Trio loci mismatching (Supplementary Table S18). Some suggested parent combinations were pairs, including genotypes differing by only one allele (for example, ‘Napoleonova’, ‘Boambe de Cotnari’, and ‘Esperens Knorpelkirsche’); these were counted only once for statistical evaluation (Supplementary Table S19). Suggested parents were compared with the parent combination for the particular variety, if known.

For the 16in1 system, 244 parent pairs were discerned for 294 analyzed samples. After removing parent pairs containing clonal varieties, 189 parent pairs remained. At least one parent combination without mismatching alleles was found for 123 out of 294 samples, i.e., 1.54 parent pairs were identified per variety on average. The whole collection contained 34 varieties with both parents known and analyzed, and 28 (82.35%) were successfully identified as the only one suggested parent combination.

In comparison, one to five parent combinations were suggested in 190 varieties using the ECPGR set of SSR markers. In total, 654 parent pairs were chosen by POLYGENE; after excluding parent pairs consisting of clones, 521 parent combinations remained, i.e., 2.74 parent pairs without mismatching alleles were found per variety on average (however, only the top five LOD score pairs with a null score of Trio loci mismatching were extracted). The correct parent combination, as the unique one, was identified in only 5 out of 34 varieties (14.71%).

Using the 16in1 SSR marker set, a single parental combination was proposed for 86/123 (69.9%) varieties for which a parental combination was suggested, and no variety had five parent pairs identified. By contrast, using the ECPGR markers, a unique parental combination was proposed for only 54/190 (28.4%) varieties, and 38/190 (20%) of varieties had five suggested parent pairs. The number of identified parent pairs for all analyzed samples is summarized in Figure 4.

Combining the results of 16in1 and ECPGR analyses, the same parent pairs were suggested for 54 varieties with unknown parents or open-pollinated mothers (Supplementary Table S19). For example, ‘Grace Star’ came from a crossing of ‘Lapins’ and ‘Burlat’; for ‘Horka’, ‘Van’ was identified as a pollinator of ‘Kordia’; parents of ‘Starking Hardy Giant’ were identified as ‘Windsor’ and ‘Schmidt’; and ‘Star’ was an offspring of ‘Deacon’ and ‘Lambert’.

2.6. Analysis of ‘Van’ Relatives: An Example of a Highly Inbred Group

The group of ‘Van’ relatives presents the most challenging task for cherry genotyping and parentage analysis using SSR markers because of a very high degree of inbreeding. Using the 16in1 system, the ‘Van’ family was defined as those from the K9.7 group based on STRUCTURE analysis that contained more than 50% of the K9.7 genome (Supplementary Table S20). This group, named K9.7-50%, contained 35 varieties; most of the members were known offspring of ‘Van’; however, some accessions had no known parentage. In the next step, the group was narrowed to those with 75% of the K9.7 genome; this K9.7-75% group comprised 14 varieties. A strong drop in heterozygosity was observed, and the PIC value decreased from 0.767 (whole collection) to 0.597 (K9.7-50%) and 0.506 (K9.7-75%), respectively. Concomitantly, Ae declined from 5.227 to 3.136 and 2.557, respectively (Table 6). Loss of alleles in all individual SSR markers was detected and may be found in Supplementary Table S20. Hence, the original PI_t_ decreased from 1.27 × 10^−19^ to 2.37 × 10^−13^ (K9.7-50%) and 1.04 × 10^−10^ (K9.7-75%), respectively.

The ECPGR system also showed a high decline in these parameters (Table 7, more details in Supplementary Table S21)—the PIC value calculated for all samples (0.564) dropped to 0.402 for K10.1-50% (the group of ‘Van’ relatives in the ECPGR system; 36 varieties) and 0.356 in the case of K10.1-75% (18 genotypes), and the original PI_t_ value decreased from 1.00 × 10^−12^ to 1.86 × 10^−8^ (K10.1-50%) and 1.27 × 10^−7^ (K10.1-75%). Loss of alleles in individual SSR markers may be found in Supplementary Table S21.

3. Discussion

Genotyping of sweet cherries is a challenging task due to low heterogeneity, especially in modern varieties, which stems from the tendency to breed self-compatible varieties that are usually derived from ‘Stella’ and its descendants [1], and/or the repeated usage of a limited number of self-incompatible varieties for breeding. For example, Choi and Kappel [31] analyzed inbreeding and founding clones for sweet cherries coming from North America, and 64 out of 66 analyzed modern varieties were found to be descendants of only five founding clones (‘Black Heart’, ‘Emperor Francis’, ‘Empress Eugenie’, ‘Napoleon’, and ‘Windsor’). A loss of heterozygosity in modern varieties was also confirmed by Campoy et al. [32], and a narrow genetic bottleneck in modern cherry varieties was already reported in 2010 [33]. The problematic differentiation of closely related varieties was also shown by the work of Wiersma et al. [34].

Both SNP and SSR markers are commonly utilized for fingerprinting plant varieties. Modern techniques frequently focus on SNP markers, and methods such as genotyping-by-sequencing and SNP chips are commonly used for their detection. Forty SNP markers screened by high-resolution melting analysis were successfully used for genotyping 143 varieties of sweet cherries by Fernandez i Marti et al. [35]. Campoy et al. [32] performed the first population genetic analysis in cultivated sweet cherries using a medium-density SNP marker array, RosBREED cherry 6 K SNP array v1, on 210 genotypes. Recently, Palasciano et al. [36] obtained 13,117 SNPs through the genotyping-by-sequencing of 143 cherry varieties, and proved that SNPs are an excellent tool for assessing genetic diversity and relationships. However, as demonstrated by Testolin et al. [6] and others, while SNPs are the markers of choice for en masse genotyping and in-depth analysis of genomes, close-kin differentiation or for constructing genetic maps, especially in large collections, SSR markers are still prioritized for routine genotyping. SSR marker analysis is relatively cheap, rapid, reliable, and simple to interpret, and using appropriate controls, the method can easily be implemented in other laboratories. A capillary genetic analyzer needed to perform fragment analysis is also available in many labs. Moreover, analysis of tens of SSR markers usually provides sufficient information for identifying varieties, i.e., only “small data” are produced that can be easily handled and shared. On the other hand, SNP genotyping is a modern trend in genetics that generates “big data” requiring a skilled bioinformatician for interpretation. In addition, expensive laboratory equipment is needed, so this type of analysis is usually outsourced, and laboratories thus lose full control over the analytical process. As there are many protocols, workflows, and pipelines for genotyping-by-sequencing, the reproducibility of results between labs may be of concern. The RosBREED cherry 6 K SNP array v1 from Thermo Fisher Scientific may represent a standardized genotyping platform, but it is no longer available. Moreover, Wiersma et al. [34] confirmed the difficulty of distinguishing closely related cultivars by both SSR and SNP markers. Unsurprisingly, informative SNP markers were located in the same regions as discriminating SSR markers. In light of these findings, the development of a cheap, fast, and easy-to-use genotyping tool based on SSRs seems ideal for routine sweet cherry genotyping.

The ECPGR recommended set of 16 SSRs [19] has long been used for cherry genotyping. However, in our study, some of these markers had very low polymorphism: 10 out of 16 SSR markers showed PIC values below 0.6, with the number of effective alleles <3. The average values for the whole ECPGR system were also rather low, with PIC of 0.564 and Ae of 3.0. It is not, therefore, surprising that some authors have tried to improve the ECPGR set [21,22]. In probably the largest genetic diversity and structure analysis of a combined sweet cherry germplasm from across Europe, Barreneche et al. [21] analyzed 314 cherry varieties obtained from 19 countries. Most of the collection represented original landraces from the country of provenance; the rest were popular old and modern varieties. Their final set of 14 SSR markers used for genotyping included nine loci from the ECPGR recommended genotyping set, two newly selected SSR markers previously used for cherry genotyping by different scientific groups, two loci associated with fruit weight/size and broadly used for molecular marker-assisted selection in breeding (BPPCT034 and CPSCT038), and one SSR marker genetically linked to fruit color (Pav-Rf-SSR) (reviewed in [37]). The PIC values of this system were 0.661 for 111 landraces and 0.639 for 104 early selections together with modern breeding varieties. It is noteworthy that among 314 varieties, 40 groups with shared genotypes were found. Some new supposed redundancies were identified, and others were confirmed; however, mislabeling and synonyms were not sufficient to explain all redundancies [21]. Furthermore, SSR markers that are widely utilized for the molecular-marker-assisted selection of seedlings with the perspective of bearing large fruit or having certain fruit color are not suitable for genotyping of modern varieties because they tend to be homozygous.

As only a limited number of highly polymorphic SSR markers suitable for the genotyping of sweet cherries have been described, we decided to identify new markers using the whole-genome sequences of 299 unique genotypes [38]. Using strict criteria, 30 SSR markers from all chromosomes were identified; most of them (20/30) were located in an intergenic region, the rest (10/30) were in non-coding sequences of predicted proteins (introns, 3’ untranslated regions). Some of the developed SSRs failed to be amplified in a multiplex PCR despite working in a simplex PCR setup. All of these contained AT-rich repetitions, underlying a poor rate of successful PCR amplification, as has also been shown for other SSRs, for example, in rice [39]. Other SSR markers were excluded due to their lower heterogeneity or close proximity to a more polymorphic SSR marker. A distance of at least 1 Mbp between SSR markers was chosen as the threshold to minimize the probability that markers showed genetic linkage. The selected distance was tenfold higher than the average distance at which intra chromosomal linkage disequilibrium (LD) declined below r^2^ = 0.2 [32], i.e. the threshold below which LD is considered low or negligible. Finally, twelve newly identified SSR markers were multiplexed with four previously known SSRs showing high polymorphism, in such a way that two SSR markers represented each of the eight Prunus avium chromosomes.

The final 16in1 genotyping set was optimized to be used in one PCR and fragment analysis reaction. It has an average PIC value of 0.767, which is higher than the PIC values of both the ECPGR system (tested on the same collection) and the set developed by Barreneche et al. [21]. The total probability of the random identity of two genetically different samples, PI_t_, was 1.27 × 10^−19^, which is seven orders of magnitude lower than the PI_t_ for the ECPGR set (1.00 × 10^−12^). This value reflects the lower heterogeneity of sweet cherries in comparison with apples and pears, where much lower values were obtained with 17 SSR markers, though based on a higher number of genotypes: PI_t_ ~ 10^−22^ for apples (880 genotypes; [25]), and PI_t_ ~ 10^−25^ for pears (322 diploid genotypes; manuscript in preparation). However, it is necessary to keep in mind that the PI_t_ value was obtained using a collection of European and North American varieties; results may differ for other collections or geographical regions, e.g., in Asia or in wild cherries.

Our newly developed 16in1 system was compared with the currently used ECPGR system, which contains the same number of SSR markers. Though both SSR sets worked similarly in genotyping tasks, their performance differed in population analyses and, in particular, in parentage assessments. The 16in1 set showed a more precise distribution into STRUCTURE groups, and critically led to an improvement in the identification of parent pairs. Whereas the ECPGR set of SSR markers provided many possible parent pairs for individual varieties, the 16in1 set usually suggested only one parent combination per variety and correctly identified the right parents in more than 82% of varieties in which both parents are known.

One of the most challenging tasks for SSR genotyping is the discrimination of genotypes within highly inbred groups with the possible existence of breeding bottlenecks and founder effects [33]. To better test our 16in1 system, cherries from the ‘Van’ family were chosen, since these varieties are repeatedly crossed with each other, and the ‘Van’ groups detected by STRUCTURE in both genotyping systems contained a sufficient number of varieties for further analysis (16in1: K9.7 consisted of 51 genotypes; ECPGR: K10.1 contained 46 varieties). As expected, losses of alleles, along with reductions in the number of effective alleles and in average PIC values, were observed, suggesting the presence of a breeding bottleneck and a founder effect in this group. However, the main goal of the ‘Van’ group analysis was to test the reduction rate of discrimination power in such a highly inbred group. When analyzed by the new 16in1 set of markers, PIC dropped from 0.767 for the whole collection to 0.597 for varieties with half of their genome corresponding to the 16in1 K9.7 group, and 0.506 for varieties whose genome harbors 75% of the K9.7 group. However, these values remain similar to the PIC values achieved by the ECPGR system across the entire collection of varieties, highlighting the improved discriminatory power of the 16in1 system. For comparison, the ECPGR system showed a decline in PIC to 0.402 for K10.1-50% and to only 0.356 for K10.1-75%.

In addition to the higher resolution capability of the 16in1 system, its formulation in a one-reaction assay makes it ideal for routine cherry genotyping since it reduces labor, time, cost per analysis, and mistakes. One-reaction cherry genotyping may now be readily available and affordable not only for scientific purposes, but also for professional cherry growers, nursery managers, cherry breeders, germplasm curators, or hobby gardeners.

4. Materials and Methods

4.1. Plant Material

A total of 235 sweet cherry accessions from the Collection of Genetic Resources of the Research and Breeding Institute of Pomology Holovousy Ltd. (RBIPH) (Holovousy, Czech Republic) and 64 high-quality breeding plant materials were used for the identification of SSR markers. These genotypes and population structure are described in Holušová et al. [38]. Principal component analysis and kinship analysis revealed no subpopulations or strong sample clustering.

For SSR marker assessments, 315 genotypes were obtained from RBIPH, the Central Institute for Supervising and Testing in Agriculture, and the Fytos company (Plzeň, Czech Republic) (Supplementary Table S1). Accessions included commercial and older historic varieties, important breeding accessions, and traditional Czech landraces. If possible, accessions were analyzed at least in biological duplicates.

4.2. DNA Extraction

Genomic DNA for all analyses was extracted from 100 mg of shoot phloem using an Exgene Plant SV kit (GeneAll Biotechnology, Seoul, Republic of Korea) following the manufacturer’s protocol. To obtain a fine phloem powder, bark was removed from 1- to 2-year-old shoots with a scalpel, phloem was scraped off, and the tissue was ground to a powder in liquid nitrogen with a mortar and pestle. DNA quality was evaluated by a NanoDrop Lite Spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA), and DNA was diluted to 10 ng/µL for further analysis.

4.3. Simple Sequence Repeat Identification

Whole-genome sequences were retrieved from the Sequence Read Archive study number PRJNA813711 [38]. Genomic variants of microsatellites identified with MISA (MIcroSAtellite Identification Tool, version 1.0) [40] from the ‘Tieton’ Genome v2.0 reference sequence [41] were detected using HipSTR (Haplotype inference and phasing for Short Tandem Repeats) [42], which determined sequence variants for individual cherry varieties (299 samples). The criteria for selecting high-quality SSR markers were as follows:

1. The sequence failed to map, at most, 2 out of 299 whole-genome sequenced samples;
2. At least three alleles exist;
3. A maximum repeat length of 60 nucleotides to limit stutter;
4. A maximum allele frequency of 0.5 for any allele;
5. A polymorphic information content (PIC_HipSTR_) of 0.8 or higher for SSR markers with a dinucleotide repeat unit, 0.7 for SSRs with a trinucleotide repetition, and 0.6 for microsatellites with a longer repeat unit.

4.4. Fragment Analysis

Primers were designed using Primer3 software (version 2.6.1) [43]. Based on the sequence identity with the ‘Tieton’ v2.0, ‘Satonishiki’ v1.0.a1, and ‘Regina’ v1.0 reference genomes (https://www.rosaceae.org/organism/24337, accessed on 2 May 2024), PCR product specificity was determined using the in_silico_PCR.pl script (https://github.com/egonozer/in_silico_pcr, accessed on 2 May 2024). Primers producing nonspecific or duplicate amplicons were discarded. Reads from all varieties were mapped to the ‘Tieton’ v2.0 reference, and primers leading to amplicons that contained large insertions or deletions substantially affecting a PCR product length were discarded.

Fragment analysis for all 30 identified SSR markers was performed in 49 varieties (Supplementary Table S1). Cultivars were chosen on both the parental combinations (aiming for the widest possible range of parents) and the inclusion of landraces and cultivars of unknown origin from various parts of the world (Belgium, Bulgaria, the Czech Republic, Italy, Moldova, Germany, Russia, Slovakia, and the USA). Statistical parameters of the collection of 49 cultivars were compared with those of the full set of 294 cultivars, and the PIC values for individual markers differed by no more than 5%. The cultivar collection used for the primary analyses therefore provided a reasonably good representation of all analyzed cultivars. PCR protocol [9] was as follows: 5 μL of Phusion Flash High-Fidelity PCR Master Mix (Thermo Fisher Scientific), 1 μL of each primer (final concentration 0.2 µM, one of them fluorescently labeled by 6-FAM, Supplementary Table S2), 2 μL DNA (10 ng/μL), and PCR grade water up to 10 μL. PCR was run on a C1000 PCR cycler (Bio-Rad, Hercules, CA, USA) with the following temperature profile: 98 °C/1 min; 23 cycles (98 °C/10 s, 60 °C/15 s, 72 °C/15 s); final extension 72 °C/30 s. Afterward, 1 μL of PCR product was mixed with 15 μL Hi-Di Formamide (Thermo Fisher Scientific) and 0.5 μL GeneScan 600 LIZ dye Size Standard v2.0 (Thermo Fisher Scientific). Samples were denatured at 95 °C for 2 min and run on a 3500 Genetic Analyzer (Thermo Fisher Scientific). Results were analyzed in v5 GeneMapper software (Thermo Fisher Scientific).”

To verify PCR specificity, all PCR fragments were sequenced. Amplicons were separated by gel electrophoresis, corresponding fragments were cut out, and they were purified by a GeneAll Expin Combo GP kit (GeneAll Biotechnology). Fragments were sequenced using a BigDye™ Terminator v3.1 Cycle Sequencing Kit (Thermo Fisher Scientific) according to the manufacturer’s instructions. Sequences were analyzed in Geneious Prime^®^ software (version 2025.0; GraphPad Software LLC d.b.a Geneious, Boston, MA, USA).

4.5. Statistical Analysis

In all analyzed populations/subpopulations, POLYGENE software, version V1.7 [44] was used to investigate general SSR marker statistical parameters, such as observed heterozygosity (Ho), expected heterozygosity (He), number of effective alleles (Ae), number of alleles (k), or polymorphic information content (PIC). In the case of PIC, two types were calculated: PIC_HipSTR_, which was calculated based on alleles differing by their sequence according to the in silico analysis of whole-genome sequences, and standard PIC, which was derived from fragment analysis and was based on the length of individual alleles. As such, PIC_HipSTR_ was expected to be the same or higher than PIC in cases that one allele is represented by more sequences with the same length (e.g., due to an SNP presence). Lengths of the shortest and longest alleles were determined for all markers. For final statistical analysis, only the 294 genotypes that were unique according to our new system were used; if two or more genotypes for samples with the same name were observed, only the prevailing genotype was used for statistical evaluation. In cases of only two available samples, one of the clones was randomly chosen. The probability of identity (PI_i_) was calculated according to the formula given in [45]:

$$[eqn]$$

where p_i_ and p_j_ are population allele proportions. The total probability of identity (PI_t_) over the 16 SSR markers used was computed as the product of the PI of each individual marker.

4.6. Development of the 16in1 SSR Marker Kit

To construct the final kit, two SSR markers per chromosome, separated by a minimally 1 Mbp region, were chosen based on statistical parameters and multiplexed into one PCR reaction. Primers with four fluorescent dyes (6-FAM, VIC, NED, and PET) were used for the PCR amplification of fragments with approximate expected lengths 100–150, 200–250, 300–350, and 400–450 nucleotides to allow unequivocal assigning of all alleles to individual markers. The multiplex PCR reaction was run under the same conditions as the simplex one, but with 25 PCR cycles. Due to the inefficient PCR amplification of some SSR markers under multiplex conditions, and in cases where no alternative high-quality SSR marker was identified in this study, a known SSR marker was used instead. The following SSR markers extracted from genetic maps on www.rosaceae.org (accessed on 30 April 2024) were tested using original primers: Chromosome 1: EMPA002, EMPA003, EMPA005, UCD-CH31, UDP98-022, Pchcms4; Chromosome 2: Pchgms1, UCD-CH11, UCD-CH12, UDAp-461; Chromosome 5: BPPCT037 (the only one tested due to very high Ho, He and PIC in [24]); Chromosome 8: CPDCT034, CPPCT006, EMPA018, PS01H03. The best SSR markers from each chromosome were chosen for the final multiplex based on good PCR amplifiability, high polymorphic content, and easily scoring of fragment analysis outputs. The final set of 16 SSR markers obtained is referred to as the 16in1 set. Concentrations of primers were fine-tuned to obtain comparable signal heights for all markers. If alleles with a low signal were identified, the corresponding primer-annealing regions were checked in whole-genome sequences, and primers were appropriately modified. The final kit was tested on the whole collection containing 315 samples to obtain detailed statistics for each marker. All varieties were tested at least twice, preferably from different sources. Finally, all alleles were verified as belonging to a particular marker by a simplex PCR and fragment analysis.

For comparison, the ECPGR recommended set of SSR markers for sweet cherry genotyping was analyzed as described by Clarke and Tobutt [19]. All SSR analyses were run at least in duplicate.

4.7. Population Analysis

POLYGENE software, version V1.7 [44], was used to analyze relationships in a collection containing only unique genotypes according to the 16in1 set (294 samples). Hierarchical clustering was performed using the unweighted pair group method with arithmetic mean (UPGMA) developed by Sokal and Michener [46] to obtain dendrograms.

The Bayesian model-based clustering method was used for the determination of the genetic structure of the analyzed population using the STRUCTURE v. 2.2.4 software [47]. Samples were ordered according to the 16in1 and ECPGR dendrograms produced by POLYGENE, respectively. Ten independent runs were conducted for each K (1-26) with a burn-in length of 200,000 and 500,000 MCMC reps after burn-in. A Web version of StructureSelector (http://lmme.qdio.ac.cn/StructureSelector/, accessed on 10 November 2024, [48]), which implements the method developed by Evanno et al. [49], was used to evaluate K values for the analyzed data. Varieties with inferred clusters ≥ 0.95 and ≥0.90, respectively, for the group with the highest and the second highest delta K values, were assumed as representative varieties for individual groups.

The tanglegram was constructed in R for the 294 genotypes unique in the 16in1 system. Genetic distance matrices were calculated based on codominant SSR data using custom functions in R (version 2.5-7) [50]. Genotypic data were converted into genind objects using the adegenet package (version 2.1.10) [51]. Hierarchical clustering was performed using the hclust() function with average linkage, and dendrograms were visualized and compared using the dendextend package (version 1.17.1) [52]. Topological similarity was assessed via cophenetic correlation and the Mantel test implemented in the vegan package [53].

Parentage analysis for 294 genotypes unique in the 16in1 system was run in POLYGENE software, version V1.7 [44] by the Phenotype method, looking for parent pairs (sex unknown); selfing was permitted. The top 5 LOD score pairs with a null score of Trio loci mismatching were extracted.

5. Conclusions

Though thousands of cherry varieties exist, only a limited number of ancestor varieties have been used for breeding. This poses a challenge for their reliable identification, as many groups of cherry varieties are genetically very similar. Therefore, a new SSR marker set was developed, validated, and formulated in a one-tube PCR assay. Our 16in1 system showed a discrimination power seven orders of magnitude better than the currently used SSR marker sets, making it especially useful for parentage and population structure analyses.