R/spp_strategy.R

In SppTrend: Analyzing Linear Trends in Species Occurrence Data

#' @title Classify species ecological strategies
#' @description This function analyses the outputs of `spp_trend()` to classify species into distinct spatial or thermal response categories based on the direction and statistical
#' significance of their species-specific trends relative to the overall trend.
#' The function incorporates hemisphere-specific logic to correctly interpret poleward shifts in latitude and can also be applied to classify elevational trends.
#'
#' @param spp_trend_result A data frame containing trend indicators per species, typically generated by the `spp_trend` function. It should include columns such as:
#'   - `species`: Name of the analyzed species.
#'    - `responses`: Name of the analyzed variable.
#'    - `trend`: Estimated slope of the linear model.
#'    - `t`: t-statistic for the species-specific trend.
#'    - `pvalue`: Statistical significance of the species-specific trend.
#'    - `dif_t`: t-statistic of the interaction term, indicating the magnitude of the difference between the species trend and the Overall Trend (OT).
#'    - `dif_pvalue`: p-values of the interaction term. A low value indicates a significant deviation from the general trend.
#'    - `n`: Total number of occurrence records (sample size) for the specific species.
#'    - `hemisphere`: Geographical subset (`North`, `South`, or `Global`) used to ensure latitudinal symmetry in the analysis.
#'
#' @param sig_level The `numeric` significance level to use for classifying trends as significant. Defaults to `0.05`. See Bonferroni correction `0.05/length(species)`.
#' @param responses A `character vector` of response variable names to analyze (`c("lat", "lon", "tme", "ele")`).
#' The function will create classification columns for responses present in this vector and in the `responses` column of `spp_trend_result`.
#'
#' @return A `data frame` summarizing the ecological strategy of each species for each analyzed response variable.
#'   The table includes:
#'    - Species name
#'    - Hemisphere
#'    - Sample size
#'    - Classification columns for:
#'        - Spatial (latitude, longitude and elevation if present) responses. Spatial Adaptation `SA`, Spatial Discordance `SD`, Spatial Conformity `SC`
#'        - Thermal (temperature if present) responses. Thermal Tolerance `TT`, Thermal Adjustment `TA`, Thermal Conformity `TC`
#'
#'   Classification for spatial responses (`lat`, `lon`, `ele`) are classified as `Spatial_lat`, `Spatial_lon` and `Spatial_ele`.
#'   Thermal responses (`tme`) are classified as `Thermal_tme`.
#'
#' @details This function takes the trend analysis results from `spp_trend` and classifies each species' response based on the
#' significance of its trend and how it differs from the general trend.
#' Applied Bonferroni correction to avoid false positives (Type I errors) due to multiple comparisons when analyzing many species.
#' The classification identifying three possible spatial responses and three thermal responses:
#'   \itemize{
#'     \item \strong{Spatial Responses:}
#'       \itemize{
#'         \item \strong{Spatial Adaptation (SA):} A significant positive temporal trend in the spatial position of species occurrences.
#'         In the context of climate change, this pattern is commonly associated with a poleward shift, corresponding to a northward displacement
#'         (towards higher latitude values) in the Northern Hemisphere and southward displacement (towards lower latitude values) in the Southern Hemisphere, as species expand into newly suitable areas.
#'         \item \strong{Spatial Discordance (SD):} A significant negative temporal trend in the spatial position of species occurrences.
#'         In the context of climate change, this pattern is often associated with an equatorward shift and may arise when other ecological and anthropogenic
#'         factors influence species distributions independently of, or in opposition to, climate-driven range shifts.
#'         \item \strong{Spatial Conformity (SC):} A spatial response pattern in which the species-specific temporal trend does not differ significantly from the overall trend.
#'         Species showing spatial conformance share the same bias structure as the complete dataset, preventing the inference of a distinct, species-specific response to climate change at the scale of analysis.
#'       }
#'     \item \strong{Thermal Responses:}
#'       \itemize{
#'         \item \strong{Thermal Tolerance (TT):} A thermal response pattern characterised by a significant positive temporal trend in the temperature conditions under which species are observed, relative to the overall trend.
#'         This pattern suggest an increased likelihood of occurrence under warmer conditions and an apparent capacity to tolerate rising temperatures through physiological, behavioural, and evolutionary mechanisms.
#'         \item \strong{Thermal Adjustment (TA):} A thermal response characterised by a significant negative temporal trend in the temperature conditions associated with species occurrences, relative to the overall trend.
#'         This indicates and increasing association with cooler temperature conditions over time, potentially reflecting microevolutionary change or phenotypic adjustment.
#'         \item \strong{Thermal Conformity (TC):} A thermal response pattern in which species-specific temperature trends do not differ significantly from the overall trend.
#'         Species showing thermal conformance share the same background thermal signal as the complete dataset, preventing the formulation of specific hypotheses regarding climate-driven thermal responses.
#'       }
#'   }
#'
#'   *Note: The interpretation of longitude trends assumes that if transformation was applied in `spp_trend`, it used the Antimeridian as 0.*
#'
#' @importFrom dplyr select mutate case_when all_of lead group_by summarise if_else ungroup relocate
#' @importFrom tidyr pivot_wider
#' @importFrom tidyselect starts_with
#' @importFrom utils head
#'
#' @examples
#'
#' # Assuming spp_trends_results is a data frame generated by spp_trend()
#'
#' spp_trends_results <- data.frame(
#'   species = paste0("spp_", 1:10),
#'   responses = rep(c("lat", "lon", "tme"), length.out = 30),
#'   trend = runif(30, -0.5, 0.5),
#'   t = runif(30, -2, 2),
#'   pvalue = runif(30, 0, 1),
#'   dif_t = runif(30, -1, 1.5),
#'   dif_pvalue = runif(30, 0.001, 0.9),
#'   n = round(runif(30, 40, 60)),
#'   hemisphere = sample(c("North", "South", "Global"), 30, replace = TRUE)
#' )
#'
#' spp <- unique(spp_trends_results$species)
#' sig_level <- 0.05 / length(spp) # Bonferroni correction
#' responses_to_analyze <- c("lat", "lon", "tme")
#'
#' spp_strategy_results <- spp_strategy(spp_trends_results,
#'                                      sig_level = sig_level,
#'                                      responses = responses_to_analyze)
#'
#' print(spp_strategy_results)
#'
#' @export
spp_strategy <- function(spp_trend_result,
                         sig_level = 0.05,
                         responses = responses) {
  classify_spatial_standard <- function(pvalue, dif_pvalue, trend) {
    dplyr::case_when(
      pvalue > sig_level ~ "SC",
      pvalue <= sig_level &
        dif_pvalue <= sig_level & trend > 0 ~ "SA",
      pvalue <= sig_level &
        dif_pvalue <= sig_level & trend < 0 ~ "SD",
      TRUE ~ "SC"
    )
  }
  classify_lat_poleward <- function(pvalue, dif_pvalue, trend, hemisphere) {
    dplyr::case_when(
      pvalue > sig_level ~ "SC",
      pvalue <= sig_level &
        dif_pvalue <= sig_level &
        trend > 0 & hemisphere == "North" ~ "SA", #SP
      pvalue <= sig_level &
        dif_pvalue <= sig_level &
        trend < 0 & hemisphere == "North" ~ "SD", #SE
      pvalue <= sig_level &
        dif_pvalue <= sig_level &
        trend > 0 & hemisphere == "South" ~ "SD", #SE
      pvalue <= sig_level &
        dif_pvalue <= sig_level &
        trend < 0 & hemisphere == "South" ~ "SA", #SP
      pvalue <= sig_level &
        dif_pvalue <= sig_level &
        trend > 0 & hemisphere == "Global" ~ "SA", #SA
      pvalue <= sig_level &
        dif_pvalue <= sig_level &
        trend < 0 & hemisphere == "Global" ~ "SD", #SD
      TRUE ~ "SC"
    )
  }
  classify_thermal <- function(pvalue, dif_pvalue, trend) {
    dplyr::case_when(
      pvalue > sig_level ~ "TC",
      pvalue <= sig_level &
        dif_pvalue <= sig_level & trend < 0 ~ "TA",
      pvalue <= sig_level &
        dif_pvalue <= sig_level & trend > 0 ~ "TT",
      TRUE ~ "TC"
    )
  }

  required_cols <- c(
    "species",
    "trend",
    "t",
    "pvalue",
    "responses",
    "dif_t",
    "dif_pvalue",
    "n",
    "hemisphere"
  )

  if (!all(required_cols %in% names(spp_trend_result))) {
    missing_cols <- setdiff(required_cols, names(spp_trend_result))
    stop(
      paste(
        "Error: The following columns were not found in 'spp_trend_result':",
        paste(missing_cols, collapse = ", "),
        ". The required columns are:",
        paste(required_cols, collapse = ", ")
      )
    )
  }
  strategies <- spp_trend_result %>%
    dplyr::select(dplyr::all_of(intersect(
      required_cols, names(spp_trend_result)
    ))) %>%
    dplyr::mutate(
      Spatial_lat_Poleward = dplyr::case_when(
        .data$responses == "lat" ~ classify_lat_poleward(
          .data$pvalue,
          .data$dif_pvalue,
          .data$trend,
          .data$hemisphere
        ),
        TRUE ~ NA_character_
      ),
      Spatial_lon = dplyr::case_when(
        .data$responses == "lon" ~ classify_spatial_standard(.data$pvalue, .data$dif_pvalue, .data$trend),
        TRUE ~ NA_character_
      ),
      Spatial_ele = dplyr::case_when(
        .data$responses == "ele" ~ classify_spatial_standard(.data$pvalue, .data$dif_pvalue, .data$trend),
        TRUE ~ NA_character_
      ),
      Thermal_tme = dplyr::case_when(
        .data$responses == "tme" ~ classify_thermal(.data$pvalue, .data$dif_pvalue, .data$trend),
        TRUE ~ NA_character_
      )
    ) %>%
    dplyr::group_by(species, hemisphere) %>%
    dplyr::summarise(
      n = sum(n),
      Spatial_lat_Poleward = dplyr::first(Spatial_lat_Poleward),
      Spatial_lon = if("lon" %in% unique(responses)){
        dplyr::first(Spatial_lon[!is.na(Spatial_lon)])
      }else{
        NA_character_
        },
      Spatial_ele = if ("ele" %in% unique(responses)){
        dplyr::first(Spatial_ele[!is.na(Spatial_ele)])
      }else{
        NA_character_
        },
      Thermal_tme = if ("tme" %in% unique(responses)){
        dplyr::first(Thermal_tme[!is.na(Thermal_tme)])
      }else{
        NA_character_
        },
      .groups = "drop"
    ) %>%
    tidyr::pivot_wider(
      names_from = hemisphere,
      values_from = Spatial_lat_Poleward,
      names_prefix = "Spatial_lat_"
    ) %>%
    dplyr::ungroup() %>%
    dplyr::relocate(tidyselect::starts_with("Spatial_lat_"), .after = Spatial_lon)
  return(strategies)
}