R/generate-cont-trend.R

In BrandonEdwards/bbsBayes: Hierarchical Bayesian Analysis of North American BBS Data

#' Generate the continental trend
#'
#' \code{generate_cont_trend} calculates the geometric mean annual changes in population size.
#'
#' @param indices Continental indices generated by \code{generate_cont_indices}
#' @param min_year Minimum year to calculate trends from
#' @param max_year Maximum year to calculate trends to
#' @param quantiles vector of quantiles to be sampled from the posterior distribution Defaults to c(0.025,0.05,0.25,0.5,0.75,0.95,0.975)
#' @param slope Logical, if TRUE, calculates an alternative trend metric, the slope of a log-linear regression through the annual indices. Default FALSE
#'
#' @return Numeric percentage of trend
#'
#' @importFrom stats lm
#'
#' @name bbsBayes-deprecated
#' @seealso \code{\link{generate_cont_trend}}
#' @keywords internal
NULL

#' @rdname bbsBayes-deprecated
#' @section \code{generate_cont_trend}:
#'   For \code{generate_cont_trend()}, use
#'   \code{generate_trends()}.
#'
#' @export
#'

generate_cont_trend <- function(indices = NULL,
                                min_year = NULL,
                                max_year = NULL,
                                quantiles = c(0.025,0.05,0.25,0.75,0.95,0.975),
                                slope = FALSE)
{

  .Deprecated(new = "generate_trends",
              msg = paste("generate_cont_trends is deprecated in favour of generate_trends()"))


  if (is.null(indices))
  {
    stop("No indices supplied to generate_cont_trend()."); return(NULL)
  }
  n = indices$samples

  if (is.null(min_year))
  {
    min_year = indices$y_min
  }
  if (is.null(max_year))
  {
    max_year = indices$y_max
  }

if(slope){
  wy = c(min_year:max_year)

  bsl = function(i){
    n = length(wy)
    sy = sum(i)
    sx = sum(wy)
    ssx = sum(wy^2)
    sxy = sum(i*wy)
    b = (n*sxy - sx*sy)/(n*ssx - sx^2)
    return(b)
  }
 ne = log(n[,wy])
  m =  t(apply(ne,1,FUN = bsl))

  sl.t = as.vector((exp(m)-1)*100)
}

  ch = n[,max_year]/n[,min_year]
  tr = 100*((ch^(1/(max_year-min_year)))-1)

  trend <- data.frame(Start_year = (indices$startyear+min_year)-1,
                              End_year = (indices$startyear+max_year)-1,
                      Region = "Continental",
                      Trend = median(tr),
                      stringsAsFactors = FALSE)
  for(qq in quantiles){
  trend[,paste0("Trend_Q",qq)] <- quantile(tr,qq,names = FALSE)
  }
  trend[,"Percent_Change"] <- median(ch)
  for(qq in quantiles){
    trend[,paste0("Percent_Change_Q",qq)] <- 100*(quantile(ch,qq,names = FALSE)-1)
  }

  if(slope){
  trend[,"Slope_Trend"] <- median(sl.t)
  for(qq in quantiles){
    trend[,paste0("Slope_Trend_Q",qq)] <- quantile(sl.t,qq,names = FALSE)
  }
  }

  return(trend)
}