Population‐level monitoring of stress in grizzly bears between 2004 and 2014

Grizzly bears (Ursus arctos) in west-central Alberta occupy an increasingly human-dominated landscape. Natural resource extraction activities are hypothesized to increase stress in animals that reside in such changing landscapes by influencing habitat and resource availability. Our study aimed to determine whether stress, represented by hair cortisol concentration (HCC), was associated with variables related to landscape conditions in a population that increased by 7% annually from 2004 to 2014. Hair samples (n = 157) were collected using barbwire hair snags placed throughout the Yellowhead bear management area in Alberta, Canada. Candidate models were developed a priori representing hypotheses related to biologically and ecologically plausible relationships between HCC and landscape variables. Generalized linear model analysis with landscape attributes representing anthropogenic disturbance, food resource availability, and terrain conditions was used to determine potential drivers of HCC. We found support (DAICc ≤ 2.00) for three models that included variables from each hypothesis. Anthropogenic variables had the greatest impact on HCC; increasing oil and gas well-site density resulted in reduced HCC, while increasing distance to coal mines resulted in elevated HCC. Hair cortisol concentration also increased as forest crown closure became more variable, while HCC decreased as the soil wetness (represented by compound topographic index) increased. Some forms of anthropogenic disturbance have been linked to increased food availability for this species. Therefore, we suggest that changes in landscape conditions from 2004 to 2014 may have indirectly increased food abundance and ultimately resulted in a reduction in HCC at a population level during this time period. Measuring HCC provides a non-invasive and important monitoring strategy to assess the impact of environmental change on residing species and should be considered in landscape management decisions.


