R/utils.R

In insurancerating: Actuarial Tools for Insurance Pricing Models

#' Set the reference level of a factor
#'
#' @description
#' Relevels a factor so that the selected category becomes the reference
#' (first) level. By default, the reference level is chosen as the level with
#' the largest total weight, for example the largest exposure in an insurance
#' portfolio. Use `method = "manual"` with `reference_level` when a specific
#' business category should be the reference level.
#'
#' Choosing a reference level does not change fitted values or the overall model
#' fit. It changes the coefficient parameterisation and therefore the level
#' against which the remaining factor relativities are expressed.
#'
#' @param x A factor (unordered). Character vectors should be converted to
#'   factor before use.
#' @param weight A numeric vector of the same length as \code{x},
#'   typically representing exposure or frequency weights. Required when
#'   `method = "largest_weight"`.
#' @param method Character. Method used to choose the reference level.
#'   Supported methods are `"largest_weight"` and `"manual"`.
#' @param reference_level Character string with the level to use as reference
#'   when `method = "manual"`.
#'
#' @details
#' `method = "largest_weight"` is useful when the reference category should
#' represent a substantial and relatively stable part of the portfolio. The
#' supplied `weight` is commonly earned exposure, but another actuarially
#' meaningful volume measure may be used.
#'
#' `method = "manual"` is appropriate when the reference category is determined
#' by tariff interpretation, governance or an established pricing convention.
#' The selected category must already be an observed factor level.
#'
#' @author Martin Haringa
#'
#' @references Kaas, Rob & Goovaerts, Marc & Dhaene, Jan & Denuit, Michel.
#' (2008). Modern Actuarial Risk Theory: Using R.
#' \doi{10.1007/978-3-540-70998-5}
#'
#' @importFrom stats relevel
#'
#' @return A factor of the same length as \code{x}, with the selected reference
#'   level set as the first level.
#'
#' @seealso [add_rebasing()] for rescaling current tariff relativities within a
#'   prepared refinement workflow after a model has been fitted.
#'
#' @examples
#' portfolio <- data.frame(
#'   region = factor(c("North", "North", "South", "West")),
#'   exposure = c(120, 80, 60, 40)
#' )
#'
#' set_reference_level(portfolio$region, portfolio$exposure)
#' set_reference_level(
#'   portfolio$region,
#'   method = "manual",
#'   reference_level = "South"
#' )
#'
#' # Apply the largest-weight reference rule to every factor in a data frame
#' library(dplyr)
#' df <- chickwts |>
#'   mutate(across(where(is.character), as.factor)) |>
#'   mutate(across(where(is.factor), ~set_reference_level(., weight)))
#'
#' @export
set_reference_level <- function(x,
                                weight = NULL,
                                method = "largest_weight",
                                reference_level = NULL) {
  if (!is.factor(x)) stop("`x` must be a factor")
  if (!is.character(method) || length(method) != 1L || is.na(method)) {
    stop("`method` must be 'largest_weight' or 'manual'.", call. = FALSE)
  }
  if (!method %in% c("largest_weight", "manual")) {
    stop("`method` must be 'largest_weight' or 'manual'.", call. = FALSE)
  }

  if (method == "manual") {
    if (!is.character(reference_level) ||
        length(reference_level) != 1L ||
        is.na(reference_level)) {
      stop("`reference_level` must be a single character string.",
           call. = FALSE)
    }
    if (!reference_level %in% levels(x)) {
      stop("`reference_level` must be one of the levels of `x`.",
           call. = FALSE)
    }

    xrelevel <- stats::relevel(x, ref = reference_level)
    attr(xrelevel, "xoriginal") <- levels(x)
    return(xrelevel)
  }

  if (!is.null(reference_level)) {
    stop("`reference_level` can only be used when `method = 'manual'`.",
         call. = FALSE)
  }
  if (!is.numeric(weight)) stop("`weight` must be numeric")
  if (length(x) != length(weight)) {
    stop("`x` and `weight` must have the same length")
  }

  counts <- sort(tapply(weight, x, FUN = sum), decreasing = TRUE)
  xrelevel <- stats::relevel(x, ref = names(counts)[1])
  attr(xrelevel, "xoriginal") <- levels(x)
  xrelevel
}


