Risk-Based Control Charts: A Complete Workflow

Purpose

A conventional control chart classifies the process from measured data. When measurement uncertainty is material, the measured decision can disagree with the state indicated by a reference or error-free value. rbcc represents the four resulting outcomes explicitly:

Reference state Measured decision Cost
In control Accept c11
In control Signal c10
Out of control Accept c01
Out of control Signal c00

The package accepts paired reference observations X and additive measurement errors UC, forms the observed values X + UC, constructs the selected chart for both series, and minimizes the total decision cost over a bounded chart coefficient.

A paired-measurement example

The example contains stable subgroups followed by a location shift. Measurement errors are generated independently and the missed-signal cost is larger than the false-signal cost.

set.seed(20261)
n <- 4L
groups <- 60L
X <- c(
  rnorm(40L * n, mean = 10.0, sd = 0.45),
  rnorm(20L * n, mean = 10.9, sd = 0.45)
)
UC <- rnorm(groups * n, mean = 0, sd = 0.08)
C <- c(c11 = 1, c10 = 8, c01 = 70, c00 = 4)

The default optimizer is the exact empirical threshold search.

xbar_fit <- rbcc_opt(
  X, UC, C,
  n = n,
  type = "xbar",
  LKL = 0,
  UKL = 5,
  optimizer = "exact"
)
xbar_fit$optimization
#> $par
#> [1] 2.847659
#> 
#> $value
#> [1] 108
#> 
#> $method
#> [1] "exact threshold search"
#> 
#> $evaluations
#> [1] 63
#> 
#> $convergence
#> [1] 0
summary(xbar_fit)
#> 
#> Risk-based XBAR control chart
#> ------------------------------------------
#>                   outcome count aggregate_cost
#>  Correct acceptance (c11)    44             44
#>        Type I error (c10)     0              0
#>       Type II error (c01)     0              0
#>     Correct control (c00)    16             64
#> Total decision cost: 108 
#> 
#> Traditional limits:  9.619397 10.949775 
#> Risk-based limits:  9.65317 10.91600 
#> Optimal coefficient: 2.847659 
#> 
#> Statistic summary:
#>              Min.     Mean     Max.
#> true     9.380060 10.28459 11.33705
#> observed 9.371395 10.28678 11.35878

The fitted object retains the four counts and their cost components, so the reported objective can be audited directly.

cbind(
  count = xbar_fit$decision_counts,
  unit_cost = C,
  aggregate_cost = xbar_fit$cost_components
)
#>     count unit_cost aggregate_cost
#> q11    44         1             44
#> q10     0         8              0
#> q01     0        70              0
#> q00    16         4             64
sum(xbar_fit$cost_components)
#> [1] 108
plot(xbar_fit)
Reference and observed X-bar statistics with conventional and optimized risk-based limits.

Reference and observed X-bar statistics with conventional and optimized risk-based limits.

The chart family

All univariate chart constructors use the same paired inputs, cost vector, optimization controls, and fitted-object conventions. Their monitored feature and memory structure differ.

fits <- list(
  Xbar = xbar_fit,
  R = rbcc_opt(X, UC, C, n = n, type = "R"),
  S = rbcc_opt(X, UC, C, n = n, type = "S"),
  MA = rbmacc_opt(X, UC, C, n = n, w = 3),
  EWMA = rbewmacc_opt(X, UC, C, n = n, lambda = 0.20),
  CUSUM = rbcusumcc_opt(X, UC, C, n = n, T = 5, se.shift = 1)
)

family_summary <- data.frame(
  chart = names(fits),
  K = vapply(fits, function(x) unname(x$Kopt), numeric(1L)),
  total_cost = vapply(fits, function(x) unname(x$cost0), numeric(1L)),
  type_I = vapply(fits, function(x) unname(x$q10), numeric(1L)),
  type_II = vapply(fits, function(x) unname(x$q01), numeric(1L)),
  evaluations = vapply(
    fits,
    function(x) as.numeric(x$optimization$evaluations),
    numeric(1L)
  )
)
knitr::kable(family_summary, digits = 3)
chart K total_cost type_I type_II evaluations
Xbar Xbar 2.848 108 0 0 63
R R 2.475 60 0 0 63
S S 1.069 60 0 0 63
MA MA 3.145 150 0 0 50
EWMA EWMA 2.736 213 0 0 36
CUSUM CUSUM 5.024 210 0 0 15

The table is not a ranking of universally superior charts. X-bar, MA, EWMA, and CUSUM emphasize location behavior with different memory; R and S monitor within-subgroup dispersion. Chart selection should therefore follow the process feature and shift pattern of interest before the cost-optimal limit is selected.

Why exact optimization is the default

