Kaggle机器学习实战（1）——泰坦尼克号生还预测

机器学习
, , , ,

泰坦尼克号生还预测是非常推荐给入门级新手的一个热身训练。这是一个二分类任务。本文主要翻译自这篇kernel，它是kaggle上点赞数最高的一篇kernel（截至2018-03-28）。采用的工具是R语言。

读取数据

1
2
3
4
5
6
7
# 调包
library('ggplot2') # 数据可视化
library('ggthemes') # 数据可视化
library('scales') # 数据可视化
library('dplyr') # 数据操作
library('mice') # 缺失值填充
library('randomForest') # 这一步就已经剧透了，我们要用随机森林

之后我们读取csv格式的数据，假设数据存放在../input文件夹下：

1
2
3
4
5
6
7
train <- read.csv('../input/train.csv', stringsAsFactors = F)
test  <- read.csv('../input/test.csv', stringsAsFactors = F)

full  <- bind_rows(train, test) # 合并训练集和测试集

# 总览数据
str(full)

经过这一步我们就对变量有了大体的印象，哪些变量是数值型的，哪些变量是字符串型的。我们总共有12个变量，1309组数据。为了更方便理解，我们再来看一遍变量代表了什么：

  • Survived：我们要预测的Y，用1表示生还，0表示死亡
  • Pclass：头等舱/二等舱/三等舱
  • Name：姓名
  • Sex：性别
  • Age：年龄
  • SibSp：船上兄弟姐妹/配偶的数量
  • Parch：船上父母/孩子的数量
  • Ticket：船票号码
  • Fare：船票费用
  • Cabin：舱位
  • Embarked：登船的港口（瑟堡/皇后镇/南安普顿）

看过电影《泰坦尼克号》或者听说过西方宣传的所谓“骑士精神”的人，都知道泰坦尼克号上妇女和小孩的生还几率比男性要高很多；此外头等舱乘客的生还率也是显然要高于二等舱、三等舱的。因此，Sex、Age、Pclass、Fare这些变量我们想当然的把它们列入考察范围。关键就在于对剩下的一些变量的运用——这就是特征工程了，在这个实战项目中，它甚至比选择哪一种分类方法更为重要。

特征工程

姓名的处理

对于“乘客的姓名”这一变量，许多人认为这是一个跟预测结果无关紧要的变量，难道叫张三还是李四会决定一个人的生死？但是，仔细观察变量会发现：这一长串的名字可不是“张三”“李四”这么简单，它们由名字、称谓（例如Mrs., Dr.）、姓氏构成。这些字符串可以拆分出一些有用的信息，例如，称谓可以体现一个人的阶级属性；有相同姓氏的人可以归为一个家庭（外国人同姓不同族的情况比中国少多了）。因此我们要做的就是尽可能利用给出的信息，来创造变量！

1
2
3
4
5
# 从姓名中提取有用的信息（Title：称谓）
full$Title <- gsub('(.*, )|(\\..*)', '', full$Name)

# 按性别显示称谓的统计信息
table(full$Sex, full$Title)

第一步是利用R语言中的gsub()函数，将不需要的字符串替换掉：

1
gsub(pattern, replacement, x)

gsub()会替换x中所有满足pattern条件的字符串转变为replacement，在此之前我们要简单了解一下R的正则表达式，.表示任意字符，\\是转义字符，因此\\.就是一般的句点，|是逻辑或，可以简单理解为前后两个条件的并集。例如Braund, Mr. Owen Harris，我们要提取Mr，就是把第一个逗号+空格之前的字符串以及第一个句点后面的字符串全部替换为空。

我们可以看到Capt、Col、Don、Dona等等这些称谓只有极个位数的人，因此要精简一下。

1
2
3
4
5
6
7
8
9
10
11
12
# 手动把一些出现频率较低的称谓归为“稀有称谓”
rare_title <- c('Dona', 'Lady', 'the Countess','Capt', 'Col', 'Don', 
                'Dr', 'Major', 'Rev', 'Sir', 'Jonkheer')

