IMPORT LIBRARY

library(dplyr)
## Warning: package 'dplyr' was built under R version 4.3.1
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(lubridate) 
## Warning: package 'lubridate' was built under R version 4.3.1
## 
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
## 
##     date, intersect, setdiff, union
library(tidyr)
## Warning: package 'tidyr' was built under R version 4.3.1
library(ggplot2)
library(TSstudio) 
## Warning: package 'TSstudio' was built under R version 4.3.1
library(forecast)
## Warning: package 'forecast' was built under R version 4.3.1
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(readr)
## Warning: package 'readr' was built under R version 4.3.1
library(padr)
## Warning: package 'padr' was built under R version 4.3.1
library(zoo)
## Warning: package 'zoo' was built under R version 4.3.1
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(ggsci)
## Warning: package 'ggsci' was built under R version 4.3.1

READ DATA

data_train <- read.csv("dataset/F&B-data-train.csv")
data_test <- read.csv("dataset/F&B-data-test.csv")

# DATA INSPECTIONS
head(data_train)
head(data_test)

DATA PREPROCESS

Di bawah ini adalah Langkah Praproses yang perlu dilakukan sebelum melakukan analisis lebih lanjut <

Bulatkan tanggalwaktu menjadi jam, dengan menggunakan tanggal_lantai, simpan hasilnya ke kolom ds <

fnb_data <- data_train %>% 
  mutate(transaction_date = ymd_hms(transaction_date),
         ds = floor_date(transaction_date, unit = "hours")) %>%
  select(transaction_date, receipt_number, ds)

Kolom ds agregat/ringkas (Date Time in hour) dan kolom nomor_tanda terima untuk mendapatkan jumlah pengunjung <

fnb_data <- fnb_data %>%
  group_by(ds) %>%
  summarise(customer_number = n_distinct(receipt_number))

Time Series Padding, untuk mengisi rangkaian data yang hilang (in hour), sehingga setiap hari data rangkaian waktu akan terdiri dari 24 record - 24 jam / per hari. <

fnb_data <- fnb_data %>% pad()
## pad applied on the interval: hour

Mengganti NA Value <

fnb_data <- fnb_data %>%
  mutate(customer_number = na.fill(customer_number, 0))

Check anyNA and tidak ditemukan adanya NA <

anyNA(fnb_data)
## [1] FALSE
colSums(is.na(data_train))
## transaction_date   receipt_number          item_id       item_group 
##                0                0                0                0 
## item_major_group         quantity        price_usd        total_usd 
##                0                0                0                0 
##     payment_type       sales_type 
##                0                0

Filter Tanggal Waktu dalam jam (ds) untuk menangkap rentang waktu Buka Outlet, dari pukul 10.00 hingga 22.00 <

fnb_data <- fnb_data %>% filter(hour(ds) >= 10 & hour(ds) <= 22)
fnb_data

Mengecheck rentan waktu awal dan akhir pada data <

range(fnb_data$ds)
## [1] "2017-12-01 13:00:00 UTC" "2018-02-18 22:00:00 UTC"
min_date <- min(fnb_data$ds)
max_date <- max(fnb_data$ds)

min_date
## [1] "2017-12-01 13:00:00 UTC"
max_date
## [1] "2018-02-18 22:00:00 UTC"

Snapshot Lihat Data Akhir <

head(fnb_data, 30)
tail(fnb_data, 30)

SEASONALITY ANALYSIS

Buat Objek Time Series menggunakan fungsi ts() dengan frekuensi 13, dimana 13 adalah rentang waktu buka outlet setiap hari, mulai pukul 10.00 hingga 22.00. <

