Набор данных «Титаник» — хорошая площадка для отработки ключевых навыков науки о данных. Здесь я хочу показать полное руководство по исследовательскому анализу данных, очистке данных, разработке функций и выбору модели с помощью python, pandas, seaborn, matplotlib и, наконец, sklearn. Затем мы подойдем к настройке гиперпараметров с помощью GridSearchCV и, наконец, объединим наши модели с ансамблем.

Полный код этой статьи доступен на github.

1. Начните импортировать библиотеки

Хорошо, главное, конечно, импортировать все библиотеки, которые нам нужны по пути: мы будем использовать pandas и seaborn для обработки, очистки и визуализации нашего набора данных. «re» будет полезен для обработки строк, а LabelEncoder будет использоваться для кодирования категориальных переменных в числовые.

import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
import seaborn as sns
from sklearn.preprocessing import LabelEncoder
import re

2. Исследовательский анализ данных

мы импортируем наш набор данных, скачанный с kaggle:

df = pd.read_csv('train.csv')

Таким образом, в столбцах указывается PassengerId, класс пассажира, его имя, пол, возраст, если он был с братом, сестрой/супругой, если у него был родитель-ребенок, оплаченный проезд и кабина и, наконец, место посадки пассажира.

Мы проверяем с помощью метода info dtypes и ненулевые счетчики для наших столбцов и видим, что Cabin практически всегда отсутствует, в то время как Age имеет некоторые пропущенные значения. Также обратите внимание, что embarked имеет два пропущенных значения.

Прежде чем устранять проблемы, мы хотим провести предварительный анализ данных. Давайте проверим, можем ли мы найти корреляции между имеющимися у нас переменными и выходными данными (выжили да/нет). Мы можем получить первое представление о наших данных, построив тепловую карту корреляции с Seaborn.

Метод корреляции по умолчанию с .corr() — метод Пирсона, но мы также можем ввести коэффициент корреляции Спирмена. Разница между ними выделена ниже, где у нас есть нелинейная связь между независимой и зависимой переменной, а корреляция Спирмена равна 1.

Мы можем использовать оба из них с методом .corr() с фреймом данных pandas, по умолчанию, если pearson.

Из наших тепловых карт мы видим, что существует отрицательная корреляция с Pclass и положительная с Fare (оба указывают на то, что чем больше вы заплатили, тем больше вероятность, что вы выжили). Однако корреляция с Age, SibSp и Parch почти нулевая. Поскольку мы еще не закодировали наши переменные здесь, у нас нет корреляции с Sex and Embark.

Мы можем делать красивые графики с помощью matplotlib и исследовать больше переменных, таких как пол и Pclass.

survivors_data = df[df.Survived==True]
non_survivors_data = df[df.Survived==False]
# calculate values for each survival status
survivors_gender = survivors_data.groupby(['Sex']).size().values
non_survivors_gender = non_survivors_data.groupby(['Sex']).size().values
# calculate totals for percentates
totals = survivors_gender + non_survivors_gender
# use calculate_percentage_function to calculate percentage of the total
data1_percentages = calculate_percentage(survivors_gender, totals)*100 
data2_percentages = calculate_percentage(non_survivors_gender, totals)*100
gender_categories = ['Female', 'Male']
f, (ax1, ax2) = plt.subplots(1, 2, figsize=(10,5))
# plot chart for count of survivors by class
ax1.bar(range(len(survivors_gender)), survivors_gender, label='Survivors', alpha=0.5, color='g')
ax1.bar(range(len(non_survivors_gender)), non_survivors_gender, bottom=survivors_gender, label='Non-Survivors', alpha=0.5, color='r')
plt.sca(ax1)
plt.xticks([0.4, 1.4], gender_categories )
ax1.set_ylabel("Count")
ax1.set_xlabel("")
ax1.set_title("Count of survivors by gender",fontsize=14)
plt.legend(loc='upper left')
# plot chart for percentage of survivors by class
ax2.bar(range(len(data1_percentages)), data1_percentages, alpha=0.5, color='g')
ax2.bar(range(len(data2_percentages)), data2_percentages, bottom=data1_percentages, alpha=0.5, color='r')
plt.sca(ax2)
plt.xticks([0.4, 1.4],  gender_categories)
ax2.set_ylabel("Percentage")
ax2.set_xlabel("")
ax2.set_title("% of survivors by gender",fontsize=14)

