bonsai and treesnip

R lightgbm, xgboost, catboost 관련 패키지

R
ligthtgbm
xgboost
catboost
Author

Don Don

Published

July 13, 2022

Library

tidymodels에서 xgboost의 경우 공식 지원한다. 다만 많이 사용되는 lightgbm, catboost의 경우 공식 지원하지 않는다. 따라서 CRAN에 등록되어있지 않은 treesnip 패키지를 설치해서 사용했어야 했다. 최근에 treesnip 패키지 개발을 멈추고 bonsai 패키지를 새롭게 만들어서 개발 중인 것 같다. bonsai 패키지 사용법에 대해 한번 정리를 해본다. bonsai 패키지의 경우 lightgbm은 지원하지만 catboost는 지원하지 않기 때문에, treesnip 패키지를 활용해서 catboost 모형을 적합하는 방법에 대해 정리해본다. bonsai, treesnip 패키지의 경우 tidymodels 워크플로우에 해당 모델이 추가된 형태이므로 tidymodels의 워크플로우를 그대로 따라간다.

#install.packages("pak")
#install.packages("lightgbm")
#pak::pak("tidymodels/bonsai")
#remotes::install_github("curso-r/treesnip@catboost")
#remotes::install_github("curso-r/treesnip")


library(tidymodels)
── Attaching packages ────────────────────────────────────── tidymodels 0.2.0 ──
✔ broom        1.0.0     ✔ recipes      1.0.1
✔ dials        1.0.0     ✔ rsample      1.0.0
✔ dplyr        1.0.9     ✔ tibble       3.1.7
✔ ggplot2      3.3.6     ✔ tidyr        1.2.0
✔ infer        1.0.2     ✔ tune         0.2.0
✔ modeldata    1.0.0     ✔ workflows    1.0.0
✔ parsnip      1.0.0     ✔ workflowsets 0.2.1
✔ purrr        0.3.4     ✔ yardstick    1.0.0
── Conflicts ───────────────────────────────────────── tidymodels_conflicts() ──
✖ purrr::discard() masks scales::discard()
✖ dplyr::filter()  masks stats::filter()
✖ dplyr::lag()     masks stats::lag()
✖ recipes::step()  masks stats::step()
• Use suppressPackageStartupMessages() to eliminate package startup messages
library(rsample)
library(bonsai)
library(skimr)
library(treesnip)
Registered S3 method overwritten by 'treesnip':
  method                     from  
  multi_predict._lgb.Booster bonsai

Attaching package: 'treesnip'
The following objects are masked from 'package:bonsai':

    predict_lightgbm_classification_class,
    predict_lightgbm_classification_prob,
    predict_lightgbm_classification_raw,
    predict_lightgbm_regression_numeric, train_lightgbm
library(catboost)

Data description

dat <- MASS::Boston
dat %>% glimpse()
Rows: 506
Columns: 14
$ crim    <dbl> 0.00632, 0.02731, 0.02729, 0.03237, 0.06905, 0.02985, 0.08829,…
$ zn      <dbl> 18.0, 0.0, 0.0, 0.0, 0.0, 0.0, 12.5, 12.5, 12.5, 12.5, 12.5, 1…
$ indus   <dbl> 2.31, 7.07, 7.07, 2.18, 2.18, 2.18, 7.87, 7.87, 7.87, 7.87, 7.…
$ chas    <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
$ nox     <dbl> 0.538, 0.469, 0.469, 0.458, 0.458, 0.458, 0.524, 0.524, 0.524,…
$ rm      <dbl> 6.575, 6.421, 7.185, 6.998, 7.147, 6.430, 6.012, 6.172, 5.631,…
$ age     <dbl> 65.2, 78.9, 61.1, 45.8, 54.2, 58.7, 66.6, 96.1, 100.0, 85.9, 9…
$ dis     <dbl> 4.0900, 4.9671, 4.9671, 6.0622, 6.0622, 6.0622, 5.5605, 5.9505…
$ rad     <int> 1, 2, 2, 3, 3, 3, 5, 5, 5, 5, 5, 5, 5, 4, 4, 4, 4, 4, 4, 4, 4,…
$ tax     <dbl> 296, 242, 242, 222, 222, 222, 311, 311, 311, 311, 311, 311, 31…
$ ptratio <dbl> 15.3, 17.8, 17.8, 18.7, 18.7, 18.7, 15.2, 15.2, 15.2, 15.2, 15…
$ black   <dbl> 396.90, 396.90, 392.83, 394.63, 396.90, 394.12, 395.60, 396.90…
$ lstat   <dbl> 4.98, 9.14, 4.03, 2.94, 5.33, 5.21, 12.43, 19.15, 29.93, 17.10…
$ medv    <dbl> 24.0, 21.6, 34.7, 33.4, 36.2, 28.7, 22.9, 27.1, 16.5, 18.9, 15…
skimr::skim(dat)
Data summary
Name dat
Number of rows 506
Number of columns 14
_______________________
Column type frequency:
numeric 14
________________________
Group variables None

