【Python】XGBoostの実装方法と特徴量重要度

Python
Image by Gerd Altmann from Pixabay
スポンサーリンク

XGBoostは簡単に実装できるわりに精度が高いので、ベースラインモデルとして優秀です。

params = {
    "objective" : "reg:squarederror", #回帰モデルの場合
    "eval_metric" : "rmse"
}

model = xgb.train(
    params = params,
    dtrain = dtrain,
    evals = [(dtrain, "train"), (dvalid, "valid")],
    num_boost_round = 100
)

似たモデルにLightGBMがあります。こっちの方がメモリ消費が小さいです。

回帰モデル

回帰とは連続した値を予測する問題のことです。

Q:この住宅の価格はいくらですか
A:〇〇万円です

load_boston

予測するデータが手元にない方もいると思います。

今回はPythonに標準搭載されているデータセットload_bostonを使いましょう。

import pandas as pd
from sklearn.datasets import load_boston

data = load_boston()
df = pd.DataFrame(data["data"], columns = data["feature_names"])
df["target"] = data["target"]
print(df.shape)
df.head()
pandas:表計算ライブラリ
load_boston:住宅価格データセット

データフレームにload_bostonのデータを入れました。

”CRIM”~”LSTAT”までの13列から”target”を予測します。

データ分割

モデルに入れる学習データと、性能を確かめる評価データに分けます。

from sklearn.model_selection import KFold
kf = KFold(n_splits = 5, shuffle = True)
df["fold"] = -1
for fold, (train_idx, valid_idx) in enumerate(kf.split(df)):
    df.loc[valid_idx, "fold"] = fold
print(df["fold"].value_counts())

KFoldはデータを等分してくれるライブラリです。

n_splits = 5にしたので5等分されます。

fold = 0
train = df.loc[df["fold"] != fold].copy()
valid = df.loc[df["fold"] == fold].copy()
print(train.shape, valid.shape)

今回は0番目のデータを評価用にしました。

次に特徴量と目的変数に分けましょう。

feat_cols = train.drop(columns = ["fold", "target"]).columns.tolist()
print(feat_cols)

X_train = train[feat_cols]
X_valid = valid[feat_cols]
y_train = train["target"]
y_valid = valid["target"]

“fold”と”target”は特徴としてふさわしくないので外しました。

Xが予測に使う特徴量で、yが予測したい目的変数です。

モデル作成

必要な手順は以下の通り。

①ライブラリのインポート
②DMatrix形式に変換
③パラメータ設定
④trainで学習
import xgboost as xgb

xgboostをインポートしましょう。xgbと略すことが多いです。

次にデータをXGBoost専用のDMatrix形式に変換します。

dtrain = xgb.DMatrix(X_train, y_train)
dvalid = xgb.DMatrix(X_valid, y_valid)

DMatrixにtrainとvalidのそれぞれの特徴量と目的変数を入れました。

これで変換は終わりです。

パラメータ(学習条件)を辞書型データで設定しましょう。

params = {
    "objective" : "reg:squarederror",
    "eval_metric" : "rmse"
}
objective:モデルの種類(回帰や分類など)
eval_metric:評価指標(rmseやaucなど)

回帰モデルを作りたいなら”objective”は”reg:squarederror”にします。

“reg”は回帰の”regression”で”squarederror”は二乗誤差の意味です。

“eval_metric”では、「何がどうなったらよいモデルなのか」を示す評価指標を指定します。

“rmse”は平均二乗誤差の平方根です。

パラメータの設定は細かいので公式ドキュメントを見るといいですよ!

とりあえず上記の設定にしておけば最低限動いてくれます。

では学習を始めましょう。↓

model = xgb.train(
    params = params,
    dtrain = dtrain,
    evals = [(dtrain, "train"), (dvalid, "valid")],
    num_boost_round = 100
)
params:さっき作ったパラメータ
dtrain:学習データ
evals:”eval_metric”で評価するデータ
num_boost_round:モデルを作り直す回数

