Conservation status of polar bears (Ursus maritimus) in relation to projected sea-ice declines

Loss of Arctic sea ice owing to climate change is the primary threat to polar bears throughout their range. We evaluated the potential response of polar bears to sea-ice declines by (i) calculating generation length (GL) for the species, which determines the timeframe for conservation assessments; (ii) developing a standardized sea-ice metric representing important habitat; and (iii) using statistical models and computer simulation to project changes in the global population under three approaches relating polar bear abundance to sea ice. Mean GL was 11.5 years. Ice-covered days declined in all subpopulation areas during 1979–2014 (median −1.26 days year−1). The estimated probabilities that reductions in the mean global population size of polar bears will be greater than 30%, 50% and 80% over three generations (35–41 years) were 0.71 (range 0.20–0.95), 0.07 (range 0–0.35) and less than 0.01 (range 0–0.02), respectively. According to IUCN Red List reduction thresholds, which provide a common measure of extinction risk across taxa, these results are consistent with listing the species as vulnerable. Our findings support the potential for large declines in polar bear numbers owing to sea-ice loss, and highlight near-term uncertainty in statistical projections as well as the sensitivity of projections to different plausible assumptions.


Introduction
Polar bears (Ursus maritimus) depend on sea ice for most aspects of their life history [1]. Anthropogenic climate change is the primary threat to the species because, over the long term, global temperatures will increase and Arctic sea ice will decrease as long as atmospheric greenhouse gas concentrations continue to rise [2,3]. The global population of approximately 26 000 polar bears [4] is divided into 19 subpopulations, which are grouped into four ecoregions reflecting sea-ice dynamics and polar bear life history (figure 1; [5]). The subpopulations currently exhibit variable status relative to climate change [6]. Two have already experienced sea-ice-related demographic declines [7,8]. Others show signs of nutritional stress [9], have been reported as stable or productive [10] or have unknown status owing to deficient data [11].
Methods to forecast the effects of continued sea-ice loss on polar bears have included structured elicitation of expert opinion [12], Bayesian network models evaluating cumulative stressors [5], demographic projections for individual subpopulations [8] and mechanistic models linking vital rates to environmental factors [13]. To date, there has been no global assessment of polar bear abundance data relative to sea ice. We explored future changes in mean global population size (MGPS) for the species, using population projections under three approaches. Approach 1 reflected the hypothesis that environmental carrying capacity (K) is directly proportional to the availability of sea ice. Approaches 2 and 3 estimated relationships between changes in sea ice and observed changes in polar bear abundance. We evaluated projection outcomes, over three polar bear generations, relative to thresholds for threatened categories under criterion A3 of the IUCN Red List of Threatened Species (hereafter Red List; [4,14]). The scientific basis of Red List categories is discussed by Mace et al. [15].

Methods
Projection timeframes can incorporate biological differences across species by referencing to generation length (GL, the average age of parents of the current cohort; [14]). We estimated GL as the mean age of adult female polar bears with new cubs based on live-capture data from 11 subpopulations. Females with 1 yearold cubs in year t þ 1 were counted as pseudo-observations in year t. Variation in GL was evaluated, using a bootstrap procedure (electronic supplementary material, table S1).
Satellite data of sea-ice concentration were collected between years 1979 -2014 to develop an index of K for polar bears (figure 2; [16]). Within each of the 19 subpopulation areas, daily sea-ice area was calculated by summing the product of ice concentration and grid cell area over all 25 Â 25 km grid cells with concentration more than 15%. We then determined the midpoint between summer-minimum and winter-maximum ice areas, and calculated the metric ice as the number of days per year that ice area was above the midpoint (i.e. the number of 'icecovered' days). Mean values of ice were projected forward, using linear models, which facilitated projections at the spatial scale of polar bear subpopulations (electronic supplementary material, table S2).
We used population projections to evaluate changes in MGPS between the years 2015 and (2015 þ 3 Â GL) based on three approaches relating ice to estimates of subpopulation abundance (N; electronic supplementary material, table S3). Approach 1 assumed a one-to-one proportional relationship between ice and N. Approaches 2 and 3 estimated linear relationships between ice and proportional changes in N, and used these relationships to predict future values of N as a function of projected ice. Approach 2 estimated a global ice-N relationship based on a maximum of two estimates of N per subpopulation, separated by at least 10 years, which were available for seven subpopulations. Approach 3 estimated a separate ice-N relationship for each polar bear ecoregion using a dataset that was similar to approach 2 but included longer time series of N available for four subpopulations. All approaches rsbl.royalsocietypublishing.org Biol. Lett. 12: 20160556 assumed that changes in N were mediated primarily through changes in K or density-independent habitat effects, and that the ratio N/K was stable relative to other factors [17]. These assumptions were established on the basis that polar bears depend fundamentally on sea ice, that sea-ice changes represent the main source of habitat modification for the species [5], and that other   rsbl.royalsocietypublishing.org Biol. Lett. 12: 20160556 potential stressors are either secondary (e.g. contaminants; [5]) or have been managed (e.g. harvest; [6]) for most subpopulations in recent decades. Projected subpopulation-specific changes in N were scaled to changes in MGPS, using the most recent estimate of N for each subpopulation. Based on 62 500 stochastic projections, we calculated the most likely change in MGPS over three generations. In addition, the probabilities of exceeding 0%, 30%, 50% and 80% reduction thresholds were generated following Red List guidelines [14]. We performed computations in R [18], using the package 'arm' [19] to simulate uncertainty in model coefficients. Data and projection methods are described fully in the electronic supplementary material.

