Distribution patterns of the native Eurasian and the non-native North American beaver in Finland—possible factors affecting the slow range expansion of the native species

Distribution patterns of species are affected by resource availability, dispersal, disturbance and population dynamics. The smaller population size and range of the native Eurasian beaver (Castor fiber) compared to the non-native North American beaver (Castor canadensis) in Finland raise questions on reasons for the slower range expansion of the native species. We compared the population growth rates and the spread of both species from their release sites. We also studied the factors possibly affecting the spread of the Eurasian beaver in South western Finland in more detail. We found that the North American beaver has spread longer distances than the Eurasian beaver, but we did not find evidence for movement barriers constraining the expansion rate of the native species. Lack of high-quality habitats does not seem to constrain the expansion to nearby areas either. Despite this, the Eurasian beaver population has grown to a high density close to its reintroduction site, and it has started to spread to novel areas only recently. We conclude that the expansion of the native beaver in Finland seems to be controlled by factors other than those related to barriers for movement: movement behavior and population dynamics, which require further investigation.


Introduction
Distribution patterns of species are largely affected by resource availability, dispersal, disturbance and population dynamics (Guisan and Thuiller 2005). Habitat fragmentation and a low proportion of suitable habitats due to e.g. anthropogenic landscape, and lack of connectivity between habitat patches can constrain the range expansion of species (Wilson et al. 2009; Barros et al. 2016). In addition, natural and anthropogenic barriers, e.g. mountains, large water bodies, dams and roads, may prevent dispersal (Shephard et al. 2008;Hapeman et al. 2011;Bracken et al. 2015;Machado et al. 2018). For a species introduced to a new area, the location of the introduction site and the size of the founder population largely determine the present distribution and population size (Stephens and Sutherland 1999). In newly formed populations, the rate of spread may also decrease due to problems related to a small population size, such as difficulty in finding mates at low densities (Allee 1931;Stephens and Sutherland 1999). Furthermore, populations can become clumped due to an aggregation of resources and high-quality habitats or due to conspecific attraction (Stephens and Sutherland 1999).
The management of non-native and native species requires information on the factors controlling the distribution patterns. For example, in Finland, there are two beaver species, the native (reintroduced) Eurasian beaver (Castor fiber) and the non-native North American beaver (C. canadensis). The native Eurasian beaver was hunted to near extinction in Europe in the late 1800s, and only eight populations and a total of 1200 individuals survived in small refuges (Nolet and Rosell 1998). In recent decades, the species has been reintroduced into several countries in Europe (e.g. Halley et al. 2021). The original Eurasian beavers were hunted to extinction in Finland in the late nineteenth century, 1 3 the last report of a hunted beaver being from 1868 (Granit 1900;Lahti 1972;Lahti and Helminen 1974). During the 1930s, beavers were reintroduced to Finland: 17-19 Eurasian beavers were brought to Finland from Norway in 1935-1936(Härkönen 1999 and seven North American beavers were introduced in 1937 ( Helminen 1969, 1980;Ermala et al. 1989). Both species were released to several places in Finland. However, according to the literature, only three individuals (one female and two males) of Eurasian beavers survived and reproduced (Lahti and Helminen 1980), resulting in low genetic diversity in the Eurasian beaver (Iso-Touru et al. 2020). The number of North American beavers was not much greater, and according to the literature, only two pairs reproduced (Härkönen 1999) but information on the genetic diversity of the species is not currently available in Finland. The Eurasian beaver population survived and started to increase slowly only in Satakunta where no North American beavers were released, and the species has not been translocated to new areas. The majority of the Finnish Eurasian beavers still occur within this area in South western Finland (in and around Satakunta; Fig. 1a). Eurasian beavers have dispersed outside Satakunta presumably at the turn of the twenty-first century but DNA-identification has been possible just in recent years. From 2010, the number and density of Eurasian beavers has increased especially in Pohjanmaa and Etelä-Pohjanmaa (Luke 2021). Smaller populations are found in western Lapland where the species has spread from Sweden (Kauhala and Timonen 2016) and close to the South eastern border of Finland where it likely dispersed from Russia (only seven observations; Luke, unpublished data). Nowadays, the population size of the Eurasian beaver is estimated to 3700-5000 individuals (Luke 2021). The Eurasian beaver is classified as near threatened in Finland (Hyvärinen et al. 2019), but hunting has been possible with a license granted by the Finnish Wildlife Agency. In 2019-2020 up to 400 hunting licenses could be permitted (The Finnish Wildlife Agency 2020), and a total of 334 Eurasian beavers were hunted during the hunting season from 20 August to 30 April (A. Impola, pers. comm.).
North American beavers flourished especially well in Sääminki in the province of Savo in South eastern Finland, and a couple of decades later, they were further translocated Fig. 1 Distribution of beavers in Finland in relation to introduction sites and main drainage basins. a Distribution of the Eurasian and North American beaver observations (lodges, dams, feeding sites, sightings and DNA-samples) with coordinates in Finland (n = 623 for the Eurasian beaver, n = 589 for the North American beaver). The original area of the Eurasian beaver in Satakunta in blue and the novel areas in Pirkanmaa, Etelä-Pohjanmaa and Pohjanmaa in light red. b The main drainage basins in Finland. Blue triangle = Eurasian beaver observation, red dot = North American beaver observation. The main drainage basins cross the Finnish border in eastern Finland and Lapland (see Fig. 1a for comparison) north and west to e.g. Lapland, Pohjois-Karjala and Häme (Ermala 1996). Beavers also dispersed naturally to the Russian side of the border where they now occur in a large area (e.g. Danilov and Fyodorov 2016). The number of North American beavers reached 10,000 in the late 1990s, but the numbers reported by hunters during monitoring counts decreased in the early twenty-first century (Ermala et al. 1999;Kauhala and Turkia 2013;Brommer et al. 2017), because since 2001 no hunting license for the North American beaver was required. The range of the North American beaver has continuously increased, and it has spread westwards close to the range of the Eurasian beaver (Kauhala and Turkia 2013; Alakoski et al. 2019). Recently, the estimated number of North American beavers has exceeded 10,000 individuals (Luke 2021). Annually approximately 2000-7000 North American beavers were hunted in Finland between 2010 and 2019 (statdb.luke.fi). The North American beaver's range covers most of eastern and central Finland as well as parts of Lapland (Fig. 1a). That is, the range is much larger than that of the Eurasian beaver, which raises the question of what constrains the population spread in the Eurasian beavers' South western range? That is, what are the underlining mechanisms behind the larger range of the North American beaver compared to that of the Eurasian beaver in Finland.
Both beaver species are semiaquatic, monogamous and territorial (Wilsson 1971;Nolet and Rosell 1994), and Parker et al. (2012) concluded that niche overlap is virtually complete between the species. They live in family groups of the reproducing female and male, offspring of the year and, often, subordinate beavers that disperse at the earliest at the age of 2 years (e.g. Müller-Schwarze 2011). The North American beaver has been reported to have larger litter sizes (Danilov 1995;Danilov et al. 2011) and family groups than the Eurasian beaver (Rosell and Parker 1995). Beavers use mainly deciduous trees for foraging, but sometimes also coniferous species are consumed (Danilov et al. 2011;Kauhala and Karvinen 2018), and mixed forests can be used as habitats (Kauhala and Turkia 2013;Kauhala and Karvinen 2018). In summer, aquatic vegetation and terrestrial herbs are also commonly utilized (Danilov et al. 2011). In addition, agricultural fields may offer extra forage (Alakoski et al. 2019). Although beavers can move on land, they mainly disperse along watercourses (e.g. Leege 1968;Novak 1987, Hartman 1994a, b, Müller-Schwarze 2011, and watershed divides are expected to decrease the expansion rate of beavers (Hartman 1995). Beavers are described as ecosystem engineers, as they can greatly modify the landscape where they live in (Brazier et al. 2020). This mainly results from foraging, and their habit of building dams to raise water levels for shelter and easier access to forage (Stringer and Gaywood 2016;Puttock et al. 2020). Simultaneously, beavers increase biodiversity by creating habitats for other species (Stringer and Gaywood 2016). Therefore, understanding the factors affecting the distribution of beavers also creates information on the location of suitable habitats of other species, and possible beaver damages to agriculture or forestry.

Objectives of the study
In this study, we (1) analyzed the distribution patterns of native and non-native beavers in Finland. For this aim, we compared (a) the distances of both beaver species from their release sites, (b) the sizes of the main drainage basins in their ranges, (c) the population growth rates and the increase in the number of occupied areas compared to the increase in population size. We hypothesized that these measures are larger for the non-native beaver with a larger present range in Finland compared to those of the native beaver.
Further, we investigated (2) the factors possibly affecting the spread of the Eurasian beaver in SW Finland. We analyzed (a) the effect of the natural barriers, i.e. watershed divides, on the occurrence of Eurasian beavers. Next, we divided the SW range of the native beaver into an original range (from now on called original area) and to newly colonized areas (from now on called novel areas). Between the original and novel areas, we analyzed (b) if artificial barriers, i.e. dams and sluices, affected the Eurasian beaver's spread along watercourses, and (c) are there differences in habitat characteristics in the sites used by the beavers. We hypothesized that (a) the occurrence of watershed divides decreases the local abundance of Eurasian beavers; (b) artificial barriers increase the distance between occupied sites; and (c) the beavers use poorer habitats in the novel areas, if habitat quality is the factor limiting spread from the original area. We discuss the possible reasons for the present range and the slow expansion rate of the Eurasian beaver population in Finland.

Study area
The landscape in Finland is dominated by bodies of water and industrial coniferous and mixed forests, where scots pine (Pinus sylvestris) is the most dominant species, along with Norway spruce (Picea abies), downy birch (Betula pubescens) and silver birch (B. pendula) as common species. Agricultural areas and denser human populations are found mainly in South western and southern Finland (for details, see Alakoski et al. 2019). In this study, we used data for the Eurasian and North American beaver observations in the whole country. In addition, we were especially interested in the SW range of the Eurasian beaver. We defined that the Eurasian beavers original reintroduction area is the region of Satakunta in SW Finland, with the novel areas being the adjacent regions of Pirkanmaa, Etelä-Pohjanmaa and Pohjanmaa in SW Finland (Fig. 1a).

Data for beavers
We used the data of the Finnish Wildlife Agency on beaver observations with exact coordinates from August 2017 to August 2018, including beaver lodges [Eurasian (E) 169; North American (NA) 179], dams (E 41; NA 53), feeding sites (E 50; NA 59) and other sightings (E 107;NA 197). The data were collected with a mobile app OmaRiista, where citizens can report their hunting bag/catch and wildlife observations directly on a digital map. Hunters are especially asked to report inhabited beaver lodges in autumn during the moose (Alces alces) hunting season every third year when Luke (Natural Resources Institute Finland) is carrying out the monitoring counts of beavers. In addition, 194 and 73 locations of DNA-samples (Iso-Touru et al. 2020 and unpublished data), and 62 and 28 locations of hunted beavers (species identified from skull morphometry; Kauhala and Timonen 2016) were included for the Eurasian beaver and the North American beaver data, respectively. Altogether, there were 623 locations for the Eurasian beaver and 589 for the North American beaver (Fig. 1a). Species determination from DNA or skulls was reliable, whereas, the beaver species of locations of lodges, dams, feeding sites and other sightings was assumed on the basis of the history of beavers in Finland (following earlier studies: Brommer et al. 2017;Alakoski et al. 2019Alakoski et al. , 2020; see also Iso-Touru et al. 2020; Fig. 1a). Although the data for individuals identified to species (267 DNA samples and 90 skulls) indicate that there are areas of sympatry, the two beaver species live mainly in separate areas in Finland ( Fig. 1a; Kauhala and Karvinen 2018;Alakoski et al. 2019;Iso-Touru et al. 2020, Sjöberg andBelova 2020). Thus, few misidentifications can occur in the data near and within the areas of sympatry.

Environmental variables
Variables used in the habitat use models were: aquatic habitat type, broadleaved, mixed, coniferous and transitional forest, forage, broadleaved trees, and agricultural and urban areas. In addition, for the analyses for the distribution patterns, data for shoreline, dams, sluices, main drainage basins and watershed areas were used. All environmental variables were analyzed in ArcMap (ESRI 2011).

Aquatic habitat
The watercourses (as polylines) and water areas (as polygons) in Finland were obtained from the data of the National Land Survey of Finland (2015; topographical map 1:100,000) and merged together as polylines, i.e. shorelines from water areas (from now on called shoreline). Therefore, the "shoreline" for watercourses includes only one line. A network was made from the shoreline in the two species' ranges in GIS. In addition, aquatic habitat type was also used as an explanatory variable. The aquatic habitat type composition was computed from a raster data (pixel size 20 × 20 m) according to the classes of the National Land Survey of Finland (excluding canals and sea water, which did not occur in the data). Also, reservoirs did not occur in the beaver habitats and they were also rare in the environment, so they were excluded from the analyses. Therefore, there were five aquatic habitat type classes: (1) lakes, (2) streams < 2 m, (3) streams 2−5 m, (4) streams 5−20 m, and (5) streams > 20 m. In addition, we used the data for dams and sluices (National Land Survey of Finland 2015; topographical map 1:100,000), and watershed areas, i.e. an area where water flows to e.g. one lake, and main drainage basins, i.e. a larger watershed area system (SYKE 2019).

Anthropogenic disturbance
Composition of agriculture was computed from the Corine land cover data for Finland for year 2018 (20 × 20 m per grid cell; SYKE 2020) classes 2.1.1−2.4.3. Similarly, the composition of urban areas (classes 1.1.1 − 1.2.2: including urban fabric, industrial, commercial and transport units) was computed.

Distribution patterns of the Eurasian and the North American beaver a. Distances to the release sites
Distances from the release sites to the present beaver observations were computed as straight-line distances and as distances along shoreline (study objective 1a in aims). There are twelve original introduction and translocation sites for the North American beaver and one reintroduction site for the Eurasian beaver (Fig. 1a). The (re)introduction and translocation sites (from now on called release sites when referring to both species) were collected from the literature (Linnamies 1956;Lahti and Helminen 1980), and only sites where the species have been continuously present since the (re)introductions were selected. Due to low reporting activity in Lapland, the presence of the North American beaver in Lapland during the past decades is uncertain (Ermala et al. 1999) but we assumed that it has occurred there since the introduction. For the Eurasian beaver, we computed only the distances from the SW range since the species has dispersed to W Lapland from Sweden and to SE Finland from Russia (only seven observations). A total of 479 and 395 routes were computed along the shorelines, and 506 and 459 as straight-line distances, for the Eurasian beaver and the North American beaver, respectively. Thus, not all observations could be connected along the shorelines.

b. Sizes of the main drainage basins
We compared the sizes of the main drainage basins occurring in the whole ranges of both species, i.e. also Lapland and SE Finland for the Eurasian beaver (Fig. 1b).

c. Population growth rates and occupied areas compared to population size
The present densities of beavers in Finland can be seen in Table 1. We used data for the estimated population sizes of the two species collected from the literature (Linnamies 1956;Lahti and Helminen 1980;Ermala et al. 1989;Ermala 1996) and from the monitoring counts of Luke (formerly the Finnish Game and Fisheries Research Institute) starting from 1995 (for details, see Brommer et al. 2017). Based on these values we computed the annual growth rate (r) relative to the previous population size of both species in different years ( where Nt 1 and Nt 2 are population sizes in consecutive estimation times t 1 and t 2 , respectively. We compared the relative population sizes (population size in year t /population size in a year with the highest population size) to the relative number of occupied areas (hunting clubs reporting beaver observations/the highest number of reporting hunting clubs; Table 2) in different years visually. This was done to see if the population sizes have increased faster than the number of occupied areas or vice versa (Hartman 1995). The number of hunting clubs that reported beaver observations in their area was collected from the data of the monitoring counts. The hunting area of a hunting club is usually 20-100 km 2 (Finnish Wildlife Agency 2020). Table 1 Number of all reported lodges in 2020 with a percentual change to the 2017 number of lodges in brackets (Luke 2021), length of shoreline computed from municipalities with coordinates for beaver observations, estimated population size and number of beavers per 100 km of shoreline in the original area of the Eurasian beaver in Satakunta, in the novel areas (here Etelä-Häme, Pohjois-Häme, Pohjanmaa and Rannikko-Pohjanmaa) and Lapland (see Fig. 1 for areas used in the study) Population size was estimated from lodge numbers and the size of family groups, which was estimated to min 2.8 and max 3.8 or 5.2, for the Eurasian and the North American beaver, respectively (Luke 2021). One family group is assumed to live in each lodge  the number of watershed divides formed by the main drainage basins along the route occurring between the watershed area and the reintroduction site. The length of shoreline within the watershed area was computed to control the effect of available aquatic habitat in the analysis.

b. The effect of artificial barriers on the spread of the Eurasian beaver
For analysis 2b and c (see aims), we divided the SW range into original (Satakunta) and novel areas (Pirkanmaa, Pohjanmaa and Etelä-Pohjanmaa; Fig. 3a). First, (2b) we computed the shortest possible routes from the original area to the novel areas ( Fig. 3b) along the shoreline between the lodges in the original area and lodges and other observations (dams, feeding sites, sightings and DNA samples) in the novel areas. This was done with the ArcMap network analyst closest facility tool. There were 116 known locations of lodges in the original area and 263 observations in the novel areas. We presumed that beavers in the novel areas have dispersed from the closest lodges in the original area. Altogether 257 routes were computed from the observations in the novel areas to the closest lodges in the original area (i.e. 6 novel locations could not be reached along shorelines). We added artificial barriers, i.e. dams and sluices, which have been suggested to inhibit the dispersal of beavers (Halley et al. 2021), to the landscape to study if they lengthened the computed routes between the novel and original areas.

c. Differences in habitat use between the original and novel areas
To compare the habitats in the original and novel areas, buffers of r = 50 m, a typical foraging distance for beavers (Müller-Schwarze 2011), were computed around the observations of the Eurasian beaver. If the buffers intersected with each other, they were dissolved together and treated as one buffer. Altogether there were 469 habitat buffers, from which 184 were in the original area and 206 in the novel areas. In addition, a total of 1939 and 3528 random points were computed along the shoreline in the municipalities with observations in the original area and the novel areas, respectively, so that there was approximately one point per 10 km of total shoreline in a municipality. To avoid overlapping buffers, the minimum distance between random points was set to 500 m. Buffers of r = 50 m were then computed around the random points in the same way as around the beaver observations. The composition of broadleaved, mixed, coniferous and transitional forest, forage combined, and urban and agricultural areas in the buffers was computed. In addition, the abundance of broadleaved forest and the composition of aquatic habitat type were computed.

Statistical analyses
Statistical analyses were computed in XLSTAT Free (Addinsoft SARL 2018) and JMP (JMP ® , Version Pro 15. SAS Institute Inc., Cary, NC, 1989. For comparisons between the beaver species or between the original and novel areas of the Eurasian beaver (study objective 1a-c, and 2b) we used the Mann-Whitney two-sample test. To study whether the watershed divides affected the occurrence of the Eurasian beavers (study objective 2a), we used a generalized linear model with negative binomial distribution. The response variable was the number of beaver observations in the watershed areas, and the length of shoreline, the distance from the release site and the number of watershed divides were used as explanatory variables.
To study the habitat use of the Eurasian beaver in the original and novel areas (study objective 2c) we used binomial logistic regressions. The response variable was the used (1) and the available habitat (0). To analyze whether the habitat use differed between the original and novel areas, we included the data of both areas and added an interaction term between "area" (1 original vs. 2 novel area) and other explanatory variables in the model (similar to Alakoski et al. 2019). The interaction terms were included separately for each explanatory variable. Thus, we could test if the habitat use differed between the areas. Because we were only interested in the interactions between the area and each explanatory variable, we included all environmental variables (composition of aquatic habitat type, broadleaved, mixed, coniferous and transitional forest, forage, agriculture and urban areas, and abundance of broadleaved trees) in the model as explanatory variables. We report only the observed statistically significant results (p < 0.05) for the interaction terms (habitat use of beavers is analyzed in more detail in Alakoski et al. 2019Alakoski et al. , 2020.

a. Distances from the release sites
The North American beaver's distances from the release sites were longer than those of the Eurasian beaver (along shoreline: Mann-Whitney two-sample test: U = 149 937, p < 0.0001; straight-line distance: U = 161 911, p < 0.0001). For the Eurasian beavers, the distances from the original reintroduction site were 109 km along shoreline (median; 24-233 km) and as a straight line 64 km (7-159 km). For the North American beaver, the distances from the closest of the 12 introduction sites were 202 km shoreline (19-558 km) and as a straight line 86 km (2-211 km). Thus, the distance which the beavers spread per year from the release site to the furthest observation point, i.e. the maximum spread rate, was on average 2.8 km per year along shoreline for the Eurasian beaver (in 82 years), and 6.9 km per year for the North American beaver (in 81 years).

b. Sizes of the main drainage basins
The sizes of the main drainage basins within the whole range did not differ between the North American beaver and the Eurasian beaver (median: E 3157 km 2 and NA 3782 km 2 , range: E 478-51 086 km 2 and NA 245-68 446 km 2 , U = 162.5, p = 0.62). The Eurasian beaver observations were from 15 different main drainage basins, four of which were in Lapland (median in the main range in SW Finland was 2147 km 2 ). The North American beaver observations were located in 24 main drainage basins but 42% of the observations were located in one main drainage basin (Fig. 1b).

c. Population growth rates and occupied areas compared to population size
The median population growth rate of the Eurasian beaver has been greater than that of the North American beaver (median: E 1.08 and NA 1.01, U = 199, p = 0.01; Table 2; Fig. 2). If the years 2004-2013 with low reporting activity were excluded for the North American beaver, the median growth rates did not differ between the species (median: E 1.08 and NA 1.02, U = 191, p = 0.43). The population of the North American beaver has grown larger, with the present estimated beaver numbers of 3700-5000 and 10,400-19,400 (Table 2; Luke 2021; for the whole study period the population growth rate was 1.09-1.10 for the Eurasian beaver and 1.10-1.11 for the North American beaver).
For both species, the relative number of occupied areas (number of reporting hunting clubs) has been smaller than the relative population size in the last decades (Fig. 4). However, the number of areas exceeded the population size for the Eurasian beaver temporarily in 1998. For the North American beaver, the number of areas has increased although the population size has remained similar in recent years.

a. The effect of the natural barriers, i.e. watershed divides, on the occurrence of Eurasian beavers
There were slightly more Eurasian beaver observations in the watershed areas that had less watershed divides between them and the release site (estimate = − 0.248, Wald ChiSq = 3.98, p = 0.046; Fig. 5). A shorter distance from the release site (estimate = − 0.021, Wald ChiSq = 59.85, p < 0.0001) and shoreline length within the watershed area (estimate = 0.004, Wald ChiSq = 20.56, p < 0.0001) had clear positive effects on beaver observations in the model. The overall model fit was R 2 = 0.19.

b. The effect of artificial barriers on the spread of the Eurasian beaver
The shortest possible routes along shoreline from the closest locations in the original area (in Satakunta) to the novel areas in SW Finland were 41 km (median, n = 257, 2-138 km). These route lengths were not related to the presence of artificial barriers within the route (U = 30 365, p = 0.11). With dams and sluices as barriers for dispersal, the routes were 42 km (n = 257, 2-164 km).

Discussion
The Eurasian beaver has spread shorter distances than the North American beaver from the original release site. We did not find, however, a reason to suspect that natural or artificial barriers would considerably constrain the range expansion of the native species. There were differences in the habitat use between the original and novel areas of the Eurasian beaver, which may relate to a higher availability

Distribution patterns of the native and non-native beaver species in Finland
The North American beaver has spread significantly further from the release sites than the Eurasian beaver, suggesting that there might be a difference in dispersal behavior. This could relate, for example, to a possibly higher fecundity of the North American beaver (Danilov 1995;Danilov et al. 2011), and to longer dispersal distances. Previous studies suggest that the Eurasian beaver may use smaller territories than the North American beaver (Alakoski et al. 2019), but more information on family groups and territory sizes of both species in Finland is needed. In Sweden, the yearly spread of the Eurasian beaver was similar to the maximum rate of spread based on our data, ca 3 km per year in one area, before a translocation of ten beavers to a new area where the beavers started to spread rapidly (Hartman 1995). In other parts of Sweden, the yearly spread was 12-20 km, i.e. faster than that of the North American beaver's maximum spread in our data (ca 7 km). In the Czech Republic, the growth rate of range diameter was reported as 15-20 km per year, except for the River Elbe where the yearly spread was only 0.8 km (Bartak et al. 2013). These reported differences in the Eurasian beavers' dispersal rates in different areas suggest that environmental factors likely have a major influence on beavers' dispersal (cf. McNew and Woolf 2005;Danilov et al. 2011). Therefore, the differences between the two species' dispersal rate in Finland may not be only due to differences between the species but might be related to the environment. For example, habitat fragmentation and distance between habitat patches, which were not studied here, often constrain the range expansion of species, but may also increase dispersal distances if the habitat patches are connected (Matthysen et al. 1995;Trenham et al. 2001;Bocedi et al. 2014). In addition, e.g. population density and habitat quality can affect the dispersal distances of beavers (e.g. Fustec et al. 2001;Danilov et al. 2011). Furthermore, watershed divides that may decrease dispersal efficiency (Hartman 1995) are more abundant in the Eurasian beaver's SW range than in the North American beaver's range in Finland (Fig. 1b). However, watershed divides were not very strongly related to the beaver abundance on our data (see below). In addition, considering the whole ranges of the two species, the sizes of the main drainage basins did not differ between the ranges of the two species. Because the present Eurasian beaver population started to expand from only one reintroduction site, and from one female, it was probably difficult for a dispersing beaver to find a mate. Thus, low genetic diversity may have affected the population dynamics of the species, but its significance requires further research (Iso-Touru et al. 2020). In addition, problems of small populations, stochasticity in deaths and births, and competition with North American beavers, may have had an impact on the failure of the other Eurasian beaver reintroductions in Finland. Mayer (2017) reported that in Norway, the expansion rate of Eurasian beavers is slow in high densities, and subordinate beavers can stay in their natal territory up to 7 years, possibly replacing their parents as habitat holders. In our data, the median population growth rate has been higher for the Eurasian beaver population that has continuously grown throughout the whole study period, whereas the North American beaver's population first increased to large size, but during the last decades has not increased continuously likely due to a higher hunting pressure compared to that of the Eurasian beaver. This suggests that lower fecundity alone might not explain the smaller range of the Eurasian beaver. However, it should be noted that the reporting activity of hunters affects the population size estimates, and it may have been smaller for the North American beaver than for the Eurasian beaver in Finland. It has also been suggested that the higher hunting pressure in Finland likely explains the slower than expected population growth rate and smaller population size of both beavers, compared to many other European countries (Ermala 1997;Lahti 1997;Hartman 1999). In addition, it cannot be ruled out that Eurasian beavers are hunted when colonizing new areas outside the areas where hunting of the beavers is restricted. However, this does not explain the high density of the species close to the reintroduction site. In any case, and against our prediction, there was no clear difference in the comparison of the relative population size versus the relative number of colonized areas between the two species (Fig. 3). We expected that for the Eurasian beaver with smaller range, the population size would have increased faster than the number of colonized areas. To highlight the comparison between the used and the available habitat, the mean (± SD) for the used habitat is represented first and the mean for the available habitat is on the second row, and the difference between the two means is on the third row. Significant differences between the original and novel areas in SW Finland are in bold font (the analysis with the combined data of areas with an interaction term "area" to study differences between areas, see the main text)

Barriers for spread
One possible reason for the slower spread of the Eurasian beaver could be watershed divides that may constrain beaver expansion (Hartman 1995). Indeed, the numbers of Eurasian beavers correlated negatively with the number of watershed divides in our analysis. However, this relationship was weak, and it seems likely that watershed divides do little to restrict the dispersal of beavers in environments where watershed areas are small and closely located, and drainage basins are not separated by e.g. mountains. In our data, the distance from the reintroduction site explained the numbers of beavers better than the number of watershed divides. In addition, artificial barriers, dams and sluices, did not lengthen the distances between beaver locations considerably, possibly because there are many alternative aquatic routes available. Therefore, it seems unlikely that natural or artificial barriers could constrain the range expansion of Eurasian beavers in SW Finland. This conclusion somewhat differs from the earlier conclusions that the range expansion of beavers mostly occurs within main drainage basins (Hartman 1995;Halley et al. 2012Halley et al. , 2020. For example, Halley et al. (2020) suggest that beaver populations should ideally be considered on a drainage basin scale. Naturally, drainage basins and larger watercourses may have had some effect to the expansion direction of beavers also in Finland. For example, the North American beaver occurs most abundantly in the lake district of Finland, and in areas with many water areas and watercourses (Brommer et al. 2017).

Differences in habitat use in the original and novel areas of the Eurasian beaver
Agriculture and wider streams and rivers were used more, whereas small streams and coniferous forest were found less in the beaver habitats in the novel areas than in the original area in SW Finland (Table 3). Based on our earlier analysis (Alakoski et al. 2019(Alakoski et al. , 2020 these patterns indicate a use of higher quality habitats in the novel areas. In a previous study with partly the same data, both beaver species seemed to prefer 5-20 m streams and avoid small streams (Alakoski et al. 2020). Wider rivers have also been preferred in other countries and they might be more utilized as dispersal routes (Ruys et al. 2011). Perhaps these preferred habitats were more available in the novel areas because of lower overall beaver densities compared to the areas near the reintroduction site with a high density of beavers. Thus, our results indicate that lack of high-quality habitats in the novel areas has not restricted the spread of Eurasian beavers from the original area. Other studies have shown that beavers usually occupy the most optimal habitats first (Halley and Rosell 2002;John et al. 2010;Halley et al. 2013), resulting to distant and irregular colonization patterns (Fustec et al. 2001), but often stay close to their natal territory if there are suitable habitats available (Hartman 1994a).

Conclusions
We conclude that the expansion of the native beaver in Finland seems to be controlled by factors other than those related to the environment. That is, we did not observe indications that habitat quality or barriers for movement could limit the spread of the Eurasian beaver in Finland. Instead, the environment in the original area may have offered enough suitable habitat patches for the Eurasian beaver, enabling the increase in the population density within this area. Perhaps only recently the population has increased close to saturation level, which may explain the recent spread of Eurasian beavers to the novel areas. Movement behavior and population dynamics are possible factors controlling the spread of beavers in Finland, which need further studies.