kaggle titanic 　xgboostを用いたデータ分析の流れ


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/863646/c070224f-27e3-8743-a5cb-9fe5171e9bbd.png)
機械学習初心者がkaggleのtitanicコンペのデータに対し、xgboostを使用したデータ分析にトライしたので備忘録として。
今回はtitanicのドメイン知識に関して深く踏み込まず、**xgboostの性能に頼った比較的シンプルな分析**を行っています。



#主なバージョン 
・python3.8.11
・xgboost 1.3.3
・scikit-learn 0.24.2

#やったこと
結果を予測する上で下記の内容を行ないました。


・特徴量に関する分析　(データの可視化)
・データ前処理（特徴量の削除、カテゴリ変数の変換、欠損値処理）
・xgboost分析  (モデル評価、パラメータ調整）

#特徴量の分析(データの可視化)
まずは、データを可視化し、自分なりのデータ分析をしてみます。その分析を基にデータ前処理やモデルによる分析を行う流れです。

データのインプット、テーブルデータの可視化

```py
import numpy as np
import pandas as pd
data =  pd.read_csv('train.csv')
data.head()
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/863646/3a6d462d-6948-74e7-a658-167274fb60bd.png)
12種類の特徴量があることが分かります。
またそれぞれの特徴量(カラム)の意味について記載しておきます。

・PassengerId：乗客のID
**・Survived：乗客の生死　（0=死亡、1=生存）　←これを予測したい**
・Name：乗客の名前
・Sex：性別
・SibSp：同乗している兄弟・配偶者の数
・Parch：同乗している親・子供の数
・Ticket：チケット番号
・Fare：料金
・Cabin：客室番号
・Embarked：出港地（タイタニック号に乗船した港）

ただの連番である乗客のIDや名前, ticketについては生存率と関係ないのでは？と推測できそうです。

また各カラムの欠損値情報等を見るため、infoメソッドで情報を可視化します。

```
data.info()
```


```
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 891 entries, 0 to 890
Data columns (total 12 columns):
 #   Column       Non-Null Count  Dtype  
---  ------       --------------  -----  
 0   PassengerId  891 non-null    int64  
 1   Survived     891 non-null    int64  
 2   Pclass       891 non-null    int64  
 3   Name         891 non-null    object 
 4   Sex          891 non-null    object 
 5   Age          714 non-null    float64
 6   SibSp        891 non-null    int64  
 7   Parch        891 non-null    int64  
 8   Ticket       891 non-null    object 
 9   Fare         891 non-null    float64
 10  Cabin        204 non-null    object 
 11  Embarked     889 non-null    object 
dtypes: float64(2), int64(5), object(5)
memory usage: 83.7+ KB
```
'Age'、'Cabin'、'Embarked'に欠損値があることが分かりました。
特に'Cabin'については欠損量の割合が多く、情報量が乏しいと判断します。
二番目に欠損値の多い'Age'については後ほど補完します。

#特徴量の削除
ここまでの簡易的な分析により、'PassengerId','Name', 'Ticket'については生存率と関係なさそうだと推定しました。また'Cabin'についても欠損値の割合が大きく、値の扱いも難しそうであるため、今回はこれらの特徴量を削除してモデルに学習させます。
（※'Name'や'Cabin'についても、上手く使えばスコアを伸ばすことができるようです。例えば'Cabin'については欠損していること自体を情報と見なし、カテゴリ変数として扱うなど。）

'PassengerId'、'Name'、'Cabin', 'Ticket'の列を削除

```
data.drop(['Cabin', 'PassengerId', 'Name', 'Ticket'], axis=1,  inplace=True)
data.head()
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/863646/bec0a525-e78c-79b7-74cd-c0384f3437c2.png)
これだけでも割と見やすくなりました。
#カテゴリ変数の変換
次にカテゴリ変数となっている特徴量の'Sex'列と'Embarked'列に対して、One-Hot-encodingを行います。


Sex'列と'Embarked'列に対するOne-Hot-Encoding

```py
# 'Sex'列を[0,1]に変換 (要素が二つの列は楽)
from sklearn.preprocessing import LabelEncoder
le = LabelEncoder()
data['Sex'] = le.fit_transform(data['Sex'])

# 'Embarked'列を[0,1]に変換(要素が三つあるため二列になる)
dummies = pd.get_dummies(data.Embarked)
data[['A','B','C']] =dummies
data.drop(['Embarked', 'C'], axis =1, inplace =True)

```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/863646/d3bf1749-feec-2c76-f06b-0092494cf21c.png)
'Sex'、'Embarked'列が二列の0,1 で構成された列に変換されました。

#欠損値補完
'Age'列の欠損値の補完するのですが、補完の方法としては、平均値や中央値などの代表値で埋める方法や、欠損値を他の変数から予測する方法などがあります(xgboostであれば、欠損値のまま扱うことも可能)。

補完方法を考えるため、先ほどまでに加工したテーブルデータをもう一度分析してみます。
今回は各行の相関係数をみてみます。

