Getting Started with scanr

The scanr package detects multiple change points in long univariate time series using window-based SCAN statistics. The method compares local windows along a series, calibrates evidence for changes using resampling, and combines information across several window sizes.

2 Load the Package

install.packages("scanr")

Alternatively, the development version of scanr can be installed directly from GitHub:

if (!requireNamespace("pak", quietly = TRUE)) { install.packages("pak") }

pak::pak("Prabashoka/scanr")

3 Example 1 - Changes in Mean

We first simulate a white-noise time series of length 20,000 with 20 change points and 21 piecewise-constant mean regimes. The noise variance is fixed, so this example focuses on changes in the mean.

set.seed(1234)

n <- 20000

change_points <- c(
  952, 1905, 2858, 3810, 4763, 5715, 6668, 7620, 8573, 9525,
  10478, 11430, 12383, 13335, 14288, 15240, 16193, 17145, 18098,
  19050
)

means <- c(
  0, 2, -1, 3, 0.5, -2, 2, 5, -0.5, 2.5, 0, -2.5, -1.5, 1.5,
  3, 1, 0, 1.25, -2, 3.5, -1.5
)

segment_starts <- c(1L, change_points + 1L)
segment_ends <- c(change_points, n)
x_mean <- numeric(n)

for (j in seq_along(means)) {
  segment_index <- segment_starts[j]:segment_ends[j]
  x_mean[segment_index] <- rnorm(length(segment_index), mean = means[j], sd = 1)
}

change_points
#>  [1]   952  1905  2858  3810  4763  5715  6668  7620  8573  9525 10478 11430
#> [13] 12383 13335 14288 15240 16193 17145 18098 19050

3.1 Run SCAN Change-Point Detection

The window_sizes argument controls the local scales used by the scan. The window sizes should be smaller than the spacing between nearby change points, while still being large enough to estimate the local distributions on each side of a candidate split. Smaller values improve localization and sensitivity to closely spaced changes, while larger values provide more stable local estimates.

default_window_sizes() randomly samples the requested number of distinct scales. By default, the largest is floor(n^(2 / 3)); it is always capped at floor(n / 2) so that both sides of a candidate split can contain a full window. Its optional seed argument makes the sampled grid reproducible without changing the caller’s random-number generator state. Supplying problem-informed bounds is useful when the approximate minimum segment length is known. Omitting window_sizes makes scan_cpd() use the same helper with its min_window, max_window, and n_windows arguments.

mean_window_sizes <- default_window_sizes(
  n = length(x_mean),
  min_window = 100,
  max_window = floor(length(x_mean)^(2 / 3)),
  n_windows = 7,
  seed = 52
)
mean_window_sizes
#> [1] 154 181 232 241 278 438 478

fit_mean <- scan_cpd(
  x_mean,
  window_sizes = mean_window_sizes,
  n_boot = 400,
  random_state = 1234,
  change_type = "mean",
  n_jobs = -1
)

fit_mean
#> scanr change-point result
#> observations: 20000 
#> change points: 952, 1905, 2859, 3809, 4760, 5715, 6668, 7620, 8573, 9525, 10478, 11437, 12383, 13337, 14288, 15241, 16193, 17145, 18098, 19050

Remark on bootstrap calibration

The adaptive threshold for declaring a change point is derived through a tapered block bootstrap across the combined windows. Increasing the number of bootstrap samples (n_boot) produces a more accurate null distribution at the cost of longer computation time; a value of 400–1,000 is typically sufficient. The method supports parallel processing via Rust’s Rayon library for computational efficiency.

fit_mean$change_points
#>  [1]   952  1905  2859  3809  4760  5715  6668  7620  8573  9525 10478 11437
#> [13] 12383 13337 14288 15241 16193 17145 18098 19050

3.2 Evaluate Detection Accuracy

For simulated data, estimated change points can be compared with the known truth. The tolerance argument controls how close an estimated change point must be to a true change point to count as a match.

