Titanic幸存預測分析(Kaggle)

分享一篇kaggle入門級案例，泰坦尼克號幸存遇難分析。

參考文章： 技術世界，原文鏈接 http://www.jasongj.com/ml/classification/

 

• 案例分析內容：

          通過訓練集分析預測什么人可能生還，並對測試集中乘客做出預測判斷

• 案例分析

• 加載包

  library(dplyr) #bind_rows()
  library(ggplot2) #繪圖
  library(ggthemes)
  library(InformationValue) #計算WOE和IV
  library(stringr) #數據處理
  library(rpart)  #預測乘客年齡
  library(scales)  #dollar_format()
  library(party)   #cforest()
  library(gbm)  #AdaBoost
  library(MLmetrics) # Mache learning metrics.e.g. Recall, Precision, Accuracy, AUC

• 加載文件

  train <- read.csv("F:\\R/泰坦尼克幸存分析/train.csv",header = T,stringsAsFactors = F) #ID 1~891乘客信息
  test <- read.csv("F:\\R/泰坦尼克幸存分析/test.csv",header = T,stringsAsFactors = F)  #ID 892~1309號乘客信息(缺少是否存活信息)
  test_survived <-read.csv("F:/R/泰坦尼克幸存分析/gender_submission.csv",header = T,stringsAsFactors = F)  #ID 892~1309號 是否存活信息
  

•  數據整理

  #合並訓練數據和測試數據
  data <- bind_rows(train,test) 
  ##Sex:性別,Age:年齡,SibSP:配偶/兄妹數,Parch:父母/子女數,Ticket:船票號
  ##Fare：費用,Cabin：艙位區域,Pclass:艙位等級,Embarked:到達碼頭,Title:頭銜
  #將是否存活設為因子
  data$Survived <- as.factor(data$Survived)
  train$Survived <-as.factor(train$Survived)
  test$Survived <- as.factor(test$Survived)

• 統計幸存和遇難人數是否與艙位等級有關

  ggplot(data = data[1:nrow(train),],aes(Pclass,..count..,fill=factor(Survived)))+                 #載入訓練數據分析
    geom_bar(stat = 'count',position = 'dodge')+        
    labs(title='艙位等級對乘客存活影響',x='艙位等級',y='存活人數',fill='Survived')+                      #fill為圖例標題屬性                                
    scale_fill_discrete(limits=c(0,1),labels=c('遇難','獲救'))+                                      #修改圖例標簽文本
    scale_x_continuous(breaks=c(1,2,3),labels=c('頭等艙','二等艙','三等艙'))+                          #修改X軸刻度文本
    geom_text(stat = "count",aes(label=..count..),position = position_dodge(width = 1),vjust=-0.3)+ #添加數據標簽
    theme(plot.title = element_text(hjust = 0.5))                                                   #修改標題位置

　　　　　可以看到,頭等艙的乘客獲救率是最高的，艙位等級越高，獲救幾率越大

 

　　

 

• 計算艙位等級(Pclass)的WOE和IV

  class(data$Pclass)  #查看變量(艙位)類型,求WOE時需要轉換為因子
  WOETable(X = factor(data$Pclass[1:nrow(train)]),Y = data$Survived[1:nrow(train)])
  IV(X = factor(data$Pclass[1:nrow(train)]),Y = data$Survived[1:nrow(train)] )

  

• 為了更為定量的計算Pclass的預測價值，可以算出Pclass的WOE和IV如下。從結果可以看出，Pclass的IV為0.5，且“Highly Predictive”。由此可以暫時將Pclass作為　預測模型的特征變量之一。

 

 

• 統計不同title(頭銜)的乘客存活率

• 訓練集中給出了乘客姓名，其中含有MR，Capt等常見稱號，這通常標志着一個人處於的社會階層，所以猜測可能與存活率存在一定聯系。接下來要進行分類整理。提取出Name中的title標簽，並進行分類。

