Mycorrhiza DOI 10.1007/s00572-013-0540-9

ORIGINAL PAPER

Population genetics of the westernmost distribution of the glaciations-surviving black truffle Tuber melanosporum Iván García-Cunchillos & Sergio Sánchez & Juan José Barriuso & Ernesto Pérez-Collazos

Received: 6 June 2013 / Accepted: 29 October 2013 # Springer-Verlag Berlin Heidelberg 2013

Abstract The black truffle (Tuber melanosporum Vittad.) is an important natural resource due to its relevance as a delicacy in gastronomy. Different aspects of this hypogeous fungus species have been studied, including population genetics of French and Italian distribution ranges. Although those studies include some Spanish populations, this is the first time that the genetic diversity and genetic structure of the wide geographical range of the natural Spanish populations have been analysed. To achieve this goal, 23 natural populations were sampled across the Spanish geographical distribution. ISSR technique demonstrated its reliability and capability to detect high levels of polymorphism in the species. Studied populations showed high levels of genetic diversity (h N =0.393, h S =0.678, Hs =0.418), indicating a non threatened genetic conservation status. These high levels may be a consequence of the wide distribution range of the species, of its spore dispersion by animals, and by its evolutionary history. AMOVA analysis showed a high degree of genetic structure among populations (47.89 %) and other partitions as geographical ranges. Bayesian genetic structure analyses

Iván García-Cunchillos and Ernesto Pérez-Collazos contributed equally to this work. Electronic supplementary material The online version of this article (doi:10.1007/s00572-013-0540-9) contains supplementary material, which is available to authorized users. I. García-Cunchillos : E. Pérez-Collazos (*) Departamento de Ciencias Agrarias y del Medio Natural, Escuela Politécnica Superior de Huesca, Universidad de Zaragoza, C/ Carretera de Cuarte km 1 s/n, 22071 Huesca, Spain e-mail: [email protected] S. Sánchez : J. J. Barriuso Centro de Investigación y Tecnología Agroalimentaria, Gobierno de Aragón, Avenida de Montañana 930, 50059 Zaragoza, Spain

differentiated two main Spanish groups separated by the Iberian Mountain System, and showed the genetic uniqueness of some populations. Our results suggest the survival of some of these populations during the last glaciation, the Spanish southern distribution range perhaps surviving as had occurred in France and Italy, but it is also likely that specific northern areas may have acted as a refugia for the later dispersion to other calcareous areas in the Iberian Peninsula and probably France. Keywords Autocorrelations . Bayesian population analysis . Biogeography . Genetic diversity and structure . Iberian refugia . ISSR . Interglacial expansion

Introduction The black truffle (Tuber melanosporum Vittad.) is a hypogeous ectomycorrhizal ascomycete that produces edible fruit bodies, which are a delicacy in increasing demand in gastronomy due to its intense flavour and taste, causing it to be an important economic resource. To cope with the high demand of the food market, the cultivation of the black truffle in manmade fields has experienced great development since the 1970s (Chevalier and Sourzat 2013). Nevertheless, due to human pressure, the natural production of black truffle decreases every year, compromising the survival of some small natural populations (Hall et al. 2007); therefore, it is important to study its genetic conservational status. The geographical distribution of T. melanosporum is restricted to the Mediterranean Basin, and its presence is closely linked to continental weather conditions, with cold winters and storm activity during the hot summer (Chevalier and Sourzat 2013). It is also directly related to the availability of limestone soils, which explains the disjointed distribution of

Mycorrhiza

