A spatially explicit analytical framework to assess wildfire risks on brown bear habitat and corridors in conservation areas

Background: Humans have altered fire regimes across ecosystems due to climate change, land use change, and increasing ignition. Unprecedented shifts in fire regimes affect animals and contribute to habitat displacement, reduced movement, and increased mortality risk. Mitigating these effects require the identification of habitats that are susceptible to wildfires. We designed an analytical framework that incorporates fire risk mapping with species distribution modeling to identify key habitats of Ursus arctos with high probability of fire in Iran. We applied the random forest algorithm for fire risk mapping. We also modeled brown bear habitats and predicted connectivity between them using species distribution models and connectivity analysis, respectively. Finally, the fire risk map, critical habitats, and corridors were overlaid to spatially identify habitats and corridors that are at high risk of fire. Results: We identified 17 critical habitats with 5245 km 2 of corridors connecting them, 40.06% and 11.34% of which are covered by conservation areas, respectively. Our analysis showed that 35.65% of key habitats and 23.56% of corridors are at high risk of fire. Conclusions: Since bears habitat in this semi-arid landscape rely on forests at higher altitudes, it is likely that shifting fire regimes due to changing climate and land use modifications reduce the extent of habitats in the future. While it is not well known how fire affects bears, identifying its key habitat where wildfires are likely to occur is the first step to manage potential impacts from increasing wildfires on this species.


Background
Natural fire regimes maintain the health of many ecosystems (Pausas and Keeley 2019) but anthropogenic climate change and land use change has altered these regimes worldwide in an accelerated rate, which poses a threat on biodiversity (Dennis et al. 2001). Forest fires can directly or indirectly impact wildlife species by increasing mortality of individuals (Clark et al. 2011), changing habitat extent and quality (Richardson et al. 2007;Wan et al. 2019b), decreasing population size (Kramp et al. 1983), changing movement patterns (Nimmo et al. 2019), causing outbreaks of diseases (Ecke et al. 2019;Albery et al. 2021), and decreasing food and shelter availability (Banks et al. 2011). The magnitude of these impacts is, however, highly dependent on the spatial and temporal scales of fires (Hobbs and Huenneke 1992;Pickett and White 2013;Wan et al. 2020), as well as the biological and ecological characteristics of species in adapting to fires (Chambers et al. 2019).
Large carnivores that depend on forest cover as their primary habitat might be particularly vulnerable to altered fire regimes. While the effect of fire on these species and their habitat is generally perceived as negative, the actual relationship is complex and ambiguous for many species. For example, prior wildfire studies on bears have produced mixed results. Wildfires can negatively impact adult sex ratio, recruitment of cubs, and population size (Cunningham and Ballard (2004). Wildfires can also destabilize and degrade habitats surrounding den sites by altering vegetation composition and increasing the depth of the active soil layer (Richardson et al. 2007). However, wildfires can sometimes lead to increased abundance and availability of and easier access to food sources (Blanchard and Knight 1990;Hamer 1999).
Because wildfires are more likely to have direct impacts on habitat than on animals, forest and species management plans can benefit from having information that identifies critical habitat susceptible to fire (e.g., Meng et al. 2016) as well as post-fire habitat change . Species distribution models have been increasingly used to provide some of this information, including the usage in predicting critical habitats, connectivity, and susceptibility of landscape to disturbance events Kaszta et al. 2020). However, one important, but rare, application is to streamline these models into a systematic spatial analytical framework for assessing fire risks on species critical habitat and corridors. Such an analytical framework would provide tremendous values in managing wildfire risks on species habitat experiencing changing fire regimes.
The brown bear (Ursus arctos; Linnaeus 1758) is the largest carnivore species of Iran, inhabiting forest and mountainous ecosystems along the Alborz and Zagros mountain ranges (Gutleb and Ziaie 1999). The Zagros Mountains forest steppe is the largest forest region in Iran, representing one of the most important habitats for the brown bears in the country. About 15.36% of this forest is located in the Fars province. The forest has been increasingly impacted by wildfires due to humans and changing climate. For example, the record for the highest number of annual forest fire incidence were broken in 2015, which burned 937 ha of forest lands (The Center for Strategic Studies, 2017). In 2020, over 1700 ha of rangeland and forests were burned by wildfires, representing one of the biggest records in burned area in this province. Recently, in 2021, 150 forest fire incidences have been reported during the hot months. Such a trend in fire activity raises concerns for natural habitats management.
While it is not well known how forest fires affect brown bears in Iran, increasing forest fires due to climate and land use changes in this region may reduce the extent of brown bear habitat and its range in the future. Therefore, identifying key habitats where forest fires are likely to occur is the first steps to mitigate the negative effects or increase the positive impacts of such events. Although there are critiques of umbrella and single-species conservation approaches (Simberloff 1998;Andelman and Fagan 2000), brown bears may serve as an "umbrella species" for structuring site-based conservation of forest habitats (Khosravi et al. 2019) because of their large area requirements, charismatic appearance, and low densities.
Using the brown bear in Iran as an example, we demonstrate a spatially explicit analytical framework that adopts a series of models for assessing wildfire risks for wildlife species. This framework contains four key components: (1) identifying critical habitats, (2) evaluating connectivity between critical habitats, (3) produce fire risk map, and (4) identifying habitats and corridors associated with high risk of wildfire.

