Kaggle Competition | Titanic Machine Learning from Disaster

The sinking of the RMS Titanic is one of the most infamous shipwrecks in history. On April 15, 1912, during her maiden voyage, the Titanic sank after colliding with an iceberg, killing 1502 out of 2224 passengers and crew. This sensational tragedy shocked the international community and led to better safety regulations for ships.

One of the reasons that the shipwreck led to such loss of life was that there were not enough lifeboats for the passengers and crew. Although there was some element of luck involved in surviving the sinking, some groups of people were more likely to survive than others, such as women, children, and the upper-class.

In this contest, we ask you to complete the analysis of what sorts of people were likely to survive. In particular, we ask you to apply the tools of machine learning to predict which passengers survived the tragedy.

This Kaggle Getting Started Competition provides an ideal starting place for people who may not have a lot of experience in data science and machine learning."

From the competition homepage.

Goal for this Notebook:

Show a simple example of an analysis of the Titanic disaster in Python using a full complement of PyData utilities. This is aimed for those looking to get into the field or those who are already in the field and looking to see an example of an analysis done with Python.

This Notebook will show basic examples of:

Data Handling

  • Importing Data with Pandas
  • Cleaning Data
  • Exploring Data through Visualizations with Matplotlib

Data Analysis

  • Supervised Machine learning Techniques:
    • Logit Regression Model
    • Plotting results
  • Unsupervised Machine learning Techniques
    • Support Vector Machine (SVM) using 3 kernels
    • Basic Random Forest
    • Plotting results

Valuation of the Analysis

  • K-folds cross validation to valuate results locally
  • Output the results from the IPython Notebook to Kaggle

Required Libraries:

To run this notebook interactively, get it from my Github here. The competition's website is located on Kaggle.com.


In [2]:
import pandas as pd
import numpy as np
import statsmodels.api as sm

from statsmodels.nonparametric import KDE, smoothers_lowess
from pandas import Series,DataFrame
from patsy import dmatrices
from sklearn import datasets, svm

In [2]:
# kaggleaux is a module that contains auxiliary functions I made for kaggle competitions. 
# the module is avaliable at github.com/agconti/AGC_KaggleAux
import sys
sys.path.append('C:/Users/Andrew/Documents/Kaggle/Kaggle Aux')
import kaggleaux as ka

Data Handling

Let's read our data in using pandas:


In [3]:
df = pd.read_csv("train.csv")

In [4]:
df # Show an overview of our data


Out[4]:
&ltclass 'pandas.core.frame.DataFrame'>
Int64Index: 891 entries, 0 to 890
Data columns (total 11 columns):
survived    891  non-null values
pclass      891  non-null values
name        891  non-null values
sex         891  non-null values
age         714  non-null values
sibsp       891  non-null values
parch       891  non-null values
ticket      891  non-null values
fare        891  non-null values
cabin       204  non-null values
embarked    889  non-null values
dtypes: float64(2), int64(4), object(5)

Let's take a look:

Above is a summary of our data contained in a Pandas DataFrame. Think of a DataFrame as a Python's super charged version of the workflow in an Excel table. As you can see the summary holds quite a bit of information. First, it lets us know we have 891 observations, or passengers, to analyze here:

Int64Index: 891 entries, 0 to 890

Next it shows us all of the columns in DataFrame. Each column tells us something about each of our observations, like their name, sex or age. These colunms are called a features of our dataset. You can think of the meaning of the words column and feature as interchangeable for this notebook.

After each feature it lets us know how many values it contains. While most of our features have complete data on every observation, like the survived feature here:

survived    891  non-null values 

some are missing information, like the age feature:

age         714  non-null values 

These missing values are represented as NaNs.

Take care of missing values:

The features ticket and cabin have many missing values and so can’t add much value to our analysis. To handle this we will drop them from the dataframe to preserve the integrity of our dataset.

To do that we'll use this line of code to drop the features entirely:

