# Whole-genome screening for near-diagnostic genetic markers for four western European white oak species identification

Antoine Kremer, Adline Delcamp, Isabelle Lesur, Stefanie Wagner, Christian Rellstab, Erwan Guichoux, Thibault Leroy

► To cite this version:

Antoine Kremer, Adline Delcamp, Isabelle Lesur, Stefanie Wagner, Christian Rellstab, et al. Whole-genome screening for near-diagnostic genetic markers for four western European white oak species identification. Annals of Forest Science, 2024, 81 (1), pp.21. 10.1186/s13595-024-01236-9. hal-04640959

HAL Id: hal-04640959

https://hal.inrae.fr/hal-04640959v1

Submitted on 9 Jul 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

<figure><img src="image_3.png" /></figure>

# Whole-genome screening for near-diagnostic genetic markers for four western European white oak species identification

Antoine Kremer¹* ID, Adline Delcamp², Isabelle Lesur³ ID, Stefanie Wagner⁴,⁵ ID, Christian Rellstab⁶ ID, Erwan Guichoux²† ID and Thibault Leroy¹,⁷† ID

**Abstract**

**Key message** Mining genome-wide DNA sequences enabled the discovery of near-diagnostic markers for species assignment in four European white oaks (Quercus petraea (Matt.) Liebl., Quercus pubescens Willd., Quercus pyrenaica Willd., and Quercus robur L.) despite their low interspecific differentiation. Near-diagnostic markers are almost fully fixed in one species and absent in the three others. As a result, only a handful of markers are needed for species identification, making this genetic assay a very promising operational taxonomic assignment procedure in research and forestry.

**Context** Identifying species in the European white oak complex has been a long-standing concern in taxonomy, evolution, forest research, and management. Quercus petraea (Matt.) Liebl., Q. robur L., Q. pubescens Willd., and Q. pyr-enaica Willd. are part of this species complex in western temperate Europe and hybridize in mixed stands, challenging species identification.

**Aims** Our aim was to identify near-diagnostic single-nucleotide polymorphisms (SNPs) for each of the four species that are suitable for routine use and rapid diagnosis in research and applied forestry.

**Methods** We first scanned existing whole-genome and target-capture data sets in a reduced number of samples (training set) to identify candidate diagnostic SNPs, i.e., genomic positions being characterized by a reference allele in one species and by the alternative allele in all other species. Allele frequencies of the candidates SNPs were then explored in a larger, range-wide sample of populations in each species (validation step).

Handling editor: Marjana Westergren.

†Erwan Guichoux and Thibault Leroy contributed equally to the research.

Reference to preprint server The preprint version of this article is available in the bioRxiv server at https://www.biorxiv.org/content/10.1101/2023.11.29.568959v1.

This paper has undergone an Optional Open Peer Review. Read the Peer review report (Kremer et al. 2024b) at https://doi.org/10.17180/73QH-D818.

<figure><img src="image_4.png" /></figure>

<footnote>*Correspondence:</footnote>

<footnote>Antoine Kremer
antoine.kremer@inrae.fr</footnote>

<footnote>Full list of author information is available at the end of the article</footnote>
**Results** We found a subset of 38 SNPs (10 for *Q. petraea*, 7 for *Q. pubescens*, 9 for *Q. pyrenaica*, and 12 for *Q. robur*) that showed near-diagnostic features across their species distribution ranges with *Q. pyrenaica* and *Q. pubescens* exhibiting the highest (0.876) and lowest (0.747) diagnosticity, respectively.

**Conclusions** We provide a new, efficient, and reliable molecular tool for the identification of the species *Q. petraea*, *Q. robur*, *Q. pubescens*, and *Q. pyrenaica*, which can be used as a routine tool in forest research and management. This study highlights the resolution offered by whole-genome sequencing data to design near-diagnostic marker sets for taxonomic assignment, even for species complexes with relatively low differentiation.

**Keywords** Quercus, Diagnosticity, Genetic differentiation, Pool-seq data, Captured sequences

# 1 Introduction

Identifying species in the European white oak complex has been a long-standing concern in evolutionary biology as well as in forest research and management. According to the latest taxonomic classification, there are about 15 oak species in Europe, which form the subsection of the Roburoids within the *Quercus* section (white oak section) (Denk et al. 2017; Hipp et al. 2020). Within the continent, however, species richness varies, with higher species diversity in the Mediterranean region and in eastern Europe compared to other areas (Camus 1938; Le Hardy de Beaulieu and Lamant, 2006). In western temperate Europe, four white oaks species occur north of the Pyrenees and Alps (Q. petraea (Matt.) Liebl., Q. robur L., Q. pubescens Willd., and Q. pyrenaica Willd.). Co-occurrence of all four species in the same forest is rare. The few reported cases indicate extensive gene flow and admixture between all four species, leading to considerable morphological variations and uncertainties when it comes to taxonomic classification based on morphological characteristics (Lepais et al. 2013; Leroy et al. 2017; Viscosi et al. 2009). The co-occurrence of the three species *Q. petraea*, *Q. robur*, and *Q. pubescens* is more common, especially in the southern parts of the temperate range, for which hybridization and morphological variation are well documented (Dupouey and Badeau 1993; Grandjean and Sigaud 1987; Macejovsky et al. 2020; Rellstab et al. 2016). Finally, forests with co-occurrences of two species and interspecific admixture have also raised questions about species classification. This is especially true for co-occurrences of *Q. petraea-Q. robur* (Bacilieri et al. 1995; Jurksiene and Baliuckas 2014; Kelleher et al. 2005; Kremer et al. 2002; Yucedag and Gailing 2013) but also for *Q. petraea-Q. pyrenaica* (Lopez de Heredia et al. 2009), and *Q. petraea-Q. pubescens* (Bruschi et al. 2000; Reutimann et al. 2020, 2023). This brief overview of species admixture and problems of taxonomic classification based on morphological characteristics highlights the pressing need for a time- and cost-efficient molecular tool for reliable species assignment within European white oaks for use in forest science and management.

In response to this challenge, molecular tools have been continuously improved, and a number of species marker kits have been developed and applied during the last decade (Guichoux et al. 2011; Neophytou 2014; Reutimann et al. 2020; Degen et al. 2021; Schroeder and Kersten 2023). These methods have set new milestones for the delimitation of oak species, but their validity has been constrained by some biological and technical limitations. From a biological point of view, the markers used in the kits are still shared between the species, although interspecific differentiation of the selected markers was higher than in earlier studies. From a technical point of view, the genomic resources explored for selecting the marker candidates were very limited until recently. Using previously published genome-wide data and genome scans targeting genomic positions that maximize differentiation between populations of *Q. robur*, *Q. petraea*, *Q. pubescens*, and *Q. pyrenaica*, we overcame these limitations and designed a new single-nucleotide (SNP) marker set for range-wide species identification in European white oaks.

Earlier genome scans for species differences showed that interspecific differentiation ($F_{ST}$) followed an L-shaped distribution suggesting that there might be highly differentiated markers at an extremely low frequency within the genome (Reutimann et al. 2020; Scotti-Saintagne et al. 2004). Recent analysis of nucleotide diversity in genes underlying species barriers between European white oaks confirmed these expectations (Leroy et al. 2020b). Our approach built on these results by launching a systematic search of so-called species "diagnostic" SNPs within existing genome-wide resources. Ideally, a diagnostic SNP contains a diagnostic allele of a given species that is fully fixed in that species and the alternate allele fixed in the other species. Earlier surveys (Scotti-Saintagne et al. 2004; Reutimann et al. 2020; Lesur et al. 2018) in European white oaks indicated that such ideal cases rarely exist. However, some markers exhibit species frequency profiles close to the ideal case (so-called near-diagnostic SNPs, for example, an allele with a frequency larger than 0.9 in the target species, and alternate allele frequency larger than 0.9 in all other species) (Schroeder and Kersten 2023). Only a few of such
markers would then be enough to assign trees to the correct species using appropriate analytical approaches. For example, Reutimann et al. (2020) showed that five SNPs were enough for correctly classifying 95% of *Q. robur* reference trees, although the single SNPs were far from being diagnostic.

In this study, we explored pool-sequenced whole-genome libraries of natural populations of four white oak species (*Q. petraea*, *Q. pubescens*, *Q. pyrenaica*, and *Q. robur*) (Leroy et al. 2020b) and genome-wide capture-based sequences of *Q. petraea* and *Q. robur* (Lesur et al. 2018) to identify near-diagnostic SNPs for each of the four species. We describe the approaches and methods used to discover near-diagnostic SNPs and explore the stability of diagnosticity over the distribution range of the four species.

Our main goal was to identify and validate a new set of near-diagnostic SNPs that can be used in the development of an efficient and cost-effective molecular tool for forest research and management. To this end, we focused on the variation of near-diagnostic SNPs across species, between populations within each species and between SNPs. We finally addressed the evolutionary drivers that may have contributed to the maintenance and/or modification of diagnosticity within the genome and throughout the distribution range of the four species.

# 2 Material and methods

## 2.1 Discovery of near-diagnostic markers

The discovery of near-diagnostic SNPs was conducted by scanning oak genomic data that have been generated in earlier studies assessing genomic diversity and differentiation in the four sympatric white oak species (*Quercus petraea*, *Q. pubescens*, *Q. pyrenaica*, *Q. robur*; Leroy et al. 2017 and 2020b, Lesur et al. 2018).

### 2.1.1 Discovery of near-diagnostic SNPs in whole genome pool-sequenced (pool-seq) resources

*Pool sequencing* In Leroy et al. (2020b), we used leaf and bud samples from up to 20 adult trees of the four species coming from four different forests located at maximum 200 km away from each other in South-West of France (Table 4 in Appendix). The sampled stands were of mixed oak composition (generally two or three species) and of natural origin. DNA extracts were pooled in equimolar amounts to obtain a single pool for each species. Libraries were then sequenced on 9 to 10 lanes for each of the four species (1 pool per species) on a Illumina HiSeq 2000 sequencing platform (Leroy et al. 2020b for details). In this study, to reduce the computation load, we only used two lanes per pool from SRA, namely ERR2215923 and ERR2215924, ERR2215937 and ERR2215938, ERR2215909 and ERR2215910, and ERR2215916 and ERR2215917 for *Q*.

*pubescens*, *Q. petraea*, *Q. pyrenaica*, and *Q. robur*, respectively. Raw reads were then trimmed using Trimmomatic (v. 0.33, Bolger et al. 2014). Bases off the start and off the end of each read were removed if their quality fell below 3 (LEADING:3, TRAILING:3). Then, a sliding window was performed using SLIDINGWINDOW:4:15. Starting at the 5' end, the scanning cut the read once the average quality within the 4-bp window falls below 15. By considering multiple bases, a single quality base below the threshold of 15 will not cause the removal of high-quality data later in the read. Following these steps, only reads longer than 50 bp were kept for the analysis (MINLEN:50).