cpd_metrics(
  true_cps = change_points,
  estimated_cps = fit_mean$change_points,
  n = length(x_mean),
  tolerance = 20
)
#> $matches
#>     true estimated distance
#> 1    952       952        0
#> 2   1905      1905        0
#> 3   5715      5715        0
#> 4   6668      6668        0
#> 5   7620      7620        0
#> 6   8573      8573        0
#> 7   9525      9525        0
#> 8  10478     10478        0
#> 9  12383     12383        0
#> 10 14288     14288        0
#> 11 16193     16193        0
#> 12 17145     17145        0
#> 13 18098     18098        0
#> 14 19050     19050        0
#> 15  2858      2859        1
#> 16  3810      3809        1
#> 17 15240     15241        1
#> 18 13335     13337        2
#> 19  4763      4760        3
#> 20 11430     11437        7
#> 
#> $precision
#> [1] 1
#> 
#> $recall
#> [1] 1
#> 
#> $f1
#> [1] 1
#> 
#> $covering
#> [1] 0.9985036

3.3 Visualize the Result

The package includes visual helpers for inspecting detected change points, window-level votes, and the overall vote scree.

vis_change_points(
  x_mean,
  fit_mean,
  true_change_points = change_points,
  x_label = "Time",
  y_label = "Value"
)

3.4 Tune the Voting Threshold

The ensemble vote threshold, vote_threshold, controls how much agreement is required across window sizes before a candidate is retained. The default value is 0.5, meaning that a candidate must be supported by at least half of the available window-level votes.

vis_vote_scree(fit_mean)

vis_window_votes(fit_mean)

4 Example 2 - Changes in Distribution

The next example uses segments from different distribution families, including normal, exponential, Poisson, gamma, uniform, Student’s t, lognormal, beta, Weibull, and chi-square distributions. Each segment is standardized before being rescaled, so the signal is not only a simple shift in the mean. This creates a more challenging distributional change-point problem.

set.seed(1234)

n <- 20000

change_points <- c(952, 1905, 2858, 3810, 4763, 5715, 6668, 7620, 8573, 9525, 10478, 11430, 12383, 13335, 14288, 15240, 16193, 17145, 18098, 19050)

segment_starts <- c(1, change_points + 1)
segment_ends <- c(change_points, n)

families <- c("normal", "exponential", "poisson", "t", "gamma", "uniform", "lognormal", "weibull", "chisq", "beta", "normal", "poisson", "exponential", "uniform", "gamma", "t", "lognormal", "beta", "weibull", "chisq", "normal")

scale_factors <- c(0.7, 1.4, 0.6, 1.8, 0.75, 2.5, 0.65, 3.0, 0.7, 2.6, 0.6, 1.7, 0.65, 2.9, 0.7, 2.8, 0.6, 2.6, 0.65, 3.0, 0.7)

simulate_base <- function(m, family) {
  z <- switch(family, normal = rnorm(m, mean = 0, sd = 1), exponential = rexp(m, rate = 1), poisson = rpois(m, lambda = 2), gamma = rgamma(m, shape = 2, rate = 1), uniform = runif(m, min = -sqrt(3), max = sqrt(3)), t = rt(m, df = 3), lognormal = rlnorm(m, meanlog = 0, sdlog = 0.7), beta = rbeta(m, shape1 = 2, shape2 = 5), weibull = rweibull(m, shape = 1.5, scale = 1), chisq = rchisq(m, df = 5))
  as.numeric((z - mean(z)) / sd(z))
}

x_dist <- numeric(n)

for (j in seq_along(families)) {
  segment_index <- segment_starts[j]:segment_ends[j]
  x_dist[segment_index] <-scale_factors[j] * simulate_base(length(segment_index), families[j])
}

change_points
#>  [1]   952  1905  2858  3810  4763  5715  6668  7620  8573  9525 10478 11430
#> [13] 12383 13335 14288 15240 16193 17145 18098 19050

4.1 Run SCAN Change-Point Detection

For distributional changes, set change_type = "distribution". This uses the distribution-sensitive local statistic rather than a mean-only statistic.

distribution_window_sizes <- default_window_sizes(
  n = length(x_dist),
  min_window = 100,
  max_window = floor(length(x_mean)^(2 / 3)),
  n_windows = 15,
  seed = 52
)