学習結果が出力されたら成功です。

XGBoostは決定木モデルを”num_boost_round”の数だけ作り直しています。

“rmse”(誤差)が徐々に減少していれば改善が進んでいるということですね。

予測結果

predictで予測結果を出力できます。

pred = model.predict(xgb.DMatrix(X_valid))
print(pred)

評価データで予測してみましょう。ここでもDMatrixの変換が必要です。

printで出力してもパッとしないので可視化しましょう。↓

import matplotlib.pyplot as plt
import numpy as np
plt.scatter(y_valid, pred, alpha = 0.5)
plt.plot(np.linspace(0, 60, 100), np.linspace(0, 60, 100), "red")
plt.show()

scatterで横軸を正解データ、縦軸を予測値として散布図を描きました。

赤線がy = xの線です。この線上にあるならピッタリ予測できていることになります。

最後に平均二乗誤差の平方根(RMSE)を計算しましょう。

from sklearn.metrics import mean_squared_error
mse = mean_squared_error(y_valid, pred)
rmse = np.sqrt(mse)
print(rmse)

“mean_squared_error”に渡すと計算してくれます。

numpyで平方根をとるとRMSEです。学習時に出力された”valid-rmse”に近い値のはず。

この結果がモデルの誤差ととらえておけばOKです。

学習状況可視化

どんな感じでモデルが成長しているか見てみたい。

というわけで学習過程を可視化してみましょう。

results_dict = {}
model = xgb.train(
    params = params,
    dtrain = dtrain,
    evals = [(dtrain, "train"), (dvalid, "valid")],
    num_boost_round = 100,
    evals_result = results_dict
)

“results_dict”という辞書型データを作って”evals_result”に渡しました。

これで学習過程が保存されます。

plt.plot(results_dict["train"]["rmse"], color = "red", label = "train")
plt.plot(results_dict["valid"]["rmse"], color = "blue", label = "valid")
plt.legend()
plt.show()

こんな感じでRMSE(誤差)が下がっていく様子が見て取れます。

赤が学習データで青が評価データですね。

EarlySroppingRounds

先ほど可視化した学習状況を見ると、10回目でほとんど頭打ちになっていました。

なので後半のほとんどは改善もなく時間の無駄になっています。

そこで”early_stopping_rounds”を設定しましょう。

results_dict = {}
model = xgb.train(
    params = params,
    dtrain = dtrain,
    evals = [(dtrain, "train"), (dvalid, "valid")],
    num_boost_round = 100,
    early_stopping_rounds = 10,
    evals_result = results_dict
)

上記コードを実行すると、評価データでの性能が10回変わらない(ほとんど同じ)場合に学習を自動で止めてくれます。

この結果を可視化すると以下の通り。

20ちょっとで止まってくれていますね。

この他にも色々なパラメータ(機能)があるので公式ドキュメントを見るといいですよ。

特徴量重要度

どの特徴が予測に効いているの?

こんな悩みに答えてくれる機能が”plot_importance”です。

xgb.plot_importance(model)
plt.show()

”CRIM”が重要そうですね。

ただしここで言う重要度とは”決定木の分岐に何回登場したか”を意味します。

つまり定量的にどれだけ予測結果に貢献したか示すわけではありません。

というわけで”importance_type”を指定しましょう!

“gain”にすると目的変数(予測したい値)への寄与度を計算してくれます。

xgb.plot_importance(model, importance_type = "gain")
plt.show()

明らかに結果が変わりましたね。どうやら”LSTAT”が重要のようです。

ちなみにデフォルトではimportance_type = “weight”になっています。

二値分類モデル

二値分類とは、二択の予測をする問題のことです。

Q:明日雨が降る?
A:降る or 降らない

パラメータ設定以外はほとんど回帰モデルと同じなので、割愛しながら解説します。

load_breast_cancer

