Package
DALEX는 순열(Permutation)을 이용하여 변수 중요도를 측정하는 함수들을 포함하고 있다. 순열에 기반한 변수 중요도는 모형을 학습한 후 특정 예측변수의 값을Permutation(= 랜덤하게 섞은)한 데이터를 이용하여 측정된 예측 성능과원래 데이터를 이용하여 측정된 예측 성능을 비교하여 해당 예측변수의 영향력을 확인한다. 만약, 특정 예측변수가 중요하다면(영향력이 크다면), 해당 예측변수의 값들을 Permutation한 후에 측정한 예측 성능은 원래 데이터를 이용한 예측 성능보다 크게 감소한다. 반면, 특정 예측변수가 중요하지 않다면, 해당 예측변수의 값들을 Permutation한 후에 측정한 예측 성능은 차이가 거의 없거나 오히려 증가한다. 이러한 변수 중요도 측정 방법은 모형 구조에 어떠한 가정도 하지 않기 때문에어떤 모형이든지쉽게 변수 중요도를 측정하여 비교할 수 있으며, 중요도를 계산하기 위해 모형을 재학습할 필요가 없다는 장점을 가지고 있다.
Package
caret을 이용하여Random Forest와XGBoost를 수행하고 순열에 기반한 변수 중요도를 적용하는 방법을 설명하기 위해 예제 데이터셋 “Heart Disease Prediction”를 사용한다. 이 데이터는 환자의 심장병을 예측하는 데이터(출처 : Package MLDataR, Gary Hutson 2021)로 총 918명의 환자에 대한 자료이며, 변수는 총 10개이다. 여기서 Target은HeartDisease이다.
1. 데이터 불러오기
pacman::p_load("MLDataR", # For Data
"data.table", "magrittr", "dplyr",
"caret", # For ML
"DALEX", # For VIM
"doParallel", "parallel")
registerDoParallel(cores=detectCores())
data(heartdisease)
data <- heartdisease
cols <- c("Sex", "RestingECG", "Angina", "HeartDisease")
data <- data %>%
mutate_at(cols, as.factor)
glimpse(data)
Rows: 918
Columns: 10
$ Age <dbl> 40, 49, 37, 48, 54, 39, 45, 54, 37, 48, 37,~
$ Sex <fct> M, F, M, F, M, M, F, M, M, F, F, M, M, M, F~
$ RestingBP <dbl> 140, 160, 130, 138, 150, 120, 130, 110, 140~
$ Cholesterol <dbl> 289, 180, 283, 214, 195, 339, 237, 208, 207~
$ FastingBS <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0~
$ RestingECG <fct> Normal, Normal, ST, Normal, Normal, Normal,~
$ MaxHR <dbl> 172, 156, 98, 108, 122, 170, 170, 142, 130,~
$ Angina <fct> N, N, N, Y, N, N, N, N, Y, N, N, Y, N, Y, N~
$ HeartPeakReading <dbl> 0.0, 1.0, 0.0, 1.5, 0.0, 0.0, 0.0, 0.0, 1.5~
$ HeartDisease <fct> 0, 1, 0, 1, 0, 0, 0, 0, 1, 0, 0, 1, 0, 1, 0~2. 데이터 분할
set.seed(100)
trainIndex <- createDataPartition(data$HeartDisease,
p = .8,
list = FALSE)
HD.train <- data[ trainIndex,]
HD.test <- data[-trainIndex,]
3. 머신러닝 분석
- 예제 데이터에 머신러닝 방법을 적용하기 위해 Package
caret을 이용한다.
3-1. Random Forest
set.seed(100)
caret.rf <- caret::train(HeartDisease~., data = HD.train,
method = "parRF",
trControl = fitControl,
tuneLength = 2,
importance = TRUE)
caret.rf
Parallel Random Forest
735 samples
9 predictor
2 classes: '0', '1'
No pre-processing
Resampling: Adaptively Cross-Validated (5 fold, repeated 5 times)
Summary of sample sizes: 588, 588, 588, 588, 588, 588, ...
Resampling results across tuning parameters:
mtry Accuracy Kappa Resamples
3 0.8247816 0.6459073 25
6 0.8217687 0.6391614 5
Accuracy was used to select the optimal model using the
largest value.
The final value used for the model was mtry = 3.caret.rf$finalModel
Call:
randomForest(x = "x", y = "y", ntree = 63, mtry = 3L, importance = TRUE)
Type of random forest: classification
Number of trees: 504
No. of variables tried at each split: 33-2. XGBoost
set.seed(100)
caret.xgb <- caret::train(HeartDisease~., data = HD.train,
method = "xgbTree",
trControl = fitControl,
tuneLength = 2)
o 1 model was skunkedcaret.xgb
eXtreme Gradient Boosting
735 samples
9 predictor
2 classes: '0', '1'
No pre-processing
Resampling: Adaptively Cross-Validated (5 fold, repeated 5 times)
Summary of sample sizes: 588, 588, 588, 588, 588, 588, ...
Resampling results across tuning parameters:
eta max_depth gamma colsample_bytree min_child_weight
0.2228220 9 1.702621 0.6528662 3
0.4876292 3 5.465586 0.5499986 5
subsample nrounds Accuracy Kappa Resamples
0.7517663 624 0.7891156 0.5733217 5
0.8219133 358 0.8168720 0.6295470 25
Accuracy was used to select the optimal model using the
largest value.
The final values used for the model were nrounds = 358, max_depth
= 3, eta = 0.4876292, gamma = 5.465586, colsample_bytree
= 0.5499986, min_child_weight = 5 and subsample = 0.8219133.caret.xgb$finalModel
##### xgb.Booster
raw: 186.8 Kb
call:
xgboost::xgb.train(params = list(eta = param$eta, max_depth = param$max_depth,
gamma = param$gamma, colsample_bytree = param$colsample_bytree,
min_child_weight = param$min_child_weight, subsample = param$subsample),
data = x, nrounds = param$nrounds, objective = "binary:logistic")
params (as set within xgb.train):
eta = "0.48762916797353", max_depth = "3", gamma = "5.46558595029637", colsample_bytree = "0.549998590815812", min_child_weight = "5", subsample = "0.82191331172362", objective = "binary:logistic", validate_parameters = "TRUE"
xgb.attributes:
niter
callbacks:
cb.print.evaluation(period = print_every_n)
# of features: 10
niter: 358
nfeatures : 10
xNames : Age SexM RestingBP Cholesterol FastingBS RestingECGNormal RestingECGST MaxHR AnginaY HeartPeakReading
problemType : Classification
tuneValue :
nrounds max_depth eta gamma colsample_bytree
2 358 3 0.4876292 5.465586 0.5499986
min_child_weight subsample
2 5 0.8219133
obsLevels : 0 1
param :
list()4. 순열 중요도
- Package
DALEX를 이용하면 순열에 기반한 변수 중요도의 결과를 그래프로 나타내어 어떤 변수가 중요한지 이해하기 쉬우며, 여러 모형의 변수 중요도 결과를 하나의 그래프로 나타낼 수 있어 모형 간에 변수 중요도 결과를 쉽게 비교할 수 있다.
4-1. 함수 explain()
- 순열에 기반한 변수 중요도를 수행하기 전에 먼저 함수
explain()을 이용하여 모형에 대한 설명(Model Explainer)을 정의한다. - 함수
explain()의 자세한 옵션은 여기를 참조한다.
explain(model, data, y, label, type)
model: 학습 모형data: Data Frame 또는 행렬로, 순열에 기반한 변수 중요도 척도값을 계산하기 위해 사용되는 데이터Target은 제외
y:data의Target에 해당하는 수치형 벡터label: 모형의 이름type: 모형의 타입(classification / regression)
Exp_rf <- DALEX::explain(model = caret.rf,
data = HD.test[, -10], # Data except for Target
y = as.numeric(HD.test$HeartDisease), # Target : Numeric Type
label = "RF")
Preparation of a new explainer is initiated
-> model label : RF
-> data : 183 rows 9 cols
-> data : tibble converted into a data.frame
-> target variable : 183 values
-> predict function : yhat.train will be used ( [33m default [39m )
-> predicted values : No value for predict function target column. ( [33m default [39m )
-> model_info : package caret , ver. 6.0.88 , task classification ( [33m default [39m )
-> predicted values : numerical, min = 0 , mean = 0.5448868 , max = 1
-> residual function : difference between y and yhat ( [33m default [39m )
-> residuals : numerical, min = 0.04761905 , mean = 1.007026 , max = 1.904762
[32m A new explainer has been created! [39m Exp_xgb <- DALEX::explain(model = caret.xgb,
data = HD.test[, -10], # Data except for Target
y = as.numeric(HD.test$HeartDisease), # Target : Numeric
label = "XGBoost")
Preparation of a new explainer is initiated
-> model label : XGBoost
-> data : 183 rows 9 cols
-> data : tibble converted into a data.frame
-> target variable : 183 values
-> predict function : yhat.train will be used ( [33m default [39m )
-> predicted values : No value for predict function target column. ( [33m default [39m )
-> model_info : package caret , ver. 6.0.88 , task classification ( [33m default [39m )
-> predicted values : numerical, min = 0.02965254 , mean = 0.5608565 , max = 0.9867967
-> residual function : difference between y and yhat ( [33m default [39m )
-> residuals : numerical, min = 0.08915342 , mean = 0.9910561 , max = 1.821626
[32m A new explainer has been created! [39m 4-2. 함수 model_parts()
- 함수
model_parts()를 이용하여 순열에 기반한 변수 중요도값을 계산할 수 있으며, 알고리즘은 다음과 같다.- 특정 예측 변수에 대해, 원래 데이터를 이용한 손실 함수값과 원래 데이터를 Permutation한 후 계산한 손실 함수값을 비교하여 특정 예측 변수의 중요도값을 얻는다.
- 함수
model_parts()의 자세한 옵션은 여기를 참조한다.
model_parts(explainer, loss_function, type, variables, B, N)
explainer: 함수explain()를 이용하여 정의한 객체loss_function: 손실 함수loss_cross_entropyloss_accuracyloss_one_minus_aucloss_root_mean_squareloss_sum_of_squares
type: 변수 중요도 척도 종류raw: 특정 예측변수의 값들을 Permutation한 후 계산한 손실 함수값difference: 특정 예측변수의 값들을 Permutation한 후 계산한 손실 함수값 \(-\) 원래 데이터를 이용하여 계산한 손실 함수값ratio: 특정 예측변수의 값들을 Permutation한 후 계산한 손실 함수값 \(/\) 원래 데이터를 이용하여 계산한 손실 함수값
variables: 변수 중요도를 계산하고자 하는 예측변수B: 변수 중요도값을 계산하기 위해 사용하는 순열의 수N: 변수 중요도값을 계산하기 위해 사용하는 데이터의 개수- 기본값은 1,000이며,
N = NULL을 하면 모든 데이터 사용
- 기본값은 1,000이며,
Caution! 함수 model_parts()의 Github에 있는 함수 feature_importance()의 자세한 구조는 여기를 참조한다. 함수 model_parts()의 옵션 type에 기반하여 계산된 변수 중요도 척도값과 그래프에 대해 자세히 설명해보자.
4-2-1. type = raw
set.seed(100)
vip.B_rf <- model_parts(explainer = Exp_rf,
# loss_function = loss_cross_entropy,
type = "raw", # raw / ratio / difference
B = 5, # B : # of permutations
N = NULL)
vip.B_rf
variable mean_dropout_loss label
1 _full_model_ 0.1281092 RF
2 RestingBP 0.1224583 RF
3 RestingECG 0.1296788 RF
4 Age 0.1313692 RF
5 Sex 0.1397006 RF
6 MaxHR 0.1489495 RF
7 FastingBS 0.1519923 RF
8 HeartPeakReading 0.1854142 RF
9 Cholesterol 0.1887708 RF
10 Angina 0.1906303 RF
11 _baseline_ 0.4968365 RFset.seed(100)
vip.B_xgb <- model_parts(explainer = Exp_xgb,
# loss_function = loss_cross_entropy,
type = "raw", # raw / ratio / difference
B = 5, # B : # of permutations
N = NULL)
vip.B_xgb
variable mean_dropout_loss label
1 _full_model_ 0.09876841 XGBoost
2 RestingBP 0.09869597 XGBoost
3 RestingECG 0.09876841 XGBoost
4 Age 0.10026564 XGBoost
5 Sex 0.11509297 XGBoost
6 FastingBS 0.11632456 XGBoost
7 MaxHR 0.13170732 XGBoost
8 HeartPeakReading 0.15667713 XGBoost
9 Cholesterol 0.16544313 XGBoost
10 Angina 0.18676648 XGBoost
11 _baseline_ 0.51137406 XGBoostCaution! Option B는 “순열의 반복의 수”이다. 즉, B = 5이면 특정 예측변수 값들을 5번 Permutation하여 5개의 손실 함수값이 계산되고 이들의 평균값을 대표값으로 사용한다. Option N은 손실 함수값을 계산하기 위해 사용할 Test Data point의 개수를 지정하는데 사용한다. 예를 들어 N = NULL이면 Test Data의 모든 데이터를 이용하여 손실 함수값을 계산하고, N = 100이면 Test Data에서 랜덤하게 선택한 100개를 이용하여 손실 함수값을 계산한다. 따라서, Option N에 특정값을 지정하면 손실 함수값을 계산할 때마다 Test Data에서 랜덤하게 그 수만큼 데이터를 선택하여 손실 함수값을 계산한다.
Result! mean_dropout_loss은 평균 손실 함수값을 의미한다. “full_model”의 mean_dropout_loss은 Permutation없이 Test Data를 이용하여 계산한 평균 손실 함수값을 의미하며, “baseline”의 mean_dropout_loss은 Test Data의 값들을 Permutation한 후 계산한 평균 손실 함수값이다.
plot(vip.B_rf, vip.B_xgb) +
ggtitle("Mean variable-importance over 5 permutations", "")
Caution! 함수 plot()은 Package ggplot을 이용하여 만들어진 함수이므로 ggplot의 다양한 옵션을 이용하여 그래프를 수정할 수 있다.
Result! X축은 평균 손실 함수값이다. 막대(Bar)의 점선 수직선은 “full_model”의 mean_dropout_loss이며, 막대의 끝은 해당 예측변수 값을 Permutation한 후 계산된 평균 손실 함수값이다. 따라서, 막대의 길이는 해당 예측변수 값을 Permutation한 후 계산된 평균 손실 함수값과 “full_model”의 mean_dropout_loss의 차이를 나타낸다. 점선 수직선을 기준으로 막대가 오른쪽으로 길수록 그 차이가 크기 때문에 해당 예측변수가 그만큼 중요하다는 것을 나타낸다.
4-2-2. type = difference
set.seed(100)
vip.B_rf <- model_parts(explainer = Exp_rf,
# loss_function = loss_cross_entropy,
type = "difference", # raw / ratio / difference
B = 5, # B : # of permutations
N = NULL)
vip.B_rf
variable mean_dropout_loss label
1 _full_model_ 0.000000000 RF
2 RestingBP -0.005650809 RF
3 RestingECG 0.001569669 RF
4 Age 0.003260082 RF
5 Sex 0.011591403 RF
6 MaxHR 0.020840377 RF
7 FastingBS 0.023883120 RF
8 HeartPeakReading 0.057304999 RF
9 Cholesterol 0.060661676 RF
10 Angina 0.062521130 RF
11 _baseline_ 0.368727361 RFset.seed(100)
vip.B_xgb <- model_parts(explainer = Exp_xgb,
# loss_function = loss_cross_entropy,
type = "difference", # raw / ratio / difference
B = 5, # B : # of permutations
N = NULL)
vip.B_xgb
variable mean_dropout_loss label
1 _full_model_ 0.000000e+00 XGBoost
2 RestingBP -7.244627e-05 XGBoost
3 RestingECG 0.000000e+00 XGBoost
4 Age 1.497223e-03 XGBoost
5 Sex 1.632456e-02 XGBoost
6 FastingBS 1.755615e-02 XGBoost
7 MaxHR 3.293890e-02 XGBoost
8 HeartPeakReading 5.790872e-02 XGBoost
9 Cholesterol 6.667472e-02 XGBoost
10 Angina 8.799807e-02 XGBoost
11 _baseline_ 4.126057e-01 XGBoostResult! mean_dropout_loss은 해당 예측변수 값을 Permutation한 후 계산된 평균 손실 함수값과 원래 Test Data의 Data Point를 통해 계산된 평균 손실 함수값의 차이를 의미한다. “full_model”의 mean_dropout_loss은 같은 값을 빼기 때문에 “0”이 된다.
Result! mean_dropout_loss은 Permutation없이 Test Data를 이용하여 계산한 평균 손실 함수값을 특정 예측변수의 값들을 Permutation한 후 계산한 평균 손실 함수값의 차이를 의미한다. “full_model”의 mean_dropout_loss은 동일한 값을 빼기 때문에 그 차이는 “0”이 된다.
plot(vip.B_rf, vip.B_xgb) +
ggtitle("Mean variable-importance over 5 permutations", "")
Result! 막대(Bar)의 점선 수직선은 “full_model”의 mean_dropout_loss으로 0의 값을 가지며, 막대의 끝은 해당 예측변수의 값들을 Permutation한 후 계산한 평균 손실 함수값과 Permutation없이 계산한 평균 손실 함수값과의 차이이다. 따라서 점선 수직선을 기준으로 오른쪽으로 막대의 길이가 길수록 그 차이가 크다는 것을 의미하며, 이는 해당 예측변수가 중요하다는 것을 나타낸다.
4-2-3. type = ratio
set.seed(100)
vip.B_rf <- model_parts(explainer = Exp_rf,
# loss_function = loss_cross_entropy,
type = "ratio", # raw / ratio / difference
B = 5, # B : # of permutations
N = NULL)
vip.B_rf
variable mean_dropout_loss label
1 _full_model_ 1.0000000 RF
2 RestingBP 0.9558907 RF
3 RestingECG 1.0122526 RF
4 Age 1.0254477 RF
5 Sex 1.0904807 RF
6 MaxHR 1.1626767 RF
7 FastingBS 1.1864279 RF
8 HeartPeakReading 1.4473139 RF
9 Cholesterol 1.4735156 RF
10 Angina 1.4880302 RF
11 _baseline_ 3.8782281 RFset.seed(100)
vip.B_xgb <- model_parts(explainer = Exp_xgb,
# loss_function = loss_cross_entropy,
type = "ratio", # raw / ratio / difference
B = 5, # B : # of permutations
N = NULL)
vip.B_xgb
variable mean_dropout_loss label
1 _full_model_ 1.0000000 XGBoost
2 RestingBP 0.9992665 XGBoost
3 RestingECG 1.0000000 XGBoost
4 Age 1.0151589 XGBoost
5 Sex 1.1652812 XGBoost
6 FastingBS 1.1777506 XGBoost
7 MaxHR 1.3334963 XGBoost
8 HeartPeakReading 1.5863081 XGBoost
9 Cholesterol 1.6750611 XGBoost
10 Angina 1.8909535 XGBoost
11 _baseline_ 5.1775061 XGBoostResult! mean_dropout_loss은 Permutation없이 계산한 평균 손실 함수값을 특정 예측변수의 값들을 Permutation한 후 계산한 평균 손실 함수값의 비율를 의미한다. “full_model”의 mean_dropout_loss은 동일한 값으로 나누기 때문에 그 비율은 “1”이 된다.
plot(vip.B_rf, vip.B_xgb) +
ggtitle("Mean variable-importance over 5 permutations", "")
Result! 막대(Bar)의 점선 수직선은 “full_model”의 mean_dropout_loss으로 1의 값을 가지며, 막대의 끝은 해당 예측변수의 값들을 Permutation한 후 계산한 평균 손실 함수값과 Permutation없이 계산한 평균 손실 함수값과의 비율이다. 따라서 점선 수직선을 기준으로 오른쪽으로 막대의 길이가 길수록 해당 예측변수가 중요하다는 것을 나타낸다.