Evaluating the efficacy of eBird data for modeling historical population trajectories of North American birds and for monitoring populations of boreal and Arctic breeding species
 
Évaluation de l'efficacité des données eBird pour la modélisation des trajectoires historiques des populations d'oiseaux nord-américains et pour la surveillance des populations d'espèces nichant dans les régions boréales et arctiques

• Abstract
• Introduction
• Methods
  • Species, areas, and dates included
  • Data filtering and statistical models
  • Assessment of population trajectories
• Results
  • Seasonal eBird trends and trajectories
  • Spruce budworm specialists
  • eBird vs BBS and CBC annual indices and trends
• Discussion
• Responses to This Article
• Acknowledgments
• Literature Cited

INTRODUCTION

Conservation and management decisions should be made using the best available knowledge regarding the population status of the species of concern. Since the mid-1960s, population sizes of many bird species in the U.S. and Canada have been reliably tracked using data from the North American Breeding Bird Survey (hereafter BBS; Sauer et al. 2013). Species not effectively monitored using that protocol include nocturnal, rare, cryptic, secretive, shorebird, water-bird, marsh-bird, and boreal and Arctic breeders (e.g., Betts et al. 2007, Harris and Haskell 2007, NABCI 2012, ECCC 2014a). For these species, other sources of data are used to model population trends, such as the Christmas Bird Count (hereafter CBC), bird observatories in the U.S. and in the Canadian Migration Monitoring Network, or various taxon-specific surveys, e.g., Waterfowl Breeding and Population Habitat Survey (hereafter WBPHS), the Mid-Winter Waterfowl Survey (hereafter MWWS), the International Shorebird Survey (hereafter ISS) and related surveys, the Nocturnal Owl Survey, the Nightjar Survey Network, the Marsh Monitoring Program, and Colonial Waterbird Counts. Broad scale syntheses of national and international bird population trends such as the Partners in Flight Landbird Conservation Plan, the State of North America’s Birds, and the Status of Birds in Canada rely upon data from the aforementioned surveys, and influence conservation and management decisions and funding strategies (ECCC 2014a, NABCI 2016, Rosenberg et al. 2016). Despite all of these diverse efforts, two major deficiencies remain: population trend data only date back to the 1970s for most species, and boreal and Arctic breeding species are poorly monitored. Citizen science data, in particular the eBird project (http://www.ebird.org), could potentially help remedy these deficiencies.

Estimating historical population change is challenging given that most long-term monitoring programs were initiated in the late 1960s or more recently. The BBS is widely considered the most reliable source of population data in the U.S. and Canada, but it was only initiated in 1966 (Sauer et al. 2013). The CBC, the longest-running bird survey in the world, was launched in 1900, but is limited to species that winter in North America (Butcher 1990a). However, the early years of the CBC may not provide sufficient data to model population trends with confidence (Bock and Smith 1971, Schreiber and Schreiber 1973, Niven et al. 2004). Numbers of CBC count circles remained relatively low until the 1950s, which is also when methodologies were standardized between circles (Butcher 1990a, Butcher et al. 2005). Consequently, most contemporary studies that use CBC data to compute population trends focus on years subsequent to 1950 (e.g., Butcher 1990b, Sauer et al. 1996, 2004, Niven et al. 2004, Butcher et al. 2005, Soykan et al. 2016). The oldest bird observatory in North America, Long Point Bird Observatory, was established in 1961, and several other major bird observatories were founded that decade (Point Reyes, Whitefish Point, and Manomet). Of the taxa specific surveys, the longest running are the Mid-winter Waterfowl Survey and WBPHS (both 1955; USFWS 2016, 2017), the American Woodcock Singing Ground Survey (1968), and the ISS (1974). Most of the other taxa specific surveys were initiated in more recent times and do not provide long-term indices of populations. In short, there is very little population trend data for bird species prior to the 1950s, and for species that winter outside of North America, data can only be analyzed reliably back to circa 1970.

We argue that to fully assess population trends of any bird species, it would be beneficial to be able to model trajectories prior to 1970. Without a historical context for bird populations and an understanding of long-term population cycles, setting target populations for species may be somewhat arbitrary. For instance, for species with large-amplitude population cycles, such as spruce budworm (Choristoneura fumiferana) specialists, a longer time series is necessary to determine their status because the BBS only dates back to the last spruce budworm outbreak in the 1970s (Patten and Burger 1998). Similarly, populations of species that were adversely affected by DDT bottomed out in the 1970s, so there is no point of reference for historical population size.

Monitoring populations of bird species that nest in the Arctic and boreal forest also remains a challenge because their breeding areas are vast and remote with limited road access. Monitoring these species on the breeding grounds is logistically and financially challenging. For boreal and Arctic species that winter largely in the U.S. and Canada, the CBC provides insight into their population status (Dunn and Sauer 1997, Niven et al. 2004, Soykan et al. 2016). For boreal and Arctic species that winter further south, trend information is more limited. Data from bird observatories that capture and record observations of these species during spring and fall migrations provide the best published trend information; however, it is still an unresolved problem whether the data from these few sites are representative of populations range-wide (Francis and Hussell 1998, Lloyd-Evans and Atwood 2004, Dunn et al. 2006, Crewe et al. 2008, 2016).

It has been shown elsewhere that data from the citizen science project eBird can be used to model population trends and that these trends broadly agree with BBS in areas and time periods with many eBird checklists (Walker and Taylor 2017, Horns et al. 2018a). The eBird project collects checklists made by birdwatchers using a suite of both generalized and specialized protocol types that describe almost any potential birdwatching scenario (Sullivan et al. 2014). The checklists are permanently archived, and the data are accessible through the eBird website in summarized and raw formats. Though launched in 2002, the flexible protocol types used by eBird allow for historical checklists to be entered, provided that at a minimum a date and location are specified. Globally there are currently (July 2020) > 45.8 million checklists in the eBird database, with ~36.9 million of these from the U.S. and Canada. Rates of checklist submission continue to grow exponentially.

Walker and Taylor (2017) estimated long-term (1970-2015) population trajectories from eBird data that produced trends comparable to those from BBS data for 22 species that reach the northern edge of their range in southern Ontario, Canada, a region with high numbers of eBird checklists. That result prompted the question of whether a similar methodology (Walker and Taylor 2017) could be applied at a broader scale and over a longer time series. Horns et al. 2018a used checklists from a broader scale, the contiguous U.S., to model trends for 574 species, but focused on a shorter timespan (1997-2016) when numbers of checklists were high. Horns et al. 2018a found trend estimates broadly agreed between eBird and BBS, but there was a large variation between species in how well the trends agreed, and many species showed trends with opposite signs (Fogarty et al. 2018, Horns et al. 2018a, b). Although both Horns et al. 2018a and Walker and Taylor 2017 found that underlying patterns of abundance were captured for many species, refinement of these techniques is necessary to produce reliable population trends using eBird data. Although numbers of historical checklists in eBird are low and geographically biased, we anticipated that range-wide aggregation of checklists might produce sample sizes sufficient to determine if there is meaningful information in the historical data. Furthermore, aggregation of lists across the continent during spring and fall migration could provide representative samples despite pronounced geographic sampling bias, because most migratory species pass through at least some areas with high eBird coverage. Additionally, during the migratory period, birds far away from their breeding ranges pass through narrower migratory pathways, facilitating the estimation of continent-wide population indices.

The objectives of this study were twofold: (1) to explore the current eBird dataset to determine how far into the past population trajectories can be estimated and (2) to produce range-wide population trajectories of boreal and Arctic breeding species for which few other sources of data exist. Because both of these objectives aim to fill information gaps in existing monitoring strategies, direct comparisons to other data sources (to assess the validity of the method) are limited or nonexistent. We therefore reasoned that we could assess such validity via indirect methods. First, we predicted that if the methods were effective at dealing with biases inherent in the eBird data set, we would (across a suite of different species) see comparable population trajectories derived from data from different seasons (spring, fall, and winter). Second, we anticipated that patterns in population change of species with well-publicized patterns of population decline and rebound associated with the use and banning of DDT, e.g., Peregrine Falcon (Falco peregrinus), Bald Eagle (Haliaeetus leucocephalus), would broadly show the known patterns. Third, we predicted that species known to be associated with spruce budworm outbreaks, i.e., Cape May Warbler (Setophaga tigrina), Bay-breasted Warbler (Setophaga castanea), Tennessee Warbler (Leiothlypis peregrina), and Evening Grosbeak (Coccothraustes vespertinus), would show patterns of increase and decrease broadly in concordance with the known timing of those outbreaks in boreal Canada (Patten and Burger 1998). To further evaluate the methodology, we also undertook comparisons with published trends from BBS and CBC data for the years and species for which trends were available. Although not unequivocal, if the results from the models produced the expected patterns, we argue that longer-term analysis of population change may be possible, particularly with the inclusion of additional historic checklists into the eBird database. Such an effort would provide additional and valuable insight into population trends for boreal (and other) species over a longer time frame than was heretofore possible.

METHODS

Species, areas, and dates included

Twenty-two species of songbirds that breed in the boreal forest or the Arctic, four songbird species that breed south of the boreal forest, and two raptor species were selected for analysis based on several criteria (Table 1). The three warbler species most associated with spruce budworm outbreaks (Cape May Warbler, Tennessee Warbler, and Bay-breasted Warbler) and Evening Grosbeak were selected initially with the expectation that populations of these four species would show similar long-term cycles, and should be positively correlated with spruce budworm outbreaks (Patten and Burger 1998, Venier et al. 2009). The additional 18 boreal or Arctic species were selected using the Status of Birds in Canada website (ECCC 2014a) and the Partners in Flight Landbird Conservation Plan 2016 Revision (Rosenberg et al. 2016). Nine of these species were selected because trend reliability scores were either low or deficient in the Status of Birds in Canada 2014. These reliability scores were based on analysis of data from various surveys that were appropriate for the species (BBS, CBC, etc.), and were scored on the precision of estimated trends, the proportion of seasonal geographic range covered by the surveys, and the reliability of the survey for the species in question (ECCC 2014a). Eight additional boreal and Arctic species were selected because they were labeled as “watch list” or “common species in steep decline” by the Partners in Flight Landbird Conservation Plan 2016 Revision (Rosenberg et al. 2016). Because the eight species selected were presumed to be declining a priori, we also selected three species of northern-breeding songbirds for which BBS data indicated survey-wide increases to ensure that results showing declines were not due to systematic bias: Magnolia Warbler (Setophaga magnolia; boreal forest), Northern Parula (Setophaga americana; eastern forest), and Black-throated Blue Warbler (Setophaga caerulescens; eastern forest; Sauer et al. 2017). Two additional eastern forest species with reliable trends in BBS were selected, of which Blue-winged Warbler (Vermivora cyanoptera) was stable and Golden-winged Warbler (Vermivora chrysoptera) was decreasing.

Finally, two raptor species were selected for which long-term population trajectories are widely known, Bald Eagle and Peregrine Falcon. Populations of both species were critically low in North America in the 1960s and 1970s due to DDT contamination, but rebounded after DDT was banned and recovery plans were implemented (USFWS 2003, 2009).

The eBird basic datasets from the August 2016 data release were downloaded for the United States (excluding Hawaii) and Canada, and filtered to include all complete checklists reporting greater than one species (Cornell Lab of Ornithology 2016). The checklists were divided into seasons and analyzed separately for each species. Spring migration and breeding (hereafter Spring) included checklists from March through June (inclusive), fall migration included lists from August through November (inclusive), and for species that winter in the United States and Canada lists from December through February (inclusive) were used. Checklists from 1928 through 2015 were included in the models. Data from 2016 were incomplete because the fall and winter seasons had not yet occurred, and checklists prior to 1928 were considered too sparse to be included.

For each species and season, unique geographic areas and date ranges were used to eliminate checklists from areas and dates from which that species rarely or never occurred. For each state/province, the total number of checklists submitted and the number of checklists reporting each species were calculated by season. If, across the entire time series, there were 10 or fewer checklists reporting a species or the frequency of checklists reporting that species was less than 0.001 in a given state/province, then checklists from that entire state or province were excluded from the analysis for that species and season. For the remaining states and provinces, 1% and 99% date quantiles were calculated for occurrence of each species by season, and all checklists from that state/province outside of the date quantiles were discarded if the date quantiles were more than 10 days from the beginning or end of the season (to avoid discarding data for species that were already or still present in the study area at the beginning or end of the season). Finally, checklists from each state/province were grouped into 40 × 40 km grid cells and cells in which a species had never been reported during a season across all years were excluded from the analysis for that species and season.

Data filtering and statistical models

The resulting data were then further filtered in a similar manner as Walker and Taylor (2017). To facilitate the inclusion of older historical and incidental checklists, all complete checklists of greater than one species were included regardless of effort information included. Shared checklists were reduced to a single checklist and assigned to the first observer associated with the list. Specialized checklist protocol types were reclassified to correspond to one of the three primary protocol types (traveling count, stationary count, or incidental) based on the effort information included, but several protocol types were excluded that did not fit within the framework of the study (see Appendix 1 for more details). Although current eBird submissions force all incidental lists to be flagged as incomplete, the original checklist submission process allowed for complete incidental lists, so these were included in the analysis.

Bayesian multilevel regression models were fit to the data for each species separately in both spring and fall, and in winter for those species which winter in the U.S. and Canada using package brms (version 1.6.1) in R Statistical programming language (version 3.2.5; R Core Team 2016, Bürkner 2017). The function brm in package brms uses syntax similar to that of package lme4 (Bates et al. 2015) in R, but uses the program Stan on the back end to specify the model (Stan Development Team 2017). Stan uses the Hamiltonian Markov chain Monte Carlo method and No-U-Turn Sampler to implement the models (Hoffman and Gelman 2014). We employed this method because models using lme4 converged inconsistently, which appeared to be an issue with using lme4 for exceptionally large datasets. Parameter estimates from models using brms were very similar to those from lme4, but the models converged for all but three of the species/season combinations.

Probability of occurrence on a checklist with a Bernoulli error distribution was used as the response variable in each model to facilitate incorporation of many checklists pre-2000 that lacked count data. The use of occurrence on a checklist as a response instead of abundance data may produce biased results for some species, especially for common and flocking species because the change in occurrence on checklists may not capture changes in abundance for those species. However, given the differences in sample size between checklists with occurrence data and checklists with count data, we opted to use occurrence for the response.

Year was fit as a factor to generate independent indices for each year, allowing for expected cyclical and other nonlinear patterns of occurrence across the time series. Years prior to 1970 were combined into bins of three years to boost sample sizes and decrease computation times. Any year or year bin with fewer than 20 checklists or zero observations of the target species was removed to facilitate model fitting. Although removing years with zero observations of the target species may introduce some negative bias in trend estimates, we believe that data from these years are simply deficient for producing annual indices. Very few years, if any, were removed for most species, and years removed were mostly from the 1920s and 1930s. The reference year for the year factor was set to 2015 because that was the year with the most data for each species. An interaction between latitude and a quadratic term for date accounted for differences in timing of peak migration based on latitude within a season because the data spanned such large geographic areas. The number of species included on a checklist was used as a surrogate for effort (as in Roberts et al. 2007, Szabo et al. 2010, Walker and Taylor 2017, and Horns et al. 2018a), and the interaction with protocol type allowed for differing rates of species accumulation by protocol type.

An observer score was calculated using similar methods as used by Kelling et al. (2015), by using a random slope term in a global mixed effects model to estimate a separate species accumulation rate for each observer. Generalized linear mixed effects models using Poisson error terms were fit within each season using package lme4 (version 1.1-11) in the R (version 3.2.5) Statistical programming language (Bates et al. 2015, R Core Team 2016). The number of species on a checklist was the response variable. Only checklists including effort were used to generate the observer scores. Species accumulation was modeled by including log(checklist duration + square root of checklist duration) as a fixed effect. A quartic polynomial for time of day was also included to account for assumed peaks in species detectability at dawn and dusk. A random intercept term was included using each unique combination of state/province and month to account for differing number of species available for detection by region and month. Finally, a random slope term was included allowing for a differing rate of species accumulation for each observer. The random slope effects (observer scores) were extracted for each observer and were included as continuous fixed effect in the occurrence models. Observers that had never entered a checklist with effort were assigned an observer score of 0, which was equivalent to the overall mean rate of species accumulation. Assigning mean observer scores for observers who had never entered checklists with effort could introduce bias because historical lists were more likely to lack effort information. Thirteen percent of observers in the pre-1970 years lacked effort information compared to five percent post-1970. However, we did not want to expunge these historical checklists from these observers because sample sizes were already low. Much of the effect of observer was likely accounted for by including the number of species on a checklist as a fixed effect already, so we were not overly concerned with setting these missing data points to the mean for all observers. To verify that the observer score models worked as intended, we compared distributions of observer scores for all eBird regional editors in the U.S. and Canada to scores from the other observers. eBird regional editors are expert bird-watchers who vet eBird data submissions in their geographic area of expertise. We found that the mean score of eBird regional editors was higher and the distribution much narrower than those of all other eBird observers, suggesting that the method was generally effective.

A random effect for geographic area (40 × 40 km cell) was included in the occurrence models to account for differences in frequency of occurrence by region and to account for geographic heterogeneity in checklist submission rates over time (Roy et al. 2012, Isaac et al. 2014). Although smaller grid cells or specific locations might be preferable to account for more habitat specific differences in occurrence, the numbers of random effects to be estimated at that scale exceeded the capacity of our computing infrastructure. Including a random effect of geographic area allows for checklists from anywhere in the species’ range to be incorporated into the models. The random intercept term can be thought of as a measure of how frequently the species of interest was recorded in each grid cell and is strongly influenced by recent years with prolific data. Although this method does not allow for estimating regional or grid cell differences in population trajectories, it facilitates the incorporation of scattered checklists from historic times. Annual indices are most reflective of areas with the highest densities of checklists, but are still influenced by checklists from areas with lower densities. An alternative methodology would be to limit the study to areas with the highest numbers of historical lists, but this would eliminate much information from areas with high numbers of contemporary lists and few historic lists. Figure 1 depicts the geographic densities of checklists in eBird over three time periods (pre-1970, 1970-2000, post-2000) and in each season as of August 2016, and Appendix 2 displays the geographic density of checklists used for each individual model based on our filtering criteria over the same time periods. Appendix 3 provides model outputs regarding the geographic random effect term for each model, including effect size, error, rhat values and the number of 40 x 40 km grid cells used.

In summary, fixed effects terms in the models included: year (fit as a factor), an interaction between latitude and a quadratic polynomial fit to date, an interaction between the log of the number of species on a checklist and protocol type, and an observer skill score. A random effect was included for location using 40 × 40 km grid cell as the grouping variable. Equation 1 provides a summary of the model in the syntax used by the brm function used in R.

(1)

The defaults of improper flat priors were used for population level effects and student t priors with three degrees of freedom used for the group level effects. The default number of iterations (1000 warm-up, 1000 sample) were used in the models. The number of chains was set at two, which ran simultaneously on two processing cores.

Models for White-crowned Sparrow (Zonotrichia leucophrys) in winter, Bald Eagle in fall, and American Tree Sparrow (Spizelloides arborea) in spring did not converge with rhat values above 1.10, but the other 65 models converged. Rhat values were 1.13 for the intercept of White-crowned Sparrow in Winter, and 1.11 for the intercept of Bald Eagle in Fall and 1.11 for the grid cell random effects for American Tree Sparrow and White-crowned Sparrow. These were some of the models with the largest sample sizes (see Appendix 4 for sample sizes for each model). Annual indices generated from these three models are included in the results, based on the close agreement between trajectories from these models and those from other seasons for the same species, but should be interpreted with caution (see Appendix 5 for plots of annual indices from each model). Smith’s Longspur (Calcarius pictus) data were much sparser than the other species included and many historical years lacked sufficient data for modeling. The species was still included in the analyses using the years with sufficient data based on the filtering criteria above, because it is considered a data deficient species based on criteria in ECCC 2014a. Posterior predictive checks (ppc) were performed to assess the fit of the models by using the pp_check function in package brms, a random subset of 3000 checklists from each model as the data, and 1000 samples from the posterior distribution for each. No ppc plots examined indicated that model fits were unreasonable.

Assessment of population trajectories

Annual indices and their estimated errors were calculated for each species in each season using fitted values from the models on the response scale of probability of occurrence on a checklist. The average latitude and average number of species per checklist were used as the new data for the fitted values, and the date was set as the peak date of occurrence. Additionally, to obtain population-level estimates, observer score was set to 0, the protocol type was set to traveling count, and random effects of geographic location were not estimated.

Several comparisons were made to evaluate the annual indices generated by the models, including: correlations between seasonal eBird trends and trajectories, correlations between spruce budworm specialists and area defoliated by spruce budworm, and correlations between eBird and BBS annual indices and trends.

Seasonal eBird trends and trajectories

Pearson correlations were calculated between seasons using annual indices from the eBird models to ascertain whether the trajectories agreed for each species. To compare magnitudes of underlying trends, models were fit through the annual indices for each species and season for the entire time series using beta regressions with logit links using the R package betareg (Cribari-Neto and Zeileis 2010). Trends were compared between seasons across species using Pearson correlations. Trends from eBird annual indices were also calculated from 1970 to present for comparison with BBS and CBC, and those trends were also compared between seasons and across species using Pearson correlations. To further assess the amount of agreement in population trajectory, beta regressions were fit for each species using the indices for all seasons combined with a fixed effect for season, to determine if residuals from these models were correlated between seasons. This allowed us to compare whether cyclical or other nonlinear patterns in the indices were correlated between seasons using Pearson correlations.

Spruce budworm specialists

Data on extent of area defoliated by spruce budworm in the boreal forest dating back to 1939 were obtained from the National Forest Pest Strategy Information System (NRCAN 2014). The number of hectares of moderate, moderate to severe, and severe defoliation were summed and used as an index of spruce budworm defoliation, and areas of light and trace defoliation were excluded. Years prior to 1970 were binned into groups of three years to match the eBird indices, and the mean defoliation index across the years in each group was used. Spearman rank correlations were calculated between annual indices of spruce budworm defoliation and annual indices from eBird models for each species in each season because the indices were on different scales.

eBird vs BBS annual indices and trends

Comparisons were made between annual indices from eBird and BBS data for those species that were detected frequently enough on BBS routes to produce indices. Survey-wide annual indices were downloaded from the U.S. Geological Survey’s BBS website (Sauer et al. 2017). Annual indices were downloaded from the Canadian Trends Website (ECCC 2014b) for three species that were not modeled by Sauer et al. (2017) because of poor coverage on BBS routes. Breeding Bird Survey indices from the years 1970-2015 were compared to annual indices from the eBird models using weighted Pearson correlations, with the inverse of the sum of the variance around each pair of BBS and eBird annual indices as the weights. Because BBS indices are shrunk toward a trend line as opposed to being estimated as individual factor levels as in the eBird models, they underestimate interannual variation in the data. As such, they, and the associated p-values, should be interpreted with caution. We provide them primarily to see whether the two types of indices broadly followed the same pattern.

Trends on a scale of annual percent change and associated confidence intervals were calculated for each species and season by fitting beta regressions with logit links using the R package betareg (Cribari-Neto and Zeileis 2010) to the annual indices from the eBird models between 1970 and 2015, for comparison with BBS and CBC data. This was done with the annual indices rather than raw eBird data because the imbalance in checklists between recent and historic times would yield trends that were almost entirely influenced by the most recent years. For comparison, trend estimates from survey-wide BBS data were downloaded from the BBS website (Sauer et al. 2017), and trends for three additional species not modeled by Sauer et al. (2017) were obtained from the Canadian Trends Website (ECCC 2014b). For species that winter in the U.S. and Canada, trends calculated from CBC data were also compared (Soykan et al. 2016). Pearson’s correlation coefficient was calculated to test if the magnitude of trend estimates from eBird data were correlated to those from BBS and CBC data. We acknowledge that annual indices and trends were calculated using different underlying model structures and response types (count vs occurrence) for the BBS and CBC than we used for the eBird data; however, we were primarily interested in the general level of agreement between datasets rather than making a formal comparison. A formal comparison of trends derived from the three datasets should focus on geographic regions with high numbers of eBird checklists and good coverage of BBS routes and CBC circles.

RESULTS

Seasonal eBird trends and trajectories

Visual inspection of plots of annual indices for each species and season indicated that there is some level of agreement between population trajectories from different seasons in eBird for most species (Appendix 5). For each species, we fit LOESS smooths with spans of 0.2 through the annual indices for each season to aid interpretation of population trajectories over time. Correlation coefficients between spring and fall annual indices were positive for all but two of the species analyzed, and were positive for all species when spring and winter, and fall and winter were compared (Table 2). The mean correlation coefficient across species was 0.36 between spring and fall, 0.68 between spring and winter, and 0.48 between fall and winter. Trends calculated from the annual indices using the entire time series were highly correlated across species between spring and fall (r = 0.70, n = 28, P < 0.0001), spring and winter (r = 0.83, n = 12, P = 0.0007), and fall and winter (r = 0.89, n = 12, P = 0.0001; Fig. 2, Table 3). Trends calculated from annual indices from 1970 to 2015 for comparison with BBS and CBC were also strongly correlated across species between spring and fall (r = 0.87, n = 28, P < 0.0001), spring and winter (r = 0.95, n = 12, P < 0.0001), and fall and winter (r = 0.89, n = 12, P = 0.0001; Fig. 3, Table 4). Correlations of residuals from trend models fit through the seasonal data were largely in the positive direction, and some species with obvious nonlinear patterns over time showed strong correlations (Table 5). Overall, the correlations between residuals suggest that for many species, the eBird dataset is able to track population trajectories beyond what was explained by the underlying trends.

Spruce budworm specialists

The spruce budworm specialists in this study (Bay-breasted Warbler, Cape May Warbler, Tennessee Warbler, and Evening Grosbeak) all showed positive correlations with the number of hectares defoliated by spruce budworm, and, of the 28 species modeled, were the species with the strongest correlations (Table 6; Patten and Burger 1998, Venier et al. 2009). Smoothed population trajectories of these three warbler species and Evening Grosbeak showed some noticeable peaks in occurrence during outbreaks of spruce budworm (Fig. 4). For Magnolia Warbler, a species that has been shown in several studies to have a negative response to spruce budworm, there was no correlation between eBird annual indices and area defoliated by spruce budworm (Patten and Burger 1998, Venier et al. 2009; Table 6). Several other species showed moderate correlations between eBird data and spruce budworm that were not as strong as the spruce budworm specialists, but could indicate a relationship (Table 6). Several species, most notably the two species of raptors, White-crowned Sparrow, Northern Parula, and Least Flycatcher (Empidonax minimus; fall model), were negatively correlated with the spruce budworm index.

eBird vs BBS and CBC annual indices and trends

Correlations between annual indices from the eBird models and BBS data were mostly positive and there were few negative correlations (Table 7). The mean correlation coefficients across species were 0.46 for spring indices, 0.32 for fall indices, and 0.55 for winter indices. eBird annual indices with smoothed trajectories plotted adjacent to those from BBS from 1970-2015 for each species are presented in Appendix 6.

Trends calculated using the annual indices from the eBird models are presented alongside those from the BBS and CBC in Table 5. Pearson’s correlation coefficients between eBird and BBS trend estimates were 0.73 (n = 25, P < 0.0001) for spring, 0.64 (n = 25, P = 0.0005) for fall, and 0.81 (n = 9, P = 0.0084) for winter (Figs. 5A and 5B). The Pearson’s correlation coefficient between winter eBird trends and CBC trends was 0.64 (n = 12, P = 0.0252), but showed obvious bias between datasets with eBird data estimating more-negative trends (Fig. 5C). Residuals were not compared between datasets because the underlying trend models were different between eBird, BBS, and CBC.

For the boreal and Arctic breeding species, population trajectories mostly corroborated the trends from the BBS and CBC. For the four species with high trend reliability from either the BBS, i.e., Olive-sided Flycatcher (Contopus cooperi), or CBC, i.e., Evening Grosbeak, Harris’s Sparrow (Zonotrichia querula), and Rusty Blackbird (Euphagus carolinus), there was agreement with eBird data, though the BBS indicated that Olive-sided Flycatchers were declining at a faster rate than was shown using eBird data. Trends for the eight boreal or Arctic breeding species with medium trend reliability also largely agreed between eBird and BBS, with overlapping trend estimates for Canada Warbler (Cardellina canadensis), Wilson’s Warbler (Cardellina pusilla), Cape May Warbler, and Connecticut Warbler (Oporornis agilis). Trends for American Tree Sparrow, Bay-breasted Warbler, and Least Flycatcher did not agree between datasets, but Least Flycatchers showed declines that disagreed only in magnitude. For the 10 species with low or deficient trend reliability from BBS and CBC, it was difficult to assess agreement of trends because of wider credible intervals that often overlapped 0. None of these 10 species showed agreement between eBird and BBS based on a combination of sign and magnitude of trend.

The four warbler species that breed just south of the boreal forest showed some agreement between eBird and BBS trends. For the species with increasing trends in BBS, Northern Parula showed similar positive trends in eBird data, and Black-throated Blue Warblers eBird showed stable trends whose credible intervals overlapped zero. The spring eBird trend for Blue-winged Warbler was similar to the BBS trend in magnitude and precision, but the credible interval overlapped zero for the BBS trend. Golden-winged Warblers showed similar negative trends in eBird and BBS data, but fall eBird data showed a stable trend.

The two raptors, Bald Eagle and Peregrine Falcon, showed similar trends between eBird and BBS, and long-term trajectories matched the known pattern of a historical decline followed by a rebound after the use of DDT was banned.

DISCUSSION

Results of the models suggest that eBird data could play an important role in estimating population trajectories for species that are poorly represented in traditional and more-structured bird monitoring efforts. Similarly, the eBird database provides a unique opportunity to aggregate and use checklists from historic times, which could be used to model population trajectories prior to the advent of current long-term surveys and provide a broader context for interpreting more recent changes in population size. As with the CBC, eBird may not currently have sufficient data to model population change reliably for most species prior to circa 1970, which was evident in the large credible intervals around annual indices from the eBird models in early years. However, the eBird dataset is continuously growing and there exists a huge repository of historical checklists in the field notebooks of living and deceased bird-watchers. We suggest that with additional eBird checklists and enhanced statistical modeling techniques that account for known sources of bias in the data, it should be possible to derive improved estimates of historical changes in population trajectories for many species.

The positive correlations shown between seasons when comparing annual indices, trends, and residuals from trend models, suggest that models of eBird data reflect underlying population trajectories, because the datasets were discrete. Although credible intervals for annual indices became much wider further back in time, we suggest this is primarily due to the smaller sample sizes in those years. Because eBird was launched in 2002, all checklists prior to that year were effectively historical checklists that were recorded before eBird protocols were in place. The correlations between seasonal trends and trajectories and between eBird, BBS, and CBC dating back to 1970 were encouraging, suggesting that the methods used for incorporating historical lists were effective at extracting signal from noise. Because the eBird dataset is continuously growing (including historical checklists) estimating historical population trajectories using these (or other) methods should only improve in the future.

Although results indicated positive relationships between seasons for most species, several species showed disagreement in sign or magnitude of trend. Focusing on years since 1970 when data were plentiful, the following seven species showed noteworthy disagreement in seasonal population trajectories: Least Flycatcher, Golden-winged Warbler, Blackpoll Warbler (Setophaga striata), Tennessee Warbler, LeConte’s Sparrow (Ammodramus leconteii), Evening Grosbeak, and Common Redpoll (Acanthis flammea). In general, when seasonal trajectories differ we believe that spring eBird trajectories are likely more representative of the population than fall trajectories for several reasons: the spring eBird dataset contains more checklists than the other seasons, birds are more likely to produce diagnostic vocalizations in spring than in other seasons, plumages of many species are more distinctive in spring, and effects of juveniles and winter mortality have been moderated. Discrepancies for three of these species (Least Flycatcher, Tennessee Warbler, and Blackpoll Warbler) may stem from difficulties in fall identification. One is an Empidonax flycatcher, a group well-known for being difficult to identify in fall when they seldom vocalize. The two warbler species are among several species often referred to as confusing fall warblers because they are relatively nondescript and similar to one another. Many fall Empidonax and confusing fall warblers are left unidentified, but advances in identification of these species may have led to increased reporting rates over time, which would introduce bias. For Golden-winged Warblers, both spring and fall trajectories showed clear but different patterns, with fall indices showing a peak in frequency in the 1990s that was not present in the spring indices. This peak may have an unknown underlying biological explanation. Evening Grosbeaks differed seasonally in the magnitude of the trend estimates, but all seasons showed steep declines. For Common Redpolls and LeConte’s Sparrows, trends agreed between two of the seasons but not the third. Common Redpolls are irruptive finches, and when irruptions occur are most likely to arrive in the northeastern U.S. in late November or December. In most years, the fall season as defined for this study would not likely overlap with the irruption from the Arctic.

Correlations between area defoliated by spruce budworm in the boreal forest and annual indices from eBird data for bird species that respond to spruce budworm outbreaks also suggest that modeling historical data is providing meaningful information on population change. Peaks of occurrence for the spruce budworm specialists largely corresponded to spruce budworm outbreaks, meaning that the models successfully captured population trajectories prior to and after 1970. Additionally, the three warbler species associated with spruce budworm showed recent increases in annual indices that coincide with the current outbreak of spruce budworm. Although peaks in occurrence of these species do not all exactly match the peaks in defoliation, we think this could have a biological explanation considering that the defoliation statistic is summed across the entire boreal forest. Regional peaks in defoliation differed by up to a decade, and so overlap with the extent of each species’ breeding range varied. An in-depth analysis of these differences may be warranted, but was beyond the scope of this study. Several of our other focal species, i.e., Yellow-bellied Flycatcher (Empidonax flaviventris), Blackpoll Warbler, and Canada Warbler, have shown positive responses to spruce budworm in some studies but not in others (Patten and Burger 1998, Venier et al. 2009). In this study, Canada Warbler trajectories from both spring and fall were moderately correlated with spruce budworm indices, but trajectories for Yellow-bellied Flycatcher and Blackpoll Warbler were not. Although many of the other focal species also showed positive correlations to spruce budworm defoliation, this could be expected for several reasons. Venier et al. 2009 noted an overall increase of territories across species during a spruce budworm outbreak, even after data from the specialist species were removed. This may be reflected in trajectories of some of the other boreal species. However, given that data from before 1970 were less influential because of the year bins, we think the long-term decline in defoliated area since 1970 might correlate well with any species that showed similar declines during that same time frame, and could explain why species that breed outside of the boreal forest were correlated with spruce budworm indices. Similarly, species that showed positive trajectories since 1970 such as Bald Eagle, Peregrine Falcon, Northern Parula, White-crowned Sparrow, and the fall model for Least Flycatcher, were all negatively correlated with the spruce budworm data. Although any species with similar declines since 1970 may be somewhat correlated with area defoliated by spruce budworm, the fact that the four species of spruce budworm specialists showed the highest correlations out of the species analyzed supports the efficacy of modeling populations using eBird data.

The annual indices and trends generated using models of eBird data agreed surprisingly well with those from the BBS, considering that many of the species have ranges that only overlap BBS coverage at their extremities. Species with higher reliability scores based on BBS, CBC, or other survey data by ECCC 2014a showed trends more similar to those in eBird than did species with low reliability (Fig. 5). Although there was no apparent systematic bias between eBird data and BBS data, eBird data produced more negative trends than CBC data.

The eBird dataset provides an opportunity to assess the status of boreal species with low or deficient trend reliability from BBS and CBC from a new perspective. Of the 10 species analyzed with low or deficient trend reliability, 8 have migratory routes that pass through areas with high eBird coverage, which highlights the utility of the eBird dataset. Some of these species, such as Gray-cheeked Thrush (Catharus minimus) and Lapland Longspur (Calcarius lapponicus), showed long-term declines in frequency of occurrence that have not been detected by other surveys. For other species, such as Blackpoll Warbler and Common Redpoll, trajectories from eBird data were more positive than those from BBS, and did not exhibit the steep declines shown by the latter. This could have important management implications. For instance, based on BBS data, Blackpoll Warbler was recently highlighted as a common species in steep decline by the Partners in Flight Landbird Conservation Plan 2016 Revision (Rosenberg et al. 2016), but long-term trajectories from eBird appeared to be stable. The BBS covers a small fraction of the breeding range of Blackpoll Warblers, whereas the majority of the population passes through areas of the eastern U.S. with high eBird coverage, suggesting the eBird trend may be more representative of the population. Thus, the use of eBird data could play an important role in assessing populations of boreal and Arctic breeding species, especially for those species that migrate through eastern North America or along the Pacific Coast where eBird coverage is the highest. For some species though, the utility of eBird data to monitor populations remains limited. For one species in this study, Smith’s Longspur, low detection rates combined with sparse eBird coverage within its range yielded small annual sample sizes and poor precision in trend estimates. The species was still included in the analyses to demonstrate the lower limits for sample sizes using this modeling technique. Additionally, very little is known about population trajectories for this species based on other surveys so the eBird data could provide the best estimates to date.

The large, nearly range-wide scale of these models during migration allowed for aggregation of checklists from a broad geographic area, which provided the sample sizes necessary to model trajectories in historic times. The 40-kilometer square grid cell that was used to group checklists geographically helped reduce bias from uneven geographic coverage, and the vast amount of data in recent years yielded estimates of frequency of occurrence for these cells. The shortcoming of this method, is that annual indices from different years may be weighted toward different regions, depending on where the checklists originated. This method essentially assumes that there are not regional differences in population trajectories for each species. Currently, there are not enough historic lists in eBird to model regional trends within a global model, though there may be several regions with enough data to model trajectories independently. The long-term trajectories we presented are largely influenced by eBird data from the eastern U.S., primarily along the east coast and in the Great Lakes area. For species whose populations primarily pass through these areas during migration, the models presented likely provide a reasonable estimate of change in frequency of occurrence. For species that pass through or winter in the center of the continent, the trajectories would be biased toward the eastern portion of their range. For the few species analyzed whose ranges reach the Pacific coast, the results would also be biased toward data from California, where there is a high density of historic lists.

The results of our indirect assessments of long-term population trajectories generated by eBird data matched our predictions: population trajectories largely agreed between seasons, Bald Eagle and Peregrine Falcon trajectories mirrored known patterns, and annual indices of species associated with spruce budworm outbreaks were strongly correlated to that of area defoliated by spruce budworm. Additionally, correlations between eBird, BBS, and CBC trajectories were positive despite differences in sampling methodology, statistical analysis, and geographic coverage, suggesting that our methodology was successful at extracting signal from noise in the eBird dataset. More work is needed to determine the reliability and precision of trajectories estimated from eBird data, and to clarify the nature of the relationship between frequency of occurrence on eBird checklists and population size for various species.

The endeavor to model population trajectories prior to 1970 was hampered by the limited number of checklists in those years, but efforts to locate and enter historic lists into eBird could remedy this. We did not attempt to assess the annual sample sizes of checklists that would be necessary to produce reliable estimates. This would not be straightforward considering that the information content of each checklist is so variable based on the amount of effort and number of species detected, and detection rates for each species are different. However, it was evident that the credible intervals around annual indices in the eBird models become much wider before 1970 for most species, so sample sizes used for 1970s in this study may serve as a guideline for the minimum at a continental scale. Regional studies may require far fewer checklists to overcome geographic bias. We believe that efforts to boost historical sample sizes are critical because examination of long-term population trajectories that date back to before 1970 is essential for assessing the status of some species. For example, the spruce budworm specialists all show long-term declines since 1970, but 1970 coincides with the peak of the last spruce budworm outbreak. Only data prior to 1970 address the long-term status of these species, and of other species that have population cycles with large amplitudes. The strength of the relationships between seasons and datasets post 1970 imply that the methodology for modeling trajectories from historic eBird checklists is sound, and that a concentrated effort to increase sample sizes from historic times by an order of magnitude would provide important context for understanding current population changes.

LITERATURE CITED

Bates, D., M. Mächler, B. Bolker, and S. Walker. 2015. Fitting linear mixed-effects models using lme4. Journal of Statistical Software 67(1):1-48. https://doi.org/10.18637/jss.v067.i01

Betts, M. G., D. Mitchell, A. W. Diamond, and J. Bêty. 2007. Uneven rates of landscape change as a source of bias in roadside wildlife surveys. Journal of Wildlife Management 71(7):2266-2273. http://dx.doi.org/10.2193/2006-004

Bock, C. E., and R. S. Smith. 1971. An analysis of Colorado Christmas counts. American Birds 25:945-947.

Bürkner, P.-C. 2017. brms: an R package for bayesian multilevel models using Stan. Journal of Statistical Software 80(1):1-28. https://doi.org/10.18637/jss.v080.i01

Butcher, G. S. 1990a. Audubon Christmas bird counts. Pages 5-13 in J. R. Sauer and S. Droege, editors. Survey designs and statistical methods for the estimation of avian population trends. Biological Report 90(1). U.S. Fish and Wildlife Service, Washington, D.C., USA. [online] URL: https://pdfs.semanticscholar.org/7d9b/959ae4b58353cd3018c9964adadc9791422c.pdf

Butcher, G. S., M. R. Fuller, L. S. McAllister, and P. H. Geissler. 1990b. An evaluation of the Christmas bird count for monitoring population trends of selected species. Wildlife Society Bulletin 18(2):129-134.

Butcher, G. S., D. K. Niven, and J. R. Sauer. 2005. Using Christmas bird count data to assess population dynamics and trends of waterbirds. American Birds 59:23-25.

Cornell Lab of Ornithology. 2016. ebird basic dataset. Version EBD_relAug-2016. Cornell Lab of Ornithology, Ithaca, New York, USA. [online] URL: https://ebird.org/science/download-ebird-data-products

Crewe, T. L., J. D. McCracken, P. D. Taylor, D. Lepage, and A. E. Heagy. 2008. The Canadian migration monitoring network - réseau canadien de surveillance des migrations: ten-year report on monitoring landbird population change. CMMN-RCSM Scientific Technical Report #1. Bird Studies Canada, Port Rowan, Ontario, Canada. [online] URL: https://www.bsc-eoc.org/download/CMMNReport2008.pdf

Crewe, T., P. D. Taylor, D. Lepage, A. C. Smith, and C. M. Francis. 2016. Quantifying regional variation in population trends using migration counts. Journal of Wildlife Management 80:245-255. https://doi.org/10.1002/jwmg.1003

Cribari-Neto, F., and A. Zeileis. 2010. Beta regression in R. Journal of Statistical Software 34(2):1-24. https://doi.org/10.18637/jss.v034.i02

Dunn, E. H., and J. R. Sauer. 1997. Monitoring Canadian bird populations with winter counts. Pages 49-55 in E. H. Dunn, M. D. Cadman, and J. B. Falls, editors. Monitoring bird populations: the Canadian experience. Occasional Paper 95. Canadian Wildlife Service, Ottawa, Ontario, Canada.

Dunn, E. H., K. A. Hobson, L. I. Wassenaar, D. J. T. Hussell, and M. L. Allen. 2006. Identification of summer origins of songbirds migrating through southern Canada in autumn. Avian Conservation and Ecology - Écologie et conservation des oiseaux 1(2):4. https://doi.org/10.5751/ACE-00048-010204

Environment and Climate Change Canada (ECCC). 2014a. Status of birds in Canada. Environment and Climate Change Canada, Ottawa, Ontario, Canada. [online] URL: https://wildlife-species.canada.ca/bird-status/index-eng.aspx?sY=2014&sL=e

Environment and Climate Change Canada (ECCC). 2014b. North American breeding bird survey - Canadian trends website. Data-version 2012. Environment and Climate Change Canada, Ottawa, Ontario, Canada.

Fogarty, F. A., M. E. Wohlfeil, and E. Fleishman. 2018. Response to Horns et al. 2018: ‘using opportunistic citizen science data to estimate avian population trends’. Biological Conservation 221, 151-159. Biological Conservation 226:329-330. https://doi.org/10.1016/j.biocon.2018.08.005

Francis, C. M., and D. J. T. Hussell. 1998. Changes in numbers of land birds counted in migration at Long Point Bird Observatory, 1961-1997. Bird Populations 4:37-66.

Harris, J. B. C., and D. G. Haskell. 2007. Land cover sampling biases associated with roadside bird surveys. Avian Conservation and Ecology - Écologie et conservation des oiseaux 2(2):12. http://dx.doi.org/10.5751/ace-00201-020212

Hoffman, M. D., and A. Gelman. 2014. The no-u-turn sampler: adaptively setting path lengths in Hamiltonian Monte Carlo. Journal of Machine Learning Research 15(1):1593-1623. [online] URL: https://dl.acm.org/doi/10.5555/2627435.2638586

Horns, J. J., F. R. Adler, and Ç. H. Şekercioğlu. 2018a. Using opportunistic citizen science data to estimate avian population trends. Biological Conservation 221:151-159. https://doi.org/10.1016/j.biocon.2018.02.027

Horns, J. J., F. R. Adler, and Ç. H. Şekercioğlu. 2018b. Rebuttal to response to Horns et al. 2018. Biological Conservation 226:331-332. https://doi.org/10.1016/j.biocon.2018.08.004

Isaac, N. J. B., A. J. van Strien, T. A. August, M. P. de Zeeuw, and D. B. Roy. 2014. Statistics for citizen science: extracting signals of change from noisy ecological data. Methods in Ecology and Evolution 5:1052-1060. https://doi.org/10.1111/2041-210X.12254

Kelling, S., A. Johnston, W. M. Hochachka, M. Iliff, D. Fink, J. Gerbracht, C. Lagoze, F. A. La Sorte, T. Moore, A. Wiggins, W.-K. Wong, C. Wood, and J. Yu. 2015. Can observation skills of citizen scientists be estimated using species accumulation curves? PLoS ONE 10(10):e0139600. https://doi.org/10.1371/journal.pone.0139600

Lloyd-Evans, T. L., and J. L. Atwood. 2004. 32 years of changes in passerine numbers during spring and fall migrations in coastal Massachusetts. Wilson Bulletin 116(1):1-16. https://doi.org/10.1676/0043-5643(2004)116[0001:yocipn]2.0.co;2

Natural Resources Canada (NRCAN). 2014. National forest pest strategy information system. Natural Resources Canada, Canadian Forest Service, Atlantic Forestry Centre, Fredericton, New Brunswick, Canada. [online] URL: https://afc-fr.cfsnet.nfis.org/NFPS-SNLRF/

North American Bird Conservation Initiative Canada (NABCI). 2016. The state of North America’s birds, 2016. Environment and Climate Change Canada, Ottawa, Ontario, Canada. [online] URL: https://www.stateofthebirds.org/2016/#_ga=2.184145958.1415599890.1599094380-928485324.1599094380

North American Bird Conservation Initiative Canada (NABCI). 2019. The state of Canada’s birds, 2019. Environment and Climate Change Canada, Ottawa, Ontario, Canada. [online] URL: http://nabci.net/wp-content/uploads/2019-State-of-Canadas-Birds-1.pdf

Niven, D. K., J. R. Sauer, G. S. Butcher, W. A. and Link. 2004. Christmas bird count provides insights into population change in land birds that breed in the boreal forest. American Birds 58:10-20.

Patten, M. A., and J. C. Burger. 1998. Spruce budworm outbreaks and the incidence of vagrancy in eastern North American wood-warblers. Canadian Journal of Zoology 76(3):433-439. https://doi.org/10.1139/z97-213

R Core Team. 2016. R: a language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. [online] URL https://www.r-project.org/

Roberts, R. L., P. F. Donald, and R. E. Green. 2007. Using simple species lists to monitor trends in animal populations: new methods and a comparison with independent data. Animal Conservation 10:332-339. https://doi.org/10.1111/j.1469-1795.2007.00117.x

Rosenberg, K. V., J. A. Kennedy, R. Dettmers, R. P. Ford, D. Reynolds, J. D. Alexander, C. J. Beardmore, P. J. Blancher, R. E. Bogart, G. S. Butcher, A. F. Camfield, A. Couturier, D. W. Demarest, W. E. Easton, R. P. Ford, J. J. Giocomo, R. H. Keller, A. E. Mini, A. O. Panjabi, D. N. Pashley, T. D. Rich, J. M. Ruth, H. Stabins, J. Stanton, and T. Will. 2016. Partners in Flight landbird conservation plan: 2016 revision for Canada and continental United States. Partners in Flight Science Committee, McAllen, Texas, USA. [online] URL: https://partnersinflight.org/wp-content/uploads/2016/08/pif-continental-plan-final-spread-single.pdf

Roy, H. E., T. Adriaens, N. J. B. Isaac, M. Kenis, T. Onkelinx, G. San Martin, P. M. J. Brown, L. Hautier, R. Poland, D. B. Roy, R. Comont, R. Eschen, R. Frost, R. Zindel, J. Van Vlaenderen, O. Nedved, H. P. Ravn, J.-C. Grégoire, J.-C. de Biseau, and D. Maes. 2012. Invasive alien predator causes rapid declines of native European ladybirds. Diversity and Distributions 18(7):717-725. https://doi.org/10.1111/j.1472-4642.2012.00883.x

Sauer, J. R., S. Schwartz, and B. Hoover. 1996. The Christmas Bird Count home page. Version 95.1. Patuxent Wildlife Research Center, Laurel, Maryland, USA. [online] URL: https://www.mbr-pwrc.usgs.gov/bbs/cbc.html

Sauer, J. S., D. K. Niven, and W. A. Link. 2004. Statistical analyses make the Christmas Bird Count relevant for conservation. American Birds 58:21-25.

Sauer, J. R., W. A. Link, J. E. Fallon, K. L. Pardieck, and D. J. Ziolkowski, Jr. 2013. The North American breeding bird survey 1966-2011: summary analysis and species accounts. North American Fauna 79:1-32. http://dx.doi.org/10.3996/nafa.79.0001

Sauer, J. R., D. K. Niven, J. E. Hines, D. J. Ziolkowski, Jr., K. L. Pardieck, J. E. Fallon, and W. A. Link. 2017. The North American breeding bird survey, results and analysis 1966-2015. Version 2.07.2017. USGS Patuxent Wildlife Research Center, Laurel, Maryland, USA. [online] URL: https://www.mbr-pwrc.usgs.gov/bbs/bbs2015.html

Schreiber, R. W., and E. A. Schreiber. 1973. Florida’s Brown Pelican population: Christmas bird count analysis. American Birds 27:711-715. [online] URL: https://sora.unm.edu/sites/default/files/journals/nab/v027n04/p00711-p00715.pdf

Soykan, C. U., J. Sauer, J. G. Schuetz G. S. LeBaron, K. Dale, and G. M. Langham. 2016. Population trends for North American winter birds based on hierarchical models. Ecosphere 7(5):e01351. https://doi.org/10.1002/ecs2.1351

Stan Development Team. 2017. Stan: a C++ library for probability and sampling. Version 2.14.0. Stan Development Team. [online] URL: http://mc-stan.org/

Sullivan, B. L., J. L. Aycrigg, J. H. Barry, R. E. Bonney, N. Bruns, C. B. Cooper, T. Damoulas, A. A. Dhondt, T. Dietterich, A. Farnsworth, D. Fink, J. W. Fitzpatrick, T. Fredericks, J. Gerbracht, C. Gomes, W. M. Hochachka, M. J. Iliff, C. Lagoze, F. A. La Sorte, M. Merrifield, W. Morris, T. B. Phillips, M. Reynolds, A. D. Rodewald, K. V. Rosenberg, N. M. Trautmann, A. Wiggins, D. W. Winkler, W.-K. Wong, C. L. Wood, J. Yu, and S. Kelling. 2014. The eBird enterprise: an integrated approach to development and application of citizen science. Biological Conservation 169:31-40. https://doi.org/10.1016/j.biocon.2013.11.003

Szabo, J. K., P. A. Vesk, P. W. J. Baxter, and H. P. Possingham. 2010. Regional avian species declines estimated from volunteer-collected long-term data using list length analysis. Ecological Applications 20:2157-2169. https://doi.org/10.1890/09-0877.1

U.S. Fish and Wildlife Service (USFWS). 2003. Monitoring plan for the American Peregrine Falcon, a species recovered under the Endangered Species Act. U.S. Fish and Wildlife Service, Divisions of Endangered Species and Migratory Birds and State Programs, Pacific Region, Portland, Oregon, USA. [online] URL: https://www.fws.gov/endangered/esa-library/pdf/Peregrineplan2003.pdf

U.S. Fish and Wildlife Service (USFWS). 2009. Post-delisting monitoring plan for the Bald Eagle (Haliaeetus leucocephalus) in the contiguous 48 States. U.S. Fish and Wildlife Service, Divisions of Endangered Species and Migratory Birds and State Programs, Midwest Regional Office, Twin Cities, Minnesota, USA. [online] URL: https://permanent.fdlp.gov/gpo58960/Postdelisting_Monitoring_Plan_for_the_Bald_Eagle_in_the_Contiguous_48_States-mar2009.pdf

U.S. Fish and Wildlife Service (USFWS). 2016. Waterfowl population status, 2016. U.S. Department of the Interior, Washington, D.C. USA. [online] URL: https://www.fws.gov/migratorybirds/pdf/surveys-and-data/Population-status/Waterfowl/WaterfowlPopulationStatusReport16.pdf

U.S. Fish and Wildlife Service (USFWS). 2017. Mid-winter waterfowl survey website. U.S. Department of the Interior, Washington, D.C. USA. [online] URL: https://migbirdapps.fws.gov/mbdc/databases/mwi/mwidb.asp

Venier, L., J. Pearce, D. Fillman, D. McNicol, and D. Welsh. 2009. Effects of spruce budworm (Choristoneura fumiferana (Clem.)) outbreaks on boreal mixed-wood bird communities. Avian Conservation and Ecology - Écologie et conservation des oiseaux 4(1):3. https://doi.org/10.5751/ace-00296-040103

Walker, J., and P. D. Taylor. 2017. Using eBird data to model population change of migratory birds. Avian Conservation and Ecology 12(1):4. https://doi.org/10.5751/ACE-00960-120104