Methods
The methodological framework Our case study area is located in the southwestern Iran, Fars province, with an extent of 125,857 km 2 . It is characterized by semi-arid climate conditions. The mean temperature and annual rainfall are 16.8°C and 315 mm, respectively. This area contains several conservation areas (a total of 23,819 km 2 ) in four main categories including national parks (NPs), wildlife refuges (WRs), protected areas (PAs), and no-hunting areas (NHAs). A number of large species of conservation importance such as brown bear, Persian leopard (Panthera pardus), caracal (caracal caracal), gray wolf (Canis lupus), mouflon (Ovis gmelini), goitered gazelle (Gazella subgutturosa), and wild goat (Capra aegagrus) are found in the region. Forests (12.49%), rangelands (55.21%), and farmlands (23.02%) cover 90.72% of the study area (Fig. 1). The status of brown bears in this part of Iran is of great concern. This is the southernmost extent of brown bear's range in the world. The bear population here is genetically distinct from other populations in Iran (Ashrafzadeh et al. 2016). They mostly rely on water sources in densely forested areas at higher altitudes where human access is limited (Ansari and Ghoddousi 2018). The forests of the Zagros Mountains are threatened by droughts, extensive land-use changes, and wildfires (Sagheb-Talebi et al. 2014). Therefore, conserving this population is a pressing priority.
We carried out a spatially explicit analytical framework that consists of four key steps, including (1) predicting the critical habitats of the species, (2) evaluating connectivity between the predicted critical habitats, (3) producing a fire risk map across the study area, and (4) synthesizing the results from those three steps to map and identify critical habitats and corridors associated with high wildfire risks (Fig. 2). The steps undertaken in the study are explained in detail below:

Species potential distribution and critical habitats
In 2015-2019, we collected opportunistic observations of brown bear signs (i.e., scats and tracks; N = 58) and direct observations of brown bears by the game guards and wildlife experts or human-caused mortality due to conflict with people (N = 113). To reduce spatial autocorrelation and sampling bias in presence points, we used a spatial filter such that each locality was at least 400-m from others (Falcucci et al. 2009). After the spatial filtering, 112 presence points remained for the modeling.
We divided the study area into two regions based on the field surveys and expert opinion: the core area covering the main current range of brown bears, with the highest population densities, and peripheral areas with occasional observations (Fig. 1). Although the peripheral regions include only sporadic occurrence of bears and might not cover the current distribution range of the species, these regions may provide suitable habitats for bears to expand their range. The model was calibrated across the core area and then projected over the peripheral region. To calibrate the species distribution model, we used seven land use/land cover types, five topographic features, two anthropogenic factors, two climatic variables, and normalized difference vegetation index (NDVI) based on literature review (Ansari and Ghoddousi 2018;Farashi 2018;Almasieh et al. 2019;Maiorano et al. 2019). All variables were calculated in a 50-m × 50-m cell size. A detailed description on the preparation of the predictor layers is given in Table S1 in Online Resource 1.
Dispersing and settled individuals may be affected by different fine-and broad-scale patterns (Mateo Sanchez et al. 2014). Therefore, we fit models for predicting critical habitats and connectivity separately. We used coarse-scale models to provide a second-order (Johnson 1980) prediction of the species' potential distribution. For this aim, we used a moving window with a radius of Fig. 2 The spatially explicit analytical framework to assess wildfire risks on brown bear habitat and corridors 450-m around each original pixel (i.e., 50-m) across the study area (Falcucci et al. 2009;Maiorano et al. 2019).
For each categorical variable (e.g., land cover classes), we calculated the percentage of cover within a moving window of 450-m using FRAGSTATS (McGarigal et al. 2012). For each continuous variable (e.g., topographic measures), the mean value inside a given radius was calculated using neighborhood statistics in ArcGIS 10.2.2 (ESRI, 2015). We applied three filtering steps to reduce the number of variables for the final species distribution model. First, we removed predictors that occurred at < 10% of presence localities because of lack of representation . Second, predictors with variance inflation factor (VIF) ≥ 2 were removed (Maiorano et al. 2019). Finally, we assessed multicollinearity between each pair of variables and for correlated pairs (R> 0.7) only the variable with the lower VIF was retained (Wasserman et al. 2012;Dormann et al. 2013). Following the reduction process, we used fourteen variables in the final model (Table S1 in Online Resource 1).
To predict critical habitats, we implemented the MaxEnt model using the "dismo" R package (Hijmans et al. 2017). The model was calibrated according to the presence localities within the core area. We ran MaxEnt using a 10-fold cross-validation and calibrated models using 75% of occurrence points as training data, and the remaining 25% of data for evaluation model predictions.
To account for sampling bias of absence localities, we randomly selected a set of 1000 background localities within the boundaries of the core area (Maiorano et al. 2019;Dondina et al. 2020). We kept other default parameters of the MaxEnt. The overall model's performance was evaluated by calculating the area under the receiver operating characteristic curve (AUC), true skill statistic (TSS; Allouche et al. 2006), and Kappa statistic (Landis and Koch 1977). Then, we projected the calibrated model in the core area over the entire landscape area (i.e., both core and peripheral area).
The final distribution map was converted into a binary layer (i.e., suitable vs. unsuitable habitats) to identify critical habitats (areas with high suitability values predicted by the model) using Maximum Training Sensitivity plus Specificity occupancy threshold (Liu et al. 2013;Poor et al. 2020). We evaluated the reliability of the projection of the model to the peripheral area using the analysis of Multivariate Environmental Similarity Surfaces (MESS; Elith et al. 2010). Finally, we used Morphological Spatial Pattern Analysis in the software GUIDOS v. 2.1 (Soille and Vogt 2009) to remove predicted critical habitats with very small extents (< 70 km 2 ) that were unlikely to serve as core habitat. This threshold corresponds to the minimum threshold to ensure the stable occurrence of at least an adult female bear (Maiorano et al. 2019).

Connectivity between critical habitats
Species distribution models trained with predictors at fine spatial resolutions can be used as a useful alternative, although weaker than models based directly on movement, to movement data to parameterize resistance surfaces in connectivity analysis (Zeller et al. 2018). Ziółkowska et al. (2016) suggested that resistance maps derived from the species distribution models with predictors at broad and medium spatial resolutions may underestimate connectivity. Therefore, we included all covariates at a fine spatial scale (50-m) into the MaxEnt model to predict the distribution map of the species for further connectivity analysis. Then, the obtained distribution map was converted to a resistance surface using an exponential decay function (Wan et al. 2019a): where R is the cost resistance value assigned to each pixel and HS shows the predicted habitat suitability value. The calculated resistance surface was rescaled to a range between 1 and 10 by linear interpolation.
We followed the method developed by Macdonald et al. (2018) and Kaszta et al. (2020) to create a set of random points for habitat connectivity. This method was implemented in the following steps: first, the predicted distribution map at fine spatial scale was rescaled between 0 and 1. In the second step, we created a raster layer in the extent same as distribution map and distributed values of pixels from 0 to 1 randomly and uniformly. Then, we subtracted the created random raster layer described above from the rescaled distribution map to determine pixels with positive values. Finally, we selected 200 points randomly from these pixels. This creates a stochastic draw of simulated occurrence points and habitat suitability weighted connectivity modeling ) which may be parallel to the idea of using density-weighted connectivity modeling (Morin et al. 2017). Potential source locations were created in such a way as we assumed that higher predicted habitat suitability means more dispersing individuals (Macdonald et al. 2019;Kaszta et al. 2020). The selected source locations were spatially rarefied to simulate the presence of the species across the study area. The selected source locations were used as starting locations on individuals for subsequent connectivity analysis.
We calculated factorial least-cost paths ) using the UNICOR software  to predict connectivity among critical habitats based on the selected source points and predicted resistance surface described above. The selected source points define starting and ending nodes of the least-cost linkages in the connectivity analysis. The factorial least-cost path method relies on Dijkstra's algorithm to predict optimal least-cost linkages from each source point to every other destination locality. The predicted least-cost paths are buffered based on kernel density estimations and Gaussian function. The buffered least-cost linkages were then combined through summation ) to produce maps of connectivity among all pairs of source points. We calculated the factorial leastcost paths without a dispersal threshold to provide a broad-scale assessment of the spatial pattern of potential corridors .

