When the Titanic sank in April 1912, more than half of the people on board died. The Kaggle platform provides a table for 891 passengers that shows who survived. This post analyses the table with R: how are sex, class and age related to survival? Lines starting with #> show the output of R.
Loading the data
Download the dataset
The data belongs to the Kaggle competition Titanic: Machine Learning from Disaster. Downloading it requires a free Kaggle account. The file train.csv contains the 891 passengers whose fate is known.
Without an account, use the R package titanic. It contains the same table under the name titanic_train:
install.packages("titanic")
train <- titanic::titanic_trainRead the dataset
Some passengers have no age, cabin or port in the file. na.strings makes R read empty fields as the missing value NA instead of as empty text.
train <- read.csv("train.csv", na.strings = c("", "NA"))Structure of the data
dim(train) # rows and columns
#> [1] 891 12
str(train)
#> 'data.frame': 891 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" ...| Column | Meaning |
|---|---|
| PassengerId | sequential number |
| Survived | 0 = died, 1 = survived |
| Pclass | class: 1, 2 or 3 |
| Name | last name, title and first name |
| Sex | female or male |
| Age | age in years |
| SibSp | number of siblings and spouses on board |
| Parch | number of parents and children on board |
| Ticket | ticket number |
| Fare | fare in British pounds |
| Cabin | cabin number |
| Embarked | port of embarkation: C = Cherbourg, Q = Queenstown, S = Southampton |
head(train[, c("Survived", "Pclass", "Name", "Sex", "Age", "Fare")])
#> Survived Pclass Name
#> 1 0 3 Braund, Mr. Owen Harris
#> 2 1 1 Cumings, Mrs. John Bradley (Florence Briggs Thayer)
#> 3 1 3 Heikkinen, Miss. Laina
#> 4 1 1 Futrelle, Mrs. Jacques Heath (Lily May Peel)
#> 5 0 3 Allen, Mr. William Henry
#> 6 0 3 Moran, Mr. James
#> Sex Age Fare
#> 1 male 22 7.2500
#> 2 female 38 71.2833
#> 3 female 26 7.9250
#> 4 female 35 53.1000
#> 5 male 35 8.0500
#> 6 male NA 8.4583
colSums(is.na(train)) # missing values per column
#> PassengerId Survived Pclass Name Sex
#> 0 0 0 0 0
#> Age SibSp Parch Ticket Fare
#> 177 0 0 0 0
#> Cabin Embarked
#> 687 2The age is missing for 177 passengers, the cabin for 687. All analyses of age therefore only cover the 714 passengers with a known age.
To make the tables readable, the column Survived gets labels:
train$Status <- factor(train$Survived, levels = c(0, 1),
labels = c("died", "survived"))Survivors overall
table(train$Status)
#>
#> died survived
#> 549 342
round(prop.table(table(train$Status)) * 100, 1) # in percent
#>
#> died survived
#> 61.6 38.4342 passengers survived, that is 38 percent. 549 passengers died.
count <- table(train$Status)
pie(count,
labels = paste0(names(count), ": ", count, " people"),
col = c("#c0392b", "#27ae60"),
main = "Passengers who survived and died")
Sex
Passengers by sex
table(train$Sex)
#>
#> female male
#> 314 577barplot(table(train$Sex),
names.arg = c("Women", "Men"),
col = c("#e377c2", "#1f77b4"),
ylab = "Number of passengers", ylim = c(0, 600),
main = "Passengers by sex")
Survival by sex
table(train$Sex, train$Status)
#>
#> died survived
#> female 81 233
#> male 468 109barplot(table(train$Sex, train$Status),
beside = TRUE,
col = c("#e377c2", "#1f77b4"),
ylab = "Number of passengers", ylim = c(0, 500),
legend.text = c("Women", "Men"),
args.legend = list(x = "topright"),
main = "Survivors and deaths by sex")
prop.table() converts the table into shares. Without further arguments, the shares refer to all 891 passengers:
round(prop.table(table(train$Sex, train$Status)) * 100, 1)
#>
#> died survived
#> female 9.1 26.2
#> male 52.5 12.2With margin = 1, the shares refer to each row, here separately to all women and to all men:
round(prop.table(table(train$Sex, train$Status), margin = 1) * 100, 1)
#>
#> died survived
#> female 25.8 74.2
#> male 81.1 18.974 percent of the women survived, but only 19 percent of the men. This matches the rule “women and children first” when the lifeboats were filled.
Class
Passengers by class
table(train$Pclass)
#>
#> 1 2 3
#> 216 184 491
table(train$Status, train$Pclass)
#>
#> 1 2 3
#> died 80 97 372
#> survived 136 87 119tab <- table(train$Status, train$Pclass)
colnames(tab) <- paste("Class", colnames(tab))
dotchart(tab, pch = 19, col = "#1f77b4", gcolor = "#c0392b",
xlab = "Number of passengers", xlim = c(0, 400),
main = "Survivors and deaths by class")
Survival rate per class
margin = 2 makes the shares refer to each column, here to each class. This replaces manual calculations such as 100 / 216 * 136.
round(prop.table(table(train$Status, train$Pclass), margin = 2) * 100, 1)
#>
#> 1 2 3
#> died 37.0 52.7 75.8
#> survived 63.0 47.3 24.2In first class 63 percent survived, in second class 47 percent and in third class 24 percent.
Sex and class together
ftable() shows a table with three features clearly in rows and columns:
ftable(table(train$Pclass, train$Sex, train$Status))
#> died survived
#>
#> 1 female 3 91
#> male 77 45
#> 2 female 6 70
#> male 91 17
#> 3 female 72 72
#> male 300 47In first and second class almost all women survived. In third class it was only half.
aggregate(Survived ~ Pclass + Sex, data = train, FUN = mean)
#> Pclass Sex Survived
#> 1 1 female 0.9680851
#> 2 2 female 0.9210526
#> 3 3 female 0.5000000
#> 4 1 male 0.3688525
#> 5 2 male 0.1574074
#> 6 3 male 0.1354467Because Survived only contains 0 and 1, the mean is the survival rate.
Age
Age distribution
library(ggplot2)
ggplot(train, aes(x = Age)) +
geom_histogram(binwidth = 5, boundary = 0, fill = "#1f77b4", colour = "white") +
labs(x = "Age in years", y = "Number of passengers", title = "Age of the passengers")
summary(train$Age)
#> Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
#> 0.42 20.12 28.00 29.70 38.00 80.00 177The average age is just under 30 years. The youngest passenger was 5 months old, the oldest 80 years.
Age of the surviving women per class
women <- subset(train, Sex == "female" & Status == "survived")
aggregate(Age ~ Pclass, data = women, FUN = median)
#> Pclass Age
#> 1 1 35
#> 2 2 28
#> 3 3 19
summary(women$Age[women$Pclass == 1])
#> Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
#> 14.00 23.25 35.00 34.94 44.00 63.00 9aggregate() with the formula Age ~ Pclass automatically leaves out passengers without an age.
Children
Here, children are all passengers under 18.
children <- subset(train, Age < 18)
nrow(children)
#> [1] 113
table(children$Pclass) # children per class
#>
#> 1 2 3
#> 12 23 78
table(children$Status, children$Pclass) # survived and died per class
#>
#> 1 2 3
#> died 1 2 49
#> survived 11 21 29In first class one child died, in second class two, in third class 49.
girls <- subset(children, Sex == "female")
table(girls$Status, girls$Pclass)
#>
#> 1 2 3
#> died 1 0 16
#> survived 7 12 19Survival rate by age group
train$AgeGroup <- cut(train$Age, breaks = c(0, 10, 18, 30, 50, 80),
labels = c("0 to 9", "10 to 17", "18 to 29",
"30 to 49", "50 to 80"),
right = FALSE, include.lowest = TRUE)
rate <- aggregate(Survived ~ AgeGroup + Sex, data = train, FUN = mean)
rate
#> AgeGroup Sex Survived
#> 1 0 to 9 female 0.6333333
#> 2 10 to 17 female 0.7600000
#> 3 18 to 29 female 0.7282609
#> 4 30 to 49 female 0.7826087
#> 5 50 to 80 female 0.9090909
#> 6 0 to 9 male 0.5937500
#> 7 10 to 17 male 0.1538462
#> 8 18 to 29 male 0.1564246
#> 9 30 to 49 male 0.2134146
#> 10 50 to 80 male 0.1346154
ggplot(rate, aes(x = AgeGroup, y = Survived * 100, fill = Sex)) +
geom_col(position = "dodge") +
scale_fill_manual(values = c(female = "#e377c2", male = "#1f77b4"),
labels = c("Women", "Men"), name = NULL) +
labs(x = "Age in years", y = "Survival rate in percent",
title = "Survival rate by age group and sex")
Boys under 10 had a much better chance than adult men. From the age of 10, the survival rate of boys drops to the level of the men.
Title in the name
Every name contains a title such as Mr, Mrs, Miss or Master. Master was the title for boys. The title sits between the comma and the first dot and can be extracted with a regular expression:
anyDuplicated(train$Name) # 0 means: no name appears twice
#> [1] 0
train$Title <- sub(".*, *([^.]+)\\..*", "\\1", train$Name)
sort(table(train$Title), decreasing = TRUE)
#>
#> Mr Miss Mrs Master Dr
#> 517 182 125 40 7
#> Rev Col Major Mlle Capt
#> 6 2 2 2 1
#> Don Jonkheer Lady Mme Ms
#> 1 1 1 1 1
#> Sir the Countess
#> 1 1ifelse() groups rare titles as “Other”:
common <- c("Mr", "Mrs", "Miss", "Master")
train$Title <- ifelse(train$Title %in% common, train$Title, "Other")
table(train$Title)
#>
#> Master Miss Mr Mrs Other
#> 40 182 517 125 27Survival by title
table(train$Title, train$Status)
#>
#> died survived
#> Master 17 23
#> Miss 55 127
#> Mr 436 81
#> Mrs 26 99
#> Other 15 12
round(prop.table(table(train$Title, train$Status), margin = 1) * 100, 1)
#>
#> died survived
#> Master 42.5 57.5
#> Miss 30.2 69.8
#> Mr 84.3 15.7
#> Mrs 20.8 79.2
#> Other 55.6 44.4Title and class
ftable(table(train$Pclass, train$Title, train$Status))
#> died survived
#>
#> 1 Master 0 3
#> Miss 2 44
#> Mr 70 37
#> Mrs 1 41
#> Other 7 11
#> 2 Master 0 9
#> Miss 2 32
#> Mr 83 8
#> Mrs 4 37
#> Other 8 1
#> 3 Master 17 11
#> Miss 51 51
#> Mr 283 36
#> Mrs 21 21
#> Other 0 0Age per title
aggregate(Age ~ Title, data = train, FUN = max)
#> Title Age
#> 1 Master 12
#> 2 Miss 63
#> 3 Mr 80
#> 4 Mrs 63
#> 5 Other 70The oldest boy with the title Master was 12 years old.
Information gain
Information gain shows how much a feature reduces the uncertainty about the outcome, here about survival. Decision trees use it to choose the feature by which they split the data first. It is based on entropy: entropy is 0 when all values are the same, and 1 when two values are equally frequent.
entropy <- function(y) {
p <- table(y) / length(y)
p <- p[p > 0]
-sum(p * log2(p))
}
information_gain <- function(y, x) {
x <- addNA(factor(x), ifany = TRUE) # missing values as a separate group
parts <- split(y, x)
conditional <- sum(sapply(parts, \(p) length(p) / length(y) * entropy(p)))
entropy(y) - conditional
}
entropy(train$Survived)
#> [1] 0.9607079The entropy of survival is 0.96, because survivors and deaths are almost equally frequent.
Features with many different values
If you apply information gain to all columns, Name, Ticket and PassengerId end up at the top:
columns <- c("PassengerId", "Name", "Ticket", "Cabin", "Sex", "Pclass", "Embarked")
sort(sapply(train[columns], \(x) information_gain(train$Survived, x)),
decreasing = TRUE)
#> PassengerId Name Ticket Cabin Sex
#> 0.96070790 0.96070790 0.82593692 0.23838554 0.21766011
#> Pclass Embarked
#> 0.08383105 0.02404709This result is misleading. Every name appears only once, so the name predicts the survival of each individual passenger perfectly, without telling anything about a new passenger. Information gain favours features with many different values.
Information gain ratio
The information gain ratio divides the information gain by the entropy of the feature itself. Features with very many different values have a high entropy and therefore get a smaller value:
gain_ratio <- function(y, x) {
x <- addNA(factor(x), ifany = TRUE)
information_gain(y, x) / entropy(x)
}
train$CabinKnown <- !is.na(train$Cabin)
train$FareGroup <- cut(train$Fare, breaks = c(0, 10, 30, 100, Inf),
right = FALSE)
features <- c("PassengerId", "Name", "Ticket", "Sex", "Title", "Pclass",
"CabinKnown", "FareGroup", "AgeGroup", "Embarked",
"SibSp", "Parch")
result <- data.frame(
feature = features,
gain = sapply(train[features], \(x) information_gain(train$Survived, x)),
ratio = sapply(train[features], \(x) gain_ratio(train$Survived, x)),
row.names = NULL
)
result[order(-result$ratio), ] |> transform(gain = round(gain, 3),
ratio = round(ratio, 3))
#> feature gain ratio
#> 4 Sex 0.218 0.232
#> 5 Title 0.239 0.143
#> 1 PassengerId 0.961 0.098
#> 2 Name 0.961 0.098
#> 7 CabinKnown 0.071 0.092
#> 3 Ticket 0.826 0.089
#> 6 Pclass 0.084 0.058
#> 8 FareGroup 0.083 0.047
#> 11 SibSp 0.033 0.025
#> 10 Embarked 0.024 0.022
#> 12 Parch 0.024 0.021
#> 9 AgeGroup 0.019 0.008According to the gain ratio, sex comes first, followed by the title. Name, Ticket and PassengerId lose most of their value, but still rank above the class. Columns that have a different value for every passenger therefore have to be removed from the table before any analysis. Neither of the two measures fully compensates for this error.
The package FSelector calculates both values with information.gain() and gain.ratio(), but needs a Java installation.
