使用LR构建评分卡模型

使用LR构建评分卡模型

第一步：获取数据和读取数据

数据为个人消费类贷款数据，目标是建成A卡的信用卡评分模型。一共有这么几个特征
SeriousDlqin2yrs(出现90天或更长时间的逾期行为)
SeriousDlqin2yrs1时，为坏客户；SeriousDlqin2yrs0时为好客户；
RevolvingUtilizationofunsecuredlines(贷款以及信用卡可用额度与总额度比例)
age（借款人借款年龄）
Numberoftime30-59dayspastduenotworse(过去两年内出现35-59天逾期但是没有发展得更坏的次数)
DebtRatio(每月偿还债务，赡养费，生活费用除以月总收入)
Monthlyincome（月收入）
Numberoftimes90dayslate（过去2年内出现90天逾期或更坏的次数）
Numberrealestateloansorlines（抵押贷款和房地产贷款数量，包括房屋净值信贷额度）
Numberoftime60-89dayspastduenotworse（过去两年内出现60-89天逾期但是没有发展得更坏的次数）
NumberofDependents（家庭中不包括自身的家属人数（配偶，子女等））

%matplotlib inline
import numpy as np
import pandas as pd
from sklearn.linear_model import LogisticRegression as LR
data = pd.read_csv(r".\rankingcard.csv",index_col=0)

第二步：数据清洗和特征工程

数据探索

#观察数据类型
data.head()#注意可以看到第一列为标签，剩下的10列为特征
 
#观察数据结构
data.shape#(150000, 11)
data.info()

#探索缺失值
data.isnull().sum()/data.shape[0]#得到缺失值的比例
#data.isnull().mean()#上一行代码的另一种形式书写

删除重复值和填补缺失值

删除重复值

data.drop_duplicates(inplace=True)
#删除后一定要恢复索引
data.index = range(data.shape[0])

对缺失比例比较少的用均值进行填补

data['NumberOfDependents'].fillna(int(data['NumberOfDependents'].mean()),inplace=True)

对缺失比例比较大，用随机森林进行填充

def fill_missing_rf(X,y,to_fill):

    """
    使用随机森林填补一个特征的缺失值的函数

    参数：
    X：要填补的特征矩阵
    y：完整的，没有缺失值的标签
    to_fill：字符串，要填补的那一列的名称
    """

    #构建我们的新特征矩阵和新标签
    df = X.copy()
    fill = df.loc[:,to_fill]
    df = pd.concat([df.loc[:,df.columns != to_fill],pd.DataFrame(y)],axis=1)

    # 找出我们的训练集和测试集
    Ytrain = fill[fill.notnull()]
    Ytest = fill[fill.isnull()]
    Xtrain = df.iloc[Ytrain.index,:]
    Xtest = df.iloc[Ytest.index,:]

    #用随机森林回归来填补缺失值
    from sklearn.ensemble import RandomForestRegressor as rfr
    rfr = rfr(n_estimators=100)
    rfr = rfr.fit(Xtrain, Ytrain)
    Ypredict = rfr.predict(Xtest)

    return Ypredict
    X = data.iloc[:,1:]
y = data["SeriousDlqin2yrs"]#y = data.iloc[:,0]
X.shape#(149391, 10)

#=====[TIME WARNING:1 min]=====#
y_pred = fill_missing_rf(X,y,"MonthlyIncome")

#注意可以通过以下代码检验数据是否数量相同
# y_pred.shape ==  data.loc[data.loc[:,"MonthlyIncome"].isnull(),"MonthlyIncome"].shape

#确认我们的结果合理之后，我们就可以将数据覆盖了
data.loc[data.loc[:,"MonthlyIncome"].isnull(),"MonthlyIncome"] = y_pred

data.info()

异常值处理

#异常值也被我们观察到，年龄的最小值居然有0，这不符合银行的业务需求，即便是儿童账户也要至少8岁，我们可以
# 查看一下年龄为0的人有多少
(data["age"] == 0).sum()
#发现只有一个人年龄为0，可以判断这肯定是录入失误造成的，可以当成是缺失值来处理，直接删除掉这个样本
data = data[data["age"] != 0]

样本不均衡的问题

通过上采样来解决样本不均衡的问题，让正负样本变成1：1

#探索标签的分布
X = data.iloc[:,1:]
y = data.iloc[:,0]
 
