天池项目——二手车预测（一）

origin from: datawhale

Baseline

import numpy as np
import pandas as pd
import matplotlib
import matplotlib.pyplot as plt
import seaborn as sns
import warnings
import time

Train_data=pd.read_csv('./used_car_train_20200313.csv', sep=' ')
TestB_data=pd.read_csv('./used_car_testB_20200421.csv', sep=' ')

Train_data.head()
Train_data.info()
Train_data.columns
TestB_data.info()

Train_data.describe()
TestB_data.describe()

#提取数值类型为数字的
numerical_cols=Train_data.select_dtypes(exclude="object").columns #选择类型object以外的，也就是int64和float64

#提取数值类型为分类的
categorical_cols=Train_data.select_dtypes(include="object").columns
print(categorical_cols)

#选择特征列
feature_cols=[col for col in numerical_cols if col not in  ['SaleID','name','regDate','creatDate','price','model','brand','regionCode','seller']]
feature_cols=[col for col in feature_cols if 'Type' not in col]

#提取特征列
X_data=Train_data[feature_cols]
Y_data=Train_data['price']

#标签列构造训练样本和测试样本
X_test=TestB_data[feature_cols]

#缺省值用-1填补
X_data=X_data.fillna(-1)
X_test=X_test.fillna(-1)

#XGB

#导入相应包

## 模型预测的
from sklearn import linear_model
from sklearn import preprocessing
from sklearn.svm import SVR
from sklearn.ensemble import RandomForestRegressor,GradientBoostingRegressor

## 数据降维处理的
from sklearn.decomposition import PCA,FastICA,FactorAnalysis,SparsePCA

import lightgbm as lgb
import xgboost as xgb

## 参数搜索和评价的
from sklearn.model_selection import GridSearchCV,cross_val_score,StratifiedKFold,train_test_split
from sklearn.metrics import mean_squared_error, mean_absolute_error

xgr=xgb.XGBRegressor(n_estimators=120, learning_rate=0.1, gamma=0, subsample=0.8, colsample_bytree=0.9, max_depth=7)

scores_train=[]
scores=[]

#5折

sk=StratifiedKFold(n_splits=5, shuffle=True, random_state=0)


for train_ind, val_ind in sk.split(X_data,Y_data):
    train_x=X_data.iloc[train_ind].values
    train_y=Y_data.iloc[train_ind]
    val_x=X_data.iloc[val_ind].values
    val_y=Y_data.iloc[val_ind]
    
    
    xgr.fit(train_x,train_y)
    pred_train_xgb=xgr.predict(train_x)
    pred_xgb=xgr.predict(val_x)
    
    score_train=mean_absolute_error(train_y,pred_train_xgb)
    scores_train.append(score_train)
    score=mean_absolute_error(val_y,pred_xgb)
    scores.append(score)
    
#定义xgb和lbg模型函数
def build_model_xgb(x_train,y_train):
    model=xgb.XGBRegressor(n_estimators=150,learning_rate=0.1, gamma=0, subsample=0.8,\
        colsample_bytree=0.9, max_depth=7)
    model.fit(x_train,y_train)
    return model

def build_model_lgb(x_train,y_train):
    estimator = lgb.LGBMRegressor(num_leaves=127,n_estimators = 150)
    param_grid = {
        'learning_rate': [0.01, 0.05, 0.1, 0.2],
    }
    gbm = GridSearchCV(estimator, param_grid)
    gbm.fit(x_train, y_train)
    return gbm

x_train,x_val,y_train,y_val=train_test_split(X_data,Y_data,test_size=0.3)

print('Train lgb...')
model_lgb = build_model_lgb(x_train,y_train)
val_lgb = model_lgb.predict(x_val)
MAE_lgb = mean_absolute_error(y_val,val_lgb)
print('MAE of val with lgb:',MAE_lgb)

print('Predict lgb...')
model_lgb_pre = build_model_lgb(X_data,Y_data)
subA_lgb = model_lgb_pre.predict(X_test)
print('Sta of Predict lgb:')
Sta_inf(subA_lgb)

print('Train xgb...')
model_xgb = build_model_xgb(x_train,y_train)
val_xgb = model_xgb.predict(x_val)
MAE_xgb = mean_absolute_error(y_val,val_xgb)
print('MAE of val with xgb:',MAE_xgb)

print('Predict xgb...')
model_xgb_pre = build_model_xgb(X_data,Y_data)
subA_xgb = model_xgb_pre.predict(X_test)
print('Sta of Predict xgb:')
Sta_inf(subA_xgb)

