本次任务为特征工程。在任务开始前，简单了解一下什么是特征工程。

一、特征工程初步了解

特征工程（FeatureEngineering）特征工程是将原始数据转化成更好的表达问题本质的特征的过程，使得将这些特征运用到预测模型中能提高对不可见数据的模型预测精度。

举一个简单的例子：我们想以某城市的区为单位来预测房价，那么原始数据中的房屋所在经纬度数据就需要经过处理，将其根据经纬度归入某一区，这就是一次简单的特征工程，从原始数据创造了新的数据，而新的数据才能更好的完成本次任务。

当然，并非所有的属性都是特征，只有对因变量有影响的才是特征。某个特征可以由一个属性衍生而来，也可以由多个属性共同衍生而来。

特征工程在机器学习里有多重要呢？从业界广泛流传的一句话便可得知：

数据和特征决定了机器学习的上限,而模型和算法只是逼近这个上限而已。

二、Code

1、异常值处理

1.1 通过箱线图分析异常值
这种方法是利用箱型图的四分位距（IQR）对异常值进行检测，也叫Tukey‘s test。四分位距(IQR)就是上四分位与下四分位的差值。而我们通过IQR的x倍为标准，规定：超过上四分位+x倍IQR距离，或者下四分位-x倍IQR距离的点为异常值。

def box_plot_outliers(data_ser, box_scale):
    """
        利用箱线图去除异常值
        :param data_ser: 接收 pandas.Series 数据格式
        :param box_scale: 箱线图尺度，
        :return:
        """
    iqr = box_scale * (data_ser.quantile(0.75) - data_ser.quantile(0.25))
    val_low =  data_ser.quantile(0.25) - iqr
    val_up = data_ser.quantile(0.75) + iqr
    rule_low = (data_ser < val_low)
    rule_up = (data_ser > val_up)
    return (rule_low, rule_up), (val_low, val_up)

1.2 通过 3-Sigma分析异常值

这个原则有个条件：数据需要服从正态分布。在3∂原则下，异常值如超过3倍标准差，那么可以将其视为异常值。正负3∂的概率是99.7%，那么距离平均值3∂之外的值出现的概率为P(|x-u| > 3∂) <= 0.003，属于极个别的小概率事件。如果数据不服从正态分布，也可以用远离平均值的多少倍标准差来描述。

def three_sigma_outliers(data_ser):
    """
      进行3sigma异常值剔除
       :param data_ser: 接收 pandas.Series 数据格式
       :return: 
    """
    val_low = data_ser.mean() - 3 * data_ser.std()
    val_up = data_ser.mean() + 3 * data_ser.std()
    rule_low = (data_ser < val_low)
    rule_up = (data_ser > val_up)
    return (rule_low, rule_up), (val_low, val_up)

1.3 删除异常值

可根据输入参数，选择是箱线图分析法还是3-sigma分析法。

def outliers_proc(data, col_name, scale=3, method='box'):
    """
    用于清洗异常值，默认用 box_plot（scale=3）进行清洗
    :param data: 接收 pandas 数据格式
    :param col_name: pandas 列名
    :param scale: 尺度
    :param method: 处理方式 box：箱线图；3sigma：3sigma
    :return:
    """
    
    data_n = data.copy()
    data_series = data_n[col_name]
    
    if method == 'box':
        rule, value = box_plot_outliers(data_series, box_scale=scale)
    if method == '3sigma':
        rule, value = three_sigma_outliers(data_series)
    
    index = np.arange(data_series.shape[0])[rule[0] | rule[1]]
    print("Delete number is: {}".format(len(index)))
    
    data_n = data_n.drop(index)
    data_n = data_n.reset_index(drop=True)
    print("Now column number is: {}".format(data_n.shape[0]))
    
    index_low = np.arange(data_series.shape[0])[rule[0]]
    outliers_low = data_series.iloc[index_low]
    print("Description of data less than the lower bound is:")
    print(outliers_low.describe())
    
    index_up = np.arange(data_series.shape[0])[rule[1]]
    outliers_up = data_series.iloc[index_up]
    print("Description of data larger than the upper bound is:")
    print(outliers_up.describe())
    
    fig, ax = plt.subplots(1, 2, figsize=(10, 10))
    sns.boxplot(y=data[col_name], data=data, ax=ax[0])
    sns.boxplot(y=data_n[col_name], data=data_n, ax=ax[1])
    return data_n