Results and discussion
The mean subpopulation-specific estimate of GL was 11.5 years (approx. 5th and 95th percentiles ¼ 9.8 and 13.6, respectively) based on 3374 observed reproductive events (electronic supplementary material, table S1). Projections were performed using GL ¼ 11.5 and 13.6 years to reflect variation in GL and approximate natural GL. We did not apply the lower fifth percentile, because harvest likely shortened several empirical estimates of GL [14]. We simulated per cent change in MGPS for six scenarios representing two values of GL and three approaches relating ice and N (table 1). Using GL ¼ 11.5 years, the most likely values for per cent change in MGPS over three generations were 230%, 24% and 243% for approaches 1, 2 and 3, respectively. Across scenarios, the estimated median probabilities of reductions greater than 30%, 50% and 80% in MGPS were 0.71 (range 0.20-0.95), 0.07 (range 0-0.35) and less than 0.01 (range 0 -0.02), respectively.
Our analyses highlight the potential for large reductions in MGPS as climate change and sea-ice loss continue [20] over the next three polar bear generations. Approach 1 was based only on projected changes in habitat, a common method when population data are lacking [14]. Approach 2 estimated a global ice-N relationship that was near 0 and not statistically significant (estimated slope coefficient [b] , 0.001, s.e. ¼ 0.005; electronic supplementary material, table S4). This finding reflects variability in current subpopulation status, uncertainty in estimates of N and the lack of empirical evidence for sea-ice mediated changes in global abundance over recent decades [6]. Approach 3 estimated a separate ice-N relationship for each ecoregion. Relationships were positive at a significance level of 0.01 for the seasonal (b ¼ 0.013, s.e. ¼ 0.002) and divergent ecoregions (b ¼ 0.032, s.e. ¼ 0.009), reflecting observed correlations between declining sea ice and declining abundance (electronic supplementary material, table S4). Relationships were not significant for the convergent (b ¼ 20.008, s.e. ¼ 0.009) and archipelago ecoregions (b ¼ 20.029, s.e. ¼ 0.030). Although approach 3 reflected regional variability in sea-ice dynamics and polar bear ecology, it was strongly influenced by several well-studied subpopulations and did not reflect finer-scale variation. For example, within the divergent ecoregion, multiple estimates of N were available for the declining Southern Beaufort sea subpopulation [7], but not for the Chukchi sea subpopulation, which inhabits a more biologically productive region and has exhibited high recruitment despite sea-ice loss [9].
Our projections (table 1) are broadly consistent with expert opinion [12] and Bayesian network model forecasts [5], although methodological differences preclude direct comparison (see electronic supplementary material). Following the Red List guidelines for risk tolerance ( [14]; §3.2.3), the high probability of reductions more than 30% in MGPS, and low probability of reductions more than 50%, were consistent with a categorization of vulnerable (i.e. facing a high risk of extinction in the wild; [4]). Our use of statistical models required estimating few parameters, consistent with sparse data available for Arctic marine mammals [21], and propagated the effects of assumptions on model outcomes in a transparent manner. Future global population assessments could explore the use of hierarchical models [22], integrate data from multiple sources [23], model population processes (e.g. density-dependent interactions between harvest and habitat loss; [17]), consider cumulative effects on polar bear health [24] or consider nonlinear or spatial responses [25].
Data accessibility. The datasets supporting this article have been uploaded as part of the electronic supplementary material. 1.00 0.88 0.35 0.02 a Approach 1 assumed a one-to-one proportional relationship between sea ice and abundance. Approaches 2 and 3 estimated global and ecoregion-specific relationships between sea ice and empirical estimates of abundance, respectively. Results from each approach are shown for the mean and 95th percentile of estimated GL. rsbl.royalsocietypublishing.org Biol. Lett. 12: 20160556