y.value_counts()#查看每一类别值得数据量，查看样本是否均衡
 
n_sample = X.shape[0]
 
n_1_sample = y.value_counts()[1]
n_0_sample = y.value_counts()[0]
 
print('样本个数：{}; 1占{:.2%}; 0占{:.2%}'.format(n_sample,n_1_sample/n_sample,n_0_sample/n_sample))
#样本个数：149165; 1占6.62%; 0占93.38%


#如果报错，就在prompt安装：pip install imblearn
import imblearn
#imblearn是专门用来处理不平衡数据集的库，在处理样本不均衡问题中性能高过sklearn很多
#imblearn里面也是一个个的类，也需要进行实例化，fit拟合，和sklearn用法相似
 
from imblearn.over_sampling import SMOTE
 
sm = SMOTE(random_state=42) #实例化
X,y = sm.fit_sample(X,y)
 
n_sample_ = X.shape[0]#278584

pd.Series(y).value_counts()
 
n_1_sample = pd.Series(y).value_counts()[1]
n_0_sample = pd.Series(y).value_counts()[0]
 
print('样本个数：{}; 1占{:.2%}; 0占{:.2%}'.format(n_sample_,n_1_sample/n_sample_,n_0_sample/n_sample_))
#样本个数：278584; 1占50.00%; 0占50.00%

划分训练集和测试集

from sklearn.model_selection import train_test_split
X = pd.DataFrame(X)
y = pd.DataFrame(y)
 
X_train, X_vali, Y_train, Y_vali = train_test_split(X,y,test_size=0.3,random_state=420)
model_data = pd.concat([Y_train, X_train], axis=1)#训练数据构建模型
model_data.index = range(model_data.shape[0])
model_data.columns = data.columns
 
vali_data = pd.concat([Y_vali, X_vali], axis=1)#验证集
vali_data.index = range(vali_data.shape[0])
vali_data.columns = data.columns
 
model_data.to_csv(r".\model_data.csv")#训练数据
vali_data.to_csv(r".\vali_data.csv")#验证数据

对变量进行分箱处理

分箱目标

分箱想要达成的效果：
我们希望不同属性的人有不同的分数，因此希望在同一个箱子内的人的属性是尽量相似的，而不同箱子的人的属性是尽量不同的，即业界常说的“组间差异大，组内差异小”。对于评分卡来说，就是说我们希望每一个箱子内的人违约概率是相似的，而不同箱子的人的违约概率差距很大，即WOE差距要大。这里可以使用卡方检验来对比两个箱子之间的相似性，如果两个箱子之间卡方减员的P值很大，则说明他们非常相似，那我们可以将这两个箱子合并为一个箱子。
基于这样的思想，我们总结出对特征进行分箱的步骤：
(1)首先将连续型变量分成一组数量较多的分类型变量（q组）
（2）确保每一组中都要包含两种类别的样本，否则IV值无法计算
（3）我们对相邻的组进行卡方检验，卡方检验的P值很大的组进行合并，直到数据中的组数小于设定的N箱为止
（4）我们让一个特征分别分成【2,3，4…20】箱，观察每个分箱个数下的IV值如何变化，找出最适合的分箱个数
（5）分箱完毕后，我们计算每个箱的WOE值，观察分箱结果

首先对各个特征做等频分箱，初始分箱个数为q，求出每个分箱中0和1的数量，可以对应求出woe和iv值，通过卡方检验