import pandas as pd
from sklearn.datasets import load_breast_cancer
data = load_breast_cancer()
df = pd.DataFrame(data["data"], columns = data["feature_names"])
df["target"] = data["target"]
print(df.shape)
df.head()

load_breast_cancerをインポートしました。

病院での検査結果を示すデータで、0が陰性、1が陽性です。

データ分割

from sklearn.model_selection import StratifiedKFold
skf = StratifiedKFold(n_splits = 5, shuffle = True)
df["fold"] = -1
for fold, (train_idx, valid_idx) in enumerate(skf.split(X = df, y = df["target"])):
    df.loc[valid_idx, "fold"] = fold
print(df["fold"].value_counts())

分類モデルでのデータ分割ではStratifiedKFoldを使いましょう。

これは予測したいラベルの割合を保ちながら分割してくれます。

例えば”学習データにラベル1が含まれない”なんてことになったらうまく学習できませんよね。

fold = 0
train = df.loc[df["fold"] != fold].copy()
valid = df.loc[df["fold"] == fold].copy()
print(train.shape, valid.shape)

feat_cols = train.drop(columns = ["fold", "target"]).columns.tolist()
print(feat_cols)

X_train = train[feat_cols]
X_valid = valid[feat_cols]
y_train = train["target"]
y_valid = valid["target"]

0番目のデータを評価用にしました。

モデル作成

回帰との違いはパラメータ設定です。

import xgboost as xgb

dtrain = xgb.DMatrix(X_train, y_train)
dvalid = xgb.DMatrix(X_valid, y_valid)

params = {
    "objective" : "binary:logistic",
    "eval_metric" : "logloss"
}

results_dict = {}
model = xgb.train(
    params = params,
    dtrain = dtrain,
    evals = [(dtrain, "train"), (dvalid, "valid")],
    num_boost_round = 100,
    early_stopping_rounds = 10,
    evals_result = results_dict
)
objective:binary:logistic
eval_metric:logloss

“objective”は二値分類の場合”binary:logistic”にします。

“binary”は2択を意味します。バイナリーオプションとか言いますよね。

“eval_metric”(評価指標)は”logloss”です。

予測結果が外れたとき、その予測確立が大きいほどペナルティを与えます。

学習過程を可視化してみましょう。↓

import matplotlib.pyplot as plt
plt.plot(results_dict["train"]["logloss"], color = "red", label = "train")
plt.plot(results_dict["valid"]["logloss"], color = "blue", label = "valid")
plt.legend()
plt.show()

logloss(誤差)が下がっていますね!

予測結果

pred = model.predict(xgb.DMatrix(X_valid))
print(pred)

predで予測結果を出力すると、0or1のラベルではなく、”ラベル1である確率”が得られます。

import numpy as np
label = np.where(pred > 0.5, 1, 0)
print(label)

numpyのwhereを使って確率0.5以上のデータをラベル1にしました。

print(label)
print("=" * 100)
print(y_valid.values)

y_validが正解データなので比較してみましょう。だいたいあってますね。

次に正解率を計算してみます。

from sklearn.metrics import accuracy_score
acc = accuracy_score(y_valid, label)
print(acc)

“accuracy_score”で正解率を計算できます。だいたい96%くらいになるはずです。

ただ、正解率はあまり良い指標ではありません。

例えばデータの99%がラベル0だった場合を考えます。

このとき”全部ラベル0にする”というふざけたモデルでも正解率99%ですよね。

なので、偽陽性や偽陰性なども考慮したAUCを計算しましょう。

from sklearn.metrics import roc_auc_score
auc = roc_auc_score(y_valid, pred)
print(auc)

”roc_auc_score”に正解データと予測確率を渡します。

AUCも0~1の値で1に近いほど優秀です。だいたい0.98くらいかと思います。

多クラス分類モデル

他クラス分類は、複数の選択肢がある問題です。

Q:明日の天気は?
A:晴れ or 曇り or 雨

これもパラメータ設定以外はこれまでと同じです。

load_iris

