Calibrate supplied predictions with binary or multiclass Inductive/Cross methods, or moderated regression, and draw their fit diagnoses and assessment graphics. These examples use base R and the installed package.
Start with
vignette("manual-calibration", package = "vennabers") for
the step-by-step guide and
vignette("diagnostics", package = "vennabers") for the
graphics. ?vennabers_examples indexes all four example
collections.
Running an example
Each heading below names an installed script. To inspect and run one:
script <- system.file("examples", "binary-options.R", package = "vennabers")
file.show(script)
source(script)Replace the filename to run another script. The listings below show the complete code without executing it during documentation builds. Sourcing a script runs its fits and checks; resampling, tuning and composition matrices can take longer than the introductory guides.
Binary precision, loss and event selection
binary-options.R: Compare Log/Brier regularization,
calibration rounding, class predictions and Cross fold diagnostics.
# After installation: source(system.file("examples", "binary-options.R", package = "vennabers"))
library(vennabers)
# Precomputed probabilities from a fixed model on its separate calibration set.
# No calibration or test truths were used to train that model or preprocessing.
truth <- c(0, 1, 0, 1)
calibration_probability <- c(.501, .502, .509, .510)
new_probability <- c(.499, .505, .511)
inductive <- venn_abers(truth, calibration_probability, precision = 2)
for (loss in c("log", "brier")) {
prediction <- predict(inductive, new_probability, loss = loss)
classes <- predict(inductive, new_probability, loss = loss, type = "class")
stopifnot(identical(classes, prediction$probability >= .5))
print(prediction)
}
# Precision affects calibration only; choosing a new precision requires a fit.
# The new values above are deliberately not rounded to two decimal places.
# To select class 0 as the event, reverse truth AND the supplied probabilities.
# Round only after selecting the event; rounded results need not be complements.
other_event <- venn_abers(truth == 0, 1 - calibration_probability, precision = 2)
print(predict(other_event, 1 - new_probability))
# Inductive bounds are VA endpoints, not confidence intervals.
# Cross uses held-out probabilities from each corresponding retained fold model.
# These illustrative folds have unequal sizes; each still receives equal weight.
# See cross-venn-abers.R for complete model fitting and held-out data generation.
cross <- cross_venn_abers(
truth = c(0, 1, 1, 1), probability = c(.501, .201, .501, .801),
fold = c("A", "B", "B", "B"), precision = 2
)
# Each column comes from the model that excluded its named calibration fold.
# All columns refer to the same new observations, in the same order.
new_probability <- cbind(B = c(.201, .501, .901), A = c(.199, .501, .801))
test_truth <- c(0, 1, 1) # Used only for assessment, after fitting.
local({
plot_file <- tempfile(fileext = ".pdf")
on.exit(unlink(plot_file), add = TRUE)
grDevices::pdf(plot_file)
on.exit(grDevices::dev.off(), add = TRUE, after = FALSE)
for (loss in c("log", "brier")) {
details <- predict(cross, new_probability, loss = loss, type = "details")
reconstructed <- vapply(seq_len(nrow(new_probability)), function(row) {
f <- details$folds[details$folds$row == row, ]
if (loss == "log") {
plogis(mean(log(f$upper)) - mean(log1p(-f$lower)))
} else mean(f$probability)
}, 0)
stopifnot(max(abs(reconstructed - details$prediction$probability)) < 1e-10,
identical(predict(cross, new_probability, loss = loss, type = "class"),
details$prediction$probability >= .5))
print(details)
print(c(brier = brier_score(test_truth, reconstructed),
log = log_loss(test_truth, reconstructed)))
calibration_plot(test_truth, reconstructed)
}
})
# Classes use the combined probability with TRUE at .5, never a fold vote.
# Fold details concern new observations; there is no aggregate Cross interval.
# No improvement in test scores or perfect Cross calibration is promised.Binary Cross calibration with retained models
cross-venn-abers.R: Fit separate GLMs, retain the
corresponding fold models and assess calibrated probabilities on
held-out observations.
# Run from the checkout root after installation:
# Rscript inst/examples/cross-venn-abers.R
library(vennabers)
source(system.file("examples", "assertions.R", package = "vennabers"))
set.seed(20260906)
observations <- data.frame(x1 = rnorm(900L), x2 = rnorm(900L))
observations$truth <- rbinom(900L, 1L, stats::plogis(
-0.4 + observations$x1 - 0.6 * observations$x2
)) # The event is 1 (equivalently TRUE).
# Set aside test observations before any model fitting or calibration.
test_rows <- sample.int(nrow(observations), 300L)
test <- observations[test_rows, ]
development <- observations[-test_rows, ]
fold_ids <- c("A", "B", "C")
fold <- sample(rep(fold_ids, length.out = nrow(development)))
models <- setNames(vector("list", length(fold_ids)), fold_ids)
out_of_fold_probability <- numeric(nrow(development))
for (id in fold_ids) {
held_out <- fold == id
# Exclusion must cover the complete learned pipeline, including any learned
# preprocessing or tuning. This example needs neither.
models[[id]] <- stats::glm(truth ~ x1 + x2,
data = development[!held_out, ], family = stats::binomial()
)
# Restore held-out probabilities to their original development positions.
out_of_fold_probability[held_out] <- predict(models[[id]],
newdata = development[held_out, ], type = "response"
)
}
# Retain these models: each predicts the same test rows in the same order.
# A single all-data refit cannot replace them. Names cannot prove model
# provenance, training exclusion, unique observations, or semantic alignment.
test_probability <- vapply(models, function(model) {
predict(model, newdata = test[c("x1", "x2")], type = "response")
}, numeric(nrow(test)))
# The package owns per-fold calibration and equal-weight log-loss aggregation.
# Model fitting and fold creation remain caller responsibilities.
calibrator <- cross_venn_abers(
development$truth, out_of_fold_probability, fold
)
prediction <- predict(calibrator, test_probability)
assert_conditions(
is.data.frame(prediction), nrow(prediction) == nrow(test),
identical(names(prediction), "probability"),
is.double(prediction$probability), all(is.finite(prediction$probability)),
all(prediction$probability >= 0 & prediction$probability <= 1)
)
# Test truth is used only for evaluation. No score improvement is promised.
scores <- c(
brier_score = brier_score(test$truth, prediction$probability),
log_loss = log_loss(test$truth, prediction$probability)
)
stopifnot(all(is.finite(scores)))
print(scores)
# A temporary PDF also works without an interactive display. Close and remove
# it even on error. For an interactive plot, call calibration_plot() directly.
local({
plot_file <- tempfile(fileext = ".pdf")
on.exit(unlink(plot_file), add = TRUE)
grDevices::pdf(plot_file)
on.exit(grDevices::dev.off(), add = TRUE, after = FALSE)
summary <- calibration_plot(test$truth, prediction$probability)
stopifnot(sum(summary$count) == nrow(test))
})
# Cross output is a regularized event probability, with no public aggregate
# interval. Individual fold diagnostics are available with type = "details".
# The combined probability inherits no IVAP perfect-calibration
# guarantee and makes no confidence-interval coverage claim.Inductive multiclass one-vs-one
inductive-ovo.R: Calibrate fixed class probabilities,
inspect pair diagnostics and score the coupled class matrix.
# After package installation: Rscript inst/examples/inductive-ovo.R
library(vennabers)
# Illustrative caller-supplied outputs from one already fixed multiclass model.
# Calibration and assessment observations are separate from that model's
# training, preprocessing, and tuning; the package does not train the model.
classes <- c("cat", "dog", "bird")
truth <- factor(rep(classes, 2), levels = classes)
p <- rbind(c(.7, .2, .1), c(.2, .6, .2), c(.1, .2, .7),
c(.5, .3, .2), c(.3, .5, .2), c(.2, .3, .5))
colnames(p) <- classes
fit <- venn_abers(truth, p, precision = 2, method = "ovo")
print(fit)
new_p <- rbind(c(.6, .3, .1), c(.1, .4, .5), c(.1, .8, .1))
colnames(new_p) <- classes
assessment <- factor(c("cat", "bird", "dog"), levels = classes)
local({
plot_file <- tempfile(fileext = ".pdf")
grDevices::pdf(plot_file)
on.exit({ grDevices::dev.off(); unlink(plot_file) })
for (loss in c("log", "brier")) {
details <- predict(fit, new_p[, rev(classes)], loss = loss, type = "details")
probability <- predict(fit, new_p, loss = loss)
predicted_class <- predict(fit, new_p, loss = loss, type = "class")
reconstructed <- probability
for (row in seq_len(nrow(new_p))) {
pairs <- details$pairs[details$pairs$row == row, ]
for (class in classes) {
q <- c(pairs$probability[pairs$event == class],
1 - pairs$probability[pairs$non_event == class])
reconstructed[row, class] <- 1 / (sum(1 / q) - (length(classes) - 2))
}
}
reconstructed <- reconstructed / rowSums(reconstructed)
stopifnot(max(abs(reconstructed - probability)) <= 1e-10,
identical(probability, details$prediction),
identical(predicted_class,
factor(classes[max.col(probability, ties.method = "first")], levels = classes)))
print(details)
print(predicted_class)
print(c(brier = brier_score(assessment, probability), log = log_loss(assessment, probability)))
groups <- calibration_plot(assessment, probability)
stopifnot(all(tapply(groups$count, groups$class, sum) == length(assessment)))
print(groups)
}
})
# Pair endpoints have only their binary interpretation. Assessment scores use
# the final matrix; Brier sums class errors without dividing by class count.Inductive multiclass one-vs-rest
inductive-ovr.R: Calibrate each class against the rest
and assess the final normalized probabilities.
# After package installation: Rscript inst/examples/inductive-ovr.R
library(vennabers)
# Illustrative probabilities from an already fixed multiclass model.
# Calibration and assessment are separate from its training and tuning.
classes <- c("cat", "dog", "bird")
truth <- factor(rep(classes, 2), levels = classes)
p <- rbind(c(.7, .2, .1), c(.2, .6, .2), c(.1, .2, .7),
c(.5, .3, .2), c(.3, .5, .2), c(.2, .3, .5))
colnames(p) <- classes
fit <- venn_abers(truth, p, method = "ovr", precision = 2)
print(fit)
new_p <- rbind(c(.6, .3, .1), c(.1, .4, .5), c(.1, .8, .1))
colnames(new_p) <- classes
assessment <- factor(c("cat", "bird", "dog"), levels = classes)
local({
plot_file <- tempfile(fileext = ".pdf")
grDevices::pdf(plot_file)
on.exit({ grDevices::dev.off(); unlink(plot_file) })
for (loss in c("log", "brier")) {
details <- predict(fit, new_p[, rev(classes)], loss = loss, type = "details")
probability <- predict(fit, new_p, loss = loss)
predicted_class <- predict(fit, new_p, loss = loss, type = "class")
points <- vapply(classes, function(class) {
parts <- details$classes[details$classes$event == class, ]
parts$probability[match(seq_len(nrow(new_p)), parts$row)]
}, numeric(nrow(new_p)))
stopifnot(max(abs(points / rowSums(points) - probability)) <= 1e-10,
identical(probability, details$prediction),
identical(predicted_class,
factor(classes[max.col(probability, ties.method = "first")], levels = classes)))
print(details)
print(predicted_class)
print(c(brier = brier_score(assessment, probability), log = log_loss(assessment, probability)))
groups <- calibration_plot(assessment, probability)
stopifnot(all(tapply(groups$count, groups$class, sum) == length(assessment)))
print(groups)
}
})
# Endpoints concern each binary class-against-rest problem. Scores and plots
# use the final normalized matrix. Manual Brier is an R extension over Python.Cross multiclass calibration
cross-multiclass.R: Supply fold-matched class matrices
and reconstruct the OvO/OvR aggregation from binary diagnostics.
# After package installation: Rscript inst/examples/cross-multiclass.R
library(vennabers)
# Precomputed out-of-fold probabilities from two retained multiclass models.
# Model A excluded rows 1:2; model B excluded rows 3:6, including preprocessing
# and tuning. Retain both models yourself. No model is trained by this package.
classes <- c("cat", "dog", "bird")
truth <- factor(c("cat", "dog", "cat", "bird", "dog", "bird"), levels = classes)
fold <- c("A", "A", "B", "B", "B", "B")
p <- rbind(c(.7, .2, .1), c(.2, .6, .2), c(.5, .3, .2),
c(.1, .2, .7), c(.3, .5, .2), c(.2, .3, .5))
colnames(p) <- classes
# The same three untouched assessment observations, in the same row order,
# predicted by each corresponding retained model. List order is immaterial.
new_p <- list(B = rbind(c(.6, .3, .1), c(.1, .4, .5), c(.1, .8, .1)),
A = rbind(c(.5, .4, .1), c(.2, .2, .6), c(.2, .7, .1)))
new_p <- lapply(new_p, function(x) { colnames(x) <- classes; x })
assessment <- factor(c("cat", "bird", "dog"), levels = classes)
local({
plot_file <- tempfile(fileext = ".pdf")
grDevices::pdf(plot_file)
on.exit({ grDevices::dev.off(); unlink(plot_file) })
for (method in c("ovo", "ovr")) {
fit <- cross_venn_abers(truth, p, fold, method = method, precision = 2)
print(fit)
for (loss in c("log", "brier")) {
details <- predict(fit, new_p, loss = loss, type = "details")
probability <- predict(fit, new_p, loss = loss)
predicted_class <- predict(fit, new_p, loss = loss, type = "class")
# Reconstruct from the fold endpoints/points, then combine subproblems.
parts <- details[[2L]]
problems <- unique(parts[c(if (method == "ovo") "non_event", "event")])
merged <- vapply(seq_len(nrow(problems)), function(i) {
selected <- parts$event == problems$event[i]
if (method == "ovo") selected <- selected & parts$non_event == problems$non_event[i]
vapply(seq_len(nrow(probability)), function(row) {
f <- parts[selected & parts$row == row, ]
if (loss == "log") {
plogis(mean(log(f$upper)) - mean(log1p(-f$lower)))
} else mean(f$probability)
}, 0)
}, numeric(nrow(probability)))
recovered <- if (method == "ovr") merged else vapply(classes, function(class) {
q <- cbind(1 - merged[, problems$non_event == class, drop = FALSE],
merged[, problems$event == class, drop = FALSE])
1 / (1 + rowSums((1 - q) / q))
}, numeric(nrow(probability)))
recovered <- recovered / rowSums(recovered)
stopifnot(max(abs(recovered - probability)) <= 1e-10,
identical(probability, details$prediction),
identical(predicted_class,
factor(classes[max.col(probability, ties.method = "first")],
levels = classes)))
print(details)
print(predicted_class)
print(c(brier = brier_score(assessment, probability),
log = log_loss(assessment, probability)))
groups <- calibration_plot(assessment, probability, bins = 1)
stopifnot(max(abs(groups$mean_probability - colMeans(probability))) <= 1e-10,
all(groups$count == length(assessment)), all(groups$event_rate == 1 / 3))
print(groups)
}
}
})
# One shared fold plan serves all pairs/classes. Fold A lacks bird but every
# OvO pair remains nonempty. Empty pairs fail with the pair and fold identified.
# Endpoints describe binary subproblems, with no aggregate confidence interval.
# Manual Cross multiclass is an R extension; Python has no matching public
# manual whole-procedure API. Numerical tests reconstruct binary components.One-hot class decisions
one-hot-classes.R: Represent binary and multiclass class
decisions as one-hot matrices, including exact ties and missing
decisions.
# After package installation: Rscript inst/examples/one-hot-classes.R
library(vennabers)
# One-hot matrices represent class decisions the predictor has already made.
# venn_abers_one_hot() never chooses a class from probabilities itself.
# Binary: the logical class decision is TRUE (the event) at probability >= 0.5,
# so an exact 0.5 tie sets the TRUE column.
binary <- venn_abers(c(FALSE, TRUE, FALSE, TRUE), c(.2, .2, .8, .8))
new_scores <- c(.2, .5, .8)
probability <- predict(binary, new_scores)$probability
binary_class <- predict(binary, new_scores, type = "class")
binary_one_hot <- venn_abers_one_hot(binary_class)
print(cbind(probability, binary_one_hot))
stopifnot(probability[[2L]] == .5,
identical(colnames(binary_one_hot), c("FALSE", "TRUE")),
identical(binary_one_hot[, "TRUE"], as.double(binary_class)))
# Multiclass: columns follow the stored factor levels, even for classes that
# no row predicts, and an exact tie keeps the first stored level.
classes <- c("cat", "dog", "bird")
truth <- factor(rep(classes, 2), levels = classes)
p <- rbind(c(.7, .2, .1), c(.2, .6, .2), c(.1, .2, .7),
c(.5, .3, .2), c(.3, .5, .2), c(.2, .3, .5))
colnames(p) <- classes
fit <- venn_abers(truth, p, method = "ovo")
new_p <- rbind(c(.6, .3, .1), c(.1, .4, .5), c(1, 1, 1) / 3)
colnames(new_p) <- classes
multiclass_class <- predict(fit, new_p, type = "class")
multiclass_one_hot <- venn_abers_one_hot(multiclass_class)
print(multiclass_one_hot)
# The third row is an exact tie, so its one sits in the first level "cat".
stopifnot(length(unique(predict(fit, new_p)[3L, ])) == 1L,
identical(multiclass_one_hot[3L, ], c(cat = 1, dog = 0, bird = 0)))
stopifnot(identical(colnames(multiclass_one_hot), classes),
all(rowSums(multiclass_one_hot) == 1),
identical(multiclass_class, factor(classes[max.col(multiclass_one_hot)], levels = classes)))
# Missing decisions give missing rows and keep the row order.
print(venn_abers_one_hot(factor(c("bird", NA, "cat"), levels = classes)))
# Workflow, Grid, parsnip and frozen predictions pass their `.pred_class` table:
# venn_abers_one_hot(predict(fitted_workflow, new_data, type = "class"))
# Schema-encoded tables keep the encoded class names; decode first for originals:
# venn_abers_one_hot(venn_abers_decode(prediction, schema)$class)Inductive regression
inductive-regression.R: Fit an independent linear model,
calibrate numeric predictions and report RMSE/MAE on separate assessment
data.
# After installation: source(system.file("examples", "inductive-regression.R",
# package = "vennabers"))
library(vennabers)
source(system.file("examples", "assertions.R", package = "vennabers"))
# Independent IID observations; roles are fixed before examining their outcomes.
set.seed(20260907)
x <- runif(240, -3, 3)
data <- data.frame(id = seq_along(x), x = x, truth = 2 + x + .4 * x^2 + rnorm(240))
training <- data[1:120, ]
calibration <- data[121:180, ]
assessment <- data[181:240, ]
stopifnot(!anyDuplicated(c(training$id, calibration$id, assessment$id)))
# Only training data determine this model; no preprocessing or tuning uses the
# calibration or assessment observations. The caller retains the fixed model.
model <- lm(truth ~ x, data = training)
calibration_prediction <- predict(model, calibration)
fit <- venn_abers_regression(calibration$truth, calibration_prediction, epsilon = .1)
uncalibrated <- predict(model, assessment)
result <- predict(fit, uncalibrated)
print(list(requested = fit$requested, m = fit$m, epsilon = fit$epsilon))
# Assessment truths enter only here. These are point-prediction scores in target
# units; the VA bounds concern regularized conditional expectations, not coverage
# of individual outcomes. Neither improvement nor a general selector guarantee
# for this particular point merger is promised.
scores <- vapply(list(uncalibrated = uncalibrated, calibrated = result$prediction),
function(prediction) {
c(RMSE = sqrt(mean((assessment$truth - prediction)^2)),
MAE = mean(abs(assessment$truth - prediction)))
}, c(RMSE = 0, MAE = 0))
assert_conditions(all(is.finite(scores)), nrow(result) == nrow(assessment),
all(result$lower <= result$prediction), all(result$prediction <= result$upper),
abs(scores["RMSE", "calibrated"]^2 -
sum((assessment$truth - result$prediction)^2) / nrow(assessment)) < 1e-12,
abs(scores["MAE", "calibrated"] -
sum(abs(assessment$truth - result$prediction)) / nrow(assessment)) < 1e-12)
print(scores)Cross regression
cross-regression.R: Retain matching fold models and
inspect the equal-weight merger of regression points and bounds.
# After installation: source(system.file("examples", "cross-regression.R",
# package = "vennabers"))
library(vennabers)
# IID observations; fold membership and the separate assessment role are fixed
# before examining outcomes. Unequal folds still receive equal aggregate weight.
set.seed(20260908)
x <- runif(240, -3, 3)
data <- data.frame(id = seq_along(x), x = x, truth = 2 + x + .4 * x^2 + rnorm(240))
development <- data[1:180, ]
assessment <- data[181:240, ]
fold <- sample(rep(c("A", "B", "C"), c(50, 60, 70)))
ids <- unique(fold)
models <- setNames(vector("list", length(ids)), ids)
out_of_fold <- rep(NA_real_, nrow(development))
newdata <- matrix(NA_real_, nrow(assessment), length(ids), dimnames = list(NULL, ids))
for (id in ids) {
held_out <- fold == id
training <- development[!held_out, ]
calibration <- development[held_out, ]
stopifnot(!any(training$id %in% c(calibration$id, assessment$id)),
!any(calibration$id %in% assessment$id))
# No learned preprocessing or tuning outside this training set. The caller
# retains this exact model for both its calibration and assessment predictions.
models[[id]] <- lm(truth ~ x, data = training)
out_of_fold[held_out] <- predict(models[[id]], calibration)
newdata[, id] <- predict(models[[id]], assessment)
}
fit <- cross_venn_abers_regression(development$truth, out_of_fold, fold, epsilon = .1)
details <- predict(fit, newdata, type = "details")
print(lapply(fit$calibrators, function(x) x[c("observations", "requested", "m", "epsilon")]))
result <- details$prediction
for (column in names(result)) {
reconstructed <- Reduce(`+`, lapply(details$folds, `[[`, column)) / length(ids)
stopifnot(all(abs(result[[column]] - reconstructed) <= 1e-10 + 1e-10 * abs(reconstructed)))
}
# Matching names cannot prove model provenance or semantic row identity. Each
# column above uses the same assessment rows in the same order. A model refitted
# on all development data cannot replace these fold models. Predict trains nothing.
# Assessment truths enter only in scoring. Aggregated bounds are descriptive
# fold summaries, with no additional coverage or selector guarantee; improvement
# in point-prediction scores is not promised.
scores <- vapply(list(uncalibrated = rowMeans(newdata), calibrated = result$prediction),
function(prediction) {
c(RMSE = sqrt(mean((assessment$truth - prediction)^2)),
MAE = mean(abs(assessment$truth - prediction)))
}, c(RMSE = 0, MAE = 0))
stopifnot(all(is.finite(scores)), nrow(result) == nrow(assessment),
all(result$lower <= result$prediction), all(result$prediction <= result$upper))
for (column in colnames(scores)) {
residual <- assessment$truth - if (column == "calibrated") result$prediction else rowMeans(newdata)
stopifnot(abs(scores["RMSE", column]^2 - sum(residual^2) / length(residual)) < 1e-12,
abs(scores["MAE", column] - sum(abs(residual)) / length(residual)) < 1e-12)
}
print(scores)Fit diagnoses and assessment graphics
calibration-graphics.R: Draw the stored calibration map
and Venn-Abers width of binary, Cross and regression fits, then assess
held-out pairs with a seeded CORP reliability diagram, its Brier
decomposition and a regression plot.
# After installation: source(system.file("examples", "calibration-graphics.R",
# package = "vennabers"))
library(vennabers)
source(system.file("examples", "assertions.R", package = "vennabers"))
# Fit diagnoses draw the stored state of a calibrator and need no new truth
# values. Assessment graphics need held-out pairs of truth and prediction.
set.seed(20261021)
observations <- data.frame(x = rnorm(900L))
observations$truth <- rbinom(900L, 1L, stats::plogis(0.3 + 1.5 * observations$x))
training <- observations[1:300, ]
calibration <- observations[301:600, ]
assessment <- observations[601:900, ]
model <- stats::glm(truth ~ x, data = training, family = stats::binomial())
calibration_probability <- predict(model, calibration, type = "response")
assessment_probability <- predict(model, assessment, type = "response")
calibrator <- venn_abers(calibration$truth, calibration_probability)
probability <- predict(calibrator, assessment_probability)$probability
# A temporary PDF also works without an interactive display. Close and remove
# it even on error. For interactive plots, call the functions directly.
local({
plot_file <- tempfile(fileext = ".pdf")
on.exit(unlink(plot_file), add = TRUE)
grDevices::pdf(plot_file)
on.exit(grDevices::dev.off(), add = TRUE, after = FALSE)
# Fit diagnosis: the calibration map and the Venn-Abers width upper - lower
# at every stored knot. The width is neither a confidence level nor a
# nominal coverage.
map <- plot(calibrator, which = c("map", "width"))
stored <- predict(calibrator, map$score)
assert_conditions(
identical(names(map), c("score", "lower", "upper", "probability", "width")),
identical(map$width, map$upper - map$lower),
identical(map$probability, stored$probability),
identical(range(map$score), c(0, 1))
)
print(head(map))
# Cross fits overlay the map of every fold calibrator; the table keeps one
# block per fold, not the aggregated Cross probability.
fold <- rep(c("A", "B", "C"), length.out = nrow(calibration))
cross <- cross_venn_abers(calibration$truth, calibration_probability, fold)
fold_maps <- plot(cross, loss = "brier")
stopifnot(identical(unique(fold_maps$fold), c("A", "B", "C")))
# Assessment graphic: the CORP reliability diagram of the held-out pairs.
# The default consistency band shows where the estimate of perfectly
# calibrated probabilities would fall; a seed makes it reproducible without
# changing the caller's random stream.
state <- .Random.seed
corp <- corp_plot(assessment$truth, probability, seed = 1)
assert_conditions(
identical(.Random.seed, state),
identical(corp, corp_plot(assessment$truth, probability, seed = 1)),
identical(corp$probability, sort(unique(probability))),
sum(corp$count) == nrow(assessment),
all(diff(corp$estimate) >= 0),
all(corp$band_lower >= 0 & corp$band_upper <= 1)
)
print(head(corp))
# The decomposition splits the Brier score of the same pairs into
# miscalibration, discrimination and uncertainty.
parts <- corp_decomposition_plot(assessment$truth, probability)
assert_conditions(
abs(parts$brier - brier_score(assessment$truth, probability)) < 1e-12,
abs(parts$brier - (parts$mcb - parts$dsc + parts$unc)) < 1e-12
)
print(parts)
# Regression: the map over the uncalibrated prediction with reference lines
# at the moderated truths a and b, then held-out pairs against the
# calibrated points.
x <- runif(600L, 0, 10)
regression_data <- data.frame(x = x, truth = x + sin(x) + rnorm(600L, sd = 0.5))
regression_model <- stats::lm(truth ~ x, data = regression_data[1:200, ])
uncalibrated <- predict(regression_model, regression_data)
regression <- venn_abers_regression(regression_data$truth[201:400],
uncalibrated[201:400])
regression_map <- plot(regression, which = c("map", "width"))
stopifnot(identical(regression_map$width,
regression_map$upper - regression_map$lower))
point <- predict(regression, uncalibrated[401:600])$prediction
held_out_truth <- regression_data$truth[401:600]
held_out_pairs <- regression_plot(held_out_truth, point,
which = c("observed", "residual"))
stopifnot(identical(held_out_pairs$residual, held_out_truth - point))
})
# CORP bands are pointwise resampling bands for the assessment pairs given.
# They are not Venn-Abers lower and upper probabilities, and neither graphic
# makes a statement about Venn-Abers validity.