「特征工程」之零基础入门数据挖掘

from scipy.stats import boxcox
boxcox_transformed_data = boxcox(original_data)

# power属性是有异常点的
num_data.boxplot(['power'])

"""这里包装了一个异常值处理的代码，可以随便调用"""
def outliers_proc(data, col_name, scale=3):
    """
        用于截尾异常值， 默认用box_plot(scale=3)进行清洗
        param:
            data：接收pandas数据格式
            col_name: pandas列名
            scale: 尺度
    """
    data_col = data[col_name]
    Q1 = data_col.quantile(0.25) # 0.25分位数
    Q3 = data_col.quantile(0.75)  # 0,75分位数
    IQR = Q3 - Q1

    data_col[data_col < Q1 - (scale * IQR)] = Q1 - (scale * IQR)
    data_col[data_col > Q3 + (scale * IQR)] = Q3 + (scale * IQR)

    return data[col_name]

num_data['power'] = outliers_proc(num_data, 'power')

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))
        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)

    data_n = data.copy()
    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]]
    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

# 删除重复值
data.drop_duplicates()
# dropna()可以直接删除缺失样本，但是有点不太好

# 填充固定值
train_data.fillna(0, inplace=True) # 填充 0
data.fillna({0:1000, 1:100, 2:0, 4:5})   # 可以使用字典的形式为不用列设定不同的填充值

train_data.fillna(train_data.mean(),inplace=True) # 填充均值
train_data.fillna(train_data.median(),inplace=True) # 填充中位数
train_data.fillna(train_data.mode(),inplace=True) # 填充众数

train_data.fillna(method='pad', inplace=True) # 填充前一条数据的值，但是前一条也不一定有值
train_data.fillna(method='bfill', inplace=True) # 填充后一条数据的值，但是后一条也不一定有值

"""插值法：用插值法拟合出缺失的数据，然后进行填充。"""
for f in features: 
    train_data[f] = train_data[f].interpolate()

train_data.dropna(inplace=True)

"""填充KNN数据：先利用knn计算临近的k个数据，然后填充他们的均值"""
from fancyimpute import KNN
train_data_x = pd.DataFrame(KNN(k=6).fit_transform(train_data_x), columns=features)

# 还可以填充模型预测的值， 这一个在我正在写的数据竞赛修炼笔记的第三篇里面可以看到，并且超级精彩，还在写

"""下面以power为例进行分桶， 当然构造一列新特征了"""
bin = [i*10 for i in range(31)]
num_data['power_bin'] = pd.cut(num_data['power'], bin, labels=False)

# 连续值经常离散化或者分离成“箱子”进行分析。
# 假设某项研究中一组人群的数据，想将他们进行分组，放入离散的年龄框中
ages = [20, 22, 25, 27, 21, 23, 37, 31, 61, 45, 41, 32]
# 如果按年龄分成18-25, 26-35, 36-60, 61以上的若干组，可以使用pandas中的cut
bins = [18, 25, 35, 60, 100]         # 定义箱子的边
cats = pd.cut(ages, bins)
print(cats)   # 这是个categories对象    通过bin分成了四个区间， 然后返回每个年龄属于哪个区间
# codes属性
print(cats.codes)    #  这里返回一个数组，指明每一个年龄属于哪个区间
print(cats.categories)
print(pd.value_counts(cats))   # 返回结果是每个区间年龄的个数

# 与区间的数学符号一致， 小括号表示开放，中括号表示封闭， 可以通过right参数改变
print(pd.cut(ages, bins, right=False))

# 可以通过labels自定义箱名或者区间名
group_names = ['Youth', 'YonngAdult', 'MiddleAged', 'Senior']
data = pd.cut(ages, bins, labels=group_names)
print(data)
print(pd.value_counts(data))

# 如果将箱子的边替代为箱子的个数，pandas将根据数据中的最小值和最大值计算出等长的箱子
data2 = np.random.rand(20)
print(pd.cut(data2, 4, precision=2))   # precision=2 将十进制精度限制在2位

# qcut是另一个分箱相关的函数， 基于样本分位数进行分箱。取决于数据的分布，使用cut不会使每个箱子具有相同数据数量的数据点，而qcut，使用
# 样本的分位数，可以获得等长的箱
data3 = np.random.randn(1000)   # 正太分布
cats = pd.qcut(data3, 4)
print(pd.value_counts(cats))

minmax = MinMaxScaler()
num_data_minmax = minmax.fit_transform(num_data)
num_data_minmax = pd.DataFrame(num_data_minmax, columns=num_data.columns, index=num_data.index)

"""类别特征某些需要独热编码一下"""
hot_features = ['bodyType', 'fuelType', 'gearbox', 'notRepairedDamage']
cat_data_hot = pd.get_dummies(cat_data, columns=hot_features)