И тот же код применяется для любой другой переменной:

Итак, у нас есть представление о том, кто, скорее всего, выжил: женщины и представители высшего класса. С Seaborn мы можем сгруппировать информацию о поле, возрасте, P-классе:

3. Очистка данных

Сделаем чистку. Начнем с возраста: поскольку мы видим, что класс и пол могут повлиять на нашу модель, мы будем вводить среднее значение категории каждой группы вместо простого среднего и медианы:

df['Filled_age'] = df['Age'].fillna(df.groupby(['Pclass','Sex'])['Age'].transform(np.mean))

Что касается Хижины, мы постараемся просто сохранить то, что у нас есть на данный момент, сохранив колоду, заполнив недостающие значения строкой, а затем выбрав первую букву строки в Хижине. Менее распространенные каюты заменены на «другие».

df['Cabin'].fillna('Missing',inplace=True)
df['Deck'] = df.Cabin.apply(lambda x: x[0])
deck_mapping= {'M':'M','C':'C','B':'B','D':'D','E':'E','A':'other','F':'other','G':'other','T':'other'}
df['Deck'] = df['Deck'].map(deck_mapping)

Что касается места посадки, то значения na заполняем наиболее частым значением.

df['Embarked'].fillna('S',inplace=True)

4. Разработка функций

Мы получаем фиктивные кодировки для местоположения Ses, Deck и Embarked, а также удаляем исходные столбцы, а также PassengerId и Ticket:

dummies = pd.get_dummies(df[['Sex','Embarked','Deck']], drop_first=True)
df2 = pd.concat([df,dummies],axis=1).drop(['Sex','Embarked','Deck'],1).drop(['PassengerId','Age','Cabin','Ticket'],1)

Обратите внимание, что мы использовали «drop_first=True» в функции get_dummies. Это удалит один столбец для каждой категории, поскольку он избыточен (для этого набора данных, если мужчина = 0, тогда женщина = 1).

Далее мы хотим обратиться к переменным Name, Fare и SibSp, Parch. Мы признаем возможность того, что у семей были разные шансы на выживание, и что у людей, путешествующих в одиночку, также были разные шансы:

df2['FamilySize'] = df2['Parch']+df2['SibSp']+1
df2['IsAlone'] = 1 #initialize to yes/1 is alone
df2['IsAlone'].loc[df2['FamilySize'] > 1] = 0 # now update to no/0 if family size is greater than 1

Затем мы экстраполируем из имени титул («мистер», «мисс» и т. д.) и фамилию, также мы хотим подсчитать общее количество фамилий, которые должны дать нам то же число, что и FamilySize, по этой причине мы также нормализуем стоимость проезда для семьи:

df2['Surname'] = df2['Name'].str.split(',').str[0]
df2['Title'] = df2['Name'].str.split(',').str[-1].str.split('.').str[0]   
df2['Surname_counts'] = df2.groupby('Surname')['Name'].transform('count')
df2['Surname_Family_counts'] = df2.groupby(['Surname','Parch','SibSp'])['Name'].transform('count')
df2['Norm_Fare'] = df2['Fare']/df2['Surname_Family_counts']

Затем мы разделяем переменные, такие как Title, сохраняя наиболее часто встречающиеся, используя метод квантильного сокращения панд для Fare и также преобразуя возраст в 5 интервалов:

top5 = df2.Title.value_counts().nlargest(5).index
df2['Less_Title'] = df2.Title.where(df2.Title.isin(top5), other='Other')
df2['FareBin'] = pd.qcut(df2['Norm_Fare'], 4)
df2['AgeBin'] = pd.cut(df2['Filled_age'].astype(int), 5)
label = LabelEncoder()
    
df2['Title_Code'] = label.fit_transform(df2['Title'])
df2['AgeBin_Code'] = label.fit_transform(df2['AgeBin'])
df2['FareBin_Code'] = label.fit_transform(df2['FareBin'])
df2.drop(['Name','Surname','Title','Ticket2','Filled_age','Ticket_num','Fare','FareBin','AgeBin','Less_Title'],1,inplace=True)

Наш набор данных не полный:

Мы можем снова проверить корреляции, чтобы увидеть, получим ли мы больше информации.

Теперь мы видим, что наиболее полезными предикторами являются Pclass, Sex_male и FareBin_code.

5. Выбор модели