fnb_ts <- ts(data = fnb_data$customer_number, frequency = 13)
fnb_ts
## Time Series:
## Start = c(1, 1) 
## End = c(80, 10) 
## Frequency = 13 
##    [1] 16 38 27 29 44 50 66 70 63 63 10 17 18 32 21 40 36 36 41 68 61 62 54  7
##   [25] 13 20 35 23 30 31 44 55 66 47 49 54  5  9 10 13 12 15 19 27 20 39 54 57
##   [49] 52  7 14 13  9 12 16 22 23 35 71 57 54 35  4 14 28 22 29 25 26 30 35 47
##   [73] 48 58 33  2 10  9 14 16 21 20 34 40 60 57 53 41  0  0  0 13 19 24 29 26
##   [97] 45 43 46 54 44 14 21 32 30 30 26 31 32 27 41 74 65 66  8 10 18 16 22 21
##  [121] 25 37 41 50 62 56 46 14 13 13 18 14 18 16 17 26 27 51 48 43  9 12 20 14
##  [145] 13 15 26 27 36 46 52 48 40 12 11 26 25 20 20 21 26 34 58 61 41 40  8  6
##  [169] 18 21 14 23 21 23 30 47 42 62 37  0  0  0 11 22 17 26 21 44 57 40 51 49
##  [193] 13 14 18 27 30 23 27 20 26 51 64 59 43  4 15 23 22 21 36 33 42 39 48 56
##  [217] 58 40  8 11 17 18 11 22 23 28 37 54 61 56 42  7 15 17 13 17 15 21 29 30
##  [241] 42 54 49 39  3 10 16 16 25 17 17 22 33 35 48 49 44 11 11 26 12 18 17 24
##  [265] 26 26 52 32 43 32  0  0  0  4 24 28 17 21 27 39 60 38 52 13 16 16 18 19
##  [289] 25 24 29 27 55 47 60 63  5 11 18 21 17 19 16 11 15 41 68 70 37  5 10 20
##  [313] 32 21 21 22 28 33 31 64 51 41  3 16 12 18 26 18 27 28 33 53 40 42 40  9
##  [337]  5  3 12 10 20 23 19 42 39 47 31 39  2 10 16 17 17 23 15 26 31 34 56 58
##  [361] 34  0  0  0  6 22 16 18 24 32 57 61 52 46  4  8 30 26 21 21 29 38 34 37
##  [385] 66 61 41  5 11 17 24 21 14 21 28 19 38 27 51 63  4 10 31 29 18 23 33 39
##  [409] 37 64 56 58 48 11 18 18 21 24 25 27 31 49 57 54 38 42 12  7 16 25 28 26
##  [433] 13 18 30 61 46 47 36  5  6 19 10 24 17  9 14 26 47 46 51 37  0  0  0  8
##  [457] 21 26 17 25 33 46 65 42 44 11 18 28 26 35 18 38 36 54 53 59 67 51 10 19
##  [481] 17 30 28 37 39 34 29 52 56 47 41  9 12  8 12 21 17 25 24 28 38 51 44 35
##  [505]  9  5  9 17 21 20 18 20 22 45 43 46 44  8 12 12 16 31 22 21 19 32 52 48
##  [529] 40 41  9  9 21 21 21 18 15 10 21 40 44 50 37  0  0  0 13 12 16 14 30 33
##  [553] 53 63 41 37  9 11 22 29 26 24 27 32 44 53 63 61 41 13 11 28 19 26 41 29
##  [577] 43 29 36 43 48 40  6  7 19 17  9 25 24 20 15 34 50 42 39  4 18 21 12 23
##  [601] 24 23 27 28 37 47 41 51  6 14 20  9 18 17 28 28 21 44 49 47 37 13  7 13
##  [625] 16 14 28 10 19 28 51 52 38 32  0  0  0  8 12 15 20 14 30 42 54 60 48  5
##  [649] 15 29 29 20 22 26 27 37 57 66 58 58  7 19 28 42 31 33 41 45 37 52 61 59
##  [673] 42  6 10 15 17 13 17 14 25 28 48 43 37 36  9 10 15 18 19 24 21 35 27 45
##  [697] 30 46 34 12 12 19 10 19 16 24 27 20 31 43 52 35  5  7 18 22 17 18 14 29
##  [721] 33 62 61 60 44  0  0  0  7 20 21 16 16 37 43 65 56 58  5 12 22 25 35 29
##  [745] 33 39 34 54 80 61 50  8 16 21 30 32 32 33 36 41 60 57 50 47 14 11 22 15
##  [769] 20 16 17 25 33 42 44 46 33  5  7 11 19 18 21 19 24 26 39 49 48 42  7 19
##  [793] 20 16 23 20 14 27 31 43 64 46 40  9 10 24 11 14 19 20 19 34 50 58 53 43
##  [817]  0  0  0  3 19 22 20 24 38 46 69 57 49 16 22 28 27 30 31 32 38 38 47 76
##  [841] 61 48  6 19 23 29 38 41 48 36 47 59 83 65 44  5  5 15 10 17 22  8 18 31
##  [865] 24 40 44 37  5 10 18 15  9 15 16 21 27 36 37 45 43  9 16 18 13 14 16 30
##  [889] 28 31 39 45 36 45  9 20 18 15 19 23 28 18 23 43 57 61 47  0  0  0  2 15
##  [913] 19 14 19 34 56 61 64 50  7 23 26 22 28 35 26 34 43 56 67 74 59  6 21 29
##  [937] 26 29 33 25 45 47 60 70 57 53  6 13 20 21 14 22 25 17 24 42 52 49 42  7
##  [961]  9 13 19 17 25 28 25 30 38 45 46 36  6  8 15 11  8 11 17 20 18 33 49 62
##  [985] 49  6  4  6 13 14 14  7 15 19 30 35 40 42  0  0  0  3 21 24 34 29 28 61
## [1009] 70 59 52  2 11 10 24 13 33 31 23 32 55 65 61 67 15 16 21 25 25 27 34 26
## [1033] 31 56 67 49 41

