一、前言
在经过一段时间的数据分析学习后，以检验学习效果及实践应用为目的，将目光投向了天池竞赛，于是乎选择了新人赛中的幸福感作为第一站，开始了寻找幸福的道路。
这里不得不提的是，由于学资尚浅，结果并不理想，所以献丑了，也希望大家可以不吝赐教，谢谢。
所有代码及相关文件请见Git仓库。
二、题目概述
题目来自天池新人赛，名为“快来一起挖掘幸福感！”（链接地址）。关于赛题总结有以下几点：
- 数据来自于《中国综合社会调查（CGSS）》项目的调查问卷结果，多为选择题的形式；
- 通过多个变量来对幸福感进行预测，变量包含有收入、健康、职业、社交等方面；
- 赛题希望预测准确不是唯一目的，更希望挖掘出变量和幸福感之间的关系；
- 测评指标为MSE，
。
三、关于数据
赛题数据主要有三部分，分别为：
- 数据字典
happiness_index.xlsx及数据来源的调查问卷happiness_survey_cgss2015.pdf - 训练集
happiness_train_complete.csv - 测试集
happiness_test_complete.csv（用于进行预测）
下面主要了解训练集。训练集中有8000行140列，测试集有2968行139列，训练集包含一列happiness为标签列，数据字典部分截图如下。
>>> df_train = pd.read_csv('../raw/happiness_train_complete.csv', encoding='gb2312')
>>> df_test = pd.read_csv('../raw/happiness_test_complete.csv', encoding='gb2312')
>>> df_train.shape
(8000, 140)
>>> df_test.shape
(2968, 139)
数据字典
通过简单的探索发现数据中存在较多的缺失值及特殊值（即数据字典中标出的-1 = 不适用; -2 = 不知道; -3 = 拒绝回答; -8 = 无法回答;）。
>>> df_train.isnull().sum().sum() # 查看数据集中缺失值总数
80805
由于变量较多，无法在输出中完整显示，所以编写了一个过程，用来逐列对数据进行描述并输出到文件，用来发现变量的特性及共性，代码及注释如下：
# 逐列对对数据进行描述的过程
def desc_bycol(df, lst_abnormal=None, output=None):
'''逐列对数据进行描述
df：需要进行描述的DataFrame
lst_abnormal：特殊值列表，默认为[-1, -2, -3, -8, None]
output：输出的文件名，默认为`describe.txt`
'''
import os
import pandas as pd
if lst_abnormal == None:
lst_abnormal = [-1, -2, -3, -8, None] # 默认的特殊值
if output == None:
output = 'describe.txt' # 默认的输出文件名
s0 = '' # 将输出的全部内容存到`s0`
for col in df.columns:
# `s1`为各列内容的输出格式模板，包括：
# - 抽取5个样例
# - pandas的describe方法，包括计数、均值、标准差、最大最小值、各分位数
# - 异常值的个数及比例
# - 各唯一值的个数（若唯一值较多则不输出）
s1 = '''Column: {}
--- Take 5 Samples ---
{}\n
--- General Describe ---
{}\n
--- Abnormal ---
{}\n
--- Unique ---
{}\n
{}\n\n'''
ser = df[col]
# 抽取5个样例
smpl = ser.sample(5)
# 计数、均值、标准差、最大最小值、各分位数[.01, .05, .1, .25, .5, .75, .9, .99, .999]
# 用于发现数据的大致分布及异常值情况
desc = ser.describe(include='all',
percentiles=[.01, .05, .1, .25, .5, .75, .9, .99, .999])
# 定位特殊值，并进行计数及计算比例
bln = ser.isnull() | ser.apply(lambda x: True if x in lst_abnormal else False)
abn = ser[bln].value_counts(dropna=False)
df_abn = abn.to_frame()
df_abn['new'] = abn / ser.shape[0]
df_abn.columns = ['Count', 'Proportion']
# 计算各唯一值的个数，若超过20个唯一值则不输出
d_unique = ser.unique()
n_unique = d_unique.shape[0]
if n_unique > 20:
unq = '{} unique values.'.format(n_unique)
else:
unq = ser.value_counts(dropna=False)
# 将计算结果填到样例模板中并返回给s0
s1 = s1.format(col, smpl, desc, df_abn, unq, '='*30)
s0 += s1
# 输出到文件
if os.path.exists(output):
f_mode = 'w'
else:
f_mode = 'a'
with open(output, f_mode) as f:
f.write(s0)
运行输出的结果如下（以happiness列为例）：
Column: happiness
--- Take 5 Samples ---
568 3
6328 2
6251 4
6578 3
226 4
Name: happiness, dtype: int64
--- General Describe ---
count 8000.000000
mean 3.850125
std 0.938228
min -8.000000
1% 1.000000
5% 2.000000
10% 3.000000
25% 4.000000
50% 4.000000
75% 4.000000
90% 5.000000
99% 5.000000
99.9% 5.000000
max 5.000000
Name: happiness, dtype: float64
--- Abnormal ---
Count Proportion
-8 12 0.0015
--- Unique ---
4 4818
5 1410
3 1159
2 497
1 104
-8 12
Name: happiness, dtype: int64
四、数据清洗
通过将输出的各列信息与数据字典相结合，得出数据清洗方案，简述如下，详情请见Git仓库中的特征梳理.xlsx文件。
-
对于有相似特性的变量，进行统一处理：
- 分类型变量：
- 对于特殊值：因为考虑到特殊值（前述的
[-1, -2, -3, -8]）表达了受访者的一种态度，所以不对其进行处理； - 对于缺失值：按照不同的缺失比例采取以下方法进行处理：
(0, 5%]：填为众数
(5%, 30%]：用随机森林分类填补
(30%, 80%]：单独作为一类，设定为-10
(80%, 100%]：删除该变量
- 对于特殊值：因为考虑到特殊值（前述的
- 数值型变量：
- 收入类变量：经检查发现收入类变量均存在异常值情况，所以这里采取盖帽法盖到99%；
- 缺失值及特殊值：按照不同的存在比例采取以下方法进行处理：
(0, 5%]：先将特殊值改为None，然后填为中位数
(5%, 30%]：先将特殊值改为None，然后用随机森林回归插补
其他情况：在统一处理后再次检查，按具体情况分别处理
- 分类型变量：
-
对于变量的个性情况，按具体情况分别处理，这里不再展开，详情请见具体文件。
-
根据变量的特性，对变量进行衍生：
根据变量province（采访的31个省份）衍生出以下12个宏观指标：（数据来自统计年鉴）衍生变量 变量含义 衍生变量 变量含义 new_urban_p_rate 城镇人口比重 new_unemp_rate 失业率 new_p_incr_rate 人口增长率 new_cpi 居民消费价格指数 new_gdp 生产总值 new_dpi 人均可支配收入 new_gdp_incr 生产总值增长 new_consp 人均消费支出 new_consum 居民消费水平 new_public_budget_ratio 地区公共预算支出/收入 new_consum_incr 居民消费水平增长 new_p_density 人口密度 根据其他变量衍生出的变量有：
衍生变量 变量含义 衍生变量 变量含义 new_age 年龄 new_inc_rate_family 个人收入占家庭比重 new_edu_age 接受教育的时长（年） new_inc_avg_family 家庭人均收入 new_party_age 党龄 new_inc_ratio_spouse 个人收入与配偶收入的比值 new_mari_age 第一次结婚年龄 new_inc_gap 实际收入与期望收入的差距 new_remari 是否有再婚
将以上步骤分过程写成Python代码，存在clean_data_happiness.py中，以方便在Jupyter中根据不同需要重复调用处理数据。
下面对变量进行筛选，这里用了两个方法：
- 变异系数
考虑到不同量级的变量间方差是不具有可比性的，所以这里没有使用scikit-learn的方差过滤，所以自行计算变异系数进行对比来筛选变量，以此过滤掉包含信息极少的变量。从输出结果，选择删除掉4个变量，invest_6、new_cpi、s_birth、birth。
>>> cv = x_train.std(axis=0) / x_train.mean(axis=0)
>>> cv = cv.apply(np.abs)
>>> cv.sort_values(ascending=True)
new_cpi 0.003562
s_birth 0.008380
birth 0.008643
new_gdp_incr 0.015037
height_cm 0.049058
...
invest_5 23.056088
trust_12 30.338438
invest_7 32.621566
invest_8 39.959353
invest_6 NaN
Length: 148, dtype: float64
>>> x_train['invest_6'].value_counts()
0 6390
Name: invest_6, dtype: int64
- 互信息法
使用互信息法将与标签Y不具有相关性的变量删除。
>>> fs_mic = mutual_info_classif(x_train, y_train)
>>> np.where(fs_mic == 0)
(array([ 10, 12, 13, 14, 16, 17, 18, 19, 20, 27, 30, 33, 35,
38, 59, 63, 64, 65, 69, 74, 76, 77, 79, 81, 82, 86,
90, 95, 97, 105, 108, 113, 116, 143], dtype=int64),)
>>> x_train.columns[np.where(fs_mic == 0)] # 与标签Y互信息值为0的变量
Index(['political', 'property_0', 'property_1', 'property_2', 'property_4',
'property_5', 'property_6', 'property_7', 'property_8', 'hukou_loc',
'media_3', 'media_6', 'leisure_2', 'leisure_5', 'work_yr', 'insur_2',
'insur_3', 'insur_4', 'house', 'invest_3', 'invest_5', 'invest_6',
'invest_8', 'daughter', 'minor_child', 's_political', 's_work_status',
'f_work_14', 'm_edu', 'trust_2', 'trust_5', 'trust_10', 'trust_13',
'new_remari'],
dtype='object')
通过上面两步，共选择出37个变量进行删除，最后留下111个变量。
>>> x_train.drop(columns=cols_drop, inplace=True)
>>> x_test.drop(columns=cols_drop, inplace=True)
>>> s = 'Shape of\n x_train: {}\t x_test: {}\n y_train: {}\t y_test: {}'
>>> print(s.format(x_train.shape, x_test.shape, y_train.shape, y_test.shape))
Shape of
x_train: (6390, 111) x_test: (1598, 111)
y_train: (6390,) y_test: (1598,)
五、选择模型
模型选择分为两个部分，具有解释性的模型和不具有解释性的模型，由于学资有限，仅从已掌握的模型中进行选择，在此之前先探索变量之间的关系。
对于线性回归、逻辑回归这类统计模型需要满足较多的假设条件，所以先画图观察变量之间的关系：
# 查看Pearson相关系数
>>> corr = df_xy.corr(method='pearson')
>>> plt.figure(figsize=(20,20))
>>> sns.heatmap(corr.apply(np.abs), vmin=0, vmax=1, cmap='Blues', annot=False, xticklabels=True, yticklabels=True)
>>> plt.show()
Pearson相关系数
# 查看Spearman相关系数
>>> corr = df_xy.corr(method='spearman')
>>> plt.figure(figsize=(20,20))
>>> sns.heatmap(corr.apply(np.abs), vmin=0, vmax=1, cmap='Blues', annot=False, xticklabels=True, yticklabels=True)
>>> plt.show()
Spearman相关系数
从图中可见，各变量和Y之间的相关性较弱，且变量之间存在较强的共线性关系。
下面通过散点图观察变量和Y之间的分布关系（因变量较多，只选取了互信息值排序前5和后5的变量），从图中看出变量的分布与Y的关联也是较弱的。
互信息排序前5的变量
互信息排序后5的变量
（一）可解释性模型
由于赛题希望可以发掘出变量和幸福感之间的关系，同时由于测评指标为MSE，所以最好从具有解释性的回归类模型中进行选择，初步将决策树回归、线性回归、逻辑回归作为备选。
从上述的探索中看出，要想拟合出较好的统计模型需要对变量做更多的操作，所以这里暂时放弃线性回归，待日后再做细致的梳理。现在只使用决策树回归、逻辑回归建立模型查看效果。（逻辑回归同样不适合现在的数据，只为确认效果）
- 尝试决策树回归
在参数保持默认的情况下尝试决策树回归，然后通过画学习曲线来调整参数，确定参数大致范围，最后在确定的小范围内用网格搜索GridSearchCV细化调参。（借鉴了贪心算法的思想，通过局部最优来接近全局最优，也是为了降低运算时间而选择这种方式，后面的模型也采用了同样的方法）
通过网格搜索调参选择出较优的参数，然后计算在验证集上的误差。
# 网格搜索调参，best_score_为-0.5264
parms = {'max_depth':np.arange(2, 10), 'min_samples_split':np.arange(195, 205)}
mdl_DTR = DecisionTreeRegressor()
mdl_DTR_gs = GridSearchCV(mdl_DTR, parms, scoring='neg_mean_squared_error', cv=10, n_jobs=-1)
mdl_DTR_gs.fit(x_train, y_train)
mdl_DTR_gs.best_score_
# 计算在验证集上的误差，返回结果为0.5741
mdl_DTR = DecisionTreeRegressor(max_depth=4, min_samples_split=199)
mdl_DTR.fit(x_train, y_train)
mdl_DTR.score(x_test, y_test)
y_pred = mdl_DTR.predict(x_test)
mse = mean_squared_error(y_test, y_pred)
效果很差，然后尝试对数据进行标准化，依据经验应该会提升模型效果及运算速度。
# 对数据进行标准化
trans_std = StandardScaler()
trans_std.fit(x_train)
x_train_std = trans_std.transform(x_train)
x_test_std = trans_std.transform(x_test)
经过尝试，网格搜索的best_score_ 为-0.5266， 验证集上的MSE为0.5741。
- 尝试逻辑回归
重复上述调参步骤，仅调整正则化系数C，使用未经标准化的数据，MSE为1.0481，标准化后为1.0212，果然现在的数据不适合这个模型，需要对变量进行更细致的处理才可以。
（二）不具解释性模型
经尝试以上两个具有解释性的模型效果较差，同时看到天池上关于赛题的讨论中，预测效果较好的都没有采用可解释性模型，所以为提升模型效果，转而尝试其他模型，如随机森林回归、支持向量回归。
- 随机森林回归
经尝试，使用标准化后的数据训练参数默认的模型，效果会有明显提升（0.5746 -> 0.5661），所以之后全部使用经标准化的数据。
重复上述调参步骤进行调参，最终得到的参数为：
-
n_estimators：300 -
max_depth：9 -
min_samples_split：94
网格搜索的best_score_ 为-0.4835， 验证集上的MSE为0.5214，有了较明显的提升。
- 支持向量回归
方法同上，经调参后的参数为：
-
C：0.88（模型倾向于稍大一点的边界） -
gamma：0.00495（可缓解过拟合情况，单个样本对整个分类超平面的影响较小）
（可见两个参数都倾向于使模型更简单）
网格搜索的best_score_ 为-0.4987， 验证集上的MSE为0.5225。
（三）模型比较及选择
尝试了以上4个模型，效果如下表所示，劣中取优，选择随机森林回归作为这次赛题的模型。
| 模型 | 网格搜索MSE | 验证集MSE |
|---|---|---|
| 决策树回归 | 0.5266 | 0.5741 |
| 逻辑回归 | 1.0481 | 1.0212 |
| 随机森林回归 | 0.4835 | 0.5214 |
| 支持向量回归 | 0.4987 | 0.5225 |
六、进行预测
在选定模型及确定参数后，便可对测试集进行预测并上传结果。
由于之前已使用scikit-learn将训练集切分为训练集及验证集，对切分后的训练集进行数据清洗，所以这里要重新读入训练数据，对训练数据整体及测试集重新执行上面的数据清洗过程，并训练模型，进行预测。
提交后，效果如预期一般：
测评结果
七、结语
此次参赛练习效果不好，总结有以下几方面原因：
-
关于特征工程
- 数据中有一些变量存在异常值的问题，通过分箱的方法可以较好的解决
- 若要选用线性回归、逻辑回归等模型，还需要进行哑变量处理，考虑到处理后会造成数据维度爆炸式增长，所以没有处理，不过可以考虑在降维后再行处理。
- 建模期间有尝试PCA降维，但验证结果发现，还是没有经过降维的数据效果要明显好于降维后的数据，所以还需要深挖变量之间的关系，组合加工，进行特征创造，从而达到降维的目的。
-
关于模型选择
经查看比赛论坛，大神们选择的都是一些比较复杂的模型，大部分在数据处理上投入的精力较少，所以我认为对于这个赛题而言，模型的选择还是很重要的，这是我今后需要努力的一个主要方向。
看到此文的读者若发现其中存在错误、不妥之处，或有更好的建议，烦请不吝指教，谢谢！
最后，欢迎批评，欢迎指导，欢迎交流！ヾ(◍°∇°◍)ﾉﾞ
网友评论