**RPI Raptor Population Trends: Interpreting Trend Graphs and Maps**

Our goal was to calculate annual indices and
population trends for all stations that had at least 10 years of
hourly data (through 2019) entered into HawkCount.org. Sites with
less than 10 years were considered to have insufficient data to
estimate a meaningful long-term trend. Species with fewer than 10
individuals per year at a site were also considered to have
insufficient data for trend analysis, and were thus excluded.

It
is important for users to recognize that the RPI partners do not have
the resources or staff to do exhaustive checking of raw data and
subsequent results. We rely on the sites contributing data to alert
us to problems and to provide accurate data in HawkCount. Any
questions on datasets or results can be sent to brown@hmana.org.
Site descriptions and protocols are posted within the site profiles
on hawkcount.org.

Results are presented in two formats:
trend graphs containing the long-term linear trend for each species
and site (see trend graphs), and trend maps which graphically show
the most recent 10-year trends across sites for a given species, in
this case for the period 2009 to 2019 (see trend maps).

**How to read the population graphs**

The following example shows the long-term trend for Northern Harrier at Hawk Mountain Sanctuary in fall. The black dots on the chart are estimated indices of annual population size as calculated by the statistical model, with error bars representing the 95% credible intervals for the annual indices. Annual indices are expressed on a linear scale and represent the average of the predicted number of individual birds detected per day of survey. If present, the blue line is a smoothed (loess) visualization of the change in annual indices over time.

The subtitle of the plot contains the estimated linear trend value for the long-term trend and for a trend based on the most recent 10-years (time periods in parentheses). These trends are expressed as annual percent change in population index size, and assume a linear change in population size among years over the time period examined. The lower and upper limits of the 95% credible interval (CI) of the trends are shown in parentheses. In the example shown, the CI do not include zero for 2009-2019, suggesting strong support for the observed decline in migration counts for this species over the past 10 years. It can be seen from the plot that this decline has extended back several decades (see species assessment). An * following the reported trend and CI provides the posterior probability of the trend: a value greater than or equal to 0.95 provides additional and strong support for the observed trend.

**Trend maps**

Trend maps display the trends calculated from all RPI stations to more easily visualize the spatial variation in trends of a species over a large geographic area. Consistent trends across sites within a region gives us confidence that the results we see are representative of trends in the underlying regional population, and not the result of site-specific biases. To ensure that all trends are directly comparable, maps show only the trend based on the most recent 10 years for each station. The example below represents the population trends of Bald Eagle based on fall count data. Upward arrows suggest a positive trend was observed at a site, and downward arrows suggest a negative trend was observed at a site. Green and red arrows suggest the 10-year trend is strongly supported (95% credible intervals exclude zero), while blue dots represent trends with little or no support, i.e., 95% credible intervals included zero, and the posterior probability was < 0.95. Note that the size of the green and red arrows is proportional to the magnitude of the increase or decline, as explained in the map legend. You can click on any of the markers to see more details about the station name, the trend, and the time period covered. You can also zoom in or out to better see the details in areas where there are several nearby stations, such as along the Appalachian Mountains. The example clearly shows that Bald Eagles are doing very well in the east but may have stabilized in the west where counts are generally much lower.

**Analysis methods**

*Data filtering*

Prior to analysis, we filtered data to exclude days of the year and hours of the day that were not typically sampled at a station by including only those days of the year and hours of the day that included the inner 95% of observations at a station. We also excluded species that were not detected during at least 50% of all years surveyed at a station. Hourly data were aggregated to daily totals, following the recommendations in Crewe et al. 2016 (accepted for publication).

*Analysis overview*