import pandas as pd
from sklearn.datasets import load_iris
data = load_iris()
df = pd.DataFrame(data["data"], columns = data["feature_names"])
df["target"] = data["target"]
print(df.shape)
df.head()

load_irisは花の種類を0,1,2の3種類で分類するデータセットです。

データ分割

Xfrom sklearn.model_selection import StratifiedKFold
skf = StratifiedKFold(n_splits = 3, shuffle = True)
df["fold"] = -1
for fold, (train_idx, valid_idx) in enumerate(skf.split(X = df, y = df["target"])):
    df.loc[valid_idx, "fold"] = fold
print(df["fold"].value_counts())

データ分割はStratifiedKFoldにしましょう。

データ数(行数)が少ないのでn_splits = 3にしておきました。

fold = 0
train = df.loc[df["fold"] != fold].copy()
valid = df.loc[df["fold"] == fold].copy()
print(train.shape, valid.shape)

feat_cols = train.drop(columns = ["fold", "target"]).columns.tolist()
print(feat_cols)

X_train = train[feat_cols]
X_valid = valid[feat_cols]
y_train = train["target"]
y_valid = valid["target"]

0番のデータを評価用にしました。

モデル作成

パラメータに注意しましょう!

import xgboost as xgb

dtrain = xgb.DMatrix(X_train, y_train)
dvalid = xgb.DMatrix(X_valid, y_valid)

params = {
    "objective" : "multi:softprob",
    "eval_metric" : "mlogloss",
    "num_class" : 3
}

results_dict = {}
model = xgb.train(
    params = params,
    dtrain = dtrain,
    evals = [(dtrain, "train"), (dvalid, "valid")],
    num_boost_round = 100,
    early_stopping_rounds = 10,
    evals_result = results_dict
)
ocjective:multi:softprob
eval_metric:mlogloss
num_class:3(ラベルの種類数)

“objective”を”multi:softprob”にしましょう。

各ラベルに該当する確率を計算してくれます。

“eval_metric”は”mlogloss”にしました。mが複数を意味します。

次に”num_class”にラベル数を渡しましょう。今回は3ですね。

import matplotlib.pyplot as plt
plt.plot(results_dict["train"]["mlogloss"], color = "red", label = "train")
plt.plot(results_dict["valid"]["mlogloss"], color = "blue", label = "valid")
plt.legend()
plt.show()

学習過程はこれまでと同じように徐々にlossが小さくなります。

予測結果

predで予測すると”各ラベルに該当する確率”を出力します。

pred = model.predict(xgb.DMatrix(X_valid))
print(pred.shape)
print(pred)

実行結果は以下の通り。

「データの行数 x 3」の配列で、左からラベル0,1,2である確率です。

label = pred.argmax(axis = 1)
print(label)
print("=" * 100)
print(y_valid.values)

argmaxを使うと”値が最も大きいデータの位置”が出力されます。

axis = 1にすると列方向において最も大きいデータの位置になります。

from sklearn.metrics import accuracy_score
acc = accuracy_score(y_valid, label)
print(acc)

正解率を計算しました。だいたい95%くらいになると思います。

他クラス分類ではAUCを計算できません。

なので適合率(Precision)と再現率(Recall)を計算しておきましょう。

from sklearn.metrics import precision_score, recall_score
precision = precision_score(y_valid.values, label, average = "micro")
recall = recall_score(y_valid.values, label, average = "micro")
print(precision, recall)
ラベル0において、
適合率:ラベル0と予測したうち、実際に0だった割合
再現率:実際にラベル0であるデータのうち、0と予測された割合

これを各ラベルにおいて計算しています。高いほど優秀です。

どちらも95%くらいになるはずです。

まとめ

今回はXGBoostの実装方法を解説しました。

パラメータの設定だけ問題によって変わるので気を付けましょう!

XGBoostが使えるようになったらLightGBMも勉強するといいですよ。

コメント

タイトルとURLをコピーしました