Information Value for Feature Selection: Corporate Bankruptcy Prediction from Kaggle 2021 Competition

R Data Science Series

Motivation

Trong post trước chúng ta đã sử dụng AUC để lựa chọn biến số cho mô hình Logistic để đạt được AUC = 0.903. Chúng ta có thể đạt được thứ hạng cao hơn nữa (cao hơn vị trí thứ hai của Team yoshida có AUC trên Test Data = 0.93293) với kĩ thuật lựa chọn biến số dựa trên Information Value (IV). Bạn đọc quan tâm có thể tham khảo thêm về IV tại trang 81 của textbook Credit Risk Scorecards: Developing and Implementing Intelligent Credit Scoring.

Vẫn như ở post trước chúng ta vẫn sử dụng bộ dữ liệu (download tại đây) đã được sử dụng cho cuộc thi Corporate Bankruptcy Prediction 2021 trên Kaggle. Trước hết load bộ dữ liệu rồi thực hiện một số thao tác xử lí sơ bộ ban đầu:

# Clear our R environment: 
rm(list = ls())

# Load tidyverse package: 
library(tidyverse)

# Load data: 
read_csv("F:/data.csv") -> data

# Rename for all columns: 

old_names <- names(data)

old_names %>% str_replace_all("[^a-z|^A-Z]", "") -> new_names

names(data) <- new_names

# Remove NetIncomeFlag column: 
data %>% select(-NetIncomeFlag)-> df

# Set response and predictors: 

response <- "Bankrupt"

predictors <- names(df %>% select(-response))

Do thời hạn cuộc thi đã hết nên chúng ta có thể phân chia dữ liệu theo tỉ lệ 70-30 như đã mô tả và sử dụng như là proxy để so sánh kết quả của những mô hình phân loại mà chúng ta xây dựng với kết quả của những đội đã tham gia cuộc thi với lưu ý rằng Validation Data có số lượng đúng bằng số quan sát ở Test Data:

# Split our data: 

library(caret)

set.seed(1)

id <- createDataPartition(y = df %>% pull(response), p = 0.7, list = FALSE)

# 70% data for training và validation: 
train_valid <- df[id, ] 

# 30% data will be used for evaluating model performance: 

df_test <- df[-id, ] # Test data.  

set.seed(1)

id_new <- createDataPartition(y = train_valid %>% pull(response), p = nrow(df_test) / nrow(train_valid), list = FALSE)

df_train <- train_valid[-id_new, ] # Train data. 

df_valid <- train_valid[id_new, ] # Validation data. 

Using IV for Feature Selection

Nhắc lại rằng WOE Binning là một kĩ thuật thuộc nhóm Data Transformation mà trong đó một biến liên lục được “rời rạc hóa” (thuật ngữ tiếng Anh là Discretization) thành các nhóm (gọi là binning). Việc chia biến liên tục thành các nhóm như vậy được thực hiện với mục đính là tính toán WoE (bằng chứng có trọng số WoE - Weight of Evidence) và giá trị thông tin IV (Information Value) của biến số. Việc lựa chọn biến số cho mô hình phân loại sẽ căn cứ vào giá trị IV của nó (xem trang 81 cuốn Credit Risk Scorecards: Developing and Implementing Intelligent Credit Scoring). Chúng ta có thể sử dụng thư viện scorecard thể thực hiện Binning dữ liệu cho 94 biến số như sau:

# Load scorecard package: 
library(scorecard)

# Binning data: 

bins <- woebin(df_train, y = response)

## [INFO] creating woe binning ... 
## [INFO] Binning on 2729 rows and 95 columns in 00:00:32

df_train_woe <- woebin_ply(df_train, bins = bins) %>% 
  mutate(Bankrupt = case_when(Bankrupt == 1 ~ "Bankrupt", TRUE ~ "NonBankrupt")) %>% 
  mutate(Bankrupt = as.factor(Bankrupt))

## [INFO] converting into woe values ...

df_valid_woe <- woebin_ply(df_valid, bins = bins) %>% 
  mutate(Bankrupt = case_when(Bankrupt == 1 ~ "Bankrupt", TRUE ~ "NonBankrupt")) %>% 
  mutate(Bankrupt = as.factor(Bankrupt))

