r语言自带数据集 r语言自带数据集Titanic

摘要

本文以R语言为基础，利用数据预览，探索式数据分析，缺失值的填补，增加新特征以及去除相关特征等方法，并通过构建随机森林模型，参数调优的方式对kaggle上的泰坦尼克项目进行了生存预测，结果是得分为0.81818，前4%。

一、项目介绍

泰坦尼克生存预测是Kaggle上参赛人数较多的竞赛之一，对于数据爱好者来说是初入机器学习领域相对比较容易的比赛，属于入门级比赛项目。

比赛的目的其实很简单，就是要求参赛选手通过训练数据集分析出什么类型的人更可能幸存，并创建模型预测出测试数据集中的所有乘客是否生还。很明显，该项目属于二分类问题。

数据来源于官方网站，包括训练数据和测试数据两部分数据。

二、加载并预览数据

在加载数据之前，先载入后续分析需要用到的R包。

library("dplyr")
library("ggplot2")
library("caret")
library("randomForest")

导入数据集，将训练集和测试集分别命名为train和test。

getwd()
setwd("../Desktop")
train <- read.csv("train.csv",stringsAsFactors = F,na.strings = "")
test <- read.csv("test.csv",stringsAsFactors = F,na.strings = "")

为了对训练集和测试集进行相同的处理和转换，避免重复操作和数据类型出现不一致的情况，这里将它们合并一起处理，并命名为data数据集。

data <- bind_rows(train,test)

预览数据。

str(data)

结果如下：

'data.frame':	1309 obs. of  12 variables:
 $ PassengerId: int  1 2 3 4 5 6 7 8 9 10 ...
 $ Survived   : int  0 1 1 1 0 0 0 0 1 1 ...
 $ Pclass     : int  3 1 3 1 3 3 1 3 3 2 ...
 $ Name       : chr  "Braund, Mr. Owen Harris" "Cumings, Mrs. John Bradley (Florence Briggs Thayer)" "Heikkinen, Miss. Laina" "Futrelle, Mrs. Jacques Heath (Lily May Peel)" ...
 $ Sex        : chr  "male" "female" "female" "female" ...
 $ Age        : num  22 38 26 35 35 NA 54 2 27 14 ...
 $ SibSp      : int  1 1 0 1 0 0 0 3 0 1 ...
 $ Parch      : int  0 0 0 0 0 0 0 1 2 0 ...
 $ Ticket     : chr  "A/5 21171" "PC 17599" "STON/O2. 3101282" "113803" ...
 $ Fare       : num  7.25 71.28 7.92 53.1 8.05 ...
 $ Cabin      : chr  NA "C85" NA "C123" ...
 $ Embarked   : chr  "S" "C" "S" "S" ...

可以看出，整个数据集包含12个变量，1309条数据，其中891条为训练数据，418条为测试数据

• PassengerId：整数型变量，表示乘客的ID。
• Survived：整数型变量，表示该乘客是否幸存。0表示遇难，1表示幸存。
• Pclass： 整数型变量，表示乘客的舱位等级，1代表1等舱，2代表2等舱，3代表3等舱。
• Name： 字符型变量，表示乘客姓名。
• Sex：字符型变量，表示乘客性别。
• Age：数值型变量，表示乘客年龄。
• SibSp：整数型变量，兄弟姐妹或配偶的人数。
• Parch： 整数型变量，代表父母或子女的人数。
• Ticket： 字符型变量，代表乘客的船票号。
• Fare： 数值型变量，代表乘客的船票价。
• Cabin： 字符型变量，代表乘客所在的舱位。
• Embarked： 字符型变量，表示乘客的登船地点，C 表示 Cherbourg（瑟堡港）, Q 表示Queenstown（皇后镇）, S 表示Southampton（南安普敦）。

为了方便后续的分析，这里先统计各列的缺失值的分布情况。

apply(data,2,function(x){sum(is.na(x))}) #汇总缺失值

结果如下：

PassengerId    Survived      Pclass        Name         Sex         Age       SibSp       Parch      Ticket        Fare       Cabin    Embarked 
          0         418           0           0           0         263           0           0           0           1        1014           2

可以看出，Age，Fare，Cabin，Embarked，以及Survived均包含缺失值，值得注意的是，虽然Survived包含418个缺失值，但它是由测试数据引入的，是需要预测的部分。

三、探索性数据分析

接下来，会逐个探索不同的特征变量与预测变量之间的关系。

（1）Pclass

ggplot(data[1:891,],aes(Pclass,fill=factor(Survived)))+
  geom_bar(position = "fill")+
  ggtitle("Pclass VS Survived")

