Forecasting Netfilx’s Stock Using Time Series

Introduction

Pada analisis kali ini, saya akan menggunakan data Netflix stock yang akan digunakan untuk melakukan forecast harga saham Netflix pada 3 bulan ke depan, dengan mengidentifikasi data Close yang telah disediakan pada netflix_stock.csv

Data Preparation

Library

library(readr)
library(dplyr) 
library(lubridate) 
library(ggplot2) 
library(TSstudio) 
library(padr) 
library(forecast) 
library(ggsci)
library(tidyr) 
library(zoo)

Import Data

stock <- read.csv("datainput/netflix_stock.csv")
head(stock)

Data Cleansing

🧮 Mengecek tipe data

glimpse(stock)

#> Rows: 1,009
#> Columns: 7
#> $ Date      <chr> "2018-02-05", "2018-02-06", "2018-02-07", "2018-02-08", "201…
#> $ Open      <dbl> 262.00, 247.70, 266.58, 267.08, 253.85, 252.14, 257.29, 260.…
#> $ High      <dbl> 267.90, 266.70, 272.45, 267.62, 255.80, 259.15, 261.41, 269.…
#> $ Low       <dbl> 250.03, 245.00, 264.33, 250.00, 236.11, 249.00, 254.70, 260.…
#> $ Close     <dbl> 254.26, 265.72, 264.56, 250.10, 249.47, 257.95, 258.27, 266.…
#> $ Adj.Close <dbl> 254.26, 265.72, 264.56, 250.10, 249.47, 257.95, 258.27, 266.…
#> $ Volume    <int> 11896100, 12595800, 8981500, 9306700, 16906900, 8534900, 685…

⚙️ Mengubah tipe data pada kolom Date

stock <- stock %>% 
  mutate(Date = ymd(Date))

Data Preprocessing

⚙️ Memfilter data yang diperlukan:

stock <- stock %>% 
  select(c(Date,Close))
stock

🧮 Mengecek range waktu data

range(stock$Date)

#> [1] "2018-02-05" "2022-02-04"

⚙️ Melakukan pengecekkan inspect periode terlewat:

stock <- stock %>% 
  arrange(Date)

complete_stock <- seq.Date(from = as.Date("2018-02-05"),
                         to = as.Date("2022-02-04"),
                         by = "day")
all(complete_stock == stock$Date)

#> [1] FALSE

️⚙️ Melakukan padding untuk memastikan interval data yang sama: pad() dari package padr dan replace nilai missing value menggunakan na.locf() agar di isi dengan nilai Close sebelum Missing:

stock <- stock %>% 
  pad(start_val = min(stock$Date), end_val = max(stock$Date)) %>% 
  mutate(Close = na.locf(Close))

stock

⚙️ Untuk memastikannya kembali dapat melakukan pengecekkan missing value

anyNA(stock)

#> [1] FALSE

Seasonality Analysis

Buat objek Time Series menggunakan fungsi ts()

stock_ts <- ts(data = stock$Close,
            frequency = 7*4*12)
#melihat plot time series
stock_dc <- stock_ts %>% decompose()
stock_dc %>% autoplot()

Pola trend diatas masih belum cukup smooth, sehingga ada kemungkinan terjadinya Multi-Seasonality. Oleh karena itu, perlu diubah ke dalam msts:

stock_msts <- msts(data = stock$Close,seasonal.periods = c(7*4, 7*4*12))

stock_msts_dc <- mstl(stock_msts)

stock_msts_dc %>% autoplot()

📈 Insight:
Berdasarkan Plot di atas, kini estimasi trend dengan Multiple Seasonality Time Series lebih smooth dan jelas.

Model Fitting and Evaluation

Cross Validation

Splitting Data :
Data Test : Data yang dijadikan test untuk validasi dengan mengambil pola 3 bulan (7 x 4 x 3)
Data Train : Total data kecuali data test

test_stock <- tail(stock_msts, 7*4*3)
train_stock <- head(stock_msts, length(stock_msts)-length(test_stock))

ARIMA MOdel

Menggunakan model ARIMA (Autoregressive Integrated Moving Average) adalah gabungan antara dua metode, yaitu Auto Regressive (AR) dan Moving Average (MA). Kemudian I nya menjelaskan Integrated. Syarat data model ARIMA adalah data stationer.

model_arima_msts <- stlm(train_stock, method = "arima")

Forecasting

arima_forecast_msts <- forecast(model_arima_msts, h = 7*4*3)
plot(arima_forecast_msts)

Model Visualization

test_stock %>%  autoplot(series = "Actual") +
  autolayer(arima_forecast_msts$mean, series = "Predict Arima")

Evaluation Performance Model

Model dievaluasi dengan melihat nilai error menggunakan MAPE (Mean Absolute Percentage Error). Untuk dapat menghitung error yang dihasilkan oleh model dengan menggunakan fungsi accuracy()

accuracy(object = arima_forecast_msts$mean, x = test_stock)

#>                 ME     RMSE      MAE       MPE   MAPE      ACF1 Theil's U
#> Test set -110.8369 149.4735 114.4361 -23.63202 24.158 0.9549567  10.34932

Holt-Winters Model

Dikarenakan model pada data memiliki trend dan seasonal, maka dapat menggunakan Triple Exponential Smoothing

model_hw_msts <- HoltWinters(train_stock)

Forecasting

