Kaggle Titanic 生存预测--特征工程

注册著名的数据科学竞赛平台 Kaggle 四个多月了,当初是为了要一份数据集,而一直没有参加比赛。这个寒假终于按耐不住,决定拿入门级的 Titanic: Machine Learning from Disaster 打响革命第一枪。

Titanic 生存预测比赛是一个二分类问题。题目提供了一份乘客名单,包含了乘客的名字、性别、年龄、船票等级等信息,以及是否成功获救的标记,最终需要提交一份对测试集中的乘客是否成功获救的 csv 文件。

经过了四次提交,最后我的 Public Score 暂时定格在 0.80861,这个成绩目前在前 7%。这篇博文主要简述一下我所做的尝试和改进,并对最后一次换用 XGBoost 所得到的最好成绩的代码进行一个详细的说明,也是对相似题目处理流程的一个总结。

Titanic-public-score

过程简述

第一次提交

在探索数据后,我决定选用以下特征进行预测'Pclass', 'Sex', 'Age', 'Embarked', 'SibSp', 'Parch', 'Fare'。其中 Age、Embarked 和 Fare 有缺失,考虑使用出现频率最高的值来填充 Embarked 特征(类别型)的缺失值,使用平均值来填充 Age 和 Fare 特征(数值型)的缺失值。而类别特征不能直接作为输入,因此采用 DictVectorizer 对特征抽取和特征向量化。

最后使用 RandomForest 分类器来进行预测。Public Score 为 0.73205。

第二次提交

Fare 只有一个缺失,而 Age 存在 86 个缺失值。因此,直接使用平均值来填充 Age 的缺失值可能对预测结果影响较大。我将 Age 从选取的有效特征中剔除,其他不变,Public Score 提高到 0.75598。

第三次提交

在对别人分享的 kernel 进行学习后,这次我做了比较详细的特征工程(具体操作在下一节),并且进行了 sklearn 中常用分类器效果的比较,最终选用了 SVC 分类器。本次的 Public Score 提高到 0.79904。

代码详述

这里针对最后一次提交对应的代码进行一个说明总结。

1
2
3
4
5
6
7
%matplotlib inline
import numpy as np
import pandas as pd
import re as re

import warnings
warnings.filterwarnings('ignore')
1
2
3
train = pd.read_csv('../input/train.csv', header=0, dtype={'Age': np.float64})
test = pd.read_csv('../input/test.csv', header=0, dtype={'Age': np.float64})
full_data = [train, test]

我们可以通过info方法来大致地了解训练集和测试集:

1
2
print(train.info())
print(test.info())

结果如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 891 entries, 0 to 890
Data columns (total 12 columns):
PassengerId 891 non-null int64
Survived 891 non-null int64
Pclass 891 non-null int64
Name 891 non-null object
Sex 891 non-null object
Age 714 non-null float64
SibSp 891 non-null int64
Parch 891 non-null int64
Ticket 891 non-null object
Fare 891 non-null float64
Cabin 204 non-null object
Embarked 889 non-null object
dtypes: float64(2), int64(5), object(5)
memory usage: 83.6+ KB
None
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 418 entries, 0 to 417
Data columns (total 11 columns):
PassengerId 418 non-null int64
Pclass 418 non-null int64
Name 418 non-null object
Sex 418 non-null object
Age 332 non-null float64
SibSp 418 non-null int64
Parch 418 non-null int64
Ticket 418 non-null object
Fare 417 non-null float64
Cabin 91 non-null object
Embarked 418 non-null object
dtypes: float64(2), int64(4), object(5)
memory usage: 36.0+ KB
None

特征工程

Pclass

Pclass 特征没有缺失值,因此可以通过groupby函数来计算船舱每一档的生还率:

1
print(train[['Pclass', 'Survived']].groupby(['Pclass'], as_index=False).mean())

可以看到结果显示船舱档位和生还率还是有较大联系的:

1
2
3
4
   Pclass  Survived
0 1 0.629630
1 2 0.472826
2 3 0.242363

Sex

1
print(train[['Sex', 'Survived']].groupby(['Sex'], as_index=False).mean())

可以看到女性的生还率更高:

1
2
3
      Sex  Survived
0 female 0.742038
1 male 0.188908

SibSp and Parch

这两个特征是船上表亲和直亲数量。通过这两个特征可以创造一个新的特征 - Family Size:

1
2
3
for dataset in full_data:
dataset['FamilySize'] = dataset['SibSp'] + dataset['Parch'] + 1
print(train[['FamilySize', 'Survived']].groupby(['FamilySize'], as_index=False).mean())

可以看到生还率和家庭成员数不是单纯的线性关系:

1
2
3
4
5
6
7
8
9
10
   FamilySize  Survived
0 1 0.303538
1 2 0.552795
2 3 0.578431
3 4 0.724138
4 5 0.200000
5 6 0.136364
6 7 0.333333
7 8 0.000000
8 11 0.000000

我们可以再创建一个新特征 IsAlone,用来表示是否是单独出行:

1
2
3
4
for dataset in full_data:
dataset['IsAlone'] = 0
dataset.loc[dataset['FamilySize'] == 1, 'IsAlone'] = 1
print(train[['IsAlone', 'Survived']].groupby(['IsAlone'], as_index=False).mean())

结果如下:

1
2
3
   IsAlone  Survived
0 0 0.505650
1 1 0.303538

现在特征值对结果的影响就比较明显了。

Embarked

注意 Embarked 特征有极少的缺失值,对于类别特征,可以考虑用最频繁的特征值进行填充:

1
train['Embarked'].value_counts()

可以看到,最频繁的特征值为S

1
2
3
4
S    644
C 168
Q 77
Name: Embarked, dtype: int64

因此,我们用S对缺失值进行填充:

1
2
3
for dataset in full_data:
dataset['Embarked'] = dataset['Embarked'].fillna('S')
print(train[['Embarked', 'Survived']].groupby(['Embarked'], as_index=False).mean())

结果如下。可以看到明显的区别:

1
2
3
4
  Embarked  Survived
0 C 0.553571
1 Q 0.389610
2 S 0.339009

Fare

测试集中 Fare 特征有一个缺失值。对于数值特征,可以用中位数(或者平均值)填充缺失值:

1
2
for dataset in full_data:
dataset['Fare'] = dataset['Fare'].fillna(train['Fare'].median())

因为每个 Fare 的值对应的样本数量太少,因此我们考虑划分区间。这里,我根据每个区间的样本数量将样本划分为四个区间,形成新特征 CategoricalFare:

1
2
train['CategoricalFare'] = pd.qcut(train['Fare'], 4)
print(train[['CategoricalFare', 'Survived']].groupby(['CategoricalFare'], as_index=False).mean())

结果如下:

1
2
3
4
5
   CategoricalFare  Survived
0 (-0.001, 7.91] 0.197309
1 (7.91, 14.454] 0.303571
2 (14.454, 31.0] 0.454955
3 (31.0, 512.329] 0.581081

Age

之前提到过,Age 特征的缺失值太多,不能简单的用平均值或者中位数进行填充。这里,我们不再简单的舍弃 Age 特征,而是换用不同的填充思路 - 根据已有数据的平均值和标准差随机生成填充数

1
2
3
4
5
6
7
8
for dataset in full_data:
age_avg = dataset['Age'].mean()
age_std = dataset['Age'].std()
age_null_count = dataset['Age'].isnull().sum()

age_null_random_list = np.random.randint(age_avg - age_std, age_avg + age_std, size=age_null_count)
dataset['Age'][np.isnan(dataset['Age'])] = age_null_random_list
dataset['Age'] = dataset['Age'].astype(int)

同样,将数据划分为区间。这里我按等区间跨度划分,生成新特征 CategoricalAge:

1
2
3
train['CategoricalAge'] = pd.cut(train['Age'], 5)

print(train[['CategoricalAge', 'Survived']].groupby(['CategoricalAge'], as_index=False).mean())

结果如下:

1
2
3
4
5
6
  CategoricalAge  Survived
0 (-0.08, 16.0] 0.504274
1 (16.0, 32.0] 0.345372
2 (32.0, 48.0] 0.394422
3 (48.0, 64.0] 0.434783
4 (64.0, 80.0] 0.090909