*Mapping and SNP calling* Data from two sequencing lanes per species (from up to 10 lanes per species in Leroy et al. 2020b) were then mapped against the oak haplome assembly ("PM1N," Plomion et al. 2018) using BWA-MEM (Li 2013). PCR duplicates were removed using Picard v. 1.140 (http://broadinstitute.github.io/picard/). Samtools v.1.1 (Li 2011) and PoPoolation2 v. 1.201 (Kofler et al. 2011) were then used to call biallelic SNPs with at least 10 reads of alternate alleles and a depth between 50 and 2000× at each position. To ensure a reasonably low rate of false positives due to Illumina sequencing errors, all SNPs with a MAF lower than 0.05 were discarded. A total of 24,345,915 SNPs were identified and then screened for their diagnostic value (see next paragraph).

*Genome scan for near-diagnostic SNPs* Allele frequencies were computed from the SNP-frequency-diff.pl script of PoPoolation2. SNPs exhibiting a high difference in allele frequency ($\Delta p$ > 0.9 between the focal species and all other species) were then selected. All candidate near-diagnostic SNPs with a coverage lower than 80 in the four populations were discarded, in order to ensure that the high $\Delta p$ was not associated with inaccurate allele frequency estimation in low coverage regions. Despite the relatively limited linkage disequilibrium in oaks (Coq-Etchegaray et al. 2023) even in species barrier regions (Leroy et al. 2020b), the relatively high nucleotide diversity in oaks (Plomion et al. 2018; Saleh et al. 2022) allows several neighboring SNPs to be identified by this screening. We therefore selected the best SNPs per identified region considering the constraints associated with the SNP design (see below).

### 2.1.2 Discovery of near-diagnostic SNPs in sequence-captured (seq-cap) genomic resources

In addition to the pool-seq resources, we mined a separate genome-wide resource that came from a sequence capture experiment of *Q. robur* and *Q*.
*petraea* aiming at calling SNPs for inferring genomic relatedness among trees (Lesur et al. 2018). Here, the discovery population consisted of a far larger panel (245 adult trees in total) equally distributed between *Q. petraea* and *Q. robur* growing in the Petite Charnie forest located in the western part of France (Table 5 in Appendix). We used the capture data in complement of the pool-seq data to ensure a higher diagnosticity of the markers for this specific pair, given the larger panel of *Q. robur* and *Q. petraea* samples available in the capture data. The capture-based assay consisted in sequencing 2.9 Mb (15,623 target regions) on an Ion Proton System (Thermo Fisher, Scientific, Waltham, MA, USA) covering both genic and intergenic regions and resulted in the calling of more than 190,000 SNPs with a coverage of more than 10× (Lesur et al. 2018). The study provided allele frequencies of each SNP, and we screened the total set of SNPs for their differentiation between *Q. petraea* and *Q. robur*, by ranking their *F*<sub>ST</sub> values to complete the discovery panel. Although limited to two species (*Q. petraea* and *Q. robur*), this data set corresponded also to a genome-wide exploration of species differentiation implemented on a larger population sample (Lesur et al. 2018). It was therefore selected for this study, pending its relevance for selecting near-diagnostic markers for the remaining two species (*Q. pubescens* and *Q. pyrenaica*), which is investigated in this study.

## 2.2 Training and validation of near-diagnostic SNPs

### 2.2.1 Training populations

The candidate diagnostic SNPs of the discovery panel were first tested on a limited number of oak individuals (training set), randomly sampled in 7 to 14 populations per species (in total 19 to 48 sampled trees per species, Table 6 in Appendix). The training experiment was conducted over two sessions that took place during two periods (training 1 and training 2, Table 6 in Appendix) with different samples (but from the same geographic range). The two sessions differed only by the samples included which was constrained by the availability of the material. The objective of the training step was to check whether the candidate SNPs exhibited near-diagnostic frequency profiles in natural populations originating mainly from the area of the discovery panel. The training step also included quality controls and repeatability assessments of the genotyping assay (see results paragraph 3.3.1).

### 2.2.2 Validation populations

Given that the discovery and training of diagnostic SNPs were implemented on limited number of trees originating mostly from the western part of the distribution of

the four species, we included a round of validation by increasing the sample sizes of the training populations and enlarging the collection of populations, studying the SNP diagnosticity across a larger part of the species' natural distribution (Fig. 1). Additionally, the validation step aimed also at reducing the number of SNPs, while still maintaining overall multilocus diagnosticity, in order to produce a low-cost and easy-to-use screening tool in operational forestry. In total, 24 populations of *Q. petraea*, 10 of *Q. pubescens*, 6 of *Q. pyrenaica*, and 19 of *Q. robur* were part of the validation set, representing in total 977 trees (Fig. 1 and Table 6 in Appendix). All samples were collected in natural populations, and their taxonomic status was assessed by the local collectors based on leaf morphology and were not checked after their receipt. A few misclassifications were detected after the molecular assay (Arbalan, Killarney, Dalkeith). Sampled populations were in most cases of mixed oak composition. Some of the populations were used in earlier large-scale genetic surveys (Gerber et al. 2014; Kremer et al. 2002); others were purposely collected for this study.

### 2.2.3 Genotyping assay

Medium-throughput SNP genotyping assays were implemented on single tree DNA extracts using the MassARRAY® technology (Agena Bioscience, San Diego, CA, USA). The assay design, using the MassARRAY Assay Designer version 4.0.0.2, was performed on candidate SNPs from pool-seq and seq-cap resources. Nine multiplexes, for a total of 359 SNP (eight 40-plex and one 39-plex), were designed for identifying the best markers. Genotyping was performed using iPLEX Gold chemistry following Ellis and Ong (2017) on a MassARRAY System CPM384 (Agena Biosciences) at the PGTB platform (https://doi.org/10.15454/1.5572396583599417E12). Data analysis was achieved using MassARRAY Typer Analyzer 4.0.4 (Agena Biosciences). After genotyping, we excluded all markers for which there was evidence that the candidate SNP identified during the discovery step was not recovered, for example, when the SNP exhibited fixation across the four species at the same allele. We also discarded loci with weak (magnitude < 5) or ambiguous signal (i.e., displaying more clusters than expected or unclear cluster delineation) and loci with more than 20% missing data. Following this selection process, 61 SNPs (in two multiplexes) were selected on the basis of their near-diagnostic value and their compatibility in one multiplex kit for subsequent genotyping on all the samples.

### 2.2.4 Diagnosticity of candidate SNPs

Standard genetic statistics (allele frequencies, diversity statistics, differentiation, and fixation indices) were
<figure><img src="image_5.png" /><figcaption>**Fig. 1** Geographic distribution of the validation populations. The green area corresponds to the distribution of the species according to Caudullo et al. (2017). a Quercus petraea. b Q. pubescens. c Q. pyrenaica. d Q. robur. Red dots correspond to the origins of the validation populations. Populations identified by their name refer to populations for which frequency profiles of near-diagnostic alleles are later illustrated and discussed (paragraphs 3.3.3 and 3.3.4)</figcaption></figure>

estimated using GENEPOP (Raymond and Rousset 1995) and ADEGENET software (Jombart 2008). Genotypic arrays of sampled trees in training and validation population are available at https://doi.org/10.57745/0JYLZU (Kremer et al. 2024a).

We defined a metric of species diagnostic accuracy, which we coined “diagnosticity” index (D) to screen SNP alleles for their ability to be close to full diagnosticity.

Full diagnosticity requires two properties: fixation of the diagnostic allele in the target species and fixation of the alternate allele in the remaining species. These two properties are included in the metric D. Considering a set of n species, diagnosticity of an allele for species x ($D_x$) regarding the remaining (n-1) species could be expressed as follows:

$$D_x = p_x - \frac{1}{n-1} \sum_{j=1}^{n-1} p_j$$

where $p_x$ is the frequency of the candidate diagnostic allele in the target species x and $p_j$ the frequency of the same allele in the alternate species j. $D_x$ amounts to the difference of allelic frequencies between species x and the remaining (n-1) species. $D_x$ is equivalent to the mean Gregorius genetic distance between species x and the three other species for a diallelic locus (Gregorius 1984).

$D_x$ has two components, which account for the two properties of diagnosticity.

* $p_x$: The higher the $p_x$, the closer the near-diagnostic allele to fixation in the target species.

* $\frac{1}{n-1} \sum_{j=1}^{n-1} p_j$: The lower the mean value of $p_j$, the closer the alternate allele to fixation in the remaining (n-1) species.

$D_x$ is more appropriate for practical diagnostic assessments than the traditional differentiation metric $F_{ST}$ when more than two species are involved (see Gregorius and Roberts (1986) for a comparison of D and $F_{ST}$). To illustrate the discrepancy between D and $F_{ST}$ regarding diagnosticity, consider the case of four species with frequency profiles ($p_1=1, p_2=1, p_3=0, p_4=0$). Addressing
diagnosticity for species 1, *F*<sub>ST</sub> would yield 1, while *D*<sub>1</sub> would yield 0.67. *D*<sub>1</sub> accounts for the lack of frequency differences between species 1 and 2, while *F*<sub>ST</sub> does not.

By extension of the definition of a diagnostic allele, a near-diagnostic SNP is a SNP bearing near-diagnostic alleles, and diagnosticity of a species (or a population of that species) refers to the mean value of all near-diagnostic SNPs assessed for that species or population. Diagnosticity of candidate SNPs is estimated in the training and validation populations.

### 2.2.5 Multilocus species clustering

To validate the selected near-diagnostic SNP for a multi-locus species assignment procedure, we implemented an empirical clustering approach using principal component analysis, free of any underlying evolutionary assumptions (ADEGENET, Jombart 2008). This method allows to check for the ability of the near-diagnostic SNPs to visually discriminate the four species.

## 3 Results

### 3.1 Discovery of near-diagnostic SNPs

All together, we recovered 61 candidate near-diagnostic alleles, 49 originating from the pool-seq study, and 12 from the seq-cap analysis (Table 7 in Appendix). The candidate SNPs are distributed over all chromosomes (except chromosome 4), and their number ranges from 1 (chromosomes 3, 9, and 12) to 17 (chromosome 2, Fig. 2).

In a few cases, near-diagnostic markers of a given species clustered in pairs in a few spots (mainly for *Q. robur* on chromosome 2, 5, 6, Fig. 2). In such cases, one marker of the pair was discarded during the validation step. Near-diagnostic markers are distributed over six chromosomes for *Q. petraea*, *Q. pubescens*, and *Q. pyrenaica* and over eight chromosomes for *Q. robur*. As indicated by their location on the chromosomes, the minimum physical distance of near-diagnostic SNPs located on the same chromosomes was 17 Kb (Table 7 in Appendix). All except two SNPs are located on scaffolds that are anchored on the pseudo-chromosome assembly of the oak genome as shown in Fig. 2.

### 3.2 Diagnosticity of candidate SNPs in the training set

The 61 candidate near-diagnostic SNPs exhibited allele frequency profiles close to the requisite properties of a diagnostic SNP but did not fulfill entirely criteria of full diagnosticity (Fig. 2, Fig. 6 in Appendix). *D* values indeed varied between 0.283 and 0.963. Most of the near-diagnostic SNPs (92%, 56/61) exhibited *D* scores greater than 0.50 (mean value 0.758). Among the 61 SNPs, 16 are candidate near-diagnostic of *Q. petraea*, 11 of *Q. pubescens*, 12 of *Q. pyrenaica*, and 22 of *Q. robur*.

Diagnosticity scores were higher in the pool-seq uncovered set (*D*=0.771) than in the seq-cap uncovered set (*D*=0.704).

<figure><img src="image_6.png" /><figcaption>**Fig. 2** Genomic location of the near-diagnostic SNPs on the 12 oak (pseudo-)chromosomes of the oak genome. The color code of the marker corresponds to the species name for which the SNP is expected to be diagnostic, with *Q. robur*, *Q. petraea*, *Q. pubescens*, and *Q. pyrenaica*, shown in pink, green, blue, and yellow, respectively. For each SNP, the diagnosticity of each marker at the training stage is indicated following the proportional and color scale shown. Thin and bold lines both indicate the location of the SNPs but separate SNPs that were excluded or included in our final set of 38 SNPs, respectively. Note that two near-diagnostic SNPs are not shown since they are located on scaffolds that are not anchored on the oak pseudo-chromosomes (see Table 7 in Appendix)</figcaption></figure>
Concerning the near-diagnostic SNPs identified with the pool-seq data, diagnosticity was highest for *Q. pyrenaica* (0.897) and *Q. robur* (0.780) and lower in *Q. petraea* (0.736) and *Q. pubescens* (0.657). Deviations to full diagnosticity in the two latter species are associated with different patterns (Fig. 6 in Appendix):

* Lower diagnosticity in *Q. petraea* was mostly related to the sharing of the near-diagnostic allele with the other species, especially with *Q. pubescens*.

* Lower diagnosticity for *Q. pubescens* was mainly due to three SNPs (Sc0000170_630013, Sc0000192_329301, and Sc0000482_334917) that showed substantial deviation from fixation within *Q. pubescens* (frequency being respectively 0.468, 0.587, 0.283), while the alternate alleles were fixed in the three other species.

Concerning the seq-cap uncovered SNPs, we selected 12 SNPs that exhibited the highest species differentiation in the Petite Charnie population. As expected, all 12 SNPs showed strong frequency differences between *Q. petraea* and *Q. robur* in our training panel. Eight out of the 12 SNPs exhibited allele frequency differences among the four species consistent with diagnosticity requirements for four species, with the near-diagnostic marker being almost fixed in the reference diagnostic species and present at very low frequencies in all the three remaining species (Fig. 6 in Appendix). The four remaining candidate SNPs exhibited near-diagnostic alleles being almost fixed, not only in one but in two species:

* Sc0000040_1694351 in *Q. petraea* and *Q. pubescens*

* Sc0000481_366275 in *Q. robur* and *Q. pyrenaica*

* Sc0000546_456229 in *Q. robur* and *Q. pyrenaica*

* Sc0000598_295142 in *Q. robur* and *Q. pyrenaica*

## 3.3 Validation of the near-diagnostic SNPs

### 3.3.1 Screening of near-diagnostic SNPs

The validation step aimed at verifying the diagnosticity of the candidate SNPs on a larger geographic scale while at the same time optimizing the assay by selecting the best SNPs according to various genetic and technical criteria. We thus attempted to optimize the MassARRAY® genotyping assays by reducing the number of near-diagnostic SNPs and combine them in one final assay, without limiting the species assignment purpose and reducing its diagnosticity. Indeed, given the frequency profiles of near-diagnostic alleles we observed in the training set (Fig. 6 in Appendix), the required number of near-diagnostic SNPs for species

assignment can be limited to a handful of markers (Reutimann et al. 2020). We aimed at selecting about 10 near-diagnostic SNPs per species for the final design of the operational assay. The following criteria were applied (Table 7 in Appendix):

* Repeatability and clarity of the cluster delimitation on the scatter plots

* Diagnosticity of SNPs

* A nearly equal number of near-diagnostic SNPs per species

Combining the remaining SNPs within one or two multiplex sets resulted in amplification incompatibilities among SNPs which lead us to discard additional SNPs. Finally, a total of 10 near-diagnostic SNPs were selected for *Q. petraea*, 7 for *Q. pubescens*, 9 for *Q. pyrenaica*, and 12 for *Q. robur* (Table 8 in Appendix).

### 3.3.2 Allele frequency profiles of near-diagnostic SNPs in the validation populations

Overall, the average diagnosticity of the 38 near-diagnostic SNPs was slightly higher in the validation than in the training populations, with the exception of *Q. pyrenaica* (Fig. 3, Fig. 6 in Appendix): 0.784 (validation) vs 0.715 (training) in *Q. petraea*, 0.747 vs 0.690 in *Q. pubescens*, 0.876 vs 0.897 in *Q. pyrenaica*, and 0.841 vs 0.758 in *Q. robur*. The lower diagnosticity of *Q. pyrenaica* in the validation set (vs the training set) was due to SNP Sc0000307_852597, which exhibited contrasting values between the training (0.753) and validation set (0.546) (Table 8 in Appendix).

However, the validation populations provided the opportunity to explore the stability of the allele frequency profiles across geographic regions and thus addressed the maintenance of diagnosticity of individual SNPs across the distribution of the four species. Most near-diagnostic SNPs exhibited larger genetic differentiation between populations within a given species than usually found (Scotti-Saintagne et al. 2004) in oak species (Tables 1, 2, 3). Mean intraspecific $F_{ST}$ values of near-diagnostic SNPs amounted to 0.104, 0.192, 0.042, and 0.104 for *Q. petraea*, *Q. pubescens*, *Q. pyrenaica*, and *Q. robur*, respectively. Furthermore, $F_{ST}$ values within a species exhibited large variation among SNPs. For example, $F_{ST}$ values of near-diagnostic SNPs of *Q. petraea* between *Q. petraea* populations varied between 0.012 and 0.252. *Quercus pyrenaica* is an exception to these general rules, as the mean $F_{ST}$ (0.042) is much lower than for the three other species and the range of variation reduced (−0.022 to 0.142, data not shown).
<figure><img src="image_7.png" /><figcaption>**Fig. 3** Heat map of frequencies of near-diagnostic alleles and diagnosticity in the validation populations of the four species. SNPs are clustered for their diagnostic value for each species (reference species). First to fourth columns correspond respectively to near-diagnostic alleles of *Q. petraea*, *Q. pubescens*, *Q. pyrenaica*, and *Q. robur*</figcaption></figure>

<table><caption>**Table 1** Frequencies and differentiation of near-diagnostic alleles of *Q. petraea* in *Q. petraea* populations</caption><thead><tr><th>SNP ID</th><th>Tronçais</th><th>Lappwald</th><th>Bezange</th><th>Pomieri</th><th>Aspromonte</th><th>Killarney</th><th>Montejo</th><th>Intraspecific F<sub>ST</sub></th><th>p-value</th></tr></thead><tbody><tr><td>Sc0000254_6223</td><td>0.867</td><td>1.000</td><td>0.883</td><td>0.975</td><td>-</td><td><strong>0.625</strong></td><td>-</td><td>0.098</td><td>0.000</td></tr><tr><td>Sc0000121_355205</td><td>0.974</td><td>1.000</td><td>0.988</td><td>0.947</td><td>-</td><td>0.906</td><td>-</td><td>0.029</td><td>0.002</td></tr><tr><td>Sc0000043_1651618</td><td>0.817</td><td>0.925</td><td>0.888</td><td><strong>0.350</strong></td><td><strong>0.344</strong></td><td><strong>0.719</strong></td><td><strong>0.763</strong></td><td>0.190</td><td>0.000</td></tr><tr><td>Sc0000083_147504</td><td>0.950</td><td>0.900</td><td>0.898</td><td>0.875</td><td>0.979</td><td><strong>0.656</strong></td><td><strong>0.776</strong></td><td>0.049</td><td>0.000</td></tr><tr><td>Sc0000118_1466708</td><td>0.933</td><td>0.950</td><td>0.929</td><td>0.800</td><td>0.917</td><td><strong>0.688</strong></td><td><strong>0.671</strong></td><td>0.052</td><td>0.000</td></tr><tr><td>Sc0000135_261350</td><td>0.948</td><td>0.800</td><td>0.929</td><td><strong>0.425</strong></td><td><strong>0.573</strong></td><td>0.875</td><td>0.974</td><td>0.133</td><td>0.000</td></tr><tr><td>Sc0000145_700044</td><td>0.983</td><td>1.000</td><td>1.000</td><td>0.850</td><td>0.927</td><td>0.875</td><td>1.000</td><td>0.039</td><td>0.000</td></tr><tr><td>Sc0000203_707735</td><td>1.000</td><td>1.000</td><td>0.981</td><td>0.875</td><td>0.979</td><td>0.800</td><td><strong>0.622</strong></td><td>0.176</td><td>0.000</td></tr><tr><td>Sc0000274_909817</td><td>0.933</td><td>0.875</td><td>0.949</td><td><strong>0.325</strong></td><td><strong>0.333</strong></td><td><strong>0.688</strong></td><td>0.816</td><td>0.259</td><td>0.000</td></tr><tr><td>Sc0000481_343721</td><td>0.983</td><td>0.900</td><td>0.949</td><td>0.816</td><td>0.958</td><td>0.900</td><td>0.973</td><td>0.012</td><td>0.063</td></tr><tr><td>Mean diagnosticity</td><td>0.852</td><td>0.848</td><td>0.852</td><td>0.637</td><td>0.664</td><td>0.686</td><td>0.737</td><td></td><td></td></tr></tbody></table>

Populations in bold characters exhibit diverging diagnosticity (DD populations), and frequencies in bold characters correspond to loci exhibiting deviations to expected frequencies of near-diagnostic alleles. Geographic locations of the populations are shown in Fig. 1a

### 3.3.3 Allele frequency profiles of near-diagnostic SNPs in *Q. petraea* populations

We examined the geographic distribution of near-diagnostic alleles between populations within a given species. To illustrate the results, we selected populations that are

representative of the variation observed among all populations. We first selected a few widely distributed populations that exhibited allele frequencies at all SNPs close to the expected diagnosticity ("ED populations": Tronçais, Lappwald, and Bézange), and added all the populations
<table><caption>**Table 2** Frequencies and differentiation of near-diagnostic alleles of *Q. pubescens* in *Q. pubescens* populations</caption><thead><tr><th>SNP ID</th><th>Auros</th><th>Pantano</th><th>Switzerland<sup>a</sup></th><th>Ventoux</th><th>Intraspecific F<sub>ST</sub></th><th>p-value</th></tr></thead><tbody><tr><td>Sc0000314_149731</td><td>0.821</td><td>-</td><td><strong>0.531</strong></td><td><strong>0.521</strong></td><td>0.113</td><td>0.003</td></tr><tr><td>Sc0000047_2398879</td><td>0.946</td><td>1.000</td><td>0.694</td><td>0.868</td><td>0.127</td><td>0.000</td></tr><tr><td>Sc0000088_1796044</td><td>1.000</td><td>1.000</td><td>0.806</td><td>0.812</td><td>0.128</td><td>0.000</td></tr><tr><td>Sc0000109_800763</td><td>0.839</td><td>0.882</td><td><strong>0.528</strong></td><td><strong>0.692</strong></td><td>0.096</td><td>0.000</td></tr><tr><td>Sc0000111_693153</td><td>1.000</td><td>0.987</td><td>1.000</td><td>1.000</td><td>-0.004</td><td>0.682</td></tr><tr><td>Sc0000170_630013</td><td>0.907</td><td>0.987</td><td><strong>0.222</strong></td><td><strong>0.400</strong></td><td>0.506</td><td>0.000</td></tr><tr><td>Sc0000192_329301</td><td>0.880</td><td>0.986</td><td><strong>0.333</strong></td><td><strong>0.487</strong></td><td>0.376</td><td>0.000</td></tr><tr><td>Mean diagnosticity</td><td>0.865</td><td>0.926</td><td><strong>0.540</strong></td><td><strong>0.635</strong></td><td></td><td></td></tr></tbody></table>

Populations in bold characters exhibit diverging diagnosticity (DD populations), and frequencies in bold characters correspond to loci exhibiting deviations to expected frequencies of near-diagnostic alleles. Geographic origins of populations are shown in Fig. 1b

<table><caption>**Table 3** Frequencies of near-diagnostic alleles and differentiation of *Q. robur* in *Q. robur* populations</caption><thead><tr><th>SNP ID</th><th>Zivinice</th><th>Sigmundsherberg</th><th>Charnie</th><th>Escherode</th><th>Pedro</th><th>Roudsea<sup>a</sup></th><th>Intraspecific F<sub>ST</sub></th><th>p-value</th></tr></thead><tbody><tr><td>Sc0000013_2578823</td><td>0.711</td><td>0.925</td><td>0.875</td><td>0.750</td><td><strong>0.273</strong></td><td><strong>0.553</strong></td><td>0.140</td><td>0.000</td></tr><tr><td>Sc0000053_1639108</td><td>1.000</td><td>1.000</td><td>1.000</td><td>0.857</td><td>0.818</td><td><strong>0.500</strong></td><td>0.098</td><td>0.009</td></tr><tr><td>Sc0000099_1839376</td><td>0.868</td><td>0.868</td><td>0.889</td><td>0.972</td><td>0.818</td><td>0.711</td><td>0.040</td><td>0.000</td></tr><tr><td>Sc0000158_462639</td><td>0.763</td><td>0.975</td><td>0.889</td><td>0.806</td><td>0.864</td><td><strong>0.658</strong></td><td>0.071</td><td>0.000</td></tr><tr><td>Sc0000203_689887</td><td>0.842</td><td>0.750</td><td>1.000</td><td>0.806</td><td>0.955</td><td><strong>0.632</strong></td><td>0.095</td><td>0.000</td></tr><tr><td>Sc0000339_4638</td><td>1.000</td><td>1.000</td><td>0.944</td><td>0.941</td><td><strong>0.227</strong></td><td>0.853</td><td>0.359</td><td>0.000</td></tr><tr><td>Sc0000381_206331</td><td>1.000</td><td>1.000</td><td>1.000</td><td>1.000</td><td><strong>0.500</strong></td><td><strong>0.737</strong></td><td>0.233</td><td>0.000</td></tr><tr><td>Sc0000447_521057</td><td>1.000</td><td>0.950</td><td>0.944</td><td>1.000</td><td>1.000</td><td>0.816</td><td>0.042</td><td>0.001</td></tr><tr><td>Sc0000517_258593</td><td>0.938</td><td>0.975</td><td>0.938</td><td>0.917</td><td>1.000</td><td>0.974</td><td>-0.003</td><td>0.166</td></tr><tr><td>Sc0000695_225347</td><td>0.947</td><td>1.000</td><td>0.889</td><td>0.917</td><td><strong>0.773</strong></td><td><strong>0.763</strong></td><td>0.050</td><td>0.000</td></tr><tr><td>Sc0000796_82698</td><td>0.842</td><td>0.925</td><td>0.944</td><td>0.972</td><td>0.864</td><td>0.789</td><td>0.021</td><td>0.016</td></tr><tr><td>Sc0000967_33996</td><td>1.000</td><td>1.000</td><td>1.000</td><td>1.000</td><td>1.000</td><td><strong>0.763</strong></td><td>0.106</td><td>0.000</td></tr><tr><td>Mean diagnosticity</td><td>0.851</td><td>0.889</td><td>0.885</td><td>0.854</td><td>0.700</td><td>0.671</td><td></td><td></td></tr></tbody></table>

Populations in bold characters exhibit diverging diagnosticity (DD populations), and frequencies in bold characters correspond to loci exhibiting deviations to expected frequencies of near-diagnostic alleles. Geographic origins of populations are indicated in Fig. 1d

that deviate from the ED population frequency profiles, which we called populations with diverging diagnosticity (“DD populations”). The DD populations included three extreme southern populations (Pomieri and Aspromonte in Italy, Montejo in Spain) and one population from the northern distribution edge (Killarney). All the remaining *Q. petraea* populations exhibited frequency profiles similar to the selected ED populations and are not shown in Table 1 and in Fig. 4. While the ED populations exhibited almost full fixation in all near-diagnostic SNPs, the DD populations showed substantial polymorphism (i.e., lower diagnosticity) at a few SNPs in Pomieri and Aspromonte (Sc0000043_1651618, Sc0000135_261350, Sc0000274_909817), and moderate polymorphism distributed among more SNPs in Killarney and Montejo.

Additionally, we examined the occurrences of near-diagnostic alleles of the other three species in *Q. petraea* populations (Fig. 4). The DD *Q. petraea* populations exhibited elevated frequencies of near-diagnostic alleles of *Q. pubescens* (Pomieri and Aspromonte) or *Q. robur* (Killarney and Montejo). Thus, they deviated also in respect to the expected frequencies of near-diagnostic alleles of the three other species, contrary to the ED populations (Fig. 4).

### 3.3.4 Allele frequency profiles of near-diagnostic SNPs in *Q. pubescens*, *Q. robur*, and *Q. pyrenaica* populations

To illustrate the intraspecific differentiation of near-diagnostic SNPs in the other three species, we followed the same procedure as for *Q. petraea*. We selected for each

<footnote><sup>a</sup> Switzerland population assembles data of populations Ayent, Cordola, Remigen, and Saillon described in Table 6 in Appendix</footnote>

<footnote><sup>a</sup> Roudsea population assembles data of populations Dalkeith and Roudsea described in Table 6 in Appendix</footnote>
<figure><img src="image_8.png" /><figcaption>**Fig. 4** Frequencies of near-diagnostic alleles of *Q. pubescens*, *Q. pyrenaica*, and *Q. robur* in *Q. petraea* populations. Populations in red and green exhibit respectively diverging diagnosticity (DD populations) and expected diagnosticity (ED populations). Shown are the mean frequencies of all near-diagnostic alleles of a given species (*Q. pubescens*, *Q. pyrenaica*, *Q. robur*) in *Q. petraea* populations. Geographic locations of *Q. petraea* populations are shown in Fig. 1a</figcaption></figure>

species two sets of populations: a subset of populations exemplifying the pattern close to full fixation of near-diagnostic loci at all SNPs (ED populations) and the set of populations that exhibited deviations to this trend (DD populations).

In the case of *Q. pubescens*, the DD populations (Switzerland and Ventoux) were located at the central northern edge of distribution. These deviations were not evenly distributed across the seven near-diagnostic SNPs of *Q. pubescens* but restricted to the same loci in the two populations (Table 2). The two populations Switzerland and Ventoux exhibited also higher frequencies of *Q. petraea* near-diagnostic alleles, in comparison to the two ED populations (Fig. 7 in Appendix).

In the case of *Q. robur*, there were also two DD populations located at the southwestern (Pedro) and north-western margin of the distribution (Roudsea) (Table 3). These two populations comprised also larger frequencies of near-diagnostic alleles of other white oak species (*Q. pubescens* and *Q. pyrenaica* in the case of Pedro; *Q. petraea* in the case of Roudsea) (Fig. 8 in Appendix). Finally, in *Q. pyrenaica*, all populations behave as ED populations (data not shown), e.g., all *Q. pyrenaica* populations exhibited frequency profiles similar to those shown for *Q. pyrenaica* in Fig. 3 and Table 8 in Appendix.

## 3.4 Multilocus structure of near-diagnostic SNPs

We used a principal component analysis (PCA) in the validation populations to assess and illustrate species differentiation (Fig. 5). We added 13 samples of known

first-generation hybrid origin to the species samples. Ten samples resulted from controlled interspecific crosses, and three came from parentage analysis conducted in a mixed *Q. petraea*-*Q. robur* stand (Truffaut et al. 2017). A combination of the three first components allowed to visually differentiate the four different species. While principal component 1 differentiated mainly *Q. petraea* and *Q. robur* (Fig. 5a), component 3 distinguished *Q. pyrenaica* from the three other species (Fig. 5b), and the biplot of components 2 and 3 provided the best visual separation between *Q. pubescens* and *Q. petraea* (Fig. 5c).

These multilocus representations showed that there is a small number of samples located at intermediate positions, especially between *Q. petraea* and *Q. robur* (Fig. 5a) and between *Q. petraea* and *Q. pubescens* (Fig. 5c). These regions of the PCA are also occupied by known interspecific hybrids, suggesting that the species samples, although identified as pure species in the field, represent either hybrids or introgressed forms. These intermediate positions are also preferentially occupied by trees belonging to diverging populations, as shown by the targeted PCA analysis on the two pairs of species sharing intermediate samples: *Q. petraea* and *Q. pubescens* (Fig. 9 in Appendix) and *Q. petraea* and *Q. robur* (Fig. 10 in Appendix).

## 4 Discussion

We explored large-scale existing genomic resources in four European white oaks of the subsection Roburoid (*Q. petraea*, *Q. pubescens*, *Q. pyrenaica*, *Q. robur*) to screen
<figure><img src="image_9.png" /><figcaption>**Fig. 5** Biplot of principal components of tree samples based on a principal component analysis (PCA) conducted in the validation populations. **a** Biplot of components 1 and 2. **b** Biplot of components 1 and 3. **c** Biplot of components 2 and 3. Numbers between brackets stand for the percentage of variation explained by the component. Red dots, *Q. petraea* samples; orange dots, *Q. pubescens* samples; green dots, *Q. pyrenaica* samples; blue dots, *Q. robur* samples; black dots, hybrids, *Pet*\*Pub, *Q. petraea*\*\*Q. pubescens* hybrids; *Pet*\*Rob, *Q. petraea*\*\*Q. robur* hybrids; *Pub*\*Pyr, *Q. pubescens*\*\*Q. pyrenaica* hybrids; *Pub*\*Rob, *Q. pubescens*\*\*Q. robur* hybrids. In Figure 5c, *Q. robur* samples are "hidden" by *Q. petraea* samples</figcaption></figure>

their genomes for near-diagnostic SNPs that could be
used for molecular fingerprinting (species and hybrid
identification) in forest research and operational forestry,
as wood or seed traceability in the wood chain and in
forest nurseries. Despite the widely reported low inter-
specific genetic differentiation among European white
oak species, we were able to identify a subset of SNPs
that exhibited near-diagnostic features across their spe-
cies' distribution ranges. Moreover, mutlivariate analysis
showed that these markers can be used for reliable hybrid
detection and accurate quantification of admixture lev-
els. However, diagnosticity varied substantially among
species, among populations within species, and among
SNPs. In the following, we discuss these variations in
relation to the known evolutionary history and genetic
interactions among and within the four species.

**4.1 Variation of diagnosticity among species**

Diagnosticity was highest in *Q. pyrenaica* (0.876) and lowest in *Q. pubescens* (0.747) with *Q. robur* and *Q. petraea* showing intermediate values. Near-diagnostic SNPs are likely located in genomic regions that exhibit larger divergence and/or regions less permeable to interspecific gene flow. The range of diagnosticity among the four species may therefore reflect the variation of divergence time and/or the variation of the intensity of gene flow during the ongoing interglacial period.

It is striking to notice that higher and lower diag-
nosticity was observed for species that showed the
older (Q. pyrenaica, Q. robur) and more recent (Q.
petraea, Q. pubescens) divergence, respectively (Leroy
et al. 2017). Fixation of near-diagnostic SNPs in spe-
cies with large population sizes as in oaks requires
long time periods. Consequently, lower diagnosticity
is likely associated with species that diverged more
recently. This is illustrated by Q. pubescens, which
shows lower diagnosticity due to the higher sharing of
near-diagnostic alleles with Q. petraea than with the
other two species (Fig. 3 and Fig. 6 in Appendix). Diag-
nosticty may in addition be dependent on the variation
of population size (Ne) among species and along diver-
gence, for which we lack any estimation today. Our
results may therefore be revisited in the light of future
evidence of Ne differences. Regarding gene flow, we
showed earlier that the four species came into contact
only recently, during the late last glacial maximum,
after being isolated for most of their earlier his-
tory (Leroy et al. 2020b, 2017), resulting in gene flow
among species. While interfertility among the four
species has been shown experimentally by controlled
crosses (Lepais et al. 2013), hybridization in natura
has also been observed among the four species in rare
mixed forests where all four species co-occur (Lepais
and Gerber 2011; Lepais et al. 2009). Interspecific mat-
ings of Q. pyrenaica in controlled crosses with the
remaining three species were quite successful; how-
ever, occurrences of natural hybridization were less
frequent due to the very late flowering of Q. pyrenaica
in comparison to the three other species (Lepais and
Gerber 2011; Lepais et al. 2013). Furthermore, Q. pyr-
enaica is mainly distributed in southwestern Europe,
where the other three species are only present in scat-
tered forests, leading, for example, to reported but rare
hybridization with Q. petraea (Valbuena-Carabana
et al. 2005) and Q. robur (Moracho et al. 2016). Alto-
gether, phenological prezygotic barriers and limited
overlapping distributions with the other three species
may have contributed to reduced genetic exchanges
between Q. pyrenaica and the other three species and
thus account for the high diagnosticity of the SNPs in
of Q. pyrenaica. In contrast to Q. pyrenaica, no repro-
ductive barriers were observed in Q. pubescens when
crosses were made with Q. petraea as female parent, as
interspecific crosses were as successful as intraspecific
crosses (Lepais et al. 2013). Reduced barriers between
these two species were corroborated by frequent
admixture detected in genetic surveys conducted in
mixed stands of Q. pubescens and Q. petraea (Alberto
et al. 2010; Neophytou 2014; Reutimann et al. 2023).
As a result, near-diagnostic SNPs of Q. pubescens and
Q. petraea were more frequently shared between the
two species (Fig. 3 and Fig. 6 in Appendix), thus con-
tributing to reduced diagnosticity. Finally, interspecific
gene exchanges involving Q. robur were mainly inves-
tigated with regard to Q. petraea. Uneven gene flow
has been repeatedly observed in mixed stands with
limited pollination from Q. robur to Q. petraea (Bacili-
eri et al. 1996; Lagache et al. 2013; Lepais et al. 2013),
with a few exceptions in stands of unbalanced mix-
tures (Gerber et al. 2014). Uneven and unidirectional
gene exchanges between these two species may have
resulted in higher diagnosticity of Q. robur in compari-
son to Q. petraea.

## 4.2 Variation of diagnosticity among populations

There are striking differences of species diagnosticity of
the markers among populations within species (Tables 1,
2, and 3). In populations of Q. petraea, Q. pubescens, and
Q. robur located in the central part of their distributions,
high levels of diagnosticity (mean values of SNP diag-
nosticity of the population) could be observed, while in
populations located at the margins of the distributions,
southern as well as northern, lower diagnosticity was
found. We further showed that populations located at
the edges of distribution are characterized by higher fre-
quencies of near-diagnostic alleles of the other three con-
generic species, suggesting extensive genetic exchanges
(Fig. 4, Figs. 7 and 8 in Appendix). More frequent inter-
specific gene flow at the northern edge of distribution
has been shown earlier in the case Q. petraea and Q.
robur (Beatty et al. 2016; Jensen et al. 2009; Gerber et al.
2014) and has been interpreted as a driver of the succes-
sion dynamics at the northern colonization front of the
two species (Kremer and Hipp 2020; Petit et al. 2003).
In our study, the sessile oak population Killarney (Fig. 4)
and the pedunculate oak population Roudsea (Fig. 8 in
Appendix) are typical examples illustrating interspecific
gene flow between the two species. Similar observations
of more frequent hybridization were made in the case of
Q. petraea and Q. pubescens at the northern edge of dis-
tribution of Q. pubescens (Neophytou et al. 2015; Reuti-
mann et al. 2020), which may have as well contributed to
the expansion of Q. pubescens.

In populations located at the southern edge of distri-
bution (Pomieri, Aspromonte, and Montejo for Q. pet-
raea, Fig. 1 and Fig. 4), the lower diagnosticity may have
resulted from more ancient genetic exchanges with Q.
pubescens and Q. robur, not excluding the potential role
of genetic drift in isolated populations. Indeed, the two
italian populations (Pomieri and Aspromonte) in Sicilia
and Calabria consist today in almost pure isolated stands,
where Q. pubescens is extremely rare, if not absent (Bag-
nato et al. 2012; Modica 2001), while our results indicated
introgression of Q. pubescens into Q. petraea (Fig. 4).
Similarly, the sessile oak population Montejo, in Central
Spain, is introgressed by Q. robur (Fig. 4), where the latter
species is absent today and where contemporary hybridi-
zation has rather been detected with Q. pyrenaica (Val-
buena-Carabana et al. 2005). Finally, a similar scenario
holds for the pedunculate oak population Pedro, which is
located at the extreme southern edge of distribution of Q.
robur (Fig. 1; Table 3, Fig. 8 in Appendix). Hybridization
has been observed with Q. pyrenaica which is today the
most frequent species in the area (Moracho et al. 2016)
and is confirmed by our results revealing the presence of
Q. pyrenaica near-diagnostic alleles in the Q. robur popu-
lation (Fig. 8 in Appendix). However, introgression by Q.
*pubescens* is even more pronounced in our data despite the today's absence of *Q. pubescens* in Extremadura (Fig. 8 in Appendix). To sum up, when comparing our results with previous investigations on interspecific gene flow, recent and/or ancient gene exchanges have faded diagnosticity in the so-called diverging populations, which are located at the northern or southern margins of the distribution.

## 4.3 Variation of diagnosticity among SNPs

Frequency profiles of near-diagnostic alleles differed markedly across SNP in diverging populations. There were cases where lack of diagnosticity affected mainly the same limited number of loci in a given species (Aspromonte and Pomieri in *Q. petraea*, Table 1; Pedro in *Q. robur*, Table 3; and to a smaller extend Switzerland and Ventoux in *Q. pubescens*, Table 2). In the remaining diverging populations (Killarney for *Q. petraea*, Table 1, and Roudsea for *Q. robur*, Table 3), reduced diagnosticity is more evenly distributed across more if not all loci. Contrasting diagnosticity distribution across loci may likely correlate to the timing of hybridization and introgression among the congeneric species. Recent gene exchanges, as first generation hybridization and subsequent backcrosses, will indistinctly impact all loci during the early phase of secondary contact among species and result in reduced diagnosticity of alleles in sympatric species. Such a scenario may hold for the two northern *Q. petraea* (Killarney) and *Q. robur* (Roudsea) populations. Continuous gene exchanges over multiple generations may ultimately result in heterogeneous genomic landscapes, shaped by variable permeability to gene flow along the chromosomes due to the presence of prezygotic or postzygotic barriers and the heterogeneous recombination landscape. This scenario leads ultimately to the maintenance of near-diagnostic loci in genomic regions impermeable to gene flow, while the remaining part of the genome will become poorly differentiated. While this scenario was supported by ABC simulations (Leroy et al. 2020b, 2017), our results further suggest that the genomic distribution of near-diagnostic loci is environment dependent. It is striking that a very limited number of near-diagnostic alleles discovered in western populations of *Q. petraea* show poor diagnosticity in the southern populations Pomieri and Aspromonte (Table 1). Our results further indicated that this low diagnosticity may be due to more interspecific gene flow with *Q. pubescens*, which suggest preferentially introgression in specific genomic regions—whether adaptive or not—resulting ultimately in heterogeneous genomic distribution of near-diagnostic SNPs especially in marginal range parts.

In a recent paper, we showed that introgressed regions between *Q. robur* and *Q. petraea* may be more frequent at higher altitudes (Leroy et al. 2020a), while in another case study in two Asian oak species, the authors found that the genomic landscape of introgression changed in different ecological settings (Fu et al. 2022). A similar picture holds for the diverging southern *Q. robur* population Pedro, where diagnosticity is substantially reduced at a few near-diagnostic SNPs in comparison to other *Q. robur* populations (Table 3), most likely due to introgression by *Q. pubescens* and *Q. pyrenaica* (Fig. 8 in Appendix). Anecdotally, the diverging status of Aspromonte, Pomieri, and Pedro echoes with the taxonomic subspecies status that has been assigned to the Sicilian and Calabrian *Q. petraea* populations (*Q. petraea* ssp. austrothyrenica, Bagnato et al. 2012; Lupini et al. 2019; Merlino et al. 2014) and to the extreme southern Spanish *Q. robur* populations (*Q. robur* ssp. estremadurensis, Vazquez-Pardo et al. 2009).

## 5 Conclusions and outlook

Here, we showed that near-diagnostic marker development for species identification is feasible despite few species barriers, extensive secondary contact, and, consequently, frequent hybridization and introgression. Recently, we demonstrated that the set of near-diagnostic markers resolved species assignment on fossil and archeological oak wood remains, where anatomical features do not allow to discriminate the four deciduous species (Wagner et al. 2024). With the steadily ongoing availability of whole genomes in non-model species including oaks (Lazic et al. 2021), the search of near-diagnostic markers could be extended to the whole Roburoid subsection facilitating white oak species assignment throughout Europe, beyond the subset of four species that we considered here. The near-diagnostic SNPs for the four white oak species could not only be used in forest research and management for reliable and affordable species assignment but also to identify admixed individuals and accurately quantify admixture levels in natural populations (Reutimann et al. 2020). Because the presented alleles are often almost fixed for the target species, these SNPs also allow the identification of hybrid state (F1, F2, backcrosses, later generation hybrids, etc.) with methods like NEWHYBRIDS (Anderson 2008) and altogether help to understand the importance of hydribization and introgression in evolutionary processes. Together with prospect of emergence of field-based genotpying techniques (Urban et al. 2021), such near-diagnostic markers would even allow fast fingerprinting in situ to make decision for forest managers and scientists.
# Appendix

<table><caption>**Table 4** Discovery samples of the whole genome pool-sequenced resources</caption><thead><tr><th>Species</th><th>Sampling site</th><th>Latitude</th><th>Longitude</th><th>Sample size</th></tr></thead><tbody><tr><td>Q. petraea</td><td>Laveyron</td><td>43.9747</td><td>0.2297</td><td>13</td></tr><tr><td rowspan="2">Q. pubescens</td><td>Branne</td><td>44.8399</td><td>-0.2049</td><td>12</td></tr><tr><td>Blaignan</td><td>45.3192</td><td>-0.8559</td><td>6</td></tr><tr><td>Q. robur</td><td>ISS Landes</td><td>44.2263</td><td>1.0112</td><td>18</td></tr><tr><td>Q. pyrenaica</td><td>ISS Landes</td><td>44.2701</td><td>1.0697</td><td>20</td></tr></tbody></table>

<table><caption>**Table 5** Discovery samples of the sequence captured genomic resources</caption><thead><tr><th>Species</th><th>Sampling site</th><th>Latitude</th><th>Longitude</th><th>Sample size</th></tr></thead><tbody><tr><td>Q. petraea</td><td>La Petite Charnie</td><td>48.086</td><td>-0.168</td><td>110</td></tr><tr><td>Q. robur</td><td>La Petite Charnie</td><td>48.086</td><td>-0.168</td><td>135</td></tr></tbody></table>

<table><caption>**Table 6** Geographic origins of training and validation samples</caption><thead><tr><th>Population</th><th>Species</th><th>Country</th><th>Latitude</th><th>Longitude</th><th>Training 1</th><th>Training 2</th><th>Validation</th></tr></thead><tbody><tr><td>Olovo</td><td>Q. petraea</td><td>Bosnia Herzegovina</td><td>44.152</td><td>18.548</td><td>11</td><td></td><td>36</td></tr><tr><td>Artouste</td><td>Q. petraea</td><td>France</td><td>42.890</td><td>-0.400</td><td></td><td></td><td>9</td></tr><tr><td>Berce</td><td>Q. petraea</td><td>France</td><td>47.813</td><td>0.391</td><td></td><td>8</td><td>20</td></tr><tr><td>Bezange</td><td>Q. petraea</td><td>France</td><td>48.759</td><td>6.493</td><td></td><td></td><td>20</td></tr><tr><td>Briouant</td><td>Q. petraea</td><td>France</td><td>43.306</td><td>1.048</td><td>2</td><td></td><td></td></tr><tr><td>Gabas</td><td>Q. petraea</td><td>France</td><td>42.880</td><td>-0.420</td><td></td><td></td><td>18</td></tr><tr><td>Gedre</td><td>Q. petraea</td><td>France</td><td>42.780</td><td>0.020</td><td></td><td></td><td>20</td></tr><tr><td>Gresigne</td><td>Q. petraea</td><td>France</td><td>44.043</td><td>1.749</td><td></td><td></td><td>20</td></tr><tr><td>Josbaig</td><td>Q. petraea</td><td>France</td><td>43.220</td><td>-0.730</td><td></td><td></td><td>19</td></tr><tr><td>Charnie</td><td>Q. petraea</td><td>France</td><td>48.086</td><td>-0.168</td><td>9</td><td></td><td>9</td></tr><tr><td>Laveyron</td><td>Q. petraea</td><td>France</td><td>43.975</td><td>-0.280</td><td>2</td><td>4</td><td>20</td></tr><tr><td>Le Hourque</td><td>Q. petraea</td><td>France</td><td>42.900</td><td>-0.430</td><td></td><td></td><td>19</td></tr><tr><td>Longchamp</td><td>Q. petraea</td><td>France</td><td>47.264</td><td>5.310</td><td></td><td>2</td><td>19</td></tr><tr><td>Papillon</td><td>Q. petraea</td><td>France</td><td>42.920</td><td>-0.030</td><td></td><td></td><td>19</td></tr><tr><td>Péguères</td><td>Q. petraea</td><td>France</td><td>42.870</td><td>-0.120</td><td></td><td>4</td><td>18</td></tr><tr><td>Saint Sauvant</td><td>Q. petraea</td><td>France</td><td>46.380</td><td>0.124</td><td></td><td></td><td>20</td></tr><tr><td>Tronçais</td><td>Q. petraea</td><td>France</td><td>46.680</td><td>2.829</td><td>11</td><td>6</td><td>30</td></tr><tr><td>Vachères</td><td>Q. petraea</td><td>France</td><td>43.983</td><td>5.633</td><td></td><td>2</td><td>20</td></tr><tr><td>Göhrde</td><td>Q. petraea</td><td>Germany</td><td>53.100</td><td>10.846</td><td></td><td>6</td><td>20</td></tr><tr><td>Lappwald</td><td>Q. petraea</td><td>Germany</td><td>52.257</td><td>10.988</td><td></td><td></td><td>20</td></tr><tr><td>Killarney</td><td>Q. petraea</td><td>Ireland</td><td>52.013</td><td>-9.504</td><td></td><td></td><td>20</td></tr><tr><td>Aspromonte</td><td>Q. petraea</td><td>Italy</td><td>38.143</td><td>15.938</td><td></td><td></td><td>48</td></tr><tr><td>Pomieri</td><td>Q. petraea</td><td>Italy</td><td>37.866</td><td>14.069</td><td></td><td></td><td>20</td></tr><tr><td>Montejo</td><td>Q. petraea</td><td>Spain</td><td>41.117</td><td>-3.500</td><td>11</td><td></td><td>38</td></tr><tr><td>Val de Seine</td><td>Q. petraea</td><td>France</td><td>48.398</td><td>3.578</td><td></td><td></td><td>7</td></tr><tr><td>Auros</td><td>Q. pubescens</td><td>France</td><td>44.492</td><td>-0.148</td><td>12</td><td>3</td><td>12</td></tr><tr><td>Briouant</td><td>Q. pubescens</td><td>France</td><td>43.306</td><td>1.048</td><td>2</td><td></td><td></td></tr><tr><td>Blaignan</td><td>Q. pubescens</td><td>France</td><td>45.319</td><td>-0.856</td><td>2</td><td></td><td>5</td></tr><tr><td>Branne</td><td>Q. pubescens</td><td>France</td><td>44.840</td><td>-0.205</td><td></td><td></td><td>11</td></tr></tbody></table>
<table><caption>**Table 7** Genetic and genomic features of near-diagnostic SNPs</caption><thead><tr><th>Population</th><th>Species</th><th>Country</th><th>Latitude</th><th>Longitude</th><th>Training 1</th><th>Training 2</th><th>Validation</th></tr></thead><tbody><tr><td>ISSVentoux</td><td>Q. pubescens</td><td>France</td><td>44.121</td><td>5.312</td><td>16</td><td>8</td><td>40</td></tr><tr><td>Pantano</td><td>Q. pubescens</td><td>Italy</td><td>40.164</td><td>16.671</td><td></td><td></td><td>38</td></tr><tr><td>Ayent</td><td>Q. pubescens</td><td>Switzerland</td><td>46.266</td><td>7.398</td><td>4</td><td>3</td><td>3</td></tr><tr><td>Cordola</td><td>Q. pubescens</td><td>Switzerland</td><td>46.195</td><td>8.863</td><td>4</td><td></td><td></td></tr><tr><td>Remigen</td><td>Q. pubescens</td><td>Switzerland</td><td>47.519</td><td>8.163</td><td>4</td><td>5</td><td>4</td></tr><tr><td>Saillon</td><td>Q. pubescens</td><td>Switzerland</td><td>46.171</td><td>7.167</td><td>4</td><td></td><td>7</td></tr><tr><td>Val de Seine</td><td>Q. pubescens</td><td>France</td><td>48.435</td><td>3.598</td><td></td><td></td><td>3</td></tr><tr><td>Briouant</td><td>Q. pyrenaica</td><td>France</td><td>43.306</td><td>1.048</td><td>2</td><td></td><td></td></tr><tr><td>ISSLandes Mont de Marsan</td><td>Q. pyrenaica</td><td>France</td><td>44.235</td><td>-1.088</td><td>10</td><td></td><td>17</td></tr><tr><td>ISSLandes</td><td>Q. pyrenaica</td><td>France</td><td>44.270</td><td>-1.070</td><td>6</td><td>8</td><td>20</td></tr><tr><td>Hoya Del Nevazo</td><td>Q. pyrenaica</td><td>Spain</td><td>36.957</td><td>-3.423</td><td>7</td><td></td><td>7</td></tr><tr><td>La Calanchera</td><td>Q. pyrenaica</td><td>Spain</td><td>39.572</td><td>-4.647</td><td>6</td><td></td><td>6</td></tr><tr><td>Pedro</td><td>Q. pyrenaica</td><td>Spain</td><td>40.079</td><td>-5.739</td><td>12</td><td>11</td><td>12</td></tr><tr><td>Rascafria</td><td>Q. pyrenaica</td><td>Spain</td><td>40.911</td><td>-3.898</td><td>3</td><td></td><td>3</td></tr><tr><td>Sigmundsherberg</td><td>Q. robur</td><td>Austria</td><td>48.683</td><td>15.750</td><td></td><td>2</td><td>20</td></tr><tr><td>Livno</td><td>Q. robur</td><td>Bosnia Herzegovina</td><td>44.015</td><td>16.630</td><td>11</td><td></td><td>45</td></tr><tr><td>Zivinice</td><td>Q. robur</td><td>Bosnia Herzegovina</td><td>44.446</td><td>18.674</td><td></td><td>5</td><td>19</td></tr><tr><td>Briouant</td><td>Q. robur</td><td>France</td><td>43.306</td><td>1.048</td><td>2</td><td></td><td></td></tr><tr><td>ISSLandes</td><td>Q. robur</td><td>France</td><td>44.226</td><td>-1.011</td><td>2</td><td></td><td>19</td></tr><tr><td>ISSLandes Mont de Marsan</td><td>Q. robur</td><td>France</td><td>44.221</td><td>-1.098</td><td></td><td></td><td>18</td></tr><tr><td>ValSeine</td><td>Q. robur</td><td>France</td><td>48.398</td><td>3.578</td><td></td><td></td><td>8</td></tr><tr><td>Charnie</td><td>Q. robur</td><td>France</td><td>48.086</td><td>-0.168</td><td>9</td><td></td><td>9</td></tr><tr><td>Escherode</td><td>Q. robur</td><td>Germany</td><td>51.333</td><td>9.400</td><td></td><td></td><td>18</td></tr><tr><td>Policoro (Pantano)</td><td>Q. robur</td><td>Italy</td><td>40.159</td><td>16.675</td><td></td><td>3</td><td>18</td></tr><tr><td>Pollutri (San Venanzio)</td><td>Q. robur</td><td>Italy</td><td>42.146</td><td>14.643</td><td></td><td>5</td><td>20</td></tr><tr><td>Arbalan</td><td>Q. robur</td><td>Spain</td><td>42.967</td><td>-2.550</td><td></td><td>6</td><td>15</td></tr><tr><td>Pedro</td><td>Q. robur</td><td>Spain</td><td>40.079</td><td>-5.739</td><td>11</td><td></td><td>11</td></tr><tr><td>Birmensdorf</td><td>Q. robur</td><td>Switzerland</td><td>47.436</td><td>8.255</td><td>3</td><td></td><td>3</td></tr><tr><td>Bonfol</td><td>Q. robur</td><td>Switzerland</td><td>47.463</td><td>7.148</td><td>3</td><td></td><td>4</td></tr><tr><td>Bueren</td><td>Q. robur</td><td>Switzerland</td><td>47.117</td><td>7.383</td><td></td><td></td><td>20</td></tr><tr><td>Cureglia</td><td>Q. robur</td><td>Switzerland</td><td>46.042</td><td>8.950</td><td>2</td><td>3</td><td>5</td></tr><tr><td>Rapperswill</td><td>Q. robur</td><td>Switzerland</td><td>47.239</td><td>8.839</td><td>3</td><td></td><td>4</td></tr><tr><td>Dalkeith</td><td>Q. robur</td><td>UK</td><td>55.917</td><td>-3.033</td><td></td><td></td><td>8</td></tr><tr><td>Roudsea</td><td>Q. robur</td><td>UK</td><td>54.232</td><td>-3.026</td><td></td><td></td><td>12</td></tr><tr><td>Total</td><td></td><td></td><td></td><td></td><td>186</td><td>94</td><td>977</td></tr></tbody></table>

<table><thead><tr><th>SNP ID<sup>f</sup></th><th>Discovery resources</th><th>Sample<sup>a</sup></th><th>Screening<sup>b</sup></th><th>Reference diagnostic species</th><th>Genotype</th><th>Diagnostic nucleotide</th><th>Expression</th><th>Gene ID</th><th>Chr<sup>c</sup></th><th>Position<sup>d</sup></th><th>Distance<sup>e</sup></th></tr></thead><tbody><tr><td>Sc0000254_6223</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. petraea</td><td>AT</td><td>T</td><td>Intergenic</td><td>NA</td><td>1</td><td>12,964,247</td><td>3,103,379</td></tr><tr><td>Sc0000121_355205</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. petraea</td><td>CT</td><td>C</td><td>Intergenic</td><td>NA</td><td>5</td><td>69,028,945</td><td>33,570,950</td></tr><tr><td>Sc0000040_1694351</td><td>Seq-cap</td><td>T</td><td>PD</td><td>Q. petraea</td><td>AT</td><td>A</td><td>Intergenic</td><td>NA</td><td>2</td><td>17,875,174</td><td>8,132,192</td></tr><tr><td>Sc0000043_1651618</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. petraea</td><td>AG</td><td>G</td><td>Intergenic</td><td>NA</td><td>2</td><td>41,927,441</td><td>2,874,167</td></tr><tr><td>Sc0000055_2262067</td><td>Pool-seq</td><td>T</td><td>VA</td><td>Q. petraea</td><td>CT</td><td>T</td><td>Intergenic</td><td>NA</td><td>5</td><td>27,820,895</td><td>3,322,468</td></tr><tr><td>Sc0000083_147504</td><td>Seq-cap</td><td>T+V</td><td></td><td>Q. petraea</td><td>AG</td><td>A</td><td>Exonic</td><td>Qrob_ G0609970.2</td><td>2</td><td>46,841,435</td><td>2,526,743</td></tr></tbody></table>
<table><thead><tr><th>SNP ID<sup>f</sup></th><th>Discovery resources</th><th>Sample<sup>a</sup></th><th>Screening<sup>b</sup></th><th>Reference diagnostic species</th><th>Genotype</th><th>Diagnostic nucleotide</th><th>Expression</th><th>Gene ID</th><th>Chr<sup>c</sup></th><th>Position<sup>d</sup></th><th>Distance<sup>e</sup></th></tr></thead><tbody><tr><td>Sc0000090_1332487</td><td>Pool-seq</td><td>T</td><td>PQ</td><td>Q. petraea</td><td>TC</td><td>T</td><td>Exonic</td><td>Qrob_ G0088630.2</td><td>7</td><td>44,619,559</td><td>28,315</td></tr><tr><td>Sc0000090_1360802</td><td>Seq-cap</td><td>T</td><td>DisC</td><td>Q. petraea</td><td>CG</td><td>C</td><td>Exonic</td><td>Qrob_ G0088640.2</td><td>7</td><td>44,647,874</td><td>28,315</td></tr><tr><td>Sc0000118_1466708</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. petraea</td><td>AG</td><td>G</td><td>Exonic</td><td>Qrob_ G0081080.2</td><td>11</td><td>4,483,502</td><td>2,258,901</td></tr><tr><td>Sc0000135_261350</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. petraea</td><td>AC</td><td>C</td><td>Exonic</td><td>Qrob_ G0222840.2</td><td>7</td><td>24,799,454</td><td>2,918,740</td></tr><tr><td>Sc0000145_700044</td><td>Seq-cap</td><td>T+V</td><td></td><td>Q. petraea</td><td>CT</td><td>T</td><td>Intergenic</td><td>NA</td><td>2</td><td>31,806,962</td><td>1,145,225</td></tr><tr><td>Sc0000203_707735</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. petraea</td><td>AT</td><td>A</td><td>Intronic</td><td>Qrob_ G0237050.2</td><td>1</td><td>51,985,124</td><td>17,848</td></tr><tr><td>Sc0000274_909817</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. petraea</td><td>AG</td><td>A</td><td>Exonic</td><td>Qrob_ G0320010.2</td><td>2</td><td>35,418,361</td><td>563,863</td></tr><tr><td>Sc0000464_236576</td><td>Pool-seq</td><td>T</td><td>PD</td><td>Q. petraea</td><td>CT</td><td>T</td><td>Exonic</td><td>Qrob_ G0523500.2</td><td>5</td><td>23,753,953</td><td>744,474</td></tr><tr><td>Sc0000481_343721</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. petraea</td><td>AG</td><td>G</td><td>Intronic</td><td>Qrob_ G0512850.2</td><td>5</td><td>3,222,156</td><td>22,554</td></tr><tr><td>Sc0000974_98303</td><td>Pool-seq</td><td>T</td><td>PD</td><td>Q. petraea</td><td>AG</td><td>A</td><td>Exonic</td><td>Qrob_ G0759540.2</td><td>3</td><td>47,813,240</td><td>NA</td></tr><tr><td>Sc0000485_93093</td><td>Pool-seq</td><td>T</td><td>VA</td><td>Q. pubescens</td><td>AG</td><td>A</td><td>Exonic</td><td>Qrob_ G0539160.2</td><td>1</td><td>16,067,626</td><td>3,103,379</td></tr><tr><td>Sc0000314_149731</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pubescens</td><td>AT</td><td>A</td><td>Intergenic</td><td>NA</td><td>2</td><td>55,831,586</td><td>5,480,810</td></tr><tr><td>Sc0000271_690874</td><td>Pool-seq</td><td>T</td><td>VA</td><td>Q. pubescens</td><td>CT</td><td>T</td><td>Intergenic</td><td>NA</td><td>7</td><td>39,945,935</td><td>897,324</td></tr><tr><td>Sc0000047_2398879</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pubescens</td><td>AG</td><td>A</td><td>Intronic</td><td>Qrob_ G0585850.2</td><td>9</td><td>16,825,011</td><td>NA</td></tr><tr><td>Sc0000062_118505</td><td>Pool-seq</td><td>T</td><td>DisC</td><td>Q. pubescens</td><td>CT</td><td>T</td><td>Intronic</td><td>Qrob_ G0091810.2</td><td>6</td><td>15,935,077</td><td>4,083,535</td></tr><tr><td>Sc0000088_1796044</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pubescens</td><td>AC</td><td>A</td><td>Intronic</td><td>Qrob_ G0328570.2</td><td>12</td><td>20,880,546</td><td>NA</td></tr><tr><td>Sc0000109_800763</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pubescens</td><td>CT</td><td>C</td><td>Exonic</td><td>Qrob_ G0124910.2</td><td>1</td><td>43,547,799</td><td>4,560,446</td></tr><tr><td>Sc0000111_693153</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pubescens</td><td>CT</td><td>T</td><td>Intergenic</td><td>NA</td><td>7</td><td>28,004,200</td><td>286,006</td></tr><tr><td>Sc0000170_630013</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pubescens</td><td>AC</td><td>A</td><td>Intronic</td><td>Qrob_ G0200770.2</td><td>6</td><td>35,588,545</td><td>14,530,676</td></tr><tr><td>Sc0000192_329301</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pubescens</td><td>AG</td><td>A</td><td>Intronic</td><td>Qrob_ G0228240.2</td><td>2</td><td>39,053,274</td><td>2,874,167</td></tr><tr><td>Sc0000482_334917</td><td>Pool-seq</td><td>T</td><td>PD</td><td>Q. pubescens</td><td>CT</td><td>C</td><td>Intronic</td><td>Qrob_ G0533790.2</td><td>NA</td><td>NA</td><td>NA</td></tr><tr><td>Sc0000403_286465</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pyrenaica</td><td>CT</td><td>T</td><td>Intergenic</td><td>NA</td><td>2</td><td>105,577,540</td><td>28,340,724</td></tr><tr><td>Sc0000006_2873224</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pyrenaica</td><td>AG</td><td>A</td><td>Intronic</td><td>Qrob_ G0005870.2</td><td>6</td><td>21,057,869</td><td>5,122,792</td></tr><tr><td>Sc0000014_2037045</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pyrenaica</td><td>AC</td><td>C</td><td>Exonic</td><td>Qrob_ G0064170.2</td><td>2</td><td>77,236,816</td><td>6,145,559</td></tr><tr><td>Sc0000053_1344456</td><td>Pool-seq</td><td>T</td><td>DisC</td><td>Q. pyrenaica</td><td>AG</td><td>A</td><td>Intergenic</td><td>NA</td><td>6</td><td>11,300,343</td><td>294,652</td></tr><tr><td>Sc0000085_73024</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pyrenaica</td><td>TG</td><td>G</td><td>Intronic</td><td>Qrob_ G0563240.2</td><td>10</td><td>14,291,154</td><td>1,973,542</td></tr><tr><td>Sc0000228_1091905</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pyrenaica</td><td>TC</td><td>C</td><td>Intergenic</td><td>NA</td><td>5</td><td>35,457,995</td><td>7,637,100</td></tr><tr><td>Sc0000269_924931</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pyrenaica</td><td>AC</td><td>A</td><td>Intronic</td><td>Qrob_ G0632320.2</td><td>2</td><td>26,007,366</td><td>561,547</td></tr><tr><td>Sc0000287_474090</td><td>Pool-seq</td><td>T</td><td>AF</td><td>Q. pyrenaica</td><td>AC</td><td>C</td><td>Exonic</td><td>Qrob_ G0459590.2</td><td>2</td><td>71,091,257</td><td>6,145,559</td></tr><tr><td>Sc0000307_852597</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pyrenaica</td><td>AG</td><td>A</td><td>Intergenic</td><td>NA</td><td>7</td><td>16,871,124</td><td>7,928,330</td></tr><tr><td>Sc0000517_383812</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pyrenaica</td><td>CG</td><td>C</td><td>Intergenic</td><td>NA</td><td>10</td><td>18,200,824</td><td>125,219</td></tr><tr><td>Sc0000695_157206</td><td>Pool-seq</td><td>T</td><td>AF</td><td>Q. pyrenaica</td><td>AT</td><td>A</td><td>Exonic</td><td>Qrob_ G0671270.2</td><td>8</td><td>55,531,952</td><td>68,141</td></tr><tr><td>Sc0000778_61930</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. pyrenaica</td><td>CT</td><td>T</td><td>Intronic</td><td>Qrob_ G0070130.2</td><td>10</td><td>12,317,612</td><td>1,973,542</td></tr></tbody></table>
<table><thead><tr><th>SNP ID<sup>f</sup></th><th>Discovery resources</th><th>Sample<sup>a</sup></th><th>Screening<sup>b</sup></th><th>Reference diagnostic species</th><th>Genotype</th><th>Diagnostic nucleotide</th><th>Expression</th><th>Gene ID</th><th>Chr<sup>c</sup></th><th>Position<sup>d</sup></th><th>Distance<sup>e</sup></th></tr></thead><tbody><tr><td>Sc0000013_2578823</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>AG</td><td>A</td><td>Intronic</td><td>Qrob_ G0010260.2</td><td>2</td><td>61,312,396</td><td>5,480,810</td></tr><tr><td>Sc0000038_794573</td><td>Pool-seq</td><td>T</td><td>PQ</td><td>Q. robur</td><td>TG</td><td>G</td><td>Intronic</td><td>Qrob_ G0701760.2</td><td>1</td><td>38,987,353</td><td>4,560,446</td></tr><tr><td>Sc0000053_1639108</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>AG</td><td>A</td><td>Intergenic</td><td>NA</td><td>6</td><td>11,594,995</td><td>256,547</td></tr><tr><td>Sc0000053_1895655</td><td>Seq-cap</td><td>T</td><td>DisC</td><td>Q. robur</td><td>CT</td><td>T</td><td>Exonic</td><td>Qrob_ G0631440.2</td><td>6</td><td>11,851,542</td><td>256,547</td></tr><tr><td>Sc0000099_1839376</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>AC</td><td>A</td><td>Intronic</td><td>Qrob_ G0084290.2</td><td>11</td><td>6,742,403</td><td>2,258,901</td></tr><tr><td>Sc0000111_979159</td><td>Pool-seq</td><td>T</td><td>AF</td><td>Q. robur</td><td>CT</td><td>C</td><td>Intronic</td><td>Qrob_ G0135420.2</td><td>7</td><td>27,718,194</td><td>286,006</td></tr><tr><td>Sc0000158_462639</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>AT</td><td>T</td><td>Exonic</td><td>Qrob_ G0304430.2</td><td>2</td><td>26,568,913</td><td>230,471</td></tr><tr><td>Sc0000158_693110</td><td>Seq-cap</td><td>T</td><td>AF</td><td>Q. robur</td><td>CG</td><td>G</td><td>Exonic</td><td>Qrob_ G0304580.2</td><td>2</td><td>26,799,384</td><td>230,471</td></tr><tr><td>Sc0000203_689887</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>TG</td><td>T</td><td>Intronic</td><td>Qrob_ G0237030.2</td><td>1</td><td>51,967,276</td><td>17,848</td></tr><tr><td>Sc0000225_507799</td><td>Seq-cap</td><td>T</td><td>AF</td><td>Q. robur</td><td>CG</td><td>G</td><td>Exonic</td><td>Qrob_ G0487320.2</td><td>5</td><td>24,498,427</td><td>744,474</td></tr><tr><td>Sc0000240_289656</td><td>Seq-cap</td><td>T</td><td>VA</td><td>Q. robur</td><td>AG</td><td>G</td><td>Exonic</td><td>Qrob_ G0318610.2</td><td>5</td><td>4,277,129</td><td>1,032,419</td></tr><tr><td>Sc0000270_806328</td><td>Pool-seq</td><td>T</td><td>AF</td><td>Q. robur</td><td>AG</td><td>G</td><td>Exonic</td><td>Qrob_ G0692980.2</td><td>7</td><td>33,545,097</td><td>5,540,897</td></tr><tr><td>Sc0000339_4638</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>CT</td><td>C</td><td>Intronic</td><td>Qrob_ G0473660.2</td><td>1</td><td>3,718,979</td><td>9,245,268</td></tr><tr><td>Sc0000381_206331</td><td>Seq-cap</td><td>T+V</td><td></td><td>Q. robur</td><td>AC</td><td>A</td><td>Intergenic</td><td>NA</td><td>7</td><td>39,048,611</td><td>897,324</td></tr><tr><td>Sc0000447_521057</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>AG</td><td>G</td><td>Exonic</td><td>Qrob_ G0543330.2</td><td>10</td><td>22,549,567</td><td>4,348,743</td></tr><tr><td>Sc0000481_366275</td><td>Seq-cap</td><td>T</td><td>PD</td><td>Q. robur</td><td>CT</td><td>T</td><td>Exonic</td><td>Qrob_ G0512860.2</td><td>5</td><td>3,244,710</td><td>22,554</td></tr><tr><td>Sc0000517_258593</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>CT</td><td>T</td><td>Intergenic</td><td>NA</td><td>10</td><td>18,075,605</td><td>125,219</td></tr><tr><td>Sc0000546_456229</td><td>Seq-cap</td><td>T</td><td>PD</td><td>Q. robur</td><td>AG</td><td>G</td><td>Exonic</td><td>Qrob_ G0761790.2</td><td>2</td><td>35,982,224</td><td>563,863</td></tr><tr><td>Sc0000598_295142</td><td>Seq-cap</td><td>T</td><td>PD</td><td>Q. robur</td><td>CT</td><td>C</td><td>Exonic</td><td>Qrob_ G0575620.2</td><td>2</td><td>32,952,187</td><td>1,145,225</td></tr><tr><td>Sc0000695_225347</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>AT</td><td>T</td><td>Intronic</td><td>Qrob_ G0671240.2</td><td>8</td><td>55,600,093</td><td>68,141</td></tr><tr><td>Sc0000796_82698</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>AT</td><td>T</td><td>Intronic</td><td>Qrob_ G0759570.2</td><td>NA</td><td>NA</td><td>NA</td></tr><tr><td>Sc0000967_33996</td><td>Pool-seq</td><td>T+V</td><td></td><td>Q. robur</td><td>AT</td><td>T</td><td>Exonic</td><td>Qrob_ G0709860.2</td><td>2</td><td>49,368,178</td><td>2,526,743</td></tr></tbody></table>

<footnote>a Study samples (T training populations, V validation population)</footnote>

<footnote>b Screening criteria from training to validation (PD poor diagnosticity, PQ poor quality of cluster delimitation, Disc genotype discrepancy between different multiplexes, VA variable success (numerous missing data), AF amplification failure after primer redesign)</footnote>

<footnote>c Chr chromosome bearing the near-diagnostic SNP</footnote>

<footnote>d Position (in bp) on the chromosome</footnote>

<footnote>e Distance (in bp) with previous near-diagnostic SNP on the same chromosome</footnote>

<footnote>f SNP identification comprises scaffold number (Sc#) and position on the scaffold (Plomion et al. 2018)</footnote>
<table><caption>**Table 8** Overall frequencies of near-diagnostic alleles in the validation populations</caption><thead><tr><th>SNP ID</th><th>Reference diagnostic species</th><th>Near-diagnostic allele</th><th>Q. petraea</th><th>Q. pubescens</th><th>Q. pyrenaica</th><th>Q. robur</th><th>Diagnosticity</th></tr></thead><tbody><tr><td>Sc0000254_6223</td><td>Q. petraea</td><td>T</td><td><strong>0.863</strong></td><td>0.116</td><td>0.074</td><td>0.023</td><td>0.792</td></tr><tr><td>Sc0000121_355205</td><td>Q. petraea</td><td>C</td><td><strong>0.966</strong></td><td>0.151</td><td>0.109</td><td>0.068</td><td>0.857</td></tr><tr><td>Sc0000043_1651618</td><td>Q. petraea</td><td>G</td><td><strong>0.738</strong></td><td>0.156</td><td>0.023</td><td>0.027</td><td>0.669</td></tr><tr><td>Sc0000083_147504</td><td>Q. petraea</td><td>A</td><td><strong>0.876</strong></td><td>0.121</td><td>0.320</td><td>0.046</td><td>0.714</td></tr><tr><td>Sc0000118_1466708</td><td>Q. petraea</td><td>G</td><td><strong>0.873</strong></td><td>0.117</td><td>0.046</td><td>0.077</td><td>0.793</td></tr><tr><td>Sc0000135_261350</td><td>Q. petraea</td><td>C</td><td><strong>0.823</strong></td><td>0.132</td><td>0.056</td><td>0.082</td><td>0.733</td></tr><tr><td>Sc0000145_700044</td><td>Q. petraea</td><td>T</td><td><strong>0.949</strong></td><td>0.074</td><td>0.015</td><td>0.035</td><td>0.908</td></tr><tr><td>Sc0000203_707735</td><td>Q. petraea</td><td>A</td><td><strong>0.903</strong></td><td>0.060</td><td>0.042</td><td>0.029</td><td>0.859</td></tr><tr><td>Sc0000274_909817</td><td>Q. petraea</td><td>A</td><td><strong>0.786</strong></td><td>0.086</td><td>0.123</td><td>0.018</td><td>0.710</td></tr><tr><td>Sc0000481_343721</td><td>Q. petraea</td><td>G</td><td><strong>0.926</strong></td><td>0.242</td><td>0.085</td><td>0.049</td><td>0.801</td></tr><tr><td>Sc0000314_149731</td><td>Q. pubescens</td><td>A</td><td>0.116</td><td><strong>0.660</strong></td><td>0.050</td><td>0.054</td><td>0.587</td></tr><tr><td>Sc0000047_2398879</td><td>Q. pubescens</td><td>A</td><td>0.043</td><td><strong>0.896</strong></td><td>0.045</td><td>0.019</td><td>0.860</td></tr><tr><td>Sc0000088_1796044</td><td>Q. pubescens</td><td>A</td><td>0.027</td><td><strong>0.906</strong></td><td>0.000</td><td>0.007</td><td>0.895</td></tr><tr><td>Sc0000109_800763</td><td>Q. pubescens</td><td>C</td><td>0.121</td><td><strong>0.752</strong></td><td>0.015</td><td>0.007</td><td>0.704</td></tr><tr><td>Sc0000111_693153</td><td>Q. pubescens</td><td>T</td><td>0.337</td><td><strong>0.996</strong></td><td>0.061</td><td>0.012</td><td>0.859</td></tr><tr><td>Sc0000170_630013</td><td>Q. pubescens</td><td>A</td><td>0.015</td><td><strong>0.657</strong></td><td>0.000</td><td>0.005</td><td>0.650</td></tr><tr><td>Sc0000192_329301</td><td>Q. pubescens</td><td>A</td><td>0.010</td><td><strong>0.697</strong></td><td>0.054</td><td>0.006</td><td>0.674</td></tr><tr><td>Sc0000403_286465</td><td>Q. pyrenaica</td><td>T</td><td>0.008</td><td>0.014</td><td><strong>0.900</strong></td><td>0.011</td><td>0.889</td></tr><tr><td>Sc0000006_2873224</td><td>Q. pyrenaica</td><td>A</td><td>0.012</td><td>0.016</td><td><strong>0.962</strong></td><td>0.009</td><td>0.950</td></tr><tr><td>Sc0000014_2037045</td><td>Q. pyrenaica</td><td>C</td><td>0.006</td><td>0.000</td><td><strong>0.908</strong></td><td>0.004</td><td>0.905</td></tr><tr><td>Sc0000085_73024</td><td>Q. pyrenaica</td><td>G</td><td>0.004</td><td>0.024</td><td><strong>0.946</strong></td><td>0.000</td><td>0.937</td></tr><tr><td>Sc0000228_1091905</td><td>Q. pyrenaica</td><td>C</td><td>0.003</td><td>0.004</td><td><strong>0.844</strong></td><td>0.002</td><td>0.841</td></tr><tr><td>Sc0000269_924931</td><td>Q. pyrenaica</td><td>A</td><td>0.003</td><td>0.000</td><td><strong>0.900</strong></td><td>0.004</td><td>0.898</td></tr><tr><td>Sc0000307_852597</td><td>Q. pyrenaica</td><td>A</td><td>0.013</td><td>0.000</td><td><strong>0.546</strong></td><td>0.022</td><td>0.534</td></tr><tr><td>Sc0000517_383812</td><td>Q. pyrenaica</td><td>C</td><td>0.001</td><td>0.022</td><td><strong>0.975</strong></td><td>0.004</td><td>0.966</td></tr><tr><td>Sc0000778_61930</td><td>Q. pyrenaica</td><td>T</td><td>0.004</td><td>0.000</td><td><strong>0.969</strong></td><td>0.000</td><td>0.968</td></tr><tr><td>Sc0000013_2578823</td><td>Q. robur</td><td>A</td><td>0.078</td><td>0.028</td><td>0.047</td><td><strong>0.814</strong></td><td>0.763</td></tr><tr><td>Sc0000053_1639108</td><td>Q. robur</td><td>A</td><td>0.094</td><td>0.015</td><td>0.071</td><td><strong>0.925</strong></td><td>0.865</td></tr><tr><td>Sc0000099_1839376</td><td>Q. robur</td><td>A</td><td>0.017</td><td>0.027</td><td>0.000</td><td><strong>0.881</strong></td><td>0.866</td></tr><tr><td>Sc0000158_462639</td><td>Q. robur</td><td>T</td><td>0.043</td><td>0.000</td><td>0.054</td><td><strong>0.799</strong></td><td>0.767</td></tr><tr><td>Sc0000203_689887</td><td>Q. robur</td><td>T</td><td>0.072</td><td>0.029</td><td>0.078</td><td><strong>0.796</strong></td><td>0.736</td></tr><tr><td>Sc0000339_4638</td><td>Q. robur</td><td>C</td><td>0.076</td><td>0.000</td><td>0.062</td><td><strong>0.922</strong></td><td>0.876</td></tr><tr><td>Sc0000381_206331</td><td>Q. robur</td><td>A</td><td>0.054</td><td>0.000</td><td>0.108</td><td><strong>0.941</strong></td><td>0.887</td></tr><tr><td>Sc0000447_521057</td><td>Q. robur</td><td>G</td><td>0.213</td><td>0.028</td><td>0.000</td><td><strong>0.956</strong></td><td>0.876</td></tr><tr><td>Sc0000517_258593</td><td>Q. robur</td><td>T</td><td>0.104</td><td>0.091</td><td>0.069</td><td><strong>0.960</strong></td><td>0.872</td></tr><tr><td>Sc0000695_225347</td><td>Q. robur</td><td>T</td><td>0.140</td><td>0.008</td><td>0.008</td><td><strong>0.921</strong></td><td>0.869</td></tr><tr><td>Sc0000796_82698</td><td>Q. robur</td><td>T</td><td>0.083</td><td>0.074</td><td>0.133</td><td><strong>0.901</strong></td><td>0.804</td></tr><tr><td>Sc0000967_33996</td><td>Q. robur</td><td>T</td><td>0.077</td><td>0.004</td><td>0.085</td><td><strong>0.965</strong></td><td>0.910</td></tr></tbody></table>

Frequencies in bold characters correspond to near-diagnostic alleles of the reference diagnostic species
<figure><img src="image_10.png" /><figcaption>Allele frequencies</figcaption></figure>

**Fig. 6** Heat map of frequencies of near-diagnostic alleles and diagnosticity in the training populations of the four species. SNPs are clustered for their diagnostic value for each species (reference species). First to fourth columns correspond respectively to near-diagnostic alleles of *Q. petraea*, *Q. pubescens*, *Q. pyrenaica*, and *Q. robur*
<figure><img src="image_11.png" /><figcaption>**Fig. 7** Frequencies of near-diagnostic alleles of *Q. petraea*, *Q. pyrenaica*, and *Q. robur* in *Q. pubescens* populations. Populations in red and green exhibit respectively diverging diagnosticity (DD populations) and expected diagnosticity (ED populations). Shown are the mean frequencies of all near-diagnostic alleles of a given species (*Q. petraea*, *Q. pyrenaica*, *Q. robur*) in *Q. pubescens* populations. Geographic locations of *Q. pubescens* populations are shown in Fig. 1b</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>**Fig. 8** Frequencies of near-diagnostic alleles of *Q. petraea*, *Q. pubescens*, and *Q. pyrenaica* in *Q. robur* populations. Populations in red and green exhibit respectively diverging diagnosticity (DD populations) and expected diagnosticity (ED populations). Shown are the mean frequencies of all near-diagnostic alleles of a given (*Q. petraea*, *Q. pubescens*, *Q. pyrenaica*) species in *Q. robur* populations. Geographic locations of *Q. robur* populations are shown in Fig. 1d</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>**Fig. 9** Biplot of principal components of tree samples based on a principal component analysis (PCA) conducted in the *Q. pubescens* and *Q. petraea* validation populations. Red dots, *Q. petraea* samples; orange dots, *Q. pubescens* samples. Blue dots, *Q. petraea* Pomieri population. Green dots, *Q. pubescens* Switzerland population</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>**Fig. 10** Biplot of principal components of tree samples based on a principal component analysis (PCA) conducted in the *Q. robur* and *Q. petraea* validation populations. Red dots, *Q. petraea* samples; blue dots, *Q. robur* samples. Green dots, *Q. petraea* Killarney population. Orange dots, *Q. robur* Roudsea population</figcaption></figure>

**Acknowledgements**

We thank colleagues that contributed to the collection of material made for this study: Dalibor Ballian (Bosnia and Herzegovina), María Valbuena Carabaña and Luis Gil (Spain), and Giovanni Giuseppe Vendramin (Italy). We extend our appreciation to partners of the former EU-supported FAIROAK and OAKFLOW projects, and of the EVOLTREE Network of Excellence, who collected material included in this study. The MassArray genotyping was performed at the PGTB (doi:10.15454/1.5572396583599417E12) with the help of Laure Dubois, Céline Lalanne, and Marie Massot. We are grateful to François Ehrenmann for his contributions to the figures of the manuscript.

**Authors' contributions**

Conception of the study, TL and AK; sampling and collection of material, AK, CR, and TL; discovery of near-diagnostic markers in pool-sequenced resources, TL; discovery of near-diagnostic markers in sequence captured resources, IL; design of multiplexes and genotyping of natural populations, EG and AD; data analysis, AK, TL, and SW; and writing of the manuscript, AK, EG, and TL. All authors reviewed and approved the final manuscript.

**Funding**

This research was supported by the European Research Council through an Advanced Grant (project TREEPEACE no. FP7-339728), by an ANR (Agence Nationale de la Recherche) Grant (project GENOAK 2022, #BSV6-009-02), and by the EVOLTREE Opportunity call (project OakID2).

**Availability of data and materials**

The data are accessible at the INRAE data repository via a DOI URL address (https://doi.org/10.57745/0JYLZU).

**Declarations**

**Ethics approval and consent to participate**

Not applicable.

**Consent for publication**

All authors gave their informed consent to this publication and its content.

**Competing interests**

The authors declare that they have no competing interests.

**Author details**

¹UMR BIOGECO, INRAE, Université de Bordeaux, 69 Route d'Arcachon, Cestas 33612, France. ²UMR BIOGECO, INRAE, PGTB, Université de Bordeaux, 69 Route d'Arcachon, Cestas 33612, France. ³Helix Venture, Merignac 33700, France. ⁴UMR CAGT, CNRS, Université Paul Sabatier, 37 Allées Jules Guesde, Toulouse 31000, France. ⁵CNRGV Plant Genomic Center, INRAE., 24 Chemin Borde Rouge, Castanet-Tolosan 31320, France. ⁶Swiss Federal Research Institute WSL, Birmensdorf 8903, Switzerland. ⁷GenPhySE, INRAE, ENVT, Université de Toulouse, Castanet-Tolosan, France.
Received: 28 November 2023 Accepted: 3 May 2024

Published online: 21 May 2024

## References

Alberto F, Niort J, Deroy J, Lepais O, Vitalis R, Galop D et al (2010) Population differentiation of sessile oak at the altitudinal front of migration in the French Pyrenees. Mol Ecol 19(13):2626-2639. https://doi.org/10.1111/j.1365-294x.2010.04631.x

Anderson EC (2008) Bayesian inference of species hybrids using multilocus dominant genetic markers. Phil Trans r Soc B 363:2841-2850. https://doi.org/10.1098/rstb.2008.0043

Bacilieri R, Ducoussou A, Petit RJ, Kremer A (1996) Mating system and asymmetric hybridization in a mixed stand of European oaks. Evolution 50(2):900-908. https://doi.org/10.1111/j.1558-5646.1996.tb03898.x

Bacilieri R, Ducoussou A, Kremer A (1995). Genetic, morphological, ecological, and phenological differentiation between *Quercus petraea* (Matt) Liebl and *Quercus robur* L in a mixed stand of northwest of France. *Silvae Genetica* 44(1): 1-10. https://www.thuenen.de/media/institute/fg/PDF/Silvae_Genetica/1995/Vol_44_Heft_1/44_1_1.pdf

Bagnato S, Merlino A, Mercurio R, Solano F, Scarfo F, Spampinato G (2012). Le basi conosctive per il restauro forestale: il caso di Bosco Pomieri (Parco Regionale delle Madonie, Sicilia). *Forest* 9: 8-19. https://doi.org/10.3832/efor0679-008

Beatty GE, Montgomery WL, Spans F, Tosh DG, Provan J (2016) Pure species in a continuum of genetic and morphological variation: sympatric oaks at the edge of their range. Ann Bot 117:541-549. https://doi.org/10.1093/aob/mcw002

Bolger AM, Lohse M, Usadel B (2014) Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics 30:2114-2120. https://doi.org/10.1093/bioinformatics/btu170

Bruschi P, Vendramin GG, Bussotti F, Grossoni P (2000). Morphological and molecular differentiation between *Quercus petraea* (Matt.) Liebl. and *Quercus pubescens* Willd. (Fagaceae) in Northern and Central Italy. *Annals of Botany* 85(3): 325-333. https://doi.org/10.1006/anbo.1999.1046

Camus A (1938). *Les chênes*. Monographie du genre Quercus. Tome II. Genre Quercus. Sous-genre Euquercus. Paul Lechevalier: Paris 830

Caudullo G, Welk E, San-Miguel-Ayanz J (2017) Chorological mapsfor the main European woody species. Data Brief 12:662-666. https://doi.org/10.1016/j.dib.2017.05.007

Coq-Etchegaray D, Bernillon S, Le-Provost G, Kremer A, Ducoussou A, Lalanne C, Bonne F, Moing A, Plomion C, Brachi B (2023) Extensive variation of leaf specialized metabolite production in sessile oak (*Quercus petraea*) populations is to a large extent genetically determined but not locally adaptive. Preprint at bioRxiv. https://doi.org/10.1101/2023.04.07.536008

Degen B, Blanc-Jolivet C, Bakhtina S, Ianbaev R, Yanbaev Y, Mader M, Nurnberg S, Schröder H (2021) Applying targeted genotyping by sequencing with a new set of nuclear and plastid SNP and indel loci for *Quercus robur* and *Quercus petraea*. Conserv Genet Resour 13:345-347. https://doi.org/10.1007/s12686-021-01207-6

Denk T, Grimm GW, Manos PS, Deng M, Hipp AL (2017) An updated infrage-
neric classification of the oaks: review of previous taxonomic schemes
and synthesis of evolutionary patterns. In: Gil-Pelegrin E, Peguero-Pina JJ,
Sancho-Knapik D (eds) Oaks physiological Ecology. Springer, Exploring
the functional diversity of the genus Quercus L, pp 13-38

Dupouey JL, Badeau V (1993) Morphological variability of oaks (*Quercus robur*
*L*, *Quercus petraea* (Matt)Liebl, *Quercus pubescens* Willd) in northeastern
France. Preliminary Results Ann for Sci 50:35s-40s. https://doi.org/10.
1051/forest:19930702

Ellis JA, Ong B (2017). The MassARRAY® System for targeted SNP genotyping. In: White S, Cantsilieris S (eds) Genotyping. Methods in Molecular Biology, vol 1492. Humana Press, New York, NY. https://doi.org/10.1007/978-1-4939-6442-0_5

Fu R, Zhu Y, Liu Y, Feng Y, Lu R, Li Y et al (2022) Genome-wide analyses of introgression bewteen two sympatric Asian oak species. Nature Ecology & Evolution 6:924-935. https://doi.org/10.1038/s41559-022-01754-7

Gerber S, Chadoeuf J, Gugerli F, Lascoux M, Buiteveld J, Cottrell J *et al* (2014).
High rates of gene flow by pollen and seed in oak populations across
Europe. Plos One 9(1). https://doi.org/10.1371/journal.pone.0085130

Grandjean G, Sigaud P (1987) Contribution à la taxonomie et à l'écologie des chênes du Berry. Ann for Sci 44:35-66. https://doi.org/10.1051/forest:19870104

Gregorius HR (1984) A Unique Genetic Distance Biometric J 26:1-14

Gregorius HR, Roberts JH (1986) Mesuring genetic differentiation in sub-
populations. Theor Appl Genet 71:826-834. https://doi.org/10.1007/
BF00276425

Guichoux E, Lagache L, Wagner S, Leger P, Petit RJ (2011). Two highly validated multiplexes (12-plex and 8-plex) for species delimitation and parent- age analysis in oaks (Quercus spp.). *Molecular Ecology Resources* **11**(3): 578-585. https://doi.org/10.1111/j.1755-0998.2011.02983.x

Hipp AL, Manos PS, Hahn M, Avishai M, Bodénés C, Cavender-Bares J et al (2020) The genomic landscape of the global oak phylogeny. New Phytol 226:1198-1212. https://doi.org/10.1111/nph.16162

Jensen J, Larsen A, Nielsen LR, Cottrell J (2009) Hybridization between *Quercus robur* and *Q. petraea* in a mixed oak stand in Denmark. *Annals of Forest Science* 66: 706. https://doi.org/10.1051/forest/2009058

Jombart T (2008) Adegenet: a R package for the multivariate analysis of genetic markers. Bioinformatics 24:1403-1405. https://doi.org/10.1093/bioinformatics/btn129

Jurksiene G, Baliuckas V (2014). Leaf morphological variation of sessile oak (*Quercus petraea* (Matt.)Liebl.) and pedunculate oak (*Quercus robur* L.) in Lithuania.. In: Treija S and Skujeniece S (eds) *Research for Rural Development 2014*, Vol 2, pp 63–69. https://hdl.handle.net/20.500.12259/87903

Kelleher CT, Hodkinson TR, Douglas GC, Kelly DL (2005). Species distinction in Irish populations of *Quercus petraea* and *Q. robur*: morphological versus molecular analyses. *Annals of Botany* 96(7): 1237–1246. https://doi.org/10.1093/aob/mci275

Kofler R, Pandey RV, Schlotterer C (2011) POPOOLATION2: identifying differen-
tiation between populations using sequencing of pooled DNA samples
(Pool-Seq). Bioinformatics 27:3435-3436. https://doi.org/10.1093/bioin
formatics/btr589

Kremer A, Hipp AL (2020) Oaks: an evolutionary success story. New Phytol 226:987-2011. https://doi.org/10.1111/nph.16274

Kremer A, Dupouey JL, Deans JD, Cottrell J, Csaikl U, Finkeldey R et al (2002) Leaf morphological differentiation between *Quercus robur* and *Quercus petraea* is stable across western European mixed oak stands. Ann for Sci 59(7):777-787. https://doi.org/10.1051/forest:2002065

Kremer A, Delcamp A, Lesur I, Wagner S, Rellstab C, Guichoux E, Leroy T (2024a) Genotypic data of near-diagnostic genetic markers for four western European white oak species identification. [dataset]. Recherche Data Gouv, V2. https://doi.org/10.57745/OJYLZU

Kremer A, Delcamp A, Lesur I, Wagner S, Rellstab C, Guichoux E, Leroy T, Anonymous A, Curtu AL, Westergren M, Dreyer E (2024b) Reviews of "Whole-genome screening for near-diagnostic genetic markers for four western European white oak species identification". [Open Peer Review Report]. HAL https://doi.org/10.17180/73QH-D818.

Lagache L, Klein EK, Guichoux E, Petit RJ (2013) Fine-scale environmental control of hybridization in oaks. Mol Ecol 22(2):423-436. https://doi.org/10.1111/mec.12121

Lazic D, Hipp AL, Carlson JE, Gailing O (2021) Use of genomic resources to assess adaptive divergence and introgression in oaks. Forests 12:690. https://doi.org/10.3390/f12060690

Le Hardy de Beaulieu A, Lamant T (2006). *Guide illustré des chênes. Tome 1*. Editions du 8ième- Paris 688

Lepais O, Gerber S (2011) Reproductive patterns shape introgression dynamics and species succession within the European white oak species complex. Evolution 65(1):156-170. https://doi.org/10.1111/j.1558-5646.2010.01101.x

Lepais O, Petit RJ, Guichoux E, Lavabre JE, Alberto F, Kremer A et al (2009) Species relative abundance and direction of introgression in oaks. Mol Ecol 18(10):2228-2242. https://doi.org/10.1111/j.1365-294x.2009.04137.x

Lepais O, Roussel G, Hubert F, Kremer A, Gerber S (2013) Strength and variability of postmating reproductive isolating barriers between four European white oak species. Tree Genet Genomes 9(3):841–853. https://doi.org/10.1007/s11295-013-0602-3

Leroy T, Roux C, Villate L, Bodenes C, Romiguier J, Paiva JAP et al (2017) Extensive recent secondary contacts between four European white oak species. New Phytol 214(2):865-878. https://doi.org/10.1111/nph.14413

Leroy T, Louvet JM, Lalanne C, Le Provost G, Labadie K, Aury JM et al (2020a) Adaptive introgression as a driver of local adaptation to climate in
European white oaks. New Phytol 226:1171–1182. https://doi.org/10.1111/nph.16095

Leroy T, Rougemont Q, Dupouey JL, Bodénès C, Lalanne C, Belser C et al (2020b) Massive postglacial gene flow between European white oaks uncovered genes underlying species barriers. New Phytol 226:1183–1197. https://doi.org/10.1111/nph.16039

Lesur I, Alexandre H, Boury C, Chancerel E, Plomion C, Kremer A (2018) Development of target sequence capture and estimation of genomic relatedness in a mixed oak stand. Front Plant Sci 9:996. https://doi.org/10.3389/fpls.2018.00996

Li H (2011) Improving SNP discovery by base alignment quality. Bioinformatics 27:1157–1158. https://doi.org/10.1093/bioinformatics/btr076

Li H (2013). Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM. *Preprint at arxiv*. https://arxiv.org/abs/1303.3997 (2013)

Lopez de Heredia U, Valbuena-Carabana M, Cordoba M, Gil L (2009). Variation components in leaf morphology of recruits of two hybridising oaks Q. petraea (Matt.) Liebl. and Q. pyrenaica Willd. at small spatial scale. Eur J For Res 128(6): 543–554. https://doi.org/10.1007/s10342-009-0302-6

Lupini A, Aci M, Mauceri A, Luzzi G, Bagnato S, Manguzzato G et al (2019) Genetic diversity in old populations of sessile oak from Calabria assessed by nuclear and chloroplast SSR. J Mt Sci 16:1111–1120. https://doi.org/10.1007/s11629-018-5335-1

Macejovsky V, Schmidtova J, Hrivnak M, Krajmerova D, Sarvasova I, Gomory D (2020). Interspecific differentiation and gene exchange among the Slovak Quercus sect. Quercus populations. *Dendrobiology* 83: 20–29. https://doi.org/10.12657/denbio.083.002

Merlino A, Baliva M, Di Filippo A, Piovesan G, Solano F. (2014). Analisi strutturali e dendroecologiche su popolamenti di *Quercus Petraea* subsp. *austrothermica* Brullo, Guarino e Siracusa nel parco regionale delle Madonie (Sicilia). *Second International Congress of Silviculture*, pp 183–189. https://doi.org/10.4129/2cis-am-ana

Modica G (2001). La rovere (Quercus petraea (Matt.) Liebl.) in Aspromonte. Monti e Boschi 3/4: 13–18. https://www.academia.edu/28363486/La_rovere_Quercus.petraea.Matt.Liebl.in_Aspromonte

Moracho E, Moreno G, Jordano P, Hampe A (2016) Unusually limited pollen dispersal and connectivity of pedunculate oak (Quercus robur) refugial populations at the species southern range margin. Mol Ecol 14:3319–3331. https://doi.org/10.1111/mec.13692

Neophytou C (2014) Bayesian clustering analyses for genetic assignment and study of hybridization in oaks: effects of asymmetric phylogenies and asymmetric sampling schemes. Tree Genet Genomes 10(2):273–285. https://doi.org/10.1007/s11295-013-0680-2

Neophytou C, Gartner SM, Vargas-Gaete R, Michiels HG (2015). Genetic variation of Central European oaks: shaped by evolutionary factors and human intervention? *Tree Genetics & Genomes* 11(4). https://doi.org/10.1007/s11295-015-0905-7

Petit RJ, Bodenes C, Ducousso A, Roussel G, Kremer A (2003) Hybridization as a mechanism of invasion in oaks. New Phytol 161(1):151–164. https://doi.org/10.1046/j.1469-8137.2003.00944.x

Plomion C, Aury JM, Amselem J, Leroy T, Murat F, Duplessis S et al (2018) Oak genome reveals facets of long lifespan. Nature Plants 4(7):440–452. https://doi.org/10.1038/s41477-018-0172-3

Raymond M, Roussel F (1995) GENEPOP(version 1.2): population genetics software for exact tests and ecumenicism. J Hered 86:248–249. https://doi.org/10.1093/oxfordjournals.jhered.a111573

Rellstab C, Buhler A, Graf R, Folly C, Gugerli F (2016). Using joint multivariate analyses of leaf morphology and molecular-genetic markers for taxon identification in three hybridizing European white oak species (Quercus spp.). *Annals of Forest Science* 73(3): 669–679. https://doi.org/10.1007/s13595-016-0552-7

Reutimann O, Gugerli F, Rellstab C (2020) A species-discrimnatory single-nucleotide polymorphism set reveals maintenance of species integrity in hybridizing European white oaks (Quercus spp.) despite high levels of admixture. Ann Bot 125:663–676. https://doi.org/10.1093/aob/mca001

Reutimann O, Dauphin B, Baltensweiler A, Gugerli F, Kremer A, Rellstab C (2023). Abiotic factors predict taxonomic composition and genetic admixture in populations of hybridizing white oak species (Quercus sect. Quercus) on a regional scale. *Tree Genetics & Genomes* 19: 22. https://doi.org/10.1007/s11295-023-01598-7

Saleh D, Chen J, Leple JC, Leroy T, Truffaut L, Dencausse B et al (2022) Genome-wide evolutionary response of European oaks during the Anthropocene. Evolution Letters 6(1):4–20. https://doi.org/10.1002/evl3.269

Schroeder H, Kersten B (2023). A small set of nuclear markers for reliable differentiation of the two closely related oak species Quercus robur and Q. petraea. Plants 12 (3), 566. https://doi.org/10.3390/plants12030566

Scotti-Saintagne C, Mariette S, Porth I, Goicoechea PG, Barreneche T, Bodenes K et al (2004). Genome scanning for interspecific differentiation between two closely related oak species Quercus robur L. and Q. petraea (Matt.) Lieb. Genetics 168(3): 1615–1626. https://doi.org/10.1534/2Fgenetics.104.026849

Truffaut L, Chancerel E, Ducousso A, Dupouey JL, Badeau V, Ehrenmann F et al (2017) Fine-scale species distribution changes in a mixed oak stand over two successive generations. New Phytol 215(1):126–139. https://doi.org/10.1111/nph.14561

Urban L, Holzer A, Jotautas Baronas J, Hall MB, Braeuninger-Weimer P, Scherm MJ, Kunz DJ, Perera SN, Martin-Herranz DE, Tipper ET, Salter SJ, Stamnitz MR (2021) Freshwater Monitoring by Nanopore Sequencing Life 10:e61504. https://doi.org/10.7554/eLife.61504

Valbuena-Carabana M, Gonzalez-Martinez S, Sork V, XCollada C, Soto A, PGoi-coechea P, et al (2005) Gene flow and hybridisation in a mixed oak forest (Quercus pyrenaica Willd. and Quercus petraea (Matt.) Liebl.) in Central Spain. Heredity 95:457–465. https://doi.org/10.1038/sj.hdy.6800752

Vazquez-Pardo F, Rincon-Hercules S, Gutierrez-Esteban M, Garcia-Alonso M, Marquez-Garcia F, Ramos -Maqueda S et al. (2009). Congreso Forestal Espanol. Montes y sociedad: Saber que hacer. Leon SECFJdCy (ed.). Socieda Espanola de Ciencias Forestales: Avila, pp 3–13

Viscosi V, Lepais O, Gerber S, Fortini P (2009) Leaf morphological analyses in four European oak species (Quercus) and their hybrids: a comparison of traditional and geometric morphometric methods. Plant Biosystems 143(3):564–574. https://doi.org/10.1080/11263500902723129

Wagner S, Seguin-Orlando A, Leplé JC, Leroy T, Lalanne C, Aury JM, Poirier S, Wincker P, Plomion C, Kremer A, Orlando L (2024) Tracking population structure and phenology through time using ancient genomes from waterlogged white oak. Mol Ecol 33(3):1–17. https://doi.org/10.1111/mec.16859

Yucedag C, Gailing O (2013). Morphological and genetic variation within and among four *Quercus petraea* and *Q. robur* natural populations. Turk J Bot 37(4): 619–629. https://doi.org/10.3906/bot-1205-18

## Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.