Variable type: numeric

skim_variable n_missing complete_rate mean sd p0 p25 p50 p75 p100 hist
crim 0 1 3.61 8.60 0.01 0.08 0.26 3.68 88.98 ▇▁▁▁▁
zn 0 1 11.36 23.32 0.00 0.00 0.00 12.50 100.00 ▇▁▁▁▁
indus 0 1 11.14 6.86 0.46 5.19 9.69 18.10 27.74 ▇▆▁▇▁
chas 0 1 0.07 0.25 0.00 0.00 0.00 0.00 1.00 ▇▁▁▁▁
nox 0 1 0.55 0.12 0.38 0.45 0.54 0.62 0.87 ▇▇▆▅▁
rm 0 1 6.28 0.70 3.56 5.89 6.21 6.62 8.78 ▁▂▇▂▁
age 0 1 68.57 28.15 2.90 45.02 77.50 94.07 100.00 ▂▂▂▃▇
dis 0 1 3.80 2.11 1.13 2.10 3.21 5.19 12.13 ▇▅▂▁▁
rad 0 1 9.55 8.71 1.00 4.00 5.00 24.00 24.00 ▇▂▁▁▃
tax 0 1 408.24 168.54 187.00 279.00 330.00 666.00 711.00 ▇▇▃▁▇
ptratio 0 1 18.46 2.16 12.60 17.40 19.05 20.20 22.00 ▁▃▅▅▇
black 0 1 356.67 91.29 0.32 375.38 391.44 396.22 396.90 ▁▁▁▁▇
lstat 0 1 12.65 7.14 1.73 6.95 11.36 16.96 37.97 ▇▇▅▂▁
medv 0 1 22.53 9.20 5.00 17.02 21.20 25.00 50.00 ▂▇▅▁▁

train/test split

splits <- initial_split(dat, prop = 0.7)

data preprocessing

rec <- recipe(medv ~ ., training(splits))

Catboost

cat_mod <- boost_tree(
        trees = 1000,
        min_n = tune(),
        tree_depth = tune(),
        learn_rate = tune()) %>%
  set_engine(engine = "catboost" 
             ) %>%
  set_mode(mode = "regression")
cat_wf <- workflow() %>% 
    add_recipe(rec) %>% 
    add_model(cat_mod)
set.seed(1)
vb_folds <- vfold_cv(training(splits), v = 5)
grid <- grid_latin_hypercube(min_n(), tree_depth(), learn_rate(),  size = 10)
res <- tune_grid(
    cat_wf, 
    resamples = vb_folds, 
    grid = grid, 
    metrics = yardstick::metric_set(rmse, rsq, mae), 
    control = control_grid(save_pred = TRUE))
show_best(res)
Warning: No value of `metric` was given; metric 'rmse' will be used.
# A tibble: 5 × 9
  min_n tree_depth learn_rate .metric .estimator  mean     n std_err .config    
  <int>      <int>      <dbl> <chr>   <chr>      <dbl> <int>   <dbl> <chr>      
