31  Bagging and Random Forests in rtemis

rtemis fits random forests with setup_Ranger(), using the ranger package (Wright and Ziegler 2017); random forests in LightGBM with setup_LightRF(); and bagged linear additive trees with setup_LINADForest(). When train() is called without hyperparameters, it fits a ranger random forest. This chapter fits forests to the three datasets, examines the out-of-bag estimate, the number of features per split, and variable importance, and returns to the grouped resampling of the Parkinson’s data.

db <- read_diabetes()
hf <- read_heart_failure()
cv <- setup_KFold(n_resamples = 10L, seed = 2026L)

31.1 A random forest

mod_rf <- train(
  db,
  hyperparameters = setup_Ranger(
    num_trees = 500L,
    importance = "permutation",
    seed = 2026L
  ),
  execution_config = setup_SerialExecution(seed = 2026L),
  question = diabetes_question,
  verbosity = 0L
)

seed fixes the random sampling of cases and features inside ranger. The fitted ranger object reports the out-of-bag mean squared error and \(R^2\):

c(
  oob_mse = mod_rf@model[["prediction.error"]],
  oob_rsq = mod_rf@model[["r.squared"]]
)
  oob_mse   oob_rsq 
3261.2501    0.4513 
mod_rf_cv <- train(
  db,
  hyperparameters = setup_Ranger(num_trees = 500L, seed = 2026L),
  outer_resampling_config = cv,
  execution_config = setup_SerialExecution(seed = 2026L),
  question = diabetes_question,
  verbosity = 0L
)

The OOB \(R^2\), computed from a single forest fitted to all 442 cases, is close to the 10-fold cross-validated \(R^2\) of 0.449, which required ten forests.

31.2 Variable importance

importance = "permutation" computes the permutation importance of Section 30.5 on the OOB cases:

plot_varimp(mod_rf)

The ranking is close to that of the linear model’s standardized coefficients (Section 15.1), with s5 and BMI first.

31.3 Features per split

rf_mtry <- rbindlist(lapply(c(1L, 3L, 5L, 10L), function(m) {
  mod <- train(
    db,
    hyperparameters = setup_Ranger(mtry = m, seed = 2026L),
    outer_resampling_config = cv,
    execution_config = setup_SerialExecution(seed = 2026L),
    verbosity = 0L
  )
  data.table(mtry = m, rsq_test = mod@metrics_test@mean_metrics[["rsq"]])
}))
rf_mtry
    mtry rsq_test
   <int>    <num>
1:     1   0.4398
2:     3   0.4489
3:     5   0.4419
4:    10   0.4242

With mtry = 10, all features are candidates at every split, and the forest is bagging of trees. Restricting the candidates to a few features per split gives a modestly higher test \(R^2\) in these data.

31.4 Classification

For the heart failure data, a probability forest with a minimum leaf size of 20 cases:

mod_rf_hf <- train(
  hf,
  hyperparameters = setup_Ranger(
    probability = TRUE,
    min_node_size = 20L,
    seed = 2026L
  ),
  outer_resampling_config = cv,
  execution_config = setup_SerialExecution(seed = 2026L),
  question = hf_question
)
<Resampled Classification Model>
Ranger (Ranger random forest)
⟳ Tested using 10 independent folds.

  <Resampled Classification Training Metrics>
  Aggregate Confusion Matrix across resamples.
                     Predicted
          Reference  Yes  No    
                Yes  182   304
                 No    2  2122

  Showing mean (sd) across resamples.
        Sensitivity: 0.374 (0.044)
        Specificity: 0.999 (2e-03)
  Balanced Accuracy: 0.687 (0.022)
                PPV: 0.990 (0.020)
                NPV: 0.875 (0.007)
                 F1: 0.542 (0.045)
           Accuracy: 0.883 (0.007)
                AUC: 0.978 (4.2e-03)
        Brier Score: 0.076 (2.1e-03)

  <Resampled Classification Test Metrics>
  Aggregate Confusion Matrix across resamples.
                     Predicted
          Reference  Yes  No   
                Yes    6   48
                 No    5  231

  Showing mean (sd) across resamples.
        Sensitivity: 0.113 (0.136)
        Specificity: 0.979 (0.036)
  Balanced Accuracy: 0.546 (0.068)
                PPV: NA (NA)
                NPV: 0.829 (0.028)
                 F1: NA (NA)
           Accuracy: 0.818 (0.041)
                AUC: 0.745 (0.118)
        Brier Score: 0.136 (0.020)

The mean test AUC is 0.745 and the Brier score 0.136, the best of the models fitted to these data so far. With 54 deaths among 290 patients, larger leaves average over more patients, and the estimated probabilities are less variable than those of leaves of one or a few patients.

31.5 Other forests

setup_LightRF() grows a random forest with LightGBM’s histogram-based tree learner. setup_LINADForest() bags linear additive trees, with each tree selecting its size on its out-of-bag cases:

forests_db <- rbindlist(lapply(
  list(
    Ranger = setup_Ranger(seed = 2026L),
    LightRF = setup_LightRF(),
    LINADForest = setup_LINADForest()
  ),
  function(hp) {
    mod <- train(
      db,
      hyperparameters = hp,
      outer_resampling_config = cv,
      execution_config = setup_SerialExecution(seed = 2026L),
      verbosity = 0L
    )
    data.table(rsq_test = mod@metrics_test@mean_metrics[["rsq"]])
  }
), idcol = "forest")
forests_db
        forest rsq_test
        <char>    <num>
1:      Ranger   0.4489
2:     LightRF   0.3948
3: LINADForest   0.4892

The forest of linear additive trees, whose leaves carry linear models, has the highest \(R^2\) of the three on the diabetes data, in which several effects are close to linear.

31.6 Grouped resampling revisited

Section 11.6 compared random and grouped assignment of the Parkinson’s recordings with a linear model. The same comparison with a random forest:

pk <- read_parkinsons()
pk_designs <- list(
  random_recordings = list(
    dat = pk[, !"subject"],
    resampler = setup_KFold(n_resamples = 5L, seed = 2026L)
  ),
  grouped_by_patient = list(
    dat = pk,
    resampler = setup_KFold(n_resamples = 5L, id_strat = "subject", seed = 2026L)
  )
)
pk_rf <- vapply(pk_designs, function(d) {
  mod <- train(
    d[["dat"]],
    hyperparameters = setup_Ranger(num_trees = 300L, seed = 2026L),
    outer_resampling_config = d[["resampler"]],
    execution_config = setup_SerialExecution(seed = 2026L),
    verbosity = 0L
  )
  mod@metrics_test@mean_metrics[["rsq"]]
}, numeric(1))
pk_rf
 random_recordings grouped_by_patient 
            0.7557            -0.0629 

With random assignment of recordings, the forest reaches a test \(R^2\) of 0.756, far above the linear model’s 0.163 under the same design. With assignment by patient, the forest’s \(R^2\) is -0.063. The forest predicts the UPDRS of recordings from patients represented in its training data far better than that of new patients, and only the second is the question. In these data, the more flexible model is the one whose performance the random design overstates more.