XGBoostでTime Series Data処理 備忘録

# 概要
仕事でtime series datasetをまた扱うことになった。LSTMに入れる前にXGBoostで味見してみようと手弁当で下調べ中、下記サイトを見つけた。
日付の処理の仕方が興味深く、だったらdatetimeデータがなくてもい行けるんじゃね？とただの数字の羅列を説明変数に入れてみた備忘録。

https://www.kaggle.com/code/robikscube/tutorial-time-series-forecasting-with-xgboost

- 実施期間: 2022年3月
- 環境：Colab
- パケージ：XGBoost

## Time series dataset (Seasonal)
はじめにColabのXGBoostが0.90とやや古いので入れ直した。

```terminal
!pip3 install xgboost==1.1.1
```

```python
import numpy as np
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
import xgboost as xgb
from xgboost import plot_importance, plot_tree
from sklearn.metrics import mean_squared_error, mean_absolute_error
from sklearn.model_selection import train_test_split
print(xgb.__version__)
```

prediction対象は、完全syntheticでseasonalな三角関数とする。そのままでは面白くないのでsin + cosとした。

```python
data_num = 15000
x = np.linspace(0, data_num * 4, data_num)
sin_arr1 = sin_arr2 = sin_arr3 = 0
sin_arr1 = np.sin(3 * x) /1.5
sin_arr2 = np.cos(1.5 * x) /2
# sin_arr3 = np.cos(8 * x + 1) /2    # 時間かかるので使わん
sin_arr = sin_arr1 + sin_arr2 + sin_arr3

plt.figure(figsize=(10,3), dpi=100)
plt.plot(x[:100], sin_arr[:100])
plt.show()
```

![Screenshot from 2022-03-27 13-07-57.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/205257/1d3c5c85-dfd7-bff6-a969-6ca5886122b6.png)

このsin_arrを予測する。


前述のKaggleの例では日付や月やwork weekなどを数値化してそれぞれをfeatureにされていた。この考え方がOKなら今回のようなdatetimeを持たいないseries dataでも、datetimeに相応するN進数データをfeatureにすればいいんじゃないかと。
ここでは8進数、9進数、...、49進数までをfeatuerとして追加する。三角関数の周期性をこのうちのいくつかが説明してくれるだろうという目論。ちなみに2~7進数はfeature importanceが小さかったので入れない。

```python
df = pd.DataFrame()
for i in range(8, 50, 1):
    x2_ = list(np.arange(0,i))
    x2_ = np.array(x2_*10000)   # 10000：スライスしてトリムするから十分大きなどうでも良い値
    x2_ = x2_[0:data_num]
    df[i] = x2_

df['val'] = sin_arr

print(df)
```

![Screenshot from 2022-03-27 13-22-29.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/205257/bbe67009-60c4-fede-bff9-05aaa9de22e3.png)


## TrainingとFeature Importance
特に説明はない。

```python
X = np.array(df.iloc[:,0:42])
Y = np.array(df.iloc[:,42])

test_size = 0.1
X_train, X_test, y_train, y_test = train_test_split(
        X, Y, 
        test_size=test_size, 
        shuffle = False, stratify = None)

reg = xgb.XGBRegressor(n_estimators=10000, 
        tree_method='gpu_hist', random_state=0, 
        objective='reg:squarederror')

reg.fit(X_train, y_train,
        eval_set=[(X_train, y_train), (X_test, y_test)],
        early_stopping_rounds = 100,
        verbose = True)
```

feature importanceを確認する。

```python
_, ax = plt.subplots(figsize=(10, 8), dpi=100)
_ = plot_importance(reg1, ax=ax, height=0.6)
```

![Screenshot from 2022-03-27 13-30-42.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/205257/ca7d57f9-ef4c-8d14-e351-3722060eb9bd.png)
上位から順に**f35**は44進数、次の**f38**は46進数、**f15**は23進数を示す。
おそらくf35,38は波高値大きなsin関数を説明し、f15は周期を倍にしたcos関数を説明しているような気がする。さながらXGBoostを使ったFFTのような感じ。

## Predictionと精度確認

Prediction結果を描画する。

```python
y_pred = reg1.predict(X_test)

plt.figure(figsize=(20,3), dpi=100)

plt.plot(x[-y_test.shape[0]:], y_test)
plt.plot(x[-y_pred.shape[0]:], y_pred)
plt.show()
```

![Screenshot from 2022-03-27 13-43-24.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/205257/ff18af28-e82d-2c7e-6143-3294cca59fec.png)

いい感じに見える。簡単にy_testとy_predの差も確認する。

```python
print(sum(y_test - y_pred))
```
-7.227987089052772



## Time series dataset (Seasonal + Trend)
前述のsin_arrに0.5/15000の勾配でtrend成分を加えてみる。
N進数はうまくsin+cosのseasonality(周期性)を説明したがこの一定の勾配を説明するfeatureがないので精度が悪くなることを期待して確認する。

```python
data_num = 15000
y_trend = np.linspace(0, 0.5, data_num)
print(y_trend)
sin_arr = sin_arr + y_trend
```

途中割愛し、feature importanceを見てみる。

![Screenshot from 2022-03-27 14-14-10.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/205257/c9a384db-c827-a991-5049-ee5558d8cd5a.png)

どのfeatureを見てよいのか迷っているみたい。
予測結果は下図のように周期はあっているが波高値は先に進むに連れ悪化している。思惑通り一定勾配を説明するfeatureがないからだろう。

![Screenshot from 2022-03-27 14-15-57.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/205257/34c3a2f8-2697-431e-197d-b0bc84cc3f1a.png)

```python
print(sum(y_test - y_pred))
```
324.42460065419976

以上