Introduction

Genetic polymorphism increases the chances of species survival in changing environments. Although a species can be very well adapted to prevailing environmental conditions, sudden environmental changes can be detrimental unless there is enough genetic potential for adaptation (Nolet and Rosell 1998; Somero 2010). Such changing environmental conditions may involve the physical properties of a species’ habitat, including climatic conditions, or changes in other biotic taxa, e.g., pathogens, food resources, or the relative abundance of interspecific competitors (Lande and Shannon 1996; Hoffmann and Sgrò 2011). High levels of genetic polymorphism provide species with a broader range of adaptive responses to such changes, through for example enhanced reproductive potential, behaviour, or immunocompetence (e.g., Hawley et al. 2005; Willi and Hoffmann 2008).

Natural- or human-mediated introductions of species into new areas usually involve only a small number of individuals. This can result in low levels of genetic polymorphism and an unsuitable genetic composition relative to the new environment because of a founder effect. For example, in the Finnish human population, a founder effect due to geographic and cultural isolation has caused enrichment of disease-causing gene variants and losses of others (Kääriäinen et al. 2017). Genetic diversity is a particularly critical determinant of the viability and conservation of rare or endangered species. (Re-)introduced populations may experience reductions in genetic polymorphism not only because of a founder effect but also since re-introductions are often followed with high initial mortality of released individuals (de Vos et al. 1956; Ebenhard 1988, McCleery et al. 2014). This hinders a population’s capacity to adapt to changes in the environment, such as climate change or an invading new competitor. Therefore, genetic monitoring in a re-introduction project is recommended, as for example, Haye et al. (2018) concluded in their study of re-introduction of highly endangered Common hamsters (Cricetus cricetus). However, the translocation of individuals into endangered populations may also serve as a powerful tool for vitalizing genetically impoverished populations, so-called genetic rescue (e.g., Frankham 2015; Whiteley et al. 2015).

The native Eurasian beaver (Castor fiber) was hunted to extinction from Finland in 1868 (Granit 1900; Lahti and Helminen 1974). The species was re-introduced to several sites in the 1930s with about 19 individuals from Norway (Lahti and Helminen 1969, 1980; Härkönen 1999). About seven North American beavers (C. canadensis) were also introduced, some of them to the same places, because their status as a distinct species was not known at the time (Ermala et al. 1989). The Eurasian beaver population became established only in the province of Satakunta where no North American beavers were introduced (Ermala et al. 1989). The founder population in Satakunta consisted, however, only of one female and two males (Linnamies 1956), suggesting that the level of genetic polymorphism in the founder population was low. This likely applies also to the original Norwegian source population which had undergone a bottleneck to only about 100 individuals in the 1880s (Collett 1897, according to Rosell and Parker 2011).

The present distribution of the Eurasian beaver in Finland covers Satakunta, parts of the provinces of Pirkanmaa, Pohjanmaa, and Etelä-Pohjanmaa, as well as western Lapland and the NE corner of Varsinais-Suomi (Fig. 1), and the population seems to be spreading in western Finland (Kauhala and Karvinen 2018; Alakoski et al. 2019). The North American beaver population increased and expanded much more rapidly than that of the native beaver and today occupies most of eastern and central Finland and some areas of northern Finland (Kauhala and Karvinen 2018).

Fig. 1
figure 1

The present distribution of the native Eurasian beaver in Finland covers the provinces of Satakunta, parts of Pirkanmaa, Pohjanmaa, and Etelä-Pohjanmaa, as well as Western Lapland and the NE corner of Varsinais-Suomi

Full size image

