Later at higher latitudes: large-scale variability in seabird breeding timing and synchronicity

. In seasonal environments, organisms are expected to optimally schedule reproduction within an annual range of environmental conditions. Latitudinal gradients generate a range of seasonality to which we can expect adaptations to have evolved, and can be used to explore drivers of timing strategies across species’ distribution ranges. This study compares the timing of egg hatching in four seabird species (Atlantic puffin Fratercula arctica, black-legged kittiwake Rissa tridactyla, common guillemot Uria aalge , and Brünnich’s guillemot U. lomvia ) covering a subarctic to Arctic latitudinal gradient along the Norwegian coast to Svalbard (65–79°N). Hatching was significantly delayed by an estimated 1.7, 2.3, and 1.9 d per latitudinal degree for puffins, kittiwakes, and common guillemots, respectively, but was not delayed for Brünnich’s guillemots. Hatching distributions revealed an increase in intra-annual breeding synchronicity along a latitudinal gradient for kittiwakes only, whereas the two guillemots exhibited high hatching syn - chronicity at all colonies. We used this large-scale, multispecies timing data series to discuss constraints, adaptations, and mechanisms affecting breeding timing, a necessary step to recognize risks to populations and predict future ecosystem change.


IntroductIon
Theory predicts that species living in seasonal systems would evolve to optimally time important life history events, such as reproduction, within a spectrum of yearly conditions (McNamara and Houston 2008).Evidence for evolutionary timing adaptations to seasonality span a wide set of taxa and ecosystems (Perrins 1965, Varpe et al. 2007, Galen and Stanton 2010).Timing of breeding is a key adaptation with large impacts on fitness.Major potential selective forces of breeding v www.esajournals.orgBURR ET AL.
timing include (1) direct physical constraints on the temporal window within which breeding can occur; for example, ice or snow cover on nesting sites limiting breeding onset, and (2) species adaptations that have arisen through other natural selection processes to increase survival or fecundity, for example, by matching breeding timing with peak prey abundances (Lack 1950), timing migration to avoid predation (Lank et al. 2003), or synchronizing fledging as a predator-swamping technique (Ims 1990).The degree of seasonality in several components of the environment, such as fluctuations in temperature and light linked to predictable within-year cycling, typically increases with latitude (Kirk 1994), often accompanied by systematic differences in the timing of seasonal events from low to high latitudes (Broms andMelle 2007, Vikebø et al. 2012).Latitude has therefore often been used as a proxy in avian phenology studies (Baker 1938, Slagsvold 1975, Olsen and Marples 1993, Wanless et al. 2008), and we choose to use the same parameter as a proxy for a general gradient of environmental conditions to address large-scale spatial comparisons and species differences.Studying variability in breeding timing along a gradient from temperate to polar ecosystems, and thus a gradient of seasonality, offers a unique opportunity to address the links between reproductive timing and seasonality.
In this study, we use data from Arctic and subarctic seabird colonies in the northeastern Atlantic to explore the variation in hatching timing along an extensive latitudinal gradient (65-79°N).First, we quantify the reproductive timing of four common species (Atlantic puffin Fratercula arctica, black-legged kittiwake Rissa tridactyla, common guillemot Uria aalge, and Brünnich's guillemot U. lomvia) across a range of Norwegian colonies.Second, we use this multispecies, multiyear, and large-scale data set to test two predictions on hatching time variability: (1) Hatching time will delay with increasing latitude.Birds in seasonal environments are expected to meet peak energy requirements by temporally matching their reproduction with the phenology and abundance of their prey in proximity of the nest site (Cushing 1975, 1990, Orians and Pearson 1979).Based on delays in timing with increasing latitude at lower trophic levels (Lie 1965, Vikebø et al. 2012), and because bottom-up control is typical of seabirds in northern Atlantic systems (Aebischer et al. 1990, Durant et al. 2003, Frederiksen et al. 2006), we hypothesize that seabird hatching timing would delay with latitude to match prey availability in local systems.Furthermore, physical conditions at higher latitudes such as snow and ice may also delay the onset of breeding (e.g., Moe et al. 2009).
(2) Intra-annual breeding synchronicity between individuals in a colony will increase with latitude.Due to a steeper seasonal gradient (e.g., shorter summers) throughout the year, it can be expected that breeding at higher latitudes must fit into a narrower temporal window.Such a trend is seen at the lowest trophic level, where phytoplankton bloom duration is shorter at higher latitudes (Zenkevitch 1963).
Finally, we also quantify inter-annual variability in hatching timing along a latitudinal gradient and compare it between species.We discuss our findings in relation to seabird life histories and behaviors as our study species differ both along the life history continuum (Gaillard et al. 1989) and with respect to feeding behaviors (Shealer 2001).

Study system and species
The 10 study colonies are situated on the coasts of mainland Norway and Svalbard in the northeastern Atlantic.Four seabird species were chosen for comparison because (1) they reflect important differences in life history patterns and feeding behaviors (Appendix S1: Table S1), (2) their breeding ranges span a substantial latitudinal range, and (3) sufficient data were available.

Data acquisition
Annual measures of breeding timing were collected for one or more of the four study species at each of the 10 study colonies (Table 1).Hatching time was used as the phenology measure, and values for a sample of individual nests were used when available.When only annual mean hatching times and standard v www.esajournals.orgBURR ET AL. deviations (SD) were available, these values were included when analyses did not require data at the individual level.Individual hatching dates were estimated (1) by nest inspections at variable intervals during hatching or egg laying and, in the latter case, adding incubation time (hatching dates with a precision level poorer than ±5 d were excluded from analyses), (2) from egg density calculated from mass and volume measurements (Hoyt 1970) and compared with local, species-specific curves for mass loss during incubation, or (3) by assessing chick age from growth functions for head-bill lengths (Anker-Nilssen and Aarvak 2006) or wing length of nestlings (Barrett 2001; see Appendix S2: Table S2 for colony-specific details).
This study used individual nest data from the years 1980-2013 and colony mean (±SD) data from years between and including 1974-2004; however, data available for each colony and spe-cies were not necessarily available for all these years (Table 1).For some years, phenology data were not available because of complete breeding failure in the colony or insufficient temporal coverage of the field season (Appendix S2: Table S2).Because of poor hatching success, we did not include the very few hatching dates for puffins at Sklinna in 2010 (n = 3) and Røst in 2011(n = 4) in the analyses.As mean hatching in both cases did not differ from the long-term average by more than 1 d, we do not expect the exclusion of these data to bias our results.
The numbers of breeding pairs and annual population trends for the past 10 yr at each colony were estimated to the best of our ability from existing sources for the year 2013; the final year included in our study (or the most applicable year for that data; Appendix S2: Table S2).Although it has been suggested that larger kittiwake colonies have a wider spread of breeding timing (Coulson Notes: Equal weights given to all years for calculation of mean of yearly mean hatching times and standard deviation of the mean.SD not applicable for colonies with data for fewer than 3 yr.v www.esajournals.orgBURR ET AL. and White 1956), we found no significant effect of colony size on breeding synchronicity (linear model with species effect, F 1,122 = 0.90, P = 0.35; Appendix S3: Fig. S1).Additionally, we found no overall latitudinal trend for colony size (linear model with species effect, F 1,122 = 2.51, P = 0.12; Appendix S3: Fig. S2).Furthermore, in this context it is difficult to choose the appropriate scale to enumerate the size of a breeding colony in cases where the colony limits are not well defined (i.e., a fjord or archipelago where many locations are breeding sites), and when population sizes vary over time.Based on this reasoning and the above statistical assessment, we did not include colony size in the following statistical models.

Statistical analysis
Latitudinal effect on hatching timing.-Totest for the effect of latitude on hatching timing, to compare hatching times between species, and to test for an interaction between latitude and species, we applied a linear mixed-effects model (LME) with hatch date at the individual level as the response variable (n = 4353), latitude and species fixed effects, and a random year effect.Though there is moderate heteroscedasticity due to larger hatch date dispersion at lower than higher latitudes, other models that are more robust to unequal variance, such as a generalized linear mixed-effects model with quasi-Poisson distribution, led to the same results.
Latitudinal effect on intra-annual hatching synchronicity.-We also applied an LME with a random year effect to measure the effect of latitude on intra-annual hatching synchronicity, to compare intra-annual variability between species and to test for an interaction between latitude and species.Hatching synchronicity within a colony was estimated for each year as the sample standard deviation (SD) from the hatching mean; a lower value describes less spread and hence greater synchronicity.For colonies without individual hatching data, any SD values that were available were used directly in this analysis (details in Table 1).
Inter-annual variability in hatching timing.-Tovisualize yearly hatching distributions, we plotted probability density functions of yearly hatching timing that were calculated using kernel density estimations.We used a standardized smoothing bandwidth for these distributions.For cases, however, when more than 50% of hatching values were reported on the same date, we increased the smoothing bandwidth.We did this because these cases correspond to field seasons with infrequent nest visits, and where we could assume that hatching occurred during the time window spanned by these visit, rather than on the day of the visit only.The two cases where this happened were (1) puffins on Anda in 2010 and (2) kittiwakes in Isfjorden 2010, where in both cases over 50% of hatching dates occurred between two consecutive nest visits.To quantify inter-annual variability, we calculated the SD from the mean of all annual mean hatching times for colonies with more than 1 yr of data.Each colony with more than 1 yr of individual hatching data thus had one value to represent the inter-annual variability.Lower SD values in this analysis indicate less inter-annual variability in timing.We then used an analysis of covariance to evaluate the effect of latitude, species, and their interaction on inter-annual variability.Only data from the years 2004-2013 were included to standardize the time frame used in this analysis, but the years of data available varied between species and colonies within this period (Table 1).R (version 3.0.2,R Core Team 2013) was used for all statistical analyses.

Timing of seabird breeding
Individual hatching dates varied between May 8 and August 3.There was, however, substantial variability between species and latitudes within this period (Table 1,Figs. 1,2,3).For puffins, with the exception of some years at lower latitudes, timing was rather constant across years and colonies.In general, Brünnich's and common guillemot hatching timing was more synchronous at all colonies compared to kittiwakes and puffins, and hatching distributions and their timing were similar in absolute timing from Bjørnøya northwards (Table 1).Conversely, kittiwake hatching timing was much more variable between colonies, both in regard to absolute timing and to the ranges of hatching dates.

Latitudinal effect on intra-annual hatching synchronicity
Hatching synchronicity increased significantly with latitude for kittiwakes only (P < 0.001; Fig. 5a).The average yearly SD for kittiwake hatching dates ranged from 7.3 to 16.6 d at the most southern colony (Røst), but only 2.4-6.7 d at the most northern site (Kongsfjorden).There were no significant changes in intra-annual hatching synchronicity with latitude for the other species (puffins: P = 0.21, common guillemots: P = 0.31, Brünnich's guillemots: P = 0.36; Fig. 5b).

Inter-annual variability in hatching timing
Probability density functions revealed notable variability in timing between years for kittiwakes at Hornøya and Bjørnøya (Fig. 2), and this pattern is more pronounced for kittiwakes than for the other species.For example at Hornøya, one of two colonies that has data for three of the four study species, density functions of hatching times for puffins (SD = 1.4 d) and common guillemots (SD = 1.4 d) were more consistent over years than for kittiwakes (SD = 2.3 d).There was no statistically significant effect of latitude on inter-annual variability in hatching timing (F 1 = 0.086, P = 0.78), or a significant difference between species (F 3 = 1.09,P = 0.41).

dIscussIon
The main findings of our study are that (1) hatching timing delayed with latitude in three of the four species considered, (2) intra-annual breeding synchronicity increased with latitude for kittiwakes but not for guillemots and puffins, and (3) the study species exhibited different levels of inter-annual variability in hatching timing with kittiwakes showing the largest variability.
The prediction that seabird hatching delays with latitude was confirmed for three of the four species in agreement with other studies on birds (Baker 1938, Slagsvold 1975, Olsen and Marples 1993, Wanless et al. 2008).Due to the large spatial and temporal scale of the data, we do not consider the noise introduced through an amalgamation of data sets that combine several methods with nonoverlapping years to be problematic.In seasonally pulsed systems, such as the northeastern Atlantic, where productivity at lower trophic levels has been shown to delay with latitude, we can assume timing has been shaped to some degree by prey phenology, although it is unknown to what extent other factors such as physical conditions (i.e., ice cover at sea and snow cover on land) constrain breeding timing in the northeastern Atlantic.Arrival dates of Brünnich's guillemots on the west coast of Greenland have been correlated with latitude, and the timing of  sea ice breakup was thought to be an important determinant of breeding initiation (Laidre et al. 2008).Overall, the cues that trigger seabirds to breed remain unclear (Frederiksen et al. 2004), yet they are important for understanding how species will respond to changes in climate and to what extent trophic links may be uncoupled due to changing phenology.
Kittiwakes exhibited an increase in synchronicity with latitude and we suggest that the latitudinal gradient of seasonality primarily drives this change.At 79°N, a sharper gradient of physical conditions within the annual cycle would select for breeding during a narrower temporal window (either directly through abiotic conditions or indirectly through prey availability), and thus more synchronously.Individuals at lower latitudes may not experience this same time constraint.In addition, we expect predation to simultaneously shape hatching timing.Higher rates of nest predation at lower latitudes (McKinnon et al. 2010) may prompt a more frequent replacement of lost clutches, and at lower latitudes where the potential window for breeding is longer, the effects of more frequent clutch replacement would manifest as wider hatching distributions.This was likely the case at Røst (T.Anker-Nilssen, personal observation), the most southern kittiwake colony in this study, where we found the widest of all the studied distributions.Some hatching dates from Røst were also estimated based on laying dates, and these values could therefore contribute to wider distributions (because a few eggs that did not hatch were also included in the analysis), yet even when Røst was removed from analysis, the trend for increasing synchronicity with latitude held for kittiwakes.In general, we used hatching date as our standardized phenology measure because it provided the most raw data; however, using laying dates for a large spatial scale would allow for the inclusion of years when eggs were laid but did not hatch, and thus the assessment of how years of low hatching success influence spatial and species-specific trends in timing.In contrast to kittiwakes, the generally high synchronicity of guillemots across their range suggests that another mechanism may be the primary driver of breeding synchronicity for these species.Flightless guillemot chicks fledge by jumping from cliffs typically within a period of only a few days as a predator-swamping strategy (Hatchwell 1991), and this breeding tactic would encourage synchronous breeding regardless of latitude.The primary driver of breeding synchronicity may thus vary depending on the species considered, and the drivers of breeding synchronicity are not necessarily the same as those of absolute timing.
In quantifying inter-annual variability in hatching timing, we have revealed that kittiwakes exhibit a higher degree of inter-annual variability than the other study species at several colonies.For example, the distributions of hatching times for kittiwakes at Hornøya are notable for their high inter-annual variability.Colony-specific findings such as this highlight the need for further investigation on local environmental drivers, but the contrast of species offers an opportunity to discuss potential mechanisms pertaining to species characteristics.For example, the effects of differences between species in regard to life history and feeding strategies may be acting in conjunction with other drivers of breeding timing.On a fast-slow continuum of life histories (Gaillard et al. 1989), kittiwakes can be described as a 'faster' species with comparatively shorter life spans and higher yearly reproductive outputs than puffins and guillemots.It is possible that life history differences along this continuum may influence reproductive timing strategies, and help explain the observed difference between species in regard to inter-annual variability.We are, however, unable to conclude on the importance of life history differences as a driving mechanism without data to evaluate the contribution of different environmental conditions.Inter-annual variability in kittiwake hatching timing may also be a response to the timing and availability of prey in surface waters.Kittiwakes are dependent on physical and biotic factors bringing their prey to the surface, whereas pursuit diving auks are able to exploit an additional dimension of the water column and may therefore be less restricted by the temporal availability of prey.With wings better adapted for flight as opposed to diving, kittiwakes can, however, forage over a greater area with lower energetic costs, potentially making them more flexible than auks in terms of feeding opportunities, and thus leading to increased variability in timing.Which component is the most important in explaining breeding timing may vary, and more geospatial data for foraging v www.esajournals.orgBURR ET AL.
birds in relation to timing are needed to disentangle this issue further.
Our large-scale comparative work sheds light on the conservation issues associated with phenology; as physical and prey-related conditions associated with breeding change, different species will face different challenges at different latitudes.A prediction and explanation of the resulting differences in phenological responses to climate change and, ultimately, their success or failure, will require a better mechanistic understanding of the main drivers and cues involved in breeding timing.This is a common challenge across taxa and environments (Winder and Schindler 2004, Visser and Both 2005, Forrest and Miller-Rushing 2010).Our study suggests that several mechanisms (e.g., to match prey availability, avoid predation, or deal with physical environmental restrictions) likely act together but differently for the various study species to shape breeding timing in highly seasonal environments.
Our large-scale spatial comparison is only one step to better understand variation in seabird phenology across space and between species.Addressing temporal trends was beyond the scope of this project, but we suggest that future studies apply time series analyses to explore how these patterns are affected by key drivers of environmental change.If the mechanisms differ from those timing their prey, seabirds may risk trophic decoupling.Although, seabird phenology has been shown to be flexible on an individual level in variable environments (Reed et al. 2009), seabirds may respond adversely if factors other than prey availability are influencing breeding timing.If seabirds are responding to cues outside of the marine environment to initiate breeding (e.g., little auks Alle alle responding to snow cover; Moe et al. 2009) or are constrained by other timing mechanisms (i.e., narrow temporal window at high latitudes), the chances of mismatch with available prey may increase.In the North Sea, changes in seabird phenology were shown to be insufficient to keep up with changes in timing of the optimal size of an important prey species (Burthe et al. 2012), and unequal rates of change in phenology between trophic levels have been documented and linked to population declines in several other taxa (Edwards and Richardson 2004, Winder and Schindler 2004, Both et al. 2009).From this perspective, our study stresses the importance of understanding to what extent different mechanisms shape breeding timing in seabirds, and thereby how phenology is a critical component for predicting future change in seabird populations.

AcknowledgMents
This research was made possible by the SEAPOP (SEAbird POPulations) program of Norway, and the many involved with the data collection.We thank the Fram Centre for their support and incentive funds to develop the later stages of the work, and an anonymous reviewer for valuable comments.Knut Helge Jensen provided helpful statistical advice, with further support from University of Bergen's Theoretical Ecology Group.Author contributions: ØV, TAN, and KEE initiated the study.ZB and ØV designed the study and analyzed the data.ZB wrote the manuscript with main contributions from ØV, TAN, KEE, SD, RTB, SCD, BM, and TKR.TAN, KEE, SD, RTB, CB, SCD, SHL, BM, TKR, and HS contributed data.lIterAture cIted

Fig. 1 .
Fig. 1.Yearly kernel density functions (nonparametric probability estimations) of puffin hatching dates at four colonies with mean of annual mean hatching times (±1 SD) indicated by point and bar above curves and corresponding values noted below the colony name.Yearly sample sizes indicated next to year legend.Notes: only a subset of the most recent years from Hornøya are portrayed in this figure for visual simplicity, and y-axes are not consistent for all graphs.The distribution of puffins at Anda in 2010 has an increased smoothing bandwidth due to sampling infrequency.

Fig. 2 .
Fig. 2. Yearly kernel density functions (nonparametric probability estimations) of kittiwakes hatching dates at six colonies with mean of annual mean hatching times (±1 SD) indicated by point and bar above curves and corresponding values noted below the colony name.Yearly sample sizes indicated next to year legend.Notes: y-axes are not consistent for all graphs.The distribution of kittiwakes at Isfjorden in 2010 has an increased smoothing bandwidth due to sampling infrequency.

Fig. 3 .
Fig. 3. Yearly kernel density functions (nonparametric probability estimations) for hatching dates of (a) Brünnich's guillemots at three colonies and (b) common guillemots at two colonies with mean of annual mean hatching times (±1 SD) indicated by point and bar above curves and corresponding values noted below the colony name.

Fig. 4 .
Fig. 4. Mean of annual mean hatching dates (± 1 SD) at each colony for four study species along a latitudinal gradient with the species-specific generalized linear model results superimposed.Only significant species-specific linear model lines are shown.Hollow triangles represent kittiwake annual mean hatching dates for which there were no individual hatching data available; these values were therefore not included in the statistical models.

Fig. 5 .
Fig. 5. Intra-annual breeding synchronicity in relation to latitude for Norwegian colonies spanning 65.2-78.9°Nfor (a) kittiwakes with significant LME model line and (b) puffins, common guillemots and Brünnich's guillemots, where points of common guillemots were offset by +0.05° to avoid overlapping.Each point represents one colony year.

Table 1 .
Colony-specific data information and hatching timing results.