def graphforbestbin(DF, X, Y, n=5,q=20,graph=True):
    '''
    自动最优分箱函数，基于卡方检验的分箱

    参数：
    DF: 需要输入的数据
    X: 需要分箱的列名
    Y: 分箱数据对应的标签 Y 列名
    n: 保留分箱个数
    q: 初始分箱的个数
    graph: 是否要画出IV图像

    区间为前开后闭 (]

    '''
    
    DF = DF[[X,Y]].copy()

    DF["qcut"],bins = pd.qcut(DF[X], retbins=True, q=q,duplicates="drop")
    coount_y0 = DF.loc[DF[Y]==0].groupby(by="qcut").count()[Y]
    coount_y1 = DF.loc[DF[Y]==1].groupby(by="qcut").count()[Y]
    num_bins = [*zip(bins,bins[1:],coount_y0,coount_y1)]
	
	#确保每个分箱中有0和1
    for i in range(q):
        if 0 in num_bins[0][2:]:
            num_bins[0:2] = [(
                num_bins[0][0],
                num_bins[1][1],
                num_bins[0][2]+num_bins[1][2],
                num_bins[0][3]+num_bins[1][3])]
            continue

        for i in range(len(num_bins)):
            if 0 in num_bins[i][2:]:
                num_bins[i-1:i+1] = [(
                    num_bins[i-1][0],
                    num_bins[i][1],
                    num_bins[i-1][2]+num_bins[i][2],
                    num_bins[i-1][3]+num_bins[i][3])]
                break
        else:
            break

    def get_woe(num_bins):
        columns = ["min","max","count_0","count_1"]
        df = pd.DataFrame(num_bins,columns=columns)
        df["total"] = df.count_0 + df.count_1
        df["percentage"] = df.total / df.total.sum()
        df["bad_rate"] = df.count_1 / df.total
        df["good%"] = df.count_0/df.count_0.sum()
        df["bad%"] = df.count_1/df.count_1.sum()
        df["woe"] = np.log(df["good%"] / df["bad%"])
        return df

    def get_iv(df):
        rate = df["good%"] - df["bad%"]
        iv = np.sum(rate * df.woe)
        return iv

    IV = []
    axisx = []
    while len(num_bins) > n:
        pvs = []
        for i in range(len(num_bins)-1):
            x1 = num_bins[i][2:]
            x2 = num_bins[i+1][2:]
            pv = scipy.stats.chi2_contingency([x1,x2])[1]
            pvs.append(pv)

        i = pvs.index(max(pvs))
        num_bins[i:i+2] = [(
            num_bins[i][0],
            num_bins[i+1][1],
            num_bins[i][2]+num_bins[i+1][2],
            num_bins[i][3]+num_bins[i+1][3])]

        bins_df = pd.DataFrame(get_woe(num_bins))
        axisx.append(len(num_bins))
        IV.append(get_iv(bins_df))
        
    if graph:
        plt.figure()
        plt.plot(axisx,IV)
        plt.xticks(axisx)
        plt.xlabel("number of box")
        plt.ylabel("IV")
        plt.show()
    return bins_df

分出来的分箱所追求的效果
1.希望每组的bad_rate相差越大越好
2.woe差异越大越好，应该具有单调性，要么由正到负，要么由负到正，只能由一个转折过程
3.如果woe值大小变化由两个转折，比如呈现w型，证明分箱过程有问题
4.num_bins保留的信息越多越好

接下来利用上述函数求出每个特征在不同的分箱下的IV值，根据IV值下降得陡峭的进行决定最后的分箱个数

model_data.columns

for i in model_data.columns[1:-1]:
    print(i)
    graphforbestbin(model_data,i,"SeriousDlqin2yrs",n=2,q=20)

分箱结果，由于有些是分不出来的，需要进行手动分箱,其他根据自动分箱

其实可以发现，不是每个特征都可以自动完成分箱这么多，像家人数量，就无法分出20组，于是可以将可以分箱的特征放出来单独分组，不能自动分箱的变量自己观察。
一样的步骤，
1.对于能够进行自动分箱的确定好其分箱个数，把其通过等频分箱和卡方检验求出的对应箱子上下限求出来，同时将手动分箱的，手动进行分箱取值，在这里都设为3个分箱，同时将np.inf替换最大值，用-np.inf替换最小值，确保模型换了个更大的值进来也能找到相应的箱子可以选择。
2.找出每个分箱中的1和0 的数量，求出每个特征对应的woe值，将每个特征对应分箱的woe值映射在特征矩阵中
3.分出训练集和测试集，准备模型拟合训练

auto_col_bins = {"RevolvingUtilizationOfUnsecuredLines":6,
                "age":5,
                "DebtRatio":4,
                "MonthlyIncome":3,
                "NumberOfOpenCreditLinesAndLoans":5}
 
#不能使用自动分箱的变量
hand_bins = {"NumberOfTime30-59DaysPastDueNotWorse":[0,1,2,13]
            ,"NumberOfTimes90DaysLate":[0,1,2,17]
            ,"NumberRealEstateLoansOrLines":[0,1,2,4,54]
            ,"NumberOfTime60-89DaysPastDueNotWorse":[0,1,2,8]
            ,"NumberOfDependents":[0,1,2,3]}
 
