1 모형 평가1

예측모형 개발에 필요한 데이터셋은 basetable 형태로 단 한벌이지만 예측모형은 서로 다른 아키텍쳐와 특성을 갖는 모형이 다수 만들어진다. 최적의 예측모형을 선택하는 것은 예측모형의 성능에 달려있다.

  • ROC_AUC
  • RMSE
  • \(R^2\)
  • Precision, Recall

상기 모형성능(Performance)의 불활실성을 고려하여 최적의 모형을 선정하고자 하는 경우 빈도주의 통계기법을 가져오는 경우와 베이지안 통계기법을 적용시키는 두가지 방법이 있다.

신뢰구간(Confidence Interval)은 다소 직관적이지 않는 추정값의 불확실성에 대한 해석을 제공하는데, “충분히 많은 실험을 반복하게 되면, 참값의 95%는 구간 [하한값, 상한값] 사이에 떨어진다.” 이와 비교하여 베이지안 방법론은 “참값이 구간 [하한값, 상한값] 사이에 떨어질 확률이 95%가 된다.”라고 직관적인 답을 제공한다. tidyposterior는 재표집(resampling) 통계량을 활용하여 예측모형을 비교하는데 필요한 팩키지다.

모형성능에 대한 추정치에 대한 불활실성(estimates of uncertainty)을 고려한 요약통계량을 구하려고 할 경우 본질적으로 반복(replicates)이 요구되고 Bootstrap, K-fold 교차검증와 같은 재표집 방법(resampling methods) 도입이 요구된다.

1.1 모형 평가 수식

베이지안 계층적 일반화선형모형(Bayesian Hierarchical Generalized Linear Model)으로 예측모형 비교 평가 후 선정 과정을 이론화시킬 수 있는데 기초 통계에서 많이 회자되는 ANOVA 모형이다.

\[\text{정확도} = b_0 + b_1 m_1 + b_2 m_2 + b_3 m_3 + b_4 m_4\]

여기서 \(m_j\)는 로지스틱 회귀모형, Random Forest, MARS, XGBoost 를 나타내는 지시변수다.
하지만 재표집에 따라 예측모형 정확도에 대한 효과가 달라지기 때문에 이를 반영하는 것이 필요한데, K-Fold 교차검증을 가정하면 다음과 같이 수식을 고쳐 작성할 수 있다.

\[\text{정확도}_i = b_{i0} + b_1 m_1 + b_2 m_2 + b_3 m_3 + b_4 m_4\]

여기서 \(i\)\(i^{th}\)번째 교차검증 재표본이 된다.

한걸음 더 들어가 정확도는 다음과 같은 정규분포를 따른다고 가정할 수 있다.

\[\text{정확도}_{ij} \sim \mathcal{N}(\beta_{i0} + \beta_{j}m_{ij},\,\sigma^{2})\,. \]

모수 \(b_{i0}\)는 평균 \(\beta_0\)와 분산을 갖는 정규분포를 갖는다. 따라서 선형혼합모형(linear mixed model)으로 흔히 불리는 Random Intercept Model이 된다. rstan 활용하여 모수를 추정한다.

2 펭귄성별 예측모형

앞서 펭귄 성별예측모형: tidymodels을 참고로 tidymodels를 활용하여 예측하는 방법을 살펴봤다. 해당 블로그 내용을 기초로 서로 다른 4가지 예측모형을 비교평가해 보자.

# 1. ingest data
library(tidyverse)
tuesdata <- tidytuesdayR::tt_load('2020-07-28')

    Downloading file 1 of 2: `penguins.csv`
    Downloading file 2 of 2: `penguins_raw.csv`
penguin <- tuesdata$penguins

# 2. preprocessing
penguin_df <- penguin %>%
  filter(!is.na(sex)) %>%
  select(-year, -island) %>% 
  mutate_if(is.character, as.factor)

# 3. split train test
library(tidymodels)

penguin_split <- initial_split(penguin_df, prop = 0.8, strata = sex)
penguin_train <- training(penguin_split)
penguin_test <- testing(penguin_split)

penguin_cv <- vfold_cv(penguin_train, v =10, repeats = 1) 

# 4. fit model
doParallel::registerDoParallel()

penguin_wf <- workflow() %>%
  add_formula(sex ~ .)

glm_spec <- logistic_reg() %>%
  set_mode("classification") %>% 
  set_engine("glm") 

rf_spec <- rand_forest() %>%
  set_mode("classification") %>%
  set_engine("ranger")

mars_spec <- mars() %>% 
  set_mode("classification") %>% 
  set_engine("earth")

xgboost_spec <- boost_tree() %>% 
  set_mode("classification") %>% 
  set_engine("xgboost")