# Internal formatter shared by discrete-axis plotting methods.
format_discrete_axis_labels <- function(labels,
                                        abbreviate_labels = TRUE,
                                        label_width = 20L,
                                        label_abbreviations = NULL) {
  if (!is.logical(abbreviate_labels) ||
      length(abbreviate_labels) != 1L ||
      is.na(abbreviate_labels)) {
    stop("`abbreviate_labels` must be TRUE or FALSE.", call. = FALSE)
  }
  if (!is.numeric(label_width) ||
      length(label_width) != 1L ||
      is.na(label_width) ||
      !is.finite(label_width) ||
      label_width < 2 ||
      label_width != floor(label_width)) {
    stop("`label_width` must be a whole number of at least 2.", call. = FALSE)
  }
  if (!is.null(label_abbreviations)) {
    if (!is.character(label_abbreviations) ||
        is.null(names(label_abbreviations)) ||
        anyNA(label_abbreviations) ||
        anyNA(names(label_abbreviations)) ||
        any(names(label_abbreviations) == "")) {
      stop(
        "`label_abbreviations` must be a named character vector.",
        call. = FALSE
      )
    }
  }

  original <- as.character(labels)
  displayed <- original
  mapped <- rep(FALSE, length(original))

  if (!is.null(label_abbreviations)) {
    matched <- match(original, names(label_abbreviations))
    mapped <- !is.na(matched)
    displayed[mapped] <- unname(label_abbreviations[matched[mapped]])
  }

  if (isTRUE(abbreviate_labels)) {
    shorten <- !mapped & !is.na(displayed) & nchar(displayed) > label_width
    displayed[shorten] <- paste0(
      substr(displayed[shorten], 1L, label_width - 1L),
      "."
    )
  }

  duplicated_labels <- duplicated(displayed) | duplicated(displayed, fromLast = TRUE)
  duplicated_labels <- duplicated_labels & !is.na(displayed)
  if (any(duplicated_labels)) {
    duplicate_groups <- split(
      which(duplicated_labels),
      displayed[duplicated_labels]
    )
    for (indices in duplicate_groups) {
      if (length(unique(original[indices])) < 2L) {
        next
      }
      sequence_number <- seq_along(indices)
      suffix_length <- nchar(sequence_number) + 1L
      stem_length <- pmax(1L, label_width - suffix_length)
      displayed[indices] <- paste0(
        substr(original[indices], 1L, stem_length),
        sequence_number,
        "."
      )
    }
  }

  displayed
}


#' Deprecated alias for `set_reference_level()`
#'
#' @description
#' `biggest_reference()` is deprecated as of version 0.9.0. Use
#' [set_reference_level()] instead.
#'
#' @param x A factor.
#' @param weight A numeric vector of the same length as \code{x}.
#' @return See [set_reference_level()].
#'
#' @export
#' @keywords internal
biggest_reference <- function(x, weight) {
  lifecycle::deprecate_warn(
    "0.9.0",
    "biggest_reference()",
    "set_reference_level()"
  )
  set_reference_level(x = x, weight = weight, method = "largest_weight")
}


#' Convert p-values into significance stars
#'
#' @param pval Numeric vector of p-values.
#'
#' @return Character vector of the same length as `pval`,
#'   containing `"***"`, `"**"`, `"*"`, `"."`, or `""`.
#'
#' @keywords internal
make_stars <- function(pval) {
  if (!is.numeric(pval)) stop("`pval` must be numeric")

  stars <- cut(
    pval,
    breaks = c(-Inf, 0.001, 0.01, 0.05, 0.1, Inf),
    labels = c("***", "**", "*", ".", ""),
    right = TRUE
  )

  stars <- as.character(stars)

  stars[is.na(pval)] <- ""
  stars
}


#' @keywords internal
elapsed_days <- function(end_date) {
  as.POSIXlt(end_date)$mday - 1
}


#' @keywords internal
matchColClasses <- function(df1, df2) {

  shared_colnames <- names(df1)[names(df1) %in% names(df2)]
  shared_coltypes <- sapply(df1[, shared_colnames, drop = FALSE], class)

  for (n in shared_colnames) {
    attributes(df2[, n]) <- attributes(df1[, n])
    class(df2[, n]) <- shared_coltypes[n]
  }

  df2
}


#' Get splits from partykit object
#' @noRd
#'
#' @param x A party object.
#'
get_splits <- function(x) {
  nodes <- partykit::nodeapply(
    x,
    ids = partykit::nodeids(x),
    FUN = function(node) {
      split <- partykit::split_node(node)
      if (is.null(split)) {
        return(numeric())
      }
      as.numeric(split$breaks)
    }
  )

  sort(unique(unlist(nodes, use.names = FALSE)))
}


#' @keywords internal
update_tickmarks_right <- function(plot_obj,
                                   cut_off,
                                   max_print) {
  ranges <- suppressMessages(
    ggplot2::ggplot_build(plot_obj)$layout$panel_params[[1]]$x
  )
  label_to_add <- sprintf("[%s, %s]", round(cut_off, 1), max_print)
  tick_positions <- ranges$get_breaks()
  tick_labels    <- ranges$get_labels()

  if (overlap_right(tick_positions, cut_off)) {
    tick_positions <- tick_positions[-length(tick_positions)]
    tick_labels    <- tick_labels[-length(tick_labels)]
  }

  return(list(tick_positions = c(tick_positions, cut_off),
              tick_labels    = c(tick_labels, label_to_add)))
}


#' @keywords internal
overlap_right <- function(positions, cut_off) {

  positions <- positions[!is.na(positions)]
  n <- length(positions)
  ticks_dif <- positions[n] - positions[n - 1]
  (cut_off - positions[n]) / ticks_dif < 0.25
}


