library(quantmod)
## Loading required package: xts
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## Loading required package: TTR
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(fpp3)
## ── Attaching packages ──────────────────────────────────────────── fpp3 0.4.0 ──
## ✔ tibble      3.1.8     ✔ tsibble     1.1.2
## ✔ dplyr       1.0.9     ✔ tsibbledata 0.4.0
## ✔ tidyr       1.2.0     ✔ feasts      0.2.2
## ✔ lubridate   1.8.0     ✔ fable       0.3.1
## ✔ ggplot2     3.3.6
## ── Conflicts ───────────────────────────────────────────────── fpp3_conflicts ──
## ✖ lubridate::date()    masks base::date()
## ✖ dplyr::filter()      masks stats::filter()
## ✖ dplyr::first()       masks xts::first()
## ✖ tsibble::index()     masks zoo::index()
## ✖ tsibble::intersect() masks base::intersect()
## ✖ tsibble::interval()  masks lubridate::interval()
## ✖ dplyr::lag()         masks stats::lag()
## ✖ dplyr::last()        masks xts::last()
## ✖ tsibble::setdiff()   masks base::setdiff()
## ✖ tsibble::union()     masks base::union()
library(readxl)
library(ggplot2)
library(tsibble)
library(lubridate)
library(tidyverse)
## ── Attaching packages
## ───────────────────────────────────────
## tidyverse 1.3.2 ──
## ✔ readr   2.1.2     ✔ stringr 1.4.1
## ✔ purrr   0.3.4     ✔ forcats 0.5.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ lubridate::as.difftime() masks base::as.difftime()
## ✖ lubridate::date()        masks base::date()
## ✖ dplyr::filter()          masks stats::filter()
## ✖ dplyr::first()           masks xts::first()
## ✖ tsibble::intersect()     masks lubridate::intersect(), base::intersect()
## ✖ tsibble::interval()      masks lubridate::interval()
## ✖ dplyr::lag()             masks stats::lag()
## ✖ dplyr::last()            masks xts::last()
## ✖ tsibble::setdiff()       masks lubridate::setdiff(), base::setdiff()
## ✖ tsibble::union()         masks lubridate::union(), base::union()
library(dplyr)
library(seasonal)
## Warning: package 'seasonal' was built under R version 4.2.2
## 
## Attaching package: 'seasonal'
## 
## The following object is masked from 'package:tibble':
## 
##     view
library(zoo)
library(tidyr)
library(tibble)
library(ggfortify)
library(moments)
library(fable)
library(brew)
library(sos)
## 
## Attaching package: 'sos'
## 
## The following object is masked from 'package:tidyr':
## 
##     matches
## 
## The following object is masked from 'package:dplyr':
## 
##     matches
## 
## The following object is masked from 'package:utils':
## 
##     ?
library(slider)
library(x13binary)
library(USgas)
library(GGally)
## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2
#importing/formating/graphing data

pbacon <- read_excel("C:/Users/natha/Desktop/Forecasting in R/Average Price of Bacon.xls")
head(pbacon)
## # A tibble: 6 × 2
##   Date                `Average Price`
##   <dttm>                        <dbl>
## 1 2000-01-01 00:00:00            2.75
## 2 2000-02-01 00:00:00            2.87
## 3 2000-03-01 00:00:00            2.93
## 4 2000-04-01 00:00:00            2.95
## 5 2000-05-01 00:00:00            3.01
## 6 2000-06-01 00:00:00            3.13
pbacont <- pbacon %>%
  mutate(Date = yearmonth(Date)) %>%
  as_tsibble(index = Date)

autoplot(pbacont) +
  labs(title = "Average Monthly Price of Bacon (USD/Pound)")
## Plot variable not specified, automatically selected `.vars = Average Price`

I do not see a cycle in the data. There doesn’t seem to be a significant seasonal component to the price changes, however, there is an upward trend in the data.

#classical additive
stladd <- pbacont %>%
  model(classical_decomposition(`Average Price`, type = "additive")) %>%
  components() %>%
  autoplot() 
stladd
## Warning: Removed 6 row(s) containing missing values (geom_path).

#classical multi
stlmulti <- pbacont %>%
  model(classical_decomposition(`Average Price`, type = "multiplicative")) %>%
  components() %>%
  autoplot() 
stlmulti
## Warning: Removed 6 row(s) containing missing values (geom_path).