1.4 Box-Cox变换

在我们进行数据分析时，遇到的数据往往不是呈正态分布的，而如果数据不是正态性的，那么在部分情况下会带来一些问题。比如某些模型的前提就是要求数据具有正态性（KNN、贝叶斯等），此外数据具有正态性可以在一定程度上提高机器学习的训练效果。解决这个问题的一个方法是使用Box-Cox转换将数据转换为正态。

from scipy import stats
convert_data = stats.boxcox(data)[0]

2、特征构造

2.1 合并训练数据和测试数据

# 训练集和测试集放在一起，方便构造特征
data_train['train'] = 1
data_test_a['train'] = 0
data = pd.concat([data_train, data_test_a], ignore_index=True)

根据常识简单判断一下，二手汽车的价格与什么有关。

• 使用时间
• 汽车行驶公里
• 汽车是否有尚未修复的损坏
• 品牌
• 销售地区
• 车型、车身类型、燃油类型、变速箱联合变量？

2.2 提取汽车使用时间

使用时间：data[‘creatDate’] - data[‘regDate’]，反应汽车使用时间，一般来说价格与使用时间成反比

# 不过要注意，数据里有时间出错的格式，所以我们需要 errors='coerce'
# errors：{'ignore'，'raise'，'coerce'}，默认为'raise'
# 如果为“ raise”，则无效的解析将引发异常
# 如果为“coerce”，则将无效解析设置为NaT
# 如果为“ ignore”，则无效的解析将返回输入

data['used_time'] = (pd.to_datetime(data['creatDate'], format='%Y%m%d', errors='coerce') -
                     pd.to_datetime(data['regDate'], format='%Y%m%d', errors='coerce')).dt.days

使用 data[‘used_time’].isnull().sum() 查看used_time空值情况，有15k个样本空数据，占比 15/200=7.5% 较大，不建议删除，直接留空值，XGBBoost之类的决策树本身可以处理空值。

2.3 提取城市信息

德国邮政编码由五个数字组成。通常，第一位数字代表大区，第二位数字代表地区，后三位数字代表具体送货地区

# 本例中地区编码为 4 位，应该省略了第一位（大区）
data['city'] = data['regionCode'].apply(lambda x : str(x)[:-3])

2.4 根据品牌统计销售情况

# 销售情况只看data_train数据集中的
data_brand = data[data['train']==1].groupby('brand')

all_info = {}
for brand, brand_data in data_brand:
    info = {}
    info['brand_price_max'] = brand_data.price.max()
    info['brand_price_median'] = brand_data.price.median()
    info['brand_price_min'] = brand_data.price.min()
    info['brand_price_sum'] = brand_data.price.sum()
    info['brand_price_std'] = brand_data.price.std()
    info['brand_price_average'] = round(brand_data.price.sum() / (len(brand_data) + 1), 2)
    all_info[brand] = info

brand_info = pd.DataFrame(all_info).T.reset_index().rename(columns={"index": "brand"})
data = pd.merge(data, brand_info, how='left', on='brand')

2.5 数据分桶

数据分桶是一种数据预处理技术，用于减少次要观察误差的影响，是一种将多个连续值分组为较少数量的“桶”的方法。

数据分桶的重要性：

• 对异常数据有比较好的鲁棒性．比如我们对年龄分桶，> 50为一个桶，那么如果出现age=200这种错误数据也不会对模型造成干扰。
• 在逻辑回归模型中，单个变量分箱之后每个箱有自己独立的权重，相当于给模型加入了非线性的能力，能够提升模型的表达能力，极大拟合．
• 缺失值也可以作为一类特殊的变量进行模型．
• 分箱之后相对于简单的one_hot编码而言能够降低模型的复杂度，提升模型运算速度，对后期生产上线较为友好．

