時系列データ分析：前処理

勉強したことをアウトプットします。
ほぼ全て北川先生の講義によるものです。
講義：https://ocwx.ocw.u-tokyo.ac.jp/course_11416/
データ：http://www.mi.u-tokyo.ac.jp/mds-oudan/lecture_document_2019_math7/time_series_analysis_2019.html


時系列データにおける前処理の２つの目的
-
１、定常性や正規性を仮定に置いているモデリングを容易に適用するため。
⇨良くわからない周期性を持つ時系列データから、なんとか周期性を見つけ出し、その周期性をモデルに落とし込み予測できるようにしたいという動機。
２、モデルのパラメータ最適化を容易にするため。

データの種類
-
４種類ある。
１、独立したデータ
各データは独立している。（例）心臓病のデータ、国勢調査、治験データなど

２、時系列データ
時間と供に生成されたデータ。データの並び順が大切になる。（例）株価、地震など

３、空間データ
mapなどの地理情報と結び付けられて生成されるデータ。

４、時空間データ
２と３を組み合わせたもの。データは時間と空間の２つに結び付けられて生成されたものとなる。（例）コロナの感染状況分析、気象情報、人口流体データなど

今回は２の分析方法を学ぶ。

変数変換
-
前処理の一種。変数変換する方法はたくさんある。
・対数
・階差
・対数階差
・ロジット
・季節変動
・前期比
・平均、メディアン移動平均
・BoxーCox変換
など

変数変換でしたい事を別観点から再説明
-
```:データをh(x)で変換
          y = h(x)
X ⇨ f(x)  -------→   Y ⇨ g(y)
```
Xはデータが持つ確率分布、変換する関数h(x)を通して、Yを正規分布などの定常性を持った確率分布に変換しモデリングできるようする（上記の前処理の目的の１つ目と同義）。


変換後データを逆関数を利用し元に戻す場合。
```:h(x)で変換したデータを元に戻す
          x = h(y)^{-1}
X ⇨ f(x)  ←-------   Y ⇨ g(y)
```

（余談：難しすぎる例）
正規分布を対数正規分布に変換する場合
![スクリーンショット 2023-01-06 16.20.51.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/4d07803b-5d18-d246-ff46-2aad7cb4d182.png)

対数変換
-
logをとる

#### 例：商品の月別売上高
下グラフを見ると、売上金額が伸びていくにつれて（月日が経つにつれて）、分散が大きくなり定常性を満たしていない。
![スクリーンショット 2023-01-06 16.24.03.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/7a795b77-e78c-ca8c-1d6f-8ddf0a5d1b44.png)
↓
対数を取ると、分散が一定になっているように見える
![スクリーンショット 2023-01-06 16.25.52.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/390313ed-ce97-712e-dcbf-9fdc62cd369a.png)

階差
-
nは時刻
```math
\nabla y_{(n)} = y_{(n)} - y_{(n-1)}
```
時刻nのデータと時刻n-1のデータの差分を取る場合、それは１階階差呼ぶ。

２階階差の定義
```math
\nabla ^ 2 y_{(n)} = \nabla y_{(n)} - \nabla y_{(n-1)}
```

k階階差のイメージ
![スクリーンショット 2023-01-06 16.35.27.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/99080487-a272-ea58-9a3e-6d852857e9a7.png)

#### 2階階差がよく使われる理由
データの確率分布を２次関数と置いた場合、２階階差は定数となり、時系列解析が容易くなるため。

$y = a + bn + cn^2$　を考える

$\nabla y = y_{(n)} - y_{(n-1)} = ・・・ = (b - c) - 2cn$

$\nabla ^ 2 y = ・・・ = 2c$

つまり非線形のデータでも周期性を見つけ出し、分析できるようにするためよく使われる。

#### 階差をとったデータを視覚化
(例):商品の月別売上高の対数を取ったものを階差をプロット
```python:コード
df = pd.read_csv('data/whard_new.csv')
df['log_whard'] = np.log(df.whard)
plt.plot(np.diff(df['log_whard'], n = kaisa)) #kaisaは自身が利用したい階差の数
plt.show()
```
１階階差
![スクリーンショット 2023-01-06 16.45.58.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/c053b654-9269-47ab-978a-2988e7133a51.png)