the black truffle across the territory, and is associated with open forests, allowing the direct insulation of its symbiotic species trees (mainly Quercus ilex L. subsp. ballota (Desf.) Samp., Quercus faginea Lam., Quercus pubescens Willd. and Quercus coccifera L; Hall et al. 2007) from the sun. The largest numbers of natural populations of T. melanosporum are found in Spain, France and Italy. Nonetheless, the presence of this species has been reported infrequently in other European countries (Riousset et al. 2001). At the Iberian Peninsula, the species is distributed mainly in the southern and middle western part, covering western Andalucía, Aragón, Cataluña, western Castilla y León, Castilla La Mancha, Comunidad Valenciana, Murcia, Navarra, La Rioja and Basque Country (Reyna 2011). Some authors have hypothesised about the origin, evolution and migration routes that shaped the present distribution of T. melanosporum (Bertault et al. 1998, 2001; Murat et al. 2004; Riccioni et al. 2008). The main hypothesis places the origin of the Melanosporum group at about 52 to 79 Ma (Jeandroz et al. 2008; Bonito et al. 2013), being followed by a speciation process that could have occurred from two independent dispersals in Asia and Europe, the first producing Tuber brumale (Europe) and Tuber pseudoexcavatum (China), and the second allowing the speciation of T. melanosporum in Europe and Tuber indicum in Asia, 19 to 25 Ma (Jeandroz et al. 2008; Bonito et al. 2013). It is possible that the populations of T. melanosporum had suffered a substantial reduction in size and distribution during glaciations, and as a consequence, only the southern distribution areas survived, as had been demonstrated in several fauna and flora species (Kaltenrieder et al. 2009; Dapporto 2010). After the last glaciations, T. melanosporum probably colonised northern areas until reaching its current geographical distribution (Bertault et al. 1998, 2001; Murat et al. 2004; Riccioni et al. 2008). Due to its economic importance and its biological relevance, studies of the species have been conducted in several different fields. Molecular techniques have been employed to improve cultivation and production, providing quick identification of T. melanosporum versus morphologically similar species (Paolocci et al. 1997; Pérez-Collazos et al. 2010; Cordero et al. 2011), detecting DNA concentration in soil (Suz et al. 2006; Parladé et al. 2013) and evaluating organoleptic properties (Talou et al. 2001; Splivallo et al. 2012). Also, reproduction and formation of mycelium and fruiting bodies have been studied and discussed (Riccioni et al. 2008; Rubini et al. 2011; Iotti et al. 2012; Linde and Selmes 2012). Population genetic analyses have been conducted over French and Italian distribution ranges, including some Spanish populations (1 to 7), evaluating their genetic diversity, structure and inferring evolutionary processes (Bertault et al. 1998, 2001; Murat et al. 2004; Riccioni et al. 2008; Martin et al. 2010; Murat et al. 2011), and recently, the spatial genetics of two Italian and French plantations were studied (Murat et al.

2013). Nevertheless, no genetic population study has been conducted over the wide distribution Iberian range of T. melanosporum, the westernmost of the species. Therefore, our aims, applied to 23 natural Spanish populations, are to (1) measure the genetic diversity to evaluate the genetic status of conservation of populations, (2) detect the genetic structure among populations and geographical ranges, (3) infer the genetic relationship among populations and individuals and (4) infer possible evolutionary processes which might have shaped the diversity and the genetic structure of the Iberian distribution range of T. melanosporum . To accomplish these aims, we used inter-simple sequence repeats (ISSR; Zietkiewicz et al. 1994) markers, a PCR-based technique with primers designed from microsatellite sequences with one to three selective bases anchored at the 3′ or 5′ ends of the primer. This method combines higher annealing temperatures and more selective primers than RAPD, which enhances reproducibility (Ge and Sun 1999; Nybom 2004). Dominant ISSR technique is quick and capable to detect high levels of polymorphism which allows easy interpretation of reliable results (Bornet and Branchard 2001; Badfar-Chaleshtori et al. 2012), and has demonstrated its usefulness in determining the influence of life history traits in the evolution and dynamics of populations (Pérez-Collazos and Catalán 2007; Coupé et al. 2011; Lurá et al. 2011; Neal et al. 2011; Shao et al. 2011).

Materials and methods Population sampling, DNA extraction and ISSR analysis One hundred and ninety individuals from 23 natural populations of T. melanosporum were sampled from throughout the main distribution range of the species in the Iberian Peninsula (Table 1, Fig. 1). The sampling included nine populations in the northern Pre-Pyrenees mountains, and small sierras across the Iberian peninsula as follows: Sierra de Izco (1 population), de Cabreja (1 population), Moncayo (1 population), Alto Tajo (1 population), Matarraña (2 populations), Javalambre (1 population), Cuenca (1 population) and six populations of the Subbetic mountains in the southern area of distribution. The number of sampled fruit bodies of each population ranged from 4 to 14, but in the PP4 and PP8, only 2 sporocarp samples were available (Table 1). Each sample was analysed macro- and microscopically to confirm its belonging to T. melanosporum, following the descriptions given by Riousset et al. (2001). Truffles were cleaned with ethanol and flamed to eliminate microorganisms present in the superficies of the peridium. Twenty grams of this material were sliced, lyophilised and pulverised. We did not apply specific methodologies to broke ascoparps (Riccioni et al. 2008). T. melanosporum is an outcrossing (heterothallic) species, but the vast majority of the tissue of the gleba belongs from one parental, therefore