#两模型结果加权融合
val_Weighted=(1-MAE_lgb/(MAE_xgb+MAE_lgb))*val_lgb+(1-MAE_xgb/(MAE_xgb+MAE_lgb))*val_xgb
val_Weighted[val_Weighted<0]=10 # 由于我们发现预测的最小值有负数，而真实情况下，price为负是不存在的，由此我们进行对应的后修正
print('MAE of val with Weighted ensemble:',mean_absolute_error(y_val,val_Weighted))

MAE of val with Weighted ensemble: 684.8596000379389

sub_Weighted = (1-MAE_lgb/(MAE_xgb+MAE_lgb))*subA_lgb+(1-MAE_xgb/(MAE_xgb+MAE_lgb))*subA_xgb

## 查看预测值的统计进行
plt.hist(Y_data)
plt.show()
plt.close()

数据分析

数据分析主要分为几个步骤：

1. 导入，读取并了解数据大致信息，判断缺失值和异常值并进行丢弃/替换补值处理

2. 查看数据大致形态如分布、峰度、偏值等，必要的话进行变换

3. 选取要进行预测的数据类型和相关列

4. 对数据进行相关性分析

导入，读取并了解数据大致信息，判断缺失值和异常值并进行丢弃/替换补值处理

import numpy as np
import pandas as pd
import matplotlib
import matplotlib.pyplot as plt
import seaborn as sns
import warnings
import missingno as msno
import time

Train_data=pd.read_csv('./used_car_train_20200313.csv', sep=' ')
TestB_data=pd.read_csv('./used_car_testB_20200421.csv', sep=' ')


#判断数据缺失和异常
Train_data.isnull().sum()
TestB_data.isnull().sum()

#缺失值可视化

missing=Train_data.isnull().sum()
missing=missing[missing>0]
missing.sort_values(inplace=True)
missing.plot.bar()

#处理异常值，除了notRepairedDamage为object类型，其他均为数字
Train_data['notRepairedDamage'].value_counts()
Train_data['notRepairedDamage'].replace('-',np.nan,inplace=True)
TestB_data['notRepairedDamage'].value_counts()
TestB_data['notRepairedDamage'].replace('-', np.nan, inplace=True)

#删掉两个类别特征严重倾斜的值
Train_data.drop(['seller','offerType'],axis=1,inplace=True)

#了解预测值的分布
Train_data['price']

查看数据大致形态如分布、峰度、偏值等，必要的话进行变换

#总体分布概况
import scipy.stats as st

y=Train_data['price']
plt.figure(1);plt.title('Johnson SU')
sns.distplot(y,kde=False, fit=st.johnsonsu)
plt.figure(2); plt.title('Normal')
sns.distplot(y, kde=False, fit=st.norm)
plt.figure(3); plt.title('Log Normal')
sns.distplot(y, kde=False, fit=st.lognorm)

#因为价格不服从正态分布，所以在进行回归之前，它必须进行转换。虽然对数变换做得很好，但最佳拟合是无界约翰逊分布

#查看偏度和峰度

sns.distplot(Train_data['price']);
print("Skewness: %f" % Train_data['price'].skew())
print("Kurtosis: %f" % Train_data['price'].kurt())

#查看预测值的具体频数
plt.hist(Train_data['price'], orientation='vertical',histtype='bar',color='red')
plt.show()

#log变换
plt.hist(np.log(Train_data['price']),orientation='vertical',histtype='bar',color='red')
plt.show()

选取要进行预测的数据类型和相关列

Y_train=Train_data['price']

#人为区分labels
numeric_features = ['power', 'kilometer', 'v_0', 'v_1', 'v_2', 'v_3', 'v_4', 'v_5', 'v_6', 'v_7', 'v_8', 'v_9', 'v_10', 'v_11', 'v_12', 'v_13','v_14' ]

categorical_features = ['name', 'model', 'brand', 'bodyType', 'fuelType', 'gearbox', 'notRepairedDamage', 'regionCode',]

#特征unique分布
for cat_fea in categorical_features:
    print(cat_fea+"的特征分布如下：")
    print("{}特征有个{}不同的值".format(cat_fea, Train_data[cat_fea].nunique()))
    print(Train_data[cat_fea].value_counts())

对数据进行相关性分析

#相关性分析
price_numeric=Train_data[numeric_features]
correlation=price_numeric.corr()
print(correlation['price'].sort_values(ascending=False),'\n')