Fire risk mapping
Based on extensive field surveys, Google Earth images, and various national reports collected by Iranian's organizations, we identified 358 fires locations in the study area from 2010 to 2018 (Fig. 1). We divided this dataset into training and validation subsets in a portion of 70/ 30, respectively. We selected the following variables for fire risk mapping based on our knowledge about the area: altitude, slope, aspect, topographic wetness index (TWI), distance from rivers, roads, urban, drainage density, annual mean rainfall, and temperature. We collected raster layers representing the above variables, and then extracted the value of each variable at each fire locations (see Text S1 in Online Resource 1). We used random forest to model the relationship between those variables and fire likelihood and then created a fire risk map by projecting that random forest model with the raster layers (Pourghasemi et al. 2020). Random forest has been shown to be one of the most powerful methods for predicting the disturbance and landscape change (e.g., Cushman et al. 2017;Abdullah et al. 2019;Murphy 2019;Milanović et al. 2021). We used mean decrease Gini and mean decrease accuracy to assess variable importance (Rahmati et al. 2019) and used OOB (Out-Of-Bag) error to assess model performance (Evans and Cushman 2009;Genuer et al. 2010). For the classification of the fire risk map, we used the natural break method (Yan et al. 2019) to divide it into four classes, including low (0-0.17), moderate (0.17-0.45), high (0.45-0.71), and very high (0.71-1.00) fire risk.

omposite map of predicted critical habitats and fire risk map
We overlaid the layers of predicted critical habitats, corridors, and conservation areas with the fire risk map to create a composite map and determine which critical habitats, corridors, and conservation areas are susceptible to fires.

Species' distribution and critical habitats
The MaxEnt model showed high performance in predicting suitable habitats for brown bear across the study area (AUC = 0.805; TSS =0.738; Kappa = 0.428). Topographic roughness (Reily et al. 1999;Evans et al. 2014), annual mean temperature, annual rainfall, and forest cover were the most important variables in predicting species' occurrence. With increasing topographic roughness and forest cover, habitat suitability increased to its maximum and then stabilized at high values of these variables. Mean annual temperature was negatively correlated with the species' occurrence and probability of brown bear occurrence was predicted to decrease with an increase in annual mean temperature. Also, with increasing rainfall habitat suitability increased to its maximum and then leveled off above 680 mm.
Using the Maximum Training Sensitivity plus Specificity occupancy threshold, 8.96% of the study area (11,276 km 2 ) was recognized as highly suitable habitat. Highly suitable habitats were more widely distributed across the northern, northeastern, and northwestern parts of the study area (Fig. 3a). According to the MSPA, we identified 17 critical habitats (larger than 70 km 2 ) with 9135 km 2 in extent across the area, 40.06% of which are covered by existing conservation areas (3,660 km 2 ). The extent and percentage of protected critical habitats differed among types of conservation areas. The maximum overlap was provided by no-hunting areas (2780 km 2 ), while no areas with critical habitats were identified inside the national park (Fig. 3b). We found some unprotected critical habitats especially north of the core areas of brown bear habitat and also northern parts of the study area (Fig. 3b). The results of Morphological Spatial Pattern Analysis showed that 18% of peripheral area were outside the range of environmental values measured in the main core region.