Mycorrhiza Table 1 Code, administrative province, sampled location and genetic diversity of the 23 populations of T. melanosporum studied Code

Province

Location

N

hN

hS

Hs

fe

fd

PLP99 % (2)

Br(2)

IZ PP1 PP2 PP3 PP4 PP5 PP6 PP7 PP8 PP9 CA MO AT MA1 MA2 JA CU

Navarra Huesca Huesca Huesca Lérida Lérida Lérida Lérida Lérida Barcelona Soria Soria Guadalajara Castellón Castellón Castellón Cuenca

Sierra de Izco Sierra de Loarre Graus1 Graus2 Sierra de Comillini Sierra de San Gervasio Sierra de Carreu Alt Urgell1 Alt Urgell2 Osona Sierra de Cabrejas Sierra de Moncayo Alto Tajo Sierra de Matarraña1 Sierra de Matarraña2 Sierra de Javalambre Serranía de Cuenca

14 4 9 10 2 8 7 12 2 7 12 7 11 7 9 6 6

0.281 0.797 0.285 0.326 0.203 0.250 0.297 0.266 0.272 0.353 0.296 0.301 0.292 0.324 0.354 0.379 0.668

0.742 1.106 0.627 0.750 0.140 0.519 0.578 0.660 0.189 0.688 0.736 0.587 0.699 0.630 0.778 0.679 1.196

0.41796 0.418061 0.417988 0.418003 0.417959 0.417870 0.418020 0.417955 0.417989 0.418214 0.417810 0.417848 0.417886 0.417977 0.418042 0.418081 0.418206

1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

0.342 0.127 0.456 0.506 0.000 0.253 0.354 0.392 0.000 0.544 0.620 0.443 0.392 0.506 0.835 0.380 0.418

1.115 1.074 1.148 1.189 1.101 1.097 1.157 1.147 1.177 1.244 1.197 1.176 1.135 1.227 1.276 1.160 1.205

SE AB1 AB2 CZ1 CZ2 MG TOTAL

Albacete Murcia Murcia Jaén Jaén Jaén

Sierra de Segura Sierra de los Álamos y del Buitre1 Sierra de los Álamos y del Buitre2 Sierra de Cazorla1 Sierra de Cazorla2 Sierra de Magina

6 10 11 8 10 13 190

0.326 0.306 0.312 0.347 0.340 0.293 0.393

0.583 0.705 0.749 0.722 0.783 0.751 0.678

0.418073 0.418136 0.418155 0.418180 0.418214 0.418006 0.418308

0 0 0 1 1 1 –

0 0 0 0 0 0 –

0.253 0.316 0.519 0.494 0.481 0.519 –

1.103 1.115 1.218 1.204 1.169 1.157 –

N sample size, h N Nei's genetic diversity index, h S normalised Shannon's diversity index, H S Bayesian estimation of the genetic diversity, f e exclusive fragments, f d diagnostic fragments, PLP99%(2) percentage of polymorphic loci at 99 % with a rarefaction index of two, Br(2) band richness with a rarefaction index of two