也有别的填充思路,例如用Sex, Title, Pclass三个特征构建随机森林模型,来生成填充值。

Name

姓名是我们一开始忽略掉的特征。实际上,通过人名前的头衔也可以进行分析:

1
2
3
4
5
6
7
8
9
10
def get_title(name):
title_search = re.search(' ([A-Za-z]+)\.', name)
if title_search:
return title_search.group(1)
return ''

for dataset in full_data:
dataset['Title'] = dataset['Name'].apply(get_title)

print(pd.crosstab(train['Title'], train['Sex']))

结果如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
Sex       female  male
Title
Capt 0 1
Col 0 2
Countess 1 0
Don 0 1
Dr 1 6
Jonkheer 0 1
Lady 1 0
Major 0 2
Master 0 40
Miss 182 0
Mlle 2 0
Mme 1 0
Mr 0 517
Mrs 125 0
Ms 1 0
Rev 0 6
Sir 0 1

由于存在一些次数较少的头衔,我们将头衔分类,把部分含义相近的头衔归在一起:

1
2
3
4
5
6
7
for dataset in full_data:
dataset['Title'] = dataset['Title'].replace(['Lady', 'Countess', 'Capt', 'Col', 'Sir', 'Don', 'Dr', 'Major', 'Rev', 'Jonkheer', 'Dona'], 'Rare')
dataset['Title'] = dataset['Title'].replace('Mlle', 'Miss')
dataset['Title'] = dataset['Title'].replace('Ms', 'Miss')
dataset['Title'] = dataset['Title'].replace('Mme', 'Mrs')

print(train[['Title', 'Survived']].groupby(['Title'], as_index=False).mean())

结果如下:

1
2
3
4
5
6
    Title  Survived
0 Master 0.575000
1 Miss 0.702703
2 Mr 0.156673
3 Mrs 0.793651
4 Rare 0.347826

数据清理

sklearn 要求数据都是数值型的,因此要进行数据的清理和转换。也可以用各种包内现成的算法,例如pd.get_dummies()。这里我们就简单的自己动手来将数据映射为数值:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
# 将数据映射为数值
for dataset in full_data:
# Mapping Sex
dataset['Sex'] = dataset['Sex'].map( {'female': 0, 'male': 1} ).astype(int)

# Mapping titles
title_mapping = {"Mr": 1, "Miss": 2, "Mrs": 3, "Master": 4, "Rare": 5}
dataset['Title'] = dataset['Title'].map(title_mapping)
dataset['Title'] = dataset['Title'].fillna(0)

# Mapping Embarked
dataset['Embarked'] = dataset['Embarked'].map( {'S': 0, 'C': 1, 'Q': 2} ).astype(int)

# Mapping Fare
dataset.loc[ dataset['Fare'] <= 7.91, 'Fare'] = 0
dataset.loc[(dataset['Fare'] > 7.91) & (dataset['Fare'] <= 14.454), 'Fare'] = 1
dataset.loc[(dataset['Fare'] > 14.454) & (dataset['Fare'] <= 31), 'Fare'] = 2
dataset.loc[ dataset['Fare'] > 31, 'Fare'] = 3
dataset['Fare'] = dataset['Fare'].astype(int)

# Mapping Age
dataset.loc[ dataset['Age'] <= 16, 'Age'] = 0
dataset.loc[(dataset['Age'] > 16) & (dataset['Age'] <= 32), 'Age'] = 1
dataset.loc[(dataset['Age'] > 32) & (dataset['Age'] <= 48), 'Age'] = 2
dataset.loc[(dataset['Age'] > 48) & (dataset['Age'] <= 64), 'Age'] = 3
dataset.loc[ dataset['Age'] > 64, 'Age'] = 4

特征选取

我们将一些无用的特征去掉,只保留部分原始特征和我们生成的新特征:

1
2
3
4
5
6
7
8
9
10
11
# Feature Selection
drop_elements = ['PassengerId', 'Name', 'Ticket', 'Cabin', 'SibSp', 'Parch', 'FamilySize']
train = train.drop(drop_elements, axis = 1)
train = train.drop(['CategoricalAge', 'CategoricalFare'], axis = 1)