Connectivity between critical habitats
Our UNICOR analysis identified a total of 5245 km 2 migration corridors among critical habitats, of which 11.34% are within the boundaries of the existing conservation areas (595 km 2 ; Fig. 3b). The most density of the predicted corridors among critical habitats concentrated in the habitat core area (40.42%; 2,120 km 2 ) and northern parts of the peripheral region. Connectivity in other parts of the peripheral area was predicted to be weak. Considering the whole study area, nearly 90% of the predicted corridors are not covered by any type of conservation areas (Fig. 3b).

Fire risk map
The result of fire risk mapping is shown in Fig. 4. The percentages of the natural breaks of the four risk classes in fire map was 74.67% (low), 8.58% (moderate), 7.86% (high), and 8.89% (very high), respectively.

Composite map of predicted critical habitats, corridors, and the fire risk map
According to the intersection maps of the predicted critical habitats, migration corridors, and conservation areas with risk map of fire, a total of 35.65% of critical habitats of brown bear (3,257 km 2 ) and 23.56% of its corridors (1,236 km 2 ) were identified as areas with high and very high risk of fires respectively. Predicted critical habitats with high risk of fire were high especially in the northern and southwestern parts of core area. According to the fire risk map, 23.69% of the existing conservation areas with 5641 km 2 in extent overlap with areas where the modeled disturbance is likely to occur. No-hunting areas (4175 km 2 ) and national parks (4.33 km 2 ) showed the maximum and minimum overlap, respectively. We also identified currently unprotected areas that showed high risk of the occurrence of fire (Fig. 4).

Discussion
Despite numerous studies mapping landscape susceptibility to forest fire, few studies have investigated the sensitivity of key habitats to these events. Such studies require an accurate identification of key habitats and connectivity between them. This study is the first to compare fire risk with key habitats, corridors, and existing conservation areas to provide a spatially explicit depiction of the brown bear's critical habitats and corridors at high risk of fire.
A systematic framework to quantitatively analyze wildfire effects on brown bears, such as the one demonstrated here (Fig. 2), is much needed in Iran. For species highly associated with mature forest, there is a strong negative relationship between large and severe fires and habitat quality (e.g., Ganey et al. 2017;Wan et al. 2018;Jones et al. 2020). However, for species that are associated with broad habitat gradients across successional stages, and especially those that forage on vegetation that is associated with early seral conditions, as is brown bear, the relationship is not expected to be usually negative. Therefore, it is necessary to systematically assess the ecological responses of brown bear to differing kinds, extents, and severity of disturbance to guide conservation and management scenarios that are effective (e.g., Chambers et al. 2019).

The effect of environmental variables on brown bear's distribution range
Our model showed that the range of brown bear was constrained mainly by topographic roughness, temperature, rainfall, and forest cover in Iran. Despite differences in environmental conditions such as rainfall and vegetation in different regions where brown bears occur, the importance of forest cover to brown bear's habitat selection has been documented across the globe, such as in Iran (Mohammadi et al. 2021), North America (Proctor et al. 2015), and Europe (Mateo Sanchez et al. 2014;Recio et al. 2021). Forest cover in southwestern Iran provides crucial food and shelter for the brown bear (Ansari and Ghoddousi 2018).
Considering the importance of forests in providing suitable habitats for bears, protecting forest land is the utmost priority for management. Habitat fragmentation

Fig. 3 Predicted habitat suitability of brown bear A and intersection map of critical habitats and corridors with existing conservation areas B
here is often associated with the conversion of forests to agricultural land, which leads to higher presence of livestock and crops and thus greater likelihood of humanbear conflict. These fragmented forests typically contain smaller core habitats with reduced natural food and water availability and access, lower overall connectivity, and greater human presences, which may alter carnivore behaviors such as increasing their aggressiveness. Therefore, the prevention of further land-use changes and forest fragmentation is urgently needed to maintain the suitability of forest habitats for bears. We encourage outreach programs to local communities to promote sustainable land use and control timber harvest and livestock grazing. We also recommend reforestation in targeted areas that can provide core areas and facilitate connectivity of brown bears.
Consistent with previous studies, we also found that rugged terrain, especially with forest cover, is a critical component of brown bear habitat (Mateo Sanchez et al. 2014;Peters et al. 2015). In general, elevation and rugged landscapes are positively related to bear habitat suitability because they provide optimal climatic and vegetation conditions and limit human presence (Almasieh et al. 2019;Zarzo-Arias et al. 2019;Mohammadi et al. 2021). However, brown bears in Spain generally avoid high-elevation meadows because the food resources are scarce compared to forests at lower elevations (Mateo- Sanchez et al. 2016), and they prefer forests with relatively dense cover further away from settlements at higher elevations (Zarzo-Arias et al. 2019). These differences were likely due to inherent differences in climate, vegetation patterns, and human density and behavior between the two study areas.
Rainfall and perennial surface water is important for providing habitat and can act as limiting factors in range expansion of the bear (Mertzanis et al. 2008;Ansari and Ghoddousi 2018). Because of the aridity in Iran, collaborative management of shared water, especially during the drought years, will be the most crucial action for habitat improvement. Examples of these actions include the installation of water troughs, as well as the restoration of depreciated local surface water resources such as springs, at strategic locations where wildlife and livestock both have access to (Ansari and Ghoddousi 2018).