from scipy.cluster.hierarchy import linkage, dendrogram
#from sklearn.cluster import AgglomerativeClustering
from sklearn.cluster import KMeans

ac = KMeans(n_clusters=3)
ac.fit(model_price_data)

model_fea = ac.predict(model_price_data)
plt.scatter(model_price_data[:,0], model_price_data[:,1], c=model_fea)

cat_data_hot['model_fea'] = model_fea
cat_data_hot = pd.get_dummies(cat_data_hot, columns=['model_fea'])

# 这里是为了标记一下哪些字段有错误
def regDateFalse(x):
    if str(x)[4:6] == '00':
        return 1
    else: 
        return 0
time_data['regDateFalse'] = time_data['regDate'].apply(lambda x: regDateFalse(x))
# 这里是改正错误字段
def changeFalse(x):
    x = str(x)
    if x[4:6] == '00':
        x = x[0:4] + '01' + x[6:]
        x = int(x)
    return x
time_data['regDate'] = time_data['regDate'].apply(lambda x: changeFalse(x))
# 使用时间：data['creatDate'] - data['regDate']，反应汽车使用时间，一般来说价格与使用时间成反比
# 不过要注意，数据里有时间出错的格式，所以我们需要 errors='coerce'
time_data['used_time'] = (pd.to_datetime(time_data['creatDate'], format='%Y%m%d') - 
                            pd.to_datetime(time_data['regDate'], format='%Y%m%d')).dt.days
# 修改错误
# 但是需要加上那一个月
time_data.loc[time_data.regDateFalse==1, 'used_time'] += 30
# 删除标记列
del time_data['regDateFalse']

# 使用时间换成年来表示
time_data['used_time'] = time_data['used_time'] / 365.0
time_data['Is_scrap'] = time_data['used_time'].apply(lambda x: 1 if x>=10 else 0)

bins = [0, 3, 7, 10, 20, 30]
time_data['estivalue'] = pd.cut(time_data['used_time'], bins, labels=False)

# 选出淡旺季
low_seasons = ['3', '6', '7', '8']
time_data['is_low_seasons'] = time_data['creatDate'].apply(lambda x: 1 if str(x)[5] in low_seasons else 0)

# 独热一下
time_data = pd.get_dummies(time_data, columns=['is_low_seasons'])

# 这样时间特征构造完毕，删除日期了
del time_data['regDate']
del time_data['creatDate']

# 构造统计特征的话需要在训练集上先计算
train_data_timestats = train_data.copy()   # 不要动train_data

train_data_timestats['estivalue'] = time_data['estivalue'][:train_data.shape[0]]
train_data_timestats['price'] = train_target

train_gt = train_data_timestats.groupby('estivalue')
all_info = {}
for kind, kind_data in train_gt:
    info = {}
    kind_data = kind_data[kind_data['price'] > 0]
    info['estivalue_count'] = len(kind_data)
    info['estivalue_price_max'] = kind_data.price.max()
    info['estivalue_price_median'] = kind_data.price.median()
    info['estivalue_price_min'] = kind_data.price.min()
    info['estivalue_price_sum'] = kind_data.price.sum()
    info['estivalueprice_std'] = kind_data.price.std()
    info['estivalue_price_average'] = round(kind_data.price.sum() / (len(kind_data) + 1), 2)
    all_info[kind] = info

estivalue_fe = pd.DataFrame(all_info).T.reset_index().rename(columns={"index": "estivalue"})
time_data = time_data.merge(estivalue_fe, how='left', on='estivalue')

cat_data['city'] = cat_data['regionCode'].apply(lambda x: str(x)[0])

com_car = [2.0, 3.0, 6.0]  # 商用车
GL_car = [0.0, 4.0, 5.0]   # 豪华系列  
self_car = [1.0, 7.0]

def class_bodyType(x):
    if x in GL_car:
        return 0
    elif x in com_car:
        return 1
    else:
        return 2

cat_data['car_class'] = cat_data['bodyType'].apply(lambda x : class_bodyType(x))

# 是否是新能源
is_fuel = [0.0, 1.0, 2.0, 3.0]
cat_data['is_fuel'] = cat_data['fuelType'].apply(lambda x: 1 if x in is_fuel else 0)

train_data_gearbox = train_data.copy()   # 不要动train_data
train_data_gearbox['gearbox'] = cat_data['gearbox'][:train_data.shape[0]]
train_data_gearbox['price'] = train_target