For a two-sided scalar chart, each observed statistic has a standardized distance from its center. Its classification changes only when the candidate coefficient crosses that distance. The empirical total cost is therefore piecewise constant, and a global minimum over the specified interval occurs at a finite set of event thresholds or interval endpoints. The exact optimizer sorts those events and updates the cost cumulatively.

Two alternatives are available for comparison or sensitivity analysis:

control <- rbcc_control(grid_size = 201L)
optimizer_fits <- list(
  exact = rbcc_opt(
    X, UC, C, n = n, type = "xbar", optimizer = "exact"
  ),
  bounded = rbcc_opt(
    X, UC, C, n = n, type = "xbar", optimizer = "optimize"
  ),
  grid = rbcc_opt(
    X, UC, C, n = n, type = "xbar", optimizer = "grid",
    control = control
  )
)
optimizer_summary <- data.frame(
  method = names(optimizer_fits),
  K = vapply(optimizer_fits, function(x) x$Kopt, numeric(1L)),
  total_cost = vapply(optimizer_fits, function(x) x$cost0, numeric(1L)),
  evaluations = vapply(
    optimizer_fits,
    function(x) as.numeric(x$optimization$evaluations),
    numeric(1L)
  )
)
knitr::kable(optimizer_summary, digits = 4)
method K total_cost evaluations
exact exact 2.8477 108 63
bounded bounded 2.8604 108 NA
grid grid 2.8500 108 201

A bounded continuous optimizer is retained for compatibility and sensitivity analysis, but the empirical objective is not smooth. A grid search is transparent and useful when users deliberately want evaluation on a prescribed grid.

Optional parallel grid evaluation

Grid points are independent. parallel = TRUE distributes their evaluation by a PSOCK cluster on Windows and by forked workers on Unix-like systems. Parallel startup and serialization overhead can dominate small examples, so parallelism is opt-in and the exact optimizer remains the recommended scalar default.

parallel_grid_fit <- rbcc_opt(
  X, UC, C,
  n = n,
  type = "xbar",
  optimizer = "grid",
  parallel = TRUE,
  workers = 2,
  control = rbcc_control(grid_size = 2001L)
)

For reproducibility studies, report the bounds, grid size, worker count, operating system, R version, and elapsed time. Compare objective values rather than expecting bit-for-bit identical timings across machines.

Multivariate monitoring

The same decision-cost logic applies to the Hotelling T-squared chart. The optimized parameter is reported as the correction subtracted from the traditional upper control limit.

data("t2uc", package = "rbcc")
X_multi <- as.matrix(t2uc[, 1:2])
UC_multi <- as.matrix(t2uc[, 5:6])
C_multi <- c(c11 = 1, c10 = 20, c01 = 160, c00 = 5)

multi_fit <- rbmcc_opt(
  X_multi, UC_multi, C_multi,
  n = 1,
  confidence_level = 0.99,
  LKL = -5,
  UKL = 5,
  optimizer = "exact"
)
summary(multi_fit)
#> 
#> Risk-based multivariate T-squared control chart
#> --------------------------------------------------
#>                   outcome count aggregate_cost
#>  Correct acceptance (c11)    49             49
#>        Type I error (c10)     0              0
#>       Type II error (c01)     0              0
#>     Correct control (c00)     1              5
#> Total decision cost: 54 
#> 
#> Traditional UCL: 8.545579 
#> Risk-based UCL: 9.551365 
#> Optimal correction: -1.005787 
#> 
#> T-squared statistic summary:
#>                Min. Mean      Max.
#> true     0.10886812 1.96 10.862370
#> observed 0.02867713 1.96  9.779256

Programmatic use

Assignments to summary() and plot() are intentionally silent. The returned objects can be inspected, stored, combined, or printed later.

s <- summary(xbar_fit)
p <- plot(xbar_fit)

names(s)
#>  [1] "call"               "chart"              "sample_size"       
#>  [4] "confidence_level"   "K"                  "LCLx"              
#>  [7] "UCLx"               "LCLy"               "UCLy"              
#> [10] "real"               "Observed"           "data_info"         
#> [13] "cost0"              "cost1"              "cost2"             
#> [16] "cost3"              "cost4"              "q11"               
#> [19] "q10"                "q01"                "q00"               
#> [22] "decision_counts"    "cost_components"    "decision"          
#> [25] "par"                "Kopt"               "optimization"      
#> [28] "statistics_summary"
class(p)
#> [1] "ggplot2::ggplot" "ggplot"          "ggplot2::gg"     "S7_object"      
#> [5] "gg"

A robust applied analysis should retain the fitted object, report the complete cost vector and decision counts, examine plausible alternative costs, and validate prospective statistical performance separately from the empirical calibration cost.