一文介绍特征工程里的卡方分箱，附代码实现

初次接触变量分箱是在做评分卡模型的时候，SAS软件里有一段宏可以直接进行连续变量的最优分箱，但如果搬到Python的话，又如何实现同样或者说类似的操作呢，今天就在这里简单介绍一个办法——卡方分箱算法。

为了让大家更好理解这个算法，我先从基础的原理开始讲起。

一、什么是卡方分布

卡方分布(chi-square distribution, χ2-distribution)是概率统计里常用的一种概率分布，也是统计推断里应用最广泛的概率分布之一，在假设检验与置信区间的计算中经常能见到卡方分布的身影。

卡方分布的定义如下：

若k个独立的随机变量Z1, Z2,..., Zk 满足标准正态分布 N(0,1) , 则这k个随机变量的平方和：

为服从自由度为k的卡方分布，记作： 或者记作 。

图1：卡方概率密度函数

图2：卡方累计分布函数

二、什么是卡方检验

χ2检验是以χ2分布为基础的一种假设检验方法，主要用于分类变量之间的独立性检验。

其基本思想是根据样本数据推断总体的分布与期望分布是否有显著性差异，或者推断两个分类变量是否相关或者独立。

一般可以设原假设为 ：观察频数与期望频数没有差异，或者两个变量相互独立不相关。

实际应用中，我们先假设原假设成立，计算出卡方的值，卡方表示观察值与理论值间的偏离程度。

卡方值的计算公式为：

其中A为实际频数，E为期望频数。卡方值用于衡量实际值与理论值的差异程度，这也是卡方检验的核心思想。

卡方值包含了以下两个信息：

1.实际值与理论值偏差的绝对大小。

2.差异程度与理论值的相对大小。

上述计算的卡方值服从卡方分布。根据卡方分布，卡方统计量以及自由度，可以确定在原假设成立的情况下获得当前统计量以及更极端情况的概率p。如果p很小，说明观察值与理论值的偏离程度大，应该拒绝原假设。否则不能拒绝原假设。

三、什么是卡方分布表

横轴为p值，纵轴为自由度。

（自由度的概念：自由度k=(行数-1)*(列数-1)，详情见实例）

四、卡方检验实例

某医院对某种病症的患者使用了A，B两种不同的疗法，结果如表1，问两种疗法有无差别？

表1 两种疗法治疗卵巢癌的疗效比较

可以计算出各格内的期望频数。

   第1行1列：43×53/87=26.2

   第1行2列：43×34/87=16.8

   第2行1列：44×53/87=26.8

   第2行2列：4×34/87=17.2

先建立原假设：A、B两种疗法没有区别。根据卡方值的计算公式，计算：

算得卡方值=10.01。

得到卡方值以后，接下来需要查询卡方分布表（见上面?）来判断p值，从而做出接受或拒绝原假设的决定。

查表自由度为1,p=0.05的卡方值为3.841，而此例卡方值10.01>3.841，因此 p < 0.05，说明原假设在0.05的显著性水平下是可以拒绝的。也就是说，原假设不成立。

五、ChiMerge分箱算法

ChiMerge卡方分箱算法由Kerber于1992提出。

它主要包括两个阶段：初始化阶段和自底向上的合并阶段。

1、初始化阶段：

首先按照属性值的大小进行排序（对于非连续特征，需要先做数值转换，比如转为坏人率，然后排序），然后每个属性值单独作为一组。

2、合并阶段：

（1）对每一对相邻的组，计算卡方值。

（2）根据计算的卡方值，对其中最小的一对邻组合并为一组。

（3）不断重复（1），（2）直到计算出的卡方值都不低于事先设定的阈值，或者分组数达到一定的条件（如最小分组数5，最大分组数8）。

值得注意的是，小编之前发现有的实现方法在合并阶段，计算的并非相邻组的卡方值（只考虑在此两组内的样本，并计算期望频数），因为他们用整体样本来计算此相邻两组的期望频数。

六、Python代码实现

1.导入相关库

import numpy as np
from scipy.stats import chi
import pandas as pd
from pandas import DataFrame,Series
import scipy

2.计算卡方值

def chi3(arr):
    '''
    计算卡方值
    arr:频数统计表,二维numpy数组。
    '''
    assert(arr.ndim==2)
    #计算每行总频数
    R_N = arr.sum(axis=1)
    #每列总频数
    C_N = arr.sum(axis=0)
    #总频数
    N = arr.sum()
    # 计算期望频数 C_i * R_j / N。
    E = np.ones(arr.shape)* C_N / N
    E = (E.T * R_N).T
    square = (arr-E)**2 / E
    #期望频数为0时，做除数没有意义，不计入卡方值
    square[E==0] = 0
    #卡方值
    v = square.sum()
    return v

3.确定卡方分箱点