individual gleba tissue can be considered as haploid (Riccioni et al. 2008). Total DNA was extracted from 20 mg of this powder using the DNeasy Plant mini kit (QIAGEN) with minor modifications to the protocol. ISSR procedure followed the PCR protocol of Wolfe et al. (1998), based on Zietkiewicz et al. (1994), with some modifications. In brief, the PCR cocktail (11 μl) included 1.25 μl Buffer (10X), 1.25 μl MgCl2 (25 mM), 0.25 μl dNTPs (10 mM), 1.25 μl primer (10 μM), 0.2 μl KapaTaq DNA polymerase (5 U/μL; Kapa Biosystems), 5.8 μl Milli-Q H 2 O, and 1 μl DNA (20–40 ng/μl). The thermocycler program for PCR was set for 7 min at 94 °C followed by 30 cycles each of 30 s at 94 °C, 45 s at 50 °C, 60 s at 72 °C and a final extension cycle of 7 min at 72 °C. Sixteen primers designed at the University of British Columbia laboratory (UBS primer set 9) were tested on five sporocarps of T. melanosporum . Four primers [(CAC)4RC, (AG)7YC, GGGC(GA)8 and CAA(GA)5] were chosen to assay the entire set of samples on the basis of the number of scoreable bands, the degree of polymorphism and the consistency of the results. Amplification products were separated on

3 % agarose gels, stained with Midori Green DNA (NIPPON Genetics EUROPE GmbH). Electrophoresis was carried out at a constant 160 V for 2 h. Data analysis ISSR bands were scored as presence (1) or absence (0). Nei (1973) and normalised Shannon genetic indices were estimated for each population and for the species. Also, Bayesian genetic diversity, defined as average panmictic heterozygosity within each population and its average value across populations (Hs) was computed with the HICKORY v.1.0.4 programme (Holsinger and Lewis 2003), assuming that the statistical inferences could be extrapolated to other potential populations not included in the present study (Holsinger 1999). Prior distributions were specified using the default values of the software (burnin set to 50,000, sampling set to 250,000 and thin set to 50). ISSR data were analysed under the free model that does not assume Hardy–Weinberg equilibrium and is therefore well suited for analysing dominant markers (Holsinger et al. 2002). As an additional measure of genetic

Mycorrhiza

NE structure

PP2 PP4 PP7

IZ MO

CA

PP1

PP6 PP5

AT CU

PP8

PP3 PP9

MA1 MA2 JA

CZ2 CZ1

MG SE

SE structure AB1 AB2

Fig. 1 Geographical distribution of T. melanosporum across the Iberian Peninsula (grey area). The map shows the location of the representative sampling of the 23 natural populations studied and the Bayesian Inference of the genetic structure of two main geographical areas (K =2; northeastern and south eastern ranges) detected in the Iberian populations of T. melanosporum. Codes IZ Sierra de Izco, PP1 Sierra de Loarre, PP2 Graus1, PP3 Graus2, PP4 Sierra de Comillini, PP5 Sierra de San

Gervasio, PP6 Sierra de Careu, PP7 Alt Urgell1, PP8 Alt Urgell2, PP9 Osona, CA Sierra de Cabrejas, MO Sierra de Moncayo, AT Alto Tajo, MA1 Sierra de Matarraña1, MA2 Sierra de Matarraña2, JA Sierra de Javalambre, CU Serranía de Cuenca, SE Sierra de Segura, AB1 Sierra de los Álamos y del Buitre1, AB2 Sierra de los Álamos y del Buitre2, CZ1 Sierra de Cazorla1, CZ2 Sierra de Cazorla2, MG Sierra de Magina. Province correspondence of codes are given in Table 1