Brown bear critical habitats and connectivity between them
We found that some of the predicted suitable habitat patches may be too small to sustain a population despite their high suitability. In a similar study, Maiorano et al. (2019) predicted critical areas for brown bear and found only 15 out of 461 predicted suitable patches met the minimum threshold to be considered as critical areas. In our study, we found 17 of 1640 patches that met our minimum size threshold, suggesting severe fragmentation of critical habitat.
Also, although the predicted critical habitats in the core region are well connected, some of the peripheral habitat patches, especially those in the western and northwestern parts of the study area, are not well connected with the core region. Other research showed that habitat connectivity in simulations strictly depends on the dispersal ability of species . Consistent with our findings, Ashrafzadeh et al. (2018) reported strong connectivity among core brown bear habitats in Iran, and that the connectivity may be constrained by several pinch points due to bottlenecks in habitat quality, or high concentrations of human activities. Our connectivity projection provides key information to assist managers in identifying and prioritizing efforts in managing key corridors for brown bear dispersal.

Overlap of key habitats and corridors by conservation areas
Conservation areas are essential to provide extensive and contiguous habitats for carnivores in Iran because they support a high density of prey species and are exposed to lower human disturbance (Ashrafzadeh et al. 2020). In contrast to conservation areas, unprotected areas in Iran are believed to exhibit high resistance to species movement because of unsuitable habitats and high density of human settlements (Moqanaki and Cushman 2017). Despite their important role in carnivore conservation, critical habitats and corridors received relatively poor protection under the current design of conservation areas. This problem has been discussed for other carnivores as well (Moqanaki and Cushman 2017;Ashrafzadeh et al. 2018;Kaszta et al. 2020). Location bias (designing new conservation areas mainly in topographically remote landscapes), conservation bias (i.e., varying conservation effort in different types of the conservation areas), and spatial bias (i.e., designing new conservation areas regardless of their spatial configuration within the compositional structure of the current conservation areas) in designing new conservation areas in Iran may lessen conservation effectiveness of the conservation network for some species, as for the brown bear. Therefore, conservation areas do not necessarily have more permeability to individuals' movement.
The highest overlap between conservation areas, critical habitats, and corridors was observed in the northwest of the core area and north, northeast, and southwest of the peripheral area, with nine conservation areas (Tang-e-Bostanak, Margoon, Male-Gale, Barm Firoz, Arjan & Parishan PAs, Basiran, Kooh-Khersi, Dare-Bagh, and Kooh-e-Khom NHAs) covering the most important identified critical habitats. Because the coverage of conservation areas is not sufficient to protect key habitats and migration corridors for brown bears, the proportion of conservation areas with highly suitable habitats should be increased to maximize the protection of habitat extent and maintain connectivity. Establishing new conservation areas in unprotected critical habitats is highly encouraged, especially in the northern parts of the study area where most of the key habitats are located.