Мы импортируем кучу модулей!

Чтобы проверить наши результаты, мы хотим проверить самую простую модель и сравнить с ней результаты: поэтому мы сохраняем нашу модель простой и выбираем только 4 переменные, мы разделяем наш набор данных на обучающий и тестовый набор данных, а затем используем прогнозы.

Точность 81% — это неплохо, но мы можем улучшить. Давайте проверим, могут ли некоторые более сложные модели работать лучше. Мы импортируем DecisionTreeClassifier, RandomForestClassifier, Bagging, AdaBoost и Gradient BoostingClassifier и KNN.

lr = LogisticRegression()
dtree = DecisionTreeClassifier(random_state=42,max_depth=3,min_samples_leaf=1)
rfc = RandomForestClassifier(random_state=42,max_depth=3)
bc = BaggingClassifier(base_estimator=lr, n_estimators=300, n_jobs=-1)
knn = KNeighborsClassifier(n_neighbors=5)
ada = AdaBoostClassifier(base_estimator=dtree,n_estimators=300,random_state=42)
bg = GradientBoostingClassifier(random_state=42, max_depth=2, subsample=0.8,n_estimators=300)
estimators = [('dtree', dtree),('bc', bc), ('rfc', rfc), ('ada', ada),('bg',bg),
              ('lr', lr),  ('knn',knn)]

теперь мы сравниваем каждую модель, для этого мы организуем наши результаты в кадре данных

results = []
for name,model in estimators:
    model.fit(X_train_std,y_train)
    preds = model.predict(X_test_std)
    rocauc = metrics.roc_auc_score(y_test, preds)
    f1 = metrics.f1_score(y_test, preds)
    prec = metrics.precision_score(y_test, preds)
    rec = metrics.recall_score(y_test, preds)
    accu = metrics.accuracy_score(y_test, preds)
     
    try:
        scores = cross_val_score(model, X, y, scoring='accuracy', cv=5)
        scores_mean = scores.mean()
        scores_std = scores.std()
        predict_prob = model.predict_proba(X_test)[:,1]
        roc_auc_model = roc_auc_score(y_test, predict_prob)
        
    except:
        scores_mean = np.nan
        scores_std = np.nan
        print('model ',name,' has not predict_proba')
        roc_auc_model = np.nan
        
    data=[name,accu,rocauc,roc_auc_model,f1,prec,rec,scores_mean, scores_std]    
    results.append(data)
mod_res = pd.DataFrame(results,columns=['model','accuracy','rocauc','roc_auc_prob','f1','prec',
                              'rec','scores_mean','scores_std']).sort_values(by='scores_mean',ascending=False)

Жесткая игра! RandomForest, DecisionTree, LogisticRegression и BaggingClassifier имеют почти одинаковую точность. RFC продемонстрировал высочайшую точность, в то время как LR и BC — более высокую полноту.

Например, мы можем также проверить важность функций RFC. Мы возвращаем некоторые знания, которые у нас были: пол был самой важной характеристикой, за ней следовали Pclass и FamilySize.

Полные переменные

Давайте проверим производительность модели со всеми переменными:

target = 'Survived'
X = df2.drop('Survived',1)
y = df2['Survived']
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.30, random_state=42)
scaler = StandardScaler()
X_train_std = scaler.fit_transform(X_train)
X_test_std = scaler.transform(X_test)

Затем мы используем тот же код, что и выше, с этими результатами:

Таким образом, DecisionTree и RandomForest работают лучше, хотя точность по-прежнему составляет около 82%. Из графика важности функций мы видим, что Title_code с FamilySize, Norm_Fare, Deck_M сыграли свою роль в модели. Как мы можем настроить модели?

6. Настройка гиперпараметров

Мы будем использовать GridSearchCV для настройки параметров наших моделей. Значение по умолчанию для переоснащения в GridSearchCV равно True, поэтому оно будет переоснащать наши модели с лучшими параметрами.

GridSearchCV(
    estimator,
    param_grid,
    scoring=None,
    n_jobs=None,
    iid='deprecated',
    refit=True,
    cv=None,
    verbose=0,
    pre_dispatch='2*n_jobs',
    error_score=nan,
    return_train_score=False,
)

Ниже код разделен для каждой модели. Что мы сделаем, так это создадим сетку параметров для каждой модели и позволим GridSearchCV найти лучшие параметры с 6 перекрестными проверками с оценкой = «точность».

