1. はじめに

機械学習では、lightGBMなどの決定木の勾配ブースティングモデルは、チューニングすべきハイパーパラメータが多いため、ベイズ最適化が良く用いられます。
以前に、Tidymodelを用いたベイズ最適化によるlightGBMのハイパーパラメータチューニングをしました。
【R】tydymodelsを用いてlightgbmをやってみる（その３）

しかしながら、現在、上記はエラーが発生しています。原因は,

'reshape' argument is no longer supported.

とでるので、どうやらパッケージ側の問題みたいです。

そこで、今回はMLR3を用いたいと思います。
MLR3とは、Tidymodelと同じくRによる機械学習フレームワークですが、文法が異なるので、好みが分かれるところですが、実は私はこちらの文法が気に入っています。

しかしながら、マイナーな存在であるので、この機会に紹介してみたいと思います。
mlr3を使ってみる（その１）
MLR3によるベイズ最適化は、以前に投稿したところです。
mlr3MBOを用いたベイズ最適化について（Bayesian Optimization）
上記は、ハイパーパラメータのチューニングではないため、今回、mlr3MBOをハイパーパラメータチューニングに適用してみます。

環境

R version 4.4.2 (2024-10-31) -- "Pile of Leaves"
RStudio 2024.12.0 Build 467
OS:Ubuntu 22.04

2.DATA

ライブラリ⁠

library(AmesHousing)
library(janitor)
library(tidyverse)
library(mlr3)
library(mlr3learners)
library(mlr3extralearners)  #lightGBMが含まれる
library(mlr3tuning)
library(mlr3pipelines)
library(mlr3mbo)
library(mlr3tuning)
library(bbotk)
library(lightgbm)

AmesHoisingパッケージで入手可能な住宅価格のデータを使用します。

# set the random seed so we can reproduce any simulated results.
set.seed(123)
# load the housing data and clean names
ames_data = make_ames() %>% clean_names()

janitorのclean_names()関数は、自動でデータ列名をクリーニングしてくれる。どのような作業が実行されたというと、
Before

"MS_SubClass"        "MS_Zoning"          "Lot_Frontage"     "Kitchen_AbvGr"

After

 "ms_sub_class"       "ms_zoning"          "lot_frontage"     "kitchen_abv_gr"

大文字を小文字で統一や、アンダーバーで区切りを自動で入れてくれたりします。地味に便利。

2.タスクの設定

多くのモデルでは、正確な予測を行うために、慎重かつ広範な変数の前処理を必要とします。
XGBoost、lightgbm、catboostのような決定木の勾配ブースティングツリーモデルは、高度に歪んだデータや相関性のあるデータに対して非常に頑健であるため、必要な前処理の量は最小限で済みます。XGBoostとは対照的に、lightgbmとcatboostはカテゴリー変数（因子）を扱うことができるので、変数をワンホットエンコードする必要はありません。
データをトレーニングデータとテストデータにわけます。
まず、データは、MLR3では最初にタスク化します。

tsk_ames = TaskRegr$new(id="ames",backend = ames_data,target = "sale_price")
tsk_ames

#> <TaskRegr:ames> (2930 x 81)
#> * Target: sale_price
#> * Properties: -
#> * Features (80):
#>   - fct (46): alley, bldg_type, bsmt_cond, bsmt_exposure,
#>     bsmt_fin_type_1, bsmt_fin_type_2, bsmt_qual, central_air,
#>     condition_1, condition_2, electrical, exter_cond, exter_qual,
#>     exterior_1st, exterior_2nd, fence, fireplace_qu, foundation,
#>     functional, garage_cond, garage_finish, garage_qual,
#>     garage_type, heating, heating_qc, house_style, kitchen_qual,
#>     land_contour, land_slope, lot_config, lot_shape, mas_vnr_type,
#>     misc_feature, ms_sub_class, ms_zoning, neighborhood,
#>     overall_cond, overall_qual, paved_drive, pool_qc, roof_matl,
#>     roof_style, sale_condition, sale_type, street, utilities
#>   - int (22): bedroom_abv_gr, enclosed_porch, fireplaces,
#>     first_flr_sf, full_bath, gr_liv_area, half_bath, kitchen_abv_gr,
#>     lot_area, low_qual_fin_sf, misc_val, mo_sold, open_porch_sf,
#>     pool_area, screen_porch, second_flr_sf, three_season_porch,
#>     tot_rms_abv_grd, wood_deck_sf, year_built, year_remod_add,
#>     year_sold
#>   - dbl (12): bsmt_fin_sf_1, bsmt_fin_sf_2, bsmt_full_bath,
#>     bsmt_half_bath, bsmt_unf_sf, garage_area, garage_cars, latitude,
#>     longitude, lot_frontage, mas_vnr_area, total_bsmt_sf