fit_dist <- scan_cpd(
  x_dist,
  window_sizes = distribution_window_sizes,
  n_boot = 1000,
  random_state = 1234,
  change_type = "distribution",
  vote_threshold = 0.5,
  n_jobs = -1
)

fit_dist$change_points
#>  [1]   954  1902  2860  3809  4764  5709  6669  7611  8573  9522 10478 11430
#> [13] 12383 13335 14289 15240 16195 17145 18100 19048

4.2 Evaluate Detection Accuracy

cpd_metrics(
  true_cps = change_points,
  estimated_cps = fit_dist$change_points,
  n = length(x_dist),
  tolerance = 20
)
#> $matches
#>     true estimated distance
#> 1   8573      8573        0
#> 2  10478     10478        0
#> 3  11430     11430        0
#> 4  12383     12383        0
#> 5  13335     13335        0
#> 6  15240     15240        0
#> 7  17145     17145        0
#> 8   3810      3809        1
#> 9   4763      4764        1
#> 10  6668      6669        1
#> 11 14288     14289        1
#> 12   952       954        2
#> 13  2858      2860        2
#> 14 16193     16195        2
#> 15 18098     18100        2
#> 16 19050     19048        2
#> 17  1905      1902        3
#> 18  9525      9522        3
#> 19  5715      5709        6
#> 20  7620      7611        9
#> 
#> $precision
#> [1] 1
#> 
#> $recall
#> [1] 1
#> 
#> $f1
#> [1] 1
#> 
#> $covering
#> [1] 0.9965096

4.3 Visualize the Result

vis_change_points(
  x_dist,
  fit_dist,
  true_change_points = change_points,
  x_label = "Time",
  y_label = "Value"
)

4.4 Real Data Example

4.5 Load the Data

Load the dataset and resample it at one-minute intervals for change-point detection.

head(swat)
#>                Timestamp   PIT502 Normal/Attack
#> 1 28/12/2015 10:00:00 AM 1.649953        Normal
#> 2 28/12/2015 10:01:00 AM 1.649953        Normal
#> 3 28/12/2015 10:02:00 AM 1.665972        Normal
#> 4 28/12/2015 10:03:00 AM 1.681991        Normal
#> 5 28/12/2015 10:04:00 AM 1.746067        Normal
#> 6 28/12/2015 10:05:00 AM 1.681991        Normal

4.6 Choose Window Sizes

The largest candidate window is ⌊n^2/3⌋. The package helper selects 15 default window sizes between 20 minutes and that upper bound.

n <- nrow(swat)
window_sizes <- default_window_sizes(
  n,
  min_window = 100,
  max_window = floor(n^(2 / 3)),
  n_windows = 10,
  seed = 400
  
)
window_sizes
#>  [1] 201 255 426 442 519 587 607 624 625 652

The PIT-502 measurements are standardized before detection. Standardization changes their units but preserves the locations of mean shifts.

x <- as.numeric(scale(swat[["PIT502"]]))

4.7 Detect change-points

The R interface uses n_boot and random_state for bootstrap calibration and reproducibility.

fit_swat <- scan_cpd(
  x,
  window_sizes = window_sizes,
  n_boot = 400,
  vote_threshold = 0.15,
  random_state = 100,
  change_type = "mean",
  n_jobs = -1
)

fit_swat
#> scanr change-point result
#> observations: 17477 
#> change points: 517, 1222, 1901, 2596, 3311, 4232, 4952, 5610, 6643, 7542, 8291, 8972, 9799, 10689, 11277, 12063, 12741, 13588, 13960, 14718, 15131, 15792, 16541, 17201

The fitted object stores change points as positions in the one-minute series. We can map those positions back to timestamps, observed pressures, and the dataset’s operating-state labels.

detected_changes <- swat[
  fit_swat$change_points,
  c("Timestamp", "PIT502", "Normal/Attack"),
  drop = FALSE
]

cat("Number of change-points:", nrow(detected_changes), "\n")
#> Number of change-points: 24