The risk of occurrence of fire within the conservation areas, critical habitats, and corridors
This study highlights the high wildfire risks that threaten some of the existing conservation areas. For example, Tang-e-Bostanak, Male-Gale, Arjan & Parishan PAs, Kooh-Khersi, Padena, and Mel-Boland NHAs had high risk of fire. No-hunting areas had the largest coverage of the predicted critical habitats than other types of conservation areas. These types of conservation areas have recently been established to increase the species populations in areas with illegal hunting. Despite the importance of no-hunting areas in the protection of key habitats, these areas are more likely to be exposed to fire due to less strict conservation, and less legal protection against anthropogenic impacts. This finding has important message for practitioners, that no-hunting areas are more in danger and therefore need further consideration in the conservation programs. An important suggestion can be to upgrade the existing no-hunting areas in the region to official protected areas (e.g., protected areas, wildlife refuge, or national parks) and to expand them to incorporate also neighboring lands where species may occur. Also, the results of species distribution modeling and habitat connectivity showed that most of the critical habitats and corridors are located outside conservation areas. For example, there is a high risk of fires in isolated and unprotected critical habitats in western and northern parts of the peripheral area. Brown bears in these regions may experience increased exposure to wildfires in the future.
Fires is an important regulator that helps maintain ecosystem health. However, fires have increasingly become a threat to many forests and the biodiversity due to altered fire regimes (Cushman et al. 2011). The results of composite maps showed that the highest fire risk in the predicted critical habitats is in the central parts of the core zone and the western parts of the peripheral areas. These areas contain the study area's largest extent of forests in the Zagros Mountains. In addition, the migration corridors between critical habitats in the central parts of the core zone and western parts of the peripheral area showed the greatest potential for fire.
While fires may have short-term positive effects on some fire-adapted species, wildfires have contributed to the reduction in populations of some large carnivores such gray wolf (Canis lupus) in Minnesota because of its fire-dependent preys (Kramp et al. 1983). Cunningham and Ballard (2004) studied short-term impacts of wildfire on black bear demographics in central Arizona. They found that wildfire had significant negative impacts on adult sex ratio, recruitment of cubs, and superpopulation size. Richardson et al. (2007) showed that forest fires decrease the stability of polar bear den sites (Ursus maritimus) in western Hudson Bay and lead to the degradation of the surrounding habitats by significantly altering vegetation composition and increasing the depth of the active layer. According to the results of Richardson et al. (2007), bears do not use burned areas for denning. The accelerated rate in forest fires as a result of changing climate may reduce the amount of suitable denning habitats in the future. Therefore, managers should be aware of possible negative short-term consequences of wildfire on cub survival. In addition, the removal of forests by increasing fires may increase forest fragmentation. As a result, bear activities that typically occur in forests will be affected. They either will have to adapt to the altered landscape, or will have to move to remnant patches of intact forests. In either case, their energy expenditure will increase and thus may lead to reduced fitness. So, improving connectivity between forest remnant patches will be useful. Additionally, regulations that limit human induced wildfires will be helpful.
Knowledge on the brown bear's key habitats with high risk of fire is necessary for wildlife managers to understand the short-term and long-term impacts of fires on bears. Our model can assist with planning the locations of future prescribed fires or early reconstruction and recovery actions following fires.

Scope and limitations
Although our framework performed well in identifying key habitats with high risk of forest fire, we recognize shortcomings of the framework. First, as we did not include GPS data to validate the predicted critical habitats and corridors, our models may not accurately predict distribution range and movement patterns of the studied species especially in the peripheral area (e.g., Zeller et al. 2018). However, we found high performance of the model to predict suitable habitats outside the calibration area, and some prior work on brown bear found that habitat suitability modeling can be an adequate surrogate for connectivity (Mateo-Sánchez et al. 2015a;Mateo-Sánchez et al. 2015b). Second, we did not include more proximate environmental variables into the models. For example, while the predicted suitable habitats are located in forest regions, we did not consider the effects of composition and configuration of different forest types in the models. Mateo Sanchez et al. (2014) showed that brown bears respond differentially to different kinds of forest and different spatial metrics measuring forest configuration in the landscape. Another factor that may pose additional limitations on the species presence in the predicted suitable habitats is human-wildlife conflict (e.g., Rostro-García et al. 2016). If the predicted critical habitats are located in close proximity to human settlements, the persecution of brown bear due to human-wildlife conflict may continue to increase and pose additional limitations on the occupation of the predicted critical habitats (e.g., Cushman et al. 2018). So, we strongly recommend incorporating human attitudes toward bears and human-wildlife conflict into models (e.g., Broekhuis et al. 2017).
Radiotelemetry or satellite tracking as a basis for testing the predicted linkages among habitats has just been started for a few species in Iran. We hope that wildlife tracking will become more feasible in the near future, but until then, we have to rely on the habitat suitabilitybased estimates of landscape resistance. Given the large assumptions made by developing a resistance surface from a habitat model based, it is important to empirically validate the predicted corridors through monitoring movements of bears using camera-trapping or satellite tacking (e.g., Proctor et al. 2015) and landscape genetic analyses (e.g., Proctor et al. 2005;Kopatz et al. 2012). For example, Proctor et al. (2015) combined resource selection functions and human presence data to predict core habitat and develop cost surfaces for connectivity network of grizzly bears (Ursus arctos) in Canada-US trans-border region, and evaluated their prediction through GPS telemetry locations. Peck et al. (2017) used step-selection functions and randomized shortest path to generate conductance surfaces and estimated the average number of net passages for grizzly bear populations in the Greater Yellowstone Ecosystem and the Northern Continental Divide Ecosystem. Their finding confirmed that randomized shortest path algorithm provides detailed, spatially explicit information for identifying and prioritizing conservation measures aiming at improving landscape connectivity.