数据分桶的分类：

常用的分桶方式：

• 等频分箱：pd.qcut()
• 等距分箱： pd.cut()
• 卡方分箱
• KS分箱：解决了卡方分箱慢的问题

以power为例，使用等距分箱（pd.cut）

pd.cut( x, bins, right=True, labels=None, retbins=False, precision=3, include_lowest=False, duplicates='raise', )

参数	含义
x	一维数组
bins	整数，标量序列或者间隔索引，是进行分组的依据。如果填入整数n，则表示将x中的数值分成等宽的n份（即每一组内的最大值与最小值之差约相等）；如果是标量序列，序列中的数值表示用来分档的分界值； 如果是间隔索引，“ bins”的间隔索引必须不重叠
right	布尔值，默认为True表示包含最右侧的数值
labels	数组或布尔值，可选指定分箱的标签。如果是数组，长度要与分箱个数一致，比如“ bins”=[1、2、3、4]表示（1,2]，（2,3],（3,4]一共3个区间，则labels的长度也就是标签的个数也要是3如果为False，则仅返回分箱的整数指示符，即x中的数据在第几个箱子里当bins是间隔索引时，将忽略此参数
retbins	是否显示分箱的分界值。默认为False，当bins取整数时可以设置retbins=True以显示分界值，得到划分后的区间
precision	整数，默认3，存储和显示分箱标签的精度
include_lowest	布尔值，表示区间的左边是开还是闭，默认为false，也就是不包含区间左边
duplicates	如果分箱临界值不唯一，则引发ValueError或丢弃非唯一

# 可以通过先看df[clo_name].describe()来看一下数据总体情况再决定分箱的间隔。
bins = [i*10 for i in range(31)]
data['power_bin'] = pd.cut(data['power'], bins, labels=False)
data[['power_bin', 'power']].head()

特征构造好后，有些原始的数据已经变成无用数据，可以删除，重新生成新的数据集。

data = data.drop(['creatDate', 'regDate', 'regionCode'], axis=1)
data.to_csv('data_for_tree.csv', index=0)

2.6 用于LR的特征构造

之前构造的数据集适应于树模型，但是不同模型对数据集要求不同，如果要适应LR（Logistic Regression）模型，需单独为其构造一份数据集。

• LR的优势：
  对观测样本的概率值输出
  实现简单高效
  多重共线性的问题可以通过L2正则化来应对
  大量的工业界解决方案
  支持online learning（个人补充）
• LR的劣势:
  特征空间太大时表现不太好
  对于大量的分类变量无能为力
  对于非线性特征需要做特征变换
  依赖所有的样本数据

• DT的优势：
  直观的决策过程
  能够处理非线性特征
  适合大量样本的数据
  对部分数据缺失不敏感
• DT的劣势：
  极易过拟合（使用RF可以一定程度防止过拟合，但是只要是模型就会过拟合！）
  无法输出score，只能给出直接的分类结果

2.6.1 异常值处理
在第1节中，我们已经对训练集的power异常值进行了处理，但是通过数据可视化可以看出，data[‘power’]的分布异常，这是因为我们并没有处理data_test_a的异常值。

data['power'].plot.hist()
data_train['power'].plot.hist()

对特征进行log变化，使数据在一定程度上可以符合正态分布。使用np.log(data+1)，+1的目的是为了防止数据等于0，而不能进行log变化。

data['power'] = np.log(data['power'] + 1) 

2.6.2 归一化

使用MinMaxScaler归一化：

from sklearn import preprocessing
min_max_scaler = preprocessing.MinMaxScaler()
data['power'] =min_max_scaler .fit_transform(data['power'].values.reshape(-1,1))
data['power'].plot.hist()

