Source code · R
crosscheck.R
# Independent check of the Python results in base R.## 1. Refit the logistic PD model with glm() on the same WoE inputs and the same# train/test split, then compare coefficients, predicted PDs and test AUC.# 2. Recompute expected loss and Basel IRB capital for the whole book from the# account-level PD, LGD and EAD.## Usage: Rscript crosscheck.R <output directory># Reads model_data.csv and portfolio_accounts.csv; writes r_crosscheck.json. args <- commandArgs(trailingOnly = TRUE)out_dir <- if (length(args) > 0) args[1] else "output" d <- read.csv(file.path(out_dir, "model_data.csv"))features <- setdiff(names(d), c("default", "test", "pd_python"))train <- d[d$test == 0, ]test <- d[d$test == 1, ] fit <- glm(reformulate(features, response = "default"), family = binomial(), data = train) # AUC through the Mann-Whitney rank statistic, the same definition as in Python.auc <- function(y, p) { r <- rank(p) n1 <- sum(y == 1) n0 <- length(y) - n1 (sum(r[y == 1]) - n1 * (n1 + 1) / 2) / (n1 * n0)} pd_r <- predict(fit, newdata = d, type = "response")test_auc <- auc(test$default, pd_r[d$test == 1])max_pd_gap <- max(abs(pd_r - d$pd_python)) # Basel IRB retail capital: K = LGD * [N((G(PD) + sqrt(R) G(0.999)) / sqrt(1 - R)) - PD]a <- read.csv(file.path(out_dir, "portfolio_accounts.csv"))rho <- 0.04k <- a$lgd * (pnorm((qnorm(a$pd) + sqrt(rho) * qnorm(0.999)) / sqrt(1 - rho)) - a$pd)capital <- sum(k * a$ead)expected_loss <- sum(a$pd * a$lgd * a$ead) # Base R has no JSON writer, so the file is assembled by hand.coefs <- coef(fit)coef_json <- paste(sprintf(' "%s": %.12g', names(coefs), coefs), collapse = ",\n")json <- sprintf( '{\n "r_version": "%s",\n "coefficients": {\n%s\n },\n "test_auc": %.12g,\n "max_pd_gap": %.6g,\n "expected_loss": %.12g,\n "irb_capital": %.12g\n}\n', R.version.string, sub("^\\(Intercept\\)", "const", coef_json), test_auc, max_pd_gap, expected_loss, capital)json <- sub('"\\(Intercept\\)"', '"const"', json)writeLines(json, file.path(out_dir, "r_crosscheck.json"))cat(json)