glm_rs <- penguin_wf  %>% 
  add_model(glm_spec) %>% 
  fit_resamples(
    resamples = penguin_cv,
    control = control_resamples(save_pred = TRUE),
    metrics = metric_set(accuracy, kap, roc_auc)
  )

rf_rs <- penguin_wf  %>% 
  add_model(rf_spec) %>% 
  fit_resamples(
    resamples = penguin_cv,
    control = control_resamples(save_pred = TRUE),
    metrics = metric_set(accuracy, kap, roc_auc)
  )

mars_rs <- penguin_wf  %>% 
  add_model(mars_spec) %>% 
  fit_resamples(
    resamples = penguin_cv,
    control = control_resamples(save_pred = TRUE),
    metrics = metric_set(accuracy, kap, roc_auc)
  )

xgboost_rs <- penguin_wf  %>% 
  add_model(xgboost_spec) %>% 
  fit_resamples(
    resamples = penguin_cv,
    control = control_resamples(save_pred = TRUE),
    metrics = metric_set(accuracy, kap, roc_auc)
  )

3 모형 평가

펭귄 성별 예측을 위한 데이터는 하나지만 모형은 4개를 비교하여 가장 성능좋은 모형을 선정해보자.

3.1 모형평가 데이터프레임

앞서 4가지 모형에 대한 모형성능 측도(Performance Metrics)를 뽑아내서 이를 데이터프레임으로 작성한다.

models_df <- tibble(
  model_name = c("GLM", "Random_Forest", "MARS", "XGBoost"),
  model_output = list(glm_rs, rf_rs, mars_rs, xgboost_rs)
)

collect_metrics_from_model <- function(model) {
  model_metrics <- model %>% 
    select(id, .metrics) %>% 
    unnest()
}

models_perf_df <- models_df %>% 
  mutate(performance = map(model_output, collect_metrics_from_model)) %>% 
  unnest(cols = performance)

models_perf_df
# A tibble: 120 x 6
   model_name model_output      id     .metric  .estimator .estimate
   <chr>      <list>            <chr>  <chr>    <chr>          <dbl>
 1 GLM        <tibble [10 × 5]> Fold01 accuracy binary         0.889
 2 GLM        <tibble [10 × 5]> Fold01 kap      binary         0.777
 3 GLM        <tibble [10 × 5]> Fold01 roc_auc  binary         0.983
 4 GLM        <tibble [10 × 5]> Fold02 accuracy binary         0.889
 5 GLM        <tibble [10 × 5]> Fold02 kap      binary         0.777
 6 GLM        <tibble [10 × 5]> Fold02 roc_auc  binary         0.967
 7 GLM        <tibble [10 × 5]> Fold03 accuracy binary         0.963
 8 GLM        <tibble [10 × 5]> Fold03 kap      binary         0.926
 9 GLM        <tibble [10 × 5]> Fold03 roc_auc  binary         0.956
10 GLM        <tibble [10 × 5]> Fold04 accuracy binary         0.741
# … with 110 more rows

3.2 성능지표별 모형평가

ROC-AUC, \(\kappa\), 정확도를 기준으로 각 모형성능을 전체적으로 평가해보자.

library(extrafont)
loadfonts()

models_perf_df %>% 
  ggplot(aes(x = model_name, y = .estimate, fill = model_name)) + 
    geom_boxplot(show.legend = FALSE) + 
    facet_wrap(~.metric, scales = "free_y") +
    theme_bw(base_family = "NanumGothic") +
    labs(x="", y="성능 추정치", title="펭귄 성별 예측모형별 성능비교")

3.3 표본내 성능평가

ROC-AUC, \(\kappa\), Precision, Recall 등 다양한 성별 예측모형 성능 비교를 위한 지표가 있지만 정확도를 기준으로 살펴보자. 우선 정확도로 필터링한 후에 K-Fold 표본내(within-sample) 각 모형별 성능을 비교해보자.

models_perf_df %>% 
  filter(.metric == "accuracy") %>% 
  select(model_name, id, .estimate)  %>% 
    ggplot(aes(x = model_name, y = .estimate, color = id, group = id)) + 
    geom_line()

4 베이지안 계층 GLM

베이지안 계층 일반환 선형모형을 적용하여 단순하게 예측모형을 비교해보자. 이를 위해서 먼저 각 재표집 K-Fold 교차검증별로 예측모형 성능(여기서는 정확도)을 추출하여 데이터프레임을 준비한다.

models_posterior_df <- models_perf_df %>% 
  filter(.metric == "accuracy") %>% 
  select(model_name, id, .estimate) %>% 
  pivot_wider(names_from = "model_name", values_from = ".estimate")