price 1.000000
v_12 0.692823
v_8 0.685798
v_0 0.628397
power 0.219834
v_5 0.164317
v_2 0.085322
v_6 0.068970
v_1 0.060914
v_14 0.035911
v_13 -0.013993
v_7 -0.053024
v_4 -0.147085
v_9 -0.206205
v_10 -0.246175
v_11 -0.275320
kilometer -0.440519
v_3 -0.730946
Name: price, dtype: float64

figsize,ax=plt.subplots(figsize=(7,7))
plt.title('Correlation of Numeric Features with Price',y =1,size=16)
sns.heatmap(correlation, square=True, vmax=0.8)

#每个数字特征分布可视化
f=pd.melt(Train_data,value_vars=numeric_features)#df.melt() 则是将宽数据集变成长数据集
g=sns.FacetGrid(f, col="variable", col_wrap=4, sharex=False, sharey=False) #分面绘制
g=g.map(sns.distplot,"value")

#数字特征相互间关系可视化
sns.set()
columns = ['price', 'v_12', 'v_8' , 'v_0', 'power', 'v_5',  'v_2', 'v_6', 'v_1', 'v_14']
sns.pairplot(Train_data[columns],size = 2 ,kind ='scatter',diag_kind='kde')
plt.show()

#箱式图
categorical_features = ['model',
 'brand',
 'bodyType',
 'fuelType',
 'gearbox',
 'notRepairedDamage']
for c in categorical_features:
    Train_data[c] = Train_data[c].astype('category')
    if Train_data[c].isnull().any():
        Train_data[c] = Train_data[c].cat.add_categories(['MISSING'])
        Train_data[c] = Train_data[c].fillna('MISSING')

def boxplot(x, y, **kwargs):
    sns.boxplot(x=x, y=y)
    x=plt.xticks(rotation=90)

f = pd.melt(Train_data, id_vars=['price'], value_vars=categorical_features)
g = sns.FacetGrid(f, col="variable",  col_wrap=2, sharex=False, sharey=False, size=5)
g = g.map(boxplot, "value", "price")

特征分析

清晰异常值

def outliers_proc(data, col_name, scale=3):
    """
    用于清洗异常值，默认用 box_plot（scale=3）进行清洗
    :param data: 接收 pandas 数据格式
    :param col_name: pandas 列名
    :param scale: 尺度
    :return:
    """

    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)) #0.75分位数的值-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)#前面返回异常的pandas.Series 数据，后面返回临界值

    data_n = data.copy() #先复制一个df
    data_series = data_n[col_name] #某一列的值
    rule, value = box_plot_outliers(data_series, box_scale=scale)
    index = np.arange(data_series.shape[0])[rule[0] | rule[1]]  #shape[0]是行数，丨是or的意思，这个就是输出有异常值的索引数
    print("Delete number is: {}".format(len(index)))#输出异常值个数
    data_n = data_n.drop(index)#删除异常值
    data_n.reset_index(drop=True, inplace=True)#重新设置索引
    print("Now column number is: {}".format(data_n.shape[0])) #删除异常值之后数值的个数
    index_low = np.arange(data_series.shape[0])[rule[0]]#低于临界值的索引数
    outliers = data_series.iloc[index_low]#低于临界值的值
    print("Description of data less than the lower bound is:")
    print(pd.Series(outliers).describe())
    index_up = np.arange(data_series.shape[0])[rule[1]]
    outliers = data_series.iloc[index_up]
    print("Description of data larger than the upper bound is:")
    print(pd.Series(outliers).describe())
    
    fig, ax = plt.subplots(1, 2, figsize=(10, 7))
    sns.boxplot(y=data[col_name], data=data, palette="Set1", ax=ax[0])
    sns.boxplot(y=data_n[col_name], data=data_n, palette="Set1", ax=ax[1])
    return data_n

Train_data = outliers_proc(Train_data, 'power', scale=3)

特征构造

#用一列做标签区分一下训练集和测试集，然后放在一起方便构造特征
Train_data['train']=1
TestB_data['train']=0
data=pd.concat([Train_data,TestB_data],ignore_index=True)

#创建汽车使用时间，因为数据里有时间出错的格式，要errors='coerce'
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')

#从邮编中提取城市信息，相当于加入了先验知识
data['city']=data['regionCode'].apply(lambda x: str(x)[:-3])
data=data

#计算某品牌的销售统计量，以train数据计算数据量
Train_gb=Train_data.groupby('brand')
all_info={}
for kind, kind_data in Train_gb:
    info={}
    kind_data=kind_data[kind_data['price']>0]
    info['brand_amount'] = len(kind_data)
    info['brand_price_max'] = kind_data.price.max()
    info['brand_price_median'] = kind_data.price.median()
    info['brand_price_min'] = kind_data.price.min()
    info['brand_price_sum'] = kind_data.price.sum()
    info['brand_price_std'] = kind_data.price.std()
    info['brand_price_average'] = round(kind_data.price.sum() / (len(kind_data) + 1), 2)
    all_info[kind]=info
    