## [INFO] converting into woe values ...

df_test_woe <- woebin_ply(df_test %>% select(-response), bins = bins) 

## [INFO] converting into woe values ...

# True label from Test Data: 

true_test_labels <- case_when(df_test$Bankrupt == 1 ~ "Bankrupt", TRUE ~ "NonBankrupt") %>% as.factor()

# Total IV for 94 features: 

do.call("bind_rows", lapply(1:length(predictors), function(j) {bins[[j]]})) %>% 
  arrange(-total_iv) %>% 
  mutate(variable_woe = str_c(variable, "_woe")) %>% 
  select(variable, variable_woe, total_iv) %>% 
  as_tibble() %>% 
  filter(!duplicated(variable)) -> df_iv

Dưới đây là IV của 6 features có Total IV cao nhất:

df_iv %>% 
  select(variable, total_iv) %>% 
  head()

## # A tibble: 6 x 2
##   variable                                        total_iv
##   <chr>                                              <dbl>
## 1 RetainedEarningstoTotalAssets                       2.82
## 2 NetprofitbeforetaxPaidincapital                     2.77
## 3 NetIncometoTotalAssets                              2.65
## 4 PerShareNetprofitbeforetaxYuan                      2.62
## 5 InterestExpenseRatio                                2.55
## 6 ROACbeforeinterestanddepreciationbeforeinterest     2.53

Dưới đây là R codes khảo sát AUC trên Validation Data trước một loạt ngưỡng IV được lựa chọn của biến số:

# Function extracts ROC/AUC on valid data: 

library(pROC)

returnROC_AUC_ValidData <- function(predictor_selected) {
  
  f <- as.formula(paste0(response, " ~ ", paste(predictor_selected, collapse = " + ")))
  
  logit <- glm(f, family = "binomial", data = df_train_woe)
  
  prob_pred <- predict(logit, df_valid_woe, type = "response")
  
  my_auc <- roc(df_valid_woe$Bankrupt, prob_pred)$auc %>% as.numeric()
  
  return(my_auc)
  
}

# Set a sequence of thresholds: 

iv_thresholds <- seq(min(df_iv$total_iv), max(df_iv$total_iv), 0.01)

auc_space <- NULL

# AUC by threshold: 

for (j in iv_thresholds) {
  
  df_iv %>% 
    filter(total_iv >= j) %>% 
    pull(variable_woe) -> predictors_for_modelling
  
  returnROC_AUC_ValidData(predictors_for_modelling) -> my_auc
  
  auc_space <- c(auc_space, my_auc)
  
}

tibble(iv_thresholds = iv_thresholds, auc = auc_space) -> df_auc_threshold

# Features create max AUC on validation data: 

auc_max <- df_auc_threshold %>% slice(which.max(auc))

auc_max

## # A tibble: 1 x 2
##   iv_thresholds   auc
##           <dbl> <dbl>
## 1          1.21 0.936

Như vậy khi chọn các biến số mà thỏa mãn IV >= 1.21 thì AUC trên Validation Data sẽ là lớn nhất và bằng 0.936. Chúng ta kì vọng rằng xu hướng này sẽ vẫn đúng cho Test Data. Figure 1 cho thấy biến động của AUC trên Validation data có dạng hình chữ U ngược trong đó điểm màu đỏ tương ứng với ngưỡng IV để AUC trên Validation Data cực đại:

# Features with IV >= 1.21: 

df_iv %>% 
  filter(total_iv >= auc_max$iv_thresholds) %>% 
  pull(variable_woe) -> var_auc_936

df_auc_threshold %>% 
  ggplot(aes(iv_thresholds, auc)) + 
  geom_line(size = 1, color = "blue") + 
  geom_point(data = auc_max, aes(iv_thresholds, auc), color = "red", size = 2) + 
  labs(x = "IV Threshold", 
       y = "AUC on Valid Data", 
       title = "Figure 1: AUC on Validation Data by Information Value Threshold")

Sử dụng các biến số mà IV >= 1.21 (có 35 biến) cho mô hình Logistic:

f936 <- as.formula(paste0(response, " ~ ", paste(var_auc_936, collapse = " + ")))

logit936 <- glm(f936, family = "binomial", data = df_train_woe)

Sử dụng mô hình Logistic này đánh giá hiệu quả dự báo - phân loại trên Test Data:

prob_pred936 <- predict(logit936, df_test_woe, type = "response")

my_auc <- roc(true_test_labels, prob_pred936)$auc %>% as.numeric()

my_auc

## [1] 0.9325832

Như vậy AUC trên Test Data = 0.93258 (Team yoshida có AUC = 0.93203, xếp thứ hai ở bảng Public Score).

Automated Machine Learning

Dựa trên kết quả thu được ở trên chúng ta có thể đạt được thứ hạng cao hơn nữa bằng sử dụng cách tiếp cận Automated Machine Learning với inputs là các biến có IV >= 1.21. Dưới đây là R codes:

# Load h2o package for Automated Machine Learning: 

library(h2o)
h2o.init(nthreads = 20, max_mem_size = "32g")

##  Connection successful!
## 
## R is connected to the H2O cluster: 
##     H2O cluster uptime:         4 hours 33 minutes 
##     H2O cluster timezone:       Asia/Bangkok 
##     H2O data parsing timezone:  UTC 
##     H2O cluster version:        3.32.1.5 
##     H2O cluster version age:    6 months and 15 days !!! 
##     H2O cluster name:           H2O_started_from_R_ADMIN_xfr920 
##     H2O cluster total nodes:    1 
##     H2O cluster total memory:   6.59 GB 
##     H2O cluster total cores:    40 
##     H2O cluster allowed cores:  20 
##     H2O cluster healthy:        TRUE 
##     H2O Connection ip:          localhost 
##     H2O Connection port:        54321 
##     H2O Connection proxy:       NA 
##     H2O Internal Security:      FALSE 
##     H2O API Extensions:         Amazon S3, Algos, AutoML, Core V3, TargetEncoder, Core V4 
##     R Version:                  R version 4.1.2 (2021-11-01)

h2o.no_progress()

# Prepare data: 

train_h2o <- as.h2o(df_train_woe) # Train data. 

valid_h2o <- as.h2o(df_valid_woe) # Validation data. 

test_h2o <- as.h2o(df_test_woe) # Convert test data to h2o frame. 


#===================================
#  Training Auto Machine Learning
#===================================

# Train Auto Machine Learning: 

autoML <- h2o.automl(x = var_auc_936, 
                     y = response, 
                     training_frame = train_h2o, 
                     leaderboard_frame = valid_h2o, 
                     stopping_metric = "AUC", 
                     stopping_rounds = 10, 
                     stopping_tolerance = 0.025, 
                     max_models = 15, 
                     max_runtime_secs = 60*60, 
                     seed = 1, 
                     sort_metric = "AUC")

## 
## 19:50:43.934: AutoML: XGBoost is not available; skipping it.

Sử dụng mô hình có AUC lớn nhất trên Validation Data và đánh giá hiệu quả phân loại - dự báo của mô hình này bằng Test Data:

prob_h20 <- h2o.predict(autoML@leader, test_h2o) %>% 
  as.data.frame() %>% 
  pull(Bankrupt)

roc(true_test_labels, prob_h20)$auc %>% as.numeric()

## [1] 0.937021

Kết quả AUC này trên Test Data này chỉ thua kém không đáng kể so với AUC = 0.93798 của Team xếp thứ nhất (Team KotaShimomura).

Summary

1. Sử dụng IV cho bước lựa chọn biến số (Feature Selection) và chỉ sử dụng mô hình thống kê truyền thống đơn giản là Logistic chúng ta cũng có thể đặt được một kết quả rất khả quan. Thực nghiệm trên bộ dữ liệu Test Data được AUC = 0.93258.

2. Sử dụng danh sách các biến số có IV >= 1.21 và sử dụng Automated Machine Learning có thể đạt được AUC cao hơn nữa trên Test Data. Kết quả thực nhiệm chỉ ra rằng AUC trên Test Data khi sử dụng Automated Machine Learning xấp xỉ với AUC của Team đang đứng vị trí thứ nhất trong cuộc thi Corporate Bankruptcy Prediction 2021.

3. Đây là dữ liệu bất cân bằng rất cao (chỉ có 3.2% các quan sát là Bankrupt) nên có thể cần xem xét đến khả năng sử dụng các giải pháp resampling dữ liệu như SMOTE, upsampling - downsampling để đạt kết quả tốt hơn nữa.

R Environment and OS

sessionInfo()

## R version 4.1.2 (2021-11-01)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 19043)
## 
## Matrix products: default
## 
## locale:
## [1] LC_COLLATE=English_United States.1252 
## [2] LC_CTYPE=English_United States.1252   
## [3] LC_MONETARY=English_United States.1252
## [4] LC_NUMERIC=C                          
## [5] LC_TIME=English_United States.1252    
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] h2o_3.32.1.5    pROC_1.17.0.1   scorecard_0.3.6 caret_6.0-88   
##  [5] lattice_0.20-45 forcats_0.5.1   stringr_1.4.0   dplyr_1.0.7    
##  [9] purrr_0.3.4     readr_2.0.0     tidyr_1.1.3     tibble_3.1.3   
## [13] ggplot2_3.3.5   tidyverse_1.3.1
## 
## loaded via a namespace (and not attached):
##  [1] nlme_3.1-153         bitops_1.0-7         fs_1.5.0            
##  [4] lubridate_1.7.10     bit64_4.0.5          doParallel_1.0.16   
##  [7] httr_1.4.2           tools_4.1.2          backports_1.2.1     
## [10] bslib_0.2.5.1        utf8_1.2.2           R6_2.5.1            
## [13] rpart_4.1-15         DBI_1.1.1            colorspace_2.0-2    
## [16] nnet_7.3-16          withr_2.4.2          tidyselect_1.1.1    
## [19] gridExtra_2.3        bit_4.0.4            compiler_4.1.2      
## [22] cli_3.0.1            rvest_1.0.1          xml2_1.3.2          
## [25] labeling_0.4.2       sass_0.4.0           scales_1.1.1        
## [28] digest_0.6.27        rmarkdown_2.9        pkgconfig_2.0.3     
## [31] htmltools_0.5.1.1    highr_0.9            dbplyr_2.1.1        
## [34] rlang_0.4.11         readxl_1.3.1         rstudioapi_0.13     
## [37] farver_2.1.0         jquerylib_0.1.4      generics_0.1.0      
## [40] jsonlite_1.7.2       vroom_1.5.4          zip_2.2.0           
## [43] ModelMetrics_1.2.2.2 RCurl_1.98-1.4       magrittr_2.0.1      
## [46] Matrix_1.3-4         Rcpp_1.0.7           munsell_0.5.0       
## [49] fansi_0.5.0          lifecycle_1.0.0      stringi_1.7.4       
## [52] yaml_2.2.1           MASS_7.3-54          plyr_1.8.6          
## [55] recipes_0.1.16       grid_4.1.2           parallel_4.1.2      
## [58] crayon_1.4.1         haven_2.4.1          splines_4.1.2       
## [61] hms_1.1.0            knitr_1.33           pillar_1.6.2        
## [64] reshape2_1.4.4       codetools_0.2-18     stats4_4.1.2        
## [67] reprex_2.0.0         glue_1.4.2           evaluate_0.14       
## [70] data.table_1.14.0    modelr_0.1.8         vctrs_0.3.8         
## [73] tzdb_0.1.2           foreach_1.5.1        cellranger_1.1.0    
## [76] gtable_0.3.0         assertthat_0.2.1     openxlsx_4.2.4      
## [79] xfun_0.25            gower_0.2.2          prodlim_2019.11.13  
## [82] broom_0.7.7          class_7.3-19         survival_3.2-13     
## [85] timeDate_3043.102    iterators_1.0.13     lava_1.6.9          
## [88] ellipsis_0.3.2       ipred_0.9-11