1 R2D3

R2D3 웹사이트에 기계학습에 대한 시각적 소개가 잘 정리되어 있다. 뉴욕과 샌프란시스코 두 도시의 부동산을 분류하는 데이터를 바탕으로 의사결정나무(Decision Tree) 기계학습 알고리즘이 동작하는 방식을 마우스 스크롤을 아래로 내리게 되면 시각적으로 확인할 수 있다.

FAQ 웹사이트에서 CSV 파일을 다운로드 받아 직접 기계학습 예측모형을 제작해보자.

2 xgBoost 예측모형 설명과 이해

2.1 데이터셋

뉴욕, 샌프란시스코 집값데이터를 다운로드 받아 작업을 하는데 xgBoost에서는 분류예측모형을 개발할 때 종속변수가 0,1로 코딩이 되어야 되서 이를 사전에 반영한다.

# 0. 환경설정 -----

library(tidyverse)
library(recipes)
library(skimr)
library(DALEX)
library(xgboost)

# 1. 데이터 -----
## 1.1. 데이터 가져오기
r2d3_df <- read_csv("https://raw.githubusercontent.com/jadeyee/r2d3-part-1-data/master/part_1_data.csv", skip=2) %>% 
    mutate(in_sf = as.integer(as.factor(in_sf))-1)

glimpse(r2d3_df)        
Observations: 492
Variables: 8
$ in_sf          <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,...
$ beds           <dbl> 2, 2, 2, 1, 0, 0, 1, 1, 1, 2, 3, 1, 1, 1, 0, 2,...
$ bath           <dbl> 1.0, 2.0, 2.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1....
$ price          <int> 999000, 2750000, 1350000, 629000, 439000, 43900...
$ year_built     <int> 1960, 2006, 1900, 1903, 1930, 1930, 1920, 1930,...
$ sqft           <int> 1000, 1418, 2150, 500, 500, 500, 500, 900, 900,...
$ price_per_sqft <int> 999, 1939, 628, 1258, 878, 878, 950, 1083, 1083...
$ elevation      <int> 10, 0, 9, 9, 10, 10, 10, 10, 12, 12, 4, 5, 5, 3...

2.2 모형

recipes 팩키지를 이용하여 요리를 한 후에 xgb.DMatrix() 함수로 데이터프레임을 xgBoost 행렬로 변환시킨다. 그리고 나서 예측모형을 생성한다.

## 1.2. X/Y 데이터 분할
x_train_tbl <- r2d3_df %>% select(-in_sf)
y_train_tbl <- r2d3_df %>% select(in_sf)   

## 1.3. X/Y 데이터 요리
rec_obj <- recipe(~ ., data = x_train_tbl) %>%
    # step_scale(all_numeric()) %>%
    prep(stringsAsFactors = FALSE)

x_train_processed_tbl <- bake(rec_obj, x_train_tbl) 

y_train_processed_tbl <- y_train_tbl

xy_train <- bind_cols(y_train_processed_tbl, x_train_processed_tbl)

# 2. 모형 -----
## 2.1. xgBoost 데이터프레임 --> 행렬 변환 -----
model_martix_train <- model.matrix(in_sf ~ . - 1, xy_train)
data_train <- xgb.DMatrix(model_martix_train, label = xy_train$in_sf)

## 2.2. CV 훈련 -----
param <- list(max_depth = 3, eta = 0.5, silent = 1, 
              nrounds =10, objective = "binary:logistic", eval_metric = "auc")

bst.cv <- xgb.cv(param = param, data_train, 
                 nfold = 5, nrounds = 10)
[1] train-auc:0.941450+0.005091 test-auc:0.917379+0.032575 
[2] train-auc:0.961154+0.005281 test-auc:0.933966+0.020488 
[3] train-auc:0.976706+0.002428 test-auc:0.946442+0.017506 
[4] train-auc:0.978945+0.002187 test-auc:0.951946+0.015630 
[5] train-auc:0.982998+0.001408 test-auc:0.954841+0.010607 
[6] train-auc:0.987774+0.001967 test-auc:0.957372+0.008616 
[7] train-auc:0.988662+0.002205 test-auc:0.958932+0.009653 
[8] train-auc:0.991129+0.002379 test-auc:0.959631+0.011810 
[9] train-auc:0.992328+0.002224 test-auc:0.961647+0.008428 
[10]    train-auc:0.993288+0.002132 test-auc:0.961771+0.010259 

## 2.3. xgBoost 모형 -----
xgb_model <- xgb.train(param, data_train, nrounds = 5)
xgb_model
##### xgb.Booster
raw: 3.6 Kb 
call:
  xgb.train(params = param, data = data_train, nrounds = 5)
params (as set within xgb.train):
  max_depth = "3", eta = "0.5", silent = "1", nrounds = "10", objective = "binary:logistic", eval_metric = "auc", silent = "1"
xgb.attributes:
  niter
callbacks:
  cb.print.evaluation(period = print_every_n)
# of features: 7 
niter: 5
nfeatures : 7

2.3 모형 이해와 설명

블랙박스 xgBoost 모형 이해와 설명을 위해서 DALEX 팩키지 explain() 함수를 사용한다. 연속형 변수 예측이 아니라 범주 예측이라 Przemyslaw Biecek(2018-04-28), “How to use DALEX with the xgboost models”을 참조하여 준용한다.

# 3. 모형 설명  -----
predict_logit <- function(model, x) {
    raw_x <- predict(model, x)
    exp(raw_x)/(1 + exp(raw_x))
}

logit <- function(x) exp(x)/(1+exp(x))

explainer_xgb <- explain(xgb_model, 
                         data = model_martix_train, 
                         y = xy_train$in_sf, 
                         predict_function = predict_logit,
                         link = logit,
                         label = "xgboost")

2.3.1 모형 이해와 설명

예측모형 성능을 이해하기 위해서 model_performance() 함수를 사용한다. 경우에 따라 차이가 있기는 하지만, 대체로 xgBoost가 회귀분석이나 분류모형 거의 모든 경우에 가장 성능이 좋게 나오거나 유사한 성능을 보이고 있다.

## 3.1. 모형성능(Performance) -----

mp_xgb <- model_performance(explainer_xgb)
plot(mp_xgb, geom = "boxplot", show_outliers = 3)

2.3.2 중요 변수

variable_importance() 함수로 중요한 변수가 무엇인지 시각적으로 파악한다.

## 3.2. 변수 중요도(Variable Importance)
vd_xgb <- variable_importance(explainer_xgb, type = "raw")
plot(vd_xgb)

2.3.3 반응 변수 연관

variable_response() 함수로 관심있는 변수가 뉴욕, 샌프란시스코 부동산 분류 예측에 어떻게 연관되는지 파악한다.

## 3.3. 변수 반응도
sv_xgb_price  <- variable_response(explainer_xgb, 
                                                variable = "price",
                                                type = "pdp")
plot(sv_xgb_price) +
    scale_x_log10(labels=scales::comma)

2.3.4 특정 관측점 예측

마지막으로 특정 관측점에 영향을 주는 요인이 무엇인지 prediction_breakdown() 함수로 파악한다. 예측에 대한 각 변수별 기여도를 분해(breakdown)하여 자세히 살펴볼 수 있다.

## 3.4. 특정 관측점 Breakdown

target_obs <- model_martix_train[1, , drop = FALSE]

sp_xgb  <- prediction_breakdown(explainer_xgb, 
                                observation = target_obs)

plot(sp_xgb)