models_posterior_df
# A tibble: 10 x 5
   id       GLM Random_Forest  MARS XGBoost
   <chr>  <dbl>         <dbl> <dbl>   <dbl>
 1 Fold01 0.889         0.926 0.926   0.926
 2 Fold02 0.889         0.852 0.889   0.889
 3 Fold03 0.963         0.926 0.963   0.889
 4 Fold04 0.741         0.926 0.704   0.852
 5 Fold05 0.963         0.889 0.963   0.889
 6 Fold06 1             0.926 1       0.926
 7 Fold07 0.889         0.963 0.889   0.926
 8 Fold08 0.926         0.889 0.852   0.889
 9 Fold09 1             0.962 1       1    
10 Fold10 0.962         0.962 0.923   0.962

4.1 베이지안 계층모형 적합

tidyposterior 팩키지에서 핵심적인 함수는 perf_mod()contrast_models()이 있다. 먼저 perf_mod()함수를 베이지안 계층 모형 적합을 한 후 예측모형 성능을 비교하여 평가해보자.

library(tidyposterior)

accuracy_model <- perf_mod(models_posterior_df, seed = 13311, verbose = FALSE)

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 1).
Chain 1: 
Chain 1: Gradient evaluation took 0.000128 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 1.28 seconds.
Chain 1: Adjust your expectations accordingly!
Chain 1: 
Chain 1: 
Chain 1: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 1: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 1: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 1: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 1: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 1: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 1: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 1: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 1: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 1: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 1: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 1: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 1: 
Chain 1:  Elapsed Time: 0.326105 seconds (Warm-up)
Chain 1:                0.162988 seconds (Sampling)
Chain 1:                0.489093 seconds (Total)
Chain 1: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 2).
Chain 2: 
Chain 2: Gradient evaluation took 3.1e-05 seconds
Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.31 seconds.
Chain 2: Adjust your expectations accordingly!
Chain 2: 
Chain 2: 
Chain 2: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 2: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 2: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 2: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 2: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 2: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 2: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 2: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 2: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 2: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 2: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 2: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 2: 
Chain 2:  Elapsed Time: 0.322286 seconds (Warm-up)
Chain 2:                0.127975 seconds (Sampling)
Chain 2:                0.450261 seconds (Total)
Chain 2: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 3).
Chain 3: 
Chain 3: Gradient evaluation took 2.6e-05 seconds
Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.26 seconds.
Chain 3: Adjust your expectations accordingly!
Chain 3: 
Chain 3: 
Chain 3: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 3: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 3: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 3: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 3: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 3: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 3: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 3: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 3: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 3: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 3: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 3: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 3: 
Chain 3:  Elapsed Time: 0.364729 seconds (Warm-up)
Chain 3:                0.200406 seconds (Sampling)
Chain 3:                0.565135 seconds (Total)
Chain 3: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 4).
Chain 4: 
Chain 4: Gradient evaluation took 4.7e-05 seconds
Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.47 seconds.
Chain 4: Adjust your expectations accordingly!
Chain 4: 
Chain 4: 
Chain 4: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 4: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 4: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 4: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 4: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 4: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 4: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 4: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 4: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 4: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 4: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 4: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 4: 
Chain 4:  Elapsed Time: 0.338557 seconds (Warm-up)
Chain 4:                0.120493 seconds (Sampling)
Chain 4:                0.45905 seconds (Total)
Chain 4: 
accuracy_model %>% tidy() %>% 
  summary() %>% 
  arrange(desc(mean))
# A tibble: 4 x 4
  model          mean lower upper
  <chr>         <dbl> <dbl> <dbl>
1 Random_Forest 0.921 0.886 0.955
2 GLM           0.921 0.885 0.956
3 XGBoost       0.914 0.879 0.948
4 MARS          0.911 0.875 0.945

4.2 모형평가 시각화

예측모형 정확도 성능에 대해서 시각화를 하게 되면 직관적으로 다수 예측모형을 시각적으로 비교할 수 있다.

accuracy_model %>% tidy() %>% 
  ggplot() +
    geom_boxplot(aes(fill = model), width=0.3, show.legend = FALSE) +
    theme_bw(base_family = "NanumGothic") +
    labs(y="사후 확률") +
    scale_y_continuous(labels = scales::percent)

4.3 모형 쌍대 평가

각 모형간 차이를 살펴보자. 많은 부분 성능차이가 겹치는 것을 확인할 수 있다. 따라서 동일한 성능을 내는 경우 단순한 모형을 펭귄 성별 예측에 사용하는 것이 권장된다.

accuracy_model %>% 
  contrast_models(.) %>% 
  ggplot() +
    theme_bw() +
    geom_vline(xintercept = 0, color = "blue") +
    theme_bw(base_family = "NanumGothic") +
    labs(y="사후 확률", x="모형간 차이")

 

데이터 과학자 이광춘 저작

kwangchun.lee.7@gmail.com