33  Boosting and Bayesian Additive Regression Trees in rtemis

rtemis fits gradient boosting with setup_LightGBM(), using LightGBM (Ke et al. 2017); Bayesian additive regression trees with setup_BART(), using the stochtree package; and rule ensembles with setup_LightRuleFit(), which extracts rules from LightGBM trees and fits a lasso with glmnet. This chapter examines the learning curve of boosting, the selection of the number of rounds, tuning, BART, and the rules of a rule ensemble.

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

33.1 The learning curve of boosting

To show the training and validation loss at every round, a model is fitted with a fixed 1000 rounds (force_nrounds), learning rate 0.05, and a validation set of 20% of the diabetes cases:

split <- resample(db, setup_StratSub(n_resamples = 1L, train_p = 0.8, seed = 2026L))
mod_lc <- train(
  db[split[[1]]],
  dat_validation = db[-split[[1]]],
  hyperparameters = setup_LightGBM(
    force_nrounds = 1000L,
    learning_rate = 0.05
  ),
  execution_config = setup_SerialExecution(seed = 2026L),
  verbosity = 0L
)
plot_learning(mod_lc)

get_learning_curve() returns the losses as data and plot_learning() draws them. The training loss, the mean squared error, falls at every round. The validation loss falls to its minimum at round 43 and then rises as the model fits the noise of the training cases (Section 32.3).

33.2 Selecting the number of rounds

By default, setup_LightGBM() leaves the number of rounds to be selected: train() runs internal resampling of the training data, stops each fit early when its validation loss has not improved for early_stopping_rounds rounds, and refits on all training cases with the selected number of rounds. Within 10-fold cross-validation, the selection is repeated in each outer training set:

mod_lgb <- train(
  db,
  hyperparameters = setup_LightGBM(),
  outer_resampling_config = cv,
  execution_config = setup_SerialExecution(seed = 2026L),
  question = diabetes_question
)

<Resampled Regression Model>
LightGBM (LightGBM gradient boosting)
⚙ Tuned using exhaustive grid search.
⟳ Tested using 10 independent folds.

  <Resampled Regression Training Metrics>
  Showing mean (sd) across resamples.
   MAE: 37.306 (1.128)
   MSE: 2038.346 (102.690)
  RMSE: 45.135 (1.139)
    R²: 0.656 (0.018)

  <Resampled Regression Test Metrics>
  Showing mean (sd) across resamples.
   MAE: 46.242 (3.170)
   MSE: 3249.394 (379.175)
  RMSE: 56.909 (3.450)
    R²: 0.449 (0.059)

vapply(mod_lgb@models, function(m) m@hyperparameters@hyperparameters[["nrounds"]], numeric(1))

 Fold_1  Fold_2  Fold_3  Fold_4  Fold_5  Fold_6  Fold_7  Fold_8  Fold_9 Fold_10 
    263     300     234     276     273     252     230     221     247     225 

The default learning rate of setup_LightGBM() is 0.01, with up to 1000 rounds and trees of at most 8 leaves.

33.3 Tuning

The learning rate and the tree size are tuned together, with the number of rounds selected by early stopping for each combination:

mod_lgb_tuned <- train(
  db,
  hyperparameters = setup_LightGBM(
    learning_rate = tune_over(0.01, 0.05, 0.1),
    num_leaves = tune_over(4L, 8L, 16L)
  ),
  tuner_config = setup_GridSearch(
    resampler_config = setup_KFold(n_resamples = 5L, seed = 2027L)
  ),
  outer_resampling_config = cv,
  execution_config = setup_SerialExecution(seed = 2026L),
  verbosity = 0L
)
rbindlist(lapply(mod_lgb_tuned@models, function(m) {
  as.list(m@hyperparameters@hyperparameters[c("learning_rate", "num_leaves", "nrounds")])
}), idcol = "fold")

       fold learning_rate num_leaves nrounds
     <char>         <num>      <int>   <int>
 1:  Fold_1          0.05          4      82
 2:  Fold_2          0.10          4      52
 3:  Fold_3          0.05          4      83
 4:  Fold_4          0.05         16      52
 5:  Fold_5          0.10          4      53
 6:  Fold_6          0.05          4      80
 7:  Fold_7          0.05          4      69
 8:  Fold_8          0.05          4      70
 9:  Fold_9          0.10          4      49
10: Fold_10          0.05          8      61

33.4 Bayesian additive regression trees

setup_BART() uses the priors of Chipman et al. (2010) by default: 200 trees, \(\alpha = 0.95\), \(\beta = 2\). The sampler starts from 5 iterations of the grow-from-root algorithm and retains 100 MCMC samples:

mod_bart <- train(
  db,
  hyperparameters = setup_BART(seed = 2026L),
  outer_resampling_config = cv,
  execution_config = setup_SerialExecution(seed = 2026L),
  question = diabetes_question,
  verbosity = 0L
)

33.5 Rule ensembles

mod_rulefit <- train(
  db,
  hyperparameters = setup_LightRuleFit(),
  execution_config = setup_SerialExecution(seed = 2026L),
  question = diabetes_question,
  verbosity = 0L
)

The boosted trees produced 288 distinct rules, of which the lasso retained 50 with nonzero coefficients. The rules with the largest coefficients:

rules <- mod_rulefit@model@rules_selected_formatted_coefs
rules[order(-abs(Coefficient))][1:8]

   Rule_ID                                                          Rule
     <int>                                                        <char>
1:      12                                    s5 > 4.7833 AND s6 > 101.5
2:      21                     bmi > 27.75 AND bp > 101.5 AND s3 <= 45.5
3:       6                  bmi > 26.25 AND s5 > 4.51635 AND bmi > 32.95
4:      42                                    bp > 101.5 AND s1 <= 180.5
5:      10    s5 > 4.7833 AND s6 <= 103.5 AND sex IN {1} AND s1 <= 199.5
6:      18                                                    bmi > 33.4
7:      30    s5 > 4.51635 AND s2 <= 147.9 AND bmi > 26.85 AND s3 > 40.5
8:      28 bmi <= 33.05 AND s5 <= 4.909 AND s5 <= 4.16665 AND s1 > 161.5
   N_Cases Coefficient
     <num>       <num>
1:      52       23.16
2:      41       22.07
3:      29       17.41
4:      35       16.35
5:      28       16.13
6:      28       14.87
7:      61       14.66
8:      37      -13.14

N_Cases is the number of training cases that satisfy the rule, and Coefficient the change in predicted progression for those cases. A rule lists the conditions along a tree path as extracted, and can contain conditions on the same feature that a simpler rule would combine. Rules overlap, so each coefficient is the contribution of its rule given all the others.

mod_rulefit_cv <- train(
  db,
  hyperparameters = setup_LightRuleFit(),
  outer_resampling_config = cv,
  execution_config = setup_SerialExecution(seed = 2026L),
  verbosity = 0L
)

33.6 Comparison

mod_glm <- train(
  db,
  hyperparameters = setup_GLM(),
  outer_resampling_config = cv,
  execution_config = setup_SerialExecution(seed = 2026L),
  verbosity = 0L
)
compare_db <- rbindlist(lapply(
  list(
    GLM = mod_glm,
    `LightGBM, default` = mod_lgb,
    `LightGBM, tuned` = mod_lgb_tuned,
    BART = mod_bart,
    LightRuleFit = mod_rulefit_cv
  ),
  function(m) {
    data.table(
      rsq = m@metrics_test@mean_metrics[["rsq"]],
      rsq_sd = m@metrics_test@sd_metrics[["rsq"]]
    )
  }
), idcol = "model")
compare_db

               model    rsq  rsq_sd
              <char>  <num>   <num>
1:               GLM 0.4935 0.08864
2: LightGBM, default 0.4487 0.05877
3:   LightGBM, tuned 0.4542 0.06918
4:              BART 0.4997 0.08151
5:      LightRuleFit 0.4460 0.05594

For the heart failure data:

compare_hf <- rbindlist(lapply(
  list(
    LightGBM = setup_LightGBM(),
    BART = setup_BART(seed = 2026L)
  ),
  function(hp) {
    mod <- train(
      hf,
      hyperparameters = hp,
      outer_resampling_config = cv,
      execution_config = setup_SerialExecution(seed = 2026L),
      verbosity = 0L
    )
    data.table(mod@metrics_test@mean_metrics[c("auc", "brier_score")])
  }
), idcol = "model")
compare_hf

      model    auc brier_score
     <char>  <num>       <num>
1: LightGBM 0.6977      0.1438
2:     BART 0.7092      0.1390

On the diabetes data, BART reaches a mean test \(R^2\) of 0.500, against 0.493 for the linear model, without tuning. The boosted models perform comparably to the random forest of Section 31.1. On these small datasets, whose structure the additive models of Chapter 23 found to be largely additive and near-linear, the boosted trees do not predict better than the linear and additive models; the differences between all of them are smaller than the variation across folds.

Chipman, Hugh A., Edward I. George, and Robert E. McCulloch. 2010. “BART: Bayesian Additive Regression Trees.” The Annals of Applied Statistics 4 (1): 266–98. https://doi.org/10.1214/09-aoas285.

Ke, Guolin, Qi Meng, Thomas Finley, et al. 2017. “LightGBM: A Highly Efficient Gradient Boosting Decision Tree.” Advances in Neural Information Processing Systems 30: 3146–54.