df = df.drop(['ticket','cabin'], axis=1) 


While this line of code removes the NaN values from every remaining column / feature:

df = df.dropna()

Now we have a clean and tidy dataset that is ready for analysis. Because .dropna() removes an observation from our data even if it only has 1 NaN in one of the features, it would have removed most of our dataset if we had not dropped the ticket and cabin features first.


In [5]:
df = df.drop(['ticket','cabin'], axis=1) 

df = df.dropna() # Remove NaN values

For a detailed look at how to use pandas for data analysis, the best resource is Wes Mckinney's book. Additional interactive tutorials that cover all of the basics can be found here (they're free). If you still need to be convinced about the power of pandas check out this wirlwhind look at all that pandas can do.

Let's take a Look at our data graphically:


In [6]:
fig = plt.figure(figsize(18,6), dpi=1600) # specifies the parameters of our graphs
a = 0.2                                   # sets the alpha level of the colors in the graph (for more attractive results)
a_bar = 0.55                              # another alpha setting

plt.subplot2grid((2,3),(0,0))             # lets us plot many diffrent shaped graphs together
df.survived.value_counts().plot(kind='bar', alpha=a_bar)# plots a bar graph of those who surived vs those who did not. 
title("Distribution of Survival, (1 = Survived)") # puts a title on our graph

plt.subplot2grid((2,3),(0,1))
plt.scatter(df.survived, df.age, alpha=a)
plt.ylabel("Age")                         # sets the y axis lable
plt.grid(b=True, which='major', axis='y') # formats the grid line style of our graphs
title("Survial by Age,  (1 = Survived)")

plt.subplot2grid((2,3),(0,2))
df.pclass.value_counts().plot(kind="barh", alpha=a_bar)
title("Class Distribution")

plt.subplot2grid((2,3),(1,0), colspan=2)
df.age[df.pclass == 1].plot(kind='kde')   # plots a kernel desnsity estimate of the subset of the 1st class passanges's age
df.age[df.pclass == 2].plot(kind='kde')
df.age[df.pclass == 3].plot(kind='kde')
plt.xlabel("Age")                         # plots an axis lable
title("Age Distribution within classes"); legend(('1st Class', '2nd Class','3rd Class'),loc='best') # sets our legend for our graph.

plt.subplot2grid((2,3),(1,2))
df.embarked.value_counts().plot(kind='bar', alpha=a_bar)
title("Passengers per boarding location");


Exploratory Visualization:

The point of this competition is to predict if an individual will survive based on the features in the data like:

  • Traveling Class (called pclass in the data)
  • Sex
  • Age
  • Fare Price

Let’s see if we can gain a better understanding of who survived and died.

First let’s plot a bar graph of those who Survived Vs. Those who did not.


In [7]:
plt.figure(figsize(6,4))
df.survived.value_counts().plot(kind='barh', color="blue", alpha=.65)
title("Survival Breakdown (1 = Survived, 0 = Died)")


Out[7]:
<matplotlib.text.Text at 0x12d00d68>

Now let’s tease more structure out of the data,

Let’s break the previous graph down by gender


In [9]:
fig = plt.figure(figsize=(18,6))

#create a plot of two subsets, male and female, of the survived variable.
#After we do that we call value_counts() so it can be easily plotted as a bar graph. 
#'barh' is just a horizontal bar graph
fig.add_subplot(121)
df.survived[df.sex == 'male'].value_counts().plot(kind='barh',label='Male')
df.survived[df.sex == 'female'].value_counts().plot(kind='barh', color='#FA2379',label='Female')
title("Who Survived? with respect to Gender, (raw value counts) "); legend(loc='best')

