Glacial allopatry vs. postglacial parapatry and peripatry: the case of hedgehogs

Abstract

Although hedgehogs are well-known examples of postglacial recolonisation, the specific processes that shape their population structures have not been examined by detailed sampling and fast-evolving genetic markers in combination with model based clustering methods. This study aims to analyse the impacts of isolation within glacial refugia and of postglacial expansion on the population structure of the Northern White-breasted hedgehog (Erinaceus roumanicus). It also discusses the role of the processes at edges of species distribution in its evolutionary history. The maternally inherited mitochondrial control region and the bi-parentally inherited nuclear microsatellites were used to examine samples within the Central Europe, Balkan Peninsula and adjacent islands. Bayesian coalescent inference and neutrality tests proposed a recent increase in the population size. The most pronounced pattern of population structure involved differentiation of the insular populations in the Mediterranean Sea and the population within the contact zone with E. europaeus in Central Europe. An interspecies hybrid was detected for the first time in Central Europe. A low genetic diversity was observed in Crete, while the highest genetic distances among individuals were found in Romania. The recent population in the post-refugial area related to the Balkan Peninsula shows a complex pattern with pronounced subpopulations located mainly in the Pannonian Basin and at the Adriatic and Pontic coasts. Detailed analyses indicate that parapatry and peripatry may not be the only factors that limit range expansion, but also strong microevolutionary forces that may change the genetic structure of the species. Here we present evidence showing that population differentiation may occur not only during the glacial restriction of the range into the refugia, but also during the interglacial range expansion. Population differentiation at the Balkan Peninsula and adjacent regions could be ascribed to diversification in steppe/forest biomes and complicated geomorphology, including pronounced geographic barriers as Carpathians.

Keywords: Balkan; Erinaceus roumanicus; Founder effect; Genetic differentiation; Interspecies interactions; Landscape genetic; Phylogeography.

Figure 1. Sampling localities for the 314 samples of the Northern White-breasted Hedgehog (Erinaceus roumanicus) in central Europe and the Balkan region.

Each sampled individual is represented by a pink dot. The distribution map of the species is shown in the smaller square. Map was done using ArcGIS (ESRI, 2014) using publicly available layers (Bathymetry-EMODnet, Administrative areas-GADM database, Elevation-WorldClim database, Rivers-EEA hydrographic data set). Inset of the figure is based on an IUCN distribution map of species (Amori et al., 2016) edited for the purposes of this article in Photoshop CS3 (Adobe).

Figure 2. Structuring of the nuclear data.

(A) Population assignment using Bayesian clustering analysis in Structure for the 260 Erinaceus roumanicus samples in Central Europe and the Balkan region. The results are shown for K = 2 to K = 7. The highest ΔK was for K = 3, and the second highest ΔK was for K = 7. (B) Geographical cluster assignments implemented by TESS. The results display five colours for Kmax = 6.

Figure 3. Median-joining network of the mitochondrial control region haplotypes for 296 Erinaceus roumanicus samples.

Haplotypes are represented by circles, which are proportional to the haplotype frequency. The hypothesized haplotypes are represented by red dots. Colour codes are used according to the population assignments of the samples using mtDNA data from the Geneland analysis (see Fig. 5C or Table S1). Lines between nodes represent a nucleotide changes.

Figure 4. Bayesian skyline plot based on mtDNA of 296 Erinaceus roumanicus samples.

The X axis represents the time in mutation units per nucleotide position. The Y axis represents a correlate of the population size (Ne). The black line shows the median of the Ne estimation, and the 95% confidence intervals are indicated by the blue areas.

Figure 5. Genetic landscape shape interpolation analysis.

Genetic landscape shape interpolation analysis based on variation of (A) the control region of mtDNA using 296 Erinaceus roumanicus samples and (B) nine nuclear microsatellite loci using 260 E. roumanicus samples. The X and Y axes reflect the geographic coordinates within the study area, while the Z axis represents the average genetic distances between the analysed samples. The Bayesian landscape genetic analysis using Geneland illustrates the spatial distribution of the populations using (C) the control region of mtDNA and (D) nine nuclear microsatellite loci. Figure was done using publicly available outlines (https://CRAN.R-project.org/ package=rworldmap) and edited in Photoshop CS3 (Adobe).