2.学習器の設定

lightGBMの学習器を設定するとともに、ハイパーパラメータのチューニング設定をします。
lightGBMでは主に次のハイパーパラメータがあります。

num_leaves
feature_fraction
num_iterations
min_data_in_leaf
max_depth
learning_rate
min_gain_to_split
bagging_fraction

今回の最適化をねらうターゲットは、num_leaves、min_data_in_leaf 、feature_fraction の３パラメータとし、learning_rate を0.1に、num_iterations を500で固定します。

lgm_lrn = lrn("regr.lightgbm",
             learning_rate = 0.1, 
             num_leaves = to_tune(10,200),
             num_iterations = 500,
             feature_fraction = to_tune(0.5,1),
             min_data_in_leaf = to_tune(1, 100)
)

3.前処理との組み合わせ

ここでは、ハイパーパラメータのチューニングが主目的なので、前処理は最小限とします。
前処理として、factor列は全てラベルエンコーディングし、ターゲットをlog変換します。
ファクター化すると自動的に水準が振られるため、ラベルエンコーディングはas.numeric（数値化）するだけで良いです。pipline処理を用いて実行します。

po("colapply")は、列に対して関数を適用するパイプラインで、as.numeric関数を、factor列にだけに適用させます。学習器とパイプラインを組み合わせて新たな学習器を作成します。

ple_lrn = as_learner(po("colapply", affect_columns=selector_type("factor"), applicator = as.numeric)　%>>% ppl("targettrafo",
  graph = lgm_lrn,
  targetmutate.trafo = function(x) log(x+1),
  targetmutate.inverter = function(x) list(response = exp(x$response)-1)) )
ple_lrn$plot()

4.ベイズ最適化の実行

サロゲートモデル（代理モデル）にガウス過程回帰を設定し、獲得関数を”EI”とし、獲得関数のオプティマイザーにDIRECT法を指定します。

bayesopt_ego = mlr_loop_functions$get("bayesopt_ego")
surrogate = srlrn(lrn("regr.km", covtype = "matern5_2",
  optim.method = "BFGS", control = list(trace = FALSE)))
acq_function = acqf("ei")
acq_optimizer = acqo(opt("nloptr", algorithm = "NLOPT_GN_ORIG_DIRECT"),
  terminator = trm("stagnation", iters = 100, threshold = 1e-5))

ベイズ最適化は上記の設定でほとんどの場合、行けそうです。
tunerにtnr("mbo")で、上記の設定でループさせます。
なお、ハイパーパラメータの場合は、オプティマイザー関数[opt("mbo")]ではなく、チューナー関数[ tnr("mbo"] を使います。tune（）関数でtunerに学習器、タスク、リサンプリング、目的関数、最適化数をセットし、ベイズ最適化のループを回します。

rsmp_cv = rsmp("cv", folds = 5)
msr_rsme = msr("regr.rmse")

tuner = tnr("mbo",
  loop_function = bayesopt_ego ,
  surrogate = surrogate,
  acq_function = acq_function,
  acq_optimizer = acq_optimizer)

instance = tune(
  tuner,
  tsk_ames,
  ple_lrn,
  rsmp_cv,
  msr_rsme,
  20
)

#> ....
#> ....
#> INFO  [20:28:47.947] [mlr3] Finished benchmark
#> INFO  [20:28:47.980] [bbotk] Result of batch 9:
#> INFO  [20:28:47.990] [bbotk]
#> INFO  [20:28:48.029] [bbotk] Finished optimizing after 20 evaluation(s)
#> INFO  [20:28:48.031] [bbotk] Result:

library(mlr3viz)
autoplot(instance)

num_leavesは明らかに小さい方が良い傾向が出ています。
では、最適なパラメータ設定を見ます。

instance$result_learner_param_vals

#> $colapply.applicator
#> function (x, ...)  .Primitive("as.double")
#> 
#> $colapply.affect_columns
#> selector_type("factor")
#> 
#> $regr.lightgbm.objective
#> [1] "regression"
#> 
#> $regr.lightgbm.verbose
#> [1] -1
#> 
#> $regr.lightgbm.learning_rate
#> [1] 0.1
#> 
#> $regr.lightgbm.num_threads
#> [1] 1
#> 
#> $regr.lightgbm.num_iterations
#> [1] 500
#> 
#> $targetmutate.trafo
#> function(x) log(x+1)
#> 
#> $targetmutate.inverter
#> function(x) list(response = exp(x$response)-1)
#> 
#> $regr.lightgbm.num_leaves
#> [1] 15
#> 
#> $regr.lightgbm.min_data_in_leaf
#> [1] 2
#> 
#> $regr.lightgbm.feature_fraction
#> [1] 0.5243773

またnum_iterations 数を500としていましたが、最適化させます。early_stopping（早期停止）を使います。

lgm_lrn$param_set$set_values(
  early_stopping_rounds = 10,
  learning_rate =0.1,
  num_iterations =500,
  num_leaves=15,
  min_data_in_leaf=2,
  feature_fraction=0.5243773
)

set_validate(lgm_lrn ,0.2)
lgm_lrn$train(tsk_ames)
lgm_lrn$internal_tuned_values[[1]]

#> [1] 116

バリデーションを20%に設定し、早期停止を使って、イテレーション数が500回から116回まで下がりました。これで必要以上な学習回数を防いで、時間の短縮とモデルの過学習を防止します。

以下、全パラメータを設定します。ここで注意する点は、早期設定時に用いた early_stoppingやバリデーション設定を元に戻（NULL)しておかないとエラーがでます。