#adjust graph to display the proportions of survival by gender
fig.add_subplot(122)
(df.survived[df.sex == 'male'].value_counts()/float(df.sex[df.sex == 'male'].size)).plot(kind='barh',label='Male')  
(df.survived[df.sex == 'female'].value_counts()/float(df.sex[df.sex == 'female'].size)).plot(kind='barh', color='#FA2379',label='Female')
title("Who Survived proportionally? with respect to Gender"); legend(loc='best')


Out[9]:
<matplotlib.legend.Legend at 0x13e65668>

Here it’s clear that although more men died and survived in raw value counts, females had a greater survival rate proportionally(~25%), than men (~20%)

Great! But let’s go down even further:

Can we capture more of the structure by using Pclass? Here we will bucket classes as lowest class or any of the high classes (classes 1 - 2). 3 is lowest class. Let’s break it down by Gender and what Class they were traveling in.


In [10]:
fig=plt.figure(figsize=(18,4), dpi=1600)
a=.65 # our alpha or opacity level.

# building on the previous code, here we create an additional subset with in the gender subset we created for the survived variable. 
# I know, thats a lot of subsets. After we do that we call value_counts() so it it can be easily plotted as a bar graph. 
# this is repeated for each gender class pair.
ax1=fig.add_subplot(141)
df.survived[df.sex == 'female'][df.pclass != 3].value_counts().plot(kind='bar', label='female highclass', color='#FA2479', alpha=a)
ax1.set_xticklabels(["Survived", "Died"], rotation=0)
title("Who Survived? with respect to Gender and Class"); legend(loc='best')

ax2=fig.add_subplot(142, sharey=ax1)
df.survived[df.sex == 'female'][df.pclass == 3].value_counts().plot(kind='bar', label='female, low class', color='pink', alpha=a)
ax2.set_xticklabels(["Died","Survived"], rotation=0)
legend(loc='best')

ax3=fig.add_subplot(143, sharey=ax1)
df.survived[df.sex == 'male'][df.pclass == 3].value_counts().plot(kind='bar', label='male, low class',color='lightblue', alpha=a)
ax3.set_xticklabels(["Died","Survived"], rotation=0)
legend(loc='best')

ax4=fig.add_subplot(144, sharey=ax1)
df.survived[df.sex == 'male'][df.pclass != 3].value_counts().plot(kind='bar', label='male highclass', alpha=a, color='steelblue')
ax4.set_xticklabels(["Died","Survived"], rotation=0)
legend(loc='best')


Out[10]:
<matplotlib.legend.Legend at 0x12eaa240>

Awesome! Now we have a lot more information on who survived and died in the tragedy. With this deeper understanding, we are better equipped to create better more insightful models. This is a typical process in interactive data analysis. First you start small and understand the most basic relationships and slowly increment the complexity of your analysis as you discover more and more about the data you’re working with. Below is the progression of process laid out together:


In [16]:
fig=plt.figure(figsize(18,12), dpi=1600)
a=0.65
# Step 1
fig.add_subplot(341)
df.survived.value_counts().plot(kind='bar', color="blue", alpha=a)
title("Step. 1")

# Step 2
fig.add_subplot(345)
df.survived[df.sex == 'male'].value_counts().plot(kind='bar',label='Male')
df.survived[df.sex == 'female'].value_counts().plot(kind='bar', color='#FA2379',label='Female')
title("Step. 2 \nWho Survied? with respect to Gender."); legend(loc='best')

fig.add_subplot(346)
(df.survived[df.sex == 'male'].value_counts()/float(df.sex[df.sex == 'male'].size)).plot(kind='bar',label='Male')
(df.survived[df.sex == 'female'].value_counts()/float(df.sex[df.sex == 'female'].size)).plot(kind='bar', color='#FA2379',label='Female')
title("Who Survied proportionally?"); legend(loc='best')


#Step 3
ax1=fig.add_subplot(349)
df.survived[df.sex == 'female'][df.pclass != 3].value_counts().plot(kind='bar', label='female highclass', color='#FA2479', alpha=a)
ax1.set_xticklabels(["Survived", "Died"], rotation=0)
title("Step. 3"); legend(loc='best')


