Marine birds in Alaska, USA, have been monitored systematically for more than 4 decades, and yet it remains unclear why some populations have increased while others have declined. We analyzed the population dynamics of 5 seabird species—Black-legged (Rissa tridactyla) and Red-legged kittiwakes (R. brevirostris), Common (Uria aalge) and Thick-billed murres (U. lomvia), and Tufted Puffins (Fratercula cirrhata)—across 4 decades in Alaska. We tested hypotheses that each species' carrying capacity varied continuously through time with climate cycles and/or in response to habitat covariates. Using an information-theoretic approach, we evaluated competing candidate stochastic growth models of each species' annual rate of change, incorporating various environmental covariates. The North Pacific Index and Pacific Decadal Oscillation were the most important climatic covariates across the whole of Alaska, where the former generally was negatively related to rates of population change, and the latter positively related. Across the 40-yr time series, we found slight decreases in zooplankton (i.e. krill) concentrations across the Gulf of Alaska, and significant increases in sea surface temperature across the Aleutian Islands. Kittiwakes showed the greatest level of sensitivity to these 2 environmental changes. Our results provide evidence that deteriorating secondary productivity (i.e. euphausiids) has contributed to declines of Black-legged Kittiwakes in the Gulf of Alaska. In contrast, the carrying capacity of murres has increased across the state, even in regions affected by warming waters and reduced productivity. These results suggest that kittiwakes act as indicators of detrimental impacts of climatic variability, whereas murres demonstrate resilience to such environmental change. Identifying the ecological factors that explain seabird population dynamics is necessary to understand the implications of climate and environmental change for long-term marine ecosystem dynamics.
INTRODUCTION
Marine birds are known to be good indicators of environmental change (Piatt et al. 2007, Sydeman et al. 2017), yet rarely are the causes of long-term population variability tested directly. This is largely due to the prohibitively expensive, but indispensable, nature of systematic long-term ecological monitoring (Lindenmayer et al. 2012). Short-term environmental perturbations can produce annual population fluctuations around a fixed or changing carrying capacity (K) that are best described by inverse density dependence (Goyert et al. 2017). Explaining changes in K requires a quantifiable assessment of its relationships with climatic and habitat variation. We examined the relationships between these 2 factors and the carrying capacity of 5 Alaskan seabird species—Black-legged (Rissa tridactyla) and Red-legged kittiwakes (R. brevirostris), Common (Uria aalge) and Thick-billed murres (U. lomvia), and Tufted Puffins (Fratercula cirrhata)—along a 4-decade time series (1974–2014). These 5 species occupy an array of foraging niches (Byrd and Williams 1993, Gaston and Hipfner 2000, Ainley et al. 2002, Piatt and Kitaysky 2002, Hatch et al. 2009): They feed on a spectrum of diverse (e.g., murres) to narrowly specific prey (e.g., Red-legged Kittiwakes), exploited at the sea surface (e.g., both species of kittiwake) or at depth (e.g., murres and puffins), in offshore (e.g., murres and kittiwakes) and/or nearshore ranges (e.g., puffins). We assessed the effects of long-term variability in the quality of foraging habitat on the number of individuals supported by this habitat (i.e. carrying capacity or quasi-equilibrium abundance, K; Garton et al. 2016).
Our objective was to determine whether changes in habitat or climate explained why murres and kittiwakes have both increased in the Arctic, but have responded differently from each other to conditions in southern Alaska, USA (Goyert et al. 2017). We also aimed to describe why puffins have declined in the Gulf of Alaska, consistent with their neighboring populations in the North Pacific (which are under review for listing pursuant to the Endangered Species Act; Hodum 2013). Links between climate and the predictability or availability of prey have been shown to affect the productivity and population trends of seabirds breeding on the Pribilof Islands of Alaska (Byrd et al. 2008, Renner et al. 2012, 2014). We evaluated how the population dynamics of our study species responded to (1) climatic indices (e.g., Pacific Decadal Oscillation and Arctic Oscillation), and (2) habitat variability, as described by sea surface temperature (SST) and secondary productivity (e.g., zooplankton concentrations; Hermann et al. 2013, M. A. Smith personal communication). Our analysis is the first attempt to explain the population dynamics of marine birds relative to these climatic indices and habitat variability across the entire extent of Alaska and its 4 large marine ecosystems: Aleutian Islands, Arctic, Eastern Bering Sea, and Gulf of Alaska (Figure 1; https://www.st.nmfs.noaa.gov/ecosystems/lme/index, M. A. Smith personal communication).
FIGURE 1.
Sites of monitored seabird colonies in Alaska, USA (yellow triangles), and their corresponding ecoregions. St. Lazaria Island is within the Gulf of Alaska, even though it is shown outside the scope of this large marine ecosystem (as delineated by https://www.st.nmfs.noaa.gov/ecosystems/lme/index; M. A. Smith personal communication). White polygons represent the spatial extent of extracted covariates: sea surface temperature (rectangles) and euphausiids (circles).
![i0004-8038-135-4-975-f01.tif](ContentImages/Journals/tauk/135/4/AUK-18-37.1/graphic/WebImages/i0004-8038-135-4-975-f01.jpg)
We tested the hypotheses that population rates of change would be positively related to the density of euphausiids (i.e. krill) and negatively related to increasing SST associated with a changing climate (Irons et al. 2008, Descamps et al. 2017). Our 5 study species feed primarily on young, cold-water forage fish, particularly ammodytids (e.g., sand lance [Ammodytes hexapterus]), osmerids (e.g., capelin [Mallotus villosus]), and gadids (e.g., pollock [Gadus chalcogrammus]; Drummond 2016). Red-legged Kittiwakes additionally specialize on several species of myctophid (lanternfish). Many seabirds, Thick-billed Murres and Black-legged Kittiwakes in particular, also feed facultatively on krill (Gaston and Hipfner 2000, Williams et al. 2008, Ito et al. 2010, Drummond 2016). While long-term data on the distribution and abundance of forage fish remain limited in scope, SST and krill biomass serve as proxies of foraging habitat quality. Therefore, secondary productivity conventionally functions as an indirect mechanistic link between primary productivity and predator abundance, where cold, nutrient-rich upwelling water positively influences the availability of prey (e.g., zooplankton and/or forage fish) via bottom-up trophic forcing (Shealer 2001, Frank et al. 2006). We expected to find variability in each species' population response to habitat and climate as a result of disparities in their foraging ecology. Thus, we provide insight into the relative vulnerability of Alaskan seabirds to long-term, ecosystem-wide change.
METHODS
Data
For our study period (1974–2014), we obtained indices of climatic (oceanic and atmospheric) variation from the National Oceanic and Atmospheric Administration's (NOAA) Bering Sea Climate database ( http://www.beringclimate.noaa.gov/data/). Seasonal and annual data were available for the Arctic Oscillation (AO), North Pacific Index (NPI, which also included the Aleutian Low Index), and Pacific Decadal Oscillation (PDO). To avoid model overparameterization (i.e. overfitting), we conducted preliminary analyses using an information-theoretic approach (see section below) to identify which of the 9 oceanic and atmospheric indices to include in a “global model,” under which we tested a nested candidate set of models (Burnham and Anderson 2002). Across the 5 species, the 4 most descriptive indices (e.g., best-supported within an information-theoretic framework, explained high amounts of variation in the data, consistent with R2 values) were: the annual AO (AOIa), the NPI spring anomaly (NPIAss), the summer PDO (PDOs), and the annual PDO (PDOa).
We limited the scope of our habitat covariates to the months when each species used foraging habitat in Alaskan waters (e.g., the breeding and fledging seasons, and post- and pre-migration, April–October). SST data layers were available from the NOAA Extended Reconstructed SST database (ERSST Version 4, based on in situ datasets; Huang et al. 2015); for our study area, we used the spatial extent of 50°–70°N and 170°–230°E. Across each ecoregion (and the entire state as a whole), we calculated annual means (°C) from monthly data spanning April–October, extracted from the 2° × 2° cell surrounding each colony (Figure 1, rectangles).
The NOAA Pacific Marine Environmental Laboratory (PMEL) provides a unique source of data on the inferred distribution and abundance of zooplankton across the geographic range of 50°–66°N and 170°–210°E (Hermann et al. 2013, M. A. Smith personal communication). These layers are available as downscaled (1) hindcasts (1974–2005) from the Coordinated Ocean-Ice Reference Experiments (CORE; Griffies et al. 2009) and (2) projections (2005–2014) from the Canadian Centre for Climate Modelling (CCCma) coupled global climate model (CGCM3; Flato et al. 2000). We selected the CCCma model over others implemented by PMEL for the Intergovernmental Panel on Climate Change because it showed intermediate levels of global warming (Hermann et al. 2013, M. A. Smith personal communication), with the others showing lower and higher warming trends. The CORE and CCCma models additionally provide information on SST, but they harbor more uncertainty than the ERSSTv4 data because they are derived from modeled output, as opposed to raw data (Huang et al. 2015). We validated our use of these sources for zooplankton data by testing the congruence among their SST layers and the ERSSTv4 data. Despite the uncertainty inherent in zooplankton data, they were extremely valuable in their ability to provide the only source of historical inference on the prey base across our entire study area. We calculated long-term trends in SST and zooplankton data using a Theil-Sen approach and Mann-Kendall nonparametric tests in package trend (Pohlert 2016) in R (R Core Team 2017).
To remain within the bounds of each species' foraging range (Smith et al. 2012), we extracted data covering a 90-km radius from each colony, at a spatial resolution of 10 × 10 km (∼0.1° × 0.1°; Figure 1, circles). We calculated annual mean euphausiid concentrations (mg C m−3) across each ecoregion (or the entire state) at the sea surface from inferred weekly data spanning April–September. While kittiwakes forage at the sea surface and murres exploit a large portion of the water column, we chose sea surface data to compare the effects of zooplankton across bird species.
The range of murres and Black-legged Kittiwakes covers all 4 ecoregions in Alaska, whereas Red-legged Kittiwakes colonize only a few islands in western Alaska, and Tufted Puffins occupy southern Alaska. We compiled a 4-decade-long time series of count data for our 5 study species from 33 monitored sites across Alaska's 4 large marine ecosystems (Goyert et al. 2017). These counts were based on breeding-age adults (murres and some kittiwake colonies), nests (2 Black-legged Kittiwake datasets), and burrow entrances (puffins; Dragoo et al. 2015). We reconstructed historical population sizes as inferred from the minimum number of nesting individuals across our study sites (Garton et al. 2011).
We used reconstructed annual population rates of change for each seabird species, rt = ln(Nt + 1) − ln(Nt), as the response variable in our modeling framework (Goyert et al. 2017). We tested for cross-correlations between these rates of change and the environmental predictors to determine whether any covariates had a lag effect (i.e. whether covariates were more likely to affect rates of change in a subsequent vs current year). To avoid collinearity among the environmental variables in our models, we removed one of any covariate pair with a Pearson's r > 0.6 or that resulted in a generalized variance inflation factor > 2 (Zuur et al. 2010). Of the 6 possible environmental predictors (AOIa, NPIAss, PDOa, PDOs, SST, euphausiids), we excluded the PDOs because it was collinear with the PDOa and NPIAss.
Information–theoretic Modeling Framework
We implemented discrete time and stochastic Ricker and Gompertz models (Dennis and Taper 1994), with density-dependent time lags of up to 2 yr (Garton et al. 2011, 2016, Goyert et al. 2017). For each species, we evaluated the global model, which included 3 climatic covariates (AOIa, NPIAss, PDOa), 2 habitat covariates (SST, euphausiids), a time trend in K (linear or quadratic), and 1 of the 6 mutually exclusive density-dependent terms (Nt, Nt − 1, Nt − 2, ln(Nt), ln(Nt − 1), ln(Nt − 2)). We conducted all our analyses in R (R Core Team 2017), and compared the relative support for models in the candidate set (n = 672 models nested within each global model) using Akaike's information criterion adjusted for small samples (AICc). We determined the most competitive models (ΔAICc < 2; Burnham and Anderson 2002) nested within the global model(s) for each species using package MuMIn (Bartoń 2015). To better represent the relative importance of covariates beyond this cutoff value, we also computed the sum of Akaike weights across all 672 candidate models for each species (in Supplementary Material Tables S1 and S2 we report the candidate set of models with cumulative Akaike weights < 0.95). We avoided model-averaging the parameter estimates, which would have misrepresented the effect size of each mutually exclusive density-dependent term; that is, we required a parameter estimate for either the Ricker (e.g., Nt), Gompertz (e.g., ln(Nt)), or time-lagged (e.g., Nt − 1) model term, as opposed to an estimate of a combination of these density-dependent terms in a single resulting averaged model.
RESULTS
Across our 4-decade-long time series, SST increased significantly in the Aleutian Islands (y = 6.6 + 0.019t, calculated using a Theil-Sen estimator) by, on average, an annual increment (time t) of 0.019°C from the mean temperature of 6.6°C in 1974 to 7.3°C in 2014 (Figure 2); the Mann-Kendall test confirmed that this increase was significant (S = 248, P = 0.006), and the residuals of this regression were not autocorrelated. The mean temperature in the first 2 decades of the Aleutian Island time series (6.8°C, including 1994) was significantly different from that in the last 2 decades (7.2°C; t38 = −2.9, P = 0.006). We found a significant positive correlation between SST from the ERSSTv4 and the derived CORE and CCCma models (Pearson's r = 0.39, P = 0.01), which supported our use of inferred zooplankton data from the latter source. This metric of euphausiid density at the sea surface decreased slightly in the Gulf of Alaska (y = 4.6 – 0.0072t; Mann-Kendall S = –174, P = 0.05; Figure 3); the residuals of this regression were not autocorrelated. Similarly, the mean euphausiid density in the first 2 decades in the Gulf of Alaska time series was significantly greater than that in the last 2 decades (t36 = 2.0, P = 0.05).
FIGURE 2.
Annual April–October mean National Oceanic and Atmospheric Administration (NOAA) extended reconstructed sea surface temperature data (ERSST Version 4) across Alaska (AK), USA. The spatial extent covered 2° cells that overlapped seabird colonies in the Aleutian Islands (AI), Arctic (ARC), Eastern Bering Sea (EBS), and Gulf of Alaska (GOA) ecoregions. The regression line indicates significant increases in sea surface temperature (SST) across the 40-yr time series in the AI.
![i0004-8038-135-4-975-f02.tif](ContentImages/Journals/tauk/135/4/AUK-18-37.1/graphic/WebImages/i0004-8038-135-4-975-f02.jpg)
FIGURE 3.
Annual April–September mean euphausiid concentration at the sea surface across Alaska (AK), USA. The spatial extent covered cells within a 90-km radius of seabird colonies in the Aleutian Islands (AI), Arctic (ARC), Eastern Bering Sea (EBS), and Gulf of Akaska (GOA) ecoregions. Data were sourced from the Coordinated Ocean-Ice Reference Experiments model (1974–2004) and from the Canadian Centre for Climate Modelling coupled global climate model (2005–2014). The regression line indicates a decrease in euphausiid concentration across the 40-yr time series in the GOA.
![i0004-8038-135-4-975-f03.tif](ContentImages/Journals/tauk/135/4/AUK-18-37.1/graphic/WebImages/i0004-8038-135-4-975-f03.jpg)
Preliminary analyses of auto- and cross-correlation functions indicated that there were no strong time lags (P > 0.05 for cross-correlation coefficients) between population rates of change and environmental covariates, and that only the PDOa was autocorrelated, at a 1-yr lag (autocorrelation coefficient = 0.482, P = 0.002). Such temporal autocorrelation in the covariates, along with the structure of the model which allowed for density dependence and time lags, helped to accommodate nonindependence in the response variable (rates of change). This resulted in successful model diagnostics: The models showed no evidence of nonindependence, as diagnosed by the lack of significant autocorrelation in the top population dynamics models and in the residuals of the aforementioned covariates regressed against time (most notably those covariates that showed significant increases or decreases through time, e.g., SST in the Aleutian Islands and euphausiids in the Gulf of Alaska). Across the entire extent of Alaska, population rates of change were generally positively related to the AOIa (Common Murres), PDOa (Tufted Puffins, Red-legged Kittiwakes, and Thick-billed Murres), and euphausiids (Black-legged Kittiwakes and Tufted Puffins). They were negatively related to the NPIAss (Common and Thick-billed murres) and SST (Common Murres; Table 1, Appendix Figure 5). These relationships varied by species at the ecoregion level (Table 2).
TABLE 1.
Parameter effects of the top kittiwake, murre, and puffin population models across Alaska, USA. We evaluated density-independent and density-dependent Ricker (linear) and Gompertz (log-linear) models, using population reconstructions for Black-legged (BK) and Red-legged kittiwakes (RK), Common (CM) and Thick-billed murres (TM), and Tufted Puffins (TP). For each species (Spp), the relative importance of covariates is shown in the first row, weighted across the candidate set of 672 models. Below this, each row represents one of the most competitive models (difference from the top model in Akaike's information criterion corrected for small sample size [ΔAICc] < 2), with its Akaike weight (wi), R2 value, intercept (μ), and covariate effects. The covariates qualified as: (1) mutually exclusive density-dependent terms, where Nt is the population size at time t and negative effects suggest density dependence, while Nt − 1 and Nt − 2 indicate a time lag of dependence on population sizes up to 2 yr past; (2) climatic covariates, including the annual Arctic Oscillation (AOIa), spring North Pacific anomaly (NPIAss), and annual Pacific Decadal Oscillation (PDOa); (3) habitat covariates, including sea surface temperature (SST) and euphausiid concentration (Euph); and (4) linear and nonlinear (quadratic) effects of time (t) on the rates of change, which suggest an increasing (+) or decreasing (−) trend in carrying capacity (K) with year. For example, the top candidate model for BK contains only the terms ln(Nt) (negative Gompertz density-dependence), t, and t2 (decreasing time trend), and the second most competitive model is the same but with a positive effect of euphausiids also included. The top null density-dependent models (with a constant time trend in K) are italicized; those with ΔAICc > 2 (i.e. significantly worse models) are listed below the dotted line. Models highlighted in bold font are depicted in Figure 4; refer to Supplementary Material Table S1 for candidate sets of models with cumulative Akaike weights < 0.95.
![i0004-8038-135-4-975-t01.tif](ContentImages/Journals/tauk/135/4/AUK-18-37.1/graphic/WebImages/i0004-8038-135-4-975-t01.jpg)
FIGURE 4.
Most competitive covariate models of population dynamics across monitored seabird colonies in Alaska, USA, for (A–B) kittiwakes, (C–D) murres, and (E) puffins. Reconstructed population sizes, N(t), are shown as black solid lines with 95% confidence intervals (gray shading). Solid red lines show model-estimated population sizes, N̂(t), and blue dashed lines indicate model-estimated carrying capacity, K̂(t). Note that 2 of the data points at either end of each time series were included in tests of time lags, and the last data point was included in rate of change calculations. The models shown here were the best of those that contained at least one environmental covariate (highlighted in bold font in Table 1), including: (A) euphausiids for Black-legged Kittiwakes, (B) PDOa for Red-legged Kittiwakes, (C) NPIAss and SST for Common Murres, (D) NPIAss for Thick-billed Murres, and (E) PDOa for Tufted Puffins.
![i0004-8038-135-4-975-f04.tif](ContentImages/Journals/tauk/135/4/AUK-18-37.1/graphic/WebImages/i0004-8038-135-4-975-f04.jpg)
TABLE 2.
Parameter effects of the top kittiwake, murre, and puffin population models across the Alaskan ecoregions (Reg) colonized by each species: Aleutian Islands (AI), Arctic (ARC), Eastern Bering Sea (EBS), and Gulf of Alaska (GOA). See Table 1 for more details; refer to Supplementary Material Table S2 for candidate sets of models with cumulative Akaike weights < 0.95.
![i0004-8038-135-4-975-t02.tif](ContentImages/Journals/tauk/135/4/AUK-18-37.1/graphic/WebImages/i0004-8038-135-4-975-t02.jpg)
Population models that assumed that the carrying capacity of Black-legged Kittiwakes varied on an annual basis with time, climate, and/or habitat performed significantly better than the null model (density dependence without any covariates, i.e. a constant trend in K; Tables 1 and 2). The time trend in K captured the most variation in their population dynamics, particularly in the Arctic, Eastern Bering Sea, and Gulf of Alaska ecoregions. K increased in the Arctic and decreased elsewhere. Population rates of change were negatively related to euphausiids in the Eastern Bering Sea and positively related elsewhere. In the Aleutian Islands, rates of change were positively related to the PDOa and negatively related to the AOIa and NPIAss.
Using covariates from the density-dependent models, including time where appropriate, explained 114% (Tufted Puffins) to 483% (Common Murres) more of the variation in growth rates for all of the species except Red-legged Kittiwakes, for which it increased the variation explained by only 13%, a negligible amount (see R2 values in Table 1 for these calculations, compared between the null model and best covariate model). For Red-legged Kittiwakes across their range in Alaska, the most descriptive models estimated either a constant K or a varying K that changed with time or the PDOa. A dynamic K was significantly better than the null model in both ecoregions occupied by this species (Aleutian Islands, Eastern Bering Sea). K decreased with time, euphausiids, and NPIAss in the Eastern Bering Sea. Increasingly positive rates of change and higher carrying capacities were associated with higher abundance of euphausiids in the Aleutian Islands. However, rates of change in this latter ecoregion were negatively related to the AOIa, NPIAss, and SST.
For Common Murres across Alaska, K varied with time, but SST, the NPIAss, and/or the AOIa described significantly more variation than the null model (constant K; except in the Eastern Bering Sea ecoregion). K decreased linearly with time in the Eastern Bering Sea, but increased elsewhere. Rates of change were inversely related to changes in SST in the Arctic and Gulf of Alaska ecoregions, and positively related to the PDOa in the Gulf of Alaska.
In the density-dependent models for Thick-billed Murres across Alaska, climate covariates were most important, where they (NPIAss, PDOa) predicted significantly more variation than the time trend alone (i.e. those covariates likely explained the perceived time trend). A dynamic K was significantly better than the null model for Thick-billed Murres in the Aleutian Islands, Eastern Bering Sea, and Gulf of Alaska ecoregions. Higher population rates of change were associated with higher estimates of euphausiids and the PDOa in the Gulf of Alaska and Aleutian Islands, and with lower values of SST in the Gulf of Alaska and Arctic. Rates of change in the Aleutian Islands and Arctic ecoregions were positively related to the NPIAss, but in the Eastern Bering Sea and Gulf of Alaska were negatively related to this index. They showed an inverse relationship with the AOIa in the Arctic and Eastern Bering Sea.
For Tufted Puffins, the most descriptive models estimated either a constant K or a K that varied linearly with time, the PDOa, or euphausiids. A dynamic K was significantly better than the null model in the Gulf of Alaska. K increased with time and the PDOa in the Aleutian Islands and decreased with time and the PDOa in the Gulf of Alaska.
In summary, a carrying capacity that varied annually with climatic and/or environmental covariates was most descriptive of the population dynamics of all 5 species (Figure 4). The order of importance of each covariate, aggregated across all 5 species and 4 ecoregions, suggested that the Gompertz model with a time trend in K was most predictive (i.e. retained the most weight), followed by models with influences from euphausiids, the PDOa, AOIa, NPIAss, and SST. Similarly, across the entire state, the most important covariates were the NPIAss and PDOa, followed by euphausiids, SST, and the AOIa.
DISCUSSION
Consistent with our hypotheses that population rates of change would be negatively related to climatic changes in SST and positively related to densities of krill (euphausiids), rates of change generally declined in response to increasing SST, and decreased in response to decreasing secondary productivity (except in the Eastern Bering Sea ecoregion, where carrying capacities were inversely related to krill). Slight decreases in inferred krill density occurred in the Gulf of Alaska, which may explain declines of Black-legged Kittiwakes in this region. Along with the PDO, krill was positively related to rates of population change of Tufted Puffins, but did not appear to contribute to puffin declines in the Gulf of Alaska. Inconsistent with our hypothesis, negative relationships between krill and kittiwake rates of change occurred in the Eastern Bering Sea. The Aleutian Island chain was the only region to show a significant increase in SST and to contribute to negative effects in Red-legged Kittiwake populations. Thus, both species of kittiwake demonstrated the greatest vulnerability to significant increases in SST (Aleutian Islands) and decreases in secondary productivity (Gulf of Alaska). Warming waters over the past 4 decades seem likely to have affected the predictability or availability of their prey base.
Inverse relationships between Common Murre population rates of change and SST occurred across Alaska, particularly in 2 of the regions (Gulf of Alaska and Arctic) where their populations grew. These relationships are consistent with other findings on the circumpolar dynamics of Common Murre populations (Irons et al. 2008). Thick-billed Murres demonstrated similar relationships (Byrd et al. 2008, Renner et al. 2012, 2014), but showed higher sensitivity to the NPIAss. Their population rates of change were negatively related to the NPIAss (analogous to the strength of the Aleutian Low Index), consistent with published effects on the productivity of other alcids (auklets [Aethia spp.]; Bond et al. 2011). Given that the carrying capacity of murres has increased over the last 4 decades, despite their inverse relationship with ocean warming, these results suggest that they may be more resilient to the observed environmental changes than the other species (Goyert et al. 2017). The high, positive rates of change that we documented in Common Murres may have represented a steady increase in K, despite drops in their population sizes as a response to environmental perturbation (Southwood 1981). This ability of murres to adapt to a volatile environment is supported by evidence of stable breeding success across several decades (Dragoo et al. 2015, Sydeman et al. 2017).
Compared with murres, the increased sensitivity of kittiwakes to degraded habitat in southern Alaska (e.g., from reduced secondary productivity and warmer waters) likely can be attributed to differences in their foraging ecology. Black-legged Kittiwake rates of population change responded positively to habitat (e.g., euphausiids) and negatively to climate (e.g., the PDO) in the Arctic and Gulf of Alaska regions. In other words, their populations decreased in the Gulf of Alaska, where secondary productivity deteriorated, suggesting that kittiwakes rely on habitat that supports healthy euphausiid populations. Meanwhile, the carrying capacity of murre populations increased in these regions and elsewhere. The ability of murres to adapt to the environmental variables examined in our study likely had to do, at least in part, with their foraging flexibility. Murres are generalists with flexible diets that exploit the entire water column, adjust their foraging time budgets, and maintain high reproductive success when food supplies shift (Burger and Piatt 1990, Harding et al. 2007). In contrast, both kittiwake species face many constraints as surface feeders with rigid diets and limited energy budgets (Furness and Tasker 2000, Einoder 2009). Forage fish, in particular the young age classes depredated by seabirds, are highly sensitive to climatic cycles and shifts in environmental regimes due to their high demands for zooplankton (Batchelder et al. 2013). Unpredictable changes in the availability of forage fish can result in spatial or temporal (i.e. phenological) mismatches with seabird distributions (Burthe et al. 2012, Staudinger et al. 2013). The resulting prey limitation likely drove the population declines that we observed in kittiwakes. Our study suggests not only that kittiwakes are susceptible to degraded foraging habitat and climate change, but also that they are useful as marine indicators (Sydeman et al. 2017).
Given the vast spatial and temporal extent of our study, the data layers available from PMEL on the inferred historical distribution of zooplankton in Alaska were a valuable and unique source of information, despite their limitations. The marine ecosystem in Alaska is highly dynamic, and using modeled data as covariates adds a degree of uncertainty to our conclusions. Furthermore, although we minimized collinearity (i.e. nonindependence) among our model covariates, we acknowledge that climatic indices are defined by several interacting factors in the environment. Long-term studies that collect raw data on oceanographic and climatological factors throughout the North Pacific are crucial to avoid further uncertainty in predicting how wildlife populations will respond to climate change (Lindenmayer et al. 2012). Our attempt to explain the population dynamics of seabirds with respect to habitat and climatic variability provides preliminary evidence of effects of multidecadal processes. However, further analysis is required to document the influence of relationships (e.g., teleconnections) between the regional and global processes driving changes in the North Pacific Ocean (Di Lorenzo and Mantua 2016). This would help to identify the nature of the trophic forcing that may serve as the link between warm waters, low availability of prey, and low breeding success in seabirds, particularly as these factors influence declines in kittiwake populations (Frank et al. 2006, Sydeman et al. 2017).
By including environmental factors in our models, we were able to estimate and explain population change in association with variation in climate and habitat quality. Further studies should examine the demographic factors that influence population change in seabirds. For example, kittiwake breeding success has been shown to vary with SST (Frederiksen et al. 2007) and prey availability (Robert et al. 2002, Frederiksen et al. 2005). Their survival, on the other hand, has been shown to be relatively constant, and likely shapes the long-term population trends of these long-lived species (as in other larids; Monaghan et al. 1989, Oro et al. 1999). While moderate changes in food supply may influence breeding success, only severe changes in prey availability are likely to affect survival, by failing to replenish body reserves (Cairns 1987, Oro and Furness 2002). Such extremely poor conditions can influence the mortality of seabirds outside the breeding season (Orben et al. 2015). Thus, breeding success most likely contributes to the annual population dynamics of our 5 study species, whereas survival probably drives their long-term population trends.
Our results provide evidence of strong bottom-up forcing in Alaskan ecosystems, likely mediated by high primary and secondary productivity during negative (i.e. cold) PDO regimes (Pinchuk et al. 2008) and low productivity during warming events associated with anomalies in the Aleution Low Index and NPI (Springer 1998). The extensive spatiotemporal coverage of our study restricted our analysis to examining bottom-up influences on our focal species. Top-down forcing also occurs, particularly at high latitudes, in low-productivity marine environments (Frank et al. 2006) or in unfavorably warm years (Hunt et al. 2002), but its source is generally specific to each colony (e.g., from predation by Larus spp., bald eagles [Haliaeetus leucocephalus], or Arctic fox [Vulpes lagopus] introductions; Regehr and Montevecchi 1997, Byrd et al. 2005). Commercial fisheries may impose some pressure on prey availability but, because they don't generally target forage fish, their impacts are more likely to be related to seabird bycatch (Žydelis et al. 2013). Competition with other top predators (e.g., marine mammals; Sinclair et al. 2008, Estes et al. 2009) could additionally affect populations of Alaskan seabirds. Kittiwakes may be more vulnerable than murres to the effects of competitors such as pink salmon (Oncorhynchus gorbuscha; Sydeman et al. 2017), particularly in warm PDO years (Mantua et al. 1997, Springer and van Vliet 2014). Future research should consider the potential of multispecies models to gauge the sensitivity of Alaskan seabird population dynamics to interspecific competition and/or facilitation.
This study supports the use of a biologically relevant scale of analysis (Scales et al. 2017) to identify environmental influences on wildlife population dynamics. Our ecoregion-level analyses were essential to revealing euphausiid concentrations as a significant predictor of rates of seabird population change. Extracting covariate data at a resolution consistent with the foraging range of our study species (e.g., ∼90 km; Smith et al. 2012) avoided diluting the effect of functional relationships. Had we analyzed these data only at the state-wide level, the climatic covariates would have emerged as the most important predictors because they did not vary among ecoregions. Understanding the indirect links between climate change and forage availability is crucial to identify the ecological drivers of population change and inform monitoring at regional scales. Isolating these factors is necessary to develop effective conservation policies, such as protections for Important Bird Areas (M. A. Smith personal communication), to enhance the success of population management and habitat mitigation under intensified environmental volatility.
ACKNOWLEDGMENTS
We are grateful for manuscript comments provided by Brie A. Drummond, Heather M. Renner, and Jeremy Littell. Courtney Conway (Idaho Cooperative Fish and Wildlife Research Unit) provided additional resources. We thank several data contributors, including: past and present Alaska Maritime National Wildlife Refuge (AMNWR) biologists (especially Marc Romano, Leslie Slater, Nora Rojek, Don Dragoo, Arthur Kettle, Jeff Williams, and Vernon Byrd); the U.S. Fish and Wildlife Service (USFWS) Office of Migratory Bird Management, particularly Dave Irons and Lisa Sheffield; Melissa Cady and the Alaska Peninsula/Becharof National Wildlife Refuge (NWR); Patrick Walsh and the Togiak NWR; Brian McCaffery and the Yukon Delta NWR; U.S. Geological Survey (USGS), especially Scott Hatch; Alaska Department of Fish and Game (ADF&G), particularly Ed Weiss; and Ed Murphy and Dean Kildaw at the University of Alaska, Fairbanks. The findings and conclusions in this article are those of the author(s) and do not necessarily represent the views of the U.S. Fish and Wildlife Service.
Funding statement: We appreciate the funding provided by the U.S. Fish and Wildlife Service Natural Resource Program Center. The Center did not have any input into the content of the manuscript and did not require approval of the manuscript before submission or publication.
Ethics statement: Our research followed the Guidelines to the Use of Wild Birds in Research.
Author contributions: H.F.G., E.O.G., and A.J.P. conceived the idea and design behind this research. H.F.G. analyzed the data and wrote the paper, with assistance from E.O.G. and A.J.P.
LITERATURE CITED
Appendices
APPENDIX FIGURE 5.
Strength of the effect from the highest-weighted climatic and/or habitat covariate(s) in the most competitive models of population dynamics across monitored seabird colonies in Alaska, USA, for (A) Black-legged Kittiwakes, (B) Red-legged Kittiwakes, (C–D) Common Murres, (E) Thick-billed Murres, and (F) Tufted Puffins. The y-axes here show residuals from a specification of the most competitive covariate model (highlighted in bold font in Table 1 and graphed in Figure 4), excluding the covariate shown on the x-axis. These residuals indicate the amount of variation captured by the specified model, to highlight the significance of the effect that each covariate had on modeling population rates of change. Covariates included: euphausiids (Euph), annual Pacific Decadal Oscillation (PDOa), North Pacific Index spring anomaly (NPIAss;), and sea surface temperature (SST).
![i0004-8038-135-4-975-f05.tif](ContentImages/Journals/tauk/135/4/AUK-18-37.1/graphic/WebImages/i0004-8038-135-4-975-f05.jpg)