test = test.drop(drop_elements, axis = 1)

print (train.head(10))

train = train.values
test = test.values

打印的结果如下:

1
2
3
4
5
6
7
8
9
10
11
   Survived  Pclass  Sex  Age  Fare  Embarked  IsAlone  Title
0 0 3 1 1 0 0 0 1
1 1 1 0 2 3 1 0 3
2 1 3 0 1 1 0 1 2
3 1 1 0 2 3 0 0 3
4 0 3 1 2 1 0 1 1
5 0 3 1 1 1 2 1 1
6 0 1 1 3 3 0 1 1
7 0 3 1 0 2 0 0 4
8 1 3 0 1 1 0 0 3
9 1 2 0 0 2 1 0 3

预测

首先生成训练集的特征和标记:

1
2
X = train[0::, 1::]
y = train[0::, 0]

引入并初始化 xgboost:

1
2
from xgboost import XGBClassifier
xgbc = XGBClassifier()

我们先使用 5 折交叉验证来看一看 xgboost 算法的预测效果如何:

1
2
3
4
from sklearn.cross_validation import cross_val_score
# sklearn.cross_validation 好像已经废除,使用 from sklearn.model_selection import cross_val_score

cross_val_score(xgbc, X, y, cv=5).mean() # 使用 5 折交叉验证

结果为 0.80595516611715701,还不错。之后就可以正式的训练和预测了:

1
2
3
4
5
6
7
8
9
xgbc.fit(X, y)
xgbc_y_predict = xgbc.predict(test)

test_data = pd.read_csv('../input/test.csv', header=0, dtype={'Age': np.float64})
xgbc_submission = pd.DataFrame({
'PassengerId': test_data['PassengerId'],
'Survived': xgbc_y_predict
})
xgbc_submission.to_csv('xgbc_submission.csv', index=False)

分类器比较

尽管这一次我直接选用了 XGBoost,但是我还是想把第三次提交时所用的分类器比较代码贴一下。首先是引入所有参与比较的分类器(全部是 sklearn 包中的),以及用于可视化的 matplotlib 和 seaborn。

1
2
3
4
5
6
7
8
9
10
11
12
import matplotlib.pyplot as plt
import seaborn as sns

from sklearn.model_selection import StratifiedShuffleSplit
from sklearn.metrics import accuracy_score, log_loss
from sklearn.neighbors import KNeighborsClassifier
from sklearn.svm import SVC
from sklearn.tree import DecisionTreeClassifier
from sklearn.ensemble import RandomForestClassifier, AdaBoostClassifier, GradientBoostingClassifier
from sklearn.naive_bayes import GaussianNB
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis, QuadraticDiscriminantAnalysis
from sklearn.linear_model import LogisticRegression

初始化各分类器,并使用数据集划分函数 StratifiedShuffleSplit:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
classifiers = [
KNeighborsClassifier(3),
SVC(probability=True),
DecisionTreeClassifier(),
RandomForestClassifier(),
AdaBoostClassifier(),
GradientBoostingClassifier(),
GaussianNB(),
LinearDiscriminantAnalysis(),
QuadraticDiscriminantAnalysis(),
LogisticRegression()
]

log_cols = ['Classifier', 'Accuracy']
log = pd.DataFrame(columns=log_cols)

sss = StratifiedShuffleSplit(n_splits=10, test_size=0.1, random_state=0)

X = train[0::, 1::]
y = train[0::, 0]

acc_dict = {}

分别计算各分类器的准确率:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
for train_index, test_index in sss.split(X, y):
X_train, X_test = X[train_index], X[test_index]
y_train, y_test = y[train_index], y[test_index]

for clf in classifiers:
name = clf.__class__.__name__
clf.fit(X_train, y_train)
train_predictions = clf.predict(X_test)
acc = accuracy_score(y_test, train_predictions)
if name in acc_dict:
acc_dict[name] += acc
else:
acc_dict[name] = acc

for clf in acc_dict:
acc_dict[clf] = acc_dict[clf] / 10.0
log_entry = pd.DataFrame([[clf, acc_dict[clf]]], columns=log_cols)
log = log.append(log_entry)