• data$Title <- sapply(data$Name,FUN=function(x){strsplit(x,split = '[,.]')[[1]][2]}) #依次提取出每行的title標簽
  #head(strsplit(data$Name,split = '[,.]')[[1]][2])
  head(data$Title)
  data$Title <- sub(pattern = ' ',replacement = '',data$Title)
  data$Title[data$Title %in%c('Mme','Mlle')] <-'Mlle'
  data$Title[data$Title %in%c('Capt','Don','Major','Sir')] <-'Sir'
  data$Title[data$Title%in%c('Dona','Lady','thhe Countess','Jonkheer')] <-'Lady'
  data$Title <- factor(data$Title)

• 抽取完乘客Title后,繪圖觀察

  ggplot(data = data[1:nrow(train),],aes(x = Title,y = ..count..,fill=factor(Survived)))+
    geom_bar(stat = 'count')+
    geom_text(stat = 'count' ,aes(label=..count..),position = position_stack(vjust = 0.85))+
    labs(title='頭銜是否影響存活率',x='尊稱/頭銜',y='人數',fill='Survived')+
    theme(plot.title = element_text(hjust =0.55))+
    scale_fill_discrete(limit=c(0,1),labels=c("遇難","獲救"))+
    theme_economist()

  

  觀察圖中不難發現，圖中Master，Miss,Mlle,Mrs,Ms獲救比例均超過50%，而Mr的獲救比例不到15.7%。接下來計算WOE和IV，

          查看Title這一變量對於最終的預測是否有用

 

• 計算頭銜(Title)的WOE和IV

 　　　　WOETable(X = factor(data$Title[1:nrow(train)]),Y = factor(data$Survived[1:nrow(train)]))

　  　　 2 IV(X = factor(data$Title[1:nrow(train)]),Y = factor(data$Survived[1:nrow(train)]) ) 

            

             IV為1.520702，且”Highly Predictive”。因此，可暫將Title作為預測模型中的一個特征變量。

• 猜測性別和存活率有關

  ggplot(data = data[1:nrow(train),],aes(x =Sex,y = ..count..,fill=factor(Survived)))+
    geom_bar(stat='count',position = 'fill')+
    geom_text(stat = 'count',aes(label=..count..),position = 'fill',vjust=1)+
    labs(title="性別是否影響存活率",fill="Survived",x='性別',y='獲救比例')+
    scale_x_discrete(breaks = c('female','male'),labels = c('女','男'))+
    scale_fill_discrete(limits=c(0,1),labels=c('遇難','獲救'))

　　　　泰坦尼克號遇難之際，船上乘客秉承‘女士優先’的原則，實際情況是，75%的女性乘客獲救，而僅有不到25%的男性乘客獲救，這也充分說明了這               一原則的真實性。

　　　　

 

• 計算性別(Sex)的WOE和IV

　　　　  WOETable(X = factor(data$Sex[1:nrow(train)]),Y = factor(data$Survived[1:nrow(train)]))

　　　　

 　　　　 IV(X =factor(data$Sex[1:nrow(train)]),Y = factor(data$Survived[1:nrow(train)]) ) 

　　　　   

　　　　為高預測變量

• 統計年齡與存活率是否有關

  summary(data$Age[1:nrow(train)])
  ggplot(data = data[!is.na(data$Age),],aes(Age,linetype=Survived,color=Survived))+
    geom_line(stat='bin',bins=10,size=0.8)+
    labs(title='年齡是否與存活率有關',x='年齡',y='人數',color="Survived",linetype="Survived")+
    scale_color_discrete(limits=c(0,1),labels=c('遇難','獲救'))+
    scale_linetype_discrete(limits=c(0,1),labels=c('遇難','獲救'))+
    theme_stata()

  除了女士優先，老弱人士可能也是優先照顧的對象，圖中顯示，20歲以下的人員獲救比例確實較高，而25歲左右的青年人士獲救人數最多，但遇難的人數也接近200人。

　　　　

