泰坦尼克號倖存預測

01-25

一、背景知識

著名的奧林匹克級郵輪Titanic在首次航行過程中，撞上冰山沉入海底，災難發生時候，船上的救生艇數量有限，只有一小部分乘客存活了下來，在救援過程中本著「婦女和兒童優先」的原則，現在就來預測一下倖存的乘客是哪些。

二、數據獲取

給定了兩個數據集，分別為訓練集（891條記錄，12個變數）和測試集（418條記錄，11個變數），測試集中缺少一個目標變數（是否存活），這需要我們去預測。通過分析訓練集構建一個機器學習模型，再將該模型用於測試集，預測測試集中的乘客是否存活。在kaggle上對比模型得到的數據與真實數據，得到該模型正確預測的比例。

數據來源：https://www.kaggle.com/c/titanic/data

library(readr)train <- read_csv("H:/R代碼/Titanic/train.csv")test <- read_csv("H:/R代碼/Titanic/test.csv")names(train) #查看變數名

變數含義：

三、數據預處理

1. 合併訓練集與測試集

為統一數據的處理過程，選擇將訓練集與數據集合併，由於測試集缺少一列Survived，先需添加一列Survived再進行合併

test$Survived <- NAcombi <- rbind(train,test)

2. 轉換數據類型

str(combi) #查看數據類型

其中的Survived、Pclass、Sex、Embarked屬於分類變數，應轉化為因子類型

combi$Survived <- as.factor(combi$Survived) combi$Pclass <- as.factor(combi$Pclass) combi$Sex <- as.factor(combi$Sex) combi$Embarked <- as.factor(combi$Embarked)

3. 查看缺失值

sapply(combi,function(x) sum(is.na(x)))

Survived列有418個缺失值，是需要預測的
Age列有263個缺失值
Fare列有1個缺失值
Cabin列有1014個缺失值
Embarked列有2個缺失值

4. 處理缺失值

Fare列只有1個缺失值，屬於數值型數據，可以利用中位數進行填補，在計算中位數時要先把缺失值去掉

combi$Fare[is.na(combi$Fare)] <- median(combi$Fare,na.rm = TRUE)

Embarked列有2個缺失值，在進行填補前先提取這兩行數據的相關信息

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

發現Pclass都為1，Fare都為80，而Embarked為C的Pclass屬於1的Fare中位數正好為80，所以將缺失值填補為C

combi$Embarked[is.na(combi$Embarked)] <- "C"

Cabin列有1014個缺失值，對預測沒有什麼意義，可以忽略
Age列有263個缺失值，由於缺失過多，需要通過其他變數來預測，填補方法會在後文專門介紹

5. 分析各變數對Survived的影響

（1）Pclass的影響

ggplot(combi[(1:891),], aes(x = Pclass, fill = Survived)) +geom_bar(position = "dodge") +labs(title = "The impacts of Pclass",x = "Pclass", y = "count") +theme(plot.title = element_text(hjust = 0.5)) #標題居中

從上圖可以知道，Pclass=1的存活率最高，Pclass=2的存活率次之，Pclass=3的存活率最低。

（2）Name的影響

Name的重複度比較低，通過提取特徵值Title來進行分類，並將人數較少的類重新進行定義：

combi$Title <- sapply(combi$Name, FUN = function(x) {strsplit(x, split = "[,.]")[[1]][2]})combi$Title <- sub(" ", "", combi$Title)table(combi$Title)combi$Title[combi$Title %in% c("Mme", "Mlle")] <- "Mlle" combi$Title[combi$Title %in% c("Capt", "Don", "Major", "Sir")] <- "Sir"combi$Title[combi$Title %in% c("Dona", "Lady", "the Countess", "Jonkheer")] <- "Lady"combi$Title <- factor(combi$Title)

繪製不同Title乘客的存活情況：

ggplot(combi[(1:891),], aes(x = Title, fill = Survived)) +geom_bar(stat = "count",position = "dodge") +labs(title = "The impacts of Title",x = "Title", y = "count") +theme(plot.title = element_text(hjust = 0.5)) #標題居中

從上圖可以知道，Title為Mr的乘客存活率非常低，而Title為Miss和Mrs的乘客存活率非常高。

（3）Sex的影響

ggplot(combi[(1:891),], aes(x = Sex, fill = Survived)) +geom_bar(stat = "count",position = "dodge") +labs(title = "The impacts of Sex",x = "Sex", y = "count") +theme(plot.title = element_text(hjust = 0.5)) #標題居中

從上圖可以知道，女性的存活率非常高，而男性的存活率非常低，很好的體現了「女士優先」的救援原則。

（4）Age的影響(注意數據是train)

train$Survived <- as.factor(train$Survived)ggplot(train[!is.na(train$Age),], aes(x = Age, color = Survived)) +geom_line(aes(label = ..count..),stat = "bin",binwidth = 5,na.rm = TRUE) +labs(title = "The impacts of Age",x = "Age", y = "count") +theme(plot.title = element_text(hjust = 0.5)) #標題居中

（5）SibSp的影響

ggplot(combi[(1:891),], aes(x = SibSp, fill = Survived)) +geom_bar(stat = "count",position = "dodge") +labs(title = "The impacts of SibSp",x = "SibSp", y = "count") +theme(plot.title = element_text(hjust = 0.5)) #標題居中

從上圖可以知道當配偶或兄弟姐妹的人數為1時，存活率很高。