hw_forecast_msts <- forecast(model_hw_msts, h = 7*4*3)
plot(hw_forecast_msts)

Model Visualization

test_stock %>% autoplot(series = "Actual") +
  autolayer(hw_forecast_msts$mean, series = "Predict Holt-Winters")

Evaluation Performance Model

accuracy(object = hw_forecast_msts$mean, x = test_stock)

#>                 ME     RMSE      MAE       MPE     MAPE      ACF1 Theil's U
#> Test set -130.4723 165.0421 130.8665 -27.23893 27.29639 0.9544257  11.29629

ETS Model

model_ets_msts <- stlm(train_stock, method = "ets")

Forecasting

ets_forecast_msts <- forecast(model_ets_msts, h = 7*4*3)
plot(ets_forecast_msts)

Model Visualization

test_stock %>% autoplot(series = "Actual") +
  autolayer(ets_forecast_msts$mean, series = "Predict ETS")

Evaluation Performance Model

accuracy(object = ets_forecast_msts$mean, x = test_stock)

#>                 ME     RMSE      MAE       MPE     MAPE      ACF1 Theil's U
#> Test set -97.03894 135.2118 101.4476 -20.92604 21.57105 0.9537757  9.442058

Visualization Actual vs Estimated

accuracydata <- data.frame(Date= stock$Date %>% tail(7*4*3),
  Actual = as.vector(test_stock) ,
  ArimaForecastmsts = as.vector(arima_forecast_msts$mean),
  HWForecastmsts = as.vector(hw_forecast_msts$mean),
  ETSForecastmsts = as.vector(ets_forecast_msts$mean))
head(accuracydata)

Visualization of Actual vs Estimated number of Close by All Models

accuracydata %>% 
 ggplot() +
 geom_line(aes(x = Date, y = Actual, colour = "Actual")) +
  geom_line(aes(x = Date, y = ArimaForecastmsts, colour = "Arima Model")) +
  geom_line(aes(x = Date, y = HWForecastmsts, colour = "Holt Winter Model")) +
  geom_line(aes(x = Date, y = ETSForecastmsts, colour = "ETS Model")) +
  labs(title = "Closing Price by Actual vs All Models",
       x = "Date",
       y = "Visitor",
       colour = "")

Prediction Performance

result <- rbind(
  accuracy(object = arima_forecast_msts$mean, x = test_stock),
  accuracy(object = hw_forecast_msts$mean, x = test_stock),
  accuracy(object = ets_forecast_msts$mean, x = test_stock))
rownames(result) <- c("ARIMA_Model", "Holt-Winters_Model", "ETS Model")
result

#>                            ME     RMSE      MAE       MPE     MAPE      ACF1
#> ARIMA_Model        -110.83687 149.4735 114.4361 -23.63202 24.15800 0.9549567
#> Holt-Winters_Model -130.47226 165.0421 130.8665 -27.23893 27.29639 0.9544257
#> ETS Model           -97.03894 135.2118 101.4476 -20.92604 21.57105 0.9537757
#>                    Theil's U
#> ARIMA_Model        10.349318
#> Holt-Winters_Model 11.296285
#> ETS Model           9.442058

📝 Insight :
Berdasarkan hasil akurasi dari ketiga model di atas, yang memiliki nilai MAE terendah yaitu Model ETS sebesar 101.44, dengan nilai MAPE 21.57%.

Assumption Test

No-autocorrelation Residual

Untuk mengecek ada atau tidaknya autokorelasi pada hasil forecasting time series bisa dengan melakukan uji Ljung-box dengan menggunakan fungsi Box.test(residual model, type = "Ljung-Box)

Box.test(model_ets_msts$residuals, type = "Ljung-Box")

#> 
#>  Box-Ljung test
#> 
#> data:  model_ets_msts$residuals
#> X-squared = 0.014894, df = 1, p-value = 0.9029

🌟 Insight :
Dari hasil di atas nilai p-value 0.9029 > alpha(0.05), maka residual pada model tidak terjadi autokorelasi.

Normality of Residual

Untuk mengecek normality residual pada hasil forecasting time series kita bisa melakukan uji normality (shapiro test) dengan menggunakan fungsi shapiro.test(residual model)

shapiro.test(model_ets_msts$residuals)

#> 
#>  Shapiro-Wilk normality test
#> 
#> data:  model_ets_msts$residuals
#> W = 0.91078, p-value < 0.00000000000000022

🌟 Insight :
Pada Uji Normalitas ini nilai p-value 0.00000000000000022 < alpha, maka residual model tidak menyebar normal.

Conclusion

Model time series dengan performa terbaik dapat dilihat dari nilai MAE (Mean Absolute Error), dan yang dapat diinterpretasikan seberapa besar penyimpangan hasil prediksi terhadap nilai aktualnya adalah Model ETS dengan nilai MAE sebesar 101.44 dan dilihat dari nilai MAPE yang dapat menunjukkan seberapa besar penyimpangannya dalam bentuk persentase sebesar 21.57%.
Pada uji asumsi autokorelasi, dari hasil Box-Ljung test Model ETS ini mendapatkan nilai p-value sebesar 0.9029 yang mengindikasikan bahwa pada model tidak terjadi autokorelasi.
Pada uji normalitas, dari hasil Saphiro test Model ETS ini mendapatkan nilai p-value sebesar < 0.00000000000000022 yang mengindikasikan bahwa residual tidak terdistribusi secara normal.