#' @importFrom ggplot2 ggplot_build
#' @keywords internal
update_tickmarks_left <- function(plot_obj,
                                  cut_off,
                                  min_print) {
  ranges <- suppressMessages(
    ggplot2::ggplot_build(plot_obj)$layout$panel_params[[1]]$x)
  label_to_add <- sprintf("[%s, %s]", min_print, round(cut_off, 1))
  tick_positions <- ranges$get_breaks()
  tick_labels    <- ranges$get_labels()

  if (overlap_left(tick_positions, cut_off)) {
    tick_positions <- tick_positions[-1]
    tick_labels    <- tick_labels[-1]
  }
  return(list(tick_positions = c(cut_off, tick_positions),
              tick_labels    = c(label_to_add, tick_labels)))
}


#' @keywords internal
overlap_left <- function(positions, cut_off) {
  positions <- positions[!is.na(positions)]
  ticks_dif <- positions[2] - positions[1]
  (positions[1] - cut_off) / ticks_dif < 0.25
}


#' @importFrom ggplot2 autoplot
#' @keywords internal
split_x_fn <- function(data, x, left = NULL, right = NULL) {

  vec <- data[[x]]
  vec_new <- vec

  if (!is.null(right)) {
    vec_new[vec_new > right] <- right
  }

  if (!is.null(left)) {
    vec_new[vec_new < left] <- left
  }

  if (!is.null(left)) {
    if (left <= min(vec, na.rm = TRUE)) {
      stop("Left should be greater than minimum value", call. = FALSE)
    }
    if (left >= max(vec, na.rm = TRUE)) {
      stop("Left should be less than maximum value", call. = FALSE)
    }
  }

  if (!is.null(right)) {
    if (right >= max(vec, na.rm = TRUE)) {
      stop("Right should be less than maximum value", call. = FALSE)
    }
    if (right <= min(vec, na.rm = TRUE)) {
      stop("Right should be greater than minimum value", call. = FALSE)
    }
  }

  if (!is.null(left) && !is.null(right)) {
    if (left >= right) {
      stop("Right should be larger than left", call. = FALSE)
    }
  }

  l1 <- split(vec_new, cut(vec_new,
                           breaks = c(min(vec, na.rm = TRUE), left,
                                      right - 1e-10, max(vec, na.rm = TRUE)),
                           include.lowest = TRUE))

  l1 <- lapply(l1, function(x) data.frame(x = x))
  if (is.null(left)) l1 <- append(list(NULL), l1)
  if (is.null(right)) l1 <- append(l1, list(NULL))
  l1
}

#' @importFrom stats as.formula
#' @keywords internal
construct_fm <- function(lhs, rhs) {
  stats::as.formula(paste0(paste0(lhs, collapse = " + "), "~ ", rhs))
}


#' @importFrom mgcv predict.gam
#' @keywords internal
confint_gam <- function(model, newdata, level = 0.95) {

  pred <- mgcv::predict.gam(
    model,
    newdata = newdata,
    type = "link",
    se.fit = TRUE
  )

  alpha <- 1 - level
  z_value <- stats::qnorm(1 - alpha / 2)

  fit_link <- pred$fit
  se_link <- pred$se.fit

  # inverse link function (family-specific)
  linkinv <- model$family$linkinv

  data.frame(
    x = newdata$x,
    fitted = linkinv(fit_link), # predicted
    conf_low = linkinv(fit_link - z_value * se_link), #lwr_95
    conf_high = linkinv(fit_link + z_value * se_link) #upr_95
  )
}

#' @keywords internal
exposure_by_factor <- function(var1, model_data, exposure) {
  x <- model_data[!is.na(model_data[[var1]]), ]
  sums <- tapply(x[[exposure]], x[[var1]], sum, na.rm = TRUE)
  df <- data.frame(
    level = names(sums),
    exposure = as.numeric(sums),
    stringsAsFactors = FALSE
  )
  df$risk_factor <- var1
  df
}


#' @keywords internal
.rating_table_model_names <- function(models, mc) {
  # mc = match.call(expand.dots = FALSE)
  dots <- mc$...
  if (is.null(dots)) return(character(0))

  # 1) Prefer explicit names in ...
  nms <- names(dots)
  if (!is.null(nms) && any(nzchar(nms))) {
    out <- nms
    if (any(!nzchar(out))) out[!nzchar(out)] <- paste0("m_", which(!nzchar(out)))
    return(out)
  }

  # 2) Otherwise: use symbol names when possible; fall back for calls (pipe etc.)
  exprs <- as.list(dots)
  out <- vapply(exprs, function(e) {
    if (is.name(e)) {
      as.character(e)
    } else {
      NA_character_
    }
  }, FUN.VALUE = character(1))

  # fallback names for non-symbols (calls, constants, etc.)
  if (length(out) == 1L && is.na(out[1L])) {
    out[1L] <- "model"
  } else {
    idx <- which(is.na(out))
    if (length(idx)) out[idx] <- paste0("m_", idx)
  }

  out
}