Dekomposisi dan Plotting untuk memeriksa Seasonal, Trend and Irregular<

fnb_ts_dc <- fnb_ts %>% decompose(type = "additive")

fnb_ts %>% tail(13*7*4) %>% decompose() %>% autoplot()

>Dari plot di atas, komponen Estimasi Tren menunjukkan sebuah pola, yang mungkin berisi ekstra musiman yang tidak tertangkap dan ini dapat dianggap sebagai data multimusim. Untuk mengatasi multi-musiman, kita perlu mengkonversi data menjadi “Multi-seasonality Time Series” yang menerima beberapa pengaturan frekuensi.

Kemudian Buat Objek “Multi-seasonality Time Series” menggunakan fungsi msts(), dengan Frequency : Harian (13) dan Mingguan (13*7), ini akan menangkap musiman dalam Harian dan Mingguan. Kemudian dekomposisi dan plot.

fnb_msts <- msts(data = fnb_data$customer_number, seasonal.periods = c(13,13*7))

fnb_msts_dc <- mstl(fnb_msts)

fnb_msts %>% tail(13*7*4) %>% stl(s.window = "periodic") %>% autoplot() 

>Berdasarkan Plot di atas, kini Estimasi Tren dengan “Multiple Seasonality Time Series” menjadi lebih halus dan jelas. Dan juga lebih jelas mengenai Musiman Harian dan Mingguan, lebih mudah dijelaskan untuk analisa lebih lanjut.

Buat Visualisasi Musiman berdasarkan jam dari Objek “Standard Time Series” object ts() dan berikan penjelasannya <

fnb_data %>% 
  mutate(Hour = hour(ds), Seasonal = fnb_ts_dc$seasonal) %>% 
  distinct(Hour, Seasonal) %>% 
  ggplot() +
  geom_bar(aes(x = Hour, y = Seasonal, fill = Seasonal), stat ="identity", width = 0.7)+
  scale_fill_gradient(low = "yellow", high = "red") +
  scale_x_continuous(breaks = seq(0,22,1)) +
  labs(title = "Seasonality Analysis  - Hourly") 

> Berdasarkan plot di atas, rentang Peak Time Outlet adalah antara pukul 19.00 - 22.00, dan pada pukul 20.00 merupakan waktu dimana Pengunjung paling banyak pengunjung datang. Hal tersebut terjadi dikarenakan pada pukul 10.00 sudah memasuki waktu makan siang. Dan pukul 07.00 merupakan waktu yang paling sedikit pengunjung yang datang ke Outlet karena Outlet baru dibuka.

Buat Visualisasi Musiman berdasarkan minggu dan jam dari Objek “Multi-Seasonal Time Series” object msts() dan berikan penjelasannya <

fnb_msts_dc_fr <- data.frame(fnb_msts_dc)