# 同义词转换
full$Title[full$Title == 'Mlle']        <- 'Miss' 
full$Title[full$Title == 'Ms']          <- 'Miss'
full$Title[full$Title == 'Mme']         <- 'Mrs' 
full$Title[full$Title %in% rare_title]  <- 'Rare Title'

# 再按性别显示称谓的统计信息
table(full$Sex, full$Title)

于是我们得到了：男性中有61人为Master称谓，25人为稀有称谓，其余757人为普通的Mr称谓；女性中有264人为未婚Miss，198人为已婚Mrs，4人为稀有称谓。

接着，我们提取名字中的姓氏，这很有可能会与ta是否生还密切相关——因为“一家人最重要的是整整齐齐”。

1
2
3
4
5
# 提取姓氏，增设一个变量：Surname
full$Surname <- sapply(full$Name,  
                      function(x) strsplit(x, split = '[,.]')[[1]][1])

nlevels(factor(full$Surname))

共有875个不重样的姓氏。（有兴趣的话，还可以从姓氏分析出国籍、种族。）

同一家族的人

现在我们有了姓氏，再结合SibSp和Parch两个反映家族特征的变量，可以拿过来制造我们所需的变量。首先，我们创造一个表示船上有多少名家族成员的变量。

1
2
3
4
5
6
7
8
9
10
11
12
# 创建一个船上同一家族人数的变量（包含自己）
full$Fsize <- full$SibSp + full$Parch + 1

# 创造一个新变量，把家族大小跟在姓氏后面更直观一些
full$Family <- paste(full$Surname, full$Fsize, sep='_')

# 把家族人数和生还情况的关系可视化一下
ggplot(full[1:891,], aes(x = Fsize, fill = factor(Survived))) +
  geom_bar(stat='count', position='dodge') +
  scale_x_continuous(breaks=c(1:11)) +
  labs(x = 'Family Size') +
  theme_few()

通过柱状图可以非常直观地看出，家庭大小在2~4人的生还率要超过50%，而相对的，单身的和家庭超过4人的则生还率有非常明显的降低。据此，我们创建一个新的变量来更好地反映家庭大小与生还情况的关系：

1
2
3
4
5
6
7
# 建立离散的家庭大小变量
full$FsizeD[full$Fsize == 1] <- 'singleton'
full$FsizeD[full$Fsize < 5 & full$Fsize > 1] <- 'small'
full$FsizeD[full$Fsize > 4] <- 'large'

# 把离散的家庭大小变量和生还情况用马赛克图表示
mosaicplot(table(full$FsizeD, full$Survived), main='Family Size by Survival', shade=TRUE)

马赛克图更为清晰地揭示了单身和大家庭的人更难以生还的事实。如果结合年龄数据，可以挖掘出更多信息——但是年龄一项有263个缺失值，因此先跳过。

最后再看一下Cabin这个变量，它有太多的缺失值了。我们能从中挖掘的信息就是它的第一位代表甲板号A-F，也就是房间所在位置，它可能会是影响能否生还的变量。

1
full$Deck<-factor(sapply(full$Cabin, function(x) strsplit(x, NULL)[[1]][1]))

另外就是有一些带有多个舱位的数据（例如“C23 C25 C27”），我们想办法在这上面做点文章。但首先，我们要做的是处理缺失值。

缺失值

接下来我们要处理最为棘手的缺失值了。也许快刀斩乱麻是一个方法，但由于这个训练样本容量较小、特征维度也较小，我们不应该选择删除某一个训练样本或者某一特征。所以我们选择用这一行或这一列的平均数、中位数、众数来填充。选择填充缺失值的方法之前，我们可以借助可视化的手段。

用合理的值填充

我们先看登船港口这一项，有62号和830号两个缺失值。他们俩有相同的票号、票价和船舱号，所以他们是一个港口上船的。让我们把票价和登船港口以及舱位等级的图画一画，关系就清晰明了。

1
2
3
4
5
6
7
8
9
10
11
# 先把缺失项排除掉
embark_fare <- full %>%
  filter(PassengerId != 62 & PassengerId != 830)

# 把登船港口跟船舱等级和票价做个可视化
ggplot(embark_fare, aes(x = Embarked, y = Fare, fill = factor(Pclass))) +
  geom_boxplot() +
  geom_hline(aes(yintercept=80), 
    colour='red', linetype='dashed', lwd=2) +
  scale_y_continuous(labels=dollar_format()) +
  theme_few()