船舱等级对生存率的影响

可以看出，船舱等级为1的乘客生存率要高于船舱等级为2和3的乘客，达到62%，而3等舱的乘客生存率仅有25%。

（2）Sex

ggplot(data[1:891,],aes(Sex,fill=factor(Survived)))+
  geom_bar(position = "fill")+
  ggtitle("Sex VS Survived")

从图中可以看出，女性的存活率要远远高于男性，达到75%，而男性的存活率还不到25%。

（3）Name

Name中人名的重复度较低，无法直接使用，但值得注意的是，人名信息中还包含了称呼，可以试着将其提取出来进行分析。

data$Title <- gsub('(.*, )|(..*)','',data$Name)

简要查看提取结果。

table(data$Title)

        Capt          Col          Don         Dona           Dr     Jonkheer         Lady        Major       Master         Miss 
           1            4            1            1            8            1            1            2           61          260 
        Mlle          Mme           Mr          Mrs           Ms          Rev          Sir the Countess 
           2            1          757          197            2            8            1            1

对其进行汇总统计：

data%>%
  group_by(Title)%>%
  count()%>%
  arrange(desc(n))

结果如下：

# A tibble: 18 x 2
# Groups:   Title [18]
   Title            n
   <chr>        <int>
 1 Mr             757
 2 Miss           260
 3 Mrs            197
 4 Master          61
 5 Dr               8
 6 Rev              8
 7 Col              4
 8 Major            2
 9 Mlle             2
10 Ms               2
11 Capt             1
12 Don              1
13 Dona             1
14 Jonkheer         1
15 Lady             1
16 Mme              1
17 Sir              1
18 the Countess     1

可以看出，Mr，Miss，Mrs，Master这4类称呼所占的数目最多，可以单独拿来进行分析，而其他称呼的数目太少，不太具备分析价值，因此将其合并为其他类，命名为“Other”。

data$Title <- ifelse(data$Title %in% c("Mr","Miss","Mrs","Master"),data$Title,"Other")

画图：

ggplot(data[1:891,],aes(Title,fill=factor(Survived)))+
  geom_bar(position = "fill")+
  ggtitle("Title VS SUrvived")

不同称呼对存活率的影响

可以看出，Mr的生存率最低，而Mrs和Miss的存活率相对较高，这和之前女性的存活率高于男性是相互吻合的。

（4）Fnumber

由于SibSp和Parch均表示乘客的亲属，因此，可以考虑将它们合起来进行分析，即为家族人数，命名为“Fnumber”。

data$Fnumber <- data$SibSp+data$Parch+1 # 家族人数=旁系亲属+直系亲属+自己

不同家族人数的存活率分布情况

从图中可以看出，家族人数在2-4人时的幸存率较高，而在1个人以及5人以上时存活率明显降低，当人数增加到8人以上时，死亡率达到100%。

（5）Fare

通过之前的统计得知，Fare列中存在1个缺失值，先定位缺失值出现的位置。

data[is.na(data$Fare),]

     PassengerId Survived Pclass               Name  Sex  Age SibSp Parch Ticket Fare Cabin Embarked Title Fnumber Fsize
1044        1044       NA      3 Storey, Mr. Thomas male 60.5     0     0   3701   NA  <NA>        S    Mr       1 Alone

结果显示缺失值出现在第1044行，其船舱等级为3，登船地点为S。在这里，可以假设相同登船地点中相同的船舱等级票价是相似的，因此，利用同为S点登船且同为3等舱的乘客票价平均值对其进行补齐。

data$Fare[1044] <- mean(data[data$Pclass=='3'&data$Embarked=='S',]$Fare,na.rm = T)

画图：

fare_cut <- cut(data[1:891,]$Fare,breaks = c(-1,25,50,100,600),
                labels = c("<25","25-50","50-100",">100"))

ggplot(data[1:891,],aes(fare_cut,fill=factor(Survived)))+
  geom_bar(position = "fill")+
  labs(x="票价")+
  ggtitle("Fare VS Survived")

不同票价的存活率分布情况

可以看出，随着票价的提升，存活率也逐渐增加。

（6）Embarked

通过之前的统计得知，Embarked列中存在2个缺失值，先定位缺失值出现的位置。

data[is.na(data$Embarked),]

    PassengerId Survived Pclass                                      Name    Sex Age SibSp Parch Ticket Fare Cabin Embarked Title Fnumber Fsize