diversity, the number of exclusive (f e ) and diagnostic (f d , frequencies over 0.99) fragments were calculated directly from the binary data matrix. Due to the limitations of having only two individuals available for sampling in the PP4 and PP8 populations, (1) the percentage of polymorphic loci (PLP ) with a significance of 1 % (p =0.99) and (2) the band richness (Br), which is the number of phenotypes expected at each locus and can be interpreted as an analogue of the allelic richness, ranging from 1 to 2 (Coart et al. 2005), were calculated. Br and PLP99 % indices were estimated according to the rarefaction method of Hurlbert (Petit et al. 1998), and restricted to the smallest studied population (2) using the software AFLPdiv v1.0 (http://www.pierroton.inra.fr/genetics/labo/Software/ Aflpdiv/). A Bayesian model-based analysis was performed to infer the genetic structure of the populations and the spatial ranges of T. melanosporum using STRUCTURE v.2.3.4 (Pritchard et al. 2000; Hubisz et al. 2009). We imposed the nonadmixture ancestry model with independent allele frequencies. We ran the analysis for K values ranging from 1 to 25, using a burnin period and a run length of the Monte Carlo Markov Chain of 50,000 and 100,000 iterations, respectively. Chain convergence was estimated through visual inspection of the posterior values excluding the burnin. Ten iterations were conducted using the ad hoc parameter ΔK of Evanno

et al. (2005) to estimate the rate of change of likelihood values between successive K values. Genetic substructuring within the inferred spatial ranges was further assessed through independent analyses of the split data matrices using the same procedures indicated above (for K groups ranging from 1 to the number of analysed populations plus two). As an additional measure of the genetic structure, analyses of molecular variance (AMOVA) were performed using Arlequin 3.11 (Excoffier et al. 2005). The AMOVAs were conducted using (1) the entire dataset considering all populations belonging to the same group (species level) and (2) using geographical division dataset to infer the partition of the variance within geographical groups. ISSR genetic relationships among all studied individuals were assessed by principal coordinates analysis (PCO) and Neighbour-Joining clustering methods (NJ; Saitou and Nei 1987). Two- and three-dimensional PCO plots were constructed based on the Nei and Li (1979) distance matrix using NTSYSpc v. 2.1 (Rohlf 2000). Relationships among individuals were also visualised through NJ clustering using MEGA 4.0.1 (Tamura et al. 2007). Numbers of migrants between pairs of populations were estimated using Arlequin following the procedure of Slatkin (1995). Isolation by distance was estimated by comparing the matrix of geographical distances between populations with an FST pairwise distance matrix obtained from Arlequin using the

Mycorrhiza

Mantel test (1967), with 1,000 non-parametric permutations in NTSYS software. To evaluate the genetic relationships between pairs of populations or regions, autocorrelograms between eight different groups were calculated (each group was defined according to the genetic structure results and the geographical distribution as follows: group 1, PP1 to PP9; group 2, IZ; group 3, CA and MO; group 4, AT; group 5, MA1-2 and JA; group 6, CU; group 7, AB1-2 and SE; group 8, CZ1-2 and MG). Briefly, the genetic ISSR distance matrices were correlated with binary matrices encoded according to the corresponding geographical relationship to be tested: 0 for the inter-populations distances within geographical ranges, and 1 for the rest of the inter-population distances (Oden and Sokal 1986; Stehlik et al. 2001; Pérez-Collazos and Catalán 2006). The results were plotted on a graph by pairs of groups in order to summarise them. Mantel tests of the autocorrelograms were computed with 1,000 nonparametric permutations in NTSYS.

Results Genetic diversity and structure of the Iberian populations of T. melanosporum The four chosen primers for the analysis provided 77 polymorphic bands. Each individual showed a single genetic fingerprint. Despite the remote risk of amplifying bacterial DNA using ISSR, due to its presence in some Tuber ascomata (Barbieri et al. 2000, 2010), the putative amount of bacterial DNA in fresh fruiting bodies may be low in respect to T. melanosporum DNA (Rivera et al. 2010) and probably it does not compete in a PCR reaction. Genetic diversity indices of Nei, Shannon, and the Bayesian heterozygosity (Table 1) showed high levels of genetic diversity throughout the Iberian populations (h N =0.393, h S =0.678 and H S =0.418). At a population level, the highest levels of genetic diversity were detected at PP1 (h N = 0.797, h S = 1.106 and H S = 0.418061) and CU (h N = 0.668, h S = 1.196 and H S =0.418206) populations. Population PP4 showed the lowest values for genetic diversity (h N =0.203; h S =0.140 and H S =0.417959). The Bayesian heterozygosity index also detected high levels of genetic diversity in the southern populations of CU, SE, AB1-2, CZ1-2, and MG (Table 1). Two exclusive fragments were detected among the studied populations. One fragment was exclusive and diagnostic of the northern Sierra de Izco population (IZ) and the other was shared by the southernmost populations of CZ1-2 and MG (Table 1). The rarefaction method applied on the band richness and percentage of polymorphic loci indices reflects high levels of genetic diversity of the MA2 (PLP99 % =0.835; Br =1.276) population, as well as for populations PP9 and AB2 (Table 1).

The AMOVA analysis revealed that a high percentage of genetic variation is found between the 23 studied populations of T. melanosporum (47.89 %), while the rest of the variance is distributed within populations (52.11 %). Genetic structure Bayesian inference analysis, among all inferred groups, detected a division in two hypothetical populations (K =2), differentiating a northeast (NE) and a southeast (SE) group (Fig. 1). The hierarchical AMOVA for this grouping showed values of 16.34, 35.87 and 47.78 % for the partition between NE and SE geographical groups, among populations and within populations, respectively. In order to discriminate the intra-genetic structure of each regional group, Bayesian analyses of partial matrices corresponding to each group were conducted. Bayesian analysis of the NE group detected two different genetic subgroups (K =2); the NE1 subgroup formed by the Pre-Pyrenees populations (PP1 to PP9), and the NE2 subgroup by the population of Sierra de Izco (IZ) and almost the whole populations of the central Iberian range (CA, MO, AT, MA1, MA2; Fig. 2a). The partition of the variance of the NE1 and NE2 subgroups was 14.69, 34.34 and 50.97 % among groups, among populations and within individuals, respectively. Bayesian genetic structure of NE1 (K = 2: NE1.1 integrated by PP4 to PP8 populations; and NE1.2 integrated by PP1 to PP3 and PP9 populations) is shown in Fig. 2c. AMOVA indicated a partition of variance among groups and populations of 23.41 and 17.48 %, respectively. The NE2 subgroup was divided into four genetic groups (K =4: NE2.1 integrated by IZ, NE2.2 integrated by CA and MO, NE2.3 integrated by MA1-2 and NE2.4 integrated by AT population; Fig. 2d), showing a partition of the variance of 28.18 and 18.73 % among groups and populations, respectively. The Bayesian analysis of the SE group differentiated the southernmost populations (SE1 group: CZ1-2 and MG; Fig. 2b) from the rest (SE2 group: AB1-2, JA, CU and SE; Fig. 2b). AMOVA analysis for this grouping was 21.72, 20.87 and 57.41 % for the three levels of partitioning. Other Bayesian analyses of partial matrices of the SE1 and SE2 groups showed the shared genetic structure of JA, CU and SE populations (SE2.1. group; Fig. 2f), the genetic identity of the Sierra de los Álamos y del Buitre populations (SE2.2 group; Fig. 2f) and the genetic substructure of the CZ1 population (Fig. 2e). AMOVA analyses showed a genetic structure among populations lower than 19 % in all subgroups (results not shown). Genetic relationship among individuals of T. melanosporum Neighbour-Joining analysis evidenced a separation of the populations of NE (Sierra de Izco, Pre-Pyrenees and central Iberian range) and SE (Serranía de Cuenca, Sub-betic Mountains, Sierra de Javalambre and one population of Sierra de Matarraña) Spain (Fig. 3). The analysis showed high

Mycorrhiza

NE2 structure

a

AT A

MA1

PP9 9

PP2 2

PP7 7

NE1.2

MG

AB1 SE AB2

AB2

SE2.2 SE

CU

SE2.1 1.00 0.80 0.60 0.40 0.20 0.00

1.00 0.80 0.60 0.40 0.20 0.00

PP4 PP8

CZ2

JA

f

c

PP5 5

MA1 MA2

PP6 6

NE1.1

CU

CZ1

CZ2 2

G MG

CZ1

SE1.1

PP9

PP1

PP6

AB1

e

PP3 PP5

JA 1.00 0.80 0.60 0.40 0 40 0.20 0.00

PP8

PP3 3

MO

PP1

AT

AB1

AB2

SE

CU

JA

CA CZ2

PP2 PP4 PP7

IZ

SE1 structure SE2 structure

MG

1.00 0.80 0.60 0.40 0 40 0.20 0.00

MA2

CA

IZ Z

b

NE2.3 NE2.4 M MO

NE2.1 NE2.2

CZ1

d 1.00 0.80 0.60 0.40 0 20 0.20 0.00

NE1 structure

Fig. 2 Bayesian inference of the genetic sub-structure of the two main geographical areas detected by STRUCTURE (Fig. 1) of the Spanish populations of T. melanosporum. a Northeastern range substructure (NE1 and NE2), b southern range substructure (SE1 and SE2), c northeastern NE1.1 and NE1.2 range substructure, d northeastern NE2.1– NE2.4 range substructure, e southern SE1.1 range substructure and f southern SE2.1 and SE2.2 range substructure. Populations are indicated on the abscise axis and the belonging proportions of each group are indicated on the ordinate axis. Codes IZ Sierra de Izco, PP1 Sierra de

Loarre, PP2 Graus1, PP3 Graus2, PP4 Sierra de Comillini, PP5 Sierra de San Gervasio, PP6 Sierra de Careu, PP7 Alt Urgell1, PP8 Alt Urgell2, PP9 Osona, CA Sierra de Cabrejas, MO Sierra de Moncayo, AT Alto Tajo, MA1 Sierra de Matarraña1, MA2 Sierra de Matarraña2, JA Sierra de Javalambre, CU Serranía de Cuenca, SE Sierra de Segura, AB1 Sierra de los Álamos y del Buitre1, AB2 Sierra de los Álamos y del Buitre2, CZ1 Sierra de Cazorla1, CZ2 Sierra de Cazorla2, MG Sierra de Magina. Province correspondence of codes are given in Table 1

genetic structure among populations due to the grouping of all individuals of one population into one cluster. In this sense, the populations that shared the same genetic structure (Fig. 2) showed some intermingling of individuals; southern populations of CZ1-2 and MG and central Iberian range populations of CA and MO. Nonetheless, only five individuals (one, one, and three individuals from PP1, JA and CZ1, respectively) intermingle with other individuals of other populations of the same geographical range. NJ analysis revealed a moderate correspondence between geography and genetic distances; nevertheless, the position of MA2 in the NJ (between the northern and central populations of IZ and AT, respectively) does not match the geographical distribution of the species. Multivariate analysis of PCO agrees with the STRUCTURE and NJ results, showing a grouping of individuals into two main geographical groups (NE and SE Spain; Supplementary data). This analysis also showed a strong clustering of individuals by populations and geographical

ranges, and the genetic uniqueness of the populations of PP4-8 and CZ1-CZ2-MG in the northern and southern distribution range, respectively. Despite the relatively wide dispersal of some individuals of the PP9 and CU populations, no intermingling among populations was detected. The first three axes of the PCO analysis explained 28.7 % of total variation.

Gene flow and isolation analyses No gene flow was detected between geographical NE and SE Spain, nor at geographical ranges inside the NE and SE, with the exception of PP1 population which presents small positive values with four populations (PP6, PP7, CA and MO; results not shown). Some nearby populations showed gene flow; 1 to 25 migrants were detected among PP1 to PP3 populations, 1 to 20 among PP4 to PP8, 1 to 14 inside the Sierra de los

Mycorrhiza

Discussion CA CZ1 2

MG

Genetic diversity of the Iberian T. melanosporum populations

MO

SE PP1 3 SE Spain

NE Spain

AB1 2

CU PP4 8

JA PP9

MA1

IZ AT

MA2

2

Fig. 3 Neighbour-Joining tree showing a high genetic structure between individuals of the different populations of T. melanosporum, and it genetic relationships in the Iberian Peninsula. Codes IZ Sierra de Izco, PP1 Sierra de Loarre, PP2 Graus1, PP3 Graus2, PP4 Sierra de Comillini, PP5 Sierra de San Gervasio, PP6 Sierra de Careu, PP7 Alt Urgell1, PP8 Alt Urgell2, PP9 Osona, CA Sierra de Cabrejas, MO Sierra de Moncayo, AT Alto Tajo, MA1 Sierra de Matarraña1, MA2 Sierra de Matarraña2, JA Sierra de Javalambre, CU Serranía de Cuenca, SE Sierra de Segura, AB1 Sierra de los Álamos y del Buitre1, AB2 Sierra de los Álamos y del Buitre2, CZ1 Sierra de Cazorla1, CZ2 Sierra de Cazorla2, MG Sierra de Magina. Province correspondence of codes are given in Table 1

Álamos y del Buitre area (AB1 and AB2), 4 inside CA and MO populations and 1 to 2 migrants among southern range of CZ1, CZ2 and MG populations (results not shown). Correlation between the genetic distance matrix and the geographical matrix was positive and statistically significant (r=0.475; p