• 統計(SibSp)配偶/兄弟姐妹人數同時在船對存活率是否有影響

  ggplot(data = train,aes(x=as.factor(train$SibSp),fill=Survived))+geom_bar(stat='count',position = 'dodge')+
    geom_text(stat = 'count',position = position_dodge(width = 1),aes(label=..count..),vjust=-0.1)+
    labs(x='親屬人數',y='人數',title='配偶/兄弟姐妹人數對存活率是否有影響')+
    scale_fill_discrete(limits=c(0,1),labels=c('遇難','獲救'))+
    theme(plot.title = element_text(hjust = 0.5))

  訓練集中提供了乘客配偶或兄弟姐妹的人數，觀察后發現沒有親屬在船上的人數較多，沉船時，獨身出行的乘客獲救幾率只有34%，有1~2名配偶或兄弟姐妹同時在船上時，該名乘客獲救幾率也較高。而人數達 4人以上時，幾乎同時遇難。

• 統計SibSp的WOE和IV

 　　　　WOETable(X = factor(train$SibSp),Y = factor(train$Survived))     

　　　　　

　　　　 IV(X = factor(train$SibSp),Y = factor(train$Survived)) 

              IV為0.1449，為高預測性變量

　　　　

• 統計Parch(父母/子女人數)對存活率影響

  ggplot(data = train,aes(Parch,fill=Survived))+
    geom_bar(stat='count',position = 'dodge')+
    labs(title='父母/子女數對存活率是否有影響',x='父母/子女數',y='人數')+
    geom_text(aes(label=..count..),stat = 'count',position = position_dodge(width = 1),vjust=-0.1)+
    scale_fill_discrete(limits=c(0,1),label=c('遇難','獲救'))+
    theme(plot.title = element_text(hjust = 0.5))

  Parch列中提供的為乘客的父母/子女人數(同時在船),探究是否該變量會影響存活率。由圖可看出，當船上沒有自己的父母或者子女時，乘客存活率與SibSp情況相仿，不足1/3。當船上Parch數為1~3人時，獲救率高於50%。

　　　　

• 計算Parch的WOE和IV

　　　　 WOETable(X = factor(train$Parch),Y = factor(train$Survived))

　　　　 IV(X = factor(train$Parch),Y = factor(train$Survived)) 

　　　　計算Parch得0.116，認為高預測變量

　　　　

• 找出Ticket與存活率關系,共享船票號的可能為一家人,單獨船票為獨身一人,分成兩組進行比較。

  ticket.count <- aggregate(data$Ticket,by=list(data$Ticket),function(x)sum(!is.na(x)))
  #整合船票號,記錄重復的次數,ticket.count記錄這兩列(有序),但data中船票號分布是無序的
  head(ticket.count)
  data$TicketCount <- apply(X = data,MARGIN = 1,FUN = function(x)ticket.count[which(ticket.count[,1]==x['Ticket']),2])
  #主體(X)為data,將ticket.count中的船票號(有序)與data$Ticket(無序)進行一一對應
  head(data$TicketCount)
  data$TicketCount <- factor(sapply(X = data$TicketCount,FUN = function(x)ifelse(x>1,'Share','Unique')))
  #重復次數>1則說明為共享船票,=1為獨自一人.比較兩組人員的存活率.

  數據集中提供了Ticket列，提供了乘客的船票號。整合船票號，發現存在重復的船票號，猜想可以與家庭共享船票號有關。前面得存活率與SibSP和Parch有關，現可將Ticket分成兩類，一類為家庭共享船票，一類為獨自乘船所用船票號。

　　 1 #重復次數>1則說明為共享船票,=1為獨自一人.比較兩組人員的存活率.
　　 2 ggplot(data,aes(TicketCount,..count..,fill=factor(Survived)))+
　　 3   geom_bar(stat = 'count',position = 'dodge')+
　　 4   labs(title='船票號與存活率聯系',x='船票號',y='人數',fill='Survived')+
　　 5   geom_text(stat = 'count',aes(label=..count..),position = position_dodge(width = 0.9),vjust=-0.1)+
　　 6   scale_fill_discrete(limits=c(0,1),labels=c('遇難','獲救'))+
　　 7   theme(plot.title = element_text(hjust = 0.5))

　　　　圖可看出，共用一張船票的家庭，存活率為50%，而單張船票(即獨自出行)的乘客，遇難的可能性高達73%

　　　　　　