（6）Parch的影響

ggplot(combi[(1:891),], aes(x = Parch, fill = Survived)) +geom_bar(stat = "count",position = "dodge") +labs(title = "The impacts of Parch",x = "Parch", y = "count") +theme(plot.title = element_text(hjust = 0.5)) #標題居中

從上圖可以知道當父母或子女的人數為1-2時，存活率相對來說較高。

（7）FamilySize的影響

由SibSp和Parch聯想到家庭規模對於生存率的影響，當乘客無親人和親人過多時，生存率都比較低，而當乘客有少數親人時，生存率明顯升高，故考慮引入新的變數FamilySize（家庭總人數）：

combi$FamilySize <- combi$Parch + combi$SibSp + 1ggplot(combi[(1:891),], aes(x = FamilySize, fill = Survived)) +geom_bar(stat = "count",position = "dodge") +labs(title = "The impacts of FamilySize",x = "FamilySize", y = "count") +theme(plot.title = element_text(hjust = 0.5)) #標題居中

從上圖可以知道，2-5口人的家庭生存情況是最好的，過大或過小家庭規模生存率都不好。

（8）Fare的影響

ggplot(combi[(1:891),], aes(x = Fare, color = Survived)) +geom_line(aes(label = ..count..),stat = "bin",binwidth = 5,na.rm = TRUE) +labs(title = "The impacts of Fare",x = "Fare", y = "count") +theme(plot.title = element_text(hjust = 0.5)) #標題居中

從上圖可以知道，票價越高，越容易存活。

（9）Embarked的影響

ggplot(combi[(1:891),], aes(x = Embarked, fill = Survived)) +geom_bar(stat = "count",position = "dodge") +labs(title = "The impacts of Embarked",x = "Embarked", y = "count") +theme(plot.title = element_text(hjust = 0.5)) #標題居中

從上圖可以知道，Embarked=C的存活率最高，Embarked=S的死亡率非常高。

6. 利用決策樹對Age的缺失值進行預測

library(rpart) Age.model <- rpart(Age~Pclass+Sex+SibSp+Parch+Fare+Embarked+Title+FamilySize,data = combi[!is.na(combi$Age),],method = "anova")combi$Age[is.na(combi$Age)] <- predict(Age.model,combi[is.na(combi$Age),])

將連續的年齡分為以下幾個等級：

combi$Age2 <- "15-"combi$Age2[combi$Age >= 15 & combi$Age < 30 ] <- "15-30"combi$Age2[combi$Age >= 30 & combi$Age < 45 ] <- "30-45"combi$Age2[combi$Age >= 45 & combi$Age < 60 ] <- "45-60"combi$Age2[combi$Age >= 60 ] <- "60+"

7.將訓練集與測試集分開：

new_train <- combi[1:891,]new_test <- combi[892:1309,]

三、樸素貝葉斯方法簡介

貝葉斯分類是一種分類演算法的總稱，以貝葉斯定理作為基礎，樸素貝葉斯分類是貝葉斯分類中最簡單的一種。

（1）分類問題：

已知類別集合C={y1,y2...yn}，其中的每一個元素都是一個類別；已知項集合I={x1,x2...xm...}，其中的每一個元素都是一個待分類項，每個元素含有某些特徵屬性。分類的目的就是構造一個分類器f，使得I中的每一個元素xi都能在C中找到一個唯一的yi滿足yi=f(xi)。

分類問題往往採用經驗性方法構造映射規則，即一般情況下的分類問題缺少足夠的信息來構造100%正確的映射規則，而是通過對經驗數據的學習從而實現一定概率意義上正確的分類，因此所訓練出的分類器並不是一定能將每個待分類項準確映射到其分類，分類器的質量與分類器構造方法、待分類數據的特性以及訓練樣本數量等諸多因素有關。

（2）貝葉斯定理

p(A|B)表示事件B已經發生的前提下，事件A發生的概率，叫做事件B發生下事件A的條件概率其基本求解公式為：p(A|B)=p(AB)/p(B)。

（3）樸素貝葉斯分類原理

樸素貝葉斯的思想：對於給出的待分類項，求解在此項出現的條件下各個類別出現的概率，哪個最大，就認為此待分類項屬於哪個類別。

樸素貝葉斯分類的正式定義如下：

設有一個待分類項x={a1,a2...am}，而每個a為x的一個特徵屬性
有類別集合C={y1,y2...yn}
計算p(y1|x)，p(y2|x)...p(yn|x)
如果p(yi|x)=max{p(y1|x)，p(y2|x)...p(yn|x)}，則x屬於yi

由上面的定義我們可以知道，只要計算出第3步中的各個條件概率，就可以判斷x屬於哪一類。計算公式為：p(yi|x)=p(x|yi)p(yi)/p(x),所有的p(x)是相同的，只需要比較p(x|yi)p(yi)的大小即可。

四、運用樸素貝葉斯方法進行預測

library(klaR) bfit <- NaiveBayes(Survived ~ Pclass + Sex + Age2 + SibSp + Parch + Embarked + Fare2 + Title + FamilySize, new_train ,na.action = na.pass)prediction <- predict(bfit,new_test)result <- data.frame(new_test$PassengerId,prediction[1])names(result) <- c("PassengerId","Survived")write.csv(result,file = "Naive_Bayesian.csv",row.names = FALSE)

最終預測模型的準確率如下：