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_train

Read 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" ...
ColumnMeaning
PassengerIdsequential number
Survived0 = died, 1 = survived
Pclassclass: 1, 2 or 3
Namelast name, title and first name
Sexfemale or male
Ageage in years
SibSpnumber of siblings and spouses on board
Parchnumber of parents and children on board
Ticketticket number
Farefare in British pounds
Cabincabin number
Embarkedport 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           2

The 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.4

342 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")

Survivors and deaths

Sex

Passengers by sex

table(train$Sex)
#> 
#> female   male 
#>    314    577
barplot(table(train$Sex),
        names.arg = c("Women", "Men"),
        col = c("#e377c2", "#1f77b4"),
        ylab = "Number of passengers", ylim = c(0, 600),
        main = "Passengers by sex")

Passengers by sex

Survival by sex

table(train$Sex, train$Status)
#>         
#>          died survived
#>   female   81      233
#>   male    468      109
barplot(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")

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.2

With 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.9

74 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 119
tab <- 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")

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.2

In 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       47

In 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.1354467

Because 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")

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     177

The 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       9

aggregate() 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 29

In 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 19

Survival 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")

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            1

ifelse() 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     27

Survival 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.4

Title 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        0

Age 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  70

The 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.9607079

The 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.02404709

This 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.008

According 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.