ax2=fig.add_subplot(3,4,10, sharey=ax1)
df.survived[df.sex=='female'][df.pclass==3].value_counts().plot(kind='bar', label='female, low class', color='pink', alpha=a)
ax2.set_xticklabels(["Died","Survived"], rotation=0)
legend(loc='best')

ax3=fig.add_subplot(3,4,11, sharey=ax1)
df.survived[df.sex == 'male'][df.pclass == 3].value_counts().plot(kind='bar', label='male, low class',color='lightblue', alpha=a)
ax3.set_xticklabels(["Died","Survived"], rotation=0)
legend(loc='best')

ax4=fig.add_subplot(3,4,12, sharey=ax1)
df.survived[df.sex == 'male'][df.pclass != 3].value_counts().plot(kind='bar', label='male highclass', alpha=a, color='steelblue')
ax4.set_xticklabels(["Died","Survived"], rotation=0)
legend(loc='best')


Out[16]:
<matplotlib.legend.Legend at 0x175f6080>

I've done my best to make the plotting code readable and intuitive, but if you’re looking for a more detailed look on how to start plotting in matplotlib, check out this beautiful notebook here.

Now that we have a basic understanding of what we are trying to predict, let’s predict it.

Supervised Machine Learning

Logistic Regression:

As explained by Wikipedia:

In statistics, logistic regression or logit regression is a type of regression analysis used for predicting the outcome of a categorical dependent variable (a dependent variable that can take on a limited number of values, whose magnitudes are not meaningful but whose ordering of magnitudes may or may not be meaningful) based on one or more predictor variables. That is, it is used in estimating empirical values of the parameters in a qualitative response model. The probabilities describing the possible outcomes of a single trial are modeled, as a function of the explanatory (predictor) variables, using a logistic function. Frequently (and subsequently in this article) "logistic regression" is used to refer specifically to the problem in which the dependent variable is binary—that is, the number of available categories is two—and problems with more than two categories are referred to as multinomial logistic regression or, if the multiple categories are ordered, as ordered logistic regression. Logistic regression measures the relationship between a categorical dependent variable and one or more independent variables, which are usually (but not necessarily) continuous, by using probability scores as the predicted values of the dependent variable.[1] As such it treats the same set of problems as does probit regression using similar techniques.

The skinny, as explained by yours truly:

Our competition wants us to predict a binary outcome. That is, it wants to know whether some will die, (represented as a 0), or survive, (represented as 1). A good place to start is to calculate the probability that an individual observation, or person, is likely to be a 0 or 1. That way we would know the chance that someone survives, and could start making somewhat informed perdictions. If we did, we'd get results like this::

(Y axis is the probability that someone survives, X axis is the passenger’s number from 1 to 891.)

While that information is useful it doesn’t let us know whether someone ended up alive or dead. It just lets us know the chance that they will survive or die. We still need to translate these probabilities into the binary decision we’re looking for. But how? We could arbitrarily say that our survival cutoff is anyone with a probability of survival over 50%. In fact, this tactic would actually perform pretty well for our data and would allow you to make decently accurate predictions. Graphically it would look something like this:

If you’re a betting man like me, you don’t like to leave everything to chance. What are the odds that setting that cutoff at 50% works? Maybe 20% or 80% would work better. Clearly we need a more exact way to make that cutoff. What can save the day? In steps the Logistic Regression.

A logistic regression follows the all steps we took above but mathematically calculates the cutoff, or decision boundary (as stats nerds call it), for you. This way it can figure out the best cut off to choose, perhaps 50% or 51.84%, that most accurately represents the training data.

The three cells below show the process of creating our Logitist regression model, training it on the data, and examining its performance.

