45  Conformal Prediction in rtemis

conformal() computes prediction intervals for regression models and prediction sets for classification models. Its configuration is one of setup_SplitConformal(), setup_CVPlus() and setup_CQR(). It returns a PredictionInterval or PredictionSet, and conformal_metrics() measures the coverage and width of the regions against outcomes not used to construct them. This chapter applies the three methods to the diabetes data and constructs prediction sets for the heart failure data.

45.1 Split conformal intervals

The diabetes data are divided into three sets of equal size: training, to fit the model; calibration, to compute the scores; and test, to measure the coverage of the resulting intervals.

db <- read_diabetes()
ex <- setup_SerialExecution(seed = 2026L)
set.seed(2026L)
idx <- sample(nrow(db))
db_training <- db[idx[1:148]]
db_calibration <- db[idx[149:295]]
db_test <- db[idx[296:442]]
mod_glm <- train(
  db_training,
  hyperparameters = setup_GLM(),
  execution_config = ex,
  question = diabetes_question,
  verbosity = 0L
)
region <- conformal(
  mod_glm,
  newdata = features(db_test),
  calibration = db_calibration,
  config = setup_SplitConformal(alpha = 0.1),
  verbosity = 0L
)
region
<PredictionInterval>
split conformal intervals for 147 cases, 90% guaranteed, mean width 188.91

alpha = 0.1 requests intervals with coverage of at least 90%. The interval is the same width for every case, twice the calibrated quantile of the absolute residuals. Coverage and width on the test cases:

conformal_metrics(region, db_test[["progression"]])
  coverage mean_width median_width   n
1   0.9388      188.9        188.9 147

The interval covers 93.9% of the test outcomes, with a width of 188.9 units of progression, against a range of 25 to 346 for the outcome. Valid coverage does not make the interval narrow: the width reflects how much of the outcome the model leaves unexplained.

45.2 Coverage across calibration sets

The coverage on one test set is one realization. Repeating the three-way division 100 times shows its distribution:

coverage_rep <- vapply(seq_len(100L), function(r) {
  set.seed(r)
  i <- sample(nrow(db))
  mod <- train(
    db[i[1:148]],
    hyperparameters = setup_GLM(),
    execution_config = ex,
    verbosity = 0L
  )
  reg <- conformal(
    mod,
    newdata = features(db[i[296:442]]),
    calibration = db[i[149:295]],
    config = setup_SplitConformal(alpha = 0.1),
    verbosity = 0L
  )
  conformal_metrics(reg, db[i[296:442], progression])[["coverage"]]
}, numeric(1))
draw_histogram(coverage_rep, xlab = "Test coverage of 90% split conformal intervals")

The mean coverage over the 100 divisions is 0.910, at least the nominal 0.9, and 90% of the divisions give a coverage between 0.850 and 0.959. With 147 calibration and 147 test cases, an individual division can fall below 0.9.

45.3 CV+

With a model trained by k-fold cross-validation, conformal() uses CV+ by default: the out-of-fold residuals serve as scores, and no cases are set aside for calibration. Here the model is fitted to the training and calibration cases together:

mod_glm_cv <- train(
  rbind(db_training, db_calibration),
  hyperparameters = setup_GLM(),
  outer_resampling_config = setup_KFold(n_resamples = 10L, seed = 2026L),
  execution_config = ex,
  verbosity = 0L
)
region_cvplus <- conformal(
  mod_glm_cv,
  newdata = features(db_test),
  config = setup_CVPlus(alpha = 0.1),
  verbosity = 0L
)
conformal_metrics(region_cvplus, db_test[["progression"]])
  coverage mean_width median_width   n
1   0.9524        192          192 147

45.4 Conformalized quantile regression

setup_CQR() requires a model of the conditional quantiles; with setup_Ranger(quantreg = TRUE), ranger grows a quantile regression forest (Meinshausen 2006):

mod_qrf <- train(
  db_training,
  hyperparameters = setup_Ranger(quantreg = TRUE, seed = 2026L),
  execution_config = ex,
  verbosity = 0L
)
region_cqr <- conformal(
  mod_qrf,
  newdata = features(db_test),
  calibration = db_calibration,
  config = setup_CQR(alpha = 0.1),
  verbosity = 0L
)
conformal_metrics(region_cqr, db_test[["progression"]])
  coverage mean_width median_width   n
1   0.9116      191.4          201 147
draw_histogram(
  region_cqr@width,
  xlab = "Width of 90% CQR interval"
)

The CQR intervals range in width from 96 to 272, following the forest’s estimate of the spread of progression for each patient, while the split conformal intervals of the linear model all have the same width.

45.5 Prediction sets

For the heart failure data, the patients are divided into training (50%), calibration (25%) and test (25%) sets, stratified on the outcome. A logistic regression is fitted, and prediction sets with a target coverage of 90% are constructed with the LAC and APS scores:

hf <- read_heart_failure()
split1 <- resample(hf, setup_StratSub(n_resamples = 1L, train_p = 0.5, seed = 2026L))
hf_training <- hf[split1[[1]]]
hf_rest <- hf[-split1[[1]]]
split2 <- resample(hf_rest, setup_StratSub(n_resamples = 1L, train_p = 0.5, seed = 2027L))
hf_calibration <- hf_rest[split2[[1]]]
hf_test <- hf_rest[-split2[[1]]]
mod_hf <- train(
  hf_training,
  hyperparameters = setup_GLM(),
  execution_config = ex,
  question = hf_question,
  verbosity = 0L
)
sets <- lapply(c("LAC", "APS"), function(score) {
  conformal(
    mod_hf,
    newdata = features(hf_test),
    calibration = hf_calibration,
    config = setup_SplitConformal(alpha = 0.1, score = score, seed = 2026L),
    verbosity = 0L
  )
})
names(sets) <- c("LAC", "APS")
set_metrics <- rbindlist(lapply(sets, conformal_metrics, hf_test[["death"]]), idcol = "score")
set_metrics
    score coverage mean_set_size singleton_rate empty_rate     n
   <char>    <num>         <num>          <num>      <num> <int>
1:    LAC   0.8356         1.096         0.9041     0.0000    73
2:    APS   0.8493         1.247         0.7260     0.0137    73

mean_set_size is the average number of classes per set, singleton_rate the proportion of sets with one class, and empty_rate the proportion of empty sets. With LAC, 9.6% of the test patients receive the set {No, Yes}, for which the model excludes neither outcome at the chosen level; with APS, 26.0%. The observed coverage on these 73 test patients, 0.836 and 0.849, is below 0.9; as Section 45.2 showed, the coverage of a single small test set varies around its expected value, and the guarantee concerns that expected value.