--- title: "Titanic Logistic Regression" author: "roosh" date: "2/24/2017" output: html_document --- ```{r setup, include=FALSE} knitr::opts_chunk$set(echo = TRUE) ``` ### Introduction This is my first submission to any Kaggle competition. Any comments/recommendations are welcome. I chose to go with the Titanic dataset as that seems to be the fundamental starting point for anyone new to Kaggle. The competition can be found on the [Kaggle Page] (https://www.kaggle.com/c/titanic/data) There are four major parts to my script: * Feature engineering * Missing value imputation * Exploratory analysis and visualizations * Prediction ### Load and check data ```{r, message = FALSE} # Load packages library(ggplot2) # visualization library(ggthemes) # visualization library(dplyr) # data manipulation library(rpart) # prediction library(rpart.plot) ``` I load the train and test files. Additionally, I combine them so I can perform feature engineering across the entire dataset. Since the test dataset has no "Survived" column I will create an empty column for now so I can merge these datasets together. ```{r, message = FALSE} train <- read.csv("../input/train.csv") test <- read.csv("../input/test.csv") test$Survived <- NA combo <- rbind(train, test) # bind training & test data ``` Lets take a look at the combined dataset ```{r} str(combo) ``` Combined, we have 1309 observations across 12 variables. Let's take a look at these visually. ```{r} summary(combo) ``` A brief description of what these variables mean: Variable Name | Description --------------|------------- Survived | Survived (1) or died (0) Pclass | Passenger's class Name | Passenger's name Sex | Passenger's sex Age | Passenger's age SibSp | Number of siblings/spouses aboard Parch | Number of parents/children aboard Ticket | Ticket number Fare | Fare Cabin | Cabin Embarked | Port of embarkation Right off the bat, we notice a significant number of NAs across the 'Age' variable. 1 NA in the 'Fare' variable, and looks like 2 passengers are missing the port of embarkment. This moves us to our next step where we derive new features that may lend predictive power to our models and then, try to fill in these blanks. ### Feature Engineering Sifting through the 'Names' column, I notice that some people have unusual titles, such as Major, Lady, Duchess, Rev., etc. My curiousity is piqued, and now I'm curious to see if having an unusual title leads to a greater probability of survival. To start with, I will create a new variable called 'Title' which I will extract from the existing 'Names' column. First, lets convert the Name type into a character string, and do a quick trial run to see how the strsplit function works in R. ```{r} combo$Name<- as.character(combo$Name) strsplit(combo$Name[1],split = '[,.]')[[1]][2] ``` Voila! We see that the title is the second word in the string. We now use the sapply function to apply this transformation to all names in our combo data frame. ```{r} combo$Title<- sapply(strsplit(combo$Name, split='[,.]'), function(x) (x[2])) combo$Title <- sub(' ','',combo$Title) table(combo$Title) ``` Using this information, I will lump the titles with few counts into one called "Unique Title". The Mlle, Mme and Ms will be combined under the 'Miss' title. ```{r} combo$Title[combo$Title %in% c('Mme','Mlle')] <- 'Miss' unique_title <- c('Dona', 'Lady', 'the Countess','Capt', 'Col', 'Don', 'Dr', 'Major', 'Rev', 'Sir', 'Jonkheer') combo$Title[combo$Title %in% unique_title] <- 'Unique Title' ``` Next, I get a breakdown of title by Sex as a sanity check ```{r} table(combo$Title, combo$Sex) ``` Next, I create a new variable called "Fsize" to measure the sizes of different families. I do so by adding the existing variables SibSp (no. of siblings) + Parch (parents or children) + 1 (to include the individual itself) ```{r} combo$Fsize<- combo$SibSp + combo$Parch + 1 table(combo$Fsize) ``` Looks like most people were travelling alone. We can now to convert this numeric variable into a categorical one using the case_when command in dplyr, this will be useful for visualizations. ```{r} combo$Famcat <- case_when( combo$Fsize >= 0 & combo$Fsize <2 ~ "Single", combo$Fsize >= 2 & combo$Fsize <4 ~ "Small Family", combo$Fsize >= 4 ~ "Large Family" ) combo$Famcat<- as.factor(combo$Famcat) table(combo$Famcat) ``` Lets also create another variable for Age category: ```{r} combo$Agecat <- case_when( combo$Age >= 0 & combo$Age <18 ~ "Child", combo$Age >= 18 & combo$Age <30 ~ "Under 30", combo$Age >= 30 & combo$Age <60 ~ "30-60", combo$Age >= 60~ "Senior" ) ``` Now let's rearrange these so they visualize correctly later: ```{r} combo$Agecat<- factor(combo$Agecat, levels=c("Child", "Under 30", "30-60", "Senior")) ``` ### Missing value imputation Let's use a simple function to see all variables missing values in the 'Combo' dataset ```{r} is.na(combo) <- combo == "" sapply(combo, function(x) sum(is.na(x))) ``` Hmmm, since there is only 1 missing Fare, let us start by fixing this first. Let's pull out the empty row using dplyr ```{r} combo %>% filter(is.na(Fare)) %>% head ``` We know this passenger departed from from Southamption (S) and was in class 3. We compute a simple mean of this group and substitute in this value for passenger 1044. ```{r} combo %>% group_by(Embarked, Pclass) %>% na.omit() %>% summarize(mean_fare = mean(Fare, na.rm = TRUE), n=n()) ``` Next, I just plug in the $11.02 for passenger 1044. ```{r} combo$Fare[1044]<- 11.02 combo[1044,] ``` Next, lets look at the 2 entries missing Embarkment ports ```{r} combo %>% filter(is.na(Embarked)) %>% head ``` The 2 interesting pieces of information here are that their fares are the same ($80) and they both travelled first class. Let's see the frequencies of all passengers by class ```{r} k <- ggplot(na.omit(combo), aes(x=factor(Embarked), fill=factor(Pclass))) k + geom_bar( position="dodge") ``` There are about 30 more S vs C passengers. But using the mean Fares derived earlier, I notice that the C mean fare is 76, closer to the fare of passengers 62 and 830 vs 106 for S passengers. I will go ahead and assign these 2 passengers to embarkment port 'C'. ```{r} combo$Embarked[c(62,830)]<- "C" table(combo$Embarked) ``` Next, I notice the number of missing values in 'Cabin'. Since 77% (1014/1309) of entries are missing these values, for now I will just ignore this variable and exclude it from any modelling. Finally, I focus on the "Age" variables. Since only 20% of passengers having missing ages, I will build a simple regression model to replace these values. ```{r} Agefit <- rpart(Age ~ Pclass + Sex + SibSp + Parch + Fare + Embarked + Title + Fsize, data = combo[!is.na(combo$Age),], method = 'anova') combo$Age[is.na(combo$Age)] <- predict(Agefit, combo[is.na(combo$Age),]) ``` Now that we have completed inputting missing values, we can start some exploratory analysis to start gleaning insights from the data. ### Exploratory analysis and visualizations First, I'm curious to see the distribution of survivors across Sexes ```{r} ggplot(na.omit(combo), aes(x=Sex, fill=factor(Survived))) +geom_bar(position="dodge") ``` Far more women survived than men! Was this different across Passenger classes? ```{r} ggplot(na.omit(combo), aes(x=Sex, fill=factor(Survived))) +geom_bar(position="dodge") + facet_grid(".~Pclass") ``` There were far more survivors in first class, its increasing looking like travelling as a first class female gave us the highest priority of survival. Next, lets do a breakdown by Age ```{r} m <- ggplot(na.omit(combo), aes(x = Agecat, fill=factor(Survived))) m + geom_bar(position = "stack") ``` Almost all children were rescued, and those Under 30 were also relatively lucky. Those over 30 have a smaller chance of survival. Finally, I combine these variables and compute a simple "Survived Ratio" using dplyr ```{r} c1<- combo %>% group_by (Pclass, Agecat) %>% summarise(survived_ratio=mean(Survived,na.rm=TRUE)) %>% arrange(desc(survived_ratio)) c1 ``` Survival rates vary greatly by class and Age category. For example, 91% of children in 1st and 2nd class survived, but only 36% of those in 3rd class made it. Survival rates of senior citizens are less than 30% across the board. ### Prediction First, we break out the "combo" data set where we have performed our feature engineering back into the training data and testing data set. ```{r} train1<- combo[1:891,] test1<- combo[892:1309,] ``` Then, I use a logistic regression model to start with the prediction. ```{r} train1$Pclass<- as.factor(train1$Pclass) log3<- glm(data=train1, Survived ~ Pclass + Sex +Fare + Embarked + Title + Fsize +Agecat + Sex*Pclass, family="binomial") summary(log3) ``` I then apply this model to the test data set and use a 50% cutoff. (ie. any survival score > 50% will be rounded to 1) ```{r} test1$Pclass<- as.factor(test1$Pclass) pred1 <- predict(log3, newdata = test1, type = "response") survival1 = ifelse(pred1 > 0.5, 1, 0) ``` Done! Now I write this file to a csv and hit submit! Not bad, top 15% with this simple model. I hope this was helpful! Hope to keep these coming in the next few months.