请在RStudio里面画一画——有了这个箱形图就直观多了。很显然，票价与登船港口和船舱等级直接相关，而80刀所对应的区间正位于头等舱-C港口的那群人的中位数附近。那么，在票价上花了80刀，并且是头等舱乘客的62号，他的登船港口我们可以放心大胆地填上C。

1
full$Embarked[c(62, 830)] <- 'C'

接着处理1044号乘客票价一栏的缺失值。很显然，我们要找和他同样的三等舱、S港口上船的乘客，比较他们的船票价格，填充一个比较合理的数值。

1
2
3
4
5
6
7
8
# 三等舱、S登船的乘客的票价，并标注出中位数
ggplot(full[full$Pclass == '3' & full$Embarked == 'S', ], 
  aes(x = Fare)) +
  geom_density(fill = '#99d6ff', alpha=0.4) + 
  geom_vline(aes(xintercept=median(Fare, na.rm=T)),
    colour='red', linetype='dashed', lwd=1) +
  scale_x_continuous(labels=dollar_format()) +
  theme_few()

好，我们完全有理由给这个缺失值填上8.05刀——三等舱S登船的票价的中位数。

预测推断

上面我们处理的都是少量的缺失值，但是看看年龄这一项，有太多的缺失值了，如果用平均数或者中位数来填充，效果肯定会不好，所以我们要建立一个回归模型来预测缺失值的年龄。

1
2
# 看一看有多少缺失值
sum(is.na(full$Age))

这里原文用了mice这个包来做预测推断。mice利用其他已知变量的数据对待填充的数据进行拟合并填充，

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
27
# 先把这些定性的变量化为factor类型
factor_vars <- c('PassengerId','Pclass','Sex','Embarked',
                 'Title','Surname','Family','FsizeD')

full[factor_vars] <- lapply(full[factor_vars], function(x) as.factor(x))

# 设置随机数种子，接下来我们用随机森林填充缺失值
set.seed(129)

# 用数据集中已知的数据拟合缺失值，首先排除一些不需要的变量
mice_mod <- mice(full[, !names(full) %in% c('PassengerId','Name','Ticket','Cabin','Family','Surname','Survived')], method='rf')

# 保存结果
mice_output <- complete(mice_mod)

# 查看补全后的和原始的年龄数值分布情况
par(mfrow=c(1,2))
hist(full$Age, freq=F, main='Age: Original Data', 
  col='darkgreen', ylim=c(0,0.04))
hist(mice_output$Age, freq=F, main='Age: MICE Output', 
  col='lightgreen', ylim=c(0,0.04))

# 结果看起来不错，所以把拟合出的结果代入到原数据集中
full$Age <- mice_output$Age

# 然后查看缺失值数量，变成0了
sum(is.na(full$Age))

这样我们就完成了年龄的缺失值填充了。前面已经说过，对于我们要预测的生还情况而言，年龄是一个相当重要的参数。所以我们可以再多做一些特征工程，创造一些更有用的变量。（下面这一步我认为是这篇kernel的精华所在！）

利用好‘年龄’变量

我们有了完整的年龄数据，再结合之前通过姓名提取的一些信息，可以创造出更具价值的变量——Child和Mother，这是两个布尔型的变量。判断一个人是否是孩子很简单，年龄小于18。判断一个人是不是妈妈需要4点：

  1. 她是女性；
  2. 她大于18岁；
  3. 她的孩子数量大于0；
  4. 她的称谓不是Miss。
1
2
3
4
5
6
7
8
9
10
11
12
# 先看一下年龄，性别和生还率的关系
ggplot(full[1:891,], aes(Age, fill = factor(Survived))) + 
  geom_histogram() + 
  facet_grid(.~Sex) + 
  theme_few()

# 创造一个变量
full$Child[full$Age < 18] <- 'Child'
full$Child[full$Age >= 18] <- 'Adult'

# 查看这个新建的变量的统计
table(full$Child, full$Survived)

从这个结果似乎看不出孩子的身份能保证一名乘客生还。接下来创造Mother这个变量，我们预想中母亲的生还率会高一些。