#保证区间覆盖使用 np.inf替换最大值，用-np.inf替换最小值 
#原因：比如一些新的值出现，例如家庭人数为30，以前没出现过，改成范围为极大值之后，这些新值就都能分到箱里边了
hand_bins = {k:[-np.inf,*v[:-1],np.inf] for k,v in hand_bins.items()}

bins_of_col = {}
 
# 生成自动分箱的分箱区间和分箱后的 IV 值
 
for col in auto_col_bins:
    bins_df = graphforbestbin(model_data,col
                             ,"SeriousDlqin2yrs"
                             ,n=auto_col_bins[col]
                             #使用字典的性质来取出每个特征所对应的箱的数量
                             ,q=20
                             ,graph=False)
    bins_list = sorted(set(bins_df["min"]).union(bins_df["max"]))
    #保证区间覆盖使用 np.inf 替换最大值 -np.inf 替换最小值
    bins_list[0],bins_list[-1] = -np.inf,np.inf
    bins_of_col[col] = bins_list
    
#合并手动分箱数据    
bins_of_col.update(hand_bins)
 
data = model_data.copy()
 
#函数pd.cut，可以根据已知的分箱间隔把数据分箱
#参数为 pd.cut(数据，以列表表示的分箱间隔)
data = data[["age","SeriousDlqin2yrs"]].copy()
 
data["cut"] = pd.cut(data["age"],[-np.inf, 48.49986200790144, 58.757170160044694, 64.0, 74.0, np.inf])
 
data.head()

#将数据按分箱结果聚合，并取出其中的标签值
data.groupby("cut")["SeriousDlqin2yrs"].value_counts()
 
#使用unstack()来将树状结构变成表状结构
data.groupby("cut")["SeriousDlqin2yrs"].value_counts().unstack()
 
bins_df = data.groupby("cut")["SeriousDlqin2yrs"].value_counts().unstack()
 
bins_df["woe"] = np.log((bins_df[0]/bins_df[0].sum())/(bins_df[1]/bins_df[1].sum()))

def get_woe(df,col,y,bins):
    df = df[[col,y]].copy()
    df["cut"] = pd.cut(df[col],bins)
    bins_df = df.groupby("cut")[y].value_counts().unstack()
    woe = bins_df["woe"] = np.log((bins_df[0]/bins_df[0].sum())/(bins_df[1]/bins_df[1].sum()))
    return woe
 
#将所有特征的WOE存储到字典当中
woeall = {}
for col in bins_of_col:
    woeall[col] = get_woe(model_data,col,"SeriousDlqin2yrs",bins_of_col[col])
    
#不希望覆盖掉原本的数据，创建一个新的DataFrame，索引和原始数据model_data一模一样
model_woe = pd.DataFrame(index=model_data.index)
 
#将原数据分箱后，按箱的结果把WOE结构用map函数映射到数据中
model_woe["age"] = pd.cut(model_data["age"],bins_of_col["age"]).map(woeall["age"])
 
#对所有特征操作可以写成：
for col in bins_of_col:
    model_woe[col] = pd.cut(model_data[col],bins_of_col[col]).map(woeall[col])
    
#将标签补充到数据中
model_woe["SeriousDlqin2yrs"] = model_data["SeriousDlqin2yrs"]
 
#这就是我们的建模数据了
model_woe.head()

vali_woe = pd.DataFrame(index=vali_data.index)
 
for col in bins_of_col:
    vali_woe[col] = pd.cut(vali_data[col],bins_of_col[col]).map(woeall[col])
vali_woe["SeriousDlqin2yrs"] = vali_data["SeriousDlqin2yrs"]
 
vali_X = vali_woe.iloc[:,:-1]
vali_y = vali_woe.iloc[:,-1]

第三步：模型开发

使用逻辑回归对模型进行拟合，先使用不加任何的参数进行训练

X = model_woe.iloc[:,:-1]
y = model_woe.iloc[:,-1]
 
from sklearn.linear_model import LogisticRegression as LR
 
lr = LR().fit(X,y)
lr.score(vali_X,vali_y)#0.8641356370249832

#结果为0.8641356370249832

对其进行学习曲线的绘制进行调参