INTRODUCTION
Over the past decade, there has been a loss of wild areas across the globe (Venter et al. 2016). The loss of these areas can lead to changes in habitat quality and food availability for wildlife, which may influence chronic stress in individuals that reside on changing landscapes (Rimbach et al. 2013, Zbyryt et al. 2018 and Barja 2019). Terrestrial animals respond to stress through the activation of the hypothalamic-pituitary-adrenal (HPA) axis (Sapolsky et al. 2000). Corticotropin-releasing hormone (CRH) is synthesized in the hypothalamus and is then followed by the release of adrenocorticotropic hormone (ACTH) from the anterior pituitary (Spencer and Deak 2017). Finally, ACTH triggers the release of glucocorticoids, such as cortisol and corticosterone, which play an important role in the processes needed for energy metabolism and biological function (Sapolsky et al. 2000, Miller andO'Callaghan 2002). Since HPA axis activity increases with exposure to stressors, cortisol concentration is typically measured as a biomarker for stress (Dallman et al. 1987, Sapolsky et al. 2000. The measurement of cortisol in hair has been recognized as a potentially important advancement in studying the role of stress in wild animals (Koren et al. 2002, Sheriff et al. 2011) and can provide a reliable biological marker of chronic stress (Macbeth et al. 2010, 2012, Russell et al. 2012, Stalder and Kirschbaum 2012, Burnard et al. 2016. During the active hair growth phase, free circulating cortisol continuously diffuses into the hair shaft from blood vessels (Meyer and Novak 2012). Therefore, it is possible to measure HPA axis activity and subsequent cortisol release over a previous time period (Russell et al. 2012), making hair cortisol concentration (HCC) a useful marker of stress. Hair cortisol concentration has been used as a biomarker of stress in a variety of animals with respect to management conditions, social behavior, nutritional status, and disease (Heimb€ urge et al. 2019). Conversely, HCC may also be influenced by acute stressors, such as capture and restraint, and this may occur when hair growth has ceased, indicating caution should be used in how hair samples are collected when investigating HCC as a biomarker of stress ). Furthermore, HCC has been shown to vary with hair type, body region, and capture method in grizzly bears (Macbeth et al. 2010).
Northwest North America is home to some of the world's remaining grizzly bears (Ursus arctos). Grizzly bears are considered a threatened species in Alberta, Canada, with fewer than 700 individuals (circa 2010; Festa-Bianchet 2010, Clark and Slocombe 2011). Grizzly bear habitat in Alberta has become highly fragmented due to resource extraction activities (i.e., forestry, oil and gas exploration, mining, and agriculture), widespread road networks, and human recreational activities (Berland et al. 2008, Proctor et al. 2012. Anthropogenic disturbed areas can cause an increase in seasonal food availability (Nielsen et al. 2016, Kearney et al. 2019, Larsen et al. 2019, such as increases in the abundance and occurrence of berries and fruits associated with more clearings and forest opening, which grizzly bears have been found to select and can result in an increase in human contact and reduced bear survival rates (McLellan and Hovey 2001, Nielsen et al. 2004a, 2004b, 2004c, Berland et al. 2008. Furthermore, this multi-use landscape has facilitated genetically fragmented sub-populations (Proctor et al. , 2012. As a result, a considerable amount of research has taken place on conservation strategies to inform policymakers and management activities in the region (Berland et al. 2008, Nielsen et al. 2009, Boulanger and Stenhouse 2014. There is still much uncertainty about current population dynamics in some areas of the province and if anthropogenic landscape change and human activities have influenced the HPA activity of grizzly bears in Alberta.
Effective monitoring systems are necessary to understand and mitigate the impacts of cumulative stressors on wildlife (Côt e et al. 2016), which has been particularly challenging for scientists and managers (Williams et al. 2002). Traditional landscape management and wildlife monitoring methods, such as counting individuals and identifying changes in habitat use, do not address the effects of landscape change on the HPA activity of individuals (Macbeth et al. 2010). Stress can influence the HPA activity and fitness in individual animals (Acevedo-Whitehouse andDuffus 2009, Schultner et al. 2013), and increased glucocorticoid concentrations may also provide early warning of declining population performance (Jachowski et al. 2012, Strasser andHeath 2013). However, the chronic elevation of glucocorticoid concentrations may also serve as an adaptive response to stressors, such as limited food availability, and subsequently facilitate fitness (Boonstra 2013). Resource availability and nutritional status influence physiological function (Busch and Hayward 2009) and therefore likely play an important role in HPA activity. Food abundance and diet type have been shown to influence glucocorticoid concentrations in grizzly bears (von der Ohe et al. 2004, Wasser et al. 2004, Bryan et al. 2013. For example, grizzly bears residing in coastal British Columbia, Canada, demonstrated an increase in HCC after a year of low dietary salmon compared to after a year of high dietary salmon (Bryan et al. 2013). Furthermore, reduced HCC values were demonstrated with improved body condition of polar bears and grizzly bears (Macbeth et al. 2012, Bourbonnais et al. 2014). It has also been shown that HCC values in female grizzly bears, and to a lower degree in males, were reduced in areas highly impacted by human disturbance (Bourbonnais et al. 2013). The authors suggest that decreased HCC values may be related to seasonal increases in food availability and distribution in areas with anthropogenic disturbance (Nielsen et al. 2016, Kearney et al. 2019, Larsen et al. 2019. This study aimed to (1) investigate how changing landscape conditions may be affecting stress in grizzly bears between 2004 and 2014 and (2) determine if monitoring HCC over time is a suitable environmental monitoring approach. We hypothesized that increased resource extraction activities have influenced glucocorticoid concentrations in grizzly bears. To test this hypothesis, the objective of this study was to determine if HCC was associated with landscape variables related to anthropogenic disturbance, food resource availability, and terrain conditions in 2004 and 2014. Two large-scale population inventory projects took place within the study area during 2004 and 2014, which demonstrated a 7% annual increase in population size from 2004 to 2014. Therefore, HCC was only measured in 2004 and 2014 and not during the time between those years. We predicted that anthropogenic disturbance would be associated with increased HCC, as our study population resides on a landscape heavily impacted by resource extraction activities. Furthermore, we predicted that increased food resource availability would be associated with decreased HCC. Terrain conditions, such as elevation, soil wetness, and ruggedness of terrain, were also expected to be related to HCC, as these variables have been shown to influence food availability for this species (Nielsen et al. 2004c). As the landscape continues to change over time, it is crucial to develop methods to monitor glucocorticoid concentrations in animals and better understand the relationships between landscape conditions and HPA activity.

Study area
This study was conducted in the Yellowhead bear management area (BMA) in Alberta, Canada (Fig. 1). The Yellowhead region includes Jasper National Park as well as provincial parks and protected areas that demonstrate little landscape change (Bourbonnais et al. 2013). Samples were collected across elevations ranging from 1100 to 2500 m. As the elevation increases from east to west, the habitat characteristics include both mixed-wood forests consisting of pine (Pinus spp.), aspen and poplar (Populus spp.), spruce (Picea spp.), and balsam fir (Abies balsamea) and alpine and sub-alpine ecosystems comprised of fir (Abies spp.), pine, and spruce, and wet meadow complexes . Grizzly bear food resources in this study area include horsetail (Equisetum spp.), sweetvetch (Hedysarum alpinum), cow-parsnip (Heracleum lanatum), buffaloberry (Shepherdia canadensis), clover (Trifolium spp.), huckleberry (Vaccinium membranaceum), insects, and ungulates (Nielsen et al. 2010).

Sample collection
Hair was collected by using barbed wire hair snags placed throughout the Yellowhead BMA from June to July 2004, and again from June to August 2014 , Nielsen et al. 2010, Bourbonnais et al. 2013. The barbed wire hair snags were randomly placed in 7 9 7 km grid cells, and hair was collected every 2 weeks . Previous studies have shown that the 7 9 7 km grid cell size is appropriate for a high enough capture-probability level and reduces sampling bias for population and density estimates in grizzly bears (Boulanger et al. 2002. Furthermore, hair was collected in each 7 9 7 km grid cell in conjunction with previous population density studies in 2004 and 2014 . Thirty-two hair samples were collected across 21 grid cells in 2004, and 125 hair samples were collected across 51 grid cells in 2014 (Appendix S1: Table S1). If wet upon collection, hair was blotted dry with paper towels and the envelope left open until the sample was dry. Following collection, hair samples were placed into envelopes and stored at room temperature. Portions of the samples were shipped to Wildlife Genetics International (Nelson, British Columbia, Canada) for DNA extraction to identify individuals and sex and to University of Saskatchewan for cortisol analysis. This research was conducted in accordance with the Alberta Environment and Parks Animal Care Committee and Parks Canada.

Extraction of cortisol
The extraction of cortisol from grizzly bear hair was completed as previously described (see Macbeth et al. 2010, 2012, Cattet et al. 2017. Briefly, external contamination of hair samples (e.g., by blood, feces, soil) was removed using three washes with methanol. Washed and dried hair was then ground to a fine powder in a mixer mill (Retsch MM400; Retsch GmbH, Haan, Germany) and placed on a spinning rotator to extract into methanol for 16-24 h. The samples were then spun in a centrifuge for 15 min (2150 g/20°C) and the supernatant was collected and dried under a gentle stream of nitrogen gas for a total of three collections. Concentrated samples were reconstituted with the appropriate volume of buffer solution (10 µL/mg). Reconstituted samples were centrifuged to remove any trace hair residue, and the supernatant was collected for storage at À80°C until analysis. ❖ www.esajournals.org

Potential drivers of HCC
In order to develop candidate predictor variables of grizzly bear HCC related to landscape conditions, we generated both annually varying and static geo-spatial variables from remotely sensed data for 2004 and 2014 (Table 1 and derived in Bourbonnais et al. 2013). Candidate variables represented terrain conditions (static), as well as anthropogenic disturbance and food resource availability (annually varying), all of which have been identified as potential drivers of glucocorticoid concentration in bears (Wasser et al. 2004, Bourbonnais et al. 2013, Bryan et al. 2013, Malcolm et al. 2014, Cristescu et al. 2016b). Candidate variables were extracted and summarized at the 7 9 7 km grid cell level using both discrete vector and continuous raster datasets. All extractions were completed using Esri's Arc-Map 10.3 GIS software.
Anthropogenic disturbance was represented by vector data for locations of roads, coal mines, rail and power lines, oil/gas well-sites, forest harvesting (cutblocks), and protected areas (Table 1). These data originated from the Human Footprint Inventory (HFI) developed by the Alberta Biodiversity Monitoring Institute (ABMI) and were mapped by the Foothills Research Institute Grizzly Bear Program (FRIGBP) based on appearance of disturbance features in Landsat Thematic Mapper imagery. For linear features (roads and rail/power lines), we calculated linear density per unit area (km/km 2 ) using a 6-km radius moving window and averaged this for each grid cell. We also calculated the mean distance to each respective feature. Well-sites were represented as point data; we calculated point density (count/km 2 ) using the same moving window and averaged this for each grid cell in addition to the mean distance across the cell. For forest harvesting, we summed the total area of all cutblocks within each grid cell and calculated the percent cutblock area, the mean distance to cutblocks, and the percent area of cutblocks less than 15 yr old, for which grizzly bears are known to select over older cutblocks (Kearney et al. 2019). Lastly, we calculated the percent protected area for each grid cell by dividing protected area by grid cell area.
Food resource availability included variables that characterize forest conditions and represent the limitation, distribution, abundance, and quality of foods, as well as the availability and distributions of foods related to disturbance features. Many bear foods are more abundant where the forest canopy is open or near forest edges (Nielsen et al. 2004c, Larsen et al. 2019, and bears using regenerating forests with a variety of age classes have shown improved body condition compared to those using older, conifer-dominated forests (Boulanger et al. 2013). Therefore, we used raster data on crown closure, percent conifer, and land cover-originally derived from Landsat 5 and 7 ETM imagery and a digital elevation model (DEM)-to represent food resource availability (McDermid 2005). Regenerated area was generated using ABMI HFI harvested area data. Mean crown closure and mean percent conifer were calculated for each grid cell to represent the degree to which a cell was covered in forest and dominated by conifers. The mean distances to upland herbaceous, upland treed, and regenerating forest land cover classes (i.e., cutblocks and wildfires ≤50 yr old) were each calculated to represent access to open and regenerating vegetation. We also calculated the standard deviation ❖ www.esajournals.org of crown closure and regenerating areas to represent the variation of forest age classes present in each grid cell.
Terrain conditions have been shown to be important predictors of HCC in grizzly bears (Bourbonnais et al. 2013), and we derived three terrain variables from a DEM: elevation, the compound topographic index (CTI), and terrain ruggedness index (TRI). The CTI represents water flow accumulation (Beven and Kirkby 1979) and soil moisture based on slope and aspect, and upstream water sources (Gessler et al. 1995). Thus, it can highlight riparian areas, which provide valuable habitat for grizzly bears, especially in multi-use landscapes (Phoebus et al. 2017). The TRI was used to quantitatively measure the complexity and variability of the terrain (Riley et al. 1999, Wilson andGallant 2000). Variation in topographic heterogeneity has been shown to influence food resource availability (Nielsen et al. 2004c(Nielsen et al. , 2010, livestock depredation events (Wells et al. 2019), den selection (Pigeon et al. 2014), and mortality risk (Nielsen et al. 2004b) in grizzly bears.

Statistical analyses
In order to evaluate if HCC data may be a useful tool to monitor bear HPA activity in response to the environment, we first compared HCC across samples with known biological attributes (sex: male and female and individual: unique bear identity). We initially analyzed these data for normality and outliers, and could not achieve normality with data transformation. Therefore, non-parametric statistical comparisons (Kruskal-Wallis rank-sum test) were used to compare samples within sex and individual. If significant (P < 0.05), post hoc multiple comparisons were completed with a Dunn test. Differences in HCC with respect to biological attributes were not expected and thus would suggest that observed differences in HCC are linked to environment, rather than demography. Each grid cell contained HCC values from different individuals. While there were multiple samples from the same individual in each year (n = 6 HCC values from 3 individuals in 2004 and n = 79 HCC values from 24 individuals), this was overcome by averaging HCC values within each cell. In 2014, For data exploration, we followed the procedure described in Zuur et al. (2010), which has been used in previous HCC studies with grizzly bears (Cattet et al. , 2017(Cattet et al. , 2018. We first logtransformed HCC values to help normalize the data and used Cleveland dot-plots to evaluate potential outliers in the HCC data relative to the order of the data (sample year; 2004 and 2014; Appendix S1: Fig. S1). Next, we used pair-plots (Pearson r ≥ 0.80) and variance inflation factors (VIFs ≥ 3.0) to identify collinear variables in each candidate model. Since highly correlated variables represent similar landscape attributes, highly correlated predictor variables were first removed from each candidate model by keeping the variable most strongly correlated with HCC (Appendix S1: Table S2; Chatterjee and Hadi 2006). Correlated variables were removed within each model (models 2-4) as well as across models (models 5-8). Covariates with the highest VIF were sequentially dropped from each candidate model, and the VIFs were recalculated. This process was repeated until all VIFs were <3.0 for each candidate model (Zuur et al. 2010).
To test the hypotheses that HCC was associated with landscape variables related to anthropogenic disturbance, food resource availability, and terrain conditions, we used generalized linear model analysis to develop candidate models representing biologically and ecologically plausible relationships between HCC and landscape variables. First, all continuous predictor variables were standardized by subtracting the mean from the observed values and dividing by the standard deviation (Cattet et al. 2018). We then developed 8 candidate models to determine the influence of predictor variables on HCC (Table 2). We used backward elimination to refine each candidate model by sequentially removing the least significant variable (P> 0.05) until all variables in the reduced model were significant (Table 3; Cattet et al. 2018). We followed an information-theoretic approach to compare the 8 refined candidate models by Akaike information criteria (AIC) weights and selected the most supported models (DAIC c ≤ 2.00) for further analysis (Burnham and Anderson 2002). Model residuals were plotted against a normal distribution, and R 2 (1-deviance/null deviance) values were calculated to assess model fit (Appendix S1: Fig. S2). Once the final candidate models were generated, the sample year (either 2004 or 2014) was included as a fixed effect and assessed only across cells that were sampled in both years (n = 6). Similarly, we compared the remaining candidate models with sample year as a fixed effect by AIC weights and selected the most supported models (DAIC c ≤ 2.00) for further analysis. We then used Tukey HSD post hoc comparison to compare HCC between 2004 and 2014 in the most supported models and with the effects of the other variables taken into account. Data were analyzed using R statistical software, and data were deemed significant when P < 0.05 and the 95% confidence interval (CI) did not include zero.

RESULTS
A total of 157 hair samples (n = 105 samples from 43 females and n = 52 samples from 50 males) were collected from grizzly bears in the Yellowhead BMA. No outliers were detected; therefore, all hair samples were used in the analysis. Furthermore, there were no missing values in the HCC or the landscape data. HCC did not vary significantly (P> 0.05) within sex (Appendix S1: Fig. S3) or individual. HCC ranged from 1.58 to 17.01 pg/mg in 2004 (4.22 AE 3.14 pg/mg; mean AE SD), while in 2014, HCC ranged from 0.10 to 15.37 pg/mg (2.53 AE 2.15 pg/mg). Fig. 2 shows the distribution of HCC values across the study area in 2004 and 2014.

Drivers of HCC
Following backward elimination, the 8 initial models were refined to 4 candidate models containing variables related to anthropogenic disturbance, terrain conditions, and food resource availability (Table 3). We found support (DAIC c ≤ 2.00) for three candidate models to predict changes in HCC: model 2: anthropogenic disturbance; model 4: terrain conditions; and model 5: anthropogenic disturbance and food resource availability (Table 4). Model 2 included one variable related to anthropogenic disturbance, oil and gas well-site density, while model 4 contained a single variable, mean compound topographic index, which represented terrain condition. Finally, model 5 was comprised of one variable related to anthropogenic disturbance, mean distance to coal mines, and one variable related to food resource availability, standard deviation of crown closure. Model 2: Anthropogenic disturbance Road density + mean distance to road + mean distance to coal mines + railway and power line density + mean distance to oil and gas well-sites + oil and gas well-site density + percent protected area + percent area of cutblock (<15 yr old) 3 Food resource availability Mean distance to upland herbaceous resources + mean distance to regenerated areas + standard deviation of distance to regenerated areas + mean crown closure + standard deviation of crown closure 4 Terrain conditions Mean compound topographic index 5 Anthropogenic disturbance and food resource availability Mean distance to road + mean distance to coal mines + railway and power line density + oil and gas well-site density + percent protected area + percent area of cutblock (<15 yr old) + mean distance to upland herbaceous resources + mean distance to regenerated areas + standard deviation of distance to regenerated areas + standard deviation of crown closure 6 Anthropogenic disturbance and terrain conditions Mean distance to cutblock + mean distance to road + mean distance to coal mines + railway and power line density + oil and gas well-site density + percent protected area + percent area of cutblock (<15 yr old) + mean compound topographic index 7 Food resource availability and terrain conditions Mean distance to upland herbaceous resources + mean distance to regenerated areas + standard deviation of distance to regenerated areas + standard deviation of crown closure + mean compound topographic index 8 Anthropogenic disturbance, food resource availability, and terrain conditions Mean distance to road + mean distance to coal mine + railway and power line density + oil and gas well-site density + percent protected area + percent area of cutblock (<15 yr old) + mean distance to upland herbaceous resources + mean distance to regenerated areas + standard deviation of distance to regenerated areas + standard deviation of crown closure + mean compound topographic index Anthropogenic disturbance and food resource availability Mean distance to coal mine + standard deviation of crown closure 6 Anthropogenic disturbance and terrain conditions Oil and gas well-site density 7 Food resource availability and terrain conditions Mean compound topographic index 8 Anthropogenic disturbance, food resource availability, and terrain conditions Mean distance to coal mine + standard deviation of crown closure anthropogenic disturbance was the most supported model (DAIC c = 0.00). HCC decreased significantly (b = À0.087; P = 0.025; 95% CI = À0.163, À0.012) with increasing oil and gas well-site density (Table 5). There was also support for model 5: anthropogenic disturbance and food resource availability (DAIC c = 0.33). HCC increased significantly (b = 0.110; P = 0.022; 95% CI = 0.017, 0.203) with increasing distance from coal mines (Table 5). Similarly, HCC increased significantly (b = 0.111; P = 0.021; 95% CI = 0.018, 0.204) when the standard deviation of crown closure was higher (Table 5). Lastly, there was support for model 4: terrain conditions (DAIC c = 1.26). HCC decreased with increasing CTI; however, this relationship was only ❖ www.esajournals.org 9 July 2020 ❖ Volume 11(7) ❖ Article e03181 trending toward significance (Table 5; b = À0.076; P = 0.052; 95% CI = À0.152, 0.000). The predicted mean values of HCC for given values of each of the significant predictor variables (effects plots) are shown in Fig. 3. All variables demonstrated a linear relationship as described above with predicted HCC when all other variables were held constant. When year was included in the final candidate models, only model 5: anthropogenic disturbance and food resource availability was supported (Appendix S1: Table S3; DAIC c = 0.00). Post hoc pairwise comparisons revealed that samples collected in 2014 had significantly lower (Tukey's HSD, P < 0.01) HCC than in 2004 when the effects of the predictor variables in model 5 were taken into account.

Predicted HCC values in 2004 and 2014
HCC values were predicted across all grid cells in the study area using the landscape variables included in model 5: anthropogenic disturbance and food resource availability (Fig. 4). Only this model was supported (DAIC c = 0.00) when year was included in the final candidate models; however, predicted hair cortisol concentrations were generated from the final model without year. Therefore, unidentified change associated with sample year is likely present. The model predicted that HCC values increased from 0.1 to 0.5 pg/mg in 43% of cells and 0.5-1.0 pg/mg in 3% of cells. The largest increases in predicted HCC were concentrated in the foothills to the east of Jasper National Park and were associated with increased distance to coal mines combined with increased variation in canopy cover. The results showed 54% of cells demonstrated either no change or a decrease in predicted HCC values. The largest decreases in HCC were associated with decreased variation in canopy cover, primarily occurring within Jasper National Park. The values of drivers of HCC included in the final model demonstrated changes from 2004 to 2014. In cells where there were HCC data for both years, the mean distance to coal mines decreased by 63.8% from 2004 to 2014, which  Table 5. Standardized coefficients (b), standard errors (SE), 95% confidence intervals (CI), and R 2 values (1-deviance/null deviance) for the parameters used in the most supported models for predicting hair cortisol concentration in samples collected from grizzly bears.

DISCUSSION
There is urgent need for biological markers of physiological function, such as cortisol, that can be used to monitor the physiological responses of individuals to their rapidly changing environments (Wikelski and Cooke 2006). West-central Alberta, Canada, is particularly susceptible to landscape change from the ongoing forestry and oil and gas industry in this area (Linke and McDermid 2012), with a total of 33.4% of forest management areas in Alberta leased to industry in 2001 (Timoney and Lee 2001). While resource extraction activities (forest harvesting, mining, and oil/gas well-sites) often result in an increase in roads and railways, and ultimately human access, these disturbances have also been linked to the presence and abundance of food resources for grizzly bears (Nielsen et al. 2004c, Roever et al. 2008, Stewart et al. 2013, resulting in individuals selecting for these habitats (Cristescu et al. 2016b). Anthropogenic disturbances on the landscape may be creating additional food resources and thus lowering HCC in individuals (Bourbonnais et al. 2013, Bryan et al. 2013, Stewart et al. 2013. Between 2004 and 2014, the majority of the study area where hair was collected demonstrated a reduction in HCC, while smaller areas showed elevated HCC. Potential drivers of changes in hair cortisol concentration were related to anthropogenic disturbance, food resource availability, and terrain conditions, with the majority of associations being related to anthropogenic variables.
In the current study, the stress response of grizzly bears was associated with variables representative of anthropogenic disturbance, with a reduction in HCC as the density of oil and gas well-sites increased and an elevation in HCC as the distance to coal mines increased. While bears have been shown to avoid active well-sites and oil and gas features during the summer and fall seasons (Laberee et al. 2014), they were also more likely to be detected near oil and gas well-sites in our study area, especially if food availability was low in the surrounding areas . Grizzly bears have also shown selection for well-sites , and this may be because these areas provide a concentrated source of bear foods at the well-site edges. Similar to road construction, forests are cleared in order to construct oil and gas well-sites, which allows for herbaceous food resources to grow (Roever et al. 2008). These food resources, including Equisetum spp., dandelions, and clover, are a significant portion of grizzly bear diet in the foothills of west-central Alberta (Munro et al. 2006). Similarly, grizzly bears have been shown to select for mining areas (Cristescu et al. 2016a(Cristescu et al. , 2016b, likely because these areas provide herbaceous material, such as forbs and graminoids, as forage (Cristescu et al. 2015). While these landscape features may provide additional food resources, these areas also increase human access and edge features, which have been shown to negatively impact grizzly bear survival rates (Nielsen et al. 2004b) with grizzly bear mortality rates linked to road and railway density in Alberta (Benn and Herrero 2002, Nielsen et al. 2004b, Boulanger and Stenhouse 2014. Over a 60-yr period, increased oil and gas extraction activities resulted in a decrease in forest crown and an increase in edge density in the boreal forests of Alberta, Canada (Pickell et al. 2015), which may be influencing food resources for grizzly bears. Crown closure is a measure of the proportion of ground (or sky, depending on view) obscured by tree canopy (McDermid et al. 2009), thus influencing forest understory vegetation abundance. Crown closure has been used to evaluate food resource availability in several grizzly bear studies (Munro et al. 2006). HCC levels have been shown to increase as forest crown closure increased (Bourbonnais et al. 2013); however, the current study found that the standard deviation of crown closure was a significant predictor of HCC and HCC increased as forest crown closure became more variable. Grizzly bear food resources tend to be on edges, which are typically areas with low canopy cover (Roever et al. 2008). Furthermore, young (≤20 yr) forest edges have been found to facilitate an increase in grizzly bear food supply (Nielsen et al. 2004c, Larsen et al. 2019. Females have Fig. 4. Predicted hair cortisol concentrations in 2004 and 2014 and over time from landscape variables in the most supported model (anthropogenic disturbance and food resource availability). Only this model was supported (DAIC c = 0.00) when year was included in the final candidate models; however, predicted hair cortisol concentrations were generated from the final model without year. Therefore, unidentified change associated with sample year is likely present. demonstrated lower HCC levels in areas that had a higher density of recent forest harvest blocks (Bourbonnais et al. 2013); therefore, the increase in grizzly bear food supply facilitated by creating new forest edges may explain the relationship between reduced glucocorticoid concentrations and anthropogenic disturbance. As the variation of crown closure changes depending on the environment type (i.e., dense forest or clearing), food resources may also be variable and thus influencing glucocorticoid concentrations in individuals.
Compound topographic index represents soil wetness as an aspect of topographic heterogeneity. In the current study, HCC was reduced in areas with increased CTI, which contrasts previous studies showing that low stress levels in male grizzly bears were associated with low soil wetness (Bourbonnais et al. 2013). Areas with higher CTI are more likely to be a riparian zone or wetland, and grizzly bears have been shown to select for such areas (Phoebus et al. 2017). Therefore, these areas are likely more productive and facilitate food resources (Nielsen et al. 2004c). Furthermore, CTI has been shown to be positively correlated with rooting resources and negatively associated with animal matter, while herbaceous matter was influenced by both high and low CTI values (Nielsen et al. 2010).
During the study period (2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014), the grizzly bear population in the study area doubled, increasing from 36 bears in 2004 to 71 bears in 2014 . Overall, HCC declined between 2004 and 2014, which may indicate that the landscape in our study area is not at carrying capacity for this species, and therefore, stress levels appear to not be density dependent. Reduced stress levels may have also positively influenced reproductive performance in this population, as the link between elevated levels of stress and reduced reproductive output has been reported in other species (Massey et al. 2016). However, the relationship between baseline glucocorticoid concentrations and population fitness in free-ranging animals remains unclear due to potential intrinsic and ecological factors that influence the physiological stress response and ultimately, confound associations between potential drivers of stress and HPA activity (Bonier et al. 2009, Dantzer et al. 2014. For example, other factors that were not measured in the current study that may influence and/or confound the stress response in grizzly bears include human recreational activity (Ladle et al. 2019), hunting suspension in 2006 (Wielgus andBunnell 1994, Festa-Bianchet 2010), and changes in snow cover (Berman et al. 2019). Moreover, interactions between predictor variables in the current study were not tested, but may further explain changes in HCC. This study was also limited by the absence of HCC values between 2004 and 2014, which could further influence the relationships between stress and landscape variables. Inter-annual variation related to anthropogenic disturbance, food resource availability, and terrain conditions may explain the differences detected in the current study. For example, elevated concentrations of hair cortisol in 2004 compared to 2014 may indicate adaptive responses by individuals to variation in natural resources, including food availability. Therefore, while individuals may be chronically stressed, indicated by increased glucocorticoid concentration, they are initiating adaptive responses that ultimately increase fitness (Boonstra 2013).
In combination with previous and ongoing research being conducted on grizzly bears in Alberta, this study provides unique insight into non-invasively measuring glucocorticoid concentrations in wildlife. The current study identified and quantified landscape changes that are associated with glucocorticoid concentrations in grizzly bears. HCC was influenced by variables representing anthropogenic disturbance (oil and gas well-sites and coal mines), food resource availability (crown closure), and terrain conditions (compound topographic index). This study demonstrates the importance of tracking HCC over time, because even where there was minimal landscape change (i.e., protected areas), there were changes in stress levels. Therefore, we suggest that hair samples should be collected annually to track these changes. Hair should be collected via barbed wire hair snags in order to monitor population-level changes in HPA activity, as individual HCC may be influenced by short-term stressors, such as capture and restraint ). This research provides further understanding of the spatial and temporal impact of landscape change on grizzly bear populations. Physiological data, such as hair cortisol concentration, can be collected non-invasively and is an important approach for monitoring wildlife. This method provides reliable information regarding potential drivers of glucocorticoid concentrations over time and at the population level. Furthermore, it is applicable to the monitoring of HPA activity in other wildlife species on changing landscapes.