Microgeographic genomic variation and connectivity in an endangered semiaquatic mammal

doi:10.21203/rs.3.rs-7018148/v1

Download PDF

Research Article

Microgeographic genomic variation and connectivity in an endangered semiaquatic mammal

https://doi.org/10.21203/rs.3.rs-7018148/v1

This work is licensed under a CC BY 4.0 License

Journal Publication

published 23 Oct, 2025

Read the published version in BMC Ecology and Evolution →

You are reading this latest preprint version

Background

The Iberian desman (Galemys pyrenaicus), a semiaquatic mammal endemic to southwestern Europe, is listed as Endangered by the IUCN due to substantial range and population decline. Its restriction to upstream locations highlights the importance of understanding its genetic variation and connectivity for effective conservation strategies. While previous studies have revealed phylogeographic structure across the range of the Iberian desman, gaps remain in our understanding of the microgeographic dynamics that shape genetic exchange within specific geographic regions.

Results

This study first combined newly generated SNP data with previously available datasets to further explore genetic structure in the Iberian desman across its entire distribution, using a set of 110 SNPs on 115 individuals. This confirmed the presence of five major phylogeographic units. Focusing on the newly generated data, we explored the microgeographic dynamics of the Occidental unit with a higher-resolution genomic dataset (7,604 SNPs, 14 individuals). This analysis provided evidence of isolation-by-distance (IBD), indicating that gene flow decreases with increasing geographic distance and that dispersal occurs primarily over short distances. Focussing on the Douro river system, our genomic clustering results showed both connectivity along the best-sampled river and between headwaters of this river and headwaters from a closely located watershed. Our IBD results were consistent with this: indicating riverine dispersal as well as a combination of riverine and overland dispersal at headwaters. These results highlight the importance of both aquatic and terrestrial corridors at upstream locations for maintaining connectivity.

Conclusion

Our findings emphasize the critical role of headwater regions in supporting gene flow and preserving genetic diversity in the Iberian desman. Conservation efforts should prioritize the protection and restoration of riparian and terrestrial corridors, particularly in fragmented landscapes, to mitigate isolation and preserve genetic diversity in the desman. This study underscores the value of genomic approaches in conservation and contributes to a deeper understanding of the ecological and evolutionary processes that maintain population connectivity in an endangered species.

population genomics

Galemys pyrenaicus

connectivity

conservation genomics

Understanding microgeographic population structure is essential for conservation biology, particularly in species inhabiting fragmented or isolated habitats (e.g. [1–3]). In freshwater ecosystems, dams and other barriers, habitat degradation and climate change can lead to genetic isolation and increase the risk of inbreeding in aquatic species by reducing gene flow and connectivity [4–6]. Similar processes can also affect semiaquatic mammals, particularly by restricting populations to headwater streams [7, 8]. However, because semiaquatic species depend on both terrestrial and aquatic habitats, overland dispersal among proximate headwaters [9] can potentially promote gene flow between populations. While semiaquatic amphibians and reptiles have been relatively well studied in this context (e.g. [7, 10]), much less is known about how semiaquatic mammals respond to habitat fragmentation, applying also to our study species, the Iberian desman [8].

The Iberian desman, Galemys pyrenaicus, is the only species of its genus and sister species to another semiaquatic species, the Russian desman, Desmana moschata. Together these two species form the tribe Desmanini within the family Talpidae (the moles and related taxa; [11, 12]). The Iberian desman is endemic to southwestern Europe, where it is restricted to the Pyrenees, northern and central Spain, and northern Portugal, and is listed as Endangered by the IUCN [11]. Over recent decades, the Iberian desman has experienced a severe decline in its range, retreating into cooler, mountainous headwater streams [8, 13]. The species is considered one of the rarest, yet least-studied European mammals [14]. The limited amount of research on such a rare species undergoing a steep decline over all its range raises concerns about undetected population isolation and inbreeding, potentially putting populations at risk.

In relation to the genetic studies that have been carried out on the Iberian desman, the first detailed phylogeographic study on the species detected low levels of genetic diversity and identified four distinct mitochondrial lineages [15]. Similar geographic structuring was found with SNP data from double digest Restriction-site Associated DNA sequencing (ddRAD), describing five distinct phylogeographic units (Occidental, Central System, Iberian Range, Pyrenees, Cantabria) [16]. These units show limited admixture between them and are largely coincident with the main mountain ranges where the Iberian desman occurs, with the exception of the Occidental phylogeographic unit (with a large range in northern Portugal and Galicia). While these broader patterns reflect major phylogeographic entities separated by low connectivity, local demography and dispersal patterns need to be defined within phylogeographic units to further understand the population structure of desmans. Recent studies suggest that the decline in Iberian desman populations may lead to high inbreeding levels and local extinctions in the Pyrenees [17] and in the Iberian Range [18]. ddRAD and relatedness estimates based on tissue samples of individuals from nearby rivers of the Iberian Range have also shown that there is lower dispersal between rivers compared to within rivers [18]. Clearly the degree of dispersal between nearby rivers is important for connectivity within geographic regions [19].

However, the existing genomic studies on the Iberian desman have provided limited information on the Occidental and Central System phylogeographic units. Populations from the Occidental unit exhibit the highest levels of heterozygosity [16, 20] and nucleotide diversity [15], whereas those from the Pyrenees and Central System show the lowest levels of heterozygosity [16, 17, 20] and nucleotide diversity [15]. This likely suggests that the Occidental unit may have served as a major glacial refugium for desman populations during the Last Glacial Maximum [15], being a key reservoir of genetic diversity for this species, underscoring its conservation value. This is particularly important considering that the Iberian desman is among the most inbred and genetically least variable mammals assessed to date [17, 20]. Maintaining the genetic variation within the Occidental phylogeographic unit may enhance the species' adaptability to changing environments, providing resilience to various threats. As such, the Occidental phylogeographic unit is not only important for conservation in its own right, but also a potential source for genetic rescue for other phylogeographic units [21].

