kaggleのtitanic xgboostを使った生存者予測 [80.1%]

#はじめに
前回はタイタニックのtrainデータを見て、生存率に関係のある特徴量とは何かを見てみました。
前回：[kaggleのtitanic 特徴量生成の前のデータ分析](https://qiita.com/ziel40/items/ae98aa1965fc29cd864d)
(こちらの記事も合わせて見ていただけると嬉しいです)

その結果を元に、今回は**欠損値の補完と必要な特徴量を追加して予測モデルに与えるデータフレームを作成します**。

なお、今回は予測モデルとして、kaggleのコンペでもよく使われている「**GBDT(勾配ブースティング木)のxgboost**」を使うのでデータをそれに適した形にしていきます。

コードはGitHubでも公開しています。
[xgboost.py](https://github.com/serikatsu28/kaggle/blob/titanic/xgboost.py)

#1. データの取得と欠損値の確認
```python
import pandas as pd
import numpy as np

train = pd.read_csv('/kaggle/input/titanic/train.csv')
test = pd.read_csv('/kaggle/input/titanic/test.csv')
#trainデータとtestデータを１つにまとめる
data = pd.concat([train,test]).reset_index(drop=True)
#欠損値が含まれる行数を確認
train.isnull().sum()
test.isnull().sum()
```
それぞれの欠損値の数は下のようになりました。

||trainデータ|testデータ|
|:---:|:---:|:---:|
PassengerId|0|0|
Survived|0|
Pclass|0|0|
Name|0|0|
Sex|0|0|
Age|177|86|
SibSp|0|0|
Parch|0|0|
Ticket|0|0|
Fare|0|1|
Cabin|687|327|
Embarked|2|0|

#2. 欠損値の補完と特徴量の作成
##2.1 Fareの補完と階級分け
欠損している行の**Pclassは3**で**EmbarkedはS**でした。
![スクリーンショット 2020-10-01 12.01.00.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/654235/d20ba268-dfdc-1a56-9c40-02c4db3e05e9.png)
ですのでこの2つの条件を満たす人の中での中央値で補完します。
その後、Fareの値による生存率の違いとPclassを考慮して階級分けを行います。


```python
#欠損値の補完
data['Fare'] = data['Fare'].fillna(data.query('Pclass==3 & Embarked=="S"')['Fare'].median())
#階級分けしたものを'Fare_bin'に入れる
data['Fare_bin'] = 0
data.loc[(data['Fare']>=10) & (data['Fare']<50), 'Fare_bin'] = 1
data.loc[(data['Fare']>=50) & (data['Fare']<100), 'Fare_bin'] = 2
data.loc[(data['Fare']>=100), 'Fare_bin'] = 3
```

##2.2 グループごとの生死の違い'Family_survival'の作成
[Titanic [0.82] - [0.83]](https://www.kaggle.com/konstantinmasich/titanic-0-82-0-83)こちらのコードで紹介されていた 'Family_survival' という特徴量を作成します。

**家族や友人同士だと船内で行動を共にしている可能性が高い**ので生存できたかどうかも**グループ内で同じ結果になりやすい**といえます。

そこで名前の名字とチケット番号でグルーピングを行い、そのグループのメンバーが生存しているかどうかで値を決めます。

```python
#名前の名字を取得して'Last_name'に入れる
data['Last_name'] = data['Name'].apply(lambda x: x.split(",")[0])


data['Family_survival'] = 0.5 #デフォルトの値
#Last_nameとFareでグルーピング
for grp, grp_df in data.groupby(['Last_name', 'Fare']):
                               
    if (len(grp_df) != 1):
        #(名字が同じ)かつ(Fareが同じ)人が2人以上いる場合
        for index, row in grp_df.iterrows():
            smax = grp_df.drop(index)['Survived'].max()
            smin = grp_df.drop(index)['Survived'].min()
            passID = row['PassengerId']
            
            if (smax == 1.0):
                data.loc[data['PassengerId'] == passID, 'Family_survival'] = 1
            elif (smin == 0.0):
                data.loc[data['PassengerId'] == passID, 'Family_survival'] = 0
            #グループ内の自身以外のメンバーについて
            #1人でも生存している → 1
            #生存者がいない(NaNも含む) → 0
            #全員NaN → 0.5

#チケット番号でグルーピング
for grp, grp_df in data.groupby('Ticket'):
    if (len(grp_df) != 1):
        #チケット番号が同じ人が2人以上いる場合
        #グループ内で1人でも生存者がいれば'Family_survival'を1にする
        for ind, row in grp_df.iterrows():
            if (row['Family_survival'] == 0) | (row['Family_survival']== 0.5):
                smax = grp_df.drop(ind)['Survived'].max()
                smin = grp_df.drop(ind)['Survived'].min()
                passID = row['PassengerId']
                if (smax == 1.0):
                    data.loc[data['PassengerId'] == passID, 'Family_survival'] = 1
                elif (smin == 0.0):
                    data.loc[data['PassengerId'] == passID, 'Family_survival'] = 0
```

##2.3 家族の人数を表す特徴量 'Family_size' の作成と階級分け
SibSpとParchの値を使って何人家族でタイタニックに乗船したかという特徴量'Family_size'をつくり、それぞれの生存率に応じた階級分けをします。

```python
#Family_sizeの作成
data['Family_size'] = data['SibSp']+data['Parch']+1
#1, 2~4, 5~の3つに分ける
data['Family_size_bin'] = 0
data.loc[(data['Family_size']>=2) & (data['Family_size']<=4),'Family_size_bin'] = 1
data.loc[(data['Family_size']>=5) & (data['Family_size']<=7),'Family_size_bin'] = 2
data.loc[(data['Family_size']>=8),'Family_size_bin'] = 3
```

##2.4 名前の敬称 'Title' の作成
Nameの列から'Mr','Miss'などの敬称を取得します。
数が少ない敬称('Mme','Mlle'など)は同じ意味を表す敬称に統合します。

```python
#名前の敬称を取得して'Title'に入れる
data['Title'] = data['Name'].map(lambda x: x.split(', ')[1].split('. ')[0])
#数の少ない敬称を統合
data['Title'].replace(['Capt', 'Col', 'Major', 'Dr', 'Rev'], 'Officer', inplace=True)
data['Title'].replace(['Don', 'Sir',  'the Countess', 'Lady', 'Dona'], 'Royalty', inplace=True)
data['Title'].replace(['Mme', 'Ms'], 'Mrs', inplace=True)
data['Title'].replace(['Mlle'], 'Miss', inplace=True)
data['Title'].replace(['Jonkheer'], 'Master', inplace=True)
```

##2.5 チケット番号のラベリング
チケット番号には数字のみのものと数字とアルファベットを含むものがあります。
その２つを分けてからチケット番号の種類ごとにラベリングを行います。

```python
#数字のみのチケットと数字とアルファベットを含むチケットに分ける
#数字のみのチケットを取得
num_ticket = data[data['Ticket'].str.match('[0-9]+')].copy()
num_ticket_index = num_ticket.index.values.tolist()
#元のdataから数字のみのチケットの行を落とした残りがアルファベットを含むチケット
num_alpha_ticket = data.drop(num_ticket_index).copy()

#数字のみのチケットの階級分け
#チケット番号は文字列になっているので数値に変換
num_ticket['Ticket'] = num_ticket['Ticket'].apply(lambda x:int(x))

num_ticket['Ticket_bin'] = 0
num_ticket.loc[(num_ticket['Ticket']>=100000) & (num_ticket['Ticket']<200000),
               'Ticket_bin'] = 1
num_ticket.loc[(num_ticket['Ticket']>=200000) & (num_ticket['Ticket']<300000),
               'Ticket_bin'] = 2
num_ticket.loc[(num_ticket['Ticket']>=300000),'Ticket_bin'] = 3

#数字とアルファベットを含むチケットの階級分け
num_alpha_ticket['Ticket_bin'] = 4
num_alpha_ticket.loc[num_alpha_ticket['Ticket'].str.match('A.+'),'Ticket_bin'] = 5
num_alpha_ticket.loc[num_alpha_ticket['Ticket'].str.match('C.+'),'Ticket_bin'] = 6
num_alpha_ticket.loc[num_alpha_ticket['Ticket'].str.match('C\.*A\.*.+'),'Ticket_bin'] = 7
num_alpha_ticket.loc[num_alpha_ticket['Ticket'].str.match('F\.C.+'),'Ticket_bin'] = 8
num_alpha_ticket.loc[num_alpha_ticket['Ticket'].str.match('PC.+'),'Ticket_bin'] = 9
num_alpha_ticket.loc[num_alpha_ticket['Ticket'].str.match('S\.+.+'),'Ticket_bin'] = 10
num_alpha_ticket.loc[num_alpha_ticket['Ticket'].str.match('SC.+'),'Ticket_bin'] = 11
num_alpha_ticket.loc[num_alpha_ticket['Ticket'].str.match('SOTON.+'),'Ticket_bin'] = 12 
num_alpha_ticket.loc[num_alpha_ticket['Ticket'].str.match('STON.+'),'Ticket_bin'] = 13
num_alpha_ticket.loc[num_alpha_ticket['Ticket'].str.match('W\.*/C.+'),'Ticket_bin'] = 14

data = pd.concat([num_ticket,num_alpha_ticket]).sort_values('PassengerId')
```

##2.6 Ageの補完と階級分け
Ageの欠損値については**中央値**や**名前の敬称ごとに年齢の平均値を求めて補完**する方法がありますが、調べてみると**ランダムフォレストを使って欠損している部分の年齢を予測させる方法**がありましたので今回はそちらを使ってみます。

```python
from sklearn.preprocessing import LabelEncoder
from sklearn.ensemble import RandomForestRegressor
#文字列になっている特徴量をラベリング
le = LabelEncoder()
data['Sex'] = le.fit_transform(data['Sex']) #生存率予測に使うのでついでにラベリング
data['Title'] = le.fit_transform(data['Title'])
#Ageの予測に使う特徴量を'age_data'にいれる
age_data = data[['Age','Pclass','Family_size',
                 'Fare_bin','Title']].copy()
#Ageが欠損している行と欠損していない行に分ける
known_age = age_data[age_data['Age'].notnull()].values  
unknown_age = age_data[age_data['Age'].isnull()].values

x = known_age[:, 1:]  
y = known_age[:, 0]
#ランダムフォレストで学習
rfr = RandomForestRegressor(random_state=0, n_estimators=100, n_jobs=-1)
rfr.fit(x, y)
#予測値を元のデータフレームに反映する
age_predict = rfr.predict(unknown_age[:, 1:])
data.loc[(data['Age'].isnull()), 'Age'] = np.round(age_predict,1)
```
これでAgeの欠損値を補完できました。
![スクリーンショット 2020-10-15 21.37.08.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/654235/316e5c81-7d0b-3d94-f1d0-bf7d08e525e0.png)

そしてAgeも階級分けしていきます。

```python
data['Age_bin'] = 0
data.loc[(data['Age']>18) & (data['Age']<=60),'Age_bin'] = 1
data.loc[(data['Age']>60),'Age_bin'] = 2
```

最後にいらない特徴量を落とします。

```python
data = data.drop(['PassengerId','Name','Age','SibSp','Parch','Ticket',
                  'Fare','Cabin','Embarked','Last_name','Family_size'], axis=1)
```

最終的にデータフレームはこのようになりました。

||Survived|Pclass|Sex|Fare_bin|Family_survival|Family_size_bin|Title|Ticket_bin|Age_bin|
|---:|---:|---:|---:|---:|---:|---:|---:|---:|---:|
|0|0.0|3|1|0|0.5|1|2|5|1|
|1|1.0|1|0|2|0.5|1|3|9|1|
|2|1.0|3|0|0|0.5|0|1|13|1|
|3|1.0|1|0|2|0.0|1|3|1|1|
|4|0.0|3|1|0|0.5|0|2|3|1|
|...|...|...|...|...|...|...|...|...|...|
|1304|NaN|3|1|0|0.5|0|2|5|1|
|1305|NaN|1|0|3|1.0|0|5|9|1|
|1306|NaN|3|1|0|0.5|0|2|12|1|
|1307|NaN|3|1|0|0.5|0|2|3|1|
|1308|NaN|3|1|1|1.0|1|0|0|0|
1309 rows × 9 columns

統合させていたデータをtrainデータとtestデータに分けて特徴量の処理は終了です。

```python
model_train = data[:891]
model_test = data[891:]

X = model_train.drop('Survived', axis=1)
Y = pd.DataFrame(model_train['Survived'])
x_test = model_test.drop('Survived', axis=1)
```

#3. xgboostでの予測
**logloss**と**accuracy**の2つを求めてモデルの性能をみてみます。


```python
from sklearn.metrics import log_loss
from sklearn.metrics import accuracy_score
from sklearn.model_selection import KFold
import xgboost as xgb
#パラメータを設定
params = {'objective':'binary:logistic',
          'max_depth':5,
          'eta': 0.1, 
          'min_child_weight':1.0,
          'gamma':0.0,
          'colsample_bytree':0.8,
          'subsample':0.8}

num_round = 1000

logloss = []
accuracy = []

kf = KFold(n_splits=4, shuffle=True, random_state=0)
for train_index, valid_index in kf.split(X):
    x_train, x_valid = X.iloc[train_index], X.iloc[valid_index] 
    y_train, y_valid = Y.iloc[train_index], Y.iloc[valid_index]
    #データフレームをxgboostに適した形に変換
    dtrain = xgb.DMatrix(x_train, label=y_train)
    dvalid = xgb.DMatrix(x_valid, label=y_valid)
    dtest = xgb.DMatrix(x_test)
    #xgboostで学習
    model = xgb.train(params, dtrain, num_round,evals=[(dtrain,'train'),(dvalid,'eval')],
                      early_stopping_rounds=50)
    
    valid_pred_proba = model.predict(dvalid)
    #loglossを求める
    score = log_loss(y_valid, valid_pred_proba)
    logloss.append(score)
    #accuracyを求める
    #valid_pred_probaは確率値なので0と1に変換
    valid_pred = np.where(valid_pred_proba >0.5,1,0)
    acc = accuracy_score(y_valid, valid_pred)
    accuracy.append(acc)
    
print(f'log_loss:{np.mean(logloss)}')
print(f'accuracy:{np.mean(accuracy)}')
```
こちらのコードで
**log_loss : 0.39114**
**accuracy : 0.8338**
という結果になりました。

モデルが出来たのでkaggleに提出する予測データを作成します。

```python
#predictで予測
y_pred_proba = model.predict(dtest)
y_pred= np.where(y_pred_proba > 0.5,1,0)
#データフレームを作成
submission = pd.DataFrame({'PassengerId':test['PassengerId'], 'Survived':y_pred})
submission.to_csv('titanic_xgboost.csv', index=False)
```
![スクリーンショット 2020-10-15 21.54.44.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/654235/59912efa-4465-eab3-c544-75e4f9e969f9.png)
正解率は**80.1%**でギリギリ8割に到達しました。

#まとめ
今回実際にkaggleに予測値を提出するところまでやってみました。

初めはCabinもラベリングして特徴量として使用していたのですが**Cabinを使わないほうが予測精度が上がりました**。やはり欠損値が多すぎてxgboostにはあまり適さなかったようです。

また、FareやAgeなどの特徴量を階級分けせずに予測を行うと**loglossやaccuracyの値は良くなる**のですが**予測値の正解率は向上しませんでした**。
特徴量をそのまま使うよりも生存率の違いなどから階級分けしたほうが過学習気味にならず適正なモデル作成ができるのではないかと感じました。

ご意見、ご指摘などがございましたらコメント、編集リクエストをしていただけるとありがたいです。

#参考にさせていただいたサイト、書籍
[KaggleチュートリアルTitanicで上位2%以内に入るノウハウ](https://qiita.com/jun40vn/items/d8a1f71fae680589e05c)
[pyhaya’s diary](https://pyhaya.hatenablog.com/entry/2018/11/26/232134)
[Titanic [0.82] - [0.83]](https://www.kaggle.com/konstantinmasich/titanic-0-82-0-83)
[Kaggleで勝つデータ分析の技術](https://www.amazon.co.jp/Kaggle%E3%81%A7%E5%8B%9D%E3%81%A4%E3%83%87%E3%83%BC%E3%82%BF%E5%88%86%E6%9E%90%E3%81%AE%E6%8A%80%E8%A1%93-%E9%96%80%E8%84%87-%E5%A4%A7%E8%BC%94-ebook/dp/B07YTDBC3Z)