Skip to contents

Calculate annual indices of relative abundance by year for different regions. These indices can then be used to plot population trajectories for the species, and to estimate trends.

Usage

generate_indices(
  model_output = NULL,
  quantiles = c(0.025, 0.05, 0.25, 0.75, 0.95, 0.975),
  regions = c("continent", "stratum"),
  regions_index = NULL,
  alternate_n = "n",
  start_year = NULL,
  max_backcast = NULL,
  drop_exclude = FALSE,
  hpdi = FALSE,
  quiet = FALSE
)

Arguments

model_output

List. Model output generated by run_model().

quantiles

Numeric. Vector of quantiles to be sampled from the posterior distribution. Default is c(0.025, 0.05, 0.25, 0.5, 0.75, 0.95, 0.975). Note that these quantiles will be used to create confidence interval bands in plot_indices() and by quantiles in generate_trends(), so make sure you specify the ones you want to use later.

regions

Character. Which region(s) to summarize and calculate indices for. Default is "continent" and "stratum". Options also include "country", "prov_state", "bcr", and "bcr_by_country". Note that some regions only apply to specific stratifications. You can also supply a custom region that exists as a column in the regions_index data frame (see examples for more details).

regions_index

Data frame. Custom regions to summarize. Data frame must include all strata in the original data in one column (strata_name), and any custom regions defined as categories in other columns. See examples for more details.

alternate_n

Character. Indicating the name of the alternative annual index parameter in a model, Default is "n" which for all models represents an index of the estimated annual relative abundance, scaled as the expected mean count averaged across all BBS routes and observers in a given stratum and year. For some of the models included in bbsBayes2, alternatives exist that provide a partial decomposition of the time-series. For the "gamye" models, the parameter "n_smooth" represents the smooth-only version of the annual index of relative abundance (i.e., the component of the annual index estimated by the spline-based smooth of the GAM). This "n_smooth" is identical to the "n" values for the same model, but excludes the annual fluctuations. For the "gamye" models, this "n_smooth" parameter is likely the most natural parameter to use in estimating trends. A similar option exists for the "slope" models, where the parameter "n_slope" represents the component of the population trajectory estimated by the log-linear regression slope parameters in the model. Users should be particularly cautious about interpreting this "n_slope" values for relatively long time-series. As a continuous regression slope, it assumes interpreting it as an estimate of population trajectory and using it to generate trend estimates assumes that there is a single continuous rate of population change across the entire time-series. Biologically, this may be reasonable for 10-20 year periods, but will be less reasonable for longer time-periods.

start_year

Numeric. Trim the data record before calculating annual indices.

max_backcast

Numeric. The number of years to back cast stratum-level estimates before the first year that species was observed on any route in that stratum. Default is NULL, which generates annual indices for the entire time series and ignores back-casting. CWS national estimates use a back cast of 5. Note that unless drop_exclude = TRUE, problematic years are only flagged, not omitted. See Details for more specifics.

drop_exclude

Logical. Whether or not strata that exceed the max_backcast threshold should be excluded from the calculations. Default is FALSE (regions are flagged and listed but not dropped).

hpdi

Logical. Should credible intervals and limits be calculated using highest posterior density intervals instead of simple quantiles of the posterior distribution. Default is FALSE. these intervals are often a better descriptor of skewed posterior distributions, such as the predicted mean counts that the indices represent. Note hpd intervals are not stable for small percentages of the posterior distribution, and so hdpi = TRUE is ignored for quantiles values between 0.33 and 0.67 (i.e., if the quantiles value defines a limit for a centered hpd interval that would include < 33% of the posterior distribution).

quiet

Logical. Suppress progress messages? Default FALSE.

Value

A list containing

  • indices - data frame of calculated regional annual indices of abundances

  • samples - array of posterior draws from the model

  • meta_data - meta data defining the analysis

  • meta_strata - data frame listing strata meta data

  • raw_data - data frame of summarized counts

indices contains the following columns:

  • year - Year of particular index

  • region - Region name

  • region_type - Type of region

  • strata_included - Strata potentially included in the annual index calculations

  • strata_excluded - Strata potentially excluded from the annual index calculations because they have no observations of the species in the first part of the time series, see arguments max_backcast and start_year

  • index - Strata-weighted count index (median)

  • index_q_XXX - Strata-weighted count index (by different quantiles)

  • obs_mean - Mean observed annual counts of birds across all routes and all years. An alternative estimate of the average relative abundance of the species in the region and year. Differences between this and the annual indices are a function of the model. For composite regions (i.e., anything other than stratum-level estimates) this average count is calculated as an area-weighted average across all strata included

  • n_routes - Number of BBS routes that contributed data for this species, region, and year

  • n_routes_total - Number of BBS routes that contributed data for this species and region for all years in the selected time-series, i.e., all years since start_year

  • n_non_zero - Number of BBS routes on which this species was observed (i.e., count is > 0) in this region and year

  • backcast_flag - Approximate annual average proportion of the covered species range that is free of extrapolated population trajectories. e.g., if 1.0, data cover full time-series; if 0.75, data cover 75 percent of time-series. Only calculated if max_backcast != NULL.

Details

max_backcast is a way to deal with the fact that the species of interest may not appear in the data until several years after the start of the time-series max_backcast specifies how many years can occur before the stratum is flagged. A max_backcast of 5 will flag any stratum without a non-zero (or non-NA) observation within the first 5 years of the time- series. Note that records are only flagged unless drop_exclude = TRUE. If you find that the early data record is sparse and results in the exclusion of many strata, consider trimming the early years by specifying a start_year.

See also

Other indices and trends functions: generate_trends(), plot_geofacet(), plot_indices(), plot_map()

Examples

# Using the example model for Pacific Wrens

# Generate the continental and stratum indices
i <- generate_indices(pacific_wren_model)
#> Processing region continent
#> Processing region stratum

# Generate the continental and stratum indices using hpdi
i <- generate_indices(pacific_wren_model, hpdi = TRUE)
#> Processing region continent
#> Processing region stratum

# Generate only country indices
i_nat <- generate_indices(pacific_wren_model, regions = "country")
#> Processing region country

# Use a custom region specification (dummy example)
library(dplyr)
#> 
#> Attaching package: ‘dplyr’
#> The following objects are masked from ‘package:stats’:
#> 
#>     filter, lag
#> The following objects are masked from ‘package:base’:
#> 
#>     intersect, setdiff, setequal, union
ri <- bbs_strata[["bbs_cws"]]
ri <- mutate(ri, my_region = if_else(prov_state %in% "ON",
                                     "Ontario", "Rest"))

# Generate indices with these custom regions
i_custom <- generate_indices(
  pacific_wren_model,
  regions = c("country", "prov_state", "my_region"),
  regions_index = ri)
#> Processing region country
#> Processing region prov_state
#> Processing region my_region