• 計算TicketCount的WOE和IV

   WOETable(X = factor(data$TicketCount),Y = factor(data$Survived))
  ##      CAT GOODS BADS TOTAL     PCT_G     PCT_B        WOE        IV
  ## 1  share   308  288   596 0.6234818 0.3533742  0.5677919 0.1533649
  ## 2 unique   186  527   713 0.3765182 0.6466258 -0.5408013 0.1460745

  > IV(X = factor(data$TicketCount),Y = factor(data$Survived))
  ## [1] 0.2994394
  ## attr(,"howgood")
  ## [1] "Highly Predictive"

  IV為0.29，且為Highly Predictive　

• 統計船費(Fare)和存活率關系

　　　　船費與艙位等級和行程距離有關，已知存活率與艙位等級(Pclass)存在一定關系,猜想船費可能也存在關系　　

　　 1 summary(data$Fare)
　　 2 class(data$Fare)
　　 3 ggplot(data[!is.na(data$Fare),],aes(x = Fare,color=factor(Survived)))+geom_line(stat = 'bin',binwidth=10,size=1)+
　　 4   labs(title='船費是否影響存活率',x='船費',y='人數',color='Survived')+
　　 5   scale_color_discrete(labels=c('遇難','獲救'))+
　　 6   theme(plot.title = element_text(hjust=0.5))

　　　　　　由圖可看出，船費超過100元的乘客幾乎都獲救　　　　

　　

• 計算Fare的WOE和IV

　　　　 WOETable(X = factor(data$Fare),Y = data$Survived) 　

　　　　IV(X = factor(data$Fare),Y = data$Survived)
　　　　[1] 0.709573
　　　　attr(,"howgood")
　　　　[1] "Highly Predictive"

　　　　　　同樣Fare為高預測變量

• 統計艙位區域(Carbin)對存活率影響

　　　　對於Cabin變量，其值以字母開始，后面伴以數字。這里有一個猜想，字母代表某個區域，數據代表該區域的序號。類似於火車票即有車箱號又有座位           號。因此，這里可嘗試將Cabin的首字母提取出來，並分別統計出不同首字母倉位對應的乘客的幸存率。

　　 1 data$Cabin_level <- substr(x = data$Cabin,start = 1,stop = 1)
　　 2 ggplot(data,aes(data$Cabin_level,fill=Survived))+geom_bar(stat = 'count',position = 'dodge')+
　　 3   geom_text(stat = 'count', aes(label=..count..),position = position_dodge(width = 1),vjust=-0.1)+
　　 4   labs(title='艙位區域對存活率影響',x='艙位',y='人數')+
　　 5   scale_fill_discrete(label=c('遇難','獲救'))

　　　　Cabin變量中存在的空字符串較多，分析其他得B,C,D,E艙的乘客幸存率遠高於50%，其他艙的乘客則低於50%。

　　　　

• 計算data$Cabin_level的WOE和IV

 　　　　WOETable(X = factor(data$Cabin_level),Y = data$Survived) 

　　　　

　　　　 IV(X = factor(data$Cabin_level),Y = data$Survived) 

　　　　

• 統計登船碼頭是否與存活率有關

  ggplot(train,aes(Embarked,fill=Survived))+geom_bar(stat = 'count',position = 'dodge')+
    geom_text(stat = 'count',aes(label=..count..),position = position_dodge(width = 1),vjust=-0.1)

　　　　到達C碼頭的乘客獲救率高於50%，而到達S碼頭的乘客遇難人數達427人，幸存率僅有29%　　　　

• 計算Embarked(登船碼頭)WOE和IV

　　　　

• 列出所有缺失數據

　　　　研究完變量后，接下來要對缺失數據進行處理

