来源:https://blog.csdn.net/LuYi_WeiLin/article/details/85060190

    贷前准入环节流程图大致如下

    为什么需要建立评分卡?
    所有的模型一定是服务于业务的,那么业务上到底出现了什么问题,需要用到评分卡模型去解决呢?我们先从金融机构传统定价模式说起。

    我们知道银行将钱借出去是要收取利息的,那么收取多少利息是合理的呢?

    利息的本质是租金,银行借钱给客户,客户获得了一定时间内这笔钱的使用权,从而需要支付租金,就好像你租房子需要付房租一个道理。但是借钱和借房子不一样的地方在于,借钱出去有信用风险,借款人有可能到期无法偿还这笔钱,所以利息包括了这两方面的价值:

    利息=时间价值+风险价值

    其中,时间价值取决于市场供需关系(银行同业拆借利率),风险价值取决于客户的违约率。

    那么核心问题就是如何量化风险价值,银行的传统模式是对同一类客户群体根据历史表现计算一个平均的违约率。这种做法会存在什么问题呢?

    现在我们假设市场上有A和B两家银行,A采取传统定价(计算的平均违约率3%),B采取风险定价(对不同违约率的客户采取不同的定价),想想会发生什么现象?

    逆向选择(劣币驱逐良币):那些违约率低于3%的优质客户会从A银行转去B银行借钱,因为B银行根据这个客户的违约率给了更低的利息。从而A银行的优质客户会流失,导致A银行的平均违约率上升,比如上升到5%,A银行被迫上调贷款利率,这时候低于5%违约率的客户又会从A银行流失…

    所以银行需要对贷款客户进行风险评级,并且对不同风险等级的客户采取不同的定价,这样就可以有效避免逆向选择问题,那么为什么由模型而不是人来完成这一任务呢?

    稳定:人工判断具有一定的随意性和波动性,没有模型稳定

    效率:模型的判断效率远远高于人工判断

    成本:模型的边际成本为0,也就是开发好模型后,再多的客户数量也不会增加成本。

    评分卡的建立步骤
    数据获取,包括获取存量客户及潜在客户的数据。存量客户是指已经在证券公司开展相关融资类业务的客户,包括个人客户和机构客户;潜在客户是指未来拟在证券公司开展相关融资类业务的客户,主要包括机构客户,这也是解决证券业样本较少的常用方法,这些潜在机构客户包括上市公司、公开发行债券的发债主体、新三板上市公司、区域股权交易中心挂牌公司、非标融资机构等。

    数据预处理,主要工作包括数据清洗、缺失值处理、异常值处理,主要是为了将获取的原始数据转化为可用作模型开发的格式化数据。
    探索性数据分析,该步骤主要是获取样本总体的大概情况,描述样本总体情况的指标主要有直方图、箱形图等。
    变量选择,该步骤主要是通过统计学的方法,筛选出对违约状态影响最显著的指标。主要有单变量特征选择方法和基于机器学习模型的方法 。
    模型开发,该步骤主要包括变量分段、变量的WOE(证据权重)变换和逻辑回归估算三部分。
    模型评估,该步骤主要是评估模型的区分能力、预测能力、稳定性,并形成模型评估报告,得出模型是否可以使用的结论。
    信用评分,根据逻辑回归的系数和WOE等确定信用评分的方法。将Logistic模型转换为标准评分的形式。
    建立评分系统,根据信用评分方法,建立自动信用评分系统。
    大家有不懂的可以去看查看我的相关博客,申请评分卡(A卡)共有六篇相关博客,代码附上,数据集在我的博客资料可以下载

    给出数据集的字典
    字段 名称
    member_id ID
    loan_amnt 申请额度
    term 产品期限
    int_rate 利率
    emp_length 工作期限
    home_ownership 是否有自有住宅
    annual_inc 年收入
    verification_status 收入核验状态
    desc 描述
    purpose 贷款目的
    title 贷款目的描述
    zip_code 联系地址邮政编码
    addr_state 联系地址所属州
    delinq_2yrs 申贷日期前2年逾期次数
    inq_last_6mths 申请日前6个月咨询次数
    mths_since_last_delinq 上次逾期距今月份数
    mths_since_last_record 上次登记公众记录距今的月份数
    open_acc 征信局中记录的信用产品数
    pub_rec 公众不良记录数
    total_acc 正在使用的信用产品数
    pub_rec_bankruptcies 公众破产记录数
    earliest_cr_line 第一次借贷时间
    loan_status 贷款状态—目标变量

    代码:

    1. # -*- coding:utf-8 -*-
    2. import pandas as pd
    3. import re
    4. import time
    5. import datetime
    6. import numpy as np
    7. import pickle
    8. from dateutil.relativedelta import relativedelta
    9. from sklearn.model_selection import train_test_split
    10. import matplotlib.pyplot as plt
    11. import seaborn as sns
    12. from statsmodels.stats.outliers_influence import variance_inflation_factor
    13. from sklearn.linear_model import LogisticRegressionCV
    14. import statsmodels.api as sm
    15. from sklearn.ensemble import RandomForestClassifier
    16. from numpy import log
    17. from sklearn.metrics import roc_auc_score
    18. '''
    19. 作者:小象学院
    20. 时间:20190224
    21. '''
    22. #当连续变量的初始取值集合太多时(>100),我们先对其进行初步划分
    23. def SplitData(df, col, numOfSplit, special_attribute=[]):
    24. '''
    25. :param df: 按照col排序后的数据集
    26. :param col: 待分箱的变量
    27. :param numOfSplit: 切分的组别数
    28. :param special_attribute: 在切分数据集的时候,某些特殊值需要排除在外
    29. :return: 在原数据集上增加一列,把原始细粒度的col重新划分成粗粒度的值,便于分箱中的合并处理
    30. '''
    31. df2 = df.copy()
    32. if special_attribute != []:
    33. df2 = df.loc[~df[col].isin(special_attribute)]
    34. N = df2.shape[0]
    35. n = N//numOfSplit #每组样本数
    36. splitPointIndex = [i*n for i in range(1,numOfSplit)] #分割点的下标
    37. rawValues = sorted(list(df2[col])) #对取值进行排序
    38. splitPoint = [rawValues[i] for i in splitPointIndex] #分割点的取值
    39. splitPoint = sorted(list(set(splitPoint)))
    40. return splitPoint
    41. def Chi2(df, total_col, bad_col, overallRate):
    42. '''
    43. :param df: 包含全部样本总计与坏样本总计的数据框
    44. :param total_col: 全部样本的个数
    45. :param bad_col: 坏样本的个数
    46. :param overallRate: 全体样本的坏样本占比
    47. :return: 卡方值
    48. '''
    49. df2 = df.copy()
    50. # 期望坏样本个数=全部样本个数*平均坏样本占比
    51. df2['expected'] = df[total_col].apply(lambda x: x*overallRate)
    52. combined = zip(df2['expected'], df2[bad_col])
    53. chi = [(i[0]-i[1])**2/i[0] for i in combined]
    54. chi2 = sum(chi)
    55. return chi2
    56. def BinBadRate(df, col, target, grantRateIndicator=0):
    57. '''
    58. :param df: 需要计算好坏比率的数据集
    59. :param col: 需要计算好坏比率的特征
    60. :param target: 好坏标签
    61. :param grantRateIndicator: 1返回总体的坏样本率,0不返回
    62. :return: 每箱的坏样本率,以及总体的坏样本率(当grantRateIndicator==1时)
    63. '''
    64. total = df.groupby([col])[target].count()
    65. total = pd.DataFrame({'total': total})
    66. bad = df.groupby([col])[target].sum()
    67. bad = pd.DataFrame({'bad': bad})
    68. regroup = total.merge(bad, left_index=True, right_index=True, how='left')
    69. regroup.reset_index(level=0, inplace=True)
    70. regroup['bad_rate'] = regroup.apply(lambda x: x.bad * 1.0 / x.total, axis=1)
    71. dicts = dict(zip(regroup[col],regroup['bad_rate']))
    72. if grantRateIndicator==0:
    73. return (dicts, regroup)
    74. N = sum(regroup['total'])
    75. B = sum(regroup['bad'])
    76. overallRate = B * 1.0 / N
    77. return (dicts, regroup, overallRate)
    78. ### ChiMerge_MaxInterval: split the continuous variable using Chi-square value by specifying the max number of intervals
    79. def ChiMerge(df, col, target, max_interval=5,special_attribute=[],minBinPcnt=0):
    80. '''
    81. :param df: 包含目标变量与分箱属性的数据框
    82. :param col: 需要分箱的属性
    83. :param target: 目标变量,取值0或1
    84. :param max_interval: 最大分箱数。如果原始属性的取值个数低于该参数,不执行这段函数
    85. :param special_attribute: 不参与分箱的属性取值
    86. :param minBinPcnt:最小箱的占比,默认为0
    87. :return: 分箱结果
    88. '''
    89. colLevels = sorted(list(set(df[col])))
    90. N_distinct = len(colLevels)#不同的取值个数
    91. if N_distinct <= max_interval: #如果原始属性的取值个数低于max_interval,不执行这段函数
    92. print ("The number of original levels for {} is less than or equal to max intervals".format(col))
    93. return colLevels[:-1]
    94. else:
    95. if len(special_attribute)>=1:
    96. df1 = df.loc[df[col].isin(special_attribute)]
    97. df2 = df.loc[~df[col].isin(special_attribute)]
    98. else:
    99. df2 = df.copy()
    100. N_distinct = len(list(set(df2[col])))#该特征不同的取值
    101. # 步骤一: 通过col对数据集进行分组,求出每组的总样本数与坏样本数
    102. if N_distinct > 100:
    103. split_x = SplitData(df2, col, 100)
    104. df2['temp'] = df2[col].map(lambda x: AssignGroup(x, split_x))
    105. else:
    106. df2['temp'] = df2[col]
    107. # 总体bad rate将被用来计算expected bad count
    108. (binBadRate, regroup, overallRate) = BinBadRate(df2, 'temp', target, grantRateIndicator=1)
    109. # 首先,每个单独的属性值将被分为单独的一组
    110. # 对属性值进行排序,然后两两组别进行合并
    111. colLevels = sorted(list(set(df2['temp'])))
    112. groupIntervals = [[i] for i in colLevels]
    113. # 步骤二:建立循环,不断合并最优的相邻两个组别,直到:
    114. # 1,最终分裂出来的分箱数<=预设的最大分箱数
    115. # 2,每箱的占比不低于预设值(可选)
    116. # 3,每箱同时包含好坏样本
    117. # 如果有特殊属性,那么最终分裂出来的分箱数=预设的最大分箱数-特殊属性的个数
    118. split_intervals = max_interval - len(special_attribute)
    119. while (len(groupIntervals) > split_intervals): # 终止条件: 当前分箱数=预设的分箱数
    120. # 每次循环时, 计算合并相邻组别后的卡方值。具有最小卡方值的合并方案,是最优方案
    121. chisqList = []
    122. for k in range(len(groupIntervals)-1):
    123. temp_group = groupIntervals[k] + groupIntervals[k+1]
    124. df2b = regroup.loc[regroup['temp'].isin(temp_group)]
    125. chisq = Chi2(df2b, 'total', 'bad', overallRate)
    126. chisqList.append(chisq)
    127. best_comnbined = chisqList.index(min(chisqList))
    128. groupIntervals[best_comnbined] = groupIntervals[best_comnbined] + groupIntervals[best_comnbined+1]
    129. # after combining two intervals, we need to remove one of them
    130. groupIntervals.remove(groupIntervals[best_comnbined+1])
    131. groupIntervals = [sorted(i) for i in groupIntervals]
    132. cutOffPoints = [max(i) for i in groupIntervals[:-1]]
    133. # 检查是否有箱没有好或者坏样本。如果有,需要跟相邻的箱进行合并,直到每箱同时包含好坏样本
    134. groupedvalues = df2['temp'].apply(lambda x: AssignBin(x, cutOffPoints))
    135. df2['temp_Bin'] = groupedvalues
    136. (binBadRate,regroup) = BinBadRate(df2, 'temp_Bin', target)
    137. [minBadRate, maxBadRate] = [min(binBadRate.values()),max(binBadRate.values())]
    138. while minBadRate ==0 or maxBadRate == 1:
    139. # 找出全部为好/坏样本的箱
    140. indexForBad01 = regroup[regroup['bad_rate'].isin([0,1])].temp_Bin.tolist()
    141. bin=indexForBad01[0]
    142. # 如果是最后一箱,则需要和上一个箱进行合并,也就意味着分裂点cutOffPoints中的最后一个需要移除
    143. if bin == max(regroup.temp_Bin):
    144. cutOffPoints = cutOffPoints[:-1]
    145. # 如果是第一箱,则需要和下一个箱进行合并,也就意味着分裂点cutOffPoints中的第一个需要移除
    146. elif bin == min(regroup.temp_Bin):
    147. cutOffPoints = cutOffPoints[1:]
    148. # 如果是中间的某一箱,则需要和前后中的一个箱进行合并,依据是较小的卡方值
    149. else:
    150. # 和前一箱进行合并,并且计算卡方值
    151. currentIndex = list(regroup.temp_Bin).index(bin)
    152. prevIndex = list(regroup.temp_Bin)[currentIndex - 1]
    153. df3 = df2.loc[df2['temp_Bin'].isin([prevIndex, bin])]
    154. (binBadRate, df2b) = BinBadRate(df3, 'temp_Bin', target)
    155. chisq1 = Chi2(df2b, 'total', 'bad', overallRate)
    156. # 和后一箱进行合并,并且计算卡方值
    157. laterIndex = list(regroup.temp_Bin)[currentIndex + 1]
    158. df3b = df2.loc[df2['temp_Bin'].isin([laterIndex, bin])]
    159. (binBadRate, df2b) = BinBadRate(df3b, 'temp_Bin', target)
    160. chisq2 = Chi2(df2b, 'total', 'bad', overallRate)
    161. if chisq1 < chisq2:
    162. cutOffPoints.remove(cutOffPoints[currentIndex - 1])
    163. else:
    164. cutOffPoints.remove(cutOffPoints[currentIndex])
    165. # 完成合并之后,需要再次计算新的分箱准则下,每箱是否同时包含好坏样本
    166. groupedvalues = df2['temp'].apply(lambda x: AssignBin(x, cutOffPoints))
    167. df2['temp_Bin'] = groupedvalues
    168. (binBadRate, regroup) = BinBadRate(df2, 'temp_Bin', target)
    169. [minBadRate, maxBadRate] = [min(binBadRate.values()), max(binBadRate.values())]
    170. # 需要检查分箱后的最小占比
    171. if minBinPcnt > 0:
    172. groupedvalues = df2['temp'].apply(lambda x: AssignBin(x, cutOffPoints))
    173. df2['temp_Bin'] = groupedvalues
    174. valueCounts = groupedvalues.value_counts().to_frame()
    175. valueCounts['pcnt'] = valueCounts['temp'].apply(lambda x: x * 1.0 / N)
    176. valueCounts = valueCounts.sort_index()
    177. minPcnt = min(valueCounts['pcnt'])
    178. while minPcnt < minBinPcnt and len(cutOffPoints) > 2:
    179. # 找出占比最小的箱
    180. indexForMinPcnt = valueCounts[valueCounts['pcnt'] == minPcnt].index.tolist()[0]
    181. # 如果占比最小的箱是最后一箱,则需要和上一个箱进行合并,也就意味着分裂点cutOffPoints中的最后一个需要移除
    182. if indexForMinPcnt == max(valueCounts.index):
    183. cutOffPoints = cutOffPoints[:-1]
    184. # 如果占比最小的箱是第一箱,则需要和下一个箱进行合并,也就意味着分裂点cutOffPoints中的第一个需要移除
    185. elif indexForMinPcnt == min(valueCounts.index):
    186. cutOffPoints = cutOffPoints[1:]
    187. # 如果占比最小的箱是中间的某一箱,则需要和前后中的一个箱进行合并,依据是较小的卡方值
    188. else:
    189. # 和前一箱进行合并,并且计算卡方值
    190. currentIndex = list(valueCounts.index).index(indexForMinPcnt)
    191. prevIndex = list(valueCounts.index)[currentIndex - 1]
    192. df3 = df2.loc[df2['temp_Bin'].isin([prevIndex, indexForMinPcnt])]
    193. (binBadRate, df2b) = BinBadRate(df3, 'temp_Bin', target)
    194. chisq1 = Chi2(df2b, 'total', 'bad', overallRate)
    195. # 和后一箱进行合并,并且计算卡方值
    196. laterIndex = list(valueCounts.index)[currentIndex + 1]
    197. df3b = df2.loc[df2['temp_Bin'].isin([laterIndex, indexForMinPcnt])]
    198. (binBadRate, df2b) = BinBadRate(df3b, 'temp_Bin', target)
    199. chisq2 = Chi2(df2b, 'total', 'bad', overallRate)
    200. if chisq1 < chisq2:
    201. cutOffPoints.remove(cutOffPoints[currentIndex - 1])
    202. else:
    203. cutOffPoints.remove(cutOffPoints[currentIndex])
    204. cutOffPoints = special_attribute + cutOffPoints
    205. return cutOffPoints
    206. def UnsupervisedSplitBin(df,var,numOfSplit = 5, method = 'equal freq'):
    207. '''
    208. :param df: 数据集
    209. :param var: 需要分箱的变量。仅限数值型。
    210. :param numOfSplit: 需要分箱个数,默认是5
    211. :param method: 分箱方法,'equal freq':,默认是等频,否则是等距
    212. :return:
    213. '''
    214. if method == 'equal freq':
    215. N = df.shape[0]
    216. n = N / numOfSplit
    217. splitPointIndex = [i * n for i in range(1, numOfSplit)]
    218. rawValues = sorted(list(df[col]))
    219. splitPoint = [rawValues[i] for i in splitPointIndex]
    220. splitPoint = sorted(list(set(splitPoint)))
    221. return splitPoint
    222. else:
    223. var_max, var_min = max(df[var]), min(df[var])
    224. interval_len = (var_max - var_min)*1.0/numOfSplit
    225. splitPoint = [var_min + i*interval_len for i in range(1,numOfSplit)]
    226. return splitPoint
    227. def AssignGroup(x, bin):
    228. '''
    229. :param x: 某个变量的某个取值
    230. :param bin: 上述变量的分箱结果
    231. :return: x在分箱结果下的映射
    232. '''
    233. N = len(bin)
    234. if x<=min(bin):
    235. return min(bin)
    236. elif x>max(bin):
    237. return 10e10
    238. else:
    239. for i in range(N-1):
    240. if bin[i] < x <= bin[i+1]:
    241. return bin[i+1]
    242. def BadRateEncoding(df, col, target):
    243. '''
    244. :param df: dataframe containing feature and target
    245. :param col: the feature that needs to be encoded with bad rate, usually categorical type
    246. :param target: good/bad indicator
    247. :return: the assigned bad rate to encode the categorical feature
    248. '''
    249. regroup = BinBadRate(df, col, target, grantRateIndicator=0)[1]
    250. br_dict = regroup[[col,'bad_rate']].set_index([col]).to_dict(orient='index')
    251. for k, v in br_dict.items():
    252. br_dict[k] = v['bad_rate']
    253. badRateEnconding = df[col].map(lambda x: br_dict[x])
    254. return {'encoding':badRateEnconding, 'bad_rate':br_dict}
    255. def AssignBin(x, cutOffPoints,special_attribute=[]):
    256. '''
    257. :param x: 某个变量的某个取值
    258. :param cutOffPoints: 上述变量的分箱结果,用切分点表示
    259. :param special_attribute: 不参与分箱的特殊取值
    260. :return: 分箱后的对应的第几个箱,从0开始
    261. for example, if cutOffPoints = [10,20,30], if x = 7, return Bin 0. If x = 35, return Bin 3
    262. '''
    263. numBin = len(cutOffPoints) + 1 + len(special_attribute)
    264. if x in special_attribute:
    265. i = special_attribute.index(x)+1
    266. return 'Bin {}'.format(0-i)
    267. if x<=cutOffPoints[0]:
    268. return 'Bin 0'
    269. elif x > cutOffPoints[-1]:
    270. return 'Bin {}'.format(numBin-1)
    271. else:
    272. for i in range(0,numBin-1):
    273. if cutOffPoints[i] < x <= cutOffPoints[i+1]:
    274. return 'Bin {}'.format(i+1)
    275. def CalcWOE(df, col, target):
    276. '''
    277. :param df: 包含需要计算WOE的变量和目标变量
    278. :param col: 需要计算WOE、IV的变量,必须是分箱后的变量,或者不需要分箱的类别型变量
    279. :param target: 目标变量,0、1表示好、坏
    280. :return: 返回WOE和IV
    281. '''
    282. total = df.groupby([col])[target].count()
    283. total = pd.DataFrame({'total': total})
    284. bad = df.groupby([col])[target].sum()
    285. bad = pd.DataFrame({'bad': bad})
    286. regroup = total.merge(bad, left_index=True, right_index=True, how='left')
    287. regroup.reset_index(level=0, inplace=True)
    288. N = sum(regroup['total'])
    289. B = sum(regroup['bad'])
    290. regroup['good'] = regroup['total'] - regroup['bad']
    291. G = N - B
    292. regroup['bad_pcnt'] = regroup['bad'].map(lambda x: x*1.0/B)
    293. regroup['good_pcnt'] = regroup['good'].map(lambda x: x * 1.0 / G)
    294. regroup['WOE'] = regroup.apply(lambda x: np.log(x.good_pcnt*1.0/x.bad_pcnt),axis = 1)
    295. WOE_dict = regroup[[col,'WOE']].set_index(col).to_dict(orient='index')
    296. for k, v in WOE_dict.items():
    297. WOE_dict[k] = v['WOE']
    298. IV = regroup.apply(lambda x: (x.good_pcnt-x.bad_pcnt)*np.log(x.good_pcnt*1.0/x.bad_pcnt),axis = 1)
    299. IV = sum(IV)
    300. return {"WOE": WOE_dict, 'IV':IV}
    301. ## 判断某变量的坏样本率是否单调
    302. def BadRateMonotone(df, sortByVar, target,special_attribute = []):
    303. '''
    304. :param df: 包含检验坏样本率的变量,和目标变量
    305. :param sortByVar: 需要检验坏样本率的变量
    306. :param target: 目标变量,0、1表示好、坏
    307. :param special_attribute: 不参与检验的特殊值
    308. :return: 坏样本率单调与否
    309. '''
    310. df2 = df.loc[~df[sortByVar].isin(special_attribute)]
    311. if len(set(df2[sortByVar])) <= 2:
    312. return True
    313. regroup = BinBadRate(df2, sortByVar, target)[1]
    314. combined = zip(regroup['total'],regroup['bad'])
    315. badRate = [x[1]*1.0/x[0] for x in combined]
    316. badRateNotMonotone = [badRate[i]<badRate[i+1] and badRate[i] < badRate[i-1] or badRate[i]>badRate[i+1] and badRate[i] > badRate[i-1]
    317. for i in range(1,len(badRate)-1)]
    318. if True in badRateNotMonotone:
    319. return False
    320. else:
    321. return True
    322. def MergeBad0(df,col,target, direction='bad'):
    323. '''
    324. :param df: 包含检验0%或者100%坏样本率
    325. :param col: 分箱后的变量或者类别型变量。检验其中是否有一组或者多组没有坏样本或者没有好样本。如果是,则需要进行合并
    326. :param target: 目标变量,0、1表示好、坏
    327. :return: 合并方案,使得每个组里同时包含好坏样本
    328. '''
    329. regroup = BinBadRate(df, col, target)[1]
    330. if direction == 'bad':
    331. # 如果是合并0坏样本率的组,则跟最小的非0坏样本率的组进行合并
    332. regroup = regroup.sort_values(by = 'bad_rate')
    333. else:
    334. # 如果是合并0好样本样本率的组,则跟最小的非0好样本率的组进行合并
    335. regroup = regroup.sort_values(by='bad_rate',ascending=False)
    336. regroup.index = range(regroup.shape[0])
    337. col_regroup = [[i] for i in regroup[col]]
    338. del_index = []
    339. for i in range(regroup.shape[0]-1):
    340. col_regroup[i+1] = col_regroup[i] + col_regroup[i+1]
    341. del_index.append(i)
    342. if direction == 'bad':
    343. if regroup['bad_rate'][i+1] > 0:
    344. break
    345. else:
    346. if regroup['bad_rate'][i+1] < 1:
    347. break
    348. col_regroup2 = [col_regroup[i] for i in range(len(col_regroup)) if i not in del_index]
    349. newGroup = {}
    350. for i in range(len(col_regroup2)):
    351. for g2 in col_regroup2[i]:
    352. newGroup[g2] = 'Bin '+str(i)
    353. return newGroup
    354. def Prob2Score(prob, basePoint, PDO):
    355. #将概率转化成分数且为正整数
    356. y = np.log(prob/(1-prob))
    357. return int(basePoint+PDO/np.log(2)*(-y))
    358. ### 计算KS值
    359. def KS(df, score, target):
    360. '''
    361. :param df: 包含目标变量与预测值的数据集,dataframe
    362. :param score: 得分或者概率,str
    363. :param target: 目标变量,str
    364. :return: KS值
    365. '''
    366. total = df.groupby([score])[target].count()
    367. bad = df.groupby([score])[target].sum()
    368. all = pd.DataFrame({'total':total, 'bad':bad})
    369. all['good'] = all['total'] - all['bad']
    370. all[score] = all.index
    371. all = all.sort_values(by=score,ascending=False)
    372. all.index = range(len(all))
    373. all['badCumRate'] = all['bad'].cumsum() / all['bad'].sum()
    374. all['goodCumRate'] = all['good'].cumsum() / all['good'].sum()
    375. KS = all.apply(lambda x: x.badCumRate - x.goodCumRate, axis=1)
    376. return max(KS)
    377. def CareerYear(x):
    378. #对工作年限进行转换
    379. if str(x).find('nan') > -1:
    380. return -1
    381. elif str(x).find("10+")>-1: #将"10+years"转换成 11
    382. return 11
    383. elif str(x).find('< 1') > -1: #将"< 1 year"转换成 0
    384. return 0
    385. else:
    386. return int(re.sub("\D", "", x)) #其余数据,去掉"years"并转换成整数
    387. def DescExisting(x):
    388. #将desc变量转换成有记录和无记录两种
    389. if type(x).__name__ == 'float':
    390. return 'no desc'
    391. else:
    392. return 'desc'
    393. def ConvertDateStr(x):
    394. mth_dict = {'Jan': 1, 'Feb': 2, 'Mar': 3, 'Apr': 4, 'May': 5, 'Jun': 6, 'Jul': 7, 'Aug': 8, 'Sep': 9, 'Oct': 10,
    395. 'Nov': 11, 'Dec': 12}
    396. if str(x) == 'nan':
    397. return datetime.datetime.fromtimestamp(time.mktime(time.strptime('9900-1','%Y-%m')))
    398. #time.mktime 不能读取1970年之前的日期
    399. else:
    400. yr = int(x[4:6])
    401. if yr <=17:
    402. yr = 2000+yr
    403. else:
    404. yr = 1900 + yr
    405. mth = mth_dict[x[:3]]
    406. return datetime.datetime(yr,mth,1)
    407. def MonthGap(earlyDate, lateDate):
    408. if lateDate > earlyDate:
    409. gap = relativedelta(lateDate,earlyDate)
    410. yr = gap.years
    411. mth = gap.months
    412. return yr*12+mth
    413. else:
    414. return 0
    415. def MakeupMissing(x):
    416. if np.isnan(x):
    417. return -1
    418. else:
    419. return x
    420. # 数据预处理
    421. # 1,读入数据
    422. # 2,选择合适的建模样本
    423. # 3,数据集划分成训练集和测试集
    424. folderOfData = 'H:/信贷风控/资料/'
    425. allData = pd.read_csv(folderOfData + '数据集.csv',header = 0, encoding = 'latin1',engine ='python')
    426. allData['term'] = allData['term'].apply(lambda x: int(x.replace(' months','')))
    427. # 处理标签:Fully Paid是正常用户;Charged Off是违约用户
    428. allData['y'] = allData['loan_status'].map(lambda x: int(x == 'Charged Off'))
    429. '''
    430. 由于存在不同的贷款期限(term),申请评分卡模型评估的违约概率必须要在统一的期限中,且不宜太长,所以选取term=36months的行本
    431. '''
    432. allData1 = allData.loc[allData.term == 36]
    433. trainData, testData = train_test_split(allData1,test_size=0.4)
    434. #固化变量
    435. trainDataFile = open(folderOfData+'trainData.pkl','wb+')
    436. pickle.dump(trainData, trainDataFile)
    437. trainDataFile.close()
    438. testDataFile = open(folderOfData+'testData.pkl','wb+')
    439. pickle.dump(testData, testDataFile)
    440. testDataFile.close()
    441. '''
    442. 第一步:数据预处理,包括
    443. (1)数据清洗
    444. (2)格式转换
    445. (3)确实值填补
    446. '''
    447. # 将带%的百分比变为浮点数
    448. trainData['int_rate_clean'] = trainData['int_rate'].map(lambda x: float(x.replace('%',''))/100)
    449. # 将工作年限进行转化,否则影响排序,CareerYear为函数前面有定义
    450. trainData['emp_length_clean'] = trainData['emp_length'].map(CareerYear)
    451. #贷款描述,有写和没写状态区分,将desc的缺失作为一种状态,非缺失作为另一种状态
    452. trainData['desc_clean'] = trainData['desc'].map(DescExisting)
    453. #处理日期。earliest_cr_line(第一次借贷时间)的格式不统一,需要统一格式且转换成python的日期
    454. trainData['app_date_clean'] = trainData['issue_d'].map(lambda x: ConvertDateStr(x))
    455. trainData['earliest_cr_line_clean'] = trainData['earliest_cr_line'].map(lambda x: ConvertDateStr(x))
    456. # 处理mths_since_last_delinq(上次逾期距今月份数)。注意原始值中有0,所以用-1代替缺失
    457. trainData['mths_since_last_delinq_clean'] = trainData['mths_since_last_delinq'].map(lambda x:MakeupMissing(x))
    458. #上次登记公众记录距今的月份数(mths_since_last_record)
    459. trainData['mths_since_last_record_clean'] = trainData['mths_since_last_record'].map(lambda x:MakeupMissing(x))
    460. #公众破产记录数(pub_rec_bankruptcies)
    461. trainData['pub_rec_bankruptcies_clean'] = trainData['pub_rec_bankruptcies'].map(lambda x:MakeupMissing(x))
    462. '''
    463. 第二步:变量衍生
    464. '''
    465. # 考虑申请额度与收入的占比(limit_income)
    466. trainData['limit_income'] = trainData.apply(lambda x: x.loan_amnt / x.annual_inc, axis = 1)
    467. #考虑earliest_cr_line到申请日期的跨度,以月份记
    468. trainData['earliest_cr_to_app'] = trainData.apply(lambda x: MonthGap(x.earliest_cr_line_clean,x.app_date_clean), axis = 1)
    469. '''
    470. 第三步:分箱,采用ChiMerge,要求分箱完之后:
    471. (1)不超过5箱
    472. (2)Bad Rate单调
    473. (3)每箱同时包含好坏样本
    474. (4)特殊值如-1,单独成一箱
    475. 连续型变量可直接分箱
    476. 类别型变量:
    477. (a)当取值较多时,先用bad rate编码,再用连续型分箱的方式进行分箱
    478. (b)当取值较少时:
    479. (b1)如果每种类别同时包含好坏样本,无需分箱
    480. (b2)如果有类别只包含好坏样本的一种,需要合并
    481. '''
    482. #对数据源变量进行划分
    483. #num_features 连续型变量
    484. num_features = ['int_rate_clean','emp_length_clean','annual_inc', 'dti', 'delinq_2yrs', 'earliest_cr_to_app','inq_last_6mths', \
    485. 'mths_since_last_record_clean', 'mths_since_last_delinq_clean','open_acc','pub_rec','total_acc','limit_income','earliest_cr_to_app']
    486. #cat_features 离散型变量
    487. cat_features = ['home_ownership', 'verification_status','desc_clean', 'purpose', 'zip_code','addr_state','pub_rec_bankruptcies_clean']
    488. #more_value_features存放类别型变量变量取值大于5的变量、less_value_features存放类别型变量变量取值小于5的变量
    489. more_value_features = []
    490. less_value_features = []
    491. # 第一步,检查类别型变量中,哪些变量取值超过5
    492. for var in cat_features:
    493. valueCounts = len(set(trainData[var]))
    494. print(valueCounts)
    495. if valueCounts > 5:
    496. more_value_features.append(var) #取值超过5的变量,需要bad rate编码,再用卡方分箱法进行分箱
    497. else:
    498. less_value_features.append(var)
    499. # (i)当取值<5时:如果每种类别同时包含好坏样本,无需分箱;如果有类别只包含好坏样本的一种,需要合并
    500. merge_bin_dict = {} #存放需要合并的变量,以及合并方法
    501. var_bin_list = [] #由于某个取值没有好或者坏样本而需要合并的变量
    502. for col in less_value_features:
    503. binBadRate = BinBadRate(trainData, col, 'y')[0]
    504. if min(binBadRate.values()) == 0 : #由于某个取值没有坏样本而进行合并
    505. print('{} need to be combined due to 0 bad rate'.format(col))
    506. combine_bin = MergeBad0(trainData, col, 'y',direction = 'bad')
    507. merge_bin_dict[col] = combine_bin
    508. newVar = col + '_Bin'
    509. trainData[newVar] = trainData[col].map(combine_bin)
    510. var_bin_list.append(newVar)
    511. if max(binBadRate.values()) == 1: #由于某个取值没有好样本而进行合并
    512. print('{} need to be combined due to 0 good rate'.format(col))
    513. combine_bin = MergeBad0(trainData, col, 'y',direction = 'good')
    514. merge_bin_dict[col] = combine_bin
    515. newVar = col + '_Bin'
    516. trainData[newVar] = trainData[col].map(combine_bin)
    517. var_bin_list.append(newVar)
    518. #保存merge_bin_dict
    519. file1 = open(folderOfData+'merge_bin_dict.pkl','wb+')
    520. pickle.dump(merge_bin_dict,file1)
    521. file1.close()
    522. #less_value_features里剩下不需要合并的变量
    523. less_value_features = [i for i in less_value_features if i + '_Bin' not in var_bin_list]
    524. # (ii)当取值>5时:用bad rate进行编码,放入连续型变量里
    525. br_encoding_dict = {} #记录按照bad rate进行编码的变量,及编码方式
    526. for col in more_value_features:
    527. br_encoding = BadRateEncoding(trainData, col, 'y')
    528. trainData[col+'_br_encoding'] = br_encoding['encoding']
    529. br_encoding_dict[col] = br_encoding['bad_rate']
    530. #bad rate进行编码后,就可以放入连续型变量中(因为使用方法都一样)
    531. num_features.append(col+'_br_encoding')
    532. file2 = open(folderOfData+'br_encoding_dict.pkl','wb+')
    533. pickle.dump(br_encoding_dict,file2)
    534. file2.close()
    535. # (iii)对连续型变量进行分箱,包括(ii)中的变量
    536. continous_merged_dict = {}
    537. for col in num_features:
    538. print("{} is in processing".format(col))
    539. if -1 not in set(trainData[col]): #-1会当成特殊值处理。如果没有-1,则所有取值都参与分箱
    540. max_interval = 5 #分箱后的最多的箱数
    541. cutOff = ChiMerge(trainData, col, 'y', max_interval=max_interval,special_attribute=[],minBinPcnt=0)
    542. trainData[col+'_Bin'] = trainData[col].map(lambda x: AssignBin(x, cutOff,special_attribute=[]))
    543. monotone = BadRateMonotone(trainData, col+'_Bin', 'y') # 检验分箱后的单调性是否满足
    544. while(not monotone):
    545. # 检验分箱后的单调性是否满足。如果不满足,则缩减分箱的个数。
    546. max_interval -= 1
    547. cutOff = ChiMerge(trainData, col, 'y', max_interval=max_interval, special_attribute=[],
    548. minBinPcnt=0)
    549. trainData[col + '_Bin'] = trainData[col].map(lambda x: AssignBin(x, cutOff, special_attribute=[]))
    550. if max_interval == 2:
    551. # 当分箱数为2时,必然单调
    552. break
    553. monotone = BadRateMonotone(trainData, col + '_Bin', 'y')
    554. newVar = col + '_Bin'
    555. trainData[newVar] = trainData[col].map(lambda x: AssignBin(x, cutOff, special_attribute=[]))
    556. var_bin_list.append(newVar)
    557. else:
    558. max_interval = 5
    559. # 如果有-1,则除去-1后,其他取值参与分箱
    560. cutOff = ChiMerge(trainData, col, 'y', max_interval=max_interval, special_attribute=[-1],
    561. minBinPcnt=0)
    562. trainData[col + '_Bin'] = trainData[col].map(lambda x: AssignBin(x, cutOff, special_attribute=[-1]))
    563. monotone = BadRateMonotone(trainData, col + '_Bin', 'y',['Bin -1'])
    564. while (not monotone):
    565. max_interval -= 1
    566. # 如果有-1,-1的bad rate不参与单调性检验
    567. cutOff = ChiMerge(trainData, col, 'y', max_interval=max_interval, special_attribute=[-1],
    568. minBinPcnt=0)
    569. trainData[col + '_Bin'] = trainData[col].map(lambda x: AssignBin(x, cutOff, special_attribute=[-1]))
    570. if max_interval == 3:
    571. # 当分箱数为3-1=2时,必然单调
    572. break
    573. monotone = BadRateMonotone(trainData, col + '_Bin', 'y',['Bin -1'])
    574. newVar = col + '_Bin'
    575. trainData[newVar] = trainData[col].map(lambda x: AssignBin(x, cutOff, special_attribute=[-1]))
    576. var_bin_list.append(newVar)
    577. continous_merged_dict[col] = cutOff
    578. file3 = open(folderOfData+'continous_merged_dict.pkl','wb+')
    579. pickle.dump(continous_merged_dict,file3)
    580. file3.close()
    581. '''
    582. 第四步:WOE编码、计算IV
    583. '''
    584. WOE_dict = {}
    585. IV_dict = {}
    586. # 分箱后的变量进行编码,包括:
    587. # 1,初始取值个数小于5,且不需要合并的类别型变量。存放在less_value_features中
    588. # 2,初始取值个数小于5,需要合并的类别型变量。合并后新的变量存放在var_bin_list中
    589. # 3,初始取值个数超过5,需要合并的类别型变量。合并后新的变量存放在var_bin_list中
    590. # 4,连续变量。分箱后新的变量存放在var_bin_list中
    591. all_var = var_bin_list + less_value_features
    592. for var in all_var:
    593. woe_iv = CalcWOE(trainData, var, 'y')
    594. WOE_dict[var] = woe_iv['WOE']
    595. IV_dict[var] = woe_iv['IV']
    596. file4 = open(folderOfData+'WOE_dict.pkl','wb+')
    597. pickle.dump(WOE_dict,file4)
    598. file4.close()
    599. #将变量IV值进行降序排列,方便后续挑选变量
    600. IV_dict_sorted = sorted(IV_dict.items(), key=lambda x: x[1], reverse=True)
    601. IV_values = [i[1] for i in IV_dict_sorted]
    602. IV_name = [i[0] for i in IV_dict_sorted]
    603. plt.title('feature IV')
    604. plt.bar(range(len(IV_values)),IV_values)
    605. '''
    606. 第五步:单变量分析和多变量分析,均基于WOE编码后的值。
    607. (1)选择IV高于0.01的变量
    608. (2)比较两两线性相关性。如果相关系数的绝对值高于阈值,剔除IV较低的一个
    609. '''
    610. #选取IV>0.01的变量
    611. high_IV = {k:v for k, v in IV_dict.items() if v >= 0.02}
    612. high_IV_sorted = sorted(high_IV.items(),key=lambda x:x[1],reverse=True)
    613. short_list = high_IV.keys()
    614. short_list_2 = []
    615. for var in short_list:
    616. newVar = var + '_WOE'
    617. trainData[newVar] = trainData[var].map(WOE_dict[var])
    618. short_list_2.append(newVar)
    619. #对于上一步的结果,计算相关系数矩阵,并画出热力图进行数据可视化
    620. trainDataWOE = trainData[short_list_2]
    621. f, ax = plt.subplots(figsize=(10, 8))
    622. corr = trainDataWOE.corr()
    623. sns.heatmap(corr, mask=np.zeros_like(corr, dtype=np.bool), cmap=sns.diverging_palette(220, 10, as_cmap=True),square=True, ax=ax)
    624. #两两间的线性相关性检验
    625. #1,将候选变量按照IV进行降序排列
    626. #2,计算第i和第i+1的变量的线性相关系数
    627. #3,对于系数超过阈值的两个变量,剔除IV较低的一个
    628. deleted_index = []
    629. cnt_vars = len(high_IV_sorted)
    630. for i in range(cnt_vars):
    631. if i in deleted_index:
    632. continue
    633. x1 = high_IV_sorted[i][0]+"_WOE"
    634. for j in range(cnt_vars):
    635. if i == j or j in deleted_index:
    636. continue
    637. y1 = high_IV_sorted[j][0]+"_WOE"
    638. roh = np.corrcoef(trainData[x1],trainData[y1])[0,1]
    639. if abs(roh)>0.7:
    640. x1_IV = high_IV_sorted[i][1]
    641. y1_IV = high_IV_sorted[j][1]
    642. if x1_IV > y1_IV:
    643. deleted_index.append(j)
    644. else:
    645. deleted_index.append(i)
    646. multi_analysis_vars_1 = [high_IV_sorted[i][0]+"_WOE" for i in range(cnt_vars) if i not in deleted_index]
    647. '''
    648. 多变量分析:VIF
    649. '''
    650. X = np.matrix(trainData[multi_analysis_vars_1])
    651. VIF_list = [variance_inflation_factor(X, i) for i in range(X.shape[1])]
    652. max_VIF = max(VIF_list)
    653. print(max_VIF)
    654. # 最大的VIF是1.32267733123(小于10),因此这一步认为没有多重共线性
    655. multi_analysis = multi_analysis_vars_1
    656. '''
    657. 第六步:逻辑回归模型。
    658. 要求:
    659. 1,变量显著
    660. 2,符号为负
    661. '''
    662. ### (1)将多变量分析的后变量带入LR模型中
    663. y = trainData['y']
    664. X = trainData[multi_analysis]
    665. X['intercept'] = [1]*X.shape[0]
    666. LR = sm.Logit(y, X).fit()
    667. #summary = LR.summary()
    668. pvals = LR.pvalues
    669. pvals = pvals.to_dict()
    670. ### 有些变量不显著,需要逐步剔除
    671. varLargeP = {k: v for k,v in pvals.items() if v >= 0.1}
    672. varLargeP = sorted(varLargeP.items(), key=lambda d:d[1], reverse = True)
    673. while(len(varLargeP) > 0 and len(multi_analysis) > 0):
    674. # 每次迭代中,剔除最不显著的变量,直到
    675. # (1) 剩余所有变量均显著
    676. # (2) 没有特征可选
    677. varMaxP = varLargeP[0][0]
    678. print(varMaxP)
    679. if varMaxP == 'intercept':
    680. print('the intercept is not significant!')
    681. break
    682. multi_analysis.remove(varMaxP)
    683. y = trainData['y']
    684. X = trainData[multi_analysis]
    685. X['intercept'] = [1] * X.shape[0]
    686. LR = sm.Logit(y, X).fit()
    687. #summary = LR.summary()
    688. pvals = LR.pvalues
    689. pvals = pvals.to_dict()
    690. varLargeP = {k: v for k, v in pvals.items() if v >= 0.1}
    691. varLargeP = sorted(varLargeP.items(), key=lambda d: d[1], reverse=True)
    692. #将模型保存
    693. saveModel =open(folderOfData+'LR_Model_Normal.pkl','wb+')
    694. pickle.dump(LR,saveModel)
    695. saveModel.close()
    696. def ModifyDf(x, new_value):
    697. if np.isnan(x):
    698. return new_value
    699. else:
    700. return x
    701. '''
    702. 将模型应用在测试数据集上
    703. '''
    704. testDataFile = open(folderOfData+'testData.pkl','rb+')
    705. testData = pickle.load((testDataFile))
    706. testDataFile.close()
    707. '''
    708. 第一步:完成数据预处理
    709. 在实际工作中,可以只清洗模型实际使用的字段
    710. '''
    711. # 将带%的百分比变为浮点数
    712. testData['int_rate_clean'] = testData['int_rate'].map(lambda x: float(x.replace('%',''))/100)
    713. # 将工作年限进行转化,否则影响排序
    714. testData['emp_length_clean'] = testData['emp_length'].map(CareerYear)
    715. # 将desc的缺失作为一种状态,非缺失作为另一种状态
    716. testData['desc_clean'] = testData['desc'].map(DescExisting)
    717. # 处理日期。earliest_cr_line的格式不统一,需要统一格式且转换成python的日期
    718. testData['app_date_clean'] = testData['issue_d'].map(lambda x: ConvertDateStr(x))
    719. testData['earliest_cr_line_clean'] = testData['earliest_cr_line'].map(lambda x: ConvertDateStr(x))
    720. # 处理mths_since_last_delinq。注意原始值中有0,所以用-1代替缺失
    721. testData['mths_since_last_delinq_clean'] = testData['mths_since_last_delinq'].map(lambda x:MakeupMissing(x))
    722. testData['mths_since_last_record_clean'] = testData['mths_since_last_record'].map(lambda x:MakeupMissing(x))
    723. testData['pub_rec_bankruptcies_clean'] = testData['pub_rec_bankruptcies'].map(lambda x:MakeupMissing(x))
    724. '''
    725. 第二步:变量衍生
    726. '''
    727. # 考虑申请额度与收入的占比
    728. testData['limit_income'] = testData.apply(lambda x: x.loan_amnt / x.annual_inc, axis = 1)
    729. # 考虑earliest_cr_line到申请日期的跨度,以月份记
    730. testData['earliest_cr_to_app'] = testData.apply(lambda x: MonthGap(x.earliest_cr_line_clean,x.app_date_clean), axis = 1)
    731. '''
    732. 第三步:分箱并代入WOE值
    733. '''
    734. modelFile =open(folderOfData+'LR_Model_Normal.pkl','rb+')
    735. LR = pickle.load(modelFile)
    736. modelFile.close()
    737. #对变量的处理只需针对入模变量即可
    738. var_in_model = list(LR.pvalues.index)
    739. var_in_model.remove('intercept')
    740. file1 = open(folderOfData+'merge_bin_dict.pkl','rb+')
    741. merge_bin_dict = pickle.load(file1)
    742. file1.close()
    743. file2 = open(folderOfData+'br_encoding_dict.pkl','rb+')
    744. br_encoding_dict = pickle.load(file2)
    745. file2.close()
    746. file3 = open(folderOfData+'continous_merged_dict.pkl','rb+')
    747. continous_merged_dict = pickle.load(file3)
    748. file3.close()
    749. file4 = open(folderOfData+'WOE_dict.pkl','rb+')
    750. WOE_dict = pickle.load(file4)
    751. file4.close()
    752. for var in var_in_model:
    753. var1 = var.replace('_Bin_WOE','')
    754. # 有些取值个数少、但是需要合并的变量
    755. if var1 in merge_bin_dict.keys():
    756. print("{} need to be regrouped".format(var1))
    757. testData[var1 + '_Bin'] = testData[var1].map(merge_bin_dict[var1])
    758. # 有些变量需要用bad rate进行编码
    759. if var1.find('_br_encoding')>-1:
    760. var2 =var1.replace('_br_encoding','')
    761. print("{} need to be encoded by bad rate".format(var2))
    762. testData[var1] = testData[var2].map(br_encoding_dict[var2])
    763. #需要注意的是,有可能在测试样中某些值没有出现在训练样本中,从而无法得出对应的bad rate是多少。故可以用最坏(即最大)的bad rate进行编码
    764. max_br = max(testData[var1])
    765. testData[var1] = testData[var1].map(lambda x: ModifyDf(x, max_br))
    766. #上述处理后,需要加上连续型变量一起进行分箱
    767. if -1 not in set(testData[var1]):
    768. testData[var1+'_Bin'] = testData[var1].map(lambda x: AssignBin(x, continous_merged_dict[var1]))
    769. else:
    770. testData[var1 + '_Bin'] = testData[var1].map(lambda x: AssignBin(x, continous_merged_dict[var1],[-1]))
    771. #WOE编码
    772. var3 = var.replace('_WOE','')
    773. testData[var] = testData[var3].map(WOE_dict[var3])
    774. '''
    775. 第四步:将WOE值代入LR模型,计算概率和分数
    776. '''
    777. testData['intercept'] = [1]*testData.shape[0]
    778. #预测数据集中,变量顺序需要和LR模型的变量顺序一致
    779. #例如在训练集里,变量在数据中的顺序是“负债比”在“借款目的”之前,对应地,在测试集里,“负债比”也要在“借款目的”之前
    780. testData2 = testData[list(LR.params.index)]
    781. testData['prob'] = LR.predict(testData2)
    782. #计算KS和AUC
    783. auc = roc_auc_score(testData['y'],testData['prob'])
    784. ks = KS(testData, 'prob', 'y')
    785. basePoint = 250
    786. PDO = 200
    787. testData['score'] = testData['prob'].map(lambda x:Prob2Score(x, basePoint, PDO))
    788. testData = testData.sort_values(by = 'score')
    789. #画出分布图
    790. plt.hist(testData['score'], 100)
    791. plt.xlabel('score')
    792. plt.ylabel('freq')
    793. plt.title('distribution')