時系列データの取り扱いを学ぶ。（Store Sales - Time Series Forecasting）

背景

Kaggleに取り組んでいて、簡単なデータの前処理やバリデーションを用いて予測をしてみたが精度が悪かった。

精度を上げるために何から始めるべきか、検討もつかなかったので
できる人のコードを模写して流れを理解しようと思った。

そのためコードは基本的には模写である。

扱うデータ

今回は下記の小売店の売上予測コンペに取り組む。
https://www.kaggle.com/competitions/store-sales-time-series-forecasting/data

コード

下記に流れを含めてコードを書いていく。

準備

ファイル名を取得する。

import os
for dirname, _, filenames in os.walk('kaggle_data'):
    for filename in filenames:
        print(os.path.join(dirname, filename))

データフレームとして読み込む。

import pandas as pd

df_train = pd.read_csv('kaggle_data/train.csv')
df_holidays_events = pd.read_csv('kaggle_data/holidays_events.csv')
df_oil = pd.read_csv('kaggle_data/oil.csv')
df_stores = pd.read_csv('kaggle_data/stores.csv')
df_transactions = pd.read_csv('kaggle_data/transactions.csv')
df_test = pd.read_csv('kaggle_data/test.csv')
df_sample_submission = pd.read_csv('kaggle_data/sample_submission.csv')