attach(data)
head(missing)
missing<- list(Pclass=nrow(data[is.na(Pclass),]))
missing$Name <- nrow(data[is.na(Name),])
missing$Sex <- nrow(data[is.na(Sex),])
missing$Age <- nrow(data[is.na(Age),])
missing$SibSp <-  nrow(data[is.na(SibSp),])
missing$Parch <-  nrow(data[is.na(Parch),])
missing$Ticket <-  nrow(data[is.na(Ticket),])
missing$Fare <-  nrow(data[is.na(Fare),])
missing$Cabin <- nrow(data[which(data$Cabin==''),])
missing$Embarked <-  nrow(data[which(data$Embarked==''),])
#names(missing)
#missing[["Cabin"]][1]
for (name in names(missing)) {
  if(missing[[name]][1]>0){
    print(paste('',name,' miss ',missing[[name]][1],' values',sep=''))
  }
}
detach(data)

　　　　　　

• 預測乘客年齡

　　　　　　 乘客年齡數據共缺失263條，缺失量較大，不適合使用中位數或均值填補，通過使用其它變量預測或者直接將缺失值設置為默認值的方法填補，這                       里通過其它變量來預測缺失的年齡信息。

　　　　1 age.model <- rpart(Age~Pclass+factor(Sex)+SibSp+Parch+Fare+factor(Embarked)+Title,data = data[!is.na(data$Age),],method = 'anova')
　　　　2 data$Age[is.na(data$Age)]
　　　　3 data$Age[is.na(data$Age)] <- predict(age.model,data[is.na(data$Age),])

• 中位數填補缺失的Embarked值

　　　　　　查看缺失碼頭,發現船費都為80,猜想船費與艙位和到達碼頭有關。繪圖查看后發現到達碼頭C的頭等艙船票為80,可以將該缺失的空值補為C

　　　　1 ggplot(data[which(data$Embarked!=''),],aes(Embarked,Fare,fill=factor(Pclass)))+
　　　　2   geom_boxplot()+
　　　　3   geom_hline(yintercept = 80,color='red',linetype=2,lwd=1)+
　　　　4   scale_y_continuous(labels = dollar_format())+
　　　　5   labs(title='船費和艙位及登船碼頭的關系',x='登船碼頭',y='船費',fill='艙位等級')+
　　　　6   theme(plot.title = element_text(hjust=0.5),panel.grid.major = element_blank())+
　　　　7   scale_fill_discrete(label=c('頭等艙','二等艙','三等艙'))

　　　　　　　　 data$Embarked[which(data$Embarked=='')] <- 'C'

　　　　　　　　　 2 data$Embarked <- as.factor(data$Embarked) 　　　　　　

　　　　　　

• 補船費的缺失值

　　　　船費和艙位等級，到達碼頭存在聯系，已知另外兩個條件，不難猜出船費為多少，將缺失的船費的數據補齊

　　　　1 data[is.na(data$Fare),c('Pclass','Embarked')]
　　　　2 summary(data[which(data$Pclass=='3'&&data$Embarked=='S'),'Fare'])
　　　　3 data[is.na(data$Fare),'Fare'] <-7.25

• 補Cabin(設為默認值)

 　　　　　　因為除去這些缺失值后,測得IV已較高,所以可直接設為一個默認值

　　　　1 summary(data$Cabin)
　　　　2 head(data$Cabin)
　　　　3 data$Cabin <- as.factor(sapply(data$Cabin,FUN = function(x) ifelse(x=='','X',str_sub(x,1,1))))

• 訓練模型

  set.seed(123)
  class(data$Embarked)
  data$Sex <- as.factor(data$Sex)
  model <- cforest(Survived~Pclass+Title+Sex+Age+SibSp+Parch+TicketCount+Fare+Cabin+Embarked,data,controls = cforest_unbiased(ntree=2000,mtry=3)　

• 交叉驗證

  cv.summarize <- function(data.true, data.predict) {
    print(paste('Recall:', Recall(data.true, data.predict)))
    print(paste('Precision:', Precision(data.true, data.predict)))
    print(paste('Accuracy:', Accuracy(data.predict, data.true)))
    print(paste('AUC:', AUC(data.predict, data.true)))

• 預測

  predict.result <-predict(model,data[(1+nrow(train)):(nrow(data)),],OOB=TRUE,type='response')
  output <- data.frame(PassengerID=test$PassengerId,Survived=predict.result)
  write.csv(output,file ='F:/R/泰坦尼克幸存分析/cit1.csv',row.names = FALSE)