brand_fe=pd.DataFrame(all_info).T.reset_index().rename(columns={"index":"brand"})
data=data.merge(brand_fe,how='left',on="brand")
data.head()

数据分箱，给树模型使用

#数据分箱的好处：

 1. 离散后稀疏向量内积乘法运算速度更快，计算结果也方便存储，容易扩展；
 2. 离散后的特征对异常值更具鲁棒性，如 age>30 为 1 否则为 0，对于年龄为 200 的也不会对模型造成很大的干扰；
 3. LR 属于广义线性模型，表达能力有限，经过离散化后，每个变量有单独的权重，这相当于引入了非线性，能够提升模型的表达能力，加大拟合；
 4. 离散后特征可以进行特征交叉，提升表达能力，由 M+N 个变量编程 M*N 个变量，进一步引入非线形，提升了表达能力；
 5. 特征离散后模型更稳定，如用户年龄区间，不会因为用户年龄长了一岁就变化
bin=[i*10 for i in range(31)]
data['power_bin']=pd.cut(data['power'], bin, labels=False)
data[['power_bin','power']].head()

# 删除不需要的数据
data = data.drop(['creatDate', 'regDate', 'regionCode'], axis=1)

# 目前的数据其实已经可以给树模型使用了，树模型不需要归一化
data.to_csv('data_for_tree.csv', index=0)

归一化，给LR使用

# 取 log做归一化
from sklearn import preprocessing
min_max_scaler=preprocessing.MinMaxScaler()
data['power']=np.log(data['power']+1)
data['power']=((data['power']-np.min(data['power']))/np.max(data['power']-np.min(data['power'])))

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

#转换其他统计量特征
def max_min(x):
    return (x - np.min(x)) / (np.max(x) - np.min(x))

data['brand_amount'] = ((data['brand_amount'] - np.min(data['brand_amount'])) / 
                        (np.max(data['brand_amount']) - np.min(data['brand_amount'])))
data['brand_price_average'] = ((data['brand_price_average'] - np.min(data['brand_price_average'])) / 
                               (np.max(data['brand_price_average']) - np.min(data['brand_price_average'])))
data['brand_price_max'] = ((data['brand_price_max'] - np.min(data['brand_price_max'])) / 
                           (np.max(data['brand_price_max']) - np.min(data['brand_price_max'])))
data['brand_price_median'] = ((data['brand_price_median'] - np.min(data['brand_price_median'])) /
                              (np.max(data['brand_price_median']) - np.min(data['brand_price_median'])))
data['brand_price_min'] = ((data['brand_price_min'] - np.min(data['brand_price_min'])) / 
                           (np.max(data['brand_price_min']) - np.min(data['brand_price_min'])))
data['brand_price_std'] = ((data['brand_price_std'] - np.min(data['brand_price_std'])) / 
                           (np.max(data['brand_price_std']) - np.min(data['brand_price_std'])))
data['brand_price_sum'] = ((data['brand_price_sum'] - np.min(data['brand_price_sum'])) / 
                           (np.max(data['brand_price_sum']) - np.min(data['brand_price_sum'])))

# 对类别特征进行 OneEncoder
data = pd.get_dummies(data, columns=['model', 'brand', 'bodyType', 'fuelType',
                                     'gearbox', 'notRepairedDamage', 'power_bin'])

相关性分析

# 相关性分析
print(data['power'].corr(data['price'], method='spearman'))
print(data['kilometer'].corr(data['price'], method='spearman'))
print(data['brand_amount'].corr(data['price'], method='spearman'))
print(data['brand_price_average'].corr(data['price'], method='spearman'))
print(data['brand_price_max'].corr(data['price'], method='spearman'))
print(data['brand_price_median'].corr(data['price'], method='spearman'))

0.5728285196051496
-0.4082569701616764
0.058156610025581514
0.3834909576057687
0.259066833880992
0.38691042393409447

data_numeric = data[['power', 'kilometer', 'brand_amount', 'brand_price_average', 
                     'brand_price_max', 'brand_price_median']]
correlation = data_numeric.corr()

f , ax = plt.subplots(figsize = (7, 7))
plt.title('Correlation of Numeric Features with Price',y=1,size=16)
sns.heatmap(correlation,square = True,  vmax=0.8)