相関係数の可視化

```py
import seaborn as sns
corr = data.corr()
sns.heatmap( corr, square =True,annot =True)
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/863646/d39b652f-81fb-b3fe-d15b-816541851ef9.png)

出力されたヒートマップの'Age'行に着目してみると、'Pclass'と負の相関が強いことが分かりました。
ついでに、箱ヒゲ図でも見てみます。

```
sns.boxplot(x='Pclass', y ='Age', data=data)
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/863646/29061c13-2fe5-3e91-3ae7-97c0c4263003.png)
Pclass1の人は年齢が高く、Pclass3の人は若い傾向がありそうです。
この情報をもとに'Age'列の補完を行います。
具体的には、各Pclass毎に所属する人の'Age'中央値を取り、'Age'の欠損値を補完するという方法です。

'Age'欠損値補完

```py
class1 = data.query('Pclass==1')['Age'].median()
class2 = data.query('Pclass==2')['Age'].median()
class3 = data.query('Pclass==3')['Age'].median()

class1_index = data[(data['Pclass']==1) & data['Age'].isnull()].index
class2_index = data[(data['Pclass']==2) & data['Age'].isnull()].index
class3_index = data[(data['Pclass']==3) & data['Age'].isnull()].index

data.loc[class1_index, 'Age']=class1
data.loc[class2_index, 'Age']= class2
data.loc[class3_index,'Age'] =class3
```


**これで前処理は終わりました。**
最後に、データを目的変数である'Survived'列とその他のデータに分けておきます。

```
X = data.drop('Survived', axis =1).values
y = data['Survived'].values
```

これで、前処理まで終わったので、xgboostモデルによる分析に入ります。
(※xgboostは決定木ベースのモデルなので数値の標準化を行っていません。)

#xgboostモデルの作成、評価、パラメータ探索

xgboostモデルで学習していきます。
(※xgboostモデルの実装方法には、Learning APIと　Scikit-Learn APIの二種類があります。自分がScikit-Learnのライブラリに慣れていることもあり、今回はScikit-Learn APIで実装します。)

モデルの実装

```py
from xgboost import XGBClassifier

classifier = XGBClassifier(use_label_encoder=False,
                           eval_metric = 'error', booster='gbtree',n_estimators=100)
```
デフォルトのままで実装してます。
use_label_encoder=Falseとeval_metric='error'については書かなかったら注意文がでたので、一応記載しました（書かなくても結果自体は変わらなかったです）。

**モデル評価にはStatified kfold(10回)を使いました。**

```py
from sklearn.model_selection import cross_val_score, StratifiedKFold
import statistics

scores = cross_val_score(estimator =classifier,X=X,y=y, cv=StratifiedKFold(n_splits=10, random_state=1, shuffle=True))
print(statistics.mean(scores))
print(scores)
```
```
0.8137078651685393
[0.8        0.76404494 0.82022472 0.8988764  0.78651685 0.79775281
 0.79775281 0.7752809  0.88764045 0.80898876]
```
結果、平均0.81の評価結果が得られました！
デフォルトでも結構高い値が出るんですね。
あと10回の評価でかなり値がバラつくものなんだなぁという印象ですね。


**GridSearchを使った最適パラメータ探索もやってみます。**
今回はn_estimatorsのみでやってみます。
n_estimatorsはブースティングを行う回数であり、少なすぎたら値の予測精度が下がってしまい、回数を増やしすぎても過学習してしまいますので、最適値を探してみます。

```
%%time
from sklearn.model_selection import GridSearchCV

estimator_range =list(range(0,100,10))
param_grid=[{'n_estimators' : estimator_range}]

gs = GridSearchCV(estimator = classifier,
                 param_grid = param_grid, scoring ='accuracy',
                 cv=StratifiedKFold(n_splits=10, random_state=1, shuffle=True), n_jobs =-1)
gs.fit(X, y)
```
これで探索結果がgs に保管されました。今回探索には、13秒くらいでしたが、探索対象のパラメータを増やしていくと計算量が倍倍に増えるのですごい時間がかかることになりますね。
探索するパラメータのチョイス仕方にも工夫が要りそうですね。

結果の可視化

```
gs.best_params_,gs.best_score_
```
```
({'n_estimators': 10}, 0.8260799001248438)
```
n_estimatorsパラメータが10の時に正答率が最大値を取り、0.826ぐらいになるという結果です。
n_estimatorのデフォルト値100と比較するとかなり少ないですね。
今回の場合は100回もブースティングすると過学習になってしまうのかもしれません。


今回の分析で、kaggleのtestデータに予測して提出してみた結果
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/863646/7e7f3eab-a511-10b0-3062-f68eb2f8231e.png)
でした。

簡易的な分析でしたが、一連の流れはなんとなく経験できました。
ここからさらに特徴量の扱いを工夫したり、他のモデルを試したり、他のパラメータを探索したりしてスコアを上げていくのだと思います。