こんにちは。CCCマーケティング株式会社 技術開発ユニットの伊藤です。前回の記事(いまさら聞けないKaggleチュートリアル(2)　登録→提出までやってみよう！)では、kaggleへ登録し、タイタニックコンペで予測の提出をしました。前回は数値の列だけをロジスティック回帰に入力するという雑なことをしたので、スコアは0.66985、順位は20000位くらいという結果でした。今回はもう少し良い結果が出るように頑張ってみます。具体的には、

  1. 特徴量の作成
  2. ロジスティック回帰とランダムフォレストのアンサンブル

を実行していきます。前回同様kaggleのノートブック上で分析をします。

ライブラリの読み込み

今回使用したライブラリです。必要なライブラリが分かるように最初のセルでライブラリのインポートを行います。

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns

from sklearn.pipeline import make_pipeline
from sklearn.impute import KNNImputer
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LogisticRegression
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import StratifiedKFold
from sklearn.metrics import accuracy_score

データの読み込み

file = '/kaggle/input/titanic/{}.csv'
train = pd.read_csv(file.format('train'))
test = pd.read_csv(file.format('test'))
sample_submission = pd.read_csv(file.format('gender_submission'))

データは以下のような列で構成されています。

  • Survival: 生存なら1、死亡なら0
  • Pclass: チケットのクラス
  • Sex: 性別
  • Age: 年齢
  • SibSp: タイタニックに乗船していた兄弟と配偶者の数
  • Parch: タイタニックに乗船していた親と子の数
  • Ticket: チケット番号
  • Cabin: 客室番号
  • Embarked: 乗船した港

プロット

グラフをプロットすることで目的変数と説明変数の関係を確認したり、外れ値や分布の偏りを確認します。

目的変数が二値なのでカテゴリ変数はカウントプロット、数値型の変数はカーネル密度推定でグラフをプロットし、データを確認しました。以下はSexとAgeをプロットしたときの例です。seabornの関数を利用してプロットしています。

sns.countplot('Survived', hue='Sex', data=train)
plt.show()

  • 性別ごとの生存(1)、死亡(0)者数

sns.kdeplot(train.query('Survived == 0')['Age'], label='0')
sns.kdeplot(train.query('Survived == 1')['Age'], label='1')
plt.show()

  • 生存(1)、死亡(1)ごとの年齢の分布

グラフから女性や子供が優先的に救命ボートに乗せられたのではないかという仮説を立てました。分布の形状から適切な変数の変換を行ったり、仮説をもとに特徴量を作ったりしていきます。

前処理

以下のような前処理の関数を作成しました。

def preprocess(train, test):
    # train, testを連結
    len_train = len(train)
    df = pd.concat([train, test])

    # 変数変換
    df['Age'] = np.log1p(df['Age'])
    df['Fare'] = np.log1p(df['Fare'])

    # 列の追加
    df['Title'] = df['Name'].str.extract('.*(Mr|Miss|Mrs|Master)\..*', expand=False).fillna('Other')
    df['Family'] = df['SibSp'] + df['Parch']
    df['Cabin_NaN'] = df['Cabin'].isna().astype('int')

    # onehot encoding
    df = pd.get_dummies(df, drop_first=True, columns=['Sex'])
    df = pd.get_dummies(df, drop_first=False, columns=['Embarked', 'Title'])

    # 列の削除
    df = df.drop(columns=[
        'PassengerId', 'Name', 'SibSp', 'Parch', 'Ticket', 'Cabin'
    ])

    # train, testの分割
    train = df.iloc[:len_train]
    test = df.iloc[len_train:].drop(columns='Survived')

    return train, test

変数変換

対数変換は決定木ベースのモデルでは関係ないですが、今回はロジスティック回帰も使用するので変数変換によって性能が変わります。グラフをプロットすることで分布に偏りがあることが分かったので、Age、Fareをlog(x + 1)で変換しました。

列の追加

与えられたデータから新たに特徴量を作成することでモデルを改善します。以下の2つの特徴量を新たに作成しました。

  • Title: Nameから抽出した敬称
  • Cabin_NaN: Cabinが欠損かどうかを表すフラグ

Title

Nameに含まれる敬称を正規表現を利用して抽出してみました。正規表現についての説明は省きますが、Mrなどの敬称の後にピリオドがついているという点に着目して敬称を抽出しています。

title, counts = np.unique(train['Name'].str.extract('.*\s(.+)\..*', expand=False), return_counts=True)
dict(zip(title, counts)))

{'Capt': 1,
 'Col': 2,
 'Countess': 1,
 'Don': 1,
 'Dr': 7,
 'Jonkheer': 1,
 'L': 1,
 'Lady': 1,
 'Major': 2,
 'Master': 40,
 'Miss': 182,
 'Mlle': 2,
 'Mme': 1,
 'Mr': 517,
 'Mrs': 124,
 'Ms': 1,
 'Rev': 6,
 'Sir': 1})

抽出した敬称から性別や年齢が推定できそうです。Mrのような見慣れた敬称が含まれていますが、私が知らない敬称も含まれていました。