def chiMerge(df,col,target,max_groups=None,threshold=None):

    '''
    卡方分箱
    df: pandas dataframe数据集
    col: 需要分箱的变量名（数值型）
    target: 类标签
    max_groups: 最大分组数。
    threshold: 卡方阈值，如果未指定max_groups，默认使用置信度95%设置threshold。
    return: 包括各组的起始值的列表.
    '''

    freq_tab = pd.crosstab(df[col],df[target])

    #转成numpy数组用于计算。
    freq = freq_tab.values

    #初始分组切分点，每个变量值都是切分点。每组中只包含一个变量值.

    #分组区间是左闭右开的，如cutoffs = [1,2,3]，则表示区间 [1,2) , [2,3) ,[3,3+)。
    cutoffs = freq_tab.index.values

    #如果没有指定最大分组
    if max_groups is None:    
        #如果没有指定卡方阈值，就以95%的置信度（自由度为类数目-1）设定阈值。
        if threshold is None:
            #类数目
            cls_num = freq.shape[-1]
            threshold = chi2.isf(0.05,df= cls_num - 1)

    while True:
        minvalue = None
        minidx = None
        #从第1组开始，依次取两组计算卡方值，并判断是否小于当前最小的卡方
        for i in range(len(freq) - 1):
            v = chi3(freq[i:i+2])
            if minvalue is None or (minvalue > v): #小于当前最小卡方，更新最小值
                minvalue = v
                minidx = i

        #如果最小卡方值小于阈值，则合并最小卡方值的相邻两组，并继续循环
        if  (max_groups is not None and  max_groups< len(freq) ) or (threshold is not None and minvalue < threshold):
            #minidx后一行合并到minidx
            tmp  = freq[minidx] + freq[minidx+1]
            freq[minidx] = tmp
            #删除minidx后一行
            freq = np.delete(freq,minidx+1,0)
            #删除对应的切分点
            cutoffs = np.delete(cutoffs,minidx+1,0)

        else: #最小卡方值不小于阈值，停止合并。
            break
    return cutoffs

4.生成分组后的新变量

def value2group(x,cutoffs):

    '''
    将变量的值转换成相应的组。
    x: 需要转换到分组的值
    cutoffs: 各组的起始值。
    return: x对应的组，如group1。从group1开始。
    '''

    #切分点从小到大排序。
    cutoffs = sorted(cutoffs)
    num_groups = len(cutoffs)

    #异常情况：小于第一组的起始值。这里直接放到第一组。
    #异常值建议在分组之前先处理妥善。
    if x < cutoffs[0]:
        return 'group1'

    for i in range(1,num_groups):
        if cutoffs[i-1] <= x < cutoffs[i]:
            return 'group{}'.format(i)

    #最后一组，也可能会包括一些非常大的异常值。
    return 'group{}'.format(num_groups)

5.实现WOE 编码

def calWOE(df ,var ,target):

    '''
    计算WOE编码
    param df：数据集pandas.dataframe
    param var：已分组的列名，无缺失值
    param target：响应变量（0,1）
    return：编码字典
    '''
    eps = 0.000001  #避免除以0
    gbi = pd.crosstab(df[var],df[target]) + eps
    gb = df[target].value_counts() + eps
    gbri = gbi/gb
    gbri['woe'] = np.log(gbri[1]/gbri[0])
    return gbri['woe'].to_dict()

6.实现IV值计算

def calIV(df,var,target):

    '''
    计算IV值
    param df：数据集pandas.dataframe
    param var：已分组的列名，无缺失值
    param target：响应变量（0,1）
    return：IV值
    '''
    eps = 0.000001  #避免除以0
    gbi = pd.crosstab(df[var],df[target]) + eps
    gb = df[target].value_counts() + eps
    gbri = gbi/gb
    gbri['woe'] = np.log(gbri[1]/gbri[0])
    gbri['iv'] = (gbri[1] - gbri[0])*gbri['woe']
    return gbri['iv'].sum()

7.Python代码使用

cutoffs = chiMerge(data,'MAX_AMOUNT','target',max_groups=8)

'''
PS：绝大数情况下，会将缺失值NaN归类到最后一组，如果不想这么简单粗暴的，需要在最开始的时候对缺失值进行填充。
'''

data['MAX_AMOUNT_chigroup'] = data['MAX_AMOUNT'].apply(value2group,args=(cutoffs,))
r = data.loc[:,['MAX_、AMOUNT','MAX_AMOUNT_chigroup','target']]
woe_map = calWOE(data,'MAX_AMOUNT_chigroup','target')
iv = calIV(data,'MAX_AMOUNT_chigroup','target')
iv

Reference

[1] Python评分卡建模—卡方分箱（1）

[2] Python评分卡建模—卡方分箱（2）之代码实现

[3] python评分卡建模—实现WOE编码及IV值计算

（以上文章均来自“风控建模”公众号，作者为东东&Monica）

[4] 维基百科-卡方分布

http://wiki.mbalib.com/wiki/%E5%8D%A1%E6%96%B9%E5%88%86%E5%B8%83