62           62        1      1                       Icard, Miss. Amelie female  38     0     0 113572   80   B28     <NA>  Miss       1 Alone
830         830        1      1 Stone, Mrs. George Nelson (Martha Evelyn) female  62     0     0 113572   80   B28     <NA>   Mrs       1 Alone

结果显示这2个缺失值分别在第62和第830行。通过比较这2行数据，发现它们的船舱等级均为1，票价均为80。和Fare的处理方式类似，可以利用同船舱等级和同票价的乘客登船点对其进行补齐。

data$Embarked <- toupper(data$Embarked) # 将小写字母转变成大写字母
pclass1 <- subset(data,Pclass==1) # 筛选出船舱等级为1的乘客
summarise(group_by(pclass1,Embarked),mfare=median(Fare),n()) # 统计不同登船点的票价中位数

结果如下：

Embarked mfare `n()`
  <chr>    <dbl> <int>
1 C         76.7   141
2 Q         90       3
3 S         52     177
4 NA        80       2

可以看出，登船点为C的票价更加接近于80，因此，用“C”对其补齐。

data[c(62,830),"Embarked"] <- "C"

画图：

ggplot(data[1:891,],aes(Embarked,fill=factor(Survived)))+
  geom_bar(position = "fill")+
  ggtitle("Embarked VS Survived")

不同登船地点的存活率分布情况

可以看出，从C地点登船的乘客存活率在50%以上，要高于Q和S，而S地点的乘客存活率较低，为33%。

（7）Age

通过之前的统计得知，Age列存在264个缺失值，缺失量较大，如果直接使用整个群体的年龄特征值进行补齐，误差会比较大。考虑到不同名头的群体可能具有不同的年龄分布，因此这里分别计算出不同名头人群年龄的中位数对缺失值进行填补。

title.age <- aggregate(data$Age,by = list(data$Title), FUN = function(x){median(x, na.rm = T)})

data[is.na(data$Age), "Age"] <- apply(data[is.na(data$Age), ] , 1, function(x) {title.age[title.age[, 1]==x["Title"], 2]})

画图：

age_cut <- cut(data[1:891,]$Age,breaks = c(0,10,20,30,40,50,60,90),
               labels = c("<10","10-20","20-30","30-40","40-50","50-60",">60"))

ggplot(data[1:891,],aes(age_cut,fill=factor(Survived)))+
  geom_bar(position = "fill")+
  ggtitle("Age VS Survived")

不同年龄区间和生存率之间的关系

可以看出，10岁以下的孩童存活率最高，其次是30-40岁之间的乘客，而60岁以上的乘客存活率最低，不到25%。

（8）Ticket & Cabin

由于Ticket列规律性不强，并且Cabin列缺失值高达1014个，因此这2个特征变量不在后续模型的考虑范围之内。

data <- data[,c(-9,-11)]

四、建模预测

重新将处理之后的data数据集拆分成训练集和测试集。

data$Sex <- factor(data$Sex)
data$Embarked <- factor(data$Embarked)
data$Title <- factor(data$Title)
data$Survived <- factor(data$Survived)

train <- data[1:891,] # 训练集
test <- data[892:1309,] # 测试集

用整个训练集训练模型，随机选取30%的训练数据用来验证模型。

set.seed(123)

par <- createDataPartition(train$Survived,times = 1,p=0.3,list = F)

test_data <- train[par,] # 30%的训练集数据，用来验证模型

特征变量的选择

在第三步中，通过创造派生特征和去掉规律性不强以及缺失值较多的变量，同时将不具备分析价值的ID列和需要进行预测的Survived列排除，最终获得了11个可供训练模型的变量。如下所示：

那么应该选择哪些特征呢？

（1）去掉关联特征。由于Fnumber是用过SibSp和Parch两列共同计算得出，因此可以尝试将SibSp和Parch两列移除。另外，Title列是由Name得出，因此将Name移除。

（2）特征离散化。由于Age和Fnumber两列属于连续型数值变量，为了降低异常值的影响并且减少模型过拟合的风险，可以考虑将这两列进行特征离散化，同时也增加了新的特征。

对于Age列，以18岁作为分界线，将其划分为未成年人和成年人两个群体。

data$Is_Child <- ifelse(data$Age<18,1,0)

对于Fnumber列，可以按照人数的不同划分成“alone”，“small”，"big"3类，取名为“家族规模”，即Fsize。

data$Fsize <- NA
data$Fsize[data$Fnumber == 1] <- 'Alone'
data$Fsize[data$Fnumber < 5 & data$Fnumber > 1] <- 'Small'
data$Fsize[data$Fnumber > 4] <- 'Big'

