R/wrap-optpart.R

Defines functions ggcpt_empty build_segments fpop_wrapper

Documented in fpop_wrapper

#' FPOP wrapper: Functional Pruning Optimal Partitioning
#'
#' Wraps the \code{fpop} package for optimal changepoint detection via
#' functional pruning.
#'
#' \pkg{fpop} was archived from CRAN in 2026 at its maintainer's request and
#' is still built on R-Forge, where it is developed. Install it with
#' \code{install.packages("fpop", repos = "https://R-Forge.R-project.org")},
#' or with \code{cpt_install_engines("core")}.
#'
#' @param x A numeric vector.
#' @param penalty Penalty value. Defaults to \code{2 * log(length(x))} (BIC).
#'   This is an \emph{absolute} penalty on the residual sum of squares, so
#'   it is only calibrated for noise of standard deviation 1: on wider data
#'   the default under-penalises badly and the segmentation shatters.
#'   Standardise the series, or scale the penalty by the noise variance
#'   (for example \code{2 * log(length(x)) * stats::var(diff(x)) / 2}).
#'   See the scale-sensitivity section of \code{\link{cpt_detect}}. This
#'   default differs from the one \code{\link{cpt_detect}} applies, which
#'   resolves its \code{"MBIC"} default to a stronger numeric value, so the
#'   two entry points need not agree unless \code{penalty} is given.
#' @param ... Additional arguments passed to \code{fpop::Fpop()}.
#' @return A \code{ggcpt} object.
#' @references
#' \insertRef{maidstone2017optimal}{ggchangepoint}
#' @export
#' @family changepoint engines
#' @examplesIf requireNamespace("fpop", quietly = TRUE)
#' set.seed(2026)
#' x <- c(rnorm(60), rnorm(60, 4))
#' fpop_wrapper(x)$changepoints
#' # the default penalty assumes sd 1; scale it for wider noise
#' y <- 5 * x
#' fpop_wrapper(y, penalty = log(length(y)) * var(diff(y)))$changepoints
fpop_wrapper <- function(x, penalty = NULL, ...) {

  need_pkg("fpop")
  reject_renamed_args(list(...), "fpop")
  validate_data(x)
  data_vec <- as_uni_vector(x, "fpop")

  # As in cpop_wrapper() and decafs_wrapper(): "Manual" is for a number the
  # caller supplied, and this wrapper's own default was reported as one.
  pen_type <- if (is.null(penalty)) "2log(n) [wrapper default]" else "Manual"
  penalty <- resolve_penalty_model(penalty, data_vec)
  if (is.null(penalty)) {
    penalty <- 2 * log(length(data_vec))
  }
  validate_scalar(penalty, "penalty", min = 0)

  fit <- fpop::Fpop(data_vec, lambda = penalty, ...)
  cp_indices <- as.integer(fit$t.est)

  ggcpt_build(
    data_vec, cp_indices,
    method = "fpop",
    change_in = "mean",
    penalty = list(type = pen_type, value = penalty),
    fit = fit,
    call = match.call()
  )
}

# Internal: build segments tibble from changepoint indices
#' @noRd
build_segments <- function(data_vec, cp_indices) {
  n <- length(data_vec)
  starts <- c(1L, cp_indices + 1L)
  ends <- c(cp_indices, n)
  n_seg <- length(starts)

  tibble::tibble(
    seg_id = seq_len(n_seg),
    start = starts,
    end = ends,
    n = ends - starts + 1L,
    # `na.rm`: under `na_action = "engine"` the series reaches here with
    # its gaps, and a segment's level is the mean of what was observed.
    param_estimate = vapply(seq_len(n_seg), function(i) {
      mean(data_vec[starts[i]:ends[i]], na.rm = TRUE)
    }, numeric(1))
  )
}

# Internal: create an empty ggcpt for no-change results
#' @noRd
ggcpt_empty <- function(data_vec, method = "unknown") {
  data_tbl <- tibble::tibble(index = seq_along(data_vec), value = as.numeric(data_vec))
  segments <- tibble::tibble(
    seg_id = 1L,
    start = 1L,
    end = length(data_vec),
    n = length(data_vec),
    param_estimate = mean(data_vec, na.rm = TRUE)
  )
  new_ggcpt(
    changepoints = tibble::tibble(cp = integer(), cp_value = numeric()),
    segments = segments,
    data = data_tbl,
    method = method,
    change_in = "mean",
    penalty = list(type = NA_character_, value = NA_real_),
    fit = NULL,
    call = NULL
  )
}

Try the ggchangepoint package in your browser

Any scripts or data that you put into this service are public.

ggchangepoint documentation built on Nov. 5, 2026, 9:06 a.m.