cv = 6
scoring='accuracy'
print('knn')
param_grid = {'n_neighbors':np.arange(2,30,2),
             'p':[1,2,3]}
gs = GridSearchCV(knn,param_grid,scoring=scoring,n_jobs=-1, cv=cv, verbose=1,)
gs.fit(X_train,y_train)
print('best score:', gs.best_score_, ' best params', gs.best_params_,' best test score', gs.score(X_test,y_test))
print('LR')
param_grid = {'penalty' : ['l1', 'l2', 'elasticnet'],
             'C':[1,2,2.2,2.5],
             'class_weight':[None,'balanced']}
gs = GridSearchCV(lr,param_grid,scoring=scoring,n_jobs=-1, cv=cv, verbose=1)
gs.fit(X_train,y_train)
print('best score:', gs.best_score_, ' best params', gs.best_params_,' best test score', gs.score(X_test,y_test))
print('DTREE')
param_grid = {'max_depth':[3,4,5,6,7],
             'class_weight':[None,'balanced'],
             'criterion':['gini'],
             'min_samples_split':[1,2,5,10],
             'min_samples_leaf':[0.005,0.01] ,
             'ccp_alpha':[0,0.005,0.007]}
gs = GridSearchCV(dtree,param_grid,scoring=scoring,n_jobs=-1, cv=cv, verbose=1)
gs.fit(X_train,y_train)
print('best score:', gs.best_score_, ' best params', gs.best_params_,' best test score', gs.score(X_test,y_test))
print('RFC')
param_grid = {'max_depth':[5,6,7],
             'class_weight':[None],
             'criterion':['entropy'],
             'min_samples_split':[1,2,4],
             'min_samples_leaf':[0.01,0.1] ,
             'ccp_alpha':[0,0.005],
             'n_estimators':[450,500,550]}
gs = GridSearchCV(rfc,param_grid,scoring=scoring,n_jobs=-1, cv=cv, verbose=1)
gs.fit(X_train,y_train)
print('best score:', gs.best_score_, ' best params', gs.best_params_,' best test score', gs.score(X_test,y_test))
print('BG')
param_grid = {'max_depth':[3,5,6,7],
             'n_estimators':[75,150,350],
             'subsample':[0.5,0.65,0.8]}
gs = GridSearchCV(bg,param_grid,scoring=scoring,n_jobs=-1, cv=cv, verbose=1)
gs.fit(X_train,y_train)
print('best score:', gs.best_score_, ' best params', gs.best_params_,' best test score', gs.score(X_test,y_test))

Это займет некоторое время. Вывод будет примерно таким:

Оценка:

Мы достигли 82% с помощью классификатора дерева решений, за которым следуют классификаторы случайного леса и бэггинга. Судя по средним баллам классификаторов мешков CV, они кажутся лучшими.

Ансамбль

Что мы можем сделать сейчас, так это сгруппировать лучшие классификаторы и попытаться увидеть, улучшим ли мы наши результаты. Создаем три ансамбля моделей. В первой модели будет использоваться голосование='жесткий', во второй и третьей модели будет использоваться голосование='мягкий'. В третьей модели мы определяем веса каждой модели, давая небольшое преимущество одной из них.

v_ests = [('rfc', rfc), ('bg', bg), ('lr', lr), ('dtree',dtree)]
eclf = VotingClassifier(estimators=v_ests, voting='hard')
v_ests = [('rfc', rfc), ('bg', bg), ('lr', lr), ('dtree',dtree)]
eclfsoft = VotingClassifier(estimators=v_ests, voting='soft')
v_ests_w = [('rfc', rfc), ('bg', bg), ('dtree',dtree)]
eclfsoft_w = VotingClassifier(estimators=v_ests, voting='soft', weights=[1,2,1])
estimators = [('Hard',eclf),   ('Soft',eclfsoft),   ('Soft_w',eclfsoft_w)]

Из приведенной ниже таблицы видно, что мы могли бы добиться большего, чем раньше, глядя на средние баллы по результатам оценки резюме. мы достигли 83%, что неплохо.

Выводы

Этот учебник охватывает большинство аспектов типичного рабочего процесса. Это далеко от совершенства, и цель состоит в том, чтобы просто пересмотреть шаги по очистке набора данных, создать новые функции, оценить их актуальность, а также внедрить и настроить модель.