逻辑回归主要的参数是C,penalty，solver和max_iter还有multi_class
C是正则化强度的倒数，C越小，损失函数会越小，模型对损失函数的惩罚越重，正则化的效力越强。防止过拟合
penalty 可以输入l1和l2来指定使用哪一种正则化方式，不填写默认L2
注意若选择L1正则化，参数solver仅能够使用求解方式‘liblinear’和‘saga’，若使用‘l2’正则化，参数solver中所有的求解方式都可以使用。
solver默认是liblinear”是二分类专用，也是现在的默认求解器
multiclass 输入ovr ，multinomial，auto来告知模型，我们要处理的分类问题的类型，默认是ovr
ovr 二分类
multinomial 多分类
auto根据标签类型自己选择

l1与l2正则化

在l1正则化在逐渐加强的过程中，携带信息量小的、对模型贡献不大的特征的参数，会比携带大量信息的、对模型由巨大贡献的特征的参数更快地变成0，所以L1正则化本质是一个特征选择的过程，掌管了参数的‘稀疏性’。一般来说，如果特征量很大，数据维度很高，我们会倾向于使用L1正则化。
相对的，L2正则化在加强的过程中，会尽量让每个特征对模型都有一些小的贡献，但携带信息少，对模型贡献不大的特征的参数非常接近于0。通常来说，如果我们的主要目的只是为了防止过拟合，选择L2正则化就足够了

c_1 = np.linspace(0.01,1,20)
c_2 = np.linspace(0.01,0.2,20)
 
score = []
for i in c_1: 
    lr = LR(solver='liblinear',C=i).fit(X,y) 
    score.append(lr.score(vali_X,vali_y))
plt.figure()
plt.plot(c_1,score)
plt.show()
 
lr.n_iter_#array([7], dtype=int32)
 
score = []
for i in [1,2,3,4,5,6]: 
    lr = LR(solver='liblinear',C=0.025,max_iter=i).fit(X,y)
    score.append(lr.score(vali_X,vali_y))
plt.figure()
plt.plot([1,2,3,4,5,6],score)
plt.show()

由于模型的特征只有10个，没有特别高维度，故采用L2正则化就行。for循环对C值和max_iter值根据C和max_iter的不同取值画学习曲线，发现准确度的上限也是逼近刚好没调参的时候，故不用进行太多的调参了。

第四步：模型检验与评估

ROC曲线

曲线落在了左上方，AUC为0.94，取得了比准确度更好的结果，对于评分卡这种更看准抓捕少数类样本的模型，取得的效果是相当不错了,召回率比较高

import scikitplot as skplt
 
#%%cmd
#pip install scikit-plot
 
vali_proba_df = pd.DataFrame(lr.predict_proba(vali_X))
skplt.metrics.plot_roc(vali_y, vali_proba_df,
                        plot_micro=False,figsize=(6,6),
                        plot_macro=False)

第五步：模型上线

建立评分卡模型

建模完毕，使用准确率和ROC曲线验证了模型的预测能力。接下来就是用逻辑回归转换为标准评分卡了，评分卡中的分数，由以下公式计算：

$$Score=A-B\ast log(odds)$$
其中A与B是常数，A叫做补偿，B叫做刻度，log（odds）代表了一个人违约的可能性。其实逻辑回归的结果取对数几率形式会得到${\theta }^{T}x$，即我们的参数*特征矩阵，所以log（odds）就是我们的参数，两个常数可以通过两个假设的分值带入公式求出，这两个假设分别是：
1、某个特定的违约概率下的预期分值
2.指定的违约概率翻倍的分数（PDO）
假如，假设对数几率为1/60时设定的特定分数为600，PDO=20，那么对数几率为1/30的分数就是620，
带上以上线性表达式，可以得到：

$$600=A-B\ast log（1/60）$$
$$620=A-B\ast log(1/30)$$
用numpy可以很容易求出A和B的值；

B = 20/np.log(2)
A = 600 + B*np.log(1/60)

有了A和B，分数就很容易得到，将截距作为log(odds)带入公式进行计算，而其他各个特征的分数，也是将系数带入进行计算

#得到基础分
base_score = A -B*lr.intercept_ 
#得到每个特征对应的分数
score_age = woeall['age'] *(-B*lr.coef_[0][1]) 

withopen(file,"w") asfdata:
    fdata.write("base_score,{}\n".format(base_score))fori,colinenumerate(X.columns):
    score = woeall[col] * (-B*lr.coef_[0][i])
    score.name = "Score"
    score.index.name = col
    score.to_csv(file,header=True,mode="a")

第六步：监测与报告

引用了菜菜老师的课程~