preprocessing.MinMaxScaler() .fit_transform() 的data为二维数组，所以要将data[‘power’]转换为矩阵。
不使用MinMaxScaler的最大最小值归一化方式：

data['power'] = ((data['power'] - np.min(data['power'])) / (np.max(data['power']) - np.min(data['power'])))

对km和刚刚构造的统计量特征做归一化

for col in ['kilometer', 'brand_amount', 'brand_price_average', 'brand_price_max','brand_price_median', 'brand_price_min', 'brand_price_std','brand_price_sum']:
    data[col] = min_max_scaler.fit_transform(data[col].values.reshape(-1,1))

2.6.3 将离散的类别特征转换成由数字代表的类别特征

categorical_features = [ 'model', 'brand', 'bodyType', 'fuelType', 'gearbox', 'notRepairedDamage', 'power_bin']
data = pd.get_dummies(data, columns=categorical_features)

pandas 下的 one hot encoder 及 pd.get_dummies() 与 sklearn.preprocessing 下的 OneHotEncoder 的区别:

• sklearn.preprocessing 下的 OneHotEncoder 不可以直接处理 string，如果数据集中的某些特征是string 类型的话，需要首先将其转换为 integers 类型
• pd.get_dummies()，则恰将 string 转换为integers 类型

对于机器学习，建议使用OneHotEncoder。举个例子：如果训练数据“颜色”这个变量有“红”，“黄”两个值，但是在测试数据的“颜色”变量除了“红”，“黄”还有“蓝”。这个时候使用get_dummies转化训练数据会生成新的column：红，黄。之后用这个数据训练模型。同时使用get_dummies转化测试数据会生成新的column：红，黄，蓝。之后应用训练好的模型进行测试数据的预测。这时就会出现错误。因为在模型训练的过程并没有“蓝”这一列。但是使用OneHotEncoder（handle_unknown=“ignore”）就可以解决上述问题。

2.6.4 生成csv文件

data.to_csv('data_for_lr.csv', index=0)

3、特征筛选

3.1 过滤式

• 原理：使用发散性或相关性指标对各个特征进行评分，选择分数大于阈值的特征或者选择前K个分数最大的特征。具体来说，计算每个特征的发散性，移除发散性小于阈值的特征/选择前k个分数最大的特征；计算每个特征与标签的相关性，移除相关性小于阈值的特征/选择前k个分数最大的特征。
• 特点：特征选择过程与学习器无关，相当于先对初始特征进行过滤，再用过滤后的特征训练学习器。过滤式特征选择法简单、易于运行、易于理解，通常对于理解数据有较好的效果，但对特征优化、提高模型泛化能力来说效果一般。
• 常用方法：方差选择法、卡方检验法、皮尔森相关系数法、互信息系数法

# 相关系数
numerical_features = ['price','power', 'kilometer', 'brand_amount', 'brand_price_average', 'brand_price_max','brand_price_median', 'brand_price_min', 'brand_price_std','brand_price_sum']
correlation = data[numerical_features].corr(method='spearman')['price']

sns.heatmap(correlation)

# 方差选择法适用于离散型特征，连续型特征需要须离散化后使用；方差较小的特征很少，方差选择法简单但不好用，一般作为特征选择的预处理步骤，先去掉方差较小的特征，然后使用其他特征选择方法选择特征
from sklearn.feature_selection import VarianceThreshold
selector= VarianceThreshold(threshold=0) # 实例化，不填参数默认方差为0
vt = selector.fit_transform(data[numerical_features]) # 获取删除不合格特征之后的新特征矩阵
vt.shape # 过滤后方差的shape

# 相关性过滤法~卡方检验 计算每个非负特征与标签之间的卡方值，每一个卡方值对应一个p值，用p值判断特征和标签之间的相关性。大众经验，p值选择0.05或0.1。如P<=0.05就说明两组数据相关。
skb = SelectKBest(chi2,k=2) # 选择k个最佳特征 
new_data = skb.fit_transform(data[data['train']==1][numerical_features],  data[data['train']==1]['price'].values)

