R/survival_table.R

In SNARL1/mrmr: Mark recapture miscellany in R

#' Generate a table of introduced survival over time
#'
#' This function generates a table of survival estimates for introduced
#' individuals after introduction. These estimates average over individuals,
#' and over primary periods within years, so that the values are comparable
#' among sites with one or multiple introduced cohorts, and
#' a variable number of primary periods per year.
#'
#' @param model A model objectgenerated by \code{mrmr:fit_model()}.
#' @param by_cohort (bool) Whether to summarize survival separately by cohort.
#' @param by_individual (bool) Whether to summarize survival separately by frog.
#' @param thin (int) optional: thinning period (e.g., 2 takes every other draw)
#' to reduce memory usage for very large datasets
#' @return A data.frame with columns year_since_introduction,
#' lo_survival, med_survival, and hi_survival,
#' where lo, med, and hi represent the 2.5%, 50%, and 97.5% posterior quantiles
#' for survival in the years following introduction.
#' @examples
#' \dontrun{
#' captures <- system.file("extdata", "capture-example.csv", package = "mrmr")
#' translocations <- system.file("extdata", "translocation-example.csv",
#'                               package = "mrmr")
#' surveys <- system.file("extdata", "survey-example.csv", package = "mrmr")
#' out <- clean_data(captures, translocations, surveys)
#' model <- fit_model(out, chains = 1, iter = 10)
#' survival_table(model)
#' }
#' @importFrom dplyr group_by summarize left_join ungroup filter distinct
#' @importFrom reshape2 melt
#' @importFrom rlang .data
#' @importFrom tibble as_tibble
#' @importFrom stats median quantile
#' @export

survival_table <- function(
  model, by_cohort = TRUE, by_individual = FALSE, thin = 1
) {

  any_translocations <- 'data.frame' %in% class(model$data$translocations)

  if (!any_translocations) {
    stop(paste("No translocation data are present, so a cohort survival",
               "table cannot be created."))
  }

   primary_period_dates <- model$data$surveys %>%
    group_by(.data$primary_period) %>%
    summarize(date = min(.data$survey_date),
              year = min(.data$year)) %>%
    ungroup

  survival_summary <- model$m_fit$draws("s", format = "draws_df") %>%
    posterior::thin_draws(thin) %>%
    tidyr::pivot_longer(tidyselect::starts_with("s")) %>%
    suppressWarnings() %>%
    mutate(
      get_numeric_indices(.data$name),
      primary_period = .data$index_2,
      pit_tag_id = dimnames(model$data$stan_d$Y)[[1]][.data$index_1]
    ) %>%
    filter(.data$pit_tag_id %in% as.character(model$data$translocations$pit_tag_id)) %>%
    left_join(distinct(model$data$translocations,
                       .data$pit_tag_id, .data$release_date)) %>%
    left_join(primary_period_dates) %>%
    filter(.data$date > .data$release_date) %>%
    dplyr::transmute(
      .data$.draw, .data$value, .data$primary_period, .data$pit_tag_id,
      .data$release_date, .data$date, .data$year,
      years_since_introduction = lubridate::year(.data$date) -
        lubridate::year(.data$release_date)
    )


  survival_summary <- survival_summary %>%
    group_by(.data$years_since_introduction, .data$.draw)

  if (by_cohort) {
    survival_summary <- survival_summary %>%
      group_by(.data$release_date, .add = TRUE)
  }

  if (by_individual) {
    survival_summary <- survival_summary %>%
      group_by(.data$pit_tag_id, .add = TRUE)
  }

  # for each group, compute fraction alive
  group_summaries <- survival_summary %>%
    summarize(fraction_alive = mean(.data$value == 2), .groups = "drop") %>%
    group_by(.data$years_since_introduction, .data$release_date)

  if (by_individual) {
    group_summaries <- group_summaries %>%
      group_by(.data$pit_tag_id, .add = TRUE)
  }

  group_summaries  %>%
    summarize(
      lo_survival = quantile(.data$fraction_alive, .025),
      median_survival = median(.data$fraction_alive),
      hi_survival = quantile(.data$fraction_alive, .975),
      .groups = "drop"
      )
}