First, we define our formula for our Logit regression. In the next cell we create a regression friendly dataframe that sets up boolean values for the categorical variables in our formula and lets our regression model know the types of inputs we're giving it. The model is then instantiated and fitted before a summary of the model's performance is printed. In the last cell we graphically compare the predictions of our model to the actual values we are trying to predict, as well as the residual errors from our model to check for any structure we may have missed.


In [17]:
# model formula
formula = 'survived ~ C(pclass) + C(sex) + age + sibsp  + C(embarked)' # here the ~ sign is an = sign, and the features of our dataset
                                                                       # are written as a formula to predict survived. The C() lets our 
                                                                       # regression know that those variables are categorical.
results = {} # create a results dictionary to hold our regression results for easy analysis later

In [18]:
# create a regression freindly dataframe using patsy's dmatrices function
y,x = dmatrices(formula, data=df, return_type='dataframe')

# instantiate our model
model = sm.Logit(y,x)

# fit our model to the training data
res = model.fit()

# save the result for outputing predictions later
results['Logit'] = [res, formula]
res.summary()


Optimization terminated successfully.
         Current function value: 316.404460
         Iterations 6
Out[18]:
Logit Regression Results
Dep. Variable: survived No. Observations: 712
Model: Logit Df Residuals: 704
Method: MLE Df Model: 7
Date: Mon, 01 Jul 2013 Pseudo R-squ.: 0.3414
Time: 14:30:56 Log-Likelihood: -316.40
converged: True LL-Null: -480.45
LLR p-value: 5.992e-67
coef std err z P>|z| [95.0% Conf. Int.]
Intercept 4.5423 0.474 9.583 0.000 3.613 5.471
C(pclass)[T.2] -1.2673 0.299 -4.245 0.000 -1.852 -0.682
C(pclass)[T.3] -2.4966 0.296 -8.422 0.000 -3.078 -1.916
C(sex)[T.male] -2.6239 0.218 -12.060 0.000 -3.050 -2.197
C(embarked)[T.Q] -0.8351 0.597 -1.398 0.162 -2.006 0.335
C(embarked)[T.S] -0.4254 0.271 -1.572 0.116 -0.956 0.105
age -0.0436 0.008 -5.264 0.000 -0.060 -0.027
sibsp -0.3697 0.123 -3.004 0.003 -0.611 -0.129

In [19]:
# Plot Predictions Vs Actual
plt.figure(figsize(18,4));
plt.subplot(121, axisbg="#DBDBDB")
# generate predictions from our fitted model
ypred = res.predict(x)
plt.plot(x.index, ypred, 'bo', x.index, y, 'mo', alpha=.25);
plt.grid(color='white', linestyle='dashed')
plt.title('Logit predictions, Blue: \nFitted/predicted values: Red');

# Residuals
plt.subplot(122, axisbg="#DBDBDB")
plt.plot(res.resid, 'r-')
plt.grid(color='white', linestyle='dashed')
plt.title('Logit Residuals');


So how well did this work?

Lets Look at the predictions we generated graphically:


In [20]:
fig = plt.figure(figsize(18,9), dpi=1600)
a = .2

# Below are examples of more advanced plotting. 
# It it looks strange check out the tutorial above.
fig.add_subplot(221, axisbg="#DBDBDB")
kde_res = KDE(res.predict())
kde_res.fit()
plt.plot(kde_res.support,kde_res.density)
plt.fill_between(kde_res.support,kde_res.density, alpha=a)
title("Distribution of our Predictions")

fig.add_subplot(222, axisbg="#DBDBDB")
plt.scatter(res.predict(),x['C(sex)[T.male]'] , alpha=a)
plt.grid(b=True, which='major', axis='x')
plt.xlabel("Predicted chance of survival")
plt.ylabel("Gender Bool")
title("The Change of Survival Probability by Gender (1 = Male)")

