申请评分卡（3）——建模（R）

理论说完了，来次实践。

数据理解与预处理

数据来自kaggle的Give Me Some Credit项目,有15万条的样本数据。
要求根据历史数据，预测申请人违约的可能性，以此作为放贷依据。

#数据加载，也可以直接放含有绝对路径的文件名
setwd('D:\\anajpnotebook\\GiveMeSomeCredit\\data')
train <- read.csv('cs-training.csv')
View(train)

首先我们看下业务含义。
数据集一共11列，分别是：
SeriousDlqin2yrs：贷款人会逾期90天以上，取值Y/N，对应1和0，这个是目标变量，也就是贷款人会不会违约。
age：年龄，一般具备独立经济能力的申请人偿还债务的可能性较高，而且因为其他因素的影响，原则上不会发放贷款给未成年人，一般国内是这样，国外还要看具体数据。
DebtRatio：负债率，如果本身要偿还的债务比较高，在月收入一定的情况下，违约的可能性就可能就比较高。
MonthlyIncome：月收入，还贷的直接来源之一，是衡量还贷能力的指标，影响是否违约的重要因素。
RevolvingUtilizationOfUnsecuredLines：信用总余额（剔除房产贷款和分期付款）/信用额度，和DebtRatio异曲同工，这个指标越低，违约的风险可能就越大。
NumberOfOpenCreditLinesAndLoans：开放式贷款数量（和信用额度），历史贷款记录，可能有正相关关系，如果历史有很多次贷款且没有违约记录，一定程度上说明该申请人信用习惯良好，申请人骗贷的情形除外。
NumberRealEstateLoansOrLines：抵押贷款和房地产贷款数量（包括房屋信贷额度）。NumberOfTime30.59DaysPastDueNotWorse：逾期30-59天的次数，近2年内没有更糟的记录，这个系列的指标对应国内的M1,M2,M3，基本类似。也就是不同程度的不良记录数。
NumberOfTime60.89DaysPastDueNotWorse：逾期60-89天的次数，近2年内没有更糟的记录。
NumberOfTimes90DaysLate：逾期90天以上的次数。
NumberOfDependents：家属个数，不含自身。

接下来看看数据的取值和分布情况：

先看缺失值情况，MonthlyIncome:29731 NumberOfDependents:3924。
以上2项含有缺失值，NumberOfDependents缺失值较少，占比约2.6%，可以考虑去掉含缺失值的样本，也可以用中位数填充。MonthlyIncome缺失值占比中等，需要填充或剔除该项属性，具体要看数据的分布和数据间的关系。
我们先来看看变量之间的相关关系：

> cor(train)

变量之间没有显著的相关关系。所以不能凭借相关关系来填充数据了。
缺失值处理比较省事的一种方式是用中位数或平均数来填充，由2个变量的统计描述指标看，用中位数比较合适。当然也可以建模来填充，比如聚类、回归或者随机森林。这里我们先用中位数填充。
由于缺乏相关业务信息，这里只对明显不符合常识的几个变量中的异常值做处理。
首先是年龄，年龄的最小值是0，样本数量只有1个，予以剔除。

> train<-train[-which(train$age==0),] 

接着是信用往来相关的变量，取值98、96的样本共269个，其中好样本和坏样本比例基本持平，且占总样本比例较小，予以剔除。

> train<-train[-which(train$NumberOfTime30.59DaysPastDueNotWorse==98),] 
> train<-train[-which(train$NumberOfTime30.59DaysPastDueNotWorse==96),] 

其他变量虽也存在异常值，不过信息不够，暂时不能剔除。
下面看看目标变量的分布情况。

> table(train$SeriousDlqin2yrs)

     0      1 
136229   9847 

可以看到，违约样本约占总样本的1/15，比例失衡，解决这个问题有2种方法：在“0”这一部分抽取与“1”这边相当的数据，也就是下采样，另外一种就是增加“1”的样本数量与“0”基本持平，也就是过采样。下采样之后数据实在有点少，这里我们采用过采样中常用的SMOTE算法。

> set.seed(1206) 
> splitIndex<-createDataPartition(train$SeriousDlqin2yrs,time=1, p=0.5,list=FALSE) 
> ttrain<-train[splitIndex,] 
> ttest<-train[-splitIndex,] 
> prop.table(table(ttrain$SeriousDlqin2yrs))

         0          1 
