tidyverts_peyton_manning_wiki.Rmd

---
title: "Peyton Manning Wiki Tidyverts"
output: html_notebook
---

Rの時系列モデリングのライブラリ `tidyverts` を用いた基本的なモデル構築と予測値算出値を行います。

### tidyverts とは？
[tidyverts](https://tidyverts.org/) は、大きく下記の3つのライブラリから構成されます。

- [tsibble](https://tsibble.tidyverts.org/) : 時系列データに特化たデータ型を提供
- [fable](https://fable.tidyverts.org/) : 時系列予測のモデルを提供
- [feasts](https://feasts.tidyverts.org/) : 時系列データの統計処理や特徴量抽出の機能を提供



### 利用するデータ
Prophet のチュートリアルで使用されている Peyton Manning の Wikipedia の閲覧者数の時系列データを利用します。下記 URL から CSV データをダウンロードしておきます。

[example_wp_log_peyton_manning.csv](https://github.com/facebook/prophet/blob/master/examples/example_wp_log_peyton_manning.csv)


### 分析環境

R Studio や Visual Studio Code で構いません。私は Azure Machine Learning で提供している PaaS の R Studio を利用しています。

[Azure Machine Learning](https://azure.microsoft.com/ja-jp/services/machine-learning/)

[R Support in Azure Machine Learning | AI Show | Channel9](https://channel9.msdn.com/Shows/AI-Show/R-in-Azure-Machine-Learning)


### ライブラリのロード

利用するライブラリをロードします。

```{r}
library(fable)
library(fable.prophet)
library(tsibble)
library(tsibbledata)
library(feasts)
library(lubridate)
library(dplyr)
library(ggplot2)
```


### データ準備と探索

CSV ファイルのパスを指定して、データをロードします。

```{r}
df <- read.csv("example_wp_log_peyton_manning.csv", colClasses=c("Date","numeric")) 
```

df に格納されているデータを `tsibble`形式に変換し可視化します。"%>%" は R 言語で利用されるパイプ処理を意味しています。autoplot を適用すると可視化まで簡単にできます。

```{r}
df %>% as_tsibble %>% autoplot
```
```{r}
df %>% as_tsibble %>% fill_gaps
```
[`fill_gaps`](https://tsibble.tidyverts.org/reference/fill_gaps.html) を利用して、データに欠損している日時項目を補完します。また `tidyr` の [`fill`](https://tidyr.tidyverse.org/reference/fill.html) で欠損値を補完します。１つ先の時系列データ値を用いて補完する場合は、`.direction  = down` と指定します。

```{r}
df_filled <- df %>% as_tsibble %>% fill_gaps %>% tidyr::fill(y, .direction = "down")
```
次に `feasts` ライブラリの [`STL`](https://feasts.tidyverts.org/reference/STL.html) を用いて、時系列データを "季節" "トレンド" "不規則成分" に分解します。

```{r}
dcmp <- df_filled %>%
    model(STL(y ~ season(window = Inf ))) %>% components
dcmp
```
可視化して特徴を確認します。Weekly の季節性トレンドは細かすぎて黒塗りに見えちゃっているようです。
```{r}
dcmp %>% autoplot
```

### 時系列モデルの構築と予測値算出

まずは非常にシンプルなモデルを試します。データには季節性などが確認されるものの、これらのモデルは考慮してくれないことが分かります。

- MEAN : 平均値
- NAIVE(=RW): ランダムウォークモデル


ここでは、ドリフトを考慮するモデルも追加します。このデータは後半にトレンドが減少している傾向にあるので、このモデルも減少方向に予測していることが確認できます。
```{r}
df_filled %>% 
    model(
        mean=MEAN(y),
        naive=NAIVE(y),
        drift=RW(y ~ drift()),  # ドリフトを考慮
        .safely=FALSE) %>% 
    forecast(h="2 years") %>% 
    autoplot(filter(as_tsibble(df_filled),year(ds) > 2007), level = NULL)
```
次に季節性を考慮します。

- SNAIVE : 季節性を考慮するランダムウォークモデル

1年周期の季節性があるデータなので、過去 1 年前のデータを予測とするモデルを採用します。

```{r}
df_filled %>% 
    model(
        mean=MEAN(y),
        naive=NAIVE(y),
        drift=RW(y ~ drift()),  # ドリフトを考慮
        snaive = SNAIVE(y ~ lag("year")), # 季節性を考慮
        .safely=FALSE) %>% 
    forecast(h="2 years") %>% 
    autoplot(filter(as_tsibble(df_filled),year(ds) > 2007), level = NULL)
```

次に指数平滑化を使用します。

- ETS : 指数平滑化

1年周期を適用してみます。
```{r}
df_filled %>% 
    model(
        est=ETS(y ~ season("A", period= "1 year")),
        .safely=FALSE) %>% 
    forecast(h="1 year") %>% 
    autoplot(filter(as_tsibble(df_filled),year(ds) > 2015), level = NULL)
```

すると、エラーが出てきました。

> Error: Seasonal periods (`period`) of length greather than 24 are not supported by ETS.

最大 24 サイクルしか考慮してくれないようです。このデータは Daily のデータなので、年周期の季節性の考慮は難しそうです。

週周期の季節性があると想定してみます。

```{r}
df_filled %>% 
    model(
        est=ETS(y ~ season("A", period= "7 days")),
        .safely=FALSE) %>% 
    forecast(h="30 days") %>% 
    autoplot(filter(as_tsibble(df_filled),year(ds) > 2014), level = NULL)
```

次に ARIMA モデルを利用します。ETS ではうまく動きませんでしたが、

- ARIMA : 自己回帰和分移動平均

```{r}
df_filled %>% 
    model(
        #arima=ARIMA(y ~ pdq(1,1,0) + PDQ(0,1,0,period="365 days")),
        arima=ARIMA(y ~ pdq(1,1,0) + PDQ(0,1,0,period="365 days")),
        .safely=FALSE) %>%
        forecast(h="2 years") %>% 
    autoplot(filter(as_tsibble(df_filled),year(ds) > 2007), level = NULL)
```

次に、Facebook の Prophet を利用します。`fable.prophet` というライブラリが公開されており、いままで使ってきた `fable`に Prophet のモデルを組み込めるようになっています。詳細は [fable.prophet](https://github.com/mitchelloharawild/fable.prophet) からご確認ください。


```{r}
df_filled %>% 
    model(
        prophet(y),
        .safely=FALSE) %>% 
    forecast(h="2 years") %>% 
    autoplot(filter(as_tsibble(df_filled),year(ds) > 2007), level = NULL)
```