R 评分卡 scorecard包

library(knitr)
library(klaR)
library(scorecard)

一、数据预处理

1.读取数据

setwd("C:\\Users\\91333\\Documents\\semiester5\\RegressionAnalysis\\hw2")
base_info <- read.csv(file = "base_info.csv", header = TRUE)
hive <- read.csv(file = "hive.csv", header = TRUE)

2.数据框的融合

mydata <- merge(base_info, hive, by.x = "ORDER_ID",by.y = "order_id")

3.去重

mydata <- mydata[!(duplicated(mydata$ORDER_ID)), ]
mydata <- mydata[, !duplicated(t(mydata))]

4.重命名：以order-id 为行名

rownames(mydata) <- mydata$ORDER_ID
mydata$ORDER_ID <- NULL

5.分离目标变量

STATUS <- mydata$STATUS

6.使用scorecard包中的函数对mydata进行初次过滤

mydata <- var_filter(mydata, y = "STATUS")
[INFO] filtering variables ...  

7.缺失值的处理

1）检查有无缺失值

sum(is.na(mydata))
[1] 18 

2）寻找缺失值的位置并做出处理

which(rowSums(is.na(mydata)) > 0)
[1]  438 1192 1229 1684 
rowSums(is.na(mydata))[rowSums(is.na(mydata)) > 0]
[1] 3 5 5 5 

只有四个样本具有缺失值，查看这四个样本，发现缺失值都为juxinliinfo相关信息，而这一信息在原数据中大范围缺失，原数据使用-99999表示数据的缺失，而且数据的缺失具有较强的聚集性，即往往出现同一样本juxinliinfo相关数据全部缺失或者全部完整的情况，所以这里直接把缺失值赋值为-99999。

mydata[c(489, 1243, 1280, 1822), ]
mydata[is.na(mydata)] <- -99999

8.数据重新编码

统一缺失数据的记号，把有一些变量的-1记号改为-99999

观察数据后发现，对于存在缺失数据的变量，缺失数据的比重要高于其他任何取值的比重，所以可以利用条件：众数为-1，来定位需要重新变慢缺失值的变量。

minus1index <-c (1 : length(mydata[1, ]))[sapply(mydata, function(x){names((table(x))[1])}) == -1]
mydata <- as.matrix(mydata)
for (i in minus1index){a = mydata[, i]
mydata[, i] <- ifelse(a == -1, -99999, a)}
mydata<-as.data.frame(mydata)

二、变量选择

1. 第一次变量选择

第一次变量选择发在数据预处理时，首先，去重的过程删去了一些变量；其次，使用var_filter函数删去了缺失率高于0.95、同一变量中数据重复率高达0.95的数据，同时也删去了IV值小于0.02的数据。

2. 第二次变量选择

观察Mydata的样本相关系数矩阵，发现存在大量变量的相关系数在0.8以上，为了提高后面的逐步回归的速度，根据相关性强弱，删去一些与其他变量高度相关的变量，降低数据维度。

mydata2 <- mydata 
mydata2$STATUS <- NULL 
while(TRUE){
   corr <- abs(cor(mydata2))
   corr_each <- colSums(corr)
      for (i in 1 : length(corr[1, ])){
          for (j in (1 : length(corr[1, ]))){
               if (i <= j){corr[i, j] <- -99999}}}
	temp <- unique(as.vector(which(corr == corr[which.max(corr)], arr.ind= T))) 
	mydata2[, temp[which.max(corr_each[temp])]] <- NULL 
	if (max(corr, na.rm = T) < 0.9) break } 

3.第三次变量选择

第三次变量选择将以IV值为依据，进一步筛选变量，虽然第一步时也用到了以IV值为依据的变量选择，但是我发现var_filter函数的IV值计算并不是传统地分箱之后再进行IV值计算，并且var_filter函数计算的IV值和分箱之后的IV值不等，要偏大一些，所以现在，我要将剩下来的变量分箱、计算IV值、计算WOE，这里用到scorecard包的bins函数。

bins <- woebin(data.frame(mydata2,STATUS), y="STATUS")
[INFO] creating woe binning ...  
[INFO] Binning on 1955 rows and 53 columns in 00:00:15 

all_iv <- sapply(bins, function(x){return(x$total_iv[1])})
mydata2[,all_iv < 0.015] <- NULL

三、建立模型

1. 划分训练集和测试集

直接利用scorecard包中的划分数据函数split_df，划分完之后计算了以下STATUS中的0和1在训练集和测试集中的比例是否大致相同，的确大致相同。

dt_list <- split_df(data.frame(mydata2, STATUS))
train <- dt_list$train
test <- dt_list$test
sum(test$STATUS)/length(test$STATUS)
[1] 0.07731959 
sum(train$STATUS)/length(train$STATUS)
[1] 0.0815732 

2. 将原变量划分区间，赋WOE值

train_woe <-  woebin_ply(train, bins)
[INFO] converting into woe values ...  

test_woe <-  woebin_ply(test, bins)
[INFO] converting into woe values ... 

3. 逐步回归选择变量

m1 <- glm( STATUS ~ ., family = binomial(), data = train_woe)
m_step <- step(m1, direction = "both", trace = FALSE)
m2 <- eval(m_step$call)
summary(m2)
Call: 
glm(formula = STATUS ~ APPLY_DT_woe + 偿债比_woe + LOAN_AMT_woe +  
    MOBILE_TIME_SPAN_woe + EDUCATION_woe + SEX_woe + UNIT_TYPE_woe +      
    APPLY_DATE_woe + MULTIPLE_LOAN_7D_woe + juxinliinfo_person3phone2_6m_ woe +      
    unionpayinfo_all_bank_cost_sum_woe + bairongmultiapplyinfo_al_m12_cel l_notbank_orgnum_woe +      
    yimeiloaninfo_bank_application_platform_num_9m_woe, family = binomial (),      
    data = train_woe) 

由上面的结果可以看出，逐步回归最终选择出13个变量，且每个变量的系数都比较显著。

四、模型性能验证

1. 绘制KS曲线和ROC曲线

train_pred <- predict(m2, train_woe, type = 'response')
test_pred <- predict(m2, test_woe, type = 'response')
train_perf <- perf_eva(train_pred,train$STATUS, title = 'train')
[INFO] The threshold of confusion matrix is 0.1643. 

test_perf <- perf_eva(test_pred,test$STATUS, title = 'test')
[INFO] The threshold of confusion matrix is 0.1315.

KS用来衡量坏客户和好客户这两个样本上的评分变量的分布之差。
$$KS=\underset{s}{\max }|B(s)-G(s)|$$
KS是最主要的模型评价指标，KS越高，模型越好,根据上图，训练集KS值为0.05327，测试集KS值为0.5011,都高于0.40，说明模型效果比较好；但是训练集KS值比测试集KS值高了0.01以上，模型拟合效果和预测效果的差距比其他模型略大。

2. 模型的稳定性度量：PSI

card <- scorecard(bins, m2)
train_score <- scorecard_ply(train, card, print_step = 0)
test_score <- scorecard_ply(test, card, print_step = 0)
psi_result <- perf_psi(
  score = list(train = train_score, test = test_score),
  label = list(train = train$STATUS, test = test$STATUS))
psi_result 

$pic
$pic$score 

$psi
    variable    dataset        
psi 1:    score train_test 0.01137085

PSI指标是为了检验建模样本和验证样本以及最新样本之间的样本变化，如果样本变化太大，可能导致评分卡失效，该模型求得的PSI = 0.0114 < 0.2
同时，从图中可以看出样本的分数接近正态，且训练集样本和测试集样本评分基本一致。

五、评分卡输出