0.93432178 0.06567822 
> prop.table(table(ttest$SeriousDlqin2yrs))

        0         1 
0.9337207 0.0662793 

训练集和测试集各一半，两者的分类结果是平衡的，坏样本比例保持在6.6%左右，我们可以用这个来建模。

逻辑回归

逻辑回归是信用评分卡的核心之一。由于其本身的特点，引入WOE后，逻辑回归的结果可以直接转换成标准评分卡格式。
首先我们用全变量来回归试试。

> fit<-glm(SeriousDlqin2yrs~.,ttrain,family = "binomial")
> summary(fit)

image.png

有3个变量未通过检验。去掉这3个变量之后再进行一次回归。

> fit2<-glm(SeriousDlqin2yrs~age+NumberOfTime30.59DaysPastDueNotWorse+MonthlyIncome+NumberOfTimes90DaysLate+NumberRealEstateLoansOrLines+NumberOfTime60.89DaysPastDueNotWorse+NumberOfDependents,ttrain,family = "binomial")
> summary(fit2)

AIC更小了，暂且就先用这些变量来建模。

模型评估

评价逻辑回归模型我们主要看ROC和AUC指标。

由图中我们可以看到FPR=0.850，TPR=0.636，AUC=0.812，准确率还不错。回到实际意义，TPR越高，本金损失风险越低；FPR越高，盈利可能性越高；具体比重要看利率和销售策略。

WOE转换

引入WOE的目的是为了精简特征。用WOE(x)替换变量x。

1. 变量处理

> bin_age<- c(-Inf,30,35,40,45,50,55,60,65,75,Inf)
> plot(cut(ttrain$age,bin_age))

bin_m2<-c(-Inf,0,1,3,5,Inf)
plot(cut(ttrain$NumberOfTime30.59DaysPastDueNotWorse,bin_m2))

> bin_m3<-c(-Inf,0,1,3,5,Inf)
> plot(cut(ttrain$NumberOfTime60.89DaysPastDueNotWorse,bin_m3))

> bin_m3_<-c(-Inf,0,1,3,5,Inf)
> plot(cut(ttrain$NumberOfTimes90DaysLate,bin_m3_))

> bin_rs<-c(-Inf,0,1,3,5,Inf)
> plot(cut(ttrain$NumberRealEstateLoansOrLines,bin_rs))

> bin_mincome<-c(-Inf,2000,4000,6000,8000,10000,12000)
> plot(cut(ttrain$MonthlyIncome,bin_mincome))

> bin_d<-c(-Inf,0,1,2,3,4,5,Inf)
> plot(cut(ttrain$NumberOfDependents,bin_d))

2. 变量变换
   首先计算WOE，WOE()=ln[(违约/总违约)/(正常/总正常)]。

#计算WOE
totalgood = as.numeric(table(ttrain$SeriousDlqin2yrs))[1]
totalbad = as.numeric(table(ttrain$SeriousDlqin2yrs))[2]
getWOE <- function(a,p,q)
{
  Good <- as.numeric(table(ttrain$SeriousDlqin2yrs[a > p & a <= q]))[1]
  Bad <- as.numeric(table(ttrain$SeriousDlqin2yrs[a > p & a <= q]))[2]
  WOE <- log((Bad/totalbad)/(Good/totalgood),base = exp(1))
  return(WOE)
}

接着进行数值替换

#age
age.WOE=c(getWOE(ttrain$age,-Inf,30),getWOE(ttrain$age,30,35),getWOE(ttrain$age,35,40),getWOE(ttrain$age,40,45),getWOE(ttrain$age,45,50),getWOE(ttrain$age,50,55),
          getWOE(ttrain$age,55,60),getWOE(ttrain$age,60,65),getWOE(ttrain$age,65,75),getWOE(ttrain$age,75,Inf))

tmp.age <- 0
for(i in 1:nrow(ttrain)) {
  if(ttrain$age[i] <= 30)
    tmp.age[i] <- Agelessthan30.WOE
  else if(ttrain$age[i] <= 35)
    tmp.age[i] <- Age30to35.WOE
  else if(ttrain$age[i] <= 40)
    tmp.age[i] <- Age35to40.WOE
  else if(ttrain$age[i] <= 45)
    tmp.age[i] <- Age40to45.WOE
  else if(ttrain$age[i] <= 50)
    tmp.age[i] <- Age45to50.WOE
  else if(ttrain$age[i] <= 55)
    tmp.age[i] <- Age50to55.WOE
  else if(ttrain$age[i] <= 60)
    tmp.age[i] <- Age55to60.WOE
  else if(ttrain$age[i] <= 65)
    tmp.age[i] <- Age60to65.WOE
  else if(ttrain$age[i] <= 75)
    tmp.age[i] <- Age65to75.WOE
  else
    tmp.age[i] <- Agemorethan.WOE
}