WikipediaによるとDonはスペイン語圏とポルトガル語圏で男性に使われる敬称らしいです。 DonをMrにまとめるなどの作業を行うのも良いかもしれませんが、面倒なので出現回数の多いMaster、Miss、Mr、Mrsを抽出し、それ以外の敬称はOtherにまとめました。

Cabin_NaN

Cabinは訓練データ中77%が欠損値でした。欠損であることにも意味があるのではないかと考え、欠損値かどうかを表す列を作成しました。

onehot encoding

pd.get_dummiesでカテゴリ変数をonehot表現に変換しました。ランダムフォレストに直接入力する場合はラベルエンコーディングを使ったほうが良いかもしれません。

モデル

ロジスティック回帰

ロジスティック回帰の訓練と予測確率を出力する関数を作成しました。前処理で欠損値の処理をしていなかったので、KNNImputerで補完することにしました。

def fit_predict_lr(train, test):
    len_train = len(train)
    df = pd.concat([train, test])

    df = df.drop(columns=['Embarked_C', 'Title_Master'])

    train = df.iloc[:len_train]
    test = df.iloc[len_train:].drop(columns='Survived')

    estimator = make_pipeline(
        StandardScaler(),
        KNNImputer(n_neighbors=5),
        LogisticRegression(C=0.1, random_state=42)
    )
    estimator.fit(train.drop(columns='Survived'), train['Survived'])
    y_proba = estimator.predict_proba(test)[:, 1]

    return y_proba

ランダムフォレスト

ランダムフォレストでも訓練と予測確率を出力する関数を作成しました。

def fit_predict_rf(train, test):
    if 'Survived' in test.columns:
        test = test.drop(columns='Survived')

    estimator = make_pipeline(
        StandardScaler(),
        KNNImputer(),
        RandomForestClassifier(n_estimators=128, max_depth=4, random_state=42)
    )

    estimator.fit(train.drop(columns='Survived'), train['Survived'])
    y_proba = estimator.predict_proba(test)[:, 1]

    return y_proba

アンサンブル

アンサンブルは複数のモデルを組み合わせる方法でバギングやブースティングなどの手法があります。今回はソフト投票で予測を作成します。

ソフト投票は複数のモデルからの予測確率の加重平均を最終的な予測確率とする手法です。scikit-learnのVotingClassifierというクラスで利用できますが、今回は使用していません。

以下のようにクロスバリデーションで正解率がどの程度か確認してみました。

fold = StratifiedKFold(n_splits=5, shuffle=True, random_state=42)
for train_idx, val_idx in fold.split(train2.drop(columns='Survived'), train2['Survived']):
    y_proba_lr = fit_predict_lr(train2.iloc[train_idx], train2.iloc[val_idx])
    y_proba_rf = fit_predict_rf(train2.iloc[train_idx], train2.iloc[val_idx])
    print('-' * 20)
    print('Accuracy(LR):', accuracy_score(train2['Survived'].iloc[val_idx], y_proba_lr >= 0.5))
    print('Accuracy(RF):', accuracy_score(train2['Survived'].iloc[val_idx], y_proba_rf >= 0.5))
    y_proba_mean = 0.5 * y_proba_lr + 0.5 * y_proba_rf
    print('Accuracy(mean):', accuracy_score(train2['Survived'].iloc[val_idx], y_proba_mean >= 0.5))
    print('-' * 20)

--------------------
Accuracy(LR): 0.8379888268156425
Accuracy(RF): 0.8379888268156425
Accuracy(mean): 0.8435754189944135
--------------------
--------------------
Accuracy(LR): 0.8202247191011236
Accuracy(RF): 0.8314606741573034
Accuracy(mean): 0.8258426966292135
--------------------
--------------------
Accuracy(LR): 0.8089887640449438
Accuracy(RF): 0.8314606741573034
Accuracy(mean): 0.8202247191011236
--------------------
--------------------
Accuracy(LR): 0.8146067415730337
Accuracy(RF): 0.8258426966292135
Accuracy(mean): 0.8258426966292135
--------------------
--------------------
Accuracy(LR): 0.8258426966292135
Accuracy(RF): 0.8370786516853933
Accuracy(mean): 0.8314606741573034
--------------------

最終的にロジスティック回帰とランダムフォレストの重みはどちらも0.5にしました。

y_proba_lr = fit_predict_lr(train2, test2)
y_proba_rf = fit_predict_rf(train2, test2)
y_proba = 0.5 * y_proba_lr + 0.5 * y_proba_rf
y_pred = (y_proba >= 0.5).astype(int)

結果

スコアは0.78708で順位は3665位でした。約21000チームいるので上位20%くらいの順位になりました。

  • リーダーボード

おわりに

今回は特徴量の作成とアンサンブルでスコアアップに挑戦しました。ほかにもk-NNなどの別のアルゴリズムを使用したり、特徴量の追加や交互作用を考慮するなど工夫できる点があると思うので今後もスコアアップに挑戦していきたいです。

本記事は、CCC MARKETING TECH BLOGを、再編集して転載したものです。