Regional environmental and climatic concerns on preserving native gene pools of a least concern species: Brown trout lineages in Mediterranean streams Manuel Vera a , Enric Aparicio b ,SandraHeras c ,AlbaAbras c ,AdriánCasanova c , Maria-Inés Roldán c , Jose-Luis García-Marin c, ⁎ a Departamento de Zoología, Xenética e Antropología Física, Campus Lugo, Universidad de Santiago de Compostela, 27002 Lugo, Spain b GRECO, Institute of Aquatic Ecology, University of Girona, 17071 Girona, Spain c Laboratori d'Ictiologia Genètica, Campus Montilivi, Universitat de Girona, 17003 Girona, Spain HIGHLIGHTS GRAPHICAL ABSTRACT •Mediterranean trout populations often reinforced with Atlantic hatchery stocks •Variable amounts of Introgressive hybridization of Atlantic genes among populations •Amount of introgression relies on low habitat quality rather than on local stocking •Warmer temperatures and lower precipitation also favored hybridization •Ongoing climatic change will promote the spread of these hybridized populations ABSTRACTARTICLE INFO Editor: Sergi Sabater Keywords: Biodiversity conservation Brown trout Climate change Genetic introgression Inland fisheries management Mediterranean streams The European brown trout, Salmo trutta,isacold-adaptedfish reported as a Least Concern species in the IUCN Red List. This species colonized new territories from southern refuges during the last glacial melting, but during the 20th century suffered from anthropic impacts on its habitats. The long-time survival of the species relies on the genetic diversity within and among populations. Brown trout is among the genetically most diverse vertebrate species; however, native populations in Mediterranean rivers have dramatically suffered of introgressive hybridization from extensive releases of evolutionary distant nonnative Atlantic stocks. In addition, in Mediterranean rivers climate change will result in unsuitable conditions for the species during the 21st century. Using brown trout populations at the headstreams of a Pyrenean river as a model, this paper revised how hatchery releases have affected the native gene pools and how environmental and climatic variables controlled the amount of local introgression at intra-basin level. Introgressive hybridization was detected in all studied sites. Ten times larger divergence was observed among populations at tributaries than among populations along the main stem. A highly impacted population distributed in a long transect in the main stem suggested that hatchery fish move towards the main stem wherever released. From already highly impacted populations and despite the cessation of hatchery releases, warmer temperatures and lower precipitation expected from climate change will extend the introgressive hybridization along the basin, contributing to the extinction of the native gene pools. Based on available morphological distinction of native, hatchery and hybrid brown trout, we advocate the involvement of regional social groups (e.g. riverside dwellers,anglers,conservationists, hikers) in citizen science programs to detect the spread of non-native phenotypes along the rivers. These are cheap and fast methods to collaborate with fishery managers in the preservation and recovery of the regional native populations. Science of the Total Environment 862 (2023) 160739 ⁎Corresponding author. E-mail addresses:
[email protected] (M. Vera),
[email protected] (S. Heras),
[email protected] (A. Abras),
[email protected] (M.-I. Roldán),
[email protected] (J.-L. García-Marin). http://dx.doi.org/10.1016/j.scitotenv.2022.160739 Received 31 August 2022; Received in revised form 10 November 2022; Accepted 3 December 2022 Available online 8 December 2022 0048-9697/© 2022 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY-NC license (http://creativecommons.org/licenses/by-nc/4.0/). Contents lists available at ScienceDirect Science of the Total Environment journal homepage: www.elsevier.com/locate/scitotenv
1. Introduction Freshwater ecosystems were the most affected by biodiversity losses during the period 1970–2012 (Almond et al., 2020). In Europe, freshwater fish species are within the groups more threatened by pollution, water extraction, habitat fragmentation, introduction of exotic species, overfishing, and releases of non-native stocks (Cowx and Gerdeaux, 2004;Freyhof and Brooks, 2011). However, the impact of climate change will become into the main threat throughout the 21st century (Sala et al., 2000). Predictions suggest the extinction of 75 % of the diversity in freshwater fish due to the reduction of rain discharge (Xenopoulos et al., 2005), and the increase of water temperature with the consequent reduction of oxygen availability (Ficke et al., 2007). Tisseuil et al. (2012) modelled the impact of climate change on freshwater biodiversity in Southern Atlantic French rivers during the 21st century and predicted the reduction in the distribution of cold-adapted fishes such as the stone loach (Barbatula barbatula) or the common minnow (Phoxinus phoxinus) requiring oxygenated waters in mountainous areas. The prediction was even more dramatic for the resident forms of brown trout (Salmo trutta) because of its restricted occurrence in headstreams. In Europe, climate change outcomes will be more accentuated in the Mediterranean region (Ormerod, 2009) and salmonid populations inhabiting headstreams in Mediterranean drainages exemplify all these threats over the freshwater native biodiversity (e.g. Almodóvar et al., 2012;Alonso et al., 2011;Freyhof et al., 2020;Nicola et al., 2009;Vera et al., 2013). Salmonids are highly appreciated by anglers in Mediterranean countries and foster regional tourism industries during summer that might be increased by promoting these fisheries for EU tourism (Torrent, 2015). However, some studies advised that ongoing climatic change might dramatically affect the sustainability of these fisheries (e.g. Almodóvar et al., 2012;Ayllon et al., 2016, 2021;Buisson and Grenouillet, 2009). The FAO European Inland Fisheries and Aquaculture Advisory Commission urged the beginning of actions focusedontheadaptation and mitigation of climatechange effects start as soon aspossible (FAO, 2017). The genetic diversity of the species maintains their ability to respond to environmental changes and to maintain abundant populations (Reed and Frankham, 2003). Diverse environmental conditions throughout the species' rank produce specific responses to distinct selective pressures. Local adaptations prospered in each territory by the selection of those genotypes giving advantages to the individuals in a specific environment, regardless of the selective value of these genotypes in other territories. Often, salmonid species distribute over wide geographical areas, environmentally diverse, determining a suitable scenario for local adaptation at different geographical scales (Fraser et al., 2011). Local selection on genes associated with variation in sexual maturation, energy homeostasis, and immune defense has been observed in Scandinavian salmon populations in the Teno River basin (Pritchard et al., 2018). Local and regional gene-environmental associations look spread throughout the genome among North-European brown trout populations (Bekkevold et al., 2020). Altogether confirming that the management of populations in the wild should be focused on the conservation of the genetic diversity at intraspecific level and, hence, local populations should be the basic units for management and conservation. In fact, by protecting local populations, the species is conserved (Laikre and Ryman, 1996). Several Salmo species colonized and diversified in Europe following the Pliocene-Pleistocene climatic cooling (Shedko et al., 2013), but likely, brown trout, Salmo trutta, accumulates one of the richest and most diverse genetic population structuring among European freshwater fish species, particularly in the Mediterranean basins (Bernatchez, 2001;García-Marín et al., 1999). This genetic diversity is often related to morphological and ecological diversity that Lelek (1987) considered could be recognized into a single polytypic species, Salmo trutta, widely distributed in Europe. Since then, it was a Least Concern (LC) species on the International Union for Conservation of Nature (IUCN) Red List. However, a more recent revision has split into many different taxa this previous single pan-European species. Currently, only Atlantic populations of brown trout are taxonomically identifiedasS. trutta, and several Mediterranean lineages are now considered asseparated species, many ofthem being listed as Vulnerable (V)or with Data Deficient (DD) in the updated Red List of European freshwater fishes (Freyhof and Brooks, 2011). Unfortunately, in some regions such as the Iberian Peninsula, brown trout lineages have not been subject to enough taxonomic reviewyet, despite they could include different endemic species (Kottelat and Freyhof, 2007), and binomials used on national legislation have not been updated. For instance, the Spanish legislation still considers a single species, S. trutta, inhabiting Atlantic and Mediterranean rivers despite four evolutionary distant mitochondrial lineages (Adriatic, AD, Mediterranean, ME, Atlantic, AT and Duero, DU) are present among these populations (Bouza et al., 2001;Cortey et al., 2004;Machordom et al., 2000;Vera et al., 2010a). Mediterranean lineages represent a longtime divergent branch from the Atlantic one distributed among Central European Atlantic rivers, that were the source of the hatchery stocks used for stocking the Mediterranean rivers. These stocking activities have resulted in introgressive hybridization not only on the native brown trout populations (Almodóvar et al., 2006;Leitwein et al., 2016;Marićet al., 2022;Splendiani et al., 2016) but also on other species of the genus Salmo (Giuffra et al., 1996;Pustovrh et al., 2012), compromising the adaptive potential and hence the long-time survival of the native Salmo biodiversity (Araguas et al., 2009;Giuffra et al., 1996;Marićet al., 2022). The historical long-time isolation of Atlantic and Mediterranean trout lineages has contributed to the accumulation of distinct genetic variants at nuclear loci useful to monitor and evaluate the stocking impact into the Mediterranean basins. In brown trout populations, microsatellite loci (e.g. Araguas et al., 2017;Sanz et al., 2009) and particularly the lactate dehydrogenase C1 locus (LDH-C1*) have been extensively used (Almodóvar et al., 2006;Splendiani et al., 2016). The LDH-C1*100 allele was fixed in native Mediterranean brown trout populations, but the releases of Atlantic stocks introduced the *90 allele, fixed in the North-European populations (Bernatchez, 2001;García-Marín et al., 1991;Hamilton et al., 1989). These stocks also have mtDNA haplotypes of the non-native Atlantic lineage (Cortey and García-Marin, 2002;Cortey et al., 2004) and phenotype traits showing a genetic basis such as spotting pattern and body pigmentation (Blanc et al., 1994;Mezzera et al., 1997;Skaala and Jorstad, 1987)contributed to the distinction between Atlantic hatchery stocks and native Mediterranean brown trout (Aparicio et al., 2005;Lorenzoni et al., 2019; Valette et al., 2022). At present, there is still limited evidence on factors promoting the introgressive hybridization of Mediterranean brown trout populations with Atlantic stocks, and on how the introgression of these gene pools originated from colderregionsmay have affectedthe adaptation of the southern populations and their resistance toongoing global warming. Splendiani et al. (2016), in a survey among Italian rivers from the Alps to Sicily, showed that the amount of introgression from non-native stocks in Italian trout populations was positively related to climatic variables such as temperature seasonality and annual precipitation. Currently, in some Mediterranean countries the reinforcement with non-native stocks is falling into disuse given the risk posed to native populations. Nevertheless, the presence of “black spots”of already naturalized or highly introgressed populations in many Mediterranean rivers is a serious threat to remnant native populations in neighboring areas within each basin (Araguas et al., 2017). Using a Pyrenean river as a modelof Mediterranean trout habitats and diagnostic mtDNA, nuclear and phenotype markers, the global aim of this study was to assess, for the first time, the factors promoting the introgressive hybridization between native and non-native brown trout at the intrabasin level. This objective includes: (i) determining the local and along basin genetic impact of hatchery releases with non-native brown trout stocks, (ii) identifying environmental and climatic variables promoting the genetic impact of the non-native stocks and (iii) indicating fishery management actions addressed to improve conservation of native trout within already impacted Mediterranean drainages. 2. Material and methods 2.1. Study region and sampling The Mediterranean region is characterized by large inter-annual variability of precipitation and the river flows show large disparities between M. Vera et al. Science of the Total Environment 862 (2023) 160739 2
wet and dry seasons, and between wet and dry years. Rivers in the Eastern Pyrenees exemplify Mediterranean streams. Suitable habitats for salmonids in this regioninclude headstreams at few coastal rivers in Spain and France. However this study analyzed the Segre River, a tributary of the Ebro River showing the main regional water discharge to the Mediterranean Sea because its upper course being partially snow-fed (pluvio-nival hydrological regime, Castelltort et al., 2018) and retains good habitat conditions for salmonids. The Segre River has a basin area of 12,880 km 2 and brown trout is present in the Segre main stem as well as in its tributaries the Noguera Pallaresa and Noguera Ribagorçana rivers, often from 500 m above sea level to the headwaters. Significant genetic differentiation has been detected among trout populations within and between the Noguera Pallaresa and Noguera Ribagorçana which confirms the restricted gene flow among salmonid populations in the Mediterranean rivers (Fernández-Cebrián et al., 2014). Ongoing conservation measures for brown trout in this basin include “genetic refuges”at several headstreams showing low impact by hatchery releases of non-native fish (Araguas et al., 2009, 2017). Additionally, brown trout populations are protected in streams and lakes into the Spanish National Park of “Aigüestortes and Estany de Sant Maurici” but the native purity of these populations into the park is questionable (García-Marin et al., 1998). A preliminary genetic study using the LDH- C1* locus genotypes in 5 localities at the Segre River suggested the existence of native populations (Aparicio et al., 2005). However, compared to the Noguera Pallaresa and Noguera Ribagorçana, the Segre River headstreams are more threatened by human-mediated pressures, with some stretches classified as “bad”or “deficient”according to the biological integrity index IBICAT2a based on fish metrics and habitat quality (Sostoa et al., 2010). From 2016 to 2018, 906 brown trout specimens from 33 localities were sampled along the Segre River basin (Fig. 1,Table 1). Fish were captured using a pulsed-DC backpack electrofisher, anesthetized with clove oil, and afterward placed in a small portable aquarium to be photographed for further classification as native Mediterranean trout, hatchery specimen or hybrid from phenotype traits using the tree model provided by Aparicio et al. (2005). Following photography, a portion of the ventral fin of each captured fish was clipped non-lethally and placed in individual vials with 95 % ethanol for further genetic analyses. To accomplish with legal regulations on wild fauna, all sampled fishes were returned alive to the river at the respective capture locality and sample sizes were limited to 30–35 (but the numbers of some collections were further restricted by the availability of fish at the time of capture). A sample of the Bagà hatchery collected in 2003 (N= 90) was also included as a reference for the non-native hatchery stock used to reinforce the regional brown trout fisheries (Araguas et al., 2017). 2.2. Environmental and climatic variables For each studied locality and using QGIS software (QGIS.org, 2022. QGIS Geographic Information System. QGIS Association. http://www. qgis.org), genetic information on hatchery impact was combined with available information on altitude and hydrology (main-stem and tributaries), management (inside or outside genetic refuge areas, GIS layer available at the Department of Climate Action, Food and Rural Agenda, of the Autonomous Government of Catalonia: http://agricultura.gencat.cat/ca/ serveis/cartografia-sig/bases-cartografiques/cacera-pesca-continental/ reserva-genetica-truites/), water quality and climatic variables. Water Fig. 1. Map of the Segre basin and sampled localities. Locality codes as in Table 1. The blue shaded areas represent genetic reservoirs. M. Vera et al. Science of the Total Environment 862 (2023) 160739 3
quality was incorporated as an indicator of human impact from a QGIS layer detailing the general status of river sections in Catalonia available at the Catalonian Water Agency (Agència Catalana de l'Aigua, ACA) (http:// aca.gencat.cat/ca/laigua/consulta-de-dades/descarrega-cartografica/). This assessment of the water quality combined the structure and functioning of the ecosystem, as well as physicochemical and hydromorphological conditions but excluding microbiological parameters. Thus, a global quality status (from 1, very good, to 4, bad) was recorded for each studied site. Finally, as in Splendiani et al. (2016) the data from 19 bioclimatic variables available in the BioClim database (Karger et al., 2017)werealsoincorporated. 2.3. DNA extraction, molecular genotyping, and phenotype analyses GenomicDNA was obtained from all individuals using the Chelex® protocol (Walsh et al., 1991). Molecular analyses included the genotyping by a single PCR multiplexof five microsatelliteloci (SsHaeIII14.20, Str591INRA, Str73, Ssa85 and SSoSL438) that provide significant divergence among wild brown trout populations in the region and among wild and nonnative hatchery fish (Araguas et al., 2017;Fernández-Cebrián et al., 2014; Sanz et al., 2009). In addition to these five loci, the genotypes at the diagnostic LDH-C1* locus obtained following the PCR Restriction Fragment Length Polymorphism (PCR-RFLP) protocol described by McMeel et al. (2001) were used as evidence for hatchery stocking in the studied area. The sequencing of the complete mtDNA control region (CR) for a subset of samples (4–6 for each site, total = 163) was used to assess the evolutionary relationships of the Segre River brown trout with already studied populations in other Mediterranean rivers and as a complementary measure of maternal genetic introgression from hatchery females. The amplification and sequencing of the complete CR (1013 bp)was carried out following the protocols described in Sanz et al. (2006). Sequences were obtained using the ABI PRISM BigDye™Terminator v3.1 Cycle Sequencing Kit protocol (Applied Biosystems, Foster City, CA) on an ABI PRISM® 3730xl Genetic Analyzer(AppliedBiosystems,FosterCity, CA). Variable sites were checked by hand with the program SEQSCAPE 2.5 (Applied Biosystems, Foster City, CA) usingas a referencethe MEcs1haplotype (GenBank Accession number: AY836350). The different haplotypes were identified using the program MEGA 7.0 (Kumar et al., 2016). Finally, a classification tree model provided by Aparicio et al. (2005) was used to identify native, hatchery and hybrid trout from phenotype traits observed in individual photographs (see Supplementary Fig. S1). 2.4. Genetic diversity and impact of hatchery releases For the five microsatellite loci and the LDH-C1* locus, allele frequencies and genotypic deviations from Hardy–Weinberg (HW) expectations in each locality were estimated using GENEPOP 4.0 (Rousset, 2008) and summarized as F IS coefficient at each one. For each locus at each studied locality, null allele frequency and its 95 % confidence interval (CI) was estimated using the EM algorithm of Dempster et al. (1977) as implemented in GENEPOP software. A lower bound of the 95 % CI >0 was considered as putative evidence for the presence of null allele. Relatedness among sampled individuals within each locality was estimated with the r xy index of Queller and Goodnight (1989), implemented in IDENTIX software (Belkhir et al., 2002). Relatedness estimators become downward biased when calculated from the same genotypes used to estimate the population allele frequencies (Wang, 2017). Thus, the null hypothesis for the random collection of individuals from a panmictic population at each study site was tested by using 5000 permutations of alleles within the respective sample. The contemporary effective population size (Ne) for each site was determined using NeEstimator v. 2.1 software (Do et al., 2014)withthe Linkage Disequilibrium (LD) method under a random mating model and the lowest limit of 0.02 for allele frequency to be used. Impact of hatchery releases into each wild site was estimated by the LDH-C1*90 allele frequency and the proportion of introgressed genome (q) calculated using microsatellite genotypes and following the incomplete baseline method assuming an admixture model with two populations (hatchery and native) as described by Sanz et al. (2009).Despiteusingareduced set of markers, at the population level these two approaches fairly agree with genomic estimates of hatchery impact using larger number of SNPs (Casanova et al., 2022). In addition, at each studied site, an estimate of the impact by released hatchery fish was also computed from the phenotype as the average impact obtained for each sampled fish (0 for native, 0.5 for hybrid and 1 for non-native hatchery fish). Using microsatellite and LDH-C1* genotypes, pairwise population differentiation (F ST ,Weir and Cockerham, 1984) and its significance (genic differentiation using the exact G test) were estimated as implemented in GENEPOP 4.0 with default parameters (10,000 dememorization steps, 100 batches and 10,000 iterations per batch). The group-level Bayesian analysis in BAPS 6.0 (Corander et al., 2008) grouped wild populations that frequently exchanged individuals. BAPS analyses were repeated 20 times, with the maximum number of clusters set to 33, the number of studied sites. In addition, the minimum number of homogeneous units (K) over sampled individuals was estimated using the Markov Chain Monte Carlo (MCMC) method in STRUCTURE 2.3.4. Runs for each possible K (1 to 15, the maximum number of units suggested by BAPS software) were repeated 20 times. Each run used a burn-in of 50,000 iterations, a run length of 200,000 iterations, and the model of correlated allele frequencies. The most likely value of K was selectedfollowing the methodology described by Evanno et al. (2005). Hierarchical analysis of the molecular variance (AMOVA; Excoffier et al., 1992) was used to estimate the percentage of genetic variation distributed within (F SC ) and among (F CT ) groups of populations, where sampled localities were hydrographically grouped (i.e. main stem and tributaries of Segre River in the study region). This test was carried out using ARLEQUIN 3.5 (Excoffier and Lischer, 2010) and the significance of the values was estimated by 10,000 permutations. Table 1 Description of studied localities. Year: sample year; Code: locality code; L: Location (M: main stem; T: tributary —number distinguish among them); Altitude (meters above sea level); GR: sampled site into a Genetic Reservoir; Status: global water quality status (1: very good, 2: good, 3: near good, 4: bad); N: sample size. Year Code Location Altitude GR Status N 2003 BA03 Bagà Hatchery 820 90 2016 QU16 M 1104 Yes 4 30 2017 SA17 M 1052 No 3 30 2016 IS16 M 1027 No 3 30 2016 PR16 M 989 No 3 30 2017 MA17 M 970 No 3 33 2017 PB17 M 823 No 3 5 2017 AR17 M 782 No 3 31 2018 SUR18 M 670 No 3 6 2018 AF18 M 640 No 3 15 2018 TS18 M 573 No 3 25 2018 OR18 M 527 No 3 9 2016 LV16 T01 1141 Yes 4 30 2016 SA16 T02 1145 Yes 4 30 2016 ME16 T03 1700 Yes 2 24 2016 OL16 T03 1134 Yes 2 30 2017 ALP17 T04 1321 No 3 33 2017 LLO17 T05 1425 No 3 33 2018 RD18 T06 1086 No 3 33 2017 VR17 T07 1342 No 3 32 2017 VI17 T08 1009 No 4 31 2018 AR18 T09 792 No 3 34 2018 BE18 T10 1777 No 2 31 2017 FM17 T11 803 No 3 33 2017 JF17 T11 1023 No 3 31 2017 CA17 T12 830 No 4 33 2017 NY17 T13 925 No 4 32 2018 PA18 T13 696 No 2 14 2018 JO18 T14 1137 No 2 31 2018 MO18 T14 1168 No 2 30 2017 LV17 T14 937 No 2 33 2018 PE18 T15 700 No 4 24 2016 R016 T16 1304 No 2 30 2016 RF16 T16 649 No 2 30 M. Vera et al. Science of the Total Environment 862 (2023) 160739 4
Phylogenetic relationships among mtDNA CR haplotypes wereanalyzed by a Median-Joining network (Bandelt et al., 1999) using the software NETWORK 4.5.1.2 (http://www.fluxus-engineering.com/sharenet.htm) and including the haplotypes previously detected in the main tributaries of the Segre River by Cortey et al. (2004): Noguera Pallaresa, Noguera Ribagorçana and Cinca rivers. In addition, haplotypes from neighboring Eastern Pyrenean drainages were added to the network (i.e. localities for the Llobregat River described in Vera et al. (2019); and for the Ter River in Cortey et al. (2004). All these additional haplotypes were available from the GenBank database (GenBank Accession Numbers for MEcs1– MEcs7: AY836350–AY836356, MEcs10: AY836359, MEcs12: AY836361, MEcs15: AY836364, MEcs23: MG970273, MEcs25–MEcs27: MG970274– MG970276, ADcs1: AY836330, ADcs4: AY836333, ADcs6–ADcs8: AY836335–AY836337, ADcs12: AY836341, ADcs13: AY836342, ADcs16: AY836345, ADcs17: AY836346). The haplotypes ATcs1 (AF273086), ATcs2 (AF273087), ATcs3 (AF274574) and ATCs4 (AF274575) reported for hatchery stocks in the Iberian Peninsula (Cortey and García-Marin, 2002) were also included. Network loops (which represent ambiguities) involving three or more haplotypes were solved using the frequency, the topological and the geographical criteria as proposed by Pfenninger and Posada (2002). 2.5. Correlates between hatchery impact and bioclimatic variables Mann-Whitney Utests were used to assess differences in the hatchery impact between sites in the main-stem and those in tributaries, as well as inside and outside genetic refuge areas. Kruskall-Wallis test checked for hatchery impact according tothe water quality status (very good, good, medium, and bad). The above analyses were performed for each of the three estimates of hatchery impact obtained from the LDH-C1* locus, the microsatellite loci and the phenotype. In addition, Spearman's rho correlations (ρ) were estimated between the proportion of the STRUCTURE clusters obtained from combined LDH-C1* and microsatellites genotypes at each site and the bioclimatic variables. Similar correlations were computed between the hatchery impact estimated from the phenotype and the bioclimatic variables. All these analyses were performed using PAST 4.0 Software (Hammer et al., 2001). 3. Results 3.1. Stocking impact and local diversity Genetic impact from non-native stocks was evidenced by the presence of the non-native LDH-C1*90 allele in all sampled localities, with frequencies ranging from 1.6 %at VI17 and JO18to 69.7 % at LLO17 (Table 1).Despite the small sample sizes used for mtDNA analyses, the presence of hatchery haplotypes of the Atlantic lineage corroborated hatchery impact in 13 localities. Estimates of hatchery impact obtained from microsatellite loci and phenotype traits were significantly correlated with the LDH- C1*90 allele frequency: Spearman ρ=0.6955(P<0.001) for microsatellites and ρ=0.7924(P<0.001) for phenotype. The estimates from the microsatellites and the phenotype traits were also highly correlated (ρ= 0.7063, P<0.001). Nevertheless, at ME16 and BE18, microsatellites clearly suggested a much higher impact than the other two markers. Based on the combined information from these three approaches to detect hatchery impact, studied localities were classified intotwo groups. A“low impact”group included those sites where the estimated impact was lower or equal to 14.4 % for the three approaches and a group of “high impact”including localities with at least one of the three impact estimates higher than 14.4 %. Seven out of 11 sites in the main stem of the Segre River were assigned to the “high impact”group, and seven out of 22 in the tributaries were also assigned to this group (Table 2). However, the Fisher exact test did not show significant differences in the abundance of low and high impact sites between the Segre main stem and its tributaries (P= 0.0856), and Mann-Whitney Utests only suggest higher impact from phenotype estimates at the main stem localities (LDH-C1* P = 0.1653, microsatellites P= 0.6435; phenotype P= 0.0048). At the regional scale, the average number of alleles (Na) per microsatellite locus was 13.6, but at local scale, the upper bound of the rank was clearly low: 8.2 at the FM17 locality (Table 2). At the hatchery collection (BA03) the mean value was 5.8. Among wild collections, allele richness (A R ) ranged from 1.697 at MO18 to 4.393 at FM17, these two localities showed also the lower and upper limit for expected heterozygosity H E (0.255 vs 0.728), and while MO18 showed the lower observed proportion of heterozygous (0.250), SUR18 had the largest one (0.833). Mann- Whitney U tests indicated that the localities in the main stem of the Segre River usually show higher diversity values (Na P= 0.0375, A R P= 0.0001, H E P= 0.0007 and observed heterozygosity (H O )P= 0.0003). Similarly, A R ,H E and H O were higher in the “high impact”localities, and only Na was similar between the high and low impact groups (P= 0.2420). Only eight out 33 studied localities showed significant departures from HW genotype expectations, being the PR16 the single one located in the main stem of the Segre River (Table 2). Null alleles were suggested for locus Str591INRA (localities VI17, AR18, and RO16), Str73 (locality AR18), Ssa85 (localities SA17, AR18, and RO16), and SSoSL438 (localities QU16, OL16 and PA18). The suggestion of null alleles at several loci in AR18 could likely reflect mating between relatives as suggested by large positive F IS values, significant relatedness, and low effective size. All Queller &Goodnight relatedness values were negative as expected when calculated from the same sample allele frequencies were estimated (Wang, 2017), but at some localities were higher than expected in an unstructured (random mating) population. In PR16, ALP17, FM17 and AR18, the significant relatedness concurred with overall departures from HW genotype expectations, and in PA18 with departures only for the LDH-C1* locus. Such results suggest recent admixtures of native and hatchery fish or some amount of inbreeding in these trout populations. Estimates on effective population size (Ne) ranged from 2 in RO16, showing also strong departures from HW genotype proportions to “very large”population size in IS16, PB17, LLO17, RD18, CA17, J018, and LV17. Twelve localities had effective population size estimates below 25, all of them but SUR18 located in tributaries. However “very large”sizes were also estimated in some tributary localities. In fact, small and large Ne values were observed among the main stem and the tributaries, and among high and low impacted populations from hatchery releases (Table 2). 3.2. Population structure An overall significant differentiation was observed among studied localities along the Segre River basin either including (F ST = 0.196; P<0.0001) or excluding the LDH-C1* locus and the hatchery sample (F ST = 0.140; P< 0.0001). The estimated overall population differentiation, F ST ,amongwild collections was 0.160(P<0.0001) at the LDH-C1* locus. Based on these observations, further analyses were performed using all loci but excluding the Bagà hatchery collection (BA03). Similar total diversity was observed among populations at the main stem (H T = 0.636) and among populations in the tributaries (H T = 0.657) that almost reach the diversity observed overall studied populations (Table 3). Permutation tests performed with FSTAT (Goudet, 1995)indicatedasignificant (P= 0.001) 10 times higher divergence among populations from the tributaries (F ST =0.183)than among those from the main stem (F ST = 0.016). The minimum number of homogeneous units (K) over sampled individuals estimated according to Evanno's method indicated two STRUCTURE clusters (Delta (K = 2) = 59.352; see Supplementary Fig. S2) named Cluster I and Cluster II, respectively. The proportion of Cluster I was significantly correlated (P<0.0001) with all estimates of hatchery impact (Fig. 2). Three localities (PR16, VR17 and RD18) showed almost identical mean proportions (48–51 %) of membership to both STRUCTURE clusters (Table 4), although most individuals at these three localities showed large (>70 %) inferred ancestry to one of the two clusters, suggesting recent admixturesbetween the two genetic clusters (Supplementary Fig. S3). Four localities from the main stem (OR18, TS18, AF18 and SUR18) and four from M. Vera et al. Science of the Total Environment 862 (2023) 160739 5
tributaries (LLO17, FM17, PE18 and BE18) showed mean proportions of memberships to Cluster I higher than 55 %. These eight samples were all from highly impacted localities by hatchery releases. All the other 22 sites had mean proportions of membership to Cluster I lower than 45 %. This later set included all the localities having low hatchery impact, and likely represent the most native gene pools for brown trout in the basin, corresponding to Cluster II. BAPS results added further structuring by suggesting up to 15 population units among the 33 sampled localities, most of them (11 out of 15) being unique to a single locality (Table 4). The largest BAPS unit grouped eight localities from the main stem and four from tributaries at sites close to main stem (Table 4,Fig. 1). The BAPS unit 4 grouped four highly impacted localities (OR18, FM17, TS18 and NY17), although the most impactedone according to the threemarkers (i.e. LDH-C1*, microsatellite loci and phenotype, Table 2) was identified as a single unit (unit 2, locality LLO17). The low impacted QU16 locality at the upper reaches of the main stem grouped together with SA16, OL16 and ME16 from nearby tributaries into the BAPS unit 12. Surprisingly, at ME16 a high estimate of hatchery impact was suggested by microsatellite loci (38.2 %) that contrasted with lower estimates obtained at the LDH-C1* locus (4.2 %) and from the phenotype (2.1 %). Finally, MO18 and JO18 from the same tributary formed another biological unit according to BAPS results (unit 15). These two sites also had the lowest representation of the STRUCTURE Cluster I. Nine different mtDNA CR haplotypes including a new one (MEcs28, GenBank Accession Number: MT457550) were detected in the Segre River basin (Fig. 3). The new MEcs28 haplotype, observed only in the LV16 locality, differed just only in one mutation (position 853) from the MEcs1 and was then incorporated into the ME lineage. Network analysis clearly distinguished among haplotypes of ME, AD and AT lineages and in addition to the new haplotype, four more haplotypes in the Segre River belonged to the native ME lineage (MEcs1, MEcs3, MEcs23 and MEcs25). The haplotype MEcs23, detected in MO18 and JO18, was recently described from a survey on brown trout diversity at the neighboring Cardener River, a tributary of the Llobregat River. This haplotype was the most abundant among the Cardener River brown trout populations, being the central one for a new haplogroup into the ME lineage (Vera et al., 2019). In addition to ME haplotypes, the ADcs1 haplotype observed in 10 localities represented the native AD lineage, while the non-native haplotypes ATcs2, ATcs3 and ATcs4 of the AT lineage confirmed the success of the hatchery releases in 13 sites. However, AT haplotypes always concurred with native ones of the ME or AD lineages (Table 2). The most frequent haplotype in the Segre basin was the MEcs1(109 out of 163 sequenced individuals)followed by the ADcs1 (21 individuals) and the hatchery ATcs4 haplotype (15 individuals). The Tamura-Nei genetic distance between the detected lineages Table 2 Genetic diversity within studied localities. Code, as in Table 1. Estimated proportion of hatchery impact according to LDH-C1*90 allele frequency (LDH), microsatellite loci (Msats),and Phenotype (Pheno.), resulting in Type L: low impact or H:High impact (seetext). Number of sequenced individuals with mtDNA haplotype of the Atlantic lineage (AT); Mediterranean lineage (ME); and Adriatic lineage (AD). Na: mean number of alleles at microsatellite locus, A R :allelerichness,H E : expected heterozygosity, H O :observed heterozygosity, F IS : Fixation index, r xy : Relatedness among individuals, N E : effective population sizes (vl: very large). (Bold values, P<0.05). Code LDH Msats Pheno. Type AT ME AD Na A R H E H O F IS r xy N E BA03 100 100 100 Hatchery 5 5.8 2.984 0.492 0.477 0.039 −0.007 77.1 QU16 5.0 2.7 6.7 L 1 4 5.6 3.077 0.523 0.528 −0.010 −0.029 40.0 SA17 8.3 6.4 8.3 L 1 3 6.0 3.658 0.629 0.617 0.021 −0.034 278.0 IS16 6.6 5.3 5.2 L 4 1 6.0 3.501 0.587 0.606 −0.032 −0.031 vl PR16 21.7 7.4 12.1 H 5 7.0 3.784 0.658 0.589 0.106 −0.024 62.7 MA17 19.6 7.1 14.1 H 4 1 7.4 3.689 0.630 0.682 −0.084 −0.036 42.5 PB17 20.0 7.2 10.0 H 3 2 4.2 3.833 0.626 0.633 −0.013 −0.271 vl AR17 8.0 5.8 13.3 L 5 7.4 3.747 0.629 0.672 −0.070 −0.029 41.5 SUR18 8.3 10.9 20.0 H 1 4 5.2 4.291 0.674 0.833 −0.266 −0.179 15.3 AF18 16.7 12.4 20.0 H 4 1 6.2 3.775 0.653 0.644 0.014 −0.069 101.2 TS18 18.0 15.4 20.0 H 2 3 6.6 3.976 0.693 0.718 −0.037 −0.051 31.7 OR18 55.6 29.3 33.3 H 1 4 5.2 3.844 0.704 0.722 −0.028 −0.113 31.5 LV16 8.3 3.7 3.6 L 1 4 5.0 2.849 0.499 0.522 −0.048 −0.064 23.9 SA16 10.0 3.1 5.4 L 5 6.4 3.326 0.592 0.609 −0.030 −0.030 30.8 ME16 4.2 38.2 2.1 H 5 5.2 2.944 0.536 0.569 −0.069 −0.027 601.8 OL16 3.3 2.9 0.0 L 5 5.8 3.087 0.526 0.522 0.007 −0.038 25.7 ALP17 19.7 12.8 18.8 H 2 3 6.6 3.409 0.585 0.556 0.052 −0.018 16.5 LLO17 69.7 53.8 46.6 H 4 1 6.0 3.436 0.618 0.631 −0.022 −0.034 vl RD18 12.1 9.6 7.4 L 2 1 1 5.8 3.398 0.591 0.581 0.018 −0.033 vl VR17 14.1 6.8 5.0 L 1 4 6.4 3.369 0.569 0.526 0.076 −0.034 10.4 VI17 1.6 3.6 3.4 L 1 4 5.6 3.266 0.543 0.559 −0.031 −0.026 15.5 AR18 7.4 8.1 4.3 L 6 5.2 3.303 0.574 0.534 0.071 −0.007 12.3 BE18 15.5 31.6 10.9 H 5 4.2 2.794 0.472 0.449 0.049 −0.037 16.2 FM17 25.8 22.2 24.2 H 1 2 2 8.2 4.393 0.728 0.682 0.064 −0.026 36.4 JF17 9.7 2.9 0.0 L 1 4 4.0 2.937 0.546 0.532 0.026 −0.039 43.3 CA17 7.6 6.1 3.1 L 4 5.8 3.367 0.582 0.551 0.053 −0.033 vl NY17 31.3 15.2 15.0 H 1 4 6.0 3.669 0.672 0.630 0.063 −0.031 3.1 PA18 7.1 7.2 6.7 L 6 5.2 3.373 0.564 0.548 0.030 −0.038 14.6 JO18 1.6 4.1 0.0 L 5 2.8 1.787 0.269 0.269 −0.001 −0.026 vl MO18 3.3 4.2 0.0 L 4 2.4 1.697 0.255 0.250 0.019 −0.011 7.2 LV17 3.0 3.5 0.0 L 5 4.4 2.880 0.513 0.455 0.111 −0.036 vl PE18 25.0 15.5 10.0 H 5 4.8 3.222 0.595 0.514 0.140 −0.048 12.6 R016 1.7 1.6 3.4 L 5 5.4 2.870 0.436 0.344 0.214 −0.024 2.0 RF16 3.3 6.3 8.0 L 5 5.6 3.397 0.563 0.628 −0.118 −0.035 92.9 Table 3 Gene diversity analysis on brown trout populations of the Segre River. Localities: number of analyzed localities, Ar: allele richness, H T : total diversity, H S :meanlocal diversity; F ST : population differentiation (in percent of H T ). Localities Ar Gene diversity H T H S F ST (%) All wild localities 33 3.332 0.663 0.570 14.1 Main stem 11 3.743 0.636 0.626 1.6 Tributaries 22 3.126 0.657 0.537 18.3 Main stem Low impact 4 3.496 0.599 0.592 1.1 High Impact 7 3.884 0.665 0.659 1.0 Tributaries Low impact 15 2.994 0.604 0.508 15.9 High impact 7 3.409 0.732 0.605 17.4 M. Vera et al. Science of the Total Environment 862 (2023) 160739 6
was 0.0040 ± 0.0017 for ME vs AD, 0.0083 ± 0.0029 for ME vs AT, and 0.0058 ± 0.0022 for AD vs AT. The distances between the MEcs1 and MEcs23 haplogroups of the ME lineage was 0.0025 ± 0.0011. Applying a divergence rate of 1–2%/My(seeCortey et al., 2004), the divergence times would be ~125,000–250,000 years between these ME haplogroups and ~200,000–400,000 years between AD and ME native lineages. 3.3. Hatchery genetic impact and environmental variables Only five out of the 33 studied localities placed inside a genetic reservoir and estimates of hatchery impact at these localities were often lower than in no refuge ones (Table 2), but these differences were not statistically significant for any marker: LDH-C1* P = 0.1904, microsatellite loci P= 0.0878, and phenotype P= 0.0700. None of the studied localities had very good water quality, only nine out of 33 were in good water quality being all of them placed in tributaries, 17 were in near to good quality and sevenin bad conditions (Table 1).The Kruskal-Wallis test indicatedsignificant differences on the amount of hatchery impact among these three groups either estimated from the LDH-C1*90 allele frequency (P= 0.0013) or the phenotype (P= 0.0033), but not from microsatellite loci (P= 0.1094). The Dunn's post hoc suggested lower impact in localities with good conditions than in the other two categories either from LDH- C1*or phenotype estimates. A positive significant correlation was indicated between the proportion of STRUCTURE Cluster I and several bioclimatic variables related with temperature: annual mean temperature (BIO 1), maximal temperature of warmest month (BIO 5), minimal temperature of coldest month (BIO 6), mean temperature of driest quarter (BIO 9), mean temperature of warmest quarter (BIO10) and mean temperature of coldest quarter (BIO 11) (Table 5). All these bioclimatic variables also showed positive significant correlation with the phenotype estimate of hatchery impact. However, the phenotype estimates of hatchery impact also correlated positively with temperature seasonality (BIO 4) and mean temperature of wettest quarter (BIO 8), but negatively with annual precipitation (BIO 12) and precipitation of wettest quarter (BIO 16). All bioclimatic variables related to temperature, but isothermally (BIO 3), showed negative highly 0 10 20 30 40 50 60 70 80 90 100 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Hatchery Impact (%) Structure Cluster I LDH-C1*90 Msats Phenotype Fig. 2. Relation between the proportion of STRUCTURE cluster I and hatchery impact analyzed from the LDH-C1*90 allele frequency, microsatellite loci (Msats) and Phenotype traits (Phenotype). Table 4 BAPS and STRUCTURE results on genetic relationships among studied localities. Cluster I refer to proportion of assignmentto cluster I from thebest STRUCTURE result of K = 2. BAPS indicated the unit assigned according to BAPS program. Localities are ordered according to BAPS unit and then according to Cluster I proportion. Locality codes as in Table 1. Code BAPS unit Cluster I BE18 1 0.95 LLO17 2 0.91 PE18 3 0.76 OR18 4 0.85 FM17 4 0.71 TS18 4 0.63 NY17 4 0.55 ALP17 5 0.44 VR17 6 0.51 SUR18 7 0.58 AF18 7 0.58 RD18 7 0.49 PR16 7 0.48 MA17 7 0.43 SA17 7 0.39 PB17 7 0.38 CA17 7 0.34 AR17 7 0.33 IS16 7 0.28 AR18 7 0.27 VI17 7 0.22 RF16 8 0.37 PA18 9 0.31 LV16 10 0.29 JF17 11 0.24 SA16 12 0.31 ME16 12 0.21 QU16 12 0.19 OL16 12 0.19 LV17 13 0.13 R016 14 0.13 MO18 15 0.11 JO18 15 0.08 M. Vera et al. Science of the Total Environment 862 (2023) 160739 7
significant correlation with altitude, while all those related to precipitation, but precipitation seasonality (BIO 15), had positive significant correlation with altitude. However, only the phenotype estimate on hatchery impact was significantly correlated with altitude. 4. Discussion In all the analyzed localities, the LDH-C1*90 allele from the non-native hatchery stock was detected, indicating that hatchery impacts are widespread throughout the main stem and tributaries of this basin. A previous survey on the regional impact of hatchery releases observed the *90 allele in four out of the five revised localities of this basin (Aparicio et al., 2005). The locality free of this allele at that time was close to our site OL16, where now the allele is present but in very low frequency. Increased impact of hatchery releases in recent times have been reported atlocallevel in other river basins in this region (Araguas et al., 2017;Vera et al., 2019). Hatchery impact was also evidenced in all revised populations inside the genetic refuges in the Segre River. Lower impacts and putative native brown trout populations has been preserved in the genetic refuges of the Noguera Pallaresa River (Araguas et al., 2017). At local level, the hatchery releases with the non-native stock have increased allele richness and heterozygosis in the highly impacted populations, as predicted by simulations (Fernández-Cebrián et al., 2014). Allele diversity in wild populations is currently higher than in the Bagà hatchery stock due to the admixture between released fish and native brown trout. In addition, the brown trout population structure inside the basin nowadays reflects the genetic impact of the hatchery releases as only two clusters (I and II) were suggested by STRUCTURE software with their abundances into the studied localities being significantly correlated with any of the estimates of hatchery impact (positively with Cluster I). However, according to BAPS, there were up to 15 reproductive units throughout the basin. Several of these units were restricted to a single studied locality and even involved highly impacted populations of small Ne in tributaries (e.g. BE18, PE18). However, other isolated and highly impacted populations in tributaries had very large Ne estimates (e.g. LLO17). Both low and highly impacted populations are contributing to the larger population divergence observed among tributaries (F ST = 0.183) that contrasted with lower divergence among localities in the main stem (F ST = 0.016). In fact, the central part of the main stem (from site SA17 to AF18, see Fig. 1) is populated by a single population unit according to BAPS results. Such pattern of population structure of low divergence in the main stem but larger among tributaries resembles observations in River Koutajoki from Finland, where the main stems was populated by migratory fish while headstreams often consist of isolated resident brown trout Fig. 3. Median-Joining network of detected haplotypes in Segre River drainage and neighbor drainages from Eastern Pyrenees. Black squares represent median vectors needed to connect all the haplotypes. Numbers besides lines indicate the number of mutational steps when more than one mutational step is involved. Dashed lines show alternative connections among haplotypes. Size of the circles is related with the haplotype abundance within Segre drainage. Rectangles indicate the different lineages. Black, light grey, white and dark grey colors in the circles represent the presence of the haplotypes in Segre River, main Segre tributaries (i.e. Cinca, Noguera Pallaresa and Noguera Ribagorzana, see Cortey et al., 2004), Llobregat River (see Vera et al., 2019) and Ter River (see Cortey et al., 2004), respectively. M. Vera et al. Science of the Total Environment 862 (2023) 160739 8
(Lemopoulos et al., 2018). Previous studies in two of the major tributaries of the Segre River, the Noguera Pallaresa and Noguera Ribagorçana rivers, also showed large population divergence among sites in tributaries (F ST > 0.160, Fernández-Cebrián et al., 2014). According to Aparicio et al. (2018), brown trout populations in these rivers are mostly composed (up to 70 %) of a stationary group, limiting dispersal activity to 20–45 m, and a more mobile group with dispersal range 10 times greater (200–500 m). Migration of a few kilometers has been also indicated for specimens of cohorts older than +2 fish in a regional stream (Vera et al., 2010b). However, none of them doing large migrations as observed in longer rivers, where downstream lakes contribute to the retention of migratory behavior in brown trout (reviewed in Ferguson et al., 2019). Experiments at the Swiss Lake Geneva suggest that native Mediterranean trout is more adapted to residency than Atlantic stocks (Caudron et al., 2009). The migratory behavior of brown trout, as in other salmonid fishes, relies on environmental parameters (Dodson et al., 2013;Ferguson et al., 2019) but also on genetic components (Lemopoulos et al., 2018). As hatchery fish originated from the AT lineage, that includes most of current anadromous trout populations in Europe, it is possible that migratory alleles have been introduced on the stocked Mediterranean populations after a long-time evolution favored residency. Substantial genomic divergence has been reported between hatchery trout of Atlantic origin and native Mediterranean populations (Casanova et al., 2022;Leitwein et al., 2016, 2018;Magris et al., 2022). Occasional gene flow between neighboring demes looks essential for preserving local diversity and limiting the effects of the genetic drift and inbreeding expected in native isolated small populations at Mediterraneanstreams (García-Marín et al., 2017). However, larger dispersal potential of released hatchery fish or descendants may promote the spread of detrimental hatchery alleles along basins. Almodóvar et al. (2006) indicated that Mediterranean mtDNA lineages of brown trout were more prone to introgression from non-native Atlantic hatchery stocks than native AT and DU lineages present in the Atlantic Iberian rivers with similar or even higher stocking efforts. As indicate above, all 33 studiedlocalities along the Segre River showedsome degree of hatchery impact, a result contrasting with native populations still present in its tributary the Noguera Pallaresa River (Araguas et al., 2017). The higher abundance of the ME lineage than the AD in the Segre River populations suggests that the ME lineage might be more sensitive to introgression from the Atlantic (AT) hatchery stocks because the AD lineage is the most abundant in the Noguera Pallaresa River (Cortey et al., 2004). In fact, at the neighboring Llobregat River, where only ME haplotypes are present, trout populations have been dramatically impacted in recent times by hatchery AT stocks (Vera et al., 2019). Considering a similar stocking history in these Pyrenean rivers as described for the regional mountains lakes by Miró and Ventura (2013), the use of non-native stocks probably took place sporadically in the first half of the 20th century, but after the Spanish Civil War (1936–1939) an intense program for the construction of hydroelectric infrastructure began in the Pyrenees (Rodriguez, 2012). Releases of non-native trout stocks and even other salmonids were used as a compensatory measure for trout populations even inside protected areas at the Aigüestortes and Estany de Sant Maurici National Park (García-Marin et al., 1998). Surprisingly, the regional hatchery operating at that time was located at the Pont de Suert town at the Noguera Ribagorçana riverside, and this river, as the closer Noguera Pallaresa, still preserve native trout populations nowadays (Aparicio et al., 2005; Araguas et al., 2004, 2017;Sanz Ball-Ilosera et al., 2002;Sanz et al., 2009). Since 1960 to the end of the 20th century a massive program of releases of non-native fish cultured at local hatcheries was used to reinforce regional fisheries accessible by forest roads or even using helicopters (Miró and Ventura, 2013). Stocking practices reduced since year 2000 in parallel with a policy of increasing genetic refuge areas, where releases of hatchery fish are banned, and promoting catch-and-release fisheries (Araguas et al., 2009). Propagule pressure combining stocking effort and distance to the stocked sites is a major factor to predict the amount of hybridization between native westslope cutthroat trout (Oncorhynchus clarki)andintroduced rainbow trout (O. mykiss) in the British Columbia (Bennett et al., 2010). However, stocking effort does not clearly correlate with the level of admixture in river stretches populated by brown trout, as other biotic and abiotic factors are often involved in the survival and reproduction of the released fishes (revised in Ferguson, 2007). In the studied basin, in addition to governmental agencies, stocking has been performed by hydroelectrical companies and angling societies and precise records of these later releases are often incomplete. Altogether, complicate to determine to what extent the observed differences in the local amount of introgression resulted from differences in the magnitude of the releases or are due to other causes. However, the 33 analyzed sites are <100 km apart along the Segre valley and all of them were into or close to fishing areas in the past. Hence, similar historical stocking efforts likely occurred among these localities because proximity to active fishing areas was a major agent to promote the introduction of hatchery stocks during the 20th century (Miró and Ventura, 2013). No significant relationships between hatchery impact and the number of trout farms or trout production were detected in a review of factors involved in the genetic impacts of hatchery releases on the Italian trout populations (Splendiani et al., 2016). This Italian study also failed to detect a relationship between the genetic impact of hatchery releases and land use variables, such as the percentage of forested, agricultural or urban territories; the density of roads or road crossing; the density of human population in the catchment area; and per capita income. To determine the status of surface water bodies, the Catalan Water Agency combines ecological and chemical indexes on quality. Reduced water quality was often the result of anthropogenic pressures, such as those that lead to water eutrophication, or high number of dams, reservoirs, or water diversions for the industry and agriculture; or pollutants among others (ACA, 2005). Using water quality as a proxy of the anthropic pressure on river segments, we detected that those in good quality conditions (all of them placed in tributaries) were less impacted. Thus, habitat recovery to reach good quality in water bodies may prevent the spread of hatchery alleles within Mediterranean basins. Our results indicated that hatchery impacts are higher in localities showing warmer temperature along the year, and estimates based on the phenotype suggest higher hatchery impact on those localities with lower precipitation particularly during the wettest quarter. Previous studies on Iberian trout populations indicated that soft slope, warm temperature in summer, and increased seasonality likely promote a higher impact of hatchery releases in Mediterranean rivers (Madeira et al., 2005). Vera et al. Table 5 Spearman correlations between genetic and environmental and bioclimatic variables. Cluster I: proportion of STRUCTURE cluster I at the studied localities; Phenotype: proportion of hatchery impact estimated at each locality from phenotype. Cluster I Phenotype Altitude Altitude −0.2982 −0.3991⁎ BIO 1 = annual mean temperature 0.3608⁎0.4295⁎−0.9806⁎⁎⁎ BIO 2 = mean diurnal range 0.2668 0.1299 −0.5845⁎⁎⁎ BIO 3 = isothermally (BIO2/BIO7) 0.1550 0.0333 0.2989 BIO 4 = temperature seasonality 0.3203 0.3794⁎−0.9758⁎⁎⁎ BIO 5 = max temperature of warmest month 0.3645⁎0.4253⁎−0.9819⁎⁎⁎ BIO 6 = min temperature of coldest month 0.3589⁎0.4354⁎−0.9785⁎⁎⁎ BIO 7 = temperature annual range (BIO5–BIO6) 0.2998 0.2976 −0.9032⁎⁎⁎ BIO 8 = mean temperature of wettest quarter 0.2703 0.4421⁎⁎ −0.8847⁎⁎⁎ BIO 9 = mean temperature of driest quarter 0.3635⁎0.4285⁎−0.9771⁎⁎⁎ BIO 10 = mean temperature of warmest quarter 0.3623⁎0.4297⁎−0.9819⁎⁎⁎ BIO 11 = mean temperature of coldest quarter 0.3598⁎0.4348⁎−0.9830⁎⁎⁎ BIO 12 = annual precipitation −0.3362 −0.3515⁎0.7658⁎⁎⁎ BIO 13 = precipitation of wettest month −0.2938 −0.3432 0.7436⁎⁎⁎ BIO 14 = precipitation of driest month −0.3363 −0.3218 0.7497⁎⁎⁎ BIO 15 = precipitation seasonality 0.1456 0.0105 −0.3126 BIO 16 = precipitation of wettest quarter −0.3352 −0.3561⁎0.7184⁎⁎⁎ BIO 17 = precipitation of driest quarter −0.3393 −0.3212 0.7907⁎⁎⁎ BIO 18 = precipitation of warmest quarter −0.3420 −0.3399 0.7831⁎⁎⁎ BIO 19 = precipitation of coldest quarter −0.3322 −0.3074 0.8075⁎⁎⁎ ⁎P<0.05. ⁎⁎ P<0.01. ⁎⁎⁎ P<0.001. M. Vera et al. Science of the Total Environment 862 (2023) 160739 9