plt.xlabel('Accuracy')
plt.title('Classifier Accuracy')

sns.set_color_codes('muted')
sns.barplot(x='Accuracy', y='Classifier', data=log, color='b')

这里的指标其实有更好的选择,因为当训练集中大量例子属于某一个类时,分类器可能通过简单地全部预测为较大类来达到更高的准确率。用 F1 Score,即精确率和召回率的调和平均数会更为科学。

比较结果如图:

根据这幅图,第三次提交时我选用了效果相对较好的 SVC 分类器。

思考与总结

特征工程

特征工程的意义是找一个更好的空间去重构表达,把原始的数据对象映射到这个空间去表达,更便于你的应用。比如分类应用,最好是找到线性可分的空间。
Stark Einstein

数据可视化

在上述代码中我们基本没有进行数据可视化,一是因为这次的数据还是比较容易去分析,不太需要数据可视化;二是我对 matplotlib、seaborn 的使用还不太熟练。

实际上,数据可视化是数据科学的核心技术之一。有效的数据可视化可以帮助深入地研究变量,因此通过数据可视化来对数据集和单独的特征进行了解应该成为机器学习项目的第一步。如果想要了解如何在这个项目中使用数据可视化进行分析,可以在最下面的参考资料中查看别人的方案;如果想全面的学习数据可视化,Kaggle 也提供了数据可视化的课程

关于缺失值

不同的机器学习模型对缺失值的敏感度不同。实际上,xgboost 对缺失值有默认的处理方法。根据作者 Tianqi Chen 在论文中的介绍,xgboost 把缺失值当做稀疏矩阵来对待,本身的在节点分裂时不考虑的缺失值的数值。缺失值数据会被分到左子树和右子树分别计算损失,选择较优的那一个。如果训练中没有数据缺失,预测时出现了数据缺失,那么默认被分类到右子树。

根据知乎上怎么理解决策树、xgboost能处理缺失值?而有的模型(svm)对缺失值比较敏感呢?问题赞同数最高的回答,有一些经验法则可供参考:

  1. 树模型对于缺失值的敏感度较低,大部分时候可以在数据有缺失时使用。
  2. 涉及到距离度量(distance measurement)时,如计算两个点之间的距离,缺失数据就变得比较重要。因为涉及到“距离”这个概念,那么缺失值处理不当就会导致效果很差,如 K 近邻算法(KNN)和支持向量机(SVM)。
  3. 线性模型的代价函数(loss function)往往涉及到距离(distance)的计算,计算预测值和真实值之间的差别,这容易导致对缺失值敏感。
  4. 神经网络的鲁棒性强,对于缺失数据不是非常敏感,但一般没有那么多数据可供使用。
  5. 贝叶斯模型对于缺失数据也比较稳定,数据量很小的时候首推贝叶斯模型。

总结来看,对于有缺失值的数据在经过缺失值处理后:

  • 数据量很小,用朴素贝叶斯
  • 数据量适中或者较大,用树模型,优先 xgboost
  • 数据量较大,也可以用神经网络
  • 避免使用距离度量相关的模型,如 KNN 和 SVM

结语

亲自参加一次比赛才发觉 Kaggle 是好文明,像我这种菜鸡可以通过高质量的比赛了解到数据分析的很多基本流程和方法,收获的经验远不是书上例题所能比的。而且社区讨论里面个个都是人才,分享的经验技巧又干货满满,我超喜欢里面的。

这篇博文主要整理一开始所用的基本方法,以及通过特征工程实现的改良。下一步我打算试着去用模型集成去进一步提高预测准确率,并且进行总结。

18.03.25 补充:尝试了模型集成,但是效果没有提升(反而降了),猜想是因为 xgboost 的效果已经足够好。系列的下篇搁置中…

参考资料

最后推荐一下这个 kernel:Data ScienceTutorial for Beginners 。包含了包括数据可视化、Python 基本语法、pandas 基本用法等 Kaggle 比赛所需要的基础知识。很全面,而且主要的数据集是精灵宝可梦数据,加分。