Introduction
Understanding large carnivore population dynamics is paramount for predicting the outcome of management interventions such as species reintroductions. Large carnivores are important regulators of ecosystems. For example, after being reintroduced to Yellowstone National Park, USA, grey wolves Canis lupus caused a marked change in the Park's ecology and triggered a trophic cascade that restored riparian vegetation and resulted in increases in two other mammal species (Ripple & Beschta, Reference Ripple and Beschta2012). Furthermore, large carnivores control populations of large herbivores and mesocarnivores through predation and intraguild competition (Prugh et al., Reference Prugh, Stoner, Epps, Bean, Ripple, Laliberte and Brashares2009). In the case of large omnivores such as the Asiatic black bear Ursus thibetanus there may be multiple effects, including large carnivore/omnivore habitat exclusion (potential avoidance of U. thibetanus by wild boars in South Korea; Kim, Reference Kim2017; Yoon, Reference Yoon2017) and aiding seed germination and dispersal (Sathyakumar & Viswanath, Reference Sathyakumar and Viswanath2003). However, many carnivores are threatened by habitat degradation and poaching, leading to population declines and range decreases (Ceballos & Ehrlich, Reference Ceballos and Ehrlich2002; Ripple et al., Reference Ripple, Estes, Beschta, Wilmers, Ritchie and Hebblewhite2014). For reintroduced populations, long-term monitoring and management are required for persistence and may also reduce negative interactions with people.
Knowledge of the population dynamics of large carnivores has been used to study carrying capacities (Hayward et al., Reference Hayward, O'Brien and Kerley2007; Kawata, Reference Kawata2008), home range and territory dynamics (Creel, Reference Creel2001), predator–prey interactions (Fryxell et al., Reference Fryxell, Falls, Falls, Brooks, Dix and Strickland1999; Herfindal et al., Reference Herfindal, Linnell, Odden, Nilsen and Andersen2005; Höner et al., Reference Höner, Wachter, East, Runyoro and Hofer2005) and impacts of human activities (Treves & Karanth, Reference Treves and Karanth2003), and these studies can inform management and conservation practices. Improved modelling techniques have provided researchers with insights into population dynamics that include spatial requirements (Gittleman & Harvey, Reference Gittleman and Harvey1982; Seaman & Powell, Reference Seaman and Powell1996), life history traits (Bekoff et al., Reference Bekoff, Daniels and Gittleman1984; Gittleman, Reference Gittleman1986), behaviour (Wrangham et al., Reference Wrangham, Gittleman and Chapman1993) and human influence (Treves & Karanth, Reference Treves and Karanth2003). Ideally, such models incorporate all these factors to the extent they influence the target species, although such studies are rare.
Ursid reintroductions have had varying levels of success (Ursus americanus: Benson & Chamberlain, Reference Benson and Chamberlain2007; Clark, Reference Clark2009; Ursus arctos: Servheen et al., Reference Servheen, Kasworm and Thier1995; Kasworm et al., Reference Kasworm, Proctor, Servheen and Paetkau2007; Zedrosser et al., Reference Zedrosser, Gerstl, Rauer and Band1999; Clark et al., Reference Clark, Huber and Servheen2002; Tremarctos ornatus, Helarctos malayanus: Clark, Reference Clark2009; Ursus arctos horribilis: Lyons et al., Reference Lyons, Gaines, Begley and Singleton2016). The most successful was the reintroduction of U. americanus to Arkansas in 1958–68 (Rogers, Reference Rogers1973; Clark, Reference Clark2009); the original 254 reintroduced bears have grown to a population c. 3,000.
To date, the Republic of Korea (known also as South Korea) is the only country to have reintroduced the Asiatic black bear (Clark Reference Clark2009; Borzée et al., Reference Borzée, Yi, Andersen, Kim, Moon and Kim2019), which is categorized as Vulnerable on the IUCN Red List (Garshelis & Steinmetz, Reference Garshelis and Steinmetz2016), with subspecies in Iran (Almasieh et al., Reference Almasieh, Kaboli and Beier2016) and Taiwan (Hwang et al., Reference Hwang, Garshelis, Wu and Wang2010) categorized as Critically Endangered. The species occurs in 20 countries across Asia, with the largest populations in China, Japan, India and Russia, but is threatened by habitat loss, hunting and bile farming practices (Garshelis & Steinmetz, Reference Garshelis and Steinmetz2016).
Following the extirpation of U. thibetanus and its subsequent reintroduction in the Republic of Korea in 2004, the government needs to determine how to continue the programme. Our study had three objectives. Firstly, to determine the population equilibrium, here defined as the naturally levelling population size (Hui, Reference Hui2006), of Jiri Mountain National Park (the site of the original reintroduction programme) and the immediate surrounding area for U. thibetanus. Secondly, to determine the species’ population trajectory throughout South Korea both with and without a second reintroduction programme, and thus whether the population could become viable and whether large-scale reintroduction would facilitate this. Thirdly, to delineate habitat patches and potential dispersal corridors for future protected area designation and habitat restoration. For the purposes of this study, we use the term population equilibrium in lieu of carrying capacity, as the former is more applicable to our model outputs. As far as we are aware, multi-method population modelling has not been conducted for any reintroduced ursid population, and such modelling is limited to a single case that assumed an existing population rather than a gradual reintroduction (Lyons et al., Reference Lyons, Gaines, Begley and Singleton2016, Reference Lyons, Gaines, Singleton, Kasworm, Proctor and Begley2018).
Study area
The study area was the U. thibetanus reintroduction site of Jiri Mountain National Park (Fig. 1), which lies at an altitude of 184–1,915 m and covers 483 km2 in three provinces in southern South Korea (Korea National Parks Service, 2018). There are two other nearby National Parks, Deogyu and Gaya Mountain National Parks (Fig. 1). These mountains collectively comprise the southern part of the Sobaek Mountain Range. Land development is not allowed within National Parks in South Korea, although existing structures can remain upon park creation. These landscapes are largely maintained for landscape conservation. There are some roads and towns between the three National Parks that could be potential barriers to dispersing individuals. Jiri Mountain National Park receives c. 3 million visitors annually (Korea National Parks Service, 2017).
The locations of mountain ranges and National Parks in the Republic of Korea. The reintroduction site is in Jiri Mountain National Park, and we propose a second reintroduction in Seorak Mountain National Park.
Methods
Species introduction
Ursus thibetanus was historically native to the Republic of Korea and the rest of the Korean peninsula. However, a combination of poaching for medicinal purposes, eradication during the Japanese occupation (1910–1945), and deforestation during the Korean War (1950–1953) and in the later 1950s led to the species’ extirpation from the Republic of Korea by the 1990s (Han, Reference Han2006; Kim et al., Reference Kim, Hong, Min, Kim, Kim and Voloshina2011). Since then, the government has sought to re-establish the species, specifically through a reintroduction programme in Jiri Mountain National Park, which is purported to have suitable habitat for the species (Lee & Jeong, Reference Lee and Jeong2009; Fig. 1). This programme started in 2004 using rehabilitated individuals from the Russian Federation, the Democratic People's Republic of Korea, and the People's Republic of China (Lee & Jeong, Reference Lee and Jeong2009; Kim et al., Reference Kim, Hong, Min, Kim, Kim and Voloshina2011). With a population of 56 individuals following the birth of 11 cubs born during the winter of 2017–2018 (Choi, Reference Choi2018), the reintroduction programme has been largely considered a success, with continued population growth and stability, and dispersal (Ok, Reference Ok2018; Borzée et al., Reference Borzée, Yi, Andersen, Kim, Moon and Kim2019). The feasibility of a second reintroduction, in either Seorak or Odae Mountain National Park, has been examined (Kim, Reference Kim2017).
Data collection
This study used data collected on U. thibetanus individuals reintroduced to Jiri Mountain National Park. All reintroduced bears (a total of 21 males and 18 females, beginning with an initial reintroduction of three juveniles of each sex in 2004) were released with radio-tracking collars (M3620 or M2320; Advanced Telemetry Systems, Isanti, USA) and/or individual ear tags (Allflex, Madison, USA). Global positioning system data were collected during 2004–2016. Locations of the bears were detected using dial tone receptors (ICOM or IC-R10; Universal Radio, Worthington, USA) with Yagi antennae (ATS, Isanti, USA). Individuals were tracked in the field and their locations were triangulated when GPS acquisition was difficult. Data were collected using triangulation, with azimuth of the strongest signal from at least two or three locations (Springer, Reference Springer1979).
Using both GPS and telemetry data, we calculated the annual 50% kernel home range (Worton, Reference Worton1989) of each individual. We grouped these ranges into classes of male or female and juvenile, subadult or adult, then compared the ranges of sex and age classes using one-way analysis of variance (ANOVA) to determine any overall significance and Tukey honest significant difference (HSD) tests to determine any significance between pairs. These home ranges and analyses were used to calculate resource targets for population simulations and habitat core sizes for connectivity mapping. Here, the definitions of resource target and habitat core follow the definitions of the software used for analyses, HexSim and Linkage Mapper (see below), respectively, where resource target refers to the sum of resource scores of occupied hexagons (Schumaker & Brookes, Reference Schumaker and Brookes2018) and habitat core refers to habitat patch.
Habitat suitability model
To model population equilibrium, we first developed a habitat suitability model of the spatial availability of resources for U. thibetanus. To do this, we created a resource selection function that models availability versus use of environmental resources and facilitates the quantification of resource availability and probability of use, which makes it particularly suited for wildlife management (Boyce & McDonald, Reference Boyce and McDonald1999). We built the model using seven remotely sensed variables that influence resource selection by U. americanus horribilis (Proctor et al., Reference Proctor, Nielsen, Kasworm, Servheen, Radandt, Machutchon and Boyce2015). The distributions and spatial predictions of habitat suitability for this subspecies are influenced by environmental variables similar to those influencing the distribution of the nominate species U. americanus (Apps et al., Reference Apps, McLellan and Woods2006), which is ecologically similar to U. thibetanus. These variables were incorporated into our model based on data availability and relevance to the study area; i.e. elevation, ruggedness, roads, alpine, riparian, normalized difference vegetation index (NDVI), greenness and wetness variables (Supplementary Table 1). To ensure NDVI and greenness accurately represented forested areas and not high productivity agricultural land, we conducted a visual analysis between these vegetation indices and high resolution satellite imagery. We removed variables with high multicollinearity (Pearson's r > 0.8).
We used a reverse stepwise approach to incorporate variables that increased the pseudo R 2 value of the model by ≥ 0.05 (Bursac et al., Reference Bursac, Gauss, Williams and Hosmer2008). For presence, we used a random subset of 5,000 GPS points (selected using the sample function in R 3.4.1; R Core Team, 2017) collected for U. thibetanus in the Park during 2004–2015. Points were randomly selected from all sex and age classes, to capture ecological variation for the population as a whole. For absence, we used the Create random points tool in ArcMap 10.5 (Esri, Redlands, USA) to create 5,000 random points within and up to 5 km from Jiri Mountain National Park (the area containing 99.9% of tracking locations), and located at least 400 m from any presence point, matching the accuracy of the GPS tracking devices.
We then projected the resulting resource selection model into ArcMap and scaled the raster values to a range of 0–1. In the final spatial model, areas with spring NDVI < 0.65 were removed as 99% of presence points were located in areas with May NDVI > 0.65. May NDVI was used as it is the peak of spring vegetation productivity and therefore represents optimal vegetation resources. The final habitat suitability model was categorized into five classes of low–high resource suitability index (Proctor et al., Reference Proctor, Nielsen, Kasworm, Servheen, Radandt, Machutchon and Boyce2015): < 0.2; 0.2 < 0.35; 0.35 < 0.5; 0.5 < 0.67; ≥ 0.67).
Simulation parameters and scenarios
We predicted population trajectories using HexSim 4.0.9, an individual-based, spatially explicit model framework for simulating populations over time steps (Schumaker & Brookes, Reference Schumaker and Brookes2018), used by Lyons et al. (Reference Lyons, Gaines, Begley and Singleton2016) to model carrying capacity for U. americanus horribilis in the North Cascades region of Washington State, USA. HexSim allows the input of a habitat suitability spatial component along with life history traits, for modelling population dynamics.
We first input our classified habitat suitability model in hexagonal cells of 500 m diameter. After preliminary simulations, we identified resource target (a parameter – the level of resource accumulation an individual of a certain class will attempt to attain; Schumaker & Brookes, Reference Schumaker and Brookes2018) as a main driver of population size. To determine resource targets (i.e. resource requirements per individual), we used home range to calculate the sum of categorized habitat suitability model values from a layer of hexagons congruent to our habitat suitability model, to represent the habitat quality for each cell. Next, we classified the population by sex and age class (juvenile: < 2 years; subadult: 2–4 years; adult: > 4 years; Yamanaka et al., Reference Yamanaka, Yamauchi, Tsujimoto, Mizoguchi, Oi and Sawada2011); i.e. a total of six classes. Time steps in HexSim are between age classes: thus for our simulations one time step was 2 years. Finally, we used the mean resource requirement for each class (Supplementary Table 2) as our computed resource target parameter in HexSim. To obtain a range of population equilibria, we used one standard error above and below the mean resource requirement parameter. In addition to resource target, we included birth, survival, grouping and dispersal parameters adapted from data collected for U. thibetanus in Jiri Mountain National Park (Supplementary Table 3). Motorways, roads and rivers were added, as barriers (Lyons et al., Reference Lyons, Gaines, Begley and Singleton2016). To mimic the Jiri Mountain reintroduction programme, we started with an initial population of six and then introduced new individuals at time steps corresponding to the reintroduction programme (i.e. 1–8 individuals released every 1–2 years; Song, Reference Song2017).
To model U. thibetanus population trajectories over time and obtain a value for the population equilibrium, we simulated three scenarios using HexSim. The first scenario restricted the population to Jiri Mountain National Park and the surrounding area within 5 km. This addressed our first objective, related to small-scale reintroduction and management that has already been conducted in the Park. The second scenario allowed the population to expand beyond Jiri Mountain National Park to other parks along the Sobaek Mountain Range (A + B in Fig. 2), with no further reintroductions. This addressed our second objective, related to the potential for future large-scale reintroduction. The third scenario, also related to our second objective, allowed the population to expand outside the Park and added a second reintroduction to Seorak Mountain National Park in the northern Taebaek Mountains (Fig. 1). This theoretical second reintroduction replicated the structure and timing of the reintroduction in Jiri Mountain National Park and occurred 12 time steps after the beginning of the initial 2004 reintroduction, approximately 2028, to allow for logistic planning.
Habitat cores (see text for details) for Ursus thibetanus, within the Baekdudaegan Mountains Reserve and other protected areas. The Southern Sobaek Mountain Range (A), comprising Jiri Mountain, Deogyu Mountain, and Gaya Mountain National Parks, is where the population is expected to persist if there are no further reintroductions. Part of the Baekdudaegan Mountains Reserve (B) will act as an important corridor for connectivity between the southern Sobaek Mountain Range and the Taebaek Mountains. The northern Taebaek Mountain Range (C), with extensive habitat cores is important potential habitat for the species. The southern Taebaek Mountain range is composed of small, fragmented habitat cores, and is probably not important habitat for the species.
Dispersal corridors
For our final objective, to determine potential dispersal corridors and connectivity barriers for U. thibetanus, we used Linkage Mapper 1.1.1 in ArcMap. Linkage Mapper requires input of habitat cores, or areas of high habitat suitability, and a landscape resistance raster. The habitat suitability model was thresholded (i.e. transformed to create a binary output of presence and absence; Liu et al., Reference Liu, White and Newell2013), at the maximum true scale statistic (Allouche et al., Reference Allouche, Tsoar and Kadmon2006) using the coords function of the pROC package (Robin et al., Reference Robin, Turck, Hainard, Tiberti, Lisacek, Sanchez and Müller2011) to create habitat cores for connectivity mapping. We then selected cores > 309 ha (3.09 km2), which is the daily foraging requirement of an adult female U. thibetanus (Sakamoto et al., Reference Sakamoto, Kunisaki, Sawaguchi, Aoi, Harashina and Deguchi2009; Proctor et al., Reference Proctor, Nielsen, Kasworm, Servheen, Radandt, Machutchon and Boyce2015). Additionally, we removed small patches (< 24 km2; mean monthly home range of an adult female, from our data) that were > 40 km away (half of maximum dispersal distance; Borzée et al., Reference Borzée, Yi, Andersen, Kim, Moon and Kim2019) from large patches (> 46 km2; mean annual home range of an adult female, from our data; methods from Lyons et al., Reference Lyons, Gaines, Begley and Singleton2016). For landscape resistance, we created a layer with 200 × 200 m resolution by combining an urbanization index (1 – NDVI) and the inverse of our habitat suitability model layer. This urbanization index represented resistance as a result of human activity and development, and the inverse of suitability represented landscape traits that are unsuitable for bear movement (Proctor et al., Reference Proctor, Nielsen, Kasworm, Servheen, Radandt, Machutchon and Boyce2015). We also used the barrier mapper tool in Linkage Mapper to identify barriers to dispersal within the landscape matrix.
As roads are also potential barriers to connectivity for ursids (Merrill et al., Reference Merrill, Mattson, Wright and Quigley1999; Doko et al., Reference Doko, Fukui, Kooiman, Toxopeus, Ichinose, Chen and Skidmore2011) but may not be well represented in the coarser resolution of our raster data (200 × 200 m for connectivity mapping in Linkage Mapper), we also analysed highways as shapefiles crossing core habitats and dispersal corridors (Proctor et al., Reference Proctor, Nielsen, Kasworm, Servheen, Radandt, Machutchon and Boyce2015). For this analysis, dispersal corridors were thresholded at the maximum value for a dispersal route used by a previously dispersed individual (Borzée et al., Reference Borzée, Yi, Andersen, Kim, Moon and Kim2019). Habitat cores and dispersal corridors were intersected separately with motorways, trunk roads, and primary and secondary highways, to identify barriers to connectivity. Additionally, road tunnels and wildlife crossings were intersected within the habitat core and dispersal corridor areas to identify any potential existing aids to wildlife dispersal. The locations of roads were obtained from OpenStreetMap (OpenStreetMap, 2019) and wildlife crossings from the National Institute of Ecology Ecological Pathway Map Service (National Institute of Ecology, 2019).
Results
Tracking and home ranges
A total of 33,924 tracking points were collected from 30 male and 26 female U. thibetanus during 2004–2016. The mean number of locations per individual was 606, with males averaging 529 locations and females 695 locations. Of these, 83.9% were within Jiri Mountain National Park and 99.9% were within 5 km of the Park. Tracking points averaged 2,827 per month, with the highest number of points during November (3,411) and the lowest number during February (2,412). The mean number of points per season was 8,481, with no significant difference in tracking points between seasons (Tukey HSD; P = 0.498).
Mean annual home range sizes as determined by a 50% kernel were 45.7 ± 37.2 km2 and 121.0 ± 79.3 km2 for adult females and males, respectively (Fig. 3). The one-way ANOVA indicated significant difference in home range size among groups (P = 0.00139) and the Tukey honest significant difference test revealed that home range size was significantly different between adult males and females (P = 0.0009) and between adult and juvenile males (P = 0.0145). Differences were not significant between all other pairs tested (P > 0.05).
Box and whisker plots (with outliers) of core annual home range sizes of U. thibetanus in Jiri Mountain National Park, Republic of Korea, grouped by sex (male or female) and age class (juvenile, subadult, or adult) for a total of six age/sex groups. Mean home range sizes (± SD) are indicated for each group.
Habitat suitability model
Our final resource selection function incorporated elevation, greenness and alpine variables (r = 0.583, P < 0.0001; Table 1). Elevation and greenness were positively correlated with U. thibetanus resource selection whereas alpine habitat was negatively correlated, meaning higher elevations and areas with more green vegetation were more likely to be selected by the bears, whereas high alpine areas were generally avoided. Both greenness and NDVI adequately differentiated between forest and agricultural areas, but NDVI was excluded in this model because of high collinearity with greenness. Presence points with normalized resource selection function values above the maximum true scale statistic threshold accounted for 88.7% of all presence points, indicating good model fit. The area of suitability above this threshold covered 6,410 km2, or 6.4% of the country, not including small or highly isolated habitat (Supplementary Fig. 1). The true scale statistic for the model was 0.7486, which is considered good for predictions of presence versus absence (Allouche et al., Reference Allouche, Tsoar and Kadmon2006; Preau et al., Reference Preau, Trochet, Bertrand and Isselin-Nondereu2018).
The final resource selection function model for Ursus thibetanus in the Republic of Korea (R 2 = 0.583, P < 0.0001).
Simulation outcomes
For Jiri Mountain National Park and surrounding suitable habitat within 5 km of the Park, our simulations resulted in a population equilibrium of 55–77 individuals at high to low resource targets (i.e. the resources required per individual), with the mean resource target yielding a population equilibrium of 64 individuals (Fig. 4a). This number was reached in 20–30 time steps, equivalent to 40–60 years, in averaged simulations. When individuals were allowed to disperse beyond Jiri Mountain National Park, the population equilibrium was 188–281, with a mean of 230 (Fig. 4b). However, the population in this scenario was largely constrained to the Sobaek Mountain Range as populations generally did not establish in the Taebaek Mountain Range. When a second reintroduction to Seorak Mountain National Park was added, the population equilibrium became 1,123–1,438, with a mean of 1,266 across all suitable habitat in the country (Fig. 4c). These population equilibria translated to population densities of 0.107–0.224 bears per km2, a value that increased with geographical range available and with a second reintroduction event. As a general trend, population densities increased along with the total area of habitat cores that we identified when conducting the analyses for connectivity mapping.
Means of 100 simulation results for the population of U. thibetanus in Jiri Mountain National Park and its immediate surrounding areas (scenario A), for the entire Republic of Korea without a second reintroduction (B), and for the entire Republic of Korea with a second reintroduction (C), with low, average and high resource requirements over 200 time steps (c. 400 years). Scenario A represents containment in Jiri Mountain National Park and the immediate vicinity and also indicates the population equilibrium for the Park. Scenario B represents the population trajectory if no further reintroductions are undertaken. Scenario C represents the population trajectory if a second reintroduction is implemented in Seorak Mountain National Park in c. 2028 (time step 12).
Habitat cores and dispersal corridors
We found 164 suitable habitat cores, with a total area of 11,038 km2 and a mean core size of 67.3 km2 (Fig. 2). These cores were mostly located in the vicinity of the Baekdudaegan Mountains Reserve and adjacent national parks. Between these cores were 328 links representing least-cost pathways averaging 11.1 km and totalling 3,649 km (Fig. 5, Table 2). There were 744 highways (motorways, trunk roads, primary and secondary roads) intersecting habitat cores and thresholded low resistance pathways (dispersal corridors), totalling 4,266 km (Fig. 6, Table 3). However, there were also 247 road tunnels along these roads with a total length of 192 km as well as 100 wildlife crossings. Within the southern Sobaek Mountain Range (Fig. 2a) 204 roads crossed habitat cores and dispersal corridors, totalling 1,192 km. Along the Baekdudaegan Mountains Reserve between the southern Sobaek Mountain Range and the Taebaek Mountain Range (Fig. 2b), there were 45 roads over 292 km, and within the Taebaek Mountain Range (Fig. 2c), there were 220 roads over 1,328 km.
Habitat cores and least-resistance pathways (i.e. dispersal corridors), for U. thibetanus in the entire Republic of Korea (left), and in Jiri, Deogyu and Gaya Mountain National Parks (right), with the dispersal route taken by a subadult male (Borzée et al., Reference Borzée, Yi, Andersen, Kim, Moon and Kim2019) providing a validation of the connectivity model.
Habitat cores and thresholded dispersal corridors for U. thibetanus overlain with motorways, trunk roads, primary and secondary roads, road tunnels and wildlife crossings (National Institute of Ecology, 2019).
Area of habitat cores and length of least-cost paths for U. thibetanus in the Republic of Korea based on data collected during 2004–2016.
Length and number of roads intersecting habitat cores and least-cost paths for U. thibetanus in the Republic of Korea based on data collected during 2004–2016.
Discussion
Our findings show that the reintroduction programme is a success, as viable populations are established in all suitable core areas in Jiri Mountain National Park, and that the population is close to our simulated equilibrium without needing further supplementation. However, our results also indicate there will be a need for a second reintroduction and continued conservation efforts for U. thibetanus to restore the species across its original range in the Republic of Korea. We have also identified core habitat areas and dispersal corridors required to preserve extant habitat and improve connectivity for the reintroduced population. This study also highlights the usefulness of integrated modelling techniques in predicting the need for further reintroduction and landscape management.
Population simulations for U. thibetanus in Jiri Mountain National Park indicate that the reintroduced population there will reach the modelled population equilibrium relatively soon (Park, Reference Park2001), with dispersals into surrounding national parks and wild, mountainous areas already occurring (Borzée et al., Reference Borzée, Yi, Andersen, Kim, Moon and Kim2019). As populations fluctuate over time and eventually stabilize around the equilibrium, we expect the population within Jiri Mountain National Park to fluctuate between low (55) and high (77) predictions and stabilize around 64 individuals. As the population in the Park was already 56 at the time of writing (Choi, Reference Choi2018), Park resources can be focused on facilitation and management rather than further population supplementation.
Simulation of our second scenario shows that if there are no further reintroductions elsewhere, the population will be restricted to the southern Sobaek Mountain Range. The mean population trajectory increases from 230 to 1,266 individuals (450% increase) if a second reintroduction is implemented in the Taebaek Mountain Range. Although it is possible that individuals or small groups from the reintroduced population in Jiri Mountain National Park may disperse via the Baekdudaegan Mountains Reserve (Kim et al., Reference Kim, Lee, Lee and Jang2016), it is unlikely that a population will become established in the northern Taebaek Mountain Range. The chance of such establishment is 0.3% within the next 200 years and 8.0% within the next 400 years (Supplementary Fig. 3), probably because of the inadequate habitat matrix in the Baekdudaegan Mountains Reserve for the long-range dispersal required for movement. Our results thus support the need for continued management and further reintroductions of U. thibetanus across the country to restore the species. We therefore recommend a second reintroduction within the Taebaek Mountain Range, ideally in Seorak Mountain National Park, for a complete restoration of U. thibetanus in the Republic of Korea. For comparison, the successful reintroduction of U. americanus to Arkansas has resulted in a population of c. 3,000 individuals and is considered viable (Rogers, Reference Rogers1973; Clark, Reference Clark2009).
For important habitat cores and dispersal corridors, we delineated three areas that will require greater protection and improved connectivity: the southern Sobaek Mountain Range (the species’ current location), the northern Taebaek Mountain Range (for a second reintroduction), and the Baekdudaegan Mountains Reserve (for connecting these two mountain ranges). Improved connectivity will benefit other species, as the Baekdudaegan Mountains Reserve is habitat for almost all terrestrial mammal species in the Republic of Korea, including the Siberian musk deer Moschus moschiferus and the long-tailed goral Naemorhedus caudatus (Jo et al., Reference Jo, Baccus and Koprowski2018). Both species are categorized as Vulnerable on the IUCN Red List (Nyambayar et al., Reference Nyambayar, Mix and Tsytsulina2015; Bragina et al., Reference Bragina, Kim, Zaumyslova, Park and Lee2020), and in the Republic of Korea M. moschiferus is categorized as Critically Endangered and N. caudatus as Endangered (Jo et al., Reference Jo, Baccus and Koprowski2018). Improvements to connectivity in the form of wildlife crossings should be made not only along motorways but also along smaller roads, as previous studies have shown roadkills occur more frequently along roads with lower traffic volume (Kim et al., Reference Kim, Serret, Clauzel, Andersen and Jang2019).
We have assumed that the environmental variables we used in the modelling are suitable for predicting habitat suitability for U. thibetanus outside Jiri Mountain National Park. The Korea National Parks Service have measured acorn production as a food resource for the species (National Parks Service, pers. comm., 2017) but it would be impossible to measure this country-wide and therefore we used greenness as a proxy for food resources, potentially overestimating the availability of suitable habitat. Furthermore, it is difficult to determine habitat quality outside protected areas, and differing land management practices or negative human–bear interactions may limit potential populations. A genetic bottleneck could also potentially influence the viability of the reintroduced population (Kim et al., Reference Kim, Hong, Min, Kim, Kim and Voloshina2011), but this variable is difficult to introduce in population simulations.
Despite these uncertainties, our estimate of the density of U. thibetanus (0.107–0.224 individuals per km2) is similar to that of populations elsewhere: 0.219 per km2 in Japan (Doko et al., Reference Doko, Fukui, Kooiman, Toxopeus, Ichinose, Chen and Skidmore2011), 0.08–0.29 per km2 in Thailand (Ngoprasert et al., Reference Ngoprasert, Reed, Steinmetz and Gale2012), and 0.17 per km2 in India (Sharma et al., Reference Sharma, Charoo and Sathyakumar2015). Additionally, our estimates fall within the range of densities reported for U. americanus in the USA (0.02–1.43 per km2; Boersen et al., Reference Boersen, Clark and King2003). The similarities between our density estimates and those of other studies suggests that our models were able to estimate the population equilibria in Jiri Mountain National Park correctly. In addition to population densities, our results also match previous studies of home range size, in that the larger home range size of adult males is similar to home range patterns of U. americanus (Pacas & Paquet, Reference Pacas and Paquet1994; Wooding & Hardisky, Reference Wooding and Hardisky1994; Koehler & Pierce, Reference Koehler and Pierce2003), which could be a result of males patrolling larger areas in search of females during mating season (Powell et al., Reference Powell, Zimmerman and Seaman1997; Koehler & Pierce, Reference Koehler and Pierce2003). Moreover, the actual dispersal pathway taken by an individual U. thibetanus that dispersed from Jiri Mountain National Park (Borzée et al., Reference Borzée, Yi, Andersen, Kim, Moon and Kim2019) provides a real-life validation of our connectivity mapping.
We posit that additional reintroductions to Jiri Mountain National Park are not necessary for that population's continued viability if current conditions are maintained and in the absence of stochastic events. Rather, efforts should be directed at continued monitoring and maintenance of the current population, improved connectivity to Baekdudaegan Mountains Reserve, and a second reintroduction programme, in the northern Taebaek Mountain Range. Although Baekdudaegan Mountains Reserve seems a sufficient habitat corridor for connecting populations in the southern Sobaek and northern Taebaek Mountain Ranges, this corridor would not easily facilitate establishment from the southern Sobaek Mountain Range to the northern Taebaek Mountain Range.
Our results indicate that a second reintroduction is required elsewhere in the country to ensure the long-term viability of this reintroduced population. Our modelling also highlights the viability of the habitat matrix around the reintroduction site and the need to improve connectivity among habitat cores in the area. Our methodology could be of value for other projects seeking to determine locations and schedules for reintroducing other large carnivore species.
Acknowledgements
This work was supported by Jirisan National Park (grant 20170453243-00), and by research grants from the Rural Development Administration of Korea (PJ01228503) and the National Research Foundation of Korea (2017R1A2B2003579).
Author contributions
Design of methods, data analysis, writing lead: DA; project conception, writing: YY, AB, KK; project conception, data collection: KM, JK, TK; project conception, project supervision: YJ; revisions: all authors.
Conflicts of interest
None.
Ethical standards
Handling of focal species was by the Korean National Parks Service in accordance with national laws and permitting requirements within the Republic of Korea, and this research otherwise abided by the Oryx guidelines on ethical standards.
Open access