トレーニングデータと各店舗のデータを'store_nbr'をkeyにマージする。
並び替える。("date"で並び替えて→"family"で並び替えて→"store_nbr"で並び替える。）
データ型を文字列に変換する。

train_merged = pd.merge(df_train, df_stores, on ='store_nbr')
train_merged = train_merged.sort_values(["store_nbr","family","date"])
train_merged = train_merged.astype({"store_nbr":'str', "family":'str', "city":'str',
                          "state":'str', "type":'str', "cluster":'str'})

テストデータも並び替える。
'onpromotion'列を削除した理由は不明。

df_test_dropped = df_test.drop(['onpromotion'], axis=1)
df_test_sorted = df_test_dropped.sort_values(by=['store_nbr','family'])

display(df_test_sorted.head())

時系列データを作成する。

Dartsという時系列に特化したライブラリを使う。
https://unit8co.github.io/darts/

調べてみるとそれなりに使われているようだが、Prophetが最も使われているようだ。
違いは今後、使っていくうちに理解しようと思う。
https://facebook.github.io/prophet/

色々データを触ってみた結果、製品群が各売上に大きな影響を及ぼすことがわかったようです。
（このあたりのデータの取っ掛かりは可視化などを通して、発見する必要がある。要勉強。）

dartsは配列型のため、numpyを使う。

import numpy as np
import darts
from darts import TimeSeries

製品群のリストの各要素ごとに辞書型でデータを保持する。
取り出し方はfor文で回す。

TimeSeries.from_group_dataframeで時系列データに変換する。
引数の役割は下記。

ーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーー

• time_col: 時系列データの時間情報を含む列を指定します。ここでは "date" が指定されており、売上データの時間情報がこの列に格納されていると想定されています。

• group_cols: 時系列データをグループ化するための列を指定します。ここでは ["store_nbr", "family"] となっており、store_nbr（店舗番号）と family（商品ファミリー）の値に基づいて時系列データがグループ化されます。つまり、各店舗と商品ファミリーの組み合わせごとに時系列データがグループ化されます。

• static_cols: 固定共変量（静的な属性情報）を指定します。ここでは ["city", "state", "type", "cluster"] が指定されており、これらの列が各時系列データの固定共変量として設定されます。これにより、各時系列データが店舗や商品ファミリーの属性情報を持つことができます。

• value_cols: 時系列データの値を含む列を指定します。ここでは "sales" が指定されており、各時系列データの売上データがこの列に格納されていると想定されています。

ーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーー

family_TS_dict = {}

for family in family_list:
  df_family = train_merged.loc[train_merged['family'] == family]

  list_of_TS_family = TimeSeries.from_group_dataframe(
                                df_family,
                                time_col="date",
                                group_cols=["store_nbr","family"],
                                static_cols=["city","state","type","cluster"],
                                value_cols="sales",
                                fill_missing_dates=True,
                                freq='D')
  for ts in list_of_TS_family:
            ts = ts.astype(np.float32)

  list_of_TS_family = sorted(list_of_TS_family, key=lambda ts: int(ts.static_covariates_values()[0,0]))
  family_TS_dict[family] = list_of_TS_family

正規化

Scaler→データのスケーリング
StaticCovariatesTransformer→ワンホットエンコーディングなど
MissingValuesFiller→欠損値処理
InvertibleMapper→対数変換など

from darts.dataprocessing import Pipeline
from darts.dataprocessing.transformers import Scaler, StaticCovariatesTransformer, MissingValuesFiller, InvertibleMapper
import sklearn

family_pipeline_dict = {}
family_TS_transformed_dict = {}

for key in family_TS_dict:
  train_filler = MissingValuesFiller(verbose=False, n_jobs=-1, name="Fill NAs")
  static_cov_transformer = StaticCovariatesTransformer(verbose=False, transformer_cat = sklearn.preprocessing.OneHotEncoder(), name="Encoder")
  log_transformer = InvertibleMapper(np.log1p, np.expm1, verbose=False, n_jobs=-1, name="Log-Transform")   
  train_scaler = Scaler(verbose=False, n_jobs=-1, name="Scaling")

  train_pipeline = Pipeline([train_filler,
                             static_cov_transformer,
                             log_transformer,
                             train_scaler])
     
  training_transformed = train_pipeline.fit_transform(family_TS_dict[key])
  family_pipeline_dict[key] = train_pipeline
  family_TS_transformed_dict[key] = training_transformed

時系列データから特徴量作成

一年の中での夏なのか、一週間の中での土日なのか、何月なのかなどで傾向が変わる可能性がある。

from darts.utils.timeseries_generation import datetime_attribute_timeseries

full_time_period = pd.date_range(start='2013-01-01', end='2017-08-31', freq='D')


year = datetime_attribute_timeseries(time_index = full_time_period, attribute="year")
month = datetime_attribute_timeseries(time_index = full_time_period, attribute="month")
day = datetime_attribute_timeseries(time_index = full_time_period, attribute="day")
dayofyear = datetime_attribute_timeseries(time_index = full_time_period, attribute="dayofyear")
weekday = datetime_attribute_timeseries(time_index = full_time_period, attribute="dayofweek")
weekofyear = datetime_attribute_timeseries(time_index = full_time_period, attribute="weekofyear")
timesteps = TimeSeries.from_times_and_values(times=full_time_period,
                                             values=np.arange(len(full_time_period)),
                                             columns=["linear_increase"])

time_cov = year.stack(month).stack(day).stack(dayofyear).stack(weekday).stack(weekofyear).stack(timesteps)
time_cov = time_cov.astype(np.float32)

もちろん正規化もする。

time_cov_scaler = Scaler(verbose=False, n_jobs=-1, name="Scaler")
time_cov_train, time_cov_val = time_cov.split_before(pd.Timestamp('20170816'))
time_cov_scaler.fit(time_cov_train)
time_cov_transformed = time_cov_scaler.transform(time_cov)

石油の特徴量

石油価格は事前にわかっている、かつ長期的な傾向を示唆するためにも
移動平均を使う。(MovingAverage)

from darts.models import MovingAverage
# Oil Price

oil = TimeSeries.from_dataframe(df_oil, 
                                time_col = 'date', 
                                value_cols = ['dcoilwtico'],
                                freq = 'D')

oil = oil.astype(np.float32)

# Transform
oil_filler = MissingValuesFiller(verbose=False, n_jobs=-1, name="Filler")
oil_scaler = Scaler(verbose=False, n_jobs=-1, name="Scaler")
oil_pipeline = Pipeline([oil_filler, oil_scaler])
oil_transformed = oil_pipeline.fit_transform(oil)

# Moving Averages for Oil Price
oil_moving_average_7 = MovingAverage(window=7)
oil_moving_average_28 = MovingAverage(window=28)

oil_moving_averages = []

ma_7 = oil_moving_average_7.filter(oil_transformed).astype(np.float32)
ma_7 = ma_7.with_columns_renamed(col_names=ma_7.components, col_names_new="oil_ma_7")
ma_28 = oil_moving_average_28.filter(oil_transformed).astype(np.float32)
ma_28 = ma_28.with_columns_renamed(col_names=ma_28.components, col_names_new="oil_ma_28")
oil_moving_averages = ma_7.stack(ma_28)

休日の特徴量

クリスマスや地震による休日、サッカーイベントなど売上に変化を起こす特徴量を入れる。

def holiday_list(df_stores):

    listofseries = []
    
    for i in range(0,len(df_stores)):
            
            df_holiday_dummies = pd.DataFrame(columns=['date'])
            df_holiday_dummies["date"] = df_holidays_events["date"]
            
            df_holiday_dummies["national_holiday"] = np.where(((df_holidays_events["type"] == "Holiday") & (df_holidays_events["locale"] == "National")), 1, 0)

            df_holiday_dummies["earthquake_relief"] = np.where(df_holidays_events['description'].str.contains('Terremoto Manabi'), 1, 0)

            df_holiday_dummies["christmas"] = np.where(df_holidays_events['description'].str.contains('Navidad'), 1, 0)

            df_holiday_dummies["football_event"] = np.where(df_holidays_events['description'].str.contains('futbol'), 1, 0)

            df_holiday_dummies["national_event"] = np.where(((df_holidays_events["type"] == "Event") & (df_holidays_events["locale"] == "National") & (~df_holidays_events['description'].str.contains('Terremoto Manabi')) & (~df_holidays_events['description'].str.contains('futbol'))), 1, 0)

            df_holiday_dummies["work_day"] = np.where((df_holidays_events["type"] == "Work Day"), 1, 0)

            df_holiday_dummies["local_holiday"] = np.where(((df_holidays_events["type"] == "Holiday") & ((df_holidays_events["locale_name"] == df_stores['state'][i]) | (df_holidays_events["locale_name"] == df_stores['city'][i]))), 1, 0)
                     
            listofseries.append(df_holiday_dummies)

    return listofseries

各ストアごとの祝日情報データフレームから、不要な0行を取り除き、各日付におけるイベントの最大値を保持する

def remove_0_and_duplicates(holiday_list):

    listofseries = []
    
    for i in range(0,len(holiday_list)):
            
            df_holiday_per_store = list_of_holidays_per_store[i].set_index('date')

            df_holiday_per_store = df_holiday_per_store.loc[~(df_holiday_per_store==0).all(axis=1)]
            
            df_holiday_per_store = df_holiday_per_store.groupby('date').agg({'national_holiday':'max', 'earthquake_relief':'max', 
                                   'christmas':'max', 'football_event':'max', 
                                   'national_event':'max', 'work_day':'max', 
                                   'local_holiday':'max'}).reset_index()

            listofseries.append(df_holiday_per_store)

    return listofseries

54店舗分を反映する。

def holiday_TS_list_54(holiday_list):

    listofseries = []
    
    for i in range(0,54):
            
            holidays_TS = TimeSeries.from_dataframe(list_of_holidays_per_store[i], 
                                        time_col = 'date',
                                        fill_missing_dates=True,
                                        fillna_value=0,
                                        freq='D')
            
            holidays_TS = holidays_TS.slice(pd.Timestamp('20130101'),pd.Timestamp('20170831'))
            holidays_TS = holidays_TS.astype(np.float32)
            listofseries.append(holidays_TS)

    return listofseries

定義した関数を呼び出す。

list_of_holidays_per_store = holiday_list(df_stores)
list_of_holidays_per_store = remove_0_and_duplicates(list_of_holidays_per_store)   
list_of_holidays_store = holiday_TS_list_54(list_of_holidays_per_store)

holidays_filler = MissingValuesFiller(verbose=False, n_jobs=-1, name="Filler")
holidays_scaler = Scaler(verbose=False, n_jobs=-1, name="Scaler")

holidays_pipeline = Pipeline([holidays_filler, holidays_scaler])
holidays_transformed = holidays_pipeline.fit_transform(list_of_holidays_store)

promotionの特徴量

製品群、店舗ごとにプロモーション情報を整理する。
時系列データの作成。

df_promotion = pd.concat([df_train, df_test], axis=0)
df_promotion = df_promotion.sort_values(["store_nbr","family","date"])
df_promotion.tail()

family_promotion_dict = {}

for family in family_list:
  df_family = df_promotion.loc[df_promotion['family'] == family]

  list_of_TS_promo = TimeSeries.from_group_dataframe(
                                df_family,
                                time_col="date",
                                group_cols=["store_nbr","family"],
                                value_cols="onpromotion",
                                fill_missing_dates=True,
                                freq='D')
  
  for ts in list_of_TS_promo:
    ts = ts.astype(np.float32)

  family_promotion_dict[family] = list_of_TS_promo

これまでのデータ処理で使ったものを同じように適用。
石油価格と同じでプロモーション費用も移動平均でトレンドを反映できるようにする。

• tqdm:プログレスバーで実行状況を可視化。
• MissingValuesFiller:欠損値処理
• Scaler:スケーリング
• Pipeline:欠損値処理とスケーリングを一元管理できる。（個別処理をしなくて良くなることが利点？　いまいち利点が理解できていない。）
• MovingAverage:移動平均でトレンドを反映する。

from tqdm import tqdm

promotion_transformed_dict = {}

for key in tqdm(family_promotion_dict):
  promo_filler = MissingValuesFiller(verbose=False, n_jobs=-1, name="Fill NAs")
  promo_scaler = Scaler(verbose=False, n_jobs=-1, name="Scaling")

  promo_pipeline = Pipeline([promo_filler,
                             promo_scaler])
  
  promotion_transformed = promo_pipeline.fit_transform(family_promotion_dict[key])
  
  # Moving Averages for Promotion Family Dictionaries
  promo_moving_average_7 = MovingAverage(window=7)
  promo_moving_average_28 = MovingAverage(window=28)

  promotion_covs = []

  for ts in promotion_transformed:
    ma_7 = promo_moving_average_7.filter(ts)
    ma_7 = TimeSeries.from_series(ma_7.pd_series())  
    ma_7 = ma_7.astype(np.float32)
    ma_7 = ma_7.with_columns_renamed(col_names=ma_7.components, col_names_new="promotion_ma_7")
    ma_28 = promo_moving_average_28.filter(ts)
    ma_28 = TimeSeries.from_series(ma_28.pd_series())  
    ma_28 = ma_28.astype(np.float32)
    ma_28 = ma_28.with_columns_renamed(col_names=ma_28.components, col_names_new="promotion_ma_28")
    promo_and_mas = ts.stack(ma_7).stack(ma_28)
    promotion_covs.append(promo_and_mas)

  promotion_transformed_dict[key] = promotion_covs

特徴量を一つの時系列データに結合

まずは石油価格と石油価格の移動平均を'general_covariates'に格納。

general_covariates = time_cov_transformed.stack(oil_transformed).stack(oil_moving_averages)

店舗ごとの祝日データと石油価格（'general_covariates'）を'store_covariates_future'に結合。

store_covariates_future = []

for store in range(0,len(store_list)):
  stacked_covariates = holidays_transformed[store].stack(general_covariates)  
  store_covariates_future.append(stacked_covariates)

各製品カテゴリ毎に店舗ごとの祝日データや石油価格（'store_covariates_future'）を、'future_covariates_dict'に辞書型で結合。

future_covariates_dict = {}

for key in tqdm(promotion_transformed_dict):

  promotion_family = promotion_transformed_dict[key]
  covariates_future = [promotion_family[i].stack(store_covariates_future[i]) for i in range(0,len(promotion_family))]

  future_covariates_dict[key] = covariates_future

transactionの特徴量

機械学習前の最後の特徴量調整になるので、流れを下記に記載する。

1. データフレームから時系列データを作成する。
   1. time_col:時間を表す列
   2. group_cols:時系列関係なく、カテゴリ毎に扱うデータ
   3. value_cols:時系列で表されているデータ
   4. fill_missing_dates:欠損している日付を埋める
   5. freq:時系列データの頻度（今回は日毎なので、'D'）
2. データを使えるようにする。
   1. floatに
   2. 店舗番号２５の20130102〜20170815ではないデータを消す。（店舗番号２５だけ20130101のデータが存在する。）
   3. datetime:日付型に変換
   4. timedelta:日付同士を計算
   5. MissingValuesFiller:欠損値処理
   6. Scaler:スケーリングで0から1の間に。
   7. Pipeline:欠損値処理とスケーリングを一元管理

df_transactions.sort_values(["store_nbr","date"], inplace=True)

TS_transactions_list = TimeSeries.from_group_dataframe(
                                df_transactions,
                                time_col="date",
                                group_cols=["store_nbr"],
                                value_cols="transactions",
                                fill_missing_dates=True,
                                freq='D')

transactions_list = []

for ts in TS_transactions_list:
            series = TimeSeries.from_series(ts.pd_series())
            series = series.astype(np.float32)
            transactions_list.append(series)

transactions_list[24] = transactions_list[24].slice(start_ts=pd.Timestamp('20130102'), end_ts=pd.Timestamp('20170815'))

from datetime import datetime, timedelta

transactions_list_full = []

for ts in transactions_list:
  if ts.start_time() > pd.Timestamp('20130101'):
    end_time = (ts.start_time() - timedelta(days=1))
    delta = end_time - pd.Timestamp('20130101')
    zero_series = TimeSeries.from_times_and_values(
                              times=pd.date_range(start=pd.Timestamp('20130101'), 
                              end=end_time, freq="D"),
                              values=np.zeros(delta.days+1))
    ts = zero_series.append(ts)
    ts = ts.with_columns_renamed(col_names=ts.components, col_names_new="transactions")
    transactions_list_full.append(ts)

transactions_filler = MissingValuesFiller(verbose=False, n_jobs=-1, name="Filler")
transactions_scaler = Scaler(verbose=False, n_jobs=-1, name="Scaler")

transactions_pipeline = Pipeline([transactions_filler, transactions_scaler])
transactions_transformed = transactions_pipeline.fit_transform(transactions_list_full)

機械学習の方針

色々試した結果、LightGBMが性能が良かったらしい。
やはり初手はこいつになるんだなと。

ちなみに試したのは下記。
後々、キャッチアップすることとする。

• NHiTSModel – score : 0.43265
• RNNModel (with LSTM layers) – score : 0.55443
• TFTModel – score : 0.43226
• ExponentialSmoothing – score : 0.37411

また製品群ごとに機械学習をするのが最も良いと考えたらしい。
（可視化やデータを弄ることを出来ていないので、ここのロジックも今後に譲る。）

＜コード手順＞
製品群ごとにデータを抽出する。
抽出したデータの特徴量を変数に入れる。

ハイパーパラメーターを設定する。
ラグ特徴量を使って予測する。

ーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーー

• lags:予測点に対して実績データを学習に利用する際の期間（今回は63で9週間）
• lags_future_covariates:将来の説明変数を学習に利用する際の期間
• lags_past_covariates:過去の説明変数を学習に利用する際の期間
• output_chunk_length:どれくらいの期間を一回で予測するか。短いと正確になるが、計算コストが高い。長いとスピードは早いが、長期的には精度が低下する。
• random_state:ランダムな操作やアルゴリズムを実行する際に使用される乱数生成器のシード値（seed）を指定する。（指定しなかったら、ランダムフォレストなどランダムなアルゴリズムの出力結果が毎回変わってしまう。）
• gpu_use_dp:GPUを使用するかどうか。

• ハイパーパラメーター:機械学習の性能を制御するための、手動で設定するもののこと。
• ラグ特徴量:過去の目的変数を特徴量として活用すること。
• GPU:同時並行に処理をすることができるプロセッサ。（CPUの一種であり、対義語。CPUは一つのタスクを超高速にできるが、同時並行が苦手。）

from darts.models import LightGBMModel

LGBM_Models_Submission = {}

display("Training...")

for family in tqdm(family_list):

  sales_family = family_TS_transformed_dict[family]
  training_data = [ts for ts in sales_family] 
  TCN_covariates = future_covariates_dict[family]
  train_sliced = [training_data[i].slice_intersect(TCN_covariates[i]) for i in range(0,len(training_data))]

  LGBM_Model_Submission = LightGBMModel(lags = 63,
                                        lags_future_covariates = (14,1),
                                        lags_past_covariates = [-16,-17,-18,-19,-20,-21,-22],
                                        output_chunk_length=1,
                                        random_state=2022,
                                        gpu_use_dp= "false",
                                        )
     
  LGBM_Model_Submission.fit(series=train_sliced, 
                        future_covariates=TCN_covariates,
                        past_covariates=transactions_transformed)

  LGBM_Models_Submission[family] = LGBM_Model_Submission

予測する

トレーニングしたものを使って、for文で商品群ごとに予測する。
辞書型で'LGBM_Forecasts_Families_Submission'に予測値を格納する。

display("Predictions...")

LGBM_Forecasts_Families_Submission = {}

for family in tqdm(family_list):

  sales_family = family_TS_transformed_dict[family]
  training_data = [ts for ts in sales_family]
  LGBM_covariates = future_covariates_dict[family]
  train_sliced = [training_data[i].slice_intersect(TCN_covariates[i]) for i in range(0,len(training_data))]

  forecast_LGBM = LGBM_Models_Submission[family].predict(n=16,
                                         series=train_sliced,
                                         future_covariates=LGBM_covariates,
                                         past_covariates=transactions_transformed)
  
  LGBM_Forecasts_Families_Submission[family] = forecast_LGBM

予測値のスケーリングを元に戻す

スケーリングをinverse_transformで元に戻す

LGBM_Forecasts_Families_back_Submission = {}

for family in tqdm(family_list):

  LGBM_Forecasts_Families_back_Submission[family] = family_pipeline_dict[family].inverse_transform(LGBM_Forecasts_Families_Submission[family], partial=True)

提出用データに整形

得られた直近21日間の予測結果が0の場合は、全てゼロにする。（予測値がマイナスになることを防ぐため。）
54の店舗と各製品群ごとの予測結果を結合して、dataframeを作成する。

ーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーー

• concat:複数のdataframeを結合する。
• merge:2つのdataframeを列を基準に結合する。
• join:複数のdataframeを行を基準に結合する。

参考サイト
図解で解決！pandasデータ結合総まとめ（concat, merge, join）

for family in tqdm(LGBM_Forecasts_Families_back_Submission):
  for n in range(0,len(LGBM_Forecasts_Families_back_Submission[family])):
    if (family_TS_dict[family][n].univariate_values()[-21:] == 0).all():
        LGBM_Forecasts_Families_back_Submission[family][n] = LGBM_Forecasts_Families_back_Submission[family][n].map(lambda x: x * 0)

listofseries = []

for store in tqdm(range(0,54)):
  for family in family_list:
      oneforecast = LGBM_Forecasts_Families_back_Submission[family][store].pd_dataframe()
      oneforecast.columns = ['fcast']
      listofseries.append(oneforecast)

df_forecasts = pd.concat(listofseries) 
df_forecasts.reset_index(drop=True, inplace=True)

# No Negative Forecasts
df_forecasts[df_forecasts < 0] = 0
forecasts_kaggle = pd.concat([df_test_sorted, df_forecasts.set_index(df_test_sorted.index)], axis=1)
forecasts_kaggle_sorted = forecasts_kaggle.sort_values(by=['id'])
forecasts_kaggle_sorted = forecasts_kaggle_sorted.drop(['date','store_nbr','family'], axis=1)
forecasts_kaggle_sorted = forecasts_kaggle_sorted.rename(columns={"fcast": "sales"})
forecasts_kaggle_sorted = forecasts_kaggle_sorted.reset_index(drop=True)

# Submission
submission_kaggle = forecasts_kaggle_sorted

ハイパーパラメータを調整

3パターンのハイパーパラメーターを設定する。
今回はラグ特徴量を1週間毎、1年毎、2年毎で設定する。

model_params = [
    {"lags" : 7, "lags_future_covariates" : (16,1), "lags_past_covariates" : [-16,-17,-18,-19,-20,-21,-22]},
    {"lags" : 365, "lags_future_covariates" : (14,1), "lags_past_covariates" : [-16,-17,-18,-19,-20,-21,-22]},
    {"lags" : 730, "lags_future_covariates" : (14,1), "lags_past_covariates" : [-16,-17,-18,-19,-20,-21,-22]}
]

設定したハイパーパラメーターでそれぞれ機械学習する

3種類のハイパーパラメーターを用いて、機械学習を進めてみる。
コードが長いが大枠では下記の流れ。

1. 前回設定したハイパーパラメーターのモデルを作成し、それぞれトレーニング
2. それぞれ予測値を出す。
3. スケーリングを元に戻す。
4. 直近21日間の売上がゼロの場合をゼロと予測するように予測値を変換する。
5. 特定の商品群と店舗の組み合わせに対する予測結果を取得して、dataframeに変換する。
6. 生成したdataframeを使うために、リストにappendしていく。
7. appendしたリストをdataframeとして結合する。
8. ゼロ以下の予測値をゼロに変換する。
9. dataframeを整形する。（列の削除、並び替え、列名の変更、インデックス番号の振り直し）
10. リストに格納する。

from sklearn.metrics import mean_squared_log_error as msle, mean_squared_error as mse
from lightgbm import early_stopping

submission_kaggle_list = []

for params in model_params:

  LGBM_Models_Submission = {}

  display("Training...")

  for family in tqdm(family_list):

    # Define Data for family
    sales_family = family_TS_transformed_dict[family]
    training_data = [ts for ts in sales_family] 
    TCN_covariates = future_covariates_dict[family]
    train_sliced = [training_data[i].slice_intersect(TCN_covariates[i]) for i in range(0,len(training_data))]

    LGBM_Model_Submission = LightGBMModel(lags = params["lags"],
                                          lags_future_covariates = params["lags_future_covariates"],
                                          lags_past_covariates = params["lags_past_covariates"],
                                          output_chunk_length=1,
                                          random_state=2022,
                                          gpu_use_dp= "false")
      
    LGBM_Model_Submission.fit(series=train_sliced, 
                          future_covariates=TCN_covariates,
                          past_covariates=transactions_transformed)

    LGBM_Models_Submission[family] = LGBM_Model_Submission
    
  display("Predictions...")


  LGBM_Forecasts_Families_Submission = {}

  for family in tqdm(family_list):

    sales_family = family_TS_transformed_dict[family]
    training_data = [ts for ts in sales_family]
    LGBM_covariates = future_covariates_dict[family]
    train_sliced = [training_data[i].slice_intersect(TCN_covariates[i]) for i in range(0,len(training_data))]

    forecast_LGBM = LGBM_Models_Submission[family].predict(n=16,
                                          series=train_sliced,
                                          future_covariates=LGBM_covariates,
                                          past_covariates=transactions_transformed)
    
    LGBM_Forecasts_Families_Submission[family] = forecast_LGBM

  # Transform Back

  LGBM_Forecasts_Families_back_Submission = {}

  for family in tqdm(family_list):

    LGBM_Forecasts_Families_back_Submission[family] = family_pipeline_dict[family].inverse_transform(LGBM_Forecasts_Families_Submission[family], partial=True)

  # Prepare Submission in Correct Format

  for family in tqdm(LGBM_Forecasts_Families_back_Submission):
    for n in range(0,len(LGBM_Forecasts_Families_back_Submission[family])):
      if (family_TS_dict[family][n].univariate_values()[-21:] == 0).all():
          LGBM_Forecasts_Families_back_Submission[family][n] = LGBM_Forecasts_Families_back_Submission[family][n].map(lambda x: x * 0)
          
  listofseries = []

  for store in tqdm(range(0,54)):
    for family in family_list:
        oneforecast = LGBM_Forecasts_Families_back_Submission[family][store].pd_dataframe()
        oneforecast.columns = ['fcast']
        listofseries.append(oneforecast)

  df_forecasts = pd.concat(listofseries) 
  df_forecasts.reset_index(drop=True, inplace=True)

  # No Negative Forecasts
  df_forecasts[df_forecasts < 0] = 0
  forecasts_kaggle = pd.concat([df_test_sorted, df_forecasts.set_index(df_test_sorted.index)], axis=1)
  forecasts_kaggle_sorted = forecasts_kaggle.sort_values(by=['id'])
  forecasts_kaggle_sorted = forecasts_kaggle_sorted.drop(['date','store_nbr','family'], axis=1)
  forecasts_kaggle_sorted = forecasts_kaggle_sorted.rename(columns={"fcast": "sales"})
  forecasts_kaggle_sorted = forecasts_kaggle_sorted.reset_index(drop=True)

  # Submission
  submission_kaggle_list.append(forecasts_kaggle_sorted)

アンサンブル学習

4つの予測結果が得られるので、それらの各予測値の平均を算出する。

df_sample_submission['sales'] = (submission_kaggle[['sales']]+submission_kaggle_list[0][['sales']]+submission_kaggle_list[1][['sales']]+submission_kaggle_list[2][['sales']])/4

csvファイルに格納

生成した予測をファイルに抽出する。

df_sample_submission.to_csv('/kaggle/working/submission.csv', index=False)

結果

結果は誤差は下記となった。
６８０人参加しているコンペで、３０位以内と上位５％の結果となった。

流石だ…

所感

今回、模写して下記がわかった。

• 時系列データの大枠の機械学習の流れ
  • 特徴量作成
    • 年、月、曜日、年の中の何周目かなどで特徴量を生成する。
  • データの正規化
    • Scaler→データのスケーリング
    • StaticCovariatesTransformer→ワンホットエンコーディングなど
    • MissingValuesFiller→欠損値処理
    • InvertibleMapper→対数変換など
  • データ結合
  • 機械学習
  • ハイパーパラメーター調整
  • アンサンブル学習
• ハイパーパラメーターの調整方法
  • 辞書型でパラメーターを保有しておくと使いやすい。
• ラグ特徴量の作成方法
  • １週間ごとや３ヶ月毎、１年毎などでラグのスパンを複数パターン試すとトレンドが見える。
• 商品、店舗ごとにデータを集計するためにはfor文が使える。
  • とにかく出てきまくったfor文

概要として触ってみた感覚が強いので、今後のKaggleや実務の中で慣れていく。

課題・今後

一方で下記は今後、理解を深める必要がある。

• Prophetなどの他の時系列データのライブラリとの違い
• ハイパーパラメーターの設定の考え方
• そもそも商品ごとで見るべきか、店舗ごとにデータを見るべきかの考え方
• データのグラフ化などによる、機械学習の方針決定

まだまだ素人だが、一歩ずつ進むことで「データ」でビジネスインパクトを起こせる人間になる。
次はProphetを触ってみようかと思う。