はじめに
SotaになっていたSignateのApple引越社 需要予測のコンペに挑戦したので、その記録を書きます。
今回はモデル作成にPyCaretというAutoMLライブラリを利用した。PyCaretは、複数のモデルの構築と評価を少ないコードで簡単に実装できるAutoMLライブラリで、個人的に最近注目しているツール。モデル学習をPyCaretのみで実施してみて、どれだけの順位に組み込めるか試してみた。
先に結果を言うと、PyCaretのみで作成したモデルの評価値は、アップロード時点では150位くらいの結果となった。しかし最近見たら200位まで下がっていたので、現時点でもかなり挑戦している人は多そう。私自身時系列データ分析の経験は少ないのだが、データ数やタスクの難易度的に時系列タスクにちょうど良い課題だったと思うのと、PyCaretの手軽さを実感できた。
実行環境
- Python 3.9.0
- PyCaret 2.3.0
実施
1. PyCaretを用いたモデル作成
PyCaretを用いたモデル作成を実施する。
SignateからダウンロードしたCSVデータを./dataフォルダに格納し、データをプログラムにインポートする。
import pandas as pd train_df = pd.read_csv("./data/train.csv") test_df = pd.read_csv("./data/test.csv") sample_submit = pd.read_csv("data/sample_submit.csv", header=None) print(train_df.shape, test_df.shape, sample_submit.shape)
(2101, 6) (365, 5) (365, 2)
学習データが2101件、テストデータが365件あることが確認できる。今回のデータは特徴量が5つしかないので特徴量を作成する。
前処理
元データに含まれる特徴量が5つしかなく非常に少ないので、特徴量を増やす目的で前処理を行う。各値を結合したり統計量を取得することで特徴量作成を行った。特徴量の作成で具体的に行った手法は以下の通り。
- 日付を月、日、曜日に分割、年は削除
- 午前の料金と午後の料金の合計(
price_am_pm) - 料金に関わる3つの特徴量(
price_am_pm、price_am、price_pm)のラグ特徴量/リード特徴量
ラグ/リード特徴量は、それぞれ過去未来の地点のデータを利用する特徴量である。本タスクにおけるラグ特徴量の例としては、1週間前の日の周辺3日間の料金区分の統計量を取得しそれをその日の特徴量とする、といったものである。過去の傾向を利用できるので非常に強力な特徴量となりるが、例えば初日のデータが過去のデータが存在しないのでラグ特徴量を取得できない、というように、ラグで遡る日数以内のデータが存在しない端のデータは特徴量が正しく取得できないため、端のデータを取り除くか、許容してデータ分析を行う、といった注意が必要となる。
今回は一括で特徴量が作成できるように、以下のようなpreprocessクラスを作成した。
def preprocess(df, column_name, is_remove=False): ''' 文字列型の日付を年月日に分割する df : 対象のデータフレーム colume_name : 対象の行 is_remove : 日付を削除するかどうか ''' # 日付から時間を計算 df[column_name + '_month'] = pd.to_datetime(df[column_name]).dt.month df[column_name + '_day'] = pd.to_datetime(df[column_name]).dt.day df[column_name + '_dayofweek'] = pd.to_datetime(df[column_name]).dt.dayofweek df[column_name] = pd.to_datetime(df[column_name]).dt.strftime('%Y-%m-%d') # その日の料金の合計 df['price_am_pm'] = df['price_am'] + df['price_am'] # Lead特徴量 df = lead_encording(df, 'price_pm',shift_num=1) df = lead_encording(df, 'price_am',shift_num=1) df = lead_encording(df, 'price_am_pm',shift_num=1) df = lead_encording(df, 'price_pm',shift_num=2) df = lead_encording(df, 'price_am',shift_num=2) df = lead_encording(df, 'price_am_pm',shift_num=2) df = lead_encording(df, 'price_pm',shift_num=7) df = lead_encording(df, 'price_am',shift_num=7) df = lead_encording(df, 'price_am_pm',shift_num=7) # Lag特徴量 df = lag_encording(df, 'price_am',shift_num=7, window_size=3) df = lag_encording(df, 'price_pm',shift_num=7, window_size=3) df = lag_encording(df, 'price_am_pm',shift_num=7, window_size=3) df = lag_encording(df, 'price_am',shift_num=7, window_size=3) df = lag_encording(df, 'price_pm',shift_num=7, window_size=3) df = lag_encording(df, 'price_am_pm',shift_num=7, window_size=3) df = lag_encording(df, 'price_am',shift_num=14, window_size=7) df = lag_encording(df, 'price_pm',shift_num=14, window_size=7) df = lag_encording(df, 'price_am_pm',shift_num=14, window_size=7) df = lag_encording(df, 'price_am',shift_num=28, window_size=7) df = lag_encording(df, 'price_pm',shift_num=28, window_size=7) df = lag_encording(df, 'price_am_pm',shift_num=28, window_size=7) if is_remove: df = df.drop(column_name, axis=1) return df def lag_encording(df, column_name, shift_num=1, window_size=3): ''' 過去地点とその周辺の統計情報を取得 ''' df[column_name + '_lag_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num) df[column_name + '_lag_mean_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).mean() df[column_name + '_lag_sum_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).sum() df[column_name + '_lag_max_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).max() df[column_name + '_lag_min_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).min() df[column_name + '_lag_median_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).median() df[column_name + '_lag_std_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).std() return df def lead_encording(df, column_name, shift_num=1, window_size=3): ''' 未来地点とその周辺の統計情報を取得 ''' df[column_name + '_lead_'+ str(shift_num) ] = df[column_name].shift(-1 * shift_num) df[column_name + '_lead_mean_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).mean() df[column_name + '_lead_sum_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).sum() df[column_name + '_lead_max_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).max() df[column_name + '_lead_min_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).min() df[column_name + '_lead_median_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).median() df[column_name + '_lead_std_'+ str(shift_num) + '_' + str(window_size)] = df[column_name].shift(shift_num).rolling(window=window_size).std() return df
上記のメソッドで、学習データとテストデータに対して前処理を実施する。また、ラグ特徴量を取得するにあたってデータが欠損した、先頭の34行を削除する。
train_x, train_y = train_df.drop("y", axis=1), train_df["y"] preprocess_all = pd.concat([train_x, test_df]) preprocessed_all_x = preprocess(preprocess_all, 'datetime') # 前処理を行なったデータを、学習データとテストデータに再分割 preprocessed_train_x, preprocessed_test_x = preprocessed_all_x[:train_x.shape[0]], preprocessed_all_x[train_x.shape[0]:] # ラグ特徴量により先頭の34行はnullのカラムが複数含まれるので、学習データから取り除く。 preprocessed_train_x = preprocessed_train_x[34:] preprocessed_train_y = train_y[34:]
前処理の結果を出力してみる。
preprocessed_train_x
学習データは2106→2067まで減少したが、136個の特徴量を生成することができた。
モデル構築
PyCaretで学習するために、追加でデータの整形を行う。PyCaretは量的データ・質的データを自動判別し、それに応じた前処理をしてくれるが、稀に量的データであっても質的データと判別されることがある。なので量的変数のカラムを明示的に渡すためにget_column_typeを用いて量的データのカラムを抽出し、それをPyCaretの引数に与える。
def get_column_type(df, type_names: []): target_column_name = [] for column_name in df: if df[column_name].dtype in type_names: target_column_name.append(column_name) del df, type_names return target_column_name # Pycaretで実行するようにデータの前処理 train_all = pd.concat([preprocessed_train_x, train_y], axis=1) numeric_columns = get_column_type(train_all, [int, float]) numeric_columns.remove('y')
PyCaretにデータをセットする。
from pycaret.regression import * # initialize setup s = setup(data = preprocessed_train_x, target = 'y', fold_strategy = 'timeseries', numeric_features = numeric_columns, silent=True, fold = 3, session_id = 123)
各引数の詳細は以下となる。
| 引数名 | 概要 |
|---|---|
| data | 対象のDataFrame |
| target | 目的変数。dataに指定したDataFrameに含まれるカラムを指定する必要がある。 |
| fold_storategy | 交差検定の手法。今回のタスクのように時系列データを利用する場合はtimeseriesを指定する |
| numeric_features | dataに与えたDataFrameの説明変数のうち、量的データのカラムを明示的に指定する場合、 量的データのカラムの配列を指定する。 |
| silent | setup() 実行前、確認するフェイズを挟む場合はFalseを指定する。 |
| fold | データの分割数 |
| session_id | 乱数シードを固定する場合に任意の数値を指定する |
上記の実行が完了すれば、PyCaretで学習するパイプラインにデータがセットされたこととなる。
はじめに、デフォルトの設定compare_models()で全てのモデルで探索をしてみる。sort='MAE'とすることで、MAEの評価値をキーに並び替えることができる。
compare_models(sort = 'MAE')
出力結果を見ると、特徴量が多いことが影響してか決定木ベースのモデルが性能がよいことが分かる。今回は、その中でも上位のモデルであるgdr(勾配ブースティング木)とrf(Random Forest)を対象にパラメータ調整してみる。
パラメータチューニング
PyCaretではパラメータチューニングも簡単に実施できる。tune_model()の引数にチューニングしたいモデルを渡し実行することで、自動でパラメータ調整とモデルの予測結果に対する平均と分散を表示してくれる。さらにoptimize=にチューニングの指標となる評価指標を指定するとその値を最小にするようにチューニングが行われるので、今回はタスクに合わせてoptimize='MAE'とする。
# GDBのパラメータチューニング gbr = create_model('gbr') gbr = tune_model(gbr, optimize='MAE')
ランダムフォレストもパラメータチューニング。
# RandomForestのパラメータチューニング rf = create_model('rf') rf = tune_model(rf, optimize='MAE')
パラメータチューニングを行うことで、ランダムフォレストでもGDBでもMAEのスコアを向上させることができた。今までOptunaとかで頑張ってパラーメター調整していた手間が、だいぶ削減されたのは素晴らしい。
ブレンド
今回はGDBとランダムフォレストのモデルの両方の結果を利用するので、モデルのブレンドを実行する。モデルのブレンドも今までと同様に、ブレンドしたいモデルと評価指標を指定するのみでブレンドが実行できる。
# 2つのモデルをブレンド blender = blend_models(estimator_list = [rf, gbr], optimize='MAE')
最終的なモデルを固定し、パラメータの詳細値を表示する。
finalized_model = finalize_model(blender) finalized_model.get_params
get_paramsにより最終的なモデルのパラメータを取得できる。
<bound method _BaseHeterogeneousEnsemble.get_params of VotingRegressor(estimators=[('rf', RandomForestRegressor(bootstrap=True, ccp_alpha=0.0, criterion='mse', max_depth=9, max_features=1.0, max_leaf_nodes=None, max_samples=None, min_impurity_decrease=0.1, min_impurity_split=None, min_samples_leaf=4, min_samples_split=7, min_weight_fraction_leaf=0.0, n_estimators=100, n_jobs=-1, oob_score=False, random_state=123, verbo... max_features=1.0, max_leaf_nodes=None, min_impurity_decrease=0.02, min_impurity_split=None, min_samples_leaf=5, min_samples_split=5, min_weight_fraction_leaf=0.0, n_estimators=230, n_iter_no_change=None, presort='deprecated', random_state=123, subsample=0.85, tol=0.0001, validation_fraction=0.1, verbose=0, warm_start=False))], n_jobs=-1, verbose=False, weights=None)>
可視化
plot_modelを実行することで、学習データと検証データに対するR2係数やその分布を取得できる。
plot_model(finalized_model)
分布を見ると学習データ検証データともにR2係数が0.97前後の値が計算されているので、非常に精度の高いモデルであることが確認できる。
ブレンドしたモデルのMAEの値と本結果から、テストデータに対して悪くても6前後のMAEを獲得できそうなので、このモデルを最終的なモデルとして提出してみる。
テストデータの予測と提出データの作成
テストデータの予測はpredict_modelで実施する。pycaretで作成したモデルと予測したいデータを入力することで簡単に予測データを作成できる。
# テストデータを予測 proba = predict_model(finalized_model, data=preprocessed_test_x) proba.tail(100)
最後にデータ整形して提出データを作成する。
# データフレームの整形 submit_df = pd.DataFrame({"Y": proba['Label']}) submit_df.index = sample_submit[0] import os from datetime import datetime # csvに保存 save_folder = "results" if not os.path.exists(save_folder): os.makedirs(save_folder) submit_df.to_csv("{}/submit_{}.csv".format(save_folder, datetime.now().strftime("%Y-%m-%d-%H%M%S")),header=False, index=True)
提出したところMAEは0.982となっており、検証した結果とは大きく異なっていた。Apple引越し社のyの値をプロットしてみると年ごとに上昇傾向にあり、一般的な決定木ベースではこの上昇の傾向を獲得できないのではと考えられる。
以下のブログにあるように、時系列的な傾向を取得するためにProphetを利用した時系列モデリングとかの出力結果を組み合わせると性能向上につながるかもしれない。
終わりに
AutoMLライブラリであるPyCaretを利用することでモデル構築の手間を削減し、特徴量作成に力を入れることができた。
本気で上位を狙うのであればそれぞれのライブラリを利用するのがベターだと思われるが、さくっと成果を見たかったり、タスクのベースラインを作成する際にはPyCaretを利用するのが便利そう。
次回はもう少し高い順位を狙いたい。