Decomposing the data confirms the fact that the data doesn’t have seasonal movement. The graph of the data resembles the graph of a well preforming stock. The lack of seasonality leads me to believe that the Holt-Winter’s method isn’t required, however, I did include it below just to prove this.

train1 <- pbacont %>%
  filter(yearmonth(Date) <= yearmonth("2015-11-01"))

test1 <- pbacont %>%
  filter(yearmonth(Date) > yearmonth("2015-11-01"))
smoothing1 <- train1 %>%
  model(
    Component = ETS(`Average Price` ~ error("A") + trend("N") + season("N")),
    Holt = ETS(`Average Price` ~ error("A") + trend("A") + season("N")),
    Damped = ETS(`Average Price` ~ error("A") + trend("Ad") + season("N")),
    Additive = ETS(`Average Price` ~ error("A") + trend("A") + season("A")),
    Multiplicative = ETS(`Average Price` ~ error("M") + trend("A") + season("M"))) %>%
  forecast(h=4)
accuracy(smoothing1, test1)
## # A tibble: 5 × 10
##   .model         .type     ME  RMSE   MAE   MPE  MAPE  MASE RMSSE  ACF1
##   <chr>          <chr>  <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Additive       Test  -0.308 0.348 0.308 -5.59  5.59   NaN   NaN 0.242
## 2 Component      Test  -0.278 0.308 0.278 -5.04  5.04   NaN   NaN 0.173
## 3 Damped         Test  -0.351 0.384 0.351 -6.37  6.37   NaN   NaN 0.222
## 4 Holt           Test  -0.329 0.362 0.329 -5.97  5.97   NaN   NaN 0.218
## 5 Multiplicative Test  -0.253 0.296 0.253 -4.61  4.61   NaN   NaN 0.178

Surprisingly, the Holt-Winter’s multiplicative method returned the lowest RMSE and MAE values suggesting it is the best method out of the five. I will be using this method to forecast the original dataset.

bacon_forecast <- pbacont %>%
  model(ETS(`Average Price` ~ error("M") + trend("A") + season("M")))
bacon_forecast%>% 
  forecast(h = 4) %>%
  autoplot(pbacont) +
  geom_line(aes(y = .fitted), color = "blue", data = augment(bacon_forecast)) +
  labs(y = "Average Price of Bacon (USD/pound)") +
  guides(color = "none")

bacon_forecast%>% gg_tsresiduals()

augment(bacon_forecast) %>%
  features(.resid, ljung_box, lag = 10, dof = 0)
## # A tibble: 1 × 3
##   .model                                                         lb_stat lb_pv…¹
##   <chr>                                                            <dbl>   <dbl>
## 1 "ETS(`Average Price` ~ error(\"M\") + trend(\"A\") + season(\…    18.5  0.0476
## # … with abbreviated variable name ¹​lb_pvalue
bacon_forecast <- pbacont %>%
  model(ETS(`Average Price` ~ error("M") + trend("A") + season("M")))%>%
  forecast(h = 4) %>%
  hilo()
bacon_forecast
## # A tsibble: 4 x 6 [1M]
## # Key:       .model [1]
##   .model                         Date Average Pri…¹ .mean                  `80%`
##   <chr>                         <mth>        <dist> <dbl>                 <hilo>
## 1 "ETS(`Average Price` ~ er… 2022 Oct N(7.3, 0.036)  7.34 [7.097398, 7.584216]80
## 2 "ETS(`Average Price` ~ er… 2022 Nov N(7.2, 0.064)  7.17 [6.840345, 7.490803]80
## 3 "ETS(`Average Price` ~ er… 2022 Dec N(7.1, 0.093)  7.10 [6.709517, 7.489646]80
## 4 "ETS(`Average Price` ~ er… 2023 Jan  N(7.1, 0.12)  7.10 [6.648480, 7.543256]80
## # … with 1 more variable: `95%` <hilo>, and abbreviated variable name
## #   ¹​`Average Price`

The point forecast for October is $7.34 with an 80% PI of [7.10, 7.58]. The point forecast for November is $7.17 with an 80% PI of [6.84, 7.49]. The point forecast for December and January is $7.10. December’s 80% PI is [6.71, 7.49] and January’s 80% PI [6.65, 7.54]. Dispite this forecast predicting the price to move against the trend of the data, I believe it to be accurate due to the low RMSE and MAE values, white noise in the residuals, and the fitting values matching the true values of the data.