Previous research, using the cytb mitochondrial marker, provided some insights into the genetic structure and dispersal patterns of the Iberian desman in the Occidental phylogeographic unit, especially across the Douro and Minho river systems [19]. The research primarily examined large-scale patterns within the phylogeographic unit. Building on this foundation, the research that we report here addresses instead the microgeographic genomic variation in the Occidental unit, focussing on populations of Iberian desman from two adjacent watersheds in Portugal (part of the Douro river system), where they occur in the upper reaches. These watersheds allowed us to test two key aspects of desman dispersal and connectivity at a fine scale: 1) Linear/river connectivity within watersheds – evaluating movement and gene flow along the river, and 2) Dispersal across headwaters – assessing movement and gene flow between adjacent watersheds.

To do this, ddRAD sequencing was applied to non-destructive tissue samples collected at the Baixo Sabor Long Term Ecological Research Site (LTER Network – http://lternet.edu) and nearby areas. The watersheds involved were the Sabor and Tua (N 41°09′–42°06′, W 7°37′–6°16′) (see Quaglietta et al. [8] for details). As a genomic methodology, ddRAD produces a substantial number of SNP markers, which permits detailed analysis of genetic diversity and population structure and has been used for other conservation genetic studies (e.g. [22]). This genomic approach uses two restriction enzymes enabling paired-end sequencing at specific loci across many samples, offering greater precision in read mapping compared to genotyping by sequencing (GBS) or single-enzyme RAD sequencing (RADseq) [23, 24].

We believe that new data on the Occidental phylogeographic unit may aid efforts to maintain the genetic diversity of the Iberian desman. Our study provides a high-resolution analysis of local dispersal patterns, offering insights into the mechanisms shaping desman genetic connectivity in fragmented riverine systems. In particular, understanding the dispersal and connectivity between watersheds may inform efforts to preserve genetic diversity and gene flow. Protecting and enhancing connectivity within the genetically diverse Occidental phylogeographic unit may help ensure the Iberian desman remains genetically viable and adaptable in the face of threats such as climate change and habitat fragmentation.

2.1. Sample collection

All Iberian desman specimens used in the present study were collected in previous conservation-related studies. These included non-destructive tissue samples collected between 2015 and 2016 in Bragança, Northern Portugal in the scope of an ecological study on the species by the LTER Network ([8]), in the Douro river system. Additionally, tissue samples were taken from animals that inadvertently died during live capture, from the Minho river system in Galicia, Spain, collected in 2012, and two from the Ulla River system in Galicia, Spain, collected in 2013 and provided by ARCEA. For the collection of tissues, sterilized tweezers were used, and samples were kept in 96% ethanol until DNA extraction. Individual desmans providing samples were georeferenced in the field using a GPS receiver. A total of 15 tissue samples from Portugal and Spain were used to extract DNA (Additional File 1 - Table S1).

2.2. DNA extraction, PCR amplification and sequencing

For each individual, DNA was extracted using the Qiagen DNeasy Blood & Tissue extraction kit, following the manufacturer’s instructions.

DNA quantification was assessed using a Qubit 2.0 Fluorometer. DNA libraries were constructed following the protocol of Peterson et al. [24] with modifications. Genomic DNA (~ 150 ng) was digested with three units of SbfI-HF and six units of MspI, while simultaneously ligating P1/SbfI (sample barcode) and P2/MspI adapters with 120 cohesive end units of T4 DNA ligase. Reaction volumes were 30 µl in 1× CutSmart buffer and 1 mM ATP. All primers were from Integrated DNA Technologies; all other reagents were from New England Biolabs. Following digestion/ligation, 1.5 µl of each sample was amplified in a 10 µl volume with 2.5 pmol Illumina Truseq P1 (AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGAT) and index (CAAGCAGAAGACGGCATACGAGATxxxxxxGTGACTGGAGTTCAGACGTGTGC, where xxxxxx is the index group sequence) primers, 0.2 mM dNTPs, and 1.25 units One Taq DNA polymerase. PCR reactions were cycled at 95°C for 40 s, 60°C for 45 s, and 68°C for 30 s, for 27 cycles. Sample aliquots were pooled and small fragments were removed with Ampure XP (0.7×). The library was diluted to 2nM and sequenced on an Illumina HiSeq 2500 platform (paired-end 155 bp reads).

2.3. Data analysis

The repetitive regions in the reference genome available in NCBI (Gpyr_1.0 - GCA_019455555.1 [20]) were predicted and masked by RepeatMasker v4.1.4 (http://www.repeatmasker.org) using the Dfam Consensus release 20220412 database [25].

2.3.1. Species-unit assignment

We first verified whether the 15 tissue samples from the western part of the Iberian Peninsula clustered as expected within the previously-identified Occidental phylogeographic unit [15, 16, 19, 20]. For that, these samples were analysed together with previous existing data from a total of 100 tissue samples of the Iberian desman (94 from ddRAD-derived SNP data and 6 from whole genomes) covering the majority of the species’ distribution range (Fig. 1, Additional File 1 - Table S1). A quality distribution plot of the read data for the 115 samples was generated using FASTQC [26]. The combined dataset was trimmed and adapters removed with AdapterRemoval v2.3.2 [27] and then aligned to the reference genome (GCA_019455555.1 - [20]) using BWA-MEM with default parameters [28]. For the whole-genome data, Picard v3.3.0 with the MarkDuplicates command (https://github.com/broadinstitute/picard) was used to identify and remove duplicate reads. To build and genotype the paired-end data and call SNPs, bcftools mpileup and bcftools call (with the --mv command) was used. After obtaining the vcf file, bfctools was used to filter with a minimum read depth of 500 and a minimum quality of 30 for each sample (with the filter -i 'QUAL > 30 || DP > 500' command). Then we kept the SNPs that were present in our samples to minimize missing data using the --positions command from vcftools. After that, vcftools was again used to filter away samples with more than 80% of missing data (--max-missing 0.8 command), kept only biallelic SNPs (--max-alleles 2 command), and retained only SNPs where the minor allele was present at a frequency of at least 5% (--maf 0.05 command). Finally, a Principal Components Analysis (PCA) was made using the --pca command from PLINK v1.9 [29]. To check phylogenetic relationships between units, the dataset was bootstrapped 100 times to ensure statistical power using a custom script (https://github.com/pdroslva84/Plink_IBS_bootstraps) and a neighbour-joining tree was constructed using PHYLIP v3.6 [30] with 100 bootstraps, using the Russian desman (Desmana moschata) as outgroup (GCA_965120235 - [31]).

2.3.2. Microgeographic genomic analysis in the Occidental phylogeographic unit

After confirming that the 15 samples clustered within the Occidental phylogeographic unit and having supported the existence of five main units in the Iberian Peninsula, we conducted a detailed analysis of our new samples. A quality distribution plot of the read data of the 15 samples was generated using FASTQC [26]. Demultiplexing by barcode was carried out using the process_radtags utility included in Stacks 2.62 package [32]. The adapters were removed with AdapterRemoval v2.3.2 [27] and then aligned to the reference genome (Gpyr_1.0 - [20]) using BWA-MEM with default parameters [28]. To build and genotype the paired-end data and call SNPs, bcftools mpileup and bcftools call (with the --mv command) was used. After obtaining the vcf file, bfctools was used to filter with a minimum read depth of 500 and a minimum quality of 30 for each sample (with the filter -i 'QUAL > 30 || DP > 500' command) and vcftools was used to remove samples with more than 90% missing data (--max-missing 0.9 command), retaining only biallelic SNPs (--max-alleles 2 command) with a minor allele frequency (MAF) of at least 5% (--maf 0.05 command). While 5% MAF filter may seem high for 15 individuals, this threshold ensures that variants are present in at least two heterozygous individuals (or one homozygote), avoiding single-individual allele occurrences and reducing noise from rare or spurious variants. Before proceeding with the analysis, we first checked if there were any related samples using the --relatedness2 command from vcftools. A Principal Components Analysis (PCA) was performed using the --pca command from PLINK v1.9 [29] and a Discriminant Analysis of Principal Components (DAPC) was made using the adegenet package from R [33]. Population structure was estimated with the program STRUCTURE 2.3.4, which implements a Bayesian model-based clustering method [34], with no prior information on population origin. A total of 100,000 generations were run after a burn-in of 10,000 generations with a number of clusters (K) ranging from 1 to 8. We implemented 3 different runs for each K value in order to evaluate the variance of the estimated posterior probability of the data, Ln P(D) [34]. In addition, the optimal K value was estimated using the ΔK method [35] using CLUMPAK [36]. Genetic distances between individuals were calculated using the --distance function of PLINK v1.9 and the --1-ibs and --square modifiers. A PERMANOVA (Permutational Multivariate Analysis of Variance; [37]) using the adonis2 from the vegan package from R [38] was used to test for significant variation in genetic distance comparisons between and within river systems. To test for isolation-by-distance (IBD), the correlation of genetic distances and Euclidean (overland) distances, Mantel tests [39] were performed in R (using the vegan package with the mantel command, applying the "spearman" method and 9999 permutations). For a subset of our data we also conducted Mantel tests using river distances rather than overland distances, making use of a river network geodatabase for Europe [40] and the sf [41] and sfnetwork [42] R packages. To account for movement between watersheds, we also carried out a third category of Mantel test, combining river distances (along rivers) and overland distances (between rivers).

An overview of the workflow used in this study is illustrated in Additional File 2 - Fig. S1. The diagram depicts the sequence of statistical and computational tests conducted to examine population structure and connectivity in the Iberian desman.

3.1 Phylogeographic subdivision in the Iberian desman

For the complete dataset of 115 Iberian desmans for which genomic data were available, a total of 110 SNPs were obtained. The PCA based on these SNPs shows five major phylogeographic units (Occidental, Central System, Cantabria, Iberian Range, Pyrenees) as described in Querejeta et al. [16] (Fig. 1). Our new samples clustered within the Occidental phylogeographic unit.

The neighbour-joining tree (Additional File 2 - Fig. S2) further supports these five phylogeographic units, emphasizing divergence of the Occidental unit from the others (although with low bootstrap support: 58). The Pyrenees unit also forms a relatively distinct cluster (bootstrap support: 78).

3.2 Microgeographic genomic analysis in the Occidental phylogeographic unit

One of our samples (from Riobó) was removed (sample 15), because of its close relatedness (r = 0.58) to another. For the 14 remaining samples, we obtained 7,604 SNPs. The genetic distances between individuals ranges from 0.2247 to 0.3196 (Additional File 1 - Table S2) and the overland distances between them ranges from 5 to 207 km (Additional File 1 - Table S3). PCA, DAPC and STRUCTURE analyses showed clear geographic clustering of the samples (Fig. 2). PCA and DAPC show similar results, with DAPC showing the best K as 4 (Additional File 1 - Table S4; Additional File 2 - Fig. S3 and S4), but also separating the northernmost samples (12, 13 and 14) along the first axis of variation (Fig. 2a and b). The same pattern was observed in STRUCTURE, despite the best K being 3 using the ΔK method (Evanno et al. 2005) (Additional File 2 - Fig. S5). PERMONOVA demonstrated that the genetic differences between the northernmost samples and the remaining ones were significant (Additional File 1 - Table S5a). The Mantel test applied to all 14 samples showed significant isolation-by-distance (IBD) with a strong correlation between genetic and overland distances (p = 0.0001; r = 0.79 - Additional File 1 - Table S6).

A hierarchical approach was used to examine genetic structure in the Douro river system, considering 11 samples from two watersheds (Tua and Sabor) and excluding the 3 Galician samples (12, 13 and 14) which are from different river systems (Minho and Ulla; Fig. 2a). A total of 8,275 SNPs were available. The genetic distances between individuals ranges from 0.2122 to 0.2744 (Additional File 1 - Table S2) and the overland and river distances between them ranges from 5 to 79 km and from 21 to 309 km, respectively (Additional File 1 - Table S3). Both STRUCTURE (Fig. 3a and b) and PCA (Fig. 3c) analyses showed geographic clustering of individuals, with STRUCTURE showing the best K as 2 using the ΔK method (Evanno et al. 2005) (Additional File 2 - Fig. S6). Moreover, there is genetic differentiation by watershed, as shown in the PERMANOVA analysis, which demonstrated significant differences between the samples from each watershed (Additional File 1 - Table S5b). Nevertheless, there are also indications of genetic exchange between the headwaters (Fig. 3a). The Mantel test for the 11 samples did not reveal IBD (p-value > 0.01) when using overland distances (p = 0.02; r = 0.34 – see Additional File 1 - Table S6), though it revealed IBD when using river distances (p = 0.003; r = 0.42). Considering the possibility of terrestrial dispersal between headwaters, another Mantel test was made using a mix of river distances for most samples but overland distances for the pairs of individuals 1–3, 9 − 7 and 10–11 where dispersal between headwaters was possible (Additional File 1 - Table S3c). This test was significant (p-value < 0.01), with the expected positive correlation for IBD (p = 0.0001; r = 0.55 - see Additional File 1 - Table S6). We also found IBD when restricting the Mantel test to river distances along the Sabor (7 samples) (p = 0.0001, r = 0.70) but not when using overland distances (p = 0.06, r = 0.45) (Additional File 1 - Table S6).

4.1 Phylogeographic subdivision in the Iberian desman

The newly generated SNP data, combined with compatible published data, distinguishes the same five major phylogeographic units of the Iberian desman as described in Querejeta et al. [16] (Fig. 1). This broad phylogeographic pattern was detected using 110 SNPs, demonstrating the utility of datasets with small numbers of genetic markers, a limitation that may often be present in low budget conservation genomic studies and when there is low genetic diversity in rare and limited range species. Reduced SNP datasets have been successfully utilised in other population studies on various taxa (e.g. [22, 43–45]).

The small SNP dataset also generated interesting phylogenetic results, although with relatively low bootstrap support (Additional File 2 - Fig. S2). Unlike previous phylogenetic analyses of the Iberian desman using mitochondrial data [15, 20], which cluster the Pyrenees, Cantabrian, and northwest Iberian Range units into one group and the southeast Iberian Range, Central System, and Occidental units into another, our SNP-based analysis does not show the Iberian Range unit split into two clades. The neighbour-joining tree separates the Occidental unit from the remaining phylogeographic units, although with bootstrap support of only 58. This separation may reflect historical isolation of the Occidental phylogeographic unit, perhaps due to its location in a glacial refugium during the Pleistocene. The Pyrenees phylogeographic unit forms a distinct cluster with stronger bootstrap support (78). Igea et al. [15] posited that this population likely originated from a distant glacial refugium in the Basque Mountains. They suggested that the population in the Pyrenees colonized this region relatively recently, following a severe bottleneck. Escoda and Castresana [20] further emphasized extremely low genomic diversity and high levels of inbreeding in Pyrenean populations, likely due to repeated bottlenecks during postglacial recolonization. While low bootstrap values limit definitive conclusions regarding finer-scale relationships, our neighbour-joining tree supports the overall phylogeographic units established by Querejeta et al. [16]. Further analysis would be desirable to determine the basis of the differences between the findings with mitochondrial and nuclear data.

4.2 Microgeographic genomic analysis in the Occidental phylogeographic unit

The 14 desman samples analysed from the Occidental phylogeographic unit (using 7,604 SNPs) showed genetic separation of the three from northernmost Galicia (12, 13 and 14) (Fig. 2; Additional File 1 - Table S5a). This likely reflects their origin from different river systems than the Douro, and therefore from a different evolutionary unit.

The Mantel test for the 14 samples that we analysed showed a significant positive correlation between genetic and overland distance, thus suggesting IBD (Additional File 1 - Table S6). This indicates that gene flow decreases with increasing geographic distance and that dispersal occurs primarily over short distances, a pattern typical of species with constrained movement capabilities (as would be expected for a small mammal). These findings at a larger scale align with earlier studies [19, 46–48], which revealed limited dispersal in the Iberian desman. Similar IBD patterns have been observed in studies of other taxa [49, 50], including other semiaquatic mammals such as otters [51], where restricted dispersal and habitat specialization contribute to localized genetic differentiation. These findings highlight the importance of local connectivity for maintaining genetic exchange within fragmented habitats.

The IBD observed for the 14 samples may substantially be the consequence of the distinctiveness of the three Galician samples, distantly located from other samples. The degree of genetic structuring within the Douro river system of the 11 samples from the Tua and Sabor watersheds (8,275 SNPs) was therefore of particular interest. The desmans from the Tua and Sabor watersheds are genetically distinct from each other (Additional File 1 - Table S5b) but our data also indicate dispersal between the watersheds across adjacent headwaters, as exemplified by the significant Mantel test using mix geographic distances (of river distances and overland distances) at the headwaters (Additional File 1 - Table S6). Similar reliance on terrestrial pathways to cross between fragmented aquatic habitats has been observed in other organisms [9]. These findings partially align and contrast with previous mitochondrial data from Querejeta et al. [19]. They analysed 157 samples (44 from Minho, 55 from Douro and 58 from other river systems) using cytb and D-loop genes (a total of 1,066 bp) and found a strong IBD signal at the regional scale (i.e., their whole study area). Similarly, we also observed strong IBD between genetic and overland distances in our analysis of 14 samples (p = 0.0001; r = 0.79). For the Douro river system, the same authors reported a weak IBD. Our findings showed a significant and strong correlation between genetic and river distances (p = 0.003; r = 0.42) and an even stronger correlation between genetic and mix (overland + river) distances (p = 0.0001; r = 0.55). This suggests that even though rivers are important for movement and gene flow, especially for semiaquatic species like the Iberian desman, both overland and riverine routes might play important roles in shaping its genetic connectivity. The overall non-significant correlation for the Douro river system from Querejeta et al. [19] could be explained by the fact that they did not capture the effects of microgeographic dispersal dynamics and instead their results reflected large-scale differentiation within this river system (i.e., without accounting for closely located headwaters).

While limited inter-river dispersal has been previously observed [18], our study provides evidence of such dispersal between the upper reaches of closely situated rivers, being the first report of this kind within the Occidental phylogeographic unit. Our findings highlight the importance of headwaters in maintaining genetic connectivity, counteracting the isolation observed downstream where larger river systems often show suboptimal conditions for dispersal (as already mentioned by Escoda et al. [18]). These headwater habitats not only maintain genetic exchange between watersheds but also provide critical refuges that buffer against environmental change, in particular in current climate change scenarios [52]. Other semiaquatic vertebrates show similar results, with the upper reaches of rivers acting as genetic reservoirs supporting metapopulations across fragmented landscapes [7].

Our microgeographic genomic study not only enhances our general understanding of habitat fragmentation and dispersal in the Iberian desman, it also provides new genetic insights into the Occidental phylogeographic unit, highlighting the role of dispersal between watersheds in maintaining genetic diversity. In particular, given its high genetic diversity compared to the other phylogeographic units [15, 16, 20], the Occidental unit could serve as an important source population for genetic rescue efforts aimed at mitigating the risks of inbreeding and genetic drift in vulnerable populations in other parts of the distribution of the Iberian desman. Examples in other species where relatively genetically diverse populations have acted as a source in genetic rescue of threatened populations include the mountain pygmy possum [53], the Allegheny woodrat [54] and the brook trout [55].

Our study highlights the genetic subdivision of the Iberian desman, confirming the existence of five major phylogeographic units, supporting the hypothesis of post-glacial isolation and then expansion although with limited recent connectivity among units. Additionally, the use of a small SNP dataset to uncover broader genetic structure demonstrates a cost-effective approach for studying endangered species in resource-limited settings.

Our genomic analysis revealed isolation-by-distance in the Occidental phylogenetic unit, as would be expected for a small mammal with limited dispersal. Furthermore, our spatial microgeographic analysis in the Douro river system showed riverine connectivity but also overland dispersal between adjacent headwaters, especially across short terrestrial gaps. This latter finding corroborates that headwater connectivity, both aquatic and terrestrial, plays a critical role in facilitating gene flow within fragmented riverine systems.

From a conservation perspective, these findings reinforce the importance of preserving both riparian and terrestrial corridors between adjacent headwaters. In particular, for the Iberian desman Occidental phylogenetic unit headwaters serve as critical pathways for gene flow, mitigating the isolation imposed by downstream habitat fragmentation. Protecting these upland areas is essential to preserving the genetic diversity of the species, which provides the basis for its adaptability and long-term survival in a rapidly changing environment. High genetic diversity may be important not only to conserve the Occidental phylogenetic unit itself, but to provide opportunities for genetic rescue of vulnerable populations in other phylogenetic units of the desman.

Future research should expand on our findings by including larger sample sizes, additional watersheds, and higher-resolution genomic datasets to better document the genetic diversity of the Iberian desman over the whole Occidental phylogenetic unit and, in particular, to understand the interplay between aquatic and terrestrial connectivity within this system. Likewise, telemetry studies would provide further complementary knowledge on fine-scale dispersal routes and behaviour within and between rivers. Such efforts will be crucial for developing targeted conservation strategies that enhance connectivity across fragmented landscapes, ensuring the Iberian desman’s persistence and the maintenance of its genetic health.

Ethics approval

There were no new collections made for the current study, and therefore no ethical approval needed to carry out our work. We made use of previously collected material. All specimens were obtained in earlier ecological and conservation-related studies requiring permitting (incorporating ethical considerations) as described below.

The Iberian desman individuals from the Douro river system were captured and small tissue biopsies taken (from the tail tip) under collection permits No. 663/2015/CAPT and No. 548/2016/CAPT issued by Instituto da Conservação da Natureza e das Florestas (ICNF), the Portuguese national authority regulating nature.

The two samples from the Minho river system correspond to tissue samples of two animals that inadvertently died on capture in the Serra do Courel with a licence no. 090/2012, requested as part of the ‘Scientific Mammalogy Work Camp in Courel - Galicia’, and issued by the Xunta de Galicia, Council of Environment, Territory and Infrastructures, the general directorate of nature conservation.

The two samples from the Ulla river system correspond to tissue samples of two animals that inadvertently died on capture collected by ARCEA Xestión de Recursos Naturais S.L. with authorisation from the Xunta de Galicia (the general directorate of natural heritage), permit No. 544/2012, within the framework of the 'Service to study the use of the Iberian desman in the Ulla river basin to determine its spatial use dynamics,' as part of the LIFE+Margall Ulla project (LIFE NAT/ES/000514), co-financed with LIFE+ NATURE & BIODIVERSITY funds, at 49.39% (No. 27/2012).

Consent for publication

Not applicable.

Availability of data and materials

The ddRAD data generated in the present study has been deposited in the European Nucleotide Archive (https://www.ebi.ac.uk/ena/browser/view/PRJEB88519). Details of all data generated or analysed during this study are included in this published article and its supplementary information files.

Competing interests

The authors declare no competing interests.

Funding

Work supported by National Funds through FCT-Fundação para a Ciência e a Tecnologia in the scope of the project UID/50027-Rede de Investigação em Biodiversidade e Biologia Evolutiva. This study was also funded by FCT - Fundação para a Ciência e a Tecnologia through the grant 2021.07072.BD (https://doi.org/10.54499/2021.07072.BD). Additional support was provided by FLAD – Luso-American Development Foundation. The collection of samples from Portugal was supported by Energias de Portugal (EDP) Biodiversity Chair and was part of the Long-Term Ecological Research (LTER) project (grant LTER/BIA-BEC/0004/2009). The collection of samples from the Ulla river system was funded by the LIFE+Margall Ulla project (LIFE NAT/ES/000514) and co-funded with LIFE+ NATURE & BIODIVERSITY (No. 27/2012).

Authors’ contributions

SS, JBS, PCA, JP, SB¹, JAA, LQ and SB² conceived and designed the study. LQ and PB provided the samples. SS, SB² and MH carried out the DNA extraction and library preparation. SS analysed the data with support from SB¹ and JAA. SS and JBS drafted the manuscript. All authors reviewed and edited the manuscript. All authors read and approved the final manuscript.

Note: SB¹ refers to Soraia Barbosa, and SB² refers to Steven Bogdanowicz.

Acknowledgements

We gratefully acknowledge the samples from Riobó, Galicia, provided by ARCEA and to Xunta de Galicia for authorising us to use these samples. We also thank Diogo Coutinho-Lima for discussions on genomic methodologies.

Author information

Authors and Affiliations

CIBIO/InBIO, Research Center in Biodiversity and Genetic Resources, University of Porto, Vairão, Portugal

Sara Sampaio, Soraia Barbosa, Paulo C Alves & Jeremy B Searle

Department of Biology, Faculty of Sciences of University of Porto, Porto, Portugal

Sara Sampaio, Paulo C Alves & Jeremy B Searle

BIOPOLIS Program in Genomics, Biodiversity and Land Planning, CIBIO, Campus de Vairão, Vairão, Portugal

Sara Sampaio, Soraia Barbosa, Lorenzo Quaglietta, Paulo C Alves & Jeremy B Searle

Department of Ecology and Evolutionary Biology, Cornell University, Ithaca, NY, USA

José A. Andrés, Steven Bogdanowicz & Jeremy B Searle

CIBIO/InBIO, Research Center in Biodiversity and Genetic Resources, Institute of Agronomy, University of Lisbon, Lisbon, Portugal

Lorenzo Quaglietta

AEPGA – Associação para o Estudo e Protecção do Gado Asinino – Largo da Igreja, Atenor (Vimioso), Portugal

Lorenzo Quaglietta

Fluvial and Terrestrial Ecology Laboratory, University of Trás-os-Montes and Alto Douro, Vila Real, Portugal

Paulo Barros

Laboratory of Molecular Ecology, Institute of Animal Physiology and Genetics, Czech Academy of Sciences, Liběchov, Czechia

Michaela Horníková

EMBL-EBI, European Molecular Biology Laboratory, European Bioinformatics Institute, Wellcome Genome Campus, Hinxton, Cambridge, UK

Joana Paupério

EBM, Biological Station of Mértola, Mértola, Portugal

Paulo C Alves

Corresponding author

Correspondence to Sara Sampaio.

Mikles CS, Aguillon SM, Chan YL, Arcese P, Benham PM, Lovette IJ, et al. Genomic differentiation and local adaptation on a microgeographic scale in a resident songbird. Mol Ecol. 2020;29:4295-307.
Gallagher BK, Fraser DJ. Microgeographic variation in demography and thermal regimes stabilize regional abundance of a widespread freshwater fish. Ecol Appl. 2024;34:e2936.
Richardson JL, Urban MC, Bolnick DI, Skelly DK. Microgeographic adaptation and the spatial scale of evolution. Trends Ecol Evol. 2014;29:165-76.
Watson RA, Culley AV, Haase CG, Thomas MR, Brandt SL, Floyd MA, et al. Instream barriers contribute to population isolation of a small‐bodied, benthic, headwater‐specialist fish (Percidae). Ecol Freshw Fish. 2024;33:e12769.
Zarri LJ, Palkovacs EP, Post DM, Therkildsen NO, Flecker AS. The evolutionary consequences of dams and other barriers for riverine fishes. BioScience. 2022;72:431-48.
Raeymaekers JA, Maes GE, Geldof S, Hontis I, Nackaerts K, Volckaert FA. Modeling genetic connectivity in sticklebacks as a guideline for river restoration. Evol Appl. 2008;1:475-88.
Lowe WH, Addis BR, Cochrane MM. Outbreeding reduces survival during metamorphosis in a headwater stream salamander. Mol Ecol. 2024;33:e17375.
Quaglietta L, Paupério J, Martins FMS, Alves PC, Beja P. Recent range contractions in the globally threatened Pyrenean desman highlight the importance of stream headwater refugia. Animal Conservation. 2018;21:515-25.
Geismar J, Haase P, Nowak C, Sauer J, Pauls SU. Local population genetic structure of the montane caddisfly Drusus discolor is driven by overland dispersal and spatial scaling. Freshw Biol. 2014;60:209-21.
Escoriza D, Hassine JB. Niche separation among north-west African semi-aquatic reptiles. Hydrobiologia. 2017;797:47-56.
Quaglietta L, Gisbert J, García-Perea R, Fernández-González A. Galemys pyrenaicus. The IUCN Red List of Threatened Species 2024: eT8826A227363991 Accessed on 18 November 2024. 2024.
Wilson DE, Mittermeier RA. Handbook of the Mammals of the World – Volume 8. Insectivores, Sloths and Colugos: Lynx Edicions, Barcelona, Spain; 2018.
Gisbert J, García-Perea R. Historia de la regresión del desmán ibérico Galemys pyrenaicus (É. Geoffroy Saint-Hilaire, 1811) en el Sistema Central (Península Ibérica). Munibe Monogr Nat Ser. 2014;3:19-35.
Charbonnel A, Buisson L, Biffi M, D’Amico F, Besnard A, Aulagnier S, et al. Integrating hydrological features and genetically validated occurrence data in occupancy modelling of an endemic and endangered semi-aquatic mammal, Galemys pyrenaicus, in a Pyrenean catchment. Biol Conserv. 2015;184:182-92.
Igea J, Aymerich P, Fernández-González A, González-Esteban J, Gómez A, Alonso R, et al. Phylogeography and postglacial expansion of the endangered semi-aquatic mammal Galemys pyrenaicus. BMC Evolutionary Biology. 2013;13:115.
Querejeta M, González-Esteban J, Gómez A, Fernández-González A, Aymerich P, Gosálbez J, et al. Genomic diversity and geographical structure of the Pyrenean desman. Conservation Genetics. 2016;17:1333-44.
Escoda L, Hawlitschek O, González-Esteban J, Castresana J. Methodological challenges in the genomic analysis of an endangered mammal population with low genetic diversity. Scientific Reports. 2022;12:21390.
Escoda L, Gonzalez-Esteban J, Gomez A, Castresana J. Using relatedness networks to infer contemporary dispersal: Application to the endangered mammal Galemys pyrenaicus. Mol Ecol. 2017;26:3343-57.
Querejeta M, Fernandez-Gonzalez A, Romero R, Castresana J. Postglacial dispersal patterns and mitochondrial genetic structure of the Pyrenean desman (Galemys pyrenaicus) in the northwestern region of the Iberian Peninsula. Ecol Evol. 2017;7:4486-95.
Escoda L, Castresana J. The genome of the Pyrenean desman and the effects of bottlenecks and inbreeding on the genomic landscape of an endangered species. Evolutionary Applications. 2021;14:1898–913.
Whiteley AR, Fitzpatrick SW, Funk WC, Tallmon DA. Genetic rescue to the rescue. Trends Ecol Evol. 2015;30:42-9.
Janjua S, Peters JL, Weckworth B, Abbas FI, Bahn V, Johansson O, et al. Improving our conservation genetic toolkit: ddRAD-seq for SNPs in snow leopards. Conserv Genet Resour. 2019;12:257-61.
Shirasawa K, Hirakawa H, Isobe S. Analytical workflow of double-digest restriction site-associated DNA sequencing based on empirical and in silico optimization in tomato. DNA Res. 2016;23:145-53.
Peterson BK, Weber JN, Kay EH, Fisher HS, Hoekstra HE. Double digest RADseq: an inexpensive method for de novo SNP discovery and genotyping in model and non-model species. PLoS One. 2012;7:e37135.
Jurka J, Kapitonov VV, Pavlicek A, Klonowski P, Kohany O, Walichiewicz J. Repbase Update, a database of eukaryotic repetitive elements. Cytogenet Genome Res. 2005;110:462-7.
Andrews S. FastQC: A Quality Control Tool for High Throughput Sequence Data. Available online at: http://www.bioinformatics.babraham.ac.uk/projects/fastqc/. 2010.
Schubert M, Lindgreen S, Orlando L. AdapterRemoval v2: rapid adapter trimming, identification, and read merging. BMC Res Notes. 2016;9:88.
Li H, Durbin R. Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics. 2009;25:1754-60.
Purcell S, Neale B, Todd-Brown K, Thomas L, Ferreira MA, Bender D, et al. PLINK: a tool set for whole-genome association and population-based linkage analyses. Am J Hum Genet. 2007;81:559-75.
Felsenstein J. PHYLIP: Phylogeny Inference Package. Cladistics 1989;5:164–6.
Sampaio S, Gomes-dos-Santos A, Froufe E, Bolotov IN, Rozhnov VV, Rutovskaya MV, et al. The complete genome sequence of Desmana moschata, the Russian desman. Biodiversity Genomes. 2025.
Rochette NC, Rivera-Colon AG, Catchen JM. Stacks 2: Analytical methods for paired-end sequencing improve RADseq-based population genomics. Mol Ecol. 2019;28:4737-54.
Jombart T. adegenet: a R package for the multivariate analysis of genetic markers. Bioinformatics. 2008;24:1403-5.
Pritchard JK, Stephens M, Donnelly P. Inference of population structure using multilocus genotype data. Genetics. 2000;155:945–59.
Evanno G, Regnaut S, Goudet J. Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study. Mol Ecol. 2005;14:2611-20.
Kopelman NM, Mayzel J, Jakobsson M, Rosenberg NA, Mayrose I. CLUMPAK: a program for identifying clustering modes and packaging population structure inferences across K. Mol Ecol Resour. 2015;15:1179-91.
Anderson MJ. Permutational Multivariate Analysis of Variance (PERMANOVA). Wiley StatsRef: Statistics Reference Online; 2017. p. 15.
Oksanen J, Simpson GL, Blanchet FG, Kindt R, Legendre P, Minchin PR, et al. vegan: Community Ecology Package (R package version 2.6-8). https://CRAN.R-project.org/package=vegan. 2024.
Legendre P, Legendre L. Numerical Ecology. 3rd English Edition. Elsevier. ed2012.
Vogt J, Soille P, De Jager A, Rimaviciute E, Mehl W, Foisneau S, et al. A Pan-European River and Catchment Database. EUR 22920 EN. Luxembourg (Luxembourg): OPOCE. 2007.
Pebesma E, Bivand R. Spatial Data Science: With applications in R. Chapman and Hall/CRC. 2023.
van der Meer L, Abad L, Gilardi A, Lovelace R. sfnetworks: Tidy Geospatial Networks. 2024.
Long JM, Taylor AT, Buonaccorsi V. A conservation-oriented SNP panel for Smallmouth Bass (Micropterus dolomieu), with emphasis on Interior Highlands lineages. Conserv Genet Resour. 2020;13:47-59.
Ekblom R, Aronsson M, Elsner-Gearing F, Johansson M, Fountain T, Persson J. Sample identification and pedigree reconstruction in Wolverine (Gulo gulo) using SNP genotyping of non-invasive samples. Conserv Genet Resour. 2021;13:261-74.
Osmond DR, King RA, Stockley B, Launey S, Stevens JR. A low-density single nucleotide polymorphism panel for brown trout (Salmo trutta L.) suitable for exploring genetic diversity at a range of spatial scales. J Fish Biol. 2023;102:258-70.
Melero Y, Aymerich P, Luque-Larena JJ, Gosálbez J. New insights into social and space use behaviour of the endangered Pyrenean desman (Galemys pyrenaicus). European Journal of Wildlife Research. 2012;58:185-93.
Melero Y, Aymerich P, Santulli G, Gosálbez J. Activity and space patterns of Pyrenean desman (Galemys pyrenaicus) suggest non-aggressive and non-territorial behaviour. European Journal of Wildlife Research. 2014;60:707-15.
Escoda L, Fernández-González Á, Castresana J. Quantitative analysis of connectivity in populations of a semi-aquatic mammal using kinship categories and network assortativity. Molecular Ecology Resources. 2019;19:310-26.
François D, Ursenbacher S, Boissinot A, Ysnel F, Lourdais O. Isolation-by-distance and male-biased dispersal at a fine spatial scale: a study of the common European adder (Vipera berus) in a rural landscape. Conserv Genet. 2021;22:823-37.
Van Buskirk J, Jansen van Rensburg A. Relative importance of isolation-by-environment and other determinants of gene flow in an alpine amphibian. Evolution. 2020;74:962-78.
Quaglietta L, Fonseca VC, Hájková P, Mira A, Boitani L. Fine-scale population genetic structure and short-range sex-biased dispersal in a solitary carnivore, Lutra lutra. J Mammal. 2013;94:561-71.
Araújo MB, Guilhaumon F, Neto DR, Pozo I, Calmaestra R. Biodiversidade e Alterações Climáticas /Biodiversidad y Alteraciones Climáticas. Ministério do Ambiente e Ordenamento do Território & Ministerio de Medio Ambiente y Medio Rural y Marino. Lisboa/Madrid; 2012. p. 656.
Warnock WG, Rasmussen JB, Taylor EB. Genetic clustering methods reveal bull trout (Salvelinus confluentus) fine-scale population structure as a spatially nested hierarchy. Conserv Genet. 2009;11:1421-33.
Smyser TJ, Johnson SA, Page LK, Hudson CM, Rhodes OE. Use of experimental translocations of Allegheny woodrat to decipher causal agents of decline. Conserv Biol. 2013;27:752-62.
Robinson ZL, Coombs JA, Hudy M, Nislow KH, Letcher BH, Whiteley AR. Experimental test of genetic rescue in isolated populations of brook trout. Mol Ecol. 2017;26:4418-33.

No competing interests reported.

AdditionalFile1.pdf
Electronic supplementary material Additional File 1 AdditionalFile1.pdf Table S1: Details of the 115 samples used in this study, including the accession number from INSDC (when available – the others were retrieved from DRYAD), locality, phylogeographic unit, latitude, longitude, sex, and year of collection. For the Occidental phylogeographic unit samples collected in this study, additional information is provided on watershed of origin. Table S2: Genetic distances between individuals (see text for further details): A. for the Occidental phylogenetic unit; B. for the Douro river system; C. for the Sabor watershed. Table S3: Matrix of overland (dark green) and/or river (dark blue) distances between samples, in meters: A. for the Occidental phylogeographic unit; B. for the Douro river system; C. for the Douro river system using a mix of overland and river distances (in bold are the sample pairs from headwaters, using the overland distances instead of the river distances); D. for the Sabor watershed. Table S4: Discriminant Analysis of Principal Components posterior probabilities for the Occidental phylogeographic unit, taking 6 principal components into account: A. for 2 groups; B. for 3 groups; C. for 4 groups. Highlighted cells indicate the samples assigned to each cluster (K) based on their highest posterior probability. Table S5: Permutational Multivariate Analysis of Variance results testing genetic differences between samples on a genetic distance matrix: A. for the Occidental phylogeographic unit (Douro vs other river systems); B. for the Douro river system (Tua vs Sabor watersheds). Df stands for degrees of freedom, SS for sum of squares, R2 for proportion of variation explained and F for F-statistic. Values with p < 0.01 are highlighted in bold. Table S6: Mantel test using genetic distances and overland, river and a mix of overland and river distances for the Occidental phylogeographic unit, Douro river system and the Sabor watershed showing the correlation coefficient (r) and p-value for each test. p < 0.01 are highlighted in bold. Missing values (-) are due to differences in the spatial context of the analyses. For the Occidental unit, river distances were not calculated because the samples are distributed across different river systems, making direct river-based distances biologically meaningless. Similarly, for the Sabor watershed, mixed distances were not applicable, as all samples are within a single watershed and no adjacent headwaters between watersheds were considered.
AdditionalFile2.pdf
Additional File 2 AdditionalFile2.pdf Fig. S1: Diagram showing the flow of analyses for this manuscript. Fig. S2: Neighbour-joining tree based on the 115 samples and 110 SNPs using PHYLIP with 100 bootstraps and the Russian desman (Desmana moschata) as outgroup. Bootstrap values from 100 replicates are displayed for nodes with support ≥ 50%. The five phylogeographic units are highlighted. Fig. S3: Bayesian Information Criterion (BIC) values for each number of clusters for the Discriminant Analysis of Principal Components analysis taking 6 principal components into account for the samples from the Occidental population, based on 7,604 SNPs. Fig. S4: Discriminant Analysis of Principal Components for the Occidental phylogeographic unit, with 6 principal components and based on 7,604 SNPs, showing 3 clusters. Fig. S5: STRUCTURE analyses for the Occidental phylogeographic unit, based on 7,604 SNPs: A. the rate of change in the likelihood of the data as K increases (Delta K) by Evanno et al. [1]; B. the probability or likelihood of the data for each K (Prob(K)) based on the Bayesian clustering algorithm by Pritchard et al. [2]. Fig. S6: STRUCTURE analyses for the Douro river system, based on 8,771 SNPs: A. the rate of change in the likelihood of the data as K increases (Delta K) by Evanno et al. [1]; B. the probability or likelihood of the data for each K (Prob(K)) based on the Bayesian clustering algorithm by Pritchard et al. [2].

Download PDF

Journal Publication

published 23 Oct, 2025

Read the published version in BMC Ecology and Evolution →

Editorial decision: Revision requested
17 Sep, 2025
Reviews received at journal
15 Sep, 2025
Reviews received at journal
12 Sep, 2025
Reviewers agreed at journal
27 Aug, 2025
Reviewers agreed at journal
21 Aug, 2025
Reviewers agreed at journal
19 Jul, 2025
Reviewers invited by journal
15 Jul, 2025
Editor assigned by journal
15 Jul, 2025
Editor invited by journal
14 Jul, 2025
Submission checks completed at journal
11 Jul, 2025
First submitted to journal
11 Jul, 2025

You are reading this latest preprint version

Microgeographic genomic variation and connectivity in an endangered semiaquatic mammal

Status:

Journal Publication

Version 1

Abstract

Figures

Background

Methods

2.1. Sample collection

2.2. DNA extraction, PCR amplification and sequencing

2.3. Data analysis

2.3.1. Species-unit assignment

2.3.2. Microgeographic genomic analysis in the Occidental phylogeographic unit

Results

3.1 Phylogeographic subdivision in the Iberian desman

3.2 Microgeographic genomic analysis in the Occidental phylogeographic unit

Discussion

4.1 Phylogeographic subdivision in the Iberian desman

4.2 Microgeographic genomic analysis in the Occidental phylogeographic unit

Conclusion

Declarations

References

Additional Declarations

Supplementary Files

Status:

Journal Publication

Version 1