fig.add_subplot(223, axisbg="#DBDBDB")
plt.scatter(res.predict(),x['C(pclass)[T.3]'] , alpha=a)
plt.xlabel("Predicted chance of survival")
plt.ylabel("Class Bool")
plt.grid(b=True, which='major', axis='x')
title("The Change of Survival Probability by Lower Class (1 = 3rd Class)")

fig.add_subplot(224, axisbg="#DBDBDB")
plt.scatter(res.predict(),x.age , alpha=a)
plt.grid(True, linewidth=0.15)
title("The Change of Survival Probability by Age")
plt.xlabel("Predicted chance of survival")
plt.ylabel("Age")


Out[20]:
<matplotlib.text.Text at 0x17b3ea90>

Now Lets Use our Model to predict the test set values and then save the results so they can be outputed to Kaggle

Read the test data


In [21]:
test_data = pd.read_csv("test.csv")

Examine our dataframe


In [22]:
test_data


Out[22]:
&ltclass 'pandas.core.frame.DataFrame'>
Int64Index: 418 entries, 0 to 417
Data columns (total 10 columns):
pclass      418  non-null values
name        418  non-null values
sex         418  non-null values
age         332  non-null values
sibsp       418  non-null values
parch       418  non-null values
ticket      418  non-null values
fare        417  non-null values
cabin       91  non-null values
embarked    418  non-null values
dtypes: float64(2), int64(3), object(5)

Add our independent variable to our test data. (Its ussually left blank by Kaggle because its the value your trying to predict.)


In [23]:
test_data['survived'] = 1.23

In [24]:
results # our binned results data


Out[24]:
{'Logit': [<statsmodels.discrete.discrete_model.BinaryResultsWrapper at 0x13072a20>,
  'survived ~ C(pclass) + C(sex) + age + sibsp  + C(embarked)']}

In [25]:
compared_resuts = ka.predict(test_data, results, 'Logit') # Use your model to make prediction on our test set. 
compared_resuts = Series(compared_resuts)                 # convert our model to a series for easy output

In [26]:
# output and submit to kaggel
compared_resuts.to_csv("logitregres.csv")

Results as scored by Kaggle: RMSE = 0.77033 That result is pretty good. ECT ECT ECT

Unsupervised Machine Learning:


In [27]:
# Create an acceptable formula for our machine learning algorithms
formula_ml = 'survived ~ C(pclass) + C(sex) + age + sibsp + parch + C(embarked)'

Support Vector Machine (SVM)

"So uhhh, what if a straight line just doesn’t cut it."

Wikipeda:

In machine learning, support vector machines (SVMs, also support vector networks[1]) are supervised learning models with associated learning algorithms that analyze data and recognize patterns, used for classification and regression analysis. The basic SVM takes a set of input data and predicts, for each given input, which of two possible classes forms the output, making it a non-probabilistic binary linear classifier. Given a set of training examples, each marked as belonging to one of two categories, an SVM training algorithm builds a model that assigns new examples into one category or the other. An SVM model is a representation of the examples as points in space, mapped so that the examples of the separate categories are divided by a clear gap that is as wide as possible. New examples are then mapped into that same space and predicted to belong to a category based on which side of the gap they fall on. In addition to performing linear classification, SVMs can efficiently perform non-linear classification using what is called the kernel trick, implicitly mapping their inputs into high-dimensional feature spaces.

From me

The logit model we just implemented was great in that it showed exactly where to draw our decision boundary or our 'survival cut off'. But if you’re like me, you could have thought, "So uhhh, what if a straight line just doesn’t cut it". A linear line is okay, but can we do better? Perhaps a more complex decision boundary like a wave, circle, or maybe some sort of strange polygon would describe the variance observed in our sample better than a line. Imagine if we were predicating survival based on age. It could be a linear decision boundary, meaning each additional time you've gone around the sun you were 1 unit more or less likely to survive. But I think it could be easy to imagine some sort of curve, where a young healthy person would have the best chance of survival, and sadly the very old and very young a like: a poor chance. Now that’s a interesting question to answer. But our logit model can only evaluate a linear decision boundary. How do we get around this? With the usual answer to life the universe and everything; $MATH$.