4.8 Visualize the Result

vis_change_points(
  x,
  fit_swat,
  index = seq_along(x),
  x_label = "Time",
  y_label = "Standardized PIT-502 pressure",
  title = "Mean changes in the SWaT PIT-502 pressure sensor"
)

4.9 Inspect the Vote Scree

The vote scree shows how much support each candidate location receives across the selected window sizes. The horizontal threshold corresponds to the vote_threshold = 0.15 used above.

vis_vote_scree(fit_swat)

5 SWAL Statistic for Single Change-Point Detection

The SCAN procedure above combines local evidence across many windows. The package also exposes lower-level localization tools for a single region that is believed to contain one change point.

For a mean change, the mean version of the SWAL statistic aligns with the CUSUM localizer. Both methods search for the split that best separates the left and right sides of the local region.

set.seed(1234)

true_single_cp <- 150
mean_region <- c(
  rnorm(true_single_cp, mean = 0, sd = 1),
  rnorm(true_single_cp, mean = 2, sd = 1)
)

c(
  truth = true_single_cp,
  cusum = ts_cusum(mean_region),
  swal_distribution = swal_statistic(mean_region, change_type = "distribution")
)
#>             truth             cusum swal_distribution 
#>               150               150               150

The distributional version is more general. It can localize changes where the mean is approximately unchanged but other aspects of the distribution, such as variance, change sharply. In the example below, both segments have mean zero, but the second segment has a much larger standard deviation.

set.seed(1234)

var_region <- c(
  rnorm(true_single_cp, mean = 0, sd = 0.5),
  rnorm(true_single_cp, mean = 0, sd = 2)
)

c(
  truth = true_single_cp,
  cusum = ts_cusum(var_region),
  swal_distribution = swal_statistic(var_region, change_type = "distribution")
)
#>             truth             cusum swal_distribution 
#>               150               164               152

The ts_wasserstein() function returns both the estimated split and the full sequence of split statistics. The plotting helper vis_swal_curve() displays that localization curve, making it easier to inspect where the statistic is maximized.

vis_swal_curve(
  var_region,
  start = 1,
  end = length(var_region),
  x_label = "Candidate split",
  y_label = "Scaled Wasserstein statistic",
  title = ""
)

6 Session Information

Recording session information helps make the vignette reproducible across R versions and operating systems.

sessionInfo()
#> R version 4.6.1 (2026-06-24 ucrt)
#> Platform: x86_64-w64-mingw32/x64
#> Running under: Windows 10 x64 (build 19045)
#> 
#> Matrix products: default
#>   LAPACK version 3.12.1
#> 
#> locale:
#> [1] LC_COLLATE=C                       LC_CTYPE=English_Australia.utf8   
#> [3] LC_MONETARY=English_Australia.utf8 LC_NUMERIC=C                      
#> [5] LC_TIME=English_Australia.utf8    
#> 
#> time zone: Australia/Sydney
#> tzcode source: internal
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] scanr_0.1.1
#> 
#> loaded via a namespace (and not attached):
#>  [1] vctrs_0.7.3        cli_3.6.6          knitr_1.51         rlang_1.2.0       
#>  [5] xfun_0.59          otel_0.2.0         generics_0.1.4     S7_0.2.2          
#>  [9] jsonlite_2.0.0     labeling_0.4.3     glue_1.8.1         htmltools_0.5.9   
#> [13] scales_1.4.0       rmarkdown_2.31     grid_4.6.1         tibble_3.3.1      
#> [17] evaluate_1.0.5     fastmap_1.2.0      yaml_2.3.12        lifecycle_1.0.5   
#> [21] compiler_4.6.1     dplyr_1.2.1        RColorBrewer_1.1-3 pkgconfig_2.0.3   
#> [25] farver_2.1.2       digest_0.6.39      R6_2.6.1           tidyselect_1.2.1  
#> [29] parallel_4.6.1     pillar_1.11.1      magrittr_2.0.5     withr_3.0.3       
#> [33] tools_4.6.1        gtable_0.3.6       ggplot2_4.0.3

1 Overview