We used hierarchical linear regression in a Bayesian framework (INLA, R-project) to estimate the linear change (%/year) in the number of individuals migrating past each count site each day over the past 10-, 20- and 30-year time periods. A trend for the full time period at each site was also estimated. The regressions for each time period assumed either a negative binomial or Poisson distribution of counts, and estimated first- and second-order effects for day of year in addition to the continuous year (trend) effect. We also assumed a random effect for year which accounted for temporal auto-correlation of counts among years, and a random effect for day of year, for which we assumed independent and identically distributed errors. Annual indices of population size (mean number of hawks per hour ± standard deviation) were estimated from the posterior distribution of each trend model.

*Regression analysis*

Trends in raptor
migration counts were estimated independently for each species, site
and season in a Bayesian framework using Integrated Nested Laplace
Approximation (R-INLA, Rue et al. 2014) in R (version 3.1.3; R Core
Team 2014; see also Ross et al. 2012). We estimated trends using
log-linear regression, which included a continuous effect for year
(*i*) to estimate log-linear change in population size over
time, first and second order effects for day of year (*j*) to
model the seasonal distribution of counts, and with hierarchical
terms to account for random variation in counts among years and among
days:

,

where
is a first-order autoregressive (AR1) random effect for year, and
is
a hierarchical term to account for random variation among days of the
year. While an AR1 random effect for day of year nested within year
might be more appropriate to properly partition variation in daily
counts among fixed and random effects, we found using simulated
migration counts that specifying the random day effect as independent
and identically distributed (IID) had no noticeable effect on trend
bias, or probability of estimating a precise trend (probability that
the simulated trend fell within the credible intervals of the
estimated trend; T. L. Crewe, unpublished data). Specifying the
random day effect as IID did, however, significantly increase the
speed at which sites could be analyzed, and given the volume of data
to analyze, an IID random effect was specified. Number of
observation hours in day *j* was included as an offset to
account for variable daily effort.

We first fit the model by assuming a negative binomial data distribution. If the estimate for the random year effect was large and spurious (>100), we instead re-fit the same model but assumed a Poisson distribution of counts. For both data distributions, year estimates and 95% credible intervals were back-transformed to annual rates of population change using 100*exp(estimate)-1. Trends were calculated using the full dataset, as well as for all 10-year subsets to estimate 10-, 20-, 30-year (etc.) trends where appropriate, such that if the final year in the dataset was 2014, a 10-year trend would include 2004-2014.

We estimated annual indices of population size for each species on the response scale (mean hawks per day), by calculating the mean and standard deviation of predicted daily count from 1000 samples from the posterior distribution of the above trend model. The trend was calculated as the % rate of change per year = 100*(exp(year coefficient)-1). This converts the year coefficient, which is an estimate of the instantaneous rate of change per year, into an estimate of the discrete rate of change per year between any year and the following year.

Crewe, T. L., P. D. Taylor, and D. Lepage. 2016. Temporal aggregation of migration counts can improve accuracy and precision of trends. Avian Conservation and Ecology (accepted for publication)11(2):8. http://dx.doi.org/10.5751/ACE-00907-110208

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

Ross, B. E., M. B. Hooten, and D. N. Koons. 2012. An accessible method for implementing hierarchical models with spatio-temporal abundance data. PLoS ONE 7:e49395.

Rue, H., S. Martino, F. Lindgren, D. Simpson, and A. Riebler. 2014. INLA: Functions which allow to perform full Bayesian analysis of latent Gaussian models using Integrated Nested Laplace Approximation. [online] URL: http://www.r-inla.org.

**Availability of data**

All graphs and maps are available from BSC's Birds
Canada’s NatureCounts web site: __http://www.naturecounts.ca/rpi/__,
as well as from __www.rpi-project.org__.
Count stations can easily copy a small line of script (available at
NatureCounts) and embed the population trends into their own web site
if they wish. To do so, simply define the population trends that you
would like to see on NatureCounts website and look for the link that
says *To insert the graphs below in your
own web site using our web service click here.*

Annual indices and trends can also be downloaded from NatureCounts in an Excel compatible format, if you wish to create your own graphs for presentation as well. Simply look for the Download button of the population trend tool on NatureCounts.