Conclusion and conservation implications
There are four major findings from this study that are especially relevant to brown bear conservation and management. First, the spatial pattern of brown bear's critical habitats showed that most of the key habitats are located in the northern parts of the study area and have moderate to high connectivity. Second, we found that 40% of the predicted critical habitats and 11% of corridors are located within the boundaries of the existing conservation areas which may increase exposure to the effects of anthropogenic disturbance due to lack of conservation programs in unprotected areas. Third, the results of forest fire risk mapping showed that 35.65% of the species' key habitats are located in regions with a high risk of fire. Fourth, although the susceptibility of the existing conservation areas (such as national parks) to fire was low, some of the conservation areas are prone to the occurrence of fires. Our results can be used to inform biodiversity managers to design conservation efforts before the occurrence of fires and early reconstruction and recovery actions following event (Chambers et al. 2019).
The accelerated rate of human-induced habitat destruction and changing climate has led to increases in the size, frequency and severity of wildfires in many parts of the world, including Iran. However, there is little evidence regarding the long-term effects of changing fire regimes on most wildlife species. Mitigation efforts should nevertheless consider the minimization of direct or indirect effects of wildfires on species when they are strongly negative and can be mitigated by practical management means (e.g., Ganey et al. 2017;Wan et al. 2018).
The habitat suitability and connectivity models we produced clarify the distribution, habitat relationships, and linkages across the brown bear population, whereas the fire risk model provides information on which habitat and corridors are most susceptible to fires. What remains to be explored is what effects of fires are on bears and if they are scale dependent in their effects (e.g., Wan et al. 2018;Wan et al. 2020). Although wildlife habitats are threatened by large and high-severity wildfires due to the warming climate, there is few data on the effect of spatial scale on the relationship between wildfires and habitat modifications. For example, Wan et al. (2020) investigated potential relationships between burn severity and changes in habitat suitability of Mexican spotted owl (Strix occidentalis) at different spatial scales though multi-scale habitat selection model. They found that fire effects on habitat suitability of the species changed at different scales and the strength and direction of relationships were scale-dependent.
Our analytical framework for predicting critical habitats and their susceptibility to fire has many management and conservation implications. First, the critical habitats with high risk of fire with negative impacts should receive research, management, and conservation priority. We recommend constructing fire trucks in areas with a high risk of fire and management according to principles of fire ecology and land use planning-LUP requirements. Also, the predicted migration corridors with low risk of forest fire should be considered in maintaining or improving connectivity between the critical conservation areas (Wasserman et al. 2013). We recommend increasing the connectivity among critical habitats associated with high forest fire risks, especially in conservation areas to allow range expansion of the species and facilitate its movement. For example, this can be accomplished through designating new conservation areas with different degrees of protection through zoninghigher restrictions in critical habitat associated with higher fire risk. Also, the creation and maintenance of intermediate stepping stone corridors might facilitate gene flow and prevent the isolation of small populations while acting as fire breaks in the event of a fire. Finally, we propose that the systematic risk assessment framework demonstrated here can be more broadly applied for the conservation of other forest species that are facing similar threats.
Additional file 1: Table S1. List of variables used for predicting critical habitats of brown bear and connectivity between key habitats. Text S1. Preprocessing of variables.