Citizen-science data provides new insight into annual and seasonal variation in migration patterns

Current rates of global environmental and climate change pose potential challenges for migratory species that must cope with or adapt to new conditions and different rates of change across broad spatial scales throughout their annual life cycle. North American migratory hummingbirds may be especially sensitive to changes in environment and climate due to their extremely small body size, high metabolic rates, and dependence on nectar as a main resource. We used occurrence information from the eBird citizen-science database to track migratory movements of five North American hummingbird species (Archilochus alexandri, A. colubris, Selasphorus calliope, S. platycercus, and S. rufus) across 6 years (2008–2013) at a daily temporal resolution to describe annual and seasonal variation in migration patterns. Our findings suggest that the timing of the onset of spring migration generally varies less than the arrival on the wintering grounds. Species follow similar routes across years, but exhibit more variation in daily longitude than latitude. Long distance migrants generally had less annual variation in geographic location and timing than shorter distance migrants. Our study is among the first to examine variation in migration routes and timing for hummingbirds, but more work is needed to understand the capacity of these species to respond to different rates of environmental change along their migratory routes.


INTRODUCTION
Global environmental changes are impacting species' populations and movement behavior (e.g., La Thompson 2007, Wilcove andWikelski 2008). Monitoring migration response to these changes is critical to predict how species will respond under future climate scenarios. Because migratory species carry out parts of their annual life cycle in different locations, they must cope with shifting conditions across regional to continental scales (Møller et al. 2008, Both et al. 2010, Klaassen et al. 2012. Across these scales, environmental change will vary in rate and magnitude (e.g., drought, extreme temperature) requiring different biological responses throughout the migratory route. Further, individuals likely have different constraints acting on them in different seasons, which may influence their migration speeds and routes (La Sorte et al. 2013). Our understanding of avian migration patterns throughout their annual life cycle (breeding, non-breeding, migration) is poor for many species (Faaborg et al. 2010). Yet fewer studies have attempted to evaluate the form and magnitude of annual variation in populationlevel migration patterns across broad scales as opposed to evaluating migration for a subset of marked individuals or at a single site (but see Marra et al. 2005, Gordo and Sanz 2006, Van Buskirk et al. 2009, Vardanis et al. 2011).
An important challenge facing ecologists is to document the consistency in seasonal migration patterns across years. Seasonal differences in migration patterns may depend on physiological or behavioral constraints (e.g., energetic demand, reproductive condition, resource acquisition; Sillett andHolmes 2002, Conklin et al. 2013). However, annual variation in migration may also depend on year-to-year differences in precipitation, temperature, or extreme weather events that in turn impact the physiological condition of migrants and may cause individuals to alter their course. For example, in birds optimal spring arrival at breeding grounds may be under stronger selection than arrival at wintering grounds (Yohannes et al. 2009, Alerstam 2011, Karlsson et al. 2012, and thus individuals may initiate spring migration with little variation across years, while adjusting their precise location, speed, and stop-over decisions en route based on daily conditions (Erni et al. 2002, Leitchi 2006, Tottrup et al. 2008, Conklin et al. 2013. Consistency in seasonal migration patterns may depend on a species' total migration distance (positively associated with consistency; Tottrup et al. 2008, Alerstam 2011 or the magnitude of environmental heterogeneity (negatively associated with consistency; La Sorte et al. 2014a) and prevailing wind patterns , Stanley et al. 2012, La Sorte et al. 2014b) along its route. Environmental shifts among seasons or years may impact individual movement decisions and thus population-level migration patterns.
Monitoring seasonal and annual variation in population migration patterns is attainable using broad-scale citizen science databases (Dickinson et al. 2010). We evaluate hummingbird migration because their high energetic demands and reliance on floral nectar resources (Lasiewski 1963, Healy and Calder 2006, Cotton 2007 suggest that they should display annual flexibility in their migration routes and timing and may be early responders to large-scale climate shifts. Movement of banded individuals can be tracked at a coarse level (USGS Bird Banding Laboratory 2104), but there are often few recaptures (Appendix A), which limits our ability to describe and predict migration patterns for hummingbirds. Technology, beyond traditional band recovery data, is rapidly advancing to monitor individual-level migratory behaviors for many species (e.g., satellite and GPS tags; for review, Robinson et al. 2009, Bridge et al. 2011), but North American hummingbirds are too small (most species ,5 g) to implement tracking methods (but see Moran et al. 2013 for stable isotope analysis). Hummingbirds appear to migrate similarly to other birds, following looped migration routes where males initiate migration earlier than females (Phillips 1975, Baltosser and Russell 2000, Camfield et al. 2013), but many fundamental questions regarding migration behavior and route remain Calder 1994, Weidensaul et al. 2013). Previous studies have evaluated hummingbird abundance, behavior and resource acquisition at specific sites during migration (e.g., Kodric-Brown and Brown 1978, Russell et al. 1994, McKinney et al. 2012), but not along the entire route. Hummingbird ability to respond to short-term annual or seasonal changes in climate and environment along their migratory pathways is not well understood, which represents a critical problem given ongoing environmental change.
Here, we use citizen-science data from eBird (Sullivan et al. 2014) to describe the routes (geographic location) and timing (onset of spring migration, end of autumn migration, and speed) of five North American migratory hummingbird species and to evaluate the form and magnitude of seasonal, annual and species-level variation in migration patterns. We expect that hummingbirds will migrate faster in the spring than in the autumn, and that species will have looped migration routes. Annually, we expect greater variation in the end of autumn migration and in autumn route (geographic location) than in the onset of spring migration or spring route, and generally greater variation in daily longitude than latitude across years. Among species, we expect birds with longer total migration distance to exhibit more consistent seasonal migration patterns, especially with respect to migration timing.

Species data
We extracted checklists of hummingbird observations from 2008 to 2013 from the eBird database (available online, www.ebird.org) for five North American latitudinal migrants: rubythroated hummingbird (Archilochus colubris), black-chinned hummingbird (Archilochus alexandri ), rufous hummingbird (Selasphorus rufus), calliope hummingbird (S. calliope), and broadtailed hummingbird (S. platycercus). We only included vetted eBird observation data collected under the ''traveling count'', ''stationary count'', and ''incidental observation'' sampling protocols. Quality control for eBird includes structured protocols for data entry, automated data filters, and a large network of regional editors (Sullivan et al. 2014). We use the eBird data as presenceonly records, with no attempt to infer abundance or absence. We initially downloaded and analyzed data for 2004-2013, but we chose to exclude data from earlier years of eBird reporting (2004)(2005)(2006)(2007) due to relatively low observer effort and poorer geographic coverage. During these years, eBird observer effort was significantly lower (Appendix B: Fig. B2) and standard error for the population daily centroids was much greater during the years, for all species (Appendix B: Fig. B4).
In addition to the population-level data, we obtained recapture data of banded individuals  from the Bird Banding Lab (BBL; USGS Bird Banding Laboratory 2014; for data: http://www.pwrc.usgs.gov/BBL/homepage/ datarequest.cfm). While there were too few recaptures within a migration season for statistical analysis, we provide these maps as a qualitative independent comparison to our population-level results (Appendix A).

Occurrence centroids for migration routes and timing
We compiled species daily observations from eBird checklists for each Julian day in each year (2008)(2009)(2010)(2011)(2012)(2013) within equal-area hexagonal cells (12,452 km 2 ) of an icosahedron map of North America (Sahr et al. 2003). The hexagon size is at a coarse enough spatial resolution to allow for detection of daily population-level movements (sensu La Sorte et al. 2013). We used the central longitude and latitude of each cell to calculate a weighted mean daily location, based on the number of checklists where the species was observed for each cell and the total number of checklists submitted for each cell on that day ( Fig  B1). In this way, weights accounted for spatial variation in the total number of checklists submitted per cell across time (Sullivan et al. 2014), which acted as a proxy for total observer effort. For each year, we took the weighted mean daily locations and used a generalized additive model (GAM; package mgcv; Wood 2006) to separately fit the latitude and longitude components of the daily locations as a function of time (Julian day of the year). Using the daily latitude and longitude predictions based on the GAM fits, we combined the predicted locations to obtain predicted population-level daily occurrence centroids (including standard error and 95% confidence intervals) for each species in each year.
We derived species-specific estimates of the onset of spring migration and end of autumn migration for each year within the study area using GAMs (Wood 2006) and segmented linear regression (Muggeo 2003). To separate seasons we used the date at which the population-level centroid was at its maximum latitude. We assumed that at this time, the entire population was fully established on the breeding ground, and southward migration had not yet begun. For each season, we applied a GAM to the predicted daily latitudes of the population-level centroids. We used the minimum upper limit from the 99% confidence band of the predicted daily latitude threshold to identify the winter season occurrences before spring migration and after autumn migration. For migratory species, the beginning of spring migration was defined as the date when the predicted latitude threshold was reached in spring (11 January-9 July) and the end of autumn migration as the date when the threshold was v www.esajournals.org reached in autumn (8 August-21 December; sensu La Sorte et al. 2013). The time within these dates contain the bulk of hummingbird observations (Appendix B: Fig. B3). We used these dates as starting points for a segmented regression on latitude by Julian date to refine the estimated breakpoints for the beginning of spring and the end of autumn migration (segmented package; Muggeo 2008). We measured the daily great circle distance traveled by the population in km/ day between all consecutive centroids from the GAM analysis (function geodist, package spaa; Zhang 2013). For spring and autumn in each year, we took the median of the five fastest speeds (to minimize the effects of potential outliers) as the maximum population-level migration speed (sensu La Sorte et al. 2013).
Total migration distance was estimated using species' breeding and winter ranges from Na-tureServe BirdLife digital range maps (Ridgely et al. 2007). We opted to use these range maps instead of calculating ranges based on eBird data because there is variation across years and across species in data quality and the geographic regions used by hummingbirds, especially in the wintering grounds. By using the NatureServe maps to estimate total migration distance, we have an independent static variable that can be used to compare results across the five species. Breeding and winter range map polygons were converted to 12,452 km 2 equal-area hexagons, and the center was calculated by averaging the latitude and longitude of the hexagon centers located in each seasonal range. We estimated total migration distance for each species using great-circle distance between the central points calculated for each breeding and winter range (La Sorte et al. 2013).

Statistical analysis
Citizen-science observations in eBird may be influenced by observations of vagrant birds, subpopulations that behave differently, or by geographic human population centers that contain large groups of active birders, such as the southeastern United States (Sullivan et al. 2014). To address these potential issues, we conducted the analyses in two ways: (1) using all eBird observations and (2) B2; Sullivan et al. 2014). To determine if we could make robust comparisons across years, we evaluated the error in the predicted daily centroids for each year.
We compared hummingbird migration routes and timing to determine the magnitude of annual and seasonal variance across species. For migration timing, we compared the estimated beginning of spring migration, peak latitude, the end of autumn migration, and population-level migration speed within seasons and across years using mean and standard deviation, and we used linear regression on the time series for each species to look for directional trends in migration timing across years. We used a generalized additive mixed model approach ( package gamm4; Wood 2006) with year as a random effect to determine the amount of variance in latitude and longitude explained by Julian date. Statistical analyses were conducted in R (R Development Core Team 2014). Code and analyses are available in the online Supplement (see Supplemental Material).

RESULTS
For the western species, population-level migration routes and timing were substantially impacted by the effect of vagrant individuals and by the high density of year-round eBird observers in the southeastern USA (Appendix B: Fig. B5). While these observations are themselves interesting, they do not represent the main population-level migratory movements that we characterize here, and require further analysis and study. Thus, we report results for the observations constrained to the western flyway for the four western species below (for results using the full dataset see Appendix B: Fig. B5; Table B1). The routes predicted using data confined to the western flyway provide a better geographic correspondence with the location of known species breeding and wintering ranges, and their likely migration routes (Fig. 1). BBL data support our decision to confine our observations of western species to the western flyway, as western birds were banded and recaptured in the east mainly during the winter season, and v www.esajournals.org The five species had a median body mass of 3.3 g (range 2.3-4.1 g) and a mean migration distance between breeding and wintering grounds of 2,716 km (range 1,721-4,103 km). Three of the five species (A. alexandri, S. calliope, and S. rufus) migrated faster in the spring than in the autumn (Table 1). The three Selasphorus species (S. platycercus, S. calliope, and S. rufus) exhibited strong looped migration patterns, but the two Archilochus species (A. alexandri and A. colubris) followed similar routes in spring and autumn (Fig. 1). All species had greater annual variation in the estimated date for end of autumn migration compared to the beginning of spring migration (Table 1; sd spring : 5.2-11.8; sd aut : 8.0-26.9, whereas less variation was displayed in the estimated date for the onset of spring migration and for the population to reach peak latitude (Table 1, Fig. 1; sd peak : 2.3-22.0). There were no strong directional trends in migration timing (beginning of spring, peak latitude or end of autumn migration) across the six years, except for A. alexandri, which demonstrated an increasingly later start of spring migration (Appendix B: Fig. B6

DISCUSSION
Our results provide among the first broadscale temporal perspective on hummingbird migration across seasons and years, and demonstrate the role of citizen-science as a powerful tool for describing population-level movements through time. We found that migration patterns  (Wethington and Russell 2003). Consistent with predictions originating from optimal migration theory that suggest birds should move more quickly in the spring to optimize arrival for breeding (Yohannes et al. 2009, Alerstam 2011, Karlsson et al. 2012), most of the species evaluated here migrated more quickly in the spring than in the autumn. In addition, longitude was less well explained by Julian date than latitude, which may indicate greater ability for populations to adjust their location east-west to track ideal resources or weather, possibly along elevation gradients (La Sorte et al. 2014a), while still moving north-south at similar rates across years. Given their small body size and high energetic requirements, hummingbirds may be early responders to climate variation and small differences in movement and timing may indicate an important biological response (Wethington and Russell 2003). Seasonal differences in migration route, timing, or variation may be explained by different selective pressures across seasons, including physiological condition post breeding, resource distribution, habitat heterogeneity, and the seasonality of weather events. For example, conditions in the wintering grounds could significantly impact spring migration timing (Robson andBarriocanal 2011, Paxton et al. 2014). Faster spring migration, observed in three of the five species, may be facilitated by taking advantage of favorable tailwinds (La Sorte et al. 2014b), while the generally slower autumn migration may depend on tracking resource distribution (Leu and Thompson 2005) or on sex/age-related differences in migration timing (Alerstam 2011). The addition of juveniles in the autumn may partially explain the more diffuse autumn migration and higher variance in the date for end of autumn migration (Phillips 1975, Kodric-Brown and Brown 1978, Carpenter et al. 1993, Wethington and Russell 2003. Annual variation in resource production, timing, and weather patterns are all likely to influence individual hummingbird migration behavior (Russell et al. 1994, Rappole andSchuchmann 2003), and may be responsible for the emergent population-level flexibility in the longitudinal component of migration and in migration timing. Studies at hummingbird migration stopover sites indicate that individuals may choose stopover locations based on a combination of nectar supply, competition, and suitable weather conditions (Kodric-Brown and Brown 1978, Russell et al. 1994, Wethington et al. 2005, McKinney et al. 2012 in addition to making adjustments to changing physiological conditions. For example, hummingbirds often use torpor as a method to conserve energy during the night (Lasiewski 1963, Carpenter and Hixon 1988, Powers et al. 2003, but this strategy has limited energetic benefits at extreme high or low temperatures (Hainsworth andWolf 1970, Carpenter 1974). Further, the temporal scale at which hummingbirds respond to changes in resource availability or weather extremes is hours (Beauchat et al. 1979, Powers andConley 1994), which requires them to have strategies to cope with high metabolic demands or to move to locations that are more suitable. Since hummingbird migration is challenging to study across large spatial scales, further research is needed to understand the scale of individual interactions and response to environmental factors versus the scale of migration movement.
Differences in total migration distance and regional heterogeneity along the migration route may also influence observed variation among species' migration patterns (Tottrup et al. 2008, La Sorte et al. 2014a). Long-distance migrants, S. rufus and S. calliope, had faster migration speeds and generally lower annual variation in migration routes and timing than shorter-distance migrants, especially during spring migration. The relatively large variance in migration patterns of shorter-distance migrants could be influenced by hummingbirds' abilities to adjust locations based upon environmental conditions coupled with their weaker constraints on timing for arrival at their breeding grounds compared to longer-distance migrants (Russell et al. 1994, Jenni and Kery 2003, Tottrup et al. 2008, Stanley et al. 2012. However, it is also possible that our observation of slower movement and greater variation in movement of the shorter-distance migrants, including A. alexandri and S. platycer-cus, may be influenced by inclusion of data from more sedentary resident populations in the southern portion of their ranges (central and southern Mexico and Guatemala; Schuchmann 1999, Baltosser and Russell 2000, Camfield et al. 2013. Species migrating in the eastern flyway of North America, including A. colubris, experience more homogeneous biotic environments and topography than in the western flyway, which allows them to more easily track productivity (La Sorte et al. 2014a), and thus exhibit greater variation in population-level migration route among years.
There are several important considerations when using citizen-science observer data, including the number of observation years and potential bias due to the behavior of bird-watchers. We were limited to six years of observation data, and low annual variation in migration patterns may indicate that there was not enough climate variation in our time series to observe large shifts in migration patterns. In terms of bias, a high number of individuals from species typically found in the west, especially S. rufus and A. alexandri, were observed in the southeastern United States. Southeast banding stations also have captured western hummingbirds (USGS Bird Banding Laboratory 2014), but generally only in the winter months and some individuals are later recaptured within the western flyway (Appendix A). This evidence suggests that S. rufus has increased its overwintering population in the southeast since the 1970s (Conway andDrennan 1979, Hill et al. 1998), which could suggest the formation of a secondary migration route, potentially in response to changing climate. Hence, while evidence is limited, and potentially confounded by the behavior of birders, citizen science data have the potential to document important distributional shifts in birds.
Our approach demonstrates that citizen-science data can be used to gain a broad-scale temporal perspective on population migration, without incurring costs associated with tracking individual birds. We were able to study migration in a way that is currently not possible using individual bird data (e.g., banding), especially for hummingbirds. We suggest that eBird and other citizen-science datasets will be useful for tracking seasonal and annual changes in migration, and can be integrated with other datasets (e.g., remote sensing of vegetation phenology) to assess mechanism and to improve predictions under global change scenarios.

ACKNOWLEDGMENTS
We would like to thank D. Fink for sharing his code for the GAMM analysis. We would like to thank the Graham lab, the hummingbird NASA project collaborators, and the Cornell Laboratory of Ornithology eBird team for continued discussions on the data, statistical analysis and citizen science, and two anonymous reviewers. We especially thank all of the volunteer birders and banders whose data was used in this study. The project was funded by NASA grant NNX11AO28G. Hummingbird banding and encounter data were downloaded from the Bird Banding Laboratory on 11 April 2014 (http://www.pwrc.usgs. gov/BBL/homepage/datarequest.cfm). We filtered the data to exclude records with incorrect months of the year (e.g., months that were .12), incorrect days of the month (e.g., days . 31), individuals that were not identified to the species-level, and that were missing latitude or longitude.

LITERATURE CITED
Several western species of hummingbirds were banded and/or recaptured in the eastern flyway (east of 1038 W longitude) by volunteers submitting data to the Bird Banding Laboratory (BBL). The number of banders in the southeastern United States has increased in recent years (Fig. A1;Bird Banding Laboratory,. Thus, it is difficult to determine if the increased number of observations of western species, particularly rufous (Selasphorus rufus) and black-chinned (Archilochus alexandri ) hummingbird, in the east are due to increased banding stations (e.g., increased detectability), formation of a secondary migration route, an ephemeral response to increased feeding stations, and/or long-term response to climate or land use changes. For S. rufus, there were 187 individuals that were captured only in the eastern flyway, and 8 individuals that were captured in the eastern flyway, but later recaptured in the western flyway. For A. alexandri, there were 80 individuals captured only in the eastern flyway and 2 individuals that were captured in the eastern flyway, but later recaptured in the western flyway (Fig. A2).
v www.esajournals.org    Table B2. Results using all the eBird observations 2008-2013 (e.g., not excluding data by flyway). Mean (and standard deviation) for migration variance explained by Julian date for latitude and longitude separately, with year as a random factor looking at the entire route, and at spring and autumn separately. Data used the western flyway data only for the western species.    (2008)(2009)(2010)(2011)(2012)(2013) across the study area at a daily timestep for each species (using all the data). Vertical lines represent the estimated start of spring migration (blue), the date maximum latitude was reached (green) and the end of fall migration (orange), using all data. Most of the observations fall between the blue and orange lines denoting migration period. Note that the green line often falls close to a dip in overall observation, which might be expected, assuming that detectability decreases during the nesting period. Fig. B4. Standard error for latitude and longitude of the daily population-level centroids, from the GAM analysis using the eBird observation data from the western flyway (and all data from ruby-throated hummingbird) and from all years. Standard error for latitude and longitude are strongly correlated, with generally higher error in predicted longitude. Before 2008 (red points), error along both axes was significantly higher and less predictable.