1    27          8 0.0703     rmse    standard    3.25     5   0.231 Preprocess…
2    31         11 0.00798    rmse    standard    4.21     5   0.403 Preprocess…
3    22          5 0.00125    rmse    standard    7.19     5   0.512 Preprocess…
4     7          9 0.0000988  rmse    standard    8.85     5   0.483 Preprocess…
5    35         12 0.00000379 rmse    standard    9.03     5   0.479 Preprocess…
best_param <- select_best(res, 'rmse')
final_model <- finalize_workflow(cat_wf, best_param)
final_model
══ Workflow ════════════════════════════════════════════════════════════════════
Preprocessor: Recipe
Model: boost_tree()

── Preprocessor ────────────────────────────────────────────────────────────────
0 Recipe Steps

── Model ───────────────────────────────────────────────────────────────────────
Boosted Tree Model Specification (regression)

Main Arguments:
  trees = 1000
  min_n = 27
  tree_depth = 8
  learn_rate = 0.0702985501331971

Computational engine: catboost 
cat_fit <- fit(final_model, data = training(splits))
pred_cat <- 
    predict(cat_fit, testing(splits)) %>% 
    mutate(modelo = "cat")


result3 <- data.frame(pred = pred_cat$.pred, 
                     real = testing(splits)$medv)
yardstick::metrics(result3, pred, real)
# A tibble: 3 × 3
  .metric .estimator .estimate
  <chr>   <chr>          <dbl>
1 rmse    standard       4.06 
2 rsq     standard       0.816
3 mae     standard       2.36

lightgbm