The answer: We could transform our logit equation from expressing a linear relationship like so:

$survived = \beta_0 + \beta_1pclass + \beta_2sex + \beta_3age + \beta_4sibsp + \beta_5parch + \beta_6embarked$

Which we'll represent for convenience as: $y = x$

to a expressing a linear expression of a non-linear relationship: $\log(y) = \log(x)$

By doing this we're not breaking the rules. Logit models are only efficient at modeling linear relationships, so we're just giving it a linear relationship of a non-linear thing.

An easy way to visualize this by looking at a graph an exponential relationship. Like the graph of $x^3$:

Here its obvious that this is not linear. If used it as an equation for our logit model, $y = x^3$; we would get bad results. But if we transformed it by taking the log of our equation, $\log(y) = \log(x^3)$. We would get a graph like this:

That looks pretty linear to me.

This process of transforming models so that they can be better expressed in a different mathematical plane is exactly what the Support Vector Machine does for us. The math behind how it does that is not trivial, so if your interested; put on your reading glasses and head over here. Below is the process of implementing a SVM model and examining the results after the SVM transforms our equation into three different mathematical plains. The first is linear, and is similar to our logic model. Next is an exponential, polynomial, transformation and finally a blank transformation.


In [40]:
#set plotting parameters
plt.figure(figsize(8,6))

# create a regression friendly data frame
y, x = dmatrices(formula_ml, data=df, return_type='matrix')

#select which features we would like to analyze
#try chaning the selection here for diffrent output.
#Choose : [2,3] - pretty sweet DBs [3,1] --standard DBs [7,3] -very cool DBs, [3,6] -- very long complex dbs, could take over an hour to calculate! 
feature_1 = 2
feature_2 = 3

X = np.asarray(x)
X = X[:,[feature_1, feature_2]]  


y = np.asarray(y)
y = y.flatten()      # needs to be 1 dimenstional so we flatten. it comes out of dmatirces with a shape. 

n_sample = len(X)

np.random.seed(0)
order = np.random.permutation(n_sample)

X = X[order]
y = y[order].astype(np.float)

# do a cross validation
X_train = X[:.9 * n_sample]
y_train = y[:.9 * n_sample]
X_test = X[.9 * n_sample:]
y_test = y[.9 * n_sample:]

#create a list of the types of kerneks we will use for your analysis
types_of_kernels = ['linear', 'rbf', 'poly']

# specify our color map for plotting the results
color_map = cm.RdBu_r

# fit the model
for fig_num, kernel in enumerate(types_of_kernels):
    clf = svm.SVC(kernel=kernel, gamma=3)
    clf.fit(X_train, y_train)

    figure(fig_num)
    #pl.clf()
    scatter(X[:, 0], X[:, 1], c=y, zorder=10, cmap=color_map)

    # Circle out the test data
    scatter(X_test[:, 0], X_test[:, 1], s=80, facecolors='none', zorder=10)
    
    axis('tight')
    x_min = X[:, 0].min()
    x_max = X[:, 0].max()
    y_min = X[:, 1].min()
    y_max = X[:, 1].max()

    XX, YY = np.mgrid[x_min:x_max:200j, y_min:y_max:200j]
    Z = clf.decision_function(np.c_[XX.ravel(), YY.ravel()])

    # Put the result into a color plot
    Z = Z.reshape(XX.shape)
    pcolormesh(XX, YY, Z > 0, cmap=color_map)
    contour(XX, YY, Z, colors=['k', 'k', 'k'], linestyles=['--', '-', '--'],
               levels=[-.5, 0, .5])

    title(kernel)
    show()


<matplotlib.figure.Figure at 0x17abe9e8>

