knitr::opts_chunk$set(
echo = TRUE,
message = FALSE,
warning = FALSE,
fig.width = 10,
fig.height = 6
)
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.2.0 ✔ readr 2.2.0
## ✔ forcats 1.0.1 ✔ stringr 1.6.0
## ✔ ggplot2 4.0.3 ✔ tibble 3.3.1
## ✔ lubridate 1.9.5 ✔ tidyr 1.3.2
## ✔ purrr 1.2.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
## Loading required package: lattice
##
## Attaching package: 'caret'
##
## The following object is masked from 'package:purrr':
##
## lift
## Metric functions loaded
## - compute_gini()
## - compute_ks()
## - compute_auc()
## - compute_calibration()
## - compute_lift()
## - compute_confusion_matrix()
## - summarise_metrics()
## Modeling & Validation
## Loaded all libraries and metric functions
This notebook implements the modeling & validation stage of the Credit Risk Intelligence Platform:
# load full dataset
df_model <- read_csv("../data/processed/modeling_dataset_woe.csv")
cat("Dataset loaded\n")## Dataset loaded
## Dimensions: 307511 × 24
## Applicants: 307,511 | Features: 22
# verify
stopifnot(nrow(df_model) == 307511, "target" %in% names(df_model))
# target distribution
cat("Target Distribution:\n")## Target Distribution:
##
## 0 1
## 282686 24825
##
## Default rate: 8.07 %
We use stratified random splitting to preserve the class distribution in both sets.
set.seed(42)
# stratified split via caret
train_idx <- caret::createDataPartition(
df_model$target,
p = 0.7,
list = FALSE
)
train <- df_model[train_idx, ]
test <- df_model[-train_idx, ]
cat("Training set:", nrow(train), "rows\n")## Training set: 215258 rows
## Test set: 92253 rows
## Total: 307511
# verify stratification
train_default_rate <- mean(train$target) * 100
test_default_rate <- mean(test$target) * 100
overall_default_rate <- mean(df_model$target) * 100
cat("Default Rates:\n")## Default Rates:
## Train: 8.061 %
## Test: 8.102 %
## Overall: 8.073 %
## → Stratification successful (rates close)
# save
saveRDS(train, "../data/processed/train_woe.rds")
saveRDS(test, "../data/processed/test_woe.rds")We fit a logistic regression on all 40 WOE features. Logistic regression is chosen for: 1. Interpretability: Coefficients map 1:1 to log-odds 2. Regulatory compliance: Explainable, monotonic 3. Scorecard conversion: Directly scalable to points
# drop 6 count features that are collinear and alias to NA coefficients
collinear_drop <- c("bureau_count_woe", "bureau_active_count_woe",
"inst_total_count_woe", "card_dpd_mean_woe",
"pos_count_woe", "prev_app_count_woe")
feature_names <- setdiff(names(train), c("target", "SK_ID_CURR", collinear_drop))
# create formula
formula_str <- paste("target ~", paste(feature_names, collapse = " + "))
# fit model
model_lr <- glm(
as.formula(formula_str),
data = train,
family = binomial(link = "logit")
)
cat("Logistic Regression Model Summary\n")## Logistic Regression Model Summary
## ─────────────────────────────────
## Observations: 215258
## Null deviance: 120650.4
## Residual deviance: 115584.8
## AIC: 115618.8
## Features: 16
coef_df <- tidy(model_lr) %>%
filter(term != "(Intercept)") %>%
arrange(desc(abs(statistic))) %>%
mutate(
estimate = round(estimate, 4),
statistic = round(statistic, 3),
p.value = round(p.value, 4),
significant = ifelse(p.value < 0.05, "***", "")
)
cat("Top 15 Features by |t-statistic|\n")## Top 15 Features by |t-statistic|
## # A tibble: 15 × 5
## term estimate statistic p.value significant
## <chr> <dbl> <dbl> <dbl> <chr>
## 1 prev_app_approval_rate_woe 0.916 29.3 0 "***"
## 2 days_employed_woe 0.686 23.4 0 "***"
## 3 days_birth_woe 0.546 16.3 0 "***"
## 4 amt_credit_woe 0.621 15.3 0 "***"
## 5 bureau_credit_sum_total_woe 0.919 14.0 0 "***"
## 6 bureau_overdue_max_ever_woe 0.883 12.3 0 "***"
## 7 bureau_months_paid_on_time_share_woe 0.860 11.0 0 "***"
## 8 inst_late_share_woe 0.520 8.14 0 "***"
## 9 inst_payment_ratio_mean_woe 0.462 6.88 0 "***"
## 10 prev_app_amount_mean_woe 0.429 6.78 0 "***"
## 11 amt_income_total_woe 0.526 6.13 0 "***"
## 12 amt_annuity_woe 0.342 5.78 0 "***"
## 13 inst_max_days_late_woe 0.190 2.63 0.0085 "***"
## 14 card_dpd_share_woe 1.67 1.00 0.317 ""
## 15 pos_dpd_mean_woe 0.0587 0.79 0.429 ""
# sign distribution
n_pos <- sum(coef_df$estimate > 0, na.rm = TRUE)
n_neg <- sum(coef_df$estimate < 0, na.rm = TRUE)
n_na <- sum(is.na(coef_df$estimate))
cat("\nPositive coefficients:", n_pos, "\n")##
## Positive coefficients: 16
## Negative coefficients: 0
## Aliased (collinear): 0
# predict probabilities
pred_prob_lr <- predict(model_lr, newdata = test, type = "response")
cat("Prediction Distribution (Test Set)\n")## Prediction Distribution (Test Set)
## ──────────────────────────────────
## Min: 0.0096
## Q1: 0.0485
## Median: 0.0713
## Mean: 0.0807
## Q3: 0.1024
## Max: 0.3869
# also on training set (to check overfitting)
pred_prob_lr_train <- predict(model_lr, newdata = train, type = "response")gini_test <- compute_gini(test$target, pred_prob_lr)
ks_test <- compute_ks(test$target, pred_prob_lr)
auc_test <- compute_auc(test$target, pred_prob_lr)
cat("Logistic Regression — Test Set Metrics\n")## Logistic Regression — Test Set Metrics
## ───────────────────────────────────────
## Gini: 0.3126
## KS: 0.2272
## AUC: 0.6563
##
## Context:
## KS 0.227 is at the 0.25 industry rule-of-thumb for a usable model.
## Gini 0.313 is modest but honest for Home Credit without leakage.
gini_train <- compute_gini(train$target, pred_prob_lr_train)
ks_train <- compute_ks(train$target, pred_prob_lr_train)
auc_train <- compute_auc(train$target, pred_prob_lr_train)
metrics_lr_summary <- tibble(
Metric = c("Gini", "KS", "AUC"),
Train = c(gini_train, ks_train, auc_train),
Test = c(gini_test, ks_test, auc_test),
Difference = Train - Test
) %>%
mutate(across(where(is.numeric), ~round(., 4)))
cat("Train vs. Test Metrics (Overfitting Check)\n")## Train vs. Test Metrics (Overfitting Check)
## ──────────────────────────────────────────
## # A tibble: 3 × 4
## Metric Train Test Difference
## <chr> <dbl> <dbl> <dbl>
## 1 Gini 0.315 0.313 0.0021
## 2 KS 0.234 0.227 0.0071
## 3 AUC 0.657 0.656 0.001
##
## Interpretation:
gini_gap <- metrics_lr_summary$Difference[1]
if (gini_gap < 0.02) {
cat("No overfitting detected (gap <0.02)\n")
} else {
cat("Possible overfitting (gap =", round(gini_gap, 4), ")\n")
}## No overfitting detected (gap <0.02)
Model calibration checks if predicted probabilities match observed default rates.
calib_df <- compute_calibration(test$target, pred_prob_lr, n_bins = 10)
cat("Calibration by Decile (Expected vs. Observed Default Rate)\n")## Calibration by Decile (Expected vs. Observed Default Rate)
print(calib_df %>%
select(bin, expected_rate, observed_rate, n) %>%
mutate(across(where(is.numeric) & !matches("^n$"),
~round(., 4))))## # A tibble: 4 × 4
## bin expected_rate observed_rate n
## <fct> <dbl> <dbl> <int>
## 1 Bin_1 0.0598 0.0591 67965
## 2 Bin_2 0.132 0.137 22586
## 3 Bin_3 0.229 0.208 1619
## 4 Bin_4 0.319 0.313 83
# calibration error
calib_error <- mean((calib_df$expected_rate - calib_df$observed_rate)^2, na.rm = TRUE)
cat("\nCalibration MSE:", round(calib_error, 6), "\n")##
## Calibration MSE: 0.000135
## (Lower is better; <0.002 is excellent)
if (calib_error < 0.002) {
cat("Model is well-calibrated\n")
} else {
cat("Model may need recalibration\n")
}## Model is well-calibrated
To justify the choice of logistic regression, we train an XGBoost model for comparison.
## Setting up XGBoost...
# use the SAME train/test split as logistic regression
# XGBoost works well with WOE features too (trees will discover splits automatically)
train_xgb <- train %>% select(all_of(feature_names))
test_xgb <- test %>% select(all_of(feature_names))
cat("Data prepared:\n")## Data prepared:
## Train shape: 215258 × 16
## Test shape: 92253 × 16
# create XGBoost matrices
dtrain <- xgb.DMatrix(
data = as.matrix(train_xgb),
label = train$target
)
dtest <- xgb.DMatrix(
data = as.matrix(test_xgb),
label = test$target
)
cat("XGBoost matrices created\n\n")## XGBoost matrices created
# set parameters
scale_pos_weight <- sum(train$target == 0) / sum(train$target == 1)
params <- list(
objective = "binary:logistic",
eval_metric = "auc",
scale_pos_weight = scale_pos_weight,
max_depth = 6,
eta = 0.1,
subsample = 0.8,
colsample_bytree = 0.8,
lambda = 1.0,
seed = 42
)
cat("XGBoost Configuration:\n")## XGBoost Configuration:
## Features: 16
## Class weight (pos_weight): 11.41
## Max depth: 6 | Learning rate: 0.1
# train model
set.seed(42)
model_xgb <- xgb.train(
params = params,
data = dtrain,
nrounds = 100,
watchlist = list(train = dtrain, test = dtest),
verbose = 1,
print_every_n = 10,
early_stopping_rounds = 20
)## Multiple eval metrics are present. Will use test_auc for early stopping.
## Will train until test_auc hasn't improved in 20 rounds.
##
## [1] train-auc:0.638962 test-auc:0.625703
## [11] train-auc:0.671921 test-auc:0.656308
## [21] train-auc:0.678520 test-auc:0.659670
## [31] train-auc:0.683407 test-auc:0.660796
## [41] train-auc:0.687392 test-auc:0.661571
## [51] train-auc:0.691835 test-auc:0.661905
## [61] train-auc:0.695652 test-auc:0.661622
## Stopping. Best iteration:
## [71] train-auc:0.699039 test-auc:0.661365
##
## [71] train-auc:0.699039 test-auc:0.661365
##
## XGBoost Training Complete
# predictions
pred_prob_xgb <- predict(model_xgb, dtest)
pred_prob_xgb_train <- predict(model_xgb, dtrain)
cat("Predictions generated\n")## Predictions generated
# quick held-out AUC sanity check (full metrics in the next section)
cat("Held-out test AUC:", round(compute_auc(test$target, pred_prob_xgb), 4), "\n")## Held-out test AUC: 0.6619
gini_test_xgb <- compute_gini(test$target, pred_prob_xgb)
ks_test_xgb <- compute_ks(test$target, pred_prob_xgb)
auc_test_xgb <- compute_auc(test$target, pred_prob_xgb)
gini_train_xgb <- compute_gini(train$target, pred_prob_xgb_train)
ks_train_xgb <- compute_ks(train$target, pred_prob_xgb_train)
auc_train_xgb <- compute_auc(train$target, pred_prob_xgb_train)
metrics_xgb_summary <- tibble(
Metric = c("Gini", "KS", "AUC"),
Train = c(gini_train_xgb, ks_train_xgb, auc_train_xgb),
Test = c(gini_test_xgb, ks_test_xgb, auc_test_xgb),
Difference = Train - Test
) %>%
mutate(across(where(is.numeric), ~round(., 4)))
cat("XGBoost - Test Set Metrics\n")## XGBoost - Test Set Metrics
## # A tibble: 3 × 4
## Metric Train Test Difference
## <chr> <dbl> <dbl> <dbl>
## 1 Gini 0.384 0.324 0.0599
## 2 KS 0.278 0.233 0.0447
## 3 AUC 0.692 0.662 0.0299
comparison <- tibble(
Metric = c("Gini (Train)", "Gini (Test)",
"KS (Train)", "KS (Test)",
"AUC (Train)", "AUC (Test)"),
Logistic = c(
gini_train, gini_test,
ks_train, ks_test,
compute_auc(train$target, pred_prob_lr_train),
compute_auc(test$target, pred_prob_lr)
),
XGBoost = c(
gini_train_xgb, gini_test_xgb,
ks_train_xgb, ks_test_xgb,
auc_train_xgb, auc_test_xgb
)
) %>%
mutate(
Difference = XGBoost - Logistic,
across(where(is.numeric), ~round(., 4))
)
cat("COMPREHENSIVE MODEL COMPARISON\n")## COMPREHENSIVE MODEL COMPARISON
## ════════════════════════════════
## # A tibble: 6 × 4
## Metric Logistic XGBoost Difference
## <chr> <dbl> <dbl> <dbl>
## 1 Gini (Train) 0.315 0.384 0.069
## 2 Gini (Test) 0.313 0.324 0.0112
## 3 KS (Train) 0.234 0.278 0.0434
## 4 KS (Test) 0.227 0.233 0.0058
## 5 AUC (Train) 0.657 0.692 0.0345
## 6 AUC (Test) 0.656 0.662 0.0056
gini_diff <- gini_test_xgb - gini_test
cat("\n\nGini Difference (XGBoost - Logistic):", round(gini_diff, 4), "\n")##
##
## Gini Difference (XGBoost - Logistic): 0.0112
Lift analysis shows what % of defaults are captured when we rank by model score.
lift <- compute_lift(test$target, pred_prob_lr)
cat("Lift-by-Decile Analysis (Logistic Regression)\n")## Lift-by-Decile Analysis (Logistic Regression)
## Ranked by model probability (highest risk first)
print(lift %>%
select(decile, n_applicants, n_defaults, observed_default_rate,
cumulative_capture_rate) %>%
mutate(across(where(is.numeric) & !contains("rate"), as.integer),
across(contains("rate"), ~round(., 4))))## # A tibble: 10 × 5
## decile n_applicants n_defaults observed_default_rate cumulative_capture_rate
## <int> <int> <int> <dbl> <dbl>
## 1 1 9225 1638 0.178 0.219
## 2 2 9225 1176 0.128 0.376
## 3 3 9225 938 0.102 0.502
## 4 4 9226 790 0.0856 0.608
## 5 5 9225 690 0.0748 0.7
## 6 6 9225 622 0.0674 0.783
## 7 7 9226 520 0.0564 0.853
## 8 8 9225 448 0.0486 0.913
## 9 9 9225 353 0.0383 0.96
## 10 10 9226 299 0.0324 1
##
## Interpretation:
cat(" Top 10% (decile 1) captures",
round(lift$cumulative_capture_rate[1] * 100, 1),
"% of all defaults\n")## Top 10% (decile 1) captures 21.9 % of all defaults
cat(" Top 30% (deciles 1–3) capture",
round(lift$cumulative_capture_rate[3] * 100, 1),
"% of all defaults\n")## Top 30% (deciles 1–3) capture 50.2 % of all defaults
# save models
saveRDS(model_lr, "../outputs/model_logistic.rds")
saveRDS(model_xgb, "../outputs/model_xgboost.rds")
# save test predictions
test_with_pred <- test %>%
mutate(
pred_lr = pred_prob_lr,
pred_xgb = pred_prob_xgb
)
saveRDS(test_with_pred, "../data/processed/test_with_predictions.rds")
# save comparison table
write_csv(comparison, "../outputs/model_comparison_metrics.csv")
cat("\nModeling complete. All artifacts saved.\n")##
## Modeling complete. All artifacts saved.
## - outputs/model_logistic.rds
## - outputs/model_xgboost.rds
## - outputs/model_comparison_metrics.csv
## - outputs/lift_by_decile.csv
## - data/processed/test_with_predictions.rds