# 皮尔森相关系数法法（相关性）
# 使用pearson系数作为特征评分标准，相关系数绝对值越大，相关性越强（相关系数越接近于1或-1时，相关性越强；相关系数越接近于0时，相关性越弱）。特点：皮尔森相关系数法能够衡量线性相关的特征集。

from scipy.stats import pearsonr                   # 计算皮尔森相关系数
from numpy import array

"""
# 保留topk特征，移除topk外特征
# 皮尔森相关系数（输入特征矩阵和目标向量，输出二元组（评分，P），二数组第i项为第i个特征的评分和p值
# topK个数
"""

skb = SelectKBest(lambda X, Y: tuple(map(tuple,array(list(map(lambda x:pearsonr(x, Y), X.T))).T)), k=3)
skb = skb.fit_transform(data[data['train']==1][numerical_features],  data[data['train']==1]['price'].values)

# 相关性过滤~互信息法
# 互信息法返回每个特征与标签之间的互信息量的估计，值越大越相关：0表现特征和标签完全独立
from sklearn.feature_selection import mutual_info_classif

# 返回每个特征与标签的互信息估计量
result=mutual_info_classif(data[data['train']==1][numerical_features],  data[data['train']==1]['price'].values) 
# 假设需要筛选出来互信息量估计量最大的前2个特征
x_new = SelectKBest(mutual_info_classif, k=2).fit_transform(x, y)

3.2 包裹式

• 包裹式从初始特征集合中不断的选择特征子集，训练学习器，根据学习器的性能来对子集进行评价，直到选择出最佳的子集。
• 包裹式特征选择直接针对给定学习器进行优化。 优点：从最终学习器的性能来看，包裹式比过滤式更好；
• 缺点：由于特征选择过程中需要多次训练学习器，因此包裹式特征选择的计算开销通常比过滤式特征选择要大得多。

sklearn.feature_selection.SequentialFeatureSelector(estimator, *, n_features_to_select=None, direction='forward', scoring=None, cv=5, n_jobs=None)

文档： https://scikit-learn.org/stable/modules/generated/sklearn.feature_selection.SequentialFeatureSelector.html

from mlxtend.feature_selection import SequentialFeatureSelector as SFS
from sklearn.linear_model import LinearRegression

sfs = SFS(LinearRegression(), k_features=10, forward=True, floating=False, scoring = 'r2', cv = 0)
x = data.drop(['price'], axis=1)
x = x.fillna(0)
y = data['price']
sfs.fit(x, y)
sfs.k_feature_names_ 

from mlxtend.plotting import plot_sequential_feature_selection as plot_sfs
fig = plot_sfs(sfs.get_metric_dict(), kind='std_dev')
plt.show()

3.3 嵌入式

下一章介绍，Lasso 回归和决策树可以完成嵌入式特征选择

三、总结

这一part下来，感觉有点懵懵的，总结一下非常有必要。

• 异常值处理：可通过箱线图或者3-sigma来分析，如果异常值占比不大可直接删除，否则填充nan，XGBBoost之类的决策树本身可以处理空值。
• 特征构造：
  • 数据分桶：pd.cut()、 pd.qcut()、卡方分桶（还未实现）
  • 特征归一化/标准化：sklearn.preprocessing提供了多种归一化函数。可以log变换等、Box-Cox变换等使数据更接近正态分布。
  • 统计量特征构造：计数、最大值、最小值、求和、比理、中位数、方差等。
  • 时间特征：从原始时间数据提取与预测值更相关的时间特征。
  • 地理信息：包括分箱、分布编码等方法。
  • 特征组合、特征交叉
• 特征筛选：过滤式、包裹式、嵌入式

难点还是在于特征的构造，在于对数据的不敏感，对业务的不熟悉，对统计模型的不了解，都会造成很大的困难。也没有其他的方法，只能多接触，多做，慢慢来。