其他变量做同样的处理。

建模并制作评分卡

替换完成后进行逻辑回归建模。

#WOE DataFrame构建
trainWOE =cbind.data.frame(tmp.age,tmp.m2,tmp.mincome,tmp.m3,tmp.m3_,tmp.rs,tmp.d)
#建模
trainWOE$SeriousDlqin2yrs = ttrain$SeriousDlqin2yrs
glm.fit = glm(SeriousDlqin2yrs~.,data = trainWOE,family = binomial(link = logit))
summary(glm.fit)
coe = (glm.fit$coefficients)

此处我们得到了计算分数需要用到的系数，下面来构建评分卡。

#构建评分卡,score=p-q*log(odds),odds=p(bad)/p(godd),比例越高，分数越低。如果是好坏比，则中间变成加号。
p <- 20/log(2)
q <- 600-20*log(15)/log(2)
getscore<-function(i,x){
  score = round(p*as.numeric(coe[i])*x,0)
  return(score)
}
#基础分计算
base <- q + p*as.numeric(coe[1])
base
#[1] 446.0211
#计算每个分箱的得分
Age.SCORE = c(getscore(2,age.WOE[1]),getscore(2,age.WOE[2]),getscore(2,age.WOE[3]),getscore(2,age.WOE[4]),getscore(2,age.WOE[5]),
              getscore(2,age.WOE[6]),getscore(2,age.WOE[7]),getscore(2,age.WOE[8]),getscore(2,age.WOE[9]),getscore(2,age.WOE[10]))
Age.SCORE
# [1]  12  11   7   6   5   2  -6 -13 -20 -28
m2.SCORE=c(getscore(3,m2.WOE[1]),getscore(3,m2.WOE[2]),getscore(3,m2.WOE[3]),getscore(3,m2.WOE[4]),getscore(3,m2.WOE[5]))
m2.SCORE
#[1] -10  17  34  45  50
m3.SCORE=c(getscore(5,m3.WOE[1]),getscore(5,m3.WOE[2]),getscore(5,m3.WOE[3]),getscore(5,m3.WOE[4]),getscore(3,m3.WOE[5]))
m3.SCORE
#[1] -4 26 38 40 55
m3_.SCORE=c(getscore(6,m3_.WOE[1]),getscore(6,m3_.WOE[2]),getscore(6,m3_.WOE[3]),getscore(6,m3_.WOE[4]),getscore(3,m3_.WOE[5]))
m3_.SCORE
#[1] -7 38 51 62 63
mincome.SCORE=c(getscore(4,mincome.WOE[1]),getscore(4,mincome.WOE[2]),getscore(4,mincome.WOE[3]),getscore(4,mincome.WOE[4]),
                getscore(4,mincome.WOE[5]),getscore(4,mincome.WOE[6]),getscore(4,mincome.WOE[7]),getscore(4,mincome.WOE[8]),
                getscore(4,mincome.WOE[9]),getscore(4,mincome.WOE[10]),getscore(4,mincome.WOE[11]))
mincome.SCORE
# [1] -6  8  7  5  2 -2 -2 -3 -5 -7 -7
rs.SCORE=c(getscore(7,rs.WOE[1]),getscore(7,rs.WOE[2]),getscore(7,rs.WOE[3]),getscore(7,rs.WOE[4]),getscore(7,rs.WOE[5]))
rs.SCORE
#[1]  5 -6 -3  9 24
d.SCORE=c(getscore(8,d.WOE[1]),getscore(8,d.WOE[2]),getscore(8,d.WOE[3]),getscore(8,d.WOE[4]),
          getscore(8,d.WOE[5]),getscore(8,d.WOE[6]),getscore(8,d.WOE[7]))
d.SCORE
#[1] -2  2  4  3  6  3  8

将这些分数合并成一张评分表，初步的评分卡就做好了。