The North American beaver is currently spreading westwards in Finland, and at present the two species meet in Pirkanmaa (Ylöjärvi) and Etelä-Pohjanmaa (Kurikka; Supplementary Fig. 1). This raises concerns about the future of the native species in Finland—the Eurasian beaver survived only in Satakunta where no North American beavers were released, which suggests that the alien species are superior competitors over the native (Ermala 1996). The possible competitive inferiority of the Eurasian beaver in Finland is surprising, as observations imply that Eurasian beavers are displacing North American beavers in Russian Karelia (Danilov et al. 2011a). The outcome of competition between two species may be affected by the environment, as well as characteristics of the species, such as their habitat use, reproductive rate, behavior, or genetic polymorphism (Ayal 1969). No genetic studies of the level of polymorphism in the two beaver species in Finland have thus far been conducted.

Native Eurasian beavers were hunted to extinction also from many other countries in the nineteenth century, and only eight small relict populations survived (Nolet and Rosell 1998; Halley et al. 2012). Eurasian beavers were, however, re-introduced to several countries in the twentieth century (Nolet and Rosell 1998 and references therein). At present, they occur in many European countries, including the neighboring countries of Finland, i.e., Sweden, Norway, Russia, and the Baltic states (Halley et al. 2012). The re-introduced Eurasian beaver populations in Finland and Sweden originate from the Norwegian refugia (C. f. fiber), whereas the populations in NW Russia and the Baltic states are derived from a mixture of different populations or subspecies which have been relocated, or migrated, from several areas (Palionene 1965; Nolet and Rosell 1998; Halley et al. 2012; Veeroja and Männil 2018; Danilov et al. 2008).

The aim of the present study was to examine genetic polymorphism of the Eurasian beaver populations in Finland, Russian Karelia, and the Baltic states, to map potential source populations for future conservation actions involving genetic rescue. We hypothesize that the Finnish Eurasian beaver population exhibits lower genetic diversity than populations in the neighboring countries. Because the beavers in the Norwegian refugia had gone through a population bottleneck (Rosell and Parker 2011) and a preliminary study from the Swedish beaver population (Hartman 1994) indicated little genetic variation, we concentrated on the beaver populations which might be genetically more diverse than these populations.

Materials and methods

Sampling and DNA extraction

Tissue samples of Eurasian beavers (n = 200) were obtained from hunters from Finland (n = 82, in the years 2013–2018), Estonia (n = 49, in the years 2016–2017), Lithuania (n = 52, in a year 2018), and Russia (n = 17, in the years 2008–2015) (Supplementary Fig. 2). The hunting season is from 20 August to 30 April in Finland, from 1 August to 15 April in Estonia and Lithuania, and from 1 October to 28 (29) February in Russian Karelia. All samples were collected during the regular hunt, or for purposes of other studies. No beavers were killed for the present study. Genomic DNA was extracted from samples using Qiagen’s blood and tissue kit following manufacturer’s instructions.

Microsatellite analysis

We analyzed 12 microsatellites identified in both C. fiber and C. canadensis (microsatellites adopted from Frosch et al. 2014). Markers were grouped and amplified in three multiplex PCR reactions using the Qiagen multiplex PCR kit (Supplementary Table 1). The PCR reaction included 10 ng of DNA, 5 μl of 2x Qiagen Multiplex PCR master mix, 0.1–0.3 μM of each primer, and RNase-free water up to 10 μl. The PCR cycling conditions were the following: (1) an initial denaturation at 95 °C for 15 min, (2) 30 cycles of 30 s denaturation (94 °C), 90 s hybridization (60 °C), 60-s elongation (72 °C), and a (3) final 30-min elongation (60 °C). One microsatellite was amplified individually because it did not successfully amplify in multiplex reaction panels (Cf32, Supplementary Table 1) using the same reaction conditions as for the others. Fragment lengths were measured using 3500xL Genetic Analyzers (Applied Biosystems, CA, US) with LIZ600 size standard. The raw data were analyzed with GeneMapper Software 5 (Applied Biosystems, CA, US). The whole dataset generated during the current study is available in Supplementary File 1.

Diversity estimations