lgm_lrn$param_set$set_values(
  learning_rate =0.1,
  num_iterations =116,
  num_leaves=15,
  min_data_in_leaf=2,
  feature_fraction=0.5243773,
  early_stopping_rounds=NULL
)
lgm_lrn$validate=NULL

前処理と組み合わせます。

ple_lrn = as_learner(po("colapply", affect_columns=selector_type("factor"), applicator = as.numeric)　%>>% ppl("targettrafo",
  graph = lgm_lrn,
  targetmutate.trafo = function(x) log(x+1),
  targetmutate.inverter = function(x) list(response = exp(x$response)-1)) )
ple_lrn 

#> <GraphLearner:colapply.targetmutate.regr.lightgbm.targetinvert>
#> * Model: -
#> * Parameters: colapply.applicator=<function>,
#>   colapply.affect_columns=<Selector>,
#>   regr.lightgbm.objective=regression, regr.lightgbm.verbose=-1,
#>   regr.lightgbm.learning_rate=0.1, regr.lightgbm.num_leaves=15,
#>   regr.lightgbm.num_threads=1, regr.lightgbm.min_data_in_leaf=2,
#>   regr.lightgbm.feature_fraction=0.5244,
#>   regr.lightgbm.num_iterations=116,
#>   regr.lightgbm.early_stopping_rounds=10,
#>   targetmutate.trafo=<function>, targetmutate.inverter=<function>
#> * Validate: NULL
#> * Packages: mlr3, mlr3pipelines, mlr3extralearners, lightgbm
#> * Predict Types:  [response], se, quantiles, distr
#> * Feature Types: logical, integer, numeric, character, factor,
#>   ordered, POSIXct
#> * Properties: featureless, hotstart_backward, hotstart_forward,
#>   importance, internal_tuning, marshal, missings,
#>   selected_features, validation, weights

5.モデルの評価

ベイズ最適化したパラメータとデフォルト設定と比較してモデルの評価をします。

lgm=lrn("regr.lightgbm")
ple_base_lrn = as_learner(po("colapply", affect_columns=selector_type("factor"), applicator = as.numeric)　%>>% ppl("targettrafo",
  graph = lgm,
  targetmutate.trafo = function(x) log(x+1),
  targetmutate.inverter = function(x) list(response = exp(x$response)-1)) )
ple_lrn$id="bays-tune"
ple_base_lrn$id="original"

デフォルトをoriginalとしました。
評価は、全データの5分割交差検証法により実施します。
ベンチマークのグリッド設計をします。学習器だけを比較しています。

design = benchmark_grid(tsk_ames,c(ple_base_lrn,ple_lrn), rsmp_cv)

ベンチマークを実行します。

bmr = benchmark(design)

目的関数により結果を出します。

bmr$aggregate(msr("regr.rmse"))[, .(task_id, learner_id, regr.rmse)]

リサンプリング結果分布を見ます。

autoplot(bmr)

ベイズ最適化でチューニングすることで、RSMEが小さくなりました。

6.Enjoy!