tidymodels에서 lightgbm을 사용할 경우 set_engine()에 ligthtgbm을 명시해주면 된다. cuda gpu를 사용할 경우 params=list(tree_method = 'gpu_hist')를 명시해줘야 한다. 또한 CRAN에 등록되어있는 lightgbm 패키지가 아닌 gpu 옵션을 지원하는 lightgbm을 따로 설치해줘야한다(참고 : https://cran.r-project.org/web/packages/lightgbm/readme/README.html#installing-a-gpu-enabled-build). 이 경우 Docker file을 다운받아서 설치하는 것이 오류가 날 가능성이 적다.

bt_mod <- 
  boost_tree(
        trees = 1000,
        min_n = tune(),
        tree_depth = tune(),
        learn_rate = tune(),
        loss_reduction = tune(),
        stop_iter = tune()) %>%
  set_engine(engine = "lightgbm" 
             #params=list(tree_method = 'gpu_hist')
             ) %>%
  set_mode(mode = "regression")

R lightgbm 공식 문서를 보면 추가적인 parameter가 많은데, 추가하고 싶은 parameter는 set_engine에 해당 parameter를 추가하면 된다.

light_wf <- workflow() %>% 
    add_recipe(rec) %>% 
    add_model(bt_mod)

워크플로우에 전처리 recipe와 정의한 모델을 업데이트한다.

vb_folds <- vfold_cv(training(splits), v = 5)

validation set을 지정한다. fold=5인 cross validation 방식으로 설정한다.

grid <- grid_latin_hypercube(min_n(), tree_depth(), learn_rate(), loss_reduction(), stop_iter(), size = 10)

tuning 방식을 정의하는데 latine hypercube 방식을 이용했다. tuning 방식은 tidymodels 공식문서를 참고하면 다양한 방법을 제공한다.

library(doFuture)
Loading required package: foreach

Attaching package: 'foreach'
The following objects are masked from 'package:purrr':

    accumulate, when
Loading required package: future
registerDoFuture()
cl <- parallel::makeCluster(6)
plan(cluster, workers = cl)

병렬 처리 옵션의 경우 doFuture 패키지를 이용했다. 해당 컴퓨터의 코어 수를 확인하고 적절하게 세팅해주면 된다.

res <- tune_grid(
    light_wf, 
    resamples = vb_folds, 
    grid = grid, 
    metrics = yardstick::metric_set(rmse, rsq, mae), 
    control = control_grid(save_pred = TRUE, allow_par = T))

tune_grid()를 이용해서, workflow, validation set, metric, 예측값 저장 옵션 등을 설정하고, 파라미터 튜닝을 진행한다.

show_best(res)
Warning: No value of `metric` was given; metric 'rmse' will be used.
# A tibble: 5 × 11
  min_n tree_depth learn_rate loss_reduction stop_iter .metric .estimator  mean
  <int>      <int>      <dbl>          <dbl>     <int> <chr>   <chr>      <dbl>
1    34          5  0.0158     0.000183              5 rmse    standard    3.55
2    27          7  0.00210    0.00000126           16 rmse    standard    4.16
3    18          3  0.000239   0.00000000667        10 rmse    standard    7.69
4    17         14  0.0000300  0.0000274             7 rmse    standard    8.86
5     6         13  0.0000199  1.05                  4 rmse    standard    8.91
# … with 3 more variables: n <int>, std_err <dbl>, .config <chr>

show_best()를 통해 tuning 후 최적의 파라미터를 확인할 수 있다.

best_param <- select_best(res, 'rmse')
final_model <- finalize_workflow(light_wf, best_param)
final_model
══ Workflow ════════════════════════════════════════════════════════════════════
Preprocessor: Recipe
Model: boost_tree()

── Preprocessor ────────────────────────────────────────────────────────────────
0 Recipe Steps

── Model ───────────────────────────────────────────────────────────────────────
Boosted Tree Model Specification (regression)

Main Arguments:
  trees = 1000
  min_n = 34
  tree_depth = 5
  learn_rate = 0.0158299029849725
  loss_reduction = 0.000183367712483135
  stop_iter = 5

Computational engine: lightgbm

최적의 파라미터를 선택한 후 finalize_workflow()를 선택된 최적의 파라미터를 workflow에 업데이트 시켜준다.

light_fit <- fit(final_model, data = training(splits))

최적의 파라미터를 업데이트한 모델을 훈련 데이터를 이용해서 학습을 진행해준다.

pred_light <- 
    predict(light_fit, testing(splits)) %>% 
    mutate(modelo = "lightgbm")

pred_light %>% head()
# A tibble: 6 × 2
  .pred modelo  
  <dbl> <chr>   
1  30.8 lightgbm
2  31.8 lightgbm
3  19.5 lightgbm
4  20.2 lightgbm
5  18.2 lightgbm
6  18.6 lightgbm
result1 <- data.frame(pred = pred_light$.pred, 
                     real = testing(splits)$medv)

마지막으로 test 데이터셋에 학습을 해서 예측값을 산출한다.

yardstick::metrics(result1, pred, real)
# A tibble: 3 × 3
  .metric .estimator .estimate
  <chr>   <chr>          <dbl>
1 rmse    standard       4.51 
2 rsq     standard       0.774
3 mae     standard       2.82

yardstick 패키지의 metrics()의 경우 회귀, 분류 세팅에 맞춰서 적절한 평가지표를 보여준다. metrics()의 경우 예측값과 실제값이 포함된 데이터셋을 넣어줘야 작동한다.

XGBOOST

xgboost의 경우 tidymodels에서 공식 지원하기 때문에, tidymodels 워크플로우에 맞춰서 그대로 사용하면 된다. 대신 gpu를 사용할 경우 lightgbm과 마찬가지로 gpu가 지원되는 xgboost 버전을 따로 설치해줘야 한다(참고 : https://github.com/dmlc/xgboost/issues/6654#issuecomment-813828671). lightgbm과 동일하게 gpu가 지원되는 xgboost를 설치했을 경우 params=list(tree_method = 'gpu_hist') 옵션을 넣어주면 된다.

xgb_mod <- boost_tree(
        trees = 1000,
        min_n = tune(),
        tree_depth = tune(),
        learn_rate = tune(),
        loss_reduction = tune(),
        stop_iter = tune()) %>%
  set_engine(engine = "xgboost" 
             #params=list(tree_method = 'gpu_hist')
             ) %>%
  set_mode(mode = "regression")
xgb_wf <- workflow() %>% 
    add_recipe(rec) %>% 
    add_model(xgb_mod)
set.seed(1)
vb_folds <- vfold_cv(training(splits), v = 5)
grid <- grid_latin_hypercube(min_n(), tree_depth(), learn_rate(), loss_reduction(), stop_iter(), size = 10)
library(doFuture)
registerDoFuture()
cl <- parallel::makeCluster(6)
plan(cluster, workers = cl)


res <- tune_grid(
    xgb_wf, 
    resamples = vb_folds, 
    grid = grid, 
    metrics = yardstick::metric_set(rmse, rsq, mae), 
    control = control_grid(save_pred = TRUE, allow_par = T))
! Fold1: internal: A correlation computation is required, but `estimate` is const...
! Fold2: internal: A correlation computation is required, but `estimate` is const...
! Fold3: internal: A correlation computation is required, but `estimate` is const...
! Fold4: internal: A correlation computation is required, but `estimate` is const...
! Fold5: internal: A correlation computation is required, but `estimate` is const...
show_best(res)
Warning: No value of `metric` was given; metric 'rmse' will be used.
# A tibble: 5 × 11
  min_n tree_depth learn_rate loss_reduction stop_iter .metric .estimator  mean
  <int>      <int>      <dbl>          <dbl>     <int> <chr>   <chr>      <dbl>
1    12          6  0.0894           6.42e-2        17 rmse    standard    3.21
2    23         11  0.00647          2.09e-3        10 rmse    standard    3.75
3     8          1  0.000643         4.19e-9         9 rmse    standard   13.7 
4    28          4  0.0000392        6.95e-6         8 rmse    standard   22.9 
5    35         11  0.0000159        1.25e+1         4 rmse    standard   23.4 
# … with 3 more variables: n <int>, std_err <dbl>, .config <chr>
best_param <- select_best(res, 'rmse')
final_model <- finalize_workflow(xgb_wf, best_param)
final_model
══ Workflow ════════════════════════════════════════════════════════════════════
Preprocessor: Recipe
Model: boost_tree()

── Preprocessor ────────────────────────────────────────────────────────────────
0 Recipe Steps

── Model ───────────────────────────────────────────────────────────────────────
Boosted Tree Model Specification (regression)

Main Arguments:
  trees = 1000
  min_n = 12
  tree_depth = 6
  learn_rate = 0.0894462983351361
  loss_reduction = 0.0642365279282407
  stop_iter = 17

Computational engine: xgboost 
xgb_fit <- fit(final_model, data = training(splits))
pred_xgb <- 
    predict(xgb_fit, testing(splits)) %>% 
    mutate(modelo = "xgb")

pred_xgb %>% head()
# A tibble: 6 × 2
  .pred modelo
  <dbl> <chr> 
1  28.5 xgb   
2  35.5 xgb   
3  22.0 xgb   
4  21.7 xgb   
5  19.3 xgb   
6  18.9 xgb   
result2 <- data.frame(pred = pred_xgb$.pred, 
                     real = testing(splits)$medv)
yardstick::metrics(result2, pred, real)
# A tibble: 3 × 3
  .metric .estimator .estimate
  <chr>   <chr>          <dbl>
1 rmse    standard       4.37 
2 rsq     standard       0.795
3 mae     standard       2.57

Citation

BibTeX citation:
@online{don2022,
  author = {Don, Don and Don, Don},
  title = {Bonsai and Treesnip},
  date = {2022-07-13},
  url = {https://dondonkim.netlify.app/posts/2022-07-12-bonsai/bonsai.html},
  langid = {en}
}
For attribution, please cite this work as:
Don, Don, and Don Don. 2022. “Bonsai and Treesnip.” July 13, 2022. https://dondonkim.netlify.app/posts/2022-07-12-bonsai/bonsai.html.