We calculated observed and expected heterozygosity at each locus and tested for deviations from Hardy-Weinberg equilibrium (HWE) for each population using ARLEQUIN ver 3.5.2.2 software (100,000 Markov chain steps; 10,000 dememorization steps) (Excoffier et al. 2005). Allelic richness was calculated with the correction of the lowest sample size 17. Allelic richness calculations were done with the HP-Rare package (Kalinowski 2005).

The within-population variation was assessed for excess or deficiency of heterozygosity relative to the allele number indicative of effective population size change using Bottleneck (Version 1.02.02; Piry et al. 1999). Expected values under mutation–drift equilibrium were based on 5000 iterations and the difference was tested using the Wilcoxon test. A two-phase mutation model was used assuming a single step mutation proportion of 90% and a variance of 12 for multistep mutations.

Population structure

Population pairwise FST using the distance method of number of different alleles was calculated with ARLEQUIN ver 3.5.2.2 software (Excoffier et al. 2005). We used a model-based Bayesian clustering method to infer the most probable number of genetic clusters (K) in the whole dataset using the structure program (Pritchard et al. 2000). K values from 2 to 10 were evaluated graphically by plotting Delta K (Evanno et al. 2005) and ln(Pr(X|K) was used to identify the k for which Pr(K = k) is the highest (calculated with Clumpak package (Kopelman et al. 2015)). Neighbor-Joining tree was calculated to describe the genetic distances between the populations using the Treefit program (Kalinowski 2009) and the similarity between the observed genetic distance between populations with the genetic distance in the tree was tested (Kalinowski 2009).

Results

All 12 microsatellites were polymorphic within our study populations. The number of alleles per locus varied from two to seven. Altogether we found 48 alleles. The two extreme populations of Eurasian beavers in our data were the Finnish and the Lithuanian populations: the Finnish population was homozygous for seven out of twelve microsatellite markers whereas the Lithuanian population was homozygous for only one microsatellite.

Observed heterozygosity (HO) varied from 0.11 in Finland to 0.42 in Lithuania (Table 1). Also, all four populations had one (not necessarily the same) marker deviating from HWE among the polymorphic loci. Considering the other diversity indices, this may be due to heterozygote deficit rather than the presence of null alleles. Allele frequencies per population are given in Supplementary Table 2.

Table 1 Summary statistics of the data. Number of samples within each population (n), number of loci, number of polymorphic loci, observed and expected mean heterozygosity from polymorphic loci, as well number of markers deviating from HWE are given s.d. = standard deviation
Full size table

Mean corrected allelic richness varied from 1.29 ± 0.39 in the Finnish population to 2.68 ± 1.17 in the Lithuanian population. Russian and Estonian beaver populations were positioned between the Finnish and Lithuanian populations based on the abovementioned genetic parameters (Table 2). The microsatellite markers used in this study partly originated from C. fiber (n = 8) and partly from C. canadensis (n = 4). This could explain the lack of heterozygosity but as seen from the allelic richness values (Table 2), markers from C. canadensis (Ccaxx) were not of poor quality in terms of diversity. In the Finnish population, the genetic diversity at the polymorphic loci deviated significantly (p = 0.06) from the expected for markers with the same number of alleles in mutation-drift equilibrium. The imbalance suggested high allele number and low genetic diversity (one-tailed P = 0.03), which is a within-population diversity pattern for a recently expanding population (Luikart and Cornuet 1998).

Table 2 Allelic richness across markers. Calculations are done for the subpopulations by country of the origin
Full size table

Based on population pairwise FST, the genetically closest populations are the Russian and Estonian populations (FST = 0.07, P < 0.001) whereas the Finnish and Russian populations seem to be genetically the most distant (FST = 0.80, P < 0.001) among the studied populations (Table 3). The most likely number of genetic clusters in the data varies according to the method used to estimate it. For our data, Delta K was not informative whereas Pr(X|K) indicated K = 5 as being the most probable number of clusters (Supplementary Fig. 3). Typically, it is better to use the value of K that captures most of the structure and is biologically sensible (Pritchard et al. 2000). At the first level (K = 2), the Bayesian clustering indicates the separation of the Finnish Eurasian beaver population from the other populations (Fig. 2). At K ≥ 4, the Lithuanian population shows some evidence of structure and supports K = 5 to be the most likely value. The neighbor-joining tree presentation of the populations shows the genetic dissimilarity of the Finnish population (Supplementary Fig. 4).

Table 3 Population pairwise FST calculated using the distance method of number of different alleles. P < 0.001 for all pairwise comparisons
Full size table
Fig. 2
figure 2

Clustering of individuals using a Bayesian clustering method from K = 2 to K = 10. Each individual is represented by a vertical bar, and populations are explained below the figure. The K is the number of genetic clusters in the Bayesian clustering method used (Pritchard et al. 2000). The plot was created with the Clumpak package (Kopelman et al. 2015)

Full size image

Discussion

Our main results support the hypothesis that the Eurasian beaver population in Finland has a lower degree of genetic diversity than populations in Russian Karelia and the Baltic states. The populations in Finland and Sweden originate from individuals re-introduced from Norway, i.e., the Norwegian refuge, whereas populations in Russian Karelia and the Baltic states are derived from a mixture of different populations originating from several areas. The Finnish founder population probably consisted only of three surviving individuals, and the source population in Norway had also gone through a bottleneck of about 100 individuals (Linnamies 1956; Nolet and Rosell 1998; Rosell and Parker 2011). The level of genetic variation seems to be low also in the Swedish population despite its larger size (ca. 4000 individuals in Finland and up to 100,000 in Sweden; Hartman 1994; Nolet and Rosell 1998; Kauhala and Karvinen 2018), which also reflects’a bottleneck in the Norwegian source population (Ellegren et al. 1993). Our results suggest that the populations which originate from several different populations are more polymorphic than the Finnish population which originates from only a few individuals from the Norwegian refuge population.

The Finnish population is genetically strikingly different from the Russian and Baltic populations (e.g., Fig. 3, Supplementary Fig. 4). The limited genetic resources in the Finnish population can be inferred from the diversity indices. For example, allelic richness, which is considered an indicator of a decrease in population size or of a past bottleneck (Nei et al. 1975; Foulley and Ollivier 2006) shows lower values for the Finnish population compared to the other populations in the present study and Eurasian beaver populations elsewhere (e.g., Frosch et al. 2014). The same trend was observed in heterozygosity, further demonstrating the reduced genetic diversity within the Finnish beaver population.

Based on the excess or deficiency of heterozygosity relative to the allele number found in the present study, and the known history of re-introductions, the Finnish population has undergone past genetic bottlenecks. While the past bottleneck is visible in the Finnish population as a low number of variable loci, the very recent population growth has given space for novel alleles in these fast mutating microsatellite markers. The limited number of loci generally reduces the statistical power to detect deviations from mutation-drift equilibrium (Piry et al. 1999). Presently, the other populations had at least twice the number of polymorphic loci available for the test, and the power to detect effective population size changes was higher for them than in the Finnish population. Results suggest the recent effective population increase has been the most dramatic in the Finnish population.

The inference is supported by the presence of rare private alleles observed only in the Finnish population (Supplementary Table 2). The same phenomenon, a slight increase of private alleles, was also seen with the recently introduced platypus (Ornithorhynchus anatinus) population (Furlan et al. 2012).

In Russian Karelia, the Eurasian beaver population originates from two localities, Archangelsk and Leningrad regions, where they were translocated mainly from Voronezh and Belarus, and to a lesser extent from other areas (Nolet and Rosell 1998; Danilov et al. 2011a, b; Danilov and Fyodorov 2016). Beavers colonized Estonia both from the southeast and from two areas of Russia: re-introduced populations from the Leningrad region and Voronezh area (Laanetu 2014). The Lithuanian population showed within population structure which indicates the presence of at least two subpopulations. At least six locations of reintroduced (78 individuals from Russia, Voronezh) or naturally immigrated animals from Belarus to Lithuania (1947–1959) have influenced the formation of the present Lithuanian population (Palionene 1965; Prūsaitė et al. 1988; Nolet and Rosell 1998; Halley et al. 2012). Phenetic (non-metric skull characters) and genetic (allozyme electrophoresis) analyses (Ulevičius and Paulauskas 2003) showed that Lithuanian beaver subpopulations of different origins have preserved their specificity even in dense populations, 30–40 years since the first releases and immigration. This specificity might be the result of the founder effect in small geographically and temporally isolated groups of reintroduced beavers (Ulevičius and Paulauskas 2003). The more diverse origins of Russian and Baltic beaver populations may well explain the higher degree of heterozygosity and allele richness in these countries, although all remnant populations before introductions were small. Our results suggest that the populations which originated from several different populations are more polymorphic than the Finnish population which originated from only a few individuals from the Norwegian refuge population.

The rapid growth of the Swedish Eurasian beaver population, relative to the Finnish population, may be due to the larger number of individuals introduced (ca. 90, Hartman 1994) and possibly better habitats. It is also possible that low genetic diversity in Eurasian beavers is associated with a reduced tolerance to interspecific competition. North American beavers do not inhabit Sweden, whereby Eurasian beavers have had the opportunity to expand in the absence of competitors. This is, however, also true for Eurasian beavers in western Finland (Satakunta) where North American beavers have never occurred. Despite the absence of North American beavers, the population of Eurasian beavers grew slowly. This might have been affected by other factors, such as poor habitat quality, as western Finland in primarily non-optimal agricultural habitat (Alakoski et al. 2019).

Although recent data indicates that in western Finland the Eurasian beaver population has been increasing in size (see also Kauhala and Karvinen 2018), interspecific competition with North American beavers may have limited the eastward dispersal of Eurasian beavers in the country (Alakoski et al. 2019). However, evidences also suggest that the expansion of Eurasian beavers in Russian Karelia led to displacement of American beaver (Danilov et al. 2011a). Data on the level of genetic diversity of North American beavers in Finland or Russia do not exist. We cannot, therefore, exclude the possibility that the difference in the competitive ability of Eurasian beavers in the two countries is due to, for example, different environments. Although the habitat requirements of the two species are similar, some differences exist (Alakoski et al. 2019; Danilov et al. 2011b).

Obligate inbreeding resulting from founder events decreases genetic diversity, and may, in turn, negatively impact the dynamics of the population. Inbreeding depression may reduce the fecundity of individuals, affect juvenile mortality, limit the breadth of their immune function, reduce their physiological capacity to adapt to environmental change, and even alter the lifespan and behaviour of individuals (Crnokrak and Roff 1999; DeRose and Roff 1999; Keller and Waller 2002; Liberg et al. 2005; Halley et al. 2012). Very little is known of how each of these effects is associated with the genetic diversity of Eurasian beavers, but they certainly warrant future research.

Genetic diversity of a population could be improved by translocating individuals from other areas, i.e., with genetic rescue (Whiteley et al. 2015). The individuals used should preferably have a higher level of heterozygosity and allele richness than found in the target population, to create hybrid vigor and thus improve the viability of the target population. Therefore, if genetic rescue were to be attempted to improve the future of the European beaver in Finland, beavers from a genetically different and more polymorphic population, e.g., from Lithuania, should be translocated to Finland.

Conclusions

The low genetic diversity of the Eurasian beaver population in Finland may be deleterious to the population especially in a changing environment. Genetic rescue could be the best solution to improve the viability of the population. The utilization of whole genome sequencing techniques (for example RADsequencing) could reveal the detailed genetic structures of populations and help understand the genetic features of adaptation in both beaver species (C. fiber and C. canadensis). Also, future research on osteoarcheological remains could be employed to resolve the origin (western or eastern) of the original Finnish beavers. Furthermore, signs of possible inbreeding depression in the Eurasian beaver population in Finland should be examined both morphologically, physiologically, and genetically.