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用来衡量坏客户和好客户这两个样本上的评分变量的分布之差。
K
S
=
max
s
∣
B
(
s
)
−
G
(
s
)
∣
KS=\max_s{|B(s)-G(s)|}
KS=smax∣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
同时,从图中可以看出样本的分数接近正态,且训练集样本和测试集样本评分基本一致。
五、评分卡输出