其他变量不做处理。

通过以上操作，最终选取了Pclass，Sex，Fare，Embarked，Is_child，Fsize，Title这7个变量。

建立模型

set.seed(1234)

rf500 <- randomForest(Survived~Pclass+Sex+Title+Fare+Fsize+Embarked+Is_Child,data = train)

利用验证集对模型进行验证，并构建混淆矩阵看看模型的效果如何。

pre_rf <- predict(rf500,test_data)

confusionMatrix(pre_rf,test_data$Survived)

结果如下：

Confusion Matrix and Statistics

          Reference
Prediction   0   1
         0 152  27
         1  13  76
                                          
               Accuracy : 0.8507          
                 95% CI : (0.8024, 0.8912)
    No Information Rate : 0.6157          
    P-Value [Acc > NIR] : < 2e-16         
                                          
                  Kappa : 0.6763          
                                          
 Mcnemar's Test P-Value : 0.03983         
                                          
            Sensitivity : 0.9212          
            Specificity : 0.7379          
         Pos Pred Value : 0.8492          
         Neg Pred Value : 0.8539          
             Prevalence : 0.6157          
         Detection Rate : 0.5672          
   Detection Prevalence : 0.6679          
      Balanced Accuracy : 0.8295          
                                          
       'Positive' Class : 0

可以看到，该模型的准确率达到85.07%，接下来对测试集进行预测。

prf <- predict(rf500,test)

该模型的预测结果在kaggle得分为0.80861。

参数调优

值得注意的是，随机森林模型中树的个数默认为500，可以尝试改变该参数，看看结果是否有提升。

当树的数量为1000时：

构建模型

set.seed(1235)

rf1000 <- randomForest(Survived~Pclass+Sex+Title+Fare+Fsize+Embarked+Is_Child,data = train,ntree=1000)

验证结果如下：

Confusion Matrix and Statistics

          Reference
Prediction   0   1
         0 152  27
         1  13  76
                                          
               Accuracy : 0.8507          
                 95% CI : (0.8024, 0.8912)
    No Information Rate : 0.6157          
    P-Value [Acc > NIR] : < 2e-16         
                                          
                  Kappa : 0.6763          
                                          
 Mcnemar's Test P-Value : 0.03983         
                                          
            Sensitivity : 0.9212          
            Specificity : 0.7379          
         Pos Pred Value : 0.8492          
         Neg Pred Value : 0.8539          
             Prevalence : 0.6157          
         Detection Rate : 0.5672          
   Detection Prevalence : 0.6679          
      Balanced Accuracy : 0.8295          
                                          
       'Positive' Class : 0

可以看出来验证的结果和之前没有什么变化。

当树的数量为1500时，验证结果如下：

Confusion Matrix and Statistics

          Reference
Prediction   0   1
         0 153  26
         1  12  77
                                          
               Accuracy : 0.8582          
                 95% CI : (0.8106, 0.8977)
    No Information Rate : 0.6157          
    P-Value [Acc > NIR] : < 2e-16         
                                          
                  Kappa : 0.6925          
                                          
 Mcnemar's Test P-Value : 0.03496         
                                          
            Sensitivity : 0.9273          
            Specificity : 0.7476          
         Pos Pred Value : 0.8547          
         Neg Pred Value : 0.8652          
             Prevalence : 0.6157          
         Detection Rate : 0.5709          
   Detection Prevalence : 0.6679          
      Balanced Accuracy : 0.8374          
                                          
       'Positive' Class : 0

准确率相比于前两种有稍许提升，达到85.82%，该模型预测结果在kaggle上得分为0.81818，排在第464位（464/12951），TOP 4%，比之前提升不少。

当树的数量为2000时，预测结果如下：

Confusion Matrix and Statistics

          Reference
Prediction   0   1
         0 153  28
         1  12  75
                                          
               Accuracy : 0.8507          
                 95% CI : (0.8024, 0.8912)
    No Information Rate : 0.6157          
    P-Value [Acc > NIR] : < 2e-16         
                                          
                  Kappa : 0.6751          
                                          
 Mcnemar's Test P-Value : 0.01771         
                                          
            Sensitivity : 0.9273          
            Specificity : 0.7282          
         Pos Pred Value : 0.8453          
         Neg Pred Value : 0.8621          
             Prevalence : 0.6157          
         Detection Rate : 0.5709          
   Detection Prevalence : 0.6754          
      Balanced Accuracy : 0.8277          
                                          
       'Positive' Class : 0

准确率有所下降，因此没有用来做测试。