1
2
3
4
5
6
7
8
9
10
# 创造一个变量
full$Mother <- 'Not Mother'
full$Mother[full$Sex == 'female' & full$Parch > 0 & full$Age > 18 & full$Title != 'Miss'] <- 'Mother'

# 查看这个新建的变量的统计
table(full$Mother, full$Survived)

# 把这两个新建的变量添加到原数据集中
full$Child  <- factor(full$Child)
full$Mother <- factor(full$Mother)

保险起见可以再复查一下有没有疏漏：

1
md.pattern(full)

填补了年龄、登船港口和船票价格的缺失值，并且创建了一些相信对预测结果能起到作用的变量，接下来便可以进行训练和预测了。

训练和预测

选择何种分类模型对最终结果的精度是至关重要的，但这篇文章的作者没有告诉我们为什么选择随机森林，以及参数如何确定。

先拆分训练集和测试集。

1
2
train <- full[1:891,]
test <- full[892:1309,]

训练模型

1
2
3
4
5
6
7
8
9
10
11
12
# 设置随机数种子
set.seed(754)

# 训练模型（只挑选了有用的变量用于训练）
rf_model <- randomForest(factor(Survived) ~ Pclass + Sex + Age + SibSp + Parch + 
                                            Fare + Embarked + Title + 
                                            FsizeD + Child + Mother,
                                            data = train)

# 可视化性能随树的数量的变化
plot(rf_model, ylim=c(0,0.36))
legend('topright', colnames(rf_model$err.rate), col=1:3, fill=1:3)

这张图上黑线表示整体的错误率（大约在20%左右），红线和绿线分别表示死亡和幸存的错误率。预测死亡比预测生还要更成功。这是为什么？

变量权重

通过绘制决策树上的基尼系数的平均值来考察变量的重要性。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
# 获取重要度
importance    <- importance(rf_model)
varImportance <- data.frame(Variables = row.names(importance), 
                            Importance = round(importance[ ,'MeanDecreaseGini'],2))

# 根据重要度创建排名变量
rankImportance <- varImportance %>%
  mutate(Rank = paste0('#',dense_rank(desc(Importance))))

# 可视化重要性排名
ggplot(rankImportance, aes(x = reorder(Variables, Importance), 
    y = Importance, fill = Importance)) +
  geom_bar(stat='identity') + 
  geom_text(aes(x = Variables, y = 0.5, label = Rank),
    hjust=0, vjust=0.55, size = 4, colour = 'red') +
  labs(x = 'Variables') +
  coord_flip() + 
  theme_few()

与乘客最终是否幸存影响最大的变量竟然是Title——也就是我们根据姓名创造出来的称谓这个变量！第二名是性别，第三名是船票费用，紧接着是年龄。而预想中关系密切的“船舱等级”竟然只排到第五位，这个结果令人惊讶。

预测

训练出了模型，对测试集进行预测。

1
2
3
4
5
6
7
prediction <- predict(rf_model, test)

# 把结果保存至两列：乘客ID和是否生还
solution <- data.frame(PassengerID = test$PassengerId, Survived = prediction)

# 导出至文件
write.csv(solution, file = 'rf_mod_Solution.csv', row.names = F)

上传，查看精度和排名。（作者的精度为0.80382）

结论

  • 这篇kernel在特征工程上面花了比较多的心思。首先是根据姓名提取出了“称谓”这个变量，最后事实证明这是一个非常有效的变量。创建“家庭大小”这个变量也很有创意。之后创建了“是母亲”“是孩子”这两个有趣的变量，这是一个非常新颖的想法，我们可以融会贯通在别的实战项目里面运用类似的思想，说不定有奇效。

  • 许多高精度的kernel证明了，对于这种变量很少的训练样本，多多进行特征工程、创造新的变量是王道。本kernel没有在选择模型以及调参方面作过多叙述，因为随机森林已经足够好。运用更加复杂的集成模型能给最终结果带来一定提升。

  • 泰坦尼克号生还预测可谓是机器学习初学者必须亲手参与的经典入门试题，可以尝试更多的特征工程、运用不同的模型多加尝试。