２階階差
![スクリーンショット 2023-01-06 16.48.14.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/fd8f1b17-26cc-5d8b-8c78-b6b3d89d8f93.png)

6階階差
![スクリーンショット 2023-01-06 16.48.25.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/54cd2d04-8975-2847-7108-e805bbbd0b21.png)

1２階階差
![スクリーンショット 2023-01-06 16.48.33.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/7ae58897-d72e-9b34-4d1b-f39e82610f15.png)

対数階差
-
対数を取ったものの階差をとる。データが一定の幅に収めたいという動機があって利用する。
↓普通の日経平均２２５とそれの対数階差変換したもの
![スクリーンショット 2023-01-06 16.53.15.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/4018ae54-93ce-4da5-e6f7-31299fd9de04.png)
なんか周期性（と異常値）を見つけられそうな気がする。

季節変動
-
sは季節（１周期）の長さ。季節変動を除去できる。
```math
\nabla y = y_{(n)} - y_{(n-s)}
```

#### グラフで確認
商品の月別売上高推移の対数変換後
![スクリーンショット 2023-01-06 16.25.52.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/193cef32-3a47-263e-9e62-85326b3f8ef4.png)
↓周期を設定したもの
・周期が12ヶ月の場合（年単位で周期があると想定）。異常が発生している時がよくわかる。
![スクリーンショット 2023-01-06 17.19.16.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/25e2d5af-32e0-824a-df20-abab0a1510ef.png)
・周期が3ヶ月の場合（四季で周期があると想定）。周期性があるだけに見える。少なくとも12ヶ月で季節変動を除去した場合に確認できた異常値は確認できなかった。
![スクリーンショット 2023-01-06 17.21.25.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/f97b0d95-079a-ed43-0c6e-3155949d9b90.png)

```python:コード
period = 12
z = np.array([])
for i in np.arange(period,len(df['log_whard'])):
    num = df.loc[i,'log_whard'] - df.loc[i - period,'log_whard']
    z = np.append(z, num)
plt.figure(figsize=(12,4))
plt.plot(z)
plt.show()
```
※pythonでstatsmodels.apiのseasonal_decomposeで季節変動を同じように除去できるらしいがよくわからなかった。

前期比
-
nは時刻。企業が決算説明の時などによく利用する前期比売上高など、時刻n-1と時刻nで比例関係がある場合などに利用できる。対数とほぼ同じ折れ線グラフになる。

```math
x(n) = \frac{t_{(n)}}{t_{(n-1)}}
```

平均、メディアン移動平均
-

大まかなトレンドを掴むために利用。

#### 移動平均
時刻nの前後k個分のデータを活用する
```math
t_{(n)} = \frac{1}{2k + 1}(y_{(n-k)} + ・・・ + y_{n} + ・・・ + y_{n+k})
```
- メリット
  - 滑らかな曲線を描ける

- デメリット
  - 外れ値に弱い
  - 急激な変化構造を捉えられない。


#### 移動median
```math
t_{(n)} = median(y_{(n-k)} + ・・・ + y_{n} + ・・・ + y_{n+k})
```

- メリット
  - 外れ値に強い
  - 急激な変化構造を捉えやすい（下グラフのn=25のとき、時系列の最後の方を確認してほしい。）

- デメリット
  - 分散が大きくなりがち

#### 日経平均２２５で確認
加工前データ
![スクリーンショット 2023-01-06 17.41.31.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/af7a9d43-a7af-8570-5460-6311dc0ac293.png)

![スクリーンショット 2023-01-06 17.42.02.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/5579fd1f-141b-2d87-0eff-85e5310ee56a.png)
![スクリーンショット 2023-01-06 17.42.07.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/3e5fa5b0-74be-747f-aab6-b12aa20c0dca.png)
![スクリーンショット 2023-01-06 17.42.12.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/eef595e5-76e3-7db1-c5fa-5f4dabacb180.png)
![スクリーンショット 2023-01-06 17.42.17.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2707093/6f1f1bc0-887d-9478-ee43-90797d1f8ae0.png)

他参考にしたもの
-
季節変換について

https://qiita.com/shinji_komine/items/9dbada69b9128f645660