Any value in the blue survived while anyone in the read did not. Checkout the graph for the linear transformation. It created its decision boundary right on 50%! That guess from earlier turned out to be pretty good. As you can see, the remaining decision boundaries are much more complex than our original linear decision boundary. These more complex boundaries may be able to capture more structure in the dataset, if that structure exists, and so might create a more powerful predictive model.

Pick a decision boundary that you like, adjust the code below, and submit the results to Kaggle to see how well it worked!


In [41]:
# Here you can output which ever result you would like by changing the Kernel and clf.predict lines

clf = svm.SVC(kernel='poly', gamma=3).fit(X_train, y_train) # Change kernel here to poly, rbf or linear
                                                            # adjusting the gamma level also changes the degree to which the model is fitted

y,x = dmatrices(formula_ml, data=test_data, return_type='dataframe')

res_svm = clf.predict(x.ix[:,[6,3]].dropna()) # Change the interger values within x.ix[:,[6,3]].dropna() explore the 
                                              # relationships between other features. the ints are column postions 
                                              # ie [6,3] 6th column and the third column are evaluated. 

res_svm = DataFrame(res_svm,columns=['Survived'])
res_svm.to_csv("svm_poly_63_g10.csv") # saves the results for you, change the name as you please.

Random Forest

"Well, What if this line / decision boundary thing doesn’t work at all."

Wikipedia, crystal clear as always:

Random forests are an ensemble learning method for classification (and regression) that operate by constructing a multitude of decision trees at training time and outputting the class that is the mode of the classes output by individual trees.

Once again, the skinny and why it matters to you:

There are always skeptics, and you just might be one about all the fancy lines we've created so far. Well for you, here’s another option; the Random Forest. This technique is a form of non-parametric modeling that does away with all those equations we created above, and uses raw computing power and a clever statistical observation to tease the structure out of the data.

An anecdote to explain how this the forest works starts with the lowly gumball jar. We've all guess how many gumballs are in that jar at one time or another, and odds are not a single one of us guessed exactly right. Interestingly though, while each of our individual guesses for probably were wrong, the average of all of the guesses, if there were enough, usually comes out to be pretty close to the actual number of gumballs in the jar. Crazy, I know. This idea is that clever statistical observation that lets random forests work.

How do they work? A random forest algorithm randomly generates many extremely simple models to explain the variance observed in random subsections of our data. These models are like our gumball guesses. They are all awful individually. Really awful. But once they are averaged, they can be powerful predictive tools. The averaging step is the secret sauce. While the vast majority of those models were extremely poor; they were all as bad as each other on average. So when their predictions are averaged together, the bad ones average their effect on our model out to zero. The thing that remains, if anything, is one or a handful of those models have stumbled upon the true structure of the data. The cell below shows the process of instantiating and fitting a random forest, generating predictions form the resulting model, and then scoring the results.


In [42]:
# import the machine learning library that holds the randomforest
import sklearn.ensemble as ske

# Create the random forest model and fit the model to our training data
y, x = dmatrices(formula_ml, data=df, return_type='dataframe')

#instantiate and fit our model
results_rf = ske.RandomForestClassifier(n_estimators=100).fit(x,y)

# Score the results
z = results_rf.score(x,y)
print "Mean accuracy of Random Forest Predictions on the data was: " + str(z)


Mean accuracy of Random Forest Predictions on the data was: 0.533802550183

Our random forest performed only slightly better than a thumb wave, meaning that if you randomly assigned 1s and 0s by waving your thumb up and down you would do almost as well on average. It seems that this time our random forest did not stumble on the true structure of the data.

These are just a few of the machine learning techniques that you can apply. Try a few for yourself and move up the leader board!

Ready to see more an example of a more advanced analysis? Check out these note books:

Follow me on github for more books to come soon!

Version Control


In [1]:
!git add Titanic.ipynb

In [2]:
!git commit


[master 3bdc76c] updated book links
 1 file changed, 5 insertions(+), 5 deletions(-)

In [ ]: