Biodiversity is currently declining at the global scale. Apart from species declines and lowered abundances, the loss of genetic diversity is equally concerning as it may undermine fitness and the potential to adapt to future environmental change. We compared genetic diversity of historical and recent Alpine populations of two butterfly species, Lycaena helle and L. hippothoe, over a period of about 10 years. Using microsatellite markers, we found no changes over time in L. helle, while genetic diversity decreased, and differentiation increased in L. hippothoe. Lycaena helle inhabits peat bogs and wetland fallows with populations being strongly isolated, while L. hippothoe used to occur in population networks on hay meadows, with the latter being strongly exposed to agricultural intensification. We conclude that currently L. hippothoe populations are strongly declining due to changes in land use, resulting in genetic erosion potentially due to the collapse of population networks.

Population sampling and study area

We studied L. helle populations in Germany (Rhineland-Palatinate, Bavaria) and L. hippothoe populations in the European Alps, including sampling sites in Austria (Tyrol) and Italy (South Tyrol) (Fig. 1; Table 1). We sampled 20 individuals each in five populations of L. helle (except for population 5 with only 15 individuals in 2005-06 and population 1 with only 14 individuals in 2017) in 2005-06 and in 2017. For L. hippothoe, eight populations with 18 individuals each were sampled in 2006-07 and in 2017-19. This design allowed us to investigate changes in the genetic constitution of five and eight populations over time in two different Lycaena species. Butterflies were caught with an insect net in the field. One leg was removed and stored in 100% ethanol until DNA extraction. Subsequently, we released the butterflies immediately. To avoid the recapture of individuals, we marked all sampled individuals with a blue dot on the left hindwing before release.

DNA extraction and microsatellite analyses

DNA was extracted with the E.Z.N.A. ® Tissue DNA Kit (Omega bio-tek, Germany) following the manufacturer’s instructions. We used five microsatellite markers for L. helle and 14 for L. hippothoe, following the polymerase chain reaction (PCR) amplifications and thermal cycling profiles given in Habel et al. (2008) and Trense et al. (2019) (Supplementary Table S1). Genotyping of microsatellites was performed on an ABI 3130XL sequencer (Applied Biosystems, Germany). The allele length of microsatellite loci was determined using the program GENEIOUS version 11.1.5 (Kearse et al., 2012; https://www.geneious.com).

Data analyses

For comparisons, we performed separate analyses for the five past (2005-06; 5_past) and recent (2017; 5_recent) population samples of L. helle, and the eight past (2006-07; 8_past) and recent (2017-19; 8_recent) population samples of L. hippothoe. The microsatellite loci were checked for potential null alleles using the program MICROCHECKER version 2.2.3 (van Oosterhout et al., 2004). Numbers of alleles, alleles per locus, polymorphic loci, gene diversity, inbreeding coefficient (F_IS), pairwise F_ST and R_ST values, observed and expected heterozygosity, and analyses of molecular variance (AMOVA) were computed with ARLEQUIN version 3.5 (Excoffier & Lischer, 2010). Numbers of private alleles, G’’_ST and Jost’s D values, Nei’s genetic distance (Nei, 1972), and unbiased Nei’s genetic distance were calculated using GENALEX version 6.51 (Peakall & Smouse, 2006, 2012). FSTAT version 2.9.3 (Goudet, 1995) was used to calculate allelic richness. FREENA was used to control for the impact of null alleles on F_ST values (Chapuis & Estoup, 2007). Hardy-Weinberg equilibrium (HWE) and linkage disequilibrium were tested with GENEPOP version 4.2 (Raymond & Rousset, 1995; Rousset, 2008). Resulting P-values for linkage disequilibrium and pairwise F_ST values were corrected for multiple testing using the Benjamini-Hochberg method (Benjamini & Hochberg, 1995). To test whether results are affected by the different numbers of microsatellite loci across species, we additionally performed analyses based on random selections of five out of 14 loci for L. hippothoe, using 15 replicates (see Supplementary Table S2 – S17). To assess statistical significance across replicates, we calculated the harmonic mean p-value (HMP) for each molecular index. The HMP significance threshold was 0.039. Variation in molecular indices across past and recent populations was analysed with Mann-Whitney U and Wilcoxon tests using STATISTICA version 12.0 (Tulsa, StatSoft). We used both approaches because the individuals tested were independent (different individuals were sampled within the different periods) while populations were sampled repeatedly.

Furthermore, the genetic population structure was examined by a Bayesian analysis, as implemented in STRUCTURE version 2.3.4 (Pritchard et al., 2000), using an admixture model to allocate individual microsatellite genotypes to K gene clusters. For the four data sets, we tested 10 iterations for each K = 1-5 for L. helle and K = 1-8 for L. hippothoe with 50 000 burn-in steps and 500 000 Markov chain Monte Carlo (MCMC) cycles. The Delta K method (Evanno et al., 2005) was used to estimate the most likely number of clusters in STRUCTURE HARVESTER version 0.6.94 (Earl & von Holdt, 2012). Afterwards, we used CLUMPAK (Kopelman et al., 2015) to combine the cluster assignments across replicate analyses and to illustrate aligned cluster assignments.