train_gb = train_data_gearbox.groupby('gearbox')
all_info = {}
for kind, kind_data in train_gb:
    info = {}
    kind_data = kind_data[kind_data['price'] > 0]
    info['gearbox_count'] = len(kind_data)
    info['gearbox_price_max'] = kind_data.price.max()
    info['gearbox_price_median'] = kind_data.price.median()
    info['gearbox_price_min'] = kind_data.price.min()
    info['gearbox_price_sum'] = kind_data.price.sum()
    info['gearbox_std'] = kind_data.price.std()
    info['gearbox_price_average'] = round(kind_data.price.sum() / (len(kind_data) + 1), 2)
    all_info[kind] = info

gearbox_fe = pd.DataFrame(all_info).T.reset_index().rename(columns={"index": "gearbox"})

cat_data = cat_data.merge(gearbox_fe, how='left', on='gearbox')

# 删掉bodyType和fuelType，然后把gearbox， car_class is_fuel独热一下， 这个不能太早，构造晚了统计特征之后再独热
del cat_data['bodyType']
del cat_data['fuelType']

cat_data = pd.get_dummies(cat_data, columns=['gearbox', 'car_class', 'is_fuel', 'notRepairedDamage'])

# 分成三段
bins = [0, 5, 10, 15]
num_data['kil_bin'] = pd.cut(num_data['kilometer'], bins, labels=False)

v_features = ['v_' + str(i) for i in range(15)]
num_data['v_sum'] = num_data[v_features].apply(lambda x: x.sum(), axis=1)
num_data['v_mean'] = num_data[v_features].apply(lambda x: x.mean(), axis=1)
num_data['v_std'] = num_data[v_features].apply(lambda x: x.std(), axis=1)

final_data = pd.concat([num_data, cat_data, time_data], axis=1)
final_data.shape  # (200000, 74)
final_data.head()

# 导入最后的数据
final_data = pd.read_csv('./pre_data/pre_data_beifen.csv')

# 对方差的大小排序
select_data.std().sort_values()    # select_data是final_data去掉了独热的那些特征

corr = select_data.corr('pearson')    # .corr('spearman')
plt.figure(figsize=(25, 15))
corr['price'].sort_values(ascending=False)[1:].plot(kind='bar')
plt.tight_layout()

# 下面看一下互相之间的关系
f, ax = plt.subplots(figsize=(50, 30))
sns.heatmap(corr, annot=True)

from sklearn.model_selection import cross_val_score, ShuffleSplit
from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import KFold

X = select_data.iloc[:, :-1]
Y = select_data['price']
names = select_data.columns

rf = RandomForestRegressor(n_estimators=20, max_depth=4)
kfold = KFold(n_splits=5, shuffle=True, random_state=7)
scores = []
for column in X.columns:
    print(column)
    tempx = X[column].values.reshape(-1, 1)
    score = cross_val_score(rf, tempx, Y, scoring="r2",
                              cv=kfold)
    scores.append((round(np.mean(score), 3), column))
print(sorted(scores, reverse=True))

# 下面再用xgboost跑一下
from xgboost import XGBRegressor
from xgboost import plot_importance

xgb = XGBRegressor()
xgb.fit(X, Y)

plt.figure(figsize=(20, 10))
plot_importance(xgb)
plt.show()

from mlxtend.feature_selection import SequentialFeatureSelector as SFS
from sklearn.linear_model import LinearRegression
from sklearn.ensemble import RandomForestRegressor
# sfs = SFS(LinearRegression(), k_features=20, forward=True, floating=False, scoring='r2', cv=0)
sfs = SFS(RandomForestRegressor(n_estimators=10, max_depth=4), k_features=20, forward=True, floating=False, scoring='r2', cv=0)

X = select_data.iloc[:, :-1]
Y = select_data['price']

sfs.fit(X, Y)
sfs.k_feature_names_    # 随机森林放这里跑太慢了，所以中断了

# 画出边际效应
from mlxtend.plotting import plot_sequential_feature_selection as plot_sfs

fig1 = plot_sfs(sfs.get_metric_dict(), kind='std_dev')
plt.grid()
plt.show()

from sklearn.linear_model import LinearRegression, Ridge,Lasso

models = [LinearRegression(), Ridge(), Lasso()]
result = dict()

for model in models:
    model_name = str(model).split('(')[0]
    scores = cross_val_score(model, X=train_X, y=train_y, verbose=0, cv=5, scoring='r2')
    result[model_name] = scores

from sklearn.decomposition import PCA
# 然后使用
pca = PCA(n_components=10)
X_new = pca.fit_transform(X)

"""查看PCA的一些属性"""
print(X_new.shape)   # （200000， 10）
print(pca.explained_variance_)    # 属性可以查看降维后的每个特征向量上所带的信息量大小（可解释性方差的大小）
print(pca.explained_variance_ratio_)  # 查看降维后的每个新特征的信息量占原始数据总信息量的百分比
print(pca.explained_variance_ratio_.sum())    # 降维后信息保留量