fnb_msts_dc_fr %>%
  mutate(ds = fnb_data$ds) %>% 
  mutate(Day_of_Week  = wday(ds, label = TRUE, abbr = FALSE), Hour = (hour(ds))) %>% 
  group_by(Day_of_Week, Hour) %>%
  summarise(Seasonal = sum(Seasonal13 + Seasonal91)) %>%
  ggplot() +
  geom_bar(aes(x = Hour, y = Seasonal, fill = Day_of_Week),stat ="identity", position = "stack", width = 0.7)+
  scale_x_continuous(breaks = seq(0,22,1)) +
  scale_fill_locuszoom()+
  labs(title = "Multi-Seasonality Analysis  - Weekly & Hourly") 
## `summarise()` has grouped output by 'Day_of_Week'. You can override using the
## `.groups` argument.

>Berdasarkan plot di atas, rentang Peak Time Outlet adalah antara pukul 19.00 - 22.00 setiap hari dalam seminggu, dan pada hari Sabtu pukul 20.00 adalah waktu dimana Pengunjung paling banyak datang. Dan pukul 10.00 - setiap hari dalam seminggu, merupakan waktu dimana pengunjung paling sedikit datang karena Outlet baru saja dibuka.

MODEL FITTING AND EVALUATION

Penyesuaian dan Evaluasi Model akan menggunakan “Multi-Seasonality Time Series”, karena dapat menangkap data multimusim, dan akan memberikan performa model/perkiraan yang lebih baik <

Membagi Data menjadi Tes (1 Minggu = 7 x 13 Jam) dan Train (Total Data tidak termasuk Data Tes) <

test_fnb_msts  <- tail(fnb_msts, 13*7)
train_fnb_msts <- head(fnb_msts, length(fnb_msts)-length(test_fnb_msts))

Karena Data Time Series mengalami Error, Trend dan Seasonal maka model akan dibangun dengan menggunakan : 1). Holt’s Winters Exponential (Triple Exponential Smoothing) 2). RIMA Model <

model_hw_msts <- HoltWinters(train_fnb_msts)
## Warning in HoltWinters(train_fnb_msts): optimization difficulties: ERROR:
## ABNORMAL_TERMINATION_IN_LNSRCH
model_arima_msts <- stlm(train_fnb_msts, method = "arima")

Perkiraan kedua model 1 minggu ke depan <

hw_forecast_msts <- forecast(model_hw_msts, h = 13*7)
arima_forecast_msts <- forecast(model_arima_msts, h=13*7)
  1. Holt Winter Model - Actual vs Estimated <
fnb_msts %>% autoplot(series = "Actual") +
  autolayer(hw_forecast_msts$mean, series = "Predict HW")

  1. Arima Model - Actual vs Estimated <
fnb_msts %>% autoplot(series = "Actual") +
  autolayer(arima_forecast_msts$mean, series = "Predict Arima")

Berikut ringkasan model yang dibuat dalam “Multi-Seasonal Time Series”: <

Evaluate Both Models Performance <

acc_hw_msts <- accuracy(object = hw_forecast_msts$mean, test_fnb_msts)
acc_arima_msts <- accuracy(object = arima_forecast_msts$mean, test_fnb_msts)

result <- rbind(acc_hw_msts, acc_arima_msts)
rownames(result) <- c("Holt Winter Model", "Arima Model")
result
##                          ME     RMSE      MAE  MPE MAPE      ACF1 Theil's U
## Holt Winter Model -2.898399 8.154527 6.451498 -Inf  Inf 0.4680130         0
## Arima Model       -2.115958 7.323714 5.694023  NaN  Inf 0.3436756         0

Berdasarkan Rangkuman Akurasi di atas, Akurasi “ARIMA Model” lebih baik, MAE 5.694023 (lebih rendah dari 6), model ini akan dipilih untuk digunakan dalam Prediksi ==> Model Terbaik “ARIMA Model”.

Data Preparation untuk Plot <

accuracyData <- data.frame(ds= fnb_data$ds %>% tail(13*7),
  actual = as.vector(test_fnb_msts) ,
  HWForecastmsts = as.vector(hw_forecast_msts$mean) ,
  ArimaForecastmsts = as.vector(arima_forecast_msts$mean))

Visualization of Actual vs Estimated number of visitors <

accuracyData %>% 
 ggplot() +
  geom_line(aes(x = ds, y = actual, colour = "Actual"),size=1)+
  geom_line(aes(x = ds, y = ArimaForecastmsts, colour = "Arima Model"),size=1)+
  labs(title = "Hourly Visitor - Actual Vs Arima Model",x = "Date",y = "Customer Number",colour = "")
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

accuracyData %>% 
 ggplot() +
  geom_line(aes(x = ds, y = actual, colour = "Actual"),size=0.5)+
  geom_line(aes(x = ds, y = HWForecastmsts, colour = "Holt Winter Model"),size=0.1)+
  geom_line(aes(x = ds, y = ArimaForecastmsts, colour = "Arima Model"),size=0.5)+
  labs(title = "Hourly Visitor - Actual Vs All Models",x = "Date",y = "Visitor",colour = "")

PREDICTION PERFORMANCE

MAE < 6 in validation dataset

#Model Terbaik adalah “Arima Model” karena memiliki Akurasi MAE 5.694023, MAE terkecil dibandingkan model lain yang dibuat <

accuracy(object = arima_forecast_msts$mean, test_fnb_msts)
##                 ME     RMSE      MAE MPE MAPE      ACF1 Theil's U
## Test set -2.115958 7.323714 5.694023 NaN  Inf 0.3436756         0

MAE < 6 in test dataset > Prediksi menggunakan “Multi-Seasonal Data” dan simpan Prediksi ke dalam File CSV “submission-fnb-test.csv” dengan menggunakan Model dan Prakiraan Terbaik : “Mode ARIMA”. Langkah : <

  1. Baca Test Data
  2. Buat Model “ARIMA” dengan Multi-Seasonal Time Series Data
  3. Perkiraan model 1 minggu ke depan
  4. Simpan hasil Prakiraan ke dalam Test Data Template
  5. Lihat hasil Test Data <
fnb_test  <- read_csv("dataset/F&B-data-test.csv")
## Rows: 91 Columns: 2
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## lgl  (1): visitor
## dttm (1): datetime
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
model_arima_msts_test <- stlm(fnb_msts, method = "arima")
arima_forecast_msts_test <- forecast(model_arima_msts_test, h=13*7)
fnb_test$visitor <- arima_forecast_msts_test$mean
fnb_test
  1. Menekstrak Templat Test Data ke dalam File CSV “submission-fnb-test.csv”<
write.csv(fnb_test,file = "submission-fnb-test.csv")

CONCLUSION

Pemeriksaan asumsi dan ringkasan <

Asumsi Autocorrelation menggunakan fungsi Ljung-box <

Box.test(model_arima_msts$residuals, type = "Ljung-Box")
## 
##  Box-Ljung test
## 
## data:  model_arima_msts$residuals
## X-squared = 7.0953e-05, df = 1, p-value = 0.9933

#Hasil : p-value = 0.9933 > 0,05 (alpha) sehingga dapat disimpulkan bahwa residual tidak terdapat autokorelasi <

Normality residu, diuji menggunakan fungsi saphiro.test <

shapiro.test(model_arima_msts$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  model_arima_msts$residuals
## W = 0.99028, p-value = 6.683e-06

Hasil : p-value = 6.452e-06 < 0.05 (alpha) sehingga dapat disimpulkan bahwa residual tidak berdistribusi normal.

Akar Penyebab dan Cara Penanganannya : Ditemukan bahwa sisa-sisa tidak terdistribusi secara normal, hal ini bisa terjadi karena ukuran data kita jika tidak cukup besar, namun bukan berarti Forecast Performance tidak baik. Kita perlu mendapatkan lebih banyak data untuk dibangun dan dimasukkan ke dalam model, sehingga sisa lebih terdistribusi secara normal dan mendapatkan hasil yang lebih baik pada Evaluasi Model dan Kinerja Perkiraan.

SUMMARY

-> Kisaran Peak Time Outlet adalah antara pukul 19.00 – 22.00 setiap hari dalam seminggu _> Pada hari Sabtu pukul 20.00 merupakan waktu dimana Pengunjung paling banyak datang. -> Pukul 10.00 - setiap hari dalam seminggu, merupakan waktu dimana pengunjung paling sedikit datang karena Outlet baru saja dibuka.