Impact assessments of wind farms on seabird populations that overlook existing drivers of demographic change should be treated with caution

Population viability analyses (PVA) are now routinely used during the consenting process for offshore wind energy developments to assess potential impacts to vulnerable species, such as seabirds. These models are typically based on mean vital rates, such as survival and fecundity, with some level of environmental stochasticity (i.e., temporal variation). However, many species of seabird are experiencing population decline due to temporal (i.e., directional) trends in their vital rates. We assess the prevalence of temporal trends in rates of fecundity for a sentinel species of seabird, the black‐legged kittiwake Rissa tridactyla, and examine how accounting for these relationships affects the predictive accuracy of PVA, as well as the projected population response to an extrinsic threat. We found that temporal trends in kittiwake rates of fecundity are widespread, and that including these trends in PVA assessments dramatically influences the projected rate of population decline. We advocate that model validation become a prerequisite step in seabird PVA assessments to identify potential biases influencing the projected population response. We also argue that environmental factors driving current population dynamics need to be incorporated in PVA impact assessments as potential “worst‐case” scenarios. These findings have immediate application for improving and reducing uncertainty in impact assessments conducted as part of the consenting process for offshore wind energy developments.


| INTRODUCTION
In many countries, it is necessary to conduct a full evaluation of potential impacts to protected wildlife populations to gain consent for certain activities, such as biological resource use, land-use change, or energy production (e.g., EU Birds Directive 2009/147/EC, EU Habitat Directive 92/43/EEC, EU EIA Directive 85/337/EEC). A commonly used approach for inferring how populations may respond to these activities is to predict how vital ratessuch as survival or fecundity-are likely to be altered and then use a population viability analysis (PVA) to project the potential population response (Reed et al., 2002). Such analyses are now widely and routinely used to conduct impact assessments on protected populations using bespoke models and "generic" computer software packages, such as VORTEX (Lacy & Pollak, 2021) and RAMAS-GIS (Akçakaya, 2000).
The data needed to inform PVA include parameters describing a species life history strategy (i.e., mean vital rates and age of first breeding). To improve realism, many PVA studies also allow vital rates to vary temporally (i.e., stochastically) (Boyce, 1992;Beissinger & Westphal, 1998;Coulson et al., 2001;Reed et al., 2002;Pe'er et al., 2013). However, population persistence depends strongly on the combination of factors that contribute to temporal variation in vital rates (Bascompte et al., 2002;Lande, 1993;Melbourne & Hastings, 2008). For example, vital rate temporal variation can be decomposed into three processes: (1) long-term directional trends, which are especially pronounced in the context of climate change and habitat degradation; (2) density dependence, which is commonly detected in wild animal populations (e.g., Saether et al., 2016); and (3) annual variation around the long-term trend, which reflects true environmental stochasticity. Separating these processes typically requires long-term datasets that are financially and logistically demanding to collect (e.g., Horswill et al., 2018), as well as advanced statistical approaches that require substantial computing power, time and technical expertise (e.g., Miller et al., 2019). Consequently, most PVA studies apply bulk temporal variation to vital rates, and overlook how extrinsic pressures or density dependence may drive long-term population change (Chaudhary & Oli, 2020).
Globally, many species of seabird are experiencing long-term population decline (Birdlife International, 2018). This trend is also mirrored in the UK, where potential impacts from offshore wind energy development threaten to further accelerate rates of change. The UK anticipates a four-fold expansion in offshore wind energy within the next decade (i.e., from 8GW installed capacity in 2020 to 40GW in 2030, BEIS, 2020). Many species of UK seabird are potentially highly sensitive to offshore wind energy developments, namely through collisions with turbine blades (Bradbury et al., 2014;Furness et al., 2013). Therefore, PVA are now routinely used during the consenting process for wind energy developments to evaluate potential impacts to seabird populations. However, as yet, current dedicated software packages for conducting this analysis (e.g., Searle et al., 2019) do not readily account for existing temporal trends in vital rates, that is, ongoing drivers of demographic change. Omitting these relationships from PVA assessments could lead to unreliable predictions of population change and misinformed conservation and management decision making (Fordham et al., 2013;Keith et al., 2008). Consequently, there is a pressing need to examine how excluding directional trends in vital rates may influence projections feeding into seabird impact assessments and the consenting process for offshore wind energy.
In this study, we assess the prevalence of temporal trends in rates of fecundity using data collected from 70 colonies of black-legged kittiwakes Rissa tridactyla (hereafter kittiwake) across the UK and Ireland (Figure 1a). Kittiwakes are a long-lived seabird that is considered globally vulnerable due to widespread declines in population size and range (Birdlife International, 2019). We constructed a series of PVAs to examine how incorporating observed temporal trends in seabird rates of fecundity can influence model validation, as well as the projected population response to an extrinsic threat that decreases rates of survival. We parameterized these models using colony-specific vital rates from kittiwakes breeding on Skomer Island, Wales (51.74 N,5.30 W, Figure 1a). This colony was selected because seabird impact assessments for offshore wind energy developments were ongoing at the time of publication. We show that including this structure in PVA assessments dramatically alters the projected population response to extrinsic threats, thus highlighting the importance of considering existing drivers of population change when conducting assessments on protected species.

| Prevalence of temporal trends in rates of fecundity
To demonstrate the potential prevalence of temporal trends in seabird rates of fecundity, we collated colonyspecific data for kittiwakes breeding in the UK and Ireland. We focused our study on rates of fecundity because mark-recapture datasets for estimating kittiwake rates of survival are considerably fewer, primarily due to the increased data-collection requirements (Horswill et al., 2018. Here, fecundity rates are defined as the proportion of apparent nesting attempts that produce a successful fledgling each year. Data were extracted from the Seabird Monitoring Programme Database (SMP, 2020), and we included colonies with at least 8 years of data between 1986 and 2019. To identify colonies with a significant temporal trend in rates of fecundity, we used Poisson generalized linear models (GLM) fitted to each colony with a log link function. In these models, the response variable was set as the number of chicks alive at fledging for each colony, and the number of nests was included on the log scale as an offset. Time was included as a fixed effect running from one to the length of each dataset. This analysis was conducted in Program R (v. 4.0.2) (R Core Team, 2020). We also examined whether kittiwakes demonstrate a temporal trend in rates of fecundity at the national level by fitting a Poisson generalized linear mixed model (GLMM) to the complete dataset using the R statistical package "lme4" (v. 1.1.23, Bates et al., 2016). 2.2 | How do temporal trends in rates of fecundity influence projected population dynamics?
We constructed three PVA models as validation simulations to examine how observed temporal trends in rates of fecundity may influence the predictive accuracy of seabird PVAs (Table 1). We also constructed two PVA projection simulations to examine how incorporating observed temporal trends in rates of fecundity may influence the projected population response to a future threat (code for all PVA models is provided as Supporting Information Material S1). We parameterized these models using vital rate data from kittiwakes breeding on Skomer Island (Figure 1a). For the three model validation simulations (Models 1-3, Table 1), we ran each PVA to match the period of observed Projection models (2020-2050) Model 4

Model 5
Note: Model simulation period detailed in brackets. Shading indicates the model structure, for example Models 3 and 5 include the observed temporal trend in rates of fecundity for kittiwakes at Skomer Island.
abundance data for kittiwakes breeding on Skomer Island . The first of these PVA models was deterministic based on mean vital rates only (Table 1; Model 1). The mean vital rates for kittiwakes breeding on Skomer Island are estimated to be 0.87 for adult survival (ϕ t ) between 1989 and 2018 (Horswill et al., unpublished) and 0.64 for fecundity (f t ) between 1989 and 2020 ( Figure 1b). For the second PVA model, we incorporated environmental stochasticity (i.e., annual variation) by assigning survival and fecundity rates from independent beta and gamma distributions, respectively (Table 1; Model 2). We assigned values from independent distributions because rates of kittiwake annual survival and fecundity at Skomer Island are not temporally correlated (Horswill et al., unpublished, spearman correlation coefficient = .08, p = .67). We selected a gamma distribution for fecundity because even though the maximum documented annual fecundity at Skomer Island is 0.98 chicks per year, kittiwakes can produce more than one offspring annually (max. clutch size for species = 3 eggs). We set the expectation of the beta and gamma distributions as the mean value of the respective vital rates and assigned the variance to reflect local annual variation: ±0.06 SD for adult survival (Horswill et al., unpublished) and ±0.21 SD for fecundity ( Figure 1b). In the third model validation PVA, we used the stochastic values of fecundity selected for Model 2 and applied a temporal trend corresponding to the relationship observed at Skomer Island for 50% kittiwake fecundity (i.e., to model female production only assuming a 1:1 sex ratio; intercept = 0.35, SE = 0.03; slope = À2.20Â10 À3 , SE = 1.85Â10 À3 , Figure 1b) (Table 1; Model 3). The two projection simulations (Models 4 and 5, Table 1) ran from 2020 to 2050 and predicted the population response to an extrinsic threat decreasing rates of survival, such as additional mortality from collisions with wind turbine blades. Models 4 and 5 followed the same structure as Model 2 (i.e., environmentally stochastic with stable fecundity) and Model 3 (i.e., environmentally stochastic with trended fecundity), respectively, with an additional 1% annual mortality applied to each age class (Table 1; Models 4 and 5). Extrinsic factors can differentially influence groups of individuals within populations of seabirds (e.g., Wood et al., 2021), however age-specific impacts to seabird populations from wind energy developments are unknown and therefore we assumed a constant level across age classes. The trended values of fecundity assumed a continued constant rate of decline to that applied in Model 3 (Figure 1b). We derived the figure of additional mortality in accordance with guidelines on acceptable thresholds set by the EU ORNIS committee (EU Birds Directive 79/409/EEC). This threshold is said to meet the condition of a negligible effect on population dynamics, although this is not a consensus view (Schippers et al., 2020). It is also often used in population assessments to determine whether additional mortality might have a significant impact on a population (Backes & Akerboom, 2018;Leopold et al., 2014).
For each of the five PVA models (Table 1), we used a Leslie matrix approach following a pre-breeding census to predict colony size (N) in year t þ 1 as a function of colony size in year t (Equation 1). We assumed a 1:1 sex ratio in the colony and halved fecundity to model female numbers only. In all models, female birds survived from egg laying to chick fledging with probability f t 2 and from fledging to age 1 year with probability ϕ j,t , such that the juvenile age class (0-1 years) has the following combined survival probability: ϕ j,t f t 2 . Birds then survived the immature age class (1-2 years) with probability ϕ i,t , and the pre-breeding class with a probability equivalent to that of breeding adults (ϕ a,t ). Previous studies have shown that kittiwakes first breed between 3 and 5 years old (Cam et al., 2002;Porter & Coulson, 1987), and we included this structure by assuming that 25% (b= 0.25) of 3 year olds, 75% of 4 year olds (1-b= 0.75) and all birds 5 years or older are breeders (Equation 1, Figure 2) (e.g. Frederiksen et al. 2004). This structure has been previously identified as the most suitable for describing the observed population dynamics of kittiwakes breeding on Skomer Island (Horswill et al., unpublished).
Like many seabirds, kittiwakes are largely unobservable during the first years of life. Therefore, colony-specific survival estimates for juveniles (i.e., during the first year post fledging ϕ j,t ) and immatures (i.e., from 1 to 2 years ϕ i,t ) are limited (Horswill & Robinson, 2015). Available estimates for black-legged kittiwakes breeding in France report fledgling survival rates to be 75% of adults, immature survival rates (i.e., from age 1 to 2 years) to be 87.5% of adults, and survival rates to be comparable to adults from age 2 onwards (Cam et al., 2005). We incorporated this structure by assuming an additive relationship between age classes (Cam et al., 2005;Horswill et al., 2014Horswill et al., , 2016, (i.e., synchronous annual fluctuations between age classes) and estimated juvenile and immature survival rates by applying the appropriate scalars to adult survival. We ran the five PVA models for 1000 iterations. Within each iteration, we held temporal residuals associated with environmental stochasticity constant across the three validation models and the two projection models by using the same stochastic draws for survival and fecundity. For the three model validation simulations, we initiated each iteration with the number of kittiwake breeding pairs (i.e., females) observed on Skomer Island in 1989 (N = 2302). For the two model projection simulations, we initiated each iteration with the number of kittiwake breeding pairs observed in 2020 (N = 1681). We estimated the initial stable age distribution for all simulations using the dominant eigenvector of the Leslie matrix based on mean vital rates (Caswell, 2001). Following other assessments of seabird colonies that use a Leslie matrix modeling approach (Miller et al., 2019), we applied a maximum colony size to prevent exponential population growth (N max = 2750). This threshold sets a carrying capacity at 448 breeding pairs above the colony size observed in 1989. All models included demographic stochasticity using binomial and Poisson distributions for survival and fecundity events, respectively. All PVA models were constructed and run using Program R (v. 4.0.2) (R Core Team, 2020).

| Prevalence of temporal trends in rates of fecundity
We identified 70 colonies of kittiwakes breeding in the UK and Ireland with more than 8 years of fecundity data between 1986 and 2019 (mean = 21 years, SD = 9 years).

| How do temporal trends in rates of fecundity influence projected population dynamics?
We ran three model validation simulations (Table 1) to examine how incorporating observed temporal trends in kittiwake rates of fecundity can influence the predictive accuracy of seabird impact assessments using PVA. The outputs of these models did not differ considerably (Figure 3). Based on the mean population trajectories, rates of decline were slightly overestimated between 1989 and 2005 by the deterministic model (Model 1) and F I G U R E 2 Lifetime cycle graph representing the Leslie matrix model describing kittiwake population dynamics at Skomer island (Equation 1). Birds reach age 1 year (N 1,t ) with probability f t 2 ϕ j,t , and age 2 years (N 2,t ) with probability ϕ i,t . At age three, 75% (1 À b) of individuals move the third class (N 3,t ) and 25% (b) of individuals enter the breeding class (N 4,t ). From the third class, 25% of individuals remain as pre-breeders and 75% become breeders. From age two, individuals have a survival probability equivalent to that of adults (ϕ a,t ) environmentally stochastic model without trended fecundity (Model 2). By contrast, the model including the observed trend in rates of fecundity (Model 3) provided a better fit to the data during this period. However, the rapid population decline observed for kittiwakes at Skomer Island between 2005 and 2020 was not accurately captured by any of the validation models, although incorporating a trend on rates of fecundity (Model 3) did result in a slightly faster rate of decline during this period (Figure 3). We ran two model projection simulations to assess how incorporating the observed temporal trend in rates of fecundity may influence the projected population response of kittiwakes at Skomer Island to an extrinsic threat decreasing rates of survival. Both projection simulations demonstrated an increased rate of population decline, compared to the model validation simulations (Figure 3). However, including the observed negative trend in rates of fecundity further accelerated population decline, such that the mean population size was approximately 46% smaller by 2050, compared to the model excluding the observed temporal trend in fecundity (Figure 3).

| DISCUSSION
Seabirds represent one of the most threatened groups of birds in the world with almost half (47%) of all species demonstrating declining population trends (Birdlife International, 2018). In the UK, PVA are now routinely used to assess the potential impacts of offshore renewable energy developments on seabird populations. However, despite available long-term data showing widespread temporal trends in seabird fecundity, these relationships are not readily incorporated into such assessments (e.g., Searle et al., 2019). In this study, we show that omitting observed trends in rates of fecundity from seabird PVAs can dramatically alter the projected population response to a perceived threat.
We show that declining temporal trends in rates of fecundity are common and widespread for kittiwakes breeding in the UK and Ireland. A recent meta-analysis demonstrates that this trend is reflected in fish-eating, surface-foraging species of seabird, like kittiwakes, across the northern hemisphere (Sydeman et al., 2021). We show that kittiwake fecundity at Skomer Island is declining at a slower rate than observed in many other kittiwake colonies across the UK and Ireland. For colonies experiencing faster rates of decline in fecundity, PVAs omitting temporal trends are likely to generate even greater biases than those demonstrated in our study. However, we also found a small number of kittiwake colonies with increasing fecundity. Omitting negative or positive temporal trends in rates of fecundity from PVA assessments risks the projected population response being under-or overestimated with subsequent consequences for the competing priorities of conservation and progressive green energy industries.
The model validation simulations demonstrated that the observed population trajectory fell within the 95% confidence intervals of both stochastic models (i.e., Models 2 and 3). The model including the observed temporal trend on rates of fecundity exhibited a slightly better fit to the observed population dynamics. However, all three validation models failed to reliably recreate the large decline in population size that occurred on Skomer Island between 2005 and 2020. Consequently, our F I G U R E 3 The three validation models (Models 1 -3) generated comparable trajectories, although the model incorporating temporal trends in rates of fecundity (Model 3 yellow) provided a slightly better fit to the data. The projection analysis showed that the population response to an extrinsic threat was considerably greater when the observed temporal trend in rates of fecundity was included (Models 4 and 5). Colored dashed lines show the 95% confidence intervals of the PVA simulations. Vertical dotted line separates the model validation and model projection periods validation models appear to be missing other demographic processes that are important for describing the population dynamics of this seabird colony.
Previous studies have shown that years with elevated breeding failure across the colony can influence dispersal dynamics in kittiwakes (Boulinier et al., 2008). Given the long-term decline in rates of fecundity observed for kittiwakes at Skomer Island, it is possible that dispersal dynamics changes in this population over the study period. This result is also supported by an integrated population modeling study of the population (Horswill et al., unpublished). Dispersal is notoriously difficult to measure empirically for seabirds due to the large sample sizes needed at multiple sites. However, the potential influence of dispersal to the population dynamics of kittiwakes breeding on Skomer Island highlights the importance of validating PVA models (i.e., simulating population dynamics for a period with available count data).
Although not a novel concept (see Beissinger & Westphal, 1998;Brook et al., 1997Brook et al., , 1999Coulson et al., 2001), tests of model validation remain largely lacking from seabird impact assessments. Consequently, we advocate that validation become an essential prerequisite step to seabird PVA studies. A mismatch between the simulated and observed trajectories will highlight that either the distribution of vital rates (mean or variance) has changed over the study period (i.e., due to environmental change) or that key demographic processes (e.g., immigration or emigration) are not considered. Being able to identify that a PVA is likely to under-or overestimate the observed dynamics also supports a more rigorous evaluation of potential biases associated with projection scenarios.
Although the validation models did not differ considerably in their estimated mean population dynamics, this result was not replicated in the projection analysis. The projection model incorporating temporal trends in rates of fecundity generated faster rates of population decline, compared to the model excluding this relationship, such that the forecasted mean population size was approximately 50% smaller by 2050. Including temporal trends in vital rates allows any long-term change in the mean value to be incorporated, such that divergent outcomes can arise given enough time, even if differences estimated during the validation period are small. This finding highlights the importance of considering a range of projection models to understand all future dynamics that may arise from conservation and management decisions, assuming the models selected for projection pass a baseline validation threshold.
To represent the hypothetical extrinsic threat in the projection simulations, we applied 1% additional annual mortality to each age class. This adjustment is considered to meet the condition of a negligible effect on population dynamics (EU Birds Directive 79/409/EEC). However, we show that this level of change can still result in severe population decline that could drive small populations to extinction. This finding supports previous work highlighting that a 1% mortality threshold is inadequate for providing safe thresholds of anthropogenic change in long-lived species, such as seabirds (Schippers et al., 2020).
Previous studies show that fecundity in kittiwakes, like other long-lived species, has a low demographic impact on population dynamics, that is, low elasticity . It would be understandable to assume that accurately quantifying vital rates with little demographic impact is of low priority when parameterising PVA models. However, environmental stochasticity can increase the demographic impact of vital rates with low elasticity (Gaillard et al., 1998), and we show that incorporating a directional trend in seabird rates of fecundity can dramatically alter the population dynamics projected by PVA. Therefore, we argue that PVA assessments of seabirds need to account for environmental factors influencing current population dynamics to reflect potential "worst-case" scenarios, acknowledging the assumption that the same drivers continue to be influential into the future.
Not all populations will have long-term monitoring data available. When using PVA to conduct impact assessments on data-limited populations, we recommend an appropriate monitoring study be undertaken to quantify the local population trajectory and stability of vital rates. Populations experiencing large annual or cyclical fluctuations in demography are likely to require longer time series of data to reliably estimate directional change through time. If data collection is not possible, an alternative approach is to reconstruct mean vital rates using predictive methods based on life-history trade-offs (Horswill et al., 2019 and consider a range of temporal trends reflective of local seabird populations. However, interpreting the results of this latter approach may be complicated by large confidence intervals that accurately reflect uncertainty in the absence of empirical evidence. We focus on the implications of overlooking temporal trends in rates of fecundity when using PVA to assess the population response of seabirds to offshore renewable energy developments. Reviews of existing PVA literature identify a range of additional factors to consider when applying this method to impact assessments and conservation decision making (e.g., Beissinger & Westphal, 1998;Reed et al., 2002;Norris, 2004;Pe'er et al., 2013;Chaudhary & Oli, 2020;Searle et al., 2020). Including density-dependent regulation in assessment focused PVA studies of seabirds remains debated Horswill et al., 2017;O'Brien et al., 2017). A lack of clear guidance on the appropriate strengths and direction of density-dependent regulation means that a densityindependent approach is generally considered precautionary, although this will not be the case for populations experiencing depensatory control (Horswill et al., 2017). As a result, there is a pressing need for more research to determine the prevalence of density-dependence feedbacks across seabird species, the relative importance of these relationships to population dynamics compared to long-term trends in the environment, and the sensitivity of results to the inclusion of different density-dependence formulations.
In this paper, we demonstrate that declining temporal trends in rates of fecundity are widespread for populations of kittiwakes breeding in the UK and Ireland. We also show that omitting this relationship from PVA models used to assess the potential impacts from offshore wind energy developments to seabirds can dramatically underestimate the projected risk. We advocate that model validation (i.e., simulating population dynamics for a period with available count data) become an essential prerequisite step for conducting PVA assessments on seabirds. If simulated dynamics do not match the observed trajectory, then likely biases (i.e., from over or underestimating population dynamics) need to be considered when interpreting the outputs from projection analyses. More broadly, this study